// JavaScript by Peter Hayes http://www.peter-hayes.freeserve.co.uk/ // Copyright 2001-2003 // This code is made freely available but please keep this notice. // The calculations are approximate but should be good enough for general use, // I accept no responsibility for errors in astronomy or coding. // WARNING moonrise code changed on 6 May 2003 to correct a systematic error // these are now local times NOT UTC as the original code did. // Meeus first edition table 45.A Longitude and distance of the moon var T45AD = new Array(0, 2, 2, 0, 0, 0, 2, 2, 2, 2, 0, 1, 0, 2, 0, 0, 4, 0, 4, 2, 2, 1, 1, 2, 2, 4, 2, 0, 2, 2, 1, 2, 0, 0, 2, 2, 2, 4, 0, 3, 2, 4, 0, 2, 2, 2, 4, 0, 4, 1, 2, 0, 1, 3, 4, 2, 0, 1, 2, 2); var T45AM = new Array(0, 0, 0, 0, 1, 0, 0, -1, 0, -1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, -1, 0, 0, 0, 1, 0, -1, 0, -2, 1, 2, -2, 0, 0, -1, 0, 0, 1, -1, 2, 2, 1, -1, 0, 0, -1, 0, 1, 0, 1, 0, 0, -1, 2, 1, 0, 0); var T45AMP = new Array( 1, -1, 0, 2, 0, 0, -2, -1, 1, 0, -1, 0, 1, 0, 1, 1, -1, 3, -2, -1, 0, -1, 0, 1, 2, 0, -3, -2, -1, -2, 1, 0, 2, 0, -1, 1, 0, -1, 2, -1, 1, -2, -1, -1, -2, 0, 1, 4, 0, -2, 0, 2, 1, -2, -3, 2, 1, -1, 3, -1); var T45AF = new Array( 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, 0, -2, 2, -2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0, 0, 0, 0, -2, 2, 0, 2, 0, 0, 0, 0, 0, 0, -2, 0, 0, 0, 0, -2, -2, 0, 0, 0, 0, 0, 0, 0, -2); var T45AL = new Array(6288774, 1274027, 658314, 213618, -185116, -114332, 58793, 57066, 53322, 45758, -40923, -34720, -30383, 15327, -12528, 10980, 10675, 10034, 8548, -7888, -6766, -5163, 4987, 4036, 3994, 3861, 3665, -2689, -2602, 2390, -2348, 2236, -2120, -2069, 2048, -1773, -1595, 1215, -1110, -892, -810, 759, -713, -700, 691, 596, 549, 537, 520, -487, -399, -381, 351, -340, 330, 327, -323, 299, 294, 0); var T45AR = new Array(-20905355, -3699111, -2955968, -569925, 48888, -3149, 246158, -152138, -170733, -204586, -129620, 108743, 104755, 10321, 0, 79661, -34782, -23210, -21636, 24208, 30824, -8379, -16675, -12831, -10445, -11650, 14403, -7003, 0, 10056, 6322, -9884, 5751, 0, -4950, 4130, 0, -3958, 0, 3258, 2616, -1897, -2117, 2354, 0, 0, -1423, -1117, -1571, -1739, 0, -4421, 0, 0, 0, 0, 1165, 0, 0, 8752); // Meeus table 45B latitude of the moon var T45BD = new Array(0, 0, 0, 2, 2, 2, 2, 0, 2, 0, 2, 2, 2, 2, 2, 2, 2, 0, 4, 0, 0, 0, 1, 0, 0, 0, 1, 0, 4, 4, 0, 4, 2, 2, 2, 2, 0, 2, 2, 2, 2, 4, 2, 2, 0, 2, 1, 1, 0, 2, 1, 2, 0, 4, 4, 1, 4, 1, 4, 2); var T45BM = new Array( 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 1, -1, -1, -1, 1, 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 0, 1, 1, 0, -1, -2, 0, 1, 1, 1, 1, 1, 0, -1, 1, 0, -1, 0, 0, 0, -1, -2); var T45BMP = new Array(0, 1, 1, 0, -1, -1, 0, 2, 1, 2, 0, -2, 1, 0, -1, 0, -1, -1, -1, 0, 0, -1, 0, 1, 1, 0, 0, 3, 0, -1, 1, -2, 0, 2, 1, -2, 3, 2, -3, -1, 0, 0, 1, 0, 1, 1, 0, 0, -2, -1, 1, -2, 2, -2, -1, 1, 1, -1, 0, 0); var T45BF = new Array( 1, 1, -1, -1, 1, -1, 1, 1, -1, -1, -1, -1, 1, -1, 1, 1, -1, -1, -1, 1, 3, 1, 1, 1, -1, -1, -1, 1, -1, 1, -3, 1, -3, -1, -1, 1, -1, 1, -1, 1, 1, 1, 1, -1, 3, -1, -1, 1, -1, -1, 1, -1, 1, -1, -1, -1, -1, -1, -1, 1); var T45BL = new Array(5128122, 280602, 277693, 173237, 55413, 46271, 32573, 17198, 9266, 8822, 8216, 4324, 4200, -3359, 2463, 2211, 2065, -1870, 1828, -1794, -1749, -1565, -1491, -1475, -1410, -1344, -1335, 1107, 1021, 833, 777, 671, 607, 596, 491, -451, 439, 422, 421, -366, -351, 331, 315, 302, -283, -229, 223, 223, -220, -220, -185, 181, -177, 176, 166, -164, 132, -119, 115, 107); // MoonPos calculates the Moon position, based on Meeus chapter 45 function MoonPos(year,month,day,hours) { // julian date var jd=julian(year,month,day,hours); var T=(jd-2451545.0)/36525; var T2=T*T; var T3=T2*T; var T4=T3*T; // Moons mean longitude L' var LP=218.3164477+481267.88123421*T-0.0015786*T2+T3/538841.0-T4/65194000.0; // Moons mean elongation var D=297.8501921+445267.1114034*T-0.0018819*T2+T3/545868.0-T4/113065000.0; // Suns mean anomaly var M=357.5291092+35999.0502909*T-0.0001536*T2+T3/24490000.0; // Moons mean anomaly M' var MP=134.9633964+477198.8675055*T+0.0087414*T2+T3/69699.0-T4/14712000.0; // Moons argument of latitude var F=93.2720950+483202.0175233*T-0.0036539*T2-T3/3526000.0+T4/863310000.0; // Additional arguments var A1=119.75+131.849*T; var A2=53.09+479264.290*T; var A3=313.45+481266.484*T; var E=1-0.002516*T-0.0000074*T2; var E2=E*E; // Sums of periodic terms from table 45.A and 45.B var Sl=0.0; var Sr=0.0; for (var i=0; i<60; i++) { var Eterm=1; if (Math.abs(T45AM[i])==1) Eterm=E; if (Math.abs(T45AM[i])==2) Eterm=E2; Sl+=T45AL[i]*Eterm*sind(rev(T45AD[i]*D+T45AM[i]*M+T45AMP[i]*MP+T45AF[i]*F)); Sr+=T45AR[i]*Eterm*cosd(rev(T45AD[i]*D+T45AM[i]*M+T45AMP[i]*MP+T45AF[i]*F)); } var Sb=0.0; for (var i=0; i<60; i++) { var Eterm=1; if (Math.abs(T45BM[i])==1) Eterm=E; if (Math.abs(T45BM[i])==2) Eterm=E2; Sb+=T45BL[i]*Eterm*sind(rev(T45BD[i]*D+T45BM[i]*M+T45BMP[i]*MP+T45BF[i]*F)); } // Additional additive terms Sl=Sl+3958*sind(rev(A1))+1962*sind(rev(LP-F))+318*sind(rev(A2)); Sb=Sb-2235*sind(rev(LP))+382*sind(rev(A3))+175*sind(rev(A1-F))+ 175*sind(rev(A1+F))+127*sind(rev(LP-MP))-115*sind(rev(LP+MP)); // geocentric longitude, latitude and distance var mglong=rev(LP+Sl/1000000.0); var mglat=rev(Sb/1000000.0); if (mglat > 180.0) mglat=mglat-360; var mr=Math.round(385000.56+Sr/1000.0); // Obliquity of Ecliptic var obl=23.4393-3.563E-9*(jd-2451543.5); // RA and dec var ra=rev(atan2d(sind(mglong)*cosd(obl)-tand(mglat)*sind(obl), cosd(mglong)))/15.0; var dec=rev(asind(sind(mglat)*cosd(obl)+cosd(mglat)*sind(obl)*sind(mglong))); if (dec > 180.0) dec=dec-360; return new Array(ra,dec,mr); } function MoonRise(year,month,day,TZ,latitude,longitude) { // returns an array containing rise and set times or one of the // following codes. // -1 rise or set event not found and moon was down at 00:00 // -2 rise or set event not found and moon was up at 00:00 // WARNING code changes on 6/7 May 2003 these are now local times // NOT UTC and rise/set not found codes changed. var hours=0; var riseset=new Array(); // elh is the elevation at the hour elhdone is true if elh calculated var elh=new Array(); var elhdone=new Array(); for (var i=0; i<=24; i++) {elhdone[i]=false;} // Compute the moon elevation at start and end of day // store elevation at the hours in an array elh to save search time var rad=MoonPos(year,month,day,hours-TZ); var altaz=m_radtoaa(rad[0],rad[1],year,month,day,hours-TZ,latitude,longitude); elh[0]=altaz[0]; elhdone[0]=true; // set the return code to allow for always up or never rises if (elh[0] > 0.0) { riseset=new Array(-2,-2); } else { riseset=new Array(-1,-1); } hours=24; rad=MoonPos(year,month,day,hours-TZ); altaz=m_radtoaa(rad[0],rad[1],year,month,day,hours-TZ,latitude,longitude); elh[24]=altaz[0]; elhdone[24]=true; // search for moonrise and set for (var rise=0; rise<2; rise++) { var found=false; var hfirst=0; var hlast=24; // Try a binary chop on the hours to speed the search while (Math.ceil((hlast-hfirst)/2) > 1) { hmid=hfirst+Math.round((hlast-hfirst)/2); if (!elhdone[hmid]) { hours=hmid; rad=MoonPos(year,month,day,hours-TZ); altaz=m_radtoaa(rad[0],rad[1],year,month,day,hours-TZ,latitude,longitude); elh[hmid]=altaz[0]; elhdone[hmid]=true; } if (((rise == 0) && (elh[hfirst] <= 0.0) && (elh[hmid] >= 0.0)) || ((rise == 1) && (elh[hfirst] >= 0.0) && (elh[hmid] <= 0.0))) { hlast=hmid; found=true; continue; } if (((rise == 0) && (elh[hmid] <= 0.0) && (elh[hlast] >= 0.0)) || ((rise == 1) && (elh[hmid] >= 0.0) && (elh[hlast] <= 0.0))) { hfirst=hmid; found=true; continue; } break; } // If the binary chop did not find a 1 hour interval if ((hlast-hfirst) > 1) { for (var i=hfirst; i= 0.0)) || ((rise == 1) && (elh[i] >= 0.0) && (elh[i+1] <= 0.0))) { hfirst=i; hlast=i+1; found=true; break; } } } // simple linear interpolation for the minutes if (found) { var elfirst=elh[hfirst]; var ellast=elh[hlast]; hours=hfirst+0.5; rad=MoonPos(year,month,day,hours-TZ); altaz=m_radtoaa(rad[0],rad[1],year,month,day,hours-TZ,latitude,longitude); // alert("day ="+day+" hour ="+hours+" altaz="+altaz[0]+" "+altaz[1]); if ((rise == 0) && (altaz[0] <= 0.0)) {hfirst=hours; elfirst=altaz[0];} if ((rise == 0) && (altaz[0] > 0.0)) {hlast=hours; ellast=altaz[0];} if ((rise == 1) && (altaz[0] <= 0.0)) {hlast=hours; ellast=altaz[0];} if ((rise == 1) && (altaz[0] > 0.0)) {hfirst=hours; elfirst=altaz[0];} var eld=Math.abs(elfirst)+Math.abs(ellast); riseset[rise]=hfirst+(hlast-hfirst)*Math.abs(elfirst)/eld; } } // End of rise/set loop return(riseset); } function m_radtoaa(ra,dec,year,month,day,hours,lat,lon) { // convert ra and dec to altitude and azimuth // year, month, day and hours are the Greenwich date and time // lat and lon are the observers latitude and longitude var lst=m_local_sidereal(year,month,day,hours,lon); var x=cosd(15.0*(lst-ra))*cosd(dec); var y=sind(15.0*(lst-ra))*cosd(dec); var z=sind(dec); // rotate so z is the local zenith var xhor=x*sind(lat)-z*cosd(lat); var yhor=y; var zhor=x*cosd(lat)+z*sind(lat); var azimuth=rev(atan2d(yhor,xhor)+180.0); // so 0 degrees is north var altitude=atan2d(zhor,Math.sqrt(xhor*xhor+yhor*yhor)); return new Array(altitude,azimuth); } function m_local_sidereal(year,month,day,hours,lon) { // Compute local siderial time in degrees // year, month, day and hours are the Greenwich date and time // lon is the observers longitude var d=m_dayno(year,month,day,hours); var lst=(98.9818+0.985647352*d+hours*15+lon); return rev(lst)/15; } function m_dayno(year,month,day,hours) { // Day number is a modified Julian date, day 0 is 2000 January 0.0 // which corresponds to a Julian date of 2451543.5 var d= 367*year-Math.floor(7*(year+Math.floor((month+9)/12))/4) +Math.floor((275*month)/9)+day-730530+hours/24; return d; } function MoonPhase(year,month,day,hours) { // the illuminated percentage from Meeus chapter 46 var j=dayno(year,month,day,hours)+2451543.5; var T=(j-2451545.0)/36525; var T2=T*T; var T3=T2*T; var T4=T3*T; // Moons mean elongation Meeus first edition // var D=297.8502042+445267.1115168*T-0.0016300*T2+T3/545868.0-T4/113065000.0; // Moons mean elongation Meeus second edition var D=297.8501921+445267.1114034*T-0.0018819*T2+T3/545868.0-T4/113065000.0; // Moons mean anomaly M' Meeus first edition // var MP=134.9634114+477198.8676313*T+0.0089970*T2+T3/69699.0-T4/14712000.0; // Moons mean anomaly M' Meeus second edition var MP=134.9633964+477198.8675055*T+0.0087414*T2+T3/69699.0-T4/14712000.0; // Suns mean anomaly var M=357.5291092+35999.0502909*T-0.0001536*T2+T3/24490000.0; // phase angle var pa=180.0-D-6.289*sind(MP)+2.1*sind(M)-1.274*sind(2*D-MP) -0.658*sind(2*D)-0.214*sind(2*MP)-0.11*sind(D); return(rev(pa)); } function MoonPhase2(jtime) { // the illuminated percentage from Meeus chapter 46 var j=jtime; var T=(j-2451545.0)/36525; var T2=T*T; var T3=T2*T; var T4=T3*T; // Moons mean elongation Meeus first edition // var D=297.8502042+445267.1115168*T-0.0016300*T2+T3/545868.0-T4/113065000.0; // Moons mean elongation Meeus second edition var D=297.8501921+445267.1114034*T-0.0018819*T2+T3/545868.0-T4/113065000.0; // Moons mean anomaly M' Meeus first edition // var MP=134.9634114+477198.8676313*T+0.0089970*T2+T3/69699.0-T4/14712000.0; // Moons mean anomaly M' Meeus second edition var MP=134.9633964+477198.8675055*T+0.0087414*T2+T3/69699.0-T4/14712000.0; // Suns mean anomaly var M=357.5291092+35999.0502909*T-0.0001536*T2+T3/24490000.0; // phase angle var pa=180.0-D-6.289*sind(MP)+2.1*sind(M)-1.274*sind(2*D-MP) -0.658*sind(2*D)-0.214*sind(2*MP)-0.11*sind(D); return(rev(pa)); } function MoonQrts(year,month) { var jstart=julian(year,month,1,0); var weekday=new Array(7); weekday[1]="Monday" ; weekday[2]="Tuesday" ; weekday[3]="Wednesday" ; weekday[4]="Thursday" ; weekday[5]="Friday" ; weekday[6]="Saturday" ; weekday[0]="Sunday" ; if ((month==4) || (month==6) || (month==9) || (month==11)) { var days=30; } else if (month==2) { var days=28; } else { var days=31; } var jstop=julian(year,month,days,24); var hours=Math.round((jstop-jstart)*24); var ev=0; var events=new Array(5); events[4]=new Array(0,0,0,0,0); // For each hour during the month, get the moon phase for (count = 0; count < hours; count++) { var time=jstart+(count/24); var phase=Math.abs(MoonPhase2(time)-360).toFixed(2); // Check if phase is Full Moon if (phase>=0 && phase<1) { var time2=jstart+((count-1)/24); var phase2=Math.abs(MoonPhase2(time2)-360).toFixed(2); var phase=(phase*1.0)+(360.0*1.0); if (phase2<300) phase2=(phase2*1.0)+(360.0*1.0); // This linear interpolation estimates the event time to the nearest minute var prectime=((360-phase2)*(time-time2)/(phase-phase2))+time2; var result=jdtocd(prectime); if (ev==0) { events[0]=new Array(0,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } else if (result[2]!=events[ev-1][2]) { events[ev]=new Array(0,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } } // Check if phase is Last Quarter if ((phase-90)>=0 && (phase-90)<1) { var time2=jstart+((count-1)/24); var phase2=Math.abs(MoonPhase2(time2)-360).toFixed(2); var prectime=((90-phase2)*(time-time2)/(phase-phase2))+time2; var result=jdtocd(prectime); if (ev==0) { events[0]=new Array(90,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } else if (result[2]!=events[ev-1][2]) { events[ev]=new Array(90,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } } // Check if phase is New Moon if ((phase-180>=0) && (phase-180)<1) { var time2=jstart+((count-1)/24); var phase2=Math.abs(MoonPhase2(time2)-360).toFixed(2); var prectime=((180-phase2)*(time-time2)/(phase-phase2))+time2; var result=jdtocd(prectime); if (ev==0) { events[0]=new Array(180,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } else if (result[2]!=events[ev-1][2]) { events[ev]=new Array(180,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } } // Check if phase is First Quarter if ((phase-270>=0) && (phase-270)<1) { var time2=jstart+((count-1)/24); var phase2=Math.abs(MoonPhase2(time2)-360).toFixed(2); var prectime=((270-phase2)*(time-time2)/(phase-phase2))+time2; var result=jdtocd(prectime); if (ev==0) { events[0]=new Array(270,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } else if (result[2]!=events[ev-1][2]) { events[ev]=new Array(270,weekday[result[3]],result[2],daytag(result[2]),(result[4]+(result[5]/60))); ev++; } } } return new Array(events,ev); }