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