program dewgeo16mm
C........................................................Dec 2024
C IRI-Plas corrected median foF2n, hmF2n: Input: Shub fs & hs => fit GTEC
C
C........................................................May 2024
C Add: input GIMs.txa ('quiet Reference') to update GIM-foF2,GIM-hmF2 real time
C
C........................................................Apr. 2024
C current day uadg
C........................................................Dec 2023
C Day-to-day running UQRG & UADG
C........................................................Dec. 2021
C TEC median for -15d
C........................................................Oct. 2021
C 1) read GIM-TEC IONEX for 00, 15, 30, 45 min (c:\web\graf\dtu\YY\MM\ut\tcMNDDUTMM.uYR)
C 2) produce WTEC from (1) in IONEX  for 00, 15, 30, 45 min =>c:\web\graf\dmu\YY\MM\ut\tmMNDDUTMM.uYR (median for 15 prec days)
C                           =>c:\web\graf\dmu\YY\MM\ut\wcMNDDUTMM.uYR (W(log deviations) for 15 prec days)
C........................................................Jan. 2019
C UPC maps IONEX glat x glon; 00, 15, 30, 45 minutes
C output : ...uYR
C+Output 
C OUTFILEI: W(TEC) index => c:\web\graf\dmu\YY\MM\ut\wcmmddutmm.uyr !+-W(TEC) index IONEX
C
C .......................................................
C
C Input: prec_mn (all days) into lines 1:14(last day of prec_mn) and current_mn (till given day=>15 line) for fixed UT
C
C
C FORMAT from 15 prec. days
C Include logarithmic scale indices -4,  -3,    -2,    -1,0, 1,    2,    3,   4
C Equivalent to log10(TEC/TECmed)= <-.301, -.155, -.045   0,  0.45, .155, .301 >
C                              
C Ref. Gulyaeva T.L. Logarithmic scale of ionospheric disturbances,
C      Geomagn.and Aeronomy, 1996, 36, 1, 160-163, 1996.
C
CC infilerz='c:/web/graf/Rz.dat'
C
      DIMENSION IM(12)
     +,DA(0:28,0:72,71),DEV(0:72) ! TEC 6 prec_days+curr_day=7,longs=0:72,lats=71
C-     +,XMED(0:72,71)
     +,xx(71,0:72)
C     +,FU1(24),FU2(24),FUA(5,24),FUB(5,24),COEFA(74),COEFB(74)	
      dimension CC(15,0:72),SMED(0:72),IMED(0:72)
	+,rzi(12)						 ! Monthly SSN1
	integer*4 IRES(0:72),IDEL(0:72)
      integer IYR,IMN,IDY,IUT,jyr_pre,jmn_pre,jdy_pre, nrun        !+++
	+,IDY1,IDY2,nndy !Maps available for idy1 to idy2; extrap. to IDY2+1,idy2+2
     +,imni,iyri,jjdy
	CHARACTER*48 OUTFILEW,infilerz,OUTFILEM ! -15d median
      CHARACTER*16 indate
	      CHARACTER*1 qs
      	CHARACTER*4 AYEAR,AYEAR_pre,CYEAR
	CHARACTER*2 AYR,AYR_pre,RDY
	+,AUT,AMN,ADY,PYR,PMN,PDY,AMN_pre,ADY_pre
     +,PDYm,PYRm,PMNm
     +,AYRm,AMNm,ADYm  ! median 
	+,AYRI,AMNI,AMM ! minute
	CHARACTER*1 ft,fht
      COMMON /SSN/RZS,dlda,dnr
         COMMON /CONST/UMR,PI 
		COMMON /BLFUR/FUA,FUB        ! Coeff. Fourie
C
      DATA IM/31,28,31,30,31,30,31,31,30,31,30,31/
C
C
      ft='t'  ! TEC-map
        PI=ATAN(1.0)*4.
        UMR=PI/180.
C++       write(*,*) ' Enter <c> or <u>:'
C++	read(*,*) qs
C	qs='c'  ! current day UADG NORM
      	qs='u'  ! TEMPPPPPPPPPPPPPPPPP any day 'DATE' UQRG
	nrun=1
C
C++
C++
C-  799	format(A4,2(1X,A2),1X,A4,2(1X,A2),2(1X,I3))				 !FEB2015 NEW COMMAND NUMBER = 799
  799	format(I4,2(1X,I2),1X,I4,2(1X,I2),2(1X,I3))				 !+++
  798	format(' YEAR, AMN, ADY = ',A4,2(1X,A2),1X,A4,2(1X,A2),2(1X,I3))
	indate='c:/web/graf/date'    !!! Tamara !!!!!!!	 !FEB2015
