C........................................................June.2015
C Interpolate maps in case of missed 1 day file, output = aver(prec_day + next_day)
C  jprg  => output file web/dtr/YY/MM/trMMDDUT.jYY
C         lati: -90.0,-88.0,...,0,...,88.0,90.0 ==> i=1,2,...91
C         longi:-180.0, -178.0,... 0...178.0 ==> j=0,...,179
C 
      PROGRAM decodjpr2
C----------------------------------------------------------------
C Select NEW+++:output 91 geographic latitudes: res1=90.0, 88.0,86,,.0..-90.0N 
C at: NEW+++: 0,1,...,91 geographic longitudes [-180,-178,...,178] = IONEX1 x 2 step
C Produce results for 0, 1,..., 23hr UT.iyr ~ 0.5, 1.5, ..., 23.5 hr UT 
C send results to d:/web/graf/dtr/YY/MM/trMMDDUT.jYY
C
C program for EXTRACTING FROM "JPRGdoy0.iYR" ionex1_TEC data 
C +++++++++
C
C*********TEST
      dimension glats(91),glons(0:179),tec_cur(91,0:179)
	+,tec_pre(91,0:179),tec_nex(91,0:179)
      integer itec(91,0:179)

	CHARACTER*2 AYR,AMN,ADY,AUT
	+,AMN_pre,ADY_pre,AYR_pre
	+,AMN_nex,ADY_nex,AYR_nex
	 INTEGER*2 iyr,imn,idy,iut
	+,imn_pre,idy_pre,iyr_pre,iyr_nex,imn_nex,idy_nex

C**********TEST
      CHARACTER*128 INFILE,OUTFILE,INFILE0
	CHARACTER*24 btxt
	CHARACTER(10) DD
      CHARACTER(10) TT
      CHARACTER(5) ZZ
	CHARACTER*4 AYEAR,NAME
	CHARACTER*20 start,enddesc,epoch,tlats
	character*3 aday,aday_pre
	dimension txx(0:179),im(12)

	DATA  IM/31,28,31,30,31,30,31,31,30,31,30,31/
C++
C++
      ICALLS=0
	epoch='EPOCH OF CURRENT MAP'
	tlats='LAT/LON1/LON2/DLON/H'
	btxt='-180.0 180.0   5.0 450.0'
	CALL DATE_AND_TIME(DATE=DD,TIME=TT,ZONE=ZZ)
      TT(5:)='      '
      WRITE(*,*) 'PC Date: Year,Month,Day = ',DD,'Time = ',TT

C cycle on year
C	AYEAR=DD(1:4)                  !Tamara               
	write(*,*) ' Enter YEAR: '
	read(*,118) AYEAR
	write(*,*) AYEAR
		ayr=ayear(3:4)
	read(ayr,*) ryr
	iyr=int(ryr)
	read(ayear,*) xyear
	IYYYY=int(xyear)
	read(ayr,*) ryr
	iyr=int(ryr)
C following comands for month
		WRITE(*,*) 'Enter AMN for input: '  !Tamara
	 	read(*,*) AMN					   !Tamara
C	AMN=DD(5:6)               !Tamara
	write(*,*) amn

	read(amn,*) xmn
	imn=int(xmn)
C
		idayy=365
	      IM(2)=28
        if(iyr/4*4.eq.iyr) then
	  idayy=366
	   	 IM(2)=29
	  endif

C  data-DAY
		WRITE(*,*) 'Enter dy for output: '  !Tamara
	 	read(*,*) idy					   !Tamara
  255	call blet2(idy,ady)
C	
  987	format(I2)
	call blet2(idy,ADY)
