先日「太陽・月・惑星視位置の海洋情報部近似式の係数の求め方」 という記事を書いたのですが、より詳細な説明・資料(Excelファイル)を孔雪さんからいただきましたので紹介します。
海王情報部から提供される近似式には水星視位置を求めるものはないのですが、水星にたいして同様の近似式を作るときの係数をどうやって求めたら(決めたら)いいか、という課題に対する解決策です。
じつは孔雪さんからは5月5日にいただいたのですが、記事を書くのが今日になってしまいました。説明・資料をこころよくご提供いただいた孔雪さんに御礼申し上げるとともに遅くなったことをお詫びします。
以下孔雪さんからいただいたメールをそのまま引用します。改行位置などをのぞいて原文のままです。
また提供していただいたファイルは
「水路部データDCT 2017略_A.xlsx」 (9,545kB)
からダウンロードできます。
注
いただいたファイルは「水路部データDCT 2017略.xlsx」という名称でした。
共有準備・問題のチェックをおこなった後保存し直したものが上のファイルです。
なお、孔雪さんから提供していただいた説明・資料(Excel)ファイルの著作権については私のブログの内容と同じに取り扱って構わない、つまり
「文章、画像、グラフ、Excelファイル、プログラムソース等すべて自由に引用あるいは再利用していただいてかまいません。 ただしお互いのオリジナリティーはたいせつにしていただきますようお願いします。 」
という扱いで差し支えないそうです。この点についても孔雪さんに感謝いたします。
------------------
以下孔雪さんからいただいたメールによるものです。
海洋情報部の係数にかなり近い係数を求める方法を探し出したわけですが、まだ、求めた係数が金星では海洋情報部のものとけっこう違うので、何とかしたいと思いました。
実は、海洋情報部のデータはN=17やN=29までとなっているので、どうもFFTでは無さそうだとは気づいていたのですが、分析ツールが最も手っ取り早いのでそうしていました。
どういうことかというと、FFTは元データも、得られる結果も、2のべき乗個でないといけないという制約が有るためで、N=32まで求めて必要な部分のみ掲載しているとは思えなかったからです。
これは、FFTではなく、高速コサイン変換を開発していたとしても同じことのようです。2のべき乗であるということは、フーリエ変換やコサイン変換を高速化する時の原理にかわってくることのようですから。
また、暦象年表からデータを得るとしたら、32個も元データを準備するより18個に減ったら大助かりということもありますので、高速化しない地道なコサイン変換をやってみようと思い、トライしてみました。
検証は、金星の2017/1~4月のデータで行なっています。
途中経過は、
・とりあえず、地道なコサイン変換バージョンを作ってみる
⇒エクセルFFTと同じ結果が得られたのでひとまず安心
・データ17個で実施
・データ18個で実施
しかし、得られた結果は、芳しくありません。
さて、ここまではデータを折り返すときに、データが有るそのポイントで折り返していました。
分析ツールは、もともとフーリエ変換をするものでコサイン変換用では無いのですが、出た結果から虚数項を無くそうとすると、そうしないとうまくいかないと思ったためです。
地道なコサイン変換ですし、ここでちょっと考えを変えて、データとデータの中間で折り返すようにデータのタイミングを1/2ずらすことにしました。
1/2ずらしたやりかたで、
・データ17個で実施 → N=17の項まで得られるがN=17の項がほぼゼロであまり意味がない
⇒
そして、
・データ18個で実施
⇒ ついにやりました!
得られたデータは金星ではほぼ誤差無し!(1ヶ所最後の桁が1違うだけ。)
あとは、検証有るのみです。
意を決して、2017年の惑星の全データを使って検証することにしました。
月は、データが多すぎて、そこまでの気合は無いし..
また、目的は水・天・海・冥なので、もし月のやり方が違っていたとしても..まぁいいか。
ということで、地道にデータを揃えて検証.....
「きっと、この方法で求めている違いにない」と言える結果です。
また、さすがに、月も1つだけは検証してみました。
これも、H.P.以外は満足いく結果でした。
H.P.ですが、暦象年表では直接地心距離を求めているからかな?と思います。
しかし、地球赤道半径=6378.137kmだと思うんだがなぁ。
ところで、さて、
暦象年表のデータを使って、海洋情報部の係数とほぼ同様の値が計算できるということは、
海洋情報部が持っている真のデータは、暦象年表のデータと同じか、極めて近い
ということだと思います。
ということは、暦象年表のデータを使って海洋情報部の式の精度を求めることが出来る
ということでしょう。
暦象年表のデータと天体の真の位置との精度はさておきですが、海洋情報部の式は近似式なので、近似式としての精度という意味でです。
水星の係数を求めるにあたり、本家海洋情報部の係数による精度を求めておきたいものです。
実施有るのみということで、これまた地道に検証しました!
フーリエ変換やコサイン変換は、変換に使ったポイントで逆変換をしたら誤差の無い値が得られてしまいます。
そこで、所定範囲内で400個以内のデータ、具体的には4ヶ月範囲なら8時間毎、3ヶ月範囲なら6時間毎の海洋情報部のデータと係数を使って計算した値を比較してみました。
幸いなことに、今回係数を求めるために使っている元データは等間隔ではないので、このような等間隔なデータで検証しても、ランダムにとったデータと同様の検証が出来ていると見ています。
結果は、どうも海洋情報部は10^-5以内程度を狙っているみたいなのですが..
そうもなっていない部分もあります。
歯切れが悪いのは、金星の1~4月のDecのデータはかなり精度が悪く、度で5.2x10^-4、秒にして1.9秒にもなる部分が有ったためです。
これを参考にして、水・天・海・冥の係数を求めることにします。
実は、天・海は、ここまでの過程で、4ヶ月使えるN=17までの係数で求めてしまいました。
というのは、検証用に海洋情報部からデータを取ってくる時に、水・金・火・木・土・天・海は同一ページから取れるので、連続でとってしまったほうが高率が良いためです。
冥は、惑星じゃないけど..過去から水路部の計算式もあることですし、歴史的背景ということで、求めることにします。
冥王星なら、たぶん1年範囲でも大丈夫だろうということで、範囲は1年としました。
水星も、ここまでの過程で、4ヶ月使えるN=17までの係数で求めていましたが、予定通り全く精度が足りません。
分析ツールでやった時にも、3ヶ月範囲でN=30以内くらいと求めましたが、N=17までの係数でならエクセルシートを直さなくても良いので、3ヶ月範囲、2ヶ月範囲と狭めてみました。
3ヶ月範囲は、予定通り全くダメでしたが、2ヶ月範囲ならOKみたいなので、1年分6シート作って出来上がり、思ったのに....
なんと、最後の11,12月用で大きな誤差が...がっかりです。
最後の期待を込めて、3ヶ月範囲でN=29までの係数で、4シート作って..
何とか、やっと2017年分全部揃いました。
..ということにしておきましょう。
3ヶ月範囲N=29までの係数で、10~12月のDecは度で1.7x10^-4、秒で0.6秒の誤差があります。
これは、目標の10^-5に届いていませんが、海洋情報部オリジナルの金星の精度よりは良いので、OKということにしておきます。
ところで、今回は3ヶ月範囲N=29までのほうが、2ヶ月範囲N=17までより、良い結果が出ましたが、これは偶然ではないかと思っています。
単純計算で、
3ヶ月範囲で30個のデータだと、平均して3.0日毎のデータ
2ヶ月範囲で18個のデータだと、平均して3.4日毎のデータ
ということで、3ヶ月範囲N=29までのほうがデータ密度が高いので当然のように
見えますが、いま求めている方法では、データの時間間隔自体が等間隔では無いので、一概には比較出来ないと思います。
未検証ですが、天体位置が急激に変化する時期が、たまたま、時間密度が高い時に来た場合に良い精度となり、時間密度が低い時に来た場合に悪い精度となるのだと思っています。
2017.5.5 孔雪
-----------------------
追記
ブログに掲載していただけるということでしたら、
文の折返し位置は、ブログで読みやすいように変更していただいて結構です。
・今回のシートは、求める係数は、76行目に出てきます。(全シート共通)
・文中では海洋情報部と書いていますが、エクセルシートは水路部で通しています。
行幅等、あらゆる場所で3文字の水路部のほうが好きですから。
だいたい「海洋情報部」って長過ぎなんだよな。文字打つのも打ちにくいし(笑)
添付したデータですが、本当は全シートに「29水星3ヶ月10月」と同様の約400個の検証データが付いています。
しかし、それだと50MB超という巨大ファイルになってしまっているので、まとめシートに結果の値をコピーして、検証データは消してあります。
コメント