Học từ Facebook #102: Vòng tròn theo mét trên bản đồ Mercator

Trang quảng cáo cho bạn nhắm "mọi người trong bán kính 10 km quanh cửa hàng". Marketplace hiện "tin đăng trong 40 km". Cả hai cần vẽ một vòng tròn bán kính r mét lên bản đồ.

Cách đầu tiên ai cũng nghĩ: đổi mét ra pixel bằng tỉ lệ ở tâm, rồi ctx.arc(x, y, r_px). Ở Hà Nội với bán kính 10 km thì ổn. Ở Oslo với 800 km thì hình tròn đó nhỏ hơn thật ở phía bắc và tâm lệch 18 px (đo ở zoom 4). Ở Svalbard với 1.000 km thì lệch gần 100 px.

Lý do là bản đồ web (Google Maps, OpenStreetMap, Leaflet, Mapbox) dùng phép chiếu Web Mercator, và phép chiếu này kéo giãn mọi thứ theo vĩ độ. Module LeafletUtils của Facebook có lời giải khá đẹp cho chuyện này, kèm một chỗ sai nhỏ ở bán kính dọc.

▶ Mở demo toàn màn hình ↗ · Sửa code trên StackBlitz ↗ · Xem source trên GitHub ↗

#Facebook làm thế nào?

#LeafletUtils (bundle 156.js)

Module là một bộ tiện ích bọc thư viện Leaflet (bản đồ mã nguồn mở). Nó dùng hệ toạ độ L.CRS.EPSG3857 (Web Mercator) và bán kính R = 6378137 m (bán trục lớn WGS84).

Đổi mét sang pixel ở mức zoom z:

metersToPixels(meters, view) {
  const circumference = 2 * Math.PI * R;            // chu vi xích đạo
  return meters / (circumference * Math.abs(Math.cos(view.center.latitude * DEG)) / 2 ** (view.zoom + 8));
}

2^(z+8) là 256·2^z: số pixel của cả thế giới ở mức zoom z (mỗi tile 256 px). Kiểm bằng Node: 10 km ở Hà Nội, zoom 12 là khoảng 280 px.

Vòng tròn thành hộp bao (đổi tên biến):

circleBounds(center, meters) {
  const δ = meters / R * RAD2DEG;                   // bán kính góc, độ
  const yN = project(center.lat + δ), yS = project(center.lat - δ);
  const h = unproject((yN + yS) / 2).lat;           // vĩ độ ở trung điểm y Mercator
  const φ1 = center.lat * DEG, φ2 = h * DEG;
  const dLng = Math.acos((Math.cos(δ * DEG) - Math.sin(φ1) * Math.sin(φ2)) / (Math.cos(φ1) * Math.cos(φ2))) * RAD2DEG;
  return bounds([h - δ, center.lng - dLng], [h + δ, center.lng + dLng]);
}

Thuật toán này giống gần như từng bước với Circle._project trong mã nguồn Leaflet 1.x (Leaflet không nằm trong bundle đã tải, nên đây là so sánh với bản công khai). Hai ý hay:

1. Tâm trên bản đồ không phải tâm trên mặt đất. Vì phía cực bị giãn nhiều hơn, tâm của hình trên bản đồ nằm lệch về phía cực. Lấy trung điểm theo trục y Mercator của hai mép bắc/nam là tâm đúng.
2. Nửa bề rộng bằng định luật cos cầu. Biểu thức acos(...) là độ chênh kinh độ giữa tâm và một điểm ở vĩ độ h cách tâm đúng δ theo vòng lớn. Khi thử bằng số, h trùng tới 9 chữ số thập phân với vĩ độ mà vòng tròn rộng nhất, asin(sin φ / cos δ). Vậy dLng chính là nửa bề rộng lớn nhất. Không có chú thích nào nói điều này; có lẽ người viết biết, có lẽ là may.

projectCircle chiếu hộp bao ra pixel, lấy tâm là bounds.getCenter() và bán kính là nửa chiều rộng, nửa chiều cao của hộp:

const c = latLngToPoint(b.getCenter(), zoom);
const sw = latLngToPoint(b.getSouthWest(), zoom), ne = latLngToPoint(b.getNorthEast(), zoom);
return { center: c, radius: point(Math.round((ne.x - sw.x) / 2), Math.round((sw.y - ne.y) / 2)) };

