getRokuyo

引数に指定した旧暦の日付から六曜を求めます。

構文
getRokuyo( year, month, day )
引数
year
month
day
戻値
六曜

プログラム

////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // 【戻値】 // 六曜 ////////////////////////////////////////////////// FUNCTION getRokuyo(year, month, day) DIM arr[6] = "大安", "赤口", "先勝", "友引", "先負", "仏滅"" DIM d = getKyureki(year, month, day) DIM n = (d[2] + d[3]) MOD 6 RESULT = arr[n] FEND ////////////////////////////////////////////////// // 【引数】 // array : 要素を追加する配列(参照引数) // str : 追加する要素 // 【戻値】 // 処理後の配列の中の要素の数 ////////////////////////////////////////////////// FUNCTION arrayPush(var arr[], str) DIM res = RESIZE(arr, UBound(arr) + 1) arr[res] = str RESULT = res + 1 FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 中気と太陽黄経を格納した配列(0 : 中気, 1 : 太陽黄経) ////////////////////////////////////////////////// FUNCTION chuki(JD) JD = JD - 9/24 DIM t = (JD + 0.5 - 2451545) / 36525 DIM λsun = longitudeSun(t) DIM λsun0 = 30 * INT(λsun/30) REPEAT t = (JD + 0.5 - 2451545) / 36525 λsun = longitudeSun(t) DIM Δλ = λsun - λsun0 SELECT TRUE CASE Δλ > 180 Δλ = Δλ - 360 CASE Δλ < -180 Δλ = Δλ + 180 SELEND DIM Δt = Δλ * 365.2 / 360 JD = JD - Δt UNTIL Δt <= 1/86400 JD = JD + 9/24 DIM arr[1] = JD, λsun RESULT = SLICE(arr) FEND ////////////////////////////////////////////////// // 【引数】 // deg : 度数法 // 【戻値】 // 度数法から弧度法に変換した値 ////////////////////////////////////////////////// FUNCTION degToRad(deg) DIM pi = 3.14159265358979 RESULT = deg * pi / 180 FEND ////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // 【戻値】 // 旧暦を格納した配列(0 : 年, 1 : 月, 2 : 日) ////////////////////////////////////////////////// FUNCTION getKyureki(year, month, day) DIM tm = YMDToJD(year, month, day, 0, 0, 0) DIM chu[-1] // n*2︰ユリウス日、n*2+1︰Δλsun0、n=0,1,2 DIM tmp tmp = nishiNibun(tm) FOR n = 0 TO UBound(tmp) arrayPush(chu, tmp[n]) NEXT // 中気の計算 3回 chu[n]︰n+1回目 FOR n = 0 TO 2 tmp = chuki(chu[n*2] + 32) FOR n = 0 TO UBound(tmp) arrayPush(chu, tmp[n]) NEXT NEXT DIM saku[5] saku[0] = saku(chu[0]) // 朔の計算 FOR n = 1 TO 4 saku[n] = saku(saku[n-1] + 30) IFB ABS(INT(saku[n-1]) - INT(saku[n])) <= 26 THEN saku[n] = saku(saku[n-1] + 35) ENDIF NEXT // 朔の修正 SELECT TRUE CASE saku[1] <= INT(chu[0*2+0]) SHIFTARRAY(saku, -1) saku[4] = saku(saku[3] + 35) CASE saku[0] > INT(chu[0*2+0]) SHIFTARRAY(saku, 1) saku[0] = saku(saku[0] - 27) SELEND DIM kyureki[3] // 0︰年、1︰閏月、2︰月、3︰日 // 閏月検索 DIM flg = FALSE IF INT(saku[4]) <= INT(chu[3*2+0]) THEN flg = TRUE DIM m[4][2] m[0][0] = INT(chu[0*2+1]/30) + 2 IF m[0][1] > 12 THEN m[0][0] = m[0][0] - 12 m[0][2] = INT(saku[0*0+0]) m[0][1] = "" FOR n = 1 TO 4 IFB flg = TRUE AND n <> 1 THEN IFB INT(chu[(n-1)*2+0]) <= INT(saku[n-1]) OR INT(chu[(n-1)*2+0]) >= INT(saku[n]) THEN m[n-1][0] = m[n-2][0] m[n-1][1] = "閏" m[n-1][2] = INT(saku[n-1]) flg = FALSE ENDIF ENDIF m[n][0] = m[n-1][0] + 1 IF m[n][0] > 12 THEN m[n][0] = m[n][0] - 12 m[n][2] = INT(saku[n]) m[n][1] = "" NEXT DIM state = 0 FOR n = 0 TO 4 IFB INT(tm) < INT(m[n][2]) THEN state = 1 BREAK ELSEIF INT(tm) = INT(m[n][2]) THEN state = 2 BREAK ENDIF NEXT IF state = 0 OR state = 1 THEN n = n - 1 DIM kyureki[3] kyureki[1] = m[n][1] kyureki[2] = m[n][0] kyureki[3] = INT(tm) - INT(m[n][2]) + 1 d = JDToYMD(tm) kyureki[0] = d[0] IF kyureki[2] > 9 AND kyureki[2] > d[1] THEN kyureki[0] = kyureki[0] - 1 RESULT = SLICE(kyureki) FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // グレゴリオ暦を格納した配列(0 : 年, 1 : 月, 2 : 日, 3 : 時, 4 : 分, 5 : 秒) ////////////////////////////////////////////////// FUNCTION JDToYMD(JD) DIM x0 = INT(JD + 68570) DIM x1 = INT(x0 / 36524.25) DIM x2 = x0 - INT(36524.25 * x1 + 0.75) DIM x3 = INT((x2 + 1) / 365.2425) DIM x4 = x2 - INT(365.25 * x3) + 31 DIM x5 = INT(INT(x4) / 30.59) DIM x6 = INT(INT(x5) / 11) DIM t2 = x4 - INT(30.59 * x5) DIM t1 = x5 - 12 * x6 + 2 DIM t0 = 100 * (x1 - 49) + x3 + x6 IFB t1 = 2 AND t2 > 28 THEN SELECT TRUE CASE t0 MOD 100 = 0 AND t0 MOD 400 = 0 t2 = 29 CASE t0 MOD 4 = 0 t2 = 29 DEFAULT t2 = 28 SELEND ENDIF DIM tm = 86400 * (JD - INT(JD)) DIM t3 = INT(tm / 3600) DIM t4 = INT((tm - 3600 * t3) / 60) DIM t5 = INT(tm - 3600 * t3 - 60 * t4) DIM t[] = t0, t1, t2, t3, t4, t5 RESULT = SLICE(t) FEND ////////////////////////////////////////////////// // 【引数】 // JC : ユリウス世紀 // 【戻値】 // 月黄経(λmoon) ////////////////////////////////////////////////// FUNCTION longitudeMoon(JC) DIM A[63] = 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0004, 0.0004, 0.0005, 0.0005, 0.0005, 0.0006, 0.0006, 0.0007, 0.0007, 0.0007, 0.0007, 0.0008, 0.0009, 0.0011, 0.0012, 0.0016, 0.0018, 0.0021, 0.0021, 0.0021, 0.0022, 0.0023, 0.0024, 0.0026, 0.0027, 0.0028, 0.0037, 0.0038, 0.004, 0.004, 0.004, 0.005, 0.0052, 0.0068, 0.0079, 0.0085, 0.01, 0.0107, 0.011, 0.0125, 0.0154, 0.0304, 0.0347, 0.0409, 0.0458, 0.0533, 0.0571, 0.0588, 0.1144, 0.1851, 0.2136, 0.6583, 1.274, 6.2888, 481267.8809 * JC, 218.3162 DIM k[63] = 2322131, 4067, 549197, 1808933, 349472, 381404, 958465, 12006, 39871, 509131, 1745069, 1908795, 2258267, 111869, 27864, 485333, 405201, 790672, 1403732, 858602, 1920802, 1267871, 1856938, 401329, 341337, 71998, 990397, 818536, 922466, 99863, 1379739, 918399, 1934, 541062, 1781068, 133, 1844932, 1331734, 481266, 31932, 926533, 449334, 826671, 1431597, 1303870, 489205, 1443603, 75870, 513197.9, 445267.1, 441199.8, 854535.2, 1367733.1, 377336.3, 63863.5, 966404, 35999, 954397.7, 890534.2, 413335.3, 477198.9, 0, 0 DIM θ0[63] = 191, 70, 220, 58, 337, 354, 340, 187, 223, 242, 24, 90, 156, 38, 127, 186, 50, 114, 98, 129, 186, 249, 152, 274, 16, 85, 357, 151, 163, 122, 17, 182, 145, 259, 21, 29, 56, 283, 205, 107, 323, 188, 111, 315, 246, 142, 52, 41, 222.5, 27.9, 47.4, 148.2, 280.7, 13.2, 124.2, 276.5, 87.53, 179.93, 145.7, 10.74, 44.963, 0, 0 DIM λmoon[63] FOR n = 0 TO 60 DIM ang = normalizeAngle(k[n] * JC + θ0[n]) λmoon[n] = A[n] * COS(degToRad(ang)) NEXT λmoon[61] = normalizeAngle(A[61]) λmoon[62] = normalizeAngle(A[62]) RESULT = normalizeAngle(CALCARRAY(λmoon, CALC_ADD)) FEND ////////////////////////////////////////////////// // 【引数】 // JC : ユリウス世紀 // 【戻値】 // 太陽黄経(λsun) ////////////////////////////////////////////////// FUNCTION longitudeSun(JC) DIM A[18] = 0.0004, 0.0004, 0.0005, 0.0005, 0.0006, 0.0007, 0.0007, 0.0007, 0.0013, 0.0015, 0.0018, 0.0018, 0.0020, 0.0200, -0.0048*JC, 1.9147, 36000.7695*JC, 280.4659 DIM k[18] = 31557.0, 29930.0, 2281.0, 155.0, 33718.0, 9038.0, 3035.0, 65929.0, 22519.0, 45038.0, 445267.0, 19.0, 32964.0, 71998.1, 35999.05, 35999.05, 0, 0 DIM θ0[18] = 161.0, 48.0, 221.0, 118.0, 316.0, 64.0, 110.0, 45.0, 352.0, 254.0, 208.0, 159.0, 158.0, 265.1, 267.52, 267.52, 0, 0 DIM λsun[18] FOR n = 0 TO 15 DIM ang = normalizeAngle(k[n] * JC + θ0[n]) λsun[n] = A[n] * COS(degToRad(ang)) NEXT λsun[16] = normalizeAngle(A[16]) λsun[17] = normalizeAngle(A[17]) RESULT = normalizeAngle(CALCARRAY(λsun, CALC_ADD)) FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 指定したユリウス日の直前の二至二分の日時 ////////////////////////////////////////////////// FUNCTION nishiNibun(JD) // ユリウス日を力学時に変換 DIM TD = JD - 9/24 DIM n = 1 REPEAT // 力学時をユリウス世紀に変換 DIM JC = (TD + 0.5 - 2451545) / 36525 DIM λsun = longitudeSun(JC) IF n = 1 THEN DIM λsun0 = INT(λsun / 90) * 90 DIM Δλ = λsun - λsun0 SELECT TRUE CASE Δλ > 180 Δλ = Δλ - 360 CASE Δλ < -180 Δλ = Δλ + 180 SELEND DIM Δt = Δλ * (365.2/360) TD = TD - Δt n = n + 1 UNTIL ABS(Δt) <= 1/86400 DIM arr[1] arr[0] = TD + 9/24 arr[1] = λsun0 RESULT = SLICE(arr) FEND ////////////////////////////////////////////////// // 【引数】 // deg : 度数法 // 【戻値】 // 度単位の角度を0~360度に正規化 ////////////////////////////////////////////////// FUNCTION normalizeAngle(deg) SELECT TRUE CASE deg >= 360 deg = deg - INT(deg / 360) * 360 CASE deg < 0 deg = deg + INT(ABS(deg / 360) + 1) * 360 SELEND RESULT = deg FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 朔の日時 ////////////////////////////////////////////////// FUNCTION saku(JD) DIM lc = 1 // loop counter DIM JD1 = INT(JD) DIM JD2 = JD - JD1 JD2 = JD2 - 9/24 DIM Δt1 = 0 DIM Δt2 = 1 WHILE ABS(Δt1+Δt2) > 1/86400 DIM JC = (JD2 + 0.5) / 36525 JC = JC + (JD1 - 2451545) / 36525 DIM λsun = longitudeSun(JC) DIM λmoon = longitudeMoon(JC) DIM Δλ = λmoon - λsun SELECT TRUE CASE lc = 1 AND Δλ < 0 Δλ = normalizeAngle(Δλ) CASE λsun >= 0 AND λsun <= 20 AND λmoon >= 300 Δλ = normalizeAngle(Δλ) Δλ = 360 - Δλ CASE ABS(Δλ) > 40 Δλ = normalizeAngle(Δλ) SELEND Δt1 = INT(Δλ * 29.530589 / 360) Δt2 = Δλ * 29.530589 / 360 Δt2 = Δt2 - Δt1 JD1 = JD1 - Δt1 JD2 = JD2 - Δt2 IFB JD2 < 0 THEN JD2 = JD2 + 1 JD1 = JD1 - 1 ENDIF IFB lc = 15 AND ABS(Δt1+Δt2) > 1/86400 THEN JD1 = INT(JD - 26) JD2 = 0 ELSEIF lc > 30 AND ABS(Δt1+Δt2) > 1/86400 JD1 = JD JD2 = 0 ENDIF lc = lc + 1 WEND RESULT = JD1 + JD2 + 9/24 FEND ////////////////////////////////////////////////// // 【引数】 // 配列 : 上限値を求める配列 // 【戻値】 // 配列の上限値 ////////////////////////////////////////////////// FUNCTION UBound(array[]) RESULT = RESIZE(array) FEND ////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // hour : 時 // minute : 分 // second : 秒 // 【戻値】 // ユリウス日 ////////////////////////////////////////////////// FUNCTION YMDToJD(year, month, day, hour = 0, minute = 0, second = 0) IFB month < 3 THEN year = year - 1 month = month + 12 ENDIF DIM JD = INT(year * 365.25) JD = JD + INT(year / 400) JD = JD - INT(year / 100) JD = JD + INT((month - 2) * 30.59) JD = JD + 1721088 JD = JD + day DIM t = second / 3600 t = t + minute / 60 t = t + hour t = t / 24 JD = JD + t RESULT = JD FEND

