//////////////////////////////////////////////////////////////////////////////////////////////////////
// This JavaScript was written by Dirk Malorny (http://www.malorny.net/)                            //
// Based on code of John Walker                                                                     //
// If you have questions or comments, please contact me via email: moonphase<at>malorny<dot>net     //
// This script may be copied if these lines remain in the script.                                   //
//////////////////////////////////////////////////////////////////////////////////////////////////////

var year, month, day, hours, minutes, seconds, today, julian, img, old, gif,
	SunDist, SunAng, MoonAge, MoonPhase, MoonDist, MoonDFrac, MoonAng, MoonPar, run,
	phases0,phases1,phases2,phases3,phases4,lastgif;

// Astronomical constants:

  epoch	      = 2444238.5;      // 1980 January 0.0

// Constants defining the Sun's apparent orbit:

  elonge      = 278.833540;     // Ecliptic longitude of the Sun at epoch 1980.0
  elongp      = 282.596403;     // Ecliptic longitude of the Sun at perigee
  eccent      = 0.016718;       // Eccentricity of Earth's orbit
  sunsmax     = 1.49598022e8;   // Semi-major axis of Earth's orbit in km
  sunangsiz   = 0.533128;       // Sun's angular size, degrees, at semi-major axis distance

// Elements of the Moon's orbit, epoch 1980.0:

  mmlong      = 64.975464;      // Moon's mean longitude at the epoch
  mmlongp     = 349.383063;     // Mean longitude of the perigee at the epoch
  mlnode      = 151.950429;     // Mean longitude of the node at the epoch
  minc        = 5.145396;       // Inclination of the Moon's orbit
  mecc        = 0.054900;       // Eccentricity of the Moon's orbit
  mangsiz     = 0.5181;         // Moon's angular size at distance a from Earth
  msmax       = 384401.0;       // Semi-major axis of Moon's orbit in km
  mparallax   = 0.9507;         // Parallax at distance a from Earth
  synmonth    = 29.53058868;    // Synodic month (new Moon to new Moon)

// Other constants:

  earthrad    = 6378.16;        // Radius of Earth in kilometres

// Useful mathematical functions:

  function abs(x)      {return ((x) < 0 ? (-(x)) : (x));}
  function torad(x)    {return ((x) * (Math.PI / 180));}
  function todeg(x)    {return ((x) * (180 / Math.PI));}
  function fixangle(x) {return (x - 360*(Math.floor(x/360)));}
  function dsin(x)     {return Math.sin(torad(x));}
  function dcos(x)     {return Math.cos(torad(x));}

// Check input fields:

  function check()
  {
    err=0;
    run=0;
    if (document.getElementById("year").value<0 || document.getElementById("year").value>2899) {error("Oooops! - Bitte Jahr korrigieren!");}
    if (document.getElementById("month").value<1 || document.getElementById("month").value>12) {error("Oooops! - Bitte Monat korrigieren!");}
    if (document.getElementById("day").value<1 || document.getElementById("day").value>31) {error("Oooops! - Bitte Tag korrigieren!");}
    if (document.getElementById("hour").value<0 || document.getElementById("hour").value>23) {error("Oooops! - Bitte Stunde korrigieren!");}
    if (document.getElementById("minute").value<0 || document.getElementById("minute").value>59) {error("Oooops! - Bitte Minute korrigieren!");}
    if (document.getElementById("second").value<0 || document.getElementById("second").value>59) {error("Oooops! - Bitte Sekunde korrigieren!");}
    if (err==0) {tick();}
  }

  function error(msg)
  {
    alert(msg);
    err=1;
  }

  function current()
  {
    run=1;
    tick();
  }

  // Convert UTC to julian date:

  function ucttoj(year, mon, mday, hour, min, sec)
  {
    var y,m,c;
    y=year;
    m=mon;
    if (m>2) m-=3;
    else {
        m+=9;
        y--;
    }
    c=y/100;
    y-=100*c;
    var jd=Math.floor((mday+(c*146097)/4+(y*1461)/4+(m*153+2)/5+1721119))-0.5;
    var jt=(sec+60*(min+60*hour))/86400;
    return jd+jt;
  }

  // Create date and time string

  function show(year, mon, mday, hour, min, sec)
  {
    var s="";
    if (mday<10) s=s+"0";
    s=s+mday+".";
    if (mon<10) s=s+"0";
    s=s+mon+"."+year+" - ";
    if (hour<10) s=s+"0";
    s=s+hour+":";
    if (min<10) s=s+"0";
    s=s+min+":";
    if (sec<10) s=s+"0";
    s=s+sec;
    return s;
  }


  function showjulian(td)
  {
    var j, ij, d, y, m, h, mi, s;
    with(Math) {
      td += 0.5;
      j = floor(td);
      j = j - 1721119.0;
      y = floor(((4 * j) - 1) / 146097.0);
      j = (j * 4.0) - (1.0 + (146097.0 * y));
      d = floor(j / 4.0);
      j = floor(((4.0 * d) + 3.0) / 1461.0);
      d = ((4.0 * d) + 3.0) - (1461.0 * j);
      d = floor((d + 4.0) / 4.0);
      m = floor(((5.0 * d) - 3) / 153.0);
      d = (5.0 * d) - (3.0 + (153.0 * m));
      d = floor((d + 5.0) / 5.0);
      y = (100.0 * y) + j;
      if (m < 10.0) {
        m = m + 3;
      } else {
        m = m - 9;
        y = y + 1;
      }
      ij = ((td - floor(td)) * 86400.0);
      h = floor(ij / 3600);
      mi = floor((ij / 60) % 60);
      s = floor(ij % 60);
    }
    return show(y,m,d,h,mi,s);
  }



