      program mapinter
C................................................Mar. 2013
C Extract JPL TEC map for maglat=-85.-80,...,80,85  (35 lines)
C                         maglon=0, 5, 10, ..., 355, 360 (0:72 columns)
C Input: <glatgrid.txt> <glongrid.txt> =geogr.coords. for mlat,mlon grids
C 
C
C
C ...............................................Jan. 2013
C extract-interpolate from JPL-TEC or UPC-TEC map 
C
C	(fm,glat,glon,syr,smn,sdy,sut,res)
C.................................................................Dec. 2012
C Gulyaeva Tamara <gulyaeva@izmiran.ru>
C
C To interpolate data from IONEX maps (TEC)
C Input : 
C		fm='t' means TEC
C         glat - geographic latitute
C         glon - geographic longitude
C         syr - YEAR = '99' for 1999,'00' for 2000, etc.
C         smn - month
C         sdy - day of month
C         sut - hour, UT
C.................................................................
C Outpur : res - value of TEC (0.1TECU)  = TEC*10
C 
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
C Use tc interpolated from directories dtc maps - update addresses for input file
C
	DIMENSION fhdat(0:72,71)	  ! maps fcF2,hcF2,TEC
	CHARACTER*1 fm
	CHARACTER*2 ayr,amn,ady,aut
crem	+,DY1,DY2
	CHARACTER*65 outfile,infile
	+,inglats,inglons
	INTEGER*2 jyr,jmn,jdy1,jdy2,jdy,jut
	DIMENSION im(12),itec(35,0:72)
	+,glats(35,0:72),glons(35,0:72)
	DATA  IM/31,28,31,30,31,30,31,31,30,31,30,31/
C+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  103	nn=1
  	fm='t'
C Read geogr. coords. for magn. grids
	inglats='glatgrid.txt'
	inglons='glongrid.txt'
  102	format(73(1X,I4))
  101	format(73(1X,F5.1))
	OPEN(10,file=inglats)
	OPEN(11,file=inglons)
C

	do i=1,35
	read(10,101) (glats(i,k),k=0,72)
	do k=0,72
	if (glats(i,k).gt.87.5) glats(i,k)=87.5
	if (glats(i,k).lt.-87.5) glats(i,k)=-87.5
	enddo
      read(11,101) (glons(i,k),k=0,72)
	enddo
	close(unit=10)
	close(unit=11)
C 
C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
	write(*,*) 'enter AYR='
	read(*,*) AYR

	read(ayr,*) ryr
	jyr=int(ryr)

	if(int(ryr/4.)*4.eq.jyr) THEN 
 	im(2)=29
	                          ELSE
	IM(2)=28
	ENDIF

crem	write(*,*) 'enter AMN='
crem	read(*,*) amn
crem	read(amn,*) rmn
crem	jmn=int(rmn)

crem	write(*,*) 'enter DY1, DY2 '
crem	read(*,*) dy1, dy2

crem	read(dy1,*) rdy1
crem	jdy1=int(rdy1)
crem	read(dy2,*) rdy2
crem	jdy2=int(rdy2)

C Start cycle on months:
	DO 163 jmn=1,12
	call blet2(jmn,AMN)
	jdy1=1
	jdy2=im(jmn)
C
	infile='d:\web\dts\YR\MN\tcMNDYUT.jYR'
	infile(12:13)=AYR
	infile(28:29)=AYR
	infile(15:16)=AMN
	infile(20:21)=AMN

C	outfile='d:\web\dtm\YR\MN\tmMNDYUT.jYR'
C Day-to-day cycle
	do 120 jdy=jdy1,jdy2
	call blet2(jdy,ADY)
	infile(22:23)=ADY
C       
C Cycle on UT
	do 115 jut=0,23
	call blet2(jut,AUT)
	infile(24:25)=AUT

C	outfile='d:\web\dtm\YR\MN\tmMNDYUT.jYR'
		outfile=infile
	outfile(10:10)='m'		!  
	outfile(19:19)='m'		!  
	OPEN(13,file=outfile)

C
 204	format(73(1X,F4.0))
C
	OPEN(12,FILE=INFILE)		

	do ilat=1,71
	read(12,204) (fhdat(k,ilat),k=0,72)	
	enddo
	close(unit=12)
C
	do 165 mlat=1,35
	  do 166 mlon=0,72
	glat=glats(mlat,mlon)
	glon=glons(mlat,mlon)

C      interpolate for glat,glon
C Linear Interpolation   
	aglat=-90.
	do 161 ilat=1,70 ! select pair of lats-lines (k) & (k+1)
c	aglat=-87.5,-85,...,87.5
	aglat=-90.+2.5*ilat
	aglat2=aglat+2.5
C
	if ((glat.ge.aglat).and.(glat.le.aglat2)) then  !   <<<<1
C Select pair of longs embracing glong
	 clon=-185.
	along=glon
	if (along.ge.180.) along=along-360.0
	do j=0,71 ! select pair of longs (j) & (j+1) >>>>>>>>>>>>>>>>
	clon=clon+5.
	clon2=clon+5.
	if ((along.ge.clon).and.(along.lt.clon2)) then  ! +++2
	res11=fhdat(j,ilat)
	res12=fhdat(j+1,ilat)
	res21=fhdat(j,ilat+1)
	res22=fhdat(j+1,ilat+1)
	res1=res11+(res12-res11)/5.*(along-clon)
	res2=res21+(res22-res21)/5.*(along-clon)
	res=res1+(res2-res1)/2.5*(glat-aglat)
	endif											 ! +++2
	enddo	! end select pair of longs (j) & (j+1)>>>>>>>>>>>>>>>
	exit
	endif											 !   <<<<1
  161 continue		   ! interpolation end
	itec(mlat,mlon)=nint(res)
  166	continue                    ! End cycle on mlong
  165	continue                    ! End cycle on mlat
	do mlat=1,35
crem 	write(*,102) (itec(mlat,k),k=0,72)
 	write(13,102) (itec(mlat,k),k=0,72)
	enddo
	close(unit=13)
 115	continue                    ! End cycle on UT
 120	continue                    ! End cycle on day-to0day
 163	continue                    ! End cycle on months
	write(*,121) ayr,amn
  121	format(1X,A2,1X,A2)
   23	pause ' '

	STOP
	END