#Chỗ sai: hộp bao h ± δ

Leaflet dừng ở pixel: bán kính dọc là p.y − top.y (từ trung điểm tới mép bắc đã chiếu), và có nhánh dự phòng khi acos ra NaN. LeafletUtils thì đổi kết quả thành một hộp bao lat/lng để dùng chung cho fitBounds, và ở bước đó sai. Mép bắc và nam thật của vòng tròn là φ ± δ. Code đã chiếu đúng hai mép này để tìm h, nhưng rồi lại dựng hộp h ± δ, tức là dịch cả hộp lên phía cực một đoạn h − φ. Phía cực giãn nhiều hơn, nên chiều cao hộp trên bản đồ lớn hơn chiều cao thật. Tâm thì vẫn đúng (getCenter() trả về h).

Tính lại bằng Node theo đúng công thức, bán kính dọc ở zoom 4 (thế giới 4.096 px):

Ở vĩ độ Việt Nam, bán kính dưới vài trăm km thì không thấy gì. Ở Bắc Âu với bán kính lớn thì vòng tròn cao hơn thật vài phần trăm. Thêm hai giới hạn:

• Vòng tròn chứa cực (φ + δ > 90°): đối số của acos vượt 1, dLng thành NaN. Trên Mercator, một chỏm cầu chứa cực không còn là hình kín: nó là một dải chạy hết chiều ngang bản đồ.
• Gần 85,05°: Web Mercator cắt ở vĩ độ này (để bản đồ vuông). Leaflet kẹp vĩ độ khi chiếu, nên mép bắc bị dồn lại và h bị kéo lệch.

#Phần còn lại của module

• wrapLongitude(lng) = mod(lng + 180, 360) − 180: đưa kinh độ về [-180, 180).
• isViewInBounds thử thêm bản dịch +360° của khung nhìn, để khung nhìn cắt kinh tuyến 180 (Fiji, Alaska) vẫn được coi là nằm trong vùng cho phép.
• getBoundsZoom duyệt zoom tăng dần, chiếu hộp ra pixel, dừng ở mức đầu tiên hộp không còn vừa khung.
• distanceFromLatLng(from, meters, bearing): công thức điểm đích (destination point) trên mặt cầu, đúng chuẩn.

flowchart TD
  A["Tâm φ, λ và bán kính r mét"] --> B["δ = r / R, đổi ra độ"]
  B --> C["Chiếu hai mép φ+δ và φ−δ ra y Mercator"]
  C --> D["h = vĩ độ ở trung điểm y"]
  D --> E["dLng bằng định luật cos cầu tại h"]
  D --> F1["Bản gốc: hộp h ± δ, bán kính dọc lớn hơn thật"]
  C --> F2["Bản sửa: bán kính dọc = nửa khoảng cách hai mép đã chiếu"]
  E --> G["Bán kính ngang = dLng × pixel mỗi độ"]

#Bên trong: Web Mercator kéo giãn theo vĩ độ

Web Mercator chiếu kinh độ thẳng ra trục x (x = (λ + 180)/360 · W) và vĩ độ ra trục y bằng y = ln(tan(π/4 + φ/2)). Phép chiếu này bảo giác (conformal): góc được giữ, nên ở một điểm bất kỳ, tỉ lệ theo x và theo y bằng nhau. Một vòng tròn rất nhỏ vẫn tròn trên bản đồ.

Nhưng tỉ lệ đó thay đổi theo vĩ độ: ở vĩ độ φ, mọi thứ bị phóng to 1 / cos φ lần. Xích đạo ×1, Hà Nội (21°) ×1,07, Oslo (60°) ×2, Svalbard (78°) ×4,8. Greenland trông to bằng châu Phi vì lẽ đó.

Với vòng tròn lớn, phía cực của nó nằm ở vĩ độ cao hơn, bị phóng to nhiều hơn phía xích đạo. Kết quả:

• Hình trên bản đồ không còn là ellipse, mà giống quả trứng: phía cực phình hơn.
• Tâm của hình trên bản đồ lệch về phía cực so với điểm chiếu của tâm thật.
• Một hình tròn vẽ bằng tỉ lệ ở tâm (r / metersPerPixel(φ)) thì nhỏ hơn ở phía cực, to hơn ở phía xích đạo, và đặt sai tâm.

