unix time

前回のunix timeのブログに書いたように, Calendricalの方法は, 確かに簡単ではあるが, 1から勘定を始めるプログラムに, 多少の違和感を持っていた私は, やはり0から始めるプログラムが書いてみたく, そのため, 366, 365, 365, 365のように, おなじものが並ぶ後続の前に, それらより1だけ多いか少ない要素が先頭にあるときの処理について考えていた.

つまり例えば1600年1月1日から1999年12月31日の400年間, 146097日には, 最初の100年は1600年がうるう年なので, 36525日あり, 後の3回の100年は, 最初の00年が平年なので, 36524日であって, 36525, 36524, 36524, 36524であり, またその36524日ある100年を4年ごと25組に分けると, 第2組以降は最初の4年が1日少なく, 1460, 1461, 1461, ..., 1461となる.

このような問題提起について, d日目を入力し, その該当する年yと, その年内の日数d'を計算する方法を探した.

何万日を扱うのも憂鬱なので, とりあえず, 先頭が1多い4日, 3日, 3日, 3日のパターンと, 先頭が1少ない, 3日, 4日, 4日, 4日のパターンを考える.

4,3,3,3のパターンでは, 総日数は13日で, 0≤d<13に対し,

d         y  d'
0,1,2,3   0  0,1,2,3
4,5,6     1  0,1,2
7,8,9     2  0,1,2
10,11,12  3  0,1,2

また, 3,4,4,4のパターンでは, 総日数は15日で, 0≤d<15に対し,

d             y  d'
0,1,2         0  0,1,2
3,4,5,6       1  0,1,2,3
7,8,9,10      2  0,1,2,3
11,12,13,14   3  0,1,2,3

となるようにしたい.

前者については

(define (foo d)
(let ((y (floor (/ d 3.25))))
(list y (+ d -13 (floor (* 3.25 (- 4 y)))))))

とする. 13は総日数, 4は年の幅, 怪しげな3.25は13/4である.

後者については

(define (bar d)
(let ((y (- 3 (floor (/ (- 14 d) 3.75)))))
(list y (- d (floor (* y 3.75))))))

とする. 14は総日数-1, 3.75は総日数の15を4で割ったものだ.

とりあえず, これでdに対して(y d')を求めてみると

(map foo (a2b 0 13))
=>((0 0.) (0. 1.) (0. 2.) (0. 3.) (1. 0.) (1. 1.) (1. 2.)
(2. 0.) (2. 1.) (2. 2.) (3. 0.) (3. 1.) (3. 2.))

(map bar (a2b 0 15))
=>((0. 0.) (0. 1.) (0. 2.) (1. 0.) (1. 1.) (1. 2.) (1. 3.)
(2. 0.) (2. 1.) (2. 2.) (2. 3.) (3. 0.) (3. 1.) (3. 2.) (3 3.))

とうまく行きそうである.

fooとbarの図を描いてみるとこうなる. 上がfoo, 下がbar.



図でははっきりしないが, dに対するfloorを取る前の値は
(0 .31 .62 .92 1.23 1.54 1.85 2.15 2.46 2.77 3.08 3.38 3.69)
と
(3.73 3.47 3.2 2.93 2.67 2.4 2.13 1.87 1.6 1.33 1.07 .8
.53 .27 0)
である.

そこで, この3とか4とかを実際の値に変えて関数を書き, クリティカルな値に対してテストすると次の通りだ.

(define (foo d)
(let ((y (floor (/ d 36524.25))))
(list y (+ d -146097 (floor (* 36524.25 (- 4 y)))))))
(map foo '(0 36524 36525 73048 73049 109572 109573 146096))
=>((0 0.) (0. 36524.) (1. 0.) (1. 36523.) (2. 0.) (2. 36523.)
  (3. 0.) (3. 36523.))

(define (bar d)
(let ((y (- 24 (floor (/ (- 36523 d) 1460.96)))))
(list y (- d (floor (* y 1460.96))))))
(map bar '(0 1459 1460 2920 2921 4381 4382 5842 33602 35062
  35063 36523))
=>((0. 0.) (0. 1459.) (1. 0.) (1. 1460.) (2. 0.) (2. 1460.)
  (3. 0.) (3. 1460.) (23. 0.) (23. 1460.) (24. 0.) (24 1460.))

(define (foo d)
(let ((y (floor (/ d 365.25))))
(list y (+ d -1461 (floor (* 365.25 (- 4 y)))))))
(map foo '(0 365 366 730 731 1095 1096 1460))
=>((0 0.) (0. 365.) (1. 0.) (1. 364.) (2. 0.) (2. 364.)
  (3. 0.) (3. 364.))

あとはこれを繋げるだけ. 注意すべきは, Calendricalの方は1年1月1日から始まるのに対し, こちらは0年から399年までの400年から始まるので, 0年分の366日を足す; Calendricalは1月1日のfixed dateが1なのに対し, こちらは, 0年1月1日が0なので1日分の修正が必要だ. この修正の結果がd0である. 400年のところは, 常に146097日なので, 計算は簡単だ. 残りは上の関数を利用している.

(define (gr rd)
(let* ((d0 (+ rd 366 -1))
       (y400 (quotient d0 146097))
       (d1 (modulo d0 146097))
       (y100 (floor (/ d1 36524.25)))
       (d2 (+ d1 -146097 (floor (* 36524.25 (- 4 y100)))))
       (y4 (if (= y100 0) (quotient d2 1461)
           (- 24 (floor (/ (- 36523 d2) 1460.96)))))
       (d3 (if (= y100 0) (modulo d2 1461)
           (- d2 (floor (* y4 1460.96)))))
       (y1 (if (and (> y100 0) (= y4 0)) (quotient d3 365)
           (floor (/ d3 365.25))))
       (d4 (if (and (> y100 0) (= y4 0)) (modulo d3 365)
           (+ d3 -1461 (floor (* 365.25 (- 4 y1))))))
       (y (+ (* y400 400) (* y100 100) (* y4 4) y1))
       (leap (if (= (modulo y 100) 0) (= (modulo y 400) 0)
           (= (modulo y 4) 0)))
       (e (if leap
        '(0 31 60 91 121 152 182 213 244 274 305 335)
        '(0 31 59 90 120 151 181 212 243 273 304 334)))
       (mon (apply + (map (lambda (f) (if (<= f d4) 1 0)) e)))
       (d (- d4 (list-ref e (- mon 1)) -1)))
    (list y mon d)))

テスト

(gr 1) => (1. 1 1.)
(gr 734274) => (2011. 5 16.)

こんな具合いだが, まぁいいか.