sys/src/cmd/astro/stars.c

   1 #include "astro.h"
   2
   3 char*   startab = "/lib/sky/estartab";
   4
   5 void
   6 stars(void)
   7 {
   8         double lomoon, himoon, sd;
   9         int wrap, f, i;
  10         char *saop;
  11         static char saoa[100];
  12
  13         sd = 1000*radsec;
  14         lomoon = omoon.point[0].ra - sd;
  15         if(lomoon < 0)
  16                 lomoon += pipi;
  17         himoon = omoon.point[NPTS+1].ra + sd;
  18         if(himoon > pipi)
  19                 himoon -= pipi;
  20         lomoon *= 12/pi;
  21         himoon *= 12/pi;
  22         wrap = 0;
  23         if(lomoon > himoon)
  24                 wrap++;
  25
  26         f = open(startab, OREAD);
  27         if(f < 0) {
  28                 fprint(2, "%s?\n", startab);
  29                 return;
  30         }
  31         epoch = 1950.0;
  32         epoch = (epoch-1900.0) * 365.24220 + 0.313;
  33         saop = saoa;
  34
  35 /*
  36  *      read mean places of stars at epoch of star table
  37  */
  38
  39 loop:
  40         if(rline(f)) {
  41                 close(f);
  42                 return;
  43         }
  44         rah = atof(line+17);
  45         ram = atof(line+20);
  46         ras = atof(line+23);
  47
  48         alpha = rah + ram/60 + ras/3600;
  49         if(wrap == 0) {
  50                 if(alpha < lomoon || alpha > himoon)
  51                         goto loop;
  52         } else
  53                 if(alpha < lomoon && alpha > himoon)
  54                         goto loop;
  55
  56         sao = atof(line+0);
  57         sprint(saop, "%ld", sao);
  58         da = atof(line+30);
  59         dday = atof(line+37);
  60         dmin = atof(line+41);
  61         dsec = atof(line+44);
  62         dd = atof(line+50);
  63         px = atof(line+57);
  64         mag = atof(line+61);
  65
  66 /*
  67  *      convert rt ascension and declination to internal format
  68  */
  69
  70         delta = fabs(dday) + dmin/60 + dsec/3600;
  71         if(dday < 0)
  72                 delta = -delta;
  73
  74         star();
  75 /*
  76  *      if(fabs(beta) > 6.55*radian)
  77  *              goto loop;
  78  */
  79         sd = .0896833e0*cos(beta)*sin(lambda-1.3820+.00092422117*eday)
  80                  + 0.99597*sin(beta);
  81         if(fabs(sd) > .0183)
  82                 goto loop;
  83
  84         for(i=0; i<=NPTS+1; i++)
  85                 setobj(&ostar.point[i]);
  86
  87         occult(&omoon, &ostar, 0);
  88         if(occ.t1 >= 0 || occ.t5 >= 0) {
  89                 i = PTIME;
  90                 if(mag > 2)
  91                         i |= DARK;
  92                 if(mag < 5)
  93                         i |= SIGNIF;
  94                 if(occ.t1 >= 0 && occ.e1 >= 0)
  95                         event("Occultation of SAO %s begins at ",
  96                                 saop, "", occ.t1, i);
  97                 if(occ.t5 >= 0 && occ.e5 >= 0)
  98                         event("Occultation of SAO %s ends at ",
  99                                 saop, "", occ.t5, i);
 100                 while(*saop++)
 101                         ;
 102         }
 103         goto loop;
 104 }