Ellipse có tâm ở trung điểm y, bán kính dọc bằng nửa khoảng cách hai mép, bán kính ngang bằng nửa bề rộng lớn nhất, là xấp xỉ tốt nhất bằng một ellipse: nó chạm đúng bốn điểm cực của hình thật. Muốn đúng tuyệt đối thì vẽ đa giác: tính 128 điểm đích theo 128 hướng, chiếu từng điểm. Demo vẽ cả hai: đường xanh lá là đa giác, đường cam là ellipse của bản viết lại, hai đường gần như trùng nhau.

#Viết lại

Phép chiếu và tỉ lệ:

export const worldSize = (zoom: number) => 256 * 2 ** zoom;

export const clampLat = (lat: number) => Math.max(-MAX_LAT, Math.min(MAX_LAT, lat));

/** y Mercator không thứ nguyên: ln(tan(π/4 + φ/2)) */
export const mercY = (lat: number) => Math.log(Math.tan(Math.PI / 4 + (clampLat(lat) * D2R) / 2));
export const latFromMercY = (y: number) => (2 * Math.atan(Math.exp(y)) - Math.PI / 2) * R2D;

export function project({ lat, lng }: LatLng, zoom: number): Pixel {
  const size = worldSize(zoom);
  return { x: ((lng + 180) / 360) * size, y: (0.5 - mercY(lat) / (2 * Math.PI)) * size };
}

export function unproject({ x, y }: Pixel, zoom: number): LatLng {
  const size = worldSize(zoom);
  return { lat: latFromMercY((0.5 - y / size) * 2 * Math.PI), lng: (x / size) * 360 - 180 };
}

/** Đưa kinh độ về [-180, 180). */
export const wrapLng = (lng: number) => ((((lng + 180) % 360) + 360) % 360) - 180;

/** Mét trên một pixel tại vĩ độ lat: chu vi xích đạo × cos φ / số pixel của thế giới. */
export const metersPerPixel = (lat: number, zoom: number) => (2 * Math.PI * EARTH_RADIUS * Math.cos(lat * D2R)) / worldSize(zoom);

Điểm đích và đa giác vòng tròn thật. Kinh độ không wrap ở đây, mà được giữ liên tục để đa giác không bị "rách" khi cắt kinh tuyến 180:

export function destination(from: LatLng, meters: number, bearingDeg: number): LatLng {
  const d = meters / EARTH_RADIUS;
  const θ = bearingDeg * D2R;
  const φ1 = from.lat * D2R;
  const sinφ2 = Math.sin(φ1) * Math.cos(d) + Math.cos(φ1) * Math.sin(d) * Math.cos(θ);
  const φ2 = Math.asin(Math.max(-1, Math.min(1, sinφ2)));
  const λ = Math.atan2(Math.sin(θ) * Math.sin(d) * Math.cos(φ1), Math.cos(d) - Math.sin(φ1) * sinφ2);
  // không wrap ở đây để đa giác liền mạch khi cắt kinh tuyến 180; người vẽ tự dời ±360
  return { lat: φ2 * R2D, lng: from.lng + λ * R2D };
}

/** Đa giác n đỉnh của vòng tròn trên mặt cầu (đường cong thật trên bản đồ). */
export function geodesicCircle(center: LatLng, meters: number, n = 128): LatLng[] {
  const out: LatLng[] = [];
  for (let i = 0; i < n; i++) out.push(destination(center, meters, (360 * i) / n));
  // kinh độ liên tục: mỗi đỉnh lệch đỉnh trước không quá 180°
  for (let i = 1; i < n; i++) {
    while (out[i].lng - out[i - 1].lng > 180) out[i].lng -= 360;
    while (out[i].lng - out[i - 1].lng < -180) out[i].lng += 360;
  }
  return out;
}

Ellipse xấp xỉ:

export function circleEllipse(center: LatLng, meters: number, zoom: number): Ellipse {
  const δ = (meters / EARTH_RADIUS) * R2D;
  const north = center.lat + δ, south = center.lat - δ;
  const coversPole = north >= 90 || south <= -90;
  const yN = project({ lat: north, lng: center.lng }, zoom).y;
  const yS = project({ lat: south, lng: center.lng }, zoom).y;
  const s = Math.sin(δ * D2R) / Math.cos(center.lat * D2R);
  const dLng = coversPole || s >= 1 ? 180 : Math.asin(s) * R2D;
  const c = project(center, zoom);
  const px = worldSize(zoom) / 360;
  return { cx: c.x, cy: (yN + yS) / 2, rx: dLng * px, ry: Math.abs(yS - yN) / 2, coversPole };
}

