混乱するので書きながらプログラム作って計算します
プログラミングには情弱初心者御用達言語なでしこを利用します
当方何の素養もない底辺高校生ですので間違いがあるかもしれない点はご容赦ください。
まずヒュベニの公式ってやつなんですが
http://www.kashmir3d.com/kash/manual/std_siki.htm
こちらを参照しますと
『 D=sqrt((M*dP)*(M*dP)+(N*cos(P)*dR)*(N*cos(P)*dR)) 』
というのが書いてあります
D: ２点間の距離(m)
P: ２点の平均緯度
dP: ２点の緯度差
dR: ２点の経度差
M: 子午線曲率半径
N: 卯酉線曲率半径
らしいですね。要するに緯度差経度差を長さに変換したあと三平方の定理使うんだと勝手に解釈してみます。曲率半径を使ってるので地球が楕円体（に近い何か）であることも考慮されてるんだろうなと思う
それで
M=6334834/sqrt((1-0.006674*sin(P)*sin(P))^3)
N=6377397/sqrt(1-0.006674*sin(P)*sin(P))
だそうです
この値の算出方法についてはまた気が向いたら調べることにしようと思う
http://www.jsurvey.jp/huroku6.pdf これの最初の方読めば分かるかも
※追記 こっちのほうがいいかな http://yamadarake.jp/trdi/report000001.html
N*cos(P)ってのはこういうことだと思う↓
分かれ！
ここも参考
http://mononofu.hatenablog.com/entry/20090324/1237894846
それでプログラム書いたよ
デフォルトは東京駅と大阪駅の距離を求める感じになってます
んで、回答が『403420.296929651』『403835.69624737』 403kmだそうです ※修正2014/03/22
googleMAPで測ると403kmくらいなので大体合ってるってかんじだと思います
ところでですね、googlemapはちょっと上の式とは違う地球楕円体を使ってるらしくてですね
そうするとMとNの値がちょっと変わるらしいんですよ
このM、Nの値を出すnazoの数値はベッセル楕円体準拠の旧日本測地系らしい
それが世界測地系になるとGRS80ってのになるらしくて
http://yamadarake.jp/trdi/report000001.html
こいつの指示通りやっていけばGRS80の定数でできて正確になるんだろうけどそれはまたこんどね～ばいばい
追記：プログラムの最後から4行目あたりが間違っていたので修正しました 2014/03/22
プログラミングには
当方何の素養もない底辺高校生ですので間違いがあるかもしれない点はご容赦ください。
まずヒュベニの公式ってやつなんですが
http://www.kashmir3d.com/kash/manual/std_siki.htm
こちらを参照しますと
『 D=sqrt((M*dP)*(M*dP)+(N*cos(P)*dR)*(N*cos(P)*dR)) 』
というのが書いてあります
D: ２点間の距離(m)
P: ２点の平均緯度
dP: ２点の緯度差
dR: ２点の経度差
M: 子午線曲率半径
N: 卯酉線曲率半径
らしいですね。要するに緯度差経度差を長さに変換したあと三平方の定理使うんだと勝手に解釈してみます。曲率半径を使ってるので地球が楕円体（に近い何か）であることも考慮されてるんだろうなと思う
それで
M=6334834/sqrt((1-0.006674*sin(P)*sin(P))^3)
N=6377397/sqrt(1-0.006674*sin(P)*sin(P))
だそうです
この値の算出方法についてはまた気が向いたら調べることにしようと思う
http://www.jsurvey.jp/huroku6.pdf これの最初の方読めば分かるかも
※追記 こっちのほうがいいかな http://yamadarake.jp/trdi/report000001.html
N*cos(P)ってのはこういうことだと思う↓
分かれ！
ここも参考
http://mononofu.hatenablog.com/entry/20090324/1237894846
それでプログラム書いたよ
| 緯度1＝34.701909 経度1＝135.494977 緯度2＝35.681382 経度2＝139.766084 緯度1をDEG2RAD 緯度3はそれ 経度1をDEG2RAD 経度3はそれ 緯度2をDEG2RAD 緯度4はそれ 経度2をDEG2RAD 経度4はそれ P=(緯度3+緯度4)/2 dP=(緯度3-緯度4)の絶対値 dR=(経度3-経度4)の絶対値 ((1-0.006674*(SIN(P))*(SIN(P)))^3)の平方根 soreはそれ M=6334834/sore (1-0.006674*(SIN(P))*(SIN(P)))の平方根 soreはそれ N=6377397/sore ((M*dP)*(M*dP)+(N*(COS(P))*dR)*(N*(COS(P))*dR))の平方根 Dはそれ Dを言う |
デフォルトは東京駅と大阪駅の距離を求める感じになってます
んで、回答が
googleMAPで測ると403kmくらいなので大体合ってるってかんじだと思います
ところでですね、googlemapはちょっと上の式とは違う地球楕円体を使ってるらしくてですね
そうするとMとNの値がちょっと変わるらしいんですよ
このM、Nの値を出すnazoの数値はベッセル楕円体準拠の旧日本測地系らしい
それが世界測地系になるとGRS80ってのになるらしくて
http://yamadarake.jp/trdi/report000001.html
こいつの指示通りやっていけばGRS80の定数でできて正確になるんだろうけどそれはまたこんどね～ばいばい
追記：プログラムの最後から4行目あたりが間違っていたので修正しました 2014/03/22