// Solve the equation of Kepler:

  function kepler(m, ecc)
  {
    var e, delta;
    EPSILON = 1e-6;
    e=m=torad(m);
    do {
       delta = e - ecc * Math.sin(e) - m;
       e -= delta / (1 - ecc * Math.cos(e));
    } while (abs(delta) > EPSILON);
    return e;
  }

// Calculating the phase of the Moon:

  function phase(pdate)
  {
    var Day, N, M, Ec, Lambdasun, ml, MM, MN, Ev, Ae, A3, MmP, F,
        mEc, A4, lP, V, lPP, NP, y, x, Lambdamoon, BetaM, mage;

    with (Math) {

// Calculation of the Sun's position:

      Day = pdate - epoch;									// Date within epoch
      N = fixangle((360 / 365.2422) * Day);					// Mean anomaly of the Sun
      M = fixangle(N + elonge - elongp);					// Convert from perigee co-ordinates to epoch 1980.0
      Ec = kepler(M, eccent);								// Solve equation of Kepler
      Ec = sqrt((1+eccent) / (1-eccent)) * tan(Ec/2);
      Ec = 2 * todeg(atan(Ec));								// True anomaly
      Lambdasun = fixangle(Ec + elongp);					// Sun's geocentric ecliptic longitude
      F = ((1 + eccent*dcos(Ec)) / (1-eccent*eccent));		// Orbital distance factor
      SunDist = sunsmax / F;								// Distance to Sun in km
      SunAng = F * sunangsiz;								// Sun's angular size in degrees

// Calculation of the Moon's position:

      ml = fixangle(13.1763966 * Day + mmlong);				// Moon's mean longitude
      MM = fixangle(ml - 0.1114041 * Day - mmlongp);		// Moon's mean anomaly
      MN = fixangle(mlnode - 0.0529539 * Day);				// Moon's ascending node mean longitude
      Ev = 1.2739 * dsin(2 * (ml - Lambdasun) - MM);		// Evection
      Ae = 0.1858 * dsin(M);								// Annual equation
      A3 = 0.37 * dsin(M);									// Correction term
      MmP = MM + Ev - Ae - A3;								// Corrected anomaly
      mEc = 6.2886 * dsin(MmP);								// Correction for the equation of the centre
      A4 = 0.214 * dsin(2 * MmP);							// Another correction term
      lP = ml + Ev + mEc - Ae + A4;							// Corrected longitude
      V = 0.6583 * dsin(2 * (lP - Lambdasun));				// Variation
      lPP = lP + V;											// True longitude
      NP = MN - 0.16 * dsin(M);								// Corrected longitude of the node
      y = dsin(lPP - NP) * dcos(minc);						// Y inclination coordinate
      x = dcos(lPP - NP);									// X inclination coordinate
      Lambdamoon = todeg(atan2(y, x))+NP;					// Ecliptic longitude
      BetaM = todeg(asin(dsin(lPP-NP)*dsin(minc))); 		// Ecliptic latitude

// Calculation of the phase of the Moon:

      mage = lPP - Lambdasun;								// Age of the Moon in degrees
      MoonPhase = (1 - dcos(mage)) / 2;						// Phase of the Moon

// Calculate distance of moon from the centre of the Earth:

      MoonDist = (msmax * (1 - mecc * mecc)) / (1 + mecc * dcos(MmP + mEc));
      MoonDFrac = MoonDist / msmax;
      MoonAng = mangsiz / MoonDFrac;						// Moon's angular diameter
      MoonPar = mparallax / MoonDFrac;						// Moon's parallax
      MoonAge = synmonth*fixangle(mage)/360;				// Age of Moon in days
    }
    return fixangle(mage) / 360;							// Position in Gif's
  }


  function meanphase(sdate, k)
  {
    var t, t2, t3, nt1;
    t = (sdate - 2415020.0) / 36525;
    t2 = t * t;
    t3 = t2 * t;
    nt1 = 2415020.75933 + synmonth * k
            + 0.0001178 * t2
            - 0.000000155 * t3
            + 0.00033 * dsin(166.56 + 132.87 * t - 0.009173 * t2);
    return nt1;
  }


  function truephase(k, phase)
  {
    var t, t2, t3, pt, m, mprime, f;
    var apcor = 0;
    with (Math) {
      k += phase;
      t = k / 1236.85;
      t2 = t * t;
      t3 = t2 * t;
      pt = 2415020.75933
         + synmonth * k
         + 0.0001178 * t2
         - 0.000000155 * t3
         + 0.00033 * dsin(166.56 + 132.87 * t - 0.009173 * t2);
      m = 359.2242
        + 29.10535608 * k
        - 0.0000333 * t2
        - 0.00000347 * t3;
      mprime = 306.0253
        + 385.81691806 * k
        + 0.0107306 * t2
        + 0.00001236 * t3;
      f = 21.2964
        + 390.67050646 * k
        - 0.0016528 * t2
        - 0.00000239 * t3;
      if ((phase < 0.01) || (abs(phase - 0.5) < 0.01)) {
         pt +=     (0.1734 - 0.000393 * t) * dsin(m)
                + 0.0021 * dsin(2 * m)
                - 0.4068 * dsin(mprime)
                + 0.0161 * dsin(2 * mprime)
                - 0.0004 * dsin(3 * mprime)
                + 0.0104 * dsin(2 * f)
                - 0.0051 * dsin(m + mprime)
                - 0.0074 * dsin(m - mprime)
                + 0.0004 * dsin(2 * f + m)
                - 0.0004 * dsin(2 * f - m)
                - 0.0006 * dsin(2 * f + mprime)
                + 0.0010 * dsin(2 * f - mprime)
                + 0.0005 * dsin(m + 2 * mprime);
         apcor = 1;
      } else if ((abs(phase - 0.25) < 0.01 || (abs(phase - 0.75) < 0.01))) {
         pt +=     (0.1721 - 0.0004 * t) * dsin(m)
                + 0.0021 * dsin(2 * m)
                - 0.6280 * dsin(mprime)
                + 0.0089 * dsin(2 * mprime)
                - 0.0004 * dsin(3 * mprime)
                + 0.0079 * dsin(2 * f)
                - 0.0119 * dsin(m + mprime)
                - 0.0047 * dsin(m - mprime)
                + 0.0003 * dsin(2 * f + m)
                - 0.0004 * dsin(2 * f - m)
                - 0.0006 * dsin(2 * f + mprime)
                + 0.0021 * dsin(2 * f - mprime)
                + 0.0003 * dsin(m + 2 * mprime)
                + 0.0004 * dsin(m - 2 * mprime)
                - 0.0003 * dsin(2 * m + mprime);
         if (phase < 0.5)
           pt += 0.0028 - 0.0004 * dcos(m) + 0.0003 * dcos(mprime);
         else
           pt += -0.0028 + 0.0004 * dcos(m) - 0.0003 * dcos(mprime);
         apcor = 1;
      }
    }
    return pt;
  }


  function phasehunt(sdate)
  {
    var adate, k1, k2, nt1, nt2, j, td, yy, mm, dd, found;
    with (Math) {
      adate = sdate - 45;
      td = adate + 0.5;
      j = floor(td);
      j = j - 1721119.0;
      yy = floor(((4 * j) - 1) / 146097.0);
      j = (j * 4.0) - (1.0 + (146097.0 * yy));
      dd = floor(j / 4.0);
      j = floor(((4.0 * dd) + 3.0) / 1461.0);
      dd = ((4.0 * dd) + 3.0) - (1461.0 * j);
      dd = floor((dd + 4.0) / 4.0);
      mm = floor(((5.0 * dd) - 3) / 153.0);
      dd = (5.0 * dd) - (3.0 + (153.0 * mm));
      dd = floor((dd + 5.0) / 5.0);
      yy = (100.0 * yy) + j;
      if (mm < 10.0) {
          mm = mm + 3;
      } else {
          mm = mm - 9;
        yy = yy + 1;
      }
      k1 = floor((yy + ((mm - 1) * (1.0 / 12.0)) - 1900) * 12.3685);
      adate = nt1 = meanphase(adate, k1);
      found = 0;
      while (found==0) {
        adate += synmonth;
        k2 = k1 + 1;
        nt2 = meanphase(adate, k2);
        if (nt1 <= sdate && nt2 > sdate) {
          found = 1;
        } else {
          nt1 = nt2;
          k1 = k2;
        }
      }
      phases0 = truephase(k1, 0.0);
      phases1 = truephase(k1, 0.25);
      phases2 = truephase(k1, 0.5);
      phases3 = truephase(k1, 0.75);
      phases4 = truephase(k1+1, 0.0);
    }
  }