/** Bán kính "tròn" ngây thơ: lấy tỉ lệ mét/pixel ở tâm rồi vẽ hình tròn. */
export function naiveCircle(center: LatLng, meters: number, zoom: number): Ellipse {
  const c = project(center, zoom);
  const r = meters / metersPerPixel(center.lat, zoom);
  return { cx: c.x, cy: c.y, rx: r, ry: r, coversPole: false };
}

/**
 * Hình học Web Mercator cho bản đồ ô vuông (tile 256 px, EPSG:3857):
 * đổi lat/lng sang pixel, mét sang pixel theo vĩ độ, điểm đích theo hướng,
 * và vẽ "vòng tròn bán kính r mét" cho đúng (nó không tròn trên bản đồ).
 * Ý tưởng từ LeafletUtils (Facebook), viết lại từ đầu.
 */

export interface LatLng {
  lat: number;
  lng: number;
}

export interface Pixel {
  x: number;
  y: number;
}

/** Bán kính Trái Đất dùng cho Web Mercator (bán trục lớn WGS84). */
export const EARTH_RADIUS = 6378137;
/** Vĩ độ lớn nhất để bản đồ thành hình vuông: atan(sinh(π)). */
export const MAX_LAT = 85.0511287798066;
const D2R = Math.PI / 180;
const R2D = 180 / Math.PI;

export const worldSize = (zoom: number) => 256 * 2 ** zoom;

export const clampLat = (lat: number) => Math.max(-MAX_LAT, Math.min(MAX_LAT, lat));

/** y Mercator không thứ nguyên: ln(tan(π/4 + φ/2)) */
export const mercY = (lat: number) => Math.log(Math.tan(Math.PI / 4 + (clampLat(lat) * D2R) / 2));
export const latFromMercY = (y: number) => (2 * Math.atan(Math.exp(y)) - Math.PI / 2) * R2D;

export function project({ lat, lng }: LatLng, zoom: number): Pixel {
  const size = worldSize(zoom);
  return { x: ((lng + 180) / 360) * size, y: (0.5 - mercY(lat) / (2 * Math.PI)) * size };
}

export function unproject({ x, y }: Pixel, zoom: number): LatLng {
  const size = worldSize(zoom);
  return { lat: latFromMercY((0.5 - y / size) * 2 * Math.PI), lng: (x / size) * 360 - 180 };
}

/** Đưa kinh độ về [-180, 180). */
export const wrapLng = (lng: number) => ((((lng + 180) % 360) + 360) % 360) - 180;

/** Mét trên một pixel tại vĩ độ lat: chu vi xích đạo × cos φ / số pixel của thế giới. */
export const metersPerPixel = (lat: number, zoom: number) => (2 * Math.PI * EARTH_RADIUS * Math.cos(lat * D2R)) / worldSize(zoom);

/** Khoảng cách vòng lớn (haversine), mét. */
export function distance(a: LatLng, b: LatLng) {
  const dLat = (b.lat - a.lat) * D2R;
  const dLng = (b.lng - a.lng) * D2R;
  const h = Math.sin(dLat / 2) ** 2 + Math.cos(a.lat * D2R) * Math.cos(b.lat * D2R) * Math.sin(dLng / 2) ** 2;
  return 2 * EARTH_RADIUS * Math.asin(Math.min(1, Math.sqrt(h)));
}

/** Điểm đích: đi `meters` theo hướng `bearingDeg` (0 = bắc, 90 = đông) trên mặt cầu. */
export function destination(from: LatLng, meters: number, bearingDeg: number): LatLng {
  const d = meters / EARTH_RADIUS;
  const θ = bearingDeg * D2R;
  const φ1 = from.lat * D2R;
  const sinφ2 = Math.sin(φ1) * Math.cos(d) + Math.cos(φ1) * Math.sin(d) * Math.cos(θ);
  const φ2 = Math.asin(Math.max(-1, Math.min(1, sinφ2)));
  const λ = Math.atan2(Math.sin(θ) * Math.sin(d) * Math.cos(φ1), Math.cos(d) - Math.sin(φ1) * sinφ2);
  // không wrap ở đây để đa giác liền mạch khi cắt kinh tuyến 180; người vẽ tự dời ±360
  return { lat: φ2 * R2D, lng: from.lng + λ * R2D };
}

