十進 BASIC で MJD
〜 修正ユリウス日を計算 〜
2023-07-07 作成 福島
TOP > asob > decbasic-mjd
[ TIPS | TOYS | OTAKU | LINK | MOVIE | CGI | AvTitle | ConfuTerm | HIST | AnSt | Asob ]

日付の加減算に関する問題

去年 (1 年前) の日付は「今年の西暦 - 1」を計算すれば求められます。
来年 (1 年後) の日付は「今年の西暦 + 1」ですね。
例: 2023-7-7 の 1 年前は 2022-7-7。

同様に、昨日の日付は「今日 - 1」で求められます。
例: 2023-7-7 の 1 日前は 2023-7-6。

では、10,950 日前の日付はどうでしょうか。(10,950 = 365 × 30)
本稿記述時の日付が 2023-7-7 なので、30 年前の日付は 1993-7-7。
···というわけには行きません。

1993 年 ~ 2023 年の間には「うるう年」が 7 回もあります。(1996, 2000, 2004, 2008, 2012, 2016, 2020)
このため、10,950 日前の日付は 1993-7-14 となります。

うるう年には、以下のような複雑なルールがあります。
  1. 4 年に一度はうるう年にする。
  2. 100 年に一度はうるう年にしない。
  3. 400 年に一度はうるう年にする。
西暦 2000 年は、上記 2 でうるう年にしないはずが、同 3 によってうるう年になります。
年月日の加減算に関する問題
時刻の計算をしてみます。

2 時間 35 分の映画が 18 時 30 分から上映されたら帰宅できるでしょうか。(例: シン・エヴァンゲリオン劇場版:||)
(途中で退席するとかはナシで)
答えは 21:05 なので、終電には余裕がありそうということが分かります。
未成年者の場合は、保護者の帯同が必要になるでしょう。

時計では 18 時 30 分 + 2 時間 35 分 = 21 時 5 分という計算をします。

同様に日付を計算して 2023-7-7 の 10 年 8 カ月 30 日後の日付はどうなるでしょう。
答えは 2033-3-4 ですが、上記うるう年のルールがあるため、即座に計算できる人はかなり少ないと思います。
(加算する 10 年 8 カ月 30 日間の中にもうるう年があることに注意)
それでも日付を計算したい
日付が(れき)のままだと計算が複雑になります。
一旦、基準日をかなり昔に設定して、そこから何日経過したか日数を揃えることにより、計算を単純にします。
その便利な基準日の定義が「修正ユリウス日」(MJD: Modified Julian Day) です。
「ユリウス通日」がオリジナルですが、そのままだと使いにくいので、これを使いやすくしたのが「修正ユリウス日」です。
「ユリウス~」という名称ですが、シーザーとは関係ありません。考案者は古典・年代学者スカリゲルで、その (天文・哲学者) の名前です。

すべての日付を「修正ユリウス日」という単純な整数に変換し、加減算が済んだのちに通常日付 (グレゴリオ暦) に戻す、という使い方をします。


十進 BASIC で実装

アルゴリズムは Wikipedia からのパクリです。
十進 BASIC の場合、INT() の挙動が他の多くのプログラミング言語と異なり、床関数となります。
Python や C 言語等では int() が小数部の切り捨てになります。

十進 BASIC の関数は、一度に複数の値を返すことができないため、
ここでは 3 つの値 y,m,d (それぞれ、年,月,日を表す) をひとつの数値に詰め込んで返しています。
このため、返戻値を受領してから 3 つの値に分解しています。

y,m,d は整数です。時刻を考慮していません。

check-mjd.BAS
DECLARE EXTERNAL FUNCTION greg2mjd
DECLARE EXTERNAL FUNCTION jd2mjd
DECLARE EXTERNAL FUNCTION mjd2jd
DECLARE EXTERNAL FUNCTION mjd2greg
DECLARE EXTERNAL FUNCTION mjd2week
DECLARE EXTERNAL FUNCTION mjd2zodiac
DECLARE EXTERNAL FUNCTION mjd2stems

