C .......................................................Mar 2014
C Extract foF2, hmF2, M3000 for => .inp to IRI-Plas at given sites
C 2 days in Sep. 2002, 2 days in Oct. 2003; 31 stations
C
C 
C Output: .inp file styrmndy.inp
C
      PROGRAM EXTRAHFM
C---------------------------------------------------------------
C data for 24 UT=0,1,2,...23hUT daily
C
C
C*********TEST
	CHARACTER*90 OUTFILE,infilest,infilef,infileh,infilem,txt
	CHARACTER*41 txt1,txt2
	character*2 ayr,amn,ST,sts(31)
	+,ADY,stinp,AUT,CYR,CMN,CDY
	INTEGER*2 IDY,iut,iut25
	CHARACTER*5 glats(31),glons(31),glat,glon
	CHARACTER(10) DD
      CHARACTER(10) TT
      CHARACTER(5) ZZ
	CHARACTER*4 STAT,STATS(31)
	CHARACTER*4 XXXX,DATF(0:23),DATH(0:23),DATM(0:23)
	dimension im(12)
	DATA  IM/31,28,31,30,31,30,31,31,30,31,30,31/
	DATA STS
	+/'at','ch','co','ds','eg','ga','gs','jr','kk','ld'
     +,'mg','mh','mm','mo','nq','ok','pa','ro','rv','sd'
     +,'sf','so','ss','tk','tq','tr','tz','vt','wk','wp'
     +,'yg'/
	data txt1/'020915  000 25.0   0  -1 -1 0  38.4  23.6'/
	data txt2/'   0.0   0.0   0.0     0  0.0     00.00  '/
C              /130313  000 25.0   0  -1 -1 0  51.6 358.7/
C              /   0.0   0.0   0.0  1HHH FF.F    TEC.XX  /
C
c
	      CALL DATE_AND_TIME(DATE=DD,TIME=TT,ZONE=ZZ)
      TT(5:)='      '
      WRITE(*,*) 'PC Date: Year,Month,Day = ',DD,'Time = ',TT
C
C read stations coordinates
	infilest='inputs1.inp'
  76	format(1X,I2,1X,A5,1X,A5,1X,A4,1X,A2)
   77	format(30X,A5,1X,A5,40X,A4,9X,A2)
	    OPEN(109,FILE=INFILEST,action='READ')
	do i=1,31
	read(109,77,end=78) glat,glon,stat,stinp
	if (stinp.eq.sts(i)) then
	write(*,76) i,glat,glon,stat,stinp
	glats(i)=glat
	glons(i)=glon
	stats(i)=stat
	 endif
	enddo
   78	close(unit=109)

C
  201	      WRITE(*,*)' ENTER yr,mn :'
		READ(*,15) AYR,AMN
   15	FORMAT(A2)
  202	      WRITE(*,*)' ENTER ADY:'
		READ(*,15) ADY
	read(ayr,*) yr
	iyr=int(yr)
	if(int(yr/4.)*4.eq.iyr) THEN 
	im(2)=29
	ndnr=366
	                          ELSE
	IM(2)=28
	ndnr=365
	ENDIF
	read(amn,*) xmn
	mn=int(xmn)
		read(ady,*) xdy
 	idy=int(xdy)

C start foF2, hmF2 or TEC

  199	nds=idy
C    
	txt(1:41)=txt1
	txt(42:82)=txt2
	txt(1:2)=AYR
	txt(3:4)=AMN
	txt(5:6)=ADY

	infilef='d:\web\graf\2010\styrf.br'
	infilef(15:16)=AYR
	infilef(20:21)=AYR
	infileh=infilef
	infileh(22:22)='h'
	infilem=infilef
	infilem(22:22)='m'

   61	format(3A2,1X,24(A4,1X))
C
	DO 88 jj=1,31
	st=sts(jj)
	stat=stats(jj)
	write(*,79) stat
	txt(42:82)=txt2
	txt(83:85)=stat
	txt(87:88)=st
	iflag=1  ! TEC input

	infilef(18:19)=st
	infileh(18:19)=st
	infilem(18:19)=st
	do i=0,23
	datf(i)='    0'
	dath(i)='    0'
	datm(i)='    0'
	dath(t)='    0'
	enddo

	    OPEN(11,FILE=INFILEF,action='READ')
   21	read(11,61,end=17) CYR,CMN,CDY,(datf(k),k=0,23)
	 if ((CYR.eq.AYR).and.(CMN.eq.AMN).and.(CDY.eq.ADY)) then
	goto 22
	else
	goto 21
	endif
   17	write(*,*) infilef
	pause ' '
   22	close(unit=11)
	    OPEN(11,FILE=INFILEH,action='READ')
   23	read(11,61,end=18) CYR,CMN,CDY,(dath(k),k=0,23)
	 if ((CYR.eq.AYR).and.(CMN.eq.AMN).and.(CDY.eq.ADY)) then
	goto 24
	else
	goto 23
	endif
  18	write(*,*) infileh
