[net.astro] shoot the moon

mom@well.UUCP (Gremlin) (03/05/86)

I am not resposible for what system's these work on
-------------------CUT---------------------------------------------------------
#! /bin/sh
# This is a shell archive, meaning:
# 1. Remove everything above the #! /bin/sh line.
# 2. Save the resulting text in a file.
# 3. Execute the file with /bin/sh (not csh) to create the files:
#	Makefile
#	moon1.c
#	moon2.c
# This archive created: Tue Mar  4 13:45:23 1986
# By:	Gremlin (After Midnight Lunch Bunch)
export PATH; PATH=/bin:$PATH
if test -f 'Makefile'
then
	echo shar: will not over-write existing file "'Makefile'"
else
cat << \SHAR_EOF > 'Makefile'

moon:
	cc -O -o moon1 moon1.c -lm
	cc -O -o moon2 moon2.c -lm
SHAR_EOF
fi # end of overwriting check
if test -f 'moon1.c'
then
	echo shar: will not over-write existing file "'moon1.c'"
else
cat << \SHAR_EOF > 'moon1.c'
/****************************************************************************
 pom.c

     Phase of the Moon. Calculates the current phase of the moon.
     Based on routines from `Practical Astronomy with Your Calculator',
        by Duffett-Smith.
     Comments give the section from the book that particular piece
        of code was adapted from.

     -- Keith E. Brandt  VIII 1984

 ****************************************************************************/

#include <stdio.h>
#include <sys/time.h>
#include <math.h>
#define PI         3.141592654
#define EPOCH   1983
#define EPSILONg 279.103035     /* solar ecliptic long at EPOCH */
#define RHOg     282.648015     /* solar ecliptic long of perigee at EPOCH */
#define e          0.01671626   /* solar orbit eccentricity */
#define lzero    106.306091     /* lunar mean long at EPOCH */
#define Pzero    111.481526     /* lunar mean long of perigee at EPOCH */
#define Nzero     93.913033     /* lunar mean long of node at EPOCH */

main()  {

double dtor();
double adj360();
double potm();

long *lo = (long *) calloc (1, sizeof(long)); /* used by time calls */
struct tm *pt; /* ptr to time structure */

double days;   /* days since EPOCH */
double phase;  /* percent of lunar surface illuminated */
double phase2; /* percent of lunar surface illuminated one day later */
int i = EPOCH;

time (lo);  /* get system time */
pt = gmtime(lo);  /* get ptr to gmt time struct */
cfree(lo);

/* calculate days since EPOCH */
days = (pt->tm_yday +1) + ((pt->tm_hour + (pt->tm_min / 60.0)
       + (pt->tm_sec / 3600.0)) / 24.0);
while (i < pt->tm_year + 1900)
   days = days + 365 + ly(i++);

phase = potm(days);
printf("The Moon is ");
if ((int)(phase + .5) == 100) {
   printf("Full\n");
   }
else if ((int)(phase + 0.5) == 0) 
   printf("New\n");
else if ((int)(phase + 0.5) == 50)  {
   phase2 = potm(++days);
   if (phase2 > phase)
      printf("at the First Quarter\n");
   else 
      printf("at the Last Quarter\n");
   }
else if ((int)(phase + 0.5) > 50) {
   phase2 = potm(++days);
   if (phase2 > phase)
      printf("Waxing ");
   else 
      printf("Waning ");
   printf("Gibbous (%1.0f%% of Full)\n", phase);
   }
else if ((int)(phase + 0.5) < 50) {
   phase2 = potm(++days);
   if (phase2 > phase)
      printf("Waxing ");
   else
      printf("Waning ");
   printf("Crescent (%1.0f%% of Full)\n", phase);
   }
}

double potm(days)
double days;
{
double N;
double Msol;
double Ec;
double LambdaSol;
double l;
double Mm;
double Ev;
double Ac;
double A3;
double Mmprime;
double A4;
double lprime;
double V;
double ldprime;
double D;
double Nm;

N = 360 * days / 365.2422;  /* sec 42 #3 */
adj360(&N);

Msol = N + EPSILONg - RHOg; /* sec 42 #4 */
adj360(&Msol);

Ec = 360 / PI * e * sin(dtor(Msol)); /* sec 42 #5 */

LambdaSol = N + Ec + EPSILONg;       /* sec 42 #6 */
adj360(&LambdaSol);

l = 13.1763966 * days + lzero;       /* sec 61 #4 */
adj360(&l);

Mm = l - (0.1114041 * days) - Pzero; /* sec 61 #5 */
adj360(&Mm);

Nm = Nzero - (0.0529539 * days);     /* sec 61 #6 */
adj360(&Nm);

Ev = 1.2739 * sin(dtor(2*(l - LambdaSol) - Mm)); /* sec 61 #7 */

Ac = 0.1858 * sin(dtor(Msol));       /* sec 61 #8 */
A3 = 0.37 * sin(dtor(Msol));

Mmprime = Mm + Ev - Ac - A3;         /* sec 61 #9 */

Ec = 6.2886 * sin(dtor(Mmprime));    /* sec 61 #10 */

A4 = 0.214 * sin(dtor(2 * Mmprime)); /* sec 61 #11 */

lprime = l + Ev + Ec - Ac + A4;      /* sec 61 #12 */

V = 0.6583 * sin(dtor(2 * (lprime - LambdaSol))); /* sec 61 #13 */

ldprime = lprime + V;                /* sec 61 #14 */

D = ldprime - LambdaSol;             /* sec 63 #2 */

return (50 * (1 - cos(dtor(D))));    /* sec 63 #3 */
}

