1. はじめに

まずは、ダウンロードしていただきましてありがとうございます。 本サンプル スクリプトは、 任意の日付（新暦）に対する旧暦と六曜を計算して、出力すると いうものです。機能としてはこれだけで、 このほかは特にこれといって何もあり ません。

このままでは、 あまりに粗末なので胸を張って『旧暦計算スクリプト』と言えま せんので、『サンプルスクリプト』と名称をつけました。 皆さんの創意工夫によ って、 より良いプログラムに育ててやって下さいませ。 本サンプルスクリプト が、その足掛かりになれば幸いと存じます...。 この説明書は、 『サンプルスクリプト』の参考資料として作成しました。 ま た、 本サンプルスクリプトの動きを理解するために必要な旧暦の計算方法や暦に 関する知識などを解説します。

本スクリプトは、MS-DOS版のjgawk 2.11.1 + 3.0 および 2.15.2 + 1.1で動作 確認しました。[*1]（jgawk は、serow氏によって配布されているGNU awk 日本語 (shift-jis)対応版です。）他の版では、 「たぶん動作するでしょう...」としか 言えません。何しろ、動作確認していませんので。（悪しからず...） なお、jgawkは本アーカイブには含まれておりません。お持ちでない方は実行形式 は以下の場所にありますので、2.11.1 + 3.0 または 2.15.2 + 1.1のいずれかを ダウンロードなさって下さいませ。 （私は、NIFTYとascii-netのIDしかありませ んので、他のNETについては解りません。） NIFTY-serve ascii-net のIDを持ちでない方は、月刊ASCII 1993月6号のお楽し みディスクを探してみてください。 もしも、見つかったならば \JGAWK\JGA_EXE.LZH を展開し、JGAWK.EXE と言う名 称の ファイルがありますのでこれ使用します。 （jgawk2.11.1+3.0(MSDOS,SJIS)が収録されている）

       NIFTY-Serve:

       jgawk 2.15.2 + 1.1 ﾌﾙｾｯﾄ ...
        FGALEL/LIB 12  #29   93/09/30  267466 B JGAWK11 .LZH

       jgawk2.11.1+3.0(MSDOS,SJIS)...
        FGALAP/LIB  7  #268  92/12/25  185212 B JGA_EXE .LZH

                            --------- + ---------

       ascii net (POOL MS-DOS):

       Japanized Gnu Awk (jgawk) 2.15.2 + 1.1 ,..
        7116 jgawk11.lzh       B pcs19614 93/09/25  267466

       jgawk 2.11.1 + 3.0 MSDOS executable ...
        6195 jga_exe.lzh       B pcs19614 92/12/24  185212

2. 起動方法

計算する日付を与えて計算させる方法と、 システム時計の日付を使って計算さ せる方法があります。以下の例を参考にして下さい。 実行時間がものすごくかかるので、（PC-9801DA（386DX）、jgawk 2.11.1+3.0 で、 約２０秒かかる）「おや？ハングアップしたのかなぁ～」と思わないで下さ い。はっきり言って効率の悪いアルゴリズムが原因です。

    1.日付を与えて計算させる方法

  たとえば、1994年11月8日の旧暦を求める場合には、以下のように入力します。

       A:\>jgawk -fqreki.awk 1994 11 8
       西暦1994年 11月  8日は、旧暦1994年 10月  6日 先負です。

    2.システム時計の日付を使って計算させる方法

       A:\>jgawk -fqreki.awk
       西暦1994年  5月  1日は、旧暦1994年  3月 21日 大安です。
       （システム時計は、5月 1日と仮定）

3. 旧暦とは...

日本で一般的に旧暦といいますと、明治 5年（1972年）に改暦される以前に使 われていた天保暦法[*2]をさします。 以降の説明では単に旧暦と言った場合に は、「天保暦法」による「こよみ」をさすと考えて下さい。 本スクリプトでは、太陽や月の位置を計算する方法や使用する定数が、 天保暦 法で規定されたものと違います。 それでは、「ニセモノ」と呼ばれそうですが旧 暦と言っているのはほとんど全て（といって間違い無い）がこのようなもので す。 （いい加減と言えばいい加減とは思いますが、 一般に広まってしまったた め、いまさら昔のやり方に戻す訳にもいかないでしょう...） 天保暦で規定された方法との相違点を列挙いたしますと...

i 中気［ちゅうき］の時刻や、朔［さく］の時刻を計算する方法の違い。

天保暦が使用されていた当時は、経験的に知られていた定数・周期に基づい て計算されていました。（中気、朔についての説明は『旧暦の計算に関連す る用語の説明』のところで行います） 本スクリプトでは天体力学に基づく略算式（計算方法については後述いたし ます）を使用して計算しています。

ii 時刻制度の違い。

天保暦が使われていた当時は、 京都における真太陽時（真太陽時と言うの は、おおざっぱに言って、日時計が示す時刻と同じです。当時、真太陽時が あったという訳ではなくて、現代流に解釈した場合にこれに該当すると言う 話です。）を使用していました。 本スクリプトでは、 日本標準時（詳しい内容については、 後述いたしま す。）を使用しています

このため、 天保暦法で規定された方法・定数に従って厳密に計算した日付と（場 合によっては）１日前後したり、月の大小が異なる可能性があります。

4. 閏月［うるうつき］と太陰太陽暦について

天保暦法は一般的に、 太陰太陽暦に分類されています。太陰太陽暦は、日付を 月の満ち欠け（朔望［さくぼう］とも言う）の周期を基準として、 月名を１年の 周期を基準として作られた暦です。

すなわち、 月の満ち欠けの周期（さくぼうげつ＝朔望月）は平均 29日半（29. 530589日）で、この12倍が 1年に近い事をもとにして組み立てられています。

           12 朔望月 =   29.530589 ［日］* 12 ［ヶ月］
                     =  354.36706  ［日］

ところが、このまま（旧暦の）1年を 354日または 355日としたのでは、 1年当 たりで10日ないし11日も日付が違ってきてしまいます。 このまま何もせずに誤差 をため込んでおくと、約 3年で 1月分に相当するため[*3]、1年を 13ヶ月として 暦と季節のずれを解消していました。 （仮に、何もせずそのものにしておくと次 第に季節と月名が合わなくなってしまい、1月だというのに真夏といった、奇妙な 現象を体験する事になります。） 暦と季節とを合わせるために、挿入した 1ヶ月を閏月［うるうつき］と呼んでい ます。 閏月の月名は、前月の月名の前に「閏」をつけて、閏XX月のように呼びま す。（閏月の置き方の規則は、後述します。）このように、 閏月によって季節と 合わせていた暦が太陰太陽暦です。

               3*（太陽年）   =  3［年］* 365.2422［日］
                              =  1095.7266 ［日］ ...... Ａ

               3*（12 朔望月）=  3［年］* 12 ［ヶ月］* 29.530589 ［日］
                              =  1063.101204 ［日］..... Ｂ

                 Ａ  －  Ｂ   =  1095.7266 ［日］－ 1063.101204 ［日］
                              =  32.625396 ［日］

5. 「こよみ」について

「こよみ」といいますと、 一般的にはカレンダーや日めくりとか運勢暦などを 連想する方が多いと思います。中国では、古くから天体の位置（月、太陽、 惑星 etc...）を掲載した書物も「こよみ」と呼んでいました。 中国から暦を輸入して 来た日本でも、同様です。（本説明書では、 天体に位置などを掲載した書物を特 に『天体暦』と言って一般の「こよみ」とは区別する事にします。）

6. 月齢について

よく旧暦の日付と月齢を混同なさる方がおられますが、 この二つは定義が違い ますので、全く別物です。 ここでは、月齢がどのようにして決められているのか、 簡単に紹介しておく事に しておきます。

月齢は新月から当日の正午（日本標準時）までの時間を、 「日」の単位で表し たものです。例えば、1994年 5月 1日の場合について計算してみますと...

       月齢＝（計算対象の時刻）－（新月の時刻）
           ＝ 1994年 05月 01日 12時 00分 －1994年 04月 11日 9時 18分
           ＝ 20日 02時間 42分
           ＝ 20.1125

といった訳で答えは、20.1日となります。

さらにもうひとつ、満月が必ず旧暦の１５日とは限りませんし、 月齢１４日と も限りません。 これらの数はあくまで平均値で、時よって±２日位の間で変化し ています。

7. 六曜について

いわゆる「先勝・友引・先負・仏滅・大安・赤口」の六つの繰り返しの事で す。『暦の百科事典』によれば... 六曜とは、六曜星の略で、六輝［ろっき］とも言う。江戸時代の暦書では、 孔 明六曜星等と記されていた。 もともとは時刻の吉凶占いに用いられて来たが、こ れが日の吉凶を占うものとなった。 孔明六曜星と呼ばれることからもわかるように、この六曜は三国史の名将 諸葛亮 孔明が発明したものと伝えられる。 孔明はこの六曜星を使って軍略をたてたとこ ろ、 こととごとく成功したという。しかし、戦国時代から六曜があったかは疑わ しい。後世のこじつけというのが定説である。 六曜の起源についてはよくわかっていない。たぶん、 もともと一ヶ月を五で割 り六日づつの小単位をつくり、 その日を数えるためのものであったのではないか と考えられる。すなわち、日の吉凶を示すものではなく、七曜と同じように、 単 に日を数えるための記号であったのではないかと考えられる。 とはいっても、 中国で発生したものには違いなく、日本へは鎌倉末期から室町 時代に伝わったと考えられる。 これがどうしたものか、 江戸時代の終わり頃か ら、民間でひそかに流行し始めた。 そして、明治の太陽暦が実施されるに及んで、皮肉なことに、 この六曜が「おば け暦」[*4]として、 人々の口々にのぼるようになり、第二次世界大戦後になって 大流行し、現在にいたっている。 六曜の吉凶を良く調べてみると、一つの特徴がある。それは、 勝負事にこだわ っている点である。 おそらく、始めのうちは、遊郭や賭場などの遊び人、勝負師 の間で使われていたのだろうと言う説もある。

六曜の選日法は、 旧暦の各月の朔日［ついたち］を「先勝・友引・先負・仏滅 ・大安・赤口」のうち何れかを配し、 翌日の２日から晦[*5]まで先の順番で繰っ て行く事になっています。従って、月が変わるとその前後で、 順番が狂う事があ ります。 また、閏月の場合には、前月と同じ規則で計算します。（例えば、閏３ 月朔日では３月と同じく、先負となります。）どの月が何にあたるのかは、 以下 の表を参照して下さい。

本スクリプトでは、［旧暦の月］＋［旧暦の日］－２ を６で割った剰余（あま り）の値をみて、零なら先勝、１なら友引.... といった具合で計算しています。

以下に六曜の吉凶の説明をします。 念のため、［］の中に一般的な読み方を記 入しておきました。

7.1 先勝［せんかち］

先づれば即ち勝つ、の意。もともとは、速喜とか、 則吉などと書かれていたの で、 「早ければ吉」・「急ぐ事良し訴訟事良し」などどと言われている。また、 午前中は吉で午後悪しとされている。

7.2 友引［ともびき］

凶事に友を引く、の意。古くは、「勝負無き日と知るべし」とあるので、 もと もとは何事も引き分けで、勝負のつかぬ日とされていた。 どの暦の解説にも、 「この日に葬儀を出す事はつつしむべし」とある。また、朝 晩は吉で正午のみ凶とされている。

7.3 先負［せんまけ］

先勝の逆で、先づれば即ち負ける、の意。従って、 勝負事や急ぐ事はなるべく 避けて、相手の出方を待つのが良いとされている。 また、朝から昼までは凶で、昼すぎから日暮れまでは吉とされている。

7.4 仏滅［ぶつめつ］

仏も滅するような最悪の日、 の意。この日は、六曜の中でも大凶とされ、祝い 事、法事など何事もうまく行かない日と言われている。 また、移転・開店なども忌み禁じられている。最近では、 葬式に関しては友引の み忌み嫌われるが、なぜか仏滅に葬式を出してもいっこう平気である人が多い。

7.5 大安［たいあん］

大いに安し、 の意。大安吉日とも称され、万事に用いても吉とか、成功せざる 事なき日とされている。 こんな所から、大変めでたい日とされ、結婚式の日となった。

7.6 赤口［しゃっく］

赤舌日［しゃくぜつにち］とも称され、 大工・板前など刃物を扱う職業の人達 は、 要注意とされていた。また、赤が火を連想される所から、「火の元に注意せ よ」とも言われていた。 この日は、午の刻（現在の時刻で12時～14時頃）だけが吉、 その外は凶とされ、 特に祝い事は大凶とされていた。

8. 旧暦の計算に関連する用語の説明

ここでは、旧暦の計算方法を説明する際に出てくる用語の説明をします。

8.1 黄道［こうどう］

地球から、 太陽・惑星を見ると決まった道筋を動いているようにみえます。こ の道筋を黄道と呼んでいます。惑星は、 太陽を焦点としてある平面（これを軌道 面という）を公転しています。 （地球の）軌道面と天球[*6]と交わる線を黄道と いいます。従いまして、太陽は黄道にそって移動ているように見えます。 また、 地球の赤道を無限に延長し、天球と交わった線を赤道と言います。

こよみを作る時には、 天が動いていると考えた方が、都合の良い事が多いため 以降の説明では、 太陽等が地球の回りを回っているような『天動説的』表現をす る場合があります。

8.2 春分点［しゅんぶんてん］

黄道と、 赤道との交点で、春分の日に太陽がこの位置にあるためこう呼んでい ます。春分点は、黄経の原点となるため、重要な意味を持っています。

8.3 黄経［こうけい］

月や太陽の位置を、表す座標の値で、 春分点を原点として黄道を基準に測った 経度を言います。東回りに測り、０度～３６０度の値をとります。 たとえで言いますと、 英国のグリニジ天文台の子午儀跡を原点として赤道を基準 に測った角度を東経（西経）ＸＸ度と言うのと同じようなものです。 また、黄道を基準として北極（または、 南極）の方向へ測った緯度を黄緯［こ うい］といいます。 北側を正、南側を負の符号をつけ０度～９０度の値をとりま す。これは、北緯（または、南緯）に対応したものです。

旧暦を計算する際に朔の時刻を求めるために、 月や太陽の黄経の値を計算しま すが、 本来のやり方では（月や地球の）軌道計算を行って求める事になります。 軌道計算は、 そう簡単にはできる代物ではありませんので本スクリプトでは、 『天体位置略算式の解説』（井上圭典、 鈴木邦裕 著 海文堂出版）と言う文献 に、軌道計算の数式解を使用した略算式があるので、 これを利用して求めていま す。

余談になりますが、 「経」、「緯」は織物の経［たていと］と緯［よこいと］ から来ているそうです。 （メルカトル図法の地図（経線と緯線が地図のどこでも 垂直に交わっている地図）を見るとその事が良く解ると思います。