C
	lda=ndoy(iyr,imn,idy)
	lda_pre=lda-1
	iyr_pre=iyr
	call submmdd(lda_pre,iyr_pre,imn_pre,idy_pre)
		lda_nex=lda
	lda_nex=lda_nex+1
 	iyr_nex=iyr
	call submmdd(lda_nex,iyr_nex,imn_nex,idy_nex)
	call blet2(idy_nex,ADY_nex)
	call blet2(imn_nex,AMN_nex)
	call blet2(iyr_nex,AYR_nex)
	if (lda_pre.le.0) then
	iyr_pre=iyr-1
	imn_pre=12
	idy_pre=31
	lda_pre=ndoy(iyr_pre,imn_pre,idy_pre)
	endif
	call blet2(iyr_pre,ayr_pre)
	call blet2(imn_pre,amn_pre)
	call blet2(idy_pre,ady_pre)
  178	call blet3(lda,ADAY)
      call blet3(lda_pre,ADAY_pre)

C
C
 1001 CONTINUE
C
	glat=92.
	DO i=1,91
	glat=glat-2.0
	glats(i)=glat
	ENDDO
C
	glon=-182.
	do i=0,179
	glon=glon+2.0
	glons(i)=glon
	enddo
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

C
	IRUN=0 ! 
	iut0=0 ! 
c-	ady0=float(idy_pre)
CEND TEMP
	enddesc='END OF HEADER       '
	start='START OF TEC MAP    '
	name='jplr' !
C--	infile='d:\utec\2012\jplrDDD0.YRi'
c		infile='d:\web\graf\TECI\JPLRDDD0.YRi'	!Tamara
c	infile(22:24)=ADAY					 !Tamara
c	infile(27:28)=ayr					 !Tamara
c	infile0=infile
c	infile0(22:24)=ADAY_pre				 !Tamara
c	infile0(27:28)=ayr_pre				 !Tamara
C
c-    outfile='d:\web\graf\TECI\YRMNDYt.txt'	 !Tamara
	outfile='d:\web\dtr\yy\mn\trMNDYUT.jYR'
	outfile(12:13)=AYR						 !Tamara
		 outfile(15:16)=AMN					 !Tamara
		 outfile(20:21)=AMN					 !Tamara
	outfile(22:23)=ADY						 !Tamara
	outfile(28:29)=AYR	
	infile0=outfile
		 infile0(15:16)=AMN_pre					 !Tamara
		 infile0(20:21)=AMN_pre					 !Tamara
	infile0(22:23)=ADY_pre						 !Tamara
	infile0(28:29)=AYR_pre	
	infile=outfile
		 infile(15:16)=AMN_nex					 !Tamara
		 infile(20:21)=AMN_nex					 !Tamara
	infile(22:23)=ADY_nex						 !Tamara
	infile(28:29)=AYR_nex	
  118	format(A4)
cc	pause  ' '
C
C =========================================
	do k=0,179
	txx(k)=0.	! output tec for -180, -178,...00,02,...,178 deg.E
	do i=1,91
	itec(i,k)=0	! output TEC lat x long (2 x 2) array 
	enddo
	enddo
c
  198 format(1X,A39,1X,A39)
	write(*,*) ' Check name of infile: ',infile
c//	pause ' '
c
CC


c-	write(*,199) ayr,amn_pre,ady_pre,amn,ady
  199	format(5(1X,A2))
   27 format(A80)
C++
	ist=0
	JOUT=0
	JSTR=0
	JUC=0
	JMAG=0
CC
  800 continue !
CC
C  +++++++++++++++++++++++++++++
  113	continue     ! 
C+++
 1100 format(180(1X,F4.0))
C 
	DO 101 iut=0,23
	call blet2(iut, AUT)
	infile(24:25)=AUT						 !Tamara
	infile0(24:25)=AUT						 !Tamara

	OPEN(11,FILE=INFILE,ACTION='READ') ! TEC input for next day
	OPEN(10,FILE=INFILE0,ACTION='READ') ! TEC input for prec day

  29	continue
C end headlines
	do 30 ii=1,91
	read (10,1100,err=33,end=30) (tec_pre(ii,k),k=0,179)
   30	continue
			close(unit=10)

