home *** CD-ROM | disk | FTP | other *** search
/ Shareware Overload / ShartewareOverload.cdr / games / viewmoon.zip / MOONBEAM.V1 < prev    next >
Text File  |  1985-03-14  |  21KB  |  570 lines
  1.  
  2. program moon;
  3.         var {global variables}
  4.           time,mm,dd,yy,jdn,jd,t,theta,lambda,beta,mmanom,msanom : real;
  5.           mra,mdec : real;
  6.           continue : boolean;
  7.           answer : char;
  8. procedure start;
  9.   var
  10.     i : integer;
  11.     stz : real;
  12. begin
  13.   clrscr;
  14.   writeln;
  15.   writeln;
  16.   writeln('               PROGRAM MOONBEAM ');
  17.   writeln('                   version 1.0');
  18.   writeln('               Copyright 1985 by');
  19.   writeln('               F. T. Mendenhall Jr.');
  20.   for i :=1 to 6 do
  21.      writeln;
  22.   writeln('Please supply the date of interest:');
  23.   write('What is the month (1-12)? ');
  24.   readln(mm);
  25.   write('What is the day (1-31)?');
  26.   readln(dd);
  27.   write('What is the year?');
  28.   readln(yy);
  29.   clrscr;
  30.   writeln('     PLEASE ENTER APPROPRIATE TIME ZONE');
  31.   writeln;
  32.   writeln(' Enter 0 for Universal Time (GMT)');
  33.   writeln;
  34.   writeln(' Enter 1 for Eastern Standard Time');
  35.   writeln;
  36.   writeln(' Enter 2 for Central Standard Time');
  37.   writeln;
  38.   writeln(' Enter 3 for Mountain Standard Time');
  39.   writeln;
  40.   writeln(' Enter 4 for Pacific Standard Time');
  41.   writeln;
  42.   write(' What time zone do you wish to use (0-4)? ');
  43.   readln(stz);
  44.   clrscr;
  45.   writeln;
  46.   writeln(' USING THE 24 HOUR CLOCK WHAT IS THE TIME OF OBSERVATION?');
  47.   writeln;
  48.   writeln('NOTE: 1pm=1300   2pm=1400   3pm=1500');
  49.   writeln('      4pm=1600   5pm=1700   6pm=1800');
  50.   writeln('      7pm=1900   8pm=2000   9pm=2100');
  51.   writeln('     10pm=2200  11pm=2300  12am=0000');
  52.   writeln;
  53.   writeln(' for example 8:45pm would entered as 2045 hours');
  54.   writeln;
  55.   write('  What is the time of your observation? ');
  56.   readln(time);
  57.   clrscr;
  58.   writeln;
  59.   write('For the date ',mm:2:0,'/',dd:2:0,'/',yy:4:0);
  60.   writeln(' At ',time:4:0,' Local Time' );
  61.   if stz <> 0 then
  62.   time := time + stz*100 +400;
  63. end;{end of procedure start}
  64. procedure julian;
  65. var
  66.   m,y,a,b : real;
  67.   begin { begin julian day procedure}
  68.     dd := dd + Int(time/100.0)/24.0 + Frac(time/100.0)/60.0;
  69.     if mm > 2 then
  70.      begin
  71.        y := yy;
  72.        m := mm;
  73.        end;{end of then clause}
  74.      if (mm=1) or (mm=2) then
  75.         begin
  76.          y := yy-1;
  77.          m := mm +12;
  78.          end; {end of then clause}
  79.      if yy < 1583 then
  80.         writeln('WARNING PROGRAM NOT VALID PRIOR TO 1583');
  81.      a := int(y/100);
  82.      b := 2-a+int(a/4);
  83.      jd := int(365.25*y)+int(30.6001*(m+1))+dd+1720994.5+b;
  84.  end; {end of procedure julian}
  85. procedure julian_new;
  86. { This procedure calculates the previous date of the new moon
  87.  from the date of interest}
  88.   var
  89.    jds,ms,ds,k,y,tt,m,m1,f,cor : real;
  90. begin { begin procedure julian_new}
  91.    ms := mm;
  92.        { store the original date and julian day}
  93.    ds := dd;
  94.    mm := 1;
  95.    dd := 1;
  96.    jds := jd;
  97.         { find the julian day of the first of the year}
  98.    julian;
  99.    tt := jds -jd;
  100.    { tt is the number of days since the first of the year}
  101.    y := yy+tt/365.;
  102.    jd := jds;
  103.     {restore the original date and jd}
  104.    mm := ms;
  105.    dd := ds;
  106.    k := (y-1900)*12.3685;
  107.    k := int(k);
  108.    tt := k/1236.85;
  109.    jdn := 2415020.75933+29.53058868*k+0.0001178*tt*tt
  110.           -0.000000155*tt*tt*tt+0.00033*sin((166.56+132.87*tt
  111.           -0.009173*tt*tt)*pi/180.0);
  112.    m := (359.2242+29.10535608*k-0.0000333*tt*tt
  113.         -0.00000347*tt*tt*tt)*pi/180.0;
  114.    m1 := (306.0256+385.81691806*k+0.0107306*tt*tt
  115.           +0.00001236*tt*tt*tt)*pi/180.0;
  116.    f := (21.2964+390.67050646*k-0.0016528*tt*tt
  117.          -0.00000239*tt*tt*tt)*pi/180.0;
  118.    cor := (0.1734-0.000393*tt)*sin(m)
  119.           +0.0021*sin(2*m)
  120.           -0.4068*sin(m1)
  121.           +0.0161*sin(2*m1)
  122.           -0.0004*sin(3*m1)
  123.           +0.0104*sin(2*f)
  124.           -0.0051*sin(m+m1)
  125.           -0.0074*sin(m-m1)
  126.           +0.0004*sin(2*f+m)
  127.           -0.0004*sin(2*f-m)
  128.           -0.0006*sin(2*f+m1)
  129.           +0.0010*sin(2*f-m1)
  130.           +0.0005*sin(m+2*m1);
  131.    jdn := jdn +cor;
  132. end;
  133. { jdn is the julian of the previous new moon }
  134. procedure phaseout;
  135. { This procedure prints out the moon phase }
  136.   var
  137.   dpn : real;
  138. begin
  139.   writeln;
  140.   dpn := jd-jdn;
  141.   writeln(' The Moon will be ',dpn:3:1,' days passed new');
  142.   writeln;
  143.   if dpn < 7.5 then
  144.     writeln(' The Moon is between New and First Quarter');
  145.   if (dpn > 7.5) and (dpn <14) then
  146.     writeln (' The Moon is between First Quarter and Full');
  147.   If (dpn >14) and (dpn < 21.5) then
  148.      writeln(' The Moon is between Full and Last Quarter');
  149.   if dpn > 21.5 then
  150.      writeln(' The Moon is between Last Quarter and New');
  151.   writeln;
  152.   writeln('STANDBY: Lunar position being calculated');
  153. end;{end of phaseout}
  154. procedure sunlat;
  155.           {This procedure finds the mean sun anomaly and the suns true
  156.           longitude theta}
  157. var
  158.   l,c : real;
  159. begin {sunlat}
  160.   t := (jd-2415020.0)/36525.0;
  161.   l := 279.69668 + 36000.76892*t + 0.0003025*t*t;
  162.   msanom := 358.475833 + 35999.04975*t - 0.000150*t*t - 0.0000033*t*t*t;
  163.   c := +(1.919460-0.004789*t-0.000014*t*t)*sin(msanom*pi/180.0)
  164.        +(0.020094 - 0.000100*t)*sin(2.0*msanom*pi/180.0)
  165.        + 0.000293*sin(3.0*msanom*pi/180.0);
  166.   theta := l + c;
  167. end; {end procedure sunlat}
  168. procedure moonpos;
  169.           {This procedure calculates the moons mean anomaly, and its
  170.            longitude, lambda, and latitude, beta}
  171. var
  172.   l,d,f,omeg,epsio,addcon,addcon2,rad : real;
  173. begin {moonpos, note that t,msanom are previously calculated in sunlat}
  174.   rad := pi/180.0;
  175.   epsio := 1.0 - 0.002495*t - 0.00000752*t*t;
  176.   addcon := 0.003964*sin((346.560+132.870*t-0.0091731*t*t)*rad);
  177.   addcon2 := sin((51.2+20.2*t)*rad);
  178.   l := 270.434164 + 481267.8831*t -0.001133*t*t + 0.0000019*t*t*t;
  179.   mmanom := 296.104608 + 477198.8491*t +0.009192*t*t +0.0000144*t*t*t;
  180.   d := 350.737486 +445267.1142*t -0.001436*t*t +0.0000019*t*t*t;
  181.   f := 11.250889 +483202.0251*t -0.003211*t*t -0.0000003*t*t*t;
  182.   omeg := 259.183275 -1934.1420*t +0.002078*t*t +0.0000022*t*t*t;
  183.   { now adding periodic variations}
  184.   l := l +addcon +0.001964*sin(omeg*rad) +0.000233*addcon2;
  185.   msanom := msanom -0.001778*addcon2;
  186.   mmanom := mmanom +addcon +0.000817*addcon2 +0.002541*sin(omeg*rad);
  187.   d := d +addcon +0.002011*addcon2 +0.001964*sin(omeg*rad);
  188.   f := f+addcon -0.024691*sin(omeg*rad)
  189.        -0.004328*sin((omeg+275.05-2.3*t)*rad);
  190.   {calculate the coordinates}
  191.   lambda := l + 6.288750*sin(mmanom*rad)
  192.               +1.274018*sin((2*d-mmanom)*rad)
  193.               +0.658309*sin(2*d*rad)
  194.               +0.213616*sin(2*mmanom*rad)
  195.               -0.185596*epsio*sin(msanom*rad)
  196.               -0.114336*sin(2*f*rad)
  197.               +0.058793*sin((2*d-2*mmanom)*rad)
  198.               +0.057212*epsio*sin((2*d-msanom-mmanom)*rad)
  199.               +0.053320*sin((2*d+mmanom)*rad)
  200.               +0.045874*epsio*sin((2*d-msanom)*rad)
  201.               +0.041024*epsio*sin((mmanom-msanom)*rad)
  202.               -0.034718*sin(d*rad)
  203.               -0.030465*epsio*sin((mmanom+msanom)*rad)
  204.               +0.015326*sin((2*d-2*f)*rad)
  205.               -0.012528*sin((2*f+mmanom)*rad)
  206.               -0.010980*sin((2*f-mmanom)*rad)
  207.               +0.010674*sin((4*d-mmanom)*rad)
  208.               +0.010034*sin(3*mmanom*rad)
  209.               +0.008548*sin((4*d-2*mmanom)*rad)
  210.               -0.007910*epsio*sin((msanom-mmanom+2*d)*rad)
  211.               -0.006783*epsio*sin((2*d+msanom)*rad)
  212.               +0.005162*sin((mmanom-d)*rad)
  213.               +0.005000*epsio*sin((msanom+d)*rad)
  214.               +0.004049*epsio*sin((mmanom-msanom+2*d)*rad)
  215.               +0.003996*sin((2*mmanom+2*d)*rad)
  216.               +0.003862*sin(4*d*rad)
  217.               +0.003665*sin((2*d-3*mmanom)*rad)
  218.               +0.002695*epsio*sin((2*mmanom-msanom)*rad)
  219.               +0.002602*sin((mmanom-2*f-2*d)*rad)
  220.               +0.002396*epsio*sin((2*d-msanom-2*mmanom)*rad)
  221.               -0.002349*sin((mmanom+d)*rad)
  222.               +0.002249*epsio*epsio*sin((2*d-2*msanom)*rad)
  223.               -0.002125*epsio*sin((2*mmanom+msanom)*rad)
  224.               -0.002079*epsio*epsio*sin(2*msanom*rad)
  225.               +0.002059*epsio*epsio*sin((2*d-mmanom-2*msanom)*rad)
  226.               -0.001773*sin((mmanom+2*d-2*f)*rad)
  227.               -0.001595*sin((2*f+2*d)*rad)
  228.               +0.001220*epsio*sin((4*d-msanom-mmanom)*rad)
  229.               -0.001110*sin((2*mmanom+2*f)*rad)
  230.               +0.000892*sin((mmanom-3*d)*rad)
  231.               -0.000811*epsio*sin((msanom+mmanom+2*d)*rad)
  232.               +0.000761*epsio*sin((4*d-msanom-2*mmanom)*rad)
  233.               +0.000717*epsio*epsio*sin((mmanom-2*msanom)*rad)
  234.               +0.000704*epsio*epsio*sin((mmanom-2*msanom-2*d)*rad)
  235.               +0.000693*epsio*sin((msanom-2*mmanom+2*d)*rad)
  236.               +0.000598*epsio*sin((2*d-msanom-2*f)*rad)
  237.               +0.000550*sin((mmanom+4*d)*rad)
  238.               +0.000538*sin(4*mmanom*rad)
  239.               +0.000521*epsio*sin((4*d-msanom)*rad)
  240.               +0.000486*sin((2*mmanom-d)*rad);
  241.               while lambda > 360.0 do
  242.                lambda := lambda -360.0;
  243.                while lambda < 0.0 do
  244.                 lambda := lambda + 360.0;
  245.           { now for Beta }
  246.  beta := +5.128189*sin(f*rad)
  247.          +0.280606*sin((mmanom+f)*rad)
  248.          +0.277693*sin((mmanom-f)*rad)
  249.          +0.173238*sin((2*d-f)*rad)
  250.          +0.055413*sin((2*d+f-mmanom)*rad)
  251.          +0.046272*sin((2*d-f-mmanom)*rad)
  252.          +0.032573*sin((2*d+f)*rad)
  253.          +0.017198*sin((2*mmanom+f)*rad)
  254.          +0.009267*sin((2*d+mmanom-f)*rad)
  255.          +0.008823*sin((2*mmanom-f)*rad)
  256.          +0.008247*epsio*sin((2*d-msanom-f)*rad)
  257.          +0.004323*sin((2*d-f-2*mmanom)*rad)
  258.          +0.004200*sin((2*d+f+mmanom)*rad)
  259.          +0.003372*epsio*sin((f-msanom-2*d)*rad)
  260.          +0.002472*epsio*sin((2*d+f-msanom-mmanom)*rad)
  261.          +0.002222*epsio*sin((2*d+f-msanom)*rad)
  262.          +0.002072*epsio*sin((2*d-f-msanom-mmanom)*rad)
  263.          +0.001877*epsio*sin((f-msanom+mmanom)*rad)
  264.          +0.001828*sin((4*d-f-mmanom)*rad)
  265.          -0.001803*epsio*sin((f+msanom)*rad)
  266.          -0.001750*sin(3*f*rad)
  267.          +0.001570*epsio*sin((mmanom-msanom-f)*rad)
  268.          -0.001487*sin((f+d)*rad)
  269.          -0.001481*epsio*sin((f+msanom+mmanom)*rad)
  270.          +0.001417*epsio*sin((f-msanom-mmanom)*rad)
  271.          +0.001350*epsio*sin((f-msanom)*rad)
  272.          +0.001330*sin((f-d)*rad)
  273.          +0.001106*sin((f+3*mmanom)*rad)
  274.          +0.001020*sin((4*d-f)*rad)
  275.          +0.000833*sin((f+4*d-mmanom)*rad)
  276.          +0.000781*sin((mmanom-3*f)*rad)
  277.          +0.000670*sin((f+4*d-2*mmanom)*rad)
  278.          +0.000606*sin((2*d-3*f)*rad)
  279.          +0.000597*sin((2*d+2*mmanom-f)*rad)
  280.          +0.000492*epsio*sin((2*d+mmanom-msanom-f)*rad)
  281.          +0.000450*sin((2*mmanom-f-2*d)*rad)
  282.          +0.000439*sin((3*mmanom-f)*rad)
  283.          +0.000423*sin((f+2*d+2*mmanom)*rad)
  284.          +0.000422*sin((2*d-f-3*mmanom)*rad)
  285.          -0.000367*epsio*sin((msanom+f+2*d-mmanom)*rad)
  286.          -0.000353*epsio*sin((msanom+f+2*d)*rad)
  287.          +0.000331*sin((f+4*d)*rad)
  288.          +0.000317*epsio*sin((2*d+f-msanom+mmanom)*rad)
  289.          +0.000306*epsio*epsio*sin((2*d-2*msanom-f)*rad)
  290.                   -0.000283*sin((mmanom+3*f)*rad);
  291.  beta := beta * (1.0 - 0.0004664*cos(omeg*rad)
  292.                    - 0.0000754*cos((omeg+275.05-2.30*t)*rad));
  293. end; {end of moonpos}
  294. procedure illum; { this procedure calculates the percent of
  295.                    the moon's disk that is illuminated}
  296. var
  297.  d,i,k,ang,y : real;
  298. begin {start of illum}
  299.   ang := cos((lambda-theta)*pi/180.0)*cos(beta*pi/180.0);
  300.   {   The line solves for the arccos. i.e
  301.           the governing equation is
  302.           cos(d)=ang    }
  303.   d := ArcTan(sqrt(1-ang*ang)/ang);
  304.   d := d*180.0/pi;
  305.   if d < 0. then
  306.      d := d + 180.0;
  307.   i := 180.0-d-0.1468*
  308.        (1-0.0549*sin(mmanom*pi/180.0))/(1-0.0167*sin(msanom*pi/180.0))
  309.        *sin(d*pi/180.0);
  310.   k := (1+cos(i*pi/180.0))/2.0;
  311.   k := k*100.0;
  312.   writeln('The Moon`s Disk will be ',k:3:0,' percent illuminated tonight.');
  313.   writeln;
  314. end; {end of illum}
  315. procedure equator; {This procedure converts the eclipetal lunar
  316.                     coordinates to equatorial coordinates
  317.                     epoch 2000}
  318. var
  319.   a,y,RA,DEC,rad,hr,min,sec,deg : real;
  320. const
  321.     epsi = 23.4392911;
  322. begin {start of equator}
  323.   rad := pi/180.0;
  324.   y := (sin(lambda*rad)*cos(epsi*rad)-(sin(beta*rad)/cos(beta*rad))
  325.          *sin(epsi*rad));
  326.   a := y/cos(lambda*rad);
  327.   RA := ArcTan(a)/rad;
  328.   if (a > 0.0) and (y > 0.0) and (RA < 0.0) then
  329.      RA := RA + 90.0;
  330.  if (a > 0.0) and (y < 0.0) then
  331.  begin
  332.     if RA > 0.0  then
  333.     RA := RA + 180.0 else
  334.     RA := RA + 270.0;
  335.   end;
  336.  if (a < 0.0) and (y > 0.0) then
  337.     begin
  338.       if RA > 0.0 then
  339.         RA := RA + 90.0 else
  340.         RA := RA + 180.0;
  341.      end;
  342.  if (a< 0.0) and (y < 0.0) then
  343.      begin
  344.        if RA > 0.0 then
  345.          RA := RA + 270.0 else
  346.          RA := RA + 360.0;
  347.      end;
  348.  RA := RA/15.0; {convert to hours from degrees}
  349.  mra := RA;
  350.  writeln(' The Moon`s Right Ascension is:');
  351.  hr := Int(RA);
  352.  min := Int(Frac(RA)*60.0);
  353.  sec := Int(Frac(Frac(RA)*60.0)*60.0);
  354.  writeln('    ',hr:2:0,' hours ',min:2:0,' min ',sec:2:0,' seconds');
  355.  writeln;
  356.  writeln(' The Moons Declination is:');
  357.  a := sin(beta*rad)*cos(epsi*rad) + cos(beta*rad)*sin(epsi*rad)*
  358.       sin(lambda*rad);
  359.  {the following equation needs to be solved:
  360.      sin(DEC)=a
  361.   the next line finds the arcsin using the ArcTan function.}
  362.  DEC := ArcTan(a/sqrt(1-a*a))/rad;
  363.  mdec := DEC;
  364.  deg := Int(DEC);
  365.  min := Int(Frac(DEC)*60.0);
  366.  sec := Int(Frac(Frac(DEC)*60.0)*60.0);
  367.  writeln('    ',deg:3:0,' degrees ',min:2:0,' min ',sec:2:0,' seconds');
  368. end; {end of equator}
  369.  
  370. {insert moonplot here}
  371. procedure moonplot ;
  372.  
  373.    TYPE
  374.       str29              = string[29];
  375.       stardata           = record
  376.                               shour     :  byte;
  377.                               smin      :  byte;
  378.                               ssec      :  byte;
  379.                               ssign     :  char;
  380.                               sdeg      :  byte;
  381.                               sdmin     :  byte;
  382.                               smsign    :  char;
  383.                               smag      :  integer;
  384.                            end;
  385.  
  386.    var
  387.       hour       :  byte;
  388.       min        :  byte;
  389.       sec        :  byte;
  390.       sign       :  string[1];
  391.       deg        :  byte;
  392.       dmin       :  byte;
  393.       msign      :  string[1];
  394.       mag        :  integer;
  395.       st         :  integer;
  396.       Irec       :  integer;
  397.       J          :  integer;
  398.       num        :  integer;
  399.       line       :  string[29];
  400.       str2       :  string[2];
  401.       str3       :  string[3];
  402.       datafile   :  text;
  403.       filename   :  string[14];
  404.       recfile    :  string[14];
  405.       starfile   :  file of stardata;
  406.       starrec    :  stardata;
  407.       RA         :  real;
  408.       DEC        :  real;
  409.       Xmag       :  real;
  410.       xpos       :  integer;
  411.       ypos       :  integer;
  412.       theta1     :  real;
  413.       mcx        :  integer;
  414.       mcy        :  integer;
  415.       radius     :  real;
  416.       minmag     :  real;
  417.       dummy      :  char;
  418.  
  419. begin
  420.    recfile := 'star.rec';
  421.    assign(starfile,recfile);
  422.    reset(starfile);
  423.    Irec := filesize(starfile);
  424.    clrscr;
  425.    writeln (' The following plots the moon`s position on a star chart.');
  426.    writeln (' R.A. is from 0 hours at the right of the screen to');
  427.    writeln (' 24 hours at the left of the screen.');
  428.    writeln (' DEC is from -50 deg at the bottom to +50 deg at the top.');
  429.    writeln (' The star data is from the Yale Bright Star data base and');
  430.    writeln (' contains stars down to 7th mag.');
  431.    writeln;
  432.    write('what is the min magnitude to be plotted (1-7)?');
  433.    readln(minmag);
  434.    minmag := minmag*100;
  435.    clrscr;
  436.    Hires;
  437.    hirescolor(15);
  438.          for J := 0 to Irec-1 do
  439.             begin
  440.                seek(starfile,J);
  441.                read(starfile,starrec);
  442.                with starrec do
  443.                begin
  444.                if (sdeg  <= 50) and (minmag > smag) then
  445.                   begin
  446.                   RA := shour + smin/60 +ssec/3600;
  447.                   DEC := sdeg + sdmin/60;
  448.                   Xmag := smag / 100.0;
  449.                   xpos := round(630 -(RA/24)*630);
  450.                   if ssign = '-' then DEC := -1*DEC;
  451.                   ypos := round((50-DEC)/50*90);
  452.                   plot(xpos,ypos,1);
  453.                   if xmag < 4 then
  454.                      begin
  455.                      xpos := xpos +1;
  456.                      plot(xpos,ypos,1);
  457.                      xpos := xpos -2;
  458.                      plot(xpos,ypos,1);
  459.                      xpos :=xpos +1;
  460.                      ypos :=ypos -1;
  461.                      plot(xpos,ypos,1);
  462.                      ypos := ypos +2;
  463.                      plot(xpos,ypos,1);
  464.                      end; {of xmag if}
  465.                   if xmag < 2 then
  466.                    begin
  467.                      xpos := xpos +1;
  468.                      plot(xpos,ypos,1);
  469.                      xpos := xpos -2;
  470.                      plot(xpos,ypos,1);
  471.                      ypos := ypos -2;
  472.                      plot (xpos,ypos,1);
  473.                      xpos := xpos +2;
  474.                      plot(xpos,ypos,1);
  475.                      end; { of mag 2 if}
  476.                   if (xmag <1) or (smsign = '-') then
  477.                      begin
  478.                         xpos := xpos -1;
  479.                         ypos := ypos -1;
  480.                         plot(xpos,ypos,1);
  481.                         ypos := ypos -1;
  482.                         plot(xpos,ypos,1);
  483.                         ypos := ypos +5;
  484.                         plot(xpos,ypos,1);
  485.                         ypos := ypos +1;
  486.                         plot(xpos,ypos,1);
  487.                         ypos := ypos -3;
  488.                         xpos := xpos - 3;
  489.                         plot(xpos,ypos,1);
  490.                         xpos := xpos +1 ;
  491.                         plot(xpos,ypos,1);
  492.                         xpos := xpos + 4;
  493.                         plot(xpos,ypos,1);
  494.                         xpos := xpos +1;
  495.                         plot(xpos,ypos,1);
  496.                       end; {of 1st mag if};
  497.                   end; {of if}
  498.                 end; {of with starrec}
  499.             end;{of for}
  500.    close(starfile);
  501. {start the moonplot code}
  502.        mcx := round(630-(mra/24)*630);
  503.        mcy := round((50-mdec)/50 *90);
  504.        theta1 := 0.;
  505.        radius := 8.0;
  506.        while radius > 1.0 do
  507.         begin
  508.            while theta1 < 2*pi do
  509.              begin
  510.              xpos := mcx+round(radius*cos(theta1));
  511.              ypos := mcy+round(radius/2.0*sin(theta1));
  512.              plot(xpos,ypos,1);
  513.              theta1 := theta1 + pi/16.0;
  514.            end;
  515.          radius := radius - 1.0;
  516.          theta1 := 0.;
  517.          plot(mcx,mcy,1);
  518.        end;
  519.        gotoxy(1,24);
  520.        write('input g to continue -');
  521.        readln(dummy);
  522.        TextMode(bw80);
  523. end;  { of procedure moonplot }
  524.  
  525. begin { start of moon}
  526. continue := true;
  527. while continue do
  528.  begin
  529.  start;
  530.  julian;
  531.  writeln ('julian date of the observation is ',jd:10:3);
  532.  julian_new;
  533.  phaseout;
  534.  sunlat;
  535.  moonpos;
  536.  writeln;
  537.  illum;
  538.  equator;
  539.  writeln;
  540.  write(' Do you want a plot of position on star chart (y/n)?');
  541.  readln(answer);
  542.  if (answer) = 'y' then moonplot;
  543.  write('Do you want to check another date (y/n)? ');
  544.  readln(answer);
  545.  if (answer) <> 'y' then continue := false;
  546.  end; {of while}
  547.  clrscr;
  548.  writeln;
  549.  writeln;
  550.  writeln('Astronomy is a beautiful adventure! I hope this program helps you');
  551.  writeln('obtain more pleasure from your investigation of the night skies.');
  552.  writeln('This program is Copyrighted. However, the author grants everyone  ');
  553.  writeln('the right to copy and distribute this program as long as they');
  554.  writeln('do not charge money or services in exchange for the program.');
  555.  writeln;
  556.  writeln('If you find this program useful and wish to encourage the ');
  557.  writeln('development of other free programs, the author requests a ');
  558.  writeln('donation of $5.00 be sent to:');
  559.  writeln;
  560.  writeln('              Fred Mendenhall');
  561.  writeln('              2209 Tam-O-Shanter Ct.');
  562.  writeln('              Carmel, Ind. 46032');
  563.  writeln;
  564.  writeln('Copies of the Pascal source file are available by');
  565.  writeln('sending a SASE and blank disk with  your donation.');
  566.  writeln;
  567.  writeln('                               Happy Observing!');
  568.  
  569. end. {end of moon}
  570.