8.4 朔［さく］

朔と言うのは新月の事です。 旧暦では朔を含む日を朔日としています。朔の時 には、太陽と同じ方向に月があるため見えません。 （これは月と太陽が重なって 見えると言う事ではない点に、 注意しなくてはいけません。仮に完全に重なって 見えた場合には、皆既日蝕または金環日蝕となります。） 天体暦の定義では、 月と太陽の黄経の差が零となる瞬間（これを黄経の合［ご う］といいます。）とされています。 本スクリプトでも天体暦と全く同じ定義を 採用しています。

8.5 二十四節気［にじゅうしせっき］、中気［ちゅうき］

もともとは、暦の日付と季節のずれを知る目安として作られ、 １年を２４等分 して約１５日毎に区切っていました。 天保暦になってからは、 太陽黄経が１５度の倍数になる日と定義が変わりまし た。１年間で太陽の黄経は３６０度変化しますので、 これを２４等分すると確か に１５度となりますが、中気と中気との時間は均等に配分されません。 太陽の黄経の変化は地球が太陽の回りを回る公転運動の結果によって起こるの ですが、 楕円軌道上を運動しているためと、他の惑星の引力によって地球の公転 運動が乱される（これを、摂動［せつどう］と言います）ため、 わずかではあり ますが運動の早さに遅速がでてきます。 （以下に示す計算例を参照して下さい） このため、太陽の黄経の変化は一様とならず、 中気と中気との時間は均等になら ない訳です。

二十四節気のうち、 太陽黄経が３０度の倍数になる日を中気［ちゅうき］とい い、その他を節気［せっき］と言います。（天保暦以前の中気は、 １年を１２等 分して約３０日毎に区切ったものを使用していた）暦と季節とを合わせるため に、 閏月を挿入しますがその目安に使用するのが中気です。天保暦以前の暦法で は、中気が無い月を閏月にしていました。 天保暦法になってからは、中気と中気の間隔が不等になったうえに、 中気と中気 との間隔が１朔望月の長さ（平均29.530589日）と比較して短く[*7]なる場合があ ります。 このため、朔と朔の間に中気が２度訪れてしまいその前後に（本来は閏 月ではないのに）中気の含まない月が出現する可能性が出て来てました。 このため、天保暦法では二分二至を含む月を必ず、 ２月・５月・８月・１１月と 決め、 この拘束条件のもとで中気の含まない月を閏月としています。従って、中 気が無いからといって直ちに閏月とはなりません。このあたりが、 閏を置く規則 を複雑にしていると言えます。 もともと中気は、 暦日と季節のずれを知る目安に過ぎませんので、１年を１２等 分して約３０日毎に区切るだけで十分に役目を果たすのに、 天保暦法の作成者は 中気の決め方を改悪し、 閏月の置き方を複雑にしてしまったのではないでしょう か？ さて愚痴はともかく、二十四節気と太陽黄経の対応表を以下にしめします。 こ の中で、 春分・秋分・夏至・冬至の４つは特に、二分二至と言います。（なお、 立春・立夏・立秋・立冬の４つを四立［よんりつ、しりゅう］と言います。）

近日点付近の中気の間隔（冬至～大寒）と、遠日点付近の中気の間隔（夏至～ 大暑）の時間を表にしました。計算は、1994年の値を本スクリプトのサブルーチ ンを使用して行いました。近日点付近では、 中気と中気の間隔が１朔望月の時間 （29.530589日 = 29日12時間44分）よりも短いことが解ると思います。

8.6 ユリウス日（ユリウス通日［つうじつ］）

暦の計算等で、 数年間にわたる２点の日数を計算するために年代学や天文学の 分野で用いられています。原点は、BC4713年 1月 1日（ユリウス暦法を使用して 遡って決めた）で、この日から（太陽暦で）連続して数えられています。

BC4713年がどうして原点になったのかと言いますと、 以下の周期の最少公倍数 （7980年）をとり各周期の第１年目にあたる年としたためです。しかし、 7980年 を経過したらまた１から数えなおすと言う規定は無いようです。

• i 太陽章（28年）ある日の七曜が同じとなる周期です。
• ii 太陰章（19年）ある日の月の欠け方（正確には月の位相）が同じとなる周 期です。（メトンと呼ばれている周期）
• iii インディクション（15年）徴税額査定更正周期で、中世の暦を研究する際 に西暦年号を確定するのに使用されている。

ユリウス通日を計算すら際に注意しなくてはならないのは、 AD1582年10月 4日 （木）までは、 ユリウス暦法で日数を数えて、 その翌日（AD1582年10月15日 （金））以降はグレゴリオ暦法によって数える事になっています。 なお、 AD1582年10月 5日から、 AD1582年10月14日の10日間は、ユリウス暦法か ら、 グレゴリオ暦法に改暦する際に日数を調整するために飛ばされましたので実 在しない日となっています。（なお、曜日については連続させています。） ただし、 実際にはグレゴリオ暦法が採用された日は、（宗教上の理由等で）国よ って異なるため、歴史上の日付に対するユリウス日を求める時には、 太陽暦で示 された暦日であっても、 ユリウス暦法によるものかグレゴリオ暦法によるものか を、常に意識しなければいけません。 また、欧州で旧暦というとユリウス暦法による暦日をさします。

両者とも同じ太陽暦なのですが、閏年の置き方が若干異なっています。 ユリウ ス暦法では、 『西暦年数を４で割り切れる年を閏年とする。』としていました が、グレゴリオ暦法では、『西暦年数を４で割り切れる年を閏年とするが、 １０ ０で割り切れる年については、 ４００で割り切れる年のみを閏年とする。』と改 正されました。 これでは、 「何の事かさっぱり解らんわ！」と言われそうなので、 1800,1900, 2000,2100年の４つの例のを挙げて考えてみましょう。 ユリウス暦法では４つの年とも（４で割り切れるため）閏年となりますが、 グレ ゴリオ暦法では、2000年だけが400で割り切れるために閏年となり、そのほか年は （400で割り切れないため）平年と言った具合になります。

もう少し詳しく説明いたしますと、1年の長さ（厳密には１太陽年と言って、太 陽が春分点に再び戻ってくるまでの時間の長さをいいます。）は、 以下のような 式であらわされます。

       S = 365.24219878-6.14e-6 * T ［日］

                                     T は y を西暦年数とすれば、
                                     T=(y-1900)/100 で与えられる。

この式は、ニューカムと言う天文学者が 19世紀末に発表した太陽の運動理論から 導かれたものです。この式によると、1太陽年の日数が年々短くなる傾向がありま すが、（その主な原因は、 潮汐摩擦によって地球の自転角速度が減少することに より１日の長さがのびるためです。）ここでは定数項の部分のみを注目し 1太陽 年の日数は一定として考えると、1年の長さは

       365.24219878 ［日］

であると言えます。 この値と、 ユリウス暦法による１年の長さを比較してみま す。ユリウス暦法では、4年に1度だけ閏日を入れるので以下のようになります。

       ( 365*4 + 1 )/4 = 1461 / 4
                       = 365.25 ［日］

一方、グレゴリオ暦法では、100で割り切れる年については、 400で割り切れる年 のみを閏年としますので、400年間で閏年は 97回となります。と言うのは、単に4 年に1度とすると、 400 年間で 100回となりますが、100で割り切れる年で、400 で割り切れない年が 400年間で 3回あるため3日少なくなると言う訳です。 したがいまして、グレゴリオ暦法による１年の長さは、以下の通りとなります。

       ( 365*400 + 100 - 3 )/400 = 146097 / 400
                                 = 365.2425 ［日］

以上の事から、 グレゴリオ暦法の方がユリウス暦法よりも真の１年の長さに近 いと言えます。 Y年M月D日（グレゴリオ暦法による日付）のユリウス日を求めるには、以下の式を 使用します。

       JD=INT(365.25*Y)+INT(Y/400)-INT(Y/100)+INT(30.59*(M-2))+D+1721088

       ※ M<3 の場合には、Mに12を加え、Yから1を引いてから計算します。
        また、INT()は () の中の整数部のみを取り出す関数です。

ユリウス日で、 日付と時刻を同時に示すため時刻を「日」の単位の小数で示す 場合があります。この場合には、 協定世界時（協定世界時の詳しい内容について は、 あとで説明します。）よりも１２時間の時差を持つ時刻（天文時）を使用す る事になっています。 このため、 ユリウス日の１日の始まりが、 （協定世界時 の）０時ではなく、１２時となっています。

例題として、 1994年5月1日3時（日本標準時）の場合について計算いたします と、以下のようになります。

       1994年5月1日 ＝ 2449473

       JD = 2449473 + (3/24) - (9/24) + (12/24)
          = 2449473 + 0.125  - 0.375  +    0.5
          = 2449473.25

なお本スクリプトでは、 天文時を使用していません。（私がデバックの目的で検 算した時に、１２時間の時差がある事をしばしば忘れてしまい、 答えを間違えて しまう事があったためと言う個人的な事情からと一般常識に合わせて、 夜の０時 に日替わりするということからそうしました。）これから示す数々の例題は、 ユ リウス日の１日の始まりが、 （協定世界時又は、力学時の）０時ですので注意し て下さい。

ユリウス日の簡単な応用として任意の日付（グレゴリオ暦法による日付）に対 する曜日を計算するというのがあります。 以下の式によって曜日番号を求めて、対応する曜日を得ます。

       曜日番号 = ( JD+2 ) % 7

例えば、1994年5月1日では、以下のように計算します。

       1994年5月1日 ＝ 2449473

       曜日番号 = ( 2449473 + 2 ) % 7
                = 0 ..... すなわち日曜日

といった具合です...。日の干支も、同様に計算できます。これは、本スクリプト の説明の範囲を越えてしまうので、皆さんの練習問題としましょう。 （といって 逃げる...）

ユリウス暦の名の由来はローマ皇帝の名Julius Caearから、 ユリウス日の方は 考案者Jeseph Justus Scaligerの父の名Julius Caesar Scaligerから由来してい ると言われ、名前の由来に関しては両者は無関係です。

8.7 日本標準時（ＪＳＴ）

放送やＮＴＴ等の時報で一般的に使用さている時刻です。 世界各地の時刻は、 協定世界時（ＵＴＣ）と呼ばれる時刻の整数の時差（地域によっては３０分の単 位の端数がつく事がある）をもつように決められています。 また、夏時間制（サ マータイム）を採用している地域では、 夏時間と冬時間で時差を変えて実施して います。 このような時刻を地方標準時と呼んでいます。日本標準時も地方標準時 の１つで、協定世界時との時差は＋９時間です。すなわち、

       ＪＳＴ＝ＵＴＣ＋９  ［時］

の関係があります。

余談ですが、理科年表などでは日本標準時の事を中央標準時と称しています。 そ の理由としては、日本が台湾を占領していた頃に出された勅令により...

• i 東経135度の子午線の時刻を中央標準時と定める
• ii 東経120度の子午線の時刻を使用し、 これを西部標準時と定める（台湾や 宮古列島などで使用されていた。）

この勅令は、 昭和12年に別の勅令で西部標準時を削除される形で改正されました が、 『中央標準時』と言う名称はそのまま残った形で現在に至っています。 な お、現在でもこの時の勅令は生きていて一応、法的な根拠があるとのことです。

8.8 協定世界時（ＵＴＣ＝Coordinated Universal Time）

定義：国際単位系（ＳＩ単位系）で定められた秒（セシュウム原子時計による １秒で原子秒と言います。）を刻み、 地球自転に基づく世界時１との差が一定の 範囲以内になるように管理された時刻。

定義にもあるように、 世界時１（ＵＴ１）と協定世界時との差（ΔＵＴ１） が、±0.9秒以内となるように管理され、これを越える場合には、１秒単位で時刻 の調整を行います。（これをステップ調整と言う。） なお、±0.9秒と言うのは天体観測で経度を測定する場合（船舶・航空機の位置を 測定する場合など）に、 協定世界時を世界時とみなしても必要な精度を維持する ことができる限界として決めらた値です。 また、 無線報時（日本の場合はＪＪ Ｙ）ではΔＵＴ１の予測値を ±0.1sec の精度で信号化されていますので結局、 ±0.1sec の精度でＵＴ１とリンクしていると言えます。

１秒のステップ調整は、12月か6月の末日（第一優先）、3月か9月の末日（第二 優先）、必要とあれば任意の月の末日の最終秒（ＵＴＣ）の後に挿入するか、 ま たは最終秒を引き抜くことによって行われます。 この際挿入される１秒を「正の 閏秒［うるうびょう］」と言い、 引き抜かれる１秒を「負の閏秒」と言います。 閏秒の前後の秒の歩みは、以下のようになります...。

       正の閏秒の場合
       23:59:58   23:59:59   23:59:60   00:00:00   00:00:01 ...
                             ~~~~~~~~

       負の閏秒の場合
       23:59:57   23:59:58   00:00:00   00:00:01   00:00:02 ...
                           ^ 23:59:59が引き抜かれます。

8.9 世界時（ＵＴ＝Universal Time）

定義：赤道上を一定の速度で運動している仮想の天体（これを平均太陽と言 う）と、 経度０度における時角（経度０度の子午線と平均太陽の角度）に１２時 間を加えた値。

平均太陽は仮想天体なので、その南中を観測する事は出来ないので、 恒星の南 中を観測して、計算によって求められています。 各国の天文台が恒星の南中や、 天頂通過を観測して求められた世界時を世界時０ （ＵＴ０）と呼んでいます。この値は、 極運動[*8]による観測地点の経度変化が 含まれているため、 （観測誤差は別として同じ結果が求められるはずですが）観 測地によってまちまちの結果となります。 極運動の補正を行って観測地に依存し ない一義的な世界時を、世界時１（ＵＴ１）と言います。 世界時１は、 瞬間自転軸に関する自転角度を基準として時刻を測っていると考え る事ができます。従って、地球の自転の変化があると、 世界時１も同じように影 響を受けます。
以前は、 世界時１に季節的変化を補正した世界時を（世界時２＝ＵＴ２）が一様 な時刻と考えられていましたが、時計の精度が向上すると、 地球の自転速度は一 定ではない事が解ったため、 時刻の基準として地球の自転を使用する事をやめま した。

変動がある時刻では、天体暦の時刻引数として使用できません。 天体暦の時刻 引数は、 一様である事が前提となっているためです。たとえば、ある時刻に月が この位置にあると言う事が、計算によって求められていても、 肝心の時刻が変化 してしまっては、正確に位置が指定できないのと同じ事になるためです。 このため、 天体暦の時刻引数として使用する時刻は、後述する力学時（ＴＤ） を使用します。
天体観測で経度を測定する場合[*9]や日常生活には、 一様である時刻よりは地球 の自転角度に忠実な時刻を使用する方が重要です。 このため、実用的には国際単 位系（ＳＩ単位系）で定められた秒を刻み、 地球自転に基づく世界時（ＵＴ１） との差が一定範囲以内になるように管理された時刻（協定世界時）を使用する事 になりました。