C#	indate='/var/www/izmiran/ionosphere/weather/graf/date' !!! LIUBA !!FEB2015
  111	open(10,file=indate)							   !FEB2015			
C	read(10,799) AYEAR,AMN,ADY,AYEAR_pre,AMN_pre,ADY_pre
	read(10,799) KYEAR,KMN,KDY,KYEAR_1,KMN_1,KDY_1  !+++
	+,lda_cur,lda_1  ! FORMAT NUMBER = 799
	close(unit=10)
c--	kyr=kyear-kyear/100*100
	kyr_1=KYEAR_1-KYEAR_1/100*100
c--      lda_nex=lda_1 !?????????????????
C
C
C
      IF (nrun.eq.1) THEN
	  KYEAR=KYEAR_1
	  KMN=KMN_1
	  KDY=KDY_1-1
	     IF (KDY.eq.0) then
		 KMN=KMN-1
		 KDY=IM(KMN)
		 ENDIF      
		     IF (KMN.eq.0) then
			 KMN=12
			 KYEAR=KYEAR-1
	          KDY=31
	         ENDIF
	KYRem=KYR_1
	KMNem=KMN_1
	KDYem=KDY_1
	ENDIF
C
C
      IF (nrun.eq.2) THEN
	  KYEAR=KYEAR_1
	  KMN=KMN_1
	  KDY=KDY_1
	KYRem=KYR
	KMNem=KMN
	KDYem=KDY
	ENDIF
C
C
      IF (nrun.eq.3) THEN
	  KYEAR=KYEAR
	  KMN=KMN
	  KDY=KDY
	ENDIF
C
		call blet2(kyrem,AYRm) !med
 	call blet2(kmnem,AMNm) !med
	call blet2(kdyem,ADYm) !med
C
      kyr=kyear-kyear/100*100
		lda_cur=ndoy(kyr,kmn,kdy)	  
	call blet2(kyr,PYR)
       call blet2(kmn,PMN)
	call blet2(kdy,PDY)
C
       kyear_pre=kyear
	kmn_pre=kmn
	kdy_pre=kdy-1
	if (kdy_pre.eq.0) then
	 kmn_pre=kmn-1
	kdy_pre=IM(kmn_pre)
	   if(kmn_pre.eq.0) then
	   kmn_pre=12
      	kdy_pre=31
      	kyear_pre=kyear-1
	endif
	                 endif
C	       
        kyr_pre=kyear_pre-kyear_pre/100*100
		lda_pre=ndoy(kyr_pre,kmn_pre,kdy_pre)	  
	call blet2(kyr_pre,PYRm)
	call blet2(kmn_pre,PMNm)
	call blet2(kdy_pre,PDYm)
C++
c-       IF (nrun.eq.3) GOTO 101
C----------------------------------------------------
c-	 IF (nrun.eq.1) THEN
c-	lda_cur=lda_pre-1      ! nrun=1
c-	if (lda_cur.eq.0) then
c-	kyr=kyr-1
c-	kyr_pre=kyr
c-	lda_cur=365
C	lda_cur=366        ! 2024
c	lda_pre=lda_cur-1
c-	endif
c-      lda_pre=lda_cur-1
c-	if (lda_pre.eq.0) then
c-	kyr_pre=kyr_pre-1
c-	lda_pre=365 ! 2025
c-	endif
c-				endif
C	 ELSE
c-        if (nrun.eq.2) then
c-	lda_cur=lda_pre	   ! nrun=2
c-      lda_pre=lda_cur-1
c-	ENDIF
c-	if (lda_cur.eq.0) then
c-	kyr=kyr_pre
c-	kyr_pre=kyr
c-	lda_cur=365
C	lda_cur=366        ! 2024
c-	lda_pre=lda_cur-1
c-	endif
c-	ndnr_pre=365  !2025
C
c-      call submmdd(lda_cur,kyr,kmn,kdy)
c-	call submmdd(lda_pre,kyr_pre,kmn_pre,kdy_pre)
C
C---------------------------------------------
  101	call blet2(kmn,AMN)
      call blet2(kmn_pre,AMN_pre)
	call blet2(kdy,ADY)
      call blet2(kdy_pre,ADY_pre)
      call blet2(kyr,AYR)
      call blet2(kyr_pre,AYR_pre)
