Ruby - 「グレゴリオ暦 -> ユリウス日」変換の検証(vol.2)！

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

以前、2種の計算式で「グレゴリオ暦 -> ユリウス日」の変換を行い、結果が同じになるかを検証しました。

今回は、別の計算式による変換も追加し、計3種で検証してみました。

0. 前提条件

  • Ruby 2.5.1-p57 での作業を想定。
  • 使用する２種の計算式は「フリーゲルの公式」とその他の計算式。
  • ここでの「ユリウス日」は JD(Julian Day) であり、JDN(Julian Day Number), CJD(Chronological Julian Day), MJD(Modified Julian Date) ではない。

1. 検証用 Ruby スクリプトの作成

以下のように作成してみた。（うるう年の２月から３月への変わり目、年末から年始への変わり目、うるう年でない年の２月から３月への変わり目をチェックするようにしている）

verify_gc2jd_2.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
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164

#! /usr/local/bin/ruby
# coding: utf-8
#=グレゴリオ暦 -> ユリウス日 変換の検証(UTC ベース)
#
# : ２種類の計算式で グレゴリオ暦 -> ユリウス日 の変換を行い、
#   相違がないかを検証する。
#
# date          name            version
# 2016.06.14    mk-mode.com     1.00 新規作成
# 2018.06.26    mk-mode.com     1.01 別の計算式(Celes Trak)による検証も追加
#
# Copyright(C) 2016-2018 mk-mode.com All Rights Reserved.
#---------------------------------------------------------------------------------