なお、 本スクリプトでは世界時１と協定世界時は等しいとみなして計算してい ます。前に説明しましたように厳密には両者は全く異なるのですが、 本スクリプ トで要求される精度では、 両者を等しいとみなしても特に問題にならない程度の 差のためです。もう一つは、ΔＵＴ１が予測不可能な値であるためです。 （従っ て、計算によって補正する事はできません。）

昔は、 世界時とは言わずに『グリニジ平均太陽時』（ＧＭＴ＝Greewich Meen solar Time）と呼んでいました。今でも航海・航空関係の方々は、 習慣でこう呼 んでいます。 国際天文学連合で1976年にだされた決議では、 「科学的用法としては、ＧＭＴの 呼称をやめ、 詳しく言う時にはＵＴ０・ＵＴ１・ＵＴ２・ＵＴＣを公式に用い区 別の必要が無い場合には単にＵＴと呼ぶ」となっています。 また、 「ＵＴは、法的・通信・常用その他の用法で、最大精度が整数秒である場 合にはＵＴＣの意味で用いられる事に注意し、 天文航法・陸地測量のための天体 暦の時刻引数としては、ＵＴ１の意味で使用する」となっています。

8.10 力学時（ＴＤ＝Dynamical Time）

天体力学理論や暦の時間引数として、 天体の位置を計算するために用いる時刻 です。本スクリプトでは、 太陽の黄経と月の黄経を計算する時に使用する時刻引 数として使用しています。 厳密には、 力学時にも太陽系重心を基準とした太陽系力学時（ＴＤＢ）と、地球 の重心を基準とした地球力学時（ＴＤＴ）の２種がありますが、 拙者は詳しい事 は解りませんし、 旧暦計算には関係無いと思われますので、説明は割愛いたしま す。

おおざっぱに言って協定世界時との関係は、 世界時＝協定世界時（ΔＵＴ１＝ ０）と仮定しますと、

       ＴＤ＝ＵＴＣ＋ΔＴ

となっています。 ΔＴは、両者の時刻の差で（1993年12月現在で）約６０［秒］ となっています。この値は、正の閏秒が入ると１秒増やし、 負の閏秒が入ると１ 秒減らす事になります。 前にも説明したとおり閏秒は地球の自転の状態によって 決るため、現在以降の値については、その時になってみないと解りません。 （今 の所は、だいたい毎年１秒づつ正の閏秒が挿入されている。） また、 協定世界時の１秒間と力学時の１秒間は全く等しいです。と言うのは、両 者とも、原子秒を採用しているからです。

本スクリプトでは、ΔＴのような不確定な値は無視（0秒として扱う）していま す。したがいまして、 ΔＴの値だけ実際の時刻とは異なった計算結果をはじきだ している訳です。 （本来は、しかるべき補正式を使用して時刻を補正するべきで しょうが、任意の時刻に対する補正値を計算できる式が見当たりませんので、 あ えて補正しないことにしました。）

この事が、 本スクリプトの計算の結果にどのように影響するかと言いますと、実 際の朔や中気の日付と異なる可能性があります。その結果、 月の大小や月名（閏 月の位置が変わってしまった場合）を誤る可能性があります。 もっとも、このよ うなケースは希ではありますが。（計算結果が、 たまたま日付の変わる時刻付近 だった場合には、 たとえ１秒の差であっても日付が違ってくる事を考えてくださ い。）

参考までに1620～1992年の期間で、5年毎（1990.0～1992.0までは毎年）のΔＴ を以下に載せておきます。出典は、 1620～1975年の期間における値は天体位置略 算式の解説 P17 から、それ以降の値は理科年表平成６年版の天83(167) からで す。なお、理科年表の値は、 年央における値が掲載されていましたので補間を行 い、年初における値を計算しました。

9. 旧暦の計算の概略

旧暦の計算の手順は、 日付を先に決めてから、月名を決めるのが本来の方法で す。具体的には、［計算対象の直前の朔からの日数＋１］を旧暦の日付とし、 月 名は、二分二至・中気と朔の関係から決定しています。 本スクリプトでは、 二分二至・中気・朔を一括して求めてから、新暦と旧暦の対 応表『朔日行列』（拙者が勝手に名前を付けた）と言う２次元配列を作成してか ら、月日を求めています。

       配列の内訳

       m[i,0] ... 旧暦の月名（1:正月 2:２月 3:３月 ....）
       m[i,1] ... 閏フラグ（1:閏月 0:閏月ではない）
       m[i,2] ... 朔日のユリウス日

本スクリプトでは、以下のような方法で行っています。

• i 計算する日付の直前の二分二至の日付を求める。
• ii 計算する日付の直前の二分二至と、中気の時刻と黄経を３回分求める。
• iii 計算する日付の直前の二分二至以前にあった朔の日付を求めます。引き続 いて、その後の朔の日付を４回分計算します。（合計で、５回だけ朔の時刻 を計算する事になります。）
• iv 中気の時刻と朔の時刻から旧暦の月名と、旧暦の朔日［ついたち］が新暦 の何月何日になるかを計算して、新暦と旧暦の対応表『朔日行列』を作成し ます。さらに、閏月の有無の判別をここで行います。閏月は、以下の条件が 同時に成り立つ時とします。（ＡＮＤ条件）

            i   ［５回目の朔の時刻］≦［３回目目の中気の時刻］
  
           ii   中気を含まない月
         

• v 計算する日付の直前の朔を（朔日行列の中から）見つけて、旧暦の日付を 求める。（朔の日付と、計算する日付が等しい場合には、朔日とする。）

i ～ iii の手順は解りづらいと思いますので、以下の図を参照して下さい。

       i と ii の手順の流れ

       ------------------ 直前の二分二至
        ^     |
        |     |
        |i    |
        |     |
        |     |
        |     v ii-1
        |  -------------- 中気１
        |     |
        |     |
        |     |
       ------ | --------- 計算対象の時刻
              |
              v ii-2
           -------------- 中気２
              |
              |
              |
              |
              |
              v ii-3
       ------------------ 中気３（直後の二分二至）



       iii の手順の流れ

       ------------------ 朔１
        ^iii-1   |
        |        |
        |        |
       --------- | ------ 直前の二分二至
                 |
                 v iii-2
       ------------------ 朔２
                 |
                 |
                 |
       --------- | ----- 中気１
                 |
                 v iii-3
          --------------- 朔３
                 |
                 |
                 |
       --------- | ------ 中気２
                 |
                 v iii-4
          --------------- 朔４
                 |
                 |
                 |
       --------- | ------ 中気３
                 |
                 v iii-5
          --------------- 朔５

この図はあくまでも説明のために作成したものですから、 必ずこのように朔と 朔の間に中気が入るとは限りません。（例題２，３を参照）

計算例をみてどんな感じかをつかんでもらった方が早いかと思います...。（と 言って逃げる...。）

  例題１
    1994年5月1日は旧暦の何月何日になるか求めなさい。

  【解】

  ・計算する日付の直前の二分二至の日付を求めます。 ５月の直前の二分二至は春
  分ですのでこの時刻を計算します。

  答え：1994年 3月 21日（ユリウス日＝2449432）

  ※ 直前の二分二至の時刻の計算方法については、後で説明します。

  ・計算する日付の直前の二分二至から、中気の時刻と黄経を３回分求めます。

  答え：

   1. 1994年 4月20日（ユリウス日＝2449462）

   2. 1994年 5月21日（ユリウス日＝2449493）

   3. 1994年 6月21日（ユリウス日＝2449524）


  ※ 中気の時刻の計算方法については、後で説明します。

  ・計算する日付の直前の二分二至を与えて、 この日付以前にあった朔の日付を求
  める。 引き続いて、その後の朔の日付を４回分計算します。（合計で、５回だけ
  朔の時刻を計算する事になります。）

  答え：

   1. 1994年 3月 12日（ユリウス日＝2449423）

   2. 1994年 4月 11日（ユリウス日＝2449453）

   3. 1994年 5月 11日（ユリウス日＝2449483）

   4. 1994年 6月  9日（ユリウス日＝2449512）

   5. 1994年 7月  9日（ユリウス日＝2449542）


  ※ 朔の時刻の計算方法については、後で説明します。

  ・中気の時刻と朔の時刻から旧暦の月名と、 旧暦の朔日［ついたち］が新暦の何
  月何日になるかを計算して、 新暦と旧暦の対応表を作成する。さらに、閏月の有
  無の判別をここで行います。条件は、

            ［５回目の朔の時刻］≦［３回目目の中気の時刻］

  の場合でかつ、中気を含まない月とする。

  答え：
    春分を含む月は旧２月と決まっていますので、1994年 3月 12日から、4月 10日
  までの間が旧暦の２月となり、 後は自動的に３月・４月・５月と決まって行きま
  す。 これをまとめて表にしますと以下のようになります。  日付のとなりにある
  （）は、ユリウス日です。

  +--------------------------------------------+
  |旧暦の月名|   朔日の日付   |   中気の日付   |
  |----------+----------------+----------------|
  |    ２    | 3/12 (2449423) |                |
  |----------+----------------+----------------|
  |          |                | 3/21 (2449432) |
  |----------+----------------+----------------|
  |    ３    | 4/11 (2449453) |                |
  |----------+----------------+----------------|
  |          |                | 4/20 (2449462) |
  |----------+----------------+----------------|
  |    ４    | 5/11 (2449483) |                |
  |----------+----------------+----------------|
  |          |                | 5/21 (2449493) |
  |----------+----------------+----------------|
  |    ５    | 6/ 9 (2449512) |                |
  |----------+----------------+----------------|
  |          |                | 6/21 (2449524) |
  |----------+----------------+----------------|
  |    ６    | 7/ 9 (2449542) |                |
  +--------------------------------------------+

    閏月の判別の条件として［５回目の朔の時刻］≦［３回目目の中気の時刻］で
  かつ、中気を含まない月としています。この問題では、

            ［５回目の朔の時刻］........ 7/ 9 (2449542)
            ［３回目目の中気の時刻］.... 6/21 (2449524)

  となって、条件が成り立ちませんので、閏月はありません。

  ・計算する日付の直前の朔を見つけて、旧暦の日付を求める。

  答え：
    問題では、 5月 1日が、旧暦の何月何日かを聞いているのでした。従って、5月
  1日の直前の朔は、4月11日となります。前に作成した表から、 月名は３月となり
  ます。日付は、［4月11日からの日数］＋１と計算します。従いまして、以下のよ
  うになります。

            日付       ユリウス日
            4/11 ..... 2449453
            5/ 1 ..... 2449473

            4/11 - 5/ 1 + 1
            = 2449473 - 2449453 +1
            = 21


    以上の事から、西暦1994年 5月 1日は、旧暦1994年 3月 21日と言う事になりま
  す。六曜は次のように計算して決めます。

            ［旧暦の月］＋［旧暦の日］－２ を６で割った剰余
            = (3 + 21 - 2 ) % 6
            = 4

  答えが 4 の場合、大安となります。と、まあこんな感じです。

  例題２
    1993年5月1日は旧暦の何月何日になるか求めなさい。

  【解】

    途中の過程を省略しまして、 いきなり旧暦の月名／朔日の日付／中気の日付の
  表を作成しますと、

  +--------------------------------------------+
  |旧暦の月名|   朔日の日付   |   中気の日付   |
  |----------+----------------+----------------|
  |    ２    | 2/21 (2449039) |                |
  |----------+----------------+----------------|
  |          |                | 3/20 (2449066) |
  |----------+----------------+----------------|
  |    ３    | 3/23 (2449069) |                |
  |----------+----------------+----------------|
  |          |                | 4/20 (2449097) |
  |----------+----------------+----------------|
  |   閏３   | 4/22 (2449099) |                |
  |----------+----------------+----------------|
  |          |                |  （中気無し）  |
  |----------+----------------+----------------|
  |    ４    | 5/21 (2449128) | 5/21 (2449128) |
  |----------+----------------+----------------|
  |          |                |                |
  |----------+----------------+----------------|
  |    ５    | 6/20 (2449158) |                |
  |----------+----------------+----------------|
  |          |                | 6/21 (2449159) |
  +--------------------------------------------+

    閏月の判別の条件として［５回目の朔の時刻］≦［３回目目の中気の時刻］で
  かつ、中気を含まない月としています。この問題では、

            ［５回目の朔の時刻］........ 6/20 (2449158)
            ［３回目目の中気の時刻］.... 6/21 (2449159)

  となっていますので、条件が成り立ちます。また、 4月22日から5月20日の間は、
  中気がありませんので、 この月は閏月になります。月名は表を見ても分かるとお
  り前月が３月でしたので、閏３月となります。閏月がある場合には、 このような
  表の形態になります。
    さて、日付の方は

            日付       ユリウス日
            4/22 ..... 2449099
            5/ 1 ..... 2449108

            4/22 - 5/ 1 + 1
            = 2449099 - 2449108 +1
            = 10

  となって、 西暦1993年  5月  1日は、旧暦1993年 閏 3月 10日と言う事になりま
  す。六曜は次のように計算して決めます。

            ［旧暦の月］＋［旧暦の日］－２ を６で割った剰余
            = (3 + 10 - 2 ) % 6
            = 5

  答えが 5 の場合、赤口となります。このように閏月の場合は、前月と同様に計算
  すれば良い訳です。

  例題３
    1985年1月1日は旧暦の何月何日になるか求めなさい。
  【解】

    途中の過程を省略しまして、 いきなり旧暦の月名／朔日の日付／中気の日付の
  表を作成しますと、

  +----------------------------------------------------+
  |旧暦の月名|     朔日の日付     |     中気の日付     |
  |----------+--------------------+--------------------|
  |  1984/11 |1984/12/22 (2446056)|1984/12/22 (2446056)|
  |----------+--------------------+--------------------|
  |          |                    |1985/ 1/20 (2446085)|
  |----------+--------------------+--------------------|
  |  1984/12 |1985/ 1/21 (2446086)|                    |
  |----------+--------------------+--------------------|
  |          |                    |1985/ 2/19 (2446115)|
  |----------+--------------------+--------------------|
  | 1985/正月|1985/ 2/20 (2446116)|                    |
  |----------+--------------------+--------------------|
  |          |                    |    （中気無し）    |
  |----------+--------------------+--------------------|
  |  1985/ 2 |1985/ 3/21 (2446145)|1985/ 3/21 (2446145)|
  |----------+--------------------+--------------------|
  |          |                    |                    |
  |----------+--------------------+--------------------|
  |  1985/ 3 |1985/ 4/20 (2446175)|                    |
  +----------------------------------------------------+


    閏月の判別の条件として［５回目の朔の時刻］≦［３回目目の中気の時刻］で
  かつ、中気を含まない月としています。この問題では、

            ［５回目の朔の時刻］........ 4/20 (2446175)
            ［３回目目の中気の時刻］.... 3/21 (2446145)

  となっていますので、条件が成り立ちません。したがって、2/20から 3/20までの
  間に中気がありませんが、この月は閏月となりません。また、1984/12/22 から、
  1985/ 1/20 の間に中気が２回もありますが、これは以下の２つの条件が重なった
  ことによって説明できます。

     i   地球が近日点付近にある時期で、  太陽の黄経の変化率が大きくなってお
       り、結果として中気の間隔が詰まった。

    ii   月の黄経の変化率が小さい時期にあった。

    このような場合に確実に閏月かどうかを判別するのは、 案外と厄介です。（実
  は、本スクリプト作成時にもこの部分に最後まで泣かされた）もっとも、 このよ
  うな月まわりになる事はそうありませんが...。

    さて、日付の方ですが、以下のように計算して求めます。

              日付        ユリウス日
            1984/12/22 .... 2446056
            1985/ 1/ 1 .... 2446066


            1985/ 1/ 1 - 1984/12/22 + 1
            = 2446056 - 2446066 + 1
            = 11

  従って、 西暦1985年 1月 1日は、旧暦1984年 11月 11日と言う事になります。本
  問のように、 太陽暦で年が改まっても旧暦ではまだ前年となる期間がありますの
  で、注意が必要になります。六曜は次のように計算して決めます。

            ［旧暦の月］＋［旧暦の日］－２ を６で割った剰余
            = (11 + 11 - 2 ) % 6
            = 2

  答えが 2 の場合、先負となります。

  例題３（1985年  1月 1日の旧暦を求める）の場合には朔日行列の各要素は以下の
  ような値になります。

       m[0,0]=11  m[0,1]=0  m[0,2]=2446056
       m[1,0]=12  m[1,1]=0  m[1,2]=2446086
       m[2,0]= 1  m[2,1]=0  m[2,2]=2446116
       m[3,0]= 2  m[3,1]=0  m[3,2]=2446145
       m[4,0]= 3  m[4,1]=0  m[4,2]=2446175

