月齢とは、地球から見て月と太陽が同じ方向になった時を0（ゼロ）として日単位で表した数字です。
月齢が、正確（ほぼ正確）に分かると二十四節気や潮の満ち引きなども求めることができます。
月齢計算は
http://park12.wakwak.com/~maki/sunmoon.htm
のプログラムを参考にさせて頂きました。
月齢計算するには修正ユリウス通日、月の視横径、太陽視横径が必要になります。
【修正ユリウス日の計算】
修正ユリウス通日とは、ユリウス通日から2400000.5を引いた数です。
ユリウス通日とは、紀元前4713年1月1日の正午を0日として数える日数のことです。
ユリウス通日や修正ユリウス通日は太陽横径や暦を求めるときによく使われます。
mjd=修正ユリウス日
求める日の月が1月、2月なら年から-1、月に+12
mjd=INT（365.25*年）+INT（年/400）-INT（年/100）+INT（30.59*月）+日-678912
MJD=修正ユリウス日+時間
MJD=mjd+時間/24+分/1440+秒/86400-0.375
【月の視横径計算】
delta=((MJD-51544.5)/365.2425+64)/86400
T=(MJD-51544.5+delte)/36525.0
DEG==1/(180/3.14159265358979)
y=218.3166+481267.8811*T-0.0015*T*T
y1=y+6.2888*COS((477198.868*T+44.963)*DEG)
y2=y1+1.274*COS((413335.35*T+10.74)*DEG)
y3=y2+0.6583*COS((890534.22*T+145.70)*DEG)
y4=y3+0.2136*COS((954397.74*T+179.93)*DEG)
y5=y4+0.1851*COS((35999.05*T+87.53)*DEG)
y6=y5+0.1144*COS((966404.0*T+276.5)*DEG)
y7=y6+0.0588*COS((63863.5*T+124.2)*DEG)
y8=y7+0.0571*COS((377336.3*T+13.2)*DEG)
y9=y8+0.0533*COS((1367733.1*T+280.7)*DEG)
y10=y9+0.0458*COS((854535.2*T+148.2)*DEG)
y11=y10+0.0409*COS((441199.8*T+47.4)*DEG)
y12=y11+0.0347*COS((445267.1*T+27.9)*DEG)
y13=y12+0.0304*COS((513197.9*T+ 222.5)*DEG)
y14=y13+0.0154*COS((75870*T+41)*DEG)
y15=y14+0.0125*COS((1443603*T+52)*DEG)
y16=y15+0.011*COS((489205*T+142)*DEG)
y17=y16+0.0107*COS((1303870*T+246)*DEG)
y18=y17+0.01*COS((1431597*T+315)*DEG)
y19=y18+0.0085*COS((826671*T+111)*DEG)
y20=y19+0.0079*COS((449334*T+188)*DEG)
y21=y20+0.0068*COS((926533*T+323)*DEG)
y22=y21+0.0052*COS((31932*T+107)*DEG)
y23=y22+0.005*COS((481266*T+205)*DEG)
y24=y23+0.004*COS((1331734*T+283)*DEG)
y25=y24+0.004*COS((1844932*T+56)*DEG)
y26=y25+0.004*COS((133*T+29)*DEG)
y27=y26+0.0038*COS((1781068*T+21)*DEG)
y28=y27+0.0037*COS((541062*T+259)*DEG)
y29=y28+0.0028*COS((1934*T+145)*DEG)
y30=y29+0.0027*COS((918399*T+182)*DEG)
y31=y30+0.0026*COS((1379739*T+17)*DEG)
y32=y31+0.0024*COS((99863*T+122)*DEG)
y33=y32+0.0023*COS((922466*T+163)*DEG)
y34=y33+0.0022*COS((818536*T+151)*DEG)
y35=y34+0.0021*COS((990397*T+357)*DEG)
y36=y35+0.0021*COS((71998*T+85)*DEG)
y37=y36+0.0021*COS((341337*T+16)*DEG)
y38=y37+0.0018*COS((401329*T+274)*DEG)
y39=y38+0.0016*COS((1856938*T+152)*DEG)
y40=y39+0.0012*COS((1267871*T+249)*DEG)
y41=y40+0.0011*COS((1920802*T+186)*DEG)
y42=y41+0.0009*COS((858602*T+129)*DEG)
y43=y42+0.0008*COS((1403732*T+98)*DEG)
y44=y43+0.0007*COS((790672*T+114)*DEG)
y45=y44+0.0007*COS((405201*T+50)*DEG)
y46=y45+0.0007*COS((485333*T+186)*DEG)
y47=y46+0.0007*COS((27864*T+127)*DEG)
y48=y47+0.0006*COS((111869*T+38)*DEG)
y49=y48+0.0006*COS((2258267*T+156)*DEG)
y50=y49+0.0005*COS((1908795*T+90)*DEG)
y51=y50+0.0005*COS((1745069*T+24 *DEG)
y52=y51+0.0005*COS((509131*T+242)*DEG)
y53=y52+0.0004*COS((39871*T+223)*DEG)
y54=y53+0.0004*COS((12006*T+187)*DEG)
Y=y54%360.0
もしYが0より小さいと(Y+360)%360
【太陽視横径計算】
Ts=(MJD-51544.5+ delta)/365.25
ys=280.4603+360.00769*Ts
ys1=ys+(1.9146-0.00005*Ts)*SIN((359.991*Ts+357.538)*DEG)
ys2=ys1+0.0200*SIN((719.981*Ts+355.05)*DEG)
ys3=ys2+0.0048*SIN((19.341*Ts+234.95)*DEG)
ys4=ys3+0.0020*SIN((329.64*Ts+247.1)*DEG)
ys5=ys4+0.0018*SIN((4452.67*Ts+297.8)*DEG)
ys6=ys5+0.0018*SIN((0.20*Ts+251.3)*DEG)
ys7=ys6+0.0015*SIN((450.37*Ts+343.2)*DEG)
ys8=ys7+0.0013*SIN((225.18*Ts+81.4)*DEG)
ys9=ys8+0.0008*SIN((659.29*Ts+132.5)*DEG)
ys10=ys9+0.0007*SIN((90.38*Ts+153.3)*DEG)
ys11=ys10+0.0007*SIN((30.35*Ts+206.8)*DEG)
ys12=ys11+0.0006*SIN((337.18*Ts+29.8)*DEG)
ys13=ys12+0.0005*SIN((1.50*Ts+207.4)*DEG)
ys14=ys13+0.0005*SIN((22.81*Ts+291.2)*DEG)
ys15=ys14+0.0004*SIN((315.56*Ts+234.9)*DEG)
ys16=ys15+0.0004*SIN((299.30*Ts+157.3)*DEG)
ys17=ys16+0.0004*SIN((720.02*Ts+21.1)*DEG)
ys18=ys17+0.0003*SIN((1079.97*Ts+352.5)*DEG)
ys19=ys18+0.0003*SIN((44.43*Ts+329.7)*DEG)
YS=ys19%360.0
もしYSが0より小さいときは(YS+360.0)%360.0
【月齢計算】
MJD=修正ユルウス通日
dif=12.1818
rdif=360/12.1818
dy=月の視横径(Y)-太陽視横径(TS)
mjdc=MJD
mjdb=mjc-(dy/dif)
月齢=MJD-mjdc
月齢が0より小さいと
(mjdb-rdif)をMJDとして月の視横径(Y)と太陽視横径(TS)も含め再計算する。
この式を19個作成してdyが0に一番近い列の月齢を月齢とする。
エクセルの図で説明します。
行2に日付と時間（○時○分○秒）を入力（リンク）します。
mjd（セルC4）の数式「=(INT(365.25*((YEAR(B2))-(IF((MONTH(B2))<3,1,0))))+INT(((YEAR(B2))-(IF((MONTH(B2))<3,1,0)))/400)-INT(((YEAR(B2))-(IF((MONTH(B2))<3,1,0)))/100)+INT(30.59*(((MONTH(B2))+(IF((MONTH(B2))<3,12,0)))-2))+(DAY(B2))-678912)」
MJD（セルC5）数式「=C4+C2/24+D2/1440+E2/86400-0.375」
dif（セルC7）「12.1818」
rdif（セルC9）「=360/12.1818」
mjdc（セルC11）数式「=C5」 セルD11の数式は「=C22」になります。
Y（セルC13」数式「=C90」
※月の視横径計算は「セルC29：C90」で上記の【月の視横径計算】計算式で計算しています。
YS（セルC15）数式「=C121」
※太陽の視横径計算は「セルC96：C121」で上記の【太陽視横径計算】計算式で計算しています。
dy（セルC17）数式「=C13-C15」
mjdb（セルC19）数式「=C11-(C17/C7)」
月齢（セルC21）数式「=C11-C19」
月齢が0（ゼロ）より小さくなるとmjdb（セルC19）からrdif（セルC9）をマイナスします。
また、月齢が29.8より大きい時はmjdbにrdifをプラスします。
（セルC22）数式「=IF(C21<0,IF(C21<0,C19-C9,C19),IF(29.8<C21,C19+C9,C19))」
この式を列CからUまで19個作成します。
この19の式の中でdy（行17）の数字が一番0（ゼロ）に近い列の月齢を求める日時の月齢とします。
http://park12.wakwak.com/~maki/sunmoon.htm
のサイトの計算と少し違っていますが、ほぼ近い値が出ているので良しとしました。
完成品はこちらでダウンロードできます。【月齢旧暦計算】（zip）クリックするとダウンロードされます。