六曜

六曜は「先勝→友引→先負→仏滅→大安→赤口」の順で繰り返すが、旧暦の毎月1日の六曜は以下のように固定されています。閏月は前月(閏月でない月)と同じになります。

六曜
1・7月 先勝せんしょう
2・8月 友引ともびき
3・9月 先負せんぶ
4・10月 仏滅ふつめつ
5・11月 大安たいあん
6・12月 赤口しゃっこう

旧暦から六曜を求める計算式

(旧暦の月 + 旧暦の日) ÷ 6の余りから求められます。余りと六曜の関係は以下のとおりです。

余り 六曜
0 大安
1 赤口
2 先勝
3 友引
4 先負
5 仏滅

プログラム実行例

六曜を求める

2020年3月14日の六曜を求める。

PRINT getRokuyo(2020, 3, 14) ////////////////////////////////////////////////// // 【引数】 // array : 要素を追加する配列(参照引数) // str : 追加する要素 // 【戻値】 // 処理後の配列の中の要素の数 ////////////////////////////////////////////////// FUNCTION arrayPush(var arr[], str) DIM res = RESIZE(arr, UBound(arr) + 1) arr[res] = str RESULT = res + 1 FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 中気と太陽黄経を格納した配列(0 : 中気, 1 : 太陽黄経) ////////////////////////////////////////////////// FUNCTION chuki(JD) JD = JD - 9/24 DIM t = (JD + 0.5 - 2451545) / 36525 DIM λsun = longitudeSun(t) DIM λsun0 = 30 * INT(λsun/30) REPEAT t = (JD + 0.5 - 2451545) / 36525 λsun = longitudeSun(t) DIM Δλ = λsun - λsun0 SELECT TRUE CASE Δλ > 180 Δλ = Δλ - 360 CASE Δλ < -180 Δλ = Δλ + 180 SELEND DIM Δt = Δλ * 365.2 / 360 JD = JD - Δt UNTIL Δt <= 1/86400 JD = JD + 9/24 DIM arr[1] = JD, λsun RESULT = SLICE(arr) FEND ////////////////////////////////////////////////// // 【引数】 // deg : 度数法 // 【戻値】 // 度数法から弧度法に変換した値 ////////////////////////////////////////////////// FUNCTION degToRad(deg) DIM pi = 3.14159265358979 RESULT = deg * pi / 180 FEND ////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // 【戻値】 // 旧暦を格納した配列(0 : 年, 1 : 月, 2 : 日) ////////////////////////////////////////////////// FUNCTION getKyureki(year, month, day) DIM tm = YMDToJD(year, month, day, 0, 0, 0) DIM chu[-1] // n*2︰ユリウス日、n*2+1︰Δλsun0、n=0,1,2 DIM tmp tmp = nishiNibun(tm) FOR n = 0 TO UBound(tmp) arrayPush(chu, tmp[n]) NEXT // 中気の計算 3回 chu[n]︰n+1回目 FOR n = 0 TO 2 tmp = chuki(chu[n*2] + 32) FOR n = 0 TO UBound(tmp) arrayPush(chu, tmp[n]) NEXT NEXT DIM saku[5] saku[0] = saku(chu[0]) // 朔の計算 FOR n = 1 TO 4 saku[n] = saku(saku[n-1] + 30) IFB ABS(INT(saku[n-1]) - INT(saku[n])) <= 26 THEN saku[n] = saku(saku[n-1] + 35) ENDIF NEXT // 朔の修正 SELECT TRUE CASE saku[1] <= INT(chu[0*2+0]) SHIFTARRAY(saku, -1) saku[4] = saku(saku[3] + 35) CASE saku[0] > INT(chu[0*2+0]) SHIFTARRAY(saku, 1) saku[0] = saku(saku[0] - 27) SELEND DIM kyureki[3] // 0︰年、1︰閏月、2︰月、3︰日 // 閏月検索 DIM flg = FALSE IF INT(saku[4]) <= INT(chu[3*2+0]) THEN flg = TRUE DIM m[4][2] m[0][0] = INT(chu[0*2+1]/30) + 2 IF m[0][1] > 12 THEN m[0][0] = m[0][0] - 12 m[0][2] = INT(saku[0*0+0]) m[0][1] = "" FOR n = 1 TO 4 IFB flg = TRUE AND n <> 1 THEN IFB INT(chu[(n-1)*2+0]) <= INT(saku[n-1]) OR INT(chu[(n-1)*2+0]) >= INT(saku[n]) THEN m[n-1][0] = m[n-2][0] m[n-1][1] = "閏" m[n-1][2] = INT(saku[n-1]) flg = FALSE ENDIF ENDIF m[n][0] = m[n-1][0] + 1 IF m[n][0] > 12 THEN m[n][0] = m[n][0] - 12 m[n][2] = INT(saku[n]) m[n][1] = "" NEXT DIM state = 0 FOR n = 0 TO 4 IFB INT(tm) < INT(m[n][2]) THEN state = 1 BREAK ELSEIF INT(tm) = INT(m[n][2]) THEN state = 2 BREAK ENDIF NEXT IF state = 0 OR state = 1 THEN n = n - 1 DIM kyureki[3] kyureki[1] = m[n][1] kyureki[2] = m[n][0] kyureki[3] = INT(tm) - INT(m[n][2]) + 1 d = JDToYMD(tm) kyureki[0] = d[0] IF kyureki[2] > 9 AND kyureki[2] > d[1] THEN kyureki[0] = kyureki[0] - 1 RESULT = SLICE(kyureki) FEND ////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // 【戻値】 // 六曜 ////////////////////////////////////////////////// FUNCTION getRokuyo(year, month, day) DIM arr[6] = "大安", "赤口", "先勝", "友引", "先負", "仏滅"" DIM d = getKyureki(year, month, day) DIM n = (d[2] + d[3]) MOD 6 RESULT = arr[n] FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // グレゴリオ暦を格納した配列(0 : 年, 1 : 月, 2 : 日, 3 : 時, 4 : 分, 5 : 秒) ////////////////////////////////////////////////// FUNCTION JDToYMD(JD) DIM x0 = INT(JD + 68570) DIM x1 = INT(x0 / 36524.25) DIM x2 = x0 - INT(36524.25 * x1 + 0.75) DIM x3 = INT((x2 + 1) / 365.2425) DIM x4 = x2 - INT(365.25 * x3) + 31 DIM x5 = INT(INT(x4) / 30.59) DIM x6 = INT(INT(x5) / 11) DIM t2 = x4 - INT(30.59 * x5) DIM t1 = x5 - 12 * x6 + 2 DIM t0 = 100 * (x1 - 49) + x3 + x6 IFB t1 = 2 AND t2 > 28 THEN SELECT TRUE CASE t0 MOD 100 = 0 AND t0 MOD 400 = 0 t2 = 29 CASE t0 MOD 4 = 0 t2 = 29 DEFAULT t2 = 28 SELEND ENDIF DIM tm = 86400 * (JD - INT(JD)) DIM t3 = INT(tm / 3600) DIM t4 = INT((tm - 3600 * t3) / 60) DIM t5 = INT(tm - 3600 * t3 - 60 * t4) DIM t[] = t0, t1, t2, t3, t4, t5 RESULT = SLICE(t) FEND ////////////////////////////////////////////////// // 【引数】 // JC : ユリウス世紀 // 【戻値】 // 月黄経(λmoon) ////////////////////////////////////////////////// FUNCTION longitudeMoon(JC) DIM A[63] = 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0003, 0.0004, 0.0004, 0.0005, 0.0005, 0.0005, 0.0006, 0.0006, 0.0007, 0.0007, 0.0007, 0.0007, 0.0008, 0.0009, 0.0011, 0.0012, 0.0016, 0.0018, 0.0021, 0.0021, 0.0021, 0.0022, 0.0023, 0.0024, 0.0026, 0.0027, 0.0028, 0.0037, 0.0038, 0.004, 0.004, 0.004, 0.005, 0.0052, 0.0068, 0.0079, 0.0085, 0.01, 0.0107, 0.011, 0.0125, 0.0154, 0.0304, 0.0347, 0.0409, 0.0458, 0.0533, 0.0571, 0.0588, 0.1144, 0.1851, 0.2136, 0.6583, 1.274, 6.2888, 481267.8809 * JC, 218.3162 DIM k[63] = 2322131, 4067, 549197, 1808933, 349472, 381404, 958465, 12006, 39871, 509131, 1745069, 1908795, 2258267, 111869, 27864, 485333, 405201, 790672, 1403732, 858602, 1920802, 1267871, 1856938, 401329, 341337, 71998, 990397, 818536, 922466, 99863, 1379739, 918399, 1934, 541062, 1781068, 133, 1844932, 1331734, 481266, 31932, 926533, 449334, 826671, 1431597, 1303870, 489205, 1443603, 75870, 513197.9, 445267.1, 441199.8, 854535.2, 1367733.1, 377336.3, 63863.5, 966404, 35999, 954397.7, 890534.2, 413335.3, 477198.9, 0, 0 DIM θ0[63] = 191, 70, 220, 58, 337, 354, 340, 187, 223, 242, 24, 90, 156, 38, 127, 186, 50, 114, 98, 129, 186, 249, 152, 274, 16, 85, 357, 151, 163, 122, 17, 182, 145, 259, 21, 29, 56, 283, 205, 107, 323, 188, 111, 315, 246, 142, 52, 41, 222.5, 27.9, 47.4, 148.2, 280.7, 13.2, 124.2, 276.5, 87.53, 179.93, 145.7, 10.74, 44.963, 0, 0 DIM λmoon[63] FOR n = 0 TO 60 DIM ang = normalizeAngle(k[n] * JC + θ0[n]) λmoon[n] = A[n] * COS(degToRad(ang)) NEXT λmoon[61] = normalizeAngle(A[61]) λmoon[62] = normalizeAngle(A[62]) RESULT = normalizeAngle(CALCARRAY(λmoon, CALC_ADD)) FEND ////////////////////////////////////////////////// // 【引数】 // JC : ユリウス世紀 // 【戻値】 // 太陽黄経(λsun) ////////////////////////////////////////////////// FUNCTION longitudeSun(JC) DIM A[18] = 0.0004, 0.0004, 0.0005, 0.0005, 0.0006, 0.0007, 0.0007, 0.0007, 0.0013, 0.0015, 0.0018, 0.0018, 0.0020, 0.0200, -0.0048*JC, 1.9147, 36000.7695*JC, 280.4659 DIM k[18] = 31557.0, 29930.0, 2281.0, 155.0, 33718.0, 9038.0, 3035.0, 65929.0, 22519.0, 45038.0, 445267.0, 19.0, 32964.0, 71998.1, 35999.05, 35999.05, 0, 0 DIM θ0[18] = 161.0, 48.0, 221.0, 118.0, 316.0, 64.0, 110.0, 45.0, 352.0, 254.0, 208.0, 159.0, 158.0, 265.1, 267.52, 267.52, 0, 0 DIM λsun[18] FOR n = 0 TO 15 DIM ang = normalizeAngle(k[n] * JC + θ0[n]) λsun[n] = A[n] * COS(degToRad(ang)) NEXT λsun[16] = normalizeAngle(A[16]) λsun[17] = normalizeAngle(A[17]) RESULT = normalizeAngle(CALCARRAY(λsun, CALC_ADD)) FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 指定したユリウス日の直前の二至二分の日時 ////////////////////////////////////////////////// FUNCTION nishiNibun(JD) // ユリウス日を力学時に変換 DIM TD = JD - 9/24 DIM n = 1 REPEAT // 力学時をユリウス世紀に変換 DIM JC = (TD + 0.5 - 2451545) / 36525 DIM λsun = longitudeSun(JC) IF n = 1 THEN DIM λsun0 = INT(λsun / 90) * 90 DIM Δλ = λsun - λsun0 SELECT TRUE CASE Δλ > 180 Δλ = Δλ - 360 CASE Δλ < -180 Δλ = Δλ + 180 SELEND DIM Δt = Δλ * (365.2/360) TD = TD - Δt n = n + 1 UNTIL ABS(Δt) <= 1/86400 DIM arr[1] arr[0] = TD + 9/24 arr[1] = λsun0 RESULT = SLICE(arr) FEND ////////////////////////////////////////////////// // 【引数】 // deg : 度数法 // 【戻値】 // 度単位の角度を0~360度に正規化 ////////////////////////////////////////////////// FUNCTION normalizeAngle(deg) SELECT TRUE CASE deg >= 360 deg = deg - INT(deg / 360) * 360 CASE deg < 0 deg = deg + INT(ABS(deg / 360) + 1) * 360 SELEND RESULT = deg FEND ////////////////////////////////////////////////// // 【引数】 // JD : ユリウス日 // 【戻値】 // 朔の日時 ////////////////////////////////////////////////// FUNCTION saku(JD) DIM lc = 1 // loop counter DIM JD1 = INT(JD) DIM JD2 = JD - JD1 JD2 = JD2 - 9/24 DIM Δt1 = 0 DIM Δt2 = 1 WHILE ABS(Δt1+Δt2) > 1/86400 DIM JC = (JD2 + 0.5) / 36525 JC = JC + (JD1 - 2451545) / 36525 DIM λsun = longitudeSun(JC) DIM λmoon = longitudeMoon(JC) DIM Δλ = λmoon - λsun SELECT TRUE CASE lc = 1 AND Δλ < 0 Δλ = normalizeAngle(Δλ) CASE λsun >= 0 AND λsun <= 20 AND λmoon >= 300 Δλ = normalizeAngle(Δλ) Δλ = 360 - Δλ CASE ABS(Δλ) > 40 Δλ = normalizeAngle(Δλ) SELEND Δt1 = INT(Δλ * 29.530589 / 360) Δt2 = Δλ * 29.530589 / 360 Δt2 = Δt2 - Δt1 JD1 = JD1 - Δt1 JD2 = JD2 - Δt2 IFB JD2 < 0 THEN JD2 = JD2 + 1 JD1 = JD1 - 1 ENDIF IFB lc = 15 AND ABS(Δt1+Δt2) > 1/86400 THEN JD1 = INT(JD - 26) JD2 = 0 ELSEIF lc > 30 AND ABS(Δt1+Δt2) > 1/86400 JD1 = JD JD2 = 0 ENDIF lc = lc + 1 WEND RESULT = JD1 + JD2 + 9/24 FEND ////////////////////////////////////////////////// // 【引数】 // 配列 : 上限値を求める配列 // 【戻値】 // 配列の上限値 ////////////////////////////////////////////////// FUNCTION UBound(array[]) RESULT = RESIZE(array) FEND ////////////////////////////////////////////////// // 【引数】 // year : 年 // month : 月 // day : 日 // hour : 時 // minute : 分 // second : 秒 // 【戻値】 // ユリウス日 ////////////////////////////////////////////////// FUNCTION YMDToJD(year, month, day, hour = 0, minute = 0, second = 0) IFB month < 3 THEN year = year - 1 month = month + 12 ENDIF DIM JD = INT(year * 365.25) JD = JD + INT(year / 400) JD = JD - INT(year / 100) JD = JD + INT((month - 2) * 30.59) JD = JD + 1721088 JD = JD + day DIM t = second / 3600 t = t + minute / 60 t = t + hour t = t / 24 JD = JD + t RESULT = JD FEND
結果
先負