10. 計算する日付の直前の二分二至の日付の求め方

旧暦の月名／朔日の日付／中気の日付の表を作成する時に、 計算する日付の直 前にあった二分二至の日付の求めますが、この方法について説明します。

• i 計算する日付における太陽黄経を計算し、90 で割って商を得る。
• ii 先に計算した商に 90 を掛けて、基準となる黄経（λsun0）を計算する。 （この値が、二分二至の黄経となる。）
• iii 仮の二分二至の時刻（ｔ）における黄経（λsun）と、 基準となる黄経 （λsun0）との差（λsun-λsun0）を計算し、 λsunの時刻に対する変化係 数（365.2/360）を掛けて、時刻引数の補正値（Δｔ）を計算する。

              Δｔ＝（λsun－λsun0）＊（365.2/360）
  

• iv 仮の二分二至の時刻を以下のように補正する。

              ｔ＝ｔ－Δｔ
  

• v 仮の二分二至の時刻の補正を行い再び太陽黄経を計算する。この値から再 び、Δｔを計算する...。 これをΔｔが１秒以下に収束するまで反復計算す る。

このようにして計算しています。 Δｔの収束が悪いと、反復する回数が多くな って、計算時間が長くなってしまいます。 例題を挙げて実際に計算してみました ので、どんな感じかつかんで下さい。

  例題
    1994年5月1日の直前の二分二至の時刻を求めなさい。

  【解】

  ・計算する日付における太陽黄経を計算し、90 で割って商を得る。

       まず、ユリウス日に変換する。時間は「日」の小数として表す。

       1994年5月1日0時 = 2449473.0 ユリウス日

       日本標準時の値であるから、力学時に変換する。ただし、 力学時－協
       定世界時は９時間とする。

       TD = JST - (9/24)
          = 2449473.0 - 0.375
          = 2449472.625

       この時の、 太陽黄経λsunは、40.034279度なので（太陽黄経の計算方
       法は、後述する。）、90で割った商は、0 となる。


  ・先に計算した商に 90 を掛けて、基準となる黄経（λsun0 ）を計算する。（こ
  の値が、二分二至の黄経となる。）

       λsun0 = 0 * 90
              = 0°（ちなみに太陽黄経 0° は、春分である。）


  ・仮の二分二至の時刻（ｔ）における黄経（λsun）と、  λsun0との（λsun-λ
  sun0）を計算し、  λ sun の時刻に対する変化係数（ 365.2/360）を掛けて、時
  刻引数の補正値（Δｔ）を計算する。

       Δｔ = （λsun－λsun0）＊（365.2/360）
            = （40.0342792282334200 - 0 ）＊（365.2/360）
            =  40.612552150419 ［日］


  ・仮の二分二至の時刻を以下のように補正する。

       ｔ = ｔ－Δｔ
          = 2449472.625 - 40.612552150419
          = 2449432.01244784958098677 ( 1994/ 3/21 0:17:55(TD) )


  ・仮の二分二至の時刻の補正を行い再び太陽黄経を計算する。 この値から再び、
  Δｔを計算する...。これをΔｔが１秒以下に収束するまで反復計算する。


       （１回目）

       ｔ = 2449432.01244784958098677
        におけるλsunは、以下の通り。

       λsun=0.1587307886278495

       この値から、Δｔを計算すると、

       Δt = （λsun - λsun0）＊（365.2/360.0）
           = ( 0.1587307886278495 - 0) * (365.2/360.0)
           = 0.1610235666858

       （Δtは約3時間51分なので、計算の打ち切り条件（1秒）よりも大きい
       ため続行する。）

       （２回目）

       ｔ = ｔ－Δｔ
          = 2449432.01244784958098677 - 0.1610235666858
          = 2449431.8514242828951795 ( 1994/ 3/20 20:26:3(TD) )

       λsun = 359.9987979993543000

       Δt = （λsun - λsun0）＊（365.2/360.0）
           = ( 359.9987979993543000 - 0 ) * ( 365.2/360 )
           = -0.0012020006456623 * ( 365.2/360 ) ※
           = -0.001219362877211

       （Δtは約105.3秒なので、計算の打ち切り条件（1秒）よりも大きいた
       め続行する。）

       ※ λsun-λsun0の結果が、180°よりも大きい場合には360から引いて
       補正する。


       （３回目）

       ｔ = ｔ－Δｔ
          = 2449431.8514242828951795 + 0.001219362877211
          = 2449431.8526436457723903 ( 1994/ 3/20  20:27:48(TD) )

       λsun = 0.0000091643490236

       Δt = （λsun - λsun0）＊（365.2/360.0）
           = ( 0.0000091643490236 - 0 ) + ( 365.2/360 )
           =  0.000009296722953991

       （約0.8秒なので、 計算の打ち切り条件（1秒）が成り立つため打ち切
       る。）

       t = ｔ－Δｔ
         = 2449431.8526436457723903 - 0.000009296722953991
         = 2449431.8526343490494362 ( 1994/ 3/20  20:27:47(TD) )

       この値は力学時(TD)による値であるから、  日本標準時(JST)に変換す
       る。ただし、力学時－協定世界時は９時間とする。

       JST = TD + (9/24)
           = 2449431.8526343490494362 + 0.375
           = 2449432.2276343490000000 <- 桁落ちによって精度が損なわ
           = 1994/ 3/21   5:27:47 (JST)  れた。

  旧暦を求める時には、 日付だけが問題となりますので時刻の部分は切り捨てしま
  います。 JSTに変換する際に桁落ちが生じていますが、 解の要求される精度（１
  秒）には影響しません。

    λsun0 = 0 として二分二至の時刻を求める場合（春分の時刻を求める場合）反
  復計算を行っている最中に、λsunが 0 と 360 の間を振動します。これは、太陽
  の黄経のとり得る値の範囲は 0 と 360 の間であるためですが、 ΔＴを計算する
  時には、補正値として扱うため -180 と 180 の間に修正しています。
  修正する方法は、θ≧180 の時に、 360 - θ としています。（例題を参照）

11. 中気の時刻の求め方

二分二至の時刻の次は、 中気の時刻を計算します。計算方法は、以下のように なっています。

• i 二分二至の時刻＋32日を初期値として、太陽黄経を求める。この値を30で 割って商を得る。（32を加えるのは、二分二至の直後の中気を求めるため）
• ii 先に計算した商に 30 を掛けて、基準となる黄経（λsun0）を計算する。 （この値が、中気の黄経となる。）
• iii 仮の中気の時刻（ｔ）における黄経（λsun）と、 λsun0との差（λsun- λsun0）を計算し、λsunの時刻に対する変化係数（365.2/360）を掛けて、 時刻引数の補正値（Δｔ）を計算する。

              Δｔ＝（λsun－λsun0）＊（365.2/360）
  

• iv 仮の中気の時刻を以下のように補正する。

              ｔ＝ｔ－Δｔ
  

• v 補正後のｔを使い太陽黄経を計算し、仮の中気の時刻の補正を行い再び太 陽黄経を計算する。これをΔｔが１秒以下に収束するまで反復計算する。
• vi このような計算を３回繰り返して、中気の時刻を３つ求める。ただし、２ 回目からは計算の初期値は前の中気の時刻＋32日とする。

初期値を計算する時に32を加える理由を説明しますと... 前に説明した方法は、 「初期値は真の値よりも先の時刻を与えなくてはならな い」という欠点があります。 このため初期値を与える時に、真の中気の時刻を確 実に飛び越す必要があります。そこで、 中気の間隔（約30日）よりもやや多めの 32日を加算することに決めた...。と言う次第です。

以下の図は初期値と真の値との関係です。

       --------------- 直前の二分二至
           |
           |
           |
           |  +32日
           |
           |
       --- | --------- 中気１
           v    |
       -------- | ---- 中気１を計算する時の初期値
                |
                |
                |  +32日
                |
                |
       -------- | ---- 中気２
           |    v
       --- | --------- 中気２を計算する時の初期値
           |
           |
           |  +32日
           |
           |
       --- | --------- 中気３
           v
       --------------- 中気３を計算する時の初期値

例題をやってみましたので、どんな感じか見て下さい。

  例題
    1994年5月1日の直前の二分二至の時刻を求め、中気の時刻を３回求めなさい。

  【解】

  先の例題の解から、 二分二至の時刻は、 1994年 3月21日 5時27分47秒(2449432.
  2276343490)だったのでこれを使い計算します。

  １回目の中気

  ・基準となる黄経を求める。

       t = t + 32
         = 2449432.2276343490 + 32
         = 2449464.2276343490 (JST)
       この時の太陽黄経は、

       λsun=31.4993009433906500 なので、基準となる黄経（λsun0）は

       λsun0=30.0


  以下、反復計算の状況を表にしました。

  +------------------------------------------------------------------+
  | N |       λsun       |        Δｔ         |          ｔ          |
  |---+------------------+--------------------+----------------------|
  | 1 | 31.49930094339065|      1.520957512573| 2449462.3316768363336|
  |---+------------------+--------------------+----------------------|
  | 2 | 30.01526745626125| 1.548798618503e-002| 2449462.3316768363336|
  |---+------------------+--------------------+----------------------|
  | 3 | 30.00014831624715| 1.504585929373e-004| 2449462.3160383915557|
  |---+------------------+--------------------+----------------------|
  | 4 | 30.00000144010887| 1.460910447454e-006| 2449462.3160369306452|
  +------------------------------------------------------------------+

       JST = TD + (9/24)
           = 2449462.3160369306452 + 0.375
           = 2449462.6910369310000
           = 1994/ 4/20  16:35:5



  ２回目の中気

  ・基準となる黄経を求める。

       t = t + 32
         = 2449462.6910369310 + 32
         = 2449494.6910369310 (JST)
       この時の太陽黄経は、

       λsun = 60.99331164862917 なので、基準となる黄経（λsun0）は

       λsun0 = 60.0


  以下、反復計算の状況を表にしました。

  +------------------------------------------------------------------+
  | N |       λsun       |        Δｔ         |          ｔ          |
  |---+------------------+--------------------+----------------------|
  | 1 | 60.99331164862917|   1.000765948355381| 2449493.3083774472196|
  |---+------------------+--------------------+----------------------|
  | 2 | 60.02436773598703| 2.471971439573e-002| 2449493.2836577328239|
  |---+------------------+--------------------+----------------------|
  | 3 | 60.00059237153853| 6.009280163095e-004| 2449493.2830568048076|
  |---+------------------+--------------------+----------------------|
  | 4 | 60.00001439707667| 1.460503444272e-005| 2449493.2830421997731|
  |---+------------------+--------------------+----------------------|
  | 5 | 60.00000034990633| 3.549605322285e-007| 2449493.2830418448126|
  +------------------------------------------------------------------+

       JST = TD + (9/24)
           = 2449493.2830418448126 + 0.375
           = 2449493.6580418450000
           = 1994/ 5/21  15:47:34



  ３回目の中気

  ・基準となる黄経を求める。

       t = t + 32
         = 2449493.6580418450000 + 32
         = 2449525.6580418450 (JST)
       この時の太陽黄経は、

       λsun = 90.63662714952599 なので、基準となる黄経（λsun0）は

       λsun0 = 90.0


  以下、反復計算の状況を表にしました。

  +------------------------------------------------------------------+
  | N |       λsun       |        Δｔ         |          ｔ          |
  |---+------------------+--------------------+----------------------|
  | 1 | 90.63662714952599|     0.6458228750191|  2449524.637218969591|
  |---+------------------+--------------------+----------------------|
  | 2 | 90.02061881486661| 2.091664219246e-002| 2449524.6163023273985|
  |---+------------------+--------------------+----------------------|
  | 3 | 90.00066674989583| 6.763807276566e-004| 2449524.6156259466709|
  |---+------------------+--------------------+----------------------|
  | 4 | 90.00002155953537| 2.187095088031e-005| 2449524.6156040757200|
  |---+------------------+--------------------+----------------------|
  | 5 | 90.00000069713205| 7.072017368317e-007| 2449524.6156033685182|
  +------------------------------------------------------------------+

       JST = TD + (9/24)
           = 2449524.6156033685182 + 0.375
           = 2449524.9906033690000
           = 1994/ 6/21  23:46:28


  結果をまとめると...

  +-------------------------------------------------------+
  | N |        ユリウス日       |        日付／時刻       |
  |---+-------------------------+-------------------------|
  | 1 |  2449462.6910369310000  |   1994/ 4/20  16:35: 5  |
  |---+-------------------------+-------------------------|
  | 2 |  2449493.6580418450000  |   1994/ 5/21  15:47:34  |
  |---+-------------------------+-------------------------|
  | 3 |  2449524.9906033690000  |   1994/ 6/21  23:46:28  |
  +-------------------------------------------------------+