/** Đa giác n đỉnh của vòng tròn trên mặt cầu (đường cong thật trên bản đồ). */
export function geodesicCircle(center: LatLng, meters: number, n = 128): LatLng[] {
  const out: LatLng[] = [];
  for (let i = 0; i < n; i++) out.push(destination(center, meters, (360 * i) / n));
  // kinh độ liên tục: mỗi đỉnh lệch đỉnh trước không quá 180°
  for (let i = 1; i < n; i++) {
    while (out[i].lng - out[i - 1].lng > 180) out[i].lng -= 360;
    while (out[i].lng - out[i - 1].lng < -180) out[i].lng += 360;
  }
  return out;
}

export interface Ellipse {
  cx: number;
  cy: number;
  rx: number;
  ry: number;
  /** vòng tròn chứa cực: không còn là hình kín trên Mercator */
  coversPole: boolean;
}

/**
 * Ellipse xấp xỉ vòng tròn r mét, tính bằng pixel ở mức zoom.
 * - Theo trục y: lấy đúng hai mép bắc/nam (φ ± δ) rồi chiếu, tâm là trung điểm y.
 * - Theo trục x: độ rộng kinh độ lớn nhất của chỏm cầu, asin(sin δ / cos φ).
 */
export function circleEllipse(center: LatLng, meters: number, zoom: number): Ellipse {
  const δ = (meters / EARTH_RADIUS) * R2D;
  const north = center.lat + δ, south = center.lat - δ;
  const coversPole = north >= 90 || south <= -90;
  const yN = project({ lat: north, lng: center.lng }, zoom).y;
  const yS = project({ lat: south, lng: center.lng }, zoom).y;
  const s = Math.sin(δ * D2R) / Math.cos(center.lat * D2R);
  const dLng = coversPole || s >= 1 ? 180 : Math.asin(s) * R2D;
  const c = project(center, zoom);
  const px = worldSize(zoom) / 360;
  return { cx: c.x, cy: (yN + yS) / 2, rx: dLng * px, ry: Math.abs(yS - yN) / 2, coversPole };
}

/** Bán kính "tròn" ngây thơ: lấy tỉ lệ mét/pixel ở tâm rồi vẽ hình tròn. */
export function naiveCircle(center: LatLng, meters: number, zoom: number): Ellipse {
  const c = project(center, zoom);
  const r = meters / metersPerPixel(center.lat, zoom);
  return { cx: c.x, cy: c.y, rx: r, ry: r, coversPole: false };
}

/**
 * Cách của LeafletUtils.circleBounds + projectCircle (viết lại để so sánh):
 * tìm vĩ độ h ở trung điểm y Mercator, nhưng rồi dựng hộp bao h ± δ (thay vì φ ± δ)
 * và lấy nửa chiều cao của hộp đó làm bán kính y.
 */
export function boundsShiftEllipse(center: LatLng, meters: number, zoom: number): Ellipse & { midLat: number } {
  const δ = (meters / EARTH_RADIUS) * R2D;
  const midLat = latFromMercY((mercY(center.lat + δ) + mercY(center.lat - δ)) / 2);
  const φ1 = center.lat * D2R, φ2 = midLat * D2R;
  const dLng = Math.acos((Math.cos(δ * D2R) - Math.sin(φ1) * Math.sin(φ2)) / (Math.cos(φ1) * Math.cos(φ2))) * R2D;
  const sw = project({ lat: midLat - δ, lng: center.lng - dLng }, zoom);
  const ne = project({ lat: midLat + δ, lng: center.lng + dLng }, zoom);
  const c = project({ lat: midLat, lng: center.lng }, zoom);
  return { cx: c.x, cy: c.y, rx: Math.round((ne.x - sw.x) / 2), ry: Math.round((sw.y - ne.y) / 2), coversPole: Number.isNaN(dLng), midLat };
}

/** Mức zoom lớn nhất mà hộp [sw, ne] vẫn vừa khung w×h (duyệt tăng dần như getBoundsZoom). */
export function fitZoom(sw: LatLng, ne: LatLng, w: number, h: number, minZoom = 0, maxZoom = 20) {
  let z = minZoom;
  for (; z <= maxZoom; z++) {
    const a = project(sw, z), b = project(ne, z);
    let dx = b.x - a.x;
    if (dx < 0) dx += worldSize(z); // hộp cắt kinh tuyến 180
    if (dx > w || a.y - b.y > h) break;
  }
  return Math.max(minZoom, z - 1);
}

