da_transform_xtogpsref_adj.inc

1 subroutine da_transform_xtogpsref_adj(grid)
2 
3    !-----------------------------------------------------------------------
4    ! Purpose: TBD
5    ! input: grid%xb%q, grid%xb%p, grid%xb%t, grid%xa%ref, and xp
6    ! output: grid%xa%q, grid%xa%p, grid%xa%t
7    !-----------------------------------------------------------------------
8 
9    implicit none
10 
11    type (domain),  intent(inout) :: grid
12 
13    integer :: i, j, k, is, ie, js, je, ks, ke
14    real    :: partone, parttwo, dividnd
15    real    :: partone9,parttwo9,dividnd9
16 
17    if (trace_use_dull) call da_trace_entry("da_transform_xtogpsref_adj")
18 
19    ! 1.0 Get the range for a tile
20 
21    is = its; ie = ite
22    js = jts; je = jte
23    ks = kts; ke = kte      
24 
25    if (test_wrfvar) then
26       is = its-1
27       js = jts-1
28 
29       ie = ite+1
30       je = jte+1
31 
32       if (is < ids) is = ids
33       if (js < jds) js = jds
34 
35       if (ie > ide) ie = ide
36       if (je > jde) je = jde
37    end if
38 
39    !------------------------------------------------------------------------------
40    !  [2.0] Calculate the adjoint for GPS refractivity:
41    !------------------------------------------------------------------------------
42 
43    do k = ks, ke
44       do j = js, je
45          do i = is, ie
46             !  Note: p in Pascal, q is the specific humidity
47 
48             !   2.1 basic state
49 
50             partone9 = 0.776*grid%xb%p(i,j,k)/grid%xb%t(i,j,k)
51             dividnd9 = grid%xb%t(i,j,k)*(0.622+0.378*grid%xb%q(i,j,k))
52             parttwo9 = 1.0+coeff*grid%xb%q(i,j,k)/dividnd9
53 
54             !   2.2 Adjoint of partone and parttwo
55 
56             partone = grid%xa%ref(i,j,k) * parttwo9
57             parttwo = grid%xa%ref(i,j,k) * partone9
58 
59             !   2.3 Adjoint of q and dividnd
60 
61             grid%xa%q(i,j,k) = grid%xa%q(i,j,k) + coeff*parttwo/dividnd9
62             dividnd=-coeff*grid%xb%q(i,j,k)*parttwo/(dividnd9*dividnd9)
63 
64             !   2.4 Adjoint of t, q, and p
65 
66             grid%xa%t(i,j,k) = grid%xa%t(i,j,k) + dividnd*(0.622+0.378*grid%xb%q(i,j,k))
67             grid%xa%q(i,j,k) = grid%xa%q(i,j,k) + grid%xb%t(i,j,k)*0.378*dividnd
68 
69             grid%xa%p(i,j,k) = grid%xa%p(i,j,k) + 0.776*partone/grid%xb%t(i,j,k)
70             grid%xa%t(i,j,k) = grid%xa%t(i,j,k) - &
71                0.776*grid%xb%p(i,j,k)*partone/(grid%xb%t(i,j,k)*grid%xb%t(i,j,k))
72 
73             !  [3.0] Clean the grid%xa%ref    
74 
75             grid%xa%ref(i,j,k) = 0.0
76          end do
77       end do
78    end do
79 
80    if (trace_use_dull) call da_trace_exit("da_transform_xtogpsref_adj")
81 
82 end subroutine da_transform_xtogpsref_adj
83 
84