といった具合です...。

12. 朔の時刻の求め方

中気の次は、 朔の時刻を計算します。 計算方法は、 以下のようになっていま す。

• i 二分二至の時刻を初期値（仮の朔の時刻）とし、月（λmoon）と太陽（λ sun）の黄経差（Δλ）を求める。ただし、 Δλ＜０の場合には正規化（値 を、0～360の範囲に制限する事）する。

              Δλ＝λmoon－λsun
  

• ii 黄経差（Δλ）の時刻に対する変化係数（29.530589/360）を掛けて、 時 刻引数の補正値（Δｔ）を計算する。

              Δｔ＝Δλ＊（29.530589/360）
  

• iii 仮の朔の時刻を以下のように補正する。

              ｔ＝ｔ－Δｔ
  

• iv 補正後のｔを使い月と太陽の黄経を計算し、Δλを計算する。ただし、春 分の近くで朔がある場合（0°≦λsun≦ 20°）で、月の黄経λmoon≧300° の場合には以下のように補正を行ってΔλとする。

              Δλ＝ 360.0 － Δλ
  

  また、Δλが引き込み範囲（±40°）を逸脱した場合には、正規化を行い Δλとする。
• v 黄経差（Δλ）の時刻に対する変化係数（29.530589/360）を掛けて、 時 刻引数の補正値（Δｔ）を計算する。

              Δｔ＝Δλ＊（29.530589/360）
  

• vi 仮の中気の時刻を以下のように補正する。

              ｔ＝ｔ－Δｔ
  

• vii 仮の朔の時刻の補正を行い再びΔλを計算する。この値から再び、Δｔを 計算する...。これをΔｔが１秒以下に収束なるまで反復計算する。 （ただ し、反復回数が15回を越える場合には、解が振動していて収束しないと考え て、i の初期値から26を引き去り計算をやり直す。それでも依然として収束 しない（反復回数30回以上）場合には、初期値を解として異常終了する。）
• viii このような計算を５回繰り返して、朔の時刻を５つ求める。ただし、２回 目からは計算の初期値を前の朔の時刻＋30日とする。

正規化処理について補足説明をしておきますと、 黄経差は、0～360の範囲で値 をとり得ますが、計算の途中では負の値になったり、360の整数倍の値が加算され ていることがあります。このような値を、0～360の範囲にする処理を言います。 例えば...

        -335°－－－＞ 360°* -1 + 25°

        1825°－－－＞ 360°*  5 + 25°

       -1775°－－－＞ 360°* -5 + 25°

となってこれらの値は、正規化後はいずれも25°となります。

反復計算の初回と、初回以降のΔλのあつかいが違う理由は...

• i 初回は、Δλの値はどのような値になるか見当もつかないので、このまま Δｔを計算すると真の朔の時刻から離れてしまう可能性があるため、Δλを 正規化して方向性を持たせます。（強制的に初期値よりも前の時間にさかの ぼる方向へもっていく。）
• ii 初回以降はΔλの絶対値は小さくなっているはずですから、Δλには方向 性を持たせません。 しかし、反復計算の最中にΔλの絶対値が小さくならない場合（Δλの範囲 が±40°を越える場合）には、強制的に収束させるため、正規化を行いひと つ前の朔の時刻を計算させる。）

最初の補正値を計算する時に、正規化をしないとΔλに方向性が無くなり、 解 が収束せず、動作が不安定になる事があります。 初期値の与え方によっては正規化を行っても目的にあった解を得る事ができず、 朔日行列が正常に組み立てられない場合があります。 本スクリプトでは、以下の ２つのケースを想定して対策を施しておきました。

     i   二分二至の前の朔の時刻を２組求めてしまう場合。以下の図ように、

            ［朔２の時刻］≦［直前の二分二至の時刻］

       の条件が成立する場合。（正しい解は、朔２の時刻が朔１となる。）

            ------------------ 朔１
             ^        |
             |        |
             |        |
             | ------ | ------ 中気（計算対象外の中気）
             |        |
             |        v
             | --------------- 朔２  <----- 飛び越してしまう
             |        |
             |        |
             |        |
            --------- | -----  直前の二分二至
                      |
                      v
               --------------- 朔３
                      |
                      |
                      |
            --------- | ------ 中気１
                      |
                      v
               --------------- 朔４
                      |
                      |
                      |
            --------- | ------ 中気２
                      |
                      v
               --------------- 朔５

    ii   二分二至の後の朔の時刻を求めてしまう場合。以下の図ように、

            ［朔１の時刻］≧［直前の二分二至の時刻］

       の条件が成立する場合。（正しい解は、朔０の時刻が朔１となる。）

               --------------- 朔０
                      |
                      |
                      |
                      |
            --------- | -----  直前の二分二至
             |        |
             v        v
            ------------------ 朔１
                      |
                      |
                      |
            --------- | ------ 中気１
                      |
                      v
               --------------- 朔２
                      |
                      |
                      |
            --------- | ------ 中気２
                      |
                      v
               --------------- 朔３


    第一点のケースは、任意の時刻に対する直前の二分二至の時刻と、 朔１が非常
  に接近している場合（１日以内）に発生する可能性があります。 このような場合
  に対しては、

       ［朔２の時刻］≦［直前の二分二至の時刻］

  の条件で検出を行い、成立した場合には、

       朔１ ← 朔２
       朔２ ← 朔３
       朔３ ← 朔４
       朔４ ← 朔５

  のようにして、朔の番号付けを繰り下げて修正します。朔５については、 （修正
  後の）朔４＋35 を初期値として再計算して補います。

    なお、 本スクリプトを作成する際にした試験では、1984年12月22日から、1985
  年3月21日の間にこの条件が成り立ちました。

    第二点のケースは、任意の時刻に対する直前の二分二至の時刻と、 朔１が非常
  に接近している場合（１日以内）に発生する可能性があります。 このような場合
  に対しては、

       ［朔１の時刻］≧［直前の二分二至の時刻］

  の条件で検出を行い、成立した場合には、

       朔５ ← 朔４
       朔４ ← 朔３
       朔３ ← 朔２
       朔２ ← 朔１

  のようにして、朔の番号付けを繰り上げて修正します。朔１については、 （修正
  後の）朔１－27 を初期値として再計算して補います。

例題をやって見ましたので、ご覧ください。

  例題
    1994年5月1日の直前の二分二至の時刻を初期値として朔の時刻を５回求めなさ
  い。

  【解】

  二分二至の時刻は、 1994年 3月21日 5時27分47秒 (2449066.22763434943)だった
  のでこれを使い計算します。

  １回目の朔

       日本標準時から力学時に変換する

       TD = JST - (9/24)
          = 2449432.22763434904943 - 0.375
          = 2449431.85263434904943

       この時の月と太陽の黄経は、

       λmoon = 359.9999999299906
       λsun  = 93.93361186916204

       なので、黄経差Δλは、

       Δλ = λmoon - λsun
            = -266.0663880608286

       Δλが負の値になったので正規化すると、

       Δλ = 93.93361193917144

       Δｔ = Δλ＊（29.530589/360）
            = 7.7053191318366

       ｔ = ｔ - Δｔ
          = 2449431.85263434943 - 7.7053191318366
          = 2449424.14731521707


  以下、反復計算の状況を表にしました。（λmoon・λsun・Δλ）
+--------------------------------------------------------------+
| N |      λmoon      |       λsun      |        Δλ        |
|---+------------------+------------------+--------------------|
| 1 | 1.776753401978283| 352.3279362096001|  9.4488171923781580|
|---+------------------+------------------+--------------------|
| 2 |  352.408282755941|  351.554257533407|   0.854025222534289|
|---+------------------+------------------+--------------------|
| 3 |  351.557149475943|  351.484313049735|   0.072836426208198|
|---+------------------+------------------+--------------------|
| 4 |  351.484525107571|  351.478347640230|   0.006177467341387|
|---+------------------+------------------+--------------------|
| 5 |  351.478365374802|  351.477841695768|   0.000523679033619|
|---+------------------+------------------+--------------------|
| 6 |  351.477843197370|  351.477798805611|   0.000044391759047|
+--------------------------------------------------------------+

（Δλ・Δｔ・ｔ）
+--------------------------------------------------------------------+
| N |        Δλ        |        Δｔ        |          ｔ          |
|---+--------------------+--------------------+----------------------|
| 1 |  9.4488171923781580|      0.775080936234| 2449423.3722342808360|
|---+--------------------+--------------------+----------------------|
| 2 |  0.8540252225342897| 7.005518845082e-002| 2449423.3021790923852|
|---+--------------------+--------------------+----------------------|
| 3 |   0.072836426208198|  5.97472935162e-003| 2449423.2962043630335|
|---+--------------------+--------------------+----------------------|
| 4 |   0.006177467341387| 5.067340253318e-004| 2449423.2956976290082|
|---+--------------------+--------------------+----------------------|
| 5 |   0.000523679033619| 4.295708419373e-005| 2449423.2956546719240|
|---+--------------------+--------------------+----------------------|
| 6 |   0.000044391759047| 3.641429976152e-006| 2449423.2956510304940|
+--------------------------------------------------------------------+

       力学時から日本標準時に変換する

       JST = TD + (9/24)
           = 2449423.2956510304940 + 0.375
           = 2449423.6706510314940
           = 1994/ 3/12  16: 5:44



  ２回目の朔

  （１回目の朔の時刻）＋３０を仮の朔の時刻とみなして計算する。

       t = t + 30
         = 2449423.670651031494 + 30
         = 2449453.670651031494

       日本標準時から力学時に変換する。

       TD = JST - (9/24)
          = 2449453.670651031494 - 0.375
          = 2449453.295651030494

       この時の月と太陽の黄経は、

       λmoon = 24.23809602459182
       λsun  = 21.16941167248130

       なので、黄経差Δλは、

       Δλ = λmoon - λsun
            = 3.06868435211052

       Δｔ = Δλ＊（29.530589/360）
            = 0.2517223788136

       ｔ = ｔ - Δｔ
          = 2449453.295651030494 - 0.2517223788136
          = 2449453.004392865191


  以下、反復計算の状況を表にしました。（λmoon・λsun・Δλ）
+--------------------------------------------------------------+
| N |      λmoon      |       λsun      |        Δλ        |
|---+------------------+------------------+--------------------|
| 1 | 21.26178396195252| 20.92228934659505|  0.3394946153574736|
|---+------------------+------------------+--------------------|
| 2 | 20.93230210025121| 20.89494753399944|  0.0373545662517785|
|---+------------------+------------------+--------------------|
| 3 | 20.89604657998604| 20.89193909000796|  0.0041074899780824|
|---+------------------+------------------+--------------------|
| 4 | 20.89205990824263| 20.89160828264627|  0.0004516255963551|
|---+------------------+------------------+--------------------|
| 5 | 20.89162156642899| 20.89157190980348|  0.0000496566255102|
+--------------------------------------------------------------+

（Δλ・Δｔ・ｔ）
+--------------------------------------------------------------------+
| N |        Δλ        |        Δｔ        |          ｔ          |
|---+--------------------+--------------------+----------------------|
| 1 |  0.3394946153574736| 2.784854431621e-002|2449453.01608010759689|
|---+--------------------+--------------------+----------------------|
| 2 |  0.0373545662517785| 3.064173175707e-003|2449453.01301593442119|
|---+--------------------+--------------------+----------------------|
| 3 |  0.0041074899780824| 3.369349954566e-004|2449453.01267899942573|
|---+--------------------+--------------------+----------------------|
| 4 |  0.0004516255963551| 3.704658296623e-005|2449453.01264195284276|
|---+--------------------+--------------------+----------------------|
| 5 |  0.0000496566255102| 4.073303886298e-006|2449453.01263787953888|
+--------------------------------------------------------------------+

       力学時から日本標準時に変換する

       JST = TD + (9/24)
           = 2449453.01263787953888 + 0.375
           = 2449453.38763788795388
           = 1994/ 4/11   9:18:11



  ３回目の朔

  （２回目の朔の時刻）＋３０を仮の朔の時刻とみなして計算する。

       t = t + 30
         = 2449453.38763788795388 + 30
         = 2449483.38763788795388

       日本標準時から力学時に変換する

       TD = JST - (9/24)
          = 2449483.38763788795388 + 0.375
          = 2449483.01263787953888

       この時の月と太陽の黄経は、

       λmoon = 53.34937446649215
       λsun  = 50.09737887498562

       なので、黄経差Δλは、

       Δλ = λmoon - λsun
            = 3.251995591506535

       Δｔ = Δλ＊（29.530589/360）
            = 0.2667592923405

       ｔ = ｔ - Δｔ
          = 2449483.01263787953888 - 0.2667592923405
          = 2449482.74587858735898


  以下、反復計算の状況を表にしました。（λmoon・λsun・Δλ）
+--------------------------------------------------------------+
| N |      λmoon      |       λsun      |        Δλ        |
|---+------------------+------------------+--------------------|
| 1 | 50.18908137261317| 49.83951459282018|  0.3495667797929940|
|---+------------------+------------------+--------------------|
| 2 | 49.84964722332275| 49.81179410307574|  0.0378531202470072|
|---+------------------+------------------+--------------------|
| 3 | 49.81289434841556| 49.80879234564254|  0.0041020027730241|
|---+------------------+------------------+--------------------|
| 4 | 49.80891161043991| 49.80846705602434|  0.0004445544155729|
|---+------------------+------------------+--------------------|
| 5 | 49.80847998175627| 49.80843180276968|  0.0000481789865887|
+--------------------------------------------------------------+

（Δλ・Δｔ・ｔ）
+--------------------------------------------------------------------+
| N |        Δλ        |        Δｔ        |          ｔ          |
|---+--------------------+--------------------+----------------------|
| 1 |  0.3495667797929940| 2.867475806145e-002| 2449482.7172038292975|
|---+--------------------+--------------------+----------------------|
| 2 |  0.0378531202470072| 3.105069267728e-003| 2449482.7140987600298|
|---+--------------------+--------------------+----------------------|
| 3 |  0.0041020027730241| 3.364848832418e-004| 2449482.7137622751465|
|---+--------------------+--------------------+----------------------|
| 4 |  0.0004445544155729| 3.646653815116e-005| 2449482.7137258086084|
|---+--------------------+--------------------+----------------------|
| 5 |  0.0000481789865887| 3.952094031632e-006| 2449482.7137218565143|
+--------------------------------------------------------------------+

       力学時から日本標準時に変換する

       JST = TD + (9/24)
           = 2449482.7137218565143 + 0.375
           = 2449483.0887218565143
           = 1994/ 5/11   2: 7:45



  ４回目の朔

  （３回目の朔の時刻）＋３０を仮の朔の時刻とみなして計算する。

       t = t + 30
         = 2449483.0887218565143 + 30
         = 2449513.0887218565143

       日本標準時から力学時に変換する

       TD = JST - (9/24)
          = 2449513.0887218565143 - 0.375
          = 2449512.7137218565143


       この時の月と太陽の黄経は、

       λmoon = 82.69589256228039
       λsun  = 78.63143984057999

       なので、黄経差Δλは、

       Δλ = λmoon - λsun
            = 4.064452721700405

       Δｔ = Δλ＊（29.530589/360）
            = 0.3334046745402

       ｔ = ｔ - Δｔ
          = 2449512.7137218565143 - 0.3334046745402
          = 2449512.3803171819347


  以下、反復計算の状況を表にしました。（λmoon・λsun・Δλ）