c-	if (nrun.eq.3) THEN 
	AYEAR='2000'
	call blet2(kyr,AYR)
	AYEAR(3:4)=AYR
	AYEAR_pre='2000'
	call blet2(kyr_pre,AYR_pre)
	AYEAR_pre(3:4)=AYR_pre
c-	GOTO 191
c-	endif
C?	if (lda.le.ndnr_pre) AYR=AYR_pre
c-	AYRm=AYR         !CHECK
c-	AYEAR='2000'
c-	AYEAR_pre='2000'
c-	AYEAR(3:4)=AYR
c-      AYEAR_pre(3:4)=AYR_pre 
  191	WRITE(*,798) AYEAR,AMN,ADY,AYEAR_pre,AMN_pre,ADY_pre
	+,lda_cur,lda_pre ! FORMAT NUMBER = 798
       dlda=float(lda_cur)
C++
  797 FORMAT(1X,A4,1X,12(1X,F5.1))
      infilerz='c:/web/graf/R12.dat'
	OPEN(10,file=infilerz,action='READ')
  796	read(10,797) CYEAR,rzi
	if (CYEAR.eq.AYEAR_pre) then
	goto 795
	                   else
	goto 796
	                   endif
  795 close(10)
C++
CC
C Start from uqrg (-2d)
CC
    1	CONTINUE		 ! cont for usrg (-1d)
C==	read(ayear,*) xyear
c-	read(ayear_pre,*) xyear
c-	IYEAR=int(xyear)
      IYEAR=KYEAR
	IYR=kyr
C==		ayr=ayear(3:4)
c-		ayr=ayear_pre(3:4)
c-	read(AYR,*) ryr
c-	IYR=int(ryr)   ! current yr
C
C-	read(AMN_pre,*) rmn
c-	imn=int(rmn)
C++
	lda=lda_cur           ! qs='q'
c-	 read(AMN,*) rmn
c-	imn=int(rmn)
        imn=kmn
	RZS=RZI(imn)            ! SSN1_12
c-		read(ayear,*) xyear
c-	IYEAR=int(xyear)
c-		ayr=ayear(3:4)
c-	read(AYR,*) ryr
c-	IYR=int(ryr)   ! current yr
C
   13 FORMAT(I2)
C>>dys=1           ! Cycle on day-to-day
C..	WRITE(*,*)' ENTER DY OF INPUT FILE=	'
c-      read(ADY,*) rdy  
c-	idy=int(rdy)
        IDY=kdy
      IDYS=IDY
      	idys2=idys
C
C
C++
	z1=iyear/4.0
      jz=int(z1)*4

      IF(jz.EQ.iyear) THEN
               IM(2)=29
	dnr=366.
	ndnr=366
        ELSE
                IM(2)=28
	  	dnr=365.
	ndnr=365
	       ENDIF
C++
 201	CONTINUE		 ! cycle on days
	AYRI=AYR_pre
   17 format (A2)
	iyri=iyr
c-	pyr=AYR ! prec_year
c-	jyr_pre=iyr
C
c-	call blet2(imn,amn)  ! current month
      IMN=kmn
	AMNI=AMN
	imni=imn
C	jmn_pre=imn
C
C
200   continue	  
	idy1=idys
C
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C
	IRUN=0 ! NORMAL COMMAND
       ndmn=im(imn)
C
       IF (nrun.eq.3) GOTO 192
	IF ((nrun.eq.2).and.(KMN_pre.eq.1).and.(KDY_pre.eq.1)) THEN
	PMN='12'
	PDY='31'
	PYR='25'  !2025
	lda=365 ! 2025
	                                    ELSE
	PYR=AYEAR_pre(3:4)
       PMN=AMN_pre
       PDY=ADY_pre
c-	lda=lda_pre
	                                ENDIF
C preceding month:
      read(PMN,*) rmn_pre
	jmn_pre=int(rmn_pre)
C        
	 read(PYR,*) ryr_pre
	jyr_pre=int(ryr_pre)
C
  192    idy1=idy
	idy2=idy
  80	call blet2(idy,ADY)
	if (qs.ne.'c') then
	iutcur=23
	immcur=45
	endif
