var/da/da_vtox_transforms/da_transform_vptox.inc

   1 subroutine da_transform_vptox(grid, vp, be, ep)
   2 !subroutine da_transform_vptox(grid, vp, be, ep, nobwin)
   3
   4    !-----------------------------------------------------------------------
   5    ! Purpose: Physical transform of analysis increment variables.
   6    !    Updated for Analysis on Arakawa-C grid
   7    !    Author: Syed RH Rizvi,  MMM/ESSL/NCAR,  Date: 10/22/2008
   8    ! Updates:
   9    !
  10    !       Implementation of multi-variate BE for cv_options=6
  11    !       Syed RH Rizvi,  MMM/NESL/NCAR,  Date: 02/01/2010
  12    !------------------------
  13    !    Zhiquan (Jake) Liu, NCAR/MMM, 2015-09
  14    !     re-order transforms to avoid local chi_u and store full variables in vp
  15    !     full vp will be written out and used as input of inverse U transform
  16    !     for multi-resolution incremental 4DVAR
  17    !   order: v4 (rh), v3 (T), v5 (Ps), v2 (Chi_u -> Chi)
  18    !-----------------------------------------------------------------------
  19
  20    implicit none
  21
  22    type (domain), intent(inout)         :: grid
  23
  24    type (vp_type), intent(inout)        :: vp  ! CV on grid structure.
  25    type (be_type), intent(in), optional :: be  ! Background errors.
  26    type (ep_type), intent(in), optional :: ep  ! Ensemble perturbations.
  27 !   integer, intent(in), optional         :: nobwin
  28
  29    integer :: i, k, j, k1, ij            ! Loop counters.
  30    !real, allocatable          :: chi_u(:,:,:)  ! Unbalanced chi
  31
  32    if (trace_use) call da_trace_entry("da_transform_vptox")
  33
  34    !---------------------------------------------------------------------------
  35    !  [1] Add flow-dependent increments in control variable space (vp):
  36    !---------------------------------------------------------------------------
  37
  38    if (be % ne > 0 .and. alphacv_method == alphacv_method_vp) then
  39       call da_add_flow_dependence_vp(be % ne, ep, vp, its,ite, jts,jte, kts,kte)
  40    end if
  41
  42    !--------------------------------------------------------------------------
  43    ! [2] Impose statistical balance constraints:
  44    !--------------------------------------------------------------------------
  45
  46 if ( jb_factor > 0.0 ) then
  47    !$OMP PARALLEL DO &
  48    !$OMP PRIVATE ( ij, k1, k, j, i)
  49    do ij = 1 , grid%num_tiles
  50
  51    ! 2.1 Pseudo rh_u to Pseudo rh (only for cv6)
  52    ! do moisture first to avoid local (chi_u,t_t,Ps_u) variables
  53    !--------------------------------------------------------------
  54    if ( cv_options == 6 ) then
  55       do k1 = kts, kte
  56          do k = kts, kte
  57             do j = grid%j_start(ij), grid%j_end(ij)
  58                do i = its, ite
  59                   vp%v4(i,j,k1) = vp%v4(i,j,k1) + be%reg_psi_rh(j,k1,k)*vp%v1(i,j,k) + &
  60                   be%reg_chi_u_rh(j,k1,k)*vp%v2(i,j,k) + be%reg_t_u_rh(j,k1,k)*vp%v3(i,j,k)
  61                end do
  62             end do
  63          end do
  64       end do
  65 !
  66       do k = kts, kte
  67          do j = grid%j_start(ij), grid%j_end(ij)
  68             do i = its, ite
  69                vp%v4(i,j,k) = vp%v4(i,j,k) + be%reg_ps_u_rh(j,k)*vp%v5(i,j,1)
  70             end do
  71          end do
  72       end do
  73    end if
  74
  75    ! 2.2 t_u --> t, do this before chi_u --> chi
  76    !----------------------------------------------
  77    if (cv_options /= 7) then
  78       do k1 = kts, kte
  79          do k = kts, kte
  80             do j = grid%j_start(ij), grid%j_end(ij)
  81                do i = its, ite
  82                   vp%v3(i,j,k) = vp%v3(i,j,k) + be%reg_psi_t(j,k,k1)*vp%v1(i,j,k1)
  83                end do
  84             end do
  85          end do
  86       end do
  87    end if
  88
  89    if ( cv_options == 6 ) then
  90       do k1 = kts, kte
  91          do k = kts, kte
  92             do j = grid%j_start(ij), grid%j_end(ij)
  93                do i = its, ite
  94                   vp%v3(i,j,k) = vp%v3(i,j,k) + be%reg_chi_u_t(j,k,k1)*vp%v2(i,j,k1)
  95                end do
  96             end do
  97          end do
  98       end do
  99    end if
 100
 101    do k = kts, kte
 102       do j = grid%j_start(ij), grid%j_end(ij)
 103          do i = its, ite
 104             grid%xa%t(i,j,k) = vp%v3(i,j,k)
 105          end do
 106       end do
 107    end do
 108
 109    ! 2.3 Ps_u --> Ps, do this before chi_u --> chi
 110    !-------------------------------------------------
 111    if (cv_options /= 7) then
 112       do k = kts,kte
 113          do j = grid%j_start(ij), grid%j_end(ij)
 114             do i = its, ite
 115                vp%v5(i,j,1) = vp%v5(i,j,1) + be%reg_psi_ps(j,k)*vp%v1(i,j,k)
 116             end do
 117          end do
 118       end do
 119    end if
 120
 121    if ( cv_options == 6 ) then
 122       do k = kts,kte
 123          do j = grid%j_start(ij), grid%j_end(ij)
 124             do i = its, ite
 125                vp%v5(i,j,1) = vp%v5(i,j,1) + be%reg_chi_u_ps(j,k)*vp%v2(i,j,k)
 126             end do
 127          end do
 128       end do
 129    end if
 130
 131    do j = grid%j_start(ij), grid%j_end(ij)
 132       do i = its, ite
 133          grid%xa%psfc(i,j) = vp%v5(i,j,1)
 134       end do
 135    end do
 136
 137    ! 2.4 Chi_u --> Chi, do this last
 138    !-----------------------------------
 139    if (cv_options /= 7) then
 140       do k = kts, kte
 141          do j = grid%j_start(ij), grid%j_end(ij)
 142             do i = its, ite
 143                vp%v2(i,j,k) = vp%v2(i,j,k) + be%reg_psi_chi(j,k)* vp%v1(i,j,k)
 144             end do
 145          end do
 146       end do
 147    end if
 148 !
 149 !   if ( cv_options == 6 ) deallocate (chi_u )
 150
 151    end do
 152    !$OMP END PARALLEL DO
 153    !--------------------------------------------------------------------------
 154    ! [3] Transform to model variable space:
 155    !--------------------------------------------------------------------------
 156
 157 #ifdef A2C
 158   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 159        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 160      ipe = ipe + 1
 161      ide = ide + 1
 162   end if
 163
 164   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 165        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 166      jpe = jpe + 1
 167      jde = jde + 1
 168   end if
 169 #endif
 170
 171 #ifdef DM_PARALLEL
 172 #include "HALO_PSICHI_UV.inc"
 173 #endif
 174
 175 #ifdef A2C
 176   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 177        fg_format==fg_format_wrf_arw_global  ) .and. ide == ipe ) then
 178      ipe = ipe - 1
 179      ide = ide - 1
 180   end if
 181
 182   if ((fg_format==fg_format_wrf_arw_regional  .or. &
 183        fg_format==fg_format_wrf_arw_global  ) .and. jde == jpe ) then
 184      jpe = jpe - 1
 185      jde = jde - 1
 186   end if
 187 #endif
 188
 189    ! Psi and chi to u and v:
 190    if ( cv_options == 5 .or. cv_options == 6 ) then
 191       call da_psichi_to_uv(vp % v1, vp % v2, grid%xb % coefx, &
 192            grid%xb % coefy , grid%xa % u, grid%xa % v)
 193    else if ( cv_options == 7 ) then
 194       grid%xa%u = vp%v1
 195       grid%xa%v = vp%v2
 196    end if
 197
 198    if ( cloud_cv_options /= 1 ) then
 199       ! Pseudo RH --> Water vapor mixing ratio:
 200       !$OMP PARALLEL DO &
 201       !$OMP PRIVATE ( ij, i, j, k )
 202       do ij = 1 , grid%num_tiles
 203          do k = kts, kte
 204             do j = grid%j_start(ij), grid%j_end(ij)
 205                do i = its, ite
 206                   grid%xa % q(i,j,k) =  vp%v4(i,j,k) * grid%xb%qs(i,j,k)
 207                enddo
 208             enddo
 209          enddo
 210       enddo
 211       !$OMP END PARALLEL DO
 212    else ! cloud_cv_options=1
 213       ! Pseudo RH --> Total water mixing ratio:
 214       !$OMP PARALLEL DO &
 215       !$OMP PRIVATE ( ij, i, j, k )
 216       do ij = 1 , grid%num_tiles
 217          do k = kts, kte
 218             do j = grid%j_start(ij), grid%j_end(ij)
 219                do i = its, ite
 220                  grid%xa % qt(i,j,k) = vp%v4(i,j,k) * grid%xb%qs(i,j,k)
 221                enddo
 222             enddo
 223          enddo
 224       enddo
 225       !$OMP END PARALLEL DO
 226    end if
 227
 228    if ( cloud_cv_options >= 2 ) then
 229       !qcloud
 230       !$OMP PARALLEL DO &
 231       !$OMP PRIVATE ( ij, i, j, k )
 232       do ij = 1 , grid%num_tiles
 233          do k = kts, kte
 234             do j = grid%j_start(ij), grid%j_end(ij)
 235                do i = its, ite
 236                   grid%xa % qcw(i,j,k) =  vp%v6(i,j,k)
 237                enddo
 238             enddo
 239          enddo
 240       enddo
 241       !$OMP END PARALLEL DO
 242       !qrain
 243       !$OMP PARALLEL DO &
 244       !$OMP PRIVATE ( ij, i, j, k )
 245       do ij = 1 , grid%num_tiles
 246          do k = kts, kte
 247             do j = grid%j_start(ij), grid%j_end(ij)
 248                do i = its, ite
 249                   grid%xa % qrn(i,j,k) =  vp%v7(i,j,k)
 250                enddo
 251             enddo
 252          enddo
 253       enddo
 254       !$OMP END PARALLEL DO
 255       !qice
 256       !$OMP PARALLEL DO &
 257       !$OMP PRIVATE ( ij, i, j, k )
 258       do ij = 1 , grid%num_tiles
 259          do k = kts, kte
 260             do j = grid%j_start(ij), grid%j_end(ij)
 261                do i = its, ite
 262                   grid%xa % qci(i,j,k) =  vp%v8(i,j,k)
 263                enddo
 264             enddo
 265          enddo
 266       enddo
 267       !$OMP END PARALLEL DO
 268       !qsnow
 269       !$OMP PARALLEL DO &
 270       !$OMP PRIVATE ( ij, i, j, k )
 271       do ij = 1 , grid%num_tiles
 272          do k = kts, kte
 273             do j = grid%j_start(ij), grid%j_end(ij)
 274                do i = its, ite
 275                   grid%xa % qsn(i,j,k) =  vp%v9(i,j,k)
 276                enddo
 277             enddo
 278          enddo
 279       enddo
 280       !$OMP END PARALLEL DO
 281       !qgraupel
 282       !$OMP PARALLEL DO &
 283       !$OMP PRIVATE ( ij, i, j, k )
 284       do ij = 1 , grid%num_tiles
 285          do k = kts, kte
 286             do j = grid%j_start(ij), grid%j_end(ij)
 287                do i = its, ite
 288                   grid%xa % qgr(i,j,k) =  vp%v10(i,j,k)
 289                enddo
 290             enddo
 291          enddo
 292       enddo
 293       !$OMP END PARALLEL DO
 294    end if ! cloud_cv_options>=2
 295
 296    if ( use_cv_w ) then
 297       !vertical velocity
 298       !$OMP PARALLEL DO &
 299       !$OMP PRIVATE ( ij, i, j, k )
 300       do ij = 1 , grid%num_tiles
 301          do k = kts, kte
 302             do j = grid%j_start(ij), grid%j_end(ij)
 303                do i = its, ite
 304                   grid%xa % w(i,j,k) =  vp%v11(i,j,k)
 305                enddo
 306             enddo
 307          enddo
 308       enddo
 309       !$OMP END PARALLEL DO
 310    end if ! use_cv_w
 311 end if ! jb_factor > 0.0
 312
 313 !   !---------------------------------------------------------------------------
 314 !   !  [4] Add flow-dependent increments in model space (grid%xa):
 315 !   !---------------------------------------------------------------------------
 316 !
 317 !!  if (be % ne > 0 .and. alphacv_method == alphacv_method_xa) then
 318 !!     call da_add_flow_dependence_xa(grid, be % ne, ep, vp)
 319 !!  end if
 320 !   if (be % ne > 0 .and. alphacv_method == alphacv_method_xa) then
 321 !      if ( anal_type_hybrid_dual_res ) then
 322 !         if( present(nobwin) ) then
 323 !            call da_add_flow_dependence_xa_dual_res(grid, be % ne, ep, vp, nobwin)
 324 !         else
 325 !            call da_add_flow_dependence_xa_dual_res(grid, be % ne, ep, vp)
 326 !         end if
 327 !      else
 328 !         if( present(nobwin) ) then
 329 !            call da_add_flow_dependence_xa(grid, be % ne, ep, vp, nobwin)
 330 !         else
 331 !            call da_add_flow_dependence_xa(grid, be % ne, ep, vp)
 332 !         end if
 333 !      endif
 334 !   end if
 335
 336    if (trace_use) call da_trace_exit("da_transform_vptox")
 337
 338 end subroutine da_transform_vptox
 339