da_cv_to_vv.inc
References to this file elsewhere.
1 subroutine da_cv_to_vv (cv_size, rcv, mzs, vv)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Fill (local) vv arrays from 1D (local) cv array.
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: cv_size ! Length of CV array.
10 real, intent(in) :: rcv(1:cv_size) ! Control variables v.
11 integer, intent(in) :: mzs(:) ! Background error structure levels.
12 type (vp_type), intent(inout) :: vv ! Grdipt/EOF cv_array.
13
14 integer :: is,ie ! Local grid range in y coordinate.
15 integer :: js,je ! Local grid range in x coordinate.
16 integer :: ix,jy ! Local grid horizontal dimensions.
17 integer :: mz ! Max vertical coordinate for v1 through v5 arrays.
18 integer :: cv_s,cv_e ! Starting and ending indices into CV array.
19 integer :: ijm ! Size of interior of v1 through v5 arrays.
20
21 if (trace_use) call da_trace_entry("da_cv_to_vv")
22
23 is = its
24 ie = ite
25 js = jts
26 je = jte
27 ix = ie-is+1
28 jy = je-js+1
29 cv_e = 0
30
31 !--------------------------------------------------------------------------
32 ! [1] Transfer components of Jb control variable:
33 !--------------------------------------------------------------------------
34
35 ! Fill v1
36 mz = mzs(1)
37 if (mz > 0) then
38 ijm = ix * jy * mz
39 cv_s = cv_e + 1
40 cv_e = cv_s + ijm - 1
41 vv % v1(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy, mz/))
42 end if
43
44 ! Fill v2
45 mz = mzs(2)
46 if (mz > 0) then
47 ijm = ix * jy * mz
48 cv_s = cv_e + 1
49 cv_e = cv_s + ijm - 1
50 vv % v2(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy, mz/))
51 end if
52
53 ! Fill v3
54 mz = mzs(3)
55 if (mz > 0) then
56 ijm = ix * jy * mz
57 cv_s = cv_e + 1
58 cv_e = cv_s + ijm - 1
59 vv % v3(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy, mz/))
60 end if
61
62 ! Fill v4
63 mz = mzs(4)
64 if (mz > 0) then
65 ijm = ix * jy * mz
66 cv_s = cv_e + 1
67 cv_e = cv_s + ijm - 1
68 vv % v4(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy, mz/))
69 end if
70
71 ! Fill v5
72 mz = mzs(5)
73 if (mz == 1) then ! Can only be 0 or 1 (2D ps_u field)
74 ijm = ix * jy * mz
75 cv_s = cv_e + 1
76 cv_e = cv_s + ijm - 1
77 vv % v5(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy,mz/))
78 end if
79
80 !--------------------------------------------------------------------------
81 ! [2] Transfer components of Je control variable:
82 !--------------------------------------------------------------------------
83
84 mz = mzs(6)
85 if (mz > 0) then
86 ijm = ix * jy * mz
87 cv_s = cv_e + 1
88 cv_e = cv_s + ijm - 1
89 vv % alpha(is:ie,js:je,1:mz) = RESHAPE(rcv(cv_s:cv_e),(/ix, jy, mz/))
90 end if
91
92 if (trace_use) call da_trace_exit("da_cv_to_vv")
93
94 end subroutine da_cv_to_vv
95
96