da_transfer_xatokma.inc
References to this file elsewhere.
1 subroutine da_transfer_xatokma(grid)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Convert analysis increments into KMA increments
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type(domain), intent(inout) :: grid
10
11 integer :: i, j, k
12 real :: PU, PD, coeff
13
14 if (trace_use) call da_trace_entry("da_transfer_xatokma")
15
16 !---------------------------------------------------------------------------
17 ! Add increment to the original guess and update xb and "grid"
18 !---------------------------------------------------------------------------
19
20 do j=jts,jte
21 do i=its,ite
22 grid%xb%w(i,j,kte+1)= grid%xb%w(i,j,kte+1) + grid%xa%w(i,j,kte+1)
23 end do
24 do i=its,ite
25 do k = kts, kte
26 grid%xb%u(i,j,k) = grid%xa%u(i,j,k) + grid%xb%u(i,j,k)
27 grid%xb%v(i,j,k) = grid%xa%v(i,j,k) + grid%xb%v(i,j,k)
28 grid%xb%t(i,j,k) = grid%xa%t(i,j,k) + grid%xb%t(i,j,k)
29 grid%xb%w(i,j,k) = grid%xa%w(i,j,k) + grid%xb%w(i,j,k)
30 grid%xb%q(i,j,k) = grid%xa%q(i,j,k) + grid%xb%q(i,j,k)
31 ! compute pressure increments at KMA full levels
32 ! Note: Psfc its in hPa in P = A + B * Psfc
33 if (k == kte) then
34 coeff = grid%xb%KMA_B(K)/ &
35 (grid%xb%KMA_A(K)+grid%xb%KMA_B(K)*grid%xb%psfc(I,J)/100.0)
36 else
37 PU = grid%xb%KMA_A(K+1) + &
38 grid%xb%KMA_B(K+1)*grid%xb%psfc(I,J)/100.0
39 PD = grid%xb%KMA_A(K ) + &
40 grid%xb%KMA_B(K )*grid%xb%psfc(I,J)/100.0
41 coeff=grid%xb%KMA_B(K) * &
42 1.0/(PD-PU)**2*(-PU*(LOG(PD)-LOG(PU))+PD-PU) &
43 + grid%xb%KMA_B(K+1)* &
44 1.0/(PD-PU)**2*(PD*(LOG(PD)-LOG(PU))-PD+PU)
45 end if
46 grid%xa%p(i,j,k) = grid%xa%psfc(i,j) * coeff
47 grid%xa%p(i,j,k) = grid%xb%psfc(i,j)*grid%xa%psfc(i,j)
48 grid%xb%p(i,j,k) = grid%xb%p(i,j,k) + grid%xa%p(I,J,k)
49 end do
50 grid%xb%psfc(i,j) = grid%xb%psfc(i,j) + grid%xa%psfc(i,j)
51 end do
52 end do
53
54 if (write_increments) call da_write_kma_increments(grid)
55
56 do j=jts,jte
57 do i=its,ite
58 grid%em_w_2(i,j,kte+1)= grid%em_w_2(i,j,kte+1) + grid%xa%w(i,j,kte+1)
59 grid%psfc(i,j) = grid%psfc(i,j) + grid%xa%psfc(i,j)
60 end do
61 end do
62
63 do k=kts,kte
64 do j=jts,jte
65 do i=its,ite
66 grid%em_u_2(i,j,k) = grid%em_u_2(i,j,k) + grid%xa%u(i,j,k)
67 grid%em_v_2(i,j,k) = grid%em_v_2(i,j,k) + grid%xa%v(i,j,k)
68 grid%em_w_2(i,j,k) = grid%em_w_2(i,j,k) + grid%xa%w(i,j,k)
69 grid%moist(i,j,k,P_QV) = grid%moist(i,j,k,P_QV) + grid%xa%q(i,j,k)
70 end do
71 end do
72 end do
73
74 if (trace_use) call da_trace_exit("da_transfer_xatokma")
75
76 end subroutine da_transfer_xatokma
77
78