var/da/da_varbc/da_varbc_pred.inc

   1 subroutine da_varbc_pred(iv)
   2
   3    !---------------------------------------------------------------------------
   4    !  PURPOSE: Calculate bias predictors
   5    !
   6    ! pred(1) - 1 (Constant)
   7    ! pred(2) - 1000-300 thickness 43 1005.43-521.46 thickness
   8    ! pred(3) - 200-50 thickness   43  194.36-56.73  thickness
   9    ! pred(4) - T_skin
  10    ! pred(5) - total column precipitable water
  11    ! pred(6) - satellite scan position
  12    ! pred(7) - satellite scan position **2
  13    ! pred(8) - satellite scan position **3
  14    ! pred(9) - CRTM Gamma Correction
  15    !
  16    !  Called from da_varbc_coldstart
  17    !
  18    !  HISTORY: 10/26/2007 - Creation                  Tom Auligne
  19    !---------------------------------------------------------------------------
  20
  21   implicit none
  22
  23   integer,parameter   :: npred_hk = 4       ! Number of H&K bias predictors
  24
  25   type (iv_type), intent(inout) :: iv       ! O-B structure.
  26
  27   integer   :: inst, nlevels, i, npredmax, n, num_rad
  28   real      :: pred_hk(npred_hk)
  29
  30   if ( iv%num_inst < 1 ) RETURN
  31
  32   if (trace_use) call da_trace_entry("da_varbc_pred")
  33
  34   do inst = 1, iv%num_inst                                  ! loop for sensor
  35     npredmax = iv%instid(inst)%varbc_info%npredmax
  36     if (npredmax <= 0) cycle                                ! VarBC instr only
  37
  38     num_rad = iv%instid(inst)%num_rad
  39     if (iv%instid(inst)%info%n2 < iv%instid(inst)%info%n1) cycle
  40     do n = iv%instid(inst)%info%n1, iv%instid(inst)%info%n2  ! loop for pixel
  41       ! get H&K predictors
  42       !--------------------
  43       if (rtm_option==rtm_option_rttov) then
  44 #ifdef RTTOV
  45         nlevels=iv%instid(inst)%nlevels
  46         call da_predictor_rttov(inst, &
  47              pred_hk(1:npred_hk), npred_hk, nlevels, &
  48              iv%instid(inst)%t(1:nlevels,n), &
  49              iv%instid(inst)%mr(1:nlevels,n)/q2ppmv, &
  50              iv%instid(inst)%ts(n))
  51 #endif
  52       else if (rtm_option==rtm_option_crtm) then
  53 #ifdef CRTM
  54         nlevels=iv%instid(inst)%nlevels-1
  55         call da_predictor_crtm( &
  56              pred_hk(1:npred_hk), npred_hk, nlevels,&
  57              iv%instid(inst)%tm(1:nlevels,n), &
  58              iv%instid(inst)%qm(1:nlevels,n), &
  59              iv%instid(inst)%ts(n), &
  60              iv%instid(inst)%pf(0:nlevels,n))
  61 #endif
  62       end if
  63
  64       ! Populate predictors
  65       !---------------------
  66       iv%instid(inst)%varbc_info%pred(1:npredmax,n) = 0.0
  67
  68       ! Constant predictor
  69       iv%instid(inst)%varbc_info%pred(1,n) = 1.0
  70
  71       ! H&K predictors
  72       if (npredmax >= 2) iv%instid(inst)%varbc_info%pred(2,n) = pred_hk(1)
  73       if (npredmax >= 3) iv%instid(inst)%varbc_info%pred(3,n) = pred_hk(2)
  74       if (npredmax >= 4) iv%instid(inst)%varbc_info%pred(4,n) = pred_hk(3)
  75       if (npredmax >= 5) iv%instid(inst)%varbc_info%pred(5,n) = pred_hk(4)
  76
  77       ! Scan predictors
  78       if (varbc_scan(inst) == 1) then ! use scanpos for polar-orbiting sensors
  79          if (npredmax >= 6) iv%instid(inst)%varbc_info%pred(6,n) = iv%instid(inst)%scanpos(n)
  80          if (npredmax >= 7) iv%instid(inst)%varbc_info%pred(7,n) = iv%instid(inst)%scanpos(n)**2
  81          if (npredmax >= 8) iv%instid(inst)%varbc_info%pred(8,n) = iv%instid(inst)%scanpos(n)**3
  82       else ! use satellite zenith angle for geostationary sensors, now turn them off in VARBC.in
  83          if (npredmax >= 6) iv%instid(inst)%varbc_info%pred(6,n) = iv%instid(inst)%satzen(n)
  84          if (npredmax >= 7) iv%instid(inst)%varbc_info%pred(7,n) = iv%instid(inst)%satzen(n)**2
  85          if (npredmax >= 8) iv%instid(inst)%varbc_info%pred(8,n) = iv%instid(inst)%satzen(n)**3
  86       end if
  87
  88     end do                       ! pixel loop
  89   end do                         ! sensor loop
  90
  91    if (trace_use) call da_trace_exit("da_varbc_pred")
  92
  93 end subroutine da_varbc_pred