var/da/da_gpspw/da_transform_xtoy_gpsztd_adj.inc

   1 subroutine da_transform_xtoy_gpsztd_adj( grid, iv, jo_grad_y, jo_grad_x )
   2 !----------------------------------------------------------------
   3 !  HISTORY
   4 !
   5 !    Purpose:  Considering the difference in the elevation
   6 !              between model surface and GPS ZTD site
   7 !
   8 !                                   Y.-R. Guo 05/21/2008
   9 !----------------------------------------------------------------
  10
  11    IMPLICIT NONE
  12
  13    type (domain),  intent(in)    :: grid
  14    type (iv_type), intent(in)    :: iv          ! obs. inc vector (o-b).
  15    type (y_type) , intent(in)    :: jo_grad_y   ! grad_y(jo)
  16    type (x_type) , intent(inout) :: jo_grad_x   ! grad_x(jo)
  17
  18    INTEGER                       :: n        ! Loop counter.
  19    INTEGER                       :: i, j     ! Index dimension.
  20    REAL                          :: dx, dxm  !
  21    REAL                          :: dy, dym  !
  22
  23 !-- 2004-04-08
  24    INTEGER                :: k        ! Index dimension.
  25    REAL                   :: dzd,ddzd     ! adjustment pw [mm]
  26    REAL                   :: obs_terr     ! real terrain height at GPS site [m]
  27    REAL                   :: model_terr   ! model terrain height at GPS site[m]
  28    REAL,DIMENSION(kts:kte):: model_z     ! model z at GPS site [m]
  29    REAL,DIMENSION(kts:kte):: model_dztd  ! model dq at GPS site [kg/kg]
  30 !--
  31    if (trace_use_dull) call da_trace_entry("da_transform_xtoy_gpspw_adj")
  32
  33    do n=iv%info(gpspw)%n1,iv%info(gpspw)%n2
  34
  35          i   = iv%info(gpspw)%i(1,n)
  36          dy  = iv%info(gpspw)%dy(1,n)
  37          dym = iv%info(gpspw)%dym(1,n)
  38          j   = iv%info(gpspw)%j(1,n)
  39          dx  = iv%info(gpspw)%dx(1,n)
  40          dxm = iv%info(gpspw)%dxm(1,n)
  41
  42 !  Initialise the varibles
  43          dzd             = 0.
  44          model_dztd(:)   = 0.
  45
  46          obs_terr   = iv%info(gpspw)%elv(n)
  47          model_terr = dym*(dxm*grid%xb%terr(i,j)   + dx*grid%xb%terr(i+1,j)) + &
  48                       dy *(dxm*grid%xb%terr(i,j+1) + dx*grid%xb%terr(i+1,j+1))
  49
  50          dzd =       1.e-4 * jo_grad_y%gpspw(n)%tpw
  51
  52          jo_grad_x%ztd(i  ,j  ) = jo_grad_x%ztd(i  ,j  ) &
  53                                 + dxm*dym*jo_grad_y%gpspw(n)%tpw
  54          jo_grad_x%ztd(i+1,j  ) = jo_grad_x%ztd(i+1,j  ) &
  55                                 + dym*dx *jo_grad_y%gpspw(n)%tpw
  56          jo_grad_x%ztd(i  ,j+1) = jo_grad_x%ztd(i  ,j+1) &
  57                                 + dxm *dy*jo_grad_y%gpspw(n)%tpw
  58          jo_grad_x%ztd(i+1,j+1) = jo_grad_x%ztd(i+1,j+1) &
  59                                 + dx *dy *jo_grad_y%gpspw(n)%tpw
  60
  61          IF ( obs_terr <= model_terr ) THEN
  62
  63             model_dztd(1)  =  (obs_terr - model_terr) * dzd
  64
  65             jo_grad_x%ref(i,j,1)     = jo_grad_x%ref(i,j,1)     + dym*dxm*model_dztd(1)
  66             jo_grad_x%ref(i+1,j,1)   = jo_grad_x%ref(i+1,j,1)   + dym*dx *model_dztd(1)
  67             jo_grad_x%ref(i,j+1,1)   = jo_grad_x%ref(i,j+1,1)   + dy *dxm*model_dztd(1)
  68             jo_grad_x%ref(i+1,j+1,1) = jo_grad_x%ref(i+1,j+1,1) + dy *dx *model_dztd(1)
  69
  70          ELSE
  71
  72             model_z(1) = dym*(dxm*grid%xb%hf(i,j,1)   + dx*grid%xb%hf(i+1,j,1)) + &
  73                          dy *(dxm*grid%xb%hf(i,j+1,1) + dx*grid%xb%hf(i+1,j+1,1))
  74
  75 ! ..............................................................................
  76 ! The following part of code is written by Y.-R. Guo             07/16/2004
  77
  78             do k = kts, kte
  79
  80                if ( model_z(k) >= obs_terr ) exit
  81
  82                model_z(k+1) = dym*(dxm*grid%xb%hf(i,j,k+1)   + dx*grid%xb%hf(i+1,j,k+1)) + &
  83                               dy *(dxm*grid%xb%hf(i,j+1,k+1) + dx*grid%xb%hf(i+1,j+1,k+1))
  84
  85                ddzd = dzd
  86
  87                if ( model_z(k+1) <= obs_terr ) then
  88                  model_dztd(k) = (model_z(k+1) - model_z(k))*ddzd
  89                else
  90                  model_dztd(k) = (obs_terr - model_z(k))*ddzd
  91                end if
  92 !
  93 ! No feedback to x%hf was considered here (Refer to comments in
  94 ! DA_Transform_XToY_Gpspw.inc).       Y.-R. Guo  04/06/2005
  95 ! .....................................................................................
  96                jo_grad_x%ref(i,j,k)     = jo_grad_x%ref(i,j,k)     + dym*dxm*model_dztd(k)
  97                jo_grad_x%ref(i+1,j,k)   = jo_grad_x%ref(i+1,j,k)   + dym*dx *model_dztd(k)
  98                jo_grad_x%ref(i,j+1,k)   = jo_grad_x%ref(i,j+1,k)   + dy *dxm*model_dztd(k)
  99                jo_grad_x%ref(i+1,j+1,k) = jo_grad_x%ref(i+1,j+1,k) + dy *dx *model_dztd(k)
 100
 101              end do
 102          end if
 103       end do
 104
 105    if (trace_use_dull) call da_trace_exit("da_transform_xtoy_gpsztd_adj")
 106
 107 end subroutine da_transform_xtoy_gpsztd_adj