var/da/da_radiance/da_predictor_crtm.inc

   1 subroutine da_predictor_crtm(pred,npred,nlevels,temp,hum,t_skin,pf)
   2
   3    implicit none
   4
   5    ! temp - model level temperatures (k)
   6    ! hum  - model level moistures    (g/kg)
   7    ! t-skin - model skin temperature (k)
   8    ! nlevels - number of model levels (0:model top)
   9    ! pm   - model level pressure (hpa)
  10    ! pf   - full level pressure  (hpa)
  11
  12    ! pred(1) - 1000-300 thickness
  13    ! pred(2) - 200-50 thickness
  14    ! pred(3) - t_skin
  15    ! pred(4) - total column precipitable water
  16
  17    integer, intent(in)  :: npred,nlevels
  18    real,    intent(in)  :: temp(nlevels), hum(nlevels), t_skin
  19    real,    intent(in)  :: pf(0:nlevels)
  20    real,    intent(out) :: pred(npred)
  21
  22    real, parameter :: kth = gas_constant/gravity
  23    real, parameter :: kpc = 100.0/gravity
  24    real :: tv(nlevels), qm(nlevels)
  25    real :: dlp(nlevels), dp(nlevels)
  26    real    :: add_thk
  27    integer :: itmp(1)
  28    integer :: index1000, index300, index200, index50
  29
  30    if (trace_use) call da_trace_entry("da_predictor_crtm")
  31
  32    qm=hum*0.001  ! g/kg to kg/kg
  33
  34    dlp(1:nlevels) = log(pf(1:nlevels)) - log(pf(0:nlevels-1))
  35     dp(1:nlevels) =     pf(1:nlevels)  -     pf(0:nlevels-1)
  36
  37    ! 0.0 find the pressure level index that is closest to
  38    ! 1000hPa, 300hPa, 200hPa, and 50hPa respectively
  39
  40    ! note: pf levels are 0:nlevels,
  41    ! return values of minloc are 1:nlevels+1
  42    itmp(1:1) = minloc(abs(pf(:)-1000.0))
  43    index1000 = itmp(1)-1
  44    itmp(1:1) = minloc(abs(pf(:)-300.0))
  45    index300  = itmp(1)-1
  46    itmp(1:1) = minloc(abs(pf(:)-200.0))
  47    index200  = itmp(1)-1
  48    itmp(1:1) = minloc(abs(pf(:)-50.0))
  49    index50   = itmp(1)-1
  50
  51    ! 1.0 convert all temperatures to virtual temperatures
  52    ! ----------------------------------------------------
  53    tv = temp*(1.0+0.608*qm)
  54
  55    if ( (1000.0 > pf(nlevels)) ) then
  56       add_thk = kth*( tv(nlevels)*(log(1000.0)-log(pf(nlevels))) )  ! approximation
  57    else
  58       add_thk = 0.0
  59    end if
  60
  61    ! 2.0 construct averages for nesdis thick layers
  62    ! ----------------------------------------------
  63
  64    pred(1) = kth*sum( tv(index300+1:index1000)*dlp(index300+1:index1000) ) + add_thk
  65    pred(2) = kth*sum( tv(index50+1:index200)*dlp(index50+1:index200) )
  66    pred(3) = t_skin
  67    pred(4) = kpc*sum( qm(1:nlevels)*dp(1:nlevels) )
  68
  69    if (trace_use) call da_trace_exit("da_predictor_crtm")
  70
  71 end subroutine da_predictor_crtm
  72