v09i040: ephem, v4.8, 3 of 5
Brandon S. Allbery - comp.sources.misc
allbery at uunet.UU.NET
Tue Nov 28 11:37:55 AEST 1989
Posting-number: Volume 9, Issue 40
Submitted-by: ecd at umn-cs.cs.umn.edu@ncs-med.UUCP (Elwood C. Downey)
Archive-name: ephem2/part03
# This is a "shell archive" file; run it with sh.
# This is file 3.
echo x io.c
cat > io.c << 'xXx'
/* this file (in principle) contains all the device-dependent code for
* handling screen movement and reading the keyboard. public routines are:
* c_pos(r,c), c_erase(), c_eol();
* chk_char(), read_char(), read_line (buf, max); and
* byetty().
* N.B. we assume output may be performed by printf(), putchar() and
* fputs(stdout). since these are buffered we flush first in read_char().
* #define UNIX or TURBO_C to give two popular versions.
* UNIX uses termcap; TURBO_C uses ANSI and the IBM-PC keyboard codes.
* TURBO_C should work for Lattice too.
*/
#define UNIX
/* #define TURBO_C */
#include <stdio.h>
#include "screen.h"
#ifdef UNIX
#include <sgtty.h>
#include <signal.h>
extern char *tgoto();
static char *cm, *ce, *cl, *kl, *kr, *ku, *kd; /* curses sequences */
static int tloaded;
static int ttysetup;
static struct sgttyb orig_sgtty;
/* move cursor to row, col, 1-based.
* we assume this also moves a visible cursor to this location.
*/
c_pos (r, c)
int r, c;
{
if (!tloaded) tload();
fputs (tgoto (cm, c-1, r-1), stdout);
}
/* erase entire screen. */
c_erase()
{
if (!tloaded) tload();
fputs (cl, stdout);
}
/* erase to end of line */
c_eol()
{
if (!tloaded) tload();
fputs (ce, stdout);
}
/* return 0 if there is a char that may be read without blocking, else -1 */
chk_char()
{
long n = 0;
if (!ttysetup) setuptty();
ioctl (0, FIONREAD, &n);
return (n > 0 ? 0 : -1);
}
/* read the next char, blocking if necessary, and return it. don't echo.
* map the arrow keys if we can too into hjkl
*/
read_char()
{
char c;
if (!ttysetup) setuptty();
fflush (stdout);
read (0, &c, 1);
return (chk_arrow (c & 0177)); /* just ASCII, please */
}
/* used to time out of a read */
static got_alrm;
static
on_alrm()
{
got_alrm = 1;
}
/* see if c is the first of any of the curses arrow key sequences.
* if it is, read the rest of the sequence, and return the hjkl code
* that corresponds.
* if no match, just return c.
*/
static
chk_arrow (c)
register char c;
{
register char *seq;
if (c == *(seq = kl) || c == *(seq = kd) || c == *(seq = ku)
|| c == *(seq = kr)) {
char seqa[10]; /* maximum arrow escape sequence ever expected */
int l = strlen(seq);
seqa[0] = c;
/* most arrow keys generate sequences starting with ESC. if so
* c might just be a lone ESC; time out if so.
*/
got_alrm=0;
if (c == '') {
signal (SIGALRM, on_alrm);
alarm(1);
}
read (0, seqa+1, l-1);
if (got_alrm == 0) {
if (c == '')
alarm(0);
seqa[l] = '\0';
if (strcmp (seqa, kl) == 0)
return ('h');
if (strcmp (seqa, kd) == 0)
return ('j');
if (strcmp (seqa, ku) == 0)
return ('k');
if (strcmp (seqa, kr) == 0)
return ('l');
}
}
return (c);
}
/* do whatever might be necessary to get the screen and/or tty back into shape.
*/
byetty()
{
ioctl (0, TIOCSETP, &orig_sgtty);
}
static
tload()
{
extern char *getenv(), *tgetstr();
extern char *UP, *BC;
char *egetstr();
static char tbuf[512];
char rawtbuf[1024];
char *tp;
char *ptr;
if (!(tp = getenv ("TERM"))) {
printf ("Must have addressable cursor\n");
exit(1);
}
if (!ttysetup) setuptty();
if (tgetent (rawtbuf, tp) != 1) {
printf ("can't find termcap for %s\n", tp);
exit (1);
}
ptr = tbuf;
ku = egetstr ("ku", &ptr);
kd = egetstr ("kd", &ptr);
kl = egetstr ("kl", &ptr);
kr = egetstr ("kr", &ptr);
cm = egetstr ("cm", &ptr);
ce = egetstr ("ce", &ptr);
cl = egetstr ("cl", &ptr);
UP = egetstr ("up", &ptr);
if (!tgetflag ("bs"))
BC = egetstr ("bc", &ptr);
tloaded = 1;
}
/* like tgetstr() but discard curses delay codes, for now anyways */
static char *
egetstr (name, sptr)
char *name;
char **sptr;
{
extern char *tgetstr();
register char c, *s;
s = tgetstr (name, sptr);
while ((c = *s) >= '0' && c <= '9')
s += 1;
return (s);
}
static
setuptty()
{
extern ospeed;
struct sgttyb sg;
ioctl (0, TIOCGETP, &orig_sgtty);
sg = orig_sgtty;
ospeed = sg.sg_ospeed;
sg.sg_flags &= ~ECHO; /* do our own echoing */
sg.sg_flags &= ~CRMOD; /* leave CR and LF unchanged */
sg.sg_flags |= XTABS; /* no tabs with termcap */
sg.sg_flags |= CBREAK; /* wake up on each char but can still kill */
ioctl (0, TIOCSETP, &sg);
ttysetup = 1;
}
#endif
#ifdef TURBO_C
/* position cursor.
* (ANSI: ESC [ r ; c f) (r/c are numbers given in ASCII digits)
*/
c_pos (r, c)
int r, c;
{
printf ("[%d;%df", r, c);
}
/* erase entire screen. (ANSI: ESC [ 2 J) */
c_erase()
{
printf ("[2J");
}
/* erase to end of line. (ANSI: ESC [ K) */
c_eol()
{
printf ("[K");
}
/* return 0 if there is a char that may be read without blocking, else -1 */
chk_char()
{
return (kbhit() == 0 ? -1 : 0);
}
/* read the next char, blocking if necessary, and return it. don't echo.
* map the arrow keys if we can too into hjkl
*/
read_char()
{
int c;
fflush (stdout);
c = getch();
if (c == 0) {
/* get scan code; convert to direction hjkl if possible */
c = getch();
switch (c) {
case 0x4b: c = 'h'; break;
case 0x50: c = 'j'; break;
case 0x48: c = 'k'; break;
case 0x4d: c = 'l'; break;
}
}
return (c);
}
/* do whatever might be necessary to get the screen and/or tty back into shape.
*/
byetty()
{
}
#endif
/* read up to max chars into buf, with cannonization.
* add trailing '\0' (buf is really max+1 chars long).
* return count of chars read (not counting '\0').
* assume cursor is already positioned as desired.
* if type END when n==0 then return -1.
*/
read_line (buf, max)
char buf[];
int max;
{
static char erase[] = "\b \b";
int n, c;
int done;
#ifdef UNIX
if (!ttysetup) setuptty();
#endif
for (done = 0, n = 0; !done; )
switch (c = read_char()) { /* does not echo */
case cntrl('h'): /* backspace or */
case 0177: /* delete are each char erase */
if (n > 0) {
fputs (erase, stdout);
n -= 1;
}
break;
case cntrl('u'): /* line erase */
while (n > 0) {
fputs (erase, stdout);
n -= 1;
}
break;
case '\r': /* EOL */
done++;
break;
default: /* echo and store, if ok */
if (n == 0 && c == END)
return (-1);
if (n >= max)
putchar (cntrl('g'));
else if (c >= ' ') {
putchar (c);
buf[n++] = c;
}
}
buf[n] = '\0';
return (n);
}
xXx
echo x main.c
cat > main.c << 'xXx'
/* main program.
* set options.
* init screen and circumstances.
* enter infinite loop updating screen and allowing operator input.
*/
#include <stdio.h>
#include <signal.h>
#include <math.h>
#include "astro.h"
#include "circum.h"
#include "screen.h"
extern char *getenv();
/* shorthands for fields of a Now structure, now.
* first undo the ones for a Now pointer from circum.h.
*/
#undef mjd
#undef lat
#undef lng
#undef tz
#undef temp
#undef pressure
#undef height
#undef epoch
#undef tznm
#define mjd now.n_mjd
#define lat now.n_lat
#define lng now.n_lng
#define tz now.n_tz
#define temp now.n_temp
#define pressure now.n_pressure
#define height now.n_height
#define epoch now.n_epoch
#define tznm now.n_tznm
static char *cfgfile = "ephem.cfg"; /* default config filename */
static Now now; /* where when and how, right now */
static double tminc; /* hrs to inc time by each loop; RTC means use clock */
static int nstep; /* steps to go before stopping */
static int optwi; /* set when want to display dawn/dusk/len-of-night */
static int oppl; /* mask of (1<<planet) bits; set when want to show it */
main (ac, av)
int ac;
char *av[];
{
void bye();
static char freerun[] =
"Running... press any key to stop to make changes.";
static char prmpt[] =
"Move to another field, RETURN to change this field, ? for help, or q to run";
static char hlp[] =
"arrow keys move to field; any key stops running; ^d exits; ^l redraws";
int curr = R_NSTEP, curc = C_NSTEPV; /* must start somewhere */
int sflag = 0; /* not silent, by default */
int one = 1; /* use a variable so optimizer doesn't get disabled */
int srchdone = 0; /* true when search funcs say so */
int newcir = 2; /* set when circumstances change - means don't tminc */
while ((--ac > 0) && (**++av == '-')) {
char *s;
for (s = *av+1; *s != '\0'; s++)
switch (*s) {
case 's': /* no credits "silent" (don't publish this) */
sflag++;
break;
case 'c': /* set name of config file to use */
if (--ac <= 0) usage("-c but no config file");
cfgfile = *++av;
break;
default:
usage("Bad - option");
}
}
if (!sflag)
credits();
/* fresh screen.
* crack config file, THEN args so args may override.
*/
c_erase();
read_cfgfile (cfgfile);
read_fieldargs (ac, av);
/* set up to clean up screen and tty if interrupted */
signal (SIGINT, bye);
/* update screen forever (until QUIT) */
while (one) {
nstep -= 1;
/* recalculate everything and update all the fields */
redraw_screen (newcir);
mm_newcir (0);
/* let searching functions change tminc and check for done */
srchdone = srch_eval (mjd, &tminc) < 0;
print_tminc(0); /* to show possibly new search increment */
/* update plot file, now that all fields are up to date and
* search function has been evaluated.
*/
plot();
/* stop loop to allow op to change parameters:
* if a search evaluation converges (or errors out),
* or if steps are done,
* or if op hits any key.
*/
newcir = 0;
if (srchdone || nstep <= 0 || (chk_char()==0 && read_char()!=0)) {
int fld;
/* update screen with the current stuff if stopped during
* unattended plotting since last redraw_screen() didn't.
*/
if (plot_ison() && nstep > 0)
redraw_screen (1);
/* return nstep to default of 1 */
if (nstep <= 0) {
nstep = 1;
print_nstep (0);
}
/* change fields until END.
* update all time fields if any are changed
* and print NEW CIRCUMSTANCES if any have changed.
* QUIT causes bye() to be called and we never return.
*/
while(fld = sel_fld(curr,curc,alt_menumask()|F_CHG,prmpt,hlp)) {
if (chg_fld ((char *)0, fld)) {
mm_now (&now, 1);
mm_newcir(1);
newcir = 1;
}
curr = unpackr (fld);
curc = unpackc (fld);
}
if (nstep > 1)
f_prompt (freerun);
}
/* increment time only if op didn't change cirumstances */
if (!newcir)
inc_mjd (&now, tminc);
}
return (0);
}
/* draw all the stuff on the screen, using the current menu.
* if how_much == 0 then just update fields that need it;
* if how_much == 1 then redraw all fields;
* if how_much == 2 then erase the screen and redraw EVERYTHING.
*/
redraw_screen (how_much)
int how_much;
{
if (how_much == 2)
c_erase();
/* print the single-step message if this is the last loop */
if (nstep < 1)
print_updating();
if (how_much == 2) {
mm_borders();
mm_labels();
srch_prstate(1);
plot_prstate(1);
alt_labels();
}
/* if just updating changed fields while plotting unattended then
* suppress most screen updates except
* always show nstep to show plot loops to go and
* always show tminc to show search convergence progress.
*/
print_nstep(how_much);
print_tminc(how_much);
if (how_much == 0 && plot_ison() && nstep > 0)
f_off();
/* print all the time-related fields */
mm_now (&now, how_much);
if (optwi)
mm_twilight (&now, how_much);
/* print solar system body info */
print_bodies (how_much);
f_on();
}
/* clean up and exit
* TODO: want to ask "are you sure?" ?
*/
void
bye()
{
c_erase();
byetty();
exit (0);
}
static
usage(why)
char *why;
{
/* don't advertise -s (silent) option */
c_erase();
f_string (1, 1, why);
f_string (2, 1, "usage: [-c <configfile>] [field=value] ...\r\n");
byetty();
exit (1);
}
/* read cfg file, fn, if present.
* if errors in file, call usage() (which exits).
* try $HOME/.ephemrc as last resort.
*/
static
read_cfgfile(fn)
char *fn;
{
char buf[128];
FILE *fp;
fp = fopen (fn, "r");
if (!fp) {
char *home = getenv ("HOME");
if (home) {
sprintf (buf, "%s/.ephemrc", home);
fp = fopen (buf, "r");
if (!fp)
return; /* oh well */
fn = buf; /* save fn for error report */
}
}
while (fgets (buf, sizeof(buf), fp)) {
buf[strlen(buf)-1] = '\0'; /* discard trailing \n */
if (crack_fieldset (buf) < 0) {
char why[128];
sprintf (why, "Unknown field spec in %s: %s\n", fn, buf);
usage (why);
}
}
fclose (fp);
}
/* process the field specs from the command line.
* if trouble call usage() (which exits).
*/
static
read_fieldargs (ac, av)
int ac; /* number of such specs */
char *av[]; /* array of strings in form <field_name value> */
{
while (--ac >= 0) {
char *fs = *av++;
if (crack_fieldset (fs) < 0) {
char why[128];
sprintf (why, "Unknown command line field spec: %s", fs);
usage (why);
}
}
}
/* process a field spec in buf, either from config file or argv.
* return 0 if recognized ok, else -1.
*/
static
crack_fieldset (buf)
char *buf;
{
if (strncmp ("LAT", buf, 3) == 0)
(void) chg_fld (buf+4, rcfpack (R_LAT,C_LATV,0));
else if (strncmp ("LONG", buf, 4) == 0)
(void) chg_fld (buf+5, rcfpack (R_LONG,C_LONGV,0));
else if (strncmp ("UT", buf, 2) == 0)
(void) chg_fld (buf+3, rcfpack (R_UT,C_UTV,0));
else if (strncmp ("UD", buf, 2) == 0)
(void) chg_fld (buf+3, rcfpack (R_UD,C_UD,0));
else if (strncmp ("TZONE", buf, 5) == 0)
(void) chg_fld (buf+6, rcfpack (R_TZONE,C_TZONEV,0));
else if (strncmp ("TZNAME", buf, 6) == 0)
(void) chg_fld (buf+7, rcfpack (R_TZN,C_TZN,0));
else if (strncmp ("HEIGHT", buf, 6) == 0)
(void) chg_fld (buf+7, rcfpack (R_HEIGHT,C_HEIGHTV,0));
else if (strncmp ("NSTEP", buf, 5) == 0)
(void) chg_fld (buf+6, rcfpack (R_NSTEP,C_NSTEPV,0));
else if (strncmp ("STPSZ", buf, 5) == 0)
(void) chg_fld (buf+6, rcfpack (R_STPSZ,C_STPSZV,0));
else if (strncmp ("TEMP", buf, 4) == 0)
(void) chg_fld (buf+5, rcfpack (R_TEMP,C_TEMPV,0));
else if (strncmp ("PRES", buf, 4) == 0)
(void) chg_fld (buf+5, rcfpack (R_PRES,C_PRESV,0));
else if (strncmp ("EPOCH", buf, 5) == 0)
(void) chg_fld (buf+6, rcfpack (R_EPOCH,C_EPOCHV,0));
else if (strncmp ("JD", buf, 2) == 0)
(void) chg_fld (buf+3, rcfpack (R_JD,C_JDV,0));
else if (strncmp ("OBJN", buf, 4) == 0)
(void) chg_fld (buf+5, rcfpack (R_OBJX,C_OBJ,0));
else if (strncmp ("OBJRA", buf, 5) == 0)
(void) chg_fld (buf+6, rcfpack (R_OBJX,C_RA,0));
else if (strncmp ("OBJDEC", buf, 6) == 0)
(void) chg_fld (buf+7, rcfpack (R_OBJX,C_DEC,0));
else if (strncmp ("PROPTS", buf, 6) == 0) {
char *bp = buf+7;
while (*bp)
switch (*bp++) {
case 'T': optwi = 1; break;
case 'S': oppl |= (1<<SUN); break;
case 'M': oppl |= (1<<MOON); break;
case 'e': oppl |= (1<<MERCURY); break;
case 'v': oppl |= (1<<VENUS); break;
case 'm': oppl |= (1<<MARS); break;
case 'j': oppl |= (1<<JUPITER); break;
case 's': oppl |= (1<<SATURN); break;
case 'u': oppl |= (1<<URANUS); break;
case 'n': oppl |= (1<<NEPTUNE); break;
case 'p': oppl |= (1<<PLUTO); break;
}
} else
return (-1);
return (0);
}
/* change the field at rcpk according to the optional string input at bp.
* if bp is != 0 use it, else issue read_line() and use buffer.
* then sscanf the buffer and update the corresponding (global) variable(s)
* or do whatever a pick at that field should do.
* return 1 if we change a field that invalidates any of the times or
* to update all related fields.
*/
static
chg_fld (bp, rcpk)
char *bp;
int rcpk;
{
char buf[NC];
int deghrs = 0, mins = 0, secs = 0;
int new = 0;
/* switch on just the row/col portion */
switch (unpackrc(rcpk)) {
case rcfpack (R_ALTM, C_ALTM, 0):
if (altmenu_setup() == 0) {
print_updating();
alt_nolabels();
clrall_bodies();
alt_labels();
print_bodies(1);
}
break;
case rcfpack (R_JD, C_JDV, 0):
if (!bp) {
static char p[] = "Julian Date (or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else
mjd = atof(bp) - 2415020L;
set_t0 (&now);
new = 1;
break;
case rcfpack (R_UD, C_UD, 0):
if (!bp) {
static char p[] = "utc date (m/d/y, or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else {
double fday, newmjd0;
int month, day, year;
mjd_cal (mjd, &month, &fday, &year); /* init with now */
day = (int)fday; /* ignore partial days */
f_sscansex (bp, &month, &day, &year);
cal_mjd (month, (double)day, year, &newmjd0);
mjd = newmjd0 + mjd_hr(mjd)/24.0;
}
set_t0 (&now);
new = 1;
break;
case rcfpack (R_UT, C_UTV, 0):
if (!bp) {
static char p[] = "utc time (h:m:s, or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else {
double newutc = (mjd-mjd_day(mjd)) * 24.0;
f_dec_sexsign (newutc, °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &newutc);
mjd = mjd_day(mjd) + newutc/24.0;
}
set_t0 (&now);
new = 1;
break;
case rcfpack (R_LD, C_LD, 0):
if (!bp) {
static char p[] = "local date (m/d/y, or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else {
double fday, newlmjd0;
int month, day, year;
mjd_cal (mjd-tz/24.0, &month, &fday, &year); /* now */
day = (int)fday; /* ignore partial days */
f_sscansex (bp, &month, &day, &year);
cal_mjd (month, (double)day, year, &newlmjd0);
mjd = newlmjd0 + mjd_hr(mjd)/24.0;
mjd += newlmjd0 - mjd_day(mjd-tz/24.0);
}
set_t0 (&now);
new = 1;
break;
case rcfpack (R_LT, C_LT, 0):
if (!bp) {
static char p[] = "local time (h:m:s, or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else {
double newlt = (mjd-mjd_day(mjd)) * 24.0 - tz;
range (&newlt, 24.0);
f_dec_sexsign (newlt, °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &newlt);
mjd = mjd_day(mjd-tz/24.0) + (newlt + tz)/24.0;
}
set_t0 (&now);
new = 1;
break;
case rcfpack (R_LST, C_LSTV, 0):
if (!bp) {
static char p[] = "local sidereal time (h:m:s, or N for Now): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'n' || bp[0] == 'N')
time_fromsys (&now);
else {
double lst, utc;
now_lst (&now, &lst);
f_dec_sexsign (lst, °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &lst);
lst -= radhr(lng); /* convert to gst */
range (&lst, 24.0);
gst_utc (mjd_day(mjd), lst, &utc);
mjd = mjd_day(mjd) + utc/24.0;
}
set_t0 (&now);
new = 1;
break;
case rcfpack (R_TZN, C_TZN, 0):
if (!bp) {
static char p[] = "timezone abbreviation (3 char max): ";
f_prompt (p);
if (read_line (buf, 3) <= 0)
break;
bp = buf;
}
strcpy (tznm, bp);
new = 1;
break;
case rcfpack (R_TZONE, C_TZONEV, 0):
if (!bp) {
static char p[] = "hours behind utc: ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
f_dec_sexsign (tz, °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &tz);
new = 1;
break;
case rcfpack (R_LONG, C_LONGV, 0):
if (!bp) {
static char p[] = "longitude (+ west) (d:m:s): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
f_dec_sexsign (-raddeg(lng), °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &lng);
lng = degrad (-lng); /* want - radians west */
new = 1;
break;
case rcfpack (R_LAT, C_LATV, 0):
if (!bp) {
static char p[] = "latitude (+ north) (d:m:s): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
f_dec_sexsign (raddeg(lat), °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &lat);
lat = degrad (lat);
new = 1;
break;
case rcfpack (R_HEIGHT, C_HEIGHTV, 0):
if (!bp) {
static char p[] = "height above sea level (ft): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
sscanf (bp, "%F", &height);
height /= 2.093e7; /* convert ft to earth radii above sea level */
new = 1;
break;
case rcfpack (R_NSTEP, C_NSTEPV, 0):
if (!bp) {
static char p[] = "number of steps to run: ";
f_prompt (p);
if (read_line (buf, 8) <= 0)
break;
bp = buf;
}
sscanf (bp, "%d", &nstep);
print_nstep (0);
break;
case rcfpack (R_TEMP, C_TEMPV, 0):
if (!bp) {
static char p[] = "temperature (deg.F): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
sscanf (bp, "%F", &temp);
temp = 5./9.*(temp - 32.0); /* want degs C */
new = 1;
break;
case rcfpack (R_PRES, C_PRESV, 0):
if (!bp) {
static char p[] =
"atmos pressure (in. Hg; 0 for no refraction correction): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
sscanf (bp, "%F", &pressure);
pressure *= 33.86; /* want mBar */
new = 1;
break;
case rcfpack (R_EPOCH, C_EPOCHV, 0):
if (!bp) {
static char p[] = "epoch (year, or E for Equinox of Date): ";
f_prompt (p);
if (read_line (buf, PW-strlen(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'e' || bp[0] == 'E')
epoch = EOD;
else {
double e;
sscanf (bp, "%F", &e);
cal_mjd (1, 1.0, (int)e, &epoch); /* convert to mjd */
epoch += (e - (int)e)*365.24;
}
new = 1;
break;
case rcfpack (R_STPSZ, C_STPSZV, 0):
if (!bp) {
static char p[] =
"step size increment (h:m:s, or <x>D, or RTC): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
if (bp[0] == 'r' || bp[0] == 'R')
tminc = RTC;
else {
int last = strlen (bp) - 1;
if (bp[last] == 'd' || bp[last] == 'D') {
/* ends in d so treat as a number of days */
double x;
sscanf (bp, "%G", &x);
tminc = x * 24.0;
} else {
if (tminc == RTC)
deghrs = mins = secs = 0;
else
f_dec_sexsign (tminc, °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &tminc);
}
}
print_tminc(0);
set_t0 (&now);
break;
case rcfpack (R_PLOT, C_PLOT, 0):
plot_setup();
if (plot_ison())
new = 1;
break;
case rcfpack (R_WATCH, C_WATCH, 0):
watch (&now, tminc, oppl);
break;
case rcfpack (R_DAWN, C_DAWN, 0):
case rcfpack (R_DUSK, C_DUSK, 0):
case rcfpack (R_LON, C_LON, 0):
if (optwi ^= 1) {
print_updating();
mm_twilight (&now, 1);
} else {
f_blanks (R_DAWN, C_DAWNV, 5);
f_blanks (R_DUSK, C_DUSKV, 5);
f_blanks (R_LON, C_LONV, 5);
}
break;
case rcfpack (R_SRCH, C_SRCH, 0):
srch_setup();
if (srch_ison())
new = 1;
break;
case rcfpack (R_SUN, C_OBJ, 0):
if ((oppl ^= (1<<SUN)) & (1<<SUN)) {
print_updating();
alt_body (SUN, 1, &now);
} else
alt_nobody (SUN);
break;
case rcfpack (R_MOON, C_OBJ, 0):
if ((oppl ^= (1<<MOON)) & (1<<MOON)) {
print_updating();
alt_body (MOON, 1, &now);
} else
alt_nobody (MOON);
break;
case rcfpack (R_MERCURY, C_OBJ, 0):
if ((oppl ^= (1<<MERCURY)) & (1<<MERCURY)) {
print_updating();
alt_body (MERCURY, 1, &now);
} else
alt_nobody (MERCURY);
break;
case rcfpack (R_VENUS, C_OBJ, 0):
if ((oppl ^= (1<<VENUS)) & (1<<VENUS)) {
print_updating();
alt_body (VENUS, 1, &now);
} else
alt_nobody (VENUS);
break;
case rcfpack (R_MARS, C_OBJ, 0):
if ((oppl ^= (1<<MARS)) & (1<<MARS)) {
print_updating();
alt_body (MARS, 1, &now);
} else
alt_nobody (MARS);
break;
case rcfpack (R_JUPITER, C_OBJ, 0):
if ((oppl ^= (1<<JUPITER)) & (1<<JUPITER)) {
print_updating();
alt_body (JUPITER, 1, &now);
} else
alt_nobody (JUPITER);
break;
case rcfpack (R_SATURN, C_OBJ, 0):
if ((oppl ^= (1<<SATURN)) & (1<<SATURN)) {
print_updating();
alt_body (SATURN, 1, &now);
} else
alt_nobody (SATURN);
break;
case rcfpack (R_URANUS, C_OBJ, 0):
if ((oppl ^= (1<<URANUS)) & (1<<URANUS)) {
print_updating();
alt_body (URANUS, 1, &now);
} else
alt_nobody (URANUS);
break;
case rcfpack (R_NEPTUNE, C_OBJ, 0):
if ((oppl ^= (1<<NEPTUNE)) & (1<<NEPTUNE)) {
print_updating();
alt_body (NEPTUNE, 1, &now);
} else
alt_nobody (NEPTUNE);
break;
case rcfpack (R_PLUTO, C_OBJ, 0):
if ((oppl ^= (1<<PLUTO)) & (1<<PLUTO)) {
print_updating();
alt_body (PLUTO, 1, &now);
} else
alt_nobody (PLUTO);
break;
case rcfpack (R_OBJX, C_OBJ, 0):
if ((oppl ^= (1<<OBJX)) & (1<<OBJX)) {
if (!bp) {
static char p[]= "object name (or RETURN for last again): ";
f_prompt (p);
if (read_line (buf, MAXOBJXNM) < 0) {
oppl ^= (1<<OBJX); /* turn back off */
break;
}
bp = buf;
}
if (bp[0] != '\0') {
/* initialize epoch to EOD to avoid initial precession. */
double e = mjd, tmp = 0.0;
objx_set (&tmp, &tmp, &e, bp);
}
objx_on();
alt_body (OBJX, 1, &now);
} else {
objx_off();
alt_nobody (OBJX);
}
break;
case rcfpack (R_OBJX, C_RA, 0):
if (oppl & (1<<OBJX)) {
double r, e;
if (!bp) {
static char p[] = "ra (current epoch, h:m:s): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
objx_get (&r, (double *)0, (double *)0, (char *)0);
f_dec_sexsign (radhr(r), °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &r);
r = hrrad(r);
e = (epoch == EOD) ? mjd : epoch;
objx_set (&r, (double *)0, &e, (char *)0);
alt_body (OBJX, 1, &now);
}
break;
case rcfpack (R_OBJX, C_DEC, 0):
if (oppl & (1<<OBJX)) {
double d, e;
if (!bp) {
static char p[] = "dec (current epoch, d:m:s): ";
f_prompt (p);
if (read_line (buf, PW-sizeof(p)) <= 0)
break;
bp = buf;
}
objx_get ((double *)0, &d, (double *)0, (char *)0);
f_dec_sexsign (raddeg(d), °hrs, &mins, &secs);
f_sscansex (bp, °hrs, &mins, &secs);
sex_dec (deghrs, mins, secs, &d);
d = degrad(d);
e = (epoch == EOD) ? mjd : epoch;
objx_set ((double *)0, &d, &e, (char *)0);
alt_body (OBJX, 1, &now);
}
break;
}
return (new);
}
static
print_tminc(force)
int force;
{
static double last;
if (force || tminc != last) {
if (tminc == RTC)
f_string (R_STPSZ, C_STPSZV, " RT CLOCK");
else if (fabs(tminc) >= 24.0)
f_double (R_STPSZ, C_STPSZV, "%6.4g dy", tminc/24.0);
else
f_signtime (R_STPSZ, C_STPSZV, tminc);
last = tminc;
}
}
static
print_bodies (force)
int force;
{
int p;
for (p = nxtbody(-1); p != -1; p = nxtbody(p))
if (oppl & (1<<p))
alt_body (p, force, &now);
}
static
clrall_bodies ()
{
int p;
for (p = nxtbody(-1); p != -1; p = nxtbody(p))
if (oppl & (1<<p))
alt_nobody (p);
}
print_updating()
{
f_prompt ("Updating...");
}
static
print_nstep(force)
int force;
{
static int last;
if (force || nstep != last) {
char buf[16];
sprintf (buf, "%8d", nstep);
f_string (R_NSTEP, C_NSTEPV, buf);
last = nstep;
}
}
xXx
echo x mainmenu.c
cat > mainmenu.c << 'xXx'
/* printing routines for the main (upper) screen.
*/
#include <stdio.h>
#include <math.h>
#include "astro.h"
#include "circum.h"
#include "screen.h"
mm_borders()
{
char line[NC+1], *lp;
register i;
lp = line;
for (i = 0; i < NC; i++)
*lp++ = '-';
*lp = '\0';
f_string (R_PLANTAB-1, 1, line);
for (i = R_TOP; i < R_PLANTAB-1; i++)
f_char (i, COL2-2, '|');
for (i = R_TOP; i < R_PLANTAB-1; i++)
f_char (i, COL3-2, '|');
for (i = R_LST; i < R_PLANTAB-1; i++)
f_char (i, COL4-2, '|');
}
mm_labels()
{
f_string (R_TZN, C_TZN, "LT");
f_string (R_UT, C_UT, "UTC");
f_string (R_JD, C_JD, "JulianDate");
f_string (R_WATCH, C_WATCH, "Watch");
f_string (R_SRCH, C_SRCH, "Search");
f_string (R_PLOT, C_PLOT, "Plot");
f_string (R_ALTM, C_ALTM, "Menu");
f_string (R_LST, C_LST, "LST");
f_string (R_DAWN, C_DAWN, "Dawn");
f_string (R_DUSK, C_DUSK, "Dusk");
f_string (R_LON, C_LON, "NiteLn");
f_string (R_NSTEP, C_NSTEP, "NStep");
f_string (R_STPSZ, C_STPSZ, "StpSz");
f_string (R_LAT, C_LAT, "Lat");
f_string (R_LONG, C_LONG, "Long");
f_string (R_HEIGHT, C_HEIGHT, "Elev");
f_string (R_TEMP, C_TEMP, "Temp");
f_string (R_PRES, C_PRES, "AtmPr");
f_string (R_TZONE, C_TZONE, "TZ");
f_string (R_EPOCH, C_EPOCH, "Epoch");
f_string (R_PLANTAB, C_OBJ, "Ob");
f_string (R_SUN, C_OBJ, "Su");
f_string (R_MOON, C_OBJ, "Mo");
f_string (R_MERCURY, C_OBJ, "Me");
f_string (R_VENUS, C_OBJ, "Ve");
f_string (R_MARS, C_OBJ, "Ma");
f_string (R_JUPITER, C_OBJ, "Ju");
f_string (R_SATURN, C_OBJ, "Sa");
f_string (R_URANUS, C_OBJ, "Ur");
f_string (R_NEPTUNE, C_OBJ, "Ne");
f_string (R_PLUTO, C_OBJ, "Pl");
}
/* print all the time/date/where related stuff: the Now structure.
* print in a nice order, based on the field locations, as much as possible.
*/
mm_now (np, all)
Now *np;
int all;
{
char buf[32];
double lmjd = mjd - tz/24.0;
double jd = mjd + 2415020L;
double tmp;
sprintf (buf, "%-3.3s", tznm);
f_string (R_TZN, C_TZN, buf);
f_time (R_LT, C_LT, mjd_hr(lmjd));
f_date (R_LD, C_LD, mjd_day(lmjd));
f_time (R_UT, C_UTV, mjd_hr(mjd));
f_date (R_UD, C_UD, mjd_day(mjd));
sprintf (buf, "%14.5f", jd);
(void) flog_log (R_JD, C_JDV, jd);
f_string (R_JD, C_JDV, buf);
now_lst (np, &tmp);
f_time (R_LST, C_LSTV, tmp);
if (all) {
f_gangle (R_LAT, C_LATV, lat);
f_gangle (R_LONG, C_LONGV, -lng); /* + west */
tmp = height * 2.093e7; /* want to see ft, not earth radii */
sprintf (buf, "%5g ft", tmp);
(void) flog_log (R_HEIGHT, C_HEIGHTV, tmp);
f_string (R_HEIGHT, C_HEIGHTV, buf);
tmp = 9./5.*temp + 32.0; /* want to see degrees F, not C */
(void) flog_log (R_TEMP, C_TEMPV, tmp);
sprintf (buf, "%6g F", tmp);
f_string (R_TEMP, C_TEMPV, buf);
tmp = pressure / 33.86; /* want to see in. Hg, not mBar */
(void) flog_log (R_PRES, C_PRESV, tmp);
sprintf (buf, "%5.2f in", tmp);
f_string (R_PRES, C_PRESV, buf);
f_signtime (R_TZONE, C_TZONEV, tz);
/* TODO: allow epoch to be plotted (?) */
if (epoch == EOD)
f_string (R_EPOCH, C_EPOCHV, "(OfDate)");
else {
mjd_year (epoch, &tmp);
f_double (R_EPOCH, C_EPOCHV, "%8.1f", tmp);
}
}
/* print the calendar for local day, if new month/year. */
mm_calendar (np, all > 1);
}
/* display dawn/dusk/length-of-night times.
*/
mm_twilight (np, force)
Now *np;
int force;
{
double dusk, dawn;
double tmp;
int status;
if (!twilight_cir (np, &dawn, &dusk, &status) && !force)
return;
switch (status) {
case -1: /* sun never sets today */
case 1: /* sun never rises today */
case 2: /* can not find where sun is! */
f_blanks (R_DAWN, C_DAWNV, 5);
f_blanks (R_DUSK, C_DUSKV, 5);
f_blanks (R_LON, C_LONV, 5);
return;
default: /* all ok */
;
}
f_mtime (R_DAWN, C_DAWNV, dawn);
f_mtime (R_DUSK, C_DUSKV, dusk);
tmp = dawn - dusk; range (&tmp, 24.0);
f_mtime (R_LON, C_LONV, tmp);
}
mm_newcir (y)
int y;
{
static char ncmsg[] = "NEW CIRCUMSTANCES";
static char nomsg[] = " ";
static int last_y = -1;
if (y != last_y) {
f_string (R_NEWCIR, C_NEWCIR, y ? ncmsg : nomsg);
last_y = y;
}
}
static
mm_calendar (np, force)
Now *np;
int force;
{
static char *mnames[] = {
"January", "February", "March", "April", "May", "June",
"July", "August", "September", "October", "November", "December"
};
static int last_m, last_y;
static double last_tz;
char str[64];
int m, y;
double d;
int f, nd;
int r;
double jd0;
/* get local m/d/y. do nothing if still same month and not forced. */
mjd_cal (mjd_day(mjd-tz/24.0), &m, &d, &y);
if (m == last_m && y == last_y && tz == last_tz && !force)
return;
last_m = m;
last_y = y;
last_tz = tz;
/* find day of week of first day of month */
cal_mjd (m, 1.0, y, &jd0);
mjd_dow (jd0, &f);
if (f < 0) {
/* can't figure it out - too hard before Gregorian */
int i;
for (i = 8; --i >= 0; )
f_string (R_CAL+i, C_CAL, " ");
return;
}
/* print header */
f_blanks (R_CAL, C_CAL, 20);
sprintf (str, "%s %4d", mnames[m-1], y);
f_string (R_CAL, C_CAL + (20 - (strlen(mnames[m-1]) + 5))/2, str);
f_string (R_CAL+1, C_CAL, "Su Mo Tu We Th Fr Sa");
/* find number of days in this month */
mjd_dpm (jd0, &nd);
/* print the calendar */
for (r = 0; r < 6; r++) {
char row[7*3+1], *rp = row;
int c;
for (c = 0; c < 7; c++) {
int i = r*7+c;
if (i < f || i >= f + nd)
sprintf (rp, " ");
else
sprintf (rp, "%2d ", i-f+1);
rp += 3;
}
row[sizeof(row)-2] = '\0'; /* don't print last blank; causes wrap*/
f_string (R_CAL+2+r, C_CAL, row);
}
/* over print the new and full moons for this month.
* TODO: don't really know which dates to use here (see moonnf())
* so try several to be fairly safe. have to go back to 4/29/1988
* to find the full moon on 5/1 for example.
*/
mm_nfmoon (jd0-3, tz, m, f);
mm_nfmoon (jd0+15, tz, m, f);
}
static
mm_nfmoon (jd, tzone, m, f)
double jd, tzone;
int m, f;
{
static char nm[] = "NM", fm[] = "FM";
double dm;
int mm, ym;
double jdn, jdf;
int di;
moonnf (jd, &jdn, &jdf);
mjd_cal (jdn-tzone/24.0, &mm, &dm, &ym);
if (m == mm) {
di = dm + f - 1;
f_string (R_CAL+2+di/7, C_CAL+3*(di%7), nm);
}
mjd_cal (jdf-tzone/24.0, &mm, &dm, &ym);
if (m == mm) {
di = dm + f - 1;
f_string (R_CAL+2+di/7, C_CAL+3*(di%7), fm);
}
}
xXx
echo x moon.c
cat > moon.c << 'xXx'
#include <stdio.h>
#include <math.h>
#include "astro.h"
/* given the mjd, find the geocentric ecliptic longitude, lam, and latitude,
* bet, and horizontal parallax, hp for the moon.
* N.B. series for long and lat are good to about 10 and 3 arcseconds. however,
* math errors cause up to 100 and 30 arcseconds error, even if use double.
* why?? suspect highly sensitive nature of difference used to get m1..6.
* N.B. still need to correct for nutation. then for topocentric location
* further correct for parallax and refraction.
*/
moon (mjd, lam, bet, hp)
double mjd;
double *lam, *bet, *hp;
{
double t, t2;
double ld;
double ms;
double md;
double de;
double f;
double n;
double a, sa, sn, b, sb, c, sc, e, e2, l, g, w1, w2;
double m1, m2, m3, m4, m5, m6;
t = mjd/36525.;
t2 = t*t;
m1 = mjd/27.32158213;
m1 = 360.0*(m1-(int)m1);
m2 = mjd/365.2596407;
m2 = 360.0*(m2-(int)m2);
m3 = mjd/27.55455094;
m3 = 360.0*(m3-(int)m3);
m4 = mjd/29.53058868;
m4 = 360.0*(m4-(int)m4);
m5 = mjd/27.21222039;
m5 = 360.0*(m5-(int)m5);
m6 = mjd/6798.363307;
m6 = 360.0*(m6-(int)m6);
ld = 270.434164+m1-(.001133-.0000019*t)*t2;
ms = 358.475833+m2-(.00015+.0000033*t)*t2;
md = 296.104608+m3+(.009192+.0000144*t)*t2;
de = 350.737486+m4-(.001436-.0000019*t)*t2;
f = 11.250889+m5-(.003211+.0000003*t)*t2;
n = 259.183275-m6+(.002078+.000022*t)*t2;
a = degrad(51.2+20.2*t);
sa = sin(a);
sn = sin(degrad(n));
b = 346.56+(132.87-.0091731*t)*t;
sb = .003964*sin(degrad(b));
c = degrad(n+275.05-2.3*t);
sc = sin(c);
ld = ld+.000233*sa+sb+.001964*sn;
ms = ms-.001778*sa;
md = md+.000817*sa+sb+.002541*sn;
f = f+sb-.024691*sn-.004328*sc;
de = de+.002011*sa+sb+.001964*sn;
e = 1-(.002495+7.52e-06*t)*t;
e2 = e*e;
ld = degrad(ld);
ms = degrad(ms);
n = degrad(n);
de = degrad(de);
f = degrad(f);
md = degrad(md);
l = 6.28875*sin(md)+1.27402*sin(2*de-md)+.658309*sin(2*de)+
.213616*sin(2*md)-e*.185596*sin(ms)-.114336*sin(2*f)+
.058793*sin(2*(de-md))+.057212*e*sin(2*de-ms-md)+
.05332*sin(2*de+md)+.045874*e*sin(2*de-ms)+.041024*e*sin(md-ms);
l = l-.034718*sin(de)-e*.030465*sin(ms+md)+.015326*sin(2*(de-f))-
.012528*sin(2*f+md)-.01098*sin(2*f-md)+.010674*sin(4*de-md)+
.010034*sin(3*md)+.008548*sin(4*de-2*md)-e*.00791*sin(ms-md+2*de)-
e*.006783*sin(2*de+ms);
l = l+.005162*sin(md-de)+e*.005*sin(ms+de)+.003862*sin(4*de)+
e*.004049*sin(md-ms+2*de)+.003996*sin(2*(md+de))+
.003665*sin(2*de-3*md)+e*.002695*sin(2*md-ms)+
.002602*sin(md-2*(f+de))+e*.002396*sin(2*(de-md)-ms)-
.002349*sin(md+de);
l = l+e2*.002249*sin(2*(de-ms))-e*.002125*sin(2*md+ms)-
e2*.002079*sin(2*ms)+e2*.002059*sin(2*(de-ms)-md)-
.001773*sin(md+2*(de-f))-.001595*sin(2*(f+de))+
e*.00122*sin(4*de-ms-md)-.00111*sin(2*(md+f))+.000892*sin(md-3*de);
l = l-e*.000811*sin(ms+md+2*de)+e*.000761*sin(4*de-ms-2*md)+
e2*.000704*sin(md-2*(ms+de))+e*.000693*sin(ms-2*(md-de))+
e*.000598*sin(2*(de-f)-ms)+.00055*sin(md+4*de)+.000538*sin(4*md)+
e*.000521*sin(4*de-ms)+.000486*sin(2*md-de);
l = l+e2*.000717*sin(md-2*ms);
*lam = ld+degrad(l);
range (lam, 2*PI);
g = 5.12819*sin(f)+.280606*sin(md+f)+.277693*sin(md-f)+
.173238*sin(2*de-f)+.055413*sin(2*de+f-md)+.046272*sin(2*de-f-md)+
.032573*sin(2*de+f)+.017198*sin(2*md+f)+.009267*sin(2*de+md-f)+
.008823*sin(2*md-f)+e*.008247*sin(2*de-ms-f);
g = g+.004323*sin(2*(de-md)-f)+.0042*sin(2*de+f+md)+
e*.003372*sin(f-ms-2*de)+e*.002472*sin(2*de+f-ms-md)+
e*.002222*sin(2*de+f-ms)+e*.002072*sin(2*de-f-ms-md)+
e*.001877*sin(f-ms+md)+.001828*sin(4*de-f-md)-e*.001803*sin(f+ms)-
.00175*sin(3*f);
g = g+e*.00157*sin(md-ms-f)-.001487*sin(f+de)-e*.001481*sin(f+ms+md)+
e*.001417*sin(f-ms-md)+e*.00135*sin(f-ms)+.00133*sin(f-de)+
.001106*sin(f+3*md)+.00102*sin(4*de-f)+.000833*sin(f+4*de-md)+
.000781*sin(md-3*f)+.00067*sin(f+4*de-2*md);
g = g+.000606*sin(2*de-3*f)+.000597*sin(2*(de+md)-f)+
e*.000492*sin(2*de+md-ms-f)+.00045*sin(2*(md-de)-f)+
.000439*sin(3*md-f)+.000423*sin(f+2*(de+md))+
.000422*sin(2*de-f-3*md)-e*.000367*sin(ms+f+2*de-md)-
e*.000353*sin(ms+f+2*de)+.000331*sin(f+4*de);
g = g+e*.000317*sin(2*de+f-ms+md)+e2*.000306*sin(2*(de-ms)-f)-
.000283*sin(md+3*f);
w1 = .0004664*cos(n);
w2 = .0000754*cos(c);
*bet = degrad(g)*(1-w1-w2);
*hp = .950724+.051818*cos(md)+.009531*cos(2*de-md)+.007843*cos(2*de)+
.002824*cos(2*md)+.000857*cos(2*de+md)+e*.000533*cos(2*de-ms)+
e*.000401*cos(2*de-md-ms)+e*.00032*cos(md-ms)-.000271*cos(de)-
e*.000264*cos(ms+md)-.000198*cos(2*f-md);
*hp = *hp+.000173*cos(3*md)+.000167*cos(4*de-md)-e*.000111*cos(ms)+
.000103*cos(4*de-2*md)-.000084*cos(2*md-2*de)-
e*.000083*cos(2*de+ms)+.000079*cos(2*de+2*md)+.000072*cos(4*de)+
e*.000064*cos(2*de-ms+md)-e*.000063*cos(2*de+ms-md)+
e*.000041*cos(ms+de);
*hp = *hp+e*.000035*cos(2*md-ms)-.000033*cos(3*md-2*de)-
.00003*cos(md+de)-.000029*cos(2*(f-de))-e*.000029*cos(2*md+ms)+
e2*.000026*cos(2*(de-ms))-.000023*cos(2*(f-de)+md)+
e*.000019*cos(4*de-ms-md);
*hp = degrad(*hp);
}
xXx
echo x moonnf.c
cat > moonnf.c << 'xXx'
#include <stdio.h>
#include <math.h>
#include "astro.h"
#define unw(w,z) ((w)-floor((w)/(z))*(z))
/* given a modified Julian date, mjd, return the mjd of the new
* and full moons about then, mjdn and mjdf.
* TODO: exactly which ones does it find? eg:
* 5/28/1988 yields 5/15 and 5/31
* 5/29 6/14 6/29
*/
moonnf (mjd, mjdn, mjdf)
double mjd;
double *mjdn, *mjdf;
{
int mo, yr;
double dy;
double mjd0;
double k, tn, tf, t;
mjd_cal (mjd, &mo, &dy, &yr);
cal_mjd (1, 0., yr, &mjd0);
k = (yr-1900+((mjd-mjd0)/365))*12.3685;
k = floor(k+0.5);
tn = k/1236.85;
tf = (k+0.5)/1236.85;
t = tn;
m (t, k, mjdn);
t = tf;
k += 0.5;
m (t, k, mjdf);
}
static
m (t, k, mjd)
double t, k;
double *mjd;
{
double t2, a, a1, b, b1, c, ms, mm, f, ddjd;
t2 = t*t;
a = 29.53*k;
c = degrad(166.56+(132.87-9.173e-3*t)*t);
b = 5.8868e-4*k+(1.178e-4-1.55e-7*t)*t2+3.3e-4*sin(c)+7.5933E-1;
ms = 359.2242+360*unw(k/1.236886e1,1)-(3.33e-5+3.47e-6*t)*t2;
mm = 306.0253+360*unw(k/9.330851e-1,1)+(1.07306e-2+1.236e-5*t)*t2;
f = 21.2964+360*unw(k/9.214926e-1,1)-(1.6528e-3+2.39e-6*t)*t2;
ms = unw(ms,360);
mm = unw(mm,360);
f = unw(f,360);
ms = degrad(ms);
mm = degrad(mm);
f = degrad(f);
ddjd = (1.734e-1-3.93e-4*t)*sin(ms)+2.1e-3*sin(2*ms)
-4.068e-1*sin(mm)+1.61e-2*sin(2*mm)-4e-4*sin(3*mm)
+1.04e-2*sin(2*f)-5.1e-3*sin(ms+mm)-7.4e-3*sin(ms-mm)
+4e-4*sin(2*f+ms)-4e-4*sin(2*f-ms)-6e-4*sin(2*f+mm)
+1e-3*sin(2*f-mm)+5e-4*sin(ms+2*mm);
a1 = (int)a;
b = b+ddjd+(a-a1);
b1 = (int)b;
a = a1+b1;
b = b-b1;
*mjd = a + b;
}
xXx
echo x nutation.c
cat > nutation.c << 'xXx'
#include <stdio.h>
#include <math.h>
#include "astro.h"
/* given the modified JD, mjd, find the nutation in obliquity, *deps, and
* the nutation in longitude, *dpsi, each in radians.
*/
nutation (mjd, deps, dpsi)
double mjd;
double *deps, *dpsi;
{
static double lastmjd, lastdeps, lastdpsi;
double ls, ld; /* sun's mean longitude, moon's mean longitude */
double ms, md; /* sun's mean anomaly, moon's mean anomaly */
double nm; /* longitude of moon's ascending node */
double t, t2; /* number of Julian centuries of 36525 days since
* Jan 0.5 1900.
*/
double tls, tnm, tld; /* twice above */
double a, b; /* temps */
if (mjd == lastmjd) {
*deps = lastdeps;
*dpsi = lastdpsi;
return;
}
t = mjd/36525.;
t2 = t*t;
a = 100.0021358*t;
b = 360.*(a-(int)a);
ls = 279.697+.000303*t2+b;
a = 1336.855231*t;
b = 360.*(a-(int)a);
ld = 270.434-.001133*t2+b;
a = 99.99736056000026*t;
b = 360.*(a-(int)a);
ms = 358.476-.00015*t2+b;
a = 13255523.59*t;
b = 360.*(a-(int)a);
md = 296.105+.009192*t2+b;
a = 5.372616667*t;
b = 360.*(a-(int)a);
nm = 259.183+.002078*t2-b;
/* convert to radian forms for use with trig functions.
*/
tls = 2*degrad(ls);
nm = degrad(nm);
tnm = 2*degrad(nm);
ms = degrad(ms);
tld = 2*degrad(ld);
md = degrad(md);
/* find delta psi and eps, in arcseconds.
*/
lastdpsi = (-17.2327-.01737*t)*sin(nm)+(-1.2729-.00013*t)*sin(tls)
+.2088*sin(tnm)-.2037*sin(tld)+(.1261-.00031*t)*sin(ms)
+.0675*sin(md)-(.0497-.00012*t)*sin(tls+ms)
-.0342*sin(tld-nm)-.0261*sin(tld+md)+.0214*sin(tls-ms)
-.0149*sin(tls-tld+md)+.0124*sin(tls-nm)+.0114*sin(tld-md);
lastdeps = (9.21+.00091*t)*cos(nm)+(.5522-.00029*t)*cos(tls)
-.0904*cos(tnm)+.0884*cos(tld)+.0216*cos(tls+ms)
+.0183*cos(tld-nm)+.0113*cos(tld+md)-.0093*cos(tls-ms)
-.0066*cos(tls-nm);
/* convert to radians.
*/
lastdpsi = degrad(lastdpsi/3600);
lastdeps = degrad(lastdeps/3600);
lastmjd = mjd;
*deps = lastdeps;
*dpsi = lastdpsi;
}
xXx
echo x objx.c
cat > objx.c << 'xXx'
/* functions to save the user-definable object.
* this way, once set, the object can be quieried for position just like the
* other bodies. also, someday we can use oribital elements to let object-x
* be any solar system body too.
*/
#include <math.h>
#include "astro.h"
#include "circum.h"
#include "screen.h"
static double objx_ra, objx_dec, objx_epoch;
static char objx_name[MAXOBJXNM+1] = "";
static int objx_onflag; /* !0 while object x is active */
/* set attributes of object x.
* use pointers so we can set just some attributes as desired.
*/
objx_set (r, d, e, name)
double *r, *d, *e;
char *name;
{
if (r) objx_ra = *r;
if (d) objx_dec = *d;
if (e) objx_epoch = *e;
if (name) strncpy (objx_name, name, MAXOBJXNM);
}
/* return those attributes of interest for object x.
* always return a 2-char name.
*/
objx_get (r, d, e, name)
double *r, *d, *e;
char *name;
{
if (r) *r = objx_ra;
if (d) *d = objx_dec;
if (e) *e = objx_epoch;
if (name) {
name[0] = name[1] = ' ';
strcpy (name, objx_name); /* includes trailing 0 */
}
}
/* turn "on" or "off" but don't forget facts about object x.
*/
objx_on ()
{
objx_onflag = 1;
}
objx_off()
{
objx_onflag = 0;
}
/* return true if object is now on, else 0.
*/
objx_ison()
{
return (objx_onflag);
}
/* fill in sp with all we can about object x. */
/* ARGSUSED */
objx_cir (as, np, sp)
double as; /* desired, accuracy, in arc seconds. ignored for now */
Now *np;
Sky *sp;
{
double xr, xd, xe;
double lst, alt, az;
double lsn, rsn, bet, lam, el;
objx_get (&xr, &xd, &xe, (char *)0);
precess (xe, epoch==EOD ? mjd : epoch, &xr, &xd);
sp->s_ra = xr;
sp->s_dec = xd;
now_lst (np, &lst);
hadec_aa (lat, hrrad(lst) - xr, xd, &alt, &az);
refract (pressure, temp, alt, &alt);
sp->s_alt = alt;
sp->s_az = az;
/* elongation */
sunpos (mjd, &lsn, &rsn);
range (&lsn, 2*PI);
eq_ecl (mjd, xr, xd, &bet, &lam);
elongation (lam, bet, lsn, &el);
sp->s_elong = raddeg (el);
/* TODO: not always new really */
return (1);
}
xXx
echo x obliq.c
cat > obliq.c << 'xXx'
#include <stdio.h>
#include "astro.h"
/* given the modified Julian date, mjd, find the obliquity of the
* ecliptic, *eps, in radians.
*/
obliquity (mjd, eps)
double mjd;
double *eps;
{
static double lastmjd, lasteps;
if (mjd != lastmjd) {
double t;
t = mjd/36525.;
lasteps = degrad(2.345229444E1 - ((((-1.81E-3*t)+5.9E-3)*t+4.6845E1)*t)/3600.0);
lastmjd = mjd;
}
*eps = lasteps;
}
xXx
echo x parallax.c
cat > parallax.c << 'xXx'
#include <stdio.h>
#include <math.h>
#include "astro.h"
/* given true ha and dec, tha and tdec, the geographical latitude, phi, the
* height above sea-level (as a fraction of the earths radius, 6378.16km),
* ht, and the equatorial horizontal parallax, ehp, find the apparent
* ha and dec, aha and adec allowing for parallax.
* all angles in radians. ehp is the angle subtended at the body by the
* earth's equator.
*/
ta_par (tha, tdec, phi, ht, ehp, aha, adec)
double tha, tdec, phi, ht, ehp;
double *aha, *adec;
{
static double last_phi, last_ht, rsp, rcp;
double rp; /* distance to object in Earth radii */
double ctha;
double stdec, ctdec;
double tdtha, dtha;
double caha;
/* avoid calcs involving the same phi and ht */
if (phi != last_phi || ht != last_ht) {
double cphi, sphi, u;
cphi = cos(phi);
sphi = sin(phi);
u = atan(9.96647e-1*sphi/cphi);
rsp = (9.96647e-1*sin(u))+(ht*sphi);
rcp = cos(u)+(ht*cphi);
last_phi = phi;
last_ht = ht;
}
rp = 1/sin(ehp);
ctha = cos(tha);
stdec = sin(tdec);
ctdec = cos(tdec);
tdtha = (rcp*sin(tha))/((rp*ctdec)-(rcp*ctha));
dtha = atan(tdtha);
*aha = tha+dtha;
caha = cos(*aha);
range (aha, 2*PI);
*adec = atan(caha*(rp*stdec-rsp)/(rp*ctdec*ctha-rcp));
}
/* given the apparent ha and dec, aha and adec, the geographical latitude, phi,
* the height above sea-level (as a fraction of the earths radius, 6378.16km),
* ht, and the equatorial horizontal parallax, ehp, find the true ha and dec,
* tha and tdec allowing for parallax.
* all angles in radians. ehp is the angle subtended at the body by the
* earth's equator.
* uses ta_par() iteratively: find a set of true ha/dec that converts back
* to the given apparent ha/dec.
*/
at_par (aha, adec, phi, ht, ehp, tha, tdec)
double aha, adec, phi, ht, ehp;
double *tha, *tdec;
{
double nha, ndec; /* ha/dec corres. to current true guesses */
double eha, edec; /* error in ha/dec */
/* first guess for true is just the apparent */
*tha = aha;
*tdec = adec;
while (1) {
ta_par (*tha, *tdec, phi, ht, ehp, &nha, &ndec);
eha = aha - nha;
edec = adec - ndec;
if (fabs(eha)<1e-6 && fabs(edec)<1e-6)
break;
*tha += eha;
*tdec += edec;
}
}
xXx
echo x pelement.c
cat > pelement.c << 'xXx'
#include <stdio.h>
#include <math.h>
#include "astro.h"
/* this array contains polynomial coefficients to find the various orbital
* elements for the mean orbit at any instant in time for each major planet.
* the first five elements are in the form a0 + a1*t + a2*t**2 + a3*t**3,
* where t is the number of Julian centuries of 36525 Julian days since 1900
* Jan 0.5. the last three elements are constants.
*
* the orbital element (column) indeces are:
* [ 0- 3]: coefficients for mean longitude, in degrees;
* [ 4- 7]: coefficients for longitude of the perihelion, in degrees;
* [ 8-11]: coefficients for eccentricity;
* [12-15]: coefficients for inclination, in degrees;
* [16-19]: coefficients for longitude of the ascending node, in degrees;
* [20]: semi-major axis, in AU;
* [21]: angular diameter at 1 AU, in arcsec;
* [22]: standard visual magnitude, ie, the visual magnitude of the planet
* when at a distance of 1 AU from both the Sun and the Earth and
* with zero phase angle.
*
* the planent (row) indeces are:
* [0]: Mercury; [1]: Venus; [2]: Mars; [3]: Jupiter; [4]: Saturn;
* [5]: Uranus; [6]: Neptune; [7]: Pluto.
*/
#define NPELE (5*4 + 3) /* 4 coeffs for ea of 5 elems, + 3 constants */
static double elements[8][NPELE] = {
{ /* mercury... */
178.179078, 415.2057519, 3.011e-4, 0.0,
75.899697, 1.5554889, 2.947e-4, 0.0,
.20561421, 2.046e-5, 3e-8, 0.0,
7.002881, 1.8608e-3, -1.83e-5, 0.0,
47.145944, 1.1852083, 1.739e-4, 0.0,
.3870986, 6.74, -0.42
},
{ /* venus... */
342.767053, 162.5533664, 3.097e-4, 0.0,
130.163833, 1.4080361, -9.764e-4, 0.0,
6.82069e-3, -4.774e-5, 9.1e-8, 0.0,
3.393631, 1.0058e-3, -1e-6, 0.0,
75.779647, .89985, 4.1e-4, 0.0,
.7233316, 16.92, -4.4
},
{ /* mars... */
293.737334, 53.17137642, 3.107e-4, 0.0,
3.34218203e2, 1.8407584, 1.299e-4, -1.19e-6,
9.33129e-2, 9.2064e-5, 7.7e-8, 0.0,
1.850333, -6.75e-4, 1.26e-5, 0.0,
48.786442, .7709917, -1.4e-6, -5.33e-6,
1.5236883, 9.36, -1.52
},
{ /* jupiter... */
238.049257, 8.434172183, 3.347e-4, -1.65e-6,
1.2720972e1, 1.6099617, 1.05627e-3, -3.43e-6,
4.833475e-2, 1.6418e-4, -4.676e-7, -1.7e-9,
1.308736, -5.6961e-3, 3.9e-6, 0.0,
99.443414, 1.01053, 3.5222e-4, -8.51e-6,
5.202561, 196.74, -9.4
},
{ /* saturn... */
266.564377, 3.398638567, 3.245e-4, -5.8e-6,
9.1098214e1, 1.9584158, 8.2636e-4, 4.61e-6,
5.589232e-2, -3.455e-4, -7.28e-7, 7.4e-10,
2.492519, -3.9189e-3, -1.549e-5, 4e-8,
112.790414, .8731951, -1.5218e-4, -5.31e-6,
9.554747, 165.6, -8.88
},
{ /* uranus... */
244.19747, 1.194065406, 3.16e-4, -6e-7,
1.71548692e2, 1.4844328, 2.372e-4, -6.1e-7,
4.63444e-2, -2.658e-5, 7.7e-8, 0.0,
.772464, 6.253e-4, 3.95e-5, 0.0,
73.477111, .4986678, 1.3117e-3, 0.0,
19.21814, 65.8, -7.19
},
{ /* neptune... */
84.457994, .6107942056, 3.205e-4, -6e-7,
4.6727364e1, 1.4245744, 3.9082e-4, -6.05e-7,
8.99704e-3, 6.33e-6, -2e-9, 0.0,
1.779242, -9.5436e-3, -9.1e-6, 0.0,
130.681389, 1.098935, 2.4987e-4, -4.718e-6,
30.10957, 62.2, -6.87
},
{ /* pluto...(osculating 1984 jan 21) */
95.3113544, .3980332167, 0.0, 0.0,
224.017, 0.0, 0.0, 0.0,
.25515, 0.0, 0.0, 0.0,
17.1329, 0.0, 0.0, 0.0,
110.191, 0.0, 0.0, 0.0,
39.8151, 8.2, -1.0
}
};
/* given a modified Julian date, mjd, return the elements for the mean orbit
* at that instant of all the major planets, together with their
* mean daily motions in longitude, angular diameter and standard visual
* magnitude.
* plan[i][j] contains all the values for all the planets at mjd, such that
* i = 0..7: mercury, venus, mars, jupiter, saturn, unranus, neptune, pluto;
* j = 0..8: mean longitude, mean daily motion in longitude, longitude of
* the perihelion, eccentricity, inclination, longitude of the ascending
* node, length of the semi-major axis, angular diameter from 1 AU, and
* the standard visual magnitude (see elements[][] comment, above).
*/
pelement (mjd, plan)
double mjd;
double plan[8][9];
{
register double *ep, *pp;
register double t = mjd/36525.;
double aa;
int planet, i;
for (planet = 0; planet < 8; planet++) {
ep = elements[planet];
pp = plan[planet];
aa = ep[1]*t;
pp[0] = ep[0] + 360.*(aa-(int)aa) + (ep[3]*t + ep[2])*t*t;
range (pp, 360.);
pp[1] = (ep[1]*9.856263e-3) + (ep[2] + ep[3])/36525;
for (i = 4; i < 20; i += 4)
pp[i/4+1] = ((ep[i+3]*t + ep[i+2])*t + ep[i+1])*t + ep[i+0];
pp[6] = ep[20];
pp[7] = ep[21];
pp[8] = ep[22];
}
}
xXx
More information about the Comp.sources.misc
mailing list