31MODULE fms_diag_bbox_mod
33 USE fms_mod,
ONLY: error_mesg, fatal, fms_error_handler, string
102 integer,
intent(in) :: subregion_start
103 integer,
intent(in) :: subregion_end
105 integer,
intent(in) :: dim
107 integer :: block_start
113 block_start = bounds%imin
114 block_end = bounds%imax
116 block_start = bounds%jmin
117 block_end = bounds%jmax
119 block_start = bounds%kmin
120 block_end = bounds%kmax
123 if (block_start < subregion_start .and. block_end < subregion_start)
then
128 if (block_start > subregion_end .and. block_end > subregion_end)
then
180 subroutine update_index(this, starting_index, ending_index, dim, ignore_halos)
182 integer,
intent(in) :: starting_index
183 integer,
intent(in) :: ending_index
184 integer,
intent(in) :: dim
185 logical,
intent(in) :: ignore_halos
194 if (ignore_halos)
then
203 this%imin = starting_index + nhalox
204 this%imax = ending_index + nhalox
206 this%jmin = starting_index + nhaloy
207 this%jmax = ending_index + nhaloy
209 this%kmin = starting_index
210 this%kmax = ending_index
216 pure integer function get_hi (this)
result(rslt)
223 pure integer function get_hj (this)
result(rslt)
230 pure integer function get_fis (this)
result(rslt)
237 pure integer function get_fie (this)
result(rslt)
244 pure integer function get_fjs (this)
result(rslt)
251 pure integer function get_fje (this)
result(rslt)
259 integer,
intent(in) :: lower_val
260 integer,
intent(in) :: upper_val
261 this%imin = lower_val
262 this%jmin = lower_val
263 this%kmin = lower_val
264 this%imax = upper_val
265 this%jmax = upper_val
266 this%kmax = upper_val
272 SUBROUTINE update_bounds(this, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k)
274 INTEGER,
INTENT(in) :: lower_i
275 INTEGER,
INTENT(in) :: upper_i
276 INTEGER,
INTENT(in) :: lower_j
277 INTEGER,
INTENT(in) :: upper_j
278 INTEGER,
INTENT(in) :: lower_k
279 INTEGER,
INTENT(in) :: upper_k
280 this%imin = min(this%imin, lower_i)
281 this%imax = max(this%imax, upper_i)
282 this%jmin = min(this%jmin, lower_j)
283 this%jmax = max(this%jmax, upper_j)
284 this%kmin = min(this%kmin, lower_k)
285 this%kmax = max(this%kmax, upper_k)
290 function set_bounds(this, field_data, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k, has_halos) &
293 class(*),
intent(in) :: field_data(:,:,:,:)
294 INTEGER,
INTENT(in) :: lower_i
295 INTEGER,
INTENT(in) :: upper_i
296 INTEGER,
INTENT(in) :: lower_j
297 INTEGER,
INTENT(in) :: upper_j
298 INTEGER,
INTENT(in) :: lower_k
299 INTEGER,
INTENT(in) :: upper_k
300 LOGICAL,
INTENT(in) :: has_halos
302 character(len=150) :: error_msg
311 this%has_halos = has_halos
317 nhalos_2 = ubound(field_data,1)-(upper_i-lower_i+1)
318 if (mod(nhalos_2, 2) .ne. 0)
then
319 error_msg =
"There are non-symmetric halos in the first dimension"
323 this%nhalo_I = nhalox
325 nhalos_2 = ubound(field_data,2)-(upper_j-lower_j + 1)
326 if (mod(nhalos_2, 2) .ne. 0)
then
327 error_msg =
"There are non-symmetric halos in the second dimension"
331 this%nhalo_J = nhaloy
333 this%imin = 1 + nhalox
334 this%imax = ubound(field_data,1) - nhalox
335 this%jmin = 1 + nhaloy
336 this%jmax = ubound(field_data,2) - nhaloy
349 class(*),
INTENT( in),
DIMENSION(:,:,:,:) :: array
350 this%imin = lbound(array,1)
351 this%imax = ubound(array,1)
352 this%jmin = lbound(array,2)
353 this%jmax = ubound(array,2)
354 this%kmin = lbound(array,3)
355 this%kmax = ubound(array,3)
357 this%has_halos = .false.
366 CLASS(*),
INTENT( in),
DIMENSION(:,:,:,:,:) :: array
367 this%imin = lbound(array,1)
368 this%imax = ubound(array,1)
369 this%jmin = lbound(array,2)
370 this%jmax = ubound(array,2)
371 this%kmin = lbound(array,3)
372 this%kmax = ubound(array,3)
380 ie_in, je_in, ke_in, err_msg)
result(ierr)
383 class(*),
intent(in) :: field(:,:,:,:)
384 integer,
intent(in),
optional :: is_in, js_in, ks_in, ie_in, je_in, ke_in
385 character(len=*),
intent(out),
optional :: err_msg
388 integer :: is, js, ks, ie, je, ke
391 integer :: twohi, twohj
392 integer :: fis, fie, fjs, fje
393 integer :: n1, n2, n3
396 if (
present(err_msg)) err_msg =
''
403 IF (
PRESENT(is_in) ) is = is_in
404 IF (
PRESENT(js_in) ) js = js_in
405 IF (
PRESENT(ks_in) ) ks = ks_in
415 IF (
PRESENT(ie_in) ) ie = ie_in
416 IF (
PRESENT(je_in) ) je = je_in
417 IF (
PRESENT(ke_in) ) ke = ke_in
419 twohi = n1 - (ie - is + 1)
420 IF ( mod(twohi, 2) /= 0 )
THEN
421 ierr = fms_error_handler(
'diag_util_mod:recondition_indices', &
422 'non-symmetric halos in first dimension', err_msg)
426 twohj = n2 - (je - js + 1)
427 IF ( mod(twohj, 2) /= 0 )
THEN
428 ierr = fms_error_handler(
'diag_util_mod:recondition_indices', &
429 'non-symmetric halos in second dimension', err_msg)
438 IF (
PRESENT(ie_in) .AND.
PRESENT(je_in) )
THEN
452 indices%bounds3D%imin = is
453 indices%bounds3D%imax = ie
454 indices%bounds3D%jmin = js
455 indices%bounds3D%jmax = je
456 indices%bounds3D%kmin = ks
457 indices%bounds3D%kmax = ke
470 integer,
intent(in) :: starting
471 integer,
intent(in) :: ending
472 integer,
intent(in) :: dim
484 bounds_out%imin = max(starting, bounds_out%imin)-starting+1
485 bounds_out%imax = min(bounds_out%imax, bounds_out%imin + ending-starting)
487 bounds_out%jmin = max(starting, bounds_out%jmin)-starting+1
488 bounds_out%jmax = min(bounds_out%jmax, bounds_out%jmin + ending-starting)
490 bounds_out%kmin =max(starting, bounds_out%kmin)-starting+1
491 bounds_out%kmax = min(bounds_out%kmax, bounds_out%kmin + ending-starting)
501 integer,
intent(in) :: starting
502 integer,
intent(in) :: ending
503 integer,
intent(in) :: dim
514 bounds_in%imin = min(abs(starting-bounds%imin+1), starting)
515 bounds_in%imax = min(bounds_in%imax, (bounds_in%imin + ending-starting))
517 bounds_in%jmin = min(abs(starting-bounds%jmin+1), starting)
518 bounds_in%jmax = min(bounds_in%jmax, (bounds_in%jmin + ending-starting))
520 bounds_in%kmin = min(abs(starting-bounds%kmin+1), starting)
521 bounds_in%kmax = min(bounds_in%kmax, (bounds_in%kmin + ending-starting))
525 END MODULE fms_diag_bbox_mod
pure integer function get_hi(this)
Gets the halo size of fmsDiagBoundsHalos_type in the I dimension.
integer, parameter ydimension
Parameter defining the y dimension.
integer, parameter xdimension
Parameter defining the x dimension.
pure integer function get_fje(this)
Gets the updated index ‘fje’ of fmsDiagBoundsHalos_type in the I dimension.
pure integer function get_fie(this)
Gets the updated index ‘fie’ of fmsDiagBoundsHalos_type in the I dimension.
pure integer function get_kmin(this)
Gets kmin of fmsDiagIbounds_type.
pure integer function get_jmin(this)
Gets jmin of fmsDiagIbounds_type.
subroutine rebase_output(bounds_out, starting, ending, dim)
Rebase the ouput bounds for a given dimension based on the starting and ending indices of a subregion...
logical function, public recondition_indices(indices, field, is_in, js_in, ks_in, ie_in, je_in, ke_in, err_msg)
Updates indices based on presence/absence of input indices is, js, ks, ie, je, and ke.
pure integer function get_kmax(this)
Gets kmax of fmsDiagIbounds_type.
pure integer function get_imin(this)
Gets imin of fmsDiagIbounds_type.
logical pure function, public determine_if_block_is_in_region(subregion_start, subregion_end, bounds, dim)
The PEs grid points are divided further into "blocks". This function determines if a block.
pure integer function get_fis(this)
Gets the updated index ‘fis’ of fmsDiagBoundsHalos_type in the I dimension.
subroutine update_index(this, starting_index, ending_index, dim, ignore_halos)
Updates the starting and ending index of a given dimension.
pure integer function get_jmax(this)
Gets jmax of fmsDiagIbounds_type.
character(len=150) function set_bounds(this, field_data, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k, has_halos)
Sets the bounds of a bounding region.
subroutine reset_bounds(this, lower_val, upper_val)
Reset the instance bounding lower and upper bounds to lower_val and upper_val, respectively.
subroutine reset_bounds_from_array_4d(this, array)
Reset the instance bounding box with the bounds determined from the first three dimensions of the 5D ...
subroutine reset_bounds_from_array_5d(this, array)
Reset the instance bounding box with the bounds determined from the first three dimensions of the 5D ...
subroutine rebase_input(bounds_in, bounds, starting, ending, dim)
Rebase the input bounds for a given dimension based on the starting and ending indices of a subregion...
pure integer function get_fjs(this)
Gets the updated index ‘fjs’ of fmsDiagBoundsHalos_type in the I dimension.
pure integer function get_hj(this)
Gets the halo size of fmsDiagBoundsHalos_type in the J dimension.
pure integer function get_imax(this)
Gets imax of fmsDiagIbounds_type.
integer, parameter zdimension
Parameter defininf the z dimension.
subroutine update_bounds(this, lower_i, upper_i, lower_j, upper_j, lower_k, upper_k)
Update the the first three (normally x, y, and z) min and max boundaries (array indices) of the insta...
Data structure holding starting and ending indices in the I, J, and K dimensions. It also has extra m...
Data structure holding a 3D bounding box. It is commonlyused to represent the interval bounds or limi...