C
	do 1130 ii=1,91
	read (11,1100,err=33,end=1130) (tec_nex(ii,k),k=0,179)
	do k=0,179
	tec_cur(ii,k)=(tec_pre(ii,k)+tec_nex(ii,k))/2.
	enddo
 1130 	continue
			close(unit=11)

C
C+++++++++
	write(*,86) ayear,amn,ady,iut,epoch
C

CADD++++++++++++++++++++++++++++++++++++++++++++++++
   86 format(2X,A4,4X,A2,4X,A2,4X,I2,2(5X,'0'),24X,A16)
   85	continue
      outfile(24:25)=AUT						 !Tamara
	OPEN(12,FILE=OUTFILE) !TEC results   !

C
C +++++++++++ 
  181 format(2X,I2,2X,F5.1)
C
C ==================================================
C OUTPUT OF RESULTS
	glat=-92.
	DO nnlat=1,91,1		!Output for glats=-90.0,88,..., 90.0
	glat=glat+2.0
c//	write(*,280) glat,btxt,tlats
c-	write(12,280) glat,btxt,tlats
CREM	DO nnlat=71,1,-1     ! Temp: output for glats=-87.5,-85,...,+87.5
	do k=0,179
	itec(nnlat,k)=nint(tec_cur(nnlat,k))
	enddo
CREM	itec(nnlat,72)=itec(nnlat,0)
crem	write(*,182) (itec(nnlat,k),k=0,179)  !Test
	write(12,182) (itec(nnlat,k),k=0,179)
	ENDDO
 182	format(180(1X,I4))							
  280 format(2X,F6.1,A24,28X,A20)
  281	format(16(1X,I4))
  282	format(9(1X,I4))	 
  775 continue   ! End calculation for given UT
C
    9	continue
C Continue to the next UT hr :
  32		irun=irun+1 
		close(unit=12)
c	if (irun.gt.iutend) goto 101 ! continue to the next UT hr =>>>>
c      if (iut.lt.22) goto 80 ! continue to the next UT hr =>>>>
C	pause ' '
C++++++++++++++++++++++++++++
C
C 
  101 continue
	irun=0
	close(unit=111)
		goto 33	   ! end of calculations
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 '
 	pause ' '
	goto 35
      endif

   33 CONTINUE
   35      write(*,*) infile,outfile
	write(*,*) irun
C temp
c-	close(unit=12)
c	idy=idy+1
c	if (idy.le.im(imn)) then
c	goto 255
c	endif
      pause ' ' 
      STOP
       END

C***********************************************************
C -----------------------------------------------------------
      subroutine blet2(ilet,alet2)
c

	CHARACTER*1 IN(0:9)
	integer*2 ilet
	DATA IN/'0','1','2','3','4','5','6','7','8','9'/
	CHARACTER*2 alet2
	a2='00'
	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
      subroutine blet3(ilet,alet3)
c     nn=1,2,3,4 

	CHARACTER*1 IN(0:9)
	DATA IN/'0','1','2','3','4','5','6','7','8','9'/
	CHARACTER*3 alet3
	a3='000'
C
	j100=ilet/100
	      j10=(ilet-j100*100)/10
	j1=ilet-j100*100-j10*10
    3	do i=0,9
		if (j100.eq.i) then
	alet3(1:1)=IN(i)
	endif
	enddo
    5	do i=0,9
	if (j10.eq.i) then
	alet3(2:2)=IN(i)
	endif
	enddo
   7	do i=0,9 
     	if (j1.eq.i) then
 	alet3(3:3)=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 ===============================================================
C
      subroutine submmdd(lda,iyr,imn,idy)
C To define imn & idy using iyr & lda=day-of-year
C
	       DIMENSION IM(0:12)
	integer*2 iyr,imn,idy,nyr,ii
	DATA IM/31,31,28,31,30,31,30,31,31,30,31,30,31/
	nyr=iyr
	imn=0
	idy=0
	z1=nyr/4.0
      jz1=int(z1)*4
      IF(jz1.EQ.nyr) THEN
               IM(2)=29
	ndyr=366
        ELSE
                IM(2)=28
	ndyr=365
	       ENDIF
