home *** CD-ROM | disk | FTP | other *** search
/ Club Amiga de Montreal - CAM / CAM_CD_1.iso / files / 383.lha / Planets / planet.c < prev    next >
Encoding:
C/C++ Source or Header  |  1990-05-25  |  15.6 KB  |  664 lines
  1. /*  Planets
  2.  *           A program to determine the position of
  3.  *         the Sun., Mer., Ven., Mars, Jup, Sat & Ura
  4.  *
  5.  *  reference Jean Meeus's book, " Astronomical Formulae for
  6.  *  Calculators
  7.  *            program by
  8.  *
  9.  *            F.T. Mendenhall
  10.  *            ihnp4!inuxe!fred
  11.  *
  12.  *        Copyright 1985 by F.T. Mendenhall
  13.  *        all rights reserved
  14.  *
  15.  * (This copyright prevents you from selling the program - the
  16.  *  author grants anyone the right to copy and install the program
  17.  *  on any machine it will run on)
  18.  !
  19.  ! Modified by Alan Paeth, awpaeth@watcgl, January 1987 for
  20.  !   (1) non-interactive mode (uses current GMT)
  21.  !   (2) output in simplified Yale format for use with star charter software
  22.  !   (3) lines 365-380 rewritten as a simplier C expression (awpaeth, Apr '87)
  23.  !
  24.  ! Modified by Petri Launiainen, pl@intrin.FI, February 1987 for
  25.  !   easier LOGFILE definition/modification in Makefile
  26.  !
  27.  ! Modified by Alan Paeth, September 1987, for command line flags.
  28.  !    The program continues to provide reduced Yale output in the old style,
  29.  !     but new versions of "starchart" accept either. The old format allows
  30.  !     representation of magnitudes below -.99 (useful for planets), but no
  31.  !     spectral or stellar identification, which is not useful for planets.
  32.  !
  33.  ! Modified by Jim Cobb, March 1988, to compute moon positions as well.
  34.  !    Also, to  change slightly the computation of epli, the obliquity
  35.  !     of the ecliptic.  Formerly epli had a term 
  36.  !        +0.00256*cos(omeg*rad).
  37.  !     This term is a compensation for computing the apparent position
  38.  !     of the sun (Chapter 18 in Meeus).  But for converting
  39.  !     ecliptical coordinates to right ascension and declination we want
  40.  !     the mean value (that is, we don't want this term).
  41.  ! 
  42.  ! Modified by Bob Leivian, Sept 1989, to compute local azimuth and altitude
  43.  !     given a lat/lon, also fixed time so you can replace a single item
  44.  !     without having to supply all (i.e. for today at 9:00 put -t 9) 
  45.  !     also ran the code through 'cb' to standardize the indention.
  46.  !     I split the code in to more managible pieces.
  47.  !     Also printed out time in local and UTC as well as JD
  48.  !
  49.  ! See "Time Notes" section below if you plan to rewrite the user interface.
  50.  */
  51.  
  52. #include <stdio.h>
  53. #include <math.h>
  54.  
  55. /* default latitude longitude (of Tempe, AZ) -- replace with your favorite place */
  56. #define DEFLONG -111.94
  57. #define DEFLAT 33.32
  58. #define DEFZONE 7 /* MST */
  59.  
  60. #ifndef PLANETFILE
  61. #define PLANETFILE "./planet.star"
  62. #endif
  63.  
  64. double    pie, rad;
  65. double    htod(), atof(), kepler(), truean();
  66. double    longi(), lati(), poly(), to_int(), range();
  67. int    cls(), helio_trans(), geo_trans(), speak();
  68.  
  69. double    planet_pos();
  70. double    sidereal_time();
  71. int    zone;
  72. double    longitude = DEFLONG;
  73. double    latitude = DEFLAT;
  74. double    GHA_Aries;
  75.  
  76.  
  77. FILE    *logfile;
  78. char    *progname;
  79.  
  80. /*
  81.  * Times Notes.
  82.  * 
  83.  * This software would ideally use a generic UNIX system call to which would
  84.  * provide day and date info as integers, plus current time west of GMT.
  85.  * These values would become the default settings for each of the -m -d -y -t
  86.  * & -z. Unfortunately, time on across many UNIX systems does not (to my
  87.  * knowledge) exist as a uniform subroutine call.
  88.  * 
  89.  * Instead, we presume that the very low level "gettimeofday" (GMT in seconds,
  90.  * from 1 Jan 1970) exists on all UNIX installations, and use this subroutine
  91.  * to use the present time of day when no command line parameters are
  92.  * present, as defaults for the current time have not been filled in.
  93.  * 
  94.  * When any flag-style command line parameters are present, hardcoded defaults
  95.  * are substitude for -m, -d, -y -t switches, which may then be overridden.
  96.  * The -z flag is filled in by a companion return value from "gettimeofday".
  97.  * This approach makes for poor (on account of hard coding) defaults.
  98.  * 
  99.  * A worthwhile change would be to simplify the command line control structure
  100.  * and default mechanism by using a higher level system call, or (safer),
  101.  * write the routine to convert Unix GMT into day, month, year, so that it
  102.  * lives within this program. This would allow for more reasonable defaults,
  103.  * and still merely one UNIX system call to get GMT time. Perhaps sources
  104.  * from a UNIX system may be studied to do this correctly.
  105.  * 
  106.  * Expertise in "time" on a specific UNIX system is *NOT* what is required --
  107.  * what *IS* required is knowledge and access to enough independent UNIX
  108.  * systems to insure (to a high probability) that the code remains portable.
  109.  * In other words, only low-level system calls known to work across a large
  110.  * range of systems should be employed.
  111.  * 
  112.  */
  113.  
  114. #include <sys/time.h>        /* for getting current GMT */
  115.  
  116. double
  117. current_time(m, d, y, t, z)
  118. int    *m;
  119. int    *d;
  120. int    *y;
  121. double    *t;
  122. int    *z;
  123. {
  124. #define SECSPERDAY (60.0 * 60.0 * 24.0)
  125.  
  126. #ifdef UNIX
  127. #define GMT1970 2440587.5
  128.     struct timeval tv;
  129.     struct timezone tz;
  130.     struct tm *p;
  131.  
  132.     gettimeofday(&tv, &tz);
  133.     p = localtime(&tv);
  134.  
  135.     *y = p->tm_year + 1900;            /* local year */
  136.     *m = p->tm_mon + 1;                /* local month */
  137.     *d = p->tm_mday;                /* local day of month */
  138.  
  139.     *t = p->tm_hour + (p->tm_min / 60.0);    /* local time (0.0-23.999) */
  140.     *z = tz.tz_minuteswest;        /* local time hrs west of Greenwich */
  141.  
  142.     return (GMT1970 + tv.tv_sec / SECSPERDAY);
  143. #endif 
  144.  
  145. #ifdef AMIGA
  146.     time_t amiga_now;
  147.     struct tm *p;
  148.  
  149. #define GMT1978 2443509.5
  150.     amiga_now = time(&amiga_now);
  151.  
  152.     p = localtime(&amiga_now);
  153.  
  154.     *t = p->tm_hour + (p->tm_min / 60.0);    /* local time (0.0-23.999) */
  155.     *z = DEFZONE * 60;
  156.     *y = p->tm_year + 1900;            /* local year */
  157.     *m = p->tm_mon + 1;                /* local month */
  158.     *d = p->tm_mday;                /* local day of month */
  159.     return GMT1978 + (amiga_now / SECSPERDAY) + (DEFZONE / 24.0);
  160. #endif
  161.  
  162. #ifdef NO_TIME_AVAIL
  163. #define MST1989 2447528.291667
  164.     /* no time is avail, pick a default i.e. Noon Jan 1 1989 */
  165.     *y = 1989;
  166.     *m = 1;
  167.     *d = 1;
  168.  
  169.     *t = 12.0;
  170.     *z = DEFZONE * 60.;
  171.     return MST1989;
  172. #endif
  173. }
  174.  
  175. double
  176. commandtime(argc, argv)
  177. char    **argv;
  178. {
  179.     int    mm, dd, yy, aa1, bb1;
  180.     double    jd, year, month, day, time, timez;
  181.     int    j;
  182.     double    now;
  183.     static char    *usage = 
  184. "\nusage:planet {current time used}\nor\tplanet juliandate\nor\tplanet [-t hh.mm -z zone -y year -m mon -d day ]";
  185.  
  186.     /* get current time for defaults */
  187.     now = current_time(&mm, &dd, &yy, &time, &zone);
  188.     timez = zone / 60.; /* time zone in hours */
  189.  
  190.     /* no args - use current time */
  191.     if (argc == 1)
  192.         return now;
  193.  
  194.     if ((argv[1][0] != '-')) {
  195.         /* one numeric arg - use it as the julian date */
  196.         if (argc == 2)
  197.             return (atof(argv[1]));
  198.         die("no switches - %s", usage);
  199.     } else {
  200.         /* flags present, set up defaults */
  201.  
  202.         /* process any user overrides */
  203.         for (j = 1; j < argc; j++) {
  204.             if (argv[j][0] != '-')
  205.                 die("unknown argument - %s", usage /*argv[j]*/);
  206.             switch (argv[j][1]) {
  207.             case 't': 
  208.                 time = htod(argv[++j]); 
  209.                 break;
  210.             case 'z': 
  211.                 timez = htod(argv[++j]);
  212.                 if (timez < 0 || timez > 24)
  213.                     die("time zone must be 0 to 24 west of UTC\n");
  214.                                 zone = timez * 60.;
  215.                 break;
  216.             case 'y': 
  217.                 yy = atoi(argv[++j]); 
  218.                 break;
  219.             case 'm': 
  220.                 mm = atoi(argv[++j]);
  221.                 if (mm < 1 || mm > 12)
  222.                     die("month must be 1..12\n");
  223.                 break;
  224.             case 'd': 
  225.                 dd = atoi(argv[++j]);
  226.                 if (dd < 1 || dd > 31)
  227.                     die("day must be 1..31\n");
  228.                 break;
  229.             case 'l':
  230.                 switch (argv[j][2]) {
  231.                 case 'o':
  232.                     longitude = atof(argv[++j]);
  233.                     if (longitude < -180 || longitude > 180)
  234.                         die("longitude must be +-180\n");
  235.                     break;
  236.                 case 'a':
  237.                     latitude = atof(argv[++j]);
  238.                     if (latitude < -90 || latitude > 90)
  239.                         die("latitude must be +-90\n");
  240.                     break;
  241.                 default:
  242.                     die("specify lat or lon\n");
  243.                 }
  244.                 break;
  245.             default:
  246.                 die("unknown switch - %s", usage /*argv[j]*/ );
  247.             }
  248.             if (j == argc)
  249.                 die("trailing command line flag - %s", argv[j - 1]);
  250.         }
  251.     }
  252.     day = ((float) dd + (time + timez) / 24.0);
  253.     if (mm > 2) {
  254.         year = yy;
  255.         month = mm;
  256.     } else {
  257.         year = yy - 1;
  258.         month = mm + 12;
  259.     }
  260.     aa1 = (int) (year / 100);
  261.     bb1 = 2 - aa1 + (int) (aa1 / 4);
  262.     jd = to_int(365.25 * year) + to_int(30.6001 * (month + 1));
  263.     jd = jd + day + 1720994.5;
  264.     if ((year + month / 100) > 1582.10)
  265.         jd = jd + bb1;
  266.     return (jd);
  267. }
  268.  
  269.  
  270. caltim(jd_frac, hour, min)
  271. double    jd_frac;
  272. int    *hour;
  273. int    *min;
  274. {
  275.     double    x;
  276.  
  277.     x = jd_frac * 24.0;
  278.     *hour = (int) floor(x);
  279.     x -= floor(x);
  280.     x *= 60.0;
  281.     *min = (int) floor(x);
  282. }
  283.  
  284.  
  285. main(argc, argv)
  286. int    argc;
  287. char    **argv;
  288. {
  289.     double    jd, T0, epli;
  290.     double    temp;
  291.         long    int_jd;
  292.     int    y, m, d, h, min, lh, lmin;
  293.     char    *str;
  294.  
  295.     progname = argv[0];
  296. #ifdef AMIGA
  297.     logfile = 0; /* don't write a log */
  298. #else
  299.     if ((logfile = fopen(PLANETFILE, "w")) == 0)
  300.         /*if ((logfile = fopen("/tmp/planet.star", "w")) == 0) {
  301.             
  302.             die("fail on opening logfile %s\n", "/tmp/planet.star");
  303.         }*/
  304. #endif
  305.  
  306.     cls();
  307.  
  308.     pie = 3.1415926535898;
  309.     rad = pie / 180.0;
  310.  
  311.     /*  calculate time T0 from 1900 January 0.5 ET */
  312.     jd = commandtime(argc, argv);
  313.  
  314.     /* jd to calendar day */
  315.         int_jd = jd +.5001;
  316.     caldat(int_jd, &m, &d, &y);
  317.  
  318.     /* julian day to UTC */
  319.     temp = jd - .4999;
  320.     caltim(temp - floor(temp), &h, &min);
  321.     printf("At %d/%d %d %d:%02d UTC (Julian: %.3f)\n", m, d, y, h, min, jd);
  322.  
  323.     /* UTC to Local */
  324.     temp -= (zone / 1440.);
  325.     caltim(temp - floor(temp), &lh, &lmin);
  326.         int_jd = jd + .5001 - (zone / 1440.);
  327.     caldat(int_jd, &m, &d, &y);
  328.     if (lh > 12) {
  329.         lh -= 12;
  330.         str = "pm";
  331.     } else {
  332.         str = "am";
  333.         if (lh == 0)
  334.             lh = 12;
  335.     }
  336.     printf("local %d/%d %d:%02d%s (", m, d, lh, lmin, str);
  337.  
  338.     /* compute exact local time for a given longitude */
  339.     temp = jd - .5 + (longitude / 360.);
  340.     caltim(temp - floor(temp), &lh, &lmin);
  341.     GHA_Aries =  15. * sidereal_time(jd);
  342.  
  343.     printf("%2d:%02d true local at lat:%1.1f lon:%1.1f, GHA Aries:%1.1f)\n",
  344.         lh, lmin, latitude, longitude, GHA_Aries);
  345.  
  346.     T0 = (jd - 2415020.0) / 36525.0;
  347.     printf("\nOBJECT      R.A.              DEC          DISTANCE   altitude   azimuth\n");
  348.  
  349.     epli = planet_pos(jd, T0);
  350.  
  351.     moon_pos(T0, epli);
  352.  
  353.     putchar('\n');
  354.     if (logfile)
  355.         close(logfile);
  356.         logfile = 0;
  357.  
  358.     /* show local position of some stars also */
  359.     speak(101.1714, -16.7013, 0., -1.46, " " , "Sirius");
  360.     speak(213.7956, 19.236, 0., -0.04, " " , "Arcturus");
  361.     speak(88.6508, 7.4057, 0., 0.5, " " , "Betelgeuse");
  362.     exit(0);
  363. }                /* end of program main */
  364.  
  365.  
  366. int
  367. cls()
  368. {
  369. #define CLEARCNT 1
  370.     int    i;
  371.     for (i = 0; i < CLEARCNT; i++)
  372.         putchar('\n');
  373. }
  374.  
  375.  
  376. double
  377. poly(a0, a1, a2, a3, t)
  378. double    a0, a1, a2, a3, t;
  379. {
  380.     return (a0 + a1 * t + a2 * t * t + a3 * t * t * t);
  381. }
  382.  
  383.  
  384. double
  385. to_int(z)
  386. double    z;
  387. {
  388.     long    trunk;
  389.  
  390.     trunk = (long) z;
  391.     z = (double) trunk;
  392.     return (z);
  393. }
  394.  
  395.  
  396. double
  397. kepler(e, M)
  398. double    e, M;
  399. {
  400.     double    corr, e0, E0, E1;
  401.     e0 = e / rad;
  402.     corr = 1;
  403.     E0 = M;
  404.     while (corr > 0.000001) {
  405.         corr = (M + e0 * sin(E0 * rad) - E0) / (1 - e * cos(E0 * rad));
  406.         E1 = E0 + corr;
  407.         if (corr < 0) {
  408.             corr = -1.0 * corr;
  409.         }
  410.         E0 = E1;
  411.     }
  412.  
  413.     return (E1);
  414. }
  415.  
  416.  
  417. double
  418. truean(e, E)
  419. double    e, E;
  420. {
  421.     double    nu;
  422.  
  423.     nu = 2.0 * atan(sqrt((1 + e) / (1 - e)) * tan(E * rad / 2));
  424.     nu = nu / rad;
  425.     if (nu < 0.0) {
  426.         nu = nu + 360.0;
  427.     }
  428.     if (fabs(nu - E) > 90.0) {
  429.         if (nu > 180.0) {
  430.             nu = nu - 180.0;
  431.         } else {
  432.             nu = nu + 180.0;
  433.         }
  434.     }
  435.     return (nu);
  436. }
  437.  
  438.  
  439. double
  440. longi(w2, i, u)
  441. double    w2, i, u;
  442. {
  443.     double    x, y, l;
  444.  
  445.     y = cos(i * rad) * sin(u * rad);
  446.     x = cos(u * rad);
  447.     l = atan2(y, x);
  448.     l = l / rad;
  449.     if (l < 0.0) {
  450.         l = l + 360.0;
  451.     }
  452.     return (l + w2);
  453. }
  454.  
  455.  
  456. double
  457. lati(u, i)
  458. double    u, i;
  459. {
  460.     double    b;
  461.     b = asin(sin(u * rad) * sin(i * rad)) / rad;
  462.     return (b);
  463. }
  464.  
  465.  
  466. int
  467. speak(ra, dec, dis, mag, sym, str)
  468. double    ra, dec, dis, mag;
  469. char    *sym, *str;
  470. {
  471.     int    h, m, s, deg;
  472.     double    lha, sha, altitude, azimuth;
  473.     ;
  474.  
  475.     if (ra < 0) {
  476.         ra = ra + 360.0;
  477.     }
  478.     sha = 360. - ra;
  479.     ra = ra / 15.0;        /* convert from degs to hours */
  480.     h = (int) to_int(ra);
  481.     m = (int) ((ra - h) * 60);
  482.     s = (int) (((ra - h) * 60 - m) * 60);
  483.     printf("%-12s%2dh %2dm %2ds ", str, h, m, s);
  484.     if (logfile) 
  485.         fprintf(logfile, "%02d%02d%02d", h, m, s);
  486.     deg = (int) to_int(dec);
  487.     m = (int) ((dec - deg) * 60);
  488.     s = (int) (((dec - deg) * 60 - m) * 60);
  489.     if (m < 0) {
  490.         m = m * -1;
  491.     }
  492.     if (s < 0) {
  493.         s = s * -1;
  494.     }
  495. #ifdef AMIGA
  496.     /* its printf doesn't know about signed */
  497.     if (deg > 0)
  498.         printf("   +%2ddeg %2dm %2ds ", deg, m, s);
  499.     else
  500.         printf("   %3ddeg %2dm %2ds ", deg, m, s);
  501. #else
  502.     printf("   %+3ddeg %2dm %2ds ", deg, m, s);
  503. #endif
  504.  
  505.     if (dis > 0.0001)
  506.         printf("  %6.3f", dis);
  507.     else
  508.         printf("        ");
  509.     if (logfile)
  510.         if (mag < 0.0)
  511.             fprintf(logfile, "%+03d%02d-%2d%s       %s\n",
  512.                 deg, m, -(int) (10 * mag), sym, str);
  513.         else
  514.             fprintf(logfile, "%+03d%02d%3d%s       %s\n",
  515.                 deg, m, (int) (100 * mag), sym, str);
  516.  
  517.     /* if we have a valid angle */
  518.     if (GHA_Aries > 0) {
  519.         /* now tell where it is at a given place on earth */
  520.         lha = range(GHA_Aries + sha + longitude);
  521.  
  522.         altitude = sin(latitude * rad) * sin(dec * rad) + 
  523.             cos(latitude * rad) * cos(dec * rad) * cos(lha * rad);
  524.         altitude = asin(altitude) / rad;
  525.  
  526.         azimuth = sin(lha * rad) / 
  527.             (cos(lha * rad) * sin(latitude * rad) - 
  528.             tan(dec * rad) * cos(latitude * rad));
  529.  
  530.         /* correct for quadrant */
  531.         if (lha <= 180.) {
  532.             if (azimuth > 0)
  533.                 azimuth = 180. + atan(azimuth) / rad;
  534.             else
  535.                 azimuth = 360. + atan(azimuth) / rad;
  536.         } else {
  537.             if (azimuth <= 0.)
  538.                 azimuth = 180. + atan(azimuth) / rad;
  539.             else
  540.                 azimuth = atan(azimuth) / rad;
  541.         }
  542.  
  543.         printf("      %5.1f     %5.1f\n", altitude, azimuth);
  544.     } else
  545.         printf("\n");
  546.     return (0);
  547. }
  548.  
  549.  
  550. double
  551. range(val)
  552. double    val;
  553. {
  554.     while (val < 0.0) {
  555.         val = val + 360.0;
  556.     }
  557.     if (val > 360.0) {
  558.         val = val - (to_int(val / 360.0) * 360.0);
  559.     }
  560.     return (val);
  561. }
  562.  
  563.  
  564. /*
  565.  * Helio_trans converts heliocentric ecliptical coordinates to geocentric
  566.  * right ascension and declination.  It also outputs the converted value.
  567.  */
  568. int
  569. helio_trans(r, b, ll, Stheta, Sr, epli, mag, sym, str)
  570. double    r, b, ll, Stheta, Sr, epli, mag;
  571. char    *sym, *str;
  572. {
  573.     double    N, D, lambda, delta, beta;
  574.  
  575.     N = r * cos(b * rad) * sin((ll - Stheta) * rad);
  576.     D = r * cos(b * rad) * cos((ll - Stheta) * rad) + Sr;
  577.     lambda = atan2(N, D) / rad;
  578.     if (lambda < 0.0) {
  579.         lambda = lambda + 360.0;
  580.     }
  581.     lambda = range(lambda + Stheta);
  582.     delta = sqrt(N * N + D * D + (r * sin(b * rad)) * (r * sin(b * rad)));
  583.     beta = asin(r * sin(b * rad) / delta) / rad;
  584.     geo_trans(lambda, beta, epli, delta, mag, sym, str);
  585.     return (0);
  586. }
  587.  
  588.  
  589. /*
  590.  * Geo_trans converts geocentric ecliptical coordinates to geocentric right
  591.  * ascension and declination.  It also outputs the converted value.  Note
  592.  * that the output coordinates are referred to the mean equinox of date.  If
  593.  * these coordinates are to be used in conjunction with a star database, use
  594.  * the epoch program to convert the planet's coordinates to the epoch for the
  595.  * star database.
  596.  */
  597. int
  598. geo_trans(lambda, beta, epli, delta, mag, sym, str)
  599. double    lambda, beta, epli, delta, mag;
  600. char    *sym, *str;
  601. {
  602.     double    N, D, RA, DEC;
  603.  
  604.     N = sin(lambda * rad) * cos(epli * rad)
  605.     -tan(beta * rad) * sin(epli * rad);
  606.     D = cos(lambda * rad);
  607.     RA = atan2(N, D) / rad;
  608.     DEC = asin(sin(beta * rad) * cos(epli * rad)
  609.          + cos(beta * rad) * sin(epli * rad) * sin(lambda * rad)) / rad;
  610.     speak(RA, DEC, delta, mag, sym, str);
  611.     return (0);
  612. }
  613.  
  614.  
  615. double
  616. htod(s)
  617. char    *s;
  618. {
  619.     /*
  620.      * htod(str) reads floating point strings of the form {+-}hh.{mm} thus
  621.      * allowing for fractional values in base 60 (ie, degrees/minutes).
  622.      */
  623.     double    x, sign;
  624.     int    full, frac;
  625.     char    *t;
  626.     t = s - 1;
  627.     while (*++t) {
  628.         if (((*t < '0') || (*t > '9')) && (*t != '.') && (*t != '+') && (*t != '-'))
  629.             die("non-digit in dd.mm style numeric argument");
  630.     }
  631.     if (s == NULL)
  632.         return 0.0;
  633.     full = frac = 0;
  634.     sign = 1.0;
  635.     if (*s == '-') {
  636.         sign = -1.0;
  637.         s++;
  638.     } else if (*s == '+')
  639.         s++;
  640.     while (*s && *s != '.')
  641.         full = 10 * full + *s++ - '0';
  642.     if (*s++ == '.') {
  643.         if (*s)
  644.             frac = 10 * (*s++ - '0');
  645.         if (*s)
  646.             frac += *s++ - '0';
  647.         if (frac > 59)
  648.             frac = 59;
  649.     }
  650.     x = (double) full + ((double) frac) / 60.0;
  651.     return sign * x;
  652. }
  653.  
  654.  
  655. die(a, b)
  656. char    *a, *b;
  657. {
  658.     fprintf(stderr, "%s: ", progname);
  659.     fprintf(stderr, a, b);
  660.     fprintf(stderr, "\n");
  661.     fflush(stderr);
  662.     exit(1);
  663. }
  664.