mk-mode BLOG

このブログは自作の自宅サーバに構築した Debian GNU/Linux で運用しています。
PC・サーバ構築等の話題を中心に公開しております。(クローンサイト: GitHub Pages

ブログ開設日2009-01-05
サーバ連続稼働時間
Reading...
Page View 合計
Reading...
今日
Reading...
昨日
Reading...

Ruby - JPL 天文暦データから地球と惑星の距離を計算!

[ プログラミング, 暦・カレンダー ] [ Ruby, カレンダー ]

こんばんは。

NASA の機関である JPL(Jet Propulsion Laboratory) が惑星探査用に編纂・発行している太陽・月・惑星の暦の最新版 DE430 には太陽・月・惑星の位置(ICRS座標系)の情報が格納されています。

それらの値を使用して、地球から太陽・月・その他の惑星との間の距離を Ruby で計算してみました。(JPL DE430 の情報を読み込むには自作の RubyGems ライブラリ “"eph_jpl” を使用する)

0. 前提条件

  • Ruby 2.3.1-p112 での作業を想定。
  • JPL DE430 のバイナリデータをこちらから取得し “JPLEPH” とリネームして適当なディレ取りに配置しておく。
  • 当方作成の RubyGems ライブラリ eph_jpl - RubyGems.org がインストール済みであること。

1. Ruby スクリプトの作成

jpl_distance_earth.rb
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
#! /usr/local/bin/ruby
# coding: utf-8
#---------------------------------------------------------------------------------
#= JPLEPH(JPL の DE430 バイナリデータ)読込後に座標(位置・速度)を計算し、さらに
#  地球から太陽・惑星・月までの距離を計算する。
#  : 自作の RubyGems ライブラリ eph_jpl を使用する
#    (JPL DE430 のバイナリデータが準備済みであること)
#
# date          name            version
# 2016.07.13    mk-mode.com     1.00 新規作成
#
# Copyright(C) 2016 mk-mode.com All Rights Reserved.
#---------------------------------------------------------------------------------
# 【引数】 ユリウス日
#          (省略可。省略時は現在日時のユリウス日を地球時TT(≒太陽系力学時TDB)とみなす)
#---------------------------------------------------------------------------------
#++
require 'eph_jpl'

class JplDistanceEarth
  FILE_BIN = "/path/to/JPLEPH"
  USAGE    = "【使用方法】 ./jpl_distance_earth.rb [ユリウス日]"
  CENTER   = 3
  KM       = false  # true: km, false: AU
  BODIES = {
     1 => "  水星 (Mercury)",  2 => "  金星 (Venus)  ",  3 => "  地球 (Earth)  ",
     4 => "  火星 (Mars)   ",  5 => "  木星 (Jupiter)",  6 => "  土星 (Saturn) ",
     7 => "天王星 (Uranus) ",  8 => "海王星 (Neptune)",  9 => "冥王星 (Pluto)  ",
    10 => "    月 (Moon)   ", 11 => "  太陽 (Sun)    "
  }

  def initialize
    @jd   = get_args
    @unit = KM ? "KM" : "AU"
  end

  def exec
    puts "JD(TT) = #{@jd}\n#{BODIES[3]}"
    BODIES.each do |target, target_name|
      next if target == 3
      eph = EphJpl.new(FILE_BIN, target, CENTER, @jd, KM)
      ps = eph.calc[0,3]
      d = calc_distance(ps)
      d = KM ? sprintf("%19.8f", d) : sprintf("%11.8f", d)
      puts "- #{target_name} : #{d} #{@unit}"
    end
  rescue => e
    msg = "[#{e.class}] #{e.message}\n"
    msg << e.backtrace.each { |tr| "\t#{tr}"}.join("\n")
    $stderr.puts msg
    exit 1
  end

  private

  #=========================================================================
  # コマンドライン引数取得
  #
  # @param:  <none>
  # @return: jd
  #=========================================================================
  def get_args
    jd = ARGV.shift
    return jd ? jd.to_f : gc2jd({
      year:  Time.now.strftime("%Y").to_i,
      month: Time.now.strftime("%m").to_i,
      day:   Time.now.strftime("%d").to_i,
      hour:  Time.now.strftime("%H").to_i,
      min:   Time.now.strftime("%M").to_i,
      sec:   Time.now.strftime("%S").to_i
    })
  rescue => e
    puts USAGE
    raise
  end

  #=========================================================================
  # 年月日(グレゴリオ暦)からユリウス日(JD)を計算する
  #
  # * フリーゲルの公式を使用する
  #   JD = int(365.25 * year)
  #      + int(year / 400)
  #      - int(year / 100)
  #      + int(30.59 (month - 2))
  #      + day
  #      + 1721088
  # * 上記の int(x) は厳密には、 x を超えない最大の整数
  # * 「ユリウス日」でなく「準ユリウス日」を求めるなら、
  #   `+ 1721088` を `- 678912` とする。
  #
  # @param:  hash = {
  #            year(int), month(int), day(int), hour(int), min(int), sec(int)
  #          }
  # @return: jd (ユリウス日)
  #=========================================================================
  def gc2jd(hash)
    year  = hash[:year]
    month = hash[:month]
    day   = hash[:day]
    hour  = hash[:hour]
    min   = hash[:min]
    sec   = hash[:sec]

    begin
      # 1月,2月は前年の13月,14月とする
      if month < 3
        year  -= 1
        month += 12
      end
      # 日付(整数)部分計算
      jd = (365.25 * year).floor \
         + (year / 400.0).floor \
         - (year / 100.0).floor \
         + (30.59 * (month - 2)).floor \
         + day \
         + 1721088
      # 時間(小数)部分計算
      t = (sec / 3600.0 + min / 60.0 + hour) / 24.0
      return jd + t
    rescue => e
      raise
    end
  end

  #=========================================================================
  # 距離計算
  #
  # @param:  ps       (x-座標, y-座標, z-座標)
  # @return: distance (距離)
  #=========================================================================
  def calc_distance(ps)
    return Math.sqrt(ps.inject(0) { |s, a| s + a * a })
  rescue => e
    raise
  end
end

exit 0 unless __FILE__ == $0
JplDistanceEarth.new.exec

2. Ruby スクリプトの実行

以下は、「地球時(TT)2016年7月13日0時0時0分0秒」で計算した例。(距離単位を AU で計算したものと km で計算したもの)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
$ ./jpl_distance_earth.rb 2457582.5
JD(TT) = 2457582.5
地球 (Earth)
-   水星 (Mercury) :  1.32217757 AU
-   金星 (Venus)   :  1.69764291 AU
-   火星 (Mars)    :  0.62571104 AU
-   木星 (Jupiter) :  5.91822023 AU
-   土星 (Saturn)  :  9.23857800 AU
- 天王星 (Uranus)  : 19.99078537 AU
- 海王星 (Neptune) : 29.30214573 AU
- 冥王星 (Pluto)   : 32.12201719 AU
-     月 (Moon)    :  0.00270216 AU
-   太陽 (Sun)     :  1.01654960 AU

$ ./jpl_distance_earth.rb 2457582.5
JD(TT) = 2457582.5
地球 (Earth)
-   水星 (Mercury) :  197794949.19472104 KM
-   金星 (Venus)   :  253963763.90310407 KM
-   火星 (Mars)    :   93605038.57629812 KM
-   木星 (Jupiter) :  885353144.47506535 KM
-   土星 (Saturn)  : 1382071597.57920527 KM
- 天王星 (Uranus)  : 2990578925.21471119 KM
- 海王星 (Neptune) : 4383538608.69842434 KM
- 冥王星 (Pluto)   : 4805385374.42336559 KM
-     月 (Moon)    :     404237.25163473 KM
-   太陽 (Sun)     :  152073656.28368542 KM

3. 計算結果の検証

国立天文台暦計算室・暦象年表のツールで確認してみる。

かなりの精度で一致することが確認できる。(但し、月は AU 単位では確認できない)

4. 参考サイト


JPL 天文暦データから太陽・月・惑星の位置情報(ICRS 座標)を取得する RubyGems ライブラリを自作していたので、距離の計算自体は中学数学レベルで済みました。

以上。

Comments