subroutine tintrp2a1(f,g,dx,dy,obstime,gridtime, &
nlevs,mype,nflds)
!$$$ subprogram documentation block
! . . . .
! subprogram: tintrp2a1
! prgmmr: parrish org: np22 date: 2013-01-26
!
! abstract: same as tintrp2a but for special case of n=1, with n argument removed.
! This has been created to solve problem of type mismatch debug
! compile error on WCOSS.
!
! program history log:
! 2013-01-26 parrish
!
! input argument list:
! f - input interpolator
! dx,dy - input x,y,z-coords of interpolation points (grid units)
! obstime - time to interpolate to
! gridtime - grid guess times to interpolate from
! nlevs - number of vertical levels over which to perform the
! 2-d intrpolation
! mype - mpi task id
! nflds - number of guess times available to interpolate from
!
! output argument list:
! g - output interpolatees
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
use kinds, only: r_kind,i_kind
use gridmod, only: istart,jstart,nlon,nlat,lon1,lon2,lat2
use constants, only: zero,one
implicit none
! Declare passed variables
integer(i_kind) ,intent(in ) :: nlevs,mype,nflds
real(r_kind),dimension(lat2,lon2,nlevs,nflds),intent(in ) :: f
real(r_kind) ,intent(in ) :: dx,dy,obstime
real(r_kind),dimension(nflds) ,intent(in ) :: gridtime
real(r_kind),dimension(nlevs) ,intent( out) :: g
! Declare local variables
integer(i_kind) m1,ix1,iy1,ix,ixp,iyp
integer(i_kind) iy,itime,itimep,j,k
real(r_kind) delx,delyp,delxp
real(r_kind) dely,delt,deltp
m1=mype+1
ix1=int(dx)
iy1=int(dy)
ix1=max(1,min(ix1,nlat))
delx=dx-float(ix1)
dely=dy-float(iy1)
delx=max(zero,min(delx,one))
ix=ix1-istart(m1)+2
iy=iy1-jstart(m1)+2
if(iy<1) then
iy1=iy1+nlon
iy=iy1-jstart(m1)+2
end if
if(iy>lon1+1) then
iy1=iy1-nlon
iy=iy1-jstart(m1)+2
end if
ixp=ix+1; iyp=iy+1
if(ix1==nlat) then
ixp=ix
end if
if(obstime > gridtime(1) .and. obstime < gridtime(nflds))then
do j=1,nflds-1
if(obstime > gridtime(j) .and. obstime <= gridtime(j+1))then
itime=j
itimep=j+1
delt=((gridtime(j+1)-obstime)/(gridtime(j+1)-gridtime(j)))
end if
end do
else if(obstime <=gridtime(1))then
itime=1
itimep=1
delt=one
else
itime=nflds
itimep=nflds
delt=one
end if
deltp=one-delt
delxp=one-delx; delyp=one-dely
do k=1,nlevs
g(k)=(f(ix,iy,k,itime)*delxp*delyp+f(ixp,iy,k,itime)*delx*delyp &
+ f(ix,iyp,k,itime)*delxp*dely+f(ixp,iyp,k,itime)*delx*dely)*delt &
+(f(ix,iy,k,itimep)*delxp*delyp+f(ixp,iy,k,itimep)*delx*delyp &
+ f(ix,iyp,k,itimep)*delxp*dely +f(ixp,iyp,k,itimep)*delx*dely)*deltp
end do ! end loop over vertical levs
return
end subroutine tintrp2a1
subroutine tintrp2a11(f,g,dx,dy,obstime,gridtime, &
mype,nflds)
!$$$ subprogram documentation block
! . . . .
! subprogram: tintrp2a11
! prgmmr: parrish org: np22 date: 2013-01-26
!
! abstract: same as intrp2a but for special case n=1 and nlevs=1, with n,nlevs arguments removed.
! This has been created to solve problem of type mismatch debug compile
! error on WCOSS.
!
! program history log:
! 2013-01-26 parrish
!
! input argument list:
! f - input interpolator
! dx,dy - input x,y,z-coords of interpolation points (grid units)
! obstime - time to interpolate to
! gridtime - grid guess times to interpolate from
! mype - mpi task id
! nflds - number of guess times available to interpolate from
!
! output argument list:
! g - output interpolatees
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
use kinds, only: r_kind,i_kind
use gridmod, only: istart,jstart,nlon,nlat,lon1,lon2,lat2
use constants, only: zero,one
implicit none
! Declare passed variables
integer(i_kind) ,intent(in ) :: mype,nflds
real(r_kind),dimension(lat2,lon2,nflds),intent(in ) :: f
real(r_kind) ,intent(in ) :: dx,dy,obstime
real(r_kind),dimension(nflds) ,intent(in ) :: gridtime
real(r_kind) ,intent( out) :: g
! Declare local variables
integer(i_kind) m1,ix1,iy1,ix,ixp,iyp
integer(i_kind) iy,itime,itimep,j
real(r_kind) delx,delyp,delxp
real(r_kind) dely,delt,deltp
m1=mype+1
ix1=int(dx)
iy1=int(dy)
ix1=max(1,min(ix1,nlat))
delx=dx-float(ix1)
dely=dy-float(iy1)
delx=max(zero,min(delx,one))
ix=ix1-istart(m1)+2
iy=iy1-jstart(m1)+2
if(iy<1) then
iy1=iy1+nlon
iy=iy1-jstart(m1)+2
end if
if(iy>lon1+1) then
iy1=iy1-nlon
iy=iy1-jstart(m1)+2
end if
ixp=ix+1; iyp=iy+1
if(ix1==nlat) then
ixp=ix
end if
if(obstime > gridtime(1) .and. obstime < gridtime(nflds))then
do j=1,nflds-1
if(obstime > gridtime(j) .and. obstime <= gridtime(j+1))then
itime=j
itimep=j+1
delt=((gridtime(j+1)-obstime)/(gridtime(j+1)-gridtime(j)))
end if
end do
else if(obstime <=gridtime(1))then
itime=1
itimep=1
delt=one
else
itime=nflds
itimep=nflds
delt=one
end if
deltp=one-delt
delxp=one-delx; delyp=one-dely
g=(f(ix,iy,itime)*delxp*delyp+f(ixp,iy,itime)*delx*delyp &
+ f(ix,iyp,itime)*delxp*dely+f(ixp,iyp,itime)*delx*dely)*delt &
+(f(ix,iy,itimep)*delxp*delyp+f(ixp,iy,itimep)*delx*delyp &
+ f(ix,iyp,itimep)*delxp*dely +f(ixp,iyp,itimep)*delx*dely)*deltp
return
end subroutine tintrp2a11
subroutine tintrp2a11_csln(f,g,gw,dx,dy,obstime,gridtime, &
mype,nflds)
!$$$ subprogram documentation block
! . . . .
! subprogram: tintrp2a11
! prgmmr: parrish org: np22 date: 2013-01-26
!
! abstract: same as intrp2a11 but calculate the background at observation location using the grid points
! with the same surface condition (land versus water).
!
! program history log:
! 2013-01-26 parrish
!
! input argument list:
! f - input interpolator
! dx,dy - input x,y,z-coords of interpolation points (grid units)
! obstime - time to interpolate to
! gridtime - grid guess times to interpolate from
! mype - mpi task id
! nflds - number of guess times available to interpolate from
!
! output argument list:
! g - output interpolatees : land or point with consistent surface
! gw - output interpolatees : water
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
use kinds, only: r_kind,i_kind
use gridmod, only: istart,jstart,nlon,nlat,lon1,lon2,lat2
use guess_grids, only: isli2
use constants, only: zero,one
implicit none
! Declare passed variables
integer(i_kind) ,intent(in ) :: mype,nflds
real(r_kind),dimension(lat2,lon2,nflds),intent(in ) :: f
real(r_kind) ,intent(in ) :: dx,dy,obstime
real(r_kind),dimension(nflds) ,intent(in ) :: gridtime
real(r_kind) ,intent( out) :: g
real(r_kind) ,intent( out) :: gw
! Declare local variables
integer(i_kind) m1,ix1,iy1,ix,ixp,iyp
integer(i_kind) iy,itime,itimep,j
real(r_kind) delx,delyp,delxp
real(r_kind) dely,delt,deltp
!
real(r_kind) gwater,gland,wt_water,wt_land,gone,wt_one
m1=mype+1
ix1=int(dx)
iy1=int(dy)
ix1=max(1,min(ix1,nlat))
delx=dx-float(ix1)
dely=dy-float(iy1)
delx=max(zero,min(delx,one))
ix=ix1-istart(m1)+2
iy=iy1-jstart(m1)+2
if(iy<1) then
iy1=iy1+nlon
iy=iy1-jstart(m1)+2
end if
if(iy>lon1+1) then
iy1=iy1-nlon
iy=iy1-jstart(m1)+2
end if
ixp=ix+1; iyp=iy+1
if(ix1==nlat) then
ixp=ix
end if
if(obstime > gridtime(1) .and. obstime < gridtime(nflds))then
do j=1,nflds-1
if(obstime > gridtime(j) .and. obstime <= gridtime(j+1))then
itime=j
itimep=j+1
delt=((gridtime(j+1)-obstime)/(gridtime(j+1)-gridtime(j)))
end if
end do
else if(obstime <=gridtime(1))then
itime=1
itimep=1
delt=one
else
itime=nflds
itimep=nflds
delt=one
end if
deltp=one-delt
delxp=one-delx; delyp=one-dely
g=(f(ix,iy,itime)*delxp*delyp+f(ixp,iy,itime)*delx*delyp &
+ f(ix,iyp,itime)*delxp*dely+f(ixp,iyp,itime)*delx*dely)*delt &
+(f(ix,iy,itimep)*delxp*delyp+f(ixp,iy,itimep)*delx*delyp &
+ f(ix,iyp,itimep)*delxp*dely +f(ixp,iyp,itimep)*delx*dely)*deltp
if((isli2(ix,iy)==0.and.isli2(ix,iyp)==0.and.isli2(ixp,iy)==0.and.isli2(ixp,iyp)==0) &
.or.(isli2(ix,iy)>=1.and.isli2(ix,iyp)>=1.and.isli2(ixp,iy)>=1.and.isli2(ixp,iyp)>=1)) then
gw=g
! obs and backgrond have consistent land/water type
else
gwater=0.0_r_kind
wt_water=0.0_r_kind
gland=0.0_r_kind
wt_land=0.0_r_kind
wt_one=delxp*delyp
gone=f(ix,iy,itime)*wt_one*delt+f(ix,iy,itimep)*wt_one*deltp
if(isli2(ix,iy)==0) then
gwater=gwater+gone
wt_water=wt_water+wt_one
else
gland=gland+gone
wt_land=wt_land+wt_one
endif
wt_one=delx*delyp
gone=f(ixp,iy,itime)*wt_one*delt+f(ixp,iy,itimep)*wt_one*deltp
if(isli2(ixp,iy)==0) then
gwater=gwater+gone
wt_water=wt_water+wt_one
else
gland=gland+gone
wt_land=wt_land+wt_one
endif
wt_one=delxp*dely
gone=f(ix,iyp,itime)*wt_one*delt+f(ix,iyp,itimep)*wt_one*deltp
if(isli2(ix,iyp)==0) then
gwater=gwater+gone
wt_water=wt_water+wt_one
else
gland=gland+gone
wt_land=wt_land+wt_one
endif
wt_one=delx*dely
gone=f(ixp,iyp,itime)*wt_one*delt+f(ixp,iyp,itimep)*wt_one*deltp
if(isli2(ixp,iyp)==0) then
gwater=gwater+gone
wt_water=wt_water+wt_one
else
gland=gland+gone
wt_land=wt_land+wt_one
endif
if(wt_water > 0.0_r_kind) then
gw=gwater/wt_water
else
gw=g
endif
if(wt_land > 0.0_r_kind) then
g=gland/wt_land
endif
endif
return
end subroutine tintrp2a11_csln
subroutine tintrp2a11_indx(dx,dy,obstime,gridtime, &
mype,nflds,ix,ixp,iy,iyp,itime,itimep)
!$$$ subprogram documentation block
! . . . .
! subprogram: tintrp2a11_indx
! prgmmr: zupanski org: CSU/CIRA/Data Assimilation group date: 2015-07-10
!
! abstract: same as tintrp2a11 but for horizontal grid indexes surrounding
! an observation point
!
! program history log:
! 2015-07-10 zupanski: add output grid indexes
! 2018-01-01 apodaca: compatibility-related updates
!
! input argument list:
! dx,dy - input x,y,z-coords of interpolation points (grid units)
! obstime - time to interpolate to
! gridtime - grid guess times to interpolate from
! mype - mpi task id
! nflds - number of guess times available to interpolate from
!
! output argument list:
! ix,iy,ixp,iyp - horizontal grid indexes
! itime,itimep - time grid indexes
!
! attributes:
! language: f90
! machine: ibm RS/6000 SP
!
!$$$
use kinds, only: r_kind,i_kind
use gridmod, only: istart,jstart,nlon,nlat,lon1
use constants, only: zero,one
implicit none
! Declare passed variables
integer(i_kind) ,intent(in ) :: mype,nflds
real(r_kind) ,intent(in ) :: dx,dy,obstime
real(r_kind),dimension(nflds) ,intent(in ) :: gridtime
integer(i_kind) ,intent(out ) :: ix,ixp
integer(i_kind) ,intent(out ) :: iy,iyp
integer(i_kind) ,intent(out ) :: itime,itimep
! Declare local variables
integer(i_kind) m1,ix1,iy1
integer(i_kind) j
real(r_kind) delx
real(r_kind) dely,delt
m1=mype+1
ix1=int(dx)
iy1=int(dy)
ix1=max(1,min(ix1,nlat))
delx=dx-float(ix1)
dely=dy-float(iy1)
delx=max(zero,min(delx,one))
ix=ix1-istart(m1)+2
iy=iy1-jstart(m1)+2
if(iy<1) then
iy1=iy1+nlon
iy=iy1-jstart(m1)+2
end if
if(iy>lon1+1) then
iy1=iy1-nlon
iy=iy1-jstart(m1)+2
end if
ixp=ix+1; iyp=iy+1
if(ix1==nlat) then
ixp=ix
end if
if(obstime > gridtime(1) .and. obstime < gridtime(nflds))then
do j=1,nflds-1
if(obstime > gridtime(j) .and. obstime <= gridtime(j+1))then
itime=j
itimep=j+1
delt=((gridtime(j+1)-obstime)/(gridtime(j+1)-gridtime(j)))
end if
end do
else if(obstime <=gridtime(1))then
itime=1
itimep=1
delt=one
else
itime=nflds
itimep=nflds
delt=one
end if
return
end subroutine tintrp2a11_indx