PHP has easter_date but I was looking for more: some code that allows to evaluate the next full- or newmoon given ANY date. The algorithm was taken from Meeus chapter 49 which is an old and very famous book. The code was on the web in delphi (www.hoerstemeier.com) and I just wrote it in php.
Enter the date as jd-format (juliantojd for example) and you get the next 'Fullmoon' or 'Newmoon' date again as jd-format. Here comes the code (first three functions necessary for evaluation):
//necessary for phase evaluation
function sin_d($x){
return sin(M_PI*($x-round($x/360)*360)/180);
}
//necessary for phase evaluation
function cos_d($x){
return cos(M_Pi*($x-round($x/360)*360)/180);
}
// based on Meeus chapter 49 (Delphi code taken from [url]www.Hoerstemeier.com[/url])
function nextphase_49($date,$phase){
//first: change date from jd to delphi
//delphi starts at 12/31/1899, php at 1/1/1
$date=$date-2415032;
if ($date<0){ //there is no date 0 in delphi
$date=$date+1;
}
// Hoerstemeier used
//kk:=int(k)+ord(phase)/8;
// with
// TMoonPhase=(Newmoon,WaxingCrescrent,FirstQuarter,WaxingGibbous,
// Fullmoon,WaningGibbous,LastQuarter,WaningCrescent);
$k=round(($date-36526)/36525*1236.85);
if ($phase=='Fullmoon') {
$kk=intval($k)+0.5;
}
elseif ($phase=='Newmoon'){
$kk=intval($k);
}
else
{
die('Phase cannot be evaluated');
}
$ts=($date-36526)/36525;
$t=$kk/1236.85;
$jde=2451550.09765+29.530588853*$kk
+$t*$t*(0.0001337-$t*(0.000000150-0.00000000073*$t));
$m=2.5534+29.10535669*$kk-$t*$t*(0.0000218+0.00000011*$t);
$ms=201.5643+385.81693528*$kk+$t*$t*(0.1017438+$t*(0.00001239-$t*0.000000058));
$f= 160.7108+390.67050274*$kk-$t*$t*(0.0016341+$t*(0.00000227-$t*0.000000011));
$o=124.7746-1.56375580*$kk+$t*$t*(0.0020691+$t*0.00000215);
$e=1-$ts*(0.002516+$ts*0.0000074);
// several corrections
if ($phase=='Newmoon'){
$korr= -0.40720*sin_d($ms)+0.17241*$e*sin_d($m)+0.01608*sin_d(2*$ms)
+0.01039*sin_d(2*$f)+0.00739*$e*sin_d($ms-$m)-0.00514*$e*sin_d($ms+$m)
+0.00208*$e*$e*sin_d(2*$m)-0.00111*sin_d($ms-2*$f)-0.00057*sin_d($ms+2*$f)
+0.00056*$e*sin_d(2*$ms+$m)-0.00042*sin_d(3*$ms)+0.00042*$e*sin_d($m+2*$f)
+0.00038*$e*sin_d($m-2*$f)-0.00024*$e*sin_d(2*$ms-$m)-0.00017*sin_d($o)
-0.00007*sin_d($ms+2*$m)+0.00004*sin_d(2*$ms-2*$f)+0.00004*sin_d(3*$m)
+0.00003*sin_d($ms+$m-2*$f)+0.00003*sin_d(2*$ms+2*$f)-0.00003*sin_d($ms+$m+2*$f)
+0.00003*sin_d($ms-$m+2*$f)-0.00002*sin_d($ms-$m-2*$f)-0.00002*sin_d(3*$ms+$m)
+0.00002*sin_d(4*$ms);
}
else {
//must be Fullmoon
$korr= -0.40614*sin_d($ms)+0.17302*$e*sin_d($m)+0.01614*sin_d(2*$ms)+0.01043*sin_d(2*$f)
+0.00734*$e*sin_d($ms-$m)-0.00515*$e*sin_d($ms+$m)+0.00209*$e*e*sin_d(2*$m)
-0.00111*sin_d($ms-2*$f)-0.00057*sin_d($ms+2*$f)+0.00056*$e*sin_d(2*$ms+$m)
-0.00042*sin_d(3*$ms)+0.00042*$e*sin_d($m+2*$f)+0.00038*$e*sin_d($m-2*$f)
-0.00024*$e*sin_d(2*$ms-$m)-0.00017*sin_d($o)-0.00007*sin_d($ms+2*$m)
+0.00004*sin_d(2*$ms-2*$f)+0.00004*sin_d(3*$m)+0.00003*sin_d($ms+$m-2*$f)
+0.00003*sin_d(2*$ms+2*$f)-0.00003*sin_d($ms+$m+2*$f)+0.00003*sin_d($ms-$m+2*$f)
-0.00002*sin_d($ms-$m-2*$f)-0.00002*sin_d(3*$ms+$m)+0.00002*sin_d(4*$ms);
}
// Additional Corrections due to planets
$korr=$korr+0.000325*sin_d(299.77+0.107408*$kk-0.009173*$t*$t)
+0.000165*sin_d(251.88+0.016321*$kk)+0.000164*sin_d(251.83+26.651886*$kk)
+0.000126*sin_d(349.42+36.412478*$kk)+0.000110*sin_d(84.66+18.206239*$kk)
+0.000062*sin_d(141.74+53.303771*$kk)+0.000060*sin_d(207.14+2.453732*$kk)
+0.000056*sin_d(154.84+7.306860*$kk)+0.000047*sin_d(34.52+27.261239*$kk)
+0.000042*sin_d(207.19+0.121824*$kk)+0.000040*sin_d(291.34+1.844379*$kk)
+0.000037*sin_d(161.72+24.198154*$kk)+0.000035*sin_d(239.56+25.513099*$kk)
+0.000023*sin_d(331.55+3.592518*$kk);
$r=$jde+$korr-2415018.5; //result as Delphi-date
//change delphi-date to jd-date
if ($r<0){ //there was no date 0 in delphi
$r=$r-1;
}
$r=$r+2415032;
return $r;
}
// evaluates next 'Newmoon' or 'Fullmoon',
// date given as jd-date
function jd_nextmoon($date,$phase){
$temp_date=$date-28;
$r=$temp_date;
while ($r<$date){
$r=nextphase_49($temp_date,$phase);
$temp_date=$temp_date+28;
}
return $r;
}
