allbery@uunet.UU.NET (Brandon S. Allbery - comp.sources.misc) (11/28/89)
Posting-number: Volume 9, Issue 40 Submitted-by: ecd@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