地図タイルにおけるピクセルの大きさをメートルに変換する

1. 概要

地図タイルにおける1ピクセルの大きさ（メートル）を計算したい。

地球は平面ではないため、ピクセルの大きさは位置に依存することになる。したがって、タイル座標が与えられた時、当該タイルにおける1ピクセルの大きさを計算できるようにしたい。

2. 計算式

Understanding Scale and Resolution - Bing Maps Articles によると、以下のように緯度とズームレベルに依存した式で計算することができるらしい。

Map resolution = 156543.04 meters/pixel * cos(latitude) / (2 ^ zoomlevel)

この式をタイル座標とピクセル座標の関係、ピクセル座標と緯度経度の関係と組み合わせれば、目的を達成することができる。

3. 実装 (JavaScript)

3.1. 緯度・ズームレベルからピクセルサイズを計算

計算式通り実装する。latitudeの単位が[rad]であることに注意。

// Calculate resolution [meter/pixel]
const resolution = (latitude, zoom) => 156543.04 * Math.cos(latitude) / (1 << zoom)

3.2. タイル座標から緯度経度を計算

ここでは、タイル座標 tileZ / tileX / tileY のタイル上の (pX, pY) 要素の緯度経度を計算する関数を実装している。

// Calculate lat/lon of pixel (pX, pY) on tile tileZ/tileX/tileY
function pixelOnTileToLatLon(pX, pY, tileZ, tileX, tileY) {
  const L = 85.05112878;
  // Pixel coordinate
  const x = 256 * tileX + pX;
  const y = 256 * tileY + pY;

  const lon = 180 * (x / (1 << (tileZ + 7)) - 1);
  const lat = (180/Math.PI) * Math.asin(Math.tanh(
    - Math.PI / (1 << (tileZ + 7)) * y + Math.atanh(Math.sin(L * Math.PI/180))
  ));
  return {lat: lat, lon: lon};
}

3.3. タイル座標からピクセルサイズを計算

上記の関数を組み合わせ、タイル座標で指定されたタイルにおけるピクセルサイズを計算する関数を実装する。
同じタイル内でも緯度によって結果が異なるため、与えられたタイル座標について北西 / 南東 / 中心の緯度経度を計算し、それぞれについて計算する。

// Calculate pixel size [m] of given tile
function pixelSize(tileZ, tileX, tileY) {
  // North west / South east / Center
  const pNW =pixelOnTileToLatLon(0, 0, tileZ, tileX, tileY);
  const pSE = pixelOnTileToLatLon(255, 255, tileZ, tileX, tileY);
  const pCenter = pixelOnTileToLatLon(128, 128, tileZ, tileX, tileY);

  const deg2rad = deg => deg / 180 * Math.PI;

  return {
    min: resolution(deg2rad(pNW.lat), tileZ),
    max: resolution(deg2rad(pSE.lat), tileZ),
    center: resolution(deg2rad(pCenter.lat), tileZ)
  };
}

4. 例

例えばタイル座標 13/7262/3232 （鍋割山付近）についてpixelSize(13, 7262, 3232)を計算すると、1ピクセルあたり15.56 - 15.57 [m] であることが分かる。

{min: 15.564756939473298, max: 15.57165416877383, center: 15.568219450043276}

以上