1. 前説

あの日打ち上がった衛星の位置を僕達はまだ知らない。

衛星データを扱う場合、ときとして観測点と衛星との位置関係（ジオメトリ）を知りたくなることがあります。衛星写真をなんとなく眺めるだけならよいのですが、画像の幾何的な歪みを補正したり、データから何か意味のある地表面物理量を推定したいとなると、観測点からのシグナル（≒電磁波）がどのような経路で衛星に届くかを考慮する必要があるためです。衛星データには多くの場合、観測時刻、観測点からセンサへの仰角・方位角、そしてもちろん観測点の地理座標などの情報がメタデータとして付帯しているので、それを参照すれば事足りることが多いです。が、プロダクトのレベルによっては（とくに提供機関がもろもろの補正を済ませ、物理量の推定まで行っている高次プロダクトの場合）、そういう細かい情報が見当たらないこともあります。あるいは衛星の時々刻々の軌道を自分で計算したい日もあると思います。人生何が起こるかわかりませんので、自分で軌道計算ができると安心です。

衛星は地球からの引力により、地球を焦点のひとつとする楕円軌道を描きます（ケプラーの第一法則）。楕円というのは円の次くらいにはシンプルです。なので観測時刻と、楕円軌道を決定する然るべきパラメータがわかれば、筆者のようなnot航空宇宙工学徒であっても衛星の位置を推定することはできます。

あらかじめ断っておくと、本当は、地表付近を周回する衛星の軌道は

地球が真球でないことによる場所ごとの引力の不均一性
大気の微妙な抵抗
他の惑星や月からの引力
太陽からの放射圧

などのせいで楕円軌道よりもっと複雑になります。高精度の位置推定が必要な場合はプロに任せることにしましょう。アプリとかも色々出ていますし。

楕円上での位置決定

右図において、赤色の楕円が衛星の軌道、バツ印が地球（焦点）だとしましょう。

楕円中心を原点とし、焦点に向かう向き（長軸方向）をU軸、それと直交する向き（短軸方向）をV軸とするU-V平面を考えます（ついでに右手系に従いUV平面に垂直なW軸も取っておきましょう）。楕円の形自体は、長径aと離心率εがあれば確定し、短軸bおよび面積Sは

となります。

衛星がある時刻に点Pにいたとしましょう。このとき、楕円を直接扱うのではなく、中心を共有し半径aをもつ補助円を考えます。点Pを補助円上に投影した点をP'とし、U軸とOP'のなす角をEとすれば、P'の座標は

とかけます。Eは離心近点角などと呼ばれます。Pの座標はこれをV軸方向にb/a倍したものなので

です。ということで離心率ε、長径a、離心近点角Eがわかれば衛星の座標がわかりますね。後の計算のためには、地球中心を原点に取り直したほうが便利です:

この座標系を改めてU-V (-W) 座標と呼ぶことにします。

座標の変換

得られた座標Pは、楕円中心を原点とする衛星の軌道面（UV平面）上での値でした。観測点との関係を調べるためには、これを地球中心を原点とする3次元座標に変換する必要があります。一般に地軸北向き方向をZ軸、春分点の方角をX軸にとるようです（True Equator, Mean Equinox (TEME) という赤道座標）。

両座標系の幾何学的関係はこのページなどを参照してください。楕円が地球の赤道面（XY平面）をつらぬくポイントは2つありますが、南→北へ通過する点を昇交点、北→南へ通過する点を降交点と言います。近地点と昇交点のなす角をω、UV平面とXY平面のなす角をi、X軸と昇交点のなす角をΩと書くことにすると、XYZ座標系をUVW座標系に一致させるためには、

Z軸まわりにXY平面をΩだけ回転させる
新たなX軸を中心にYZ平面をiだけ引き起こす（この時点で旧Z軸がW軸に一致）
新たなZ軸（W軸）を中心にXY平面をωだけ回転させる

の順で操作すればよいです。それを変換行列で表しますと、

となります。ωを近地点引数、iを軌道傾斜角、Ωを昇交点赤経といいます。勢いで書いてしまいましたがW成分はないのでw=0ですね。

さて、TEMEは天文学などでは便利なようですが、我々がほしいのは緯度経度座標です。X軸をグリニッジ子午線の方向に一致させましょう。地球中心・地球固定直交座標系（Earth Centered, Earth Fixed: ECEF）と言うようです。春分点の方向から測ったグリニッジ子午線方向の角度（赤経）をθGとし、

によって座標系を回転させます。

ECEF座標系において、X'軸を始線とするZ'軸まわりの回転角が経度λとなります。そして赤道面からの仰角が緯度Φ…と言いたいところなのですがこれは地心緯度と呼ばれる量で、地球が真球でないせいで、地表でみた天頂方向を基準に考えた地理緯度（普通に「緯度」というときにはこちらを指す）とは微妙に異なります。詳細は省きますがx', y', z'座標と地理緯度Φ, 経度λとの関係は

となります。Nは卯酉線曲率半径といい、

で定まる量です。Hは衛星の楕円体高、bE, aE, εEはそれぞれ回転楕円体（地球）の短径、長径、離心率です。この関係式からΦとλを逆にとけば、衛星の緯度経度座標が求まります。

観測点からみた衛星の方位角・天頂角も求めておきましょう。観測点の緯度経度を既知とすれば（ふつう既知だと思います）、先ほどの変換式を使って観測点のECEF座標が求まります。衛星座標（位置ベクトル）から観測点座標（位置ベクトル）を引けば、観測点から衛星への視線ベクトルになります。これを観測者を基準にした地平直交座標 (East-North-Up: ENU座標) に変換すればよいです。ENU座標は、

E軸回りに -(π/2-Φs) 回転させる（これでEN平面が赤道面と平行になる）
新たなU軸まわりに -(π/2+λs) 回転させる

ことによりECEF座標を平行移動させた状態になります。添字sは観測点を表します。右手系で、反時計回りを回転の正方向とすることに注意しましょう。平行移動（i.e., 座標原点の違い）は視差ベクトルを考える限り関係ありませんので、これでENU座標でみた視差ベクトルが得られたことになります。変換行列は

わかると思いますがここでいうuは最初のuvw座標のuとは違うので念のため。方位角はふつう北向きから東回りに測られることに注意して、方位角θaz・天頂角θznはそれぞれ

大きな流れとしては以上です。まあまあ長かったですね…では次章で各変数を具体的に計算していってみるとしましょう。

次章：楕円上での位置決定

Google Sites

Report abuse