#Đi qua code

• mercY kẹp vĩ độ về ±85,0511° (bằng atan(sinh π)), giống Leaflet và Google Maps. Không kẹp thì tan(π/2) ra vô cực.
• project trả pixel thế giới ở mức zoom z: góc trên trái là (0, 0), xích đạo ở giữa. Muốn toạ độ trên màn hình thì trừ pixel của góc trên trái khung nhìn.
• destination là công thức điểm đích chuẩn: φ2 = asin(sin φ1 cos d + cos φ1 sin d cos θ). Math.max(-1, Math.min(1, …)) chặn sai số làm tròn vượt miền của asin.
• geodesicCircle sửa kinh độ để hai đỉnh liền nhau không lệch quá 180°. Người vẽ tự vẽ thêm bản sao dời ±1 thế giới (demo làm vậy với Fiji).
• circleEllipse: ry là nửa khoảng cách y của hai mép thật φ ± δ; cy là trung điểm của chúng. rx dùng asin(sin δ / cos φ), công thức bề rộng lớn nhất của chỏm cầu. Khi vòng tròn chứa cực (sin δ / cos φ ≥ 1) thì trả coversPole: true và bề rộng 180° để người gọi biết đây là một dải.
• naiveCircle giữ lại để so sánh trong demo.

Khác bản gốc:

• Bán kính dọc lấy từ hai mép φ ± δ đã chiếu, không dựng lại hộp h ± δ.
• Bán kính ngang dùng asin(sin δ / cos φ): cùng kết quả với acos(...) của bản gốc, nhưng không cần h và không ra NaN khi chạm cực (trả về cờ thay vì NaN).
• Không làm tròn bán kính ra số nguyên (bản gốc Math.round), để vẽ trên canvas có khử răng cưa không bị nhảy 1 px khi kéo.

#Tại sao không làm cách khác?

• Hình tròn theo tỉ lệ ở tâm: đủ tốt khi bán kính nhỏ (dưới vài chục km) hoặc ở gần xích đạo. Bảng sai số trong demo cho thấy lúc nào không đủ.
• Đa giác điểm đích: đúng nhất, cũng là cách Turf.js (turf.circle) làm. Đổi lại phải chiếu 64 đến 128 điểm mỗi lần vẽ và phải tự xử lý kinh tuyến 180. Facebook chọn ellipse có lẽ vì vẽ bằng L.Circle/SVG ellipse rẻ và dễ chỉnh khi kéo.
• Dùng thẳng L.Circle của Leaflet: nếu chỉ cần vẽ thì đây là lựa chọn tốt, nó làm đúng như bản viết lại (tâm ở trung điểm y, bán kính dọc từ mép đã chiếu). Facebook cần thêm hộp bao lat/lng để tính zoom vừa khung, nên tự viết lại.
• Tính khoảng cách trên mặt phẳng Mercator: sai nhanh khi đi xa xích đạo. Đo khoảng cách luôn dùng haversine hoặc Vincenty.

#Edge case và cạm bẫy

