44use fms_mod,
only: check_nml_error, stdlog, &
45 mpp_pe, mpp_root_pe, write_version_number, &
46 error_mesg, fatal, note, &
52use constants_mod,
only: pi
53use mpp_mod,
only: input_nml_file
54use platform_mod,
only: r4_kind, r8_kind, fms_path_len
59public :: topography_init, &
85 module procedure get_topog_mean_1d_r4, get_topog_mean_1d_r8
86 module procedure get_topog_mean_2d_r4, get_topog_mean_2d_r8
110 module procedure get_topog_stdev_1d_r4, get_topog_stdev_1d_r8
111 module procedure get_topog_stdev_2d_r4, get_topog_stdev_2d_r8
130 module procedure get_ocean_frac_1d_r4, get_ocean_frac_1d_r8
131 module procedure get_ocean_frac_2d_r4, get_ocean_frac_2d_r8
151 module procedure get_ocean_mask_1d_r4, get_ocean_mask_1d_r8
152 module procedure get_ocean_mask_2d_r4, get_ocean_mask_2d_r8
171 module procedure get_water_frac_1d_r4, get_water_frac_1d_r8
172 module procedure get_water_frac_2d_r4, get_water_frac_2d_r8
191 module procedure get_water_mask_1d_r4, get_water_mask_1d_r8
192 module procedure get_water_mask_2d_r4, get_water_mask_2d_r8
196 module procedure interp_topog_1d_r4, interp_topog_1d_r8
197 module procedure interp_topog_2d_r4, interp_topog_2d_r8
201 module procedure find_indices_r4, find_indices_r8
205 module procedure input_data_r4, input_data_r8
209 module procedure interp_water_1d_r4, interp_water_1d_r8
210 module procedure interp_water_2d_r4, interp_water_2d_r8
214 module procedure determine_ocean_points_r4, determine_ocean_points_r8
220character(len=FMS_PATH_LEN) :: topog_file =
'DATA/navy_topography.data', &
221 water_file =
'DATA/navy_pctwater.data'
222namelist /topography_nml/ topog_file, water_file
224integer,
parameter :: TOPOG_INDEX = 1
225integer,
parameter :: WATER_INDEX = 2
226logical :: file_is_opened(2) = .false.
257integer,
parameter :: compute_stdev = 123
263#include<file_version.h>
265logical :: module_is_initialized = .false.
273subroutine topography_init ()
274 if ( module_is_initialized )
return
275 call write_version_number(
"TOPOGRAPHY_MOD", version)
277 module_is_initialized = .true.
278end subroutine topography_init
285function open_topog_file ( )
286logical :: open_topog_file
287real(kind=r4_kind) :: r_ipts, r_jpts
290namelen = len(trim(topog_file))
291 if ( file_exists(topog_file) .AND. topog_file(namelen-2:namelen) ==
'.nc')
then
292 if (mpp_pe() == mpp_root_pe())
call mpp_error (
'topography_mod', &
293 'Reading NetCDF formatted input data file: '//trim(topog_file), note)
294 if(.not. file_is_opened(topog_index) )
then
296 call mpp_error(fatal,
'topography_mod: Error in opening file '//trim(topog_file))
304 open_topog_file = .true.
305 file_is_opened(topog_index) = .true.
307 open_topog_file = .false.
310end function open_topog_file
312function open_water_file ( )
313logical :: open_water_file
314real(kind=r4_kind) :: r_ipts, r_jpts
317namelen = len(trim(water_file))
318if ( file_exists(water_file) .AND. water_file(namelen-2:namelen) ==
'.nc')
then
319 if (mpp_pe() == mpp_root_pe())
call mpp_error (
'topography_mod', &
320 'Reading NetCDF formatted input data file: '//trim(water_file), note)
321 if(.not. file_is_opened(water_index) )
then
323 call mpp_error(fatal,
'topography_mod: Error in opening file '//trim(water_file))
331 open_water_file = .true.
332 file_is_opened(water_index) = .true.
334 open_water_file = .false.
337end function open_water_file
346integer :: iunit, ierr, io
350read (input_nml_file, topography_nml, iostat=io)
351ierr = check_nml_error(io,
'topography_nml')
355if (mpp_pe() == mpp_root_pe())
then
357 write (iunit, nml=topography_nml)
362#include "topography_r4.fh"
363#include "topography_r8.fh"
365end module topography_mod
Opens a given netcdf or domain file.
Read data from a defined field in a file.
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.
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.