c-	pause ' '
   24	close(unit=11)
C M3000 files reading/input of M3000
	    OPEN(11,FILE=INFILEM,action='READ')
   25	read(11,61,end=19) CYR,CMN,CDY,(datm(k),k=0,23)
	 if ((CYR.eq.AYR).and.(CMN.eq.AMN).and.(CDY.eq.ADY)) then
	goto 26
	else
	goto 25
	endif
   19	write(*,*) infilem
c-	pause ' '
   26	close(unit=11)
C +++++++++++++++++++++++++++++++++++++++++++++++++++++++

	outfile='styrmndy.ins'	! Create .inp file for IRI-Plas
	outfile(3:4)=ayr
	outfile(5:6)=amn
	outfile(7:8)=ady

C

C start cycle on stations
   79 format(1X,A4)

c	if (jj.lt.65) goto 88   ! last station only
	outfile(1:2)=st
	
      OPEN(110,FILE=outfile)

	glat=glats(jj)
	glon=glons(jj)
	txt(31:35)=glat
	txt(37:41)=glon
c
c	
C
c ***************************************************************
C
  107	FORMAT(3A2,2X,24(A3,2X))
  108	FORMAT(3A2,2X,24(A3,A2))
  109	FORMAT(A90)
C
C              /130313  000 25.0   0  -1 -1 0  51.6 358.7/
C              /   0.0   0.0   0.0  1HHH FF.F    TEC.XX  /
C start cycle on UT
	do iut=0,23
      iut25=iut+25
	call blet2(iut25,AUT)
	txt(13:14)=AUT
	txt(61:70)='    0  0.0'
	XXXX=dath(iut)
	if ((XXXX.eq.'   0').or.(XXXX(4:4).eq.' ')) goto 41 ! No hmF2 for UT
C input hmF2
	txt(63:65)=XXXX(2:4)
	goto 42  ! avoid input of M3000

C input of M3000
   41	XXXX=datm(iut) 
	if ((XXXX.eq.'   0').or.(XXXX(4:4).eq.' ')) goto 42 ! No M3000 for UT
	txt(62:62)='1'
	txt(63:65)=XXXX(2:4)

C input of foF2
   42	XXXX=datf(iut) 
	if ((XXXX.eq.'   0').or.(XXXX(4:4).eq.' ')) goto 44 ! No foF2 for UT
		txt(67:68)=XXXX(2:3)
		txt(70:70)=XXXX(4:4)
C
C output 
   44	write(*,109) txt
      write(110,109) txt
	enddo                   ! end cycle on UT
   80	continue
C
   81	close(unit=110)
  88	continue              ! end cycle on stations
 111	format(A6,1X,24(I4,1X))
	GOTO 33
C
    2 CONTINUE
       write(*,*) 'INPUT FILE IS NOT IN YOUR DIRECTORY '
 	pause ' '
	goto 34

   33	continue
      write(*,*) outfile
      pause ' ' 
   34   STOP
       END

C***********************************************************
      subroutine blet2(ilet,alet2)
c     nn=1,2,3,4 

	CHARACTER*1 IN(0:9)
	INTEGER*2 ilet
	DATA IN/'0','1','2','3','4','5','6','7','8','9'/
	CHARACTER*2 alet2
CLjuba	a2='00'
	alet2='00'             !NEW+7
	j10=ilet/10
	j1=ilet-j10*10
    5	do i=0,9
	if (j10.eq.i) then
	alet2(1:1)=IN(i)
	endif
	enddo
   7	do i=0,9 
     	if (j1.eq.i) then
 	alet2(2:2)=IN(i)
	endif
	enddo 
	return
	end
C
      integer function ndoy(iyr,imn,idy)
C To define day-of-year
C
	       DIMENSION IM(12)
	integer*2 iyr,imn,idy
	DATA IM/31,28,31,30,31,30,31,31,30,31,30,31/
	z1=iyr/4.0
      jz1=int(z1)*4
      IF(jz1.EQ.iyr) THEN
               IM(2)=29
        ELSE
                IM(2)=28
	       ENDIF
C Day-of-year LDA1
		mosum=0
      if(imn.gt.1) then
         do 1234 i=1,imn-1
1234    mosum=mosum+im(i)
         endif
      ndoy=mosum+idy	
	END FUNCTION
C ===============================================================