ly(yr)
int yr;
{
/* returns 1 if leapyear, 0 otherwise */
return (yr % 4 == 0 && yr % 100 != 0 || yr % 400 == 0);
}

double dtor(deg)
double deg;
{
/* convert degrees to radians */
return (deg * PI / 180);
}

double adj360(deg)
double *deg;
{
/* adjust value so 0 <= deg <= 360 */
do if (*deg < 0)
   *deg += 360;
else if (*deg > 360)
   *deg -= 360;
while (*deg < 0 || *deg > 360);
}


SHAR_EOF
fi # end of overwriting check
if test -f 'moon2.c'
then
	echo shar: will not over-write existing file "'moon2.c'"
else
cat << \SHAR_EOF > 'moon2.c'
/*$The phase of the moon, for your safety and convenience.
**
**  Stolen from ArchMach & converted from PL/I by Brian Hess.
**  Extensively cleaned up by Rich $alz.
**
**  If you can figure out how this program works, then YOU deserve
**  to be working on it, sucker!  Here's a hint:  The epoch is 13:23,
**  10/1/78.
*/


#include <math.h>
#include <sys/types.h>
#include <time.h>


/* Globals. */
long		 Day;
long		 Hour;
long		 Minute;
double		 Fraction;
char		*Moon[] =
{
    "new",
    "first quarter of the",
    "full moon",
    "last quarter of the"
};


/* Linked in later. */
char		*rindex();
time_t		 time();
struct tm	*localtime();


#define LINES		24
#define WIDTH		80
#define CENTER		((WIDTH - 2 * LINES) / 2)
#define BRIGHT		'@'
#define LEDGE		'('
#define REDGE		')'
#define FULL		0.5
#define TwoPi		(2 * 3.14159)
#define ZERO		0.03

#define Plural(X)	if (X != 1) printf("s");
/*!*/
long
Calculate()
{
    register struct tm	*tm;
    register long	 Length;
    register long	 Phase;
    register long	 DeltaH;
    register long	 Delta;
    register long	 offset;
    time_t		 tick;
    long		 julian;
    long		 year;
    long		 hour;
    long		 minute;

    tick	= time((time_t *)0);
    tm		= localtime(&tick);
    julian	= tm->tm_yday + 1;
    year	= tm->tm_year - 78;
    hour	= tm->tm_hour;
    minute	= tm->tm_min;

    Length	= (double)2551 / 60 * 1000 + (double)443 / 60;
    offset	= ((year * 365L + julian) * 24L + hour) * 60L + minute;
    Delta	= offset - (273L * 24L + 13L) * 60L + 23L;
    Phase	= Delta - (Delta / Length) * Length;

    Fraction	= Phase / Length;

    Length	= Length / 4;
    Phase	= Phase / Length;
    Delta	= Delta - (Delta / Length) * Length;
    DeltaH	= Delta / 60;
    Minute	= Delta - DeltaH * 60;
    Day		= DeltaH / 24;
    Hour	= DeltaH - Day * 24;
    return(Phase);
}
/*!*/
int
CharPos(x)
    double		x;
{
    register int	i;

    i = x * LINES + 0.5;
    if ((i += LINES + CENTER) < 1)
	i = 1;
    return(i);
}


void
Draw()
{
    register char	*p;
    register int	 i;
    register int	 end;
    register double	 y;
    register double	 cht;
    register double	 squisher;
    register double	 horizon;
    register double	 terminator;
    char		 Buffer[256];

    /* Clear screen? */

    if (Fraction < FULL)
	squisher = cos(TwoPi * Fraction);
    else
	squisher = cos(TwoPi * (Fraction - FULL));

    cht = (double)2.0 / (LINES - 6.0);
    for (y = 0.93; y > -1.0; y -= cht)
    {
	for (i = sizeof Buffer, p = Buffer; --i >= 0; )
	    *p++ = ' ';
	horizon = sqrt(1.0 - y * y);
	Buffer[    CharPos(-horizon)]	= LEDGE;
	Buffer[i = CharPos( horizon)]	= REDGE;
	Buffer[++i]			= '\0';
	terminator = horizon * squisher;
	if (Fraction > ZERO && Fraction < (1.0 - ZERO))
	{
	    if (Fraction < FULL)
	    {
		i   = CharPos( terminator);
		end = CharPos( horizon);
	    }
	    else
	    {
		i   = CharPos(-horizon);
		end = CharPos( terminator);
	    }
	    while (i <= end)
		Buffer[i++] = BRIGHT;
	}
/*	printf(" %s\n", Buffer);*/
    }
}
/*!*/
main(ac)
    int			 ac;
{
    register long	 Phase;

    Phase = Calculate();
    if (ac == 1)
	Draw();

    printf("It is ");
    if (Day == 0 && Hour == 0 && Minute == 0)
	printf("exactly");
    else
    {
	if (Day)
	{
	    printf("%ld day", Day);
	    Plural(Day);
	    if (Hour | Day)
		printf(", ");
	}
	if (Hour)
	{
	    printf("%ld hour", Hour);
	    Plural(Hour);
	    if (Minute)
		printf(", ");
	}
	if (Minute)
	{
	    printf("%ld minute", Minute);
	    Plural(Minute);
	}
	printf(" since");
    }
    printf(" the %s moon.\n", Moon[Phase]);

    exit(0);
}
SHAR_EOF
fi # end of overwriting check
#	End of shell archive
exit 0