測地線航海算法(Geodesic Sailing)
 球面や回転楕円体面等の曲面上における任意の2地点間の最短経路を測地線
(Geodesic)と言う。ここでは地球を回転楕円体と考え、測地線に沿って航行する
場合の航海算法について紹介する。

1.測地線長の計算(Calculation of Geodetic Distance)
  回転楕円体上の2地点間の最短距離(測地線長:Geodetic Distance or Geodesic 
Distance)を計算する方法は何種類か存在するが、一般にLambertとAndoyerにより
考案された公式がよく用いられている。ここでは、測地線長計算公式について紹介する。
[1] Lambert-Andoyer Method for Computation of Geodetic Distance
 A点(lA,LA)と B点(lB,LB)の測地線長(ρ)を求める。
(lA,LA): A点の測地緯度,測地経度  (lB,LB): B点の測地緯度,測地経度

(1) 測地緯度を化成緯度に変換
 φ=tan-1(B/A・tanl)                                                 
 φ:化成緯度(Parametric Latitude), l:測地緯度(Geodetic Latitude)
  A:地球赤道半径(Semi-major Axis of the Earth)
  B:地球極半径(Semi-minor Axis of the Earth)

(2) 球面上の距離(Spherical Distance:X)の計算
 X=cos-1[sinφA・sinφB+cosφA・cosφB・cos(LA-LB)]
  φA: A点の化成緯度  φA: B点の化成緯度

(3) Lambert-Andoyer Correction(Δρ)
    Δρ=F/8*{(sinX-X)*(sinφA+sinφB)2/cos2(X/2)
             -(sinX+X)*(sinφA-sinφB)2/sin2(X/2)}
      F=(A-B)/A    F:偏平率(Flattening of the Earth)

(4) 測地線長(ρ:Geodetic Distance)
   ρ=A・(X+Δρ)    ρ:測地線長

[2] 小野の公式
Lambert-Andoyerの式を若干変形したものであり,基本的にはLambert-Andoyerの式と
同じである。
 A点(lA,LA)と B点(lB,LB)の測地線長(ρ)を求める。
(lA,LA): A点の測地緯度,測地経度  (lB,LB): B点の測地緯度,測地経度

(1) 測地緯度を化成緯度に変換
 φ=tan-1(B/A・tanl)                                                  
 φ:化成緯度(Parametric Latitude), l:測地緯度(Geodetic Latitude)
  A:地球赤道半径(Semi-major Axis of the Earth)
  B:地球極半径(Semi-minor Axis of the Earth)

(2) 球面上の距離(Spherical Distance:X)の計算
 X=cos-1[sinφA・sinφB+cosφA・cosφB・cos(LA-LB)]
  φA: A点の化成緯度  φB: B点の化成緯度

(3) 測地線長(ρ:Geodetic Distance)
   ρ=A・X-C・P-D・Q     ρ:測地線長
   C=(sinφA+sinφB)2, D=(sinφA-sinφB)2
   P=(A-B)・(X-sinX)/[4(1+cosX)], Q=(A-B)・(X+sinX)/[4(1-cosX)]


<参考文献> 
(1)Lambert,W.D.,The distance between two widely separated points on the 
surface of the earth, Jour. Washington Acad. Sci.,32,125-130,1942.
(2) 小野房吉:電波航法の新しい測位原理と測位精度の評価,航海,第79号,P.35-40,
1984.(この文献の方位角計算式は誤りである.)
(3) 平岩節:航海算法の修正について−U.,航海,第81号,P.8-14,1984
(4) 石田正一:大圏航海算法の精度について,航海,第123号,P.35-40,1995

2.方位角の計算(Calculation of Azimuth)
 地球楕円体上において A点(lA,LA)から B点(lB,LB)を見たときの方位角の
計算方法について次に述べる。

(1) 2地点の測地緯度を化成緯度に変換
 φ=tan-1(BE/AE・tanl)
 φ:化成緯度(Reduced Latitude), l:測地緯度(Geodetic Latitude)
  AE:地球赤道半径(Semi-major Axis of the Earth)
  BE:地球極半径(Semi-minor Axis of the Earth)

(2) 方位角の計算
 Zc=tan-1[sinh/(cosφAtanφB-sinφAcosh)] 
 h=LA-LB (0≦h<360゚とする.もし h<0 の時は h=h+360゚)
φA:A点の化成緯度,φB:B点の化成緯度,LA,LB:A点,B点の測地経度
 [Zcの前後に付ける符号]
  a)h<180゚ならW かつ Zc>0 ならN,Zc<0 ならS
  b)h>180゚ならE かつ Zc>0 ならS,Zc<0 ならN

 Zcの前に N又は S,後に E又は Wの符号を付けたものが, A点(lA,LA)から
 B点(lB,LB)を見たときの方位角Z である。


[問題] A点(43゚35'55.075"N, 142゚26'58.607"E)と B点(43゚03'51.223"N,
 144゚47'40.571"E)の測地線長と A点から B点を見たときの方位角を求めよ。
(地球赤道半径 A=6377397.155m, 地球偏平率 F=1/299.152813)
    答 測地線長 199.202 Km,出発針路 106.6゚, 到着針路 108.2゚
      (大圏距離 198.618 Km,出発針路 106.6゚, 到着針路 108.2゚)

☆補足
 一般的に、楕円体上の2地点間の測地線長と方位角は、Jordan法を用いて
計算されている。