+--------------------------------------------------------------+
| N |      λmoon      |       λsun      |        Δλ        |
|---+------------------+------------------+--------------------|
| 1 | 78.62336688640019| 78.31247708690155|  0.3108897994986393|
|---+------------------+------------------+--------------------|
| 2 | 78.31267542440926| 78.28807857854559|  0.0245968458636696|
|---+------------------+------------------+--------------------|
| 3 | 78.28809913764832| 78.28614822140861|  0.0019509162397071|
|---+------------------+------------------+--------------------|
| 4 | 78.28614988257434| 78.28599511372572|  0.0001547688486170|
|---+------------------+------------------+--------------------|
| 5 | 78.28599524569950| 78.28598296748378|  0.0000122782157206|
+--------------------------------------------------------------+

（Δλ・Δｔ・ｔ）
+--------------------------------------------------------------------+
| N |        Δλ        |        Δｔ        |          ｔ          |
|---+--------------------+--------------------+----------------------|
| 1 |  0.3108897994986393| 2.550210803691e-002| 2449512.3548150738978|
|---+--------------------+--------------------+----------------------|
| 2 |  0.0245968458636696| 2.017664849712e-003| 2449512.3527974090481|
|---+--------------------+--------------------+----------------------|
| 3 |  0.0019509162397071| 1.600325156895e-004| 2449512.3526373765324|
|---+--------------------+--------------------+----------------------|
| 4 |  0.0001547688486170| 1.269559794031e-005| 2449512.3526246809344|
|---+--------------------+--------------------+----------------------|
| 5 |  0.0000122782157206| 1.007174839161e-006| 2449512.3526236737596|
+--------------------------------------------------------------------+

       力学時から日本標準時に変換する

       JST = TD + (9/24)
           = 2449512.3526236737596 + 0.375
           = 2449512.7276236737596
           = 1994/ 6/ 9  17:27:46



  ５回目の朔

  （４回目の朔の時刻）＋３０を仮の朔の時刻とみなして計算する。

       t = t + 30
         = 2449512.7276236737596 + 30
         = 2449542.7276236737596

       日本標準時から力学時に変換する

       TD = JST - (9/24)
          = 2449542.7276236737596 - 0.375
          = 2449542.3526236737596


       この時の月と太陽の黄経は、

       λmoon = 112.2766488159473
       λsun  = 106.9141295248953

       なので、黄経差Δλは、

       Δλ = λmoon - λsun
            = 5.362519291051981

       Δｔ = Δλ＊（29.530589/360）
            = 0.4398843144129

       ｔ = ｔ - Δｔ
          = 2449542.3526236737596 - 0.4398843144129
          = 2449541.9127393592087


  以下、反復計算の状況を表にしました。（λmoon・λsun・Δλ）
+--------------------------------------------------------------+
| N |      λmoon      |       λsun      |        Δλ        |
|---+------------------+------------------+--------------------|
| 1 | 106.6230744129181| 106.4945366712398|  0.1285377416783149|
|---+------------------+------------------+--------------------|
| 2 | 106.4881685146168| 106.4844792357378|  0.0036892788798184|
|---+------------------+------------------+--------------------|
| 3 | 106.4842968794353| 106.4841905681576|  0.0001063112777189|
+--------------------------------------------------------------+

（Δλ・Δｔ・ｔ）
+--------------------------------------------------------------------+
| N |        Δλ        |        Δｔ        |          ｔ          |
|---+--------------------+--------------------+----------------------|
| 1 |  0.1285377416783149| 1.054387561247e-002| 2449541.9021954835962|
|---+--------------------+--------------------+----------------------|
| 2 |  0.0036892788798184| 3.026293841842e-004| 2449541.9018928542121|
|---+--------------------+--------------------+----------------------|
| 3 |  0.0001063112777189| 8.720651801058e-006| 2449541.9018841335603|
+--------------------------------------------------------------------+

       力学時から日本標準時に変換する

       JST = TD + (9/24)
           = 2449541.9018841335603 + 0.375
           = 2449542.2768841335603
           = 1994/ 7/ 9   6:38:42



  結果をまとめると...

  +-------------------------------------------------------+
  | N |        ユリウス日       |        日付／時刻       |
  |---+-------------------------+-------------------------|
  | 1 |  2449423.6706510314940  |   1994/ 3/12  16: 5:44  |
  |---+-------------------------+-------------------------|
  | 2 |  2449453.3876378879538  |   1994/ 4/11   9:18:11  |
  |---+-------------------------+-------------------------|
  | 3 |  2449483.0887218565143  |   1994/ 5/11   2: 7:45  |
  |---+-------------------------+-------------------------|
  | 4 |  2449512.7276236737596  |   1994/ 6/ 9  17:27:46  |
  |---+-------------------------+-------------------------|
  | 5 |  2449542.2768841335603  |   1994/ 7/ 9   6:38:42  |
  +-------------------------------------------------------+

と言った次第です。

13. 太陽の黄経の求め方

今まで、黄経の計算方法については後で説明すると言っていましたが、 本章で 取り上げる事にします。

本計算式は、 『天体位置略算式の解説』という文献の「太陽の位置」の項にあ る太陽の黄経の計算式を引用しました。 太陽の黄経は（ユリウス世紀を単位とした時刻の）cos関数の和になっており、各 係数は表の形式に編集されています。すなわち、

              18
       λsun＝Σ Ａ*ｃｏｓ（ｋ*ｔ+θ0）
              n=1

で表される式の係数 Ａ・ｋ・ｔ・θ0 が表になっています。

       ※各変数の意味は以下の通りです...
       Ａ：振幅(Amplitude Of Vibration)
       ｋ：角速度(Angular Velocity)
       ｔ：時刻(Time)
       θ0：初期位相角(Initial Phase)

       Σ（シグマ）をつかわないで（係数をあてはめて）表記すると、

       λsun＝ 0.0004     * cos( 31557.0  * t + 161.0 )
             + 0.0004     * cos( 29930.0  * t +  48.0 )
             + 0.0005     * cos(  2281.0  * t + 221.0 )
                             | 途中省略
             + 0.0048 * t * cos( 35999.05 * t + 267.52 )
             + 1.9147     * cos( 35999.05 * t + 267.52 )
             + 36000.7695 * t
             + 280.4659

といった式で表されます。途中省略された部分の係数は表を参照して下さい。 な お、 第1項から第16項が摂動項、第17,18項は比例項となっています。（比例項は 上のように、cosの関数ではなく単なる一次関数の形で表される。） ｔは時刻引数で、西暦2000年1月1日12時（力学時）（2451545.0JD）からの経過時 間をユリウス世紀（36525日）を単位として測った時間です。なお、この時刻は力 学時（DT）によって表します。任意の時刻（時刻は力学時）JDに対するｔは、 以 下のように計算します。

       ｔ＝（JD+0.5-2451545.0）／36525

三角関数の引数は、ラジアンが単位ですから、ＣＯＳの引数に換算係数 k をかけ てあげます。（案外と忘れやすいので注意しましょう。）換算係数 k は、

       ｋ=π／180.0
         =0.01745329251994298

と言う値です。本スクリプトでは、BEGIN アクションの冒頭で定義しています。 1994年11月 8日 16:00(JST)の場合には（UTCとDTの差（ΔＴ）は0とする。）以 下の要領で計算します。

       1994年11月 8日  = 2449664
       16:00(JST) = 7:00(DT) ...... DT=JST-9:00
                  = 0.2916666665（日）

       ｔ＝（24496640.2916666665+0.5-2451545）／36525
         ＝ -0.05147729865

と言ったぐあいです。

補足説明しておきますと、（）の中は元式では、JD-2451545 となっています。 JDは日付と時刻を同時に示すため時刻を「日」の単位の小数で表しますが、 この 時に使用する時刻が天文時を使用する事になっています。 本スクリプトでは、天文時を使用しないようにしたために 0.5 日だけ時刻が違っ てきます。そこで、不足している 0.5日を補って JD+0.5-2451545 となる訳で す。 （詳しくは、 ユリウス日（ユリウス通日［つうじつ］）の項を参照下さ い。）

さらに、第15項の振幅は定数ではなくて t に関する一次関数の形（-0.0048*t） になっていますので注意して下さい。 また、元の表では、第18項から順に掲載さ れていますが、 実際に適用する際にアンダーフローを避けるためこのような順番 になっています。

    太陽黄経の計算式 係数表
  +--------------------------------------+
  | n |     Ａ     |    ｋ    |    θ0   |
  |---+------------+----------+----------|
  |  1|      0.0004|   31557.0|     161.0|
  |---+------------+----------+----------|
  |  2|      0.0004|   29930.0|      48.0|
  |---+------------+----------+----------|
  |  3|      0.0005|    2281.0|     221.0|
  |---+------------+----------+----------|
  |  4|      0.0005|     155.0|     118.0|
  |---+------------+----------+----------|
  |  5|      0.0006|   33718.0|     316.0|
  |---+------------+----------+----------|
  |  6|      0.0007|    9038.0|      64.0|
  |---+------------+----------+----------|
  |  7|      0.0007|    3035.0|     110.0|
  |---+------------+----------+----------|
  |  8|      0.0007|   65929.0|      45.0|
  |---+------------+----------+----------|
  |  9|      0.0013|   22519.0|     352.0|
  |---+------------+----------+----------|
  | 10|      0.0015|   45038.0|     254.0|
  |---+------------+----------+----------|
  | 11|      0.0018|  445267.0|     208.0|
  |---+------------+----------+----------|
  | 12|      0.0018|      19.0|     159.0|
  |---+------------+----------+----------|
  | 13|      0.0020|   32964.0|     158.0|
  |---+------------+----------+----------|
  | 14|      0.0200|   71998.1|     265.1|
  |---+------------+----------+----------|
  | 15|   -0.0048*t|  35999.05|    267.52|
  |---+------------+----------+----------|
  | 16|      1.9147|  35999.05|    267.52|
  |---+------------+----------+----------|
  | 17| 36000.7695t|         0|         0|
  |---+------------+----------+----------|
  | 18|    280.4659|         0|         0|
  +--------------------------------------+

本計算式を実際に適用する際の注意事項として、 「アンダーフローを避ける」 の一言に尽きます。特に、ユリウス日の小数で時刻を示す部分では、 （アンダー フローの余裕を考えて）15桁以上確保したいところです。 （アンダーフローのチ ェックスクリプト参照）また、略算式と精密な天体暦との差は、 1970年から2030 年の間0.1分角（約0.0017°）以内ということです。（もちろん、これ以外の期間 ではもっと大きくなりますが、旧暦計算には影響しない精度と思います。）

jgawkではどの位の桁数でアンダーフローが発生するのか簡単なスクリプトで調 べてみましたので参考にして下さい。 （jgawk 2.11.1 + 3.0 を使用した結果で す。）

            アンダーフローのチェックスクリプト

          1 BEGIN{
          2 printf("16桁 --> %20.20f\n", (1+1e-16) );
          3 printf("15桁 --> %20.20f\n", (1+1e-15) );
          4 printf("14桁 --> %20.20f\n", (1+1e-14) );
          5 }

                            --------- + ---------

       実行結果

       a:\jgawk -ftest.awk
       16桁 --> 1.00000000000000000000
       15桁 --> 1.00000000000000100000
       14桁 --> 1.00000000000001000000

このスクリプトの結果をみる限りでは、 15桁程度の精度を持っている（らしい） 事が解りました。

アンダーフローを避けるために本スクリプトでは、 以下のように行っています。 （このほかにもアンダーフローを避ける方法がありますが、 凝り過ぎると計算が 遅くなってしまいますのでこの辺でやめておくことにしました。）

• i ｃｏｓの引数（すなわち、ｋ＊ｔ＋θ0の計算結果）を正規化（値を0°～ 360°の範囲に制限する事）してから加算する。 （比例項の計算も同様に行 います。）
• ii 振幅（Ａのこと）の値が小さいものから順に加算する。
• iii ユリウス日を整数部と、 小数部に分割しておき、 時刻引数を計算する際 に、別々計算してから、合成する。

では、 例によって例題をやってみます。 計算する時刻は、 1994年11月 8日 16:00(JST)とします（ｔは、 前の例で求めてあるのでこれを使用する。前にスク ロールして戻ってもらうのも何ですので、 もう一度のせましょう。 ｔ＝ -0. 05147729865です。）項番号nで、「比」とあるのは、比例項の合計です。

    1994年11月 8日 16:00(JST)の太陽黄経の計算過程
  +------------------------------------------------------------+
  | n |     ｋ*ｔθ0     | A*cos（k*t+θ0） |       合計       |
  |---+------------------+------------------+------------------|
  |  1|    336.5308863791|      0.0003669100| ---------------- |
  |---+------------------+------------------+------------------|
  |  2|    307.2844512889|      0.0002423090|      0.0006092190|
  |---+------------------+------------------+------------------|
  |  3|    103.5802817705|     -0.0001174038|      0.0004918152|
  |---+------------------+------------------+------------------|
  |  4|    110.0210187086|     -0.0001711824|      0.0003206327|
  |---+------------------+------------------+------------------|
  |  5|     20.2884439880|      0.0005627753|      0.0008834081|
  |---+------------------+------------------+------------------|
  |  6|    318.7481747661|      0.0005262732|      0.0014096812|
  |---+------------------+------------------+------------------|
  |  7|    313.7663985854|      0.0004842038|      0.0018938851|
  |---+------------------+------------------+------------------|
  |  8|    251.1531770474|     -0.0002261274|      0.0016677576|
  |---+------------------+------------------+------------------|
  |  9|    272.7827116130|      0.0000631129|      0.0017308705|
  |---+------------------+------------------+------------------|
  | 10|     95.5654232259|     -0.0001454734|      0.0015853971|
  |---+------------------+------------------+------------------|
  | 11|    326.8576602765|      0.0015071669|      0.0030925640|
  |---+------------------+------------------+------------------|
  | 12|    158.0219313256|     -0.0016691889|      0.0014233751|
  |---+------------------+------------------+------------------|
  | 13|    261.1023271730|     -0.0003093405|      0.0011140346|
  |---+------------------+------------------+------------------|
  | 14|    158.8323037871|     -0.0186505508|     -0.0175365162|
  |---+------------------+------------------+------------------|
  | 15|    214.3861518935|      0.0002039119|     -0.0177404281|
  |---+------------------+------------------+------------------|
  | 16|    214.3861518935|     -1.5801062208|     -1.5978466489|
  |---+------------------+------------------+------------------|
  | 比| ---------------- |    227.2435366785|    225.6456900296|
  +------------------------------------------------------------+

  答え  λsun=225.6456900296°

