N88-BASICで天体の軌道 (2回目)
2021/10/19(火)
N88-BASICで天体の軌道 (2回目)
軌道要素、ユリウス日、回転行列を使用して
実際の天体の軌道を描画します
使用している式などの詳しい説明は
https://ulprojectmail.blogspot.com/2021/06/n88-basicmatrix-1.html
N88-BASICで行列(matrix) (1回目)～
https://ulprojectmail.blogspot.com/2021/09/n88-basicjd-1.html
N88-BASICでユリウス日(JD) (1回目)～
https://ulprojectmail.blogspot.com/2021/10/kepler-8.html
天体の軌道(Kepler) (8回目)
を参照して下さい
平均運動nは
n = 2π/P
で求めています
周期Pが分からない天体は
地球の軌道長半径a=1AU,周期P=1年より
M = nt(t=1年で1周なのでM = 2π)
2π = n = √(μ/a3) = √(μ)
√(μ) = 2πとして下記式でnを求めています
n = √(μ/a3) (if e < 1)(楕)円
n = √(μ) (if e = 1)放物線
n = √(μ/|a|3) (if e > 1)双曲線
このnを1年の日数で割って1日毎の
平均運動として以下の様に使用しています
平均近点角M = nt (日数t)
行列Mを*IDENTITYで単位行列にし、
dに回転角、x, y, zに回転軸のベクトルの
成分を入れて*ROTATEを呼ぶと、
行列Mに回転行列が入ります。
*IDENTITYを呼ぶまで、
つづけて移動や回転する事ができますが、
初めに呼出した方が、あとから計算されます。
例、
GOSUB *IDENTITY ’Mを単位行列にする
d = PI: x = 0: y = 0: z = 1
GOSUB *ROTATE ’次にz軸に対して180°回転
x = 100: y = 0: z = 0
GOSUB *TRANSLATE ’初めx軸方向に100移動
なら、
p’= Mp、行列と点p(位置ベクトル)の積で
p’は点pを、
x軸方向に100移動させて、
z軸に対して180°回転させた点になります
Rを回転行列、Tを移動行列とすると、
M = RTより、Mp = RTpなので
Tpの計算結果にRが掛けられるので、
この様な順で作用します
VIEWMODE=1とすると、TRpとなり、
Rが先に作用します。
VIEWMODEは*ROTETEや*TRANSLATEを
呼ぶ前に設定しておきます
物体を動かすときはRTp (VIEWMODE=0)
見ている人(視点や視線)が動く場合は
pRT (VIEWMODE=1)を使用すると考えると
良いと思います
理由は、VIEWMODE=0では座標軸が物体に固定
され(物体とともに動き)、VIEWMODE=1では
座標軸は見ている人に固定されます
昇交点黄経Ωで軌道傾斜軸ベクトル
(x, y, z)を求めます。
(x, y, z) = (cosΩ, sinΩ, 0)
これと、軌道傾斜角iで*ROTATEを
呼出し、回転行列をMに追加しています
また、さらに回転させて黄道(地球の軌道)面
が画面と一致するようにしています
表示はx, yのみで行っているため遠くも
同じ大きさになります
(z軸によって、遠くならx, yを小さくすれば
遠近感を表現できますが、今回あまり効果が
ないと思い遠近法は無視しました。)
軌道傾斜の結果は、傾斜角が小さすぎて
見た目はほとんど効果がありませんでしので
横(斜め)から見た軌道の表示を追加しました
ついでに、小惑星2つ(竜宮と糸川)を追加しました
軌道周期Pが載っていなかったので、P=0とすると、
ケプラーの第3法則T2 / a3 = const.(一定)
(T軌道周期、a軌道長半径)より、
地球のデータ(周期P3, 長半径a3)より
k = P32 / a33とし、P = √(k･a3)
で計算するようにしました
P=0指定した小惑星の軌道を消す色は暗い紫
近日点方向は青にしています
小惑星と惑星の元期(軌道要素を決めた日)が
違うので、DATA文の構造を変更し、
複数の元期と、表示開始日、日数、経過日数
を指定できるようにしました
地球の自転軸の傾き方向が一定ではないため
秋分の日などは毎年時刻が変化します
このプログラムでは平均値のみの表示です
黄道を真上から見るときは、0(ﾟ)
黄道を手前真横から見るときは、90(ﾟ)です
次に、曜日、惑星(木星、土星、天王星、海王星)
と準惑星(冥王星)を追加します
ユリウス日は、
紀元前4713(-4712)年1月1日(月)12:00(世界時) = 0日
とした日数なので、0.5を足せば、
紀元前4713(-4712)年1月1日(月)0:00(世界時) = 0日
とした日数になります。
ユリウス日+0.5をdとし、
d / 7の余りをwとすると、
w = d - 7[d / 7] … []はガウス記号
[]はBASICではINT()、C言語ではfloor()です。
[-1.2] = -2、[1.2] = 1
N88-BASICやNL-BASICの場合MODはdが整数型の
範囲を超えるので使用できません
w=0～6が曜日(月～日)に対応します
全軌道が入るよう1AU=4dotにしたので
火星以内の惑星の軌道はほぼ分からなく
なりました
その他の小惑星や準惑星も追加できますが、
惑星数(最後の惑星番号)を変更し忘れると
追加分が表示されません
軌道データは、理科年表2020(国立天文台)
を使用しました
2021/9/23(木)の秋分の日から
2021/10/31(日)までの
不具合があるかもしれませんがご了承下さい
NL-BASICとblg～.zip(plan002.bas)は
このブログ(以下のリンク)から
ダウンロードできます
軌道要素、ユリウス日、回転行列を使用して
実際の天体の軌道を描画します
使用している式などの詳しい説明は
https://ulprojectmail.blogspot.com/2021/06/n88-basicmatrix-1.html
N88-BASICで行列(matrix) (1回目)～
https://ulprojectmail.blogspot.com/2021/09/n88-basicjd-1.html
N88-BASICでユリウス日(JD) (1回目)～
https://ulprojectmail.blogspot.com/2021/10/kepler-8.html
天体の軌道(Kepler) (8回目)
を参照して下さい
平均運動nは
n = 2π/P
で求めています
周期Pが分からない天体は
地球の軌道長半径a=1AU,周期P=1年より
M = nt(t=1年で1周なのでM = 2π)
2π = n = √(μ/a3) = √(μ)
√(μ) = 2πとして下記式でnを求めています
n = √(μ/a3) (if e < 1)(楕)円
n = √(μ) (if e = 1)放物線
n = √(μ/|a|3) (if e > 1)双曲線
このnを1年の日数で割って1日毎の
平均運動として以下の様に使用しています
平均近点角M = nt (日数t)
行列Mを*IDENTITYで単位行列にし、
dに回転角、x, y, zに回転軸のベクトルの
成分を入れて*ROTATEを呼ぶと、
行列Mに回転行列が入ります。
*IDENTITYを呼ぶまで、
つづけて移動や回転する事ができますが、
初めに呼出した方が、あとから計算されます。
例、
GOSUB *IDENTITY ’Mを単位行列にする
d = PI: x = 0: y = 0: z = 1
GOSUB *ROTATE ’次にz軸に対して180°回転
x = 100: y = 0: z = 0
GOSUB *TRANSLATE ’初めx軸方向に100移動
なら、
p’= Mp、行列と点p(位置ベクトル)の積で
p’は点pを、
x軸方向に100移動させて、
z軸に対して180°回転させた点になります
Rを回転行列、Tを移動行列とすると、
M = RTより、Mp = RTpなので
Tpの計算結果にRが掛けられるので、
この様な順で作用します
VIEWMODE=1とすると、TRpとなり、
Rが先に作用します。
VIEWMODEは*ROTETEや*TRANSLATEを
呼ぶ前に設定しておきます
物体を動かすときはRTp (VIEWMODE=0)
見ている人(視点や視線)が動く場合は
pRT (VIEWMODE=1)を使用すると考えると
良いと思います
理由は、VIEWMODE=0では座標軸が物体に固定
され(物体とともに動き)、VIEWMODE=1では
座標軸は見ている人に固定されます
昇交点黄経Ωで軌道傾斜軸ベクトル
(x, y, z)を求めます。
(x, y, z) = (cosΩ, sinΩ, 0)
これと、軌道傾斜角iで*ROTATEを
呼出し、回転行列をMに追加しています
また、さらに回転させて黄道(地球の軌道)面
が画面と一致するようにしています
表示はx, yのみで行っているため遠くも
同じ大きさになります
(z軸によって、遠くならx, yを小さくすれば
遠近感を表現できますが、今回あまり効果が
ないと思い遠近法は無視しました。)
軌道傾斜の結果は、傾斜角が小さすぎて
見た目はほとんど効果がありませんでしので
横(斜め)から見た軌道の表示を追加しました
ついでに、小惑星2つ(竜宮と糸川)を追加しました
軌道周期Pが載っていなかったので、P=0とすると、
ケプラーの第3法則T2 / a3 = const.(一定)
(T軌道周期、a軌道長半径)より、
地球のデータ(周期P3, 長半径a3)より
k = P32 / a33とし、P = √(k･a3)
で計算するようにしました
P=0指定した小惑星の軌道を消す色は暗い紫
近日点方向は青にしています
小惑星と惑星の元期(軌道要素を決めた日)が
違うので、DATA文の構造を変更し、
複数の元期と、表示開始日、日数、経過日数
を指定できるようにしました
地球の自転軸の傾き方向が一定ではないため
秋分の日などは毎年時刻が変化します
このプログラムでは平均値のみの表示です
黄道を真上から見るときは、0(ﾟ)
黄道を手前真横から見るときは、90(ﾟ)です
次に、曜日、惑星(木星、土星、天王星、海王星)
と準惑星(冥王星)を追加します
ユリウス日は、
紀元前4713(-4712)年1月1日(月)12:00(世界時) = 0日
とした日数なので、0.5を足せば、
紀元前4713(-4712)年1月1日(月)0:00(世界時) = 0日
とした日数になります。
ユリウス日+0.5をdとし、
d / 7の余りをwとすると、
w = d - 7[d / 7] … []はガウス記号
[]はBASICではINT()、C言語ではfloor()です。
[-1.2] = -2、[1.2] = 1
N88-BASICやNL-BASICの場合MODはdが整数型の
範囲を超えるので使用できません
w=0～6が曜日(月～日)に対応します
全軌道が入るよう1AU=4dotにしたので
火星以内の惑星の軌道はほぼ分からなく
なりました
その他の小惑星や準惑星も追加できますが、
惑星数(最後の惑星番号)を変更し忘れると
追加分が表示されません
軌道データは、理科年表2020(国立天文台)
を使用しました
2021/9/23(木)の秋分の日から
2021/10/31(日)までの
軌道を表示しています
不具合があるかもしれませんがご了承下さい
このブログ(以下のリンク)から
ダウンロードできます
