c irirtam.for c --------------------------------------------------------------- c IRI subroutines READIRTAMCOF, FOUT1, and GAMMA2 that read the c IRTAM coefficients and calculate NmF2, hmF2, and B0 for the c Real-Time IRI. The process of assimilating ionsonde and c Digisonde data into a Real-Time IRI is described in: Galkin et c al.,Radio Sci., 47, RS0L07, doi:10.1029/2011RS004952, 2012. c --------------------------------------------------------------- C 2016.01 06/17/16 First release C 2016.02 11/03/17 Added B1 input C 2016.03 12/20/19 Correct directory and name for IRTAM coeffs C c --------------------------------------------------------------- c SUBROUTINE READIRTAMCOF(ISEL,IDATE,IHHMM,MFF,FF) C ----------------------------------------------------------- C Finds and reads IRTAM coefficients for foF2, hmF2, B0, and C B1 for date IDATE (yyyymmdd) and time HOURUT (decimal hours, C Universal Time): C ISEL parameter filename C 0 foF2 IRTAM_foF2_COEFFS_yyyymmdd_hhmm.asc C 1 hmF2 IRTAM_hmF2_COEFFS_yyyymmdd_hhmm.asc C 2 B0 IRTAM_B0in_COEFFS_yyyymmdd_hhmm.asc C 3 B1 IRTAM_B1in_COEFFS_yyyymmdd_hhmm.asc c The coefficient array is stored in FF(MFF). C D. Bilitza, Jun 17.2016 C Added B1 input D. Bilitza, Nov 3, 2017 C ----------------------------------------------------------- c CHARACTER FILNAM*28 integer iuccir CHARACTER*100 FILNAM CHARACTER*120 LINE CHARACTER*12 INAME DIMENSION FF(MFF) konsol=6 IUCCIR=10 iname='foF2_COEFFS_' if(isel.gt.0) iname='hmF2_COEFFS_' if(isel.gt.1) iname='B0in_COEFFS_' if(isel.gt.2) iname='B1in_COEFFS_' c read foF2 coefficients WRITE(FILNAM,104) iname,idate,ihhmm ccc104 FORMAT('/home/bilitza/tango_home/IRI_data/iri/', ccc & 'iri_dev/IRTAM_dir/IRTAM_',A12,I8,'_',I4.4,'.ASC') 104 FORMAT('IRTAM_',A12,I8,'_',I4.4,'.ASC') OPEN(IUCCIR,FILE=FILNAM,STATUS='OLD',ERR=201, & FORM='FORMATTED') print*,mff,filnam C skip header with comments c do 1,13 READ(iuccir,1289) LINE 4686 READ(iuccir,1289) LINE 1289 Format(A120) if(LINE(1:12).NE."# END_HEADER") goto 4686 c read coefficients READ(iuccir,4689) FF 4699 FORMAT(E16.8) 4689 FORMAT(4E16.8) 4690 CLOSE(10) goto 300 201 WRITE(6,203) FILNAM 203 FORMAT(1X////, & ' The file ',A100,' is not in your directory.') 300 CONTINUE RETURN END C C real function FOUT1(XMODIP,XLATI,XLONGI,UT,TOV,FF0) c-------------------------------------------------------------- C CALCULATES CRITICAL FREQUENCY FOF2/MHZ USING SUBROUTINE GAMMA2. C XMODIP = MODIFIED DIP LATITUDE, XLATI = GEOG. LATITUDE, XLONGI= C LONGITUDE (ALL IN DEG.), MONTH = MONTH, UT = UNIVERSAL TIME C (DEC. HOURS), FF0 = ARRAY WITH RZ12-ADJUSTED CCIR/URSI COEFF. C D.BILITZA,JULY 85. C Modified to accept IRTAM coefficients FFO and TOV C D. Bilitza, Jun 17.2016 c-------------------------------------------------------------- DIMENSION FF0(1064) INTEGER QF(9) DATA QF/11,11,8,4,1,0,0,0,0/ FOUT1=GAMMA2(XMODIP,XLATI,XLONGI,UT,TOV,6,QF,9,76,13,1064,FF0) RETURN END C C REAL FUNCTION GAMMA2(SMODIP,SLAT,SLONG,HOUR,TOV, & IHARM,NQ,K1,M,MM,M3,SFE) C--------------------------------------------------------------- C CALCULATES GAMMA2=FOF2 OR M3000 USING CCIR NUMERICAL MAP C COEFFICIENTS SFE(M3) FOR MODIFIED DIP LATITUDE (SMODIP/DEG) C GEOGRAPHIC LATITUDE (SLAT/DEG) AND LONGITUDE (SLONG/DEG) C AND UNIVERSIAL TIME (HOUR/DECIMAL HOURS). IHARM IS THE MAXIMUM C NUMBER OF HARMONICS USED FOR DESCRIBING DIURNAL VARIATION. C NQ(K1) IS AN INTEGER ARRAY GIVING THE HIGHEST DEGREES IN C LATITUDE FOR EACH LONGITUDE HARMONIC WHERE K1 GIVES THE NUMBER C OF LONGITUDE HARMONICS. M IS THE NUMBER OF COEFFICIENTS FOR C DESCRIBING VARIATIONS WITH SMODIP, SLAT, AND SLONG. MM IS THE C NUMBER OF COEFFICIENTS FOR THE FOURIER TIME SERIES DESCRIBING C VARIATIONS WITH UT. C M=1+NQ(1)+2*[NQ(2)+1]+2*[NQ(3)+1]+... , MM=2*IHARM+1, M3=M*MM C SHEIKH,4.3.77. C Modified to accept IRTAM coefficients SFE(M3) and TOV C D. Bilitza, Jun 17.2016 C--------------------------------------------------------------- REAL*8 C(12),S(12),COEF(100),SUM DIMENSION NQ(K1),XSINX(13),SFE(M3) COMMON/CONST/UMR,PI XMIN=(HOUR-TOV)*60.+720 HOU=(15.0*HOUR-180.0)*UMR S(1)=SIN(HOU) C(1)=COS(HOU) DO 250 I=2,IHARM C(I)=C(1)*C(I-1)-S(1)*S(I-1) S(I)=C(1)*S(I-1)+S(1)*C(I-1) 250 CONTINUE MMM=M*MM DO 300 I=1,M MI=(I-1)*MM COEF(I)=SFE(MI+1) DO 301 J=1,IHARM COEF(I)=COEF(I)+SFE(MI+2*J)*S(J)+SFE(MI+2*J+1)*C(J) 301 CONTINUE COEF(I)=COEF(I)+SFE(MMM+I)*XMIN 300 CONTINUE SUM=COEF(1) SS=SIN(SMODIP*UMR) S3=SS XSINX(1)=1.0 INDEX=NQ(1) DO 350 J=1,INDEX SUM=SUM+COEF(1+J)*SS XSINX(J+1)=SS SS=SS*S3 350 CONTINUE XSINX(NQ(1)+2)=SS NP=NQ(1)+1 SS=COS(SLAT*UMR) S3=SS DO 400 J=2,K1 S0=SLONG*(J-1.)*UMR S1=COS(S0) S2=SIN(S0) INDEX=NQ(J)+1 DO 450 L=1,INDEX NP=NP+1 SUM=SUM+COEF(NP)*XSINX(L)*SS*S1 NP=NP+1 SUM=SUM+COEF(NP)*XSINX(L)*SS*S2 450 CONTINUE SS=SS*S3 400 CONTINUE GAMMA2=SUM RETURN END C C