Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Ruby script to verify of converting Gregorian Calendar to Julian Day.(v2)
#! /usr/local/bin/ruby
# coding: utf-8
#=グレゴリオ暦 -> ユリウス日 変換の検証(UTC ベース)
#
# : 2種類の計算式で グレゴリオ暦 -> ユリウス日 の変換を行い、
# 相違がないかを検証する。
#
# 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
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.