緯度経度からメートル座標系へ変換（ヒュベニの公式）

メートル変換

方法は色々ありますが、今回は計算量の少ないヒュベニの公式を使います。

ヒュベニの公式

ヒュベニの公式は、２つの緯度経度の差をメートル座標系で求める公式です。
そのため、絶対座標は求まりませんが、
基準点を一点に定めることにより、１つの座標系内での位置を求めることが可能です。

公式

\begin{align}
a:&=長半径 &
F:&=逆扁平率 \\
e:&=第一離心率
=\frac{\sqrt{2F-1}}{F} \\
\\
\phi:&=緯度&
\lambda:&=経度 \\
\phi0:&=基準緯度&
\lambda0:&=基準経度 \\
\\
\Delta \phi&=\phi-\phi0 &
\Delta \lambda&=\lambda-\lambda0 \\
\\
\mu&=\frac{\lambda+\lambda0}{2} \\
\\
W&=\sqrt{1-e^2 \sin^2 \mu} \\
N&=\frac{a}{W} &
M&=a \left( \frac{1-e^2}{W^3} \right) \\
\\
\Delta X&=\Delta \lambda \cdot N \cdot \cos \mu &
\Delta Y&=\Delta \phi \cdot M
\end{align}

実装

/// <summary>ラジアン変換</summary>
const double DEG2RAD = Math.PI / 180;

/// <summary>長半径[GRS80,WGS-84]</summary>
const double a = 6378137;

/// <summary>逆扁平率[GRS80]</summary>
const double F = 298.257222101;

/// <summary>第一離心率の二乗</summary>
const double e2 = (2 * F - 1) / F / F;

/// <summary>メートル座標を求める</summary>
/// <param name="lat">緯度（Degree）</param>
/// <param name="lon">経度（Degree）</param>
/// <param name="lat0">基準緯度（Degree）</param>
/// <param name="lon0">基準経度（Degree）</param>
public static (double X, double Y) Convert(double lat, double lon, double lat0, double lon0)
{
    var nlat = (lat + lat0) / 2 * DEG2RAD;
    var dlat = (lat - lat0) * DEG2RAD;
    var dlon = (lon - lon0) * DEG2RAD;

    var (s, c) = Math.SinCos(nlat);

    var w = Math.Sqrt(1 - e2 * s * s);
    var n = a / w;
    var m = a * (1 - e2) / (w * w * w);

    return (dlon * n * c, dlat * m);
}

まとめ

公式と定数さえ分かっていれば結構簡単です。
長半径と扁平率は使う座標系によって変わりますが、他の計算式はそのまま使えます。