1 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
4 ! This module provides routines for smoothing.
5 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
14 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
17 ! Purpose: Apply the 1-2-1 smoother from the MM5 program TERRAIN
18 ! (found in smth121.F) to array.
19 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
20 subroutine one_two_one(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
21 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
)
26 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
27 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
28 integer, intent(in
) :: npass
29 real, intent(in
) :: msgval
30 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), intent(inout
) :: array
33 integer :: ix
, iy
, iz
, ipass
34 real, pointer, dimension(:,:,:) :: scratch
36 allocate(scratch(start_x
+1:end_x
-1, start_y
:end_y
, start_z
:end_z
))
41 do ix
=start_x
+1,end_x
-1
43 scratch(ix
,iy
,iz
) = 0.50*array(ix
,iy
,iz
)+0.25*(array(ix
-1,iy
,iz
)+array(ix
+1,iy
,iz
))
48 do iy
=start_y
+1,end_y
-1
49 do ix
=start_x
+1,end_x
-1
51 array(ix
,iy
,iz
) = 0.50*scratch(ix
,iy
,iz
)+0.25*(scratch(ix
,iy
-1,iz
)+scratch(ix
,iy
+1,iz
))
56 call exchange_halo_r(array
, &
57 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
58 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
64 end subroutine one_two_one
67 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
70 ! Purpose: Apply the smoother-desmoother from the MM5 program TERRAIN
71 ! (found in smther.F) to array.
72 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
73 subroutine smth_desmth(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
74 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
)
79 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
80 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
81 integer, intent(in
) :: npass
82 real, intent(in
) :: msgval
83 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), intent(inout
) :: array
86 integer :: ix
, iy
, iz
, ipass
87 real, pointer, dimension(:,:,:) :: scratch
89 allocate(scratch(start_x
+1:end_x
-1, start_y
:end_y
, start_z
:end_z
))
97 do ix
=start_x
+1,end_x
-1
99 scratch(ix
,iy
,iz
) = 0.5*array(ix
,iy
,iz
) + 0.25*(array(ix
-1,iy
,iz
)+array(ix
+1,iy
,iz
))
104 do iy
=start_y
+1,end_y
-1
105 do ix
=start_x
+1,end_x
-1
107 array(ix
,iy
,iz
) = 0.5*scratch(ix
,iy
,iz
) + 0.25*(scratch(ix
,iy
-1,iz
)+scratch(ix
,iy
+1,iz
))
112 call exchange_halo_r(array
, &
113 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
114 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
120 do ix
=start_x
+1,end_x
-1
122 scratch(ix
,iy
,iz
) = 1.52*array(ix
,iy
,iz
) - 0.26*(array(ix
-1,iy
,iz
)+array(ix
+1,iy
,iz
))
127 do iy
=start_y
+1,end_y
-1
128 do ix
=start_x
+1,end_x
-1
130 array(ix
,iy
,iz
) = 1.52*scratch(ix
,iy
,iz
) - 0.26*(scratch(ix
,iy
-1,iz
)+scratch(ix
,iy
+1,iz
))
135 call exchange_halo_r(array
, &
136 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
137 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
143 end subroutine smth_desmth
146 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
147 ! Name: smth_desmth_special
149 ! Purpose: Apply the smoother-desmoother from the MM5 program TERRAIN
150 ! (found in smther.F) to array; however, any grid points that were not
151 ! originally negative but which have been smoothed to a negative value
152 ! will be restored to their original values.
153 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
154 subroutine smth_desmth_special(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
155 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
)
160 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
161 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
162 integer, intent(in
) :: npass
163 real, intent(in
) :: msgval
164 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), intent(inout
) :: array
167 integer :: ix
, iy
, iz
, ipass
168 real, pointer, dimension(:,:,:) :: scratch
, orig_array
170 allocate(scratch(start_x
+1:end_x
-1, start_y
:end_y
, start_z
:end_z
))
171 allocate(orig_array(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
))
181 do ix
=start_x
+1,end_x
-1
183 scratch(ix
,iy
,iz
) = 0.5*array(ix
,iy
,iz
) + 0.25*(array(ix
-1,iy
,iz
)+array(ix
+1,iy
,iz
))
188 do iy
=start_y
+1,end_y
-1
189 do ix
=start_x
+1,end_x
-1
191 array(ix
,iy
,iz
) = 0.5*scratch(ix
,iy
,iz
) + 0.25*(scratch(ix
,iy
-1,iz
)+scratch(ix
,iy
+1,iz
))
196 call exchange_halo_r(array
, &
197 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
198 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
204 do ix
=start_x
+1,end_x
-1
206 scratch(ix
,iy
,iz
) = 1.52*array(ix
,iy
,iz
) - 0.26*(array(ix
-1,iy
,iz
)+array(ix
+1,iy
,iz
))
211 do iy
=start_y
+1,end_y
-1
212 do ix
=start_x
+1,end_x
-1
214 array(ix
,iy
,iz
) = 1.52*scratch(ix
,iy
,iz
) - 0.26*(scratch(ix
,iy
-1,iz
)+scratch(ix
,iy
+1,iz
))
219 call exchange_halo_r(array
, &
220 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
221 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
225 ! Remove artificially negative values
229 if (array(ix
,iy
,iz
) < 0. .and
. orig_array(ix
,iy
,iz
) >= 0.) then
230 array(ix
,iy
,iz
) = orig_array(ix
,iy
,iz
)
237 deallocate(orig_array
)
239 end subroutine smth_desmth_special
243 ! Smoothing routines for E-grid, contributed by Matthew Pyle
246 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
247 ! Name: one_two_one_egrid
249 ! Purpose: Apply the 1-2-1 smoother from the MM5 program TERRAIN
250 ! (found in smth121.F) to array.
251 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
252 subroutine one_two_one_egrid(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
253 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
, hflag
)
258 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
259 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
260 integer, intent(in
) :: npass
261 real, intent(in
) :: msgval
, hflag
262 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), intent(inout
) :: array
265 integer :: ix
, iy
, iz
, ipass
266 real, pointer, dimension(:,:,:) :: scratch
267 integer, dimension(start_y
:end_y
) :: ihe
, ihw
, istart
269 allocate(scratch(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
))
272 if (hflag
== 1.0) then
273 ihe(iy
) = abs(mod(iy
+1,2))
276 ! assign ive,ivw equivs to ihe,ihw
277 ihe(iy
) = abs(mod(iy
,2))
283 if (hflag
== 1.0) then
284 if (mod(iy
,2) == 0) then
287 istart(iy
) = start_x
+1
290 if (abs(mod(iy
,2)) == 1) then
293 istart(iy
) = start_x
+1
302 scratch(ix
,iy
,1) = array(ix
,iy
,1) ! for points used in 2nd computation but not defined in 1st computation
307 do iy
=start_y
+1,end_y
-1
308 do ix
=istart(iy
),end_x
-1
310 if ( (msgval
== 1.0 .and
. array(ix
,iy
,iz
) /= 0.) .or
. msgval
/= 1.0) then
311 scratch(ix
,iy
,iz
) = 0.50*array(ix
,iy
,iz
)+ &
312 0.25*(array(ix
+ihw(iy
),iy
-1,iz
)+array(ix
+ihe(iy
),iy
+1,iz
))
319 do iy
=start_y
+1,end_y
-1
320 do ix
=istart(iy
),end_x
-1
322 if ( (msgval
== 1.0 .and
. array(ix
,iy
,iz
) /= 0.) .or
. msgval
/= 1.0) then
323 array(ix
,iy
,iz
) = 0.50*scratch(ix
,iy
,iz
)+ &
324 0.25*(scratch(ix
+ihe(iy
),iy
-1,iz
)+scratch(ix
+ihw(iy
),iy
+1,iz
))
330 call exchange_halo_r(array
, &
331 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
332 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
338 end subroutine one_two_one_egrid
341 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
342 ! Name: smth_desmth_egrid_old
344 ! Purpose: Apply the smoother-desmoother for E grid
345 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
346 subroutine smth_desmth_egrid_old(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
347 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
, hflag
)
352 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
353 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
354 integer, intent(in
) :: npass
355 real, intent(in
) :: msgval
, hflag
356 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), &
357 intent(inout
) :: array
360 integer :: ix
, iy
, iz
, ipass
361 real, pointer, dimension(:,:,:) :: scratch
362 integer, dimension(start_y
:end_y
) :: ihe
, ihw
, istart
363 real, parameter:: cenwgt
= 1.52
364 real, parameter:: endwgt
= 0.13
366 allocate(scratch(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
))
369 if (hflag
== 1.0) then
370 ihe(iy
) = abs(mod(iy
+1,2))
373 ! assign ive,ivw equivs to ihe,ihw
374 ihe(iy
) = abs(mod(iy
,2))
380 if (hflag
== 1.0) then
381 if (mod(iy
,2) == 0) then
384 istart(iy
) = start_x
+1
387 if (abs(mod(iy
,2)) == 1) then
390 istart(iy
) = start_x
+1
403 scratch(ix
,iy
,1) = array(ix
,iy
,1)
407 do iy
=start_y
+1,end_y
-1
408 do ix
=istart(iy
),end_x
-1
410 if ( (msgval
== 1.0 .and
. array(ix
,iy
,iz
) /= 0.) .or
. msgval
/= 1.0) then
411 scratch(ix
,iy
,iz
) = 0.50*array(ix
,iy
,iz
)+ &
412 0.125*(array(ix
+ihw(iy
),iy
-1,iz
)+array(ix
+ihe(iy
),iy
+1,iz
)+ &
413 array(ix
+ihw(iy
),iy
+1,iz
)+array(ix
+ihe(iy
),iy
-1,iz
))
424 do iy
=start_y
+2,end_y
-2
425 do ix
=istart(iy
),end_x
-1
427 if ( (msgval
== 1.0 .and
. scratch(ix
,iy
,iz
) /= 0.) .or
. msgval
/= 1.0) then
428 array(ix
,iy
,iz
) = cenwgt
*scratch(ix
,iy
,iz
) - &
429 endwgt
*(scratch(ix
+ihw(iy
),iy
-1,iz
)+scratch(ix
+ihe(iy
),iy
+1,iz
) + &
430 scratch(ix
+ihw(iy
),iy
+1,iz
)+scratch(ix
+ihe(iy
),iy
-1,iz
))
440 end subroutine smth_desmth_egrid_old
443 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
444 ! Name: smth_desmth_egrid
446 ! Purpose: Apply the smoother-desmoother for E grid
447 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
448 subroutine smth_desmth_egrid(array
, start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, &
449 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, npass
, msgval
, hflag
)
454 integer, intent(in
) :: start_dom_x
, start_dom_y
, start_x
, start_y
, start_z
455 integer, intent(in
) :: end_dom_x
, end_dom_y
, end_x
, end_y
, end_z
456 integer, intent(in
) :: npass
457 real, intent(in
) :: msgval
, hflag
458 real, dimension(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
), &
459 intent(inout
) :: array
462 integer :: ix
, iy
, iz
, ipass
463 real, pointer, dimension(:,:,:) :: scratch
464 integer, dimension(start_y
:end_y
) :: ihe
, ihw
, istart
465 real, parameter :: cenwgt
= 1.52
466 real, parameter :: endwgt
= 0.26
468 allocate(scratch(start_x
:end_x
, start_y
:end_y
, start_z
:end_z
))
472 if (hflag
.eq
. 1.0) then
473 ihe(iy
)=abs(mod(iy
+1,2))
476 ! assign ive,ivw equivs to ihe,ihw
478 ihe(iy
)=abs(mod(iy
,2))
487 if (hflag
.eq
. 1.0) then
488 if (mod(iy
,2) .eq
. 0) then
495 if (abs(mod(iy
,2)) .eq
. 1) then
514 scratch(ix
,iy
,1)=array(ix
,iy
,1) ! for points used in 2nd computation but
520 do iy
=start_y
+1,end_y
-1
521 do ix
=istart(iy
),end_x
-1
523 if ( (msgval
.eq
. 1.0 .and
. array(ix
,iy
,iz
) .ne
. 0.) .or
. msgval
.ne
. 1.0) then
524 scratch(ix
,iy
,iz
) = 0.50*array(ix
,iy
,iz
)+ &
525 0.25*(array(ix
+ihw(iy
),iy
-1,iz
)+array(ix
+ihe(iy
),iy
+1,iz
))
532 do iy
=start_y
+1,end_y
-1
533 do ix
=istart(iy
),end_x
-1
535 if ( (msgval
.eq
. 1.0 .and
. array(ix
,iy
,iz
) .ne
. 0.) .or
. msgval
.ne
. 1.0) then
536 array(ix
,iy
,iz
) = 0.50*scratch(ix
,iy
,iz
)+ &
537 0.25*(scratch(ix
+ihe(iy
),iy
-1,iz
)+scratch(ix
+ihw(iy
),iy
+1,iz
))
543 call exchange_halo_r(array
, &
544 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
545 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
554 do iy
=start_y
+2,end_y
-2
555 do ix
=istart(iy
),end_x
-1
557 if ( (msgval
.eq
. 1.0 .and
. array(ix
,iy
,iz
) .ne
. 0.) .or
. msgval
.ne
. 1.0) then
558 scratch(ix
,iy
,iz
) = cenwgt
*array(ix
,iy
,iz
) - &
559 endwgt
*(array(ix
+ihw(iy
),iy
-1,iz
)+array(ix
+ihe(iy
),iy
+1,iz
))
566 do iy
=start_y
+2,end_y
-2
567 do ix
=istart(iy
),end_x
-1
569 if ( (msgval
.eq
. 1.0 .and
. array(ix
,iy
,iz
) .ne
. 0.) .or
. msgval
.ne
. 1.0) then
570 array(ix
,iy
,iz
) = cenwgt
*scratch(ix
,iy
,iz
) - &
571 endwgt
*(scratch(ix
+ihe(iy
),iy
-1,iz
)+scratch(ix
+ihw(iy
),iy
+1,iz
))
577 call exchange_halo_r(array
, &
578 start_x
, end_x
, start_y
, end_y
, start_z
, end_z
, &
579 start_dom_x
, end_dom_x
, start_dom_y
, end_dom_y
, start_z
, end_z
)
585 end subroutine smth_desmth_egrid
587 end module smooth_module