class VerifyGc2Jd
  def exec
    55.upto(59) { |sec| generate(2016,  2, 28, 23, 59, sec) }
     0.upto( 5) { |sec| generate(2016,  2, 29,  0,  0, sec) }
    puts "..."
    55.upto(59) { |sec| generate(2016,  2, 29, 23, 59, sec) }
     0.upto( 5) { |sec| generate(2016,  3,  1,  0,  0, sec) }
    puts "..."
    55.upto(59) { |sec| generate(2016, 12, 31, 23, 59, sec) }
     0.upto( 5) { |sec| generate(2017,  1,  1,  0,  0, sec) }
    puts "..."
    55.upto(59) { |sec| generate(2017,  2, 28, 23, 59, sec) }
     0.upto( 5) { |sec| generate(2017,  3,  1,  0,  0, sec) }
  rescue => e
    msg = "[#{e.class}] #{e.message}\n"
    msg << e.backtrace.each { |tr| "\t#{tr}" }.join("\n")
    $stderr.puts msg
  end

  private

  def generate(*dt)
    str =  Time.new(*dt).strftime("%Y-%m-%d %H:%M:%SUTC")
    jd_1, jd_2, jd_3 = gc2jd1(*dt), gc2jd2(*dt), gc2jd3(*dt)
    str << " -> JD: %16.8f, %16.8f, %16.8f" % [jd_1, jd_2, jd_3]
    str << " (NOT MATCH!)" if jd_1 != jd_2 || jd_1 != jd_3
    puts str
  rescue => e
    raise
  end

  #=========================================================================
  # グレゴリオ暦 -> ユリウス日(JD)
  #
  # * フリーゲルの公式を使用する
  #     month = 1, 2 の場合、year = year - 1, month = month + 12 とする
  #     JD = int(365.25 * year)
  #        + int(year / 400)
  #        - int(year / 100)
  #        + int(30.59 * (month - 2))
  #        + day
  #        + 1721088.5
  #   ※上記の int(x) は厳密には、x を超えない最大の整数
  #   ※正子を整数とする JDN を求めるなら、 1721088.5 を 1721089 に置き換える
  #   ※修正ユリウス日 MJDを求めるなら + 1721088.5 を - 678912 に置き換える。
  #
  # @param:  year
  # @param:  month
  # @param:  day
  # @param:  hour
  # @param:  minute
  # @param:  second
  # @return: JD
  #=========================================================================
  def gc2jd1(year, month, day, hour, min, sec)
    # 1月,2月は前年の13月,14月とする
    if month < 3
      year  -= 1
      month += 12
    end
    # 日付(整数)部分計算
    jd  = (365.25 * year).truncate \
        + (year / 400.0).truncate \
        - (year / 100.0).truncate \
        + (30.59 * (month - 2)).truncate \
        + day \
        + 1721088.5
    # 時間(小数)部分計算
    f  = (hour + min / 60.0 + sec / 3600.0) / 24.0
    return jd + f
  rescue => e
    raise
  end

  #=========================================================================
  # グレゴリオ暦 -> ユリウス日(JD)
  #
  # * フリーゲルの公式ではない計算式を使用する
  #   [The Jurian Pieiod](http://www.tondering.dk/claus/cal/julperiod.php#formula)
  #     a = int((14 - month) / 12)
  #     y = year + 4800 - a
  #     m = month + 12a - 3
  #    JD = day + int((153m + 2) / 5) + 365y
  #       + int(y / 4) - int(y / 100) - int (y / 400) - 32045.5
  #   ※上記の int(x) は厳密には、x を超えない最大の整数
  #   ※正子を整数とする JDN を求めるなら、 32045.5 を 32045 に置き換える
  #   ※修正ユリウス日 MJD を求めるなら、 32045.5 を 2432046 に置き換える
  #
  # @param:  year
  # @param:  month
  # @param:  day
  # @param:  hour
  # @param:  minute
  # @param:  second
  # @return: JD
  #=========================================================================
  def gc2jd2(year, month, day, hour, min, sec)
    # 日付(整数)部分計算
    a = ((14 - month) / 12.0).truncate
    y = year + 4800 - a
    m = month + 12 * a - 3
    jd = day \
       + ((153 * m + 2) / 5.0).truncate \
       + 365 * y \
       + (y / 4.0).truncate \
       - (y / 100.0).truncate \
       + (y / 400.0).truncate \
       - 32045.5
    # 時間(小数)部分計算
    f = (hour + min / 60.0 + sec / 3600.0) / 24.0
    return jd + f
  rescue => e
    raise
  end

  #=========================================================================
  # グレゴリオ暦 -> ユリウス日(JD)
  #
  # * フリーゲルの公式ではない計算式を使用する(Celes Trak)
  #   JD = 367.0 * year
  #      - int((7 * (year + int((month + 9) / 12.0))) * 0.25)
  #      + int(275 * month / 9.0)
  #      + day + 1721013.5
  #   ※修正ユリウス日 MJD を求めるなら、 JD -= 678987.0 とする
  #   ※上記の int(x) は厳密には、x を超えない最大の整数
  #   ※修正ユリウス日 MJD を求めるなら、 32045.5 を 2432046 に置き換える
  #
  # @param:  year
  # @param:  month
  # @param:  day
  # @param:  hour
  # @param:  minute
  # @param:  second
  # @return: JD
  #=========================================================================
  def gc2jd3(year, month, day, hour, min, sec)
    # 日付(整数)部分計算
    jd = 367.0 * year \
       - ((7 * (year + ((month + 9) / 12.0).truncate)) * 0.25).truncate \
       + (275 * month / 9.0).truncate \
       + day + 1721013.5
    # 時間(小数)部分計算
    f = (hour + min / 60.0 + sec / 3600.0) / 24.0
    return jd + f
  rescue => e
    raise
  end
end

VerifyGc2Jd.new.exec if __FILE__ == $0

2. 検証用 Ruby スクリプトの実行

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

