ここは何処？　わたしは誰？？？

　「わたしは誰？？？」は、ご自分で考えていただくとして、「ここは何処？」について考えてみましょう。

　一般的に、場所は○○県△△市□□町××番地と表現され、郵便や宅配便も届きます。　しかし、幾何学的な計算を行う衛星の位置計算では「ここは何処」を数値で表現しなければ計算ができ ません。

経度、緯度、高さ

　地面の上で生活しているとき便利な場所を表す方法としては、北緯や東経などの経緯度があります。　地上ではあまり意識しませんが、海上の船の位置などは近くに島などの目標物がない場合 でも経緯度で表現すれば１点に定まります。　また、海面や平地では気にしませんが、山の上や飛行機の位置を示すためには高さの情報も必要になります。

　一方、衛星の位置を表す軌道要素で使われている表現は、地心赤道直交座標（空間に固定された座標）と呼ばれるもので、地球の重心を原点に赤道平面で春分点方向をＸｓ、Ｘｓから東に９０度方向をＹｓ、北極方向をＺｓで表します。

　したがって、自分がいる経度λ、緯度φ、高さｈ と Ｘｓ、Ｙｓ、Ｚｓ を関係づける必要があります。

ベッセル回転楕円体

　日本で作られている地図のほとんどはベッセル回転楕円体を採用してい ます。　この楕円体にはつぎの数値が与えられています。

　　赤道半径　　　　　　ａ＝６３７７．３９７１５５００ｋｍ

　　離心率の２乗　　ｅ＾２＝０．００６６７４３７２２３０６１４

　この形は地心から赤道までの距離に比べて、地心から北極（南極）までの距離が短くなります。　

　それではベッセル回転楕円体で、自分がいる経緯度、北極、南極、地心を通る子午面を考えると、その面は楕円となり、中心が地球の重心になります。　楕円の性質として楕円の接線に垂直な線（法線）は赤道上と両極以外では中心を通らないと言うことです。　

　一方、緯度は法線と赤道面のなす角であり、中心を通る角度と違いが生じてしまいます。　また、自分のいる場所と中心までの距離に比べて、法線に沿って北極と南極を結ぶ線と交わる点までの距離（東西線曲率半径： Ｎ と呼ぶ）に違いが生じことになります。　東西線曲率半径 Ｎ は 赤道半径 ａ 、離心率 ｅ 、緯度 φ からつぎの式で求めることができます。

　　Ｎ＝ａ／（１−ｅ＾２ ＊ sin＾２φ）＾０．５

　したがって、経度、緯度、高さから、経度０度方向をＸｊ、東経９０方向をＹｊ、北極方向をＺｊとする直交座標系に変換するときには、以上のことを考慮しなければなりません。　これらの関係はつぎのように表すことができます。

　　Ｘｊ＝（Ｎ＋ｈ） cosφ cosλ

　　Ｙｊ＝（Ｎ＋ｈ） cosφ sinλ

　　Ｚｊ＝（Ｎ（１−ｅ＾２）＋ｈ） sinφ

となる。　このＸｊ、Ｙｊ、Ｚｊは日本 （東京）測地系と呼ばれ ていますが、この中心は地球の重心を通っていません。　地球の重心に一致させるには、各座標をつぎの値だけ平行移動する必要があります。　ここで、平行移動した後の座標をＸ ｍｇ、Ｙｍｇ、Ｚｍｇとするとつぎのような関係があ ります。

　　Ｘｍｇ　＝　Ｘｊ　−　１３６　ｍ

　　Ｙｍｇ　＝　Ｙｊ　＋　５２１　ｍ

　　Ｚｍｇ　＝　Ｚｊ　＋　６８１　ｍ

以上の手順で座標変換を行うと、経緯度と高さを地球重心を中心としたＧ系地心直交座標（地球に固定した座標）に変換することができます。

　それでは、地球の重心から見た「ここは何処」を計算して、フォームに表示するプログラムを使ってみましょう。

Calsat32Keiido2Xyz.zipをダウンロードして、解凍してください。

　解凍したフォルダに含まれる"Calsat32Keiido2Xyz.exe"を起動して、「 経緯度と高さを地心直交座標に変換ボタン」を押すとつぎのようなフォームが表示されます。

　また、"Calsat32Keiido2Xyz.vbp"をクリックすると Visual BASIC が起動して、つぎのコードを見ることができます。

Private Sub Command1_Click()

Dim dblPI As Double
Dim dblA As Double, dblE2 As Double, dblN As Double
Dim dblKeido As Double, dblIdo As Double, dblTakasa As Double
Dim dblXj As Double, dblYj As Double, dblZj As Double
Dim dblXmg As Double, dblYmg As Double, dblZmg As Double


'定数の設定
    dblPI = 3.14159265358979
    dblA = 6377397.155
    dblE2 = 0.006674372230614

'観測点の経度、緯度、高さをＪ系直交座標に変換
    dblKeido = Val(Text1.Text) * 2 * dblPI / 360#
    dblIdo = Val(Text2.Text) * 2 * dblPI / 360#
    dblTakasa = Val(Text3.Text)
    dblN = dblA / ((1# - dblE2 * Sin(dblIdo) ^ 2#) ^ 0.5)
    dblXj = (dblN + dblTakasa) * Cos(dblIdo) * Cos(dblKeido)
    dblYj = (dblN + dblTakasa) * Cos(dblIdo) * Sin(dblKeido)
    dblZj = (dblN * (1# - dblE2) + dblTakasa) * Sin(dblIdo)

'観測点のＪ系直交座標をＧ系地心直交座標に変換
    dblXmg = dblXj - 136#
    dblYmg = dblYj + 521#
    dblZmg = dblZj + 681#

'結果の表示
    Cls
    Print "地心直交座標では、"
    Print "X:"; dblXmg / 1000#; "(km)"
    Print "Y:"; dblYmg / 1000#; "(km)"
    Print "Z:"; dblZmg / 1000#; "(km)"

End Sub

　

　テキストボックスには東経１３５度、北緯３５度、高さ５０ｍがセットされていますが、場所を変更して計算してみましょう。

　なお、２００２年４月からは測量法が改正され、楕円体としてはＧＲＳ−８０、測地系としてはＩＴＲＦ−９４系が採用されることになります。　このため、ＣＡＬＳＡＴ３２　Ｖｅｒ１．１．４以降では新しい測地系を採用することとしてプログラムの改訂を行いますので、新しい経緯度で自局の位置を設定する必要があります。　しかし、実用的には新旧の差は他の誤差要因に比べて小さなものであると考えられ、問題になることはありません。

つぎのＵＲＬに、新経緯度に変換するためのプログラムと概要があります。

    http://vldb.gsi.go.jp/sokuchi/tky2jgd/about.html   

　これで、自分の 経緯度からＧ系地心直交座標への変換手順については理解いただけたのではないかと思います。　残っているのは、Ｇ系地心直交座標（地球に固定した座標）と地心赤道直交座標（ 空間に固定した座標）を結びつける春分点の方向、すなわちグリニッジ平均恒星時を求めることです。　グリニッジ平均恒星時と時刻 は深く結びついているので改めて別の項でお話ししたいと思います。