da_transform_xtovp.inc
References to this file elsewhere.
1 subroutine da_transform_xtovp(xb, xbx, xa, xp, vp, be)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Transforms analysis to control variables (Vp-type)
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (xb_type), intent(in) :: xb ! First guess structure.
10 type (xbx_type),intent(in) :: xbx ! Header/non-gridded vars.
11 type (x_type), intent(inout) :: xa ! Analysis increments.
12 type (xpose_type), intent(inout) :: xp ! Dimensions and xpose buffers.
13 type (vp_type), intent(out) :: vp ! CV on grid structure.
14 type (be_type), intent(in), optional :: be ! Background errors.
15
16 real, dimension(ims:ime,jms:jme,kms:kme) :: vor, & ! Vorticity.
17 div ! Divergence.
18
19 real, dimension(ims:ime,jms:jme) :: one_over_m2 ! Multiplicative coeff.
20
21 integer :: i, j, k ! Loop counters.
22
23 !----------------------------------------------------------------
24 ! [1.0] Perform transform v = U^{-1} x
25 !----------------------------------------------------------------
26
27 call da_zero_vp_type (vp)
28
29 ! [2.2] Transform u, v to streamfunction via vorticity:
30
31 ! Communicate halo region.
32 call wrf_dm_halo(xp%domdesc,xp%comms,xp%halo_id3)
33
34 call da_uv_to_vorticity(xb, xa % u, xa % v, vor)
35
36 ! Convert vorticity to Del**2 psi:
37 if (.not. global) then
38 if (fg_format == fg_format_wrf) then
39 one_over_m2(its:ite,jts:jte) = 1.0 / (xb % map_factor(its:ite,jts:jte) * &
40 xb % map_factor(its:ite,jts:jte))
41 do k = kts, kte
42 vor(its:ite,jts:jte,k) = &
43 one_over_m2(its:ite,jts:jte)*vor(its:ite,jts:jte,k)
44 end do
45 else
46 write(unit=message(1),fmt='(A,I5,A,L10)') &
47 ' Wrong choice of fg_format= ',fg_format,' with global = ',global
48 call da_error(__FILE__,__LINE__,message(1:1))
49 end if
50 end if
51
52 ! Calculate psi:
53 write (unit=stdout,fmt=*) ' calling Solve_PoissonEquation for Psi'
54 call da_solve_poissoneqn_fct(xbx, vor, vp%v1, xp)
55
56 ! [2.3] Transform u, v to velocity potential via divergence:
57
58 call da_message((/'calling UV_To_Divergence'/))
59 call da_uv_to_divergence(xb, xa % u, xa % v, div)
60
61 ! Convert divergence to Del**2 chi:
62 if (.not. global) then
63 if (fg_format == fg_format_wrf) then
64 do k = kts, kte
65 div(its:ite,jts:jte,k) = &
66 one_over_m2(its:ite,jts:jte) * div(its:ite,jts:jte,k)
67 end do
68 else
69 write(unit=message(1),fmt='(A,I5,A,L10)') &
70 ' Wrong choice of fg_format= ',fg_format,' with global = ',global
71 call da_error(__FILE__,__LINE__,message(1:1))
72 end if
73 end if
74
75 ! Calculate chi:
76
77 call da_message((/' calling Solve_PoissonEquation for Chi'/))
78 call da_solve_poissoneqn_fct(xbx, div, vp%v2, xp)
79
80 ! [2.4] Transform chi to chi_u:
81 call da_message((/' calculating chi_u'/))
82 do k=kts,kte
83 do j=jts,jte
84 vp%v2(its:ite,j,k) = vp%v2(its:ite,j,k) - &
85 be%reg_chi(j,k)*vp%v1(its:ite,j,k)
86 end do
87 end do
88
89 call da_message((/' calculating t_u'/))
90 ! [2.4] Compute t_u:
91 do k=kts,kte
92 do j=jts,jte
93 do i=its,ite
94 vp%v3(i,j,k) = xa%t(i,j,k) - &
95 sum(be%reg_t(j,k,kts:kte)*vp%v1(i,j,kts:kte))
96 end do
97 end do
98 end do
99
100 ! [2.6] Choice of moisture control variable:
101
102 call da_message((/' calculating psudo rh'/))
103 vp % v4(its:ite,jts:jte,kts:kte) = xa % q (its:ite,jts:jte,kts:kte) / &
104 xb % qs (its:ite,jts:jte,kts:kte)
105
106 ! [2.7] Choice of surface pressure control variable:
107
108 ! [2.7] compute psfc_u:
109 call da_message((/' calculating psfc_u '/))
110 do j=jts,jte
111 do i=its,ite
112 vp % v5(i,j,1) = xa%psfc(i,j) - &
113 sum(be%reg_ps(j,kts:kte)*vp%v1(i,j,kts:kte))
114 end do
115 end do
116
117 end subroutine da_transform_xtovp
118
119