Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Blame
Revision Log

interp_fcn.F @ 354

Last change on this file since 354 was 1, checked in by lfita, 10 years ago

-- --- Opening of the WRF+LMDZ coupling repository --- -- -

WRF: version v3.3
LMDZ: version v1818

More details in:

File size: 201.1 KB

Line
1	!WRF:MEDIATION_LAYER:INTERPOLATIONFUNCTION
2	!
3
4	#if (DA_CORE != 1)
5	#define MM5_SINT
6	#endif
7	!#define DUMBCOPY
8
9	! Note, NMM-specific routines moved to end. 20080612. JM
10
11	SUBROUTINE interp_fcn ( cfld, & ! CD field
12	cids, cide, ckds, ckde, cjds, cjde, &
13	cims, cime, ckms, ckme, cjms, cjme, &
14	cits, cite, ckts, ckte, cjts, cjte, &
15	nfld, & ! ND field
16	nids, nide, nkds, nkde, njds, njde, &
17	nims, nime, nkms, nkme, njms, njme, &
18	nits, nite, nkts, nkte, njts, njte, &
19	shw, & ! stencil half width for interp
20	imask, & ! interpolation mask
21	xstag, ystag, & ! staggering of field
22	ipos, jpos, & ! Position of lower left of nest in CD
23	nri, nrj ) ! nest ratios
24	USE module_timing
25	USE module_configure
26	IMPLICIT NONE
27
28
29	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
30	cims, cime, ckms, ckme, cjms, cjme, &
31	cits, cite, ckts, ckte, cjts, cjte, &
32	nids, nide, nkds, nkde, njds, njde, &
33	nims, nime, nkms, nkme, njms, njme, &
34	nits, nite, nkts, nkte, njts, njte, &
35	shw, &
36	ipos, jpos, &
37	nri, nrj
38	LOGICAL, INTENT(IN) :: xstag, ystag
39
40	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
41	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
42	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
43
44	! Local
45
46	!logical first
47
48	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, nioff, njoff
49	#ifdef MM5_SINT
50	INTEGER nfx, ior
51	PARAMETER (ior=2)
52	INTEGER nf
53	REAL psca(cims:cime,cjms:cjme,nri*nrj)
54	LOGICAL icmask( cims:cime, cjms:cjme )
55	INTEGER i,j,k
56	INTEGER nrio2, nrjo2
57	#endif
58
59	! Iterate over the ND tile and compute the values
60	! from the CD tile.
61
62	#ifdef MM5_SINT
63
64	ioff = 0 ; joff = 0
65	nioff = 0 ; njoff = 0
66	IF ( xstag ) THEN
67	ioff = (nri-1)/2
68	nioff = nri
69	ENDIF
70	IF ( ystag ) THEN
71	joff = (nrj-1)/2
72	njoff = nrj
73	ENDIF
74
75	nrio2 = nri/2
76	nrjo2 = nrj/2
77
78	nfx = nri * nrj
79	!$OMP PARALLEL DO &
80	!$OMP PRIVATE ( i,j,k,ni,nj,ci,cj,ip,jp,nk,ck,nf,icmask,psca )
81	DO k = ckts, ckte
82	icmask = .FALSE.
83	DO nf = 1,nfx
84	DO j = cjms,cjme
85	nj = (j-jpos) * nrj + ( nrjo2 + 1 ) ! j point on nest
86	DO i = cims,cime
87	ni = (i-ipos) * nri + ( nrio2 + 1 ) ! i point on nest
88	if ( ni .ge. nits-nioff-nrio2 .and. &
89	ni .le. nite+nioff+nrio2 .and. &
90	nj .ge. njts-njoff-nrjo2 .and. &
91	nj .le. njte+njoff+nrjo2 ) then
92	! if ( imask(ni,nj) .eq. 1 .or. imask(ni-nioff,nj-njoff) .eq. 1 ) then
93	! icmask( i, j ) = .TRUE.
94	! endif
95	if ( ni.ge.nims.and.ni.le.nime.and.nj.ge.njms.and.nj.le.njme) then
96	if ( imask(ni,nj) .eq. 1 ) then
97	icmask( i, j ) = .TRUE.
98	endif
99	endif
100	if ( ni-nioff.ge.nims.and.ni.le.nime.and.nj-njoff.ge.njms.and.nj.le.njme) then
101	if (ni .ge. nits-nioff .and. nj .ge. njts-njoff ) then
102	if ( imask(ni-nioff,nj-njoff) .eq. 1) then
103	icmask( i, j ) = .TRUE.
104	endif
105	endif
106	endif
107	endif
108	psca(i,j,nf) = cfld(i,k,j)
109	ENDDO
110	ENDDO
111	ENDDO
112
113	! tile dims in this call to sint are 1-over to account for the fact
114	! that the number of cells on the nest local subdomain is not
115	! necessarily a multiple of the nest ratio in a given dim.
116	! this could be a little less ham-handed.
117
118	!call start_timing
119
120	CALL sint( psca, &
121	cims, cime, cjms, cjme, icmask, &
122	cits-1, cite+1, cjts-1, cjte+1, nrj*nri, xstag, ystag )
123
124	!call end_timing( ' sint ' )
125
126	DO nj = njts, njte+joff
127	cj = jpos + (nj-1) / nrj ! j coord of CD point
128	jp = mod ( nj-1 , nrj ) ! coord of ND w/i CD point
129	nk = k
130	ck = nk
131	DO ni = nits, nite+ioff
132	ci = ipos + (ni-1) / nri ! i coord of CD point
133	ip = mod ( ni-1 , nri ) ! coord of ND w/i CD point
134	if ( imask ( ni, nj ) .eq. 1 .or. imask ( ni-ioff, nj-joff ) .eq. 1 ) then
135	nfld( ni-ioff, nk, nj-joff ) = psca( ci , cj, ip+1 + (jp)*nri )
136	endif
137	ENDDO
138	ENDDO
139	ENDDO
140	!$OMP END PARALLEL DO
141	#endif
142
143	#ifdef DUMBCOPY
144	!write(0,'(") cims:cime, ckms:ckme, cjms:cjme ",6i4)')cims,cime, ckms,ckme, cjms,cjme
145	!write(0,'(") nims:nime, nkms:nkme, njms:njme ",6i4)')nims,nime, nkms,nkme, njms,njme
146	!write(0,'(") cits:cite, ckts:ckte, cjts:cjte ",6i4)')cits,cite, ckts,ckte, cjts,cjte
147	!write(0,'(") nits:nite, nkts:nkte, njts:njte ",6i4)')nits,nite, nkts,nkte, njts,njte
148
149	DO nj = njts, njte
150	cj = jpos + (nj-1) / nrj ! j coord of CD point
151	jp = mod ( nj , nrj ) ! coord of ND w/i CD point
152	DO nk = nkts, nkte
153	ck = nk
154	DO ni = nits, nite
155	ci = ipos + (ni-1) / nri ! j coord of CD point
156	ip = mod ( ni , nri ) ! coord of ND w/i CD point
157	! This is a trivial implementation of the interp_fcn; just copies
158	! the values from the CD into the ND
159	if ( imask ( ni, nj ) .eq. 1 ) then
160	nfld( ni, nk, nj ) = cfld( ci , ck , cj )
161	endif
162	ENDDO
163	ENDDO
164	ENDDO
165	#endif
166
167	RETURN
168
169	END SUBROUTINE interp_fcn
170
171	!==================================
172	! this is the default function used in feedback.
173
174	SUBROUTINE copy_fcn ( cfld, & ! CD field
175	cids, cide, ckds, ckde, cjds, cjde, &
176	cims, cime, ckms, ckme, cjms, cjme, &
177	cits, cite, ckts, ckte, cjts, cjte, &
178	nfld, & ! ND field
179	nids, nide, nkds, nkde, njds, njde, &
180	nims, nime, nkms, nkme, njms, njme, &
181	nits, nite, nkts, nkte, njts, njte, &
182	shw, & ! stencil half width for interp
183	imask, & ! interpolation mask
184	xstag, ystag, & ! staggering of field
185	ipos, jpos, & ! Position of lower left of nest in CD
186	nri, nrj ) ! nest ratios
187	USE module_configure
188	IMPLICIT NONE
189
190
191	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
192	cims, cime, ckms, ckme, cjms, cjme, &
193	cits, cite, ckts, ckte, cjts, cjte, &
194	nids, nide, nkds, nkde, njds, njde, &
195	nims, nime, nkms, nkme, njms, njme, &
196	nits, nite, nkts, nkte, njts, njte, &
197	shw, &
198	ipos, jpos, &
199	nri, nrj
200	LOGICAL, INTENT(IN) :: xstag, ystag
201
202	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(OUT) :: cfld
203	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ),INTENT(IN) :: nfld
204	INTEGER, DIMENSION ( nims:nime, njms:njme ),INTENT(IN) :: imask
205
206	! Local
207
208	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
209	INTEGER :: icmin,icmax,jcmin,jcmax
210	INTEGER :: istag,jstag, ipoints,jpoints,ijpoints
211	INTEGER , PARAMETER :: passes = 2
212	INTEGER spec_zone
213
214	! Loop over the coarse grid in the area of the fine mesh. Do not
215	! process the coarse grid values that are along the lateral BC
216	! provided to the fine grid. Since that is in the specified zone
217	! for the fine grid, it should not be used in any feedback to the
218	! coarse grid as it should not have changed.
219
220	! Due to peculiarities of staggering, it is simpler to handle the feedback
221	! for the staggerings based upon whether it is a even ratio (2::1, 4::1, etc.) or
222	! an odd staggering ratio (3::1, 5::1, etc.).
223
224	! Though there are separate grid ratios for the i and j directions, this code
225	! is not general enough to handle aspect ratios .NE. 1 for the fine grid cell.
226
227	! These are local integer increments in the looping. Basically, istag=1 means
228	! that we will assume one less point in the i direction. Note that ci and cj
229	! have a maximum value that is decreased by istag and jstag, respectively.
230
231	! Horizontal momentum feedback is along the face, not within the cell. For a
232	! 3::1 ratio, temperature would use 9 pts for feedback, while u and v use
233	! only 3 points for feedback from the nest to the parent.
234
235	CALL nl_get_spec_zone( 1 , spec_zone )
236	istag = 1 ; jstag = 1
237	IF ( xstag ) istag = 0
238	IF ( ystag ) jstag = 0
239
240	IF( MOD(nrj,2) .NE. 0) THEN ! odd refinement ratio
241
242	IF ( ( .NOT. xstag ) .AND. ( .NOT. ystag ) ) THEN
243	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
244	nj = (cj-jpos)*nrj + jstag + 1
245	DO ck = ckts, ckte
246	nk = ck
247	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
248	ni = (ci-ipos)*nri + istag + 1
249	cfld( ci, ck, cj ) = 0.
250	DO ijpoints = 1 , nri * nrj
251	ipoints = MOD((ijpoints-1),nri) + 1 - nri/2 - 1
252	jpoints = (ijpoints-1)/nri + 1 - nrj/2 - 1
253	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
254	1./REAL(nrinrj) nfld( ni+ipoints , nk , nj+jpoints )
255	END DO
256	! cfld( ci, ck, cj ) = 1./9. * &
257	! ( nfld( ni-1, nk , nj-1) + &
258	! nfld( ni , nk , nj-1) + &
259	! nfld( ni+1, nk , nj-1) + &
260	! nfld( ni-1, nk , nj ) + &
261	! nfld( ni , nk , nj ) + &
262	! nfld( ni+1, nk , nj ) + &
263	! nfld( ni-1, nk , nj+1) + &
264	! nfld( ni , nk , nj+1) + &
265	! nfld( ni+1, nk , nj+1) )
266	ENDDO
267	ENDDO
268	ENDDO
269
270	ELSE IF ( ( xstag ) .AND. ( .NOT. ystag ) ) THEN
271	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
272	nj = (cj-jpos)*nrj + jstag + 1
273	DO ck = ckts, ckte
274	nk = ck
275	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
276	ni = (ci-ipos)*nri + istag + 1
277	cfld( ci, ck, cj ) = 0.
278	DO ijpoints = (nri+1)/2 , (nri+1)/2 + nri*(nri-1) , nri
279	ipoints = MOD((ijpoints-1),nri) + 1 - nri/2 - 1
280	jpoints = (ijpoints-1)/nri + 1 - nrj/2 - 1
281	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
282	1./REAL(nri ) * nfld( ni+ipoints , nk , nj+jpoints )
283	END DO
284	! cfld( ci, ck, cj ) = 1./3. * &
285	! ( nfld( ni , nk , nj-1) + &
286	! nfld( ni , nk , nj ) + &
287	! nfld( ni , nk , nj+1) )
288	ENDDO
289	ENDDO
290	ENDDO
291
292	ELSE IF ( ( .NOT. xstag ) .AND. ( ystag ) ) THEN
293	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
294	nj = (cj-jpos)*nrj + jstag + 1
295	DO ck = ckts, ckte
296	nk = ck
297	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
298	ni = (ci-ipos)*nri + istag + 1
299	cfld( ci, ck, cj ) = 0.
300	DO ijpoints = ( nrjnrj +1 )/2 - nrj/2 , ( nrjnrj +1 )/2 - nrj/2 + nrj-1
301	ipoints = MOD((ijpoints-1),nri) + 1 - nri/2 - 1
302	jpoints = (ijpoints-1)/nri + 1 - nrj/2 - 1
303	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
304	1./REAL( nrj) * nfld( ni+ipoints , nk , nj+jpoints )
305	END DO
306	! cfld( ci, ck, cj ) = 1./3. * &
307	! ( nfld( ni-1, nk , nj ) + &
308	! nfld( ni , nk , nj ) + &
309	! nfld( ni+1, nk , nj ) )
310	ENDDO
311	ENDDO
312	ENDDO
313
314	END IF
315
316	! Even refinement ratio
317
318	ELSE IF ( MOD(nrj,2) .EQ. 0) THEN
319	IF ( ( .NOT. xstag ) .AND. ( .NOT. ystag ) ) THEN
320
321	! This is a simple schematic of the feedback indexing used in the even
322	! ratio nests. For simplicity, a 2::1 ratio is depicted. Only the
323	! mass variable staggering is shown.
324	! Each of
325	! the boxes with a "T" and four small "t" represents a coarse grid (CG)
326	! cell, that is composed of four (2::1 ratio) fine grid (FG) cells.
327
328	! Shown below is the area of the CG that is in the area of the FG. The
329	! first grid point of the depicted CG is the starting location of the nest
330	! in the parent domain (ipos,jpos - i_parent_start and j_parent_start from
331	! the namelist).
332
333	! For each of the CG points, the feedback loop is over each of the FG points
334	! within the CG cell. For a 2::1 ratio, there are four total points (this is
335	! the ijpoints loop). The feedback value to the CG is the arithmetic mean of
336	! all of the FG values within each CG cell.
337
338	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
339	! \| t t \|\| t t \| \| t t \|\| t t \|
340	! jpos+ \| \|\| \| \| \|\| \|
341	! (njde-njds)- \| T \|\| T \| \| T \|\| T \|
342	! jstag \| \|\| \| \| \|\| \|
343	! \| t t \|\| t t \| \| t t \|\| t t \|
344	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
345	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
346	! \| t t \|\| t t \| \| t t \|\| t t \|
347	! \| \|\| \| \| \|\| \|
348	! \| T \|\| T \| \| T \|\| T \|
349	! \| \|\| \| \| \|\| \|
350	! \| t t \|\| t t \| \| t t \|\| t t \|
351	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
352	!
353	! ...
354	! ...
355	! ...
356	! ...
357	! ...
358
359	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
360	! jpoints = 1 \| t t \|\| t t \| \| t t \|\| t t \|
361	! \| \|\| \| \| \|\| \|
362	! \| T \|\| T \| \| T \|\| T \|
363	! \| \|\| \| \| \|\| \|
364	! jpoints = 0, \| t t \|\| t t \| \| t t \|\| t t \|
365	! nj=3 \|-------------\|\|-------------\| \|-------------\|\|-------------\|
366	! \|-------------\|\|-------------\| \|-------------\|\|-------------\|
367	! jpoints = 1 \| t t \|\| t t \| \| t t \|\| t t \|
368	! \| \|\| \| \| \|\| \|
369	! jpos \| T \|\| T \| ... \| T \|\| T \|
370	! \| \|\| \| ... \| \|\| \|
371	! jpoints = 0, \| t t \|\| t t \| ... \| t t \|\| t t \|
372	! nj=1 \|-------------\|\|-------------\| \|-------------\|\|-------------\|
373	! ^ ^
374	! \| \|
375	! \| \|
376	! ipos ipos+
377	! ni = 1 3 (nide-nids)/nri
378	! ipoints= 0 1 0 1 -istag
379	!
380
381	! For performance benefits, users can comment out the inner most loop (and cfld=0) and
382	! hardcode the loop feedback. For example, it is set up to run a 2::1 ratio
383	! if uncommented. This lacks generality, but is likely to gain timing benefits
384	! with compilers unable to unroll inner loops that do not have parameterized sizes.
385
386	! The extra +1 ---------/ and the extra -1 ----\ (both for ci and cj)
387	! / \ keeps the feedback out of the
388	! / \ outer row/col, since that CG data
389	! / \ specified the nest boundary originally
390	! / \ This
391	! / \ is just
392	! / \ a sentence to not end a line
393	! / \ with a stupid backslash
394	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
395	nj = (cj-jpos)*nrj + jstag
396	DO ck = ckts, ckte
397	nk = ck
398	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
399	ni = (ci-ipos)*nri + istag
400	cfld( ci, ck, cj ) = 0.
401	DO ijpoints = 1 , nri * nrj
402	ipoints = MOD((ijpoints-1),nri)
403	jpoints = (ijpoints-1)/nri
404	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
405	1./REAL(nrinrj) nfld( ni+ipoints , nk , nj+jpoints )
406	END DO
407	! cfld( ci, ck, cj ) = 1./4. * &
408	! ( nfld( ni , nk , nj ) + &
409	! nfld( ni+1, nk , nj ) + &
410	! nfld( ni , nk , nj+1) + &
411	! nfld( ni+1, nk , nj+1) )
412	END DO
413	END DO
414	END DO
415
416	! U
417
418	ELSE IF ( ( xstag ) .AND. ( .NOT. ystag ) ) THEN
419	! \|---------------\|
420	! \| \|
421	! jpoints = 1 u u \|
422	! \| \|
423	! U \|
424	! \| \|
425	! jpoints = 0, u u \|
426	! nj=3 \| \|
427	! \|---------------\|
428	! \|---------------\|
429	! \| \|
430	! jpoints = 1 u u \|
431	! \| \|
432	! jpos U \|
433	! \| \|
434	! jpoints = 0, u u \|
435	! nj=1 \| \|
436	! \|---------------\|
437	!
438	! ^
439	! \|
440	! \|
441	! ipos
442	! ni = 1 3
443	! ipoints= 0 1 0
444	!
445
446	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
447	nj = (cj-jpos)*nrj + 1
448	DO ck = ckts, ckte
449	nk = ck
450	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
451	ni = (ci-ipos)*nri + 1
452	cfld( ci, ck, cj ) = 0.
453	DO ijpoints = 1 , nri*nrj , nri
454	ipoints = MOD((ijpoints-1),nri)
455	jpoints = (ijpoints-1)/nri
456	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
457	1./REAL(nri ) * nfld( ni+ipoints , nk , nj+jpoints )
458	END DO
459	! cfld( ci, ck, cj ) = 1./2. * &
460	! ( nfld( ni , nk , nj ) + &
461	! nfld( ni , nk , nj+1) )
462	ENDDO
463	ENDDO
464	ENDDO
465
466	! V
467
468	ELSE IF ( ( .NOT. xstag ) .AND. ( ystag ) ) THEN
469	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
470	nj = (cj-jpos)*nrj + 1
471	DO ck = ckts, ckte
472	nk = ck
473	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
474	ni = (ci-ipos)*nri + 1
475	cfld( ci, ck, cj ) = 0.
476	DO ijpoints = 1 , nri
477	ipoints = MOD((ijpoints-1),nri)
478	jpoints = (ijpoints-1)/nri
479	cfld( ci, ck, cj ) = cfld( ci, ck, cj ) + &
480	1./REAL(nri ) * nfld( ni+ipoints , nk , nj+jpoints )
481	END DO
482	! cfld( ci, ck, cj ) = 1./2. * &
483	! ( nfld( ni , nk , nj ) + &
484	! nfld( ni+1, nk , nj ) )
485	ENDDO
486	ENDDO
487	ENDDO
488	END IF
489	END IF
490
491	RETURN
492
493	END SUBROUTINE copy_fcn
494
495	!==================================
496	! this is the 1pt function used in feedback.
497
498	SUBROUTINE copy_fcnm ( cfld, & ! CD field
499	cids, cide, ckds, ckde, cjds, cjde, &
500	cims, cime, ckms, ckme, cjms, cjme, &
501	cits, cite, ckts, ckte, cjts, cjte, &
502	nfld, & ! ND field
503	nids, nide, nkds, nkde, njds, njde, &
504	nims, nime, nkms, nkme, njms, njme, &
505	nits, nite, nkts, nkte, njts, njte, &
506	shw, & ! stencil half width for interp
507	imask, & ! interpolation mask
508	xstag, ystag, & ! staggering of field
509	ipos, jpos, & ! Position of lower left of nest in CD
510	nri, nrj ) ! nest ratios
511	USE module_configure
512	USE module_wrf_error
513	IMPLICIT NONE
514
515
516	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
517	cims, cime, ckms, ckme, cjms, cjme, &
518	cits, cite, ckts, ckte, cjts, cjte, &
519	nids, nide, nkds, nkde, njds, njde, &
520	nims, nime, nkms, nkme, njms, njme, &
521	nits, nite, nkts, nkte, njts, njte, &
522	shw, &
523	ipos, jpos, &
524	nri, nrj
525	LOGICAL, INTENT(IN) :: xstag, ystag
526
527	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(OUT) :: cfld
528	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ), INTENT(IN) :: nfld
529	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: imask
530
531	! Local
532
533	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
534	INTEGER :: icmin,icmax,jcmin,jcmax
535	INTEGER :: istag,jstag, ipoints,jpoints,ijpoints
536	INTEGER , PARAMETER :: passes = 2
537	INTEGER spec_zone
538
539	CALL nl_get_spec_zone( 1, spec_zone )
540	istag = 1 ; jstag = 1
541	IF ( xstag ) istag = 0
542	IF ( ystag ) jstag = 0
543
544	IF( MOD(nrj,2) .NE. 0) THEN ! odd refinement ratio
545
546	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
547	nj = (cj-jpos)*nrj + jstag + 1
548	DO ck = ckts, ckte
549	nk = ck
550	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
551	ni = (ci-ipos)*nri + istag + 1
552	cfld( ci, ck, cj ) = nfld( ni , nk , nj )
553	ENDDO
554	ENDDO
555	ENDDO
556
557	ELSE ! even refinement ratio, pick nearest neighbor on SW corner
558	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
559	nj = (cj-jpos)*nrj + 1
560	DO ck = ckts, ckte
561	nk = ck
562	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
563	ni = (ci-ipos)*nri + 1
564	ipoints = nri/2 -1
565	jpoints = nrj/2 -1
566	cfld( ci, ck, cj ) = nfld( ni+ipoints , nk , nj+jpoints )
567	END DO
568	END DO
569	END DO
570
571	END IF
572
573	RETURN
574
575	END SUBROUTINE copy_fcnm
576
577	!==================================
578	! this is the 1pt function used in feedback for integers
579
580	SUBROUTINE copy_fcni ( cfld, & ! CD field
581	cids, cide, ckds, ckde, cjds, cjde, &
582	cims, cime, ckms, ckme, cjms, cjme, &
583	cits, cite, ckts, ckte, cjts, cjte, &
584	nfld, & ! ND field
585	nids, nide, nkds, nkde, njds, njde, &
586	nims, nime, nkms, nkme, njms, njme, &
587	nits, nite, nkts, nkte, njts, njte, &
588	shw, & ! stencil half width for interp
589	imask, & ! interpolation mask
590	xstag, ystag, & ! staggering of field
591	ipos, jpos, & ! Position of lower left of nest in CD
592	nri, nrj ) ! nest ratios
593	USE module_configure
594	USE module_wrf_error
595	IMPLICIT NONE
596
597
598	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
599	cims, cime, ckms, ckme, cjms, cjme, &
600	cits, cite, ckts, ckte, cjts, cjte, &
601	nids, nide, nkds, nkde, njds, njde, &
602	nims, nime, nkms, nkme, njms, njme, &
603	nits, nite, nkts, nkte, njts, njte, &
604	shw, &
605	ipos, jpos, &
606	nri, nrj
607	LOGICAL, INTENT(IN) :: xstag, ystag
608
609	INTEGER, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(OUT) :: cfld
610	INTEGER, DIMENSION ( nims:nime, nkms:nkme, njms:njme ), INTENT(IN) :: nfld
611	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: imask
612
613	! Local
614
615	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
616	INTEGER :: icmin,icmax,jcmin,jcmax
617	INTEGER :: istag,jstag, ipoints,jpoints,ijpoints
618	INTEGER , PARAMETER :: passes = 2
619	INTEGER spec_zone
620
621	CALL nl_get_spec_zone( 1, spec_zone )
622	istag = 1 ; jstag = 1
623	IF ( xstag ) istag = 0
624	IF ( ystag ) jstag = 0
625
626	IF( MOD(nrj,2) .NE. 0) THEN ! odd refinement ratio
627
628	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
629	nj = (cj-jpos)*nrj + jstag + 1
630	DO ck = ckts, ckte
631	nk = ck
632	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
633	ni = (ci-ipos)*nri + istag + 1
634	cfld( ci, ck, cj ) = nfld( ni , nk , nj )
635	ENDDO
636	ENDDO
637	ENDDO
638
639	ELSE ! even refinement ratio
640	DO cj = MAX(jpos+spec_zone,cjts),MIN(jpos+(njde-njds)/nrj-jstag-spec_zone,cjte)
641	nj = (cj-jpos)*nrj + 1
642	DO ck = ckts, ckte
643	nk = ck
644	DO ci = MAX(ipos+spec_zone,cits),MIN(ipos+(nide-nids)/nri-istag-spec_zone,cite)
645	ni = (ci-ipos)*nri + 1
646	ipoints = nri/2 -1
647	jpoints = nrj/2 -1
648	cfld( ci, ck, cj ) = nfld( ni+ipoints , nk , nj+jpoints )
649	END DO
650	END DO
651	END DO
652
653	END IF
654
655	RETURN
656
657	END SUBROUTINE copy_fcni
658
659	!==================================
660
661	SUBROUTINE p2c ( cfld, & ! CD field
662	cids, cide, ckds, ckde, cjds, cjde, &
663	cims, cime, ckms, ckme, cjms, cjme, &
664	cits, cite, ckts, ckte, cjts, cjte, &
665	nfld, & ! ND field
666	nids, nide, nkds, nkde, njds, njde, &
667	nims, nime, nkms, nkme, njms, njme, &
668	nits, nite, nkts, nkte, njts, njte, &
669	shw, & ! stencil half width
670	imask, & ! interpolation mask
671	xstag, ystag, & ! staggering of field
672	ipos, jpos, & ! Position of lower left of nest in CD
673	nri, nrj, & ! nest ratios
674	cbdy_xs, nbdy_xs, & ! boundary arrays
675	cbdy_xe, nbdy_xe, & ! no needed, but registry builds interface
676	cbdy_ys, nbdy_ys, & ! for lateral bc forcing this way
677	cbdy_ye, nbdy_ye, &
678	cbdy_txs, nbdy_txs, &
679	cbdy_txe, nbdy_txe, &
680	cbdy_tys, nbdy_tys, &
681	cbdy_tye, nbdy_tye &
682	)
683	USE module_configure
684	IMPLICIT NONE
685
686	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
687	cims, cime, ckms, ckme, cjms, cjme, &
688	cits, cite, ckts, ckte, cjts, cjte, &
689	nids, nide, nkds, nkde, njds, njde, &
690	nims, nime, nkms, nkme, njms, njme, &
691	nits, nite, nkts, nkte, njts, njte, &
692	shw, &
693	ipos, jpos, &
694	nri, nrj
695
696	LOGICAL, INTENT(IN) :: xstag, ystag
697
698	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
699	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
700	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
701	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_xs, cbdy_txs, nbdy_xs, nbdy_txs
702	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_xe, cbdy_txe, nbdy_xe, nbdy_txe
703	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_ys, cbdy_tys, nbdy_ys, nbdy_tys
704	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_ye, cbdy_tye, nbdy_ye, nbdy_tye
705
706	CALL interp_fcn (cfld, & ! CD field
707	cids, cide, ckds, ckde, cjds, cjde, &
708	cims, cime, ckms, ckme, cjms, cjme, &
709	cits, cite, ckts, ckte, cjts, cjte, &
710	nfld, & ! ND field
711	nids, nide, nkds, nkde, njds, njde, &
712	nims, nime, nkms, nkme, njms, njme, &
713	nits, nite, nkts, nkte, njts, njte, &
714	shw, & ! stencil half width for interp
715	imask, & ! interpolation mask
716	xstag, ystag, & ! staggering of field
717	ipos, jpos, & ! Position of lower left of nest in CD
718	nri, nrj ) ! nest ratios
719
720	END SUBROUTINE p2c
721
722	!==================================
723
724	SUBROUTINE bdy_interp ( cfld, & ! CD field
725	cids, cide, ckds, ckde, cjds, cjde, &
726	cims, cime, ckms, ckme, cjms, cjme, &
727	cits, cite, ckts, ckte, cjts, cjte, &
728	nfld, & ! ND field
729	nids, nide, nkds, nkde, njds, njde, &
730	nims, nime, nkms, nkme, njms, njme, &
731	nits, nite, nkts, nkte, njts, njte, &
732	shw, & ! stencil half width
733	imask, & ! interpolation mask
734	xstag, ystag, & ! staggering of field
735	ipos, jpos, & ! Position of lower left of nest in CD
736	nri, nrj, & ! nest ratios
737	cbdy_xs, nbdy_xs, &
738	cbdy_xe, nbdy_xe, &
739	cbdy_ys, nbdy_ys, &
740	cbdy_ye, nbdy_ye, &
741	cbdy_txs, nbdy_txs, &
742	cbdy_txe, nbdy_txe, &
743	cbdy_tys, nbdy_tys, &
744	cbdy_tye, nbdy_tye, &
745	cdt, ndt &
746	) ! boundary arrays
747	USE module_configure
748	IMPLICIT NONE
749
750	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
751	cims, cime, ckms, ckme, cjms, cjme, &
752	cits, cite, ckts, ckte, cjts, cjte, &
753	nids, nide, nkds, nkde, njds, njde, &
754	nims, nime, nkms, nkme, njms, njme, &
755	nits, nite, nkts, nkte, njts, njte, &
756	shw, &
757	ipos, jpos, &
758	nri, nrj
759
760	LOGICAL, INTENT(IN) :: xstag, ystag
761
762	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
763	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
764	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
765	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_xs, cbdy_txs, nbdy_xs, nbdy_txs
766	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_xe, cbdy_txe, nbdy_xe, nbdy_txe
767	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_ys, cbdy_tys, nbdy_ys, nbdy_tys
768	REAL, DIMENSION( * ), INTENT(INOUT) :: cbdy_ye, cbdy_tye, nbdy_ye, nbdy_tye
769	REAL cdt, ndt
770
771	! Local
772
773	INTEGER nijds, nijde, spec_bdy_width
774
775	nijds = min(nids, njds)
776	nijde = max(nide, njde)
777	CALL nl_get_spec_bdy_width( 1, spec_bdy_width )
778
779	CALL bdy_interp1( cfld, & ! CD field
780	cids, cide, ckds, ckde, cjds, cjde, &
781	cims, cime, ckms, ckme, cjms, cjme, &
782	cits, cite, ckts, ckte, cjts, cjte, &
783	nfld, & ! ND field
784	nijds, nijde , spec_bdy_width , &
785	nids, nide, nkds, nkde, njds, njde, &
786	nims, nime, nkms, nkme, njms, njme, &
787	nits, nite, nkts, nkte, njts, njte, &
788	shw, imask, &
789	xstag, ystag, & ! staggering of field
790	ipos, jpos, & ! Position of lower left of nest in CD
791	nri, nrj, &
792	cbdy_xs, nbdy_xs, &
793	cbdy_xe, nbdy_xe, &
794	cbdy_ys, nbdy_ys, &
795	cbdy_ye, nbdy_ye, &
796	cbdy_txs, nbdy_txs, &
797	cbdy_txe, nbdy_txe, &
798	cbdy_tys, nbdy_tys, &
799	cbdy_tye, nbdy_tye, &
800	cdt, ndt &
801	)
802
803	RETURN
804
805	END SUBROUTINE bdy_interp
806
807	SUBROUTINE bdy_interp1( cfld, & ! CD field
808	cids, cide, ckds, ckde, cjds, cjde, &
809	cims, cime, ckms, ckme, cjms, cjme, &
810	cits, cite, ckts, ckte, cjts, cjte, &
811	nfld, & ! ND field
812	nijds, nijde, spec_bdy_width , &
813	nids, nide, nkds, nkde, njds, njde, &
814	nims, nime, nkms, nkme, njms, njme, &
815	nits, nite, nkts, nkte, njts, njte, &
816	shw1, &
817	imask, & ! interpolation mask
818	xstag, ystag, & ! staggering of field
819	ipos, jpos, & ! Position of lower left of nest in CD
820	nri, nrj, &
821	cbdy_xs, bdy_xs, &
822	cbdy_xe, bdy_xe, &
823	cbdy_ys, bdy_ys, &
824	cbdy_ye, bdy_ye, &
825	cbdy_txs, bdy_txs, &
826	cbdy_txe, bdy_txe, &
827	cbdy_tys, bdy_tys, &
828	cbdy_tye, bdy_tye, &
829	cdt, ndt &
830	)
831
832	USE module_configure
833	use module_state_description
834	IMPLICIT NONE
835
836	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
837	cims, cime, ckms, ckme, cjms, cjme, &
838	cits, cite, ckts, ckte, cjts, cjte, &
839	nids, nide, nkds, nkde, njds, njde, &
840	nims, nime, nkms, nkme, njms, njme, &
841	nits, nite, nkts, nkte, njts, njte, &
842	shw1, & ! ignore
843	ipos, jpos, &
844	nri, nrj
845	INTEGER, INTENT(IN) :: nijds, nijde, spec_bdy_width
846	LOGICAL, INTENT(IN) :: xstag, ystag
847
848	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(INOUT) :: cfld
849	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ), INTENT(INOUT) :: nfld
850	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
851	REAL, DIMENSION ( * ), INTENT(INOUT) :: cbdy_xs, cbdy_txs ! not used
852	REAL, DIMENSION ( * ), INTENT(INOUT) :: cbdy_xe, cbdy_txe ! not used
853	REAL, DIMENSION ( * ), INTENT(INOUT) :: cbdy_ys, cbdy_tys ! not used
854	REAL, DIMENSION ( * ), INTENT(INOUT) :: cbdy_ye, cbdy_tye ! not used
855	REAL :: cdt, ndt
856	REAL, DIMENSION ( njms:njme, nkms:nkme, spec_bdy_width ), INTENT(INOUT) :: bdy_xs, bdy_txs
857	REAL, DIMENSION ( njms:njme, nkms:nkme, spec_bdy_width ), INTENT(INOUT) :: bdy_xe, bdy_txe
858	REAL, DIMENSION ( nims:nime, nkms:nkme, spec_bdy_width ), INTENT(INOUT) :: bdy_ys, bdy_tys
859	REAL, DIMENSION ( nims:nime, nkms:nkme, spec_bdy_width ), INTENT(INOUT) :: bdy_ye, bdy_tye
860
861	! Local
862
863	REAL*8 rdt
864	INTEGER ci, cj, ck, ni, nj, nk, ni1, nj1, nk1, ip, jp, ioff, joff
865	#ifdef MM5_SINT
866	INTEGER nfx, ior
867	PARAMETER (ior=2)
868	INTEGER nf
869	REAL psca1(cims:cime,cjms:cjme,nri*nrj)
870	REAL psca(cims:cime,cjms:cjme,nri*nrj)
871	LOGICAL icmask( cims:cime, cjms:cjme )
872	INTEGER i,j,k
873	#endif
874	INTEGER shw
875	INTEGER spec_zone
876	INTEGER relax_zone
877	INTEGER sz
878	INTEGER n2ci,n
879	INTEGER n2cj
880
881	! statement functions for converting a nest index to coarse
882	n2ci(n) = (n+ipos*nri-1)/nri
883	n2cj(n) = (n+jpos*nrj-1)/nrj
884
885	rdt = 1.D0/cdt
886
887	shw = 0
888
889	ioff = 0 ; joff = 0
890	IF ( xstag ) THEN
891	ioff = (nri-1)/2
892	ENDIF
893	IF ( ystag ) THEN
894	joff = (nrj-1)/2
895	ENDIF
896
897	! Iterate over the ND tile and compute the values
898	! from the CD tile.
899
900	#ifdef MM5_SINT
901	CALL nl_get_spec_zone( 1, spec_zone )
902	CALL nl_get_relax_zone( 1, relax_zone )
903	sz = MIN(MAX( spec_zone, relax_zone + 1 ),spec_bdy_width)
904
905	nfx = nri * nrj
906
907	!$OMP PARALLEL DO &
908	!$OMP PRIVATE ( i,j,k,ni,nj,ni1,nj1,ci,cj,ip,jp,nk,ck,nf,icmask,psca,psca1 )
909	DO k = ckts, ckte
910
911	DO nf = 1,nfx
912	DO j = cjms,cjme
913	nj = (j-jpos) * nrj + ( nrj / 2 + 1 ) ! j point on nest
914	DO i = cims,cime
915	ni = (i-ipos) * nri + ( nri / 2 + 1 ) ! i point on nest
916	psca1(i,j,nf) = cfld(i,k,j)
917	ENDDO
918	ENDDO
919	ENDDO
920	! hopefully less ham handed but still correct and more efficient
921	! sintb ignores icmask so it does not matter that icmask is not set
922	!
923	! SOUTH BDY
924	IF ( njts .ge. njds .and. njts .le. njds + sz + joff ) THEN
925	CALL sintb( psca1, psca, &
926	cims, cime, cjms, cjme, icmask, &
927	n2ci(nits)-1, n2ci(nite)+1, n2cj(MAX(njts,njds)), n2cj(MIN(njte,njds+sz+joff)), nrj*nri, xstag, ystag )
928	ENDIF
929	! NORTH BDY
930	IF ( njte .le. njde .and. njte .ge. njde - sz - joff ) THEN
931	CALL sintb( psca1, psca, &
932	cims, cime, cjms, cjme, icmask, &
933	n2ci(nits)-1, n2ci(nite)+1, n2cj(MAX(njts,njde-sz-joff)), n2cj(MIN(njte,njde-1+joff)), nrj*nri, xstag, ystag )
934	ENDIF
935	! WEST BDY
936	IF ( nits .ge. nids .and. nits .le. nids + sz + ioff ) THEN
937	CALL sintb( psca1, psca, &
938	cims, cime, cjms, cjme, icmask, &
939	n2ci(MAX(nits,nids)), n2ci(MIN(nite,nids+sz+ioff)), n2cj(njts)-1, n2cj(njte)+1, nrj*nri, xstag, ystag )
940	ENDIF
941	! EAST BDY
942	IF ( nite .le. nide .and. nite .ge. nide - sz - ioff ) THEN
943	CALL sintb( psca1, psca, &
944	cims, cime, cjms, cjme, icmask, &
945	n2ci(MAX(nits,nide-sz-ioff)), n2ci(MIN(nite,nide-1+ioff)), n2cj(njts)-1, n2cj(njte)+1, nrj*nri, xstag, ystag )
946	ENDIF
947
948	DO nj1 = MAX(njds,njts-1), MIN(njde+joff,njte+joff+1)
949	cj = jpos + (nj1-1) / nrj ! j coord of CD point
950	jp = mod ( nj1-1 , nrj ) ! coord of ND w/i CD point
951	nk = k
952	ck = nk
953	DO ni1 = MAX(nids,nits-1), MIN(nide+ioff,nite+ioff+1)
954	ci = ipos + (ni1-1) / nri ! j coord of CD point
955	ip = mod ( ni1-1 , nri ) ! coord of ND w/i CD point
956
957	ni = ni1-ioff
958	nj = nj1-joff
959
960	IF ( ( ni.LT.nids) .OR. (nj.LT.njds) ) THEN
961	CYCLE
962	END IF
963
964	!bdy contains the value at t-dt. psca contains the value at t
965	!compute dv/dt and store in bdy_t
966	!afterwards store the new value of v at t into bdy
967	! WEST
968	IF ( ni .ge. nids .and. ni .lt. nids + sz ) THEN
969	bdy_txs( nj,k,ni ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
970	bdy_xs( nj,k,ni ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
971	ENDIF
972
973	! SOUTH
974	IF ( nj .ge. njds .and. nj .lt. njds + sz ) THEN
975	bdy_tys( ni,k,nj ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
976	bdy_ys( ni,k,nj ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
977	ENDIF
978
979	! EAST
980	IF ( xstag ) THEN
981	IF ( ni .ge. nide - sz + 1 .AND. ni .le. nide ) THEN
982	bdy_txe( nj,k,nide-ni+1 ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
983	bdy_xe( nj,k,nide-ni+1 ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
984	ENDIF
985	ELSE
986	IF ( ni .ge. nide - sz .AND. ni .le. nide-1 ) THEN
987	bdy_txe( nj,k,nide-ni ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
988	bdy_xe( nj,k,nide-ni ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
989	ENDIF
990	ENDIF
991
992	! NORTH
993	IF ( ystag ) THEN
994	IF ( nj .ge. njde - sz + 1 .AND. nj .le. njde ) THEN
995	bdy_tye( ni,k,njde-nj+1 ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
996	bdy_ye( ni,k,njde-nj+1 ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
997	ENDIF
998	ELSE
999	IF ( nj .ge. njde - sz .AND. nj .le. njde-1 ) THEN
1000	bdy_tye(ni,k,njde-nj ) = rdt(psca(ci+shw,cj+shw,ip+1+(jp)nri)-nfld(ni,k,nj))
1001	bdy_ye( ni,k,njde-nj ) = psca(ci+shw,cj+shw,ip+1+(jp)*nri )
1002	ENDIF
1003	ENDIF
1004
1005	ENDDO
1006	ENDDO
1007	ENDDO
1008	!$OMP END PARALLEL DO
1009	#endif
1010
1011	RETURN
1012
1013	END SUBROUTINE bdy_interp1
1014
1015
1016
1017	SUBROUTINE interp_fcni( cfld, & ! CD field
1018	cids, cide, ckds, ckde, cjds, cjde, &
1019	cims, cime, ckms, ckme, cjms, cjme, &
1020	cits, cite, ckts, ckte, cjts, cjte, &
1021	nfld, & ! ND field
1022	nids, nide, nkds, nkde, njds, njde, &
1023	nims, nime, nkms, nkme, njms, njme, &
1024	nits, nite, nkts, nkte, njts, njte, &
1025	shw, & ! stencil half width
1026	imask, & ! interpolation mask
1027	xstag, ystag, & ! staggering of field
1028	ipos, jpos, & ! Position of lower left of nest in CD
1029	nri, nrj ) ! nest ratios
1030	USE module_configure
1031	IMPLICIT NONE
1032
1033
1034	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1035	cims, cime, ckms, ckme, cjms, cjme, &
1036	cits, cite, ckts, ckte, cjts, cjte, &
1037	nids, nide, nkds, nkde, njds, njde, &
1038	nims, nime, nkms, nkme, njms, njme, &
1039	nits, nite, nkts, nkte, njts, njte, &
1040	shw, &
1041	ipos, jpos, &
1042	nri, nrj
1043	LOGICAL, INTENT(IN) :: xstag, ystag
1044
1045	INTEGER, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1046	INTEGER, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
1047	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
1048
1049	! Local
1050
1051	INTEGER ci, cj, ck, ni, nj, nk, ip, jp
1052
1053	! Iterate over the ND tile and compute the values
1054	! from the CD tile.
1055
1056	!write(0,'("cits:cite, ckts:ckte, cjts:cjte ",6i4)')cits,cite, ckts,ckte, cjts,cjte
1057	!write(0,'("nits:nite, nkts:nkte, njts:njte ",6i4)')nits,nite, nkts,nkte, njts,njte
1058
1059	DO nj = njts, njte
1060	cj = jpos + (nj-1) / nrj ! j coord of CD point
1061	jp = mod ( nj , nrj ) ! coord of ND w/i CD point
1062	DO nk = nkts, nkte
1063	ck = nk
1064	DO ni = nits, nite
1065	ci = ipos + (ni-1) / nri ! j coord of CD point
1066	ip = mod ( ni , nri ) ! coord of ND w/i CD point
1067	! This is a trivial implementation of the interp_fcn; just copies
1068	! the values from the CD into the ND
1069	nfld( ni, nk, nj ) = cfld( ci , ck , cj )
1070	ENDDO
1071	ENDDO
1072	ENDDO
1073
1074	RETURN
1075
1076	END SUBROUTINE interp_fcni
1077
1078	SUBROUTINE interp_fcnm( cfld, & ! CD field
1079	cids, cide, ckds, ckde, cjds, cjde, &
1080	cims, cime, ckms, ckme, cjms, cjme, &
1081	cits, cite, ckts, ckte, cjts, cjte, &
1082	nfld, & ! ND field
1083	nids, nide, nkds, nkde, njds, njde, &
1084	nims, nime, nkms, nkme, njms, njme, &
1085	nits, nite, nkts, nkte, njts, njte, &
1086	shw, & ! stencil half width
1087	imask, & ! interpolation mask
1088	xstag, ystag, & ! staggering of field
1089	ipos, jpos, & ! Position of lower left of nest in CD
1090	nri, nrj ) ! nest ratios
1091	USE module_configure
1092	IMPLICIT NONE
1093
1094
1095	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1096	cims, cime, ckms, ckme, cjms, cjme, &
1097	cits, cite, ckts, ckte, cjts, cjte, &
1098	nids, nide, nkds, nkde, njds, njde, &
1099	nims, nime, nkms, nkme, njms, njme, &
1100	nits, nite, nkts, nkte, njts, njte, &
1101	shw, &
1102	ipos, jpos, &
1103	nri, nrj
1104	LOGICAL, INTENT(IN) :: xstag, ystag
1105
1106	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1107	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
1108	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
1109
1110	! Local
1111
1112	INTEGER ci, cj, ck, ni, nj, nk, ip, jp
1113
1114	! Iterate over the ND tile and compute the values
1115	! from the CD tile.
1116
1117	!write(0,'("mask cits:cite, ckts:ckte, cjts:cjte ",6i4)')cits,cite, ckts,ckte, cjts,cjte
1118	!write(0,'("mask nits:nite, nkts:nkte, njts:njte ",6i4)')nits,nite, nkts,nkte, njts,njte
1119
1120	DO nj = njts, njte
1121	cj = jpos + (nj-1) / nrj ! j coord of CD point
1122	jp = mod ( nj , nrj ) ! coord of ND w/i CD point
1123	DO nk = nkts, nkte
1124	ck = nk
1125	DO ni = nits, nite
1126	ci = ipos + (ni-1) / nri ! j coord of CD point
1127	ip = mod ( ni , nri ) ! coord of ND w/i CD point
1128	! This is a trivial implementation of the interp_fcn; just copies
1129	! the values from the CD into the ND
1130	nfld( ni, nk, nj ) = cfld( ci , ck , cj )
1131	ENDDO
1132	ENDDO
1133	ENDDO
1134
1135	RETURN
1136
1137	END SUBROUTINE interp_fcnm
1138
1139	SUBROUTINE interp_mask_land_field ( enable, & ! says whether to allow interpolation or just the bcasts
1140	cfld, & ! CD field
1141	cids, cide, ckds, ckde, cjds, cjde, &
1142	cims, cime, ckms, ckme, cjms, cjme, &
1143	cits, cite, ckts, ckte, cjts, cjte, &
1144	nfld, & ! ND field
1145	nids, nide, nkds, nkde, njds, njde, &
1146	nims, nime, nkms, nkme, njms, njme, &
1147	nits, nite, nkts, nkte, njts, njte, &
1148	shw, & ! stencil half width
1149	imask, & ! interpolation mask
1150	xstag, ystag, & ! staggering of field
1151	ipos, jpos, & ! Position of lower left of nest in CD
1152	nri, nrj, & ! nest ratios
1153	clu, nlu )
1154
1155	USE module_configure
1156	USE module_wrf_error
1157
1158	IMPLICIT NONE
1159
1160
1161	LOGICAL, INTENT(IN) :: enable
1162	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1163	cims, cime, ckms, ckme, cjms, cjme, &
1164	cits, cite, ckts, ckte, cjts, cjte, &
1165	nids, nide, nkds, nkde, njds, njde, &
1166	nims, nime, nkms, nkme, njms, njme, &
1167	nits, nite, nkts, nkte, njts, njte, &
1168	shw, &
1169	ipos, jpos, &
1170	nri, nrj
1171	LOGICAL, INTENT(IN) :: xstag, ystag
1172
1173	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1174	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
1175	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
1176
1177	REAL, DIMENSION ( cims:cime, cjms:cjme ) :: clu
1178	REAL, DIMENSION ( nims:nime, njms:njme ) :: nlu
1179
1180	! Local
1181
1182	INTEGER ci, cj, ck, ni, nj, nk, ip, jp
1183	INTEGER :: icount , ii , jj , ist , ien , jst , jen , iswater
1184	REAL :: avg , sum , dx , dy
1185	INTEGER , PARAMETER :: max_search = 5
1186	CHARACTER*120 message
1187
1188	! Find out what the water value is.
1189
1190	CALL nl_get_iswater(1,iswater)
1191
1192	! Right now, only mass point locations permitted.
1193
1194	IF ( ( .NOT. xstag) .AND. ( .NOT. ystag ) ) THEN
1195
1196	! Loop over each i,k,j in the nested domain.
1197
1198	IF ( enable ) THEN
1199
1200	DO nj = njts, njte
1201	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1202	cj = ( nj + (nrj/2)-1 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1203	ELSE
1204	cj = ( nj + (nrj-1)/2 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1205	END IF
1206	DO nk = nkts, nkte
1207	ck = nk
1208	DO ni = nits, nite
1209	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1210	ci = ( ni + (nri/2)-1 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1211	ELSE
1212	ci = ( ni + (nri-1)/2 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1213	END IF
1214
1215
1216
1217
1218	!
1219	! (ci,cj+1) (ci+1,cj+1)
1220	! - -------------
1221	! 1-dy \| \| \|
1222	! \| \| \|
1223	! - \| * \|
1224	! dy \| \| (ni,nj) \|
1225	! \| \| \|
1226	! - -------------
1227	! (ci,cj) (ci+1,cj)
1228	!
1229	! \|--\|--------\|
1230	! dx 1-dx
1231
1232
1233	! For odd ratios, at (nri+1)/2, we are on the coarse grid point, so dx = 0
1234
1235	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1236	dx = ( REAL ( MOD ( ni+(nri-1)/2 , nri ) ) + 0.5 ) / REAL ( nri )
1237	ELSE
1238	dx = REAL ( MOD ( ni+(nri-1)/2 , nri ) ) / REAL ( nri )
1239	END IF
1240	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1241	dy = ( REAL ( MOD ( nj+(nrj-1)/2 , nrj ) ) + 0.5 ) / REAL ( nrj )
1242	ELSE
1243	dy = REAL ( MOD ( nj+(nrj-1)/2 , nrj ) ) / REAL ( nrj )
1244	END IF
1245
1246	! This is a "land only" field. If this is a water point, no operations required.
1247
1248	IF ( ( NINT(nlu(ni ,nj )) .EQ. iswater ) ) THEN
1249	! noop
1250	! nfld(ni,nk,nj) = 1.e20
1251	nfld(ni,nk,nj) = -1
1252
1253	! If this is a nested land point, and the surrounding coarse values are all land points,
1254	! then this is a simple 4-pt interpolation.
1255
1256	ELSE IF ( ( NINT(nlu(ni ,nj )) .NE. iswater ) .AND. &
1257	( NINT(clu(ci ,cj )) .NE. iswater ) .AND. &
1258	( NINT(clu(ci+1,cj )) .NE. iswater ) .AND. &
1259	( NINT(clu(ci ,cj+1)) .NE. iswater ) .AND. &
1260	( NINT(clu(ci+1,cj+1)) .NE. iswater ) ) THEN
1261	nfld(ni,nk,nj) = ( 1. - dx ) * ( ( 1. - dy ) * cfld(ci ,ck,cj ) + &
1262	dy * cfld(ci ,ck,cj+1) ) + &
1263	dx * ( ( 1. - dy ) * cfld(ci+1,ck,cj ) + &
1264	dy * cfld(ci+1,ck,cj+1) )
1265
1266	! If this is a nested land point and there are NO coarse land values surrounding,
1267	! we temporarily punt.
1268
1269	ELSE IF ( ( NINT(nlu(ni ,nj )) .NE. iswater ) .AND. &
1270	( NINT(clu(ci ,cj )) .EQ. iswater ) .AND. &
1271	( NINT(clu(ci+1,cj )) .EQ. iswater ) .AND. &
1272	( NINT(clu(ci ,cj+1)) .EQ. iswater ) .AND. &
1273	( NINT(clu(ci+1,cj+1)) .EQ. iswater ) ) THEN
1274	! nfld(ni,nk,nj) = -1.e20
1275	nfld(ni,nk,nj) = -1
1276
1277	! If there are some water points and some land points, take an average.
1278
1279	ELSE IF ( NINT(nlu(ni ,nj )) .NE. iswater ) THEN
1280	icount = 0
1281	sum = 0
1282	IF ( NINT(clu(ci ,cj )) .NE. iswater ) THEN
1283	icount = icount + 1
1284	sum = sum + cfld(ci ,ck,cj )
1285	END IF
1286	IF ( NINT(clu(ci+1,cj )) .NE. iswater ) THEN
1287	icount = icount + 1
1288	sum = sum + cfld(ci+1,ck,cj )
1289	END IF
1290	IF ( NINT(clu(ci ,cj+1)) .NE. iswater ) THEN
1291	icount = icount + 1
1292	sum = sum + cfld(ci ,ck,cj+1)
1293	END IF
1294	IF ( NINT(clu(ci+1,cj+1)) .NE. iswater ) THEN
1295	icount = icount + 1
1296	sum = sum + cfld(ci+1,ck,cj+1)
1297	END IF
1298	nfld(ni,nk,nj) = sum / REAL ( icount )
1299	END IF
1300	END DO
1301	END DO
1302	END DO
1303
1304
1305	! Get an average of the whole domain for problem locations.
1306
1307	sum = 0
1308	icount = 0
1309	DO nj = njts, njte
1310	DO nk = nkts, nkte
1311	DO ni = nits, nite
1312	IF ( ( nfld(ni,nk,nj) .GT. -1.e19 ) .AND. ( nfld(ni,nk,nj) .LT. 1.e19 ) ) THEN
1313	icount = icount + 1
1314	sum = sum + nfld(ni,nk,nj)
1315	END IF
1316	END DO
1317	END DO
1318	END DO
1319	ELSE
1320	sum = 0.
1321	icount = 0
1322	ENDIF
1323	CALL wrf_dm_bcast_real( sum, 1 )
1324	CALL wrf_dm_bcast_integer( icount, 1 )
1325	IF ( enable ) THEN
1326	IF ( icount .GT. 0 ) THEN
1327	avg = sum / REAL ( icount )
1328
1329	! OK, if there were any of those island situations, we try to search a bit broader
1330	! of an area in the coarse grid.
1331
1332	DO nj = njts, njte
1333	DO nk = nkts, nkte
1334	DO ni = nits, nite
1335	IF ( nfld(ni,nk,nj) .LT. -1.e19 ) THEN
1336	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1337	cj = ( nj + (nrj/2)-1 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1338	ELSE
1339	cj = ( nj + (nrj-1)/2 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1340	END IF
1341	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1342	ci = ( ni + (nri/2)-1 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1343	ELSE
1344	ci = ( ni + (nri-1)/2 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1345	END IF
1346	ist = MAX (ci-max_search,cits)
1347	ien = MIN (ci+max_search,cite,cide-1)
1348	jst = MAX (cj-max_search,cjts)
1349	jen = MIN (cj+max_search,cjte,cjde-1)
1350	icount = 0
1351	sum = 0
1352	DO jj = jst,jen
1353	DO ii = ist,ien
1354	IF ( NINT(clu(ii,jj)) .NE. iswater ) THEN
1355	icount = icount + 1
1356	sum = sum + cfld(ii,nk,jj)
1357	END IF
1358	END DO
1359	END DO
1360	IF ( icount .GT. 0 ) THEN
1361	nfld(ni,nk,nj) = sum / REAL ( icount )
1362	ELSE
1363	! CALL wrf_error_fatal ( "horizontal interp error - island" )
1364	write(message,*) 'horizontal interp error - island, using average ', avg
1365	CALL wrf_message ( message )
1366	nfld(ni,nk,nj) = avg
1367	END IF
1368	END IF
1369	END DO
1370	END DO
1371	END DO
1372	ENDIF
1373	ENDIF
1374	ELSE
1375	CALL wrf_error_fatal ( "only unstaggered fields right now" )
1376	END IF
1377
1378	END SUBROUTINE interp_mask_land_field
1379
1380	SUBROUTINE interp_mask_water_field ( enable, & ! says whether to allow interpolation or just the bcasts
1381	cfld, & ! CD field
1382	cids, cide, ckds, ckde, cjds, cjde, &
1383	cims, cime, ckms, ckme, cjms, cjme, &
1384	cits, cite, ckts, ckte, cjts, cjte, &
1385	nfld, & ! ND field
1386	nids, nide, nkds, nkde, njds, njde, &
1387	nims, nime, nkms, nkme, njms, njme, &
1388	nits, nite, nkts, nkte, njts, njte, &
1389	shw, & ! stencil half width
1390	imask, & ! interpolation mask
1391	xstag, ystag, & ! staggering of field
1392	ipos, jpos, & ! Position of lower left of nest in CD
1393	nri, nrj, & ! nest ratios
1394	clu, nlu, cflag, nflag )
1395
1396	USE module_configure
1397	USE module_wrf_error
1398
1399	IMPLICIT NONE
1400
1401
1402	LOGICAL, INTENT(IN) :: enable
1403	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1404	cims, cime, ckms, ckme, cjms, cjme, &
1405	cits, cite, ckts, ckte, cjts, cjte, &
1406	nids, nide, nkds, nkde, njds, njde, &
1407	nims, nime, nkms, nkme, njms, njme, &
1408	nits, nite, nkts, nkte, njts, njte, &
1409	shw, &
1410	ipos, jpos, &
1411	nri, nrj, cflag, nflag
1412	LOGICAL, INTENT(IN) :: xstag, ystag
1413
1414	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1415	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
1416	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
1417
1418	REAL, DIMENSION ( cims:cime, cjms:cjme ) :: clu
1419	REAL, DIMENSION ( nims:nime, njms:njme ) :: nlu
1420
1421	! Local
1422
1423	INTEGER ci, cj, ck, ni, nj, nk, ip, jp
1424	INTEGER :: icount , ii , jj , ist , ien , jst , jen
1425	REAL :: avg , sum , dx , dy
1426	INTEGER , PARAMETER :: max_search = 5
1427
1428	! Right now, only mass point locations permitted.
1429
1430	IF ( ( .NOT. xstag) .AND. ( .NOT. ystag ) ) THEN
1431
1432	IF ( enable ) THEN
1433	! Loop over each i,k,j in the nested domain.
1434
1435	DO nj = njts, njte
1436	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1437	cj = ( nj + (nrj/2)-1 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1438	ELSE
1439	cj = ( nj + (nrj-1)/2 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1440	END IF
1441	DO nk = nkts, nkte
1442	ck = nk
1443	DO ni = nits, nite
1444	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1445	ci = ( ni + (nri/2)-1 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1446	ELSE
1447	ci = ( ni + (nri-1)/2 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1448	END IF
1449
1450
1451
1452
1453	!
1454	! (ci,cj+1) (ci+1,cj+1)
1455	! - -------------
1456	! 1-dy \| \| \|
1457	! \| \| \|
1458	! - \| * \|
1459	! dy \| \| (ni,nj) \|
1460	! \| \| \|
1461	! - -------------
1462	! (ci,cj) (ci+1,cj)
1463	!
1464	! \|--\|--------\|
1465	! dx 1-dx
1466
1467
1468	! At ni=2, we are on the coarse grid point, so dx = 0
1469
1470	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1471	dx = ( REAL ( MOD ( ni+(nri-1)/2 , nri ) ) + 0.5 ) / REAL ( nri )
1472	ELSE
1473	dx = REAL ( MOD ( ni+(nri-1)/2 , nri ) ) / REAL ( nri )
1474	END IF
1475	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1476	dy = ( REAL ( MOD ( nj+(nrj-1)/2 , nrj ) ) + 0.5 ) / REAL ( nrj )
1477	ELSE
1478	dy = REAL ( MOD ( nj+(nrj-1)/2 , nrj ) ) / REAL ( nrj )
1479	END IF
1480
1481	! This is a "water only" field. If this is a land point, no operations required.
1482
1483	IF ( ( NINT(nlu(ni ,nj )) .NE. nflag ) ) THEN
1484	! noop
1485	! nfld(ni,nk,nj) = 1.e20
1486	nfld(ni,nk,nj) = 0
1487
1488	! If this is a nested water point, and the surrounding coarse values are all water points,
1489	! then this is a simple 4-pt interpolation.
1490
1491	ELSE IF ( ( NINT(nlu(ni ,nj )) .EQ. nflag ) .AND. &
1492	( NINT(clu(ci ,cj )) .EQ. nflag ) .AND. &
1493	( NINT(clu(ci+1,cj )) .EQ. nflag ) .AND. &
1494	( NINT(clu(ci ,cj+1)) .EQ. nflag ) .AND. &
1495	( NINT(clu(ci+1,cj+1)) .EQ. nflag ) ) THEN
1496	nfld(ni,nk,nj) = ( 1. - dx ) * ( ( 1. - dy ) * cfld(ci ,ck,cj ) + &
1497	dy * cfld(ci ,ck,cj+1) ) + &
1498	dx * ( ( 1. - dy ) * cfld(ci+1,ck,cj ) + &
1499	dy * cfld(ci+1,ck,cj+1) )
1500
1501	! If this is a nested water point and there are NO coarse water values surrounding,
1502	! we temporarily punt.
1503
1504	ELSE IF ( ( NINT(nlu(ni ,nj )) .EQ. nflag ) .AND. &
1505	( NINT(clu(ci ,cj )) .NE. nflag ) .AND. &
1506	( NINT(clu(ci+1,cj )) .NE. nflag ) .AND. &
1507	( NINT(clu(ci ,cj+1)) .NE. nflag ) .AND. &
1508	( NINT(clu(ci+1,cj+1)) .NE. nflag ) ) THEN
1509	! nfld(ni,nk,nj) = -1.e20
1510	nfld(ni,nk,nj) = -1
1511
1512	! If there are some land points and some water points, take an average.
1513
1514	ELSE IF ( NINT(nlu(ni ,nj )) .EQ. nflag ) THEN
1515	icount = 0
1516	sum = 0
1517	IF ( NINT(clu(ci ,cj )) .EQ. nflag ) THEN
1518	icount = icount + 1
1519	sum = sum + cfld(ci ,ck,cj )
1520	END IF
1521	IF ( NINT(clu(ci+1,cj )) .EQ. nflag ) THEN
1522	icount = icount + 1
1523	sum = sum + cfld(ci+1,ck,cj )
1524	END IF
1525	IF ( NINT(clu(ci ,cj+1)) .EQ. nflag ) THEN
1526	icount = icount + 1
1527	sum = sum + cfld(ci ,ck,cj+1)
1528	END IF
1529	IF ( NINT(clu(ci+1,cj+1)) .EQ. nflag ) THEN
1530	icount = icount + 1
1531	sum = sum + cfld(ci+1,ck,cj+1)
1532	END IF
1533	nfld(ni,nk,nj) = sum / REAL ( icount )
1534	END IF
1535	END DO
1536	END DO
1537	END DO
1538
1539	! Get an average of the whole domain for problem locations.
1540
1541	sum = 0
1542	icount = 0
1543	DO nj = njts, njte
1544	DO nk = nkts, nkte
1545	DO ni = nits, nite
1546	IF ( ( nfld(ni,nk,nj) .GT. -1.e19 ) .AND. ( nfld(ni,nk,nj) .LT. 1.e19 ) ) THEN
1547	icount = icount + 1
1548	sum = sum + nfld(ni,nk,nj)
1549	END IF
1550	END DO
1551	END DO
1552	END DO
1553	ELSE
1554	sum = 0.
1555	icount = 0
1556	ENDIF
1557	CALL wrf_dm_bcast_real( sum, 1 )
1558	CALL wrf_dm_bcast_integer( icount, 1 )
1559	IF ( enable ) THEN
1560	IF ( icount .NE. 0 ) THEN
1561	avg = sum / REAL ( icount )
1562
1563
1564	! OK, if there were any of those lake situations, we try to search a bit broader
1565	! of an area in the coarse grid.
1566
1567	DO nj = njts, njte
1568	DO nk = nkts, nkte
1569	DO ni = nits, nite
1570	IF ( nfld(ni,nk,nj) .LT. -1.e19 ) THEN
1571	IF ( MOD ( nrj , 2 ) .EQ. 0 ) THEN
1572	cj = ( nj + (nrj/2)-1 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1573	ELSE
1574	cj = ( nj + (nrj-1)/2 ) / nrj + jpos -1 ! first coarse position equal to or below nest point
1575	END IF
1576	IF ( MOD ( nri , 2 ) .EQ. 0 ) THEN
1577	ci = ( ni + (nri/2)-1 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1578	ELSE
1579	ci = ( ni + (nri-1)/2 ) / nri + ipos -1 ! first coarse position equal to or to the left of nest point
1580	END IF
1581	ist = MAX (ci-max_search,cits)
1582	ien = MIN (ci+max_search,cite,cide-1)
1583	jst = MAX (cj-max_search,cjts)
1584	jen = MIN (cj+max_search,cjte,cjde-1)
1585	icount = 0
1586	sum = 0
1587	DO jj = jst,jen
1588	DO ii = ist,ien
1589	IF ( NINT(clu(ii,jj)) .EQ. nflag ) THEN
1590	icount = icount + 1
1591	sum = sum + cfld(ii,nk,jj)
1592	END IF
1593	END DO
1594	END DO
1595	IF ( icount .GT. 0 ) THEN
1596	nfld(ni,nk,nj) = sum / REAL ( icount )
1597	ELSE
1598	! CALL wrf_error_fatal ( "horizontal interp error - lake" )
1599	print *,'horizontal interp error - lake, using average ',avg
1600	nfld(ni,nk,nj) = avg
1601	END IF
1602	END IF
1603	END DO
1604	END DO
1605	END DO
1606	ENDIF
1607	ENDIF
1608	ELSE
1609	CALL wrf_error_fatal ( "only unstaggered fields right now" )
1610	END IF
1611
1612	END SUBROUTINE interp_mask_water_field
1613
1614	SUBROUTINE none
1615	END SUBROUTINE none
1616
1617	SUBROUTINE smoother ( cfld , &
1618	cids, cide, ckds, ckde, cjds, cjde, &
1619	cims, cime, ckms, ckme, cjms, cjme, &
1620	cits, cite, ckts, ckte, cjts, cjte, &
1621	nids, nide, nkds, nkde, njds, njde, &
1622	nims, nime, nkms, nkme, njms, njme, &
1623	nits, nite, nkts, nkte, njts, njte, &
1624	xstag, ystag, & ! staggering of field
1625	ipos, jpos, & ! Position of lower left of nest in
1626	nri, nrj &
1627	)
1628
1629	USE module_configure
1630	IMPLICIT NONE
1631
1632	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1633	cims, cime, ckms, ckme, cjms, cjme, &
1634	cits, cite, ckts, ckte, cjts, cjte, &
1635	nids, nide, nkds, nkde, njds, njde, &
1636	nims, nime, nkms, nkme, njms, njme, &
1637	nits, nite, nkts, nkte, njts, njte, &
1638	nri, nrj, &
1639	ipos, jpos
1640	LOGICAL, INTENT(IN) :: xstag, ystag
1641	INTEGER :: smooth_option, feedback , spec_zone
1642
1643	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1644
1645	! If there is no feedback, there can be no smoothing.
1646
1647	CALL nl_get_feedback ( 1, feedback )
1648	IF ( feedback == 0 ) RETURN
1649	CALL nl_get_spec_zone ( 1, spec_zone )
1650
1651	! These are the 2d smoothers used on the fedback data. These filters
1652	! are run on the coarse grid data (after the nested info has been
1653	! fedback). Only the area of the nest in the coarse grid is filtered.
1654
1655	CALL nl_get_smooth_option ( 1, smooth_option )
1656
1657	IF ( smooth_option == 0 ) THEN
1658	! no op
1659	ELSE IF ( smooth_option == 1 ) THEN
1660	CALL sm121 ( cfld , &
1661	cids, cide, ckds, ckde, cjds, cjde, &
1662	cims, cime, ckms, ckme, cjms, cjme, &
1663	cits, cite, ckts, ckte, cjts, cjte, &
1664	xstag, ystag, & ! staggering of field
1665	nids, nide, nkds, nkde, njds, njde, &
1666	nims, nime, nkms, nkme, njms, njme, &
1667	nits, nite, nkts, nkte, njts, njte, &
1668	nri, nrj, &
1669	ipos, jpos & ! Position of lower left of nest in
1670	)
1671	ELSE IF ( smooth_option == 2 ) THEN
1672	CALL smdsm ( cfld , &
1673	cids, cide, ckds, ckde, cjds, cjde, &
1674	cims, cime, ckms, ckme, cjms, cjme, &
1675	cits, cite, ckts, ckte, cjts, cjte, &
1676	xstag, ystag, & ! staggering of field
1677	nids, nide, nkds, nkde, njds, njde, &
1678	nims, nime, nkms, nkme, njms, njme, &
1679	nits, nite, nkts, nkte, njts, njte, &
1680	nri, nrj, &
1681	ipos, jpos & ! Position of lower left of nest in
1682	)
1683	END IF
1684
1685	END SUBROUTINE smoother
1686
1687	SUBROUTINE sm121 ( cfld , &
1688	cids, cide, ckds, ckde, cjds, cjde, &
1689	cims, cime, ckms, ckme, cjms, cjme, &
1690	cits, cite, ckts, ckte, cjts, cjte, &
1691	xstag, ystag, & ! staggering of field
1692	nids, nide, nkds, nkde, njds, njde, &
1693	nims, nime, nkms, nkme, njms, njme, &
1694	nits, nite, nkts, nkte, njts, njte, &
1695	nri, nrj, &
1696	ipos, jpos & ! Position of lower left of nest in
1697	)
1698
1699	USE module_configure
1700	IMPLICIT NONE
1701
1702	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1703	cims, cime, ckms, ckme, cjms, cjme, &
1704	cits, cite, ckts, ckte, cjts, cjte, &
1705	nids, nide, nkds, nkde, njds, njde, &
1706	nims, nime, nkms, nkme, njms, njme, &
1707	nits, nite, nkts, nkte, njts, njte, &
1708	nri, nrj, &
1709	ipos, jpos
1710	LOGICAL, INTENT(IN) :: xstag, ystag
1711
1712	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1713	REAL, DIMENSION ( cims:cime, cjms:cjme ) :: cfldnew
1714
1715	INTEGER :: i , j , k , loop
1716	INTEGER :: istag,jstag
1717
1718	INTEGER, PARAMETER :: smooth_passes = 1 ! More passes requires a larger stencil (currently 48 pt)
1719
1720	istag = 1 ; jstag = 1
1721	IF ( xstag ) istag = 0
1722	IF ( ystag ) jstag = 0
1723
1724	! Simple 1-2-1 smoother.
1725
1726	smoothing_passes : DO loop = 1 , smooth_passes
1727
1728	DO k = ckts , ckte
1729
1730	! Initialize dummy cfldnew
1731
1732	DO i = MAX(ipos,cits-3) , MIN(ipos+(nide-nids)/nri,cite+3)
1733	DO j = MAX(jpos,cjts-3) , MIN(jpos+(njde-njds)/nrj,cjte+3)
1734	cfldnew(i,j) = cfld(i,k,j)
1735	END DO
1736	END DO
1737
1738	! 1-2-1 smoothing in the j direction first,
1739
1740	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1741	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1742	cfldnew(i,j) = 0.25 * ( cfld(i,k,j+1) + 2.*cfld(i,k,j) + cfld(i,k,j-1) )
1743	END DO
1744	END DO
1745
1746	! then 1-2-1 smoothing in the i direction last
1747
1748	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1749	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1750	cfld(i,k,j) = 0.25 * ( cfldnew(i+1,j) + 2.*cfldnew(i,j) + cfldnew(i-1,j) )
1751	END DO
1752	END DO
1753
1754	END DO
1755
1756	END DO smoothing_passes
1757
1758	END SUBROUTINE sm121
1759
1760	SUBROUTINE smdsm ( cfld , &
1761	cids, cide, ckds, ckde, cjds, cjde, &
1762	cims, cime, ckms, ckme, cjms, cjme, &
1763	cits, cite, ckts, ckte, cjts, cjte, &
1764	xstag, ystag, & ! staggering of field
1765	nids, nide, nkds, nkde, njds, njde, &
1766	nims, nime, nkms, nkme, njms, njme, &
1767	nits, nite, nkts, nkte, njts, njte, &
1768	nri, nrj, &
1769	ipos, jpos & ! Position of lower left of nest in
1770	)
1771
1772	USE module_configure
1773	IMPLICIT NONE
1774
1775	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1776	cims, cime, ckms, ckme, cjms, cjme, &
1777	cits, cite, ckts, ckte, cjts, cjte, &
1778	nids, nide, nkds, nkde, njds, njde, &
1779	nims, nime, nkms, nkme, njms, njme, &
1780	nits, nite, nkts, nkte, njts, njte, &
1781	nri, nrj, &
1782	ipos, jpos
1783	LOGICAL, INTENT(IN) :: xstag, ystag
1784
1785	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
1786	REAL, DIMENSION ( cims:cime, cjms:cjme ) :: cfldnew
1787
1788	REAL , DIMENSION ( 2 ) :: xnu
1789	INTEGER :: i , j , k , loop , n
1790	INTEGER :: istag,jstag
1791
1792	INTEGER, PARAMETER :: smooth_passes = 1 ! More passes requires a larger stencil (currently 48 pt)
1793
1794	xnu = (/ 0.50 , -0.52 /)
1795
1796	istag = 1 ; jstag = 1
1797	IF ( xstag ) istag = 0
1798	IF ( ystag ) jstag = 0
1799
1800	! The odd number passes of this are the "smoother", the even
1801	! number passes are the "de-smoother" (note the different signs on xnu).
1802
1803	smoothing_passes : DO loop = 1 , smooth_passes * 2
1804
1805	n = 2 - MOD ( loop , 2 )
1806
1807	DO k = ckts , ckte
1808
1809	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1810	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1811	cfldnew(i,j) = cfld(i,k,j) + xnu(n) * ((cfld(i,k,j+1) + cfld(i,k,j-1)) * 0.5-cfld(i,k,j))
1812	END DO
1813	END DO
1814
1815	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1816	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1817	cfld(i,k,j) = cfldnew(i,j)
1818	END DO
1819	END DO
1820
1821	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1822	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1823	cfldnew(i,j) = cfld(i,k,j) + xnu(n) * ((cfld(i+1,k,j) + cfld(i-1,k,j)) * 0.5-cfld(i,k,j))
1824	END DO
1825	END DO
1826
1827	DO j = MAX(jpos+2,cjts-2) , MIN(jpos+(njde-njds)/nrj-2-jstag,cjte+2)
1828	DO i = MAX(ipos+2,cits-2) , MIN(ipos+(nide-nids)/nri-2-istag,cite+2)
1829	cfld(i,k,j) = cfldnew(i,j)
1830	END DO
1831	END DO
1832
1833	END DO
1834
1835	END DO smoothing_passes
1836
1837	END SUBROUTINE smdsm
1838
1839	!==================================
1840	! this is used to modify a field over the nest so we can see where the nest is
1841
1842	SUBROUTINE mark_domain ( cfld, & ! CD field
1843	cids, cide, ckds, ckde, cjds, cjde, &
1844	cims, cime, ckms, ckme, cjms, cjme, &
1845	cits, cite, ckts, ckte, cjts, cjte, &
1846	nfld, & ! ND field
1847	nids, nide, nkds, nkde, njds, njde, &
1848	nims, nime, nkms, nkme, njms, njme, &
1849	nits, nite, nkts, nkte, njts, njte, &
1850	shw, & ! stencil half width for interp
1851	imask, & ! interpolation mask
1852	xstag, ystag, & ! staggering of field
1853	ipos, jpos, & ! Position of lower left of nest in CD
1854	nri, nrj ) ! nest ratios
1855	USE module_configure
1856	USE module_wrf_error
1857	IMPLICIT NONE
1858
1859
1860	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1861	cims, cime, ckms, ckme, cjms, cjme, &
1862	cits, cite, ckts, ckte, cjts, cjte, &
1863	nids, nide, nkds, nkde, njds, njde, &
1864	nims, nime, nkms, nkme, njms, njme, &
1865	nits, nite, nkts, nkte, njts, njte, &
1866	shw, &
1867	ipos, jpos, &
1868	nri, nrj
1869	LOGICAL, INTENT(IN) :: xstag, ystag
1870
1871	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(OUT) :: cfld
1872	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ), INTENT(IN) :: nfld
1873	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: imask
1874
1875	! Local
1876
1877	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
1878	INTEGER :: icmin,icmax,jcmin,jcmax
1879	INTEGER :: istag,jstag, ipoints,jpoints,ijpoints
1880
1881	istag = 1 ; jstag = 1
1882	IF ( xstag ) istag = 0
1883	IF ( ystag ) jstag = 0
1884
1885	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-jstag-1,cjte)
1886	nj = (cj-jpos)*nrj + jstag + 1
1887	DO ck = ckts, ckte
1888	nk = ck
1889	DO ci = MAX(ipos+1,cits),MIN(ipos+(nide-nids)/nri-istag-1,cite)
1890	ni = (ci-ipos)*nri + istag + 1
1891	cfld( ci, ck, cj ) = 9021000. !magic number: Beverly Hills * 100.
1892	ENDDO
1893	ENDDO
1894	ENDDO
1895
1896	END SUBROUTINE mark_domain
1897
1898	#if ( NMM_CORE == 1 )
1899	!=======================================================================================
1900	! E grid interpolation for mass with addition of terrain adjustments. First routine
1901	! pertains to initial conditions and the next one corresponds to boundary conditions
1902	! This is gopal's doing
1903	!=======================================================================================
1904
1905	SUBROUTINE interp_mass_nmm (cfld, & ! CD field
1906	cids, cide, ckds, ckde, cjds, cjde, &
1907	cims, cime, ckms, ckme, cjms, cjme, &
1908	cits, cite, ckts, ckte, cjts, cjte, &
1909	nfld, & ! ND field
1910	nids, nide, nkds, nkde, njds, njde, &
1911	nims, nime, nkms, nkme, njms, njme, &
1912	nits, nite, nkts, nkte, njts, njte, &
1913	shw, & ! stencil half width for interp
1914	imask, & ! interpolation mask
1915	xstag, ystag, & ! staggering of field
1916	ipos, jpos, & ! Position of lower left of nest in CD
1917	nri, nrj, & ! nest ratios
1918	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
1919	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
1920	CBWGT4, HBWGT4, & ! dummys for weights
1921	CZ3d, Z3d, & ! Z3d interpolated from CZ3d
1922	CFIS,FIS, & ! CFIS dummy on fine domain
1923	CSM,SM, & ! CSM is dummy
1924	CPDTOP,PDTOP, &
1925	CPTOP,PTOP, &
1926	CPSTD,PSTD, &
1927	CKZMAX,KZMAX )
1928
1929	USE MODULE_MODEL_CONSTANTS
1930	USE module_timing
1931	IMPLICIT NONE
1932
1933	LOGICAL,INTENT(IN) :: xstag, ystag
1934	INTEGER,INTENT(IN) :: ckzmax,kzmax
1935	INTEGER,INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
1936	cims, cime, ckms, ckme, cjms, cjme, &
1937	cits, cite, ckts, ckte, cjts, cjte, &
1938	nids, nide, nkds, nkde, njds, njde, &
1939	nims, nime, nkms, nkme, njms, njme, &
1940	nits, nite, nkts, nkte, njts, njte, &
1941	shw,ipos,jpos,nri,nrj
1942
1943	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IMASK
1944
1945	! parent domain
1946
1947	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
1948	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3
1949	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT4,CFIS,CSM
1950	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CFLD
1951	REAL,DIMENSION(cims:cime,cjms:cjme,1:KZMAX), INTENT(IN) :: CZ3d
1952	REAL,DIMENSION(1:KZMAX), INTENT(IN) :: CPSTD
1953	REAL,INTENT(IN) :: CPDTOP,CPTOP
1954
1955	! nested domain
1956
1957	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIH,JJH
1958	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3
1959	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT4
1960	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: FIS,SM
1961	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(INOUT) :: NFLD
1962	REAL,DIMENSION(1:KZMAX), INTENT(IN) :: PSTD
1963	REAL,DIMENSION(nims:nime,njms:njme,1:KZMAX), INTENT(OUT) :: Z3d
1964	REAL,INTENT(IN) :: PDTOP,PTOP
1965
1966	! local
1967
1968	INTEGER,PARAMETER :: JTB=134
1969	REAL, PARAMETER :: LAPSR=6.5E-3,GI=1./G, D608=0.608
1970	REAL, PARAMETER :: COEF3=R_DGILAPSR
1971	INTEGER :: I,J,K,IDUM
1972	REAL :: dlnpdz,tvout,pmo
1973	REAL,DIMENSION(nims:nime,njms:njme) :: ZS,DUM2d
1974	REAL,DIMENSION(JTB) :: PIN,ZIN,Y2,PIO,ZOUT,DUM1,DUM2
1975	!-----------------------------------------------------------------------------------------------------
1976	!
1977	!*** CHECK DOMAIN BOUNDS BEFORE INTERPOLATION
1978	!
1979	DO J=NJTS,MIN(NJTE,NJDE-1)
1980	DO I=NITS,MIN(NITE,NIDE-1)
1981	IF(IIH(i,j).LT.(CIDS-shw) .OR. IIH(i,j).GT.(CIDE+shw)) &
1982	CALL wrf_error_fatal ('mass points:check domain bounds along x' )
1983	IF(JJH(i,j).LT.(CJDS-shw) .OR. JJH(i,j).GT.(CJDE+shw)) &
1984	CALL wrf_error_fatal ('mass points:check domain bounds along y' )
1985	ENDDO
1986	ENDDO
1987
1988	IF(KZMAX .GT. (JTB-10)) &
1989	CALL wrf_error_fatal ('mass points: increase JTB in interp_mass_nmm')
1990
1991	! WRITE(21,*)'------------- MED NEST INITIAL 1 ----------------'
1992	! DO J=NJTS,MIN(NJTE,NJDE-1)
1993	! DO I=NITS,MIN(NITE,NIDE-1)
1994	! WRITE(21,*)I,J,IMASK(I,J),NFLD(I,1,J)
1995	! ENDDO
1996	! ENDDO
1997	! WRITE(21,*)
1998
1999	!
2000	!*** DEFINE LOCAL TOPOGRAPHY ON THE BASIS OF FIS. ALSO CHECK IF SM IS LAND (SM=0) OVER TOPO
2001	!*** YOU DON'T WANT MOUNTAINS INSIDE WATER BODIES!
2002	!
2003
2004	DO J=NJTS,MIN(NJTE,NJDE-1)
2005	DO I=NITS,MIN(NITE,NIDE-1)
2006	ZS(I,J)=FIS(I,J)/G
2007	ENDDO
2008	ENDDO
2009
2010	!
2011	!*** Interpolate GPMs DERIVED FROM STANDARD ATMOSPHERIC LAPSE RATE FROM THE PARENT TO
2012	!*** THE NESTED DOMAIN
2013	!
2014	!*** INDEX CONVENTIONS
2015	!*** HBWGT4
2016	!*** 4
2017	!***
2018	!***
2019	!***
2020	!*** h
2021	!*** 1 2
2022	!*** HBWGT1 HBWGT2
2023	!***
2024	!***
2025	!*** 3
2026	!*** HBWGT3
2027
2028	Z3d=0.0
2029	DO K=NKTS,KZMAX ! Please note that we are still in isobaric surfaces
2030	DO J=NJTS,MIN(NJTE,NJDE-1)
2031	DO I=NITS,MIN(NITE,NIDE-1)
2032	!
2033	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
2034	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2035	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2036	+ HBWGT3(I,J)*CZ3d(IIH(I,J), JJH(I,J)-1,K) &
2037	+ HBWGT4(I,J)*CZ3d(IIH(I,J), JJH(I,J)+1,K)
2038	ELSE
2039	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2040	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2041	+ HBWGT3(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2042	+ HBWGT4(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)+1,K)
2043
2044	ENDIF
2045	!
2046	ENDDO
2047	ENDDO
2048	ENDDO
2049
2050	! RECONSTRUCT PDs ON THE BASIS OF TOPOGRAPHY AND THE INTERPOLATED HEIGHTS
2051
2052	DO J=NJTS,MIN(NJTE,NJDE-1)
2053	DO I=NITS,MIN(NITE,NIDE-1)
2054	!
2055	IF (ZS(I,J) .LT. Z3d(I,J,1)) THEN
2056	dlnpdz = (log(PSTD(1))-log(PSTD(2)) )/(Z3d(i,j,1)-Z3d(i,j,2))
2057	dum2d(i,j) = exp(log(PSTD(1)) + dlnpdz*(ZS(I,J) - Z3d(i,j,1)))
2058	dum2d(i,j) = dum2d(i,j) - PDTOP -PTOP
2059	ELSE ! target level bounded by input levels
2060	DO K =NKTS,KZMAX-1 ! still in the isobaric surfaces
2061	IF(ZS(I,J) .GE. Z3d(I,J,K) .AND. ZS(I,J) .LT. Z3d(I,J,K+1))THEN
2062	dlnpdz = (log(PSTD(K))-log(PSTD(K+1)) ) /(Z3d(I,J,K)-Z3d(I,J,K+1))
2063	dum2d(i,j) = exp(log(PSTD(K)) + dlnpdz*(ZS(I,J)- Z3d(I,J,K)))
2064	dum2d(i,j) = dum2d(i,j) - PDTOP -PTOP
2065	ENDIF
2066	ENDDO
2067	ENDIF
2068	IF(ZS(I,J) .GE. Z3d(I,J,KZMAX))THEN
2069	WRITE(0,*)'I=',I,'J=',J,'K=',KZMAX,'TERRAIN HEIGHT',ZS(I,J),'Z3d',Z3d(I,J,KZMAX)
2070	CALL wrf_error_fatal3 ( "interp_fcn.b" , 176 , "MOUNTAIN TOO HIGH TO FIT THE MODEL DEPTH")
2071	ENDIF
2072	!
2073	ENDDO
2074	ENDDO
2075
2076	DO K=NKDS,NKDE ! NKTE is 1, nevertheless let us pretend religious
2077	DO J=NJTS,MIN(NJTE,NJDE-1)
2078	DO I=NITS,MIN(NITE,NIDE-1)
2079	IF(IMASK(I,J) .NE. 1)THEN
2080	NFLD(I,J,K)= dum2d(i,j) ! PD defined in the nested domain
2081	ENDIF
2082	ENDDO
2083	ENDDO
2084	ENDDO
2085
2086	!
2087	END SUBROUTINE interp_mass_nmm
2088	!
2089	!--------------------------------------------------------------------------------------
2090
2091	SUBROUTINE nmm_bdymass_hinterp ( cfld, & ! CD field
2092	cids, cide, ckds, ckde, cjds, cjde, &
2093	cims, cime, ckms, ckme, cjms, cjme, &
2094	cits, cite, ckts, ckte, cjts, cjte, &
2095	nfld, & ! ND field
2096	nids, nide, nkds, nkde, njds, njde, &
2097	nims, nime, nkms, nkme, njms, njme, &
2098	nits, nite, nkts, nkte, njts, njte, &
2099	shw, & ! stencil half width
2100	imask, & ! interpolation mask
2101	xstag, ystag, & ! staggering of field
2102	ipos, jpos, & ! Position of lower left of nest in CD
2103	nri, nrj, & ! nest ratios
2104	c_bxs,n_bxs, &
2105	c_bxe,n_bxe, &
2106	c_bys,n_bys, &
2107	c_bye,n_bye, &
2108	c_btxs,n_btxs, &
2109	c_btxe,n_btxe, &
2110	c_btys,n_btys, &
2111	c_btye,n_btye, &
2112	CTEMP_B,NTEMP_B, & ! These temp arrays should be removed
2113	CTEMP_BT,NTEMP_BT, & ! later on
2114	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
2115	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
2116	CBWGT4, HBWGT4, & ! dummys
2117	CZ3d, Z3d, & ! Z3d dummy on nested domain
2118	CFIS,FIS, & ! CFIS dummy on fine domain
2119	CSM,SM, & ! CSM is dummy
2120	CPDTOP,PDTOP, &
2121	CPTOP,PTOP, &
2122	CPSTD,PSTD, &
2123	CKZMAX,KZMAX )
2124
2125
2126	USE MODULE_MODEL_CONSTANTS
2127	USE module_configure
2128	USE module_wrf_error
2129
2130	IMPLICIT NONE
2131
2132
2133	INTEGER, INTENT(IN) :: ckzmax,kzmax
2134	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
2135	cims, cime, ckms, ckme, cjms, cjme, &
2136	cits, cite, ckts, ckte, cjts, cjte, &
2137	nids, nide, nkds, nkde, njds, njde, &
2138	nims, nime, nkms, nkme, njms, njme, &
2139	nits, nite, nkts, nkte, njts, njte, &
2140	shw, &
2141	ipos, jpos, &
2142	nri, nrj
2143
2144
2145	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ), INTENT(OUT) :: ctemp_b,ctemp_bt
2146	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ), INTENT(OUT) :: ntemp_b,ntemp_bt
2147	LOGICAL, INTENT(IN) :: xstag, ystag
2148	REAL, DIMENSION( * ), INTENT(INOUT) :: c_bxs,n_bxs,c_bxe,n_bxe,c_bys,n_bys,c_bye,n_bye
2149	REAL, DIMENSION( * ), INTENT(INOUT) :: c_btxs,n_btxs,c_btxe,n_btxe,c_btys,n_btys,c_btye,n_btye
2150
2151	! parent domain
2152
2153	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IMASK
2154	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
2155	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3
2156	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT4,CFIS,CSM
2157	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CFLD
2158	REAL,DIMENSION(cims:cime,cjms:cjme,1:KZMAX), INTENT(IN) :: CZ3d
2159	REAL,DIMENSION(1:KZMAX), INTENT(IN) :: CPSTD
2160	REAL,INTENT(IN) :: CPDTOP,CPTOP
2161
2162	! nested domain
2163
2164	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIH,JJH
2165	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3
2166	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT4
2167	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: FIS,SM
2168	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(INOUT) :: NFLD
2169	REAL,DIMENSION(1:KZMAX), INTENT(IN) :: PSTD
2170	REAL,DIMENSION(nims:nime,njms:njme,1:KZMAX), INTENT(OUT) :: Z3d
2171	REAL,INTENT(IN) :: PDTOP,PTOP
2172
2173	! Local
2174
2175	INTEGER :: nijds, nijde, spec_bdy_width,i,j,k
2176	REAL :: dlnpdz,dum2d
2177	REAL,DIMENSION(nims:nime,njms:njme) :: zs
2178
2179	INTEGER,PARAMETER :: JTB=134
2180	INTEGER :: ii,jj
2181	REAL, DIMENSION (nims:nime,njms:njme) :: CWK1,CWK2,CWK3,CWK4
2182
2183	nijds = min(nids, njds)
2184	nijde = max(nide, njde)
2185	CALL nl_get_spec_bdy_width( 1, spec_bdy_width )
2186
2187	!
2188	!*** DEFINE LOCAL TOPOGRAPHY ON THE BASIS OF FIS. ASLO CHECK IF SM IS LAND (SM=0) OVER TOPO
2189	!*** YOU DON'T WANT MOUNTAINS INSIDE WATER BODIES!
2190	!
2191	DO J=NJTS,MIN(NJTE,NJDE-1)
2192	DO I=NITS,MIN(NITE,NIDE-1)
2193	ZS(I,J)=FIS(I,J)/G
2194	ENDDO
2195	ENDDO
2196
2197	! X start boundary
2198
2199	NMM_XS: IF(NITS .EQ. NIDS)THEN
2200	! WRITE(0,*)'ENTERING X1 START BOUNDARY AT MASS POINTS',NJTS,MIN(NJTE,NJDE-1)
2201	I = NIDS
2202
2203	DO K=NKTS,KZMAX
2204	DO J = NJTS,MIN(NJTE,NJDE-1)
2205	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
2206	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2207	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2208	+ HBWGT3(I,J)*CZ3d(IIH(I,J), JJH(I,J)-1,K) &
2209	+ HBWGT4(I,J)*CZ3d(IIH(I,J), JJH(I,J)+1,K)
2210	ELSE
2211	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2212	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2213	+ HBWGT3(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2214	+ HBWGT4(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)+1,K)
2215	ENDIF
2216	END DO
2217	END DO
2218
2219	DO J = NJTS,MIN(NJTE,NJDE-1)
2220	IF(MOD(J,2) .NE. 0)THEN
2221	IF (ZS(I,J) .LT. Z3d(I,J,2)) THEN ! level 2 has to be changed
2222	dlnpdz = (log(PSTD(1))-log(PSTD(2)) )/(Z3d(I,J,1)-Z3d(I,J,2))
2223	dum2d = exp(log(PSTD(1)) + dlnpdz*(ZS(I,J) - Z3d(I,J,1)))
2224	CWK1(I,J) = dum2d -PDTOP -PTOP
2225	ELSE ! target level bounded by input levels
2226	DO K =NKTS,KZMAX-1
2227	IF(ZS(I,J) .GE. Z3d(I,J,K) .AND. ZS(I,J) .LT. Z3d(I,J,K+1))THEN
2228	dlnpdz = (log(PSTD(K))-log(PSTD(K+1)) ) /(Z3d(I,J,K)-Z3d(I,J,K+1))
2229	dum2d = exp(log(PSTD(K)) + dlnpdz*(ZS(I,J)- Z3d(I,J,K)))
2230	CWK1(I,J) = dum2d -PDTOP -PTOP
2231	ENDIF
2232	ENDDO
2233	ENDIF
2234	IF(ZS(I,J) .GE. Z3d(I,J,KZMAX))THEN
2235	WRITE(0,*)'I=',I,'J=',J,'K=',K,'TERRAIN HEIGHT',ZS(I,J),'Z3d',Z3d(I,J,KZMAX)
2236	CALL wrf_error_fatal("BC:MOUNTAIN TOO HIGH TO FIT THE MODEL DEPTH")
2237	ENDIF
2238	ELSE
2239	CWK1(I,J)=0.
2240	ENDIF
2241	ENDDO
2242
2243	DO J = NJTS,MIN(NJTE,NJDE-1)
2244	DO K = NKDS,NKDE
2245	ntemp_b(i,j,k) = CWK1(I,J)
2246	ntemp_bt(i,j,k) = 0.0
2247	END DO
2248	END DO
2249	ENDIF NMM_XS
2250
2251	! X end boundary
2252
2253	NMM_XE: IF(NITE-1 .EQ. NIDE-1)THEN
2254	! WRITE(0,*)'ENTERING X END BOUNDARY AT MASS POINTS',NJTS,MIN(NJTE,NJDE-1)
2255	I = NIDE-1
2256	II = NIDE - I
2257
2258	DO K=NKTS,KZMAX
2259	DO J=NJTS,MIN(NJTE,NJDE-1)
2260	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
2261	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2262	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2263	+ HBWGT3(I,J)*CZ3d(IIH(I,J), JJH(I,J)-1,K) &
2264	+ HBWGT4(I,J)*CZ3d(IIH(I,J), JJH(I,J)+1,K)
2265	ELSE
2266	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2267	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2268	+ HBWGT3(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2269	+ HBWGT4(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)+1,K)
2270	ENDIF
2271	ENDDO
2272	ENDDO
2273
2274	DO J = NJTS,MIN(NJTE,NJDE-1)
2275	IF(MOD(J,2) .NE.0)THEN ! 1,3,5,7 of nested domain
2276	IF (ZS(I,J) .LT. Z3d(I,J,2)) THEN ! level 2 has to be changed
2277	dlnpdz = (log(PSTD(1))-log(PSTD(2)) )/(Z3d(I,J,1)-Z3d(I,J,2))
2278	dum2d = exp(log(PSTD(1)) + dlnpdz*(ZS(I,J) - Z3d(I,J,1)))
2279	CWK2(I,J) = dum2d -PDTOP -PTOP
2280	ELSE ! target level bounded by input levels
2281	DO K =NKTS,KZMAX-1
2282	IF(ZS(I,J) .GE. Z3d(I,J,K) .AND. ZS(I,J) .LT. Z3d(I,J,K+1))THEN
2283	dlnpdz = (log(PSTD(K))-log(PSTD(K+1)) ) /(Z3d(I,J,K)-Z3d(I,J,K+1))
2284	dum2d = exp(log(PSTD(K)) + dlnpdz*(ZS(I,J)- Z3d(I,J,K)))
2285	CWK2(I,J) = dum2d -PDTOP -PTOP
2286	ENDIF
2287	ENDDO
2288	ENDIF
2289	IF(ZS(I,J) .GE. Z3d(I,J,KZMAX))THEN
2290	WRITE(0,*)'I=',I,'J=',J,'K=',K,'TERRAIN HEIGHT',ZS(I,J),'Z3d',Z3d(I,J,KZMAX)
2291	CALL wrf_error_fatal("BC:MOUNTAIN TOO HIGH TO FIT THE MODEL DEPTH")
2292	ENDIF
2293	ELSE
2294	CWK2(I,J) = 0.0
2295	ENDIF
2296	ENDDO
2297
2298	DO J = NJTS,MIN(NJTE,NJDE-1)
2299	DO K = NKDS,NKDE
2300	ntemp_b(i,j,k) = CWK2(I,J)
2301	ntemp_bt(i,j,k) = 0.0
2302	END DO
2303	END DO
2304	ENDIF NMM_XE
2305
2306	! Y start boundary
2307
2308	NMM_YS: IF(NJTS .EQ. NJDS)THEN
2309	! WRITE(20,*)'ENTERING Y START BOUNDARY AT MASS POINTS',NITS,MIN(NITE,NIDE-1)
2310	J = NJDS
2311	DO K=NKTS,KZMAX
2312	DO I = NITS,MIN(NITE,NIDE-1)
2313	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
2314	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2315	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2316	+ HBWGT3(I,J)*CZ3d(IIH(I,J), JJH(I,J)-1,K) &
2317	+ HBWGT4(I,J)*CZ3d(IIH(I,J), JJH(I,J)+1,K)
2318	ELSE
2319	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2320	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2321	+ HBWGT3(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2322	+ HBWGT4(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)+1,K)
2323	ENDIF
2324	END DO
2325	END DO
2326
2327	DO I = NITS,MIN(NITE,NIDE-1)
2328	IF (ZS(I,J) .LT. Z3d(I,J,2)) THEN ! level 2 has to be changed
2329	dlnpdz = (log(PSTD(1))-log(PSTD(2)) )/(Z3d(I,J,1)-Z3d(I,J,2))
2330	dum2d = exp(log(PSTD(1)) + dlnpdz*(ZS(I,J) - Z3d(I,J,1)))
2331	CWK3(I,J) = dum2d -PDTOP -PTOP
2332	ELSE ! target level bounded by input levels
2333	DO K =NKTS,KZMAX-1
2334	IF(ZS(I,J) .GE. Z3d(I,J,K) .AND. ZS(I,J) .LT. Z3d(I,J,K+1))THEN
2335	dlnpdz = (log(PSTD(K))-log(PSTD(K+1)) ) /(Z3d(I,J,K)-Z3d(I,J,K+1))
2336	dum2d = exp(log(PSTD(K)) + dlnpdz*(ZS(I,J)- Z3d(I,J,K)))
2337	CWK3(I,J) = dum2d -PDTOP -PTOP
2338	ENDIF
2339	ENDDO
2340	ENDIF
2341	IF(ZS(I,J) .GE. Z3d(I,J,KZMAX))THEN
2342	WRITE(0,*)'I=',I,'J=',J,'K=',K,'TERRAIN HEIGHT',ZS(I,J),'Z3d',Z3d(I,J,KZMAX)
2343	CALL wrf_error_fatal("BC:MOUNTAIN TOO HIGH TO FIT THE MODEL DEPTH")
2344	ENDIF
2345	ENDDO
2346
2347	DO K = NKDS, NKDE
2348	DO I = NITS,MIN(NITE,NIDE-1)
2349	ntemp_b(i,j,k) = CWK3(I,J)
2350	ntemp_bt(i,j,k) = 0.0
2351	END DO
2352	END DO
2353	END IF NMM_YS
2354
2355	! Y end boundary
2356
2357	NMM_YE: IF(NJTE-1 .EQ. NJDE-1)THEN
2358	! WRITE(20,*)'ENTERING Y END BOUNDARY AT MASS POINTS',NITS,MIN(NITE,NIDE-1)
2359	J = NJDE-1
2360	JJ = NJDE - J
2361	DO K=NKTS,KZMAX
2362	DO I = NITS,MIN(NITE,NIDE-1)
2363	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
2364	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2365	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2366	+ HBWGT3(I,J)*CZ3d(IIH(I,J), JJH(I,J)-1,K) &
2367	+ HBWGT4(I,J)*CZ3d(IIH(I,J), JJH(I,J)+1,K)
2368	ELSE
2369	Z3d(I,J,K) = HBWGT1(I,J)*CZ3d(IIH(I,J), JJH(I,J) ,K) &
2370	+ HBWGT2(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J) ,K) &
2371	+ HBWGT3(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2372	+ HBWGT4(I,J)*CZ3d(IIH(I,J)+1,JJH(I,J)+1,K)
2373	ENDIF
2374	END DO
2375	END DO
2376
2377	DO I = NITS,MIN(NITE,NIDE-1)
2378	IF (ZS(I,J) .LT. Z3d(I,J,2)) THEN ! level 2 has to be changed
2379	dlnpdz = (log(PSTD(1))-log(PSTD(2)) )/(Z3d(I,J,1)-Z3d(I,J,2))
2380	dum2d = exp(log(PSTD(1)) + dlnpdz*(ZS(I,J) - Z3d(I,J,1)))
2381	CWK4(I,J) = dum2d -PDTOP -PTOP
2382	ELSE ! target level bounded by input levels
2383	DO K =NKTS,KZMAX-1
2384	IF(ZS(I,J) .GE. Z3d(I,J,K) .AND. ZS(I,J) .LT. Z3d(I,J,K+1))THEN
2385	dlnpdz = (log(PSTD(K))-log(PSTD(K+1)) ) /(Z3d(I,J,K)-Z3d(I,J,K+1))
2386	dum2d = exp(log(PSTD(K)) + dlnpdz*(ZS(I,J)- Z3d(I,J,K)))
2387	CWK4(I,J) = dum2d -PDTOP -PTOP
2388	ENDIF
2389	ENDDO
2390	ENDIF
2391	IF(ZS(I,J) .GE. Z3d(I,J,KZMAX))THEN
2392	WRITE(0,*)'I=',I,'J=',J,'K=',K,'TERRAIN HEIGHT',ZS(I,J),'Z3d',Z3d(I,J,KZMAX)
2393	CALL wrf_error_fatal("BC:MOUNTAIN TOO HIGH TO FIT THE MODEL DEPTH")
2394	ENDIF
2395	ENDDO
2396
2397	DO K = NKDS,NKDE
2398	DO I = NITS,MIN(NITE,NIDE-1)
2399	ntemp_b(i,j,k) = CWK4(I,J)
2400	ntemp_bt(i,j,k) = 0.0
2401	END DO
2402	END DO
2403	END IF NMM_YE
2404
2405	RETURN
2406
2407	END SUBROUTINE nmm_bdymass_hinterp
2408	!
2409	!=======================================================================================
2410	!
2411	! ADDED FOR INCLUDING MOISTURE AND THERMODYNAMIC ENERGY BALANCE
2412	!
2413	!=======================================================================================
2414
2415	SUBROUTINE interp_scalar_nmm (cfld, & ! CD field
2416	cids,cide,ckds,ckde,cjds,cjde, &
2417	cims,cime,ckms,ckme,cjms,cjme, &
2418	cits,cite,ckts,ckte,cjts,cjte, &
2419	nfld, & ! ND field
2420	nids,nide,nkds,nkde,njds,njde, &
2421	nims,nime,nkms,nkme,njms,njme, &
2422	nits,nite,nkts,nkte,njts,njte, &
2423	shw, & ! stencil half width for interp
2424	imask, & ! interpolation mask
2425	xstag,ystag, & ! staggering of field
2426	ipos,jpos, & ! Position of lower left of nest in CD
2427	nri,nrj, & ! nest ratios
2428	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
2429	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
2430	CBWGT4, HBWGT4, & ! dummys for weights
2431	CC3d,C3d, &
2432	CPD,PD, &
2433	CPSTD,PSTD, &
2434	CPDTOP,PDTOP, &
2435	CPTOP,PTOP, &
2436	CETA1,ETA1,CETA2,ETA2 )
2437
2438	USE MODULE_MODEL_CONSTANTS
2439	USE module_timing
2440	IMPLICIT NONE
2441
2442	LOGICAL,INTENT(IN) :: xstag, ystag
2443	INTEGER,INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
2444	cims, cime, ckms, ckme, cjms, cjme, &
2445	cits, cite, ckts, ckte, cjts, cjte, &
2446	nids, nide, nkds, nkde, njds, njde, &
2447	nims, nime, nkms, nkme, njms, njme, &
2448	nits, nite, nkts, nkte, njts, njte, &
2449	shw,ipos,jpos,nri,nrj
2450
2451	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IMASK
2452
2453	! parent domain
2454
2455	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
2456	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2
2457	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT3,CBWGT4
2458
2459	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CFLD
2460	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CC3d ! scalar input on constant pressure levels
2461	REAL,DIMENSION(ckms:ckme), INTENT(IN) :: CPSTD
2462	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CPD
2463	REAL,DIMENSION(ckms:ckme), INTENT(IN) :: CETA1,CETA2
2464	REAL, INTENT(IN) :: CPDTOP,CPTOP
2465
2466	! nested domain
2467
2468	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIH,JJH
2469	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT1,HBWGT2
2470	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT3,HBWGT4
2471
2472	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(OUT):: NFLD ! This is scalar on hybrid levels
2473	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(OUT):: C3d ! Scalar on constant pressure levels
2474	REAL,DIMENSION(nkms:nkme), INTENT(IN) :: PSTD
2475	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: PD
2476	REAL,DIMENSION(nkms:nkme), INTENT(IN) :: ETA1,ETA2
2477	REAL,INTENT(IN) :: PDTOP,PTOP
2478
2479	! local
2480
2481	INTEGER,PARAMETER :: JTB=134
2482	INTEGER :: I,J,K
2483	REAL,DIMENSION(JTB) :: PIN,CIN,Y2,PIO,PTMP,COUT,DUM1,DUM2
2484
2485	!-----------------------------------------------------------------------------------------------------
2486	!
2487	!
2488	! *** CHECK VERTICAL BOUNDS BEFORE USING SPLINE OR LINEAR INTERPOLATION
2489	!
2490	IF(nkme .GT. (JTB-10) .OR. NKDE .GT. (JTB-10)) &
2491	CALL wrf_error_fatal ('mass points: increase JTB in interp_mass_nmm')
2492
2493	!
2494	! FIRST, HORIZONTALLY INTERPOLATE MOISTURE NOW AVAILABLE ON CONSTANT PRESSURE SURFACE (LEVELS) FROM THE
2495	! PARENT TO THE NESTED DOMAIN
2496	!
2497	!*** INDEX CONVENTIONS
2498	!*** HBWGT4
2499	!*** 4
2500	!***
2501	!***
2502	!***
2503	!*** h
2504	!*** 1 2
2505	!*** HBWGT1 HBWGT2
2506	!***
2507	!***
2508	!*** 3
2509	!*** HBWGT3
2510
2511	C3d=0.0
2512	DO K=NKDS,NKDE-1 ! Please note that we are still in isobaric surfaces
2513	DO J=NJTS,MIN(NJTE,NJDE-1)
2514	DO I=NITS,MIN(NITE,NIDE-1)
2515	IF(IMASK(I,J) .NE. 1)THEN
2516	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
2517	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2518	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2519	+ HBWGT3(I,J)*CC3d(IIH(I,J), JJH(I,J)-1,K) &
2520	+ HBWGT4(I,J)*CC3d(IIH(I,J), JJH(I,J)+1,K)
2521
2522	ELSE
2523	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2524	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2525	+ HBWGT3(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2526	+ HBWGT4(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)+1,K)
2527
2528	ENDIF
2529	ENDIF
2530	ENDDO
2531	ENDDO
2532	ENDDO
2533
2534	!
2535	! RECOVER THE SCALARS FROM CONSTANT PRESSURE SURFACES (LEVELS) ON TO HYBRID SURFACES
2536	!
2537	DO J=NJTS,MIN(NJTE,NJDE-1)
2538	DO I=NITS,MIN(NITE,NIDE-1)
2539	IF(IMASK(I,J) .NE. 1)THEN
2540
2541	! clean local array before use of spline or linear interpolation
2542
2543	CIN=0.;PIN=0.;Y2=0;PIO=0.;PTMP=0.;COUT=0.;DUM1=0.;DUM2=0.
2544
2545	DO K=NKDS+1,NKDE ! inputs at standard levels
2546	PIN(K-1) = EXP((ALOG(PSTD(NKDE-K+1))+ALOG(PSTD(NKDE-K+2)))*0.5)
2547	CIN(K-1) = C3d(I,J,NKDE-K+1)
2548	ENDDO
2549
2550	Y2(1 )=0.
2551	Y2(NKDE-1)=0.
2552
2553	DO K=NKDS,NKDE ! target points in model interface levels (pint)
2554	PTMP(K) = ETA1(K)PDTOP + ETA2(K)PD(I,J) + PTOP
2555	ENDDO
2556
2557	DO K=NKDS,NKDE-1 ! target points in model levels
2558	PIO(K) = EXP((ALOG(PTMP(K))+ALOG(PTMP(K+1)))*0.5)
2559	ENDDO
2560
2561	IF(PTMP(1) .GE. PSTD(1))THEN ! if lower boundary is higher than PMSL(1) re-set lower boundary
2562	PIN(NKDE-1) = PIO(1) ! be consistent with target. This may not happen at all
2563	WRITE(0,*)'WARNING: NESTED DOMAIN PRESSURE AT LOWEST LEVEL HIGHER THAN PSTD'
2564	WRITE(0,*)'I,J,PIO(1),PSTD(1)',I,J,PIO(1),PSTD(1)
2565	ENDIF
2566
2567	CALL SPLINE2(I,J,JTB,NKDE-1,PIN,CIN,Y2,NKDE-1,PIO,COUT,DUM1,DUM2) ! interpolate
2568
2569	DO K=1,NKDE-1
2570	NFLD(I,J,K)= COUT(K) ! scalar in the nested domain
2571	ENDDO
2572
2573	ENDIF
2574	ENDDO
2575	ENDDO
2576
2577	END SUBROUTINE interp_scalar_nmm
2578	!
2579	!===========================================================================================
2580	!
2581	SUBROUTINE nmm_bdy_scalar (cfld, & ! CD field
2582	cids,cide,ckds,ckde,cjds,cjde, &
2583	cims,cime,ckms,ckme,cjms,cjme, &
2584	cits,cite,ckts,ckte,cjts,cjte, &
2585	nfld, & ! ND field
2586	nids,nide,nkds,nkde,njds,njde, &
2587	nims,nime,nkms,nkme,njms,njme, &
2588	nits,nite,nkts,nkte,njts,njte, &
2589	shw, & ! stencil half width for interp
2590	imask, & ! interpolation mask
2591	xstag,ystag, & ! staggering of field
2592	ipos,jpos, & ! Position of lower left of nest in CD
2593	nri,nrj, & ! nest ratios
2594	c_bxs,n_bxs, &
2595	c_bxe,n_bxe, &
2596	c_bys,n_bys, &
2597	c_bye,n_bye, &
2598	c_btxs,n_btxs, &
2599	c_btxe,n_btxe, &
2600	c_btys,n_btys, &
2601	c_btye,n_btye, &
2602	cdt, ndt, &
2603	CTEMP_B,NTEMP_B, & ! to be removed
2604	CTEMP_BT,NTEMP_BT, &
2605	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
2606	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
2607	CBWGT4, HBWGT4, & ! dummys for weights
2608	CC3d,C3d, &
2609	CPD,PD, &
2610	CPSTD,PSTD, &
2611	CPDTOP,PDTOP, &
2612	CPTOP,PTOP, &
2613	CETA1,ETA1,CETA2,ETA2 )
2614	USE MODULE_MODEL_CONSTANTS
2615	USE module_timing
2616	IMPLICIT NONE
2617
2618	LOGICAL,INTENT(IN) :: xstag, ystag
2619	REAL, INTENT(INOUT) :: cdt, ndt
2620	INTEGER,INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
2621	cims, cime, ckms, ckme, cjms, cjme, &
2622	cits, cite, ckts, ckte, cjts, cjte, &
2623	nids, nide, nkds, nkde, njds, njde, &
2624	nims, nime, nkms, nkme, njms, njme, &
2625	nits, nite, nkts, nkte, njts, njte, &
2626	shw,ipos,jpos,nri,nrj
2627	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ), INTENT(OUT) :: ctemp_b,ctemp_bt
2628	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ), INTENT(OUT) :: ntemp_b,ntemp_bt
2629	REAL, DIMENSION( * ), INTENT(INOUT) :: c_bxs,n_bxs,c_bxe,n_bxe,c_bys,n_bys,c_bye,n_bye
2630	REAL, DIMENSION( * ), INTENT(INOUT) :: c_btxs,n_btxs,c_btxe,n_btxe,c_btys,n_btys,c_btye,n_btye
2631
2632
2633	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IMASK
2634
2635	! parent domain
2636
2637	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
2638	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2
2639	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT3,CBWGT4
2640	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CFLD
2641	REAL,DIMENSION(cims:cime,cjms:cjme,ckms:ckme), INTENT(IN) :: CC3d ! scalar input on constant pressure levels
2642	REAL,DIMENSION(ckms:ckme), INTENT(IN) :: CPSTD
2643	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CPD
2644	REAL,DIMENSION(ckms:ckme), INTENT(IN) :: CETA1,CETA2
2645	REAL, INTENT(IN) :: CPDTOP,CPTOP
2646
2647	! nested domain
2648
2649	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIH,JJH
2650	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT1,HBWGT2
2651	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT3,HBWGT4
2652	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(OUT):: NFLD
2653	REAL,DIMENSION(nims:nime,njms:njme,nkms:nkme), INTENT(OUT):: C3d !Scalar on constant pressure levels
2654	REAL,DIMENSION(nkms:nkme), INTENT(IN) :: PSTD
2655	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: PD
2656	REAL,DIMENSION(nkms:nkme), INTENT(IN) :: ETA1,ETA2
2657	REAL,INTENT(IN) :: PDTOP,PTOP
2658
2659	! local
2660
2661	INTEGER,PARAMETER :: JTB=134
2662	INTEGER :: I,J,K,II,JJ
2663	REAL,DIMENSION(JTB) :: PIN,CIN,Y2,PIO,PTMP,COUT,DUM1,DUM2
2664	REAL, DIMENSION (nims:nime,njms:njme,nkms:nkme) :: CWK1,CWK2,CWK3,CWK4
2665	!-----------------------------------------------------------------------------------------------------
2666	!
2667	!
2668	! *** CHECK VERTICAL BOUNDS BEFORE USING SPLINE INTERPOLATION
2669	!
2670	IF(nkme .GT. (JTB-10) .OR. NKDE .GT. (JTB-10)) &
2671	CALL wrf_error_fatal ('mass points: increase JTB in interp_mass_nmm')
2672
2673	! X start boundary
2674
2675	NMM_XS: IF(NITS .EQ. NIDS)THEN
2676	! WRITE(0,*)'ENTERING X1 START BOUNDARY AT T POINTS',NJTS,MIN(NJTE,NJDE-1)
2677	I = NIDS
2678	DO K=NKDS,NKDE-1 ! Please note that we are still in isobaric surfaces
2679	DO J = NJTS,MIN(NJTE,NJDE-1)
2680	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
2681	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2682	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2683	+ HBWGT3(I,J)*CC3d(IIH(I,J), JJH(I,J)-1,K) &
2684	+ HBWGT4(I,J)*CC3d(IIH(I,J), JJH(I,J)+1,K)
2685	ELSE
2686	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2687	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2688	+ HBWGT3(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2689	+ HBWGT4(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)+1,K)
2690	ENDIF
2691	ENDDO
2692	ENDDO
2693	!
2694	DO J=NJTS,MIN(NJTE,NJDE-1)
2695	IF(MOD(J,2) .NE. 0)THEN
2696	CIN=0.;PIN=0.;Y2=0;PIO=0.;PTMP=0.;COUT=0.;DUM1=0.;DUM2=0. ! clean up local array
2697	DO K=NKDS+1,NKDE ! inputs at standard levels
2698	PIN(K-1) = EXP((ALOG(PSTD(NKDE-K+1))+ALOG(PSTD(NKDE-K+2)))*0.5)
2699	CIN(K-1) = C3d(I,J,NKDE-K+1)
2700	ENDDO
2701	Y2(1 )=0.
2702	Y2(NKDE-1)=0.
2703	DO K=NKDS,NKDE ! target points in model interface levels (pint)
2704	PTMP(K) = ETA1(K)PDTOP + ETA2(K)PD(I,J) + PTOP
2705	ENDDO
2706	DO K=NKDS,NKDE-1 ! target points in model levels
2707	PIO(K) = EXP((ALOG(PTMP(K))+ALOG(PTMP(K+1)))*0.5)
2708	ENDDO
2709	IF(PTMP(1) .GE. PSTD(1))THEN ! if lower boundary is higher than PMSL(1) re-set lower boundary
2710	PIN(NKDE-1) = PIO(1) ! be consistent with target. This may not happen at all
2711	WRITE(0,*)'WARNING: NESTED DOMAIN PRESSURE AT LOWEST LEVEL HIGHER THAN PSTD'
2712	WRITE(0,*)'I,J,PIO(1),PSTD(1)',I,J,PIO(1),PSTD(1)
2713	ENDIF
2714
2715	CALL SPLINE2(I,J,JTB,NKDE-1,PIN,CIN,Y2,NKDE-1,PIO,COUT,DUM1,DUM2) ! interpolate
2716
2717	DO K=1,NKDE-1
2718	CWK1(I,J,K)= COUT(K) ! scalar in the nested domain
2719	ENDDO
2720	ELSE
2721	DO K=NKDS,NKDE-1
2722	CWK1(I,J,K)=0.0
2723	ENDDO
2724	ENDIF
2725	ENDDO
2726
2727	DO J = NJTS,MIN(NJTE,NJDE-1)
2728	DO K = NKDS,NKDE-1
2729	ntemp_b(i,j,k) = CWK1(I,J,K)
2730	ntemp_bt(i,j,k) = 0.0
2731	END DO
2732	END DO
2733
2734	ENDIF NMM_XS
2735
2736
2737	! X end boundary
2738
2739	NMM_XE: IF(NITE-1 .EQ. NIDE-1)THEN
2740	! WRITE(0,*)'ENTERING X END BOUNDARY AT T POINTS',NJTS,MIN(NJTE,NJDE-1)
2741	I = NIDE-1
2742	DO K=NKDS,NKDE-1 ! Please note that we are still in isobaric surfaces
2743	DO J = NJTS,MIN(NJTE,NJDE-1)
2744	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
2745	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2746	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2747	+ HBWGT3(I,J)*CC3d(IIH(I,J), JJH(I,J)-1,K) &
2748	+ HBWGT4(I,J)*CC3d(IIH(I,J), JJH(I,J)+1,K)
2749	ELSE
2750	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2751	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2752	+ HBWGT3(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2753	+ HBWGT4(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)+1,K)
2754	ENDIF
2755	ENDDO
2756	ENDDO
2757
2758	DO J=NJTS,MIN(NJTE,NJDE-1)
2759	IF(MOD(J,2) .NE. 0)THEN
2760	CIN=0.;PIN=0.;Y2=0;PIO=0.;PTMP=0.;COUT=0.;DUM1=0.;DUM2=0. ! clean up local array
2761	DO K=NKDS+1,NKDE ! inputs at standard levels
2762	PIN(K-1) = EXP((ALOG(PSTD(NKDE-K+1))+ALOG(PSTD(NKDE-K+2)))*0.5)
2763	CIN(K-1) = C3d(I,J,NKDE-K+1)
2764	ENDDO
2765	Y2(1 )=0.
2766	Y2(NKDE-1)=0.
2767	DO K=NKDS,NKDE ! target points in model interface levels (pint)
2768	PTMP(K) = ETA1(K)PDTOP + ETA2(K)PD(I,J) + PTOP
2769	ENDDO
2770	DO K=NKDS,NKDE-1 ! target points in model levels
2771	PIO(K) = EXP((ALOG(PTMP(K))+ALOG(PTMP(K+1)))*0.5)
2772	ENDDO
2773	IF(PTMP(1) .GE. PSTD(1))THEN ! if lower boundary is higher than PMSL(1) re-set lower boundary
2774	PIN(NKDE-1) = PIO(1) ! be consistent with target. This may not happen at all
2775	WRITE(0,*)'WARNING: NESTED DOMAIN PRESSURE AT LOWEST LEVEL HIGHER THAN PSTD'
2776	WRITE(0,*)'I,J,PIO(1),PSTD(1)',I,J,PIO(1),PSTD(1)
2777	ENDIF
2778
2779	CALL SPLINE2(I,J,JTB,NKDE-1,PIN,CIN,Y2,NKDE-1,PIO,COUT,DUM1,DUM2) ! interpolate
2780
2781	DO K=1,NKDE-1
2782	CWK2(I,J,K)= COUT(K) ! scalar in the nested domain
2783	ENDDO
2784	ELSE
2785	DO K=NKDS,NKDE-1
2786	CWK2(I,J,K)=0.0
2787	ENDDO
2788	ENDIF
2789	ENDDO
2790
2791	DO J = NJTS,MIN(NJTE,NJDE-1)
2792	DO K = NKDS,MIN(NKTE,NKDE-1)
2793	ntemp_b(i,j,k) = CWK2(I,J,K)
2794	ntemp_bt(i,j,k) = 0.0
2795	END DO
2796	END DO
2797
2798	ENDIF NMM_XE
2799
2800	! Y start boundary
2801
2802	NMM_YS: IF(NJTS .EQ. NJDS)THEN
2803	! WRITE(0,*)'ENTERING Y START BOUNDARY AT T POINTS',NITS,MIN(NITE,NIDE-1)
2804	J = NJDS
2805	DO K=NKDS,NKDE-1
2806	DO I = NITS,MIN(NITE,NIDE-1)
2807	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
2808	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2809	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2810	+ HBWGT3(I,J)*CC3d(IIH(I,J), JJH(I,J)-1,K) &
2811	+ HBWGT4(I,J)*CC3d(IIH(I,J), JJH(I,J)+1,K)
2812	ELSE
2813	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2814	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2815	+ HBWGT3(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2816	+ HBWGT4(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)+1,K)
2817
2818	ENDIF
2819	ENDDO
2820	ENDDO
2821	!
2822	DO I=NITS,MIN(NITE,NIDE-1)
2823	CIN=0.;PIN=0.;Y2=0;PIO=0.;PTMP=0.;COUT=0.;DUM1=0.;DUM2=0. ! clean up local array
2824	DO K=NKDS+1,NKDE ! inputs at standard levels
2825	PIN(K-1) = EXP((ALOG(PSTD(NKDE-K+1))+ALOG(PSTD(NKDE-K+2)))*0.5)
2826	CIN(K-1) = C3d(I,J,NKDE-K+1)
2827	ENDDO
2828	Y2(1 )=0.
2829	Y2(NKDE-1)=0.
2830	DO K=NKDS,NKDE ! target points in model interface levels (pint)
2831	PTMP(K) = ETA1(K)PDTOP + ETA2(K)PD(I,J) + PTOP
2832	ENDDO
2833	DO K=NKDS,NKDE-1 ! target points in model levels
2834	PIO(K) = EXP((ALOG(PTMP(K))+ALOG(PTMP(K+1)))*0.5)
2835	ENDDO
2836	IF(PTMP(1) .GE. PSTD(1))THEN ! if lower boundary is higher than PMSL(1) re-set lower boundary
2837	PIN(NKDE-1) = PIO(1) ! be consistent with target. This may not happen at all
2838	WRITE(0,*)'WARNING: NESTED DOMAIN PRESSURE AT LOWEST LEVEL HIGHER THAN PSTD'
2839	WRITE(0,*)'I,J,PIO(1),PSTD(1)',I,J,PIO(1),PSTD(1)
2840	ENDIF
2841
2842	CALL SPLINE2(I,J,JTB,NKDE-1,PIN,CIN,Y2,NKDE-1,PIO,COUT,DUM1,DUM2) ! interpolate
2843
2844	DO K=1,NKDE-1
2845	CWK3(I,J,K)= COUT(K) ! scalar in the nested domain
2846	ENDDO
2847	ENDDO
2848
2849	DO K = NKDS,NKDE-1
2850	DO I = NITS,MIN(NITE,NIDE-1)
2851	ntemp_b(i,J,K) = CWK3(I,J,K)
2852	ntemp_bt(i,J,K) = 0.0
2853	ENDDO
2854	ENDDO
2855
2856	ENDIF NMM_YS
2857
2858	! Y end boundary
2859
2860	NMM_YE: IF(NJTE-1 .EQ. NJDE-1)THEN
2861	! WRITE(0,*)'ENTERING Y END BOUNDARY AT T POINTS',NITS,MIN(NITE,NIDE-1)
2862	J = NJDE-1
2863	DO K=NKDS,NKDE-1
2864	DO I = NITS,MIN(NITE,NIDE-1)
2865	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
2866	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2867	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2868	+ HBWGT3(I,J)*CC3d(IIH(I,J), JJH(I,J)-1,K) &
2869	+ HBWGT4(I,J)*CC3d(IIH(I,J), JJH(I,J)+1,K)
2870	ELSE
2871	C3d(I,J,K) = HBWGT1(I,J)*CC3d(IIH(I,J), JJH(I,J) ,K) &
2872	+ HBWGT2(I,J)*CC3d(IIH(I,J)+1,JJH(I,J) ,K) &
2873	+ HBWGT3(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)-1,K) &
2874	+ HBWGT4(I,J)*CC3d(IIH(I,J)+1,JJH(I,J)+1,K)
2875
2876	ENDIF
2877	ENDDO
2878	ENDDO
2879
2880	DO I=NITS,MIN(NITE,NIDE-1)
2881	CIN=0.;PIN=0.;Y2=0;PIO=0.;PTMP=0.;COUT=0.;DUM1=0.;DUM2=0. ! clean up local array
2882	DO K=NKDS+1,NKDE ! inputs at standard levels
2883	PIN(K-1) = EXP((ALOG(PSTD(NKDE-K+1))+ALOG(PSTD(NKDE-K+2)))*0.5)
2884	CIN(K-1) = C3d(I,J,NKDE-K+1)
2885	ENDDO
2886	Y2(1 )=0.
2887	Y2(NKDE-1)=0.
2888	DO K=NKDS,NKDE ! target points in model interface levels (pint)
2889	PTMP(K) = ETA1(K)PDTOP + ETA2(K)PD(I,J) + PTOP
2890	ENDDO
2891	DO K=NKDS,NKDE-1 ! target points in model levels
2892	PIO(K) = EXP((ALOG(PTMP(K))+ALOG(PTMP(K+1)))*0.5)
2893	ENDDO
2894	IF(PTMP(1) .GE. PSTD(1))THEN ! if lower boundary is higher than PMSL(1) re-set lower boundary
2895	PIN(NKDE-1) = PIO(1) ! be consistent with target. This may not happen at all
2896	WRITE(0,*)'WARNING: NESTED DOMAIN PRESSURE AT LOWEST LEVEL HIGHER THAN PSTD'
2897	WRITE(0,*)'I,J,PIO(1),PSTD(1)',I,J,PIO(1),PSTD(1)
2898	ENDIF
2899
2900	CALL SPLINE2(I,J,JTB,NKDE-1,PIN,CIN,Y2,NKDE-1,PIO,COUT,DUM1,DUM2) ! interpolate
2901
2902	DO K=1,NKDE-1
2903	CWK4(I,J,K)= COUT(K) ! scalar in the nested domain
2904	ENDDO
2905	ENDDO
2906
2907	DO K = NKDS,NKDE-1
2908	DO I = NITS,MIN(NITE,NIDE-1)
2909	ntemp_b(i,J,K) = CWK4(I,J,K)
2910	ntemp_bt(i,J,K) = 0.0
2911	END DO
2912	END DO
2913
2914	ENDIF NMM_YE
2915
2916	END SUBROUTINE nmm_bdy_scalar
2917	!
2918	!
2919	!=======================================================================================
2920	SUBROUTINE SPLINE2(I,J,JTBX,NOLD,XOLD,YOLD,Y2,NNEW,XNEW,YNEW,P,Q)
2921	!
2922	! ******************************************************************
2923	! * *
2924	! * THIS IS A ONE-DIMENSIONAL CUBIC SPLINE FITTING ROUTINE *
2925	! * PROGRAMED FOR A SMALL SCALAR MACHINE. *
2926	! * *
2927	! * PROGRAMER Z. JANJIC *
2928	! * *
2929	! * NOLD - NUMBER OF GIVEN VALUES OF THE FUNCTION. MUST BE GE 3. *
2930	! * XOLD - LOCATIONS OF THE POINTS AT WHICH THE VALUES OF THE *
2931	! * FUNCTION ARE GIVEN. MUST BE IN ASCENDING ORDER. *
2932	! * YOLD - THE GIVEN VALUES OF THE FUNCTION AT THE POINTS XOLD. *
2933	! * Y2 - THE SECOND DERIVATIVES AT THE POINTS XOLD. IF NATURAL *
2934	! * SPLINE IS FITTED Y2(1)=0. AND Y2(NOLD)=0. MUST BE *
2935	! * SPECIFIED. *
2936	! * NNEW - NUMBER OF VALUES OF THE FUNCTION TO BE CALCULATED. *
2937	! * XNEW - LOCATIONS OF THE POINTS AT WHICH THE VALUES OF THE *
2938	! * FUNCTION ARE CALCULATED. XNEW(K) MUST BE GE XOLD(1) *
2939	! * AND LE XOLD(NOLD). *
2940	! * YNEW - THE VALUES OF THE FUNCTION TO BE CALCULATED. *
2941	! * P, Q - AUXILIARY VECTORS OF THE LENGTH NOLD-2. *
2942	! * *
2943	! ******************************************************************
2944	!---------------------------------------------------------------------
2945	IMPLICIT NONE
2946	!---------------------------------------------------------------------
2947	INTEGER,INTENT(IN) :: I,J,JTBX,NNEW,NOLD
2948	REAL,DIMENSION(JTBX),INTENT(IN) :: XNEW,XOLD,YOLD
2949	REAL,DIMENSION(JTBX),INTENT(INOUT) :: P,Q,Y2
2950	REAL,DIMENSION(JTBX),INTENT(OUT) :: YNEW
2951	!
2952	INTEGER :: II,JJ,K,K1,K2,KOLD,NOLDM1
2953	REAL :: AK,BK,CK,DEN,DX,DXC,DXL,DXR,DYDXL,DYDXR &
2954	,RDX,RTDXC,X,XK,XSQ,Y2K,Y2KP1
2955	!---------------------------------------------------------------------
2956
2957	! debug
2958
2959	II=9999
2960	JJ=9999
2961	IF(I.eq.II.and.J.eq.JJ)THEN
2962	WRITE(0,*)'DEBUG in SPLINE2: I,J',I,J
2963	WRITE(0,*)'DEBUG in SPLINE2:HSO= ',xnew(1:nold)
2964	DO K=1,NOLD
2965	WRITE(0,*)'DEBUG in SPLINE2:L,ZETAI,PINTI= ' &
2966	,K,YOLD(K),XOLD(K)
2967	ENDDO
2968	ENDIF
2969	!
2970	NOLDM1=NOLD-1
2971	!
2972	DXL=XOLD(2)-XOLD(1)
2973	DXR=XOLD(3)-XOLD(2)
2974	DYDXL=(YOLD(2)-YOLD(1))/DXL
2975	DYDXR=(YOLD(3)-YOLD(2))/DXR
2976	RTDXC=0.5/(DXL+DXR)
2977	!
2978	P(1)= RTDXC(6.(DYDXR-DYDXL)-DXL*Y2(1))
2979	Q(1)=-RTDXC*DXR
2980	!
2981	IF(NOLD.EQ.3)GO TO 150
2982	!---------------------------------------------------------------------
2983	K=3
2984	!
2985	100 DXL=DXR
2986	DYDXL=DYDXR
2987	DXR=XOLD(K+1)-XOLD(K)
2988	DYDXR=(YOLD(K+1)-YOLD(K))/DXR
2989	DXC=DXL+DXR
2990	DEN=1./(DXL*Q(K-2)+DXC+DXC)
2991	!
2992	P(K-1)= DEN(6.(DYDXR-DYDXL)-DXL*P(K-2))
2993	Q(K-1)=-DEN*DXR
2994	!
2995	K=K+1
2996	IF(K.LT.NOLD)GO TO 100
2997	!-----------------------------------------------------------------------
2998	150 K=NOLDM1
2999	!
3000	200 Y2(K)=P(K-1)+Q(K-1)*Y2(K+1)
3001	!
3002	K=K-1
3003	IF(K.GT.1)GO TO 200
3004	!-----------------------------------------------------------------------
3005	K1=1
3006	!
3007	300 XK=XNEW(K1)
3008	!
3009	DO 400 K2=2,NOLD
3010	!
3011	IF(XOLD(K2).GT.XK)THEN
3012	KOLD=K2-1
3013	GO TO 450
3014	ENDIF
3015	!
3016	400 CONTINUE
3017	!
3018	YNEW(K1)=YOLD(NOLD)
3019	GO TO 600
3020	!
3021	450 IF(K1.EQ.1)GO TO 500
3022	IF(K.EQ.KOLD)GO TO 550
3023	!
3024	500 K=KOLD
3025	!
3026	Y2K=Y2(K)
3027	Y2KP1=Y2(K+1)
3028	DX=XOLD(K+1)-XOLD(K)
3029	RDX=1./DX
3030	!
3031	AK=.1666667RDX(Y2KP1-Y2K)
3032	BK=0.5*Y2K
3033	CK=RDX(YOLD(K+1)-YOLD(K))-.1666667DX*(Y2KP1+Y2K+Y2K)
3034	!
3035	550 X=XK-XOLD(K)
3036	XSQ=X*X
3037	!
3038	YNEW(K1)=AKXSQX+BKXSQ+CKX+YOLD(K)
3039
3040	! debug
3041
3042	IF(I.eq.II.and.J.eq.JJ)THEN
3043	WRITE(0,*) 'DEBUG:: k1,xnew(k1),ynew(k1): ', K1,XNEW(k1),YNEW(k1)
3044	ENDIF
3045
3046	!
3047	600 K1=K1+1
3048	IF(K1.LE.NNEW)GO TO 300
3049
3050	RETURN
3051
3052	END SUBROUTINE SPLINE2
3053
3054	!=======================================================================================
3055	! E grid interpolation for H and V points
3056	!=======================================================================================
3057
3058	SUBROUTINE interp_h_nmm (cfld, & ! CD field
3059	cids, cide, ckds, ckde, cjds, cjde, &
3060	cims, cime, ckms, ckme, cjms, cjme, &
3061	cits, cite, ckts, ckte, cjts, cjte, &
3062	nfld, & ! ND field
3063	nids, nide, nkds, nkde, njds, njde, &
3064	nims, nime, nkms, nkme, njms, njme, &
3065	nits, nite, nkts, nkte, njts, njte, &
3066	shw, & ! stencil half width for interp
3067	imask, & ! interpolation mask
3068	xstag, ystag, & ! staggering of field
3069	ipos, jpos, & ! Position of lower left of nest in CD
3070	nri, nrj, & ! nest ratios
3071	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
3072	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
3073	CBWGT4, HBWGT4 ) ! dummys for weights
3074	USE module_timing
3075	IMPLICIT NONE
3076
3077	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3078	cims, cime, ckms, ckme, cjms, cjme, &
3079	cits, cite, ckts, ckte, cjts, cjte, &
3080	nids, nide, nkds, nkde, njds, njde, &
3081	nims, nime, nkms, nkme, njms, njme, &
3082	nits, nite, nkts, nkte, njts, njte, &
3083	shw, &
3084	ipos, jpos, &
3085	nri, nrj
3086	LOGICAL, INTENT(IN) :: xstag, ystag
3087
3088	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3089	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3090	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3091	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
3092	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3093	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
3094	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3095
3096	! local
3097	INTEGER i,j,k
3098	!
3099	!*** CHECK DOMAIN BOUNDS BEFORE INTERPOLATION
3100	!
3101	DO J=NJTS,MIN(NJTE,NJDE-1)
3102	DO I=NITS,MIN(NITE,NIDE-1)
3103	IF(IIH(i,j).LT.(CIDS-shw) .OR. IIH(i,j).GT.(CIDE+shw)) &
3104	CALL wrf_error_fatal ('hpoints:check domain bounds along x' )
3105	IF(JJH(i,j).LT.(CJDS-shw) .OR. JJH(i,j).GT.(CJDE+shw)) &
3106	CALL wrf_error_fatal ('hpoints:check domain bounds along y' )
3107	ENDDO
3108	ENDDO
3109	!
3110	!*** INDEX CONVENTIONS
3111	!*** HBWGT4
3112	!*** 4
3113	!***
3114	!***
3115	!***
3116	!*** h
3117	!*** 1 2
3118	!*** HBWGT1 HBWGT2
3119	!***
3120	!***
3121	!*** 3
3122	!*** HBWGT3
3123
3124	DO K=NKDS,NKDE
3125	DO J=NJTS,MIN(NJTE,NJDE-1)
3126	DO I=NITS,MIN(NITE,NIDE-1)
3127	IF(IMASK(I,J) .NE. 1)THEN
3128	!
3129	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3130	NFLD(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3131	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3132	+ HBWGT3(I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3133	+ HBWGT4(I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3134	ELSE
3135	NFLD(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3136	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3137	+ HBWGT3(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3138	+ HBWGT4(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3139	ENDIF
3140	!
3141	ENDIF
3142	ENDDO
3143	ENDDO
3144	ENDDO
3145
3146	END SUBROUTINE interp_h_nmm
3147	!
3148	SUBROUTINE interp_v_nmm (cfld, & ! CD field
3149	cids, cide, ckds, ckde, cjds, cjde, &
3150	cims, cime, ckms, ckme, cjms, cjme, &
3151	cits, cite, ckts, ckte, cjts, cjte, &
3152	nfld, & ! ND field
3153	nids, nide, nkds, nkde, njds, njde, &
3154	nims, nime, nkms, nkme, njms, njme, &
3155	nits, nite, nkts, nkte, njts, njte, &
3156	shw, & ! stencil half width for interp
3157	imask, & ! interpolation mask
3158	xstag, ystag, & ! staggering of field
3159	ipos, jpos, & ! Position of lower left of nest in CD
3160	nri, nrj, & ! nest ratios
3161	CII, IIV, CJJ, JJV, CBWGT1, VBWGT1, & ! south-western grid locs and weights
3162	CBWGT2, VBWGT2, CBWGT3, VBWGT3, & ! note that "C"ourse grid ones are
3163	CBWGT4, VBWGT4 ) ! dummys
3164	USE module_timing
3165	IMPLICIT NONE
3166
3167	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3168	cims, cime, ckms, ckme, cjms, cjme, &
3169	cits, cite, ckts, ckte, cjts, cjte, &
3170	nids, nide, nkds, nkde, njds, njde, &
3171	nims, nime, nkms, nkme, njms, njme, &
3172	nits, nite, nkts, nkte, njts, njte, &
3173	shw, &
3174	ipos, jpos, &
3175	nri, nrj
3176	LOGICAL, INTENT(IN) :: xstag, ystag
3177
3178	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3179	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3180	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3181	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: VBWGT1,VBWGT2,VBWGT3,VBWGT4
3182	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3183	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIV,JJV
3184	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3185
3186	! local
3187	INTEGER i,j,k
3188
3189
3190	!
3191	!*** CHECK DOMAIN BOUNDS BEFORE INTERPOLATION
3192	!
3193	DO J=NJTS,MIN(NJTE,NJDE-1)
3194	DO I=NITS,MIN(NITE,NIDE-1)
3195	IF(IIV(i,j).LT.(CIDS-shw) .OR. IIV(i,j).GT.(CIDE+shw)) &
3196	CALL wrf_error_fatal ('vpoints:check domain bounds along x' )
3197	IF(JJV(i,j).LT.(CJDS-shw) .OR. JJV(i,j).GT.(CJDE+shw)) &
3198	CALL wrf_error_fatal ('vpoints:check domain bounds along y' )
3199	ENDDO
3200	ENDDO
3201	!
3202	!*** INDEX CONVENTIONS
3203	!*** VBWGT4
3204	!*** 4
3205	!***
3206	!***
3207	!***
3208	!*** h
3209	!*** 1 2
3210	!*** VBWGT1 VBWGT2
3211	!***
3212	!***
3213	!*** 3
3214	!*** VBWGT3
3215
3216	DO K=NKDS,NKDE
3217	DO J=NJTS,MIN(NJTE,NJDE-1)
3218	DO I=NITS,MIN(NITE,NIDE-1)
3219	IF(IMASK(I,J) .NE. 1)THEN
3220
3221	IF(MOD(JJV(I,J),2) .NE. 0)THEN ! 1,3,5,7
3222	NFLD(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3223	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3224	+ VBWGT3(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)-1,K) &
3225	+ VBWGT4(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)+1,K)
3226	ELSE
3227	NFLD(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3228	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3229	+ VBWGT3(I,J)*CFLD(IIV(I,J), JJV(I,J)-1,K) &
3230	+ VBWGT4(I,J)*CFLD(IIV(I,J), JJV(I,J)+1,K)
3231	ENDIF
3232
3233	ENDIF
3234	ENDDO
3235	ENDDO
3236	ENDDO
3237
3238	END SUBROUTINE interp_v_nmm
3239	!
3240	!=======================================================================================
3241	! E grid nearest neighbour interpolation for H points.
3242	! This routine assumes cfld and nfld are in IJK
3243	!=======================================================================================
3244	!
3245	SUBROUTINE interp_hnear_nmm (cfld, & ! CD field
3246	cids, cide, ckds, ckde, cjds, cjde, &
3247	cims, cime, ckms, ckme, cjms, cjme, &
3248	cits, cite, ckts, ckte, cjts, cjte, &
3249	nfld, & ! ND field
3250	nids, nide, nkds, nkde, njds, njde, &
3251	nims, nime, nkms, nkme, njms, njme, &
3252	nits, nite, nkts, nkte, njts, njte, &
3253	shw, & ! stencil half width for interp
3254	imask, & ! interpolation mask
3255	xstag, ystag, & ! staggering of field
3256	ipos, jpos, & ! Position of lower left of nest in CD
3257	nri, nrj, & ! nest ratios
3258	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
3259	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
3260	CBWGT4, HBWGT4 ) ! just dummys
3261	USE module_timing
3262	IMPLICIT NONE
3263
3264	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3265	cims, cime, ckms, ckme, cjms, cjme, &
3266	cits, cite, ckts, ckte, cjts, cjte, &
3267	nids, nide, nkds, nkde, njds, njde, &
3268	nims, nime, nkms, nkme, njms, njme, &
3269	nits, nite, nkts, nkte, njts, njte, &
3270	shw, &
3271	ipos, jpos, &
3272	nri, nrj
3273	LOGICAL, INTENT(IN) :: xstag, ystag
3274
3275	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3276	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3277	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3278	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
3279	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3280	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
3281	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3282
3283	! local
3284
3285	LOGICAL FLIP
3286	INTEGER i,j,k,n
3287	REAL SUM,AMAXVAL
3288	REAL, DIMENSION (4, nims:nime, njms:njme ) :: NBWGT
3289
3290	!
3291	!*** INDEX CONVENTIONS
3292	!*** NBWGT4=0
3293	!*** 4
3294	!***
3295	!***
3296	!***
3297	!*** h
3298	!*** 1 2
3299	!*** NBWGT1=1 NBWGT2=0
3300	!***
3301	!***
3302	!*** 3
3303	!*** NBWGT3=0
3304
3305	DO J=NJTS,MIN(NJTE,NJDE-1)
3306	DO I=NITS,MIN(NITE,NIDE-1)
3307	IF(IMASK(I,J) .NE. 1)THEN
3308	NBWGT(1,I,J)=HBWGT1(I,J)
3309	NBWGT(2,I,J)=HBWGT2(I,J)
3310	NBWGT(3,I,J)=HBWGT3(I,J)
3311	NBWGT(4,I,J)=HBWGT4(I,J)
3312	ENDIF
3313	ENDDO
3314	ENDDO
3315
3316	DO J=NJTS,MIN(NJTE,NJDE-1)
3317	DO I=NITS,MIN(NITE,NIDE-1)
3318	IF(IMASK(I,J) .NE. 1)THEN
3319
3320	AMAXVAL=0.
3321	DO N=1,4
3322	AMAXVAL=amax1(NBWGT(N,I,J),AMAXVAL)
3323	ENDDO
3324
3325	FLIP=.TRUE.
3326	SUM=0.0
3327	DO N=1,4
3328	IF(AMAXVAL .EQ. NBWGT(N,I,J) .AND. FLIP)THEN
3329	NBWGT(N,I,J)=1.0
3330	FLIP=.FALSE.
3331	ELSE
3332	NBWGT(N,I,J)=0.0
3333	ENDIF
3334	SUM=SUM+NBWGT(N,I,J)
3335	IF(SUM .GT. 1.0)CALL wrf_error_fatal ( "horizontal interp error - interp_hnear_nmm" )
3336	ENDDO
3337
3338	ENDIF
3339	ENDDO
3340	ENDDO
3341
3342	DO K=NKDS,NKDE
3343	DO J=NJTS,MIN(NJTE,NJDE-1)
3344	DO I=NITS,MIN(NITE,NIDE-1)
3345	IF(IMASK(I,J) .NE. 1)THEN
3346	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3347	NFLD(I,J,K) = NBWGT(1,I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3348	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3349	+ NBWGT(3,I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3350	+ NBWGT(4,I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3351	ELSE
3352	NFLD(I,J,K) = NBWGT(1,I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3353	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3354	+ NBWGT(3,I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3355	+ NBWGT(4,I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3356	ENDIF
3357	ENDIF
3358	ENDDO
3359	ENDDO
3360	ENDDO
3361
3362	END SUBROUTINE interp_hnear_nmm
3363	!
3364	!=======================================================================================
3365	! E grid nearest neighbour interpolation for H points.
3366	! This routine assumes cfld and nfld are in IKJ or ILJ
3367	!=======================================================================================
3368	!
3369	SUBROUTINE interp_hnear_ikj_nmm (cfld, & ! CD field
3370	cids, cide, ckds, ckde, cjds, cjde, &
3371	cims, cime, ckms, ckme, cjms, cjme, &
3372	cits, cite, ckts, ckte, cjts, cjte, &
3373	nfld, & ! ND field
3374	nids, nide, nkds, nkde, njds, njde, &
3375	nims, nime, nkms, nkme, njms, njme, &
3376	nits, nite, nkts, nkte, njts, njte, &
3377	shw, & ! stencil half width for interp
3378	imask, & ! interpolation mask
3379	xstag, ystag, & ! staggering of field
3380	ipos, jpos, & ! Position of lower left of nest in CD
3381	nri, nrj, & ! nest ratios
3382	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
3383	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
3384	CBWGT4, HBWGT4 ) ! just dummys
3385	USE module_timing
3386	IMPLICIT NONE
3387
3388	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3389	cims, cime, ckms, ckme, cjms, cjme, &
3390	cits, cite, ckts, ckte, cjts, cjte, &
3391	nids, nide, nkds, nkde, njds, njde, &
3392	nims, nime, nkms, nkme, njms, njme, &
3393	nits, nite, nkts, nkte, njts, njte, &
3394	shw, &
3395	ipos, jpos, &
3396	nri, nrj
3397	LOGICAL, INTENT(IN) :: xstag, ystag
3398
3399	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfld
3400	REAL, DIMENSION ( nims:nime, nkms:nkme, njms:njme ) :: nfld
3401	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3402	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
3403	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3404	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
3405	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3406
3407	! local
3408
3409	LOGICAL FLIP
3410	INTEGER i,j,k,n
3411	REAL SUM,AMAXVAL
3412	REAL, DIMENSION (4, nims:nime, njms:njme ) :: NBWGT
3413
3414	!
3415	!*** INDEX CONVENTIONS
3416	!*** NBWGT4=0
3417	!*** 4
3418	!***
3419	!***
3420	!***
3421	!*** h
3422	!*** 1 2
3423	!*** NBWGT1=1 NBWGT2=0
3424	!***
3425	!***
3426	!*** 3
3427	!*** NBWGT3=0
3428
3429	DO J=NJTS,MIN(NJTE,NJDE-1)
3430	DO I=NITS,MIN(NITE,NIDE-1)
3431	IF(IMASK(I,J) .NE. 1)THEN
3432	NBWGT(1,I,J)=HBWGT1(I,J)
3433	NBWGT(2,I,J)=HBWGT2(I,J)
3434	NBWGT(3,I,J)=HBWGT3(I,J)
3435	NBWGT(4,I,J)=HBWGT4(I,J)
3436	ENDIF
3437	ENDDO
3438	ENDDO
3439
3440	DO J=NJTS,MIN(NJTE,NJDE-1)
3441	DO I=NITS,MIN(NITE,NIDE-1)
3442	IF(IMASK(I,J) .NE. 1)THEN
3443
3444	AMAXVAL=0.
3445	DO N=1,4
3446	AMAXVAL=amax1(NBWGT(N,I,J),AMAXVAL)
3447	ENDDO
3448
3449	FLIP=.TRUE.
3450	SUM=0.0
3451	DO N=1,4
3452	IF(AMAXVAL .EQ. NBWGT(N,I,J) .AND. FLIP)THEN
3453	NBWGT(N,I,J)=1.0
3454	FLIP=.FALSE.
3455	ELSE
3456	NBWGT(N,I,J)=0.0
3457	ENDIF
3458	SUM=SUM+NBWGT(N,I,J)
3459	IF(SUM .GT. 1.0)CALL wrf_error_fatal ( "horizontal interp error - interp_hnear_nmm" )
3460	ENDDO
3461
3462	ENDIF
3463	ENDDO
3464	ENDDO
3465
3466	DO J=NJTS,MIN(NJTE,NJDE-1)
3467	DO K=NKDS,NKDE
3468	DO I=NITS,MIN(NITE,NIDE-1)
3469	IF(IMASK(I,J) .NE. 1)THEN
3470	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3471	NFLD(I,K,J) = NBWGT(1,I,J)*CFLD(IIH(I,J), K,JJH(I,J) ) &
3472	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,K,JJH(I,J) ) &
3473	+ NBWGT(3,I,J)*CFLD(IIH(I,J), K,JJH(I,J)-1) &
3474	+ NBWGT(4,I,J)*CFLD(IIH(I,J), K,JJH(I,J)+1)
3475	ELSE
3476	NFLD(I,K,J) = NBWGT(1,I,J)*CFLD(IIH(I,J), K,JJH(I,J) ) &
3477	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,K,JJH(I,J) ) &
3478	+ NBWGT(3,I,J)*CFLD(IIH(I,J)+1,K,JJH(I,J)-1) &
3479	+ NBWGT(4,I,J)*CFLD(IIH(I,J)+1,K,JJH(I,J)+1)
3480	ENDIF
3481	ENDIF
3482	ENDDO
3483	ENDDO
3484	ENDDO
3485
3486	END SUBROUTINE interp_hnear_ikj_nmm
3487	!
3488	!=======================================================================================
3489	! E grid nearest neighbour interpolation for integer H points
3490	!=======================================================================================
3491	!
3492	SUBROUTINE interp_int_hnear_nmm (cfld, & ! CD field; integers
3493	cids, cide, ckds, ckde, cjds, cjde, &
3494	cims, cime, ckms, ckme, cjms, cjme, &
3495	cits, cite, ckts, ckte, cjts, cjte, &
3496	nfld, & ! ND field; integers
3497	nids, nide, nkds, nkde, njds, njde, &
3498	nims, nime, nkms, nkme, njms, njme, &
3499	nits, nite, nkts, nkte, njts, njte, &
3500	shw, & ! stencil half width for interp
3501	imask, & ! interpolation mask
3502	xstag, ystag, & ! staggering of field
3503	ipos, jpos, & ! lower left of nest in CD
3504	nri, nrj, & ! nest ratios
3505	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! s-w grid locs and weights
3506	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
3507	CBWGT4, HBWGT4 ) ! just dummys
3508	USE module_timing
3509	IMPLICIT NONE
3510
3511	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3512	cims, cime, ckms, ckme, cjms, cjme, &
3513	cits, cite, ckts, ckte, cjts, cjte, &
3514	nids, nide, nkds, nkde, njds, njde, &
3515	nims, nime, nkms, nkme, njms, njme, &
3516	nits, nite, nkts, nkte, njts, njte, &
3517	shw, &
3518	ipos, jpos, &
3519	nri, nrj
3520	LOGICAL, INTENT(IN) :: xstag, ystag
3521
3522	INTEGER, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3523	INTEGER, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3524	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3525	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
3526	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3527	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
3528	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3529
3530	! local
3531
3532	LOGICAL FLIP
3533	INTEGER i,j,k,n
3534	REAL SUM,AMAXVAL
3535	REAL, DIMENSION (4, nims:nime, njms:njme ) :: NBWGT
3536
3537	!
3538	!*** INDEX CONVENTIONS
3539	!*** NBWGT4=0
3540	!*** 4
3541	!***
3542	!***
3543	!***
3544	!*** h
3545	!*** 1 2
3546	!*** NBWGT1=1 NBWGT2=0
3547	!***
3548	!***
3549	!*** 3
3550	!*** NBWGT3=0
3551
3552	DO J=NJTS,MIN(NJTE,NJDE-1)
3553	DO I=NITS,MIN(NITE,NIDE-1)
3554	IF(IMASK(I,J) .NE. 1)THEN
3555	NBWGT(1,I,J)=HBWGT1(I,J)
3556	NBWGT(2,I,J)=HBWGT2(I,J)
3557	NBWGT(3,I,J)=HBWGT3(I,J)
3558	NBWGT(4,I,J)=HBWGT4(I,J)
3559	ENDIF
3560	ENDDO
3561	ENDDO
3562
3563	DO J=NJTS,MIN(NJTE,NJDE-1)
3564	DO I=NITS,MIN(NITE,NIDE-1)
3565	IF(IMASK(I,J) .NE. 1)THEN
3566
3567	AMAXVAL=0.
3568	DO N=1,4
3569	AMAXVAL=amax1(NBWGT(N,I,J),AMAXVAL)
3570	ENDDO
3571
3572	FLIP=.TRUE.
3573	SUM=0.0
3574	DO N=1,4
3575	IF(AMAXVAL .EQ. NBWGT(N,I,J) .AND. FLIP)THEN
3576	NBWGT(N,I,J)=1.0
3577	FLIP=.FALSE.
3578	ELSE
3579	NBWGT(N,I,J)=0.0
3580	ENDIF
3581	SUM=SUM+NBWGT(N,I,J)
3582	IF(SUM .GT. 1.0)CALL wrf_error_fatal ( "horizontal interp error - interp_hnear_nmm" )
3583	ENDDO
3584	!
3585	ENDIF
3586	ENDDO
3587	ENDDO
3588
3589	DO J=NJTS,MIN(NJTE,NJDE-1)
3590	DO K=NKTS,NKTS
3591	DO I=NITS,MIN(NITE,NIDE-1)
3592	IF(IMASK(I,J) .NE. 1)THEN
3593	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3594	NFLD(I,J,K) = NBWGT(1,I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3595	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3596	+ NBWGT(3,I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3597	+ NBWGT(4,I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3598	ELSE
3599	NFLD(I,J,K) = NBWGT(1,I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3600	+ NBWGT(2,I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3601	+ NBWGT(3,I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3602	+ NBWGT(4,I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3603	ENDIF
3604	ENDIF
3605	ENDDO
3606	ENDDO
3607	ENDDO
3608
3609	END SUBROUTINE interp_int_hnear_nmm
3610	!
3611	!--------------------------------------------------------------------------------------
3612	!
3613	SUBROUTINE nmm_bdy_hinterp (cfld, & ! CD field
3614	cids, cide, ckds, ckde, cjds, cjde, &
3615	cims, cime, ckms, ckme, cjms, cjme, &
3616	cits, cite, ckts, ckte, cjts, cjte, &
3617	nfld, & ! ND field
3618	nids, nide, nkds, nkde, njds, njde, &
3619	nims, nime, nkms, nkme, njms, njme, &
3620	nits, nite, nkts, nkte, njts, njte, &
3621	shw, & ! stencil half width
3622	imask, & ! interpolation mask
3623	xstag, ystag, & ! staggering of field
3624	ipos, jpos, & ! Position of lower left of nest in CD
3625	nri, nrj, & ! nest ratios
3626	c_bxs,n_bxs, &
3627	c_bxe,n_bxe, &
3628	c_bys,n_bys, &
3629	c_bye,n_bye, &
3630	c_btxs,n_btxs, &
3631	c_btxe,n_btxe, &
3632	c_btys,n_btys, &
3633	c_btye,n_btye, &
3634	CTEMP_B,NTEMP_B, & ! These temp arrays should be removed
3635	CTEMP_BT,NTEMP_BT, & ! later on
3636	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
3637	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
3638	CBWGT4, HBWGT4 ) ! dummys
3639
3640	! use module_state_description
3641	USE module_configure
3642	USE module_wrf_error
3643
3644	IMPLICIT NONE
3645
3646
3647	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3648	cims, cime, ckms, ckme, cjms, cjme, &
3649	cits, cite, ckts, ckte, cjts, cjte, &
3650	nids, nide, nkds, nkde, njds, njde, &
3651	nims, nime, nkms, nkme, njms, njme, &
3652	nits, nite, nkts, nkte, njts, njte, &
3653	shw, &
3654	ipos, jpos, &
3655	nri, nrj
3656
3657	LOGICAL, INTENT(IN) :: xstag, ystag
3658
3659	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3660	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3661	!
3662	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: ctemp_b,ctemp_bt
3663	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: ntemp_b,ntemp_bt
3664	!
3665	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3666	REAL, DIMENSION( * ), INTENT(INOUT) :: c_bxs,n_bxs,c_bxe,n_bxe,c_bys,n_bys,c_bye,n_bye
3667	REAL, DIMENSION( * ), INTENT(INOUT) :: c_btxs,n_btxs,c_btxe,n_btxe,c_btys,n_btys,c_btye,n_btye
3668	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3669	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
3670	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3671	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
3672
3673	! Local
3674
3675	INTEGER :: i,j,k
3676	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: cwk1,cwk2,cwk3,cwk4
3677
3678	! X start boundary
3679
3680	NMM_XS: IF(NITS .EQ. NIDS)THEN
3681	! WRITE(0,*)'ENTERING X1 START BOUNDARY AT MASS POINTS',NJTS,MIN(NJTE,NJDE-1)
3682	I = NIDS
3683	DO K = NKDS,NKDE
3684	DO J = NJTS,MIN(NJTE,NJDE-1)
3685	IF(MOD(J,2) .NE.0)THEN ! 1,3,5,7 of nested domain
3686	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
3687	CWK1(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3688	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3689	+ HBWGT3(I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3690	+ HBWGT4(I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3691
3692
3693	ELSE
3694	CWK1(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3695	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3696	+ HBWGT3(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3697	+ HBWGT4(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3698	ENDIF
3699	ELSE
3700	CWK1(I,J,K) = 0.0 ! even rows at mass points of the nested domain
3701	ENDIF
3702	ntemp_b(i,J,K) = CWK1(I,J,K)
3703	ntemp_bt(i,J,K) = 0.0
3704	END DO
3705	END DO
3706	ENDIF NMM_XS
3707
3708	! X end boundary
3709
3710	NMM_XE: IF(NITE-1 .EQ. NIDE-1)THEN
3711	! WRITE(0,*)'ENTERING X END BOUNDARY AT MASS POINTS',NJTS,MIN(NJTE,NJDE-1)
3712	I = NIDE-1
3713	DO K = NKDS,NKDE
3714	DO J = NJTS,MIN(NJTE,NJDE-1)
3715	IF(MOD(J,2) .NE.0)THEN ! 1,3,5,7 of the nested domain
3716	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
3717	CWK2(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3718	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3719	+ HBWGT3(I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3720	+ HBWGT4(I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3721	ELSE
3722	CWK2(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3723	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3724	+ HBWGT3(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3725	+ HBWGT4(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3726	ENDIF
3727	ELSE
3728	CWK2(I,J,K) = 0.0 ! even rows at mass points
3729	ENDIF
3730	ntemp_b(i,J,K) = CWK2(I,J,K)
3731	ntemp_bt(i,J,K) = 0.0
3732	END DO
3733	END DO
3734	ENDIF NMM_XE
3735
3736	! Y start boundary
3737
3738	NMM_YS: IF(NJTS .EQ. NJDS)THEN
3739	! WRITE(0,*)'ENTERING Y START BOUNDARY AT MASS POINTS',NITS,MIN(NITE,NIDE-1)
3740	J = NJDS
3741	DO K = NKDS, NKDE
3742	DO I = NITS,MIN(NITE,NIDE-1)
3743	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3744	CWK3(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3745	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3746	+ HBWGT3(I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3747	+ HBWGT4(I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3748	ELSE
3749	CWK3(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3750	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3751	+ HBWGT3(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3752	+ HBWGT4(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3753	ENDIF
3754	ntemp_b(i,J,K) = CWK3(I,J,K)
3755	ntemp_bt(i,J,K) = 0.0
3756	END DO
3757	END DO
3758	END IF NMM_YS
3759
3760	! Y end boundary
3761
3762	NMM_YE: IF(NJTE-1 .EQ. NJDE-1)THEN
3763	! WRITE(0,*)'ENTERING Y END BOUNDARY AT MASS POINTS',NITS,MIN(NITE,NIDE-1)
3764	J = NJDE-1
3765	DO K = NKDS,NKDE
3766	DO I = NITS,MIN(NITE,NIDE-1)
3767	IF(MOD(JJH(I,J),2) .NE. 0)THEN ! 1,3,5,7
3768	CWK4(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3769	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3770	+ HBWGT3(I,J)*CFLD(IIH(I,J), JJH(I,J)-1,K) &
3771	+ HBWGT4(I,J)*CFLD(IIH(I,J), JJH(I,J)+1,K)
3772	ELSE
3773	CWK4(I,J,K) = HBWGT1(I,J)*CFLD(IIH(I,J), JJH(I,J) ,K) &
3774	+ HBWGT2(I,J)*CFLD(IIH(I,J)+1,JJH(I,J) ,K) &
3775	+ HBWGT3(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)-1,K) &
3776	+ HBWGT4(I,J)*CFLD(IIH(I,J)+1,JJH(I,J)+1,K)
3777
3778	ENDIF
3779	ntemp_b(i,J,K) = CWK4(I,J,K)
3780	ntemp_bt(i,J,K) = 0.0
3781	END DO
3782	END DO
3783	END IF NMM_YE
3784
3785	RETURN
3786
3787	END SUBROUTINE nmm_bdy_hinterp
3788
3789	!--------------------------------------------------------------------------------------
3790
3791	SUBROUTINE nmm_bdy_vinterp ( cfld, & ! CD field
3792	cids, cide, ckds, ckde, cjds, cjde, &
3793	cims, cime, ckms, ckme, cjms, cjme, &
3794	cits, cite, ckts, ckte, cjts, cjte, &
3795	nfld, & ! ND field
3796	nids, nide, nkds, nkde, njds, njde, &
3797	nims, nime, nkms, nkme, njms, njme, &
3798	nits, nite, nkts, nkte, njts, njte, &
3799	shw, & ! stencil half width
3800	imask, & ! interpolation mask
3801	xstag, ystag, & ! staggering of field
3802	ipos, jpos, & ! Position of lower left of nest in CD
3803	nri, nrj, & ! nest ratios
3804	c_bxs,n_bxs, &
3805	c_bxe,n_bxe, &
3806	c_bys,n_bys, &
3807	c_bye,n_bye, &
3808	c_btxs,n_btxs, &
3809	c_btxe,n_btxe, &
3810	c_btys,n_btys, &
3811	c_btye,n_btye, &
3812	CTEMP_B,NTEMP_B, & ! These temp arrays should be removed
3813	CTEMP_BT,NTEMP_BT, & ! later on
3814	CII, IIV, CJJ, JJV, CBWGT1, VBWGT1, & ! south-western grid locs and weights
3815	CBWGT2, VBWGT2, CBWGT3, VBWGT3, & ! note that "C"ourse grid ones are
3816	CBWGT4, VBWGT4 ) ! dummys
3817
3818	! use module_state_description
3819	USE module_configure
3820	USE module_wrf_error
3821
3822	IMPLICIT NONE
3823
3824
3825	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3826	cims, cime, ckms, ckme, cjms, cjme, &
3827	cits, cite, ckts, ckte, cjts, cjte, &
3828	nids, nide, nkds, nkde, njds, njde, &
3829	nims, nime, nkms, nkme, njms, njme, &
3830	nits, nite, nkts, nkte, njts, njte, &
3831	shw, &
3832	ipos, jpos, &
3833	nri, nrj
3834
3835	LOGICAL, INTENT(IN) :: xstag, ystag
3836
3837	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
3838	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
3839	!
3840	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: ctemp_b,ctemp_bt
3841	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: ntemp_b,ntemp_bt
3842	!
3843	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
3844	REAL, DIMENSION( * ), INTENT(INOUT) :: c_bxs,n_bxs,c_bxe,n_bxe,c_bys,n_bys,c_bye,n_bye
3845	REAL, DIMENSION( * ), INTENT(INOUT) :: c_btxs,n_btxs,c_btxe,n_btxe,c_btys,n_btys,c_btye,n_btye
3846	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
3847	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: VBWGT1,VBWGT2,VBWGT3,VBWGT4
3848	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
3849	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIV,JJV
3850
3851	! Local
3852
3853	INTEGER :: i,j,k
3854	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: cwk1,cwk2,cwk3,cwk4
3855
3856	! X start boundary
3857
3858	NMM_XS: IF(NITS .EQ. NIDS)THEN
3859	! WRITE(0,*)'ENTERING X START BOUNDARY AT VELOCITY POINTS',NITS,NIDS,NJTS,MIN(NJTE,NJDE-1)
3860	I = NIDS
3861	DO K = NKDS,NKDE
3862	DO J = NJTS,MIN(NJTE,NJDE-1)
3863	IF(MOD(J,2) .EQ.0)THEN ! 1,3,5,7 of nested domain
3864	IF(MOD(JJV(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
3865	CWK1(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3866	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3867	+ VBWGT3(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)-1,K) &
3868	+ VBWGT4(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)+1,K)
3869	ELSE
3870	CWK1(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3871	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3872	+ VBWGT3(I,J)*CFLD(IIV(I,J), JJV(I,J)-1,K) &
3873	+ VBWGT4(I,J)*CFLD(IIV(I,J), JJV(I,J)+1,K)
3874	ENDIF
3875	ELSE
3876	CWK1(I,J,K) = 0.0 ! odd rows along J, at mass points have zero velocity
3877	ENDIF
3878	ntemp_b(i,J,K) = CWK1(I,J,K)
3879	ntemp_bt(i,J,K) = 0.0
3880	END DO
3881	END DO
3882	ENDIF NMM_XS
3883
3884	! X end boundary
3885
3886	NMM_XE: IF(NITE-1 .EQ. NIDE-1)THEN
3887	! WRITE(0,*)'ENTERING X END BOUNDARY AT VELOCITY POINTS',NITE-1,NIDE-1,NJTS,MIN(NJTE,NJDE-1)
3888	I = NIDE-1
3889	DO K = NKDS,NKDE
3890	DO J = NJTS,MIN(NJTE,NJDE-1)
3891	IF(MOD(J,2) .EQ.0)THEN ! 1,3,5,7 of the nested domain
3892	IF(MOD(JJV(I,J),2) .NE. 0)THEN ! 1,3,5,7 of the parent domain
3893	CWK2(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3894	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3895	+ VBWGT3(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)-1,K) &
3896	+ VBWGT4(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)+1,K)
3897	ELSE
3898	CWK2(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3899	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3900	+ VBWGT3(I,J)*CFLD(IIV(I,J), JJV(I,J)-1,K) &
3901	+ VBWGT4(I,J)*CFLD(IIV(I,J), JJV(I,J)+1,K)
3902	ENDIF
3903	ELSE
3904	CWK2(I,J,K) = 0.0 ! odd rows at mass points
3905	ENDIF
3906	ntemp_b(i,J,K) = CWK2(I,J,K)
3907	ntemp_bt(i,J,K) = 0.0
3908	END DO
3909	END DO
3910	ENDIF NMM_XE
3911
3912	! Y start boundary
3913
3914	NMM_YS: IF(NJTS .EQ. NJDS)THEN
3915	! WRITE(0,*)'ENTERING Y START BOUNDARY AT VELOCITY POINTS',NJTS,NJDS,NITS,MIN(NITE,NIDE-1)
3916	J = NJDS
3917	DO K = NKDS, NKDE
3918	DO I = NITS,MIN(NITE,NIDE-2) ! NIDE-1 SHOULD NOT MATTER IF WE FILL UP PHANTOM CELL
3919	IF(MOD(JJV(I,J),2) .NE. 0)THEN ! 1,3,5,7
3920	CWK3(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3921	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3922	+ VBWGT3(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)-1,K) &
3923	+ VBWGT4(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)+1,K)
3924	ELSE
3925	CWK3(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3926	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3927	+ VBWGT3(I,J)*CFLD(IIV(I,J), JJV(I,J)-1,K) &
3928	+ VBWGT4(I,J)*CFLD(IIV(I,J), JJV(I,J)+1,K)
3929	ENDIF
3930	ntemp_b(i,J,K) = CWK3(I,J,K)
3931	ntemp_bt(i,J,K) = 0.0
3932	END DO
3933	END DO
3934	END IF NMM_YS
3935
3936	! Y end boundary
3937
3938	NMM_YE: IF(NJTE-1 .EQ. NJDE-1)THEN
3939	! WRITE(0,*)'ENTERING Y END BOUNDARY AT VELOCITY POINTS',NJTE-1,NJDE-1,NITS,MIN(NITE,NIDE-1)
3940	J = NJDE-1
3941	DO K = NKDS,NKDE
3942	DO I = NITS,MIN(NITE,NIDE-2) ! NIDE-1 SHOULD NOT MATTER IF WE FILL UP PHANTOM CELL
3943	IF(MOD(JJV(I,J),2) .NE. 0)THEN ! 1,3,5,7
3944	CWK4(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3945	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3946	+ VBWGT3(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)-1,K) &
3947	+ VBWGT4(I,J)*CFLD(IIV(I,J)+1,JJV(I,J)+1,K)
3948	ELSE
3949	CWK4(I,J,K) = VBWGT1(I,J)*CFLD(IIV(I,J), JJV(I,J) ,K) &
3950	+ VBWGT2(I,J)*CFLD(IIV(I,J)+1,JJV(I,J) ,K) &
3951	+ VBWGT3(I,J)*CFLD(IIV(I,J), JJV(I,J)-1,K) &
3952	+ VBWGT4(I,J)*CFLD(IIV(I,J), JJV(I,J)+1,K)
3953	ENDIF
3954	ntemp_b(i,J,K) = CWK4(I,J,K)
3955	ntemp_bt(i,J,K) = 0.0
3956	END DO
3957	END DO
3958	END IF NMM_YE
3959
3960	RETURN
3961
3962	END SUBROUTINE nmm_bdy_vinterp
3963
3964	!
3965	!=======================================================================================
3966	! E grid interpolation: simple copy from parent to mother domain
3967	!=======================================================================================
3968	!
3969
3970	SUBROUTINE nmm_copy ( cfld, & ! CD field
3971	cids, cide, ckds, ckde, cjds, cjde, &
3972	cims, cime, ckms, ckme, cjms, cjme, &
3973	cits, cite, ckts, ckte, cjts, cjte, &
3974	nfld, & ! ND field
3975	nids, nide, nkds, nkde, njds, njde, &
3976	nims, nime, nkms, nkme, njms, njme, &
3977	nits, nite, nkts, nkte, njts, njte, &
3978	shw, & ! stencil half width
3979	imask, & ! interpolation mask
3980	xstag, ystag, & ! staggering of field
3981	ipos, jpos, & ! Position of lower left of nest in CD
3982	nri, nrj, & ! nest ratios
3983	CII, IIH, CJJ, JJH )
3984
3985	USE module_timing
3986	IMPLICIT NONE
3987
3988	LOGICAL, INTENT(IN) :: xstag, ystag
3989	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
3990	cims, cime, ckms, ckme, cjms, cjme, &
3991	cits, cite, ckts, ckte, cjts, cjte, &
3992	nids, nide, nkds, nkde, njds, njde, &
3993	nims, nime, nkms, nkme, njms, njme, &
3994	nits, nite, nkts, nkte, njts, njte, &
3995	shw, &
3996	ipos, jpos, &
3997	nri, nrj
3998	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ), INTENT(IN) :: cfld
3999	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ), INTENT(INOUT) :: nfld
4000	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: imask
4001	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
4002	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
4003
4004	! local
4005	INTEGER i,j,k
4006
4007	DO J=NJTS,MIN(NJTE,NJDE-1)
4008	DO K=NKTS,NKTE
4009	DO I=NITS,MIN(NITE,NIDE-1)
4010	NFLD(I,J,K) = CFLD(IIH(I,J),JJH(I,J),K)
4011	ENDDO
4012	ENDDO
4013	ENDDO
4014
4015	RETURN
4016
4017	END SUBROUTINE nmm_copy
4018	!
4019	!=======================================================================================
4020	! E grid test for mass point coincidence
4021	!=======================================================================================
4022	!
4023	SUBROUTINE test_nmm (cfld, & ! CD field
4024	cids, cide, ckds, ckde, cjds, cjde, &
4025	cims, cime, ckms, ckme, cjms, cjme, &
4026	cits, cite, ckts, ckte, cjts, cjte, &
4027	nfld, & ! ND field
4028	nids, nide, nkds, nkde, njds, njde, &
4029	nims, nime, nkms, nkme, njms, njme, &
4030	nits, nite, nkts, nkte, njts, njte, &
4031	shw, & ! stencil half width for interp
4032	imask, & ! interpolation mask
4033	xstag, ystag, & ! staggering of field
4034	ipos, jpos, & ! Position of lower left of nest in CD
4035	nri, nrj, & ! nest ratios
4036	CII, IIH, CJJ, JJH, CBWGT1, HBWGT1, & ! south-western grid locs and weights
4037	CBWGT2, HBWGT2, CBWGT3, HBWGT3, & ! note that "C"ourse grid ones are
4038	CBWGT4, HBWGT4 ) ! dummys for weights
4039	USE module_timing
4040	IMPLICIT NONE
4041
4042	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
4043	cims, cime, ckms, ckme, cjms, cjme, &
4044	cits, cite, ckts, ckte, cjts, cjte, &
4045	nids, nide, nkds, nkde, njds, njde, &
4046	nims, nime, nkms, nkme, njms, njme, &
4047	nits, nite, nkts, nkte, njts, njte, &
4048	shw, &
4049	ipos, jpos, &
4050	nri, nrj
4051	LOGICAL, INTENT(IN) :: xstag, ystag
4052
4053	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ) :: cfld
4054	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ) :: nfld
4055	REAL, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4 ! dummy
4056	REAL, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
4057	INTEGER, DIMENSION ( cims:cime, cjms:cjme ), INTENT(IN) :: CII,CJJ ! dummy
4058	INTEGER, DIMENSION ( nims:nime, njms:njme ), INTENT(IN) :: IIH,JJH
4059	INTEGER, DIMENSION ( nims:nime, njms:njme ) :: imask
4060
4061	! local
4062	INTEGER i,j,k
4063	REAL,PARAMETER :: error=0.0001,error1=1.0
4064	REAL :: diff
4065	!
4066	!*** CHECK DOMAIN BOUNDS BEFORE INTERPOLATION
4067	!
4068	DO J=NJTS,MIN(NJTE,NJDE-1)
4069	DO I=NITS,MIN(NITE,NIDE-1)
4070	IF(IIH(i,j).LT.(CIDS-shw) .OR. IIH(i,j).GT.(CIDE+shw)) &
4071	CALL wrf_error_fatal ('hpoints:check domain bounds along x' )
4072	IF(JJH(i,j).LT.(CJDS-shw) .OR. JJH(i,j).GT.(CJDE+shw)) &
4073	CALL wrf_error_fatal ('hpoints:check domain bounds along y' )
4074	ENDDO
4075	ENDDO
4076
4077	!
4078	!*** INDEX CONVENTIONS
4079	!*** HBWGT4
4080	!*** 4
4081	!***
4082	!***
4083	!***
4084	!*** h
4085	!*** 1 2
4086	!*** HBWGT1 HBWGT2
4087	!***
4088	!***
4089	!*** 3
4090	!*** HBWGT3
4091
4092
4093	! WRITE(0,*)NITS,MIN(NITE,NIDE-1),CITS,CITE
4094	DO J=NJTS,MIN(NJTE,NJDE-1)
4095	DO K=NKDS,NKDE
4096	DO I=NITS,MIN(NITE,NIDE-1)
4097	IF(ABS(1.0-HBWGT1(I,J)) .LE. ERROR)THEN
4098	DIFF=ABS(NFLD(I,J,K)-CFLD(IIH(I,J),JJH(I,J),K))
4099	IF(DIFF .GT. ERROR)THEN
4100	CALL wrf_debug(1,"dyn_nmm: NON-COINCIDENT, NESTED MASS POINT")
4101	WRITE(0,*)I,IIH(I,J),J,JJH(I,J),HBWGT1(I,J),NFLD(I,J,K),CFLD(IIH(I,J),JJH(I,J),K),DIFF
4102	ENDIF
4103	IF(DIFF .GT. ERROR1)THEN
4104	WRITE(0,*)I,IIH(I,J),J,JJH(I,J),HBWGT1(I,J),NFLD(I,J,K),CFLD(IIH(I,J),JJH(I,J),K),DIFF
4105	CALL wrf_error_fatal ('dyn_nmm: NON-COINCIDENT, NESTED MASS POINT')
4106	ENDIF
4107	ENDIF
4108	ENDDO
4109	ENDDO
4110	ENDDO
4111
4112	END SUBROUTINE test_nmm
4113
4114	!==================================
4115	! this is the default function used in nmm feedback at mass points.
4116
4117	SUBROUTINE nmm_feedback ( cfld, & ! CD field
4118	cids, cide, ckds, ckde, cjds, cjde, &
4119	cims, cime, ckms, ckme, cjms, cjme, &
4120	cits, cite, ckts, ckte, cjts, cjte, &
4121	nfld, & ! ND field
4122	nids, nide, nkds, nkde, njds, njde, &
4123	nims, nime, nkms, nkme, njms, njme, &
4124	nits, nite, nkts, nkte, njts, njte, &
4125	shw, & ! stencil half width for interp
4126	imask, & ! interpolation mask
4127	xstag, ystag, & ! staggering of field
4128	ipos, jpos, & ! Position of lower left of nest in CD
4129	nri, nrj, & ! nest ratios
4130	CII, IIH, CJJ, JJH, &
4131	CBWGT1, HBWGT1, CBWGT2, HBWGT2, &
4132	CBWGT3, HBWGT3, CBWGT4, HBWGT4 )
4133	USE module_configure
4134	IMPLICIT NONE
4135
4136
4137	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
4138	cims, cime, ckms, ckme, cjms, cjme, &
4139	cits, cite, ckts, ckte, cjts, cjte, &
4140	nids, nide, nkds, nkde, njds, njde, &
4141	nims, nime, nkms, nkme, njms, njme, &
4142	nits, nite, nkts, nkte, njts, njte, &
4143	shw, &
4144	ipos, jpos, &
4145	nri, nrj
4146	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
4147	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIH,JJH
4148	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4
4149	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: HBWGT1,HBWGT2,HBWGT3,HBWGT4
4150	LOGICAL, INTENT(IN) :: xstag, ystag
4151
4152	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ), INTENT(OUT) :: cfld
4153	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ), INTENT(IN) :: nfld
4154	INTEGER, DIMENSION ( nims:nime, njms:njme ),INTENT(IN) :: imask
4155
4156	! Local
4157
4158	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
4159	INTEGER :: icmin,icmax,jcmin,jcmax
4160	INTEGER :: is, ipoints,jpoints,ijpoints
4161	INTEGER , PARAMETER :: passes = 2
4162	REAL :: AVGH
4163
4164	!=====================================================================================
4165	!
4166
4167	IF(nri .ne. 3 .OR. nrj .ne. 3) &
4168	CALL wrf_error_fatal ('Feedback works for only 1:3 ratios, currently. Modify the namelist' )
4169
4170	! WRITE(0,*)'SIMPLE FEED BACK IS SWITCHED ON FOR MASS'
4171
4172	CFLD = 9999.0
4173
4174	DO ck = ckts, ckte
4175	nk = ck
4176	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4177	nj = (cj-jpos)*nrj + 1
4178	if(mod(cj,2) .eq. 0)THEN
4179	is=0 ! even rows for mass points (2,4,6,8)
4180	else
4181	is=1 ! odd rows for mass points (1,3,5,7)
4182	endif
4183	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4184	ni = (ci-ipos)*nri + 2 -is
4185	IF(IS==0)THEN ! (2,4,6,8)
4186	! AVGH = NFLD(NI,NJ+1,NK) + NFLD(NI,NJ-1,NK) + NFLD(NI+1,NJ+1,NK)+ NFLD(NI+1,NJ-1,NK) &
4187	! + NFLD(NI+1,NJ,NK) + NFLD(NI-1,NJ,NK) + NFLD(NI,NJ+2,NK) + NFLD(NI,NJ-2,NK) &
4188	! + NFLD(NI+1,NJ+2,NK)+ NFLD(NI-1,NJ+2,NK)+ NFLD(NI+1,NJ-2,NK)+ NFLD(NI-1,NJ-2,NK)
4189
4190	AVGH = NFLD(NI,NJ+2,NK) &
4191	+ NFLD(NI ,NJ+1,NK) + NFLD(NI+1,NJ+1,NK) &
4192	+ NFLD(NI-1,NJ ,NK) + NFLD(NI,NJ ,NK) + NFLD(NI+1,NJ ,NK) &
4193	+ NFLD(NI ,NJ-1,NK) + NFLD(NI+1,NJ-1,NK) &
4194	+ NFLD(NI,NJ-2,NK)
4195
4196	ELSE
4197	! AVGH = NFLD(NI,NJ+1,NK) + NFLD(NI,NJ-1,NK) + NFLD(NI-1,NJ+1,NK)+ NFLD(NI-1,NJ-1,NK) &
4198	! + NFLD(NI+1,NJ,NK) + NFLD(NI-1,NJ,NK) + NFLD(NI,NJ+2,NK) + NFLD(NI,NJ-2,NK) &
4199	! + NFLD(NI+1,NJ+2,NK)+ NFLD(NI-1,NJ+2,NK)+ NFLD(NI+1,NJ-2,NK)+ NFLD(NI-1,NJ-2,NK)
4200
4201	AVGH = NFLD(NI,NJ+2,NK) &
4202	+ NFLD(NI-1,NJ+1,NK) + NFLD(NI,NJ+1,NK) &
4203	+ NFLD(NI-1,NJ ,NK) + NFLD(NI,NJ ,NK) + NFLD(NI+1,NJ ,NK) &
4204	+ NFLD(NI-1,NJ-1,NK) + NFLD(NI,NJ-1,NK) &
4205	+ NFLD(NI,NJ-2,NK)
4206
4207	ENDIF
4208	!dusan CFLD(CI,CK,CJ) = 0.5CFLD(CI,CK,CJ) + 0.5(NFLD(NI,NK,NJ)+AVGH)/13.0
4209	! CFLD(CI,CJ,CK) = (NFLD(NI,NJ,NK)+AVGH)/13.0
4210	CFLD(CI,CJ,CK) = AVGH/9.0
4211	ENDDO
4212	ENDDO
4213	ENDDO
4214
4215	END SUBROUTINE nmm_feedback
4216
4217	!===========================================================================================
4218
4219	SUBROUTINE nmm_vfeedback ( cfld, & ! CD field
4220	cids, cide, ckds, ckde, cjds, cjde, &
4221	cims, cime, ckms, ckme, cjms, cjme, &
4222	cits, cite, ckts, ckte, cjts, cjte, &
4223	nfld, & ! ND field
4224	nids, nide, nkds, nkde, njds, njde, &
4225	nims, nime, nkms, nkme, njms, njme, &
4226	nits, nite, nkts, nkte, njts, njte, &
4227	shw, & ! stencil half width for interp
4228	imask, & ! interpolation mask
4229	xstag, ystag, & ! staggering of field
4230	ipos, jpos, & ! Position of lower left of nest in CD
4231	nri, nrj, & ! nest ratios
4232	CII, IIV, CJJ, JJV, &
4233	CBWGT1, VBWGT1, CBWGT2, VBWGT2, &
4234	CBWGT3, VBWGT3, CBWGT4, VBWGT4 )
4235	USE module_configure
4236	IMPLICIT NONE
4237
4238
4239	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
4240	cims, cime, ckms, ckme, cjms, cjme, &
4241	cits, cite, ckts, ckte, cjts, cjte, &
4242	nids, nide, nkds, nkde, njds, njde, &
4243	nims, nime, nkms, nkme, njms, njme, &
4244	nits, nite, nkts, nkte, njts, njte, &
4245	shw, &
4246	ipos, jpos, &
4247	nri, nrj
4248	INTEGER,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CII,CJJ ! dummy
4249	INTEGER,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: IIV,JJV
4250	REAL,DIMENSION(cims:cime,cjms:cjme), INTENT(IN) :: CBWGT1,CBWGT2,CBWGT3,CBWGT4
4251	REAL,DIMENSION(nims:nime,njms:njme), INTENT(IN) :: VBWGT1,VBWGT2,VBWGT3,VBWGT4
4252	LOGICAL, INTENT(IN) :: xstag, ystag
4253
4254	REAL, DIMENSION ( cims:cime, cjms:cjme, ckms:ckme ), INTENT(OUT) :: cfld
4255	REAL, DIMENSION ( nims:nime, njms:njme, nkms:nkme ), INTENT(IN) :: nfld
4256	INTEGER, DIMENSION ( nims:nime, njms:njme ),INTENT(IN) :: imask
4257
4258	! Local
4259
4260	INTEGER ci, cj, ck, ni, nj, nk, ip, jp, ioff, joff, ioffa, joffa
4261	INTEGER :: icmin,icmax,jcmin,jcmax
4262	INTEGER :: is, ipoints,jpoints,ijpoints
4263	INTEGER , PARAMETER :: passes = 2
4264	REAL :: AVGV
4265
4266	!=====================================================================================
4267	!
4268
4269	IF(nri .ne. 3 .OR. nrj .ne. 3) &
4270	CALL wrf_error_fatal ('Feedback works for only 1:3 ratios, currently. Modify the namelist')
4271
4272	! WRITE(0,*)'SIMPLE FEED BACK IS SWITCHED ON FOR VELOCITY'
4273
4274	CFLD = 9999.0
4275
4276	DO ck = ckts, ckte
4277	nk = ck
4278	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4279	nj = (cj-jpos)*nrj + 1
4280	if(mod(cj,2) .eq. 0)THEN
4281	is=1 ! even rows for velocity points (2,4,6,8)
4282	else
4283	is=0 ! odd rows for velocity points (1,3,5,7)
4284	endif
4285	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4286	ni = (ci-ipos)*nri + 2 -is
4287	IF(IS==0)THEN ! (1,3,5,7)
4288	! AVGV = NFLD(NI,NK,NJ+1) + NFLD(NI,NK,NJ-1) + NFLD(NI+1,NK,NJ+1)+ NFLD(NI+1,NK,NJ-1) &
4289	! + NFLD(NI+1,NK,NJ) + NFLD(NI-1,NK,NJ) + NFLD(NI,NK,NJ+2) + NFLD(NI,NK,NJ-2) &
4290	! + NFLD(NI+1,NK,NJ+2)+ NFLD(NI-1,NK,NJ+2)+ NFLD(NI+1,NK,NJ-2)+ NFLD(NI-1,NK,NJ-2)
4291
4292	AVGV = NFLD(NI,NJ+2,NK) &
4293	+ NFLD(NI ,NJ+1,NK) + NFLD(NI+1,NJ+1,NK) &
4294	+ NFLD(NI-1,NJ ,NK) + NFLD(NI,NJ ,NK) + NFLD(NI+1,NJ ,NK) &
4295	+ NFLD(NI ,NJ-1,NK) + NFLD(NI+1,NJ-1,NK) &
4296	+ NFLD(NI,NJ-2,NK)
4297
4298	ELSE
4299	! AVGV = NFLD(NI,NK,NJ+1) + NFLD(NI,NK,NJ-1) + NFLD(NI-1,NK,NJ+1)+ NFLD(NI-1,NK,NJ-1) &
4300	! + NFLD(NI+1,NK,NJ) + NFLD(NI-1,NK,NJ) + NFLD(NI,NK,NJ+2) + NFLD(NI,NK,NJ-2) &
4301	! + NFLD(NI+1,NK,NJ+2)+ NFLD(NI-1,NK,NJ+2)+ NFLD(NI+1,NK,NJ-2)+ NFLD(NI-1,NK,NJ-2)
4302
4303	AVGV = NFLD(NI,NJ+2,NK) &
4304	+ NFLD(NI-1,NJ+1,NK) + NFLD(NI,NJ+1,NK) &
4305	+ NFLD(NI-1,NJ ,NK) + NFLD(NI,NJ ,NK) + NFLD(NI+1,NJ ,NK) &
4306	+ NFLD(NI-1,NJ-1,NK) + NFLD(NI,NJ-1,NK) &
4307	+ NFLD(NI,NJ-2,NK)
4308
4309	ENDIF
4310	!dusan CFLD(CI,CK,CJ) = 0.5CFLD(CI,CK,CJ) + 0.5(NFLD(NI,NK,NJ)+AVGV)/13.0
4311	! CFLD(CI,CK,CJ) = (NFLD(NI,NK,NJ)+AVGV)/13.0
4312	CFLD(CI,CJ,CK) = AVGV/9.0
4313	ENDDO
4314	ENDDO
4315	ENDDO
4316
4317	END SUBROUTINE nmm_vfeedback
4318
4319
4320	SUBROUTINE nmm_smoother ( cfld , &
4321	cids, cide, ckds, ckde, cjds, cjde, &
4322	cims, cime, ckms, ckme, cjms, cjme, &
4323	cits, cite, ckts, ckte, cjts, cjte, &
4324	nids, nide, nkds, nkde, njds, njde, &
4325	nims, nime, nkms, nkme, njms, njme, &
4326	nits, nite, nkts, nkte, njts, njte, &
4327	xstag, ystag, &
4328	ipos, jpos, &
4329	nri, nrj &
4330	)
4331
4332	USE module_configure
4333	IMPLICIT NONE
4334
4335	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
4336	cims, cime, ckms, ckme, cjms, cjme, &
4337	cits, cite, ckts, ckte, cjts, cjte, &
4338	nids, nide, nkds, nkde, njds, njde, &
4339	nims, nime, nkms, nkme, njms, njme, &
4340	nits, nite, nkts, nkte, njts, njte, &
4341	nri, nrj, &
4342	ipos, jpos
4343	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(INOUT) :: cfld
4344	LOGICAL, INTENT(IN) :: xstag, ystag
4345
4346
4347	! Local
4348
4349	INTEGER :: feedback
4350	INTEGER, PARAMETER :: smooth_passes = 5
4351
4352	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfldnew
4353	INTEGER :: ci, cj, ck
4354	INTEGER :: is, npass
4355	REAL :: AVGH
4356
4357	RETURN
4358	! If there is no feedback, there can be no smoothing.
4359
4360	CALL nl_get_feedback ( 1, feedback )
4361	IF ( feedback == 0 ) RETURN
4362
4363	WRITE(0,*)'SIMPLE SMOOTHER IS SWITCHED ON FOR HEIGHT'
4364
4365	DO npass = 1, smooth_passes
4366
4367	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4368	if(mod(cj,2) .eq. 0)THEN
4369	is=0 ! even rows for mass points (2,4,6,8)
4370	else
4371	is=1 ! odd rows for mass points (1,3,5,7)
4372	endif
4373	DO ck = ckts, ckte
4374	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4375	IF(IS==0)THEN ! (2,4,6,8)
4376	AVGH = CFLD(CI,CK,CJ+1) + CFLD(CI,CK,CJ-1) + CFLD(CI+1,CK,CJ+1) + CFLD(CI+1,CK,CJ-1)
4377	ELSE
4378	AVGH = CFLD(CI,CK,CJ+1) + CFLD(CI,CK,CJ-1) + CFLD(CI-1,CK,CJ+1) + CFLD(CI-1,CK,CJ-1)
4379	ENDIF
4380	CFLDNEW(CI,CK,CJ) = (AVGH + 4*CFLD(CI,CK,CJ)) / 8.0
4381	ENDDO
4382	ENDDO
4383	ENDDO
4384
4385	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4386	if(mod(cj,2) .eq. 0)THEN
4387	is=0 ! even rows for mass points (2,4,6,8)
4388	else
4389	is=1 ! odd rows for mass points (1,3,5,7)
4390	endif
4391	DO ck = ckts, ckte
4392	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4393	CFLD(CI,CK,CJ) = CFLDNEW(CI,CK,CJ)
4394	ENDDO
4395	ENDDO
4396	ENDDO
4397
4398	ENDDO ! do npass
4399
4400	END SUBROUTINE nmm_smoother
4401
4402
4403	SUBROUTINE nmm_vsmoother ( cfld , &
4404	cids, cide, ckds, ckde, cjds, cjde, &
4405	cims, cime, ckms, ckme, cjms, cjme, &
4406	cits, cite, ckts, ckte, cjts, cjte, &
4407	nids, nide, nkds, nkde, njds, njde, &
4408	nims, nime, nkms, nkme, njms, njme, &
4409	nits, nite, nkts, nkte, njts, njte, &
4410	xstag, ystag, &
4411	ipos, jpos, &
4412	nri, nrj &
4413	)
4414
4415	USE module_configure
4416	IMPLICIT NONE
4417
4418	INTEGER, INTENT(IN) :: cids, cide, ckds, ckde, cjds, cjde, &
4419	cims, cime, ckms, ckme, cjms, cjme, &
4420	cits, cite, ckts, ckte, cjts, cjte, &
4421	nids, nide, nkds, nkde, njds, njde, &
4422	nims, nime, nkms, nkme, njms, njme, &
4423	nits, nite, nkts, nkte, njts, njte, &
4424	nri, nrj, &
4425	ipos, jpos
4426	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ), INTENT(INOUT) :: cfld
4427	LOGICAL, INTENT(IN) :: xstag, ystag
4428
4429
4430	! Local
4431
4432	INTEGER :: feedback
4433	INTEGER, PARAMETER :: smooth_passes = 5
4434
4435	REAL, DIMENSION ( cims:cime, ckms:ckme, cjms:cjme ) :: cfldnew
4436	INTEGER :: ci, cj, ck
4437	INTEGER :: is, npass
4438	REAL :: AVGV
4439
4440	RETURN
4441	! If there is no feedback, there can be no smoothing.
4442
4443	CALL nl_get_feedback ( 1, feedback )
4444	IF ( feedback == 0 ) RETURN
4445
4446	WRITE(0,*)'SIMPLE SMOOTHER IS SWITCHED ON FOR VELOCITY'
4447
4448	DO npass = 1, smooth_passes
4449
4450	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4451	if(mod(cj,2) .eq. 0)THEN
4452	is=1 ! even rows for mass points (2,4,6,8)
4453	else
4454	is=0 ! odd rows for mass points (1,3,5,7)
4455	endif
4456	DO ck = ckts, ckte
4457	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4458	IF(IS==0)THEN ! (2,4,6,8)
4459	AVGV = CFLD(CI,CK,CJ+1) + CFLD(CI,CK,CJ-1) + CFLD(CI+1,CK,CJ+1) + CFLD(CI+1,CK,CJ-1)
4460	ELSE
4461	AVGV = CFLD(CI,CK,CJ+1) + CFLD(CI,CK,CJ-1) + CFLD(CI-1,CK,CJ+1) + CFLD(CI-1,CK,CJ-1)
4462	ENDIF
4463	CFLDNEW(CI,CK,CJ) = (AVGV + 4*CFLD(CI,CK,CJ)) / 8.0
4464	ENDDO
4465	ENDDO
4466	ENDDO
4467
4468	DO cj = MAX(jpos+1,cjts),MIN(jpos+(njde-njds)/nrj-1,cjte) ! exclude top and bottom BCs
4469	if(mod(cj,2) .eq. 0)THEN
4470	is=1 ! even rows for mass points (2,4,6,8)
4471	else
4472	is=0 ! odd rows for mass points (1,3,5,7)
4473	endif
4474	DO ck = ckts, ckte
4475	DO ci = MAX(ipos+is,cits),MIN(ipos+(nide-nids)/nri-1,cite) ! excludes LBCs
4476	CFLD(CI,CK,CJ) = CFLDNEW(CI,CK,CJ)
4477	ENDDO
4478	ENDDO
4479	ENDDO
4480
4481	ENDDO
4482
4483	END SUBROUTINE nmm_vsmoother
4484	!======================================================================================
4485	! End of gopal's doing
4486	!======================================================================================
4487	#endif

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format