LET y = 2012
LET m=1
LET d=1
PRINT "Gregorian y,m,d"; y; m; d

LET mjd = greg2mjd(y,m,d)
LET w = mjd2week(mjd)
PRINT "mjd";mjd;mid$("WedThuFriSatSunMonTue",1+w*3,3)

LET ymd = mjd2greg(mjd)
! 返戻値を 3 つに分解する。
LET y = INT(ymd / 10000)
LET m = MOD(INT(ymd/100),100)
LET d = MOD(ymd,100)
PRINT "Gregorian y,m,d"; y; m; d

PRINT "---"

LET y = 1582
LET m = 2
LET d = 1
PRINT "Julian y,m,d"; y; m; d

LET mjd = jd2mjd(y,m,d)
PRINT "mjd"; mjd

LET ymd = mjd2jd(mjd)
! 返戻値を 3 つに分解する。
LET y = INT(ymd / 10000)
LET m = MOD(INT(ymd/100),100)
LET d = MOD(ymd,100)
PRINT "Julian y,m,d"; y; m; d

END


MERGE "mjd.BAS"

mjd.BAS (上記 check-mjd.BAS と同じ場所に保存する)
! mjd.BAS version 1.0 written by fuku@rouge.gr.jp
!
! ユリウス日 / 修正ユリウス日 / グレゴリオ暦 の相互変換を行う。
!
! https://ja.wikipedia.org/wiki/ユリウス通日#修正ユリウス日(MJD)
! より。

! Modified Julian Day (MJD)
!
!Mutual conversion:
!    Julian Date / Modified Julian Day / Gregorian Calendar


! Convert Gregorian-calendar to MJD. EXTERNAL FUNCTION greg2mjd(y,m,d) IF m <= 2 THEN LET m = m + 12 LET y = y - 1 END IF LET greg2mjd = INT(365.25*y) + INT(y/400) - INT(y/ 100) & & + INT(30.59*(m-2)) + d - 678912 END FUNCTION
! Convert Julian-date to MJD. EXTERNAL FUNCTION jd2mjd(y,m,d) IF m <= 2 THEN LET m = m + 12 LET y = y - 1 END IF LET jd2mjd = INT(365.25*y) + INT(30.59*(m-2)) + d - 678914 END FUNCTION
! Convert MJD to Julian-date. EXTERNAL FUNCTION mjd2jd(days) LET n = days + 678883 LET a = 4 * n + 3 LET b = 5 * INT(MOD(a,1461)/4) + 2 LET y = INT(a/1461) LET m = INT(b/153) + 3 LET d = INT(MOD(b,153)/5) + 1 IF m > 12 THEN LET m = m - 12 LET y = y + 1 END IF LET mjd2jd = y * 10000 + m * 100 + d END FUNCTION
! Convert MJD to Gregorian-date. EXTERNAL FUNCTION mjd2greg(days) LET n = days + 678881 LET a = 4 * n + 3 + 4 * INT(3 / 4 * INT(4*(n+1)/146097+1)) LET b = 5 * INT(MOD(a,1461)/4) + 2 LET y = INT(a/1461) LET m = INT(b/153) + 3 LET d = INT(MOD(b,153)/5) + 1 IF m > 12 THEN LET m = m - 12 LET y = y + 1 END IF LET mjd2greg = y * 10000 + m * 100 + d END FUNCTION
! Get day of week from MJD. EXTERNAL FUNCTION mjd2week(days) LET mjd2week = MOD(days, 7) ! 0:Wed, 1:Thu, .., 6:Tue END FUNCTION
! Get Oriental-Zodiac from MJD. EXTERNAL FUNCTION mjd2zodiac(days) LET mjd2zodiac = MOD(days, 12) ! 0:Tora, 1:U, 2:Tatsu, .., 11:Ushi END FUNCTION
! Get Heavenly-Stems from MJD. EXTERNAL FUNCTION mjd2stems(days) LET mjd2stems = MOD(days, 10) ! 0:Kinoe, 1:Kinoto, .., 9:Mizunoto END FUNCTION