C Day-of-year = LDA
c//	write(*,*) lda
c//	pause ' '
C Define imn & idy
C NEW++++++++++++++++++++
	mosum=0
         do 1234 ii=1,12
	ndsum=mosum
	mosum=mosum+im(ii)
	imn=ii
	if (lda.le.mosum) then
	idy=lda-ndsum
	 goto 33
	exit
	endif
1234    continue
 33     continue
	RETURN
	END 
C ===============================================================
C ===========================================================
	subroutine SUBTREAD1(nn,yy)
C                                                             ! June 2015   
C
C   EXTRACT TEC FROM IONEX1 format to array yy(glat,glon) relevant to IONEX format
c
C
C
c
      DIMENSION txx(0:179),yy(91,0:179),xxx(16)
	character*60 txt

  90	format(3X,F5.1,A60)
180	format(16(1X,F4.0))
179	format(4(1X,F4.0))

c-	do j=1,3
	do k=0,179
	txx(k)=0.
	enddo
c-	enddo
	alat=92.	! resuts alat
	glat1=92.
	trem=20.

C===

	DO 20 ii=1,91 ! cycly on output lats
	alat=alat-2.0
   10	read(nn,90,err=2,end=9) glat,txt
	glat2=glat
C
	do j=1,11
	j1=(j-1)*16
    1	READ (nn,180,END=9,ERR=2) (xxx(i),i=1,16)
	do k=1,16
	if ((xxx(k).le.0.).or.(xxx(k).gt.9990.)) then
	if (txx(j1+k-1).gt.0.) then
	xxx(k)=txx(j1+k-1)	! replace wrong TEC
	                          else
	xxx(k)=trem ! replace wrong TEC
	endif
	                                       endif
	txx(j1+k-1)=xxx(k)
	trem=xxx(k)
	enddo
	enddo
	 READ (nn,179,END=9,ERR=2) (xxx(i),i=1,4)
		do k=1,4
	if ((xxx(k).le.0.).or.(xxx(k).gt.9990.)) then
	if (txx(j1+k-1).gt.0.) then
	xxx(k)=txx(175+k)	! replace wrong TEC
	                          else
	xxx(k)=trem ! replace wrong TEC
	endif
	                                       endif
	txx(175+k)=xxx(k)
	trem=xxx(k)
	enddo
c-	if ((alat.le.glat1).and.(alat.gt.glat2)) then
c-	rat=(glat1-alat)/2.
c-      do k=0,179
c-	txx(3,k)=txx(1,k)+(txx(2,k)-txx(1,k))*rat
c-	txx(1,k)=txx(2,k)
c-	enddo
	glat1=glat2
C Interpolate in longi
	alon=-182.
c-      alon=-180.
c-	glon1=-182.
c-	jout=0
crem	DO jj=0,71 ! cycle on output longs
	DO jlon=0,179 ! cycle on intput longs
	alon=alon+2.
c-   27	glon2=glon1+2.
	yy(ii,jlon)=txx(jlon)
crem	jg1=jg1+1
c-	if ((alon.ge.glon1).and.(alon.lt.glon2)) then
c-	yy(ii,jout)=txx(3,jlon)+(txx(3,jlon+1)-txx(3,jlon))*(alon-glon1)/2.
c-	glon1=glon2
c-	alon=alon+5.
c-	jout=jout+1
c-	endif 
c-	 glon1=glon2
   	ENDDO
c-	yy(ii,72)=yy(ii,0)
c-	                        else
c-	glat1=glat2
c-	do k=0,179
c-	txx(1,k)=txx(2,k)
c-	enddo
c-	goto 10
c-	endif
   20	continue
	goto 22
C ============================================	
    9	write(*,*) ' end of file'
	pause ' '
	goto 22
    2	write(*,*) ' error in input file'
 	pause ' '
C
   22 continue
	return
	end
C
C ===============================================================