14. 月の黄経の求め方

太陽の黄経の計算方法の次は、月の黄経の計算方法の説明をします。 太陽と同 じく（ユリウス世紀を単位とした時刻の）cos関数の和になっており、各係数は表 の形式に編集されています。 項数は63項と大変多くなっていますが計算の要領は 太陽の場合と全く同様です。

              63
       λsun＝Σ Ａ*ｃｏｓ（ｋ*ｔ+θ0）
              n=1


       Σ（シグマ）をつかわないで（係数をあてはめて）表記すると、

       λsun＝ 0.0003 * cos(  2322131.0  * t + 191.0 )
             + 0.0003 * cos(     4067.0  * t +  70.0 )
             + 0.0003 * cos(   549197.0  * t + 220.0 )
                         | 以下省略
             + 1.2740 * cos(  413335.35 * t +  10.74 )
             + 6.2888 * cos( 477198.868 * t + 44.963 )
             + 481267.8809 * t +  218.3162

ｔは、 （太陽の黄経の時と同じく）時刻引数で西暦2000年1月1日12時（力学 時）（2451545.0JD）からの経過時間をユリウス世紀（36525日）を単位として測 った時間です。なお、この時刻は力学時（DT）によって表します。 三角関数の引数は、ラジアンが単位ですから、ＣＯＳの引数に換算係数 k をかけ てあげます。（案外と忘れやすいので注意しましょう。）換算係数 k は、

       ｋ=π／180.0
         =0.01745329251994298...

と言う値です。本スクリプトでは、BEGIN アクションの冒頭で定義しています。 元の表では、 第63項から順に掲載されていますが、実際に適用する際にアンダー フローを避けるためこのような順番になっています。

    月の黄経の計算式 係数表
  +-------------------------------------------+
  | n |        Ａ       |    ｋ    |    θ0   |
  |---+-----------------+----------+----------|
  |  1|           0.0003| 2322131.0|     191.0|
  |---+-----------------+----------+----------|
  |  2|           0.0003|    4067.0|      70.0|
  |---+-----------------+----------+----------|
  |  3|           0.0003|  549197.0|     220.0|
  |---+-----------------+----------+----------|
  |  4|           0.0003| 1808933.0|      58.0|
  |---+-----------------+----------+----------|
  |  5|           0.0003|  349472.0|     337.0|
  |---+-----------------+----------+----------|
  |  6|           0.0003|  381404.0|     354.0|
  |---+-----------------+----------+----------|
  |  7|           0.0003|  958465.0|     340.0|
  |---+-----------------+----------+----------|
  |  8|           0.0004|   12006.0|     187.0|
  |---+-----------------+----------+----------|
  |  9|           0.0004|   39871.0|     223.0|
  |---+-----------------+----------+----------|
  | 10|           0.0005|  509131.0|     242.0|
  |---+-----------------+----------+----------|
  | 11|           0.0005| 1745069.0|      24.0|
  |---+-----------------+----------+----------|
  | 12|           0.0005| 1908795.0|      90.0|
  |---+-----------------+----------+----------|
  | 13|           0.0006| 2258267.0|     156.0|
  |---+-----------------+----------+----------|
  | 14|           0.0006|  111869.0|      38.0|
  |---+-----------------+----------+----------|
  | 15|           0.0007|   27864.0|     127.0|
  |---+-----------------+----------+----------|
  | 16|           0.0007|  485333.0|     186.0|
  |---+-----------------+----------+----------|
  | 17|           0.0007|  405201.0|      50.0|
  |---+-----------------+----------+----------|
  | 18|           0.0007|  790672.0|     114.0|
  |---+-----------------+----------+----------|
  | 19|           0.0008| 1403732.0|      98.0|
  |---+-----------------+----------+----------|
  | 20|           0.0009|  858602.0|     129.0|
  |---+-----------------+----------+----------|
  | 21|           0.0011| 1920802.0|     186.0|
  |---+-----------------+----------+----------|
  | 22|           0.0012| 1267871.0|     249.0|
  |---+-----------------+----------+----------|
  | 23|           0.0016| 1856938.0|     152.0|
  |---+-----------------+----------+----------|
  | 24|           0.0018|  401329.0|     274.0|
  |---+-----------------+----------+----------|
  | 25|           0.0021|  341337.0|      16.0|
  |---+-----------------+----------+----------|
  | 26|           0.0021|   71998.0|      85.0|
  |---+-----------------+----------+----------|
  | 27|           0.0021|  990397.0|     357.0|
  |---+-----------------+----------+----------|
  | 28|           0.0022|  818536.0|     151.0|
  |---+-----------------+----------+----------|
  | 29|           0.0023|  922466.0|     163.0|
  |---+-----------------+----------+----------|
  | 30|           0.0024|   99863.0|     122.0|
  |---+-----------------+----------+----------|
  | 31|           0.0026| 1379739.0|      17.0|
  |---+-----------------+----------+----------|
  | 32|           0.0027|  918399.0|     182.0|
  |---+-----------------+----------+----------|
  | 33|           0.0028|    1934.0|     145.0|
  |---+-----------------+----------+----------|
  | 34|           0.0037|  541062.0|     259.0|
  |---+-----------------+----------+----------|
  | 35|           0.0038| 1781068.0|      21.0|
  |---+-----------------+----------+----------|
  | 36|           0.0040|     133.0|      29.0|
  |---+-----------------+----------+----------|
  | 37|           0.0040| 1844932.0|      56.0|
  |---+-----------------+----------+----------|
  | 38|           0.0040| 1331734.0|     283.0|
  |---+-----------------+----------+----------|
  | 39|           0.0050|  481266.0|     205.0|
  |---+-----------------+----------+----------|
  | 40|           0.0052|   31932.0|     107.0|
  |---+-----------------+----------+----------|
  | 41|           0.0068|  926533.0|     323.0|
  |---+-----------------+----------+----------|
  | 42|           0.0079|  449334.0|     188.0|
  |---+-----------------+----------+----------|
  | 43|           0.0085|  826671.0|     111.0|
  |---+-----------------+----------+----------|
  | 44|           0.0100| 1431597.0|     315.0|
  |---+-----------------+----------+----------|
  | 45|           0.0107| 1303870.0|     246.0|
  |---+-----------------+----------+----------|
  | 46|           0.0110|  489205.0|     142.0|
  |---+-----------------+----------+----------|
  | 47|           0.0125| 1443603.0|      52.0|
  |---+-----------------+----------+----------|
  | 48|           0.0154|   75870.0|      41.0|
  |---+-----------------+----------+----------|
  | 49|           0.0304|  513197.9|     222.5|
  |---+-----------------+----------+----------|
  | 50|           0.0347|  445267.1|      27.9|
  |---+-----------------+----------+----------|
  | 51|           0.0409|  441199.8|      47.4|
  |---+-----------------+----------+----------|
  | 52|           0.0458|  854535.2|     148.2|
  |---+-----------------+----------+----------|
  | 53|           0.0533| 1367733.1|     280.7|
  |---+-----------------+----------+----------|
  | 54|           0.0571|  377336.3|      13.2|
  |---+-----------------+----------+----------|
  | 55|           0.0588|   63863.5|     124.2|
  |---+-----------------+----------+----------|
  | 56|           0.1144|  966404.0|     276.5|
  |---+-----------------+----------+----------|
  | 57|           0.1851|   35999.0|     87.53|
  |---+-----------------+----------+----------|
  | 58|           0.2136|  954397.7|    179.93|
  |---+-----------------+----------+----------|
  | 59|           0.6583|  890534.2|     145.7|
  |---+-----------------+----------+----------|
  | 60|           1.2740|  413335.3|     10.74|
  |---+-----------------+----------+----------|
  | 61|           6.2888| 477198.86|    44.963|
  |---+-----------------+----------+----------|
  | 62|         218.3162|         0|         0|
  |---+-----------------+----------+----------|
  | 63|  481267.8809 * t|         0|         0|
  +-------------------------------------------+

では、 例によって例題をやってみます。 計算する時刻は、 1994年11月 8日 16:00(JST)とします（ｔは、前の例で求めてあるのでこれを使用する。 ちなみに ｔ＝ -0.05147729865です。）項番号nで、 「比」とあるのは、 比例項の合計で す。

    1994年11月 8日 16:00(JST)の太陽黄経の計算過程
  +------------------------------------------------------------+
  | n |     ｋ*ｔθ0     | A*cos（k*t+θ0） |       合計       |
  |---+------------------+------------------+------------------|
  |  1|    173.9689995338|     -0.0002983396| ---------------- |
  |---+------------------+------------------+------------------|
  |  2|    220.6418263746|     -0.0005259784|     -0.0002276388|
  |---+------------------+------------------+------------------|
  |  3|     28.8220111772|     -0.0002631419|      0.0002628365|
  |---+------------------+------------------+------------------|
  |  4|    179.0157141151|     -0.0005630976|     -0.0002999557|
  |---+------------------+------------------+------------------|
  |  5|    347.1254848263|     -0.0002706395|      0.0002924581|
  |---+------------------+------------------+------------------|
  |  6|    160.3523842101|     -0.0005531730|     -0.0002825335|
  |---+------------------+------------------+------------------|
  |  7|    320.8109456952|     -0.0003206535|      0.0002325196|
  |---+------------------+------------------+------------------|
  |  8|    288.9635523613|     -0.0001906668|      0.0001299866|
  |---+------------------+------------------+------------------|
  |  9|    330.5486253706|      0.0001576425|      0.0003483093|
  |---+------------------+------------------+------------------|
  | 10|    313.3114590442|      0.0005006244|      0.0003429819|
  |---+------------------+------------------+------------------|
  | 11|    192.5619153474|      0.0000125937|     -0.0004880308|
  |---+------------------+------------------+------------------|
  | 12|    110.3897159399|     -0.0001616082|     -0.0001742019|
  |---+------------------+------------------+------------------|
  | 13|    186.5152007661|     -0.0007577333|     -0.0005961251|
  |---+------------------+------------------+------------------|
  | 14|     39.2860768875|     -0.0002933368|      0.0004643965|
  |---+------------------+------------------+------------------|
  | 15|    132.6365503079|     -0.0007674786|     -0.0004741418|
  |---+------------------+------------------+------------------|
  | 16|     42.3682124095|     -0.0002502981|      0.0005171805|
  |---+------------------+------------------+------------------|
  | 17|     71.3471081434|     -0.0000264142|      0.0002238839|
  |---+------------------+------------------+------------------|
  | 18|     92.3413187281|     -0.0000550109|     -0.0000285966|
  |---+------------------+------------------+------------------|
  | 19|    197.6686059717|     -0.0008172732|     -0.0007622623|
  |---+------------------+------------------+------------------|
  | 20|    210.4884211691|     -0.0015928318|     -0.0007755585|
  |---+------------------+------------------+------------------|
  | 21|    308.3017910026|     -0.0009110478|      0.0006817839|
  |---+------------------+------------------+------------------|
  | 22|    142.4258783884|     -0.0018621260|     -0.0009510782|
  |---+------------------+------------------+------------------|
  | 23|    321.8479922349|     -0.0006039266|      0.0012581994|
  |---+------------------+------------------+------------------|
  | 24|    134.6672085313|     -0.0018693046|     -0.0012653780|
  |---+------------------+------------------+------------------|
  | 25|     84.8933093757|     -0.0016823824|      0.0001869223|
  |---+------------------+------------------+------------------|
  | 26|    338.7374515169|      0.0002746674|      0.0019570498|
  |---+------------------+------------------+------------------|
  | 27|    134.0378450791|     -0.0011851126|     -0.0014597801|
  |---+------------------+------------------+------------------|
  | 28|    134.9778690360|     -0.0027401466|     -0.0015550339|
  |---+------------------+------------------+------------------|
  | 29|    196.9422199368|     -0.0049403245|     -0.0022001780|
  |---+------------------+------------------+------------------|
  | 30|     21.3225245262|     -0.0027046085|      0.0022357160|
  |---+------------------+------------------+------------------|
  | 31|    271.7634325746|     -0.0026245991|      0.0000800094|
  |---+------------------+------------------+------------------|
  | 32|     65.3003935621|     -0.0014963748|      0.0011282242|
  |---+------------------+------------------+------------------|
  | 33|     45.4429044034|      0.0004681603|      0.0019645351|
  |---+------------------+------------------+------------------|
  | 34|    126.5898357267|     -0.0017373448|     -0.0022055050|
  |---+------------------+------------------+------------------|
  | 35|    136.4306411058|     -0.0044905989|     -0.0027532541|
  |---+------------------+------------------+------------------|
  | 36|     22.1535192790|     -0.0007858917|      0.0037047072|
  |---+------------------+------------------+------------------|
  | 37|    123.8844398735|     -0.0030159704|     -0.0022300787|
  |---+------------------+------------------+------------------|
  | 38|    128.9311544548|     -0.0055295149|     -0.0025135445|
  |---+------------------+------------------+------------------|
  | 39|    270.7263860349|     -0.0054661275|      0.0000633874|
  |---+------------------+------------------+------------------|
  | 40|    263.2268993838|     -0.0060794039|     -0.0006132764|
  |---+------------------+------------------+------------------|
  | 41|    147.5840463114|     -0.0118198191|     -0.0057404151|
  |---+------------------+------------------+------------------|
  | 42|     97.4994866511|     -0.0128509058|     -0.0010310867|
  |---+------------------+------------------+------------------|
  | 43|     36.2100444866|     -0.0059926233|      0.0068582825|
  |---+------------------+------------------+------------------|
  | 44|     60.2536789809|      0.0068582825|      0.0068582825|
  |---+------------------+------------------+------------------|
  | 45|     86.2946041469|      0.0075497837|      0.0006915013|
  |---+------------------+------------------+------------------|
  | 46|    159.0481120216|     -0.0027229075|     -0.0102726913|
  |---+------------------+------------------+------------------|
  | 47|    259.2172313422|     -0.0050614811|     -0.0023385736|
  |---+------------------+------------------+------------------|
  | 48|     95.4173511290|     -0.0074000548|     -0.0023385736|
  |---+------------------+------------------+------------------|
  | 49|     84.4584331486|     -0.0044643916|      0.0029356631|
  |---+------------------+------------------+------------------|
  | 50|    146.7525125466|     -0.0334843555|     -0.0290199639|
  |---+------------------+------------------+------------------|
  | 51|     15.6261293616|      0.0059039734|      0.0393883289|
  |---+------------------+------------------+------------------|
  | 52|     79.0362993347|      0.0146145405|      0.0087105671|
  |---+------------------+------------------+------------------|
  | 53|     73.4947324833|      0.0087105671|      0.0087105671|
  |---+------------------+------------------+------------------|
  | 54|     28.9465919446|      0.0586771349|      0.0499665678|
  |---+------------------+------------------+------------------|
  | 55|     76.6795374170|      0.0722244952|      0.0135473602|
  |---+------------------+------------------+------------------|
  | 56|    208.6326716819|     -0.0281855138|     -0.1004100090|
  |---+------------------+------------------+------------------|
  | 57|     34.3961518935|      0.1245500178|      0.1527355317|
  |---+------------------+------------------+------------------|
  | 58|     10.1125034182|      0.2772855495|      0.1527355317|
  |---+------------------+------------------+------------------|
  | 59|     23.4039955472|      0.8814251863|      0.6041396368|
  |---+------------------+------------------+------------------|
  | 60|    333.3527438381|      2.0201068075|      1.1386816211|
  |---+------------------+------------------+------------------|
  | 61|    320.0543546637|      6.8414397925|      4.8213329850|
  |---+------------------+------------------+------------------|
  | 比| ---------------- |    283.9457623839|    290.7872021764|
  +------------------------------------------------------------+

  λmoon= 290.7872021764

