distance between two points

hubeny_distance.js

/* ヒュベニの公式より世界測地系1984(WGS84)での二点間距離を算出します */
/* http://yamadarake.jp/trdi/report000001.htmlを基に作成 */
function calc_len(lat1,lon1,lat2,lon2){
	// 角度をラジアンに
	lat1 = lat1/180*Math.PI;
	lat2 = lat2/180*Math.PI;
	lon1 = lon1/180*Math.PI;
	lon2 = lon2/180*Math.PI;
	
	// 長半径(赤道半径)
	var a = 6378137.0;
	// 短半径（極半径）
	var b = 6356752.314245;

	// 緯度経度の差
	var dx = Math.abs(lon1 - lon2);
	var dy = Math.abs(lat1 - lat2);

	// 緯度の平均
	var uy = (lat1 + lat2)/2;

	// 第一離心率
	var e = Math.sqrt((Math.pow(a, 2) - Math.pow(b, 2))/Math.pow(a, 2));

	var W = Math.sqrt(1 - Math.pow(e*Math.sin(uy), 2));

	// 卯酉線曲率半径
	var N = a/W;

	// 子午線極率半径
	var M = a*((1 - Math.pow(e, 2))/Math.pow(W, 3));

	// 二点間の距離
	var d = Math.sqrt(Math.pow(dy*M, 2) + Math.pow(dx*N*Math.cos(uy), 2));

	return d;
}

// 追記:これらをまとめて使いやすく見にくくしたらこうなった
function calc_len2(lat1,lon1,lat2,lon2){
	return  Math.sqrt(Math.pow(Math.abs(lat1/180*Math.PI - lat2/180*Math.PI)*(6378137.0*(0.9933056200098024/Math.pow(Math.sqrt(1 - Math.pow( 0.08181919084296535*Math.sin( (lat1/180*Math.PI + lat2/180*Math.PI)/2), 2)), 3))), 2) + Math.pow(Math.abs(lon1/180*Math.PI - lon2/180*Math.PI)*( 6378137.0/Math.sqrt(1 - Math.pow( 0.08181919084296535*Math.sin( (lat1/180*Math.PI + lat2/180*Math.PI)/2), 2)))*Math.cos( (lat1/180*Math.PI + lat2/180*Math.PI)/2), 2));
}

// 関数使用例
var distance = calc_len(34.992831, 135.714261, 34.993661, 135.714672);
console.log(String(distance)," m");