// We have to update the Text permanently:
  function tick()
  {

    if (run==1) {
		today=new Date();
    } else {
		today=new Date(document.getElementById("year").value,document.getElementById("month").value-1,document.getElementById("day").value,document.getElementById("hour").value,document.getElementById("minute").value,document.getElementById("second").value);
    }
    year    = today.getYear();
    if (year<1000) {year+=1900;}
    month   = today.getMonth()+1;
    day     = today.getDate();
    hours   = today.getHours();
    minutes = today.getMinutes();
    seconds = today.getSeconds();

    hours  += today.getTimezoneOffset() / 60;			// Time zone correction
    if (hours>23) {
      hours-=24;
      day++;
    }
    if (hours<0) {
      hours+=24;
      day--;
    }

    with (Math) {
      julian  = ucttoj(year,month,day,hours,minutes,seconds);
      img = round(99*phase(julian)+1);

      gif = 'moon000'+img+'.gif';
      if (img>9)  { gif = 'moon00'+img+'.gif'; }
      if (img>99) { gif = 'moon0'+img+'.gif'; }
      mad = floor(MoonAge);
      mah = floor(24*(MoonAge-mad));
      mam = floor(60*(24*(MoonAge-mad)-mah));
      mas = floor(60*(60*(24*(MoonAge-mad)-mah)-mam));

      phasehunt(julian+0.5);

      var gifStream =   "<img src='images/"+gif+"' alt='"+round(MoonPhase*10000)/100+"%' />";
      
      if (gifStream != lastgif) {
		document.getElementById("image").innerHTML=gifStream;
		lastgif=gifStream;
	  }
		
		
		document.getElementById("local_time").innerHTML=show(year,month,today.getDate(),today.getHours(),minutes,seconds);
		document.getElementById("utc").innerHTML=show(year,month,day,hours,minutes,seconds);
		document.getElementById("julian_date").innerHTML=round(julian*100000)/100000;
		document.getElementById("moon_age").innerHTML=mad+" Tage, "+mah+" Stunden, "+mam+" Minuten, "+mas+" Sekunden";
		document.getElementById("moon_phase").innerHTML=round(MoonPhase*10000)/100+"%";
		document.getElementById("moon_dist").innerHTML=round(MoonDist)+" km = "+round(MoonDist/earthrad*100)/100;
		document.getElementById("moon_diameter").innerHTML=round(MoonAng*10000)/10000+" Grad";
		document.getElementById("sun_dist").innerHTML=round(SunDist)+" km = "+round(SunDist/sunsmax*1000)/1000;
		document.getElementById("sun_diameter").innerHTML=round(SunAng*10000)/10000+" Grad";
		document.getElementById("phases0").innerHTML=showjulian(phases0);
		document.getElementById("phases1").innerHTML=showjulian(phases1);
		document.getElementById("phases2").innerHTML=showjulian(phases2);
		document.getElementById("phases3").innerHTML=showjulian(phases3);
		document.getElementById("phases4").innerHTML=showjulian(phases4);
    }

    if (run==1) window.setTimeout("tick();", 100);
  }

// main
run=1;
window.onload = tick;