15. 計算を早くするために（高速化へのアプローチ）

二十四節気の説明の項で述べましたが、 天保暦法では二分二至を含む月を必 ず、 ２月・５月・８月・１１月と決め、この拘束条件のもとで中気の含まない月 を閏月としています。従って、 朔の時刻を５個・中気の時刻を４個（計算対象日 の直前にある二分二至の時刻を含める）を計算しないと閏月が決定できないの で、このあたりで時間を食ってしまうようです。

この事から、 朔の時刻の計算や中の時刻の計算で収束の早いアルゴリズムを採用 すれば良い訳です。（それが解っていながらできなかったのが、 悲しいような口 惜しいような...。）

高速化で思い付くのは、 黄経計算の項数を減らすのと、収束精度を甘くすること が考えられますが、 どちらにしても精度を落とした割に高速化するとは考えられ ず良い方法と言えないと思います。

そのほか、 初期値の与え方を工夫すると収束が早くなるかも知れません。（拙者 は未確認です。）収束が早いアルゴリズムを開発された場合には、 御教示下れば 幸いです。

最も簡単でかつ確実に高速化する方法は、 『表引き』です。つまり、まじめに シコシコと計算するのをサッとあきらめて、 朔日が何年何月何日（年月日でなく ユリウス日で定義しても良い）になるかと旧暦の月名と閏月かどうかの判別フラ グを予めに表にしておき、 実行時にこれを検索して新暦に対する旧暦を求めよう と言うわけです。（本スクリプトの『朔日行列』の形式で、 ファイル化するのも 良いかもしれない。）

この方法では、 表を検索するだけで良い訳ですから計算するよりも早いはずで す。 ただし、 表に載っていない日付の計算は（残念ですが）できませんが...。 （適用期間を限定して、 逃げる手もあるが。）もう一つの欠点は、適用期間を長 くするとその分だけ表が巨大になりディスクスペース上の制約がでてきます。 こ のような事から、、 任意の日付に対する旧暦を求めると言う用途には向いていな い方法と言えます。

16. 本スクリプトを過去の暦日を知るために使用する場合の注意事項

過去の暦日を求める場合に注意しなくてはいけないのが、 計算で求めた日付で はなくて、 あくまで歴史上の実施が基準となります。（過去に何度となく暦日変 更を行っていた（らしい）事が確認されています。）そのため、 色々な史料を見 ないと正確な暦日がはっきりしない場合があります。

また、 本スクリプトで計算された暦日は、『天保暦法に準拠』して求めたもので あることも忘れてはいけません。（「旧暦とは...」の項を参照して下さい。）

17. 無保証

本スクリプトによって得られた結果に対しては、 作者は一切責任を負いません し、損害を補償する義務も負いません。

18. 著作権について

著作権は（一応）保留いたしますが、 後述する再配布規定だけ厳守して下され ば、『煮るなと焼くなと勝手』とさせて頂きます。

冒頭でも述べましたが、『皆さんの創意工夫によって、 より良いプログラムに 育ててやって下さいませ。 本サンプルスクリプトが、その足掛かりになれば幸い と存じます...。』と言うのが、そのココロです。したがいまして、その範囲を逸 脱しない限り、『お咎めなし』と言う事にいたしょう？。

19. 配布規定について

本当は、「一般常識の許す範囲で...」としたい所ですが、世の中には色々な人 がいますので、ある程度の枠組みが必要なのは仕方の無いことなのでしょう。 と もあれ拙者が希望するのは以下のようなものです。

• i 本スクリプト・説明書を再配布する場合はオリジナルのアーカイブファイ ルに含まれるファイルを全て含み、オリジナルのスクリプト・説明書を改変 しないで下さい。 また、 頒布の都合でアーカイブ形式（.arc .zip .zoo 等）を変更する場合も同様に扱って下さい。
• ii 再配布する際の媒体に要するコストを除いて一切の金銭等の授受は禁止い たします。 （他の言語に移植したものを配布する場合にも適用いたしま す。）
• iii 他の言語に移植したり、 本スクリプトの一部または全部を引用して作成 し、これを配布する場合には、オリジナルのスクリプトと本説明書を必ず同 梱して下さい。 その場合の著作権につきましては、引用した部分のみ著者に帰属し、その外 は製作者に帰属します。
• iv 配布によって著作者が一切の制限を受ける可能性がないようにして下さ い。

なお、 いづれの場合におきましても、承諾・転載の報告など一切不要といたしま す。

20. 連絡先

何か気付いた点や、 質問・意見等がありましたら、 連絡いただければ幸いで す。当方の都合により、 月３～４回程度の頻度でアクセスしておりますのでせっ かくメールをいただいても、 返事が遅れるまたは、できない場合があるかと思い ますが、その際にはなにとぞ御容赦願います。悪しからず。また、 頂いたメール は会議室で公開することもありますので、予めご了承下さい。 できましたら、 電子会議室等にアップロードしていただければ幸いとと存じま す。

       各ネットのＩＤ

       NIFTY-Serve    HCB03121
       ASCII-NET      net06887

                            --------- + ---------

       電子会議室

       NIFTY-Serve    FGALEL 会議室 #12  ｽｸﾘﾌﾟﾄ天国ＡＷＫﾏﾙ秘物語
       ASCII-NET      salon.pool

21. 謝辞

tadfさんには、計算の結果が間違っている部分と、 πの値が間違っていたのを 指摘して頂きありがとうございました。

本説明書は ntf （ぬてふ）を使用して作成いたしました。拙者がめちゃくちゃ に書いた文章も（中身はともかくとして）体裁を整え見栄え良くしていただきと ても助かりました。作者の PoorTom さんに厚く御礼申し上げます。

22. 最後に...

当初は、today と言うソフトのソースの中にある inreki.c / inreki.h を移植 しようと考えておりましたが、 いくらソースを眺めても中で何をやっているのか さっぱり理解できなかったので、 移植は諦め（１からではなくて）零から製作す る事にしました。

『TODAY』について簡単な説明を加えておきますと、 TODAY は、森 佳史氏が製 作されたもので、今日が“何々の日”なのかを教えてくれるツールです。 実行形 式 or ソースのある場所は、以下の通りです。

       NIFTY-Serve

       FGALAP lib 4  688 GFA02107 931016  13K TD210MPB.LZH
              Disc: TODAY v2.10A → 2.10B 差分
       FGALAP lib 4  671 GFA02107 930915  14K TD210MPA.LZH
              Disc: TODAY Ver.2.10→2.10A 差分
       FGALAP lib 4  643 GFA02107 930730 104K TD210SRC.LZH
              Disc: TODAY Ver 2.10 ソース一式
       FGALAP lib 4  642 GFA02107 930730 213K TD210M  .LZH
              Disc: 今日は何の日 TODAY Ver 2.10
       FGALAP lib 4   91 PBA02244 910706 107K TODAY29M.LZH
              Disc: 今日は何の日 TODAY Ver 2.9


       ascii-net(POOL MS-DOS)

        No. File Name             Type Author   YY/MM/DD Bytes
       ---- ------------------------ - -------- -------- -------
       7202 TD210MPB.ISH             T pcs30807 93/10/16   15513
       Desc: TODAY Ver.2.10A から Ver.2.10B への差分データです。

       7071 TD210MPA.ISH             T pcs30807 93/09/15   16413
       Desc: TODAY Ver.2.10A への差分データです。

       6922 TD210SRC.ISH             T pcs30807 93/07/29  124025
       Desc: TODAY Ver.2.10 のソースプログラム一式です。

       6921 TD210M.ISH               T pcs30807 93/07/29  254573
       Desc: TODAY Ver.2.10 です。日めくりカレンダーのようなものです。

       6220 today29p.ish             T pcs30807 93/01/02   25481
       Desc: TODAY Ver.2.9 のバグを取る差分です。要BUPDATE

       4068 TDV29SRC.ISH             T pcs30807 91/06/19   57422
       Desc: TODAY Ver.2.9 のソースファイル一式です。

       4067 TODAY29M.ISH             T pcs30807 91/06/19  128405
       Desc: TODAY Ver.2.9 です。日めくりカレンダーのようなものです。

拙者の言葉の足らない所があり、 適切な説明ができたかどうかあやしい所も数 々あると思います。 構成を全く考えず、思い付くままに書きましたので、 最初から読んでも話がつな がらない所があると思いますが、あまり気にしないでください。 また、 誤って理解し説明してしまった点もあるかもしれません。バージョンアッ プについては今の所は全く考えておりませんが、 御指摘事項の内容によっては、 修正を行う用意はあります...。 ほとんど校正をしていませんので、誤字・脱字 等々結構あるかもしれませんが、 なにとぞ勘弁してください。 こんな「ごみファイル」をダウンロードして頂いた 上に、ここまで目を通して下さいましたあなたがたに深く感謝いたします。 最後 まで、付き合っていただきまして誠にありがとうございました。

23. 参考図書

本スクリプトを作成するにあたり、 以下の参考書を使用いたしました。コメン トは、拙者の個人的なものですので、あまり役に立たないかと思いますが、 参考 ということで理解して下さいませ。

23.1 旧暦計算のアルゴリズムについて

日経バイト: NIKKEI BYTE/JANUARY 1988 hobby p.239-243

大安,仏滅が一発でわかる六曜プログラム 鈴木 隆 著

本スクリプトで閏月の判別を行う部分のアルゴリズムを参考にさせていただき ました。 記事を見る限りでは、プログラムソースは、ＦＤサービスで入手可能だ そうです。しかし、だいぶ昔の記事ですので現在も入手ができるかどうかは、 定 かではありません。（なお、拙者はプログラムソースは持っておりません。）

23.2 暦の一般的な説明

暦の科学: 山崎 昭,久保良雄 著 ブルーバックス B-583 講談社

暦に関する基礎知識と、天体暦について述べられています。

理科年表読本 こよみと天文・今昔:内田正男 著 丸善株式会社

理科年表読本と言う事で、 理科年表の暦部を理解できるような構成となってい ます。 旧暦のしくみなども述べられておりスクリプトを作成する際にいろいろ と、 参考にさせていただきました。本書の著者は、『日本暦日原典』と言う新旧 暦対照表の著者でもあるため、専門家としての立場からの意見が目を引きます。 ちょっとお高くなりますが、 なるべくでしたら理科年表とセットで入手されると 良いかと思います。 （拙者は、丸善さんの回し者ではないのは言うまでもありま せんが...。）

時と暦:青木信仰 著 （財）東京大学出版会

時刻と暦にまつわる各種の問題を科学的な見地から述べられています。 また、 歴史・法律の分野からの問題も述べられており、 興味のある内容になっていま す。

暦の百科事典:暦の会 著 新人物往来社

百科事典と銘打っているだけあって、 暦に関する事柄はほとんど総ての方面で 網羅されており、 ほかに類を見ないものとなっています。暦や歴史に興味のある 方は、一読をお薦めいたします。 また、巻末の資料編にある新旧暦対照表（西暦593年～2000年）は目を引きます。

日本暦日原典:内田正男 著 雄山閣

445年から1872年までの朔・中気・節気の時刻を掲載したものです。 445年から 1684年までは、 計算値（暦日変更があった場合には、 脚注にその旨が記してあ る。）を掲載し、 1685年から1872年までは一般に使われた暦（頒暦）が連続して 残っているため、これをそのまま掲載しています。 新旧暦対照表の中でも最も信 頼性が高いものです。 巻末に暦を理解する目的で「暦法編」が付けられています。 ただ、やや難解のよ うに感じましたので、各種の入門書を読んでからの方が良いと思います。 もっとも、表を参照するだけであれば、 「暦法編」は理解しなくても一向に差し 支えありません。また過去に記録されている日蝕（日食）について論じられ

23.3 天体の位置計算など

天体位置略算式の解説: 暦計算研究会（井上圭典・鈴木邦裕）著 海文堂出版

太陽と惑星（水・金・火・木・土・天・海・冥）と月の位置の略算式とＢＡＳ ＩＣ（N88 BASIC）によるプログラム例が掲載されております。なおこのプログラ ムは、 『暦泉』と言う別売りのＦＤでも入手できます。本スクリプトでは、太陽 と月の黄経の略算式を本書のものを使用しております。

新こよみ便利帳: 暦計算研究会 編 恒星社厚生閣

西暦2000年までの二十四節気・雑節・月の朔望など暦の計算に必要な表はほと んど一通り揃っています。また、日食表・月食表なども掲載されています。 この ほか、 太陽と惑星（水・金・火・木・土・天・海・冥）と月の位置の略算式が掲 載されております。

天文計算入門: 長谷川一郎 著 恒星社厚生閣

球面三角法・各種座標系・二体問題・太陽系内天体の位置推算・軌道決定など 天文計算の入門書となっています。 本格的に天文計算をやってみたいとお考えで したら、お薦めいたします。

23.4 天体暦など

天文年鑑: 1994年版 誠文堂新光社

天文観測年表: 1994年版 地人書館

この２つの年鑑（年表）は、 ほとんど同じ構成となっています。天体観測用に 作られていますので、暦以外に各種の数表などが豊富に掲載されております。

理科年表: 平成６年版 第67冊 国立天文台 編 丸善株式会社

暦部・天文部・気象部・物理／化学部・地学部・生物部から構成され、 暦部お よび天文部は、 国立天文台・海上保安庁水路部の計算値を、気象部では、気象庁 のデータが掲載されています。

23.5 運勢暦など

神宮宝暦: 平成６年版 神宮館

「暦の予備知識」として、各種暦注の説明が詳しい。 九星占星術よる毎日の運 勢が掲載されています。