SUBROUTINE PRINTR(NWHY,CARD) PRIN001 c$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ c note format 1200 changed 17-dec-86 3P2f6.0 --> 2P2f6.0 c$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ c there is an inconsistency in the way POLAR is treated c This gives 4 values in the print-out e.g PRIN166. c The group path is the last column printed (no header) c$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$ implicit real*8 (a-h,o-z) character*8 id,ndate,modx,mody,modz,modrin,koll,modsav character*8 headr1,headr2,head1,head2,unit,units character*1 type C PRINTS OUTPUT AND PUNCHES RAYSETS WHEN REQUESTED PRIN002 DIMENSION G(3,3),G1(3,3),TYPE(3),HEADR1(20),HEADR2(20),UNITS(20), PRIN003 1 HEAD1(20),HEAD2(20),UNIT(20),RPRINT(20),NPR(20) PRIN004 COMMON /CONST/ PI,PIT2,PID2,DEGS,RAD,DUM(3) PRIN005 COMMON /FLG/ NTYP,NEWWR,NEWWP,PENET,LINES,IHOP,HPUNCH PRIN006 COMMON /RIN/ MODRIN(3),COLL,FIELD,SPACE,N2,N2I,PNP(10),POLAR(2), PRIN007 1 LPOLAR(2) PRIN008 COMMON R(20),T /WW/ ID(10),W0,W(400) PRIN009 common / leo / apogee,range EQUIVALENCE (THETA,R(2)),(PHI,R(3)) PRIN010 EQUIVALENCE (EARTHR,W(2)),(XMTRH,W(3)),(TLAT,W(4)),(TLON,W(5)), PRIN011 1 (F,W(6)),(AZ1,W(10)),(BETA,W(14)),(RCVRH,W(20)),(HOP,W(22)), PRIN012 2 (PLAT,W(24)),(PLON,W(25)),(RAYSET,W(72)) PRIN013 LOGICAL SPACE,NEWWR,NEWWP,PENET PRIN014 REAL*8 N2,N2I PRIN015 COMPLEX*16 PNP,polar,lpolar PRIN016 DATA TYPE / 'X','N','O'/ PRIN017 2,HEADR1(7)/' PHAS'/,HEADR2(7)/'E PATH'/,UNITS(7)/' KM' /, PRIN018 3 HEADR1(8)/' ABSO'/,HEADR2(8)/'RPTION'/,UNITS(8)/' DB' /, PRIN019 4 HEADR1(9)/' DOP'/,HEADR2(9)/'PLER '/,UNITS(9)/' C/S' /, PRIN020 5 HEADR1(10)/' PATH' /,HEADR2(10)/'LENGTH'/,UNITS(10)/' KM' / PRIN023 c data cos_ell_save / 1.d0/ c magnetic dip angle - r(2) = colatitude. Dipole field tan_dip = 2.d0 / dtan(r(2)) dipole_dip = datan(tan_dip) c CALL RINDEX PRIN024 IF (.NOT.NEWWP) GO TO 10 PRIN025 C********* NEW W ARRAY -- REINITIALIZE PRIN026 NEWWP=.FALSE. PRIN027 SPL=dSIN(PLON-TLON) PRIN028 CPL=dSIN(PID2-(PLON-TLON)) PRIN029 SP=dSIN(PLAT) PRIN030 CP=dSIN(PID2-PLAT) PRIN031 SL=dSIN(TLAT) PRIN032 CL=dSIN(PID2-TLAT) PRIN033 C********* MATRIX TO ROTATE COORDINATES PRIN034 G(1,1)=CPL*SP*CL-CP*SL PRIN035 G(1,2)=SPL*SP PRIN036 G(1,3)=-SL*SP*CPL-CL*CP PRIN037 G(2,1)=-SPL*CL PRIN038 G(2,2)=CPL PRIN039 G(2,3)=SL*SPL PRIN040 G(3,1)=CL*CP*CPL+SP*SL PRIN041 G(3,2)=CP*SPL PRIN042 G(3,3)=-SL*CP*CPL+SP*CL PRIN043 DENM=G(1,1)*G(2,2)*G(3,3)+G(1,2)*G(3,1)*G(2,3)+G(2,1)*G(3,2)*G(1,3PRIN044 1)-G(2,2)*G(3,1)*G(1,3)-G(1,2)*G(2,1)*G(3,3)-G(1,1)*G(3,2)*G(2,3) PRIN045 C********* THE MATRIX G1 IS THE INVERSE OF THE MATRIX G PRIN046 G1(1,1)=(G(2,2)*G(3,3)-G(3,2)*G(2,3))/DENM PRIN047 G1(1,2)=(G(3,2)*G(1,3)-G(1,2)*G(3,3))/DENM PRIN048 G1(1,3)=(G(1,2)*G(2,3)-G(2,2)*G(1,3))/DENM PRIN049 G1(2,1)=(G(3,1)*G(2,3)-G(2,1)*G(3,3))/DENM PRIN050 G1(2,2)=(G(1,1)*G(3,3)-G(3,1)*G(1,3))/DENM PRIN051 G1(2,3)=(G(2,1)*G(1,3)-G(1,1)*G(2,3))/DENM PRIN052 G1(3,1)=(G(2,1)*G(3,2)-G(3,1)*G(2,2))/DENM PRIN053 G1(3,2)=(G(3,1)*G(1,2)-G(1,1)*G(3,2))/DENM PRIN054 G1(3,3)=(G(1,1)*G(2,2)-G(2,1)*G(1,2))/DENM PRIN055 R0=EARTHR+XMTRH PRIN056 C********* CARTESIAN COORDINATES OF TRANSMITTER PRIN057 XR=R0*G(1,1) PRIN058 YR=R0*G(2,1) PRIN059 ZR=R0*G(3,1) PRIN060 CTHR=G(3,1) PRIN061 STHR=dSIN(dACOS (CTHR)) PRIN062 PHIR=dATAN2(YR,XR) PRIN063 ALPH=dATAN2(G(3,2),G(3,3)) PRIN064 C********* PRIN065 NR=6 PRIN066 NP=0 PRIN067 DO 7 NN=7,20 PRIN068 IF (W(NN+50).EQ.0.) GO TO 7 PRIN069 C********* DEPENDENT VARIABLE NUMBER NN IS BEING INTEGRATED PRIN070 C********* NR IS THE NUMBER OF DEPENDENT VARIABLES BEING INTEGRATED PRIN071 NR=NR+1 PRIN072 IF (W(NN+50).NE.2.) GO TO 7 PRIN073 C********* DEPENDENT VARIABLE NUMBER NN IS BEING INTEGRATED AND PRINTED.PRIN074 C********* NP IS THE NUMBER OF DEPENDENT VARIABLES BEING INTEGRATED AND PRIN075 C********* PRINTED PRIN076 NP=NP+1 PRIN077 C********* SAVE THE INDEX OF THE DEPENDENT VARIABLE TO PRINT PRIN078 NPR(NP)=NR PRIN079 HEAD1(NP)=HEADR1(NN) PRIN080 HEAD2(NP)=HEADR2(NN) PRIN081 UNIT(NP)=UNITS(NN) PRIN082 7 CONTINUE PRIN083 NP1=MIN0(NP,3) PRIN084 ccc PDEV=ABSORB=DOPP=0. PRIN085 pdev = 0.d0 absorb = 0.d0 dopp = 0.d0 C********* PRINT COLUMN HEADINGS AT THE BEGINNING OF EACH RAY PRIN086 10 IF (IHOP.NE.0.d0) GO TO 12 PRIN087 write(2,1100)(HEAD1(NN),HEAD2(NN),NN=1,NP1) PRIN088 1100 FORMAT (44X,7HAZIMUTH/43X,9HDEVIATION,8X,9HELEVATION/ PRIN089 1 19X,16HHEIGHT RANGE,1X,2(5X,12HXMTR LOCAL),5X,26HPOLARIZATIPRIN090 2ON GROUP PATH,5A6,A5) PRIN091 write(2,1150) (UNIT(NN),NN=1,NP1) PRIN092 1150 FORMAT (13X,2(8X,2HKM),2X,2(6X,3HDEG,5X,3HDEG),6X,12HREAL IMAG,PRIN093 1 7X,2HKM,4X,3(4X,A6,2X)) PRIN094 IF (RAYSET.EQ.0.d0) GO TO 12 PRIN095 C********* PUNCH A TRANSMITTER RAYSET PRIN096 TLOND=TLON*DEGS PRIN097 IF (TLOND.LT.0.d0) TLOND=TLOND+360.d0 PRIN098 TLATD=TLAT*DEGS PRIN099 IF (TLATD.LT.0.d0) TLATD=TLATD+360.d0 PRIN100 AZ=AZ1*DEGS PRIN101 EL=BETA*DEGS PRIN102 NHOP=HOP PRIN103 ccccc write(7,1200) c a ID(1),TYPE(NTYP),XMTRH,TLATD,TLOND,RCVRH,F,AZ,EL,POLAR PRIN104 c 1,NHOP,1HT PRIN105 1200 FORMAT (A3,A1,4PF9.0,2P2F6.0,4P2F9.0,5P2F10.0,5X,2P2F5.0,I1,A1) PRIN106 C********* PRIN107 12 V=0.d0 PRIN108 IF (N2.NE.0.d0) V=(R(4)**2+R(5)**2+R(6)**2)/N2-1.d0 PRIN109 H=R(1)-EARTHR PRIN110 STH=dSIN(THETA) PRIN111 CTH=dSIN(PID2-THETA) PRIN112 C********* CARTESIAN COORDINATES OF RAY POINT, ORIGIN AT TRANSMITTER PRIN113 XP=R(1)*STH*dSIN(PID2-PHI)-XR PRIN114 YP=R(1)*STH*dSIN(PHI)-YR PRIN115 ZP=R(1)*CTH-ZR PRIN116 C********* CARTESIAN COORDINATES OF RAY POINT, ORIGIN AT TRANSMITTER ANDPRIN117 C********* ROTATED PRIN118 EPS=XP*G1(1,1)+YP*G1(1,2)+ZP*G1(1,3) PRIN119 ETA=XP*G1(2,1)+YP*G1(2,2)+ZP*G1(2,3) PRIN120 ZETA=XP*G1(3,1)+YP*G1(3,2)+ZP*G1(3,3) PRIN121 RCE2=ETA**2+ZETA**2 PRIN122 RCE=dSQRT(RCE2) PRIN123 C********* GROUND RANGE PRIN124 RANGE=EARTHR*dATAN2(RCE,EARTHR+EPS+XMTRH) PRIN125 C********* ANGLE OF WAVE NORMAL WITH LOCAL HORIZONTAL PRIN126 ELL=dATAN2(R(4),dSQRT(R(5)**2+R(6)**2))*DEGS PRIN127 c c see if ray was perpendicular to dipole field line between c successive prints - check for sign change in cosine ell_dip = ell/degs + dipole_dip cos_ell = cos(ell_dip) ccc print*,' dip, local elev ', dipole_dip*degs,ell if(cos_ell_save / cos_ell .lt. 0.) print*,' Ray perp to field' if(cos_ell_save / cos_ell .lt. 0.) a write(2,2222) dipole_dip*degs,ell 2222 format(5x,'Ray was perp to the field', a ' -- DIP & ELEV =',2f8.3) cos_ell_save = cos_ell c C********* STRAIGHT LINE DISTANCE FROM TRANSMITTER TO RAY POINT PRIN128 SR=dSQRT(RCE2+EPS**2) PRIN129 IF (NP.LT.1) GO TO 16 PRIN130 DO 15 I=1,NP PRIN131 NN=NPR(I) PRIN132 15 RPRINT(I)=R(NN) PRIN133 16 IF (SR.GE.1.d-6) GO TO 20 PRIN134 C********* TOO CLOSE TO TRANSMITTER TO CALCULATE DIRECTION FROM PRIN135 C********* TRANSMITTER PRIN136 write(2,1500) V,NWHY,H,RANGE,ELL,POLAR,T,(RPRINT(NN),NN=1,NP1) PRIN137 1500 FORMAT (1X,1pE7.0,A8,0pF10.4,F11.4,26X,F8.3,F9.3,F8.3,8F12.4) PRIN138 GO TO 40 PRIN139 C********* ELEVATION ANGLE OF RAY POINT FROM TRANSMITTER PRIN140 20 EL=dATAN2(EPS,RCE)*DEGS PRIN141 IF (RCE.GE.1.d-6) GO TO 30 PRIN142 C********* NEARLY DIRECTLY ABOVE OR BELOW TRANSMITTER. CAN NOT CALCULATEPRIN143 C********* AZIMUTH DIRECTION FROM TRANSMITTER ACCURATELY PRIN144 write(2,2500) V,NWHY,H,RANGE,EL,ELL,POLAR,T,(RPRINT(NN),NN=1,NP1) PRIN145 2500 FORMAT (1X,1pE7.0,A8,0pF10.4,F11.4,17X,F9.3,F8.3,F9.3,F8.3, PRIN146 1 8F12.4) PRIN147 GO TO 40 PRIN148 C********* AZIMUTH ANGLE OF RAY POINT FROM TRANSMITTER PRIN149 30 ANGA=dATAN2(ETA,ZETA) PRIN150 AZDEV=180.d0-dMOD(540.d0-(AZ1-ANGA)*DEGS,360.d0) PRIN151 IF (R(5).NE.0.d0.OR.R(6).NE.0.d0) GO TO 34 PRIN152 C********* WAVE NORMAL IS VERTICAL, SO AZIMUTH DIRECTION CANNOT BE PRIN153 C********* CALCULATED PRIN154 write(2,3000)V,NWHY,H,RANGE,AZDEV,EL,ELL,POLAR,T,(RPRINT(NN),NN=1,PRIN155 1 NP1) PRIN156 3000 FORMAT (1X,1pE7.0,A8,0pF10.4,F11.4,F9.3,8X,F9.3,F8.3,F9.3,F8.3, PRIN157 1 8F12.4) PRIN158 GO TO 40 PRIN159 34 ANA=ANGA-ALPH PRIN160 SANA=dSIN(ANA) PRIN161 SPHI=SANA*STHR/STH PRIN162 CPHI=-dSIN(PID2-ANA)*dSIN(PID2-(PHI-PHIR))+SANA*dSIN(PHI-PHIR) PRIN163 1 *CTHR PRIN164 AZA=180.d0-dMOD (540.d0 a -(dATAN2(SPHI,CPHI)-dATAN2(R(6),R(5)))*DEGS,360.d0) PRIN165 cc type*,' PRINTR t ',t write(2,3500) a V,NWHY,H,RANGE,AZDEV,AZA,EL,ELL,POLAR,T,(RPRINT(NN),NN PRIN166 1 =1,NP1) PRIN167 3500 FORMAT (1X,1pE7.0,A8,0pF10.4,F11.4,2(F9.3,F8.3),F9.3,F8.3, PRIN168 1 8F12.4) PRIN169 C********* PRIN170 40 LINES=LINES+1 PRIN171 IF (NP.LE.3) GO TO 45 PRIN172 C********* ADDITIONAL LINE TO PRINT REMAINING DEPENDENT INTEGRATION PRIN173 C********* VARIABLES PRIN174 write(2,4000) (RPRINT(NN),NN=4,NP) PRIN175 4000 FORMAT (99X,3F12.4) PRIN176 LINES=LINES+1 PRIN177 45 IF (CARD.EQ.0.d0) RETURN PRIN178 C PRIN179 C********* PUNCH A RAYSET PRIN180 IF (AZDEV.LT.-90.d0) AZDEV=AZDEV+360.d0 PRIN181 IF (AZA.LT.-90.d0) AZA=AZA+360.d0 PRIN182 TDEV=T-SR PRIN183 NR=6 PRIN184 IF (W(57).EQ.0.d0) GO TO 47 PRIN185 C********* PHASE PATH PRIN186 NR=NR+1 PRIN187 PDEV=R(NR)-SR PRIN188 print*, ' Phase path ', pdev 47 IF (W(58).EQ.0.d0) GO TO 48 PRIN189 C********* ABSORPTION PRIN190 NR=NR+1 PRIN191 ABSORB=R(NR) PRIN192 C********* DOPPLER SHIFT PRIN193 48 IF (W(59).NE.0.d0) DOPP=R(NR+1) PRIN194 c write(7,4500)HPUNCH,RANGE,AZDEV,AZA,ELL,SR,TDEV,PDEV,ABSORB,DOPP, PRIN195 c 1 POLAR,IHOP,NWHY PRIN196 4500 FORMAT (4P2F9.0,3P3F6.0,3PF8.0,3P4F6.0,2P2F5.0,I1,A1) PRIN197 RETURN PRIN198 END PRIN199-