var/da/da_vtox_transforms/da_transform_xtoxa.inc

   1 subroutine da_transform_xtoxa(grid, wpec_on)
   2
   3    !--------------------------------------------------------------------------------------------------
   4    ! Purpose: Transfers fields from WRF TL fields to diagnostic fields needed by observational operators
   5    !    Updated for Analysis on Arakawa-C grid
   6    !    Author: Syed RH Rizvi,  MMM/ESSL/NCAR,  Date: 10/22/2008
   7    !--------------------------------------------------------------------------------------------------
   8
   9    implicit none
  10
  11    type(domain),   intent(inout) :: grid
  12    logical, intent(in), optional :: wpec_on
  13
  14    integer :: i, j, k
  15    real    :: sdmd, s1md, mu
  16    real    :: p(kms:kme)
  17    real    :: mr_a(kms:kme)
  18    real    :: mr_b(kms:kme)
  19    real    :: PU, PD, coeff
  20    real    :: geoh(kms:kme)
  21    logical :: calc_wpec_term
  22
  23    if (trace_use) call da_trace_entry("da_transform_xtoxa")
  24
  25    !----------------------------------------------------------------------
  26    ! [4.0]: Move the following:
  27    !----------------------------------------------------------------------
  28
  29    calc_wpec_term = .false.
  30    if ( present(wpec_on) ) then
  31       if ( wpec_on ) calc_wpec_term =.true.
  32    end if
  33
  34    if ( .not. var4d ) then ! for 4dvar, xa%p comes from the model
  35
  36    do j=jts,jte
  37       do i=its,ite
  38         if ((fg_format==fg_format_wrf_arw_regional) .or. &
  39             (fg_format==fg_format_wrf_arw_global  ) ) then
  40
  41             sdmd=0.0
  42             s1md=0.0
  43             do k=kts,kte
  44                mr_a(k) = grid%xa%q(i,j,k)/(1.0 - grid%xb%q(i,j,k))**2
  45                mr_b(k) = grid%xb%q(i,j,k)/(1.0 - grid%xb%q(i,j,k))
  46
  47                sdmd=sdmd+mr_a(k)*grid%xb%dnw(k)
  48                s1md=s1md+(1.0+mr_b(k))*grid%xb%dnw(k)
  49             end do
  50
  51             mu=-(grid%xa%psfc(i,j)+grid%xb%psac(i,j)*sdmd)/s1md
  52             if ( calc_wpec_term ) then
  53                grid%xa%mu(i,j) = mu
  54             end if
  55
  56             p(kte+1)=0.0
  57
  58             do k=kte,kts,-1
  59                p(k)=p(k+1)-(c1h(k)*mu*(1.0+mr_b(k)) &
  60                        + (c1h(k)*grid%xb%psac(i,j)+c2h(k))*mr_a(k))*grid%xb%dnw(k)
  61
  62                grid%xa%p(i,j,k)=0.5*(p(k)+p(k+1))
  63             end do
  64          else if (fg_format == fg_format_kma_global) then
  65             do k=kts,kte
  66                if (k == kte) then
  67                   coeff=grid%xb%KMA_B(K)/(grid%xb%KMA_A(K)+grid%xb%KMA_B(K)*grid%xb%psfc(I,J)/100.0)
  68                else
  69                   PU = grid%xb%KMA_A(K+1) + grid%xb%KMA_B(K+1)*grid%xb%psfc(I,J)/100.0
  70                   PD = grid%xb%KMA_A(K ) + grid%xb%KMA_B(K )*grid%xb%psfc(I,J)/100.0
  71                   coeff=grid%xb%KMA_B(K)  *1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU)) &
  72                     + PD-PU)&
  73                     + grid%xb%KMA_B(K+1)*1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
  74                end if
  75                ! Here since grid%xa%psfc holds value in Pa. dlnp -> dp
  76                grid%xa%p(i,j,k) =  grid%xb%p(i,j,k) * grid%xa%psfc(I,J)/100.0 * coeff
  77
  78             end do
  79          end if
  80       end do
  81    end do
  82
  83    endif ! only for 3dvar
  84
  85    call da_pt_to_rho_lin(grid)
  86
  87    if ( .not. var4d .and. calc_wpec_term ) then
  88       do j=jts,jte
  89          do i=its,ite
  90             if ((fg_format==fg_format_wrf_arw_regional) .or. &
  91                (fg_format==fg_format_wrf_arw_global  ) ) then
  92                geoh(kts)=0.0
  93                do k=kts,kte
  94                   geoh(k+1)=geoh(k)+(-c1h(k)*grid%xa%mu(i,j) &
  95                           + grid%xa%rho(i,j,k)*(c1h(k)*grid%xb%psac(i,j)+c2h(k)) &
  96                           /grid%xb%rho(i,j,k))* grid%xb%dnw(k)/grid%xb%rho(i,j,k)
  97                   grid%xa%geoh(i,j,k)=0.5*(geoh(k)+geoh(k+1))
  98                end do
  99             end if
 100          end do
 101       end do
 102    end if
 103
 104 #ifdef A2C
 105   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 106        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 107      ipe = ipe + 1
 108      ide = ide + 1
 109   end if
 110
 111   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 112        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 113      jpe = jpe + 1
 114      jde = jde + 1
 115   end if
 116 #endif
 117 #ifdef DM_PARALLEL
 118 #include "HALO_XA.inc"
 119 #endif
 120
 121   if ( .not. var4d .and. calc_wpec_term ) then
 122 #ifdef DM_PARALLEL
 123 #include "HALO_XA_WPEC.inc"
 124 #include "HALO_XB_WPEC.inc"
 125 #endif
 126   end if
 127
 128 #ifdef A2C
 129   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 130        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 131      ipe = ipe - 1
 132      ide = ide - 1
 133   end if
 134
 135   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 136        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 137      jpe = jpe - 1
 138      jde = jde - 1
 139   end if
 140 #endif
 141
 142    ! If test_transforms = .true., not "XToY" transform needed to do here:
 143
 144    if (.not.test_transforms) then
 145       ! Exchange grid%xa halo region.
 146
 147       if (sfc_assi_options == 2) then
 148          call da_transform_xtowtq (grid)
 149          ! Exchange grid%xa (surface variable) halo region.
 150 #ifdef DM_PARALLEL
 151 #include "HALO_SFC_XA.inc"
 152 #endif
 153       end if
 154
 155       if (use_ssmt1obs .or. use_ssmt2obs .or. use_gpspwobs .or. use_gpsztdobs .or. &
 156           use_ssmitbobs .or. use_ssmiretrievalobs .or. use_gpsrefobs .or. use_gpsephobs) then
 157
 158          ! Now do something for PW
 159          call da_transform_xtotpw(grid)
 160
 161          ! Space-based GPS Refractivity, GPS Excess Phase, and Ground-based GPS ZTD:
 162          if ( use_gpsrefObs .or. use_gpsephObs .or. use_gpsztdObs) then
 163             call da_transform_xtogpsref_lin(grid)
 164             if (use_GpsztdObs) call da_transform_xtoztd_lin(grid)
 165          endif
 166
 167          if (use_ssmt1obs .or. use_ssmt2obs .or. &
 168               use_ssmitbobs .or. use_ssmiretrievalobs) then
 169             if (global) then
 170                call da_error(__FILE__,__LINE__, &
 171                   (/"grid%xb%speed is not available, see da_transfer_kmatoxb.inc"/))
 172             end if
 173             call da_transform_xtoseasfcwind_lin(grid)
 174          end if
 175          if (use_ssmitbobs) call da_transform_xtotb_lin (grid)
 176
 177          ! Exchange grid%xa halo region.
 178 #ifdef DM_PARALLEL
 179 #include "HALO_SSMI_XA.inc"
 180 #endif
 181       end if
 182    end if
 183
 184    ! Compute w increments using Richardson's eqn.
 185
 186    if ( Use_RadarObs ) then
 187       if ( .not. var4d .and. cloud_cv_options == 1) call da_uvprho_to_w_lin(grid)
 188
 189       do k=kts,kte
 190          do j=jts,jte
 191             do i=its,ite
 192                grid%xa%wh(i,j,k)=0.5*(grid%xa%w(i,j,k)+grid%xa%w(i,j,k+1))
 193             end do
 194          end do
 195       end do
 196
 197 #ifdef DM_PARALLEL
 198 #include "HALO_RADAR_XA_W.inc"
 199 #endif
 200    end if
 201
 202    if ( cloud_cv_options == 1 )then
 203       ! Partition of hydrometeor increments via warm rain process
 204       call da_moist_phys_lin(grid)
 205    end if
 206
 207    !---------------------------------------------------------------
 208    ! Polar treatment for Global
 209    !---------------------------------------------------------------
 210
 211    if (global)  then
 212       call da_get_vpoles(grid%xa%u,grid%xa%v, &
 213          ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 214       call da_get_spoles(grid%xa%t, &
 215          ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 216       call da_get_spoles(grid%xa%p, &
 217          ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 218       call da_get_spoles(grid%xa%q, &
 219          ids, ide, jds, jde, ims, ime, jms, jme, kms, kme, its, ite, jts, jte, kts, kte)
 220       call da_get_spoles(grid%xa%psfc, &
 221          ids, ide, jds, jde, ims, ime, jms, jme,   1,   1, its, ite, jts, jte,   1,   1)
 222       call da_set_boundary_xa(grid)
 223    end if
 224
 225    if (trace_use) call da_trace_exit("da_transform_xtoxa")
 226
 227 end subroutine da_transform_xtoxa
 228
 229