da_vv_to_cv.inc
References to this file elsewhere.
1 subroutine da_vv_to_cv(vv, xp, mzs, cv_size, rcv)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Fill (local) 1D cv array from 3D (local) vv arrays.
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (vp_type), intent(in) :: vv ! Grdipt/EOF rcv.
10 type (xpose_type), intent(in) :: xp ! Dimensions and xpose buffers.
11 integer, intent(in) :: mzs(:) ! Background error structure levels.
12 integer, intent(in) :: cv_size ! Length of CV array.
13 real, intent(out) :: rcv(1:cv_size) ! Control variables v.
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 #ifdef NORESHAPE
23 integer :: i,j,k,ijk
24 #endif
25
26 if (trace_use) call da_trace_entry("da_vv_to_cv")
27
28 is = xp % its
29 ie = xp % ite
30 js = xp % jts
31 je = xp % jte
32 ix = ie-is+1
33 jy = je-js+1
34 cv_e = 0
35
36 ! Store v1
37 mz = mzs(1)
38 if (mz > 0) then
39 ijm = ix * jy * mz
40 cv_s = cv_e + 1
41 cv_e = cv_s + ijm - 1
42 #ifdef NORESHAPE
43 ijk = 0
44 do k=1,mz
45 do j=js,je
46 do i=is,ie
47 ijk = ijk + 1
48 rcv(ijk) = vv%v1(i,j,k)
49 end do
50 end do
51 end do
52 #else
53 rcv(cv_s:cv_e) = RESHAPE(vv % v1(is:ie,js:je,1:mz), (/ijm/))
54 #endif
55 end if
56
57 ! Store v2
58 mz = mzs(2)
59 if (mz > 0) then
60 ijm = ix * jy * mz
61 cv_s = cv_e + 1
62 cv_e = cv_s + ijm - 1
63 #ifdef NORESHAPE
64 do k=1,mz
65 do j=js,je
66 do i=is,ie
67 ijk = ijk + 1
68 rcv(ijk) = vv%v2(i,j,k)
69 end do
70 end do
71 end do
72 #else
73 rcv(cv_s:cv_e) = RESHAPE(vv % v2(is:ie,js:je,1:mz), (/ijm/))
74 #endif
75 end if
76
77 ! Store v3
78 mz = mzs(3)
79 if (mz > 0) then
80 ijm = ix * jy * mz
81 cv_s = cv_e + 1
82 cv_e = cv_s + ijm - 1
83
84 #ifdef NORESHAPE
85 do k=1,mz
86 do j=js,je
87 do i=is,ie
88 ijk = ijk + 1
89 rcv(ijk) = vv%v3(i,j,k)
90 end do
91 end do
92 end do
93 #else
94 rcv(cv_s:cv_e) = RESHAPE(vv % v3(is:ie,js:je,1:mz), (/ijm/))
95 #endif
96 end if
97
98 ! Store v4
99 mz = mzs(4)
100 if (mz > 0) then
101 ijm = ix * jy * mz
102 cv_s = cv_e + 1
103 cv_e = cv_s + ijm - 1
104 #ifdef NORESHAPE
105 do k=1,mz
106 do j=js,je
107 do i=is,ie
108 ijk = ijk + 1
109 rcv(ijk) = vv%v4(i,j,k)
110 end do
111 end do
112 end do
113 #else
114 rcv(cv_s:cv_e) = RESHAPE(vv % v4(is:ie,js:je,1:mz), (/ijm/))
115 #endif
116 end if
117
118 ! Store v5
119 mz = mzs(5)
120 if (mz > 0) then
121 ijm = ix * jy * mz
122 cv_s = cv_e + 1
123 cv_e = cv_s + ijm - 1
124 #ifdef NORESHAPE
125 do k=1,mz
126 do j=js,je
127 do i=is,ie
128 ijk = ijk + 1
129 rcv(ijk) = vv%v5(i,j,k)
130 end do
131 end do
132 end do
133 #else
134 rcv(cv_s:cv_e) = RESHAPE(vv % v5(is:ie,js:je,1:mz), (/ijm/))
135 #endif
136 end if
137
138 ! Store alpha:
139 mz = mzs(6)
140 if (mz > 0) then
141 ijm = ix * jy * mz
142 cv_s = cv_e + 1
143 cv_e = cv_s + ijm - 1
144 #ifdef NORESHAPE
145 do k=1,mz
146 do j=js,je
147 do i=is,ie
148 ijk = ijk + 1
149 rcv(ijk) = vv%alpha(i,j,k)
150 end do
151 end do
152 end do
153 #else
154 rcv(cv_s:cv_e) = RESHAPE(vv % alpha(is:ie,js:je,1:mz), (/ijm/))
155 #endif
156 end if
157
158 if (trace_use) call da_trace_exit("da_vv_to_cv")
159
160 end subroutine da_vv_to_cv
161
162