Skip to content

Instantly share code, notes, and snippets.

@komasaru
Last active April 19, 2018 06:09
Show Gist options
  • Save komasaru/cfdab2ecfe998661f771643504ac875f to your computer and use it in GitHub Desktop.
Save komasaru/cfdab2ecfe998661f771643504ac875f to your computer and use it in GitHub Desktop.
Ruby script to calc a Mean Obliquity of the Ecliptic.
#! /usr/local/bin/ruby
# coding: utf-8
#---------------------------------------------------------------------------------
#= 平均黄道傾斜角の計算
#
# date name version
# 2016.05.26 mk-mode.com 1.00 新規作成
#
# Copyright(C) 2016 mk-mode.com All Rights Reserved.
#---------------------------------------------------------------------------------
# 引数 : 日時(TT(地球時))
# 書式:YYYYMMDD or YYYYMMDDHHMMSS
# 無指定なら現在(システム日時)を地球時とみなす。
#---------------------------------------------------------------------------------
#++
require 'date'
class MeanObliquityEcliptic
def initialize
get_now
get_arg
end
def calc
jd = calc_jd(@tt)
t = calc_t(jd)
moe = calc_eps_a(t)
puts " 地球時(TT): #{@tt.strftime("%Y-%m-%d %H:%M:%S")}"
puts " ユリウス日(JD): #{jd}"
puts " ユリウス世紀数(JD): #{t}"
puts "平均黄道傾斜角(eps_a): #{moe} °"
rescue => e
$stderr.puts "[#{e.class}] #{e.message}"
e.backtrace.each { |tr| $stderr.puts "\t#{tr}"}
exit 1
end
private
# 現在日時の取得
def get_now
t = Time.now
@tt = Time.new(t.year, t.month, t.day, t.hour, t.min, t.sec)
rescue => e
raise
end
# コマンドライン引数の取得
def get_arg
return unless arg = ARGV.shift
exit 0 unless arg =~ /^\d{8}$|^\d{14}$/
year, month, day = arg[ 0, 4].to_i, arg[ 4, 2].to_i, arg[ 6, 2].to_i
hour, min, sec = arg[ 8, 2].to_i, arg[10, 2].to_i, arg[12, 2].to_i
exit 0 unless Date.valid_date?(year, month, day)
exit 0 if hour > 23 || min > 59 || sec > 59
@tt = Time.new(year, month, day, hour, min, sec)
rescue => e
raise
end
# ユリウス日の計算
def calc_jd(tt)
year, month, day = tt.year, tt.month, tt.day
hour, min, sec = tt.hour, tt.min, tt.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.5
# 時間(小数)部分計算
t = (sec / 3600.0 + min / 60.0 + hour) / 24.0
return jd + t
rescue => e
raise
end
end
# ユリウス世紀数の計算
def calc_t(jd)
return (jd - 2451545) / 36525.0
rescue => e
raise
end
# 黄道傾斜角(εA)の計算
def calc_eps_a(t)
return (84381.406 + \
( -46.836769 + \
( -0.0001831 + \
( 0.00200340 + \
( -5.76 * 10 ** -7 + \
( -4.34 * 10 ** -8) \
* t) * t) * t) * t) * t) / 3600.0
rescue => e
raise
end
end
MeanObliquityEcliptic.new.calc if __FILE__ == $0
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment