ここから以下の3回に分けて数値積分は使わない軌道計算を行います。計算例は先に結果を示したサイディング・スプリング(C/2013A1)彗星です。
@地心太陽座標の計算
A軌道要素の計算
B要素変化法による軌道要素の最適化
計算内容の詳細は以下の本を参照下さい。
・「天体の軌道計算」中野主一著、誠文堂新光社(1992)
・「マイコン天文学 I」中野主一、恒星社厚生閣(1983)
・「天体軌道論」長谷川一郎、恒星社厚生閣(1983)
@観測データのダウンロード(今回は以下のファイルに添付済み)
まずは軌道要素を計算するためにMPCのWEBからC2013A1の観測データをダウンロードしたものがC2013A1.txtです。
具体的には以下のWEBページにアクセスします。
http://www.minorplanetcenter.net/db_search
C2013A1の場合は一番上のBOXに"C/2013 A1"と入力し、[SHOW]のボタンを押すと、計算された軌道要素とそれのもととなった観測データが表示されます。ここで[Observations]の項目の2行下にある[Dowinload]を押すと、テキスト形式で表示されるので、ファイルにコピーして[shift jis]モードでセーブします。ブラウザメニューのツールにある[ファイルで保存]で保存してもいいですが必ず[shift jis]モードで保存ください。
ファイル名の拡張子は[xxxxx.car]としてください。
Aプログラムの説明
今回のプログラムでは以下のファイルを使用/出力します。
1.TEXT "*******.car"観測データの読み込み(@のファイル)
2.TEXT "*******.dat"観測データ/視差補正太陽座標データの書き込み
3.BINARY "obsdata.dat"観測データ/視差補正太陽座標データの書き込み(元データを含む)
4.TEXT "ObsCodes.txt"観測地(天文台)情報 (MPC WEBより)
5.BINARY "lnxp1600p2200.405" JPL/NASA天体暦(今回は期間が短いのでDE405を使用)
JPL/NASAのwebよりダウンロードください。(ftp://ssd.jpl.nasa.gov/pub/eph/planets/Linux/de405/)
・[******]には最初に定義した天体名が入ります。今回の場合[C2013A1]
・[obsdata.dat]はプログラム内部で使うデータやり取り用のバイナリファイルです。
計算手順としては
1)"*******.car"から「観測時刻(UT)」「赤経・赤緯」「天文台コード」を読み込む。
(現時点でフォーマットが違うデータははじいています。)
2)「観測時刻(UT)」にΔTを加えてJD(TT)形式に変換
3)「赤経・赤緯」をラジアンに変換
4)「天文台コード」をもとに"ObsCodes.txt"にある対応する固有情報を読み出し、
JD(TT)時の地心太陽座標を補正する。
5)計算結果を"obsdata.dat"と"*******.dat"に書き込む。
・あらかじめ観測時刻に対応した「地心太陽座標」を計算しておくと、後に「日心彗星座標」を
計算した時に「地心太陽座標」を足すことで、すぐに「地心彗星座標」が計算できる。
「天体の軌道計算」や「マイコン天文学 I」でとられている手法です。
B計算結果
上記手順で計算したときのディスプレイ表示結果がresults1.txtです。
軌道計算(1) 地心太陽座標の計算のC言語ソース:
(いずれUPする最終形を使用してもらった方が良いとは思いますが・・・・)
comet_v1.0.zip
注:windows上のCコンパイラ(VC等)でコンパイルください。
JPL/NASAのDE405のファイル(lnxp1600p2200.405)がほかに必要です。(上記URL参照)
"jpl_eph.c"の217行目の天体暦ファイルのディレクトリ/ファイル名は環境に合わせて修整ください。
|