緯度経度間の直線距離は、地球が楕円型の所為でそれほど平面上の座標の距離ほど簡単に算出できない。
そのため、緯度軽度間の距離を計算するにはいくつかの方法があるらしいが、
最初に見つけた以下のJavaのサンプルコードが私の用途に削ぐわない結果が出たので、
別の方法を調べた。
それがヒュベニの公式(Hybeny's Distance Formula)を使った方法である。
こちらの実装例を紹介するページは多かったが、Javaはぱっと見つからなかったので、
書いてみた。
尚、ヒュベニの公式を使った距離の精度は近い位置だとそれらしい結果だったが、
距離が離れていると誤差が大きそうだった。
とはいえ、私の用途だと大きな問題なさそうだ。
double distance(final double latitude1,
final double longitude1,
final double latitude2,
final double longitude2) {
// このあたりはstatic化できる
final double POLE_RADIUS = 6356752.314245; // 極半径(短半径)
final double EQUATORIAL_RADIUS = 6378137.0; // 赤道半径(長半径)
final double ECCENTRICITY2 =
(Math.pow(EQUATORIAL_RADIUS, 2) - Math.pow(POLE_RADIUS, 2))
/ Math.pow(EQUATORIAL_RADIUS, 2); // 離心率の2乗
// ラジアンへ変換
final double rlat1 = Math.toRadians(latitude1);
final double rlng1 = Math.toRadians(longitude1);
final double rlat2 = Math.toRadians(latitude2);
final double rlng2 = Math.toRadians(longitude2);
final double latDiff = rlat1 - rlat2; // 緯度差
final double lngDiff = rlng1 - rlng2; // 経度差
final double latAvg = (rlat1 + rlat2) / 2; // 緯度経度平均
final double w = Math.sqrt(1 - ECCENTRICITY2 * Math.pow(Math.sin(latAvg), 2));
final double m = EQUATORIAL_RADIUS * (1 - ECCENTRICITY2) / Math.pow(w, 3); // 子午線曲率半径
final double n = EQUATORIAL_RADIUS / w; // 卯酉線曲率半径
return Math.sqrt(Math.pow(m * latDiff, 2)
+ Math.pow(n * lngDiff * Math.cos(latAvg), 2));
}