c	if ((AMN.eq.'12').and.(AYR.ne.PYR)) then
c	call subtecut00c(qs,PYR,AMN,ADY,iutcur,immcur)
c	                    else
	call subtecut00c(qs,AYR,AMN,ADY,iutcur,immcur)
C?	call subtecut00c(qs,PYR,PMN,PDY,iutcur,immcur)
c      endif
C
C
	fht='t'				! only TEC source map
      lda1=lda_cur
C
C
	outfilew='c:\web\graf\dwa\YY\MM\ut\wcMMDDUTMM.uYY'   ! current mn W(TEC)
CC
	outfilem='c:\web\graf\dtu\YR\MN\ti\tmMNDDUTMM.uYR' !! uadg med -15d
		outfilem(17:18)=AYRm
c-		outfilem(17:18)=PYRm
	outfilem(20:21)=AMNm
c-      	outfilem(20:21)=PMNm
	outfilem(28:29)=AMNm
c-       	outfilem(28:29)=PMNm
	outfilem(30:31)=ADYm
c-      	outfilem(30:31)=PDYm
      	outfilem(38:39)=AYRm 
c-             	outfilem(38:39)=PYRm
CC
       if (qs.eq.'u') then
      outfilew(15:15)='u'
	else
      outfilew(15:15)='a'
	endif
	outfilew(17:18)=ayr
	outfilew(20:21)=amn
	outfilew(38:39)=ayr
	outfilew(28:29)=amn
CD	outfiled=outfilew
CD	outfiled(9:10)='el'
CD	outfiled(18:19)='dt'
C	 
C
C  Start cycle on UT		=======================================================================
C
C     	 infile1='c:\web\dtm\YY\MM\tcMMDDUTMM.jYY'   ! preceding mn
         iutend=23
	if (qs.ne.'c') then
	iutcur=23
	immcur=45
	endif
      if (qs.eq.'c') iutend=iutcur
C
      DO 300 iut=0,iutend   ! norm
CREM	 iut=23   ! temp
	call blet2(iut,AUT)
	outfilew(32:33)=aut
CD	outfiled(24:25)=aut
	outfilem(32:33)=AUT
C add 
		ut=float(iut)
C Start cycle on minutes
C
      immend=45
	if ((qs.eq.'c').and.(iut.eq.iutcur)) immend=immcur
      DO 301 imm=0,immend,15	  ! norm
CREM	 imm=45  ! temp
	
      call blet2(imm,AMM)
C
	outfilew(34:35)=amm
CD	outfiled(26:27)=amm
	outfilem(34:35)=AMM
CC
	IF (nrun.eq.1) THEN
	OPEN(113,FILE=OUTFILEM)
	 ENDIF

C	DA(0:28,0:72,71)
	 DO n=1,71				! lats
	      DO I=0,28			! 27 prec days
       DO K=0,72				! lons
	      DA(I,K,n)=0.
	 enddo
	      ENDDO
	ENDDO

c	 ,XMED(0:72,71)
 	do k=0,72
	ires(k)=0
	idel(k)=0
	imed(k)=0
c-	xmed(k,i)=0.
       DEV(k)=0.
	enddo
C  Avoid input of prec. month:
      ncnt=0        ! count of days 1,2,...,7   
Cpre      if (idy1.ge.15) goto 230	  ! goto input of current month
      if (idy1.gt.15) goto 230	  ! goto input of current month
C  1st input of preceding mon-data

C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
	nndy=im(jmn_pre)	  ! prec. mn
C	nnd1=nndy-5
Cpre	nnd1=nndy-13		  !-15d
	nnd1=nndy-14		  !-15d
	if (idys.gt.1) nnd1=nnd1+idys-1
C
C  day-to-day input for prec mon
C
	DO 400 jdy_pre=nnd1,nndy  !=15 days for fixed UT
	call blet2(jdy_pre,PDY)
c-	if (nrun.eq.3) then
c-	CALL subreadmmu2c(qs,pyrm,pmnm,pdy,aut,amm,xx)
c-	               else
	CALL subreadmmu2c(qs,pyr,pmn,pdy,aut,amm,xx)
c-	                endif
C
C DA(0:28,0:72,71)
c.		! TEC at glon=-180,-175,...,180
	 	DO ilat=1,71	  ! < LAT-BY-LAT
	   do k=0,72		 !<<<<<<<<<<<<<<<<<<
	da(28,k,ilat)=xx(ilat,k) ! TEC 
	enddo
C
	ENDDO				  !< LAT-BY-LAT
