C........................................................Dec 2023 C C........................................................Nov. 2021 C C........................................................Feb. 2016 C input file toMNDYUT.jYR fro magn. coords. C CC.......................................................Mar. 2015 C C........................................................Nov. 2014 C Produce Volume map (VEC) from TECeff array => /dev/evMMDDUT.jYR C C........................................................Jan. 2013 C (1/0)=North/South magnetic hemisphere C........................................................Nov.2012 C JPL hourly results C........................................................Apr.2012 C C Input: TECeff calculated with Cv along the Ne(h) profile C including ratio ov volume in shell of altitude = 200 km in [0 20,200]km C Ne coefficients Pi * dlat, dlat= 2.5 deg (const km size) C Pi(lat) = (dlon1+dlon2)/2, lat = -90.0,-87.5,...,90.0 [0,1,...,71] C C Dec 2011 C C Add GEC calculation for ut=24h C Add daily rate-of-change, daily mean GEC -> OUTFILE2 C C T.L.Gulyaeva...........................................July 2011 C subroutine SUBGEC8oc(AYEAR,AMN,ADY,gecdy) C---------------------------------------------------------------- C Select all latitudes: res1=-87.5N,-85N,,,,-5N,00,5,...,87.5 C C each line in input file for TEC at:[-180, -165,...0,15E,30E,...,165E] C separate files for UT=0,1,2,3,...23hUT C Produce GYRMNDY.res files (UT GEC) C output array: res(72) C program for EXTRACTING FROM "uhrgdoy0.iYR" ionex_TEC data C Output: GYRMNDY.res data for day into annual file gecyr.br C all results for given year, includint YRMNDY GEC(t), t=0,1,...,23 UT + daily mean, GEC(t)=GEC/1.E+32, GECU C outfile='gecyr.br' C if input file jplg (2h UT step), linear interpolation is applied for 1h UT results C C*********TEST CHARACTER(10) DD CHARACTER(10) TT CHARACTER(5) ZZ CHARACTER*80 INFILE,OUTFILE,OUTFILE2 +,OUTFILE1,OUTFILE0 CHARACTER*4 AYEAR CHARACTER*131 txt character*2 ayr,amn,ady,aut c +,dy1,dy2 character*3 aday integer*2 iidy,iyr,imn,iut dimension itxx(0:72),im(12),itxx0(0:72),rtav(72) +,gecoef(71),iamp(0:7) C- +,vec(0:72),ivec(0:72) dimension res(72),rgec(0:24),igec(0:24),iroc(0:24) c +,pgec(0:24) ! prec.day +,res1(72),res0(72),rgec1(0:24),rgec0(0:24),igec1(0:24),igec0(0:24) COMMON /CONST/UMR DATA IM/31,28,31,30,31,30,31,31,30,31,30,31/ C ! 87.5,85,...,2.5,0 DATA GECOEF/ + .67443E+10, .134533E+11, .201528E+11, .268159E+11, .334247E+11, + .399717E+11, .464414E+11, .528289E+11, .591099E+11, .652797E+11, + .713261E+11, .772395E+11, .830030E+11, .886092E+11, .940534E+11, + .993114E+11,.1043856E+12,.1092541E+12,.1139244E+12,.1183770E+12, +.1225998E+12,.1265905E+12,.1303368E+12,.1338486E+12,.1370943E+12, +.1400788E+12,.1428022E+12,.1452546E+12,.1474168E+12,.1493105E+12, +.1509128E+12,.1522346E+12,.1532637E+12,.1539977E+12,.1544452E+12, +.1545900E+12,.1544452E+12,.1539977E+12,.1532637E+12,.1522346E+12, +.1509128E+12,.1493105E+12,.1474168E+12,.1452546E+12,.1428022E+12, +.1400788E+12,.1370943E+12,.1338486E+12,.1303368E+12,.1265905E+12, +.1225998E+12,.1183770E+12,.1139244E+12,.1092541E+12,.1043856E+12, + .993114E+11, .940534E+11, .886092E+11, .830030E+11, .772395E+11, + .713261E+11, .652797E+11, .591099E+11, .528289E+11, .464414E+11, + .399717E+11, .334247E+11, .268159E+11, .201528E+11, .134533E+11, + .67443E+10/ C c comp=' ' CALL DATE_AND_TIME(DATE=DD,TIME=TT,ZONE=ZZ) TT(5:)=' ' c- WRITE(*,*) 'PC Date: Year,Month,Day = ',DD,'Time = ',TT c 118 format(A4) c pause ' ' c- WRITE(*,'(\A\)') ' ENTER YEAR OF INPUT FILE=' c- READ(*,'(A)') AYEAR ayr=ayear(3:4) read(ayr,*) ryr iyr=int(ryr) read(ayear,*) year 1 IRUN=0 c c- WRITE(*,'(\A\)') ' ENTER MN OF INPUT FILE=' c- READ(*,'(A)') AMN read(amn,*) xmn imn=int(xmn) c- WRITE(*,'(\A\)') ' ENTER , OF INPUT FILE=' c- READ(*,'(A)') DY1,DY2 CTemp DY1='01' C? read(dy1,*) xdy C? idy1=int(xdy) C? read(dy2,*) xdy2 C? idy2=int(xdy2) CC?? nrdy2=idy2 read(ADY,*) xdy c outfiles ================================= outfile2='gecoYR.jtc' ! 3h Amplitude of GEC C++ outfile2='c:\MSDEV\Projects\GEC\gecoYR.jtc' ! 3h Amplitude of GEC C outfile='etYRMNj.txt' outfile='eoYRMNj.txt' outfile(3:4)=ayr outfile(5:6)=amn outfile2(5:6)=ayr C++ outfile2(27:28)=ayr outfile1=outfile outfile0=outfile outfile1(10:10)='1' ! North outfile0(10:10)='0' ! South OPEN(12,FILE=OUTFILE,ACCESS='APPEND') OPEN(11,FILE=OUTFILE1,ACCESS='APPEND') OPEN(10,FILE=OUTFILE0,ACCESS='APPEND') OPEN(13,FILE=OUTFILE2,ACCESS='APPEND') CC CALCULATION OF DAY OF YEAR (LDA) iiyr=int(year) if(int(year/4.)*4.eq.iiyr) THEN im(2)=29 ELSE IM(2)=28 ENDIF c xdy=1. ! start day-of-month C C Start input for given UT:++++++++++++++++++++ C 273 FORMAT(1X,I2) C Cycle day-to-day +++++++++++++++++++++++++++++ 700 iidy=int(xdy) write(*,273) iidy do n=0,24 rgec(n)=0. enddo call blet2(iidy,ADY) C C lda=ndoy(iyr,imn,iidy) daynr=float(lda) 178 call blet3(lda,ADAY) c// infile='jplg0010.05i' !!!!!!!! crem infile='c:\web\des\YR\MN\ecMNDYUT.jYR' cold infile='c:\web\des\YR\MN\vcMNDYUT.jYR' Crem infile='c:\web\des\YR\MN\etMNDYUT.jYR' ! Plas1 results NEW geogr.coord!!! infile='c:\web\des\YR\MN\eoMNDYUT.jYR' ! Plas1 results NEW magn.coord!!! C infile(12:13)=ayr infile(15:16)=amn infile(20:21)=amn infile(22:23)=ady infile(28:29)=ayr C c 11 continue ! c 27 format(6X,F3.1) c *************************************************************** D start cycle on hr-to-hr input: DO 775 iut=0,23 rut=float(iut) call blet2(iut,AUT) infile(24:25)=aUT C- outvec=infile C- outvec(10:10)='v' C- outvec(19:19)='e' OPEN(111,FILE=INFILE) C- OPEN(14,FILE=OUTVEC) 80 GEC=0. GEC1=0. GEC0=0. cntav1=0. cntav0=0. do k=0,72 itxx0(k)=0. ! prec. lati TECeff enddo C D start cycle on glatitudes: glati=-90.0 DO 777 lat=1,71 ! Start cycle on glati glati=glati+2.5 res(lat)=0. ! total gec res1(lat)=0. ! GEC1 for each glat res0(lat)=0. ! GEC0 for each glat do k=0,72 itxx(k)=0. ! curr. lati TECeff enddo C Select coefficient: c>> if (glati.le.0) then C_ acoef=gecoef(lat)/1.0E+11 !scale acoef=gecoef(lat) c>> else c>> acoef=gecoef(73-lat) c>> endif C++++++++++++++++++ if(lat.eq.72) then goto 171 ! for glat=90.0 endif C read TECeff for current lati C READ (111,187,END=33,ERR=2) (itxx(k),k=0,72) READ (111,187,END=35,ERR=2) (itxx(k),k=0,72) 187 format(73(1X,I4)) c 88 continue c 85 continue C++++++++++++++++++ C C Select values for Long=-180,-165,...0,15,30,..,165...!!!!!!!!!!!!: C 171 nlon=-1 C sum=0. ! total GEC cnt=0. C prepare calculation of GEC1 (North) & GEC0 (South): C 1-North, 0- South Magnetic Hemisphere sum1=0. ! GEC1 cnt1=0. sum0=0. ! GEC0 cnt0=0. do 776 il=-180,180,5 ! for long-cycle nlon=nlon+1 glong=float(il) cc= if (glong.lt.0.) glong=glong+360. C cc= ires(nlati,iidy,ii)=itxx(nlon) C++++++++ c> if(lat.eq.1) then c> rtav(nlon+1)=(itxx(nlon)+itxx(nlon+1))/2. c> goto 172 c> endif c> if(lat.eq.72) then c> rtav(nlon+1)=(itxx0(nlon)+itxx0(nlon+1))/2. c> goto 172 c> endif C c> rtav(nlon+1)=(itxx0(nlon)+itxx0(nlon+1)+itxx(nlon)+itxx(nlon+1))/4. c// isum=isum+itxx(nlon) rtav(nlon)=float(itxx(nlon)) C> vec(nlon)=rtav(nlon)*acoef/1.0E+11 ! scale by E+11 => VECU*10 above grid point C- vec(nlon)=rtav(nlon)/acoef ! scale by E+11 => VECU*10 above grid point C- ivec(nlon)=nint(vec(nlon)) C- 172 if (nlon.lt.72) sum=sum+rtav(nlon) 172 sum=sum+rtav(nlon) cnt=cnt+1. C++ C Check position North/South of Magnetic equator call submageq(glati,glong,ns) C if (ns.eq.1) then cnt1=cnt1+1 ! North Hemisphere sum1=sum1+rtav(nlon) else cnt0=cnt0+1 !South hemisphere sum0=sum0+rtav(nlon) endif C 776 continue ! long cycle res(lat)=sum*1.0E+15*acoef c// res1(lat)=sum1*1.0E+15*acoef/cnt*cnt1 c// res0(lat)=sum0*1.0E+15*acoef/cnt*cnt0 res1(lat)=sum1*1.0E+15*acoef res0(lat)=sum0*1.0E+15*acoef 182 format(3A2,24(1X,I4)) 121 GEC=GEC+res(lat) GEC1=GEC1+res1(lat) GEC0=GEC0+res0(lat) do i=0,72 itxx0(i)=itxx(i) ! prec. lati TEC itxx(i)=0. enddo cntav1=cntav1+cnt1 cntav0=cntav0+cnt0 C C Output of VEC map for lati: C C> write(*,187) (ivec(k),k=0,72) C- write(14,187) (ivec(k),k=0,72) 777 continue ! lati cycle C --------------- C OUTPUT OF RESULTS 9 continue CC 176 FORMAT(1X,I2,2X,E12.6) c do kk=1,72 c1 write(*,176) iut,GEC c1 write(12,176) iut,GEC c enddo rgec(iut)=GEC/1.0E+32 rgec1(iut)=GEC1/1.0E+32 rgec0(iut)=GEC0/1.0E+32 775 continue ! End calculation for given UT C Continue to the next UT hr : 32 irun=irun+1 c?? iutend=iruns c?? if (irun.lt.iutend) goto 80 ! continue to the next UT hr =>>>> C CC C- GOTO 1277 C Output of daily results C C sum=0. sum1=0. sum0=0. do i=0,23 sum=sum+rgec(i) ! aver using UT=0,...,23 h sum1=sum1+rgec1(i) ! aver using UT=0,...,23 h sum0=sum0+rgec0(i) ! aver using UT=0,...,23 h igec(i)=nint(rgec(i)*100.) igec1(i)=nint(rgec1(i)*100.) igec0(i)=nint(rgec0(i)*100.) c// pgec(i)=rgec(i) ! keep results for prec.day+++++++++++++++++++++ c// pgec1(i)=rgec1(i) ! keep results for prec.day+++++++++++++++++++++ c// pgec0(i)=rgec0(i) ! keep results for prec.day+++++++++++++++++++++ c- if (i.lt.24) iroc(i)=nint((rgec(i+1)-rgec(i))*100.) ! rate-of-change enddo sum=sum/24. ! aver using UT=0,...,23 h sum1=sum1/24. ! aver using UT=0,...,23 h sum0=sum0/24. ! aver using UT=0,...,23 h c- iroc(24)=nint(sum*100) iaver=nint(sum*100.) iaver1=nint(sum1*100.) iaver0=nint(sum0*100.) gecdy=sum ! dayly GECaver using UT=0,...,23 h C C Add output of 4 title lines IF (iidy.eq.1) THEN txt(1:3)='GEC' write(12,275) txt txt(5:9)='South' write(10,275) txt txt(5:9)='North' write(11,275) txt write(*,275) txt 275 format(A10) WRITE(12,*) 'PC Date: Year,Month,Day = ',ADY,'Time = ',TT WRITE(11,*) 'PC Date: Year,Month,Day = ',ADY,'Time = ',TT WRITE(10,*) 'PC Date: Year,Month,Day = ',ADY,'Time = ',TT txt(1:9)='YYMMDD|UT' do k=0,24 iroc(k)=k enddo 276 format(A9,24(I2,3X)) write(*,276) txt,(iroc(k),k=0,23) write(12,276) txt,(iroc(k),k=0,23) write(11,276) txt,(iroc(k),k=0,23) write(10,276) txt,(iroc(k),k=0,23) txt(1:50)='--------------------------------------------------' txt(51:100)='--------------------------------------------------' txt(101:131)='-------------------------- ' 274 format(A131) write(12,274) txt write(11,274) txt write(10,274) txt write(*,274) txt ENDIF C write(*,177) ayr,amn,ady,(igec(i),i=0,23) ! GEC write(12,177) ayr,amn,ady,(igec(i),i=0,23) write(11,177) ayr,amn,ady,(igec1(i),i=0,23) write(10,177) ayr,amn,ady,(igec0(i),i=0,23) 177 format(3A2,24(1X,I4)) C+++++++++++++++++++++++++++++++++++++++++++++++ Ctemp Produce 3-h UT amplitude of GEC C do i=0,7 i1=i*3 xmin=9999. xmax=0. C do k=1,3 if (rgec(i1+k-1).lt.xmin) xmin=rgec(i1+k-1) ! 3h Amplitude if (rgec(i1+k-1).gt.xmax) xmax=rgec(i1+k-1) enddo iamp(i)=nint((xmax-xmin)*100.) enddo write(*,277) ayr,amn,ady,(iamp(i),i=0,7),iaver,iaver1,iaver0 ! 3h Amplitude write(13,277) ayr,amn,ady,(iamp(i),i=0,7),iaver,iaver1,iaver0 277 format(3A2,8(1X,I4),3(2X,I4)) C pause ' ' C 1277 CONTINUE C xdy=xdy+1. irun=0 cc= irun=irun+1 c pause ' ' close(unit=111) C iidy=int(xdy) C? if ((iidy.le.idy2).and.(idy1.lt.idy2)) then C? goto 700 ! continue for the next day of month C? else C- goto 33 ! end of calculations goto 35 ! end of calculations C? endif C c write (*,888) alat,alon,aday c 888 format(1X,' alat=',A3,' alon=',A3,' day=',A3) 2 CONTINUE if (irun.eq.0) then write(*,*) 'INPUT FILE IS NOT IN YOUR DIRECTORY ',infile pause ' ' goto 35 endif C 33 if ((iidy.le.idy2).and.(idy1.lt.idy2)) goto 32 35 write(*,*) infile,outfile write(*,*) irun close (unit=12) close (unit=13) C if ((imn.lt.12).and.(iyr.ne.15)) GOTO 1 C pause ' ' C STOP return END C***********************************************************