program irpred
c
c Version 2.3.0
c
c This program calculates the predictors from IR data
c using the radial profiles in IRRP1.dat (first choice)
c or IRRP2.dat (second choice).
c
c The IR predictor data is added to the SHIPS
c predictors on the file irpred.inp and the combined
c predictor set is written to the file irpred.dat.
c The format of the irpred.inp and .dat files is the same as the
c input for the lsdiag program.
c
c Modified May 2005 to print IR00 and IRM3 lines (MD)
c Modified May 2016 by MD for ver2.0.0
c IRM1 line added, minutes added to time variable
c Modified Nov 2016 by MD for ver2.1.0
c ircal3.f routine separated, ircal2.f removed (not needed)
c jday.f, jdayi.f, tdiff.f routines removed (libs linked instead)
c Modified Mar 2017 by KM for ver2.2.0
c ircal4.f routine added
c Modified Apr 2020 by SS for ver2.3.0
c updated for 7-day forecast
c
parameter (mxl=300)
character *180 iline(mxl)
character *180 irline1,irline2,irline3
character *180 hline,tline,tirline
c
character *25 fninf
c
c* -Changed mxr: 130-->150 and numpred:16-->19 for ircal3-->ircal4
parameter (mxr=150,mxv=10)
parameter (numpred=19)
dimension rad(mxr),vrad(mxr,mxv)
dimension idat(0:numpred)
c
data luinp,ludat,luout /10,11,12/
c
c Variables for ircal4 (slon read but not used)
real slat,slon,vmax
integer ivmax
c
c Initialize the IR profile variables to missing
do i=1,mxr
rad(i) = -99.
do j=1,mxv
vrad(i,j) = -99.
enddo
enddo
c
c Make the default IR data lines
m=9999
write(irline1,200) m,m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m
write(irline2,200) m,m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m
write(irline3,200) m,m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m,m,m,
+ m,m,m,m,m,m,m,m
200 format(11x,29(i4,1x))
irline1(157:160) = 'IR00'
irline2(157:160) = 'IRM1'
irline3(157:160) = 'IRM3'
c
c Open and read the input file
open(unit=luinp,file='irpred.inp',form='formatted',
+ status='old',err=900)
c
icount=0
1000 continue
ii = icount + 1
read(luinp,100,end=900) iline(ii)
100 format(a180)
c
icount = icount + 1
if (icount .gt. mxl) go to 900
c
if (iline(ii)(157:160) .eq. 'LAST') go to 2000
go to 1000
c
2000 continue
c
c Make sure the first line of the input file is the header file
hline = iline(1)
if (hline(157:160) .ne. 'HEAD') go to 900
c
c Read the date/time information from the header file
c* read(hline,400) iyr2,imon,iday,itime
c* 400 format(6x,3(i2),1x,i2)
read(hline,400) iyr2,imon,iday,itime,ivmax,
+ slat,slon
400 format(6x,3(i2),1x,i2,2x,i3,1x,f6.1,1x,f6.1)
vmax=float(ivmax)
slat=abs(slat)
c
c Calculate 4 digit year
if (iyr2 .gt. 50) then
iyear=iyr2+1900
else
iyear=iyr2+2000
endif
c
itime = 100*itime
c
c
c Open the t=0 file with the IR radial profile information
open(unit=ludat,file='IRRP1.dat',form='formatted',
+ status='old',err=800)
c
c Check to make sure the 1st GOES file is current
fninf='IRRP1.inf'
call irdtc(fninf,iyear ,imon ,iday ,itime ,
+ iyeari,imoni,idayi,itimei,idelhr,lrtime)
c
c write(6,*) iyear,imon,iday,itime
c write(6,*) iyeari,imoni,idayi,itimei
c write(6,*) 'idelhr=',idelhr
c
adt = abs(float(idelhr))
if (adt .gt. 3.5) go to 800
c
c irtime=idelhr
irtime=lrtime
c
c Read the profile file
do i=1,mxr
read(ludat,*,err=800,end=800) rad(i),(vrad(i,k),k=1,mxv)
enddo
c
c Calcualte the radial profiles
c* call ircal3(irtime,rad,vrad,mxr,mxv,idat,ierr)
call ircal4(slat,vmax,irtime,rad,vrad,mxr,mxv,idat,ierr)
if (ierr .ne. 0) go to 800
c
c Write IR predictors for output file
write(irline1(1:110),205) (idat(k),k=0,numpred)
205 format(10x,20(1x,i4))
c
800 continue
close(ludat)
c
c Open the t= -1.5 hr file with the IR radial profile information
open(unit=ludat,file='IRRP2.dat',form='formatted',
+ status='old',err=810)
c
c Check to make sure the 2nd GOES file is current
fninf='IRRP2.inf'
call irdtc(fninf,iyear ,imon ,iday ,itime ,
+ iyeari,imoni,idayi,itimei,idelhr,lrtime)
c
adt = abs(float(idelhr))
if (adt .gt. 4.0) go to 810
c
c write(6,*) iyear,imon,iday,itime
c write(6,*) iyeari,imoni,idayi,itimei
c write(6,*) 'idelhr=',idelhr
c
c irtime=idelhr
irtime=lrtime
c
c Read the profile file
do i=1,mxr
read(ludat,*,err=820,end=820) rad(i),(vrad(i,k),k=1,mxv)
enddo
c
c Calcualte the radial profiles
c* call ircal3(irtime,rad,vrad,mxr,mxv,idat,ierr)
call ircal4(slat,vmax,irtime,rad,vrad,mxr,mxv,idat,ierr)
if (ierr .ne. 0) go to 810
c
c Write IR predictors for output file
write(irline2(1:110),205) (idat(k),k=0,numpred)
c
810 continue
close(ludat)
c
c Open the t= -3 hr file with the IR radial profile information
open(unit=ludat,file='IRRP3.dat',form='formatted',
+ status='old',err=820)
c
c Check to make sure the 2nd GOES file is current
fninf='IRRP3.inf'
call irdtc(fninf,iyear ,imon ,iday ,itime ,
+ iyeari,imoni,idayi,itimei,idelhr,lrtime)
c
adt = abs(float(idelhr))
if (adt .gt. 6.5) go to 820
c
c write(6,*) iyear,imon,iday,itime
c write(6,*) iyeari,imoni,idayi,itimei
c write(6,*) 'idelhr=',idelhr
c
c irtime=idelhr
irtime=lrtime
c
c Read the profile file
do i=1,mxr
read(ludat,*,err=820,end=820) rad(i),(vrad(i,k),k=1,mxv)
enddo
c
c Calcualte the radial profiles
c* call ircal3(irtime,rad,vrad,mxr,mxv,idat,ierr)
call ircal4(slat,vmax,irtime,rad,vrad,mxr,mxv,idat,ierr)
if (ierr .ne. 0) go to 820
c
c Write IR predictors for output file
write(irline3(1:110),205) (idat(k),k=0,numpred)
c
820 continue
close(ludat)
c
c Open and write the output file
open(unit=luout,file='irpred.dat',form='formatted',
+ status='replace',err=900)
c
do i=1,icount-1
write(luout,100) iline(i)
enddo
write(luout,100) irline1
write(luout,100) irline2
write(luout,100) irline3
write(luout,100) iline(icount)
write(6,100) irline1
c
900 continue
c
close(luinp)
close(luout)
c
stop
end
subroutine irdtc(fninf,iyearn,imonn,idayn,itimen,
+ iyeari,imoni,idayi,itimei,idelhr,lrtime)
c This routine opens the file fninp that contains the date/time information
c in an IR file obtained from the NCEP IBM. The number of hours between the
c expected date/time (iyearn (i4), imon (i2), idayn (i2), itimen (UTC, i4))
c is calculated (idelhr) to make sure they are consistent with the date/time
c in the .inf file. The year,month,day and time on the .inf file are returned.
c If there is any problem opening or reading the .inf file, idelhr is set to -99.
c
c Modified 5/2/2016 to also provide relative time variable in hhmm
c format (lrtime) for lsdiag file.
c
character *3 cmoni
character *25 fninf,cdum
c
c Set the output variables to default values
idelhr=-99
iyeari=-99
imoni =-99
idayi =-99
itimei=-99
c
c Open and read the input file
luinf=75
open(file=fninf,unit=luinf,form='formatted',status='old',
+ err=9000)
c
c Skip 4 lines
read(luinf,100,end=9000,err=9000) cdum
read(luinf,100,end=9000,err=9000) cdum
read(luinf,100,end=9000,err=9000) cdum
read(luinf,100,end=9000,err=9000) cdum
100 format(a25)
c
c Read the day, month, year and time
read(luinf,110,err=9000,end=9000)
+ idayi,cmoni,iyeari,itimei
110 format(1x,i2,2x,a3,2x,i4,1x,i4)
c
close(luinf)
c
c Find the numerical month
call nummon(cmoni,imoni)
c
if (imoni .lt. 1 .or. imoni .gt. 12) go to 9000
9000 continue
c
c Calcuate the number of hours between the .dat and .inf date/times
itimeih = itimei/100
itimenh = itimen/100
itimeim = itimei - 100*itimeih
itimenm = itimen - 100*itimenh
c
call tdiff(iyeari,imoni,idayi,itimeih,
+ iyearn,imonn,idayn,itimenh,idelhr)
idelmin = itimeim-itimenm
c
if (idelmin .eq. 0) then
lrtime = 100*idelhr
else
if (idelhr .ge. 0) then
lrtime = 100*idelhr + idelmin
else
lrtime = -1*(100*iabs(idelhr+1) + (60-idelmin))
endif
endif
c
return
end
subroutine nummon(cmon,imon)
c This routine finds the numerical month for the 3 character cmon
c where nmon is the list of 3 character month abbreviations.
c If the cmon is not found in the list, imon is set to zero.
c
character *3 nmon(12),cmon
c
nmon( 1) = 'JAN'
nmon( 2) = 'FEB'
nmon( 3) = 'MAR'
nmon( 4) = 'APR'
nmon( 5) = 'MAY'
nmon( 6) = 'JUN'
nmon( 7) = 'JUL'
nmon( 8) = 'AUG'
nmon( 9) = 'SEP'
nmon(10) = 'OCT'
nmon(11) = 'NOV'
nmon(12) = 'DEC'
c
imon=0
c
do i=1,12
if (cmon .eq. nmon(i)) then
imon=i
go to 1000
endif
enddo
c
1000 continue
c
return
end