C  107	FORMAT(180(1X,F4.1))		   
C Move data day-by-day up so that last day of month nndy=>day_27:
C
 	do ilat=1,71		  ! < LAT-BY-LAT
	do j=1,28			  ! <day-by-day
	do k=0,72			  ! all longi
	da(j-1,k,ilat)=da(j,k,ilat) ! move data one day up  
	enddo				  ! all longi
	enddo				  ! all days
	enddo				   ! all lati
	ncnt=ncnt+1
  400	CONTINUE		  !<day-by-day
C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
C
C start cycle on day-to-day for current mon  ======================================
C
  230	continue
	nndy=im(imn)	 ! end of current month

  490 format(1X,' year=',I2,' month=',I2,' day=',I2,' ut=',I2)
C	idy=idy1-1
C	lday2=ndoy(iyr,imn,idy2) ! function to define day-of-year
C
C Current month ONLY FOR ONE DAY
	call blet2(iyr,AYR)
	call blet2(imn,AMN)
C
C  499	idy=idy+1	 ! day-by-day input for idy=nndy,...,idy3, current month
  499	CONTINUE
  	   DO 500 jjdy=1,idy2
	call blet2(jjdy,RDY)
C 
c-		if ((AMN.eq.'12').and.(AYR.ne.PYR)) then
c-      CALL subreadmmu2c(qs,pyr,amn,r dy,aut,amm,xx)
c-	                                     else
      CALL subreadmmu2c(qs,ayr,amn,rdy,aut,amm,xx)
c-	                               endif
C DA(0:28,0:72,71)
	 	DO ilat=1,71	  ! < LAT-BY-LAT
	   do k=0,72		 !<<<<<<<<<<<<<<<<<<
	da(28,k,ilat)=xx(ilat,k) ! TEC 
	enddo
C
	ENDDO				  !< LAT-BY-LAT
C++
  700 format(A128)
C++
C
C   reading map for current month
C
cr	glat=-90.
 704 	 continue
C
C Move data one day up: (1)	if(idy.lt.idy1)  (2) a
C
 	do ilat=1,71		  ! < LAT-BY-LAT
 	do j=1,28			  ! <day-by-day
	do k=0,72			  ! all longi
	da(j-1,k,ilat)=da(j,k,ilat) ! move data one day up 
	enddo				  ! all longi
	enddo				  !<day-by-day
      enddo                 ! all lati 
	ncnt=ncnt+1
		if ((ncnt.lt.15).or.(jjdy.lt.idy2)) goto 500
C
C
C Calculations for  day=idy: da(28,long,lat)
C	
  240	CONTINUE
	outfilew(30:31)=ady
CD	outfiled(22:23)=ady
C
		OPEN(112,FILE=OUTFILEW)
CD	OPEN(113,FILE=OUTFILED)
C
C  VERIFICATION OF MEDIAN	FOR TEC_MAP(long,lat) for one day (28)
  115 CONTINUE
c
      DO 150 ilat=1,71			  ! Cycle on lati
C
C Median calculation:
C
	DO ilon=0,72	
C Prepare common array for -7 prec days, CC(7,0:72), for median calculation:
C     C DA(0:28,0:72,71) ! TEC -15 prec_days+curr_day=15,longs=73,lats=71

C__	DO kk=21,27
	DO kk=13,27
C__	kmed=kk-20
	kmed=kk-12
	CC(kmed,ilon)=da(kk,ilon,ilat)
	ENDDO	   ! cycle kk
	  
	  ENDDO      !cycle ilon
C__	ndays=7
  	ndays=15
	nlons=73
C-	CALL smedstd3mm(CC,SMED,STD,STD1,STD2)
	CALL smed15mm1(CC,SMED)
        DO 144 M=0,72				  ! cycle on longi
CC      xpar=da(28,m,ilat) ??
	apar=da(27,m,ilat)
	amed=smed(m)
	IF (nrun.eq.1) THEN
	imed(m)=nint(smed(m)*10.) !
	ENDIF
	delog=0.
C-	if ((apar.le.0.).or.(amed.le.0.)) then
C-	write(*,*) ilat,ilon
C-	continue 
C-	endif
	if (apar.le.0.) then
	apar=1.
	endif
	call windex(apar,amed,delog,iw)
C
CD	DEV(m)=delog
	ires(m)=iw  
CD	IDEL(m)=nint(delog*1000.)
 144  CONTINUE        ! cycle M ~ glong

