moon phase

here is my easy way to calculate moon phase

int moon_phase(int y, int m, int d)
{
    /*
      calculates the moon phase (0-7), accurate to 1 segment.
      0 = > new moon.
      4 => full moon.
      */

    int c,e;
    double jd;
    int b;

    if (m < 3) {
        y--;
        m += 12;
    }
    ++m;
    c = 365.25*y;
    e = 30.6*m;
    jd = c+e+d-694039.09;  /* jd is total days elapsed */
    jd /= 29.53;           /* divide by the moon cycle (29.53 days) */
    b = jd;		   /* int(jd) -> b, take integer part of jd */
    jd -= b;		   /* subtract integer part to leave fractional part of original jd */
    b = jd*8 + 0.5;	   /* scale fraction from 0-8 and round by adding 0.5 */
    b = b & 7;		   /* 0 and 8 are the same so turn 8 into 0 */
    return b;
}
i developed this algorithm using floating point rather than integer for use in pocket calculators and only just managed to fit inside the fx-201p. an older, but integer formula is the following:

int Moon_phase(int year,int month,int day)
{
    /*k
      Calculates the moon phase (0-7), accurate to 1 segment.
      0 = > new moon.
      4 => Full moon.
    */
   
    int g, e;

    if (month == 1) --day;
    else if (month == 2) day += 30;
    else // m >= 3
    {
        day += 28 + (month-2)*3059/100;

        // adjust for leap years
        if (!(year & 3)) ++day;
        if ((year%100) == 0) --day;
    }
   
    g = (year-1900)%19 + 1;
    e = (11*g + 18) % 30;
    if ((e == 25 && g > 11) || e == 24) e++;
    return ((((e + day)*6+11)%177)/22 & 7);
}
both formulae are simplified to work from 1900 to 2199 inclusive. however, ive discovered that they disagree. consider september 23, 2002. the second formula claims this a full moon, the first does not. the moon is not full on this night so the first seems more accurate.

here is a command line pc program that is more accurate than the above for reference.

phase.exe    phase.c

and here's another program to calculate the instant of full moon

fullmoon.exe