29 module fms_diag_axis_object_mod
34 use platform_mod,
only: r8_kind, r4_kind, i4_kind, i8_kind
41 use fms2_io_mod,
only: fmsnetcdffile_t, fmsnetcdfdomainfile_t, fmsnetcdfunstructureddomainfile_t, &
43 use fms_diag_yaml_mod,
only:
subregion_type, diag_yaml, max_subaxes
67 procedure :: get_ntiles
88 INTEGER ,
private :: axis_id
106 class(fmsDiagAxis_type),
allocatable :: axis
112 CHARACTER(len=:),
ALLOCATABLE ,
private :: subaxis_name
113 INTEGER ,
private :: starting_index
115 INTEGER ,
private :: ending_index
117 INTEGER ,
private :: parent_axis_id
118 INTEGER ,
private :: compute_idx(2)
119 INTEGER,
allocatable,
private :: global_idx(:)
120 real(kind=r4_kind),
allocatable,
private :: zbounds(:)
132 INTEGER ,
private :: ndiurnal_samples
133 CHARACTER(len=:),
ALLOCATABLE,
private :: axis_name
134 CHARACTER(len=:),
ALLOCATABLE,
private :: long_name
135 CHARACTER(len=:),
ALLOCATABLE,
private :: units
136 INTEGER ,
private :: edges_id
137 CHARACTER(len=:),
ALLOCATABLE,
private :: edges_name
138 CLASS(*),
ALLOCATABLE,
private :: diurnal_data(:)
148 CHARACTER(len=:),
ALLOCATABLE,
private :: axis_name
149 CHARACTER(len=:),
ALLOCATABLE,
private :: units
150 CHARACTER(len=:),
ALLOCATABLE,
private :: long_name
151 CHARACTER(len=1) ,
private :: cart_name
152 CLASS(*),
ALLOCATABLE,
private :: axis_data(:)
153 CHARACTER(len=:),
ALLOCATABLE,
private :: type_of_data
155 integer,
ALLOCATABLE,
private :: subaxis(:)
156 integer ,
private :: nsubaxis
158 INTEGER ,
private :: type_of_domain
160 INTEGER ,
private :: length
161 INTEGER ,
private :: direction
162 INTEGER,
ALLOCATABLE,
private :: edges_id
164 CHARACTER(len=:),
ALLOCATABLE,
private :: edges_name
166 CHARACTER(len=:),
ALLOCATABLE,
private :: aux
168 CHARACTER(len=128) ,
private :: req
169 INTEGER ,
private :: tile_count
171 INTEGER ,
private :: num_attributes
172 INTEGER ,
private :: domain_position
173 integer,
allocatable ,
private :: structured_ids(:)
175 CHARACTER(len=:),
ALLOCATABLE,
private :: set_name
205 & set_name, Domain, Domain2, DomainU, aux, req, tile_count, domain_position, axis_length )
207 CHARACTER(len=*),
INTENT(in) :: axis_name
208 class(*),
INTENT(in) :: axis_data(:)
209 CHARACTER(len=*),
INTENT(in) :: units
210 CHARACTER(len=1),
INTENT(in) :: cart_name
211 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: long_name
212 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: set_name
213 INTEGER,
INTENT(in),
OPTIONAL :: direction
214 TYPE(
domain1d),
INTENT(in),
OPTIONAL :: Domain
215 TYPE(
domain2d),
INTENT(in),
OPTIONAL :: Domain2
216 TYPE(
domainug),
INTENT(in),
OPTIONAL :: DomainU
217 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: aux
219 CHARACTER(len=*),
INTENT(in),
OPTIONAL :: req
220 INTEGER,
INTENT(in),
OPTIONAL :: tile_count
221 INTEGER,
INTENT(in),
OPTIONAL :: domain_position
222 integer,
intent(in),
optional :: axis_length
224 this%axis_name = trim(axis_name)
225 this%units = trim(units)
226 this%cart_name = uppercase(cart_name)
229 if (
present(long_name)) this%long_name = trim(long_name)
231 select type (axis_data)
232 type is (real(kind=r8_kind))
233 allocate(real(kind=r8_kind) :: this%axis_data(axis_length))
234 this%axis_data = axis_data
235 this%length = axis_length
236 this%type_of_data =
"double"
237 type is (real(kind=r4_kind))
238 allocate(real(kind=r4_kind) :: this%axis_data(axis_length))
239 this%axis_data = axis_data
240 this%length = axis_length
241 this%type_of_data =
"float"
243 call mpp_error(fatal,
"The axis_data in your diag_axis_init call is not a supported type. &
244 & Currently only r4 and r8 data is supported.")
248 if (
present(domain))
then
249 if (
present(domain2) .or.
present(domainu))
call mpp_error(fatal, &
250 "The presence of Domain with any other domain type is prohibited. "//&
251 "Check you diag_axis_init call for axis_name:"//trim(axis_name))
253 call this%axis_domain%set(domain=domain)
254 else if (
present(domain2))
then
255 if (
present(domainu))
call mpp_error(fatal, &
256 "The presence of Domain2 with any other domain type is prohibited. "//&
257 "Check you diag_axis_init call for axis_name:"//trim(axis_name))
259 call this%axis_domain%set(domain2=domain2)
261 else if (
present(domainu))
then
263 call this%axis_domain%set(domainu=domainu)
268 if (
present(tile_count)) this%tile_count = tile_count
270 this%domain_position = center
271 if (
present(domain_position)) this%domain_position = domain_position
275 if (
present(direction)) this%direction = direction
278 if (
present(aux)) this%aux = trim(aux)
279 if (
present(req)) this%req = trim(req)
281 if (
present(set_name)) this%set_name = trim(set_name)
283 if (max_subaxes .gt. 0)
then
284 allocate(this%subaxis(max_subaxes))
285 this%subaxis = diag_null
289 this%num_attributes = 0
295 character(len=*),
intent(in) :: att_name
296 class(*),
intent(in) :: att_value(:)
300 if (.not.
allocated(this%attributes)) &
303 this%num_attributes = this%num_attributes + 1
305 j = this%num_attributes
306 call this%attributes(j)%add(att_name, att_value)
312 class(fmsnetcdffile_t),
INTENT(INOUT) :: fms2io_fileobj
313 logical,
INTENT(IN) :: edges_in_file
319 character(len=:),
ALLOCATABLE :: axis_edges_name
320 character(len=:),
pointer :: axis_name
321 integer :: axis_length
325 integer :: type_of_domain
326 logical :: is_subaxis
327 logical :: needs_domain_decomposition
329 integer :: domain_decomposition(4)
332 needs_domain_decomposition = .false.
336 axis_name => this%axis_name
337 axis_length = this%length
339 type_of_domain = this%type_of_domain
342 axis_name => this%subaxis_name
343 axis_length = this%ending_index - this%starting_index + 1
344 if (
allocated(this%global_idx))
then
345 needs_domain_decomposition = .true.
346 domain_decomposition(1:2) = this%global_idx
347 domain_decomposition(3) = this%starting_index
348 domain_decomposition(4) = this%ending_index
351 if (
present(parent_axis))
then
352 select type(parent_axis)
354 diag_axis => parent_axis
359 call this%write_diurnal_metadata(fms2io_fileobj)
364 select type (fms2io_fileobj)
367 type is (fmsnetcdffile_t)
370 call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
371 if (needs_domain_decomposition)
then
372 call register_variable_attribute(fms2io_fileobj, axis_name,
"domain_decomposition", &
373 domain_decomposition)
375 type is (fmsnetcdfdomainfile_t)
376 select case (type_of_domain)
381 call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
384 call register_axis(fms2io_fileobj, axis_name, diag_axis%cart_name, domain_position=diag_axis%domain_position)
385 call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
387 type is (fmsnetcdfunstructureddomainfile_t)
388 select case (type_of_domain)
392 call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
397 call register_field(fms2io_fileobj, axis_name, diag_axis%type_of_data, (/axis_name/))
402 if(
allocated(diag_axis%long_name)) &
403 call register_variable_attribute(fms2io_fileobj, axis_name,
"long_name", diag_axis%long_name, &
404 str_len=len_trim(diag_axis%long_name))
406 if (diag_axis%cart_name .NE.
"N") &
407 call register_variable_attribute(fms2io_fileobj, axis_name,
"axis", diag_axis%cart_name, str_len=1)
409 if (trim(diag_axis%units) .NE.
"none") &
410 call register_variable_attribute(fms2io_fileobj, axis_name,
"units", diag_axis%units, &
411 str_len=len_trim(diag_axis%units))
413 select case (diag_axis%direction)
415 call register_variable_attribute(fms2io_fileobj, axis_name,
"positive",
"up", str_len=2)
417 call register_variable_attribute(fms2io_fileobj, axis_name,
"positive",
"down", str_len=4)
421 if (.not. edges_in_file .and.
allocated(diag_axis%edges_name) .and. .not. is_subaxis)
then
422 call register_variable_attribute(fms2io_fileobj, axis_name,
"edges", diag_axis%edges_name, &
423 str_len=len_trim(diag_axis%edges_name))
426 if(
allocated(diag_axis%attributes))
then
427 do i = 1, diag_axis%num_attributes
428 select type (att_value => diag_axis%attributes(i)%att_value)
429 type is (
character(len=*))
430 call register_variable_attribute(fms2io_fileobj, axis_name, diag_axis%attributes(i)%att_name, &
431 trim(att_value(1)), str_len=len_trim(att_value(1)))
433 call register_variable_attribute(fms2io_fileobj, axis_name, diag_axis%attributes(i)%att_name, att_value)
443 class(fmsnetcdffile_t),
INTENT(INOUT) :: fms2io_fileobj
448 integer :: global_io_index(2)
451 call this%get_global_io_domain(global_io_index, fms2io_fileobj%is_file_using_netcdf_mpi())
452 call write_data(fms2io_fileobj, this%axis_name, this%axis_data(global_io_index(1):global_io_index(2)))
454 i = this%starting_index
455 j = this%ending_index
457 if (
present(parent_axis))
then
458 select type(parent_axis)
460 call write_data(fms2io_fileobj, this%subaxis_name, parent_axis%axis_data(i:j))
464 call write_data(fms2io_fileobj, this%axis_name, this%diurnal_data)
472 integer,
intent(inout) :: naxis
474 integer,
intent(in) :: n_diurnal_samples
476 logical,
intent(in) :: is_edges
479 CHARACTER(32) :: axis_name
480 CHARACTER(32) :: long_name
481 CHARACTER(32) :: edges_name
482 CHARACTER(128) :: units
483 real(kind=r8_kind),
allocatable :: diurnal_data(:)
492 WRITE (axis_name,
'(a,i2.2)')
'time_of_day_edges_', n_diurnal_samples
493 long_name =
"time of day edges"
494 allocate(diurnal_data(n_diurnal_samples + 1))
495 diurnal_data(1) = 0.0
497 do i = 1, n_diurnal_samples
498 diurnal_data(i+1) = 24.0* real(i)/n_diurnal_samples
501 WRITE (axis_name,
'(a,i2.2)')
'time_of_day_', n_diurnal_samples
502 long_name =
"time of day"
503 allocate(diurnal_data(n_diurnal_samples))
505 do i = 1, n_diurnal_samples
506 diurnal_data(i) = 24.0*(real(i)-0.5)/n_diurnal_samples
508 WRITE (edges_name,
'(a,i2.2)')
'time_of_day_edges_', n_diurnal_samples
513 11
FORMAT(a,
' since ',i4.4,
'-',i2.2,
'-',i2.2,
' ',i2.2,
':',i2.2,
':',i2.2)
516 select type (diurnal_axis => diag_axis(naxis)%axis)
518 diurnal_axis%axis_id = naxis
519 diurnal_axis%ndiurnal_samples = n_diurnal_samples
520 diurnal_axis%axis_name = trim(axis_name)
521 diurnal_axis%long_name = trim(long_name)
522 diurnal_axis%units = trim(units)
523 diurnal_axis%diurnal_data = diurnal_data
524 diurnal_axis%edges_id = edges_id
526 WRITE (edges_name,
'(a,i2.2)')
'time_of_day_edges_', n_diurnal_samples
527 diurnal_axis%edges_name = trim(edges_name)
546 integer,
intent(in) :: axis_ids(2)
550 allocate(this%structured_ids(2))
551 this%structured_ids = axis_ids
565 rslt = this%structured_ids
578 if (
allocated(this%edges_id))
get_edges_id = this%edges_id
585 integer,
intent(out) :: global_io_index(2)
586 logical,
intent(in) :: use_collective_writes
588 type(
domain2d),
pointer :: io_domain
590 global_io_index(1) = 1
591 global_io_index(2) = this%length
593 if (
allocated(this%axis_domain))
then
594 select type(domain => this%axis_domain)
596 if (use_collective_writes)
then
597 io_domain => domain%domain2
602 if (this%cart_name .eq.
"X")
then
604 position=this%domain_position)
605 elseif (this%cart_name .eq.
"Y")
then
607 position=this%domain_position)
621 if (
allocated(this%axis_domain))
then
622 axis_length = this%axis_domain%length(this%cart_name, this%domain_position, this%length)
638 if (
allocated(this%aux)) rslt = trim(this%aux) .ne.
""
649 if (
allocated(this%set_name)) rslt = trim(this%set_name) .ne.
""
657 integer,
optional,
intent(inout) :: x_or_y
661 select case (trim(this%cart_name))
670 if (
present(x_or_y)) x_or_y = diag_null
679 integer,
intent(out) :: dim_size
680 integer,
intent(out) :: layout
683 integer :: layout_xy(2)
685 select type (domain => this%axis_domain)
690 if (this%cart_name .eq.
"X")
then
692 layout = layout_xy(1)
693 else if (this%cart_name .eq.
"Y")
then
695 layout = layout_xy(2)
705 character(len=:),
allocatable :: rslt
715 character(len=:),
allocatable :: rslt
723 integer,
intent(in) :: axis_id
725 this%axis_id = axis_id
732 CHARACTER(len=*),
intent(in) :: edges_name
733 integer,
intent(in) :: edges_id
739 this%edges_name = edges_name
740 this%edges_id = edges_id
745 subroutine get_indices(this, compute_idx, corners_indices, starting_index, ending_index, need_to_define_axis)
747 integer,
intent(in) :: compute_idx(:)
748 class(*),
intent(in) :: corners_indices(:)
749 integer,
intent(out) :: starting_index
751 integer,
intent(out) :: ending_index
753 logical,
intent(out) :: need_to_define_axis
756 integer :: subregion_start
757 integer :: subregion_end
762 select type (corners_indices)
763 type is (
integer(kind=i4_kind))
764 subregion_start = minval(corners_indices)
765 subregion_end = maxval(corners_indices)
769 need_to_define_axis = .false.
770 starting_index = diag_null
771 ending_index = diag_null
774 if (compute_idx(1) < subregion_start .and. compute_idx(2) < subregion_start)
return
775 if (compute_idx(1) > subregion_end .and. compute_idx(2) > subregion_end)
return
777 need_to_define_axis = .true.
778 if (compute_idx(1) >= subregion_start .and. compute_idx(2) >= subregion_end)
then
780 starting_index = compute_idx(1)
781 ending_index = subregion_end
782 else if (compute_idx(1) >= subregion_start .and. compute_idx(2) <= subregion_end)
then
784 starting_index = compute_idx(1)
785 ending_index = compute_idx(2)
786 else if (compute_idx(1) <= subregion_start .and. compute_idx(2) <= subregion_end)
then
788 starting_index = subregion_start
789 ending_index = compute_idx(2)
790 else if (compute_idx(1) <= subregion_start .and. compute_idx(2) >= subregion_end)
then
792 starting_index = subregion_start
793 ending_index = subregion_end
796 if (this%domain_position .ne. center)
then
797 if (subregion_end - subregion_start + 1 .eq. 1)
then
799 if (ending_index .eq. compute_idx(2)) need_to_define_axis = .false.
801 if (ending_index - starting_index + 1 .eq. 1)
then
803 if (starting_index .eq. compute_idx(2) .or. ending_index .eq. compute_idx(1)) &
804 need_to_define_axis = .false.
814 integer,
intent(inout) :: compute_idx(:)
815 logical,
intent(out) :: need_to_define_axis
816 integer,
optional,
intent(in) :: tile_number
819 need_to_define_axis = .false.
820 compute_idx = diag_null
822 if (.not.
allocated(this%axis_domain))
then
824 if (this%cart_name .eq.
"X" .or. this%cart_name .eq.
"Y")
then
826 compute_idx(2) =
size(this%axis_data)
827 need_to_define_axis = .true.
832 select type(domain => this%axis_domain)
834 if (
present(tile_number))
then
838 need_to_define_axis = .false.
844 select case (this%cart_name)
847 & position=this%domain_position)
848 need_to_define_axis = .true.
851 & position=this%domain_position)
852 need_to_define_axis = .true.
860 subroutine fill_subaxis(this, starting_index, ending_index, axis_id, parent_id, parent_axis_name, compute_idx, &
861 global_idx, zbounds, nz_subaxis)
863 integer ,
intent(in) :: starting_index
864 integer ,
intent(in) :: ending_index
865 integer ,
intent(in) :: axis_id
866 integer ,
intent(in) :: parent_id
867 character(len=*) ,
intent(in) :: parent_axis_name
868 integer ,
intent(in) :: compute_idx(2)
870 integer,
optional,
intent(in) :: global_idx(2)
872 real(kind=r4_kind),
optional,
intent(in) :: zbounds(2)
873 integer,
optional,
intent(in) :: nz_subaxis
877 character(len=2) :: nsubaxis_char
880 if (
present(nz_subaxis)) nsubaxis = nz_subaxis
882 this%axis_id = axis_id
883 this%starting_index = starting_index
884 this%ending_index = ending_index
885 this%parent_axis_id = parent_id
886 write(nsubaxis_char,
'(i2.2)') nsubaxis
887 this%subaxis_name = trim(parent_axis_name)//
"_sub"//nsubaxis_char
888 this%compute_idx = compute_idx
890 if (
present(zbounds))
then
892 allocate(this%zbounds(2))
893 this%zbounds = zbounds
896 if (
present(global_idx))
then
898 allocate(this%global_idx(2))
899 this%global_idx = global_idx
910 res = this%ending_index - this%starting_index + 1
918 indx = this%starting_index
926 indx = this%ending_index
934 indx = this%compute_idx
953 function get_length(this, cart_axis, domain_position, global_length) &
956 character(len=*),
INTENT(IN) :: cart_axis
957 integer,
INTENT(IN) :: domain_position
958 integer,
INTENT(IN) :: global_length
964 if (trim(cart_axis) ==
"X")
call mpp_get_compute_domain(this%Domain2, xsize=length, position=domain_position)
965 if (trim(cart_axis) ==
"Y")
call mpp_get_compute_domain(this%Domain2, ysize=length, position=domain_position)
968 length = global_length
977 TYPE(
domain1d),
INTENT(in),
OPTIONAL :: Domain
978 TYPE(
domain2d),
INTENT(in),
OPTIONAL :: Domain2
979 TYPE(
domainug),
INTENT(in),
OPTIONAL :: DomainU
985 this%Domain2 = domain2
987 this%DomainUG = domainu
996 if (
allocated(axis_array))
call mpp_error(fatal,
"The diag_axis containers is already allocated")
1008 if (
allocated(axis_array))
deallocate(axis_array)
1021 logical,
intent(in),
optional :: is_regional
1023 character(len=:),
allocatable :: axis_name
1027 axis_name = this%axis_name
1028 if (
present(is_regional))
then
1029 if (is_regional)
then
1030 if (this%cart_name .eq.
"X" .or. this%cart_name .eq.
"Y") axis_name = axis_name//
"_sub01"
1034 axis_name = this%subaxis_name
1045 if (this%cart_name .eq.
"Z")
is_z_axis = .true.
1051 character(len=*),
intent(in) :: cart_name
1053 select case (cart_name)
1054 case (
"X",
"Y",
"Z",
"T",
"U",
"N")
1056 call mpp_error(fatal,
"diag_axit_init: Invalid cart_name: "//cart_name//&
1057 "The acceptable values are X, Y, Z, T, U, N.")
1063 integer,
INTENT(IN) :: domain_position
1065 select case (domain_position)
1066 case (center, north, east)
1068 call mpp_error(fatal,
"diag_axit_init: Invalid domain_positon. &
1069 &The acceptable values are NORTH, EAST, CENTER")
1075 integer,
INTENT(IN) :: direction
1077 select case(direction)
1080 call mpp_error(fatal,
"diag_axit_init: Invalid direction. &
1081 &The acceptable values are-1 0 1")
1088 integer,
INTENT(IN) :: axis_id(:)
1089 integer,
INTENT(OUT) :: domain_type
1091 character(len=*),
INTENT(IN) :: var_name
1099 do i = 1,
size(axis_id)
1101 select type (axis => diag_axis(j)%axis)
1104 if (domain_type .ne. axis%type_of_domain)
then
1109 if ((axis%type_of_domain .eq.
two_d_domain .and.
size(axis_id) > 1) &
1110 & .or. axis%type_of_domain .eq.
ug_domain)
then
1111 domain_type = axis%type_of_domain
1112 domain => axis%axis_domain
1115 call mpp_error(fatal,
"The variable:"//trim(var_name)//
" has axis that are not in the same domain")
1126 integer,
intent(inout) :: naxis
1128 integer,
intent(in) :: is_x_or_y
1130 logical,
intent(out) :: write_on_this_pe
1132 integer :: compute_idx(2)
1133 integer :: global_idx(2)
1134 integer :: starting_index
1135 integer :: ending_index
1137 call parent_axis%get_compute_domain(compute_idx, write_on_this_pe, tile_number=subregion%tile)
1138 if (.not. write_on_this_pe)
return
1142 call parent_axis%get_indices(compute_idx, subregion%corners(:,is_x_or_y), starting_index, ending_index, &
1145 if (.not. write_on_this_pe)
return
1147 select type(corners=> subregion%corners)
1148 type is (
integer(kind=i4_kind))
1149 global_idx(1) = minval(corners(:,is_x_or_y))
1150 global_idx(2) = maxval(corners(:,is_x_or_y))
1154 call define_new_axis(diag_axis, parent_axis, naxis, parent_axis%axis_id, &
1155 starting_index, ending_index, compute_idx, global_idx)
1162 integer,
INTENT(in) :: axis_ids(:)
1163 integer,
intent(inout) :: naxis
1165 logical,
intent(in) :: is_cube_sphere
1166 logical,
intent(out) :: write_on_this_pe
1171 integer :: lat_indices(2)
1172 integer :: lon_indices(2)
1173 integer :: compute_idx(2)
1174 integer :: starting_index(2)
1176 integer :: ending_index(2)
1177 logical :: need_to_define_axis(2)
1179 integer :: parent_axis_ids(2)
1180 logical :: is_x_y_axis
1181 integer :: compute_idx_2(2, 2)
1182 integer :: global_idx (2, 2)
1184 write_on_this_pe = .false.
1185 need_to_define_axis = .true.
1186 parent_axis_ids = diag_null
1191 select type (corners => subregion%corners)
1192 type is (real(kind=r4_kind))
1193 lon(1) = minval(corners(:,1))
1194 lon(2) = maxval(corners(:,1))
1195 lat(1) = minval(corners(:,2))
1196 lat(2) = maxval(corners(:,2))
1199 if_is_cube_sphere:
if (is_cube_sphere)
then
1201 call get_local_indices_cubesphere(lat(1), lat(2), lon(1), lon(2),&
1202 & lon_indices(1), lon_indices(2), lat_indices(1), lat_indices(2))
1203 loop_over_axis_ids:
do i = 1,
size(axis_ids)
1204 select_axis_type:
select type (parent_axis => diag_axis(axis_ids(i))%axis)
1207 call parent_axis%get_compute_domain(compute_idx, is_x_y_axis)
1210 if (.not. is_x_y_axis) cycle
1214 if (parent_axis%cart_name .eq.
"X")
then
1215 call parent_axis%get_indices(compute_idx, lon_indices, starting_index(1), ending_index(1), &
1216 need_to_define_axis(1))
1217 parent_axis_ids(1) = axis_ids(i)
1218 compute_idx_2(1,:) = compute_idx
1219 global_idx(1,:) = lon_indices
1220 else if (parent_axis%cart_name .eq.
"Y")
then
1221 call parent_axis%get_indices(compute_idx, lat_indices, starting_index(2), ending_index(2), &
1222 need_to_define_axis(2))
1223 parent_axis_ids(2) = axis_ids(i)
1224 compute_idx_2(2,:) = compute_idx
1225 global_idx(2,:) = lat_indices
1227 end select select_axis_type
1228 enddo loop_over_axis_ids
1229 else if_is_cube_sphere
1230 loop_over_axis_ids2:
do i = 1,
size(axis_ids)
1231 select type (parent_axis => diag_axis(axis_ids(i))%axis)
1234 call parent_axis%get_compute_domain(compute_idx, is_x_y_axis)
1237 if (.not. is_x_y_axis) cycle
1240 if (parent_axis%cart_name .eq.
"X")
then
1241 select type(adata=>parent_axis%axis_data)
1242 type is (real(kind=r8_kind))
1243 lon_indices(1) =
nearest_index(real(lon(1), kind=r8_kind), adata)
1244 lon_indices(2) =
nearest_index(real(lon(2), kind=r8_kind), adata)
1245 type is (real(kind=r4_kind))
1246 lon_indices(1) =
nearest_index(real(lon(1), kind=r4_kind), adata)
1247 lon_indices(2) =
nearest_index(real(lon(2), kind=r4_kind), adata)
1249 call parent_axis%get_indices(compute_idx, lon_indices, starting_index(1), ending_index(1), &
1250 need_to_define_axis(1))
1251 parent_axis_ids(1) = axis_ids(i)
1252 compute_idx_2(1,:) = compute_idx
1253 global_idx(1,:) = lon_indices
1254 else if (parent_axis%cart_name .eq.
"Y")
then
1255 select type(adata=>parent_axis%axis_data)
1256 type is (real(kind=r8_kind))
1257 lat_indices(1) =
nearest_index(real(lat(1), kind=r8_kind), adata)
1258 lat_indices(2) =
nearest_index(real(lat(2), kind=r8_kind), adata)
1259 type is (real(kind=r4_kind))
1260 lat_indices(1) =
nearest_index(real(lat(1), kind=r4_kind), adata)
1261 lat_indices(2) =
nearest_index(real(lat(2), kind=r4_kind), adata)
1263 call parent_axis%get_indices(compute_idx, lat_indices, starting_index(2), ending_index(2), &
1264 need_to_define_axis(2))
1265 parent_axis_ids(2) = axis_ids(i)
1266 compute_idx_2(2,:) = compute_idx
1267 global_idx(2,:) = lat_indices
1270 enddo loop_over_axis_ids2
1271 endif if_is_cube_sphere
1274 if (any(.not. need_to_define_axis ))
return
1277 write_on_this_pe = .true.
1279 do i = 1,
size(parent_axis_ids)
1280 if (parent_axis_ids(i) .eq. diag_null) cycle
1281 select type (parent_axis => diag_axis(parent_axis_ids(i))%axis)
1283 call define_new_axis(diag_axis, parent_axis, naxis, parent_axis_ids(i), &
1284 starting_index(i), ending_index(i), compute_idx_2(i,:), global_idx(i,:))
1292 starting_index, ending_index, compute_idx, global_idx, new_axis_id, zbounds, &
1297 integer,
intent(inout) :: naxis
1299 integer,
intent(in) :: parent_id
1300 integer,
intent(in) :: starting_index
1301 integer,
intent(in) :: ending_index
1302 integer,
intent(in) :: compute_idx(2)
1304 integer,
optional,
intent(in) :: global_idx(2)
1306 integer,
optional,
intent(out) :: new_axis_id
1307 real(kind=r4_kind),
optional,
intent(in) :: zbounds(2)
1308 integer,
optional,
intent(in) :: nz_subaxis
1314 parent_axis%nsubaxis = parent_axis%nsubaxis + 1
1315 parent_axis%subaxis(parent_axis%nsubaxis) = naxis
1319 diag_axis(naxis)%axis%axis_id = naxis
1320 if (
present(new_axis_id)) new_axis_id = naxis
1322 select type (sub_axis => diag_axis(naxis)%axis)
1324 call sub_axis%fill_subaxis(starting_index, ending_index, naxis, parent_id, &
1325 parent_axis%axis_name, compute_idx, global_idx=global_idx, zbounds=zbounds, nz_subaxis=nz_subaxis)
1332 result(parent_axis_id)
1335 integer :: parent_axis_id
1339 parent_axis_id = diag_null
1341 parent_axis_id = this%parent_axis_id
1343 parent_axis_id = diag_null
1354 integer :: sub_axis_id
1356 sub_axis_id = this%axis_id
1359 if (this%cart_name .ne.
"Z") sub_axis_id = this%subaxis(this%nsubaxis)
1368 class(*),
intent(in) :: compress_att(:)
1369 character(len=120) :: axis_names(2)
1373 select type (compress_att)
1374 type is (
character(len=*))
1375 read(compress_att(1),*, iostat=ios) axis_names
1376 if (ios .ne. 0) axis_names =
""
1387 character(len=*),
intent(in) :: axis_name
1388 integer,
intent(in) :: naxis
1389 character(len=*),
intent(in) :: set_name
1396 select type(axis => diag_axis(i)%axis)
1398 if (trim(axis%axis_name) .eq. trim(axis_name))
then
1399 if (trim(axis%set_name) .eq. trim(set_name))
then
1412 result(n_diurnal_samples)
1415 integer :: n_diurnal_samples
1417 n_diurnal_samples = this%ndiurnal_samples
1423 class(fmsnetcdffile_t),
intent(inout) :: fms2io_fileobj
1425 call register_axis(fms2io_fileobj, this%axis_name,
size(this%diurnal_data))
1427 call register_variable_attribute(fms2io_fileobj, this%axis_name,
"units", &
1428 &trim(this%units), str_len=len_trim(this%units))
1429 call register_variable_attribute(fms2io_fileobj, this%axis_name,
"long_name", &
1430 &trim(this%long_name), str_len=len_trim(this%long_name))
1431 if (this%edges_id .ne. diag_null) &
1432 call register_variable_attribute(fms2io_fileobj, this%axis_name,
"edges", &
1433 &trim(this%edges_name), str_len=len_trim(this%edges_name))
1438 real(kind=r4_kind),
intent(in) :: zbounds(2)
1439 integer,
intent(inout) :: var_axis_ids(:)
1441 integer,
intent(inout) :: naxis
1443 integer,
intent(inout) :: file_axis_id(:)
1444 integer,
intent(inout) :: nfile_axis
1446 integer,
intent(inout) :: nz_subaxis
1449 class(*),
pointer :: zaxis_data(:)
1450 integer :: subaxis_indices(2)
1453 integer :: subaxis_id
1454 logical :: axis_found
1457 axis_found = .false.
1458 do i = 1, nfile_axis
1459 select type (axis => diag_axis(file_axis_id(i))%axis)
1461 if (axis%zbounds(1) .eq. zbounds(1) .and. axis%zbounds(2) .eq. zbounds(2))
then
1463 subaxis_id = file_axis_id(i)
1470 do i = 1,
size(var_axis_ids)
1471 select type (parent_axis => diag_axis(var_axis_ids(i))%axis)
1473 if (parent_axis%cart_name .eq.
"Z")
then
1475 if (axis_found)
then
1476 var_axis_ids(i) = subaxis_id
1479 zaxis_data => parent_axis%axis_data
1481 select type(zaxis_data)
1482 type is (real(kind=r4_kind))
1484 subaxis_indices(1) =
nearest_index(real(zbounds(1)), real(zaxis_data))
1485 subaxis_indices(2) =
nearest_index(real(zbounds(2)), real(zaxis_data))
1486 type is (real(kind=r8_kind))
1487 subaxis_indices(1) =
nearest_index(real(zbounds(1)), real(zaxis_data))
1488 subaxis_indices(2) =
nearest_index(real(zbounds(2)), real(zaxis_data))
1491 nz_subaxis = nz_subaxis + 1
1492 call define_new_axis(diag_axis, parent_axis, naxis, parent_axis%axis_id, &
1493 &subaxis_indices(1), subaxis_indices(2), (/lbound(zaxis_data,1), ubound(zaxis_data,1)/), &
1494 &new_axis_id=subaxis_id, zbounds=zbounds, nz_subaxis=nz_subaxis)
1495 var_axis_ids(i) = subaxis_id
1507 integer,
intent(in) :: axis_id
1508 integer,
intent(in) :: parent_axis_id
1514 select type(axis => diag_axis(axis_id)%axis)
1516 if (axis%parent_axis_id .eq. parent_axis_id) rslt = .true.
1521 end module fms_diag_axis_object_mod
integer, parameter direction_down
The axis points down if positive.
integer function get_base_minute()
gets the module variable base_minute
integer function get_base_year()
gets the module variable base_year
integer function get_base_hour()
gets the module variable base_hour
integer, parameter no_domain
Use the FmsNetcdfFile_t fileobj.
integer max_axis_attributes
Maximum number of user definable attributes per axis.
character(len=6) pack_size_str
Pack size as a string to be used in fms2_io register call set to "double" or "float".
integer max_axes
Maximum number of independent axes.
integer, parameter is_x_axis
integer indicating that it is a x axis
integer, parameter is_y_axis
integer indicating that it is a y axis
integer function get_base_day()
gets the module variable base_day
integer, parameter ug_domain
Use the FmsNetcdfUnstructuredDomainFile_t fileobj.
integer, parameter direction_up
The axis points up if positive.
integer function get_base_month()
gets the module variable base_month
integer function get_base_second()
gets the module variable base_second
integer, parameter two_d_domain
Use the FmsNetcdfDomainFile_t fileobj.
Attribute type for diagnostic fields.
Type to hold the attributes of the field/axis/file.
subroutine, public get_local_indexes(latStart, latEnd, lonStart, lonEnd, istart, iend, jstart, jend)
Find the local start and local end indexes on the local PE for regional output.
Add a dimension to a given file.
Defines a new field within the given file Example usage:
Write data to a defined field within a file Example usage:
subroutine, public define_new_subaxis_latlon(diag_axis, axis_ids, naxis, subRegion, is_cube_sphere, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by lat lon.
subroutine get_dim_size_layout(this, dim_size, layout)
pure character(len=:) function, allocatable get_axis_name(this, is_regional)
integer function get_length(this, cart_axis, domain_position, global_length)
Get the length of a 2D domain.
integer function get_ending_index(this)
Accesses its member ending_index.
logical function, public fms_diag_axis_object_end(axis_array)
subroutine check_if_valid_domain_position(domain_position)
Check if a domain_position is valid and crashes if it isn't.
logical function is_x_or_y_axis(this, x_or_y)
Determine if an axis object is an x or y axis.
pure logical function is_unstructured_grid(this)
subroutine set_axis_id(this, axis_id)
Set the axis_id.
pure logical function is_z_axis(this)
integer function axis_length(this)
Get the axis length of a subaxis.
pure integer function get_diurnal_axis_samples(this)
integer function get_ntiles(this)
Get the ntiles in a domain.
subroutine, public define_new_subaxis_index(parent_axis, subRegion, diag_axis, naxis, is_x_or_y, write_on_this_pe)
Fill in the subaxis object for a subRegion defined by index.
pure integer function, public get_axis_id_from_name(axis_name, diag_axis, naxis, set_name)
subroutine get_compute_domain(this, compute_idx, need_to_define_axis, tile_number)
subroutine, public define_diurnal_axis(diag_axis, naxis, n_diurnal_samples, is_edges)
Defined a new diurnal axis.
pure integer function, dimension(2) get_structured_axis(this)
pure logical function has_aux(this)
Determine if an axis object has an auxiliary name.
subroutine add_structured_axis_ids(this, axis_ids)
subroutine check_if_valid_cart_name(cart_name)
Check if a cart_name is valid and crashes if it isn't.
subroutine write_axis_metadata(this, fms2io_fileobj, edges_in_file, parent_axis)
Write the axis meta data to an open fileobj.
subroutine fill_subaxis(this, starting_index, ending_index, axis_id, parent_id, parent_axis_name, compute_idx, global_idx, zbounds, nz_subaxis)
Fills in the information needed to define a subaxis.
subroutine get_global_io_domain(this, global_io_index, use_collective_writes)
Get the starting and ending indices of the global io domain of the axis.
logical function, public fms_diag_axis_object_init(axis_array)
pure integer function get_subaxes_id(this)
integer function get_starting_index(this)
Accesses its member starting_index.
subroutine get_indices(this, compute_idx, corners_indices, starting_index, ending_index, need_to_define_axis)
Determine if the subRegion is in the current PE. If it is, determine the starting and ending indices ...
pure integer function get_parent_axis_id(this)
subroutine check_if_valid_direction(direction)
Check if a direction is valid and crashes if it isn't.
subroutine, public create_new_z_subaxis(zbounds, var_axis_ids, diag_axis, naxis, file_axis_id, nfile_axis, nz_subaxis)
Creates a new z subaxis to use.
subroutine register_diag_axis_obj(this, axis_name, axis_data, units, cart_name, long_name, direction, set_name, Domain, Domain2, DomainU, aux, req, tile_count, domain_position, axis_length)
Initialize the axis.
subroutine write_axis_data(this, fms2io_fileobj, parent_axis)
Write the axis data to an open fms2io_fileobj.
subroutine, public define_new_axis(diag_axis, parent_axis, naxis, parent_id, starting_index, ending_index, compute_idx, global_idx, new_axis_id, zbounds, nz_subaxis)
Creates a new subaxis and fills it will all the information it needs.
integer function get_axis_length(this)
Get the length of the axis.
pure character(len=:) function, allocatable get_set_name(this)
Get the set name of an axis object.
pure character(len=120) function, dimension(2), public parse_compress_att(compress_att)
subroutine set_edges(this, edges_name, edges_id)
Set the name and ids of the edges.
logical function, public is_parent_axis(axis_id, parent_axis_id, diag_axis)
Determine if the diag_axis(parent_axis_id) is the parent of diag_axis(axis_id)
subroutine add_axis_attribute(this, att_name, att_value)
Add an attribute to an axis.
subroutine write_diurnal_metadata(this, fms2io_fileobj)
pure logical function has_set_name(this)
Determine if an axis object has a set_name.
pure integer function get_edges_id(this)
integer function, dimension(2) get_compute_indices(this)
Accesses its member compute_indices.
subroutine, public get_domain_and_domain_type(diag_axis, axis_id, domain_type, domain, var_name)
Loop through a variable's axis_id to determine and return the domain type and domain to use.
subroutine set_axis_domain(this, Domain, Domain2, DomainU)
Set the axis domain.
pure character(len=:) function, allocatable get_aux(this)
Get the auxiliary name of an axis object.
integer function, dimension(size(domain%tile_id(:))) mpp_get_tile_id(domain)
Returns the tile_id on current pe.
integer function mpp_get_ntile_count(domain)
Returns number of tiles in mosaic.
type(domain2d) function, pointer mpp_get_io_domain(domain)
Set user stack size.
These routines retrieve the axis specifications associated with the compute domains....
These routines retrieve the axis specifications associated with the global domains....
Retrieve layout associated with a domain decomposition The 1D version of this call returns the number...
One dimensional domain used to manage shared data access between pes.
The domain2D type contains all the necessary information to define the global, compute and data domai...
Domain information for managing data on unstructured grids.
integer function stdout()
This function returns the current standard fortran unit numbers for output.
integer function mpp_pe()
Returns processor ID.
Type to hold the 1d domain.
Type to hold the 2d domain.
Type to hold the domain info for an axis This type was created to avoid having to send in "Domain",...
Type to hold the unstructured domain.
Type to hold the diagnostic axis description.
Type to hold the diag_axis (either subaxis or a full axis)
Type to hold the diurnal axis.
Type to hold the diagnostic axis description.
Type to hold the subaxis.
type to hold the sub region information about a file