diff --git a/config_src/infra/FMS1/MOM_interp_infra.F90 b/config_src/infra/FMS1/MOM_interp_infra.F90 index ca5b2b8516..170573f7ec 100644 --- a/config_src/infra/FMS1/MOM_interp_infra.F90 +++ b/config_src/infra/FMS1/MOM_interp_infra.F90 @@ -4,9 +4,9 @@ module MOM_interp_infra ! This file is part of MOM6. See LICENSE.md for the license. use MOM_domain_infra, only : MOM_domain_type, domain2d -use MOM_io_infra, only : axistype use MOM_time_manager, only : time_type use horiz_interp_mod, only : horiz_interp_new, horiz_interp, horiz_interp_init, horiz_interp_type +use mpp_io_mod, only : axistype, mpp_get_axis_data use time_interp_external_mod, only : time_interp_external use time_interp_external_mod, only : init_external_field, time_interp_external_init use time_interp_external_mod, only : get_external_field_size @@ -16,7 +16,7 @@ module MOM_interp_infra public :: horiz_interp_type, horiz_interp_init public :: time_interp_extern, init_extern_field, time_interp_external_init -public :: get_external_field_info +public :: get_external_field_info, axistype, get_axis_data public :: run_horiz_interp, build_horiz_interp_weights !> Read a field based on model time, and rotate to the model domain. @@ -114,6 +114,15 @@ subroutine build_horiz_interp_weights_2d_to_2d(Interp, lon_in, lat_in, lon_out, end subroutine build_horiz_interp_weights_2d_to_2d +!> Extracts and returns the axis data stored in an axistype. +subroutine get_axis_data( axis, dat ) + type(axistype), intent(in) :: axis !< An axis type + real, dimension(:), intent(out) :: dat !< The data in the axis variable + + call mpp_get_axis_data( axis, dat ) +end subroutine get_axis_data + + !> get size of an external field from field index function get_extern_field_size(index) diff --git a/config_src/infra/FMS1/MOM_io_infra.F90 b/config_src/infra/FMS1/MOM_io_infra.F90 index 3ea201235a..0d4cc0deb5 100644 --- a/config_src/infra/FMS1/MOM_io_infra.F90 +++ b/config_src/infra/FMS1/MOM_io_infra.F90 @@ -312,13 +312,12 @@ subroutine get_filename_suffix(suffix) end subroutine get_filename_suffix -!> Get information about the number of dimensions, variables, global attributes and time levels +!> Get information about the number of dimensions, variables and time levels !! in the file associated with an open file unit -subroutine get_file_info(IO_handle, ndim, nvar, natt, ntime) +subroutine get_file_info(IO_handle, ndim, nvar, ntime) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for I/O integer, optional, intent(out) :: ndim !< The number of dimensions in the file integer, optional, intent(out) :: nvar !< The number of variables in the file - integer, optional, intent(out) :: natt !< The number of global attributes in the file integer, optional, intent(out) :: ntime !< The number of time levels in the file ! Local variables @@ -328,7 +327,6 @@ subroutine get_file_info(IO_handle, ndim, nvar, natt, ntime) if (present(ndim)) ndim = ndims if (present(nvar)) nvar = nvars - if (present(natt)) natt = natts if (present(ntime)) ntime = ntimes end subroutine get_file_info @@ -683,7 +681,7 @@ subroutine write_field_4d(IO_handle, field_md, MOM_domain, field, tstamp, tile_c type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value @@ -697,7 +695,7 @@ subroutine write_field_3d(IO_handle, field_md, MOM_domain, field, tstamp, tile_c type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value @@ -711,7 +709,7 @@ subroutine write_field_2d(IO_handle, field_md, MOM_domain, field, tstamp, tile_c type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value @@ -724,7 +722,7 @@ subroutine write_field_1d(IO_handle, field_md, field, tstamp) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, dimension(:), intent(in) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field call mpp_write(IO_handle%unit, field_md, field, tstamp=tstamp) end subroutine write_field_1d @@ -734,7 +732,7 @@ subroutine write_field_0d(IO_handle, field_md, field, tstamp) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, intent(in) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field call mpp_write(IO_handle%unit, field_md, field, tstamp=tstamp) end subroutine write_field_0d @@ -750,7 +748,8 @@ end subroutine MOM_write_axis !> Store information about an axis in a previously defined axistype and write this !! information to the file indicated by unit. -subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian, sense, domain, data, calendar) +subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian, sense, domain, & + data, edge_axis, calendar) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(axistype), intent(inout) :: axis !< The axistype where this information is stored. character(len=*), intent(in) :: name !< The name in the file of this axis @@ -763,6 +762,7 @@ subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian !! -1 if they increase downward. type(domain1D), optional, intent(in) :: domain !< The domain decomposion for this axis real, dimension(:), optional, intent(in) :: data !< The coordinate values of the points on this axis + logical, optional, intent(in) :: edge_axis !< If true, this axis marks an edge of the tracer cells character(len=*), optional, intent(in) :: calendar !< The name of the calendar used with a time axis call mpp_write_meta(IO_handle%unit, axis, name, units, longname, cartesian=cartesian, sense=sense, & @@ -772,19 +772,13 @@ end subroutine write_metadata_axis !> Store information about an output variable in a previously defined fieldtype and write this !! information to the file indicated by unit. subroutine write_metadata_field(IO_handle, field, axes, name, units, longname, & - min, max, fill, scale, add, pack, standard_name, checksum) + pack, standard_name, checksum) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(inout) :: field !< The fieldtype where this information is stored type(axistype), dimension(:), intent(in) :: axes !< Handles for the axis used for this variable character(len=*), intent(in) :: name !< The name in the file of this variable character(len=*), intent(in) :: units !< The units of this variable character(len=*), intent(in) :: longname !< The long description of this variable - real, optional, intent(in) :: min !< The minimum valid value for this variable - real, optional, intent(in) :: max !< The maximum valid value for this variable - real, optional, intent(in) :: fill !< Missing data fill value - real, optional, intent(in) :: scale !< An multiplicative factor by which to scale - !! the variable before output - real, optional, intent(in) :: add !< An offset to add to the variable before output integer, optional, intent(in) :: pack !< A precision reduction factor with which the !! variable. The default, 1, has no reduction, !! but 2 is not uncommon. @@ -793,8 +787,9 @@ subroutine write_metadata_field(IO_handle, field, axes, name, units, longname, & optional, intent(in) :: checksum !< Checksum values that can be used to verify reads. - call mpp_write_meta(IO_handle%unit, field, axes, name, units, longname, min=min, max=max, & - fill=fill, scale=scale, add=add, pack=pack, standard_name=standard_name, checksum=checksum) + call mpp_write_meta(IO_handle%unit, field, axes, name, units, longname, & + pack=pack, standard_name=standard_name, checksum=checksum) + ! unused opt. args: min=min, max=max, fill=fill, scale=scale, add=add, & end subroutine write_metadata_field diff --git a/config_src/infra/FMS2/MOM_axis.F90 b/config_src/infra/FMS2/MOM_axis.F90 deleted file mode 100644 index b5d2b3ed88..0000000000 --- a/config_src/infra/FMS2/MOM_axis.F90 +++ /dev/null @@ -1,616 +0,0 @@ -!> This module contains routines that define and register axes to files -module MOM_axis - -! This file is part of MOM6. See LICENSE.md for the license. -use MOM_domains, only : MOM_domain_type -use MOM_error_handler, only : MOM_error, NOTE, FATAL, WARNING -use MOM_grid, only : ocean_grid_type -use MOM_dyn_horgrid, only : dyn_horgrid_type -use MOM_string_functions, only : lowercase -use MOM_verticalGrid, only : verticalGrid_type -use fms2_io_mod, only : is_dimension_registered, register_axis, is_dimension_unlimited -use fms2_io_mod, only : FmsNetcdfDomainFile_t, FmsNetcdfFile_t, unlimited -use fms2_io_mod, only : get_variable_size, get_variable_num_dimensions, check_if_open -use fms2_io_mod, only : fms2_open_file=>open_file, fms2_close_file=>close_file -use fms2_io_mod, only : get_variable_dimension_names, read_data, get_unlimited_dimension_name -use fms2_io_mod, only : get_dimension_size -use mpp_domains_mod, only : domain2d, CENTER, CORNER, NORTH_FACE=>NORTH, EAST_FACE=>EAST -use mpp_domains_mod, only : mpp_get_compute_domain -use netcdf -implicit none ; private - -public MOM_register_diagnostic_axis, get_var_dimension_metadata, get_time_units -public MOM_get_diagnostic_axis_data, MOM_register_variable_axes, get_time_index -public convert_checksum_to_string -!> A type for making arrays of pointers to real 1-d arrays -type p1d - real, dimension(:), pointer :: p => NULL() !< A pointer to a 1d array -end type p1d - -!> A structure with information about a single axis variable -type axis_atts - character(len=64) :: name !< Names of the axis - character(len=48) :: units !< Physical dimensions of the axis - character(len=240) :: longname !< Long name of the axis - character(len=8) :: positive !< Positive-definite direction: up, down, east, west, north, south - integer :: horgrid_position !< Horizontal grid position - logical :: is_domain_decomposed !< if .true. the axis data are domain-decomposed - !! and need to be indexed by the compute domain - !! before passing to write_data -end type axis_atts - -!> Type for describing an axis variable (e.g., lath, lonh, Time) -type, public :: axis_data_type - !> An array of descriptions of the registered axes - type(axis_atts), pointer :: axis(:) => NULL() !< structure with axis attributes - type(p1d), pointer :: data(:) => NULL() !< pointer to the axis data -end type axis_data_type - -!> interface for registering axes associated with a variable to a netCDF file object -interface MOM_register_variable_axes - module procedure MOM_register_variable_axes_subdomain - module procedure MOM_register_variable_axes_full -end interface MOM_register_variable_axes - -contains - -!> register a MOM diagnostic axis to a domain-decomposed file -subroutine MOM_register_diagnostic_axis(fileObj, axisName, axisLength) - type(FmsNetcdfDomainFile_t), intent(inout) :: fileObj !< netCDF file object returned by call to open_file - character(len=*), intent(in) :: axisName !< name of the axis to register to file - integer, intent(in), optional :: axisLength !< length of axis/dimension ;only needed for Layer, Interface, Time, - !! Period - select case (trim(lowercase(axisName))) - case ('latq'); call register_axis(fileObj,'latq','y', domain_position=NORTH_FACE) - case ('lath'); call register_axis(fileObj,'lath','y', domain_position=CENTER) - case ('lonq'); call register_axis(fileObj,'lonq','x', domain_position=EAST_FACE) - case ('lonh'); call register_axis(fileObj,'lonh','x', domain_position=CENTER) - case default - if (.not. present(axisLength)) call MOM_error(FATAL,"MOM_io:register_diagnostic_axis: "//& - "An axis_length argument is required to register the axis "//trim(axisName)) - call register_axis(fileObj, trim(axisName), axisLength) - end select -end subroutine MOM_register_diagnostic_axis - - -!> Get the horizontal grid, vertical grid, and/or time dimension names and lengths -!! for a single variable from the hor_grid, t_grid, and z_grid values returned by a prior call to query_vardesc -subroutine get_var_dimension_metadata(hor_grid, z_grid, t_grid_in, & - dim_names, dim_lengths, num_dims, G, dG, GV) - - character(len=*), intent(in) :: hor_grid !< horizontal grid - character(len=*), intent(in) :: z_grid !< vertical grid - character(len=*), intent(in) :: t_grid_in !< time grid - character(len=*), dimension(:), intent(inout) :: dim_names !< array of dimension names - integer, dimension(:), intent(inout) :: dim_lengths !< array of dimension sizes - integer, intent(inout) :: num_dims !< number of axes to register in the restart file - type(ocean_grid_type), optional, intent(in) :: G !< The ocean's grid structure - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure - - ! local - logical :: use_lath - logical :: use_lonh - logical :: use_latq - logical :: use_lonq - character(len=8) :: t_grid - character(len=8) :: t_grid_read - integer :: isg, ieg, jsg, jeg, IsgB, IegB, JsgB, JegB - !integer :: npes - real, pointer, dimension(:) :: gridLatT => NULL(), & ! The latitude or longitude of T or B points for - gridLatB => NULL(), & ! the purpose of labeling the output axes. - gridLonT => NULL(), & - gridLonB => NULL() - type(MOM_domain_type), pointer :: domain => NULL() ! Domain used to get the pe count - - use_lath = .false. - use_lonh = .false. - use_latq = .false. - use_lonq = .false. - - ! set the ocean grid coordinates - - if (present(G)) then - gridLatT => G%gridLatT ; gridLatB => G%gridLatB - gridLonT => G%gridLonT ; gridLonB => G%gridLonB - isg = G%isg ; ieg = G%ieg ; jsg = G%jsg ; jeg = G%jeg - IsgB = G%IsgB ; IegB = G%IegB ; JsgB = G%JsgB ; JegB = G%JegB - - call get_horizontal_grid_logic(hor_grid, use_lath, use_lonh, use_latq, use_lonq) - elseif (present(dG)) then - gridLatT => dG%gridLatT ; gridLatB => dG%gridLatB - gridLonT => dG%gridLonT ; gridLonB => dG%gridLonB - isg = dG%isg ; ieg = dG%ieg ; jsg = dG%jsg ; jeg = dG%jeg - IsgB = dG%IsgB ; IegB = dG%IegB ; JsgB = dG%JsgB ; JegB = dG%JegB - - call get_horizontal_grid_logic(hor_grid, use_lath, use_lonh, use_latq, use_lonq) - endif - - ! add longitude name to dimension name array - if (use_lonh) then - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("lonh")) = "lonh" - dim_lengths(num_dims) = size(gridLonT(isg:ieg)) - elseif (use_lonq) then - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("lonq")) = "lonq" - dim_lengths(num_dims) = size(gridLonB(IsgB:IegB)) - endif - ! add latitude name to dimension name array - if (use_lath) then - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("lath")) = "lath" - dim_lengths(num_dims) = size(gridLatT(jsg:jeg)) - elseif (use_latq) then - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("latq")) = "latq" - dim_lengths(num_dims) = size(gridLatB(JsgB:JegB)) - endif - - if (present(GV)) then - ! vertical grid - select case (trim(z_grid)) - case ('L') - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("Layer")) = "Layer" - dim_lengths(num_dims) = GV%ke - case ('i') - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("Interface")) = "Interface" - dim_lengths(num_dims) = GV%ke+1 - case ('1') ! Do nothing. - case default - call MOM_error(FATAL, "MOM_io: get_var_dimension_features: "//& - " has an unrecognized z_grid argument"//trim(z_grid)) - end select - endif - ! time - t_grid = adjustl(t_grid_in) - select case (t_grid(1:1)) - case ('s', 'a', 'm') - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("Time")) = "Time" - dim_lengths(num_dims) = unlimited - case ('p') - if (len_trim(t_grid(2:8)) <= 0) then - call MOM_error(FATAL,"MOM_io:get_var_dimension_features: "//& - "No periodic axis length was specified in "//trim(t_grid)) - endif - num_dims = num_dims+1 - dim_names(num_dims) = "" - dim_names(num_dims)(1:len_trim("Period")) = "Period" - dim_lengths(num_dims) = unlimited - case ('1') ! Do nothing. - case default - call MOM_error(WARNING, "MOM_io: get_var_dimension_metadata: "//& - "Unrecognized t_grid "//trim(t_grid)) - end select -end subroutine get_var_dimension_metadata - - -!> Populate the axis_data structure with axis data and attributes for diagnostic and restart files -subroutine MOM_get_diagnostic_axis_data(axis_data_CS, axis_name, axis_number, G, dG, GV, time_val, time_units) - - type(axis_data_type), intent(inout) :: axis_data_CS !< structure containing the axis data and metadata - character(len=*), intent(in) :: axis_name !< name of the axis - integer, intent(in) :: axis_number !< positional value (wrt to file) of the axis to register - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the file uses any - !! horizontal grid axes. - type(verticalGrid_type), target, optional, intent(in) :: GV !< ocean vertical grid structure - real,dimension(:), target, optional, intent(in) :: time_val !< time value - character(len=*), optional,intent(in) :: time_units!< units for non-periodic time axis - ! local - character(len=40) :: x_axis_units='', y_axis_units='' - integer :: isg, ieg, jsg, jeg, IsgB, IegB, JsgB, JegB - real, pointer, dimension(:) :: gridLatT => NULL(), & ! The latitude or longitude of T or B points for - gridLatB => NULL(), & ! the purpose of labeling the output axes. - gridLonT => NULL(), & - gridLonB => NULL() - - ! initialize axis_data_CS elements - axis_data_CS%axis(axis_number)%name = '' - axis_data_CS%axis(axis_number)%longname = '' - axis_data_CS%axis(axis_number)%units = '' - axis_data_CS%axis(axis_number)%horgrid_position = 0 - axis_data_CS%axis(axis_number)%is_domain_decomposed = .false. - axis_data_CS%axis(axis_number)%positive = '' - axis_data_CS%data(axis_number)%p => NULL() - - ! set the ocean grid coordinates and metadata - if (present(G)) then - gridLatT => G%gridLatT ; gridLatB => G%gridLatB - gridLonT => G%gridLonT ; gridLonB => G%gridLonB - x_axis_units = G%x_axis_units ; y_axis_units = G%y_axis_units - isg = G%isg ; ieg = G%ieg ; jsg = G%jsg ; jeg = G%jeg - IsgB = G%IsgB ; IegB = G%IegB ; JsgB = G%JsgB ; JegB = G%JegB - elseif (present(dG)) then - gridLatT => dG%gridLatT ; gridLatB => dG%gridLatB - gridLonT => dG%gridLonT ; gridLonB => dG%gridLonB - x_axis_units = dG%x_axis_units ; y_axis_units = dG%y_axis_units - isg = dG%isg ; ieg = dG%ieg ; jsg = dG%jsg ; jeg = dG%jeg - IsgB = dG%IsgB ; IegB = dG%IegB ; JsgB = dG%JsgB ; JegB = dG%JegB - endif - - select case(trim(lowercase(axis_name))) - case('lath') - if (associated(gridLatT)) & - axis_data_CS%data(axis_number)%p=>gridLatT(jsg:jeg) - - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Latitude' - axis_data_CS%axis(axis_number)%units = y_axis_units - axis_data_CS%axis(axis_number)%horgrid_position = CENTER - axis_data_CS%axis(axis_number)%is_domain_decomposed = .true. - case('lonh') - if (associated(gridLonT)) & - axis_data_CS%data(axis_number)%p=>gridLonT(isg:ieg) - - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%horgrid_position = CENTER - axis_data_CS%axis(axis_number)%longname = 'Longitude' - axis_data_CS%axis(axis_number)%units = x_axis_units - axis_data_CS%axis(axis_number)%is_domain_decomposed = .true. - case('latq') - if (associated(gridLatB)) & - axis_data_CS%data(axis_number)%p=>gridLatB(JsgB:JegB) - - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Latitude' - axis_data_CS%axis(axis_number)%units = y_axis_units - axis_data_CS%axis(axis_number)%horgrid_position = NORTH_FACE - axis_data_CS%axis(axis_number)%is_domain_decomposed = .true. - case('lonq') - if (associated(gridLonB)) & - axis_data_CS%data(axis_number)%p=>gridLonB(IsgB:IegB) - - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Longitude' - axis_data_CS%axis(axis_number)%units = x_axis_units - axis_data_CS%axis(axis_number)%horgrid_position = EAST_FACE - axis_data_CS%axis(axis_number)%is_domain_decomposed = .true. - case('layer') - if (present(GV)) then - axis_data_CS%data(axis_number)%p=>GV%sLayer(1:GV%ke) - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Layer pseudo-depth, -z*' - axis_data_CS%axis(axis_number)%units = GV%zAxisUnits - axis_data_CS%axis(axis_number)%positive = 'up' - endif - case('interface') - if (present(GV)) then - axis_data_CS%data(axis_number)%p=>GV%sInterface(1:GV%ke+1) - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Interface pseudo-depth, -z*' - axis_data_CS%axis(axis_number)%units = GV%zAxisUnits - axis_data_CS%axis(axis_number)%positive = 'up' - endif - case('time') - if (.not.(present(time_val))) & - call MOM_error(FATAL, "MOM_io::get_diagnostic_axis_data: requires time_val"//& - " and time_units arguments for "//trim(axis_name)) - - axis_data_CS%data(axis_number)%p=>time_val - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Time' - - if (present(time_units)) then - axis_data_CS%axis(axis_number)%units = time_units - else - axis_data_CS%axis(axis_number)%units = 'days' - endif - case('period') - if (.not.(present(time_val))) & - call MOM_error(FATAL, "MOM_axis::get_diagnostic_axis_data: requires a time_val argument "// & - "for "//trim(axis_name)) - axis_data_CS%data(axis_number)%p=>time_val - axis_data_CS%axis(axis_number)%name = trim(axis_name) - axis_data_CS%axis(axis_number)%longname = 'Periods for cyclical variables' - case default - call MOM_error(WARNING, "MOM_axis::get_diagnostic_axis_data:"//trim(axis_name)//" is an unrecognized axis") - end select - -end subroutine MOM_get_diagnostic_axis_data - - -!> set the logical variables that determine which diagnositic axes to use -subroutine get_horizontal_grid_logic(grid_string_id, use_lath, use_lonh, use_latq, use_lonq) - character(len=*), intent(in) :: grid_string_id !< horizontal grid string - logical, intent(out) :: use_lath !< if .true., y-axis is oriented in CENTER position - logical, intent(out) :: use_lonh !< if .true., x-axis is oriented in CENTER position - logical, intent(out) :: use_latq !< if .true., y-axis is oriented in NORTH_FACE position - logical, intent(out) :: use_lonq !< if .true., x-axis is oriented in EAST_FACE position - - use_lath = .false. - use_lonh = .false. - use_latq = .false. - use_lonq = .false. - select case (trim(grid_string_id)) - case ('h') ; use_lath = .true. ; use_lonh = .true. ! x=CENTER, y=CENTER - case ('q') ; use_latq = .true. ; use_lonq = .true. ! x=EAST_FACE, y=NORTH_FACE - case ('u') ; use_lath = .true. ; use_lonq = .true. ! x=EAST_FACE, y=CENTER - case ('v') ; use_latq = .true. ; use_lonh = .true. ! x=CENTER, y=NORTH_FACE - case ('T') ; use_lath = .true. ; use_lonh = .true. ! x=CENTER, y=CENTER - case ('Bu') ; use_latq = .true. ; use_lonq = .true. ! x=EAST_FACE, y=NORTH_FACE - case ('Cu') ; use_lath = .true. ; use_lonq = .true. ! x=EAST_FACE, y=CENTER - case ('Cv') ; use_latq = .true. ; use_lonh = .true. ! x=CENTER, y=NORTH_FACE - case ('1') ; ! x=0, y=0 - case default - call MOM_error(FATAL, "MOM_axis:get_var_dimension_features "//& - "Unrecognized hor_grid argument "//trim(grid_string_id)) - end select -end subroutine get_horizontal_grid_logic - -!> Define the time units for the input time value -function get_time_units(time_value) result(time_units_out) - real, intent(in) :: time_value !< numerical time value in seconds - !! i.e., before dividing by 86400. - ! local - character(len=10) :: time_units ! time units - character(len=10) :: time_units_out ! time units trimmed - time_units = '' - time_units_out = '' - if (time_value < 0.0) then - time_units = "days" ! The default value. - elseif (mod(time_value,86400.0)==0.0) then - time_units = "days" - elseif ((time_value >= 0.99) .and. (time_value < 1.01)) then - time_units = "seconds" - elseif ((time_value >= 3599.0) .and. (time_value < 3601.0)) then - time_units = "hours" - elseif ((time_value >= 86399.0) .and. (time_value < 86401.0)) then - time_units = "days" - elseif ((time_value >= 3.0e7) .and. (time_value < 3.2e7)) then - time_units = "years" - else - write(time_units,'(es8.2," s")') time_value - endif - time_units_out = trim(time_units) -end function get_time_units - -!> function to get the index of a time_value from a netCDF file -function get_time_index(filename, time_to_find) result (time_index) - character(len=*) :: filename ! name of the file to read in - real, intent(in) :: time_to_find ! time value to search for in file - ! local - type(fmsNetcdfFile_t) :: fileobj ! netCDF file object returned by open_file - real, allocatable, dimension(:) :: file_times ! array of time values read from file - integer :: dim_unlim_size, i, time_index - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - logical :: file_open_success - - time_index = 1 - dim_unlim_size = 0 - dim_unlim_name = "" - file_open_success = .false. - - if (.not. check_if_open(fileobj)) & - !call MOM_error(FATAL, "get_time_index_nodd: netcdf file object must be open.") - file_open_success=fms2_open_file(fileobj, trim(filename), "read", is_restart=.false.) - - call get_unlimited_dimension_name(fileobj, dim_unlim_name) - call get_dimension_size(fileObj, trim(dim_unlim_name), dim_unlim_size) - ! time index will be one more than the unlimited dimension size if the time_to_find is not in the file - if (dim_unlim_size .gt. 0) then - time_index = dim_unlim_size+1 - allocate(file_times(dim_unlim_size)) - call read_data(fileobj,trim(dim_unlim_name), file_times) - - do i=1,dim_unlim_size - if (ABS(file_times(i)-time_to_find) .gt. TINY(time_to_find)) then - continue - else - time_index = i - exit - endif - enddo - deallocate(file_times) - endif - if (check_if_open(fileobj)) call fms2_close_file(fileobj) -end function get_time_index - -!> register axes associated with a variable from a domain-decomposed netCDF file that are mapped to -!! a sub-domain (e.g., a supergrid). -!> \note The user must specify units for variables with longitude/x-axis and/or latitude/y-axis axes to obtain -!! the correct domain decomposition for the data buffer. -subroutine MOM_register_variable_axes_subdomain(fileObj, variableName, io_domain, position) - type(FmsNetcdfFile_t), intent(inout) :: fileObj !< netCDF file object returned by call to open_file - character(len=*), intent(in) :: variableName !< name of the variable - type(domain2d), intent(in) :: io_domain !< type that contains the mpp io domain - integer, optional, intent(in) :: position !< A flag indicating where this data is discretized - - ! Local variables - character(len=40) :: units ! units corresponding to a specific variable dimension - character(len=40), allocatable, dimension(:) :: dim_names ! variable dimension names - integer :: i, isg, ieg, isc, iec, jsg, jeg, jsc, jec, xlen, ylen - integer :: ndims ! number of dimensions - integer :: pos ! Discrete variable position. Default is CENTER - integer, allocatable, dimension(:) :: dimSizes ! variable dimension sizes - - if (.not. check_if_open(fileObj)) call MOM_error(FATAL,"MOM_axis:register_variable_axes_subdomain: The fileObj "// & - " has not been opened. Call fms2_open_file(fileObj,...) "// & - "before passing the fileObj argument to this function.") - - ! get variable dimension names and lengths - ndims = get_variable_num_dimensions(fileObj, trim(variableName)) - allocate(dimSizes(ndims)) - allocate(dim_names(ndims)) - call get_variable_size(fileObj, trim(variableName), dimSizes, broadcast=.true.) - call get_variable_dimension_names(fileObj, trim(variableName), dim_names) - - ! Get the lengths of the global indicies, using the discrete position of this variable - pos = CORNER ; if (present(position)) pos = position - call mpp_get_compute_domain(io_domain, xsize=xlen, ysize=ylen, position=pos) - ! register the axes - !>\note: This is not a comprehensive check for all possible supported horizontal axes associated with variables - !! read from netCDF files. Developers should add/remove cases as needed. - do i=1,ndims - !if (.not.(is_dimension_registered(fileObj, trim(dim_names(i))))) then - select case(trim(lowercase(dim_names(i)))) - case ("grid_x_t") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case ("nx") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("nxp") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("longitude") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("long") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("lon") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("lonh") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("lonq") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case("xh") - call register_axis(fileObj, trim(dim_names(i)), xlen) - case ("grid_y_t") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case ("ny") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("nyp") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("latitude") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("lat") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("lath") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("latq") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case("yh") - call register_axis(fileObj, trim(dim_names(i)), ylen) - case default ! assumes that the axis is not domain-decomposed - if (.not. is_dimension_unlimited(fileObj, trim(dim_names(i)))) & - call MOM_error(WARNING,"MOM_register_variable_axes_subdomain: the axis "//trim(dim_names(i))//& - "is not included in the valid x and y dimension cases. If the code hangs, check the whether "//& - "an x or y axis is being registered as a non-domain-decomposed variable, "//& - "and add it to the accepted cases if necessary.") - call register_axis(fileObj, trim(dim_names(i)), dimSizes(i)) - end select - ! endif - enddo - - if (allocated(dimSizes)) deallocate(dimSizes) - if (allocated(dim_names)) deallocate(dim_names) -end subroutine MOM_register_variable_axes_subdomain - -!> register axes associated with a variable from a domain-decomposed netCDF file -!> @note The user must specify units for variables with longitude/x-axis and/or latitude/y-axis axes -!! to obtain the correct domain decomposition for the data buffer. -subroutine MOM_register_variable_axes_full(fileObj, variableName, position) - type(FmsNetcdfDomainFile_t), intent(inout) :: fileObj !< netCDF file object returned by call to open_file - character(len=*), intent(in) :: variableName !< name of the variable - integer, optional, intent(in) :: position !< A flag indicating where this data is discretized - - ! Local variables - character(len=40) :: units ! units corresponding to a specific variable dimension - character(len=40), allocatable, dimension(:) :: dim_names ! variable dimension names - integer :: i - integer :: ndims ! number of dimensions - integer :: xPos, yPos ! domain positions for x and y axes. Default is CENTER - integer, allocatable, dimension(:) :: dimSizes ! variable dimension sizes - - if (.not. check_if_open(fileObj)) call MOM_error(FATAL,"MOM_axis:register_variable_axes: The fileObj has "// & - "not been opened. Call fms2_open_file(fileObj,...) before "// & - "passing the fileObj argument to this function.") - xpos = CENTER ; ypos = CENTER - if (present(position)) then - if ((position == CORNER) .or. (position == EAST_FACE)) xpos = EAST_FACE - if ((position == CORNER) .or. (position == NORTH_FACE)) ypos = NORTH_FACE - endif - - ! get variable dimension names and lengths - ndims = get_variable_num_dimensions(fileObj, trim(variableName)) - allocate(dimSizes(ndims)) - allocate(dim_names(ndims)) - call get_variable_size(fileObj, trim(variableName), dimSizes) - call get_variable_dimension_names(fileObj, trim(variableName), dim_names) - ! register the axes - !>@note: This is not a comprehensive check for all possible supported horizontal axes associated with variables - !! read from netCDF files. Developers should add/remove cases as needed. - do i=1,ndims - if (.not.(is_dimension_registered(fileobj, trim(dim_names(i))))) then - select case(trim(lowercase(dim_names(i)))) - case ("grid_x_t") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case ("nx") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("nxp") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("longitude") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("long") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("lon") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("lonh") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("lonq") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("xh") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case("i") - call register_axis(fileObj, trim(dim_names(i)),"x", domain_position=xPos) - case ("grid_y_t") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case ("ny") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("nyp") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("latitude") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("lat") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("lath") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("latq") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("yh") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case("j") - call register_axis(fileObj, trim(dim_names(i)),"y", domain_position=yPos) - case default ! assumes that the axis is not domain-decomposed - if (.not. is_dimension_unlimited(fileObj, trim(dim_names(i)))) & - call MOM_error(WARNING,"MOM_register_variable_axes_full: the axis "//trim(dim_names(i))//" is not "//& - "included in the valid x and y dimension cases. If the code hangs, check the whether "//& - "an x or y axis is being registered as a non-domain-decomposed variable, "//& - "and add it to the accepted cases if necessary.") - call register_axis(fileObj, trim(dim_names(i)), dimSizes(i)) - end select - endif - enddo - - deallocate(dimSizes) - deallocate(dim_names) -end subroutine MOM_register_variable_axes_full - - -!> convert the variable checksum integer(s) to a single string -!! If there is more than 1 checksum, commas are inserted between -!! each checksum value in the output string -function convert_checksum_to_string(checksum_int) result (checksum_string) - integer(kind=8), intent(in) :: checksum_int !< checksum integer values -! local - character(len=64) :: checksum_string - integer :: i - - checksum_string = '' - - write (checksum_string,'(Z16)') checksum_int ! Z16 is the hexadecimal format code - -end function convert_checksum_to_string - - -end module MOM_axis diff --git a/config_src/infra/FMS2/MOM_interp_infra.F90 b/config_src/infra/FMS2/MOM_interp_infra.F90 index ca5b2b8516..170573f7ec 100644 --- a/config_src/infra/FMS2/MOM_interp_infra.F90 +++ b/config_src/infra/FMS2/MOM_interp_infra.F90 @@ -4,9 +4,9 @@ module MOM_interp_infra ! This file is part of MOM6. See LICENSE.md for the license. use MOM_domain_infra, only : MOM_domain_type, domain2d -use MOM_io_infra, only : axistype use MOM_time_manager, only : time_type use horiz_interp_mod, only : horiz_interp_new, horiz_interp, horiz_interp_init, horiz_interp_type +use mpp_io_mod, only : axistype, mpp_get_axis_data use time_interp_external_mod, only : time_interp_external use time_interp_external_mod, only : init_external_field, time_interp_external_init use time_interp_external_mod, only : get_external_field_size @@ -16,7 +16,7 @@ module MOM_interp_infra public :: horiz_interp_type, horiz_interp_init public :: time_interp_extern, init_extern_field, time_interp_external_init -public :: get_external_field_info +public :: get_external_field_info, axistype, get_axis_data public :: run_horiz_interp, build_horiz_interp_weights !> Read a field based on model time, and rotate to the model domain. @@ -114,6 +114,15 @@ subroutine build_horiz_interp_weights_2d_to_2d(Interp, lon_in, lat_in, lon_out, end subroutine build_horiz_interp_weights_2d_to_2d +!> Extracts and returns the axis data stored in an axistype. +subroutine get_axis_data( axis, dat ) + type(axistype), intent(in) :: axis !< An axis type + real, dimension(:), intent(out) :: dat !< The data in the axis variable + + call mpp_get_axis_data( axis, dat ) +end subroutine get_axis_data + + !> get size of an external field from field index function get_extern_field_size(index) diff --git a/config_src/infra/FMS2/MOM_io_infra.F90 b/config_src/infra/FMS2/MOM_io_infra.F90 index 22548218d1..df9d6dc7ca 100644 --- a/config_src/infra/FMS2/MOM_io_infra.F90 +++ b/config_src/infra/FMS2/MOM_io_infra.F90 @@ -5,20 +5,28 @@ module MOM_io_infra use MOM_domain_infra, only : MOM_domain_type, rescale_comp_data, AGRID, BGRID_NE, CGRID_NE use MOM_domain_infra, only : domain2d, domain1d, CENTER, CORNER, NORTH_FACE, EAST_FACE -use MOM_error_infra, only : MOM_error=>MOM_err, NOTE, FATAL, WARNING - -use MOM_read_data_fms2, only : prepare_to_read_var -use fms2_io_mod, only : fms2_open_file => open_file, fms2_close_file => close_file -use fms2_io_mod, only : FmsNetcdfDomainFile_t, fms2_read_data => read_data, check_if_open +use MOM_error_infra, only : MOM_error=>MOM_err, NOTE, FATAL, WARNING, is_root_PE +use MOM_string_functions, only : lowercase + +use fms2_io_mod, only : fms2_open_file => open_file, check_if_open, fms2_close_file => close_file +use fms2_io_mod, only : FmsNetcdfDomainFile_t, FmsNetcdfFile_t, fms2_read_data => read_data +use fms2_io_mod, only : get_unlimited_dimension_name, get_num_dimensions, get_num_variables +use fms2_io_mod, only : get_variable_names, variable_exists, get_variable_size, get_variable_units +use fms2_io_mod, only : register_field, write_data, register_variable_attribute +use fms2_io_mod, only : variable_att_exists, get_variable_attribute, get_variable_num_dimensions +use fms2_io_mod, only : get_variable_dimension_names, is_dimension_registered, get_dimension_size +use fms2_io_mod, only : is_dimension_unlimited, register_axis, unlimited +use fms2_io_mod, only : get_global_io_domain_indices use fms_mod, only : write_version_number, open_namelist_file, check_nml_error use fms_io_mod, only : file_exist, field_exist, field_size, read_data use fms_io_mod, only : fms_io_exit, get_filename_appendix +use mpp_domains_mod, only : mpp_get_compute_domain, mpp_get_global_domain use mpp_io_mod, only : mpp_open, mpp_close, mpp_flush use mpp_io_mod, only : mpp_write_meta, mpp_write use mpp_io_mod, only : mpp_get_atts, mpp_attribute_exist -use mpp_io_mod, only : mpp_get_axes, axistype, mpp_get_axis_data -use mpp_io_mod, only : mpp_get_fields, fieldtype +use mpp_io_mod, only : mpp_get_axes, mpp_axistype=>axistype, mpp_get_axis_data +use mpp_io_mod, only : mpp_get_fields, mpp_fieldtype=>fieldtype use mpp_io_mod, only : mpp_get_info, mpp_get_times use mpp_io_mod, only : mpp_io_init ! These are encoding constants. @@ -36,9 +44,9 @@ module MOM_io_infra public :: MOM_read_data, MOM_read_vector, write_metadata, write_field public :: field_exists, get_field_atts, get_field_size, get_axis_data, read_field_chksum public :: io_infra_init, io_infra_end, MOM_namelist_file, check_namelist_error, write_version -! These types are inherited from underlying infrastructure code, to act as containers for -! information about fields and axes, respectively, and are opaque to this module. -public :: fieldtype, axistype +! These types act as containers for information about files, fields and axes, respectively, +! and may also wrap opaque types from the underlying infrastructure. +public :: file_type, fieldtype, axistype ! These are encoding constant parmeters. public :: ASCII_FILE, NETCDF_FILE, SINGLE_FILE, MULTIPLE public :: APPEND_FILE, READONLY_FILE, OVERWRITE_FILE, WRITEONLY_FILE @@ -97,15 +105,40 @@ module MOM_io_infra end interface flush_file !> Type for holding a handle to an open file and related information -type, public :: file_type ; private +type :: file_type ; private integer :: unit = -1 !< The framework identfier or netCDF unit number of an output file + type(FmsNetcdfDomainFile_t), pointer :: fileobj => NULL() !< A domain-decomposed + !! file object that is open for writing character(len=:), allocatable :: filename !< The path to this file, if it is open logical :: open_to_read = .false. !< If true, this file or fileset can be read logical :: open_to_write = .false. !< If true, this file or fileset can be written to + integer :: num_times !< The number of time levels in this file + real :: file_time !< The time of the latest entry in the file. + logical :: FMS2_file !< If true, this file-type is to be used with FMS2 interfaces. end type file_type -!> For now, this is hard-coded to exercise the new FMS2 interfaces. -logical :: FMS2_reads = .true. +!> This type is a container for information about a variable in a file. +type :: fieldtype ; private + character(len=256) :: name !< The name of this field in the files. + type(mpp_fieldtype) :: FT !< The FMS1 field-type that this type wraps + character(len=:), allocatable :: longname !< The long name for this field + character(len=:), allocatable :: units !< The units for this field + integer(kind=int64) :: chksum_read !< A checksum that has been read from a file + logical :: valid_chksum !< If true, this field has a valid checksum value. + logical :: FMS2_field !< If true, this field-type should be used with FMS2 interfaces. +end type fieldtype + +!> This type is a container for information about an axis in a file. +type :: axistype ; private + character(len=256) :: name !< The name of this axis in the files. + type(mpp_axistype) :: AT !< The FMS1 axis-type that this type wraps + real, allocatable, dimension(:) :: ax_data !< The values of the data on the axis. + logical :: domain_decomposed = .false. !< True if axis is domain-decomposed +end type axistype + +!> For now, these module-variables are hard-coded to exercise the new FMS2 interfaces. +logical :: FMS2_reads = .true. +logical :: FMS2_writes = .true. contains @@ -115,17 +148,11 @@ subroutine read_field_chksum(field, chksum, valid_chksum) type(fieldtype), intent(in) :: field !< The field whose checksum attribute is to be read. integer(kind=int64), intent(out) :: chksum !< The checksum for the field. logical, intent(out) :: valid_chksum !< If true, chksum has been successfully read. - ! Local variables - integer(kind=int64), dimension(3) :: checksum_file - checksum_file(:) = -1 - valid_chksum = mpp_attribute_exist(field, "checksum") - if (valid_chksum) then - call get_field_atts(field, checksum=checksum_file) - chksum = checksum_file(1) - else - chksum = -1 - endif + chksum = -1 + valid_chksum = field%valid_chksum + if (valid_chksum) chksum = field%chksum_read + end subroutine read_field_chksum !> Returns true if the named file or its domain-decomposed variant exists. @@ -156,7 +183,7 @@ end function FMS_file_exists logical function file_is_open(IO_handle) type(file_type), intent(in) :: IO_handle !< Handle to a file to inquire about - file_is_open = (IO_handle%unit >= 0) + file_is_open = ((IO_handle%unit >= 0) .or. associated(IO_handle%fileobj)) end function file_is_open !> closes a file (or fileset). If the file handle does not point to an open file, @@ -164,9 +191,16 @@ end function file_is_open subroutine close_file_type(IO_handle) type(file_type), intent(inout) :: IO_handle !< The I/O handle for the file to be closed - call mpp_close(IO_handle%unit) + if (associated(IO_handle%fileobj)) then + call fms2_close_file(IO_handle%fileobj) + deallocate(IO_handle%fileobj) + else + call mpp_close(IO_handle%unit) + endif if (allocated(IO_handle%filename)) deallocate(IO_handle%filename) IO_handle%open_to_read = .false. ; IO_handle%open_to_write = .false. + IO_handle%num_times = 0 ; IO_handle%file_time = 0.0 + IO_handle%FMS2_file = .false. end subroutine close_file_type !> closes a file. If the unit does not point to an open file, @@ -178,10 +212,14 @@ subroutine close_file_unit(unit) end subroutine close_file_unit !> Ensure that the output stream associated with a file handle is fully sent to disk. -subroutine flush_file_type(file) - type(file_type), intent(in) :: file !< The I/O handle for the file to flush +subroutine flush_file_type(IO_handle) + type(file_type), intent(in) :: IO_handle !< The I/O handle for the file to flush - call mpp_flush(file%unit) + if (associated(IO_handle%fileobj)) then + ! There does not appear to be an fms2 flush call. + else + call mpp_flush(IO_handle%unit) + endif end subroutine flush_file_type !> Ensure that the output stream associated with a unit is fully sent to disk. @@ -272,20 +310,80 @@ subroutine open_file_type(IO_handle, filename, action, MOM_domain, threading, fi !! to threading=MULTIPLE write to the same file (SINGLE_FILE) !! or to one file per PE (MULTIPLE, the default). - if (present(MOM_Domain)) then - call mpp_open(IO_handle%unit, filename, action=action, form=NETCDF_FILE, threading=threading, & + ! Local variables + type(FmsNetcdfDomainFile_t) :: fileobj_read ! A handle to a domain-decomposed file for obtaining information + ! about the exiting time axis entries in append mode. + logical :: success ! If true, the file was opened successfully + integer :: file_mode ! An integer that encodes whether the file is to be opened for + ! reading, writing or appending + character(len=40) :: mode ! A character string that encodes whether the file is to be opened for + ! reading, writing or appending + character(len=:), allocatable :: filename_tmp ! A copy of filename with .nc appended if necessary. + character(len=256) :: dim_unlim_name ! name of the unlimited dimension in the file + integer :: index_nc + + if (IO_handle%open_to_write) then + call MOM_error(WARNING, "open_file_type called for file "//trim(filename)//& + " with an IO_handle that is already open to to write.") + return + endif + if (IO_handle%open_to_read) then + call MOM_error(FATAL, "open_file_type called for file "//trim(filename)//& + " with an IO_handle that is already open to to read.") + endif + + file_mode = WRITEONLY_FILE ; if (present(action)) file_mode = action + + if (FMS2_writes .and. present(MOM_Domain)) then + if (.not.associated(IO_handle%fileobj)) allocate (IO_handle%fileobj) + + ! The FMS1 interface automatically appends .nc if necessary, but FMS2 interface does not. + index_nc = index(trim(filename), ".nc") + if (index_nc > 0) then + filename_tmp = trim(filename) + else + filename_tmp = trim(filename)//".nc" + if (is_root_PE()) call MOM_error(WARNING, "Open_file is appending .nc to the filename "//trim(filename)) + endif + + if (file_mode == WRITEONLY_FILE) then ; mode = "write" + elseif (file_mode == APPEND_FILE) then ; mode = "append" + elseif (file_mode == OVERWRITE_FILE) then ; mode = "overwrite" + elseif (file_mode == READONLY_FILE) then ; mode = "read" + else + call MOM_error(FATAL, "open_file_type called with unrecognized action.") + endif + + IO_handle%num_times = 0 + IO_handle%file_time = 0.0 + if ((file_mode == APPEND_FILE) .and. file_exists(filename_tmp, MOM_Domain)) then + ! Determine the latest file time and number of records so far. + success = fms2_open_file(fileObj_read, trim(filename_tmp), "read", MOM_domain%mpp_domain) + call get_unlimited_dimension_name(fileObj_read, dim_unlim_name) + if (len_trim(dim_unlim_name) > 0) & + call get_dimension_size(fileObj_read, trim(dim_unlim_name), IO_handle%num_times) + if (IO_handle%num_times > 0) & + call fms2_read_data(fileObj_read, trim(dim_unlim_name), IO_handle%file_time, & + unlim_dim_level=IO_handle%num_times) + call fms2_close_file(fileObj_read) + endif + + success = fms2_open_file(IO_handle%fileobj, trim(filename_tmp), trim(mode), MOM_domain%mpp_domain) + if (.not.success) call MOM_error(FATAL, "Unable to open file "//trim(filename_tmp)) + IO_handle%FMS2_file = .true. + elseif (present(MOM_Domain)) then + call mpp_open(IO_handle%unit, filename, action=file_mode, form=NETCDF_FILE, threading=threading, & fileset=fileset, domain=MOM_Domain%mpp_domain) + IO_handle%FMS2_file = .false. else - call mpp_open(IO_handle%unit, filename, action=action, form=NETCDF_FILE, threading=threading, & + call mpp_open(IO_handle%unit, filename, action=file_mode, form=NETCDF_FILE, threading=threading, & fileset=fileset) + IO_handle%FMS2_file = .false. endif IO_handle%filename = trim(filename) - if (present(action)) then - if (action == READONLY_FILE) then - IO_handle%open_to_read = .true. ; IO_handle%open_to_write = .false. - else - IO_handle%open_to_read = .false. ; IO_handle%open_to_write = .true. - endif + + if (file_mode == READONLY_FILE) then + IO_handle%open_to_read = .true. ; IO_handle%open_to_write = .false. else IO_handle%open_to_read = .false. ; IO_handle%open_to_write = .true. endif @@ -319,43 +417,60 @@ subroutine get_filename_suffix(suffix) end subroutine get_filename_suffix -!> Get information about the number of dimensions, variables, global attributes and time levels +!> Get information about the number of dimensions, variables and time levels !! in the file associated with an open file unit -subroutine get_file_info(IO_handle, ndim, nvar, natt, ntime) +subroutine get_file_info(IO_handle, ndim, nvar, ntime) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for I/O integer, optional, intent(out) :: ndim !< The number of dimensions in the file integer, optional, intent(out) :: nvar !< The number of variables in the file - integer, optional, intent(out) :: natt !< The number of global attributes in the file integer, optional, intent(out) :: ntime !< The number of time levels in the file ! Local variables + character(len=256) :: dim_unlim_name ! name of the unlimited dimension in the file integer :: ndims, nvars, natts, ntimes - call mpp_get_info(IO_handle%unit, ndims, nvars, natts, ntimes ) + if (IO_handle%FMS2_file) then + if (present(ndim)) ndim = get_num_dimensions(IO_handle%fileobj) + if (present(nvar)) nvar = get_num_variables(IO_handle%fileobj) + if (present(ntime)) then + ntime = 0 + call get_unlimited_dimension_name(IO_handle%fileobj, dim_unlim_name) + if (len_trim(dim_unlim_name) > 0) & + call get_dimension_size(IO_handle%fileobj, trim(dim_unlim_name), ntime) + endif + else + call mpp_get_info(IO_handle%unit, ndims, nvars, natts, ntimes ) - if (present(ndim)) ndim = ndims - if (present(nvar)) nvar = nvars - if (present(natt)) natt = natts - if (present(ntime)) ntime = ntimes + if (present(ndim)) ndim = ndims + if (present(nvar)) nvar = nvars + if (present(ntime)) ntime = ntimes + endif end subroutine get_file_info !> Get the times of records from a file - !### Modify this to also convert to time_type, using information about the dimensions? subroutine get_file_times(IO_handle, time_values, ntime) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for I/O real, allocatable, dimension(:), intent(inout) :: time_values !< The real times for the records in file. integer, optional, intent(out) :: ntime !< The number of time levels in the file - integer :: ntimes + character(len=256) :: dim_unlim_name ! name of the unlimited dimension in the file + integer :: ntimes ! The number of time levels in the file + + !### Modify this routine to optionally convert to time_type, using information about the dimensions? if (allocated(time_values)) deallocate(time_values) call get_file_info(IO_handle, ntime=ntimes) if (present(ntime)) ntime = ntimes if (ntimes > 0) then allocate(time_values(ntimes)) - call mpp_get_times(IO_handle%unit, time_values) + if (IO_handle%FMS2_file) then + call get_unlimited_dimension_name(IO_handle%fileobj, dim_unlim_name) + call fms2_read_data(IO_handle%fileobj, trim(dim_unlim_name), time_values) + else + call mpp_get_times(IO_handle%unit, time_values) + endif endif end subroutine get_file_times @@ -365,7 +480,48 @@ subroutine get_file_fields(IO_handle, fields) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for I/O type(fieldtype), dimension(:), intent(inout) :: fields !< Field-type descriptions of all of !! the variables in a file. - call mpp_get_fields(IO_handle%unit, fields) + type(mpp_fieldtype), dimension(size(fields)) :: mpp_fields ! Fieldtype structures for the variables + character(len=256), dimension(size(fields)) :: var_names ! The names of all variables + character(len=256) :: units ! The units of a variable as recorded in the file + character(len=2048) :: longname ! The long-name of a variable as recorded in the file + character(len=64) :: checksum_char ! The hexadecimal checksum read from the file + integer(kind=int64), dimension(3) :: checksum_file ! The checksums for a variable in the file + integer :: nvar ! The number of variables in the file + integer :: i + + nvar = size(fields) + ! Local variables + if (IO_handle%FMS2_file) then + call get_variable_names(IO_handle%fileobj, var_names) + do i=1,nvar + fields(i)%name = trim(var_names(i)) + longname = "" + call get_variable_attribute(IO_handle%fileobj, var_names(i), 'long_name', longname) + fields(i)%longname = trim(longname) + units = "" + call get_variable_attribute(IO_handle%fileobj, var_names(i), 'units', units) + fields(i)%units = trim(units) + + fields(i)%valid_chksum = variable_att_exists(IO_handle%fileobj, var_names(i), "checksum") + if (fields(i)%valid_chksum) then + call get_variable_attribute(IO_handle%fileobj, var_names(i), 'checksum', checksum_char) + ! If there are problems, there might need to be code added to handle commas. + read (checksum_char(1:16), '(Z16)') fields(i)%chksum_read + endif + enddo + else + call mpp_get_fields(IO_handle%unit, mpp_fields) + do i=1,nvar + fields(i)%FT = mpp_fields(i) + call mpp_get_atts(fields(i)%FT, name=fields(i)%name, units=units, longname=longname, & + checksum=checksum_file) + fields(i)%longname = trim(longname) + fields(i)%units = trim(units) + fields(i)%valid_chksum = mpp_attribute_exist(fields(i)%FT, "checksum") + if (fields(i)%valid_chksum) fields(i)%chksum_read = checksum_file(1) + enddo + endif + end subroutine get_file_fields !> Extract information from a field type, as stored or as found in a file @@ -376,7 +532,12 @@ subroutine get_field_atts(field, name, units, longname, checksum) character(len=*), optional, intent(out) :: longname !< The long name of the variable integer(kind=int64), dimension(:), & optional, intent(out) :: checksum !< The checksums of the variable in a file - call mpp_get_atts(field, name=name, units=units, longname=longname, checksum=checksum) + + if (present(name)) name = trim(field%name) + if (present(units)) units = trim(field%units) + if (present(longname)) longname = trim(field%longname) + if (present(checksum)) checksum = field%chksum_read + end subroutine get_field_atts !> Field_exists returns true if the field indicated by field_name is present in the @@ -389,7 +550,44 @@ function field_exists(filename, field_name, domain, no_domain, MOM_domain) type(MOM_domain_type), optional, intent(in) :: MOM_Domain !< A MOM_Domain that describes the decomposition logical :: field_exists !< True if filename exists and field_name is in filename - if (present(MOM_domain)) then + ! Local variables + type(FmsNetcdfDomainFile_t) :: fileObj_dd ! A handle to a domain-decomposed file for obtaining information + ! about the exiting time axis entries in append mode. + type(FmsNetcdfFile_t) :: fileObj_simple ! A handle to a non-domain-decomposed file for obtaining information + ! about the exiting time axis entries in append mode. + logical :: success ! If true, the file was opened successfully + logical :: domainless ! If true, this file does not use a domain-decomposed file. + + domainless = .not.(present(MOM_domain) .or. present(domain)) + if (present(no_domain)) then + if (domainless .and. .not.no_domain) call MOM_error(FATAL, & + "field_exists: When no_domain is present and false, a domain must be supplied in query about "//& + trim(field_name)//" in file "//trim(filename)) + domainless = no_domain + endif + + if (FMS2_reads) then + field_exists = .false. + if (file_exists(filename)) then + if (domainless) then + success = fms2_open_file(fileObj_simple, trim(filename), "read") + if (success) then + field_exists = variable_exists(fileObj_simple, field_name) + call fms2_close_file(fileObj_simple) + endif + else + if (present(MOM_domain)) then + success = fms2_open_file(fileObj_dd, trim(filename), "read", MOM_domain%mpp_domain) + else + success = fms2_open_file(fileObj_dd, trim(filename), "read", domain) + endif + if (success) then + field_exists = variable_exists(fileobj_dd, field_name) + call fms2_close_file(fileObj_dd) + endif + endif + endif + elseif (present(MOM_domain)) then field_exists = field_exist(filename, field_name, domain=MOM_domain%mpp_domain, no_domain=no_domain) else field_exists = field_exist(filename, field_name, domain=domain, no_domain=no_domain) @@ -408,7 +606,32 @@ subroutine get_field_size(filename, fieldname, sizes, field_found, no_domain) logical, optional, intent(in) :: no_domain !< If present and true, do not check for file !! names with an appended tile number - call field_size(filename, fieldname, sizes, field_found=field_found, no_domain=no_domain) + ! Local variables + type(FmsNetcdfFile_t) :: fileobj_read ! A handle to a non-domain-decomposed file for obtaining information + ! about the exiting time axis entries in append mode. + logical :: success ! If true, the file was opened successfully + logical :: field_exists ! True if filename exists and field_name is in filename + integer :: i, ndims + + if (FMS2_reads) then + field_exists = .false. + if (file_exists(filename)) then + success = fms2_open_file(fileObj_read, trim(filename), "read") + if (success) then + field_exists = variable_exists(fileobj_read, fieldname) + if (field_exists) then + ndims = get_variable_num_dimensions(fileobj_read, fieldname) + if (ndims > size(sizes)) call MOM_error(FATAL, & + "get_field_size called with too few sizes for "//trim(fieldname)//" in "//trim(filename)) + call get_variable_size(fileobj_read, fieldname, sizes(1:ndims)) + do i=ndims+1,size(sizes) ; sizes(i) = 0 ; enddo + endif + endif + endif + if (present(field_found)) field_found = field_exists + else + call field_size(filename, fieldname, sizes, field_found=field_found, no_domain=no_domain) + endif end subroutine get_field_size @@ -417,7 +640,17 @@ subroutine get_axis_data( axis, dat ) type(axistype), intent(in) :: axis !< An axis type real, dimension(:), intent(out) :: dat !< The data in the axis variable - call mpp_get_axis_data( axis, dat ) + integer :: i + + ! This routine might not be needed for MOM6. + if (allocated(axis%ax_data)) then + if (size(axis%ax_data) > size(dat)) call MOM_error(FATAL, & + "get_axis_data called with too small of an output data array for "//trim(axis%name)) + do i=1,size(axis%ax_data) ; dat(i) = axis%ax_data(i) ; enddo + elseif (.not.FMS2_writes) then + call mpp_get_axis_data( axis%AT, dat ) + endif + end subroutine get_axis_data !> This routine uses the fms_io subroutine read_data to read a scalar named @@ -433,20 +666,40 @@ subroutine MOM_read_data_0d(filename, fieldname, data, timelevel, scale, MOM_Dom optional, intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition ! Local variables - type(FmsNetcdfDomainFile_t) :: fileobj ! A handle to a domain-decomposed file object - logical :: has_time_dim ! True if the variable has an unlimited time axis. + type(FmsNetcdfFile_t) :: fileObj ! A handle to a non-domain-decomposed file + type(FmsNetcdfDomainFile_t) :: fileobj_DD ! A handle to a domain-decomposed file object character(len=96) :: var_to_read ! Name of variable to read from the netcdf file - logical :: success ! True if the file was successfully opened + logical :: has_time_dim ! True if the variable has an unlimited time axis. + logical :: success ! True if the file was successfully opened if (present(MOM_Domain) .and. FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj_DD, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive variable name in the file and prepare to read it. - call prepare_to_read_var(fileobj, fieldname, MOM_domain, "MOM_read_data_1d: ", filename, & + call prepare_to_read_var(fileobj_DD, fieldname, "MOM_read_data_0d: ", filename, & var_to_read, has_time_dim, timelevel) + ! Read the data. + if (present(timelevel) .and. has_time_dim) then + call fms2_read_data(fileobj_DD, var_to_read, data, unlim_dim_level=timelevel) + else + call fms2_read_data(fileobj_DD, var_to_read, data) + endif + + ! Close the file-set. + if (check_if_open(fileobj_DD)) call fms2_close_file(fileobj_DD) + elseif (FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileObj, trim(filename), "read") + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file, and determine whether it + ! has a time dimension. + call find_varname_in_file(fileObj, fieldname, "MOM_read_data_0d: ", filename, & + var_to_read, has_time_dim, timelevel) + ! Read the data. if (present(timelevel) .and. has_time_dim) then call fms2_read_data(fileobj, var_to_read, data, unlim_dim_level=timelevel) @@ -481,20 +734,40 @@ subroutine MOM_read_data_1d(filename, fieldname, data, timelevel, scale, MOM_Dom optional, intent(in) :: MOM_Domain !< The MOM_Domain that describes the decomposition ! Local variables - type(FmsNetcdfDomainFile_t) :: fileobj ! A handle to a domain-decomposed file object - logical :: has_time_dim ! True if the variable has an unlimited time axis. + type(FmsNetcdfFile_t) :: fileObj ! A handle to a non-domain-decomposed file + type(FmsNetcdfDomainFile_t) :: fileobj_DD ! A handle to a domain-decomposed file object character(len=96) :: var_to_read ! Name of variable to read from the netcdf file - logical :: success ! True if the file was successfully opened + logical :: has_time_dim ! True if the variable has an unlimited time axis. + logical :: success ! True if the file was successfully opened if (present(MOM_Domain) .and. FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj_DD, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive variable name in the file and prepare to read it. - call prepare_to_read_var(fileobj, fieldname, MOM_domain, "MOM_read_data_1d: ", filename, & + call prepare_to_read_var(fileobj_DD, fieldname, "MOM_read_data_1d: ", filename, & var_to_read, has_time_dim, timelevel) + ! Read the data. + if (present(timelevel) .and. has_time_dim) then + call fms2_read_data(fileobj_DD, var_to_read, data, unlim_dim_level=timelevel) + else + call fms2_read_data(fileobj_DD, var_to_read, data) + endif + + ! Close the file-set. + if (check_if_open(fileobj_DD)) call fms2_close_file(fileobj_DD) + elseif (FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileObj, trim(filename), "read") + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file, and determine whether it + ! has a time dimension. + call find_varname_in_file(fileObj, fieldname, "MOM_read_data_1d: ", filename, & + var_to_read, has_time_dim, timelevel) + ! Read the data. if (present(timelevel) .and. has_time_dim) then call fms2_read_data(fileobj, var_to_read, data, unlim_dim_level=timelevel) @@ -533,17 +806,17 @@ subroutine MOM_read_data_2d(filename, fieldname, data, MOM_Domain, & ! Local variables type(FmsNetcdfDomainFile_t) :: fileobj ! A handle to a domain-decomposed file object - logical :: has_time_dim ! True if the variable has an unlimited time axis. character(len=96) :: var_to_read ! Name of variable to read from the netcdf file - logical :: success ! True if the file was successfully opened + logical :: has_time_dim ! True if the variable has an unlimited time axis. + logical :: success ! True if the file was successfully opened if (FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive variable name in the file and prepare to read it. - call prepare_to_read_var(fileobj, fieldname, MOM_domain, "MOM_read_data_2d: ", filename, & + call prepare_to_read_var(fileobj, fieldname, "MOM_read_data_2d: ", filename, & var_to_read, has_time_dim, timelevel, position) ! Read the data. @@ -587,9 +860,41 @@ subroutine MOM_read_data_2d_region(filename, fieldname, data, start, nread, MOM_ real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied !! by before it is returned. - ! This subroutine does not have an FMS-2 variant yet. + ! Local variables + type(FmsNetcdfFile_t) :: fileObj ! A handle to a non-domain-decomposed file + type(FmsNetcdfDomainFile_t) :: fileobj_DD ! A handle to a domain-decomposed file object + character(len=96) :: var_to_read ! Name of variable to read from the netcdf file + logical :: success ! True if the file was successfully opened - if (present(MOM_Domain)) then + if (present(MOM_Domain) .and. FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileobj_DD, filename, "read", MOM_domain%mpp_domain) + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file and prepare to read it. + call prepare_to_read_var(fileobj_DD, fieldname, "MOM_read_data_2d_region: ", & + filename, var_to_read) + + ! Read the data. + call fms2_read_data(fileobj_DD, var_to_read, data, corner=start(1:2), edge_lengths=nread(1:2)) + + ! Close the file-set. + if (check_if_open(fileobj_DD)) call fms2_close_file(fileobj_DD) + elseif (FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileObj, trim(filename), "read") + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file, and determine whether it + ! has a time dimension. + call find_varname_in_file(fileObj, fieldname, "MOM_read_data_2d_region: ", filename, var_to_read) + + ! Read the data. + call fms2_read_data(fileobj, var_to_read, data, corner=start(1:2), edge_lengths=nread(1:2)) + + ! Close the file-set. + if (check_if_open(fileobj)) call fms2_close_file(fileobj) + elseif (present(MOM_Domain)) then ! Read the variable using the FMS-1 interface. call read_data(filename, fieldname, data, start, nread, domain=MOM_Domain%mpp_domain, & no_domain=no_domain) else @@ -624,17 +929,17 @@ subroutine MOM_read_data_3d(filename, fieldname, data, MOM_Domain, & ! Local variables type(FmsNetcdfDomainFile_t) :: fileobj ! A handle to a domain-decomposed file object - logical :: has_time_dim ! True if the variable has an unlimited time axis. character(len=96) :: var_to_read ! Name of variable to read from the netcdf file - logical :: success ! True if the file was successfully opened + logical :: has_time_dim ! True if the variable has an unlimited time axis. + logical :: success ! True if the file was successfully opened if (FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive variable name in the file and prepare to read it. - call prepare_to_read_var(fileobj, fieldname, MOM_domain, "MOM_read_data_3d: ", filename, & + call prepare_to_read_var(fileobj, fieldname, "MOM_read_data_3d: ", filename, & var_to_read, has_time_dim, timelevel, position) ! Read the data. @@ -681,11 +986,11 @@ subroutine MOM_read_data_4d(filename, fieldname, data, MOM_Domain, & if (FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive variable name in the file and prepare to read it. - call prepare_to_read_var(fileobj, fieldname, MOM_domain, "MOM_read_data_4d: ", filename, & + call prepare_to_read_var(fileobj, fieldname, "MOM_read_data_4d: ", filename, & var_to_read, has_time_dim, timelevel, position) ! Read the data. @@ -716,7 +1021,35 @@ subroutine MOM_read_data_0d_int(filename, fieldname, data, timelevel) integer, intent(inout) :: data !< The 1-dimensional array into which the data integer, optional, intent(in) :: timelevel !< The time level in the file to read - call read_data(filename, fieldname, data, timelevel=timelevel, no_domain=.true.) + ! Local variables + type(FmsNetcdfFile_t) :: fileObj ! A handle to a non-domain-decomposed file + logical :: has_time_dim ! True if the variable has an unlimited time axis. + character(len=96) :: var_to_read ! Name of variable to read from the netcdf file + logical :: success ! If true, the file was opened successfully + + ! This routine might not be needed for MOM6. + if (FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileObj, trim(filename), "read") + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file, and determine whether it + ! has a time dimension. + call find_varname_in_file(fileObj, fieldname, "MOM_read_data_0d_int: ", filename, & + var_to_read, has_time_dim, timelevel) + + ! Read the data. + if (present(timelevel) .and. has_time_dim) then + call fms2_read_data(fileobj, var_to_read, data, unlim_dim_level=timelevel) + else + call fms2_read_data(fileobj, var_to_read, data) + endif + + ! Close the file-set. + if (check_if_open(fileobj)) call fms2_close_file(fileobj) + else + call read_data(filename, fieldname, data, timelevel=timelevel, no_domain=.true.) + endif end subroutine MOM_read_data_0d_int @@ -728,7 +1061,36 @@ subroutine MOM_read_data_1d_int(filename, fieldname, data, timelevel) integer, dimension(:), intent(inout) :: data !< The 1-dimensional array into which the data integer, optional, intent(in) :: timelevel !< The time level in the file to read - call read_data(filename, fieldname, data, timelevel=timelevel, no_domain=.true.) + ! Local variables + type(FmsNetcdfFile_t) :: fileObj ! A handle to a non-domain-decomposed file for obtaining information + ! about the exiting time axis entries in append mode. + logical :: has_time_dim ! True if the variable has an unlimited time axis. + character(len=96) :: var_to_read ! Name of variable to read from the netcdf file + logical :: success ! If true, the file was opened successfully + + ! This routine might not be needed for MOM6. + if (FMS2_reads) then + ! Open the FMS2 file-set. + success = fms2_open_file(fileObj, trim(filename), "read") + if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) + + ! Find the matching case-insensitive variable name in the file, and determine whether it + ! has a time dimension. + call find_varname_in_file(fileObj, fieldname, "MOM_read_data_1d_int: ", filename, & + var_to_read, has_time_dim, timelevel) + + ! Read the data. + if (present(timelevel) .and. has_time_dim) then + call fms2_read_data(fileobj, var_to_read, data, unlim_dim_level=timelevel) + else + call fms2_read_data(fileobj, var_to_read, data) + endif + + ! Close the file-set. + if (check_if_open(fileobj)) call fms2_close_file(fileobj) + else + call read_data(filename, fieldname, data, timelevel=timelevel, no_domain=.true.) + endif end subroutine MOM_read_data_1d_int @@ -767,13 +1129,13 @@ subroutine MOM_read_vector_2d(filename, u_fieldname, v_fieldname, u_data, v_data if (FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive u- and v-variable names in the file and prepare to read them. - call prepare_to_read_var(fileobj, u_fieldname, MOM_domain, "MOM_read_vector_2d: ", filename, & + call prepare_to_read_var(fileobj, u_fieldname, "MOM_read_vector_2d: ", filename, & u_var, has_time_dim, timelevel, position=u_pos) - call prepare_to_read_var(fileobj, v_fieldname, MOM_domain, "MOM_read_vector_2d: ", filename, & + call prepare_to_read_var(fileobj, v_fieldname, "MOM_read_vector_2d: ", filename, & v_var, has_time_dim, timelevel, position=v_pos) ! Read the u-data and v-data. There would already been an error message for one @@ -837,13 +1199,13 @@ subroutine MOM_read_vector_3d(filename, u_fieldname, v_fieldname, u_data, v_data if (FMS2_reads) then ! Open the FMS2 file-set. - success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain, is_restart=.false.) + success = fms2_open_file(fileobj, filename, "read", MOM_domain%mpp_domain) if (.not.success) call MOM_error(FATAL, "Failed to open "//trim(filename)) ! Find the matching case-insensitive u- and v-variable names in the file and prepare to read them. - call prepare_to_read_var(fileobj, u_fieldname, MOM_domain, "MOM_read_vector_3d: ", filename, & + call prepare_to_read_var(fileobj, u_fieldname, "MOM_read_vector_3d: ", filename, & u_var, has_time_dim, timelevel, position=u_pos) - call prepare_to_read_var(fileobj, v_fieldname, MOM_domain, "MOM_read_vector_3d: ", filename, & + call prepare_to_read_var(fileobj, v_fieldname, "MOM_read_vector_3d: ", filename, & v_var, has_time_dim, timelevel, position=v_pos) ! Read the u-data and v-data, dangerously assuming either both or neither have time dimensions. @@ -873,80 +1235,521 @@ subroutine MOM_read_vector_3d(filename, u_fieldname, v_fieldname, u_data, v_data end subroutine MOM_read_vector_3d +!> Find the case-sensitive name of the variable in a netCDF file with a case-insensitive name match. +!! Optionally also determine whether this variable has an unlimited time dimension. +subroutine find_varname_in_file(fileobj, fieldname, err_header, filename, var_to_read, has_time_dim, timelevel) + type(FmsNetcdfFile_t), intent(inout) :: fileobj !< An FMS2 handle to an open NetCDF file + character(len=*), intent(in) :: fieldname !< The variable name to seek in the file + character(len=*), intent(in) :: err_header !< A descriptive prefix for error messages + character(len=*), intent(in) :: filename !< The name of the file to read + character(len=*), intent(out) :: var_to_read !< The variable name to read from the file + logical, optional, intent(out) :: has_time_dim !< Indicates whether fieldname has a time dimension + integer, optional, intent(in) :: timelevel !< A time level to read + + ! Local variables + logical :: variable_found ! Is a case-insensitive version of the variable found in the netCDF file? + character(len=256), allocatable, dimension(:) :: var_names ! The names of all the variables in the netCDF file + character(len=256), allocatable :: dim_names(:) ! The names of a variable's dimensions + integer :: nvars ! The number of variables in the file + integer :: dim_unlim_size ! The current size of the unlimited (time) dimension in the file. + integer :: num_var_dims ! The number of dimensions a variable has in the file. + integer :: time_dim ! The position of the unlimited (time) dimension for a variable, or -1 + ! if it has no unlimited dimension. + integer :: i + + ! Open the file if necessary + if (.not.check_if_open(fileobj)) & + call MOM_error(FATAL, trim(err_header)//trim(filename)//" was not open in call to find_varname_in_file.") + + ! Search for the variable in the file, looking for the case-sensitive name first. + if (variable_exists(fileobj, trim(fieldname))) then + var_to_read = trim(fieldname) + else ! Look for case-insensitive variable name matches. + nvars = get_num_variables(fileobj) + if (nvars < 1) call MOM_error(FATAL, "nvars is less than 1 for file "//trim(filename)) + allocate(var_names(nvars)) + call get_variable_names(fileobj, var_names) + + ! search for the variable in the file + variable_found = .false. + do i=1,nvars + if (lowercase(trim(var_names(i))) == lowercase(trim(fieldname))) then + variable_found = .true. + var_to_read = trim(var_names(i)) + exit + endif + enddo + if (.not.(variable_found)) & + call MOM_error(FATAL, trim(err_header)//trim(fieldname)//" not found in "//trim(filename)) + deallocate(var_names) + endif + + ! FMS2 can not handle a timelevel argument if the variable does not have one in the file, + ! so some error checking and logic are required. + if (present(has_time_dim) .or. present(timelevel)) then + time_dim = -1 + + num_var_dims = get_variable_num_dimensions(fileobj, trim(var_to_read)) + allocate(dim_names(num_var_dims)) ; dim_names(:) = "" + call get_variable_dimension_names(fileobj, trim(var_to_read), dim_names) + + do i=1,num_var_dims + if (is_dimension_unlimited(fileobj, dim_names(i))) then + time_dim = i + if (present(timelevel)) then + call get_dimension_size(fileobj, dim_names(i), dim_unlim_size) + if ((timelevel > dim_unlim_size) .and. is_root_PE()) call MOM_error(FATAL, & + trim(err_header)//"Attempting to read a time level of "//trim(var_to_read)//& + " that exceeds the size of the time dimension in "//trim(filename)) + endif + exit + endif + enddo + deallocate(dim_names) + + if (present(timelevel) .and. (time_dim < 0) .and. is_root_PE()) & + call MOM_error(WARNING, trim(err_header)//"time level specified, but the variable "//& + trim(var_to_read)//" does not have an unlimited dimension in "//trim(filename)) + if ((.not.present(timelevel)) .and. (time_dim > 0) .and. is_root_PE()) & + call MOM_error(WARNING, trim(err_header)//"The variable "//trim(var_to_read)//& + " has an unlimited dimension in "//trim(filename)//" but no time level is specified.") + if (present(has_time_dim)) has_time_dim = (time_dim > 0) + endif + +end subroutine find_varname_in_file + + +!> Find the case-insensitive name match with a variable in an open domain-decomposed file-set, +!! prepare FMS2 to read this variable, and return some information needed to call fms2_read_data +!! correctly for this variable and file. +subroutine prepare_to_read_var(fileobj, fieldname, err_header, filename, var_to_read, & + has_time_dim, timelevel, position) + type(FmsNetcdfDomainFile_t), intent(inout) :: fileobj !< An FMS2 handle to an open domain-decomposed file + character(len=*), intent(in) :: fieldname !< The variable name to seek in the file + character(len=*), intent(in) :: err_header !< A descriptive prefix for error messages + character(len=*), intent(in) :: filename !< The name of the file to read + character(len=*), intent(out) :: var_to_read !< The variable name to read from the file + logical, optional, intent(out) :: has_time_dim !< Indicates whether fieldname has a time dimension + integer, optional, intent(in) :: timelevel !< A time level to read + integer, optional, intent(in) :: position !< A flag indicating where this variable is discretized + + ! Local variables + logical :: variable_found ! Is a case-insensitive version of the variable found in the netCDF file? + character(len=256), allocatable, dimension(:) :: var_names ! The names of all the variables in the netCDF file + character(len=256), allocatable :: dim_names(:) ! The names of a variable's dimensions + integer :: nvars ! The number of variables in the file. + integer :: dim_unlim_size ! The current size of the unlimited (time) dimension in the file. + integer :: num_var_dims ! The number of dimensions a variable has in the file. + integer :: time_dim ! The position of the unlimited (time) dimension for a variable, or -1 + ! if it has no unlimited dimension. + integer :: i + + ! Open the file if necessary + if (.not.check_if_open(fileobj)) & + call MOM_error(FATAL, trim(err_header)//trim(filename)//" was not open in call to prepare_to_read_var.") + + ! Search for the variable in the file, looking for the case-sensitive name first. + if (variable_exists(fileobj, trim(fieldname))) then + var_to_read = trim(fieldname) + else ! Look for case-insensitive variable name matches. + nvars = get_num_variables(fileobj) + if (nvars < 1) call MOM_error(FATAL, "nvars is less than 1 for file "//trim(filename)) + allocate(var_names(nvars)) + call get_variable_names(fileobj, var_names) + + variable_found = .false. + do i=1,nvars + if (lowercase(trim(var_names(i))) == lowercase(trim(fieldname))) then + variable_found = .true. + var_to_read = trim(var_names(i)) + exit + endif + enddo + if (.not.(variable_found)) & + call MOM_error(FATAL, trim(err_header)//trim(fieldname)//" not found in "//trim(filename)) + deallocate(var_names) + endif + + ! FMS2 can not handle a timelevel argument if the variable does not have one in the file, + ! so some error checking and logic are required. + if (present(has_time_dim) .or. present(timelevel)) then + time_dim = -1 + + num_var_dims = get_variable_num_dimensions(fileobj, trim(var_to_read)) + allocate(dim_names(num_var_dims)) ; dim_names(:) = "" + call get_variable_dimension_names(fileobj, trim(var_to_read), dim_names) + + do i=1,num_var_dims + if (is_dimension_unlimited(fileobj, dim_names(i))) then + time_dim = i + if (present(timelevel)) then + call get_dimension_size(fileobj, dim_names(i), dim_unlim_size) + if ((timelevel > dim_unlim_size) .and. is_root_PE()) call MOM_error(FATAL, & + trim(err_header)//"Attempting to read a time level of "//trim(var_to_read)//& + " that exceeds the size of the time dimension in "//trim(filename)) + endif + exit + endif + enddo + deallocate(dim_names) + + if (present(timelevel) .and. (time_dim < 0) .and. is_root_PE()) & + call MOM_error(WARNING, trim(err_header)//"time level specified, but the variable "//& + trim(var_to_read)//" does not have an unlimited dimension in "//trim(filename)) + if ((.not.present(timelevel)) .and. (time_dim > 0) .and. is_root_PE()) & + call MOM_error(WARNING, trim(err_header)//"The variable "//trim(var_to_read)//& + " has an unlimited dimension in "//trim(filename)//" but no time level is specified.") + if (present(has_time_dim)) has_time_dim = (time_dim > 0) + endif + + ! Registering the variable axes essentially just specifies the discrete position of this variable. + call MOM_register_variable_axes(fileobj, var_to_read, filename, position) + +end subroutine prepare_to_read_var + +!> register axes associated with a variable from a domain-decomposed netCDF file +subroutine MOM_register_variable_axes(fileObj, variableName, filename, position) + type(FmsNetcdfDomainFile_t), intent(inout) :: fileObj !< Handle to an open FMS2 netCDF file object + character(len=*), intent(in) :: variableName !< name of the variable + character(len=*), intent(in) :: filename !< The name of the file to read + integer, optional, intent(in) :: position !< A flag indicating where this data is discretized + + ! Local variables + character(len=256), allocatable, dimension(:) :: dim_names ! variable dimension names + integer, allocatable, dimension(:) :: dimSizes ! variable dimension sizes + logical, allocatable, dimension(:) :: is_x ! Is this a (likely domain-decomposed) x-axis + logical, allocatable, dimension(:) :: is_y ! Is this a (likely domain-decomposed) y-axis + logical, allocatable, dimension(:) :: is_t ! Is this a time axis or another unlimited axis + integer :: ndims ! number of dimensions + integer :: xPos, yPos ! Discrete positions for x and y axes. Default is CENTER + integer :: i + + xPos = CENTER ; yPos = CENTER + if (present(position)) then + if ((position == CORNER) .or. (position == EAST_FACE)) xPos = EAST_FACE + if ((position == CORNER) .or. (position == NORTH_FACE)) yPos = NORTH_FACE + endif + + ! get variable dimension names and lengths + ndims = get_variable_num_dimensions(fileObj, trim(variableName)) + allocate(dimSizes(ndims)) + allocate(dim_names(ndims)) + allocate(is_x(ndims)) ; is_x(:) = .false. + allocate(is_y(ndims)) ; is_y(:) = .false. + allocate(is_t(ndims)) ; is_t(:) = .false. + call get_variable_size(fileObj, trim(variableName), dimSizes) + call get_variable_dimension_names(fileObj, trim(variableName), dim_names) + call categorize_axes(fileObj, filename, ndims, dim_names, is_x, is_y, is_t) + + ! register the axes + do i=1,ndims + if ( .not.is_dimension_registered(fileobj, trim(dim_names(i))) ) then + if (is_x(i)) then + call register_axis(fileObj, trim(dim_names(i)), "x", domain_position=xPos) + elseif (is_y(i)) then + call register_axis(fileObj, trim(dim_names(i)), "y", domain_position=yPos) + else + call register_axis(fileObj, trim(dim_names(i)), dimSizes(i)) + endif + endif + enddo + + deallocate(dimSizes, dim_names, is_x, is_y, is_t) +end subroutine MOM_register_variable_axes + +!> Determine whether a variable's axes are associated with x-, y- or time-dimensions. Other +!! unlimited dimensions are also labeled as time axes for these purposes. +subroutine categorize_axes(fileObj, filename, ndims, dim_names, is_x, is_y, is_t) + type(FmsNetcdfDomainFile_t), intent(in) :: fileObj !< Handle to an open FMS2 netCDF file object + character(len=*), intent(in) :: filename !< The name of the file to read + integer, intent(in) :: ndims !< The number of dimensions associated with a variable + character(len=*), dimension(ndims), intent(in) :: dim_names !< Names of the dimensions associated with a variable + logical, dimension(ndims), intent(out) :: is_x !< Indicates if each dimension a (likely decomposed) x-axis + logical, dimension(ndims), intent(out) :: is_y !< Indicates if each dimension a (likely decomposed) y-axis + logical, dimension(ndims), intent(out) :: is_t !< Indicates if each dimension unlimited (usually time) axis + + ! Local variables + character(len=128) :: cartesian ! A flag indicating a Cartesian direction - usually a single character. + character(len=512) :: dim_list ! A concatenated list of dimension names. + character(len=128) :: units ! units corresponding to a specific variable dimension + logical :: x_found, y_found ! Indicate whether an x- or y- dimension have been found. + integer :: i + + x_found = .false. ; y_found = .false. + is_x(:) = .false. ; is_y(:) = .false. + do i=1,ndims + is_t(i) = is_dimension_unlimited(fileObj, trim(dim_names(i))) + ! First look for indicative variable attributes + if (.not.is_t(i)) then + if (variable_exists(fileobj, trim(dim_names(i)))) then + if (variable_att_exists(fileobj, trim(dim_names(i)), "cartesian_axis")) then + call get_variable_attribute(fileobj, trim(dim_names(i)), "cartesian_axis", cartesian) + cartesian = adjustl(cartesian) + if ((index(cartesian, "X") == 1) .or. (index(cartesian, "x") == 1)) is_x(i) = .true. + if ((index(cartesian, "Y") == 1) .or. (index(cartesian, "y") == 1)) is_y(i) = .true. + if ((index(cartesian, "T") == 1) .or. (index(cartesian, "t") == 1)) is_t(i) = .true. + endif + endif + endif + if (is_x(i)) x_found = .true. + if (is_y(i)) y_found = .true. + enddo + + if (.not.(x_found .and. y_found)) then + ! Next look for hints from axis names for uncharacterized axes + do i=1,ndims ; if (.not.(is_x(i) .or. is_y(i) .or. is_t(i))) then + call categorize_axis_from_name(dim_names(i), is_x(i), is_y(i)) + if (is_x(i)) x_found = .true. + if (is_y(i)) y_found = .true. + endif ; enddo + endif + + if (.not.(x_found .and. y_found)) then + ! Look for hints from CF-compliant axis units for uncharacterized axes + do i=1,ndims ; if (.not.(is_x(i) .or. is_y(i) .or. is_t(i))) then + call get_variable_units(fileobj, trim(dim_names(i)), units) + call categorize_axis_from_units(units, is_x(i), is_y(i)) + if (is_x(i)) x_found = .true. + if (is_y(i)) y_found = .true. + endif ; enddo + endif + + if (.not.(x_found .and. y_found) .and. ((ndims>2) .or. ((ndims==2) .and. .not.is_t(ndims)))) then + ! This is a case where one would expect to find x-and y-dimensions, but none have been found. + if (is_root_pe()) then + dim_list = trim(dim_names(1))//", "//trim(dim_names(2)) + do i=3,ndims ; dim_list = trim(dim_list)//", "//trim(dim_names(i)) ; enddo + call MOM_error(WARNING, "categorize_axes: Failed to identify x- and y- axes in the axis list ("//& + trim(dim_list)//") of a variable being read from "//trim(filename)) + endif + endif + +end subroutine categorize_axes + +!> Determine whether an axis is associated with the x- or y-directions based on a comparison of +!! its units with CF-compliant variants of latitude or longitude units. +subroutine categorize_axis_from_units(unit_string, is_x, is_y) + character(len=*), intent(in) :: unit_string !< string of units + logical, intent(out) :: is_x !< Indicates if the axis units are associated with an x-direction axis + logical, intent(out) :: is_y !< Indicates if the axis units are associated with an y-direction axis + + is_x = .false. ; is_y = .false. + select case (lowercase(trim(unit_string))) + case ("degrees_north"); is_y = .true. + case ("degree_north") ; is_y = .true. + case ("degrees_n") ; is_y = .true. + case ("degree_n") ; is_y = .true. + case ("degreen") ; is_y = .true. + case ("degreesn") ; is_y = .true. + case ("degrees_east") ; is_x = .true. + case ("degree_east") ; is_x = .true. + case ("degreese") ; is_x = .true. + case ("degreee") ; is_x = .true. + case ("degree_e") ; is_x = .true. + case ("degrees_e") ; is_x = .true. + case default ; is_x = .false. ; is_y = .false. + end select + +end subroutine categorize_axis_from_units + +!> Tries to determine whether the axis name is commonly associated with an x- or y- axis. This +!! approach is fragile and unreliable, but it a backup to reading a CARTESIAN file attribute. +subroutine categorize_axis_from_name(dimname, is_x, is_y) + character(len=*), intent(in) :: dimname !< A dimension name + logical, intent(out) :: is_x !< Indicates if the axis name is associated with an x-direction axis + logical, intent(out) :: is_y !< Indicates if the axis name is associated with an y-direction axis + + is_x = .false. ; is_y = .false. + select case(trim(lowercase(dimname))) + case ("grid_x_t") ; is_x = .true. + case ("nx") ; is_x = .true. + case ("nxp") ; is_x = .true. + case ("longitude") ; is_x = .true. + case ("long") ; is_x = .true. + case ("lon") ; is_x = .true. + case ("lonh") ; is_x = .true. + case ("lonq") ; is_x = .true. + case ("xh") ; is_x = .true. + case ("xq") ; is_x = .true. + case ("i") ; is_x = .true. + + case ("grid_y_t") ; is_y = .true. + case ("ny") ; is_y = .true. + case ("nyp") ; is_y = .true. + case ("latitude") ; is_y = .true. + case ("lat") ; is_y = .true. + case ("lath") ; is_y = .true. + case ("latq") ; is_y = .true. + case ("yh") ; is_y = .true. + case ("yq") ; is_y = .true. + case ("j") ; is_y = .true. + + case default ; is_x = .false. ; is_y = .false. + end select + +end subroutine categorize_axis_from_name + + !> Write a 4d field to an output file. subroutine write_field_4d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, fill_value) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value - call mpp_write(IO_handle%unit, field_md, MOM_domain%mpp_domain, field, tstamp=tstamp, & - tile_count=tile_count, default_data=fill_value) + + ! Local variables + integer :: time_index + + if (IO_handle%FMS2_file .and. present(tstamp)) then + time_index = write_time_if_later(IO_handle, tstamp) + call write_data(IO_handle%fileobj, trim(field_md%name), field, unlim_dim_level=time_index) + elseif (IO_handle%FMS2_file) then + call write_data(IO_handle%fileobj, trim(field_md%name), field) + else + call mpp_write(IO_handle%unit, field_md%FT, MOM_domain%mpp_domain, field, tstamp=tstamp, & + tile_count=tile_count, default_data=fill_value) + endif end subroutine write_field_4d !> Write a 3d field to an output file. subroutine write_field_3d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, fill_value) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value - call mpp_write(IO_handle%unit, field_md, MOM_domain%mpp_domain, field, tstamp=tstamp, & + ! Local variables + integer :: time_index + + if (IO_handle%FMS2_file .and. present(tstamp)) then + time_index = write_time_if_later(IO_handle, tstamp) + call write_data(IO_handle%fileobj, trim(field_md%name), field, unlim_dim_level=time_index) + elseif (IO_handle%FMS2_file) then + call write_data(IO_handle%fileobj, trim(field_md%name), field) + else + call mpp_write(IO_handle%unit, field_md%FT, MOM_domain%mpp_domain, field, tstamp=tstamp, & tile_count=tile_count, default_data=fill_value) + endif end subroutine write_field_3d !> Write a 2d field to an output file. subroutine write_field_2d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, fill_value) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:), intent(inout) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field integer, optional, intent(in) :: tile_count !< PEs per tile (default: 1) real, optional, intent(in) :: fill_value !< Missing data fill value - call mpp_write(IO_handle%unit, field_md, MOM_domain%mpp_domain, field, tstamp=tstamp, & + ! Local variables + integer :: time_index + + if (IO_handle%FMS2_file .and. present(tstamp)) then + time_index = write_time_if_later(IO_handle, tstamp) + call write_data(IO_handle%fileobj, trim(field_md%name), field, unlim_dim_level=time_index) + elseif (IO_handle%FMS2_file) then + call write_data(IO_handle%fileobj, trim(field_md%name), field) + else + call mpp_write(IO_handle%unit, field_md%FT, MOM_domain%mpp_domain, field, tstamp=tstamp, & tile_count=tile_count, default_data=fill_value) + endif end subroutine write_field_2d !> Write a 1d field to an output file. subroutine write_field_1d(IO_handle, field_md, field, tstamp) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, dimension(:), intent(in) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field + + ! Local variables + integer :: time_index - call mpp_write(IO_handle%unit, field_md, field, tstamp=tstamp) + if (IO_handle%FMS2_file .and. present(tstamp)) then + time_index = write_time_if_later(IO_handle, tstamp) + call write_data(IO_handle%fileobj, trim(field_md%name), field, unlim_dim_level=time_index) + elseif (IO_handle%FMS2_file) then + call write_data(IO_handle%fileobj, trim(field_md%name), field) + else + call mpp_write(IO_handle%unit, field_md%FT, field, tstamp=tstamp) + endif end subroutine write_field_1d !> Write a 0d field to an output file. subroutine write_field_0d(IO_handle, field_md, field, tstamp) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, intent(in) :: field !< Field to write - real, optional, intent(in) :: tstamp !< Model timestamp + real, optional, intent(in) :: tstamp !< Model time of this field + + ! Local variables + integer :: time_index - call mpp_write(IO_handle%unit, field_md, field, tstamp=tstamp) + if (IO_handle%FMS2_file .and. present(tstamp)) then + time_index = write_time_if_later(IO_handle, tstamp) + call write_data(IO_handle%fileobj, trim(field_md%name), field, unlim_dim_level=time_index) + elseif (IO_handle%FMS2_file) then + call write_data(IO_handle%fileobj, trim(field_md%name), field) + else + call mpp_write(IO_handle%unit, field_md%FT, field, tstamp=tstamp) + endif end subroutine write_field_0d +!> Returns the integer time index for a write in this file, also writing the time variable to +!! the file if this time is later than what is already in the file. +integer function write_time_if_later(IO_handle, field_time) + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing + real, intent(in) :: field_time !< Model time of this field + + ! Local variables + character(len=256) :: dim_unlim_name ! name of the unlimited dimension in the file + + if ((field_time > IO_handle%file_time) .or. (IO_handle%num_times == 0)) then + IO_handle%file_time = field_time + IO_handle%num_times = IO_handle%num_times + 1 + if (IO_handle%FMS2_file) then + call get_unlimited_dimension_name(IO_handle%fileobj, dim_unlim_name) + call write_data(IO_handle%fileobj, trim(dim_unlim_name), (/field_time/), & + corner=(/IO_handle%num_times/), edge_lengths=(/1/)) + endif + endif + + write_time_if_later = IO_handle%num_times +end function write_time_if_later + !> Write the data for an axis subroutine MOM_write_axis(IO_handle, axis) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(axistype), intent(in) :: axis !< An axis type variable with information to write - call mpp_write(IO_handle%unit, axis) + integer :: is, ie + + if (IO_handle%FMS2_file) then + if (axis%domain_decomposed) then + ! FMS2 does not domain-decompose 1d arrays, so we explicitly slice it + call get_global_io_domain_indices(IO_handle%fileobj, trim(axis%name), is, ie) + call write_data(IO_handle%fileobj, trim(axis%name), axis%ax_data(is:ie)) + else + call write_data(IO_handle%fileobj, trim(axis%name), axis%ax_data) + endif + else + call mpp_write(IO_handle%unit, axis%AT) + endif end subroutine MOM_write_axis !> Store information about an axis in a previously defined axistype and write this !! information to the file indicated by unit. -subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian, sense, domain, data, calendar) +subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian, sense, domain, data, & + edge_axis, calendar) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(axistype), intent(inout) :: axis !< The axistype where this information is stored. character(len=*), intent(in) :: name !< The name in the file of this axis @@ -959,28 +1762,115 @@ subroutine write_metadata_axis(IO_handle, axis, name, units, longname, cartesian !! -1 if they increase downward. type(domain1D), optional, intent(in) :: domain !< The domain decomposion for this axis real, dimension(:), optional, intent(in) :: data !< The coordinate values of the points on this axis + logical, optional, intent(in) :: edge_axis !< If true, this axis marks an edge of the tracer cells character(len=*), optional, intent(in) :: calendar !< The name of the calendar used with a time axis - call mpp_write_meta(IO_handle%unit, axis, name, units, longname, cartesian=cartesian, sense=sense, & - domain=domain, data=data, calendar=calendar) + character(len=:), allocatable :: cart ! A left-adjusted and trimmed copy of cartesian + logical :: is_x, is_y, is_t ! If true, this is a domain-decomposed axis in one of the directions. + integer :: position ! A flag indicating the axis staggering position. + integer :: i, isc, iec, global_size + + if (IO_handle%FMS2_file) then + if (is_dimension_registered(IO_handle%fileobj, trim(name))) then + call MOM_error(FATAL, "write_metadata_axis was called more than once for axis "//trim(name)//& + " in file "//trim(IO_handle%filename)) + return + endif + endif + + axis%name = trim(name) + if (present(data) .and. allocated(axis%ax_data)) call MOM_error(FATAL, & + "Data is already allocated in a call to write_metadata_axis for axis "//& + trim(name)//" in file "//trim(IO_handle%filename)) + + if (IO_handle%FMS2_file) then + is_x = .false. ; is_y = .false. ; is_t = .false. + position = CENTER + if (present(cartesian)) then + cart = trim(adjustl(cartesian)) + if ((index(cart, "X") == 1) .or. (index(cart, "x") == 1)) is_x = .true. + if ((index(cart, "Y") == 1) .or. (index(cart, "y") == 1)) is_y = .true. + if ((index(cart, "T") == 1) .or. (index(cart, "t") == 1)) is_t = .true. + endif + + ! For now, we assume that all horizontal axes are domain-decomposed. + if (is_x .or. is_y) & + axis%domain_decomposed = .true. + + if (is_x) then + if (present(edge_axis)) then ; if (edge_axis) position = EAST_FACE ; endif + call register_axis(IO_handle%fileobj, trim(name), 'x', domain_position=position) + elseif (is_y) then + if (present(edge_axis)) then ; if (edge_axis) position = NORTH_FACE ; endif + call register_axis(IO_handle%fileobj, trim(name), 'y', domain_position=position) + elseif (is_t .and. .not.present(data)) then + ! This is the unlimited (time) dimension. + call register_axis(IO_handle%fileobj, trim(name), unlimited) + else + if (.not.present(data)) call MOM_error(FATAL,"MOM_io:register_diagnostic_axis: "//& + "An axis_length argument is required to register the axis "//trim(name)) + call register_axis(IO_handle%fileobj, trim(name), size(data)) + endif + + if (present(data)) then + ! With FMS2, the data for the axis labels has to match the computational domain on this PE. + if (present(domain)) then + ! The commented-out code on the next ~11 lines runs but there is missing data in the output file + ! call mpp_get_compute_domain(domain, isc, iec) + ! call mpp_get_global_domain(domain, size=global_size) + ! if (size(data) == global_size) then + ! allocate(axis%ax_data(iec+1-isc)) ; axis%ax_data(:) = data(isc:iec) + ! ! A simpler set of labels: do i=1,iec-isc ; axis%ax_data(i) = real(isc + i) - 1.0 ; enddo + ! elseif (size(data) == global_size+1) then + ! ! This is an edge axis. Note the effective SW indexing convention here. + ! allocate(axis%ax_data(iec+2-isc)) ; axis%ax_data(:) = data(isc:iec+1) + ! ! A simpler set of labels: do i=1,iec+1-isc ; axis%ax_data(i) = real(isc + i) - 1.5 ; enddo + ! else + ! call MOM_error(FATAL, "Unexpected size of data for "//trim(name)//" in write_metadata_axis.") + ! endif + + ! This works for a simple 1x1 IO layout, but gives errors for nontrivial IO layouts + allocate(axis%ax_data(size(data))) ; axis%ax_data(:) = data(:) + + else ! Store the entire array of axis labels. + allocate(axis%ax_data(size(data))) ; axis%ax_data(:) = data(:) + endif + endif + + + ! Now create the variable that describes this axis. + call register_field(IO_handle%fileobj, trim(name), "double", dimensions=(/name/)) + if (len_trim(longname) > 0) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'long_name', & + trim(longname), len_trim(longname)) + if (len_trim(units) > 0) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'units', & + trim(units), len_trim(units)) + if (present(cartesian)) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'cartesian_axis', & + trim(cartesian), len_trim(cartesian)) + if (present(sense)) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'sense', sense) + else + if (present(data)) then + allocate(axis%ax_data(size(data))) ; axis%ax_data(:) = data(:) + endif + + call mpp_write_meta(IO_handle%unit, axis%AT, name, units, longname, cartesian=cartesian, sense=sense, & + domain=domain, data=data, calendar=calendar) + endif end subroutine write_metadata_axis !> Store information about an output variable in a previously defined fieldtype and write this !! information to the file indicated by unit. subroutine write_metadata_field(IO_handle, field, axes, name, units, longname, & - min, max, fill, scale, add, pack, standard_name, checksum) + pack, standard_name, checksum) type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(inout) :: field !< The fieldtype where this information is stored type(axistype), dimension(:), intent(in) :: axes !< Handles for the axis used for this variable character(len=*), intent(in) :: name !< The name in the file of this variable character(len=*), intent(in) :: units !< The units of this variable character(len=*), intent(in) :: longname !< The long description of this variable - real, optional, intent(in) :: min !< The minimum valid value for this variable - real, optional, intent(in) :: max !< The maximum valid value for this variable - real, optional, intent(in) :: fill !< Missing data fill value - real, optional, intent(in) :: scale !< An multiplicative factor by which to scale - !! the variable before output - real, optional, intent(in) :: add !< An offset to add to the variable before output integer, optional, intent(in) :: pack !< A precision reduction factor with which the !! variable. The default, 1, has no reduction, !! but 2 is not uncommon. @@ -988,9 +1878,45 @@ subroutine write_metadata_field(IO_handle, field, axes, name, units, longname, & integer(kind=int64), dimension(:), & optional, intent(in) :: checksum !< Checksum values that can be used to verify reads. + ! Local variables + character(len=256), dimension(size(axes)) :: dim_names ! The names of the dimensions + type(mpp_axistype), dimension(size(axes)) :: mpp_axes ! The array of mpp_axistypes for this variable + character(len=16) :: prec_string ! A string specifying the precision with which to save this variable + character(len=64) :: checksum_string ! checksum character array created from checksum argument + integer :: i, ndims + + ndims = size(axes) + if (IO_handle%FMS2_file) then + do i=1,ndims ; dim_names(i) = trim(axes(i)%name) ; enddo + prec_string = "double" ; if (present(pack)) then ; if (pack > 1) prec_string = "float" ; endif + call register_field(IO_handle%fileobj, trim(name), trim(prec_string), dimensions=dim_names) + if (len_trim(longname) > 0) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'long_name', & + trim(longname), len_trim(longname)) + if (len_trim(units) > 0) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'units', & + trim(units), len_trim(units)) + if (present(standard_name)) & + call register_variable_attribute(IO_handle%fileobj, trim(name), 'standard_name', & + trim(standard_name), len_trim(standard_name)) + if (present(checksum)) then + write (checksum_string,'(Z16)') checksum(1) ! Z16 is the hexadecimal format code + call register_variable_attribute(IO_handle%fileobj, trim(name), "checksum", & + trim(checksum_string), len_trim(checksum_string)) + endif + else + do i=1,ndims ; mpp_axes(i) = axes(i)%AT ; enddo + call mpp_write_meta(IO_handle%unit, field%FT, mpp_axes, name, units, longname, & + pack=pack, standard_name=standard_name, checksum=checksum) + ! unused opt. args: min=min, max=max, fill=fill, scale=scale, add=add, & + endif - call mpp_write_meta(IO_handle%unit, field, axes, name, units, longname, min=min, max=max, & - fill=fill, scale=scale, add=add, pack=pack, standard_name=standard_name, checksum=checksum) + ! Store information in the field-type, regardless of which interfaces are used. + field%name = trim(name) + field%longname = trim(longname) + field%units = trim(units) + field%chksum_read = -1 + field%valid_chksum = .false. end subroutine write_metadata_field diff --git a/config_src/infra/FMS2/MOM_read_data_fms2.F90 b/config_src/infra/FMS2/MOM_read_data_fms2.F90 deleted file mode 100644 index 83a10e7e30..0000000000 --- a/config_src/infra/FMS2/MOM_read_data_fms2.F90 +++ /dev/null @@ -1,555 +0,0 @@ -!> This module contains routines that wrap the fms2 read_data calls -module MOM_read_data_fms2 - -! This file is part of MOM6. See LICENSE.md for the license. -use MOM_error_infra, only : MOM_error=>MOM_err, NOTE, FATAL, WARNING, is_root_PE -use MOM_domain_infra, only : MOM_domain_type, AGRID, BGRID_NE, CGRID_NE -use MOM_domain_infra, only : domain2d, CENTER, CORNER, NORTH_FACE, EAST_FACE -use MOM_string_functions, only : lowercase -use fms2_io_mod, only : FmsNetcdfDomainFile_t, FmsNetcdfFile_t -use fms2_io_mod, only : fms2_open_file => open_file, fms2_close_file => close_file -use fms2_io_mod, only : get_num_variables, get_variable_names, check_if_open -use fms2_io_mod, only : read_data, variable_exists, get_variable_size, get_variable_units -use fms2_io_mod, only : get_variable_attribute, attribute_exists => variable_att_exists -use fms2_io_mod, only : get_variable_num_dimensions, get_variable_dimension_names -use fms2_io_mod, only : is_dimension_unlimited, get_dimension_size -use fms2_io_mod, only : is_dimension_registered, register_axis - -implicit none ; private - -public prepare_to_read_var -! public MOM_read_data_scalar, MOM_read_data_2d_noDD, MOM_read_data_1d_noDD - -contains - -!> Find the case-insensitive name match with a variable in a domain-decomposed file-set -!! opening the file(s) as necessary, prepare FMS2 to read this variable, and return some -!! information needed to call read_data correctly for this variable and file. -subroutine prepare_to_read_var(fileobj, fieldname, domain, err_header, filename, var_to_read, & - has_time_dim, timelevel, position) - type(FmsNetcdfDomainFile_t), intent(inout) :: fileobj !< A handle to an FMS2 file object, that - !! will be opened if necessary - character(len=*), intent(in) :: fieldname !< The variable name to seek in the file - type(MOM_domain_type), intent(in) :: domain !< MOM domain attribute with the mpp_domain decomposition - character(len=*), intent(in) :: err_header !< A descriptive prefix for error messages - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(out) :: var_to_read !< The variable name to read from the file - logical, optional, intent(out) :: has_time_dim !< Indicates whether fieldname has a time dimension - integer, optional, intent(in) :: timelevel !< A time level to read - integer, optional, intent(in) :: position !< A flag indicating where this variable is discretized - - ! Local variables - logical :: file_open_success !.true. if call to open_file is successful - logical :: variable_found ! Is a case-insensitive version of the variable found in the netCDF file? - character(len=96), allocatable, dimension(:) :: var_names !< array for names of variables in a netCDF - !! file opened to read - character(len=96), allocatable :: dim_names(:) ! variable dimension names - integer :: nvars ! The number of variables in the file. - integer :: i, dim_unlim_size, num_var_dims, time_dim - - ! Open the file if necessary - if (.not.(check_if_open(fileobj))) then - file_open_success = fms2_open_file(fileobj, filename, "read", domain%mpp_domain, is_restart=.false.) - if (.not.file_open_success) call MOM_error(FATAL, trim(err_header)//" failed to open "//trim(filename)) - endif - - ! Search for the variable in the file, looking for the case-sensitive name first. - if (variable_exists(fileobj, trim(fieldname))) then - var_to_read = trim(fieldname) - variable_found = .true. - else ! Look for case-insensitive variable name matches. - var_to_read = "" - variable_found = .false. - - nvars = get_num_variables(fileobj) - if (nvars < 1) call MOM_error(FATAL, "nvars is less than 1 for file "//trim(filename)) - allocate(var_names(nvars)) - call get_variable_names(fileobj, var_names) - - do i=1,nvars - if (lowercase(trim(var_names(i))) == lowercase(trim(fieldname))) then - variable_found = .true. - var_to_read = trim(var_names(i)) - exit - endif - enddo - if (.not.(variable_found)) & - call MOM_error(FATAL, trim(err_header)//trim(fieldname)//" not found in "//trim(filename)) - deallocate(var_names) - endif - - ! FMS2 can not handle a timelevel argument if the variable does not have one in the file, - ! so some error checking and logic are required. - if (present(has_time_dim) .or. present(timelevel)) then - time_dim = -1 - - num_var_dims = get_variable_num_dimensions(fileobj, trim(var_to_read)) - allocate(dim_names(num_var_dims)) ; dim_names(:) = "" - call get_variable_dimension_names(fileobj, trim(var_to_read), dim_names) - - do i=1,num_var_dims - if (is_dimension_unlimited(fileobj, dim_names(i))) then - time_dim = i - if (present(timelevel)) then - call get_dimension_size(fileobj, dim_names(i), dim_unlim_size) - if ((timelevel > dim_unlim_size) .and. is_root_PE()) call MOM_error(FATAL, & - trim(err_header)//"Attempting to read a time level of "//trim(var_to_read)//& - " that exceeds the size of the time dimension in "//trim(filename)) - endif - exit - endif - enddo - deallocate(dim_names) - - if (present(timelevel) .and. (time_dim < 0) .and. is_root_PE()) & - call MOM_error(WARNING, trim(err_header)//"time level specified, but the variable "//& - trim(var_to_read)//" does not have an unlimited dimension in "//trim(filename)) - if ((.not.present(timelevel)) .and. (time_dim > 0) .and. is_root_PE()) & - call MOM_error(WARNING, trim(err_header)//"The variable "//trim(var_to_read)//& - " has an unlimited dimension in "//trim(filename)//" but no time level is specified.") - if (present(has_time_dim)) has_time_dim = (time_dim > 0) - endif - - ! Registering the variable axes essentially just specifies the discrete position of this variable. - call MOM_register_variable_axes(fileobj, var_to_read, filename, position) - -end subroutine prepare_to_read_var - -!> register axes associated with a variable from a domain-decomposed netCDF file -!> @note The user must specify units for variables with longitude/x-axis and/or latitude/y-axis axes -!! to obtain the correct domain decomposition for the data buffer. -subroutine MOM_register_variable_axes(fileObj, variableName, filename, position) - type(FmsNetcdfDomainFile_t), intent(inout) :: fileObj !< Handle to an open FMS2 netCDF file object - character(len=*), intent(in) :: variableName !< name of the variable - character(len=*), intent(in) :: filename !< The name of the file to read - integer, optional, intent(in) :: position !< A flag indicating where this data is discretized - - ! Local variables - character(len=40) :: units ! units corresponding to a specific variable dimension - character(len=40), allocatable, dimension(:) :: dim_names ! variable dimension names - integer, allocatable, dimension(:) :: dimSizes ! variable dimension sizes - logical, allocatable, dimension(:) :: is_x ! Is this a (likely domain-decomposed) x-axis - logical, allocatable, dimension(:) :: is_y ! Is this a (likely domain-decomposed) y-axis - logical, allocatable, dimension(:) :: is_t ! Is this a time axis or another unlimited axis - integer :: ndims ! number of dimensions - integer :: i - integer :: xPos, yPos ! domain positions for x and y axes. Default is CENTER - - if (.not. check_if_open(fileObj)) call MOM_error(FATAL,"MOM_axis:register_variable_axes: The fileObj has "// & - "not been opened. Call fms2_open_file(fileObj,...) before "// & - "passing the fileObj argument to this function.") - xPos = CENTER ; yPos = CENTER - if (present(position)) then - if ((position == CORNER) .or. (position == EAST_FACE)) xPos = EAST_FACE - if ((position == CORNER) .or. (position == NORTH_FACE)) yPos = NORTH_FACE - endif - - ! get variable dimension names and lengths - ndims = get_variable_num_dimensions(fileObj, trim(variableName)) - allocate(dimSizes(ndims)) - allocate(dim_names(ndims)) - allocate(is_x(ndims)) ; is_x(:) = .false. - allocate(is_y(ndims)) ; is_y(:) = .false. - allocate(is_t(ndims)) ; is_t(:) = .false. - call get_variable_size(fileObj, trim(variableName), dimSizes) - call get_variable_dimension_names(fileObj, trim(variableName), dim_names) - call categorize_axes(fileObj, filename, ndims, dim_names, is_x, is_y, is_t) - - ! register the axes - do i=1,ndims - if ( .not.is_dimension_registered(fileobj, trim(dim_names(i))) ) then - if (is_x(i)) then - call register_axis(fileObj, trim(dim_names(i)), "x", domain_position=xPos) - elseif (is_y(i)) then - call register_axis(fileObj, trim(dim_names(i)), "y", domain_position=yPos) - else - call register_axis(fileObj, trim(dim_names(i)), dimSizes(i)) - endif - endif - enddo - - deallocate(dimSizes) - deallocate(dim_names) - deallocate(is_x, is_y, is_t) -end subroutine MOM_register_variable_axes - -!> Determine whether a variable's axes are associated with x-, y- or time-dimensions. Other -!! unlimited dimensions are also labeled as time axes for these purposes. -subroutine categorize_axes(fileObj, filename, ndims, dim_names, is_x, is_y, is_t) - type(FmsNetcdfDomainFile_t), intent(in) :: fileObj !< Handle to an open FMS2 netCDF file object - character(len=*), intent(in) :: filename !< The name of the file to read - integer, intent(in) :: ndims !< The number of dimensions associated with a variable - character(len=*), dimension(ndims), intent(in) :: dim_names !< Names of the dimensions associated with a variable - logical, dimension(ndims), intent(out) :: is_x !< Indicates if each dimension a (likely decomposed) x-axis - logical, dimension(ndims), intent(out) :: is_y !< Indicates if each dimension a (likely decomposed) y-axis - logical, dimension(ndims), intent(out) :: is_t !< Indicates if each dimension unlimited (usually time) axis - - integer :: i - character(len=256) :: cartesian ! A flag indicating a Cartesian direction - usually a single character. - character(len=512) :: dim_list ! A concatenated list of dimension names. - character(len=40) :: units ! units corresponding to a specific variable dimension - logical :: x_found, y_found ! Indicate whether an x- or y- dimension have been found. - - x_found = .false. ; y_found = .false. - is_x(:) = .false. ; is_y(:) = .false. - do i=1,ndims - is_t(i) = is_dimension_unlimited(fileObj, trim(dim_names(i))) - ! First look for indicative variable attributes - if (.not.is_t(i)) then - if (variable_exists(fileobj, trim(dim_names(i)))) then - if (attribute_exists(fileobj, trim(dim_names(i)), "cartesian_axis")) then - call get_variable_attribute(fileobj, trim(dim_names(i)), "cartesian_axis", cartesian) - cartesian = adjustl(cartesian) - if ((index(cartesian, "X") == 1) .or. (index(cartesian, "x") == 1)) is_x(i) = .true. - if ((index(cartesian, "Y") == 1) .or. (index(cartesian, "y") == 1)) is_y(i) = .true. - if ((index(cartesian, "T") == 1) .or. (index(cartesian, "t") == 1)) is_t(i) = .true. - ! if (is_root_pe() .and. is_x(i)) & - ! call MOM_error(NOTE, "X-dimension determined from cartesian_axis for "//trim(dim_names(i))) - ! if (is_root_pe() .and. is_y(i)) & - ! call MOM_error(NOTE, "Y-dimension determined from cartesian_axis for "//trim(dim_names(i))) - endif - endif - endif - if (is_x(i)) x_found = .true. - if (is_y(i)) y_found = .true. - enddo - - if (.not.(x_found .and. y_found)) then - ! Next look for hints from axis names for uncharacterized axes - do i=1,ndims ; if (.not.(is_x(i) .or. is_y(i) .or. is_t(i))) then - call categorize_axis_from_name(dim_names(i), is_x(i), is_y(i)) - if (is_x(i)) x_found = .true. - if (is_y(i)) y_found = .true. - endif ; enddo - endif - - if (.not.(x_found .and. y_found)) then - ! Look for hints from CF-compliant axis units for uncharacterized axes - do i=1,ndims ; if (.not.(is_x(i) .or. is_y(i) .or. is_t(i))) then - call get_variable_units(fileobj, trim(dim_names(i)), units) - call categorize_axis_from_units(units, is_x(i), is_y(i)) - if (is_x(i)) x_found = .true. - if (is_y(i)) y_found = .true. - endif ; enddo - endif - - if (.not.(x_found .and. y_found) .and. (ndims>2) .or. ((ndims==2) .and. .not.is_t(ndims))) then - ! This is a case where one would expect to find x-and y-dimensions, but none have been found. - if (is_root_pe()) then - dim_list = trim(dim_names(1))//", "//trim(dim_names(2)) - do i=3,ndims ; dim_list = trim(dim_list)//", "//trim(dim_names(i)) ; enddo - call MOM_error(WARNING, "categorize_axes: Failed to identify x- and y- axes in the axis list ("//& - trim(dim_list)//") of a variable being read from "//trim(filename)) - endif - endif - -end subroutine categorize_axes - -!> Determine whether an axis is associated with the x- or y-directions based on a comparison of -!! its units with CF-compliant variants of latitude or longitude units. -subroutine categorize_axis_from_units(unit_string, is_x, is_y) - character(len=*), intent(in) :: unit_string !< string of units - logical, intent(out) :: is_x !< Indicates if the axis units are associated with an x-direction axis - logical, intent(out) :: is_y !< Indicates if the axis units are associated with an y-direction axis - - is_x = .false. ; is_y = .false. - select case (lowercase(trim(unit_string))) - case ("degrees_north"); is_y = .true. - case ("degree_north") ; is_y = .true. - case ("degrees_n") ; is_y = .true. - case ("degree_n") ; is_y = .true. - case ("degreen") ; is_y = .true. - case ("degreesn") ; is_y = .true. - case ("degrees_east") ; is_x = .true. - case ("degree_east") ; is_x = .true. - case ("degreese") ; is_x = .true. - case ("degreee") ; is_x = .true. - case ("degree_e") ; is_x = .true. - case ("degrees_e") ; is_x = .true. - case default ; is_x = .false. ; is_y = .false. - end select - -end subroutine categorize_axis_from_units - -!> Tries to determine whether the axis name is commonly associated with an x- or y- axis. This -!! approach is fragile and unreliable, but it a backup to reading a CARTESIAN file attribute. -subroutine categorize_axis_from_name(dimname, is_x, is_y) - character(len=*), intent(in) :: dimname !< A dimension name - logical, intent(out) :: is_x !< Indicates if the axis name is associated with an x-direction axis - logical, intent(out) :: is_y !< Indicates if the axis name is associated with an y-direction axis - - is_x = .false. ; is_y = .false. - select case(trim(lowercase(dimname))) - case ("grid_x_t") ; is_x = .true. - case ("nx") ; is_x = .true. - case ("nxp") ; is_x = .true. - case ("longitude") ; is_x = .true. - case ("long") ; is_x = .true. - case ("lon") ; is_x = .true. - case ("lonh") ; is_x = .true. - case ("lonq") ; is_x = .true. - case ("xh") ; is_x = .true. - case ("xq") ; is_x = .true. - case ("i") ; is_x = .true. - - case ("grid_y_t") ; is_y = .true. - case ("ny") ; is_y = .true. - case ("nyp") ; is_y = .true. - case ("latitude") ; is_y = .true. - case ("lat") ; is_y = .true. - case ("lath") ; is_y = .true. - case ("latq") ; is_y = .true. - case ("yh") ; is_y = .true. - case ("yq") ; is_y = .true. - case ("j") ; is_y = .true. - - case default ; is_x = .false. ; is_y = .false. - end select - -end subroutine categorize_axis_from_name - -!===== Everything below this pertains to reading non-decomposed variables ===! -!===== using FMS2 interfaces will probably be discarded eventually. =========! - -!!> This routine calls the fms_io read_data subroutine to read a scalar (0-D) field named "fieldname" -!! from file "filename". -subroutine MOM_read_data_scalar(filename, fieldname, data, timelevel, scale, leave_file_open) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, intent(inout) :: data !< The variable to read from read_data - integer, optional, intent(in) :: timelevel !< time level to read - real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied by - logical, optional, intent(in) :: leave_file_open !< if .true., leave file open - - ! Local variables - type(FmsNetcdfFile_t) :: fileobj ! A handle to a simple netCDF file - logical :: close_the_file ! indicates whether to close the file after read_data is called. - character(len=96) :: var_to_read ! variable to read from the netcdf file - character(len=48) :: err_header ! A preamble for error messages - - err_header = "MOM_read_data_fms2:MOM_read_data_scalar: " - - ! Find the matching variable name in the file, opening it and reading metadata if necessary. - call find_varname_in_file(fileobj, fieldname, err_header, filename, var_to_read) - - ! read the data - if (present(timelevel)) then - call read_data(fileobj, trim(var_to_read), data, unlim_dim_level=timelevel) - else - call read_data(fileobj, trim(var_to_read), data) - endif - - ! Close the file, if necessary - close_the_file = .true. ; if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - if (close_the_file .and. check_if_open(fileobj)) call fms2_close_file(fileobj) - - ! Rescale the data that was read if necessary. - if (present(scale)) then ; if (scale /= 1.0) then - data = scale*data - endif ; endif - -end subroutine MOM_read_data_scalar - -!> This routine calls the fms_io read_data subroutine to read 1-D non-domain-decomposed data field named "fieldname" -!! from file "filename". The routine multiplies the data by "scale" if the optional argument is included in the call. -subroutine MOM_read_data_1d_noDD(filename, fieldname, data, start_index, & - edge_lengths, timelevel, scale, leave_file_open) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, dimension(:), intent(inout) :: data !< The 1-dimensional data array to pass to read_data - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is - !! the variable size - integer, optional, intent(in) :: timelevel !< time level to read - real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied by - logical, optional, intent(in) :: leave_file_open !< if .true., leave file open - - ! Local variables - type(FmsNetcdfFile_t) :: fileobj ! A handle to a simple netCDF file - logical :: close_the_file ! indicates whether to close the file after read_data is called. - integer :: time_dim ! The dimension position of a variables unlimited time axis, or -1 if it has none. - integer, parameter :: ndim = 1 ! The dimensionality of the array being read - integer, dimension(ndim) :: start, nread ! indices for first data value and number of values to read - character(len=96) :: var_to_read ! variable to read from the netcdf file - character(len=48) :: err_header ! A preamble for error messages - - err_header = "MOM_read_data_fms2:MOM_read_data_1d_noDD: " - - ! Find the matching case-insensitive variable name in the file, opening the file if necessary. - call find_varname_in_file(fileobj, fieldname, err_header, filename, var_to_read) - - ! set the start and nread values that will be passed as the read_data corner and edge_lengths arguments - start(:) = 1 ; if (present(start_index)) start(:) = start_index(:) - nread(:) = shape(data) ; if (present(edge_lengths)) nread(:) = edge_lengths(:) - - time_dim = -1 - if (present(timelevel)) then - time_dim = get_time_dim(fileobj, var_to_read, err_header, filename, timelevel) - if (time_dim == ndim) then ; nread(ndim) = 1 ; start(ndim) = timelevel ; endif - endif - - ! read the data - if (time_dim > 0) then - call read_data(fileobj, trim(var_to_read), data, corner=start, edge_lengths=nread, & - unlim_dim_level=timelevel) - else - call read_data(fileobj, trim(var_to_read), data, corner=start, edge_lengths=nread) - endif - - ! Close the file, if necessary - close_the_file = .true. ; if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - if (close_the_file .and. check_if_open(fileobj)) call fms2_close_file(fileobj) - - ! Rescale the data that was read if necessary. - if (present(scale)) then ; if (scale /= 1.0) then - data(:) = scale*data(:) - endif ; endif - -end subroutine MOM_read_data_1d_noDD - -!> This routine calls the fms_io read_data subroutine to read a 2-D non-domain-decomposed data field named "fieldname" -!! from file "filename". The routine multiplies the data by "scale" if the optional argument is included in the call. -subroutine MOM_read_data_2d_noDD(filename, fieldname, data, start_index, & - edge_lengths, timelevel, position, scale, leave_file_open) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, dimension(:,:), intent(inout) :: data !< The 2-dimensional data array to pass to read_data - integer, dimension(2), optional, intent(in) :: start_index !< starting indices of data buffer. Default is 1 - integer, dimension(2), optional, intent(in) :: edge_lengths !< number of data values to read in. - !! Default values are the variable dimension sizes - integer, optional, intent(in) :: timelevel !< time level to read - integer, optional, intent(in) :: position !< A flag indicating where this data is located - real, optional, intent(in) :: scale !< A scaling factor that the field is multiplied by - logical, optional, intent(in) :: leave_file_open !< if .true., leave file open - - ! Local variables - type(FmsNetcdfFile_t) :: fileobj ! A handle to a simple netCDF file - logical :: close_the_file ! indicates whether to close the file after read_data is called. - integer :: time_dim ! The dimension position of a variables unlimited time axis, or -1 if it has none. - integer, parameter :: ndim = 2 ! The dimensionality of the array being read - integer, dimension(ndim) :: start, nread ! indices for first data value and number of values to read - character(len=96) :: var_to_read ! variable to read from the netcdf file - character(len=48) :: err_header ! A preamble for error messages - - err_header = "MOM_read_data_fms2:MOM_read_data_2d_DD: " - - ! Find the matching case-insensitive variable name in the file, opening the file if necessary. - call find_varname_in_file(fileobj, fieldname, err_header, filename, var_to_read) - - ! set the start and nread values that will be passed as the read_data corner and edge_lengths arguments - start(:) = 1 ; if (present(start_index)) start(:) = start_index(:) - nread(:) = shape(data) ; if (present(edge_lengths)) nread(:) = edge_lengths(:) - - time_dim = -1 - if (present(timelevel)) then - time_dim = get_time_dim(fileobj, var_to_read, err_header, filename, timelevel) - if (time_dim == ndim) then ; nread(ndim) = 1 ; start(ndim) = timelevel ; endif - endif - - ! read the data - if (time_dim > 0) then - call read_data(fileobj, trim(var_to_read), data, corner=start, edge_lengths=nread, & - unlim_dim_level=timelevel) - else - call read_data(fileobj, trim(var_to_read), data, corner=start, edge_lengths=nread) - endif - - ! Close the file, if necessary - close_the_file = .true. ; if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - if (close_the_file .and. check_if_open(fileobj)) call fms2_close_file(fileobj) - - ! Rescale the data that was read if necessary. - if (present(scale)) then ; if (scale /= 1.0) then - data(:,:) = scale*data(:,:) - endif ; endif - -end subroutine MOM_read_data_2d_noDD - - -!> Find the case-sensitive name of the variable in a netCDF file with a case-insensitive name match. -subroutine find_varname_in_file(fileobj, fieldname, err_header, filename, var_to_read) - type(FmsNetcdfFile_t), intent(inout) :: fileobj !< A handle to a file object, that - !! will be opened if necessary - character(len=*), intent(in) :: fieldname !< The variable name to seek in the file - character(len=*), intent(in) :: err_header !< A descriptive prefix for error messages - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(out) :: var_to_read !< The variable name to read from the file - - ! Local variables - logical :: file_open_success !.true. if call to open_file is successful - logical :: variable_found ! Is a case-insensitive version of the variable found in the netCDF file? - character(len=96), allocatable, dimension(:) :: var_names !< array for names of variables in a netCDF - !! file opened to read - integer :: nvars ! The number of variables in the file. - integer :: i - - var_to_read = "" - - ! Open the file if necessary - if (.not.(check_if_open(fileobj))) then - file_open_success = fms2_open_file(fileobj, filename, "read", is_restart=.false.) - if (.not.file_open_success) call MOM_error(FATAL, trim(err_header)//" failed to open "//trim(filename)) - endif - - if (variable_exists(fileobj, fieldname)) then - var_to_read = fieldname - else - variable_found = .false. - nvars = get_num_variables(fileobj) - if (nvars < 1) call MOM_error(FATAL, "nvars is less than 1 for file "//trim(filename)) - allocate(var_names(nvars)) - call get_variable_names(fileobj, var_names) - - ! search for the variable in the file - do i=1,nvars - if (lowercase(trim(var_names(i))) == lowercase(trim(fieldname))) then - variable_found = .true. - var_to_read = trim(var_names(i)) - exit - endif - enddo - if (.not.(variable_found)) & - call MOM_error(FATAL, trim(err_header)//trim(fieldname)//" not found in "//trim(filename)) - deallocate(var_names) - endif - -end subroutine find_varname_in_file - -!> Return the number of the time dimension for a variable in an open non-domain-decomposed file, -!! or -1 if it has no time (or other unlimited) dimension. -integer function get_time_dim(fileobj, var_to_read, err_header, filename, timelevel) - type(FmsNetcdfFile_t), intent(in) :: fileobj !< A handle to an open file object - character(len=*), intent(in) :: var_to_read !< The variable name to read from the file - character(len=*), intent(in) :: err_header !< A descriptive prefix for error messages - character(len=*), intent(in) :: filename !< The name of the file to read - integer, optional, intent(in) :: timelevel !< A time level to read - - ! Local variables - integer :: i, dim_unlim_size, num_var_dims - character(len=96), allocatable :: dim_names(:) ! variable dimension names - - num_var_dims = get_variable_num_dimensions(fileobj, trim(var_to_read)) - allocate(dim_names(num_var_dims)) ; dim_names(:) = "" - call get_variable_dimension_names(fileobj, trim(var_to_read), dim_names) - - get_time_dim = -1 - do i=1,num_var_dims - if (is_dimension_unlimited(fileobj, dim_names(i))) then - get_time_dim = i - if (present(timelevel)) then - call get_dimension_size(fileobj, dim_names(i), dim_unlim_size) - if (timelevel > dim_unlim_size) call MOM_error(FATAL, trim(err_header)//& - "Attempting to read a time level of "//trim(var_to_read)//& - " that exceeds the size of "//trim(filename)) - endif - exit - endif - enddo - if (get_time_dim < 0) & - call MOM_error(WARNING, trim(err_header)//"time level specified, but the variable "//& - trim(var_to_read)//" does not have an unlimited dimension in "//trim(filename)) - deallocate(dim_names) - -end function get_time_dim - -end module MOM_read_data_fms2 diff --git a/config_src/infra/FMS2/MOM_write_field_fms2.F90 b/config_src/infra/FMS2/MOM_write_field_fms2.F90 deleted file mode 100644 index 24ba5ebb50..0000000000 --- a/config_src/infra/FMS2/MOM_write_field_fms2.F90 +++ /dev/null @@ -1,1503 +0,0 @@ -!> This module contains wrapper functions to write data to netcdf files -module MOM_write_field_fms2 - -! This file is part of MOM6. See LICENSE.md for the license. - -use MOM_axis, only : MOM_get_diagnostic_axis_data, MOM_register_diagnostic_axis -use MOM_axis, only : axis_data_type, get_time_index, get_var_dimension_metadata -use MOM_axis, only : get_time_units, convert_checksum_to_string -use MOM_coms_infra, only : PE_here, root_PE, num_PEs -use MOM_domain_infra, only : MOM_domain_type -use MOM_domain_infra, only : domain2d, CENTER, CORNER, NORTH_FACE, EAST_FACE -use MOM_error_handler, only : MOM_error, NOTE, FATAL, WARNING -use MOM_grid, only : ocean_grid_type -use MOM_dyn_horgrid, only : dyn_horgrid_type -use MOM_string_functions, only : lowercase, append_substring -use MOM_verticalGrid, only : verticalGrid_type - -use netcdf, only : nf90_max_name -! fms2_io -use fms2_io_mod, only : check_if_open, get_dimension_size -use fms2_io_mod, only : get_num_dimensions, get_num_variables, get_variable_names -use fms2_io_mod, only : get_unlimited_dimension_name, get_variable_dimension_names -use fms2_io_mod, only : get_variable_num_dimensions, get_variable_size, get_variable_units -use fms2_io_mod, only : get_variable_unlimited_dimension_index, is_dimension_unlimited -use fms2_io_mod, only : is_dimension_registered, register_axis -use fms2_io_mod, only : register_field, register_variable_attribute, fms2_open_file => open_file -use fms2_io_mod, only : fms2_close_file => close_file, write_data, variable_exists -use fms2_io_mod, only : FmsNetcdfDomainFile_t, FmsNetcdfFile_t, unlimited - -implicit none; private - -public write_field - -! CAUTION: The following variables are saved by default, and are only necessary for consecutive calls to -! write_field with the same file name. The user should ensure that fms2_close_file on -! the fileobj_write_field structures are called at every requisite time step at after the last -! variable is written to the file by omitting the optional leave_file_open argument, or setting it to .false. - -!> netCDF non-domain-decomposed file object returned by call to open_file in write_field calls -type(FmsNetcdfFile_t), private :: fileobj_write_field - -!> netCDF domain-decomposed file object returned by call to open_file in write_field calls -type(FmsNetcdfDomainFile_t), private :: fileobj_write_field_dd - -!> index of the time_level value that is written to netCDF file by the write_field routines -integer, private :: write_field_time_index - -!> interface to write data to a netcdf file generated by create_file -interface write_field - module procedure write_field_4d_DD - module procedure write_field_3d_DD - module procedure write_field_2d_DD - module procedure write_field_1d_DD - module procedure write_scalar - module procedure write_field_4d_noDD - module procedure write_field_3d_noDD - module procedure write_field_2d_noDD - module procedure write_field_1d_noDD -end interface - -contains -!> This function uses the fms_io function write_data to write a 1-D domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_1d_DD(filename, fieldname, data, mode, domain, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - type(MOM_domain_type),intent(in) :: domain !< MOM domain attribute with the mpp_domain decomposition - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is - !! the variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< if .true., leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: num_dims, substring_index - integer :: dim_unlim_size! size of the unlimited dimension - integer, dimension(1) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=1024) :: filename_temp - character(len=48), dimension(2) :: dim_names !< variable dimension names (or name, in the 1-D case); 1 extra - !! dimension in case appending along the time axis - integer, dimension(2) :: dim_lengths !< variable dimension lengths (or length, in the 1-D case) - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size=0 - dim_unlim_name="" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension is user-specified - ! and not assumed from the t_grid value - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - ! define the start and edge_length arguments - start(:) = 1 - nwrite(:) = dim_lengths(1) - if (present(start_index)) then - start(1) = max(1, start_index(1)) - endif - - if (present(edge_lengths)) then - nwrite(1) = max(dim_lengths(1),edge_lengths(1)) - endif - - if (.not.(check_if_open(fileobj_write_field_dd))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_1d_DD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field_dd, trim(filename_temp), lowercase(trim(mode)), & - domain%mpp_domain, is_restart=.false.) - ! register the diagnostic axis associated with the variable - call MOM_register_diagnostic_axis(fileobj_write_field_dd, trim(dim_names(1)), dim_lengths(1)) - endif - ! register and write the time_level - if (present(time_level)) then - if (.not. (variable_exists(fileobj_write_field_dd, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field_dd, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field_dd, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time_level if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), & - corner=(/write_field_time_index/), edge_lengths=(/1/)) - endif - endif - ! register the field if it is not already in the file - if (.not.(variable_exists(fileobj_write_field_dd, trim(fieldname)))) then - call register_field(fileobj_write_field_dd, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the variable - if (present(time_level)) then - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field_dd)) call fms2_close_file(fileobj_write_field_dd) - write_field_time_index = 0 - endif - -end subroutine write_field_1d_DD - -!> This function uses the fms_io function write_data to write a 2-D domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_2d_DD(filename, fieldname, data, mode, domain, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - type(MOM_domain_type),intent(in) :: domain !< MOM domain attribute with the mpp_domain decomposition - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is - !! the variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: i, is, ie, js, je, j, ndims, num_dims, substring_index - integer, allocatable, dimension(:) :: x_inds, y_inds - integer :: dim_unlim_size ! size of the unlimited dimension - integer :: file_dim_length - integer, dimension(2) :: start, nwrite ! indices for starting points and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(3) :: dim_names ! variable dimension names; 1 extra dimension in case appending - ! along the time axis - character(len=48), allocatable, dimension(:) :: file_dim_names - integer, dimension(3) :: dim_lengths ! variable dimension lengths - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_lengths(:) = 0 - dim_names(:) = "" - dim_unlim_size = 0 - dim_unlim_name = "" - ndims = 2 - num_dims = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension - ! is user-specified rather than derived from the t_grid value - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - start(:) = 1 - nwrite(:) = dim_lengths(1:ndims) - - if (present(start_index)) then - do i=1,ndims - start(i) = max(1,start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i),edge_lengths(i)) - enddo - endif - - if (.not.(check_if_open(fileobj_write_field_dd))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_2d_DD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field_dd, trim(filename_temp), lowercase(trim(mode)), & - domain%mpp_domain, is_restart=.false.) - endif - ! register the horizontal diagnostic axes associated with the variable - do i=1,num_dims - if (.not.(is_dimension_registered(fileobj_write_field_dd, trim(dim_names(i))))) & - call MOM_register_diagnostic_axis(fileobj_write_field_dd, trim(dim_names(i)), dim_lengths(i)) - enddo - ! register and write the time_level - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field_dd, dim_unlim_name) - call get_dimension_size(fileobj_write_field_dd, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field_dd, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field_dd, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field_dd, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time_level if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), & - corner=(/write_field_time_index/), edge_lengths=(/1/)) - endif - endif - ! register the variable if it is not already in the file - if (.not.(variable_exists(fileobj_write_field_dd, trim(fieldname)))) then - call register_field(fileobj_write_field_dd, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the variable - if (present(time_level)) then - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field_dd)) call fms2_close_file(fileobj_write_field_dd) - write_field_time_index=0 - if (allocated(file_dim_names)) deallocate(file_dim_names) - endif - -end subroutine write_field_2d_DD - -!> This function uses the fms_io function write_data to write a 3-D domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_3d_DD(filename, fieldname, data, mode, domain, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:,:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - type(MOM_domain_type),intent(in) :: domain !< MOM domain attribute with the mpp_domain decomposition - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is - !! the variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: i, is, ie, js, je, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(3) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(4) :: dim_names !< variable dimension names; 1 extra dimension in case appending - !! along the time axis - integer, dimension(4) :: dim_lengths !< variable dimension lengths - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size = 0 - dim_unlim_name = "" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension is user-specified - ! and not assumed from the t_grid value - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - ndims = 3 - start(:) = 1 - nwrite(:) = dim_lengths(1:3) - if (present(start_index)) then - do i=1,ndims - start(i) = max(1,start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i), edge_lengths(i)) - enddo - endif - - ! open the file - if (.not.(check_if_open(fileobj_write_field_dd))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_3d_DD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field_dd, trim(filename_temp), lowercase(trim(mode)), & - domain%mpp_domain, is_restart=.false.) - ! register the horizontal and vertical diagnostic axes associated with the variable - do i=1,ndims - call MOM_register_diagnostic_axis(fileobj_write_field_dd, trim(dim_names(i)), dim_lengths(i)) - enddo - endif - ! register and write the time_level - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field_dd ,dim_unlim_name) - call get_dimension_size(fileobj_write_field_dd, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field_dd, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field_dd, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field_dd, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time_level if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size ) & - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), & - corner=(/write_field_time_index/), edge_lengths=(/1/)) - endif - endif - ! register the field if it is not already in the file - if (.not.(variable_exists(fileobj_write_field_dd, trim(fieldname)))) then - call register_field(fileobj_write_field_dd, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the data - if (present(time_level)) then - call get_dimension_size(fileobj_write_field_dd, trim(dim_unlim_name), dim_unlim_size) - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field_dd)) call fms2_close_file(fileobj_write_field_dd) - write_field_time_index=0 - endif - - -end subroutine write_field_3d_DD - -!> This function uses the fms_io function write_data to write a 4-D domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_4d_DD(filename, fieldname, data, mode, domain, hor_grid, z_grid, t_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:,:,:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - type(MOM_domain_type),intent(in) :: domain !< MOM domain attribute with the mpp_domain decomposition - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - character(len=*), intent(in) :: t_grid !< time descriptor - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is - !! the variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - real :: file_time ! most recent time currently written to file - integer :: i, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(4) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(4) :: dim_names ! variable dimension names - integer, dimension(4) :: dim_lengths ! variable dimension lengths - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - num_dims = 0 - dim_unlim_size = 0 - dim_unlim_name = "" - dim_names(:) = "" - dim_lengths(:) = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - ndims = 4 - start(:) = 1 - nwrite(:) = dim_lengths(:) - if (present(start_index)) then - do i=1,ndims - start(i) = max(1,start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i), edge_lengths(i)) - enddo - endif - - ! open the file - if (.not.(check_if_open(fileobj_write_field_dd))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_4d_DD:mode argument must be write, overwrite, or append") - ! get the index of the corresponding time_level the first time the file is opened - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field_dd, trim(filename_temp), lowercase(trim(mode)), & - domain%mpp_domain, is_restart=.false.) - ! register the horizontal and vertical diagnostic axes associated with the variable - do i=1,ndims - call MOM_register_diagnostic_axis(fileobj_write_field_dd, trim(dim_names(i)), dim_lengths(i)) - enddo - endif - ! register the time dimension and write the time_level - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field_dd, dim_unlim_name) - call get_dimension_size(fileobj_write_field_dd, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field_dd, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field_dd, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field_dd, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time_level if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field_dd, trim(dim_unlim_name), (/time_level/), & - corner=(/write_field_time_index/), edge_lengths=(/1/)) - endif - endif - ! register the variable if it is not already in the file - if (.not.(variable_exists(fileobj_write_field_dd, trim(fieldname)))) then - call register_field(fileobj_write_field_dd, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field_dd, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the data - if (present(time_level)) then - call get_dimension_size(fileobj_write_field_dd, trim(dim_unlim_name), dim_unlim_size) - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field_dd, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field_dd)) call fms2_close_file(fileobj_write_field_dd) - write_field_time_index=0 - endif - -end subroutine write_field_4d_DD - -!> This routine uses the fms_io function write_data to write a scalar variable named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_scalar(filename, fieldname, data, mode, time_level, time_units, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=48), dimension(1) :: dim_names ! variable dimension names - integer :: i, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - real, allocatable, dimension(:) :: file_times - integer, dimension(1) :: dim_lengths ! variable dimension lengths - integer, allocatable, dimension(:) :: pelist ! list of pes associated with the netCDF file - - dim_unlim_size = 0 - dim_unlim_name= "" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - - if (.not.(check_if_open(fileobj_write_field))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_scaler:mode argument must be write, overwrite, or append") - ! get the index of the corresponding time_level the first time the file is opened - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! get the pes associated with the file. - !>\note this is required so that only pe(1) is identified as the root pe to create the file - !! Otherwise, multiple pes may try to open the file in write (NC_NOCLOBBER) mode, leading to failure - if (.not.(allocated(pelist))) then - allocate(pelist(num_PEs())) - pelist(:) = 0 - do i=1,size(pelist) - pelist(i) = i-1 - enddo - endif - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field, trim(filename_temp), trim(mode), is_restart=.false., & - pelist=pelist) - endif - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field, dim_unlim_name) - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - ! write the time value if it is not already written to the file - if (.not.(variable_exists(fileobj_write_field, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/)) - else - ! write the next time value if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/), & - edge_lengths=(/1/)) - endif - endif - ! register the variable - if (.not.(variable_exists(fileobj_write_field, trim(fieldname)))) then - if (present(time_level)) then - call register_field(fileobj_write_field, trim(fieldname), "double", dimensions=(/trim(dim_unlim_name)/)) - else - call register_field(fileobj_write_field, trim(fieldname), "double") - endif - if (present(units)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'long_name', trim(longname)) - endif - ! write the data - if (present(time_level)) then - call write_data(fileobj_write_field, trim(fieldname), data, unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field, trim(fieldname), data) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field)) call fms2_close_file(fileobj_write_field) - if (allocated(pelist)) deallocate(pelist) - write_field_time_index=0 - endif -end subroutine write_scalar - -!> This function uses the fms_io function write_data to write a 1-D non-domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_1d_noDD(filename, fieldname, data, mode, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(1), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(1), optional, intent(in) :: edge_lengths !< number of data values to read in; default is the - !! variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< if .true., leave file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - integer :: i, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(1) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(2) :: dim_names ! variable dimension names (up to 2 if appended at time level) - integer, dimension(2) :: dim_lengths ! variable dimension lengths - integer, allocatable, dimension(:) :: pelist ! list of pes associated with the netCDF file - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size = 0 - dim_unlim_name= "Time" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension is user-specified - ! and not assumed from the t_grid value. - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - start(:) = 1 - nwrite(:) = dim_lengths(1) - if (present(start_index)) then - start(1) = max(1,start_index(1)) - endif - - if (present(edge_lengths)) then - nwrite(1) = max(dim_lengths(1),edge_lengths(1)) - endif - - if (.not.(check_if_open(fileobj_write_field))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_1d_noDD:mode argument must be write, overwrite, or append") - ! get the index of the corresponding time_level the first time the file is opened - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! get the pes associated with the file. - !>\note this is required so that only pe(1) is identified as the root pe to create the file - !! Otherwise, multiple pes may try to open the file in write (NC_NOCLOBBER) mode, leading to failure - if (.not.(allocated(pelist))) then - allocate(pelist(num_PEs())) - pelist(:) = 0 - do i=1,size(pelist) - pelist(i) = i-1 - enddo - endif - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field, trim(filename_temp), lowercase(trim(mode)), & - is_restart=.false., pelist=pelist) - endif - ! write the data, and the time value if it is not already written to the file - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field,dim_unlim_name) - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time value if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/), & - edge_lengths=(/1/)) - endif - endif - ! register the field if it is not already in the file - if (.not.(variable_exists(fileobj_write_field, trim(fieldname)))) then - call register_field(fileobj_write_field, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = '' - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the variable to the file - if (present(time_level)) then - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field)) call fms2_close_file(fileobj_write_field) - if (allocated(pelist)) deallocate(pelist) - write_field_time_index = 0 - endif - - -end subroutine write_field_1d_noDD - -!> This function uses the fms_io function write_data to write a scalar variable named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_2d_noDD(filename, fieldname, data, mode, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:), intent(in) :: data !< The 2-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(2), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(2), optional, intent(in) :: edge_lengths !< number of data values to read in; default is the - !! variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success ! .true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: i, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(2) :: start, nwrite ! indices for starting points and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(3) :: dim_names ! variable dimension names - integer, dimension(3) :: dim_lengths ! variable dimension lengths - integer, allocatable, dimension(:) :: pelist ! list of pes associated with the netCDF file - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size = 0 - dim_unlim_name = "" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension is user-specified - ! and not assumed from the t_grid value - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - - ! set the start (start_index) and nwrite (edge_lengths) values - ndims=2 - start(:) = 1 - nwrite(:) = dim_lengths(1:2) - if (present(start_index)) then - do i=1,ndims - start(i) = max(1,start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i),edge_lengths(i)) - enddo - endif - - if (.not.(check_if_open(fileobj_write_field))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_2d_noDD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! get the pes associated with the file. - !>\note this is required so that only pe(1) is identified as the root pe to create the file - !! Otherwise, multiple pes may try to open the file in write (NC_NOCLOBBER) mode, leading to failure - if(.not.(allocated(pelist))) then - allocate(pelist(num_PEs())) - pelist(:) = 0 - do i=1,size(pelist) - pelist(i) = i-1 - enddo - endif - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field, trim(filename_temp), lowercase(trim(mode)), & - is_restart=.false., pelist=pelist) - endif - - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field,dim_unlim_name) - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time value if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/), & - edge_lengths=(/1/)) - endif - endif - - ! register the variable to the file - if (.not.(variable_exists(fileobj_write_field, trim(fieldname)))) then - call register_field(fileobj_write_field, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the variable to the file - if (present(time_level)) then - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field)) call fms2_close_file(fileobj_write_field) - if (allocated(pelist)) deallocate(pelist) - write_field_time_index=0 - endif - -end subroutine write_field_2d_noDD - -!> This function uses the fms_io function write_data to write a 3-D non-domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_3d_noDD(filename, fieldname, data, mode, hor_grid, z_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:,:), intent(in) :: data !< The 3-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - integer, dimension(4), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(4), optional, intent(in) :: edge_lengths !< number of data values to read in; default is the - !! variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: i, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(3) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time_units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(4) :: dim_names ! variable dimension names - integer, dimension(4) :: dim_lengths ! variable dimension lengths - integer, allocatable, dimension(:) :: pelist ! list of pes associated with the netCDF file - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size = 0 - dim_unlim_name = "" - dim_names(:) = "" - dim_lengths(:) = 0 - num_dims = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - ! get the dimension names and lengths - ! NOTE: the t_grid argument is set to '1' (do nothing) because the presence of a time dimension is user-specified - ! and not assumed from the t_grid value - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, dG=dG) - endif - - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, '1', dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - ndims = 3 - start(:) = 1 - nwrite(:) = dim_lengths(1:3) - if (present(start_index)) then - do i=1,ndims - start(i) = max(1,start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i), edge_lengths(i)) - enddo - endif - - ! open the file - if (.not.(check_if_open(fileobj_write_field))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_io:write_3d_noDD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! get the pes associated with the file. - !>\note this is required so that only pe(1) is identified as the root pe to create the file - !! Otherwise, multiple pes may try to open the file in write (NC_NOCLOBBER) mode, leading to failure - if (.not.(allocated(pelist))) then - allocate(pelist(num_PEs())) - pelist(:) = 0 - do i=1,size(pelist) - pelist(i) = i-1 - enddo - endif - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field, trim(filename_temp), lowercase(trim(mode)), & - is_restart=.false., pelist=pelist) - endif - ! register and write the time_level - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field,dim_unlim_name) - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - - if (.not. (variable_exists(fileobj_write_field, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - ! write the time_level if it is larger than the most recent file time - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/), & - edge_lengths=(/1/)) - endif - endif - ! register the field if it is not already in the file - if (.not.(variable_exists(fileobj_write_field, trim(fieldname)))) then - call register_field(fileobj_write_field, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field, trim(fieldname), "checksum", checksum_char) - endif - endif - - if (present(time_level)) then - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field)) call fms2_close_file(fileobj_write_field) - if (allocated(pelist)) deallocate(pelist) - write_field_time_index=0 - endif - -end subroutine write_field_3d_noDD - -!> This function uses the fms_io function write_data to write a 4-D non-domain-decomposed data field named "fieldname" -!! to the file "filename" in "write", "overwrite", or "append" mode. It should be called after create_file in the MOM -!! file write procedure. -subroutine write_field_4d_noDD(filename, fieldname, data, mode, hor_grid, z_grid, t_grid, & - start_index, edge_lengths, time_level, time_units, & - checksums, G, dG, GV, leave_file_open, units, longname) - character(len=*), intent(in) :: filename !< The name of the file to read - character(len=*), intent(in) :: fieldname !< The variable name of the data in the file - real, target, dimension(:,:,:,:), intent(in) :: data !< The 1-dimensional data array to pass to read_data - character(len=*), intent(in) :: mode !< "write", "overwrite", or "append" - character(len=*), intent(in) :: hor_grid !< horizontal grid descriptor - character(len=*), intent(in) :: z_grid !< vertical grid descriptor - character(len=*), intent(in) :: t_grid !< time descriptor - integer, dimension(4), optional, intent(in) :: start_index !< starting index of data buffer. Default is 1 - integer, dimension(4), optional, intent(in) :: edge_lengths !< number of data values to read in; default is the - !! variable size - real, optional, intent(in) :: time_level !< time value to write - real, optional, intent(in) :: time_units !< length of the units for time [s]. The - !! default value is 86400.0, for 1 day. - integer(kind=8), dimension(:,:), optional, intent(in) :: checksums !< variable checksum - type(ocean_grid_type), optional, intent(in) :: G !< ocean horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(dyn_horgrid_type), optional, intent(in) :: dG !< dynamic horizontal grid structure; G or dG - !! is required if the new file uses any - !! horizontal grid axes. - type(verticalGrid_type), optional, intent(in) :: GV !< ocean vertical grid structure, which is - !! required if the new file uses any - !! vertical grid axes. - logical, optional, intent(in) :: leave_file_open !< flag indicating whether to leave the file open - character(len=*), optional, intent(in) :: units !< variable units - character(len=*), optional, intent(in) :: longname !< long name variable attribute - ! local - logical :: file_open_success !.true. if call to open_file is successful - logical :: close_the_file ! indicates whether to close the file after write_data is called; default is .true. - integer :: i, ndims, num_dims, substring_index - integer :: dim_unlim_size ! size of the unlimited dimension - integer, dimension(4) :: start, nwrite ! indices for first data value and number of values to write - character(len=20) :: t_units ! time units - character(len=nf90_max_name) :: dim_unlim_name ! name of the unlimited dimension in the file - character(len=1024) :: filename_temp - character(len=64) :: checksum_char ! checksum character array created from checksum argument - character(len=48), dimension(4) :: dim_names ! variable dimension names - integer, dimension(4) :: dim_lengths ! variable dimension lengths - integer, allocatable, dimension(:) :: pelist ! list of pes associated with the netCDF file - - close_the_file = .true. - if (present(leave_file_open)) close_the_file = .not.(leave_file_open) - - dim_unlim_size = 0 - dim_unlim_name = "" - dim_names(:) = "" - dim_lengths(:) = 0 - ndims = 4 - num_dims = 0 - ! append '.nc' to the file name if it is missing - filename_temp = "" - substring_index = 0 - substring_index = index(trim(filename), ".nc") - if (substring_index <= 0) then - filename_temp = append_substring(filename,".nc") - else - filename_temp = filename - endif - - ! get the dimension names and lengths - if (present(G)) then - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, G=G) - elseif(present(dG)) then - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, dG=dG) - endif - if (present(GV)) & - call get_var_dimension_metadata(hor_grid, z_grid, t_grid, dim_names, & - dim_lengths, num_dims, GV=GV) - ! set the start (start_index) and nwrite (edge_lengths) values - start(:) = 1 - nwrite(:) = dim_lengths(:) - if (present(start_index)) then - do i=1,ndims - start(i) = max(1, start_index(i)) - enddo - endif - - if (present(edge_lengths)) then - do i=1,ndims - nwrite(i) = max(dim_lengths(i), edge_lengths(i)) - enddo - endif - - ! open the file - if (.not.(check_if_open(fileobj_write_field))) then - if ((lowercase(trim(mode)) .ne. "write") .and. (lowercase(trim(mode)) .ne. "append") .and. & - (lowercase(trim(mode)) .ne. "overwrite")) & - call MOM_error(FATAL,"MOM_write_field_fms2:write_4d_noDD:mode argument must be write, overwrite, or append") - ! get the time_level index - if (present(time_level)) write_field_time_index = get_time_index(trim(filename_temp), time_level) - ! get the pes associated with the file. - !>\note this is required so that only pe(1) is identified as the root pe to create the file - !! Otherwise, multiple pes may try to open the file in write (NC_NOCLOBBER) mode, leading to failure - if (.not.(allocated(pelist))) then - allocate(pelist(num_PEs())) - pelist(:) = 0 - do i=1,size(pelist) - pelist(i) = i-1 - enddo - endif - ! open the file in write or append mode - file_open_success = fms2_open_file(fileobj_write_field, trim(filename_temp), lowercase(trim(mode)), & - is_restart=.false., pelist=pelist) - endif - ! register and write the time_level - if (present(time_level)) then - call get_unlimited_dimension_name(fileobj_write_field,dim_unlim_name) - call get_dimension_size(fileobj_write_field, trim(dim_unlim_name), dim_unlim_size) - num_dims=num_dims+1 - dim_names(num_dims) = trim(dim_unlim_name) - ! write the time value if it is not already written to the file - if (.not. (variable_exists(fileobj_write_field, trim(dim_unlim_name)))) then - ! set the time units - t_units = "" - if (present(time_units)) then - t_units = get_time_units(time_units) - else - t_units = "days" - endif - - call register_field(fileobj_write_field, trim(dim_unlim_name), "double", dimensions=(/trim(dim_unlim_name)/)) - call register_variable_attribute(fileobj_write_field, trim(dim_unlim_name), 'units', trim(t_units)) - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/)) - else - if (write_field_time_index .gt. dim_unlim_size) & - call write_data(fileobj_write_field, trim(dim_unlim_name), (/time_level/), corner=(/write_field_time_index/), & - edge_lengths=(/1/)) - endif - endif - ! register the variable - if (.not.(variable_exists(fileobj_write_field, trim(fieldname)))) then - call register_field(fileobj_write_field, trim(fieldname), "double", dimensions=dim_names(1:num_dims)) - if (present(units)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'units', trim(units)) - if (present(longname)) & - call register_variable_attribute(fileobj_write_field, trim(fieldname), 'long_name', trim(longname)) - ! write the checksum attribute - if (present(checksums)) then - ! convert the checksum to a string - checksum_char = "" - checksum_char = convert_checksum_to_string(checksums(1,1)) - call register_variable_attribute(fileobj_write_field, trim(fieldname), "checksum", checksum_char) - endif - endif - ! write the variable to the file - if (present(time_level)) then - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite, & - unlim_dim_level=write_field_time_index) - else - call write_data(fileobj_write_field, trim(fieldname), data, corner=start, edge_lengths=nwrite) - endif - ! close the file - if (close_the_file) then - if (check_if_open(fileobj_write_field)) call fms2_close_file(fileobj_write_field) - deallocate(pelist) - write_field_time_index=0 - endif -end subroutine write_field_4d_nodd - -end module MOM_write_field_fms2 diff --git a/src/framework/MOM_horizontal_regridding.F90 b/src/framework/MOM_horizontal_regridding.F90 index 8d02dfdf3f..73fb1f0a41 100644 --- a/src/framework/MOM_horizontal_regridding.F90 +++ b/src/framework/MOM_horizontal_regridding.F90 @@ -11,9 +11,9 @@ module MOM_horizontal_regridding use MOM_error_handler, only : callTree_enter, callTree_leave, callTree_waypoint use MOM_file_parser, only : get_param, log_param, log_version, param_file_type use MOM_grid, only : ocean_grid_type -use MOM_interpolate, only : time_interp_external, get_external_field_info, horiz_interp_init +use MOM_interpolate, only : time_interp_external, horiz_interp_init use MOM_interpolate, only : build_horiz_interp_weights, run_horiz_interp, horiz_interp_type -use MOM_io_infra, only : axistype, get_axis_data +use MOM_interp_infra, only : axistype, get_external_field_info, get_axis_data use MOM_time_manager, only : time_type use netcdf, only : NF90_OPEN, NF90_NOWRITE, NF90_GET_ATT, NF90_GET_VAR diff --git a/src/framework/MOM_interpolate.F90 b/src/framework/MOM_interpolate.F90 index f282d03ff6..4a931d0bf3 100644 --- a/src/framework/MOM_interpolate.F90 +++ b/src/framework/MOM_interpolate.F90 @@ -9,7 +9,6 @@ module MOM_interpolate use MOM_interp_infra, only : time_interp_external_init, get_external_field_info use MOM_interp_infra, only : horiz_interp_type, horiz_interp_init use MOM_interp_infra, only : run_horiz_interp, build_horiz_interp_weights -use MOM_io_infra, only : axistype use MOM_time_manager, only : time_type implicit none ; private diff --git a/src/framework/MOM_io.F90 b/src/framework/MOM_io.F90 index 1159ac87e1..247a0a9678 100644 --- a/src/framework/MOM_io.F90 +++ b/src/framework/MOM_io.F90 @@ -254,11 +254,11 @@ subroutine create_file(IO_handle, filename, vars, novars, fields, threading, tim if (use_latq) & call write_metadata(IO_handle, axis_latq, name="latq", units=y_axis_units, longname="Latitude", & - cartesian='Y', domain=y_domain, data=gridLatB(JsgB:JegB)) + cartesian='Y', domain=y_domain, data=gridLatB(JsgB:JegB), edge_axis=.true.) if (use_lonq) & call write_metadata(IO_handle, axis_lonq, name="lonq", units=x_axis_units, longname="Longitude", & - cartesian='X', domain=x_domain, data=gridLonB(IsgB:IegB)) + cartesian='X', domain=x_domain, data=gridLonB(IsgB:IegB), edge_axis=.true.) if (use_layer) & call write_metadata(IO_handle, axis_layer, name="Layer", units=trim(GV%zAxisUnits), & @@ -383,7 +383,7 @@ subroutine reopen_file(IO_handle, filename, vars, novars, fields, threading, tim type(MOM_domain_type), pointer :: Domain => NULL() character(len=200) :: check_name, mesg - integer :: length, ndim, nvar, natt, ntime, thread + integer :: length, nvar, thread logical :: exists, one_file, domain_set thread = SINGLE_FILE @@ -418,7 +418,7 @@ subroutine reopen_file(IO_handle, filename, vars, novars, fields, threading, tim endif if (.not.file_is_open(IO_handle)) return - call get_file_info(IO_handle, ndim, nvar, natt, ntime) + call get_file_info(IO_handle, nvar=nvar) if (nvar == -1) then write (mesg,*) "Reopening file ",trim(filename)," apparently had ",nvar,& @@ -1343,7 +1343,7 @@ end subroutine query_vardesc !> Write a 4d field to an output file, potentially with rotation subroutine MOM_write_field_4d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, & fill_value, turns, scale) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:,:), intent(inout) :: field !< Unrotated field to write @@ -1378,7 +1378,7 @@ end subroutine MOM_write_field_4d !> Write a 3d field to an output file, potentially with rotation subroutine MOM_write_field_3d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, & fill_value, turns, scale) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:,:), intent(inout) :: field !< Unrotated field to write @@ -1413,7 +1413,7 @@ end subroutine MOM_write_field_3d !> Write a 2d field to an output file, potentially with rotation subroutine MOM_write_field_2d(IO_handle, field_md, MOM_domain, field, tstamp, tile_count, & fill_value, turns, scale) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata type(MOM_domain_type), intent(in) :: MOM_domain !< The MOM_Domain that describes the decomposition real, dimension(:,:), intent(inout) :: field !< Unrotated field to write @@ -1447,7 +1447,7 @@ end subroutine MOM_write_field_2d !> Write a 1d field to an output file subroutine MOM_write_field_1d(IO_handle, field_md, field, tstamp, fill_value, scale) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, dimension(:), intent(in) :: field !< Field to write real, optional, intent(in) :: tstamp !< Model timestamp @@ -1476,7 +1476,7 @@ end subroutine MOM_write_field_1d !> Write a 0d field to an output file subroutine MOM_write_field_0d(IO_handle, field_md, field, tstamp, fill_value, scale) - type(file_type), intent(in) :: IO_handle !< Handle for a file that is open for writing + type(file_type), intent(inout) :: IO_handle !< Handle for a file that is open for writing type(fieldtype), intent(in) :: field_md !< Field type with metadata real, intent(in) :: field !< Field to write real, optional, intent(in) :: tstamp !< Model timestamp