Usenet 1994 January/usenetsourcesnewsgroupsinfomagicjanuary1994.iso/sources/misc/volume28/ephem/part05

Text File | 1992-03-15 | 55KB | 2,001 lines

Newsgroups: comp.sources.misc
From: [e-mail] (Elwood C. Downey)
Subject:  v28i088:  ephem - an interactive astronomical ephemeris, v4.28, Part05/09
Message-ID: <[e-mail]>
X-Md4-Signature: 4ed4c4e87680e897b120410d7e931830
Date: Tue, 10 Mar 1992 21:58:23 GMT
Approved: [e-mail]

Submitted-by: [e-mail] (Elwood C. Downey)
Posting-number: Volume 28, Issue 88
Archive-name: ephem/part05
Environment: UNIX, VMS, DOS, MAC
Supersedes: ephem-4.21: Volume 14, Issue 76-81

#! /bin/sh
# into a shell via "sh file" or similar.  To overwrite existing files,
# type "sh file -c".
# The tool that generated this appeared in the comp.sources.unix newsgroup;
# send mail to [e-mail] if you want that tool.
# Contents:  ephem.cfg main.c mainmenu.c plans.c
# Wrapped by kent@sparky on Tue Mar 10 14:34:07 1992
PATH=/bin:/usr/bin:/usr/ucb ; export PATH
echo If this archive is complete, you will see the following message:
echo '          "shar: End of archive 5 (of 9)."'
if test -f 'ephem.cfg' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'ephem.cfg'\"
else
  echo shar: Extracting \"'ephem.cfg'\" \(168 characters\)
  sed "s/^X//" >'ephem.cfg' <<'END_OF_FILE'
XUT=0:0:0
XUD=5/1/1990
XTZNAME=CDT
XTZONE=5
XLONG=93:42:8
XLAT=44:50:37
XHEIGHT=800
XTEMP=40
XPRES=29.5
XSTPSZ=RTC
XPROPTS=TSMevmjsunpxy
XEPOCH=2000
XNSTEP=1
X
XOBJX=Austin
XOBJY=Juno
END_OF_FILE
  if test 168 -ne `wc -c <'ephem.cfg'`; then
    echo shar: \"'ephem.cfg'\" unpacked with wrong size!
  fi
  # end of 'ephem.cfg'
fi
if test -f 'main.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'main.c'\"
else
  echo shar: Extracting \"'main.c'\" \(26117 characters\)
  sed "s/^X//" >'main.c' <<'END_OF_FILE'
