1 subroutine da_vv_to_cv(vv, xp, mzs, nmzs, cv_size, rcv &
6 !---------------------------------------------------------------------------
7 ! Purpose: Fill (local) 1D cv array from 3D (local) vv arrays.
8 !---------------------------------------------------------------------------
12 type (vp_type), intent(in) :: vv ! Grdipt/EOF rcv.
13 type (xpose_type), intent(in) :: xp ! Dimensions and xpose buffers.
14 integer, intent(in) :: nmzs ! number of cv for mz
15 integer, intent(in) :: mzs(nmzs) ! Background error structure levels.
16 integer, intent(in) :: cv_size ! Length of CV array.
17 real, intent(inout) :: rcv(1:cv_size) ! Control variables v.
19 integer, optional, intent(in) :: mzs_chem(:) ! chem Background error structure levels.
20 type (xchem_type), optional, intent(inout) :: vchem ! Grdipt/EOF rcv for chem
22 integer :: is,ie ! Local grid range in y coordinate.
23 integer :: js,je ! Local grid range in x coordinate.
24 integer :: ix,jy ! Local grid horizontal dimensions.
25 integer :: mz ! Max vertical coordinate for v1 through v5 arrays.
26 integer :: ne ! Ensemble size.
27 integer :: cv_s,cv_e ! Starting and ending indices into CV array.
28 integer :: ijm ! Size of interior of v1 through v5 arrays.
29 integer :: ijmn ! Size of interior of alpha cv arrays.
31 integer :: ic, offset,tt
35 integer :: i,j,k,ijk,m,n
38 if (trace_use) call da_trace_entry("da_vv_to_cv")
46 call da_error(__FILE__,__LINE__,(/"This subroutine should not be called for use_rf = .false."/))
66 rcv(ijk) = vv%v1(i,j,k)
71 rcv(cv_s:cv_e) = RESHAPE(vv % v1(is:ie,js:je,1:mz), (/ijm/))
85 rcv(ijk) = vv%v2(i,j,k)
90 rcv(cv_s:cv_e) = RESHAPE(vv % v2(is:ie,js:je,1:mz), (/ijm/))
105 rcv(ijk) = vv%v3(i,j,k)
110 rcv(cv_s:cv_e) = RESHAPE(vv % v3(is:ie,js:je,1:mz), (/ijm/))
118 cv_e = cv_s + ijm - 1
124 rcv(ijk) = vv%v4(i,j,k)
129 rcv(cv_s:cv_e) = RESHAPE(vv % v4(is:ie,js:je,1:mz), (/ijm/))
137 cv_e = cv_s + ijm - 1
143 rcv(ijk) = vv%v5(i,j,k)
148 rcv(cv_s:cv_e) = RESHAPE(vv % v5(is:ie,js:je,1:mz), (/ijm/))
151 if ( nmzs > 10 ) then
157 cv_e = cv_s + ijm - 1
163 rcv(ijk) = vv%v6(i,j,k)
168 rcv(cv_s:cv_e) = RESHAPE(vv % v6(is:ie,js:je,1:mz), (/ijm/))
176 cv_e = cv_s + ijm - 1
182 rcv(ijk) = vv%v7(i,j,k)
187 rcv(cv_s:cv_e) = RESHAPE(vv % v7(is:ie,js:je,1:mz), (/ijm/))
195 cv_e = cv_s + ijm - 1
201 rcv(ijk) = vv%v8(i,j,k)
206 rcv(cv_s:cv_e) = RESHAPE(vv % v8(is:ie,js:je,1:mz), (/ijm/))
214 cv_e = cv_s + ijm - 1
220 rcv(ijk) = vv%v9(i,j,k)
225 rcv(cv_s:cv_e) = RESHAPE(vv % v9(is:ie,js:je,1:mz), (/ijm/))
233 cv_e = cv_s + ijm - 1
239 rcv(ijk) = vv%v10(i,j,k)
244 rcv(cv_s:cv_e) = RESHAPE(vv % v10(is:ie,js:je,1:mz), (/ijm/))
248 if ( nmzs == 8 .or. nmzs == 13 ) then
254 cv_e = cv_s + ijm - 1
260 rcv(ijk) = vv%v11(i,j,k)
265 rcv(cv_s:cv_e) = RESHAPE(vv % v11(is:ie,js:je,1:mz), (/ijm/))
274 ix = ite_int - its_int + 1
275 jy = jte_int - jts_int + 1
276 ijmn = ix * jy * mz * ne
278 cv_e = cv_s + ijmn - 1
282 do j = jts_int, jte_int !js, je
283 do i = its_int, ite_int ! is, ie
285 rcv(ijk) = vv%alpha(i,j,m,n)
291 ! rcv(cv_s:cv_e) = RESHAPE(vv % alpha(is:ie,js:je,1:mz,1:ne), (/ijmn/))
292 rcv(cv_s:cv_e) = RESHAPE(vv % alpha(its_int:ite_int,jts_int:jte_int,1:mz,1:ne), (/ijmn/))
297 if (present(mzs_chem) .and. present(vchem)) then
302 if ( num_chem >= PARAM_FIRST_SCALAR ) then
303 do ic = PARAM_FIRST_SCALAR, num_chem
307 cv_e = cv_s + ijm - 1
313 rcv(ijk) = vchem % chem_ic(i,j,k,ic)
318 rcv(cv_s:cv_e) = RESHAPE(vchem % chem_ic(is:ie,js:je,1:mz,ic))
327 if (trace_use) call da_trace_exit("da_vv_to_cv")
329 end subroutine da_vv_to_cv