var/da/da_vtox_transforms/da_transform_xtoxa_adj.inc

   1 subroutine da_transform_xtoxa_adj(grid, wpec_on)
   2
   3    !--------------------------------------------------------------------------------------------------
   4    ! Purpose: Transfers fields from WRF TL fields to diagnostic fields needed by observational operators
   5    !          Adjoint
   6    !--------------------------------------------------------------------------------------------------
   7
   8    implicit none
   9
  10    type(domain),   intent(inout) :: grid
  11    logical, intent(in), optional :: wpec_on
  12
  13    integer :: i, j, k
  14    real    :: sdmd, s1md, mu
  15    real    :: p(kms:kme), mr_a(kms:kme), mr_b(kms:kme)
  16    real    :: PU, PD, coeff
  17    real    :: geoh(kms:kme)
  18    logical :: calc_wpec_term
  19
  20    if (trace_use) call da_trace_entry("da_transform_xtoxa_adj")
  21
  22    !-------------------------------------------------------------------------
  23    ! Polar treatment for Global
  24    !-------------------------------------------------------------------------
  25
  26    if (global) then
  27       ! Poles treatment for global WRFVAR
  28       call da_get_avpoles(grid%xa%u,grid%xa%v, &
  29          ids,ide,jds,jde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  30       call da_get_aspoles(grid%xa%t, &
  31          ids,ide,jds,jde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  32       call da_get_aspoles(grid%xa%p, &
  33          ids,ide,jds,jde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  34       call da_get_aspoles(grid%xa%q, &
  35          ids,ide,jds,jde,ims,ime,jms,jme,kms,kme,its,ite,jts,jte,kts,kte)
  36       call da_get_aspoles(grid%xa%psfc, &
  37          ids,ide,jds,jde,ims,ime,jms,jme,1,1,its,ite,jts,jte,1,1)
  38    end if
  39
  40    ! Compute w increments using Richardson's eqn.
  41    if ( Use_RadarObs)  then
  42       do k=kts,kte
  43          do j=jts,jte
  44             do i=its,ite
  45                grid%xa%w(i,j,k)=grid%xa%w(i,j,k)+0.5*grid%xa%wh(i,j,k)
  46                grid%xa%w(i,j,k+1)=grid%xa%w(i,j,k+1)+0.5*grid%xa%wh(i,j,k)
  47                grid%xa%wh(i,j,k)=0.0
  48             end do
  49          end do
  50       end do
  51       if ( .not. var4d .and. cloud_cv_options == 1 ) call da_uvprho_to_w_adj(grid)
  52    end if
  53
  54    if ( cloud_cv_options == 1 ) then
  55       ! Partition of hydrometeor increments via warm rain process
  56       call da_moist_phys_adj(grid)
  57    end if
  58
  59    !-------------------------------------------------------------------------
  60    ! If test_transforms = .true., not "XToY" transform needed to do here (YRG):
  61
  62    if (.not.test_transforms) then
  63       if (use_ssmt1obs .or. use_ssmt2obs .or. use_gpspwobs .or. use_gpsztdobs .or. &
  64           use_ssmitbobs .or. use_ssmiretrievalobs .or. use_gpsrefobs .or. use_gpsephobs) then
  65          if (use_ssmitbobs) call da_transform_xtotb_adj(grid)
  66          if (use_ssmt1obs .or. use_ssmt2obs .or. &
  67              use_ssmitbobs .or. use_ssmiretrievalobs) then
  68             if (global) then
  69                call da_error(__FILE__,__LINE__, &
  70                   (/"grid%xb%speed is not available, see da_transfer_kmatoxb.inc"/))
  71             end if
  72             call da_transform_xtoseasfcwind_adj(grid)
  73          end if
  74
  75          ! GPS Refractivity:
  76          if ( use_gpsrefObs .or. use_gpsephObs .or. use_gpsztdObs) then
  77             if (use_gpsztdObs) call da_transform_xtoztd_adj(grid)
  78             call da_transform_XTogpsref_adj(grid)
  79          endif
  80
  81          ! Now for PW.
  82          call da_transform_xtotpw_adj(grid)
  83       end if
  84
  85       if (sfc_assi_options == 2) call da_transform_xtowtq_adj(grid)
  86    end if
  87
  88    calc_wpec_term = .false.
  89    if ( present(wpec_on) ) then
  90       if ( wpec_on ) calc_wpec_term = .true.
  91    end if
  92
  93    if ( .not. var4d .and. calc_wpec_term ) then
  94
  95       grid%xa%mu=0.0
  96       grid%xa%rho=0.0
  97
  98       do j=jts,jte
  99          do i=its,ite
 100             if ((fg_format==fg_format_wrf_arw_regional) .or. &
 101                 (fg_format==fg_format_wrf_arw_global  )  ) then
 102                geoh(:)=0.0
 103                do k=kte,kts,-1
 104                   geoh(k) = geoh(k) + 0.5*grid%xa%geoh(i,j,k)
 105                   geoh(k+1) = geoh(k+1) + 0.5*grid%xa%geoh(i,j,k)
 106                   grid%xa%mu(i,j)=grid%xa%mu(i,j)-c1h(k)*geoh(k+1)*grid%xb%dnw(k)/grid%xb%rho(i,j,k)
 107                   grid%xa%rho(i,j,k)=grid%xa%rho(i,j,k)+ &
 108                      geoh(k+1)*(c1h(k)*grid%xb%psac(i,j)+c2h(k))*grid%xb%dnw(k)/(grid%xb%rho(i,j,k)**2)
 109                   geoh(k) = geoh(k) + geoh(k+1)
 110                end do
 111             end if
 112          end do
 113       end do
 114    end if
 115
 116    call da_pt_to_rho_adj(grid)
 117
 118    if ( .not. var4d ) then
 119
 120    do j=jts,jte
 121       do i=its,ite
 122          if ((fg_format==fg_format_wrf_arw_regional) .or. &
 123              (fg_format==fg_format_wrf_arw_global  )  ) then
 124             if ( .not. calc_wpec_term ) mu=0.0
 125             s1md=0.0
 126
 127             p(:)=0.0
 128
 129             do k=kts,kte
 130                mr_b(k) = grid%xb%q(i,j,k)/(1.0 - grid%xb%q(i,j,k))
 131                s1md=s1md+(1.0+mr_b(k))*grid%xb%dnw(k)
 132
 133                p(k) = p(k) + 0.5*grid%xa%p(i,j,k)
 134                p(k+1) = p(k+1) + 0.5*grid%xa%p(i,j,k)
 135
 136                if ( calc_wpec_term ) then
 137                   grid%xa%mu(i,j) = grid%xa%mu(i,j) - c1h(k)*p(k)*(1.0+mr_b(k))*grid%xb%dnw(k)
 138                   mu = grid%xa%mu(i,j)
 139                else
 140                   mu = mu - c1h(k)*p(k)*(1.0+mr_b(k))*grid%xb%dnw(k)
 141                end if
 142
 143                mr_a(k) = - p(k)*(c1h(k)*grid%xb%psac(i,j)+c2h(k))*grid%xb%dnw(k)
 144                p(k+1) = p(k+1) + p(k)
 145             end do
 146
 147             grid%xa%psfc(i,j) = grid%xa%psfc(i,j) - mu/s1md
 148             sdmd=-mu*grid%xb%psac(i,j)/s1md
 149
 150             do k=kts,kte
 151                mr_a(k) = mr_a(k) + sdmd*grid%xb%dnw(k)
 152                grid%xa%q(i,j,k) = grid%xa%q(i,j,k) + mr_a(k)/(1.0 - grid%xb%q(i,j,k))**2
 153             end do
 154          else if (fg_format == fg_format_kma_global)then
 155             do k=kts,kte
 156                if (k == kte) then
 157                   coeff = grid%xb%KMA_B(K)/(grid%xb%KMA_A(K)+grid%xb%KMA_B(K)* &
 158                      grid%xb%psfc(I,J)/100.0)
 159                else
 160                   PU = grid%xb%KMA_A(K+1) + grid%xb%KMA_B(K+1)*grid%xb%psfc(I,J)/100.0
 161                   PD = grid%xb%KMA_A(K ) + grid%xb%KMA_B(K )*grid%xb%psfc(I,J)/100.0
 162                   coeff=grid%xb%KMA_B(K)*1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU))+PD-PU)&
 163                      + grid%xb%KMA_B(K+1)*1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
 164                end if
 165
 166                grid%xa%psfc(i,j) = grid%xa % psfc(i,j) + &
 167                   grid%xb%p(i,j,k) * grid%xa % p(i,j,k)/100.0 * coeff
 168             end do
 169          end if
 170       end do
 171    end do
 172
 173    endif ! only for 3dvar
 174
 175    if (global) then
 176       call da_set_boundary_xa(grid)
 177    end if
 178
 179    if (trace_use) call da_trace_exit("da_transform_xtoxa_adj")
 180
 181 end subroutine da_transform_xtoxa_adj
 182
 183