Ruby script to calc a Mean Obliquity of the Ecliptic.

mean_obliquity_ecliptic.rb

#! /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