関数xy2latlon マニュアル

(The documentation of function xy2latlon)
Last Update: 2023/3/30


◆機能・用途(Purpose)

平面直角座標で与えられた座標を緯度経度に変換する。
Convert a coordinate given by a cartesian coordinate to the latitude and longitude.

◆形式(Format)

#include <coordinate.h>
inline void xy2latlon
(const double x,const double y,const double refN,const double refE,
 double ∗N,double ∗E)


◆引数(Arguments)

x 変換したい地点の\(x\)座標(基準点から東方向、m)。
The \(x\)-coordinate of the point to convert; east (m) from the reference point.
y 変換したい地点の\(y\)座標(基準点から北方向、m)。
The \(y\)-coordinate of the point to convert; north (m) from the reference point.
refN 投影の基準点の緯度(\(^{\circ}\))。
The latitude (\(^{\circ}\)) of the reference point of the projection.
refE 投影の基準点の経度(\(^{\circ}\))。
The longitude (\(^{\circ}\)) of the reference point of the projection.
N 変換後の緯度(\(^{\circ}\))の代入先。 宣言しただけのdouble型変数を&を付けて与える。
The memory into which the latitude (\(^{\circ}\)) computed by the conversion will be inserted. Give an empty double-type variable with &.
E 変換後の\(y\)座標(基準点から北方向、m)の代入先。 宣言しただけのdouble型変数を&を付けて与える。
The memory into which the longitude (\(^{\circ}\)) computed by the conversion will be inserted. Give an empty double-type variable with &.


◆使用例(Example)

double N,E;
xy2latlon(100.0,200.0,35.89278,137.48028,&N,&E);


◆用いている計算式(Formula used)

計算式は河瀬(2011)に基づく。 なお元の式は\(x\)が北方向、\(y\)が東方向になっているが、 本プログラムでは\(x\)を東方向、\(y\)を北方向に取る。
The formula is based on Kawase (2011). The original formula defines \(x\) and \(y\) directions to be north and east, respectively, whereas this program uses \(x\) to be east and \(y\) to be north.
このプログラムでは表1に示す定数と表2に示す変数を用いる。 なお楕円体としてはWGS84で用いられているITRF座標系GRS80の値を使用した。
This program uses constants listed in Table 1 and variables listed in Table 2. Here the GRS80 ellipsoid in the ITRF coordinate for WGS84 is used.

表1. このプログラムで使用する定数 (Constants used in this program)
記号
Notation 		意味
Description 		値
Value
\(a\) 楕円体の長半径(m)
The longer radius (m) of the ellipsoid 		6378137.0
\(F\) 楕円体の扁平率の逆数
The inverse of the ellipticity of the ellipsoid 		298.257222101
\(m_0\) 縮尺係数
The scale factor 		0.9999


表2. このプログラムで使用する変数 (Variables used in this program)
記号
Notation 		意味
Description
\(x\) 変換したい地点の\(x\)座標(基準点から東方向, m)
The \(x\)-coordinate of the point to convert; east (m) from the reference point
\(y\) 変換したい地点の\(y\)座標(基準点から北方向, m)
The \(y\)-coordinate of the point to convert; north (m) from the reference point
\(\phi_0\) 投影の基準点の緯度(rad)
The latitude (rad) of the reference point of the projection
\(\lambda_0\) 投影の基準点の経度(rad)
The longitude (rad) of the reference point of the projection


これらの定数及び変数を用いて \((x,y)\)に対応する地点の緯度\(\phi\), 経度\(\lambda\)は 以下の式で計算できる。 式の中での角度は全てrad単位とする。
Using these constants and variables, the latitude \(\phi\) and longitude \(\lambda\) corresponding to the point \((x,y)\) are computed by the formula below. Note that all angles in the formula are in the radian unit.

●\((x,y)\)によらない定数 (Constants independent of the latitude and longitude)

\[\begin{equation} n=\frac{1}{2F-1} \end{equation}\] \[\begin{equation} \beta_1= \frac{1}{2}n-\frac{2}{3}n^2+\frac{37}{96}n^3 -\frac{1}{360}n^4-\frac{81}{512}n^5 \end{equation}\] \[\begin{equation} \beta_2= \frac{1}{48}n^2+\frac{1}{15}n^3-\frac{437}{1440}n^4+\frac{46}{105}n^5 \end{equation}\] \[\begin{equation} \beta_3=\frac{17}{480}n^3-\frac{37}{840}n^4-\frac{209}{4480}n^5 \end{equation}\] \[\begin{equation} \beta_4=\frac{4397}{161280}n^4-\frac{11}{504}n^5 \end{equation}\] \[\begin{equation} \beta_5=\frac{4583}{161280}n^5 \end{equation}\] \[\begin{equation} \delta_1= 2n-\frac{2}{3}n^2-2n^3+\frac{116}{45}n^4 +\frac{26}{45}n^5-\frac{2854}{675}n^6 \end{equation}\] \[\begin{equation} \delta_2= \frac{7}{3}n^2-\frac{8}{5}n^3-\frac{227}{45}n^4 +\frac{2704}{315}n^5+\frac{2323}{945}n^6 \end{equation}\] \[\begin{equation} \delta_3= \frac{56}{15}n^3-\frac{136}{35}n^4 -\frac{1262}{105}n^5+\frac{73814}{2835}n^6 \end{equation}\] \[\begin{equation} \delta_4= \frac{4279}{630}n^4-\frac{332}{35}n^5-\frac{399572}{14175}n^6 \end{equation}\] \[\begin{equation} \delta_5=\frac{4174}{315}n^5-\frac{144838}{6237}n^6 \end{equation}\] \[\begin{equation} \delta_6=\frac{601676}{22275}n^6 \end{equation}\]

●関数 (Function)

\[\begin{equation} S(\phi)=\frac{a}{1+n}\sum_{j=0}^5S_{1j}^2S_{2j}(\phi) \end{equation}\] \[\begin{equation} S_{1j}=\prod_{k=1}^j\left(\frac{3n}{2k}-n\right) \end{equation}\] \[\begin{equation} S_{2j}(\phi)=\phi+\sum_{l=1}^{2j}S_{3jl} \end{equation}\] \[\begin{equation} S_{3jl}=\left(\frac{1}{l}-4l\right)\sin(2l\phi)\prod_{m=0}^lS_{4jm} \end{equation}\] \[\begin{equation} S_{4jm}= \begin{cases} \frac{3n}{2j+m}-n & (m\mbox{が偶数}\textcolor[named]{Gray}{\mbox{ / for even $m$}}) \\ \frac{1}{\frac{3n}{2j-(m-1)}-n} & (m\mbox{が奇数}\textcolor[named]{Gray}{\mbox{ / for odd $m$}}) \end{cases} \end{equation}\] \[\begin{equation} S_p=S\left(\frac{\pi}{2}\right) \end{equation}\] \[\begin{equation} S_{\phi_0}=S\left(\phi_0\right) \end{equation}\]

●緯度経度を用いて計算される値 (Values computed from the latitude and longitude)

\[\begin{equation} \eta=\frac{\pi}{2m_0S_p}x \end{equation}\] \[\begin{equation} \xi=\frac{\pi}{2m_0S_p}\left(y+m_0S_{\phi_0}\right) \end{equation}\] \[\begin{equation} \eta’= \eta-\sum_{j=1}^5 \beta_j\cos(2j\xi)\sinh(2j\eta) \end{equation}\] \[\begin{equation} \xi’= \xi-\sum_{j=1}^5 \beta_j\sin(2j\xi)\cosh(2j\eta) \end{equation}\] \[\begin{equation} \chi=\sin^{-1}\left(\frac{\sin\xi’}{\cosh\eta’}\right) \end{equation}\]

● 緯度\(\phi\)の計算式 (Formula for the latetide \(\phi\))

\[\begin{equation} \phi= \chi+\sum_{j=1}^6\delta_j\sin(2j\chi) \end{equation}\]

● 経度\(\lambda\)の計算式 (Formula for the longitude \(\lambda\))

\[\begin{equation} \lambda= \lambda_0+\tan^{-1}\left(\frac{\sinh \eta’}{\cos\xi’}\right) \end{equation}\]

◆検証(Validation)

御嶽山山頂の緯度経度 \((35^{\circ}53’34”N=35.89278^{\circ}N\), \(137^{\circ}28’49”E=137.48028^{\circ}E\)) をまず関数latlon2xyを用いて直交座標に変換する。 世界測地系8系(基準点:36\(^{\circ}\)N, 138.5\(^{\circ}\)E)を用いた場合、 結果は \(x\)(東)\(=-92058.3366\)(m)、\(y\)(北)\(=-11415.4240\)(m)となる。 次にこの座標を用いて関数xy2latlonにより緯度経度を計算すると 35.89278N, 137.48028Eとなって元の座標が復元される。
Let us first convert the latitude and longitude of the summit of Mt. Ontake, Japan (\(35^{\circ}53’34”N=35.89278^{\circ}N\), \(137^{\circ}28’49”E=137.48028^{\circ}E\)) to the cartesian coordinate using function latlon2xy. The reference point of zone 8 (36\(^{\circ}\)N, 138.5\(^{\circ}\)E) is used. The result is \(x\) (east) \(=-92058.3366\) (m) and \(y\) (north) \(=-11415.4240\) (m). Next, using these coordinates, compute the latitude and longitude by function xy2latlon. The result is 35.89278N and 137.48028E, consistent with the original coordinate.

◆引用文献(References)