C- 	WRITE(*,498) (IRES(K),K=0,72)		!Temp
 	WRITE(112,498) (IRES(K),K=0,72)  !W(TEC)
C- 	WRITE(*,498) (IDEL(K),K=0,72)      ! TEMP
CD 	WRITE(113,498) (IDEL(K),K=0,72)  ! log(TEC/TECmed)*1000
	IF (nrun.eq.1) THEN
C- 	WRITE(*,498) (IMED(K),K=0,72)		!Temp
 	WRITE(113,498) (IMED(K),K=0,72)  !W(TEC)
	ENDIF
C		 
 150	continue         ! cycle ilat
C
C END OF MEDIAN_TECxhi_MAP FOR DAY idy(28)
C ------------------------------------------------------
C?? 498  FORMAT(360(1X,I4))
 498  FORMAT(73(1X,I4))
	CLOSE(unit=112)
	IF (nrun.eq.1) THEN
	CLOSE(unit=113)
	ENDIF
C

 202  CONTINUE
C
C__	 ncnt=6
C	 ncnt=14
C	  						  !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
C
C =========================================
  500 continue   ! end day-by-day cycle     !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
c-        call blet2(iut,AUT)
C-  	call subratmapfh(qs,AYR,AMN,ADY,AUT,AMM)

  301 CONTINUE		! end MIN cycle
	 
C	idy=indy
  300 CONTINUE		! end UT cycle
C
C  	call subratmapfh(qs,AYR,AMN,ADY)
Cpre 	call subratmapfh_15(qs,AYR,AMN,ADY,AYRm,AMNm,ADYm)
C-- 	call subratmapfhn_15(qs,AYR,AMN,ADY,PYRm,AMNm,ADYm) ! Shub corr foF2n, hmF2n
c-		if ((AMN.eq.'12').and.(AYR.ne.PYR)) then
c-  	call subratmapfhs_15(qs,PYR,AMN,ADY,PYRm,AMNm,ADYm) ! Shub foF2s, hmF2s
c-	                               ELSE
 	call subratmapfhs_15(qs,AYR,AMN,ADY,AYRm,AMNm,ADYm) ! Shub foF2s, hmF2s
c-	                                ENDIF
c       idy=idy2
c	call blet2(idy,ADY)
c-                 	if ((AMN.eq.'12').and.(AYR.ne.PYR)) then
c-      call mapionexmm2c('w',qs,ayear_pre,amn,ady,iutend,immend)  !CC18????????????
c-	call mapionexmm2c('f',qs,ayear_pre,amn,ady,iutend,immend)  !CC18????????????
c-	call mapionexmm2c('h',qs,ayear_pre,amn,ady,iutend,immend)  !CC18????????????
c-	                            ELSE
C	                   
      call mapionexmm2c('w',qs,ayear,amn,ady,iutend,immend)  !CC18????????????
	call mapionexmm2c('f',qs,ayear,amn,ady,iutend,immend)  !CC18????????????
	call mapionexmm2c('h',qs,ayear,amn,ady,iutend,immend)  !CC18????????????
c-	                           ENDIF
C	call subpianiaWmm1(qs,AYEAR_pre,AMN,ADY)
CC
CC
	if (nrun.eq.2) THEN
	DO 505 imm=0,45,15
	           	if ((AMN.eq.'12').and.(AYR.ne.PYR)) then
	fht='f'
      call subextr2upc(fht,kyr_pre,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	fht='h'
      call subextr2upc(fht,kyr_pre,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	fht='t'
      call subextr2upc(fht,kyr_pre,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	              ELSE
	fht='f'
      call subextr2upc(fht,iyr,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	fht='h'
      call subextr2upc(fht,iyr,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	fht='t'
      call subextr2upc(fht,iyr,imn,idy,iyr_nex,imn_nex,idy_nex,imm)
	                ENDIF
  505	CONTINUE
	ENDIF
C
      GOTO 102
    2 CONTINUE
	       write(*,*) 'INPUT FILE IS NOT IN YOUR DIRECTORY '
	   	pause ' '
    5	iend=1
  102	continue
            write(*,*) ' yr=',AYR,'  mn=',AMN,' dy1=',idy1,' dy2=',idy2 
 1333	continue
C	if (imn.lt.12) GOTO 1333
       nrun=nrun+1
	qs='c'
        if (nrun.le.3) goto 111
	pause ' ' 
      STOP
	END