X/* main "ephem" program.
X * -------------------------------------------------------------------
X * Copyright (c) 1990,1991,1992 by Elwood Charles Downey
X *
X * Permission is granted to make and distribute copies of this program
X * free of charge, provided the copyright notice and this permission
X * notice are preserved on all copies.  All other rights reserved.
X * -------------------------------------------------------------------
X * set options.
X * init screen and circumstances.
X * enter infinite loop updating screen and allowing operator input.
X */
X
X#include <stdio.h>
X#include <ctype.h>
X#include <signal.h>
X#include <setjmp.h>
X#include <math.h>
X#ifdef VMS
X#include <stdlib.h>
X#endif
X#include "astro.h"
X#include "circum.h"
X#include "screen.h"
X
Xextern char *strncpy();
Xextern char *getenv();
X
X/* shorthands for fields of a Now structure, now.
X * first undo the ones for a Now pointer from circum.h.
X */
X#undef mjd
X#undef lat
X#undef lng
X#undef tz
X#undef temp
X#undef pressure
X#undef height
X#undef epoch
X#undef tznm
X
X#define mjd    now.n_mjd
X#define lat    now.n_lat
X#define lng    now.n_lng
X#define tz    now.n_tz
X#define temp    now.n_temp
X#define pressure now.n_pressure
X#define height    now.n_height
X#define epoch    now.n_epoch
X#define tznm    now.n_tznm
X
Xstatic jmp_buf fpe_err_jmp;    /* used to recover from SIGFPE */
Xstatic char *cfgfile;        /* !0 if -c used */
Xstatic char cfgdef[] = "ephem.cfg"; /* default configuration file name */
Xstatic Now now;        /* where when and how, right now */
Xstatic double tminc;    /* hrs to inc time by each loop; RTC means use clock */
Xstatic int nstep;    /* steps to go before stopping */
Xstatic int spause;    /* secs to pause between steps */
Xstatic int optwi;    /* set when want to display dawn/dusk/len-of-night */
Xstatic int oppl;    /* mask of (1<<planet) bits; set when want to show it */
X
Xmain (ac, av)
Xint ac;
Xchar *av[];
X{
X    void bye();
X    void on_fpe();
X    static char freerun[] =
X        "Running... press any key to stop to make changes.";
X    static char prmpt[] =
X"Move to another field, RETURN to change this field, ? for help, or q to run";
X    static char hlp[] =
X    "arrow keys move to field; any key stops running; ^d exits; ^l redraws";
X    int fld = rcfpack(R_NSTEP, C_NSTEPV, 0); /* initial cursor loc */
X    int sflag = 0;    /* not silent, by default */
X    int one = 1;    /* use a variable so optimizer doesn't get disabled */
X    int srchdone = 0; /* true when search funcs say so */
X    int newcir = 2;    /* set when circumstances change - means don't tminc */
X
X    while ((--ac > 0) && (**++av == '-')) {
X        char *s;
X        for (s = *av+1; *s != '\0'; s++)
X        switch (*s) {
X        case 's': /* no credits "silent" (don't publish this) */
X            sflag++;
X            break;
X        case 'c': /* set name of config file to use */
X            if (--ac <= 0) usage("-c but no config file");
X            cfgfile = *++av;
X            break;
X        case 'd': /* set alternate database file name */
X            if (--ac <= 0) usage("-d but no database file");
X            obj_setdbfilename (*++av);
X            break;
X        default:
X            usage("Bad - option");
X        }
X    }
X
X    if (!sflag)
X        credits();
X
X    /* fresh screen.
X     * crack config file, THEN args so args may override.
X     */
X    c_erase();
X    read_cfgfile ();
X    read_fieldargs (ac, av);
X
X    /* set up to clean up screen and tty if interrupted.
X     * also set up to stop if get floating error.
X     */
X    (void) signal (SIGINT, bye);
X    (void) signal (SIGFPE, on_fpe);
X
X    /* update screen forever (until QUIT) */
X    while (one) {
X
X        /* if get a floating error, longjmp() here and stop looping */
X        if (setjmp (fpe_err_jmp))
X        nstep = 0;
X        else {
X        nstep -= 1;
X
X        /* recalculate everything and update all the fields */
X        redraw_screen (newcir);
X        mm_newcir (0);
X
X        /* let searching functions change tminc and check for done */
X        srchdone = srch_eval (mjd, &tminc) < 0;
X        print_tminc(0);    /* to show possibly new search increment */
X
X        /* update plot and listing files, now that all fields are up
X         * to date and search function has been evaluated.
X         */
X        plot();
X        listing();
X
X        /* handle spause if we are really looping */
X        if (nstep > 0)
X            slp_sync();
X        }
X
X        /* stop loop to allow op to change parameters:
X         * if a search evaluation converges (or errors out),
X         * or if steps are done,
X         * or if op hits any key.
X         */
X        newcir = 0;
X        if (srchdone || nstep <= 0 || (chk_char()==0 && read_char()!=0)) {
X        int nfld;
X
X        /* update screen with the current stuff if stopped during
X         * unattended plotting or listing since last redraw_screen()
X         * didn't.
X         */
X        if ((plot_ison() || listing_ison()) && nstep > 0)
X            redraw_screen (1);
X
X        /* return nstep to default of 1 */
X        if (nstep <= 0) {
X            nstep = 1;
X            print_nstep (0);
X        }
X
X        /* change fields until END.
X         * update all time fields if any are changed
X         * and print NEW CIRCUMSTANCES if any have changed.
X         * QUIT causes bye() to be called and we never return.
X         */
X        while(nfld = sel_fld(fld,alt_menumask()|F_CHG,prmpt,hlp)) {
X            if (chg_fld ((char *)0, &nfld)) {
X            mm_now (&now, 1);
X            mm_newcir(1);
X            newcir = 1;
X            }
X            fld = nfld;
X        }
X        if (nstep > 1)
X            f_prompt (freerun);
X        }
X
X        /* increment time only if op didn't change cirumstances */
X        if (!newcir)
X        inc_mjd (&now, tminc);
X    }
X
X    return (0);
X}
X
X/* read in ephem's configuration file, if any.
X * if errors in file, call usage() (which exits).
X * if use -d, require it; else try $EPHEMCFG and ephem.cfg but don't
X *   complain if can't find these since, after all, one is not required.
X * skip all lines that doesn't begin with an alpha char.
X */
Xstatic
Xread_cfgfile()
X{
X    char buf[128];
X    FILE *fp;
X    char *fn;
X
X    /* open the config file.
X     * only REQUIRED if used -d option.
X     * if succcessful, fn points to file name.
X     */
X    if (cfgfile) {
X        fn = cfgfile;
X        fp = fopen (fn, "r");
X        if (!fp) {
X        (void) sprintf (buf, "Can not open %s", fn);
X        usage (buf);    /* does not return */
X        }
X    } else {
X        fn = getenv ("EPHEMCFG");
X        if (!fn)
X        fn = cfgdef;
X    }
X    fp = fopen (fn, "r");
X    if (!fp)
X        return;    /* oh well; after all, it's not required */
X
X    while (fgets (buf, sizeof(buf), fp)) {
X        if (!isalpha(buf[0]))
X        continue;
X        buf[strlen(buf)-1] = '\0';        /* discard trailing \n */
X        if (crack_fieldset (buf) < 0) {
X        char why[NC];
X        (void) sprintf (why, "Bad field spec in %s: %s\n", fn, buf);
X        usage (why);
X        }
X    }
X    (void) fclose (fp);
X}
X
X
X/* draw all the stuff on the screen, using the current menu.
X * if how_much == 0 then just update fields that need it;
X * if how_much == 1 then redraw all fields;
X * if how_much == 2 then erase the screen and redraw EVERYTHING.
X */
Xredraw_screen (how_much)
Xint how_much;
X{
X    if (how_much == 2)
X        c_erase();
X
X    /* print the single-step message if this is the last loop */
X    if (nstep < 1)
X        print_updating();
X
X    if (how_much == 2) {
X        mm_borders();
X        mm_labels();
X        srch_prstate(1);
X        plot_prstate(1);
X        listing_prstate(1);
X    }
X
X    /* if just updating changed fields while plotting or listing
X     * unattended then suppress most screen updates except
X     * always show nstep to show plot loops to go and
X     * always show tminc to show search convergence progress.
X     */
X    print_nstep(how_much);
X    print_tminc(how_much);
X    print_spause(how_much);
X    if (how_much == 0 && (plot_ison() || listing_ison()) && nstep > 0)
X        f_off();
X
X    /* print all the time-related fields */
X    mm_now (&now, how_much);
X
X    if (optwi)
X        mm_twilight (&now, how_much);
X
X    /* print stuff on bottom menu */
X    print_alt (how_much);
X
X    f_on();
X}
X
X/* clean up and exit.
X */
Xvoid
Xbye()
X{
X    c_erase();
X    byetty();
X    exit (0);
X}
X
X/* this gets called when a floating point error occurs.
X * we force a jump back into main() with looping terminated.
X */
Xstatic
Xvoid
Xon_fpe()
X{
X    extern void longjmp();
X
X    (void) signal (SIGFPE, on_fpe);
X    f_msg ("Floating point error has occurred - computations aborted.");
X    longjmp (fpe_err_jmp, 1);
X}
X
Xusage(why)
Xchar *why;
X{
X    /* don't advertise -s (silent) option */
X    c_erase();
X    f_string (1, 1, why);
X    f_string (2, 1,
X        "usage: [-c <configfile>] [-d <database>] [field=value ...]\r\n");
X    byetty();
X    exit (1);
X}
X
X/* process the field specs from the command line.
X * if trouble call usage() (which exits).
X */
Xstatic
Xread_fieldargs (ac, av)
Xint ac;        /* number of such specs */
Xchar *av[];    /* array of strings in form <field_name value> */
X{
X    while (--ac >= 0) {
X        char *fs = *av++;
X        if (crack_fieldset (fs) < 0) {
X        char why[NC];
X        (void) sprintf (why, "Bad command line spec: %.*s",
X                            sizeof(why)-26, fs);
X        usage (why);
X        }
X    }
X}
X
X/* process a field spec in buf, either from config file or argv.
X * return 0 if recognized ok, else -1.
X */
Xstatic
Xcrack_fieldset (buf)
Xchar *buf;
X{
X#define    ARRAY_SIZ(a)    (sizeof(a)/sizeof((a)[0]))
X#define    MAXKW        6    /* longest keyword, not counting trailing 0 */
X    /* N.B. index of item is its case value, below.
X     * N.B. if add an item, keep it no longer than MAXKW chars.
X     */
X    static char keywords[][MAXKW+1] = {
X        /*  0 */    "LAT",
X        /*  1 */    "LONG",
X        /*  2 */    "UT",
X        /*  3 */    "UD",
X        /*  4 */    "TZONE",
X        /*  5 */    "TZNAME",
X        /*  6 */    "HEIGHT",
X        /*  7 */    "NSTEP",
X        /*  8 */    "PAUSE",
X        /*  9 */    "STPSZ",
X        /* 10 */    "TEMP",
X        /* 11 */    "PRES",
X        /* 12 */    "EPOCH",
X        /* 13 */    "JD",
X        /* 14 */    "OBJX",
X        /* 15 */    "OBJY",
X        /* 16 */    "PROPTS",
X        /* 17 */    "MENU"
X    };
X    int i;
X    int l;
X    int f;
X
X    for (i = 0; i < ARRAY_SIZ(keywords); i++)
X        if (strncmp (keywords[i], buf, l = strlen(keywords[i])) == 0) {
X        buf += l+1;    /* skip keyword and its subsequent delimiter */
X        break;
X        }
X
X    switch (i) {
X    case 0: f = rcfpack (R_LAT,C_LATV,0); (void) chg_fld (buf, &f);
X        break;
X    case 1: f = rcfpack (R_LONG,C_LONGV,0), (void) chg_fld (buf, &f);
X        break;
X    case 2: f = rcfpack (R_UT,C_UTV,0), (void) chg_fld (buf, &f);
X        break;
X    case 3: f = rcfpack (R_UD,C_UD,0), (void) chg_fld (buf, &f);
X        break;
X    case 4: f = rcfpack (R_TZONE,C_TZONEV,0), (void) chg_fld (buf, &f);
X        break;
X    case 5: f = rcfpack (R_TZN,C_TZN,0), (void) chg_fld (buf, &f);
X        break;
X    case 6: f = rcfpack (R_HEIGHT,C_HEIGHTV,0), (void) chg_fld (buf, &f);
X        break;
X    case 7: f = rcfpack (R_NSTEP,C_NSTEPV,0), (void) chg_fld (buf, &f);
X        break;
X    case 8: f = rcfpack (R_PAUSE,C_PAUSEV,0), (void) chg_fld (buf, &f);
X        break;
X    case 9: f = rcfpack (R_STPSZ,C_STPSZV,0), (void) chg_fld (buf, &f);
X        break;
X    case 10: f = rcfpack (R_TEMP,C_TEMPV,0), (void) chg_fld (buf, &f);
X        break;
X    case 11: f = rcfpack (R_PRES,C_PRESV,0), (void) chg_fld (buf, &f);
X        break;
X    case 12: f = rcfpack (R_EPOCH,C_EPOCHV,0), (void) chg_fld (buf, &f);
X        break;
X    case 13: f = rcfpack (R_JD,C_JDV,0), (void) chg_fld (buf, &f);
X        break;
X    case 14: (void) obj_filelookup (OBJX, buf);
X        break;
X    case 15: (void) obj_filelookup (OBJY, buf);
X        break;
X    case 16:
X        if (buf[-1] != '+')
X        optwi = oppl = 0;
X        while (*buf)
X        switch (*buf++) {
X        case 'T': optwi = 1; break;
X        case 'S': oppl |= (1<<SUN); break;
X        case 'M': oppl |= (1<<MOON); break;
X        case 'e': oppl |= (1<<MERCURY); break;
X        case 'v': oppl |= (1<<VENUS); break;
X        case 'm': oppl |= (1<<MARS); break;
X        case 'j': case 'J': oppl |= (1<<JUPITER); break;
X        case 's': oppl |= (1<<SATURN); break;
X        case 'u': oppl |= (1<<URANUS); break;
X        case 'n': oppl |= (1<<NEPTUNE); break;
X        case 'p': oppl |= (1<<PLUTO); break;
X        case 'x': oppl |= (1<<OBJX); obj_on(OBJX); break;
X        case 'y': oppl |= (1<<OBJY); obj_on(OBJY); break;
X        }
X        break;
X    case 17:
X        if (strncmp (buf, "DATA", 4) == 0)
X        altmenu_init (F_MNU1);
X        else if (strncmp (buf, "RISET", 5) == 0)
X        altmenu_init (F_MNU2);
X        else if (strncmp (buf, "SEP", 3) == 0)
X        altmenu_init (F_MNU3);
X        else if (strncmp (buf, "JUP", 3) == 0)
X        altmenu_init (F_MNUJ);
X        break;
X    default:
X        return (-1);
X    }
X    return (0);
X}
X
X/* react to the field at *fld according to the optional string input at bp.
X * if bp is != 0 use it, else issue read_line() and use buffer.
X * then sscanf the buffer and update the corresponding (global) variable(s)
X * or do whatever a pick at that field should do.
X * we might also change *fld if we want to change the current cursor location.
X * return 1 if we change a field that invalidates any of the times or
X * to update all related fields.
X */
Xstatic
Xchg_fld (bp, fld)
Xchar *bp;
Xint *fld;
X{
X    char buf[NC];
X    int deghrs = 0, mins = 0, secs = 0;
X    int new = 0;
X
X    /* switch on just the row/col portion */
X    switch (unpackrc(*fld)) {
X    case rcfpack (R_ALTM, C_ALTM, 0):
X        if (altmenu_setup() == 0) {
X        print_updating();
X        alt_erase();
X        print_alt(2);
X        }
X        break;
X    case rcfpack (R_JD, C_JDV, 0):
X        if (!bp) {
X        static char p[] = "Julian Date (or n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else
X        mjd = atof(bp) - 2415020L;
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_UD, C_UD, 0):
X        if (!bp) {
X        static char p[] = "utc date (m/d/y, or year.d, or n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else {
X        if (decimal_year(bp)) {
X            double y = atof (bp);
X            year_mjd (y, &mjd);
X        } else {
X            double day, newmjd0;
X            int month, year;
X            mjd_cal (mjd, &month, &day, &year); /* init with now */
X            f_sscandate (bp, &month, &day, &year);
X            cal_mjd (month, day, year, &newmjd0);
X            /* if don't give a fractional part to days
X             * then retain current hours.
X             */
X            if ((long)day == day)
X            mjd = newmjd0 + mjd_hr(mjd)/24.0;
X            else
X            mjd = newmjd0;
X        }
X        }
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_UT, C_UTV, 0):
X        if (!bp) {
X        static char p[] = "utc time (h:m:s, or n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else {
X        double newutc = (mjd-mjd_day(mjd)) * 24.0;
X        f_dec_sexsign (newutc, °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &newutc);
X        mjd = mjd_day(mjd) + newutc/24.0;
X        }
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_LD, C_LD, 0):
X        if (!bp) {
X        static char p[] = "local date (m/d/y, or year.d, n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else {
X        if (decimal_year(bp)) {
X            double y = atof (bp);
X            year_mjd (y, &mjd);
X            mjd += tz/24.0;
X        } else {
X            double day, newlmjd0;
X            int month, year;
X            mjd_cal (mjd-tz/24.0, &month, &day, &year); /* now */
X            f_sscandate (bp, &month, &day, &year);
X            cal_mjd (month, day, year, &newlmjd0);
X            /* if don't give a fractional part to days
X             * then retain current hours.
X             */
X            if ((long)day == day)
X            mjd = newlmjd0 + mjd_hr(mjd-tz/24.0)/24.0;
X            else
X            mjd = newlmjd0;
X            mjd += tz/24.0;
X        }
X        }
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_LT, C_LT, 0):
X        if (!bp) {
X        static char p[] = "local time (h:m:s, or n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else {
X        double newlt = (mjd-mjd_day(mjd)) * 24.0 - tz;
X        range (&newlt, 24.0);
X        f_dec_sexsign (newlt, °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &newlt);
X        mjd = mjd_day(mjd-tz/24.0) + (newlt + tz)/24.0;
X        }
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_LST, C_LSTV, 0):
X        if (!bp) {
X        static char p[] = "local sidereal time (h:m:s, or n for Now): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'n' || bp[0] == 'N')
X        time_fromsys (&now);
X        else {
X        double lst, utc;
X        now_lst (&now, &lst);
X        f_dec_sexsign (lst, °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &lst);
X        lst -= radhr(lng); /* convert to gst */
X        range (&lst, 24.0);
X        gst_utc (mjd_day(mjd), lst, &utc);
X        mjd = mjd_day(mjd) + utc/24.0;
X        }
X        set_t0 (&now);
X        new = 1;
X        break;
X    case rcfpack (R_TZN, C_TZN, 0):
X        if (!bp) {
X        static char p[] = "timezone abbreviation (3 char max): ";
X        f_prompt (p);
X        if (read_line (buf, 3) <= 0)
X            break;
X        bp = buf;
X        }
X        (void) strncpy (tznm, bp, sizeof(tznm)-1);
X        new = 1;
X        break;
X    case rcfpack (R_TZONE, C_TZONEV, 0):
X        if (!bp) {
X        static char p[] = "hours behind utc: ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        f_dec_sexsign (tz, °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &tz);
X        new = 1;
X        break;
X    case rcfpack (R_LONG, C_LONGV, 0):
X        if (!bp) {
X        static char p[] = "longitude (+ west) (d:m:s): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        f_dec_sexsign (-raddeg(lng), °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &lng);
X        lng = degrad (-lng);         /* want - radians west */
X        new = 1;
X        break;
X    case rcfpack (R_LAT, C_LATV, 0):
X        if (!bp) {
X        static char p[] = "latitude (+ north) (d:m:s): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        f_dec_sexsign (raddeg(lat), °hrs, &mins, &secs);
X        f_sscansex (bp, °hrs, &mins, &secs);
X        sex_dec (deghrs, mins, secs, &lat);
X        lat = degrad (lat);
X        new = 1;
X        break;
X    case rcfpack (R_HEIGHT, C_HEIGHTV, 0):
X        if (!bp) {
X        static char p[] = "height above sea level (ft): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (sscanf (bp, "%lf", &height) == 1) {
X        height /= 2.093e7; /*convert ft to earth radii above sea level*/
X        new = 1;
X        }
X        break;
X    case rcfpack (R_NSTEP, C_NSTEPV, 0):
X        if (!bp) {
X        static char p[] = "number of steps to run: ";
X        f_prompt (p);
X        if (read_line (buf, 8) <= 0)
X            break;
X        bp = buf;
X        }
X        (void) sscanf (bp, "%d", &nstep);
X        print_nstep (0);
X        break;
X    case rcfpack (R_PAUSE, C_PAUSEV, 0):
X        if (!bp) {
X        static char p[] = "seconds to pause between steps: ";
X        f_prompt (p);
X        if (read_line (buf, 8) <= 0)
X            break;
X        bp = buf;
X        }
X        (void) sscanf (bp, "%d", &spause);
X        print_spause (0);
X        break;
X    case rcfpack (R_TEMP, C_TEMPV, 0):
X        if (!bp) {
X        static char p[] = "temperature (deg.F): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (sscanf (bp, "%lf", &temp) == 1) {
X        temp = 5./9.*(temp - 32.0);    /* want degs C */
X        new = 1;
X        }
X        break;
X    case rcfpack (R_PRES, C_PRESV, 0):
X        if (!bp) {
X        static char p[] =
X            "atmos pressure (in. Hg; 0 for no refraction correction): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (sscanf (bp, "%lf", &pressure) == 1) {
X        pressure *= 33.86;        /* want mBar */
X        new = 1;
X        }
X        break;
X    case rcfpack (R_EPOCH, C_EPOCHV, 0):
X        if (!bp) {
X        static char p[] = "epoch (year, or e for Equinox of Date): ";
X        f_prompt (p);
X        if (read_line (buf, PW-strlen(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'e' || bp[0] == 'E')
X        epoch = EOD;
X        else {
X        double e;
X        e = atof(bp);
X        year_mjd (e, &epoch);
X        }
X        new = 1;
X        break;
X    case rcfpack (R_STPSZ, C_STPSZV, 0):
X        if (!bp) {
X        static char p[] =
X            "step size increment (h:m:s, or <x>d for x days, or r for RTC): ";
X        f_prompt (p);
X        if (read_line (buf, PW-sizeof(p)) <= 0)
X            break;
X        bp = buf;
X        }
X        if (bp[0] == 'r' || bp[0] == 'R')
X        tminc = RTC;
X        else {
X        int last = strlen (bp) - 1;
X        if (bp[last] == 'd') {
X            /* ends in d so treat as a number of days */
X            double x;
X            if (sscanf (bp, "%lf", &x) == 1)
X            tminc = x * 24.0;
X        } else {
X            if (tminc == RTC)
X            deghrs = mins = secs = 0;
X            else
X            f_dec_sexsign (tminc, °hrs, &mins, &secs);
X            f_sscansex (bp, °hrs, &mins, &secs);
X            sex_dec (deghrs, mins, secs, &tminc);
X        }
X        }
X        print_tminc(0);
X        set_t0 (&now);
X        break;
X    case rcfpack (R_PLOT, C_PLOT, 0):
X        plot_setup();
X        if (plot_ison())
X        new = 1;
X        break;
X    case rcfpack (R_LISTING, C_LISTING, 0):
X        listing_setup();
X        if (listing_ison())
X        new = 1;
X        break;
X    case rcfpack (R_WATCH, C_WATCH, 0):
X        watch (&now, tminc, oppl);
X        /* set new reference time to what watch left it.
X         * no need to set new since watch just did a redraw.
X         */
X        set_t0 (&now);
X        break;
X    case rcfpack (R_DAWN, C_DAWN, 0):
X    case rcfpack (R_DUSK, C_DUSK, 0):
X    case rcfpack (R_LON, C_LON, 0):
X        if (optwi ^= 1) {
X        print_updating();
X        mm_twilight (&now, 1);
X        } else {
X        f_blanks (R_DAWN, C_DAWNV, 5);
X        f_blanks (R_DUSK, C_DUSKV, 5);
X        f_blanks (R_LON, C_LONV, 5);
X        }
X        break;
X    case rcfpack (R_SRCH, C_SRCH, 0):
X        srch_setup();
X        if (srch_ison())
X        new = 1;
X        break;
X    case rcfpack (R_SUN, C_OBJ, 0):
X        toggle_body (SUN);
X        break;
X    case rcfpack (R_SUN, C_CONSTEL, 0):
X        if (oppl & (1<<SUN))
X        constellation_msg (SUN, &now);
X        break;
X    case rcfpack (R_MOON, C_OBJ, 0):
X        toggle_body (MOON);
X        break;
X    case rcfpack (R_MOON, C_CONSTEL, 0):
X        if (oppl & (1<<MOON))
X        constellation_msg (MOON, &now);
X        break;
X    case rcfpack (R_MERCURY, C_OBJ, 0):
X        toggle_body (MERCURY);
X        break;
X    case rcfpack (R_MERCURY, C_CONSTEL, 0):
X        if (oppl & (1<<MERCURY))
X        constellation_msg (MERCURY, &now);
X        break;
X    case rcfpack (R_VENUS, C_OBJ, 0):
X        toggle_body (VENUS);
X        break;
X    case rcfpack (R_VENUS, C_CONSTEL, 0):
X        if (oppl & (1<<VENUS))
X        constellation_msg (VENUS, &now);
X        break;
X    case rcfpack (R_MARS, C_OBJ, 0):
X        toggle_body (MARS);
X        break;
X    case rcfpack (R_MARS, C_CONSTEL, 0):
X        if (oppl & (1<<MARS))
X        constellation_msg (MARS, &now);
X        break;
X    case rcfpack (R_JUPITER, C_OBJ, 0):
X        toggle_body (JUPITER);
X        break;
X    case rcfpack (R_JUPITER, C_CONSTEL, 0):
X        if (oppl & (1<<JUPITER))
X        constellation_msg (JUPITER, &now);
X        break;
X    case rcfpack (R_JUPITER, C_XTRA, 0):
X        if (oppl & (1<<JUPITER)) {
X        print_updating();
X        alt_erase();
X        altmenu_init (F_MNUJ);
X        print_alt (2);
X        *fld = rcfpack(R_NSTEP, C_NSTEPV, 0);
X        }
X        break;
X    case rcfpack (R_SATURN, C_OBJ, 0):
X        toggle_body (SATURN);
X        break;
X    case rcfpack (R_SATURN, C_CONSTEL, 0):
X        if (oppl & (1<<SATURN))
X        constellation_msg (SATURN, &now);
X        break;
X    case rcfpack (R_URANUS, C_OBJ, 0):
X        toggle_body (URANUS);
X        break;
X    case rcfpack (R_URANUS, C_CONSTEL, 0):
X        if (oppl & (1<<URANUS))
X        constellation_msg (URANUS, &now);
X        break;
X    case rcfpack (R_NEPTUNE, C_OBJ, 0):
X        toggle_body (NEPTUNE);
X        break;
X    case rcfpack (R_NEPTUNE, C_CONSTEL, 0):
X        if (oppl & (1<<NEPTUNE))
X        constellation_msg (NEPTUNE, &now);
X        break;
X    case rcfpack (R_PLUTO, C_OBJ, 0):
X        toggle_body (PLUTO);
X        break;
X    case rcfpack (R_PLUTO, C_CONSTEL, 0):
X        if (oppl & (1<<PLUTO))
X        constellation_msg (PLUTO, &now);
X        break;
X    case rcfpack (R_OBJX, C_OBJ, 0):
X        /* this might change which columns are used so erase all when
X         * returns and redraw if still on.
X         */
X        obj_setup (OBJX);
X        alt_nobody (OBJX);
X        if (obj_ison (OBJX)) {
X        oppl |= 1 << OBJX;
X        print_updating();
X        alt_body (OBJX, 1, &now);
X        } else
X        oppl &= ~(1 << OBJX);    /* already erased; just clear flag */
X        break;
X    case rcfpack (R_OBJX, C_CONSTEL, 0):
X        if (oppl & (1<<OBJX))
X        constellation_msg (OBJX, &now);
X        break;
X    case rcfpack (R_OBJY, C_OBJ, 0):
X        /* this might change which columns are used so erase all when
X         * returns and redraw if still on.
X         */
X        obj_setup (OBJY);
X        alt_nobody (OBJY);
X        if (obj_ison (OBJY)) {
X        oppl |= 1 << OBJY;
X        print_updating();
X        alt_body (OBJY, 1, &now);
X        } else
X        oppl &= ~(1 << OBJY);    /* already erased; just clear flag */
X        break;
X    case rcfpack (R_OBJY, C_CONSTEL, 0):
X        if (oppl & (1<<OBJY))
X        constellation_msg (OBJY, &now);
X        break;
X    }
X
X    return (new);
X}
X
Xstatic
Xprint_tminc(force)
Xint force;
X{
X    static double last = -123.456;    /* anything unlikely */
X
X    if (force || tminc != last) {
X        if (tminc == RTC)
X        f_string (R_STPSZ, C_STPSZV, " RT CLOCK");
X        else if (fabs(tminc) >= 24.0)
X        f_double (R_STPSZ, C_STPSZV, "%6.4g dy", tminc/24.0);
X        else
X        f_signtime (R_STPSZ, C_STPSZV, tminc);
X        last = tminc;
X    }
X}
X
X/* print stuff on bottom menu */
Xstatic
Xprint_alt (howmuch)
Xint howmuch;
X{
X    if (howmuch == 2)
X        alt_labels();
X    if (alt_menumask() == F_MNUJ)
X        altj_display (howmuch, &now);
X    else {
X        int p;
X        for (p = nxtbody(-1); p != -1; p = nxtbody(p))
X        if (oppl & (1<<p))
X            alt_body (p, howmuch, &now);
X    }
X}
X
Xprint_updating()
X{
X    f_prompt ("Updating...");
X}
X
Xstatic
Xprint_nstep(force)
Xint force;
X{
X    static int last;
X
X    if (force || nstep != last) {
X        char buf[16];
X        (void) sprintf (buf, "%8d", nstep);
X        f_string (R_NSTEP, C_NSTEPV, buf);
X        last = nstep;
X    }
X}
X
Xstatic
Xprint_spause(force)
Xint force;
X{
X    static int last;
X
X    if (force || spause != last) {
X        char buf[16];
X        (void) sprintf (buf, "%8d", spause);
X        f_string (R_PAUSE, C_PAUSEV, buf);
X        last = spause;
X    }
X}
X
X/* if not plotting/listing/searching then sleep spause seconds.
X * if time is being based on the real-time clock, sync on the next
X *   integral multiple of spause seconds after the minute.
X * check for keyboard action once each second to let it break out early.
X */
Xslp_sync()
X{
X    extern long time();
X
X    if (spause > 0 && !plot_ison() && !srch_ison() && !listing_ison()) {
X        int n;
X        if (tminc == RTC) {
X        long t;
X        (void) time (&t);
X        n = spause - (t % spause);
X        } else
X        n = spause;
X        while (--n >= 0)
X        if (chk_char() == 0)
X            break;
X        else
X            (void) sleep (1);
X    }
X}
X
Xstatic
Xtoggle_body (p)
Xint p;
X{
X    if ((oppl ^= (1<<p)) & (1<<p)) {
X        print_updating();
X        alt_body (p, 1, &now);
X    } else
X        alt_nobody (p);
X}
END_OF_FILE
  if test 26117 -ne `wc -c <'main.c'`; then
    echo shar: \"'main.c'\" unpacked with wrong size!
  fi
  # end of 'main.c'
fi
if test -f 'mainmenu.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'mainmenu.c'\"
else
  echo shar: Extracting \"'mainmenu.c'\" \(6349 characters\)
  sed "s/^X//" >'mainmenu.c' <<'END_OF_FILE'
X/* printing routines for the main (upper) screen.
X */
X
X#include <stdio.h>
X#include <math.h>
X#include "astro.h"
X#include "circum.h"
X#include "screen.h"
X
X/* #define PC_GRAPHICS */
X#ifdef PC_GRAPHICS
X#define    JOINT    207
X#define    VERT    179
X#define    HORIZ    205
X#else
X#define    JOINT    '-'
X#define    VERT    '|'
X#define    HORIZ    '-'
X#endif
X
Xmm_borders()
X{
X    char line[NC+1], *lp;
X    register i;
X
X    lp = line;
X    for (i = 0; i < NC; i++)
X        *lp++ = HORIZ;
X    *lp = '\0';
X    f_string (R_PLANTAB-1, 1, line);
X    for (i = R_TOP; i < R_PLANTAB-1; i++)
X        f_char (i, COL2-2, VERT);
X    f_char (R_PLANTAB-1, COL2-2, JOINT);
X    for (i = R_TOP; i < R_PLANTAB-1; i++)
X        f_char (i, COL3-2, VERT);
X    f_char (R_PLANTAB-1, COL3-2, JOINT);
X    for (i = R_LST; i < R_PLANTAB-1; i++)
X        f_char (i, COL4-2, VERT);
X    f_char (R_PLANTAB-1, COL4-2, JOINT);
X}
X
X/* print the permanent labels on the top menu */
Xmm_labels()
X{
X    f_string (R_TZN,    C_TZN,        "LT");
X    f_string (R_UT,        C_UT,        "UTC");
X    f_string (R_JD,        C_JD,        "JulianDate");
X    f_string (R_LISTING,    C_LISTING,    "Listing");
X    f_string (R_WATCH,    C_WATCH,    "Watch");
X    f_string (R_SRCH,    C_SRCH,        "Search");
X    f_string (R_PLOT,    C_PLOT,        "Plot");
X    f_string (R_ALTM,    C_ALTM,        "Menu");
X
X    f_string (R_LST,    C_LST,        "LST");
X    f_string (R_DAWN,    C_DAWN,        "Dawn");
X    f_string (R_DUSK,    C_DUSK,        "Dusk");
X    f_string (R_LON,    C_LON,        "NiteLn");
X    f_string (R_PAUSE,    C_PAUSE,    "Pause");
X    f_string (R_NSTEP,    C_NSTEP,    "NStep");
X    f_string (R_STPSZ,    C_STPSZ,    "StpSz");
X
X    f_string (R_LAT,    C_LAT,        "Lat");
X    f_string (R_LONG,    C_LONG,        "Long");
X    f_string (R_HEIGHT,    C_HEIGHT,    "Elev");
X    f_string (R_TEMP,    C_TEMP,        "Temp");
X    f_string (R_PRES,    C_PRES,        "AtmPr");
X    f_string (R_TZONE,    C_TZONE,    "TZ");
X    f_string (R_EPOCH,    C_EPOCH,    "Epoch");
X}
X
X/* print all the time/date/where related stuff: the Now structure.
X * print in a nice order, based on the field locations, as much as possible.
X */
Xmm_now (np, all)
XNow *np;
Xint all;
X{
X    char buf[32];
X    double lmjd = mjd - tz/24.0;
X    double jd = mjd + 2415020L;
X    double tmp;
X
X    (void) sprintf (buf, "%-3.3s", tznm);
X    f_string (R_TZN, C_TZN, buf);
X    f_time (R_LT, C_LT, mjd_hr(lmjd));
X    f_date (R_LD, C_LD, lmjd);
X
X    f_time (R_UT, C_UTV, mjd_hr(mjd));
X    f_date (R_UD, C_UD, mjd);
X
X    (void) sprintf (buf, "%14.5f", jd);
X    (void) flog_log (R_JD, C_JDV, jd, buf);
X    f_string (R_JD, C_JDV, buf);
X
X    now_lst (np, &tmp);
X    f_time (R_LST, C_LSTV, tmp);
X
X    if (all) {
X        f_gangle (R_LAT, C_LATV, lat);
X        f_gangle (R_LONG, C_LONGV, -lng);    /* + west */
X
X        tmp = height * 2.093e7;    /* want to see ft, not earth radii */
X        (void) sprintf (buf, "%5g ft", tmp);
X        (void) flog_log (R_HEIGHT, C_HEIGHTV, tmp, buf);
X        f_string (R_HEIGHT, C_HEIGHTV, buf);
X
X        tmp = 9./5.*temp + 32.0;     /* want to see degrees F, not C */
X        (void) sprintf (buf, "%6g F", tmp);
X        (void) flog_log (R_TEMP, C_TEMPV, tmp, buf);
X        f_string (R_TEMP, C_TEMPV, buf);
X
X        tmp = pressure / 33.86;    /* want to see in. Hg, not mBar */
X        (void) sprintf (buf, "%5.2f in", tmp);
X        (void) flog_log (R_PRES, C_PRESV, tmp, buf);
X        f_string (R_PRES, C_PRESV, buf);
X
X        f_signtime (R_TZONE, C_TZONEV, tz);
X
X        if (epoch == EOD)
X        f_string (R_EPOCH, C_EPOCHV, "(OfDate)");
X        else {
X        mjd_year (epoch, &tmp);
X        f_double (R_EPOCH, C_EPOCHV, "%8.1f", tmp);
X        }
X    }
X
X    /* print the calendar for local day, if new month/year.  */
X    mm_calendar (np, all > 1);
X}
X
X/* display dawn/dusk/length-of-night times.
X */
Xmm_twilight (np, force)
XNow *np;
Xint force;
X{
X    double dusk, dawn;
X    double tmp;
X    int status;
X
X    if (!twilight_cir (np, &dawn, &dusk, &status) && !force)
X        return;
X
X    if (status != 0) {
X        f_blanks (R_DAWN, C_DAWNV, 5);
X        f_blanks (R_DUSK, C_DUSKV, 5);
X        f_string (R_LON, C_LONV, "-----");
X        return;
X    }
X
X    f_mtime (R_DAWN, C_DAWNV, dawn);
X    f_mtime (R_DUSK, C_DUSKV, dusk);
X    tmp = dawn - dusk; range (&tmp, 24.0);
X    f_mtime (R_LON, C_LONV, tmp);
X}
X
Xmm_newcir (y)
Xint y;
X{
X    static char ncmsg[] = "NEW CIRCUMSTANCES";
X    static char nomsg[] = "                 ";
X    static int last_y = -1;
X
X    if (y != last_y) {
X        f_string (R_NEWCIR, C_NEWCIR, y ? ncmsg : nomsg);
X        last_y = y;
X    }
X}
X
Xstatic
Xmm_calendar (np, force)
XNow *np;
Xint force;
X{
X    static char *mnames[] = {
X        "January", "February", "March", "April", "May", "June",
X        "July", "August", "September", "October", "November", "December"
X    };
X    static int last_m, last_y;
X    static double last_tz = -100;
X    char str[64];
X    int m, y;
X    double d;
X    int f, nd;
X    int r;
X    double jd0;
X
X    /* get local m/d/y. do nothing if still same month and not forced. */
X    mjd_cal (mjd_day(mjd-tz/24.0), &m, &d, &y);
X    if (m == last_m && y == last_y && tz == last_tz && !force)
X        return;
X    last_m = m;
X    last_y = y;
X    last_tz = tz;
X
X    /* find day of week of first day of month */
X    cal_mjd (m, 1.0, y, &jd0);
X    mjd_dow (jd0, &f);
X    if (f < 0) {
X        /* can't figure it out - too hard before Gregorian */
X        int i;
X        for (i = 8; --i >= 0; )
X        f_string (R_CAL+i, C_CAL, "                    ");
X        return;
X    }
X
X    /* print header */
X    f_blanks (R_CAL, C_CAL, 20);
X    (void) sprintf (str, "%s %4d", mnames[m-1], y);
X    f_string (R_CAL, C_CAL + (20 - (strlen(mnames[m-1]) + 5))/2, str);
X    f_string (R_CAL+1, C_CAL, "Su Mo Tu We Th Fr Sa");
X
X    /* find number of days in this month */
X    mjd_dpm (jd0, &nd);
X
X    /* print the calendar */
X    for (r = 0; r < 6; r++) {
X        char row[7*3+1], *rp = row;
X        int c;
X        for (c = 0; c < 7; c++) {
X        int i = r*7+c;
X        if (i < f || i >= f + nd)
X            (void) sprintf (rp, "   ");
X        else
X            (void) sprintf (rp, "%2d ", i-f+1);
X        rp += 3;
X        }
X        row[sizeof(row)-2] = '\0';    /* don't print last blank; causes wrap*/
X        f_string (R_CAL+2+r, C_CAL, row);
X    }
X
X    /* over print the new and full moons for this month.
X     * TODO: don't really know which dates to use here (see moonnf())
X     *   so try several to be fairly safe. have to go back to 4/29/1988
X     *   to find the full moon on 5/1 for example.
X     */
X    mm_nfmoon (jd0-3, tz, m, f);
X    mm_nfmoon (jd0+15, tz, m, f);
X}
X
Xstatic
Xmm_nfmoon (jd, tzone, m, f)
Xdouble jd, tzone;
Xint m, f;
X{
X    static char nm[] = "NM", fm[] = "FM";
X    double dm;
X    int mm, ym;
X    double jdn, jdf;
X    int di;
X
X    moonnf (jd, &jdn, &jdf);
X    mjd_cal (jdn-tzone/24.0, &mm, &dm, &ym);
X    if (m == mm) {
X        di = dm + f - 1;
X        f_string (R_CAL+2+di/7, C_CAL+3*(di%7), nm);
X    }
X    mjd_cal (jdf-tzone/24.0, &mm, &dm, &ym);
X    if (m == mm) {
X        di = dm + f - 1;
X        f_string (R_CAL+2+di/7, C_CAL+3*(di%7), fm);
X    }
X}
END_OF_FILE
  if test 6349 -ne `wc -c <'mainmenu.c'`; then
    echo shar: \"'mainmenu.c'\" unpacked with wrong size!
  fi
  # end of 'mainmenu.c'
fi
if test -f 'plans.c' -a "${1}" != "-c" ; then
  echo shar: Will not clobber existing file \"'plans.c'\"
else
  echo shar: Extracting \"'plans.c'\" \(17647 characters\)
  sed "s/^X//" >'plans.c' <<'END_OF_FILE'
X#include <stdio.h>
X#include <math.h>
X#include "astro.h"
X
X#define    TWOPI        (2*PI)
X#define    mod2PI(x)    ((x) - (long)((x)/TWOPI)*TWOPI)
X
X/* given a modified Julian date, mjd, and a planet, p, find:
X *   lpd0: heliocentric longitude,
X *   psi0: heliocentric latitude,
X *   rp0:  distance from the sun to the planet,
X *   rho0: distance from the Earth to the planet,
X *         none corrected for light time, ie, they are the true values for the
X *         given instant.
X *   lam:  geocentric ecliptic longitude,
X *   bet:  geocentric ecliptic latitude,
X *         each corrected for light time, ie, they are the apparent values as
X *       seen from the center of the Earth for the given instant.
X *   dia:  angular diameter in arcsec at 1 AU,
X *   mag:  visual magnitude when 1 AU from sun and earth at 0 phase angle.
X *
X * all angles are in radians, all distances in AU.
X * the mean orbital elements are found by calling pelement(), then mutual
X *   perturbation corrections are applied as necessary.
X *
X * corrections for nutation and abberation must be made by the caller. The RA
X *   and DEC calculated from the fully-corrected ecliptic coordinates are then
X *   the apparent geocentric coordinates. Further corrections can be made, if
X *   required, for atmospheric refraction and geocentric parallax although the
X *   intrinsic error herein of about 10 arcseconds is usually the dominant
X *   error at this stage.
X * TODO: combine the several intermediate expressions when get a good compiler.
X */
Xplans (mjd, p, lpd0, psi0, rp0, rho0, lam, bet, dia, mag)
Xdouble mjd;
Xint p;
Xdouble *lpd0, *psi0, *rp0, *rho0, *lam, *bet, *dia, *mag;
X{
X    static double plan[8][9];
X    static double lastmjd = -10000;
X    double dl;    /* perturbation correction for longitude */
X    double dr;    /*  "   orbital radius */
X    double dml;    /*  "   mean longitude */
X    double ds;    /*  "   eccentricity */
X    double dm;    /*  "   mean anomaly */
X    double da;    /*  "   semi-major axis */
X    double dhl;    /*  "   heliocentric longitude */
X    double lsn, rsn;/* true geocentric longitude of sun and sun-earth rad */
X    double mas;    /* mean anomaly of the sun */
X    double re;    /* radius of earth's orbit */
X    double lg;    /* longitude of earth */
X    double map[8];    /* array of mean anomalies for each planet */
X    double lpd, psi, rp, rho;
X    double ll, sll, cll;
X    double t;
X    double dt;
X    int pass;
X    int j;
X    double s, ma;
X    double nu, ea;
X    double lp, om;
X    double lo, slo, clo;
X    double inc, y;
X    double spsi, cpsi;
X    double rpd;
X
X    /* only need to fill in plan[] once for a given mjd */
X    if (mjd != lastmjd) {
X        pelement (mjd, plan);
X        lastmjd = mjd;
X    }
X
X    dt = 0;
X    t = mjd/36525.;
X    sunpos (mjd, &lsn, &rsn);
X    masun (mjd, &mas);
X        re = rsn;
X    lg = lsn+PI;
X
X    /* first find the true position of the planet at mjd.
X     * then repeat a second time for a slightly different time based
X     * on the position found in the first pass to account for light-travel
X     * time.
X     */
X    for (pass = 0; pass < 2; pass++) {
X
X        for (j = 0; j < 8; j++)
X        map[j] = degrad(plan[j][0]-plan[j][2]-dt*plan[j][1]);
X
X        /* set initial corrections to 0.
X         * then modify as necessary for the planet of interest.
X         */
X        dl = 0;
X        dr = 0;
X        dml = 0;
X        ds = 0;
X        dm = 0;
X        da = 0;
X        dhl = 0;
X
X        switch (p) {
X
X        case MERCURY:
X        p_mercury (map, &dl, &dr);
X        break;
X
X        case VENUS:
X        p_venus (t, mas, map, &dl, &dr, &dml, &dm);
X        break;
X
X        case MARS:
X        p_mars (mas, map, &dl, &dr, &dml, &dm);
X        break;
X
X        case JUPITER:
X        p_jupiter (t, plan[p][3], &dml, &ds, &dm, &da);
X        break;
X
X        case SATURN:
X        p_saturn (t, plan[p][3], &dml, &ds, &dm, &da, &dhl);
X        break;
X
X        case URANUS:
X        p_uranus (t, plan[p][3], &dl, &dr, &dml, &ds, &dm, &da, &dhl);
X        break;
X
X        case NEPTUNE:
X        p_neptune (t, plan[p][3], &dl, &dr, &dml, &ds, &dm, &da, &dhl);
X        break;
X
X        case PLUTO:
X        /* no perturbation theory for pluto */
X        break;
X        }
X
X        s = plan[p][3]+ds;
X        ma = map[p]+dm;
X        anomaly (ma, s, &nu, &ea);
X        rp = (plan[p][6]+da)*(1-s*s)/(1+s*cos(nu));
X        lp = raddeg(nu)+plan[p][2]+raddeg(dml-dm);
X        lp = degrad(lp);
X        om = degrad(plan[p][5]);
X        lo = lp-om;
X        slo = sin(lo);
X        clo = cos(lo);
X        inc = degrad(plan[p][4]);
X        rp = rp+dr;
X        spsi = slo*sin(inc);
X        y = slo*cos(inc);
X        psi = asin(spsi)+dhl;
X        spsi = sin(psi);
X        lpd = atan(y/clo)+om+degrad(dl);
X        if (clo<0) lpd += PI;
X        range (&lpd, TWOPI);
X        cpsi = cos(psi);
X        rpd = rp*cpsi;
X        ll = lpd-lg;
X        rho = sqrt(re*re+rp*rp-2*re*rp*cpsi*cos(ll));
X
X        /* when we view a planet we see it in the position it occupied
X         * dt days ago, where rho is the distance between it and earth,
X         * in AU. use this as the new time for the next pass.
X         */
X        dt = rho*5.775518e-3;
X
X        if (pass == 0) {
X        /* save heliocentric coordinates after first pass since, being
X         * true, they are NOT to be corrected for light-travel time.
X         */
X        *lpd0 = lpd;
X        range (lpd0, TWOPI);
X        *psi0 = psi;
X        *rp0 = rp;
X        *rho0 = rho;
X        }
X    }
X
X        sll = sin(ll);
X    cll = cos(ll);
X        if (p < MARS)
X        *lam = atan(-1*rpd*sll/(re-rpd*cll))+lg+PI;
X    else
X        *lam = atan(re*sll/(rpd-re*cll))+lpd;
X    range (lam, TWOPI);
X        *bet = atan(rpd*spsi*sin(*lam-lpd)/(cpsi*re*sll));
X    *dia = plan[p][7];
X    *mag = plan[p][8];
X}
X
X/* set auxilliary variables used for jupiter, saturn, uranus, and neptune */
Xstatic
Xaux_jsun (t, x1, x2, x3, x4, x5, x6)
Xdouble t;
Xdouble *x1, *x2, *x3, *x4, *x5, *x6;
X{
X        *x1 = t/5+0.1;
X        *x2 = mod2PI(4.14473+5.29691e1*t);
X        *x3 = mod2PI(4.641118+2.132991e1*t);
X        *x4 = mod2PI(4.250177+7.478172*t);
X        *x5 = 5 * *x3 - 2 * *x2;
X    *x6 = 2 * *x2 - 6 * *x3 + 3 * *x4;
X}
X
X/* find the mean anomaly of the sun at mjd.
X * this is the same as that used in sun() but when it was converted to C it
X * was not known it would be required outside that routine.
X * TODO: add an argument to sun() to return mas and eliminate this routine.
X */
Xstatic
Xmasun (mjd, mas)
Xdouble mjd;
Xdouble *mas;
X{
X    double t, t2;
X    double a, b;
X
X    t = mjd/36525;
X    t2 = t*t;
X    a = 9.999736042e1*t;
X    b = 360.*(a-(long)a);
X    *mas = degrad (3.5847583e2-(1.5e-4+3.3e-6*t)*t2+b);
X}
X
X/* perturbations for mercury */
Xstatic
Xp_mercury (map, dl, dr)
Xdouble map[];
Xdouble *dl, *dr;
X{
X    *dl = 2.04e-3*cos(5*map[2-1]-2*map[1-1]+2.1328e-1)+
X         1.03e-3*cos(2*map[2-1]-map[1-1]-2.8046)+
X         9.1e-4*cos(2*map[3]-map[1-1]-6.4582e-1)+
X         7.8e-4*cos(5*map[2-1]-3*map[1-1]+1.7692e-1);
X
X    *dr = 7.525e-6*cos(2*map[3]-map[1-1]+9.25251e-1)+
X         6.802e-6*cos(5*map[2-1]-3*map[1-1]-4.53642)+
X         5.457e-6*cos(2*map[2-1]-2*map[1-1]-1.24246)+
X         3.569e-6*cos(5*map[2-1]-map[1-1]-1.35699);
X}
X
X/* ....venus */
Xstatic
Xp_venus (t, mas, map, dl, dr, dml, dm)
Xdouble t, mas, map[];
Xdouble *dl, *dr, *dml, *dm;
X{
X    *dml = degrad (7.7e-4*sin(4.1406+t*2.6227));
X    *dm = *dml;
X
X    *dl = 3.13e-3*cos(2*mas-2*map[2-1]-2.587)+
X         1.98e-3*cos(3*mas-3*map[2-1]+4.4768e-2)+
X         1.36e-3*cos(mas-map[2-1]-2.0788)+
X         9.6e-4*cos(3*mas-2*map[2-1]-2.3721)+
X         8.2e-4*cos(map[3]-map[2-1]-3.6318);
X
X    *dr = 2.2501e-5*cos(2*mas-2*map[2-1]-1.01592)+
X         1.9045e-5*cos(3*mas-3*map[2-1]+1.61577)+
X         6.887e-6*cos(map[3]-map[2-1]-2.06106)+
X         5.172e-6*cos(mas-map[2-1]-5.08065e-1)+
X         3.62e-6*cos(5*mas-4*map[2-1]-1.81877)+
X         3.283e-6*cos(4*mas-4*map[2-1]+1.10851)+
X         3.074e-6*cos(2*map[3]-2*map[2-1]-9.62846e-1);
X}
X
X/* ....mars */
Xstatic
Xp_mars (mas, map, dl, dr, dml, dm)
Xdouble mas, map[];
Xdouble *dl, *dr, *dml, *dm;
X{
X    double a;
X
X    a = 3*map[3]-8*map[2]+4*mas;
X    *dml = degrad (-1*(1.133e-2*sin(a)+9.33e-3*cos(a)));
X    *dm = *dml;
X
X    *dl = 7.05e-3*cos(map[3]-map[2]-8.5448e-1)+
X         6.07e-3*cos(2*map[3]-map[2]-3.2873)+
X         4.45e-3*cos(2*map[3]-2*map[2]-3.3492)+
X         3.88e-3*cos(mas-2*map[2]+3.5771e-1)+
X         2.38e-3*cos(mas-map[2]+6.1256e-1)+
X         2.04e-3*cos(2*mas-3*map[2]+2.7688)+
X         1.77e-3*cos(3*map[2]-map[2-1]-1.0053)+
X         1.36e-3*cos(2*mas-4*map[2]+2.6894)+
X         1.04e-3*cos(map[3]+3.0749e-1);
X
X    *dr = 5.3227e-5*cos(map[3]-map[2]+7.17864e-1)+
X         5.0989e-5*cos(2*map[3]-2*map[2]-1.77997)+
X         3.8278e-5*cos(2*map[3]-map[2]-1.71617)+
X         1.5996e-5*cos(mas-map[2]-9.69618e-1)+
X         1.4764e-5*cos(2*mas-3*map[2]+1.19768)+
X         8.966e-6*cos(map[3]-2*map[2]+7.61225e-1);
X     *dr += 7.914e-6*cos(3*map[3]-2*map[2]-2.43887)+
X         7.004e-6*cos(2*map[3]-3*map[2]-1.79573)+
X         6.62e-6*cos(mas-2*map[2]+1.97575)+
X         4.93e-6*cos(3*map[3]-3*map[2]-1.33069)+
X         4.693e-6*cos(3*mas-5*map[2]+3.32665)+
X         4.571e-6*cos(2*mas-4*map[2]+4.27086)+
X         4.409e-6*cos(3*map[3]-map[2]-2.02158);
X}
X
X/* ....jupiter */
Xstatic
Xp_jupiter (t, s, dml, ds, dm, da)
Xdouble t, s;
Xdouble *dml, *ds, *dm, *da;
X{
X    double dp;
X    double x1, x2, x3, x4, x5, x6, x7;
X    double sx3, cx3, s2x3, c2x3;
X        double sx5, cx5, s2x5;
X    double sx6;
X        double sx7, cx7, s2x7, c2x7, s3x7, c3x7, s4x7, c4x7, c5x7;
X
X    aux_jsun (t, &x1, &x2, &x3, &x4, &x5, &x6);
X        x7 = x3-x2;
X    sx3 = sin(x3);
X    cx3 = cos(x3);
X        s2x3 = sin(2*x3);
X    c2x3 = cos(2*x3);
X        sx5 = sin(x5);
X    cx5 = cos(x5);
X        s2x5 = sin(2*x5);
X    sx6 = sin(x6);
X        sx7 = sin(x7);
X    cx7 = cos(x7);
X        s2x7 = sin(2*x7);
X    c2x7 = cos(2*x7);
X        s3x7 = sin(3*x7);
X    c3x7 = cos(3*x7);
X        s4x7 = sin(4*x7);
X    c4x7 = cos(4*x7);
X        c5x7 = cos(5*x7);
X
X    *dml = (3.31364e-1-(1.0281e-2+4.692e-3*x1)*x1)*sx5+
X          (3.228e-3-(6.4436e-2-2.075e-3*x1)*x1)*cx5-
X          (3.083e-3+(2.75e-4-4.89e-4*x1)*x1)*s2x5+
X          2.472e-3*sx6+1.3619e-2*sx7+1.8472e-2*s2x7+6.717e-3*s3x7+
X          2.775e-3*s4x7+6.417e-3*s2x7*sx3+
X          (7.275e-3-1.253e-3*x1)*sx7*sx3+
X          2.439e-3*s3x7*sx3-(3.5681e-2+1.208e-3*x1)*sx7*cx3;
X        *dml += -3.767e-3*c2x7*sx3-(3.3839e-2+1.125e-3*x1)*cx7*sx3-
X          4.261e-3*s2x7*cx3+
X          (1.161e-3*x1-6.333e-3)*cx7*cx3+
X          2.178e-3*cx3-6.675e-3*c2x7*cx3-2.664e-3*c3x7*cx3-
X          2.572e-3*sx7*s2x3-3.567e-3*s2x7*s2x3+2.094e-3*cx7*c2x3+
X          3.342e-3*c2x7*c2x3;
X    *dml = degrad(*dml);
X
X    *ds = (3606+(130-43*x1)*x1)*sx5+(1289-580*x1)*cx5-6764*sx7*sx3-
X         1110*s2x7*sx3-224*s3x7*sx3-204*sx3+(1284+116*x1)*cx7*sx3+
X         188*c2x7*sx3+(1460+130*x1)*sx7*cx3+224*s2x7*cx3-817*cx3+
X         6074*cx3*cx7+992*c2x7*cx3+
X         508*c3x7*cx3+230*c4x7*cx3+108*c5x7*cx3;
X    *ds += -(956+73*x1)*sx7*s2x3+448*s2x7*s2x3+137*s3x7*s2x3+
X         (108*x1-997)*cx7*s2x3+480*c2x7*s2x3+148*c3x7*s2x3+
X         (99*x1-956)*sx7*c2x3+490*s2x7*c2x3+
X         158*s3x7*c2x3+179*c2x3+(1024+75*x1)*cx7*c2x3-
X         437*c2x7*c2x3-132*c3x7*c2x3;
X    *ds *= 1e-7;
X
X    dp = (7.192e-3-3.147e-3*x1)*sx5-4.344e-3*sx3+
X         (x1*(1.97e-4*x1-6.75e-4)-2.0428e-2)*cx5+
X         3.4036e-2*cx7*sx3+(7.269e-3+6.72e-4*x1)*sx7*sx3+
X         5.614e-3*c2x7*sx3+2.964e-3*c3x7*sx3+3.7761e-2*sx7*cx3+
X         6.158e-3*s2x7*cx3-
X         6.603e-3*cx7*cx3-5.356e-3*sx7*s2x3+2.722e-3*s2x7*s2x3+
X         4.483e-3*cx7*s2x3-2.642e-3*c2x7*s2x3+4.403e-3*sx7*c2x3-
X         2.536e-3*s2x7*c2x3+5.547e-3*cx7*c2x3-2.689e-3*c2x7*c2x3;
X
X    *dm = *dml-(degrad(dp)/s);
X
X    *da = 205*cx7-263*cx5+693*c2x7+312*c3x7+147*c4x7+299*sx7*sx3+
X         181*c2x7*sx3+204*s2x7*cx3+111*s3x7*cx3-337*cx7*cx3-
X         111*c2x7*cx3;
X    *da *= 1e-6;
X}
X
X/* ....saturn */
Xstatic
Xp_saturn (t, s, dml, ds, dm, da, dhl)
Xdouble t, s;
Xdouble *dml, *ds, *dm, *da, *dhl;
X{
X    double dp;
X    double x1, x2, x3, x4, x5, x6, x7, x8;
X    double sx3, cx3, s2x3, c2x3, s3x3, c3x3, s4x3, c4x3;
X        double sx5, cx5, s2x5, c2x5;
X    double sx6;
X        double sx7, cx7, s2x7, c2x7, s3x7, c3x7, s4x7, c4x7, c5x7, s5x7;
X    double s2x8, c2x8, s3x8, c3x8;
X
X    aux_jsun (t, &x1, &x2, &x3, &x4, &x5, &x6);
X        x7 = x3-x2;
X    sx3 = sin(x3);
X    cx3 = cos(x3);
X        s2x3 = sin(2*x3);
X    c2x3 = cos(2*x3);
X        sx5 = sin(x5);
X    cx5 = cos(x5);
X        s2x5 = sin(2*x5);
X    sx6 = sin(x6);
X        sx7 = sin(x7);
X    cx7 = cos(x7);
X        s2x7 = sin(2*x7);
X    c2x7 = cos(2*x7);
X        s3x7 = sin(3*x7);
X    c3x7 = cos(3*x7);
X        s4x7 = sin(4*x7);
X    c4x7 = cos(4*x7);
X        c5x7 = cos(5*x7);
X
X    s3x3 = sin(3*x3);
X    c3x3 = cos(3*x3);
X    s4x3 = sin(4*x3);
X    c4x3 = cos(4*x3);
X    c2x5 = cos(2*x5);
X    s5x7 = sin(5*x7);
X    x8 = x4-x3;
X    s2x8 = sin(2*x8);
X    c2x8 = cos(2*x8);
X    s3x8 = sin(3*x8);
X    c3x8 = cos(3*x8);
X
X    *dml = 7.581e-3*s2x5-7.986e-3*sx6-1.48811e-1*sx7-4.0786e-2*s2x7-
X          (8.14181e-1-(1.815e-2-1.6714e-2*x1)*x1)*sx5-
X          (1.0497e-2-(1.60906e-1-4.1e-3*x1)*x1)*cx5-1.5208e-2*s3x7-
X          6.339e-3*s4x7-6.244e-3*sx3-1.65e-2*s2x7*sx3+
X          (8.931e-3+2.728e-3*x1)*sx7*sx3-5.775e-3*s3x7*sx3+
X          (8.1344e-2+3.206e-3*x1)*cx7*sx3+1.5019e-2*c2x7*sx3;
X    *dml += (8.5581e-2+2.494e-3*x1)*sx7*cx3+1.4394e-2*c2x7*cx3+
X          (2.5328e-2-3.117e-3*x1)*cx7*cx3+
X          6.319e-3*c3x7*cx3+6.369e-3*sx7*s2x3+9.156e-3*s2x7*s2x3+
X          7.525e-3*s3x8*s2x3-5.236e-3*cx7*c2x3-7.736e-3*c2x7*c2x3-
X          7.528e-3*c3x8*c2x3;
X    *dml = degrad(*dml);
X
X    *ds = (-7927+(2548+91*x1)*x1)*sx5+(13381+(1226-253*x1)*x1)*cx5+
X         (248-121*x1)*s2x5-(305+91*x1)*c2x5+412*s2x7+12415*sx3+
X         (390-617*x1)*sx7*sx3+(165-204*x1)*s2x7*sx3+26599*cx7*sx3-
X         4687*c2x7*sx3-1870*c3x7*sx3-821*c4x7*sx3-
X         377*c5x7*sx3+497*c2x8*sx3+(163-611*x1)*cx3;
X    *ds += -12696*sx7*cx3-4200*s2x7*cx3-1503*s3x7*cx3-619*s4x7*cx3-
X         268*s5x7*cx3-(282+1306*x1)*cx7*cx3+(-86+230*x1)*c2x7*cx3+
X         461*s2x8*cx3-350*s2x3+(2211-286*x1)*sx7*s2x3-
X         2208*s2x7*s2x3-568*s3x7*s2x3-346*s4x7*s2x3-
X         (2780+222*x1)*cx7*s2x3+(2022+263*x1)*c2x7*s2x3+248*c3x7*s2x3+
X         242*s3x8*s2x3+467*c3x8*s2x3-490*c2x3-(2842+279*x1)*sx7*c2x3;
X    *ds += (128+226*x1)*s2x7*c2x3+224*s3x7*c2x3+
X         (-1594+282*x1)*cx7*c2x3+(2162-207*x1)*c2x7*c2x3+
X         561*c3x7*c2x3+343*c4x7*c2x3+469*s3x8*c2x3-242*c3x8*c2x3-
X         205*sx7*s3x3+262*s3x7*s3x3+208*cx7*c3x3-271*c3x7*c3x3-
X         382*c3x7*s4x3-376*s3x7*c4x3;
X    *ds *= 1e-7;
X
X    dp = (7.7108e-2+(7.186e-3-1.533e-3*x1)*x1)*sx5-7.075e-3*sx7+
X         (4.5803e-2-(1.4766e-2+5.36e-4*x1)*x1)*cx5-7.2586e-2*cx3-
X         7.5825e-2*sx7*sx3-2.4839e-2*s2x7*sx3-8.631e-3*s3x7*sx3-
X         1.50383e-1*cx7*cx3+2.6897e-2*c2x7*cx3+1.0053e-2*c3x7*cx3-
X         (1.3597e-2+1.719e-3*x1)*sx7*s2x3+1.1981e-2*s2x7*c2x3;
X    dp += -(7.742e-3-1.517e-3*x1)*cx7*s2x3+
X         (1.3586e-2-1.375e-3*x1)*c2x7*c2x3-
X         (1.3667e-2-1.239e-3*x1)*sx7*c2x3+
X         (1.4861e-2+1.136e-3*x1)*cx7*c2x3-
X         (1.3064e-2+1.628e-3*x1)*c2x7*c2x3;
X
X    *dm = *dml-(degrad(dp)/s);
X
X    *da = 572*sx5-1590*s2x7*cx3+2933*cx5-647*s3x7*cx3+33629*cx7-
X         344*s4x7*cx3-3081*c2x7+2885*cx7*cx3-1423*c3x7+
X         (2172+102*x1)*c2x7*cx3-671*c4x7+296*c3x7*cx3-320*c5x7-
X         267*s2x7*s2x3+1098*sx3-778*cx7*s2x3-2812*sx7*sx3;
X    *da += 495*c2x7*s2x3+688*s2x7*sx3+250*c3x7*s2x3-393*s3x7*sx3-
X         856*sx7*c2x3-228*s4x7*sx3+441*s2x7*c2x3+2138*cx7*sx3+
X         296*c2x7*c2x3-999*c2x7*sx3+211*c3x7*c2x3-642*c3x7*sx3-
X         427*sx7*s3x3-325*c4x7*sx3+398*s3x7*s3x3-890*cx3+
X         344*cx7*c3x3+2206*sx7*cx3-427*c3x7*c3x3;
X    *da *= 1e-6;
X
X    *dhl = 7.47e-4*cx7*sx3+1.069e-3*cx7*cx3+2.108e-3*s2x7*s2x3+
X          1.261e-3*c2x7*s2x3+1.236e-3*s2x7*c2x3-2.075e-3*c2x7*c2x3;
X    *dhl = degrad(*dhl);
X}
X
X/* ....uranus */
Xstatic
Xp_uranus (t, s, dl, dr, dml, ds, dm, da, dhl)
Xdouble t, s;
Xdouble *dl, *dr, *dml, *ds, *dm, *da, *dhl;
X{
X    double dp;
X    double x1, x2, x3, x4, x5, x6;
X    double x8, x9, x10, x11, x12;
X    double sx4, cx4, s2x4, c2x4;
X    double sx9, cx9, s2x9, c2x9;
X    double sx11, cx11;
X
X    aux_jsun (t, &x1, &x2, &x3, &x4, &x5, &x6);
X
X        x8 = mod2PI(1.46205+3.81337*t);
X        x9 = 2*x8-x4;
X    sx9 = sin(x9);
X    cx9 = cos(x9);
X        s2x9 = sin(2*x9);
X    c2x9 = cos(2*x9);
X
X    x10 = x4-x2;
X    x11 = x4-x3;
X    x12 = x8-x4;
X
X    *dml = (8.64319e-1-1.583e-3*x1)*sx9+(8.2222e-2-6.833e-3*x1)*cx9+
X          3.6017e-2*s2x9-3.019e-3*c2x9+8.122e-3*sin(x6);
X    *dml = degrad(*dml);
X
X    dp = 1.20303e-1*sx9+6.197e-3*s2x9+(1.9472e-2-9.47e-4*x1)*cx9;
X    *dm = *dml-(degrad(dp)/s);
X
X    *ds = (163*x1-3349)*sx9+20981*cx9+1311*c2x9;
X    *ds *= 1e-7;
X
X    *da = -3.825e-3*cx9;
X
X    *dl = (1.0122e-2-9.88e-4*x1)*sin(x4+x11)+
X         (-3.8581e-2+(2.031e-3-1.91e-3*x1)*x1)*cos(x4+x11)+
X         (3.4964e-2-(1.038e-3-8.68e-4*x1)*x1)*cos(2*x4+x11)+
X         5.594e-3*sin(x4+3*x12)-1.4808e-2*sin(x10)-
X         5.794e-3*sin(x11)+2.347e-3*cos(x11)+9.872e-3*sin(x12)+
X         8.803e-3*sin(2*x12)-4.308e-3*sin(3*x12);
X
X    sx11 = sin(x11);
X    cx11 = cos(x11);
X    sx4 = sin(x4);
X    cx4 = cos(x4);
X    s2x4 = sin(2*x4);
X    c2x4 = cos(2*x4);
X    *dhl = (4.58e-4*sx11-6.42e-4*cx11-5.17e-4*cos(4*x12))*sx4-
X          (3.47e-4*sx11+8.53e-4*cx11+5.17e-4*sin(4*x11))*cx4+
X          4.03e-4*(cos(2*x12)*s2x4+sin(2*x12)*c2x4);
X    *dhl = degrad(*dhl);
X
X    *dr = -25948+4985*cos(x10)-1230*cx4+3354*cos(x11)+904*cos(2*x12)+
X         894*(cos(x12)-cos(3*x12))+(5795*cx4-1165*sx4+1388*c2x4)*sx11+
X         (1351*cx4+5702*sx4+1388*s2x4)*cos(x11);
X    *dr *= 1e-6;
X}
X
X/* ....neptune */
Xstatic
Xp_neptune (t, s, dl, dr, dml, ds, dm, da, dhl)
Xdouble t, s;
Xdouble *dl, *dr, *dml, *ds, *dm, *da, *dhl;
X{
X    double dp;
X    double x1, x2, x3, x4, x5, x6;
X    double x8, x9, x10, x11, x12;
X    double sx8, cx8;
X    double sx9, cx9, s2x9, c2x9;
X    double s2x12, c2x12;
X
X    aux_jsun (t, &x1, &x2, &x3, &x4, &x5, &x6);
X
X        x8 = mod2PI(1.46205+3.81337*t);
X        x9 = 2*x8-x4;
X    sx9 = sin(x9);
X    cx9 = cos(x9);
X        s2x9 = sin(2*x9);
X    c2x9 = cos(2*x9);
X
X    x10 = x8-x2;
X    x11 = x8-x3;
X    x12 = x8-x4;
X
X    *dml = (1.089e-3*x1-5.89833e-1)*sx9+(4.658e-3*x1-5.6094e-2)*cx9-
X          2.4286e-2*s2x9;
X    *dml = degrad(*dml);
X
X    dp = 2.4039e-2*sx9-2.5303e-2*cx9+6.206e-3*s2x9-5.992e-3*c2x9;
X
X    *dm = *dml-(degrad(dp)/s);
X
X    *ds = 4389*sx9+1129*s2x9+4262*cx9+1089*c2x9;
X    *ds *= 1e-7;
X
X    *da = 8189*cx9-817*sx9+781*c2x9;
X    *da *= 1e-6;
X
X    s2x12 = sin(2*x12);
X    c2x12 = cos(2*x12);
X    sx8 = sin(x8);
X    cx8 = cos(x8);
X    *dl = -9.556e-3*sin(x10)-5.178e-3*sin(x11)+2.572e-3*s2x12-
X         2.972e-3*c2x12*sx8-2.833e-3*s2x12*cx8;
X
X    *dhl = 3.36e-4*c2x12*sx8+3.64e-4*s2x12*cx8;
X    *dhl = degrad(*dhl);
X
X    *dr = -40596+4992*cos(x10)+2744*cos(x11)+2044*cos(x12)+1051*c2x12;
X    *dr *= 1e-6;
X}
X
END_OF_FILE
  if test 17647 -ne `wc -c <'plans.c'`; then
    echo shar: \"'plans.c'\" unpacked with wrong size!
  fi
  # end of 'plans.c'
fi
echo shar: End of archive 5 \(of 9\).
cp /dev/null ark5isdone
MISSING=""
for I in 1 2 3 4 5 6 7 8 9 ; do
    if test ! -f ark${I}isdone ; then
    MISSING="${MISSING} ${I}"
    fi
done
if test "${MISSING}" = "" ; then
    echo You have unpacked all 9 archives.
    rm -f ark[1-9]isdone ark[1-9][0-9]isdone
else
    echo You still must unpack the following archives:
    echo "        " ${MISSING}
fi
exit 0
exit 0 # Just in case...