home *** CD-ROM | disk | FTP | other *** search
/ ftp.barnyard.co.uk / 2015.02.ftp.barnyard.co.uk.tar / ftp.barnyard.co.uk / cpm / walnut-creek-CDROM / ENTERPRS / CPM / UTILS / F / F80.LBR / PLANET.FOR < prev    next >
Text File  |  2000-06-30  |  17KB  |  507 lines
  1.  
  2.  
  3. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  4. C
  5. C Program PLANETF
  6. C
  7. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  8.         PROGRAM XPLANET
  9.         BYTE OPTION, ZOD1(12),CALEND,
  10.      1  RETRO(8),BLAH,GLYPH(21,14)
  11.         LOGICAL HLP
  12.         INTEGER ID(21,14),IM(21,14),IS(21,14),IH(4,3)
  13.         INTEGER*4 JD1
  14.         DOUBLE PRECISION T,AN1(9),XJD,CONV,GMT,XL,YL,XP,YP,XN,
  15.      1  DAY,PI2,WEEK(7),WEEKR(7),TP,
  16.      2  YBES,TT
  17.         DIMENSION EC0(9),EC1(9),AN0(9),P0(9),P1(9),EC(9),
  18.      1  TL(9),TB(9),TH0(9),TH1(9),XI(9),R(9),A(9),G(9),GE(9),
  19.      2  X(9),Y(9),Z(9),XX(9),YY(9),ZZ(9),H(4),
  20.      3  P(21,14),
  21.      8  MONTH(12),MOND(12),XNODE(2,10),YNODE(2,10)
  22. CCCCC8CCCC
  23. C
  24. CCCCCCCCC
  25.         DATA MOND/31,0,31,30,31,30,31,31,30,31,30,31/
  26.         DATA MONTH/'JA','FE','MR','AP','MY','JN','JL','AU',
  27.      1  'SE','OC','NO','DE'/
  28.         DATA WEEK/'Monday  ','Tuesday ','Wednsday',
  29.      1  'Thursday','Friday  ','Saturday','Sunday  '/
  30.         DATA WEEKR/'Moon)   ','Mars)   ','Mercury)','Jupiter)',
  31.      1  'Venus)  ','Saturn) ','Sun)    '/
  32.         DATA ZOD1/'A','T','G','C','L','V','=','S','/','K','Q','P'/
  33. CCCCCCCCC
  34. C
  35. C     PLUTO ELEMENTS FROM SHARAF (1964)
  36. C
  37. CCCCCCCCC
  38.         DATA A/1.00000023,.387098599,.723331619,1.523688395,
  39.      1  5.202802875,9.53884320,19.19097811,30.0706724,39.672599/
  40.         DATA AN0/358.475833,102.279381,212.603222,319.529425,
  41.      1  225.444651,175.758444,74.313628,41.269550,231.002308/
  42.         DATA AN1/35999.049750D0,149472.515289D0,58517.803875D0,
  43.      1  19139.858500D0,3034.906654D0,1222.116782D0,428.502578D0,
  44.      2  218.466783D0,144.072477D0/
  45.         DATA EC0/.01675104,.20561421,.00682069,.09331290,
  46.      1  .04825382,.05606075,.04704433,.00853341,.24706226/
  47.         DATA EC1/-.00004180,.00002046,-.00004774,.00009206,
  48.      1  .0,.0,.0,.0,.0/
  49.         DATA P0/101.220833,75.899697,130.163833,334.218203,
  50.      1  11.907422,90.110981,169.048778,43.755611,221.592475/
  51.         DATA P1/1.719175,1.555489,1.408036,1.840758,
  52.      1  .0,.0,.0,.0,1.388888/
  53.         DATA TH0/0.,47.145944,75.779647,48.786442,
  54.      1  98.932822,112.347606,73.490250,130.678889,108.937165/
  55.         DATA TH1/0.,1.185208,.899850,.770992,
  56.      1  .0,.0,.509667,1.100972,1.358056/
  57.         DATA XI/0.,7.002881,3.393630,1.850333,
  58.      1  1.311614,2.494239,.7726658,1.779256,17.109816/
  59. CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
  60. C
  61. C Code
  62. C
  63. CCCCCCCCC
  64.         ASIN(DUM)=ATAN(DUM/SQRT(-DUM*DUM+1))
  65.         TAN(DUMMY)=SIN(DUMMY)/COS(DUMMY)
  66.         ARG(P,Q)=P-Q*INT(P/Q)-Q*INT(SIGN(.5,P-Q*INT(P/Q))-.5)
  67.         HLP=.FALSE.
  68.         LINEP=0
  69.         PI=ATAN(1.0)
  70.         PI2=6.283185307179586476D0
  71.         CONV=PI2/360.D0
  72.         CONVS=CONV
  73.         DO 130 I=1,21
  74.         DO 130 J=1,14
  75. 130     P(I,J)=0.
  76.         WRITE(5,10000)
  77. 10000   FORMAT('$Sample input (Gregorian date,time',
  78.      1  'conversion to GMT,latitude,longitude,OPTION)'//
  79.      2  ' 1899/12/31 23:59 +08:00 34N59 122W59) : ')
  80. CCCCCCCCC
  81. C
  82. C     HOROSCOPE LOOP STARTS HERE
  83. C
  84. CCCCCCCCC
  85.         READ(5,11000)NY,CALEND,NM,ND,NH,NMI,NCH,NCM,LATD,LATS ,
  86.      1  LATM,LONGD,LONGS,LONGM,OPTION
  87. 11000   FORMAT(I4,A1,4(I2,1X),I3,2(1X,I2),A1,I2,1X,I3,A1,I2,1X,A1)
  88.         WRITE(5,11005)
  89. 11005   FORMAT('$Do you care for a hardcopy ? ')
  90.         READ(5, 11006) BLAH
  91. 11006   FORMAT(A1)
  92. C       IF (BLAH.NE.'Y'.AND.BLAH.NE.'y') GO TO 11035
  93. C       CALL OPEN(7,'B:PLANET.DAT')
  94. C       IOUT=7
  95. C       GO TO 11008
  96. 11035   IOUT=5
  97. 11008   LINEP=1
  98.         WRITE(IOUT,11500)NY,CALEND,NM,ND,NH,NMI,NCH,NCM,LATD,
  99.      1  LATS,LATM,LONGD,LONGS,LONGM,OPTION
  100. 11500   FORMAT('1'//' ---Date--- -Time -Zone- -Lat- -Long- ',
  101.      1  'Comments--------'/' ',I4,A1,I2,'/',I2,I3,':',I2,
  102.      1  I4,':',I2,I3,A1,I2,I4,A1,I2,1X,A1)
  103. 144     GMT=((NH+NCH-12)*60+NMI+SIGN(NCM*1.,NCH*1.))/1440.D0
  104.         XJD=(ND-32075+1461.D0*(NY+4800+(NM-14)/12)/4+367.D0*(NM-2-
  105.      1  (NM-14)/12*12)/12-3.D0*((NY+4900+(NM-14)/12)/100)/4)+GMT
  106.         IF(CALEND.NE.'B')GO TO 150
  107.         NY=1-NY
  108. 150     IF(CALEND.NE.'J'.AND.CALEND.NE.'B')GO TO 160
  109.         XJD=(ND-32075+1461.D0*(NY+4800+(NM-14)/12)/4
  110.      1  +367.D0*(NM-2-(NM-14)/12*12)/12-38)+GMT
  111. 160     T=(XJD-2415020.D0)/36525.D0
  112. C
  113.         JFINAL=2
  114. CCCCCCCCC
  115. C
  116. C     Planetary Ephemerides comxzted below (Newcomb/Hill,1898)
  117. C
  118. CCCCCCC
  119. 170     DO 270 J=1,JFINAL
  120.         TT=T+(J-1)/36525.D0
  121.         TP=TT+18262.D0/36525.D0
  122.         W=TT+0.0
  123.         DO 174 I=1,9
  124.         EC(I)=EC0(I)+EC1(I)*W
  125.         N=AN1(I)*TT/360.D0
  126. 174     G(I)=(AN1(I)*TT-N*360.D0+DBLE(AN0(I)))*CONV
  127.         G5=G(5)
  128.         G6=G(6)
  129.         G7=(220.169542+428.49311*TP)*CONVS
  130.         G(5)=G(5)+(.6506*SIN(2*G6-2*G5+336.9*CONVS)
  131.      1  +(3.9987-.002213*36525./4332.58*TP)
  132.      2  *SIN(5*G6-2*G5+(67.15-8197.0/3600.*TP)*CONVS)
  133.      3  +.5380*SIN(5*G6-3*G5+176.5*CONVS)
  134.      4  +.4112*SIN(2*G6-G5+1.4*CONVS)
  135.      4  +.0399278*36525./4332.58*TP*SIN(-G5+227.46*CONVS)
  136.      5  +.2763*SIN(3*G6-2*G5+127.4*CONVS)
  137.      6  +.2669*SIN(G6-G5+79.2*CONVS)  )*299.12837/3600.*CONVS
  138. CCCCCCCCC
  139. C
  140. CCCCCCCCC
  141.         G(6)=G(6)+(24.153*SIN(5*G6-2*G5+(247.11-2.277*TP)*CONVS)
  142.      1  +5.679*SIN(4*G6-2*G5+277.39*CONVS)
  143.      2  +3.505*SIN(2*G6-G5+181.43*CONVS)
  144.      3  +.657765*36525./10759.20*TP*SIN(G6+238.0*CONVS)
  145.      4  +.278*SIN(3*G6-G5+121.2*CONVS)
  146.      5  +.266*SIN(2*G6-2*G5+157.0*CONVS)
  147.      6  +.238*SIN(6*G6-2*G5-3*G7+6.9*CONVS)
  148.      7  +.223*SIN(10*G6-4*G5+(133.6-14814.5/3600.*W)*CONVS)
  149.      8  +.234*SIN(3*G7-G6+321.7*CONVS)  )*120.455/3600.*CONVS
  150. CCCCCCCCC
  151. C
  152. C
  153. CCCCCCCCC
  154.         G(7)=G(7)+33.086*W*W/3600.*CONVS
  155.         G(8)=G(8)-22.401/3600.*W*W*CONVS
  156.         DO 177 I=1,9
  157.         GE(I)=G(I)+EC(I)*SIN(G(I)+EC(I)*SIN(G(I)+EC(I)*SIN(
  158.      1  G(I)+EC(I)*SIN(G(I)+EC(I)*SIN(G(I))))))
  159. 177     X(I)=2.*ATAN(SQRT((1.+EC(I))/(1.-EC(I)))*TAN(.5*GE(I)))
  160.      1  +(P0(I)+P1(I)*W)*CONVS + (I/7-I/9)*W*5025.3/3600.*CONVS
  161.         GT=(358.415+35998.928*W)*CONVS
  162.         GJ=(225.209+3034.462*W+.332*SIN((134.4+38.5*W)*CONVS))
  163.      1  *CONVS
  164.         X(4)=X(4)+(25.384*COS(GJ-G(4)-48.9*CONVS)
  165.      1  +52.490*SIN((47.48+19.771*W)*CONVS)-37.05-13.50*W
  166.      2  +21.869*COS(2*GJ-G(4)-188.3*CONVS)
  167.      3  +16.035*COS(2*GJ-2*G(4)-191.9*CONVS)
  168.      4  +13.966*COS(-GT+2*G(4)-20.5*CONVS)
  169.      5  +8.559*COS(-GT+G(4)-35.1*CONVS)  )/3600.*CONVS
  170. CCCCCCCCC
  171. C
  172. C
  173. CCCCCCCCC
  174.         GN=(225.417+3034.904*W)*CONVS
  175.         GG=(175.753+1222.113*W)*CONVS
  176.         G1=(74.412+428.498*W)*CONVS
  177.         GP=(74.320+428.498*W)*CONVS
  178.         GPP=(41.339+218.467*W)*CONVS
  179.         X(7)=X(7)+(142.938*SIN(GG-2*G1)+19.508*COS(GG-2*G1)
  180.      1  +75.70*COS(3*G1-GG)-102.30*SIN(3*G1-GG)
  181.      2  -48.623*SIN(GN-G1)-21.320*COS(GN-G1)
  182.      3  -27.871*COS(GP-GPP)+19.869*SIN(GP-GPP)
  183.      4  +28.793*COS(2*GP-2*GPP)+10.035*SIN(2*GP-2*GPP)
  184.      5  +(18.37*SIN(3*GP-3*GPP)+8.91*COS(3*GP-3*GPP))*COS(G(7))
  185.      6  +(8.35*SIN(3*GP-3*GPP)-16.44*COS(3*GP-3*GPP))*SIN(G(7))
  186.      7  -18.585*COS(GG-G1)+12.603*SIN(GG-G1)
  187.      8  +4.327*COS(3*GP-3*GPP)+14.280*SIN(3*GP-3*GPP)
  188.      9  )/3600.*CONVS
  189. CCCCCCCCC
  190. C
  191. C
  192. CCCCCCCCC
  193.         X(7)=X(7)+((112.317*W-1.551*W*W-.516*W*W*W)*SIN(G(7))
  194.      1  -(68.339*W+9.721*W*W)*COS(G(7))+6.605*W*SIN(2*G(7))
  195.      2  -(29.44-.410*W)*SIN((20.45-22.61*W)*CONVS) )/3600.*CONVS
  196. CCCCCCCCC
  197. C
  198. C     JUP/SAT TERM & QUADRATURES OF LONG PERIOD UR`WNEP TERMS
  199. C
  200. CCCCCCCCC
  201.         X(8)=X(8)+(  33.972*W*COS(G(8))
  202.      1  +18.553*SIN((180.966+1004.034*W+.1403*W*W)*CONVS)
  203.      2  +34.138*SIN((153.267+2816.296*W-.0573*W*W)*CONVS)
  204.      3  +(26.50*W+3.92*W*W)*SIN(G(8))   )/3600.*CONVS
  205. CCCCCCCCC
  206. C
  207. C
  208. CCCCCCCCC
  209.         DO 185 I=1,9
  210.         TH=(TH0(I)+TH1(I)*W)*CONVS
  211.         DO 180 K=1,2
  212. CCCCCCCCC
  213. C
  214. C     MOTION OF JUP/SAT NODES FROM LEVERRIER
  215. C
  216. CCCCCCCCC
  217.         TH2=TH+(I/5*3636.6-I/6*493.1-I/7*3143.5)/3600.*TP*CONVS
  218.         R(I)=A(I)*(1.-EC(I)*COS(2.*ATAN(TAN((TH2+(K-1)*PI-(P0(I)
  219.      1  +P1(I)*W)*CONVS)*.5)*SQRT((1.-EC(I))/(1.+EC(I))))))
  220.         XNODE(K,I)=R(I)*COS(TH2+(K-1)*PI)
  221. 180     YNODE(K,I)=R(I)*SIN(TH2+(K-1)*PI)
  222.         TL(I)=TH+ATAN2(SIN(X(I)-TH)*COS(XI(I)*CONVS),
  223.      1  COS(X(I)-TH))+(I/5-I/7)*TP*5026.1/3600.*CONVS
  224.         TB(I)=ASIN(SIN(XI(I)*CONVS)*SIN(X(I)-TH))
  225. 185     R(I)=A(I)*(1.-EC(I)*COS(GE(I)))
  226.         TB(5)=TB(5)+4.37431*36525./4332.58*SIN(X(5)+23.62*CONVS)
  227.      1  *TP/3600.*CONVS
  228.         TB(6)=TB(6)+24.266*36525./10759.2*SIN(X(6)-13.05*CONVS)
  229.      1  *TP/3600.*CONVS
  230. CCCCCCCCC
  231. C
  232. C     TB(5) EFFECTS LESS THAN 35"T, TB(6) LESS THAN 84"T
  233. C
  234. CCCCCCCCC
  235.         R(5)=R(5)*10.**(.0002303*COS(2*G6-2*G5+336.9*CONVS)
  236.      1  +.0001679*COS(5*G6-3*G5+176.4*CONVS)
  237.      2  +.0000125634*36525./4332.58*TP*COS(-G5+227.4*CONVS)  )
  238. CCCCCCCCC
  239. C
  240. C       22",16",10"T,(7",5",3")
  241. C
  242. CCCCCCCCC
  243.         R(6)=R(6)*10.**(.0007005*COS(4*G6-2*G5+277.3*CONVS)
  244.      1  +.0003783*COS(G6-G5+79.8*CONVS)
  245.      2  +.000083491*36525./10759.20*TP*COS(G6+58.0*CONVS)
  246.      3  +.0002443*COS(2*G6-G5+176.0*CONVS)  )
  247. CCCCCCCCC
  248. C
  249. C
  250. CCCCCCCCC
  251.         DO 190 I=1,9
  252.         X(I)=R(I)*COS(TB(I))*COS(TL(I))
  253.         Y(I)=R(I)*COS(TB(I))*SIN(TL(I))
  254. 190     Z(I)=R(I)*SIN(TB(I))
  255.         DO 210 I=2,9
  256.         DO 200 K=1,2
  257.         XNODE(K,I)=XNODE(K,I)-X(1)
  258.         YNODE(K,I)=YNODE(K,I)-Y(1)
  259. 200     P(K+17,I)=ARG(ATAN2(YNODE(K,I),XNODE(K,I))/6.2831853,1.)
  260.         XX(I)=X(I)-X(1)
  261.         YY(I)=Y(I)-Y(1)
  262. 210     ZZ(I)=Z(I)-Z(1)
  263.         XX(1)=-X(1)
  264.         YY(1)=-Y(1)
  265.         ZZ(1)=-Z(1)
  266.         XNUT=-17.2327/1296000.*SIN((259.18-1934.142*W)*CONVS)
  267.         DO 230 I=1,9
  268.         P(J+19,I)=ARG(ATAN2(YY(I),XX(I))/6.28318531+XNUT,1.)
  269.         IF(J-1)230,220,230
  270. 220     DIST=SQRT(XX(I)*XX(I)+YY(I)*YY(I)+ZZ(I)*ZZ(I))
  271.         P(1,I)=ASIN(ZZ(I)/DIST)/CONVS
  272.         E=(23.4522944-.0130125*W+.002558*COS((259.18-1934.142
  273.      1  *W)*CONVS))*CONVS
  274.         P(2,I)=ASIN(COS(P(1,I)*CONVS)*SIN(P(20,I)*6.28318531)
  275.      1  *SIN(E)+SIN(P(1,I)*CONVS)*COS(E))/CONVS
  276.         P(6,I)=ARG(TL(I)/6.28318531,1.)
  277.         P(3,I)=TB(I)/CONVS
  278.         P(4,I)=ASIN(COS(TB(I))*SIN(TL(I))*SIN(E)
  279.      1  +SIN(TB(I))*COS(E))/CONVS
  280. 230     CONTINUE
  281. CCCCCCCCC
  282. C
  283. C     LUNAR EPHEMERIDES COMPUTED TO WITHIN 1' (BROWN/I.L.E.)
  284. C
  285. CCCCCCCCC
  286.         DAY=XJD-2415020.D0+(J-1.D0)
  287.         XL=.751206D0+DAY*.0366011014634D0
  288.         YL=.776935D0+DAY*.0027379092649D0
  289.         XP=.928693D0+DAY*.0003094557786D0
  290.         YP=.781169D0+DAY*.0000001307457D0
  291.         XN=.719954D0-DAY*.0001470942283D0
  292.         AL=(XL-IDINT(XL))*PI2
  293.         BL=(YL-IDINT(YL))*PI2
  294.         AP=(XP-IDINT(XP))*PI2
  295.         BP=(YP-IDINT(YP))*PI2
  296.         AN=(XN-IDINT(XN))*PI2
  297.         U=AL-AP
  298.         V=BL-BP
  299.         F=AL-AN
  300.         D=AL-BL
  301.         DL=22639*SIN(U)-4586*SIN(U-D-D)+2370*SIN(D+D)+769*SIN(U+U)
  302.      1  -668*SIN(V)-412*SIN(F+F)-212*SIN(U+U-D-D)-206*SIN(U+V-D-D)
  303.      2  +192*SIN(U+D+D)-165*SIN(V-D-D)+148*SIN(U-V)-125*SIN(D)
  304.      3  -110*SIN(U+V)-55*SIN(F+F-D-D)-45*SIN(U+F+F)+40*SIN(U-F-F)
  305.      4  -38*SIN(U-4*D)+36*SIN(3*U)-31*SIN(U+U-4*D)+28*SIN(U-V-D-D)
  306.      5  -24*SIN(V+D+D)+19*SIN(U-D)+18*SIN(V+D)+15*SIN(U-V+D+D)
  307.      6  +14*SIN(4*D)+14*SIN(U+U+D+D)-13*SIN(3*U-D-D)
  308.         P(J+19,10)=ARG(XL-IDINT(XL)+DL/1296000.+XNUT,1.)
  309.         P(J+19,11)=ARG(AN/6.28318531,1.)
  310.         IF(J-1)270,260,270
  311. 260     P(1,10)=(18461*SIN(F)+1010*SIN(U+F)-1000*SIN(F-U)
  312.      1  -624*SIN(F-2*D)-167*SIN(U+F-D-D)+199*SIN(F-U+D+D)
  313.      2  +117*SIN(F+D+D)+62*SIN(U+U+F)-33*SIN(F-U-D-D)
  314.      3  -32*SIN(F-U-U)-30*SIN(V+F-D-D)-16*SIN(U+U+F-D-D)
  315.      4  +15*SIN(U+F+D+D)+12*SIN(F-V-D-D)+9*SIN(F-U-V+D+D))/3600.
  316.         P(2,10)=ASIN(COS(P(1,10)*CONVS)*SIN(P(20,10)*6.28318531)
  317.      1  *SIN(E)+SIN(P(1,10)*CONVS)*COS(E))/CONVS
  318. 270     CONTINUE
  319.         WRITE(5,271)
  320. 271     FORMAT(' ','EPH.READY'/)
  321. CCCCCCCCC
  322. C
  323. CCCCCCCCC
  324.         E=(23.4522944-.0130125*T+.002558*DCOS((259.18-1934.142
  325.      1  *T)*CONV))*CONVS
  326.         DO 277 I=2,9
  327.         RETRO(I-1)='D'
  328.         IF(ARG(P(21,I)-P(20,I),1.)-.5)277,277,275
  329. 275     RETRO(I-1)='R'
  330. 277     CONTINUE
  331. CCCCCCCCC
  332. C
  333. CCCCCCCCC
  334.         XLW=(LONGD+LONGM/60.)/360.
  335.         XLN=(LATD+LATM/60.)*CONVS
  336.         IF(LATS.EQ.'N')GO TO 280
  337.         XLN=-XLN
  338. 280     IF(LONGS.EQ.'W')GO TO 290
  339.         XLW=-XLW
  340. CCCCCCCCC
  341. C
  342. C     H(1)=GMT  H(2)=GST  H(3)=LMT  H(4)=LST
  343. C
  344. CCCCCCCCC
  345. 290     H(1)=ARG(SNGL(GMT)+.5,1.)
  346.         W=100.00213590D0*T-IDINT(100.D0*T)
  347.         H(2)=ARG(H(1)+.27691940+W,1.)
  348.         H(3)=ARG(H(1)-XLW,1.)
  349.         H(4)=ARG(H(2)-XLW,1.)
  350.         DO 300 I=1,4
  351.         W=24.*H(I)+.00013888888
  352.         IH(I,1)=W+0
  353.         W=60.*(W-IH(I,1))
  354.         IH(I,2)=W+0
  355. 300     IH(I,3)=60.*(W-IH(I,2))
  356.         WRITE(IOUT,12000),((IH(I,J),J=1,3),I=1,4)
  357. 12000   FORMAT(/' GMT=',I2,2(':',I2),'    GST=',I2,2(':',I2),
  358.      1  '    LMT=',I2,2(':',I2),'    LST=',I2,2(':',I2))
  359.         JD1=IDINT(XJD)
  360.         XJD1=XJD-JD1
  361.         I=JD1+INT(XJD1+1.5-XLW)-(JD1-1+INT(XJD1+1.5-XLW))/7*7
  362.         YBES=1900.D0+(XJD-2415020.31351528D0)/365.24219878125D0
  363.      1  -.00000107523D0*T*T
  364. C       IYBES=YBES+0
  365. C       YBES1=YBES-IYBES
  366.         WRITE(IOUT,13000)WEEK(I),WEEKR(I),XJD,YBES
  367. 13000   FORMAT(' Day of the week (from LMT) IS ',A9,
  368.      1  ' (ruled by ',A8/' JD=',D14.7,7X,' Besselian year = ',
  369.      2  D14.7)
  370. CCCCCCCCC
  371. C
  372. C     HOUSE CUSP COMPUTATIONS (ACCURATE TO LAST DIGIT)
  373. C
  374. CCCCCCCCC
  375. 310   ITER=1+LATD/45
  376. CCCCCCCCC
  377. C
  378. C     "ITER" IS THE NUMBER OF PLACIDEAN ITERATIONS
  379. C@@@
  380. C@@@
  381. C
  382. CCCCCCCCC
  383.         RA=6.28318531*H(4)
  384.         DO 340 N=1,6
  385.         AA=6.28318531*(N+2)/12.
  386.         P(7,N)=ARG(ATAN2(SIN(RA)*COS(AA)+COS(RA)*COS(XLN)*
  387.      1  SIN(AA),COS(E)*(COS(RA)*COS(AA)-SIN(RA)*COS(XLN)*
  388.      2  SIN(AA))-SIN(E)*SIN(XLN)*SIN(AA))/6.28318531,1.)
  389.         P(8,N)=ARG(ATAN2(SIN(RA+AA)*COS(E),COS(RA+AA))
  390.      1  /6.28318531,1.)
  391.         P(9,N)=ARG(ATAN2(SIN(RA+AA),COS(E)*COS(RA+AA)-SIN((3.-
  392.      1  ABS(1.-N))/3.*ASIN(TAN(E)*TAN(XLN)))*COS(E))/PI*.5,1.)
  393.         DO 330 J=1,ITER
  394.         X(1)=P(9,N)*6.28318531
  395.         X(2)=X(1)+CONVS*.5
  396. CCCCCCCCC
  397. C
  398. C     INCREMENT "CONVS*.5" IS ARBITRARILY CHOSEN (]W5 DEG)
  399. C
  400. CCCCCCCCC
  401.         DO 320 K=1,2
  402.         AD=ASIN(TAN(XLN)*TAN(ASIN(SIN(X(K))*SIN(E))))
  403.         DRA=ARG(ATAN2(SIN(X(K))*COS(E),COS(X(K)))-RA-(N-1)*PI/6.
  404.      1  ,PI)-PI*.5
  405. 320     Y(K)=DRA-(4.-N)*AD/3.
  406. 330     P(9,N)=P(9,N)-Y(1)/(Y(2)-Y(1))*.5/360.
  407. CCCCCCCCC
  408. C
  409. C     THE ".5" CONSTANT HERE MUST AGREE WITH THE ONE ABOVE
  410. C
  411. CCCCCCCCC
  412.         P(10,N)=ARG(ATAN2(SIN(RA+AA),COS(E)*COS(RA+AA)-SIN(E)
  413.      1  *TAN(XLN)*SIN(AA))/6.28318531,1.)
  414. 340     P(12,N)=ARG(P(7,1)+(N-1)/12.,1.)
  415.         DO 350 N=1,6
  416.         P(11,N)=ARG(P(7,1)+(N-1)/12.+(ARG(P(7,4)-P(7,1),1.)
  417.      1  -.25)*(3.-ABS(3.-ABS(N-7.)))/3.,1.)
  418. 350     P(13,N)=ARG(P(7,4)+(N-4)/12.,1.)
  419.         W=ARG(ATAN2(SQRT((COS(XLN)/SIN(E))**2-SIN(P(7,4)
  420.      1  *6.28318531)**2),SIN(XLN)*SIN(P(7,4)*6.28318531)),PI)/3.
  421.         DO 351 N=1,3
  422. 351     P(14,N)=ARG(ATAN2(COS(RA+(N-1)*W),-SIN(RA+(N-1)*W)
  423.      1  *COS(E)-SIN(E)*TAN(XLN))/6.28318531,1.)
  424.         W=ARG(ATAN2(SIN(P(7,1)*6.28318531)*COS(E),COS(P(7,1)
  425.      1  *6.28318531))-RA,PI)/3.
  426.         DO 352 N=4,6
  427. 352     P(14,N)=ARG(ATAN2(SIN(RA+(N-4)*W),COS(RA+(N-4)*W)*COS(E)
  428.      1  )/6.28318531+.5,1.)
  429.         DO 355 J=7,14
  430.         DO 355 N=1,6
  431. 355     P(J,N+6)=ARG(P(J,N)+.5,1.)
  432. CCCCCCCCC
  433. C
  434. C     PARS FORTUNA & PLANETARY NODES FOLLOW
  435. C
  436. CCCCCCCCC
  437.         P(20,12)=ARG(P(20,10)-P(20,1)+P(11,1),1.)
  438.         DO 360 I=2,9
  439. 360     P(17,I)=ARG(P(20,I)-P(20,1)+P(11,1),1.)
  440. CCCCCCCCC
  441. C
  442. C     ARGUMENTS TRANSFORMED TO SIGNS/DEGREES/MINUTES
  443. C
  444. CCCCCCCCC
  445.         DO 390 I=1,4
  446.         DO 390 J=1,12
  447.         GLYPH(I,J)='N'
  448.         IF(P(I,J))370,380,380
  449. 370     GLYPH(I,J)='S'
  450. 380     DEG=ABS(P(I,J))+.0083333333
  451.         ID(I,J)=DEG+0
  452. 390     IM(I,J)=60*(DEG-ID(I,J))
  453.         IFINAL=JFINAL+19
  454.         DO 400 I=6,IFINAL
  455.         DO 400 J=1,12
  456.         W=12.*P(I,J)+.00027777777
  457.         IS(I,J)=W+0
  458.         ISS=W+1
  459.         GLYPH(I,J)=ZOD1(ISS)
  460.         W=30.*(W-IS(I,J))
  461.         ID(I,J)=W+0
  462. 400     IM(I,J)=60.*(W-ID(I,J))
  463. CCCCCCCCC
  464. C
  465. C     PROGRESSION & EPHEMERIS OUTPUT ROUTINE FOLLOWS
  466. C
  467. CCCCCCCCC
  468.         IF(OPTION.NE.'P')GO TO 420
  469.         WRITE(IOUT,14000)
  470. 14000   FORMAT(/19X,'Progressed (geocentric) longitudes'//
  471.      1  ' YEAR  SUN   MERC VENUS  MARS JUPIT SATUR ',
  472.      2  ' URAN  NEPT PLUTO  MOON  NODE')
  473.         DO 410 I=1,21
  474.         II=I
  475. 410     WRITE(IOUT,15000)II,(ID(I,J),GLYPH(I,J),IM(I,J),J=1,11)
  476. 15000   FORMAT(' ',I6,1X,11(1X,I6,A1,I6))
  477.         GO TO 9999
  478. 420     IF(OPTION.NE.'E')GO TO 460
  479.         WRITE(IOUT,16000)IH(1,1),IH(1,2),NY
  480. 16000   FORMAT(//19X,'Daily ephemeris for ',I2,':',I2,' GMT'//
  481.      1  ' DATE  SUN   MERC VENUS  MARS JUPIT SATUR  URAN',
  482.      2  '  NEPT PLUTO  MOON  NODE'/I6)
  483.         MY=NY+(NM+7)/10
  484.         MOND(2)=1461*(MY+4800)/4-1461*(MY+4799)/4-337
  485.         IF(CALEND.EQ.'J')GO TO 430
  486.         MOND(2)=MOND(2)-(MY+4900)/100*3/4+(MY+4899)/100*3/4
  487. 430     DO 460 I=1,21
  488.         WRITE(IOUT,17000)MONTH(NM),ND,(ID(I,J),GLYPH(I,J),
  489.      1  IM(I,J),J=1,11)
  490. 17000   FORMAT(' ',A2,I2,11(1X,I2,A1,I2))
  491.         ND=ND+1
  492.         IF(ND-MOND(NM))460,460,440
  493. 440     ND=ND-MOND(NM)
  494.         NM=NM+1
  495.         IF(NM-12)460,460,450
  496. 450     NM=1
  497.         NY=NY+1
  498.         WRITE(IOUT,18000)NY
  499. 18000   FORMAT(I6)
  500. 460     CONTINUE
  501.         GO TO 9999
  502. CCCCCCCCC
  503. 9999    IF(LINEP)9997,9997,9996
  504. 9996    WRITE(5,96000)
  505. 96000   FORMAT(' $DONE')
  506. 9997    END
  507.