Context Navigation

← Previous Revision
Latest Revision
Next Revision →
Normal
Revision Log

interp_fcn.F @ 3026

Last change on this file since 3026 was 2759, checked in by aslmd, 3 years ago
adding unmodified code from WRFV3.0.1.1, expurged from useless data +1M size
File size: 197.1 KB

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

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

Download in other formats:

Original Format