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