wrf svn trunk commit r4103
[wrffire.git] / wrfv2_fire / var / da / da_vtox_transforms / da_transform_bal_adj.inc
blob5b0d04ef7da15a86043748be77d567f1f8ae4eb0
1 SUBROUTINE da_transform_bal_adj( vp, be, grid )
3    IMPLICIT NONE
5    TYPE (vp_type), INTENT(INOUT)        :: vp ! CV on grid structure.out
6    TYPE (be_type), INTENT(IN)           :: be ! Background errors.
7    type (domain) , intent(inout) :: grid   ! Dimensions and xpose buffers.
9    INTEGER                              :: i, j, k, kk, ij, iunit  ! Loop counters.
10 !-------------------------------------------------------------------
11 !  [1.0] Initialise:
12 !------------------------------------------------------------------- 
13   
14 !  linear balance btw psi and t-b, Psfc_b and chi_b 
17 !  [3.4]  adj of Transform psi and chi to u and v:        
19 #ifdef DM_PARALLEL
20 #include "HALO_PSICHI_UV_ADJ.inc"
21 #endif
23    !$OMP PARALLEL DO &
24    !$OMP PRIVATE ( i, j, k, ij )
25    DO ij = 1 , grid%num_tiles
26       DO k = kts,kte
27          DO j = grid%j_start(ij), grid%j_end(ij)
28             DO i= its,ite
29                vp%v1(i,j,k)=0.
30                vp%v2(i,j,k)=0.
31             END DO
32          END DO
33       END DO
34    ENDDO
35    !$OMP END PARALLEL DO
37    call da_psichi_to_uv_adj( grid%xa % u, grid%xa % v, grid%xb % coefx,   &
38                             grid%xb % coefy,  vp % v1, vp % v2  )
40 !  [3.3]  adj Calculate Psfc_b from psi
42 !--convert from delt.ps to delt.ln(ps)
43    !$OMP PARALLEL DO &
44    !$OMP PRIVATE ( i, j, ij )
45    DO ij = 1 , grid%num_tiles
46       DO j = grid%j_start(ij), grid%j_end(ij)
47          DO i= its,ite
48             grid%xa%psfc(i,j) = grid%xa%psfc(i,j) * grid%xb%psfc(i,j)
49          END DO
50       END DO
51    END DO
52    !$OMP END PARALLEL DO
54    !$OMP PARALLEL DO &
55    !$OMP PRIVATE ( i, j, k, ij )
56    DO ij = 1 , grid%num_tiles
57       DO k = kts,kte
58          DO j = grid%j_start(ij), grid%j_end(ij)
59             DO i= its,ite
60                vp % v1(i,j,k)= vp % v1(i,j,k)+ &
61                                 be % wgvz(i,j,k) * grid%xa % psfc(i,j)
62             END DO
63          END DO
64       END DO
65    END DO
66    !$OMP END PARALLEL DO
68    !$OMP PARALLEL DO &
69    !$OMP PRIVATE ( i, j, ij )
70    DO ij = 1 , grid%num_tiles
71       DO j = grid%j_start(ij), grid%j_end(ij)
72          DO i= its,ite
73             vp % v5(i,j,1)=grid%xa % psfc(i,j)
74          END DO
75       END DO
76    END DO
77    !$OMP END PARALLEL DO
79 !  [3.2] adj Calculate chi_b from psi
81    !$OMP PARALLEL DO &
82    !$OMP PRIVATE ( i, j, k, ij )
83    DO ij = 1 , grid%num_tiles
84       DO k = kts,kte
85          DO j = grid%j_start(ij), grid%j_end(ij)
86             DO i= its,ite
87                vp % v1(i,j,k)= vp % v1(i,j,k)+ &
88                          be % bvz(i,j,k) * vp % v2(i,j,k) 
89             END DO
90          END DO
91       END DO
92    END DO
93    !$OMP END PARALLEL DO
95 !  [3.1] Calculate t_b from psi
97    !$OMP PARALLEL DO &
98    !$OMP PRIVATE ( i, j, k, kk, ij )
99    DO ij = 1 , grid%num_tiles
100       DO kk = kts,kte
101          DO k = kts,kte
102             DO j = grid%j_start(ij), grid%j_end(ij)
103                DO i= its,ite
104                   vp % v1(i,j,kk)= vp % v1(i,j,kk)+ &
105                           be % agvz(i,j,k,kk) * grid%xa % t(i,j,k)     
106                END DO
107             END DO
108          END DO
109       END DO
110    END DO
111    !$OMP END PARALLEL DO
113    !$OMP PARALLEL DO &
114    !$OMP PRIVATE ( i, j, k, ij )
115    DO ij = 1 , grid%num_tiles
116       DO k = kts,kte
117          DO j = grid%j_start(ij), grid%j_end(ij)
118             DO i= its,ite
119                vp % v3(i,j,k)= grid%xa % t(i,j,k)
120             END DO
121          END DO
122       END DO
123    END DO
124    !$OMP END PARALLEL DO
126 !  [3.5] treat humidity                         
128    IF ( cv_options == 3 ) THEN
129    !$OMP PARALLEL DO &
130    !$OMP PRIVATE ( i, j, k, ij )
131    DO ij = 1 , grid%num_tiles
132       DO k = kts,kte
133          DO j = grid%j_start(ij), grid%j_end(ij)
134             DO i= its,ite
135                vp % v4(i,j,k) = grid%xa % q(i,j,k) * grid%xb % qs(i,j,k)
136             END DO
137          END DO
138       END DO
139    END DO
140    !$OMP END PARALLEL DO
142    ELSE IF ( cv_options_hum == 1 ) THEN
144       vp % v4(its:ite,jts:jte,kts:kte) = vp % v4(its:ite,jts:jte,kts:kte)+&
145                                          grid%xa % q(its:ite,jts:jte,kts:kte)
147    ELSE IF ( cv_options_hum == 2 ) THEN
149       CALL DA_TPRH_To_Q_Adj( grid )
151       vp % v4(its:ite,jts:jte,kts:kte) = vp % v4(its:ite,jts:jte,kts:kte)+&
152                                          grid%xa % rh(its:ite,jts:jte,kts:kte)
154    END IF
156 END SUBROUTINE da_transform_bal_adj