$ ./verify_gc2jd_2.rb
2016-02-28 23:59:55UTC -> JD: 2457447.49994213, 2457447.49994213, 2457447.49994213
2016-02-28 23:59:56UTC -> JD: 2457447.49995370, 2457447.49995370, 2457447.49995370
2016-02-28 23:59:57UTC -> JD: 2457447.49996528, 2457447.49996528, 2457447.49996528
2016-02-28 23:59:58UTC -> JD: 2457447.49997685, 2457447.49997685, 2457447.49997685
2016-02-28 23:59:59UTC -> JD: 2457447.49998843, 2457447.49998843, 2457447.49998843
2016-02-29 00:00:00UTC -> JD: 2457447.50000000, 2457447.50000000, 2457447.50000000
2016-02-29 00:00:01UTC -> JD: 2457447.50001157, 2457447.50001157, 2457447.50001157
2016-02-29 00:00:02UTC -> JD: 2457447.50002315, 2457447.50002315, 2457447.50002315
2016-02-29 00:00:03UTC -> JD: 2457447.50003472, 2457447.50003472, 2457447.50003472
2016-02-29 00:00:04UTC -> JD: 2457447.50004630, 2457447.50004630, 2457447.50004630
2016-02-29 00:00:05UTC -> JD: 2457447.50005787, 2457447.50005787, 2457447.50005787
...
2016-02-29 23:59:55UTC -> JD: 2457448.49994213, 2457448.49994213, 2457448.49994213
2016-02-29 23:59:56UTC -> JD: 2457448.49995370, 2457448.49995370, 2457448.49995370
2016-02-29 23:59:57UTC -> JD: 2457448.49996528, 2457448.49996528, 2457448.49996528
2016-02-29 23:59:58UTC -> JD: 2457448.49997685, 2457448.49997685, 2457448.49997685
2016-02-29 23:59:59UTC -> JD: 2457448.49998843, 2457448.49998843, 2457448.49998843
2016-03-01 00:00:00UTC -> JD: 2457448.50000000, 2457448.50000000, 2457448.50000000
2016-03-01 00:00:01UTC -> JD: 2457448.50001157, 2457448.50001157, 2457448.50001157
2016-03-01 00:00:02UTC -> JD: 2457448.50002315, 2457448.50002315, 2457448.50002315
2016-03-01 00:00:03UTC -> JD: 2457448.50003472, 2457448.50003472, 2457448.50003472
2016-03-01 00:00:04UTC -> JD: 2457448.50004630, 2457448.50004630, 2457448.50004630
2016-03-01 00:00:05UTC -> JD: 2457448.50005787, 2457448.50005787, 2457448.50005787
...
2016-12-31 23:59:55UTC -> JD: 2457754.49994213, 2457754.49994213, 2457754.49994213
2016-12-31 23:59:56UTC -> JD: 2457754.49995370, 2457754.49995370, 2457754.49995370
2016-12-31 23:59:57UTC -> JD: 2457754.49996528, 2457754.49996528, 2457754.49996528
2016-12-31 23:59:58UTC -> JD: 2457754.49997685, 2457754.49997685, 2457754.49997685
2016-12-31 23:59:59UTC -> JD: 2457754.49998843, 2457754.49998843, 2457754.49998843
2017-01-01 00:00:00UTC -> JD: 2457754.50000000, 2457754.50000000, 2457754.50000000
2017-01-01 00:00:01UTC -> JD: 2457754.50001157, 2457754.50001157, 2457754.50001157
2017-01-01 00:00:02UTC -> JD: 2457754.50002315, 2457754.50002315, 2457754.50002315
2017-01-01 00:00:03UTC -> JD: 2457754.50003472, 2457754.50003472, 2457754.50003472
2017-01-01 00:00:04UTC -> JD: 2457754.50004630, 2457754.50004630, 2457754.50004630
2017-01-01 00:00:05UTC -> JD: 2457754.50005787, 2457754.50005787, 2457754.50005787
...
2017-02-28 23:59:55UTC -> JD: 2457813.49994213, 2457813.49994213, 2457813.49994213
2017-02-28 23:59:56UTC -> JD: 2457813.49995370, 2457813.49995370, 2457813.49995370
2017-02-28 23:59:57UTC -> JD: 2457813.49996528, 2457813.49996528, 2457813.49996528
2017-02-28 23:59:58UTC -> JD: 2457813.49997685, 2457813.49997685, 2457813.49997685
2017-02-28 23:59:59UTC -> JD: 2457813.49998843, 2457813.49998843, 2457813.49998843
2017-03-01 00:00:00UTC -> JD: 2457813.50000000, 2457813.50000000, 2457813.50000000
2017-03-01 00:00:01UTC -> JD: 2457813.50001157, 2457813.50001157, 2457813.50001157
2017-03-01 00:00:02UTC -> JD: 2457813.50002315, 2457813.50002315, 2457813.50002315
2017-03-01 00:00:03UTC -> JD: 2457813.50003472, 2457813.50003472, 2457813.50003472
2017-03-01 00:00:04UTC -> JD: 2457813.50004630, 2457813.50004630, 2457813.50004630
2017-03-01 00:00:05UTC -> JD: 2457813.50005787, 2457813.50005787, 2457813.50005787

やはり、結果として一致した。

グレゴリオ暦からユリウス日を計算するのによく「フリーゲルの公式」を使用するが、別の計算式でも問題ないことが分かりました。

Web 上で検索すると他にも計算式がヒットしますが、どれも同じでしょうね。（単に計算式をどう変形させるかだけの違い）

以上。