51 use constants_mod,
only: pi
52 use mpp_mod,
only: input_nml_file
53 use platform_mod,
only: r4_kind, r8_kind, fms_path_len
58 public :: topography_init, &
84 module procedure get_topog_mean_1d_r4, get_topog_mean_1d_r8
85 module procedure get_topog_mean_2d_r4, get_topog_mean_2d_r8
109 module procedure get_topog_stdev_1d_r4, get_topog_stdev_1d_r8
110 module procedure get_topog_stdev_2d_r4, get_topog_stdev_2d_r8
129 module procedure get_ocean_frac_1d_r4, get_ocean_frac_1d_r8
130 module procedure get_ocean_frac_2d_r4, get_ocean_frac_2d_r8
150 module procedure get_ocean_mask_1d_r4, get_ocean_mask_1d_r8
151 module procedure get_ocean_mask_2d_r4, get_ocean_mask_2d_r8
170 module procedure get_water_frac_1d_r4, get_water_frac_1d_r8
171 module procedure get_water_frac_2d_r4, get_water_frac_2d_r8
190 module procedure get_water_mask_1d_r4, get_water_mask_1d_r8
191 module procedure get_water_mask_2d_r4, get_water_mask_2d_r8
195 module procedure interp_topog_1d_r4, interp_topog_1d_r8
196 module procedure interp_topog_2d_r4, interp_topog_2d_r8
200 module procedure find_indices_r4, find_indices_r8
204 module procedure input_data_r4, input_data_r8
208 module procedure interp_water_1d_r4, interp_water_1d_r8
209 module procedure interp_water_2d_r4, interp_water_2d_r8
213 module procedure determine_ocean_points_r4, determine_ocean_points_r8
219 character(len=FMS_PATH_LEN) :: topog_file =
'DATA/navy_topography.data', &
220 water_file =
'DATA/navy_pctwater.data'
221 namelist /topography_nml/ topog_file, water_file
223 integer,
parameter :: TOPOG_INDEX = 1
224 integer,
parameter :: WATER_INDEX = 2
225 logical :: file_is_opened(2) = .false.
255 integer :: ipts, jpts
256 integer,
parameter :: compute_stdev = 123
262 #include<file_version.h>
264 logical :: module_is_initialized = .false.
272 subroutine topography_init ()
273 if ( module_is_initialized )
return
276 module_is_initialized = .true.
277 end subroutine topography_init
284 function open_topog_file ( )
285 logical :: open_topog_file
286 real(kind=r4_kind) :: r_ipts, r_jpts
289 namelen = len(trim(topog_file))
290 if ( file_exists(topog_file) .AND. topog_file(namelen-2:namelen) ==
'.nc')
then
292 'Reading NetCDF formatted input data file: '//trim(topog_file), note)
293 if(.not. file_is_opened(topog_index) )
then
295 call mpp_error(fatal,
'topography_mod: Error in opening file '//trim(topog_file))
303 open_topog_file = .true.
304 file_is_opened(topog_index) = .true.
306 open_topog_file = .false.
309 end function open_topog_file
311 function open_water_file ( )
312 logical :: open_water_file
313 real(kind=r4_kind) :: r_ipts, r_jpts
316 namelen = len(trim(water_file))
317 if ( file_exists(water_file) .AND. water_file(namelen-2:namelen) ==
'.nc')
then
319 'Reading NetCDF formatted input data file: '//trim(water_file), note)
320 if(.not. file_is_opened(water_index) )
then
322 call mpp_error(fatal,
'topography_mod: Error in opening file '//trim(water_file))
330 open_water_file = .true.
331 file_is_opened(water_index) = .true.
333 open_water_file = .false.
336 end function open_water_file
345 integer :: iunit, ierr, io
349 read (input_nml_file, topography_nml, iostat=io)
354 if (
mpp_pe() == mpp_root_pe())
then
356 write (iunit, nml=topography_nml)
361 #include "topography_r4.fh"
362 #include "topography_r8.fh"
364 end module topography_mod
Opens a given netcdf or domain file.
Read data from a defined field in a file.
integer function, public check_nml_error(IOSTAT, NML_NAME)
Checks the iostat argument that is returned after reading a namelist and determines if the error code...
subroutine, public write_version_number(version, tag, unit)
Prints to the log file (or a specified unit) the version id string and tag name.
subroutine, public error_mesg(routine, message, level)
Print notes, warnings and error messages; terminates program for warning and error messages....
subroutine, public horiz_interp_del(Interp)
Deallocates memory used by "horiz_interp_type" variables. Must be called before reinitializing with h...
Subroutine for performing the horizontal interpolation between two grids.
integer function stdlog()
This function returns the current standard fortran unit numbers for log messages. Log messages,...
integer function mpp_pe()
Returns processor ID.
subroutine read_namelist
Reads the namelist file, write namelist to log file, and initializes constants.
type(fmsnetcdffile_t), dimension(2) fileobj
needed for fms2_io
Returns fractional area covered by ocean in a grid box. Returns fractional area covered by ocean in t...
Returns a land-ocean mask in a grid box.
Returns a "realistic" mean surface height field.
Returns fractional area covered by water.
Returns a land-water mask in a grid box.
Allocates space and initializes a derived-type variable that contains pre-computed interpolation indi...
Returns a standard deviation of higher resolution topography with the given model grid boxes.