• Kinh tuyến 180: tâm ở Fiji (178°E) với bán kính 500 km thì nửa vòng tròn nằm ở 180°W. Phải vẽ hai lần (dời ±1 thế giới) hoặc dùng kinh độ "không wrap" và để bản đồ tự lặp.
• Cực và vĩ độ trên 85°: Mercator không hiển thị được. Nhắm quanh Svalbard với bán kính lớn thì hình bị cắt ở mép trên.
• Tỉ lệ lấy theo tâm khung nhìn: metersToPixels của bản gốc lấy vĩ độ của tâm khung nhìn, không phải của vòng tròn. Vòng tròn ở mép trên của một khung nhìn rộng sẽ bị tính sai tỉ lệ. Luôn dùng vĩ độ của chính đối tượng.
• Trái Đất không tròn: WGS84 là ellipsoid, bán kính xích đạo và bán kính cực lệch nhau 21 km. Dùng hình cầu bán kính 6.378.137 m sai tới khoảng 0,5% khoảng cách: đủ cho hiển thị, không đủ cho đo đạc.
• Ngưỡng pháp lý: nếu bán kính dùng để kiểm luật (như ngưỡng 15 dặm cho quảng cáo nhà ở trong bài #104), hãy kiểm trên số mét, đừng trên pixel.

#Dùng trong dự án

import { circleEllipse, project, worldSize } from './mercator';

/** Vẽ vòng tròn r mét lên canvas khớp với tile bản đồ đang hiện. */
function drawRadius(ctx: CanvasRenderingContext2D, view: { center: LatLng; zoom: number; width: number; height: number }, c: LatLng, meters: number) {
  const origin = project(view.center, view.zoom);
  const left = origin.x - view.width / 2, top = origin.y - view.height / 2;
  const e = circleEllipse(c, meters, view.zoom);
  for (const shift of [-worldSize(view.zoom), 0, worldSize(view.zoom)]) {
    ctx.beginPath();
    ctx.ellipse(e.cx - left + shift, e.cy - top, e.rx, e.ry, 0, 0, Math.PI * 2);
    ctx.stroke();
  }
}

Cần đúng hình (vòng tròn rất lớn, gần cực) thì thay ellipse bằng geodesicCircle(c, meters).map((p) => project(p, zoom)).

#Kiểm thử

import { describe, expect, it } from 'vitest';
import { boundsShiftEllipse, circleEllipse, destination, distance, geodesicCircle, metersPerPixel, project, wrapLng } from './mercator';

const truth = (c: { lat: number; lng: number }, m: number, z: number) => {
  const pts = geodesicCircle(c, m, 720).map((p) => project(p, z));
  const xs = pts.map((p) => p.x), ys = pts.map((p) => p.y);
  return { rx: (Math.max(...xs) - Math.min(...xs)) / 2, ry: (Math.max(...ys) - Math.min(...ys)) / 2 };
};

describe('mercator', () => {
  it('10 km ở Hà Nội, zoom 12 ≈ 280 px', () => {
    expect(10000 / metersPerPixel(21.03, 12)).toBeCloseTo(280, 0);
  });

  it('điểm đích đi đúng khoảng cách', () => {
    const a = { lat: 59.9, lng: 10.7 };
    expect(distance(a, destination(a, 500_000, 37))).toBeCloseTo(500_000, -1);
  });

  it('ellipse khớp đa giác thật, bản gốc lệch ở vĩ độ cao', () => {
    const c = { lat: 70, lng: 0 };
    const t = truth(c, 1_000_000, 4);
    const e = circleEllipse(c, 1_000_000, 4);
    expect(e.rx).toBeCloseTo(t.rx, 0);
    expect(e.ry).toBeCloseTo(t.ry, 0);
    expect(boundsShiftEllipse(c, 1_000_000, 4).ry / t.ry).toBeGreaterThan(1.1);
  });

  it('chứa cực: bản gốc NaN, bản sửa báo cờ', () => {
    expect(circleEllipse({ lat: 88, lng: 0 }, 300_000, 3).coversPole).toBe(true);
    expect(Number.isNaN(boundsShiftEllipse({ lat: 88, lng: 0 }, 300_000, 3).rx)).toBe(true);
  });

  it('wrap kinh độ', () => {
    expect(wrapLng(190)).toBe(-170);
    expect(wrapLng(-180)).toBe(-180);
  });
});

#So với thư viện và chuẩn web

• Leaflet L.Circle: bán kính theo mét, Leaflet 1.x tự đổi ra ellipse pixel bằng chính thuật toán trung điểm y + định luật cos cầu, có dự phòng khi NaN.
• Turf.js circle(center, radius, { steps }): trả về đa giác GeoJSON bằng điểm đích, đúng hình nhưng phải tự xử lý kinh tuyến 180.
• Mapbox GL / MapLibre: không có "circle theo mét" sẵn; layer circle có bán kính theo pixel. Cách phổ biến là dùng đa giác Turf.
• Google Maps google.maps.Circle: vẽ đúng hình trên mặt cầu (bản đồ tự chia đa giác).
• Geolocation API trả về accuracy theo mét: vòng tròn "vị trí của bạn" chính là bài toán này, thường nhỏ nên hình tròn ngây thơ là đủ.

Bài này thuộc series Học từ mã nguồn Facebook (106 bài, đọc theo thứ tự).

Code + demo: examples/facebook-patterns · Mở trên StackBlitz