summaryrefslogtreecommitdiff
path: root/info/digests/texhax/txh/ramsdell.txh
diff options
context:
space:
mode:
Diffstat (limited to 'info/digests/texhax/txh/ramsdell.txh')
-rw-r--r--info/digests/texhax/txh/ramsdell.txh477
1 files changed, 477 insertions, 0 deletions
diff --git a/info/digests/texhax/txh/ramsdell.txh b/info/digests/texhax/txh/ramsdell.txh
new file mode 100644
index 0000000000..2bad22f1ad
--- /dev/null
+++ b/info/digests/texhax/txh/ramsdell.txh
@@ -0,0 +1,477 @@
+Date: Thu, 5 Nov 87 09:45:56 EST
+From: ramsdell%linus@mitre-bedford.ARPA
+To: texhax@score.stanford.edu, rs@uunet.UU.NET
+Subject: Lunisolar calendar generator for LaTeX
+
+Here is a program that prints the current phase of the moon. It also
+generates the LaTeX commands that produce lunisolar calendars. It is
+a handy program for astronomers and sailors, as a lunisolar calendar
+tells when the moon will be bright.
+John
+
+
+#! /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:
+# lunisolar.c
+# This archive created: Thu Nov 5 08:53:57 1987
+export PATH; PATH=/bin:$PATH
+if test -f 'lunisolar.c'
+then
+ echo shar: will not over-write existing file "'lunisolar.c'"
+else
+cat << \SHAR_EOF > 'lunisolar.c'
+/* Prints the phase of the moon and generates LaTeX commands */
+/* that produce lunisolar calendars. */
+/* Usage: lunisolar
+ gives the phase of the moon,
+ and: lunisolar <year> <time_zone>
+ generates a lunisolar calendar for LaTeX.
+ Time zone may be one of:
+ GMT NST AST EST CST
+ MST PST YST HST BST
+ JST
+*/
+/* Construct with the command "cc -O -o lunisolar lunisolar.c -lm". */
+/* John D. Ramsdell - November 1987 */
+static char copyright[] =
+ "Copyright 1987 by The MITRE Corporation. All rights reserved.";
+/*
+ * Permission to use, copy, modify, and distribute this
+ * software and its documentation for any purpose and without
+ * fee is hereby granted, provided that the above copyright
+ * notice appear in all copies. The MITRE Corporation
+ * makes no representations about the suitability of this
+ * software for any purpose. It is provided "as is" without
+ * express or implied warranty.
+ */
+
+#include <stdio.h>
+#include <math.h>
+
+#define PI 3.141592653589793
+#define NEGATIVE_PI (-1.0 * PI)
+#define TWO_PI (2.0 * PI)
+#define RADIANS_PER_DEGREE (PI / 180.0)
+
+struct {
+ char *name; /* Name of time zone. */
+ int offset; /* Offset in minutes. */
+} tz_map[] =
+{
+ { "GMT", 0*60 }, /* Greenwich Mean Time */
+ { "NST", 7*30 }, /* Newfoundland is 3.5 hours */
+ /* different from GMT. */
+ { "AST", 4*60 }, /* Alantic Standard Time. */
+ { "EST", 5*60 }, /* Eastern Standard Time. */
+ { "CST", 6*60 }, /* Central Standard Time. */
+ { "MST", 7*60 }, /* Mountain Standard Time. */
+ { "PST", 8*60 }, /* Pacific Standard Time. */
+ { "YST", 9*60 }, /* Yukon Standard Time. */
+ { "HST", 10*60 }, /* Hawaiian Standard Time. */
+ { "BST", 11*60 }, /* Bering Standard Time. */
+ { "JST", -9*60 }, /* Japan Standard Time. */
+ { "" , 0 } /* Mark end of list with "". */
+};
+
+char *time_zone_name; /* Selected time zone name. */
+double time_zone_offset; /* Selected offset in minutes. */
+
+int select_time_zone (name)
+ char *name;
+{
+ int i;
+ if (strcmp (name, "") == 0) return 1; /* Return error. */
+ for (i = 0; strcmp (tz_map[i].name, "") != 0; i++)
+ if (strcmp (name, tz_map[i].name) == 0) {
+ time_zone_name = name;
+ time_zone_offset = (double) tz_map[i].offset;
+ return 0; /* Found match. */
+ }
+ return 1; /* No match found. */
+}
+
+int leap_year (year) /* True if year is a leap_year. */
+ int year;
+{
+ return year % 4 == 0 && year % 100 != 0 || year % 400 == 0;
+}
+
+/* Time is most often represented as a double precision number */
+/* in units of days. Angles are in radians. */
+
+/* J2000 is the number of days between January 1, 2000; 12:00:00 GMT */
+/* and the local origin of time. */
+double J2000; /* This date is really called J2000.0. */
+
+void make_J2000 () /* Assumes time origin of */
+{ /* January 1, 1970; 00:00:00 GMT. */
+ int year;
+ J2000 = (2000 - 1970) * 365 + 0.5;
+ for (year = 1970; year < 2000; year++) /* Account for leap years. */
+ if (leap_year (year)) J2000 += 1.0;
+}
+
+double days_after_J2000 () /* Returns the current time, */
+{ /* in units of days, after J2000.0. */
+ long seconds = time((long *) 0);
+ double seconds_per_day = 24 * 60 * 60;
+ return seconds / seconds_per_day - J2000;
+}
+
+double normalize_angle (angle) /* Returns the angle between */
+ double angle; /* -PI < angle <= PI. */
+{
+ if (angle > PI)
+ do angle -= TWO_PI; while (angle > PI);
+ else
+ while (angle <= NEGATIVE_PI) angle += TWO_PI;
+ return angle;
+}
+
+/*******************************************************************/
+
+/* Astronomical almanac */
+
+/*
+ * All formulas are from:
+ * The Astronomical Almanac for the Year 1984,
+ * US Naval Observatory and Royal Greenwich Observatory,
+ * US Government Printing Office, Washington DC, 1984.
+ */
+
+/* Angular position of the sun to a */
+/* precision of 0.01 degrees. (Page C24). */
+
+#define SUN0 (RADIANS_PER_DEGREE * 280.460)
+#define SUN1 (RADIANS_PER_DEGREE * 0.9856474)
+#define SUN2 (RADIANS_PER_DEGREE * 357.528)
+#define SUN3 (RADIANS_PER_DEGREE * 0.9856003)
+#define SUN4 (RADIANS_PER_DEGREE * 1.915)
+#define SUN5 (RADIANS_PER_DEGREE * 0.020)
+
+double sun_position (days)
+ double days;
+{
+ double mean_longitude_of_sun, mean_anomaly, ecliptic_longitude;
+ mean_longitude_of_sun =
+ normalize_angle (SUN0 + SUN1 * days);
+ mean_anomaly =
+ normalize_angle (SUN2 + SUN3 * days);
+ ecliptic_longitude =
+ normalize_angle (mean_longitude_of_sun
+ + SUN4 * sin (mean_anomaly)
+ + SUN5 * sin (2.0 * mean_anomaly));
+ return ecliptic_longitude;
+}
+
+/* Angular velocity of the sun. Derivative of sun_position. */
+
+double sun_velocity (days)
+ double days;
+{
+ double mean_anomaly =
+ normalize_angle (SUN2 + SUN3 * days);
+ return SUN1 + SUN4 * SUN3 * cos (mean_anomaly)
+ + SUN5 * 2.0 * SUN3 * cos (2.0 * mean_anomaly);
+}
+
+/* Angular position of the moon to a */
+/* precision of 0.3 degrees. (Page D46). */
+
+#define DAYS_PER_JULEAN_CENTURY 36525.0
+#define RADIAN_CENTURY (RADIANS_PER_DEGREE / DAYS_PER_JULEAN_CENTURY)
+
+#define MOON0 (RADIANS_PER_DEGREE * 218.32)
+#define MOON1 (RADIAN_CENTURY * 481267.883)
+#define MOON2A (RADIANS_PER_DEGREE * 6.29)
+#define MOON2B (RADIANS_PER_DEGREE * 134.9)
+#define MOON2C (RADIAN_CENTURY * 477198.85)
+#define MOON3A (RADIANS_PER_DEGREE * -1.27)
+#define MOON3B (RADIANS_PER_DEGREE * 259.2)
+#define MOON3C (RADIAN_CENTURY * -413335.38)
+#define MOON4A (RADIANS_PER_DEGREE * 0.66)
+#define MOON4B (RADIANS_PER_DEGREE * 235.7)
+#define MOON4C (RADIAN_CENTURY * 890534.23)
+#define MOON5A (RADIANS_PER_DEGREE * 0.21)
+#define MOON5B (RADIANS_PER_DEGREE * 269.9)
+#define MOON5C (RADIAN_CENTURY * 954397.70)
+#define MOON6A (RADIANS_PER_DEGREE * -0.19)
+#define MOON6B (RADIANS_PER_DEGREE * 357.5)
+#define MOON6C (RADIAN_CENTURY * 035999.05)
+#define MOON7A (RADIANS_PER_DEGREE * -0.11)
+#define MOON7B (RADIANS_PER_DEGREE * 186.6)
+#define MOON7C (RADIAN_CENTURY * 966404.05)
+
+double moon_position (days)
+ double days;
+{
+ return normalize_angle (MOON0
+ + MOON1 * days
+ + MOON2A * sin (MOON2B + MOON2C * days)
+ + MOON3A * sin (MOON3B + MOON3C * days)
+ + MOON4A * sin (MOON4B + MOON4C * days)
+ + MOON5A * sin (MOON5B + MOON5C * days)
+ + MOON6A * sin (MOON6B + MOON6C * days)
+ + MOON7A * sin (MOON7B + MOON7C * days));
+}
+
+/****************************************************************/
+
+/* Prints an English sentence giving the current phase of the moon. */
+#define PHASE_LIMIT MOON1
+void print_moon ()
+{
+ double days, phase; /* Computes the moon's phase by */
+ int percent; /* computing the difference between */
+ make_J2000 (); /* the sun and moon's */
+ days = days_after_J2000 (); /* ecliptic longitude. */
+ phase = sun_position (days);
+ phase = normalize_angle (moon_position (days) - phase);
+ percent = 50.0 * (1.0 - cos (phase)) + 0.5; /* Visable fraction. */
+ printf("The moon is ");
+ if (fabs (phase) < PHASE_LIMIT)
+ printf ("new");
+ else if (fabs (normalize_angle (phase + PI)) < PHASE_LIMIT)
+ printf ("full");
+ else if (fabs (phase - PI/2.0) < PHASE_LIMIT)
+ printf ("first quarter (%d%% of full)", percent);
+ else if (fabs (phase + PI/2.0) < PHASE_LIMIT)
+ printf ("last quarter (%d%% of full)", percent);
+ else if (phase > PI/2.0)
+ printf ("waxing and gibbous (%d%% of full)", percent);
+ else if (phase > 0.0)
+ printf ("a waxing crescent (%d%% of full)", percent);
+ else if (phase > PI/-2.0)
+ printf ("a waning crescent (%d%% of full)", percent);
+ else
+ printf ("waning and gibbous (%d%% of full)", percent);
+ printf (".\n");
+}
+
+/**********************************************************/
+
+/* lunisolar calendar routines. */
+
+int first_day_of_year (year) /* Returns the integer number of days */
+ int year; /* between the start of year and */
+{ /* J2000.0. */
+ int days;
+ days = 365 * (year - 2000);
+ if (year > 2000)
+ for (; year > 2000; year--)
+ if (leap_year (year)) days += 1;
+ else; /* Needed else! */
+ else
+ for (; year < 2000; year++)
+ if (leap_year (year)) days -= 1;
+ return days;
+}
+
+/* Routines that find the seasons. */
+
+#define DIGITS 15
+int zero (x, f, fp) /* Root finder using */
+ int x; /* Newton's method. */
+ double (*f) ();
+ double (*fp) ();
+{
+ int i;
+ double y, midnite, noon;
+ y = x;
+ for (i = 0; i < DIGITS; i++)
+ y = y - f(y)/fp(y);
+ noon = 0.5 + (time_zone_offset / 24.0 / 60.0);
+ midnite = floor (y - noon) + noon;
+ if (f (midnite) * f (midnite + 1.0) <= 0.0)
+ return midnite;
+ x = midnite;
+ printf ("%%Not sure about the season change for day %d.\n", x);
+ return x;
+}
+
+double phase; /* sun_zero has a root at the */
+double sun_zero (days) /* desired day. Used with zero */
+ double days; /* to find the seasons. */
+{
+ return normalize_angle (sun_position (days) - phase);
+}
+
+void find_seasons (first_day, seasons)
+ int first_day, *seasons;
+{ /* Remember Spring is the */
+ int i; /* first season of a year. */
+ phase = PI/-2.0; /* Find start of winters. */
+ seasons[0] = zero (first_day - 11, sun_zero, sun_velocity);
+ seasons[4] = zero (seasons[0] + 365, sun_zero, sun_velocity);
+ phase = 0.0; /* Find start of other seasons. */
+ for (i = 1; i < 4; i++, phase += PI/2.0)
+ seasons[i] = zero (seasons[i-1] + 91, sun_zero, sun_velocity);
+ printf ("%% Seasons relative to January 1:");
+ for (i = 0; i < 5; i++)
+ printf (" %d", seasons[i] - first_day);
+ printf (".\n");
+}
+
+void make_moon_table (seasons, moon) /* Computes the position of */
+ int *seasons; /* the moon for each day at */
+ float *moon; /* noon local time. */
+{
+ int i, day;
+ for (i = 0, day = seasons[0]; day < seasons[4]; i++, day++) {
+ double dday = day + time_zone_offset / (24.0 * 60.0);
+ moon[i] = normalize_angle (moon_position (dday) - sun_position (dday));
+ }
+}
+
+/* Routines that output LaTeX commands. */
+
+/* Dates spiral inward by an amount DELTA_RADIUS. */
+#define START_RADIUS 1.0
+#define DELTA_RADIUS 0.005
+float radius;
+
+int month, day, moon_index;
+int days_per_month[12] = {31, 28, 31, 30, 31, 30, 31, 31, 30, 31, 30, 31};
+
+void mark_dates (year, from, to, moon) /* Makes LaTeX statements */
+ int year, from, to; /* that place the dates. */
+ float *moon;
+{
+ radius = START_RADIUS;
+ for (; from < to; moon_index++, from++) {
+ printf ("\\put(%1.6f,%1.6f){\\makebox(0,0){%d/%d}}\n",
+ radius * sin (moon[moon_index]),
+ radius * cos (moon[moon_index]),
+ month, day);
+ radius -= DELTA_RADIUS;
+ day++;
+ if (day > days_per_month[month-1])
+ if (month == 2 && leap_year (year) && day == 29);
+ else {
+ day = 1;
+ month++;
+ if (month > 12) month = 1;
+ }
+ }
+}
+
+void header (season, year) /* Start of each season. */
+ char *season;
+ int year;
+{
+ printf ("\\begin{figure}\n");
+ printf ("\\begin{center}\n");
+ printf ("\\begin{picture}(2.0,2.0)(-1.0,-1.0)\n");
+ printf ("\\tiny\n");
+ printf ("\\put(0,0){\\makebox(0,0){\\Huge %s %d}}\n",
+ season, year);
+}
+
+void trailer () /* End of each season. */
+{
+ printf ("\\put(-1.0,0.0){\\line(1,0){0.5}}\n");
+ printf ("\\put(0.5,0.0){\\line(1,0){0.5}}\n");
+ printf ("\\put(0.0,-1.0){\\line(0,1){0.5}}\n");
+ printf ("\\put(0.0,-0.4){\\circle{0.1}}\n");
+ printf ("\\put(0.0,-0.3){\\makebox(0,0)[b]{\\large Full Moon}}\n");
+ printf ("\\put(0.0,0.5){\\line(0,1){0.5}}\n");
+ printf ("\\put(0.0,0.4){\\circle*{0.1}}\n");
+ printf ("\\put(0.0,0.3){\\makebox(0,0)[t]{\\large New Moon}}\n");
+ printf ("\\end{picture}\n");
+ printf ("\\\\ {\\Large Lunisolar Calendar}\n");
+ printf ("\\\\ {\\large Dates mark the lunar phase at noon %s.}\n",
+ time_zone_name);
+ printf ("\\end{center}\n");
+ printf ("\\end{figure}\n");
+}
+
+char *season_titles[4] =
+{ "Winter", "Spring", "Summer", "Fall"};
+
+void LaTeXize_tables (year, first_day, seasons, moon)
+ int year, first_day, *seasons;
+ float *moon;
+{
+ int a_season;
+ printf ("\\documentstyle{article}\n");
+ printf ("\\pagestyle{empty}\n");
+ printf ("\\begin{document}\n");
+ printf ("\\Large\n");
+ printf ("\\setlength{\\unitlength}{60mm}\n");
+ month = 12; /* December */
+ day = 32 - first_day + seasons[0];
+ moon_index = 0;
+ for (a_season = 0; a_season < 4; a_season++) {
+ header (season_titles[a_season], a_season == 0 ? year - 1 : year);
+ mark_dates (year, seasons[a_season], seasons[a_season+1], moon);
+ trailer ();
+ }
+ printf ("\\end{document}\n");
+}
+
+/* Lunisolar master routine. */
+
+int seasons[5]; /* Stores days that mark season changes. */
+float moon[370]; /* Stores moon phases for each day. */
+
+int lunisolar (year, tz) /* Constructs a LaTeX file that */
+ int year; char *tz; /* generates a lunisolar calendar */
+{ /* for the year year and time zone tz. */
+ if (year < 1950 || year > 2050) {
+ printf ("Program useful between the years 1950 and 2050.\n");
+ return 1; /* error return. */
+ }
+ else if (select_time_zone (tz) != 0)
+ return 1;
+ else {
+ int day_of_Jan1 = first_day_of_year (year);
+ printf ("%% Lunisolar calendar for %d.\n", year);
+ printf ("%% Constructed for %s, %1.2f hours %s of Greenwich.\n",
+ time_zone_name, fabs (time_zone_offset) / 60.0,
+ (time_zone_offset >= 0.0 ? "west" : "east"));
+ find_seasons (day_of_Jan1, seasons);
+ make_moon_table (seasons, moon);
+ LaTeXize_tables (year, day_of_Jan1, seasons, moon);
+ return 0;
+ }
+}
+
+main (argc, argv) /* Invokes print_moon with */
+ int argc; /* no arguments, and */
+ char **argv; /* lunisolar with one. */
+{
+ int i;
+ if (argc == 1) print_moon ();
+ else {
+ if (argc == 3) {
+ int year;
+ if (sscanf (argv[1], "%d", &year) == 1 &&
+ lunisolar (year, argv[2]) == 0)
+ exit (0);
+ }
+ fprintf (stderr, "Bad args:");
+ for (i = 0; i < argc; i++)
+ fprintf (stderr, " %s", argv[i]);
+ fprintf (stderr,
+ "\nUsage: %s\ngives the phase of the moon,\n",
+ argv[0]);
+ fprintf (stderr, "and: %s <year> <time_zone>\n", argv[0]);
+ fprintf (stderr, "generates a lunisolar calendar for LaTeX.\n");
+ fprintf (stderr, "Time zone may be one of:");
+ for (i = 0; strcmp(tz_map[i].name, "") != 0; i++) {
+ if (i % 5 == 0) fprintf (stderr, "\n");
+ fprintf (stderr, "%s\t", tz_map[i].name);
+ }
+ fprintf (stderr, "\n");
+ exit (1);
+ }
+}
+SHAR_EOF
+fi # end of overwriting check
+# End of shell archive
+exit 0
+-------
+-------