CrocoDash.rm6.regional_mom6 package#
Submodules#
CrocoDash.rm6.regional_mom6.MOM_parameter_tools module#
- CrocoDash.rm6.regional_mom6.MOM_parameter_tools.change_MOM_parameter(directory, param_name, param_value=None, comment=None, override=True, delete=False)#
Changes a parameter in the
MOM_input
orMOM_override
file. Returns original value, if there was one. If delete keyword argument is True the parameter is removed. But note, that only parameters from theMOM_override
are be deleted, since deleting parameters fromMOM_input
is not safe and can lead to errors. If the parameter does not exist, it will be added to the MOM_override file.MOM parameter files need to be present in the run directory
- Parameters:
param_name (str) – Parameter name to modify
directory (str) – Directory where the MOM_input and MOM_override files are located
param_value (Optional[str]) – New assigned Value
comment (Optional[str]) – Any comment to add
delete (Optional[bool]) – Whether to delete the specified
param_name
- CrocoDash.rm6.regional_mom6.MOM_parameter_tools.read_MOM_file_as_dict(filename, directory)#
Read the MOM_input file from directory and return a dictionary of the variables and their values.
- CrocoDash.rm6.regional_mom6.MOM_parameter_tools.write_MOM_file(MOM_file_dict, directory)#
Write the MOM_input file from a dictionary of variables and their values to directory. Does not support removing fields.
CrocoDash.rm6.regional_mom6.config module#
- class CrocoDash.rm6.regional_mom6.config.Config#
Bases:
object
Handles saving and loading experiment class configurations, including a function.
- static load_from_json(filename, mom_input_dir='mom_input/from_config', mom_run_dir='mom_run/from_config', create_hgrid_and_vgrid=True)#
Load experiment variables from a configuration file and generate the horizontal and vertical grids (
hgrid
/vgrid
).(This is basically another way to initialize an experiment.)
Arguments: config_file_path (str): Path to the config file. mom_input_dir (str): Path to the MOM6 input directory. Default:
"mom_input/from_config"
. mom_run_dir (str): Path to the MOM6 run directory. Default:"mom_run/from_config"
. create_hgrid_and_vgrid (bool): Whether to create the hgrid and the vgrid. Default is True.Returns: An experiment object with the fields from the configuration at
config_file_path
.
- static save_to_json(obj, path=None, export=True)#
Write a
json
configuration file for the experiment. Thejson
file contains the experiment variable information to allow for easy pass off to other users, with a strict computer independence restriction. It also makes information about the expirement readable, and can be also useful for simply printing out information about the experiment.- Parameters:
obj (rm6.Experiment) – The experiment object to save.
path (str) – Path to write the config file to. If not provided, the file is written to the
mom_run_dir
directory.export (bool) – If
True
(default), the configuration file is written to disk on the givenpath
- Returns:
A dictionary containing the configuration information.
- Return type:
Dict
- CrocoDash.rm6.regional_mom6.config.convert_string_to_value(value, value_type: str)#
Helper function to handle deserialization of different types.
- CrocoDash.rm6.regional_mom6.config.convert_value_to_string(value)#
Helper function to handle serialization of different types.
CrocoDash.rm6.regional_mom6.regional_mom6 module#
- class CrocoDash.rm6.regional_mom6.regional_mom6.experiment(*, date_range, resolution, number_vertical_layers, layer_thickness_ratio, depth, mom_run_dir, mom_input_dir, fre_tools_dir=None, longitude_extent=None, latitude_extent=None, hgrid_type='even_spacing', hgrid_path=None, vgrid_type='hyperbolic_tangent', vgrid_path=None, repeat_year_forcing=False, minimum_depth=4, tidal_constituents=['M2', 'S2', 'N2', 'K2', 'K1', 'O1', 'P1', 'Q1', 'MM', 'MF'], create_empty=False, expt_name=None, boundaries=['south', 'north', 'west', 'east'], regridding_method='bilinear', fill_method=<function fill_missing_data>)#
Bases:
object
The main class for setting up a regional experiment.
Everything about the regional experiment.
Methods in this class generate the various input files needed for a MOM6 experiment forced with open boundary conditions (OBCs). The code is agnostic to the user’s choice of boundary forcing, bathymetry, and surface forcing; users need to prescribe what variables are all called via mapping dictionaries from MOM6 variable/coordinate name to the name in the input dataset.
The class can be used to generate the grids for a new experiment, or to read in an existing one (see argument description below).
- Parameters:
date_range (Tuple[str]) – Start and end dates of the boundary forcing window. For example:
("2003-01-01", "2003-01-31")
.resolution (float) – Lateral resolution of the domain (in degrees).
number_vertical_layers (int) – Number of vertical layers.
layer_thickness_ratio (float) – Ratio of largest to smallest layer thickness; used as input in
hyperbolictan_thickness_profile()
.depth (float) – Depth of the domain.
mom_run_dir (str) – Path of the MOM6 control directory.
mom_input_dir (str) – Path of the MOM6 input directory, to receive the forcing files.
fre_tools_dir (str) – Path of GFDL’s FRE tools (NOAA-GFDL/FRE-NCtools) binaries.
longitude_extent (Tuple[float]) – Extent of the region in longitude (in degrees). For example:
(40.5, 50.0)
.latitude_extent (Tuple[float]) – Extent of the region in latitude (in degrees). For example:
(-20.0, 30.0)
.hgrid_type (str) – Type of horizontal grid to generate. Currently, only
'even_spacing'
is supported. Setting this argument to'from_file'
requires the additional hgrid_path argumenthgrid_path (str) – Path to the horizontal grid file if the hgrid_type is
'from_file'
.vgrid_type (str) – Type of vertical grid to generate. Currently, only
'hyperbolic_tangent'
is supported. Setting this argument to'from_file'
requires the additional vgrid_path argumentvgrid_path (str) – Path to the vertical grid file if the vgrid_type is
'from_file'
.repeat_year_forcing (bool) – When
True
the experiment runs with repeat-year forcing. WhenFalse
(default) then inter-annual forcing is used.minimum_depth (int) – The minimum depth in meters of a grid cell allowed before it is masked out and treated as land.
tidal_constituents (List[str]) – List of tidal constituents to be used in the experiment. Default is
["M2", "S2", "N2", "K2", "K1", "O1", "P1", "Q1", "MM", "MF"]
.create_empty (bool) – If
True
, the experiment object is initialized empty. This is used for testing and experienced user manipulation.expt_name (str) – The name of the experiment (for config file use)
boundaries (List[str]) – List of (rectangular) boundaries to be set. Default is
["south", "north", "west", "east"]
. The boundaries are set as (list index + 1) in MOM_override in the order of the list, and less than 4 boundaries can be set.regridding_method (str) – regridding method to use throughout the entire experiment. Default is
'bilinear'
. Any other xesmf regridding method can be used.fill_method (Function) – The fill function to be used after regridding datasets. it takes a xarray DataArray and returns a filled DataArray. Default is
rgd.fill_missing_data
.
- property bathymetry#
- property bathymetry_path#
Finds the bathymetry file from disk, and returns the file path.
- configure_cpu_layout(layout)#
Wrapper for the
check_mask
function of GFDL’s FRE Tools. User provides processorlayout
tuple of processing units.
- classmethod create_empty(longitude_extent=None, latitude_extent=None, date_range=None, resolution=None, number_vertical_layers=None, layer_thickness_ratio=None, depth=None, mom_run_dir=None, mom_input_dir=None, fre_tools_dir=None, hgrid_type='even_spacing', repeat_year_forcing=False, minimum_depth=4, tidal_constituents=['M2', 'S2', 'N2', 'K2', 'K1', 'O1', 'P1', 'Q1', 'MM', 'MF'], expt_name=None, boundaries=['south', 'north', 'west', 'east'], regridding_method='bilinear', fill_method=<function fill_missing_data>)#
Note: This method is unsafe; only experience users are urged to use it!
Alternative to the initialisation method to create an empty expirement object, with the opportunity to override whatever values wanted.
This method allows developers and experienced users to set specific variables for specific function requirements, like just regridding the initial condition or subsetting bathymetry, instead of having to set so many other variables that aren’t needed.
- property era5_paths#
Finds the ERA5 files from disk, and prints the file paths
- find_MOM6_rectangular_orientation(input)#
Convert between MOM6 boundary and the specific segment number needed, or the inverse.
- get_glorys(raw_boundaries_path)#
This is a wrapper that calls
get_glorys_data()
once for each of the rectangular boundary segments and the initial condition. For more complex boundary shapes, callget_glorys_data()
directly for each of your boundaries that aren’t parallel to lines of constant latitude or longitude. For example, for an angled Northern boundary that spans multiple latitudes, we need to download a wider rectangle containing the entire boundary.- Parameters:
raw_boundaries_path (str) – Path to the directory containing the raw boundary forcing files.
boundaries (List[str]) – List of cardinal directions for which to create boundary forcing files. Default is
["south", "north", "west", "east"]
.
- property init_tracers#
- property init_velocities#
- property initial_condition_paths#
Finds the initial condition files from disk, and prints the file paths
- property ocean_state_boundary_paths#
Finds the ocean state files from disk, and prints the file paths
- run_FRE_tools(layout=None)#
A wrapper for FRE Tools
check_mask
,make_solo_mosaic
, andmake_quick_mosaic
. User provides processorlayout
tuple of processing units.
- setup_bathymetry(*, bathymetry_path, longitude_coordinate_name='lon', latitude_coordinate_name='lat', vertical_coordinate_name='elevation', fill_channels=False, positive_down=False, write_to_file=True, regridding_method=None)#
Cut out and interpolate the chosen bathymetry and then fill inland lakes.
Users can optionally fill narrow channels (see
fill_channels
keyword argument below). Note, however, that narrow channels are less of an issue for models that are discretized on an Arakawa C grid, like MOM6.Output is saved in the input directory of the experiment.
- Parameters:
bathymetry_path (str) – Path to the netCDF file with the bathymetry.
longitude_coordinate_name (Optional[str]) – The name of the longitude coordinate in the bathymetry dataset at
bathymetry_path
. For example, for GEBCO bathymetry:'lon'
(default).latitude_coordinate_name (Optional[str]) – The name of the latitude coordinate in the bathymetry dataset at
bathymetry_path
. For example, for GEBCO bathymetry:'lat'
(default).vertical_coordinate_name (Optional[str]) – The name of the vertical coordinate in the bathymetry dataset at
bathymetry_path
. For example, for GEBCO bathymetry:'elevation'
(default).fill_channels (Optional[bool]) – Whether or not to fill in diagonal channels. This removes more narrow inlets, but can also connect extra islands to land. Default:
False
.positive_down (Optional[bool]) – If
True
, it assumes that the bathymetry vertical coordinate is positive downwards. Default:False
.write_to_file (Optional[bool]) – Whether to write the bathymetry to a file. Default:
True
.regridding_method (Optional[str]) – The type of regridding method to use. Defaults to self.regridding_method
- setup_boundary_tides(tpxo_elevation_filepath, tpxo_velocity_filepath, tidal_constituents=None, bathymetry_path=None, rotational_method=RotationMethod.EXPAND_GRID, regridding_method=None, fill_method=None)#
Subset the tidal data and generate more boundary files.
- Parameters:
path_to_td (str) – Path to boundary tidal file.
tpxo_elevation_filepath – Filepath to the TPXO elevation product. Generally of the form
h_tidalversion.nc
tpxo_velocity_filepath – Filepath to the TPXO velocity product. Generally of the form
u_tidalversion.nc
tidal_constituents – List of tidal constituents to include in the regridding. Default is set in the experiment constructor (See
Experiment
)bathymetry_path (str) – Path to the bathymetry file. Default is
None
, in which case the boundary condition is not maskedrotational_method (str) – Method to use for rotating the tidal velocities. Default is ‘EXPAND_GRID’.
regridding_method (Optional[str]) – The type of regridding method to use. Defaults to self.regridding_method
fill_method (Function) – Fill method to use throughout the function. Default is
self.fill_method
- Returns:
Regridded tidal velocity and elevation files in ‘inputdir/forcing’
- Return type:
netCDF files
The tidal data functions are sourced from the GFDL NWA25 and modified so that:
Converted code for regional-mom6
segment()
classImplemented horizontal subsetting.
Combined all functions of NWA25 into a four function process (in the style of regional-mom6), i.e.,
setup_boundary_tides()
,coords()
,segment.regrid_tides()
, andsegment.encode_tidal_files_and_output()
.
Code credit:
Author(s): GFDL, James Simkins, Rob Cermak, and contributors Year: 2022 Title: "NWA25: Northwest Atlantic 1/25th Degree MOM6 Simulation" Version: N/A Type: Python Functions, Source Code Web Address: https://github.com/jsimkins2/nwa25
- setup_era5(era5_path)#
Setup the ERA5 forcing files for the experiment. This assumes that all of the ERA5 data in the prescribed date range are downloaded. We need the following fields: “2t”, “10u”, “10v”, “sp”, “2d”, “msdwswrf”, “msdwlwrf”, “lsrr”, and “crr”.
- Parameters:
era5_path (str) – Path to the ERA5 forcing files. Specifically, the single-level reanalysis product. For example,
'SOMEPATH/era5/single-levels/reanalysis'
- setup_initial_condition(raw_ic_path, varnames, arakawa_grid='A', vcoord_type='height', rotational_method=RotationMethod.EXPAND_GRID, regridding_method=None)#
Reads the initial condition from files in
raw_ic_path
, interpolates to the model grid, fixes up metadata, and saves back to the input directory.- Parameters:
raw_ic_path (Union[str, Path, list[str]]) – Path(s) to raw initial condition file(s) to read in.
varnames (Dict[str, str]) – Mapping from MOM6 variable/coordinate names to the names in the input dataset. For example,
{'xq': 'lonq', 'yh': 'lath', 'salt': 'so', ...}
.arakawa_grid (Optional[str]) – Arakawa grid staggering type of the initial condition. Either
'A'
(default),'B'
, or'C'
.vcoord_type (Optional[str]) – The type of vertical coordinate used in the forcing files. Either
'height'
or'thickness'
.rotational_method (Optional[RotationMethod]) – The method used to rotate the velocities.
regridding_method (Optional[str]) – The type of regridding method to use. Defaults to self.regridding_method
- setup_ocean_state_boundaries(raw_boundaries_path, varnames, arakawa_grid='A', bathymetry_path=None, rotational_method=RotationMethod.EXPAND_GRID, regridding_method=None, fill_method=None)#
A wrapper for
setup_single_boundary()
. Given a list of up to four cardinal directions, it creates a boundary forcing file for each one. Ensure that the raw boundaries are all saved in the same directory, and that they are named using the formateast_unprocessed.nc
.- Parameters:
raw_boundaries_path (str) – Path to the directory containing the raw boundary forcing files.
varnames (Dict[str, str]) – Mapping from MOM6 variable/coordinate names to the name in the input dataset.
boundaries (List[str]) – List of cardinal directions for which to create boundary forcing files. Default is
["south", "north", "west", "east"]
.arakawa_grid (Optional[str]) – Arakawa grid staggering type of the boundary forcing. Either
'A'
(default),'B'
, or'C'
.bathymetry_path (Optional[str]) – Path to the bathymetry file. Default is
None
, in which case the boundary condition is not masked.rotational_method (Optional[str]) – Method to use for rotating the boundary velocities. Default is
EXPAND_GRID
.regridding_method (Optional[str]) – The type of regridding method to use. Defaults to self.regridding_method
fill_method (Function) – Fill method to use throughout the function. Default is
self.fill_method
- setup_run_directory(surface_forcing=None, using_payu=False, overwrite=False, with_tides=False)#
Set up the run directory for MOM6. Either copy a pre-made set of files, or modify existing files in the ‘rundir’ directory for the experiment.
- Parameters:
surface_forcing (Optional[str]) – Specify the choice of surface forcing, one of:
'jra'
or'era5'
. If not prescribed then constant fluxes are used.using_payu (Optional[bool]) – Whether or not to use payu (payu-org/payu) to run the model. If
True
, a payu configuration file will be created. Default:False
.overwrite (Optional[bool]) – Whether or not to overwrite existing files in the run directory. If
False
(default), will only modify theMOM_layout
file and will not re-copy across the rest of the default files.
- setup_single_boundary(path_to_bc, varnames, orientation, segment_number, arakawa_grid='A', bathymetry_path=None, rotational_method=RotationMethod.EXPAND_GRID, regridding_method=None, fill_method=None)#
Set up a boundary forcing file for a given
orientation
.- Parameters:
path_to_bc (str) – Path to boundary forcing file. Ideally this should be a pre cut-out netCDF file containing only the boundary region and 3 extra boundary points on either side. Users can also provide a large dataset containing their entire domain but this will be slower.
varnames (Dict[str, str]) – Mapping from MOM6 variable/coordinate names to the name in the input dataset.
orientation (str) – Orientation of boundary forcing file, i.e.,
'east'
,'west'
,'north'
, or'south'
.segment_number (int) – Number the segments according to how they’ll be specified in the
MOM_input
.arakawa_grid (Optional[str]) – Arakawa grid staggering type of the boundary forcing. Either
'A'
(default),'B'
, or'C'
.bathymetry_path (str) – Path to the bathymetry file. Default is
None
, in which case the boundary condition is not masked.rotational_method (Optional[str]) – Method to use for rotating the boundary velocities. Default is ‘EXPAND_GRID’.
regridding_method (Optional[str]) – The type of regridding method to use. Defaults to self.regridding_method
fill_method (Function) – Fill method to use throughout the function. Default is
rgd.fill_missing_data
- property tides_boundary_paths#
Finds the tides files from disk, and prints the file paths
- tidy_bathymetry(fill_channels=False, positive_down=False, vertical_coordinate_name='depth', bathymetry=None, write_to_file=True, longitude_coordinate_name='lon', latitude_coordinate_name='lat')#
An auxiliary method for bathymetry used to fix up the metadata and remove inland lakes after regridding the bathymetry. Having
tidy_bathymetry()
as a separate method fromsetup_bathymetry()
allows for the regridding to be done separately, since regridding can be really expensive for large domains.If the bathymetry is already regridded and what is left to be done is fixing the metadata or fill in some channels, then
tidy_bathymetry()
directly can read the existingbathymetry_unfinished.nc
file that should be in the input directory.- Parameters:
fill_channels (Optional[bool]) – Whether to fill in diagonal channels. This removes more narrow inlets, but can also connect extra islands to land. Default:
False
.positive_down (Optional[bool]) – If
False
(default), assume that bathymetry vertical coordinate is positive down, as is the case in GEBCO for example.bathymetry (Optional[xr.Dataset]) – The bathymetry dataset to tidy up. If not provided, it will read the bathymetry from the file
bathymetry_unfinished.nc
in the input directory that was created bysetup_bathymetry()
.
- CrocoDash.rm6.regional_mom6.regional_mom6.generate_rectangular_hgrid(lons, lats)#
Construct a horizontal grid with all the metadata required by MOM6, based on arrays of longitudes (
lons
) and latitudes (lats
) on the supergrid. Here, ‘supergrid’ refers to both cell edges and centres, meaning that there are twice as many points along each axis than for any individual field.Caution
It is assumed the grid’s boundaries are lines of constant latitude and longitude. Rotated grids need to be handled differently.
It is also assumed here that the longitude array values are uniformly spaced.
Ensure both
lons
andlats
are monotonically increasing.- Parameters:
lons (numpy.array) – All longitude points on the supergrid. Must be uniformly spaced.
lats (numpy.array) – All latitude points on the supergrid.
- Returns:
An FMS-compatible horizontal grid (
hgrid
) that includes all required attributes.- Return type:
xarray.Dataset
- CrocoDash.rm6.regional_mom6.regional_mom6.get_glorys_data(longitude_extent, latitude_extent, timerange, segment_name, download_path, modify_existing=True)#
Generates a bash script to download all of the required ocean forcing data.
- Parameters:
longitude_extent (tuple of floats) – Westward and Eastward extents of the segment
latitude_extent (tuple of floats) – Southward and Northward extents of the segment
timerange (tuple of datetime strings) – Start and end of the segment, each in format %Y-%m-%d %H:%M:%S
segment_range (str) – name of the segment (without the
.nc
extension, e.g.,east_unprocessed
)download_path (str) – Location of where the script is saved
modify_existing (bool) – Whether to add to an existing script or start a new one
- Returns:
file path
- CrocoDash.rm6.regional_mom6.regional_mom6.hyperbolictan_thickness_profile(nlayers, ratio, total_depth)#
Generate a hyperbolic tangent thickness profile with
nlayers
vertical layers and total depth oftotal_depth
whose bottom layer is (about)ratio
times larger than the top layer.The thickness profile transitions from the top-layer thickness to the bottom-layer thickness via a hyperbolic tangent proportional to
tanh(2π * (k / (nlayers - 1) - 1 / 2))
, wherek = 0, 1, ..., nlayers - 1
is the layer index withk = 0
corresponding to the top-most layer.The sum of all layer thicknesses is
total_depth
.Positive parameter
ratio
prescribes (approximately) the ratio of the thickness of the bottom-most layer to the top-most layer. The final ratio of the bottom-most layer to the top-most layer ends up a bit different from the prescribedratio
. In particular, the final ratio of the bottom over the top-most layer thickness is(1 + ratio * exp(2π)) / (ratio + exp(2π))
. This slight departure comes about because of the value of the hyperbolic tangent profile at the end-pointstanh(π)
, which is approximately 0.9963 and not 1. Note that becauseexp(2π)
is much greater than 1, the value of the actual ratio is not that different from the prescribed valueratio
, e.g., forratio
values between 1/100 and 100 the final ratio of the bottom-most layer to the top-most layer only departs from the prescribedratio
by ±20%.- Parameters:
nlayers (int) – Number of vertical layers.
ratio (float) – The desired value of the ratio of bottom-most to the top-most layer thickness. Note that the final value of the ratio of bottom-most to the top-most layer thickness ends up
(1 + ratio * exp(2π)) / (ratio + exp(2π))
. Must be positive.total_depth (float) – The total depth of grid, i.e., the sum of all thicknesses.
- Returns:
An array containing the layer thicknesses.
- Return type:
numpy.array
Examples
The spacings for a vertical grid with 20 layers, with maximum depth 1000 meters, and for which the top-most layer is about 4 times thinner than the bottom-most one.
>>> from regional_mom6 import hyperbolictan_thickness_profile >>> nlayers, total_depth = 20, 1000 >>> ratio = 4 >>> dz = hyperbolictan_thickness_profile(nlayers, ratio, total_depth) >>> dz array([20.11183771, 20.2163053 , 20.41767549, 20.80399084, 21.53839043, 22.91063751, 25.3939941 , 29.6384327 , 36.23006369, 45.08430684, 54.91569316, 63.76993631, 70.3615673 , 74.6060059 , 77.08936249, 78.46160957, 79.19600916, 79.58232451, 79.7836947 , 79.88816229]) >>> dz.sum() 1000.0 >>> dz[-1] / dz[0] 3.9721960481753706
If we want the top layer to be thicker then we need to prescribe
ratio < 1
.>>> from regional_mom6 import hyperbolictan_thickness_profile >>> nlayers, total_depth = 20, 1000 >>> ratio = 1/4 >>> dz = hyperbolictan_thickness_profile(nlayers, ratio, total_depth) >>> dz array([79.88816229, 79.7836947 , 79.58232451, 79.19600916, 78.46160957, 77.08936249, 74.6060059 , 70.3615673 , 63.76993631, 54.91569316, 45.08430684, 36.23006369, 29.6384327 , 25.3939941 , 22.91063751, 21.53839043, 20.80399084, 20.41767549, 20.2163053 , 20.11183771]) >>> dz.sum() 1000.0 >>> dz[-1] / dz[0] 0.25174991059652
Now how about a grid with the same total depth as above but with equally-spaced layers.
>>> from regional_mom6 import hyperbolictan_thickness_profile >>> nlayers, total_depth = 20, 1000 >>> ratio = 1 >>> dz = hyperbolictan_thickness_profile(nlayers, ratio, total_depth) >>> dz array([50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50., 50.])
- CrocoDash.rm6.regional_mom6.regional_mom6.longitude_slicer(data, longitude_extent, longitude_coords)#
Slices longitudes while handling periodicity and the ‘seams’, that is the longitude values where the data wraps around in a global domain (for example, longitudes are defined, usually, within domain [0, 360] or [-180, 180]).
The algorithm works in five steps:
Determine whether we need to add or subtract 360 to get the middle of the
longitude_extent
to lie withindata
’s longitude range (herebyold_lon
).Shift the dataset so that its midpoint matches the midpoint of
longitude_extent
(up to a multiple of 360). Now, the modifiedold_lon
does not increase monotonically from West to East since the ‘seam’ has moved.Fix
old_lon
to make it monotonically increasing again. This uses the information we have about the way the dataset was shifted/rolled.Slice the
data
index-wise. We know that|longitude_extent[1] - longitude_extent[0]| / 360
multiplied by the number of discrete longitude points in the global input data gives the number of longitude points in our slice, and we’ve already set the midpoint to be the middle of the target domain.Add back the correct multiple of 360 so the whole domain matches the target.
- Parameters:
data (xarray.Dataset) – The global data you want to slice in longitude.
longitude_extent (Tuple[float, float]) – The target longitudes (in degrees) we want to slice to. Must be in increasing order.
longitude_coords (Union[str, list[str]) – The name or list of names of the longitude coordinates(s) in
data
.
- Returns:
The sliced
data
.- Return type:
xarray.Dataset
- class CrocoDash.rm6.regional_mom6.regional_mom6.segment(*, hgrid, infile, outfolder, varnames, segment_name, orientation, startdate, bathymetry_path=None, arakawa_grid='A', time_units='days', repeat_year_forcing=False)#
Bases:
object
Class to turn raw boundary and tidal segment data into MOM6 boundary and tidal segments.
Boundary segments should only contain the necessary data for that segment. No horizontal chunking is done here, so big fat segments will process slowly.
Data should be at daily temporal resolution, iterating upwards from the provided startdate. Function ignores the time metadata and puts it on Julian calendar.
Note
Only supports z-star (z*) vertical coordinate.
- Parameters:
hgrid (xarray.Dataset) – The horizontal grid used for domain.
infile (Union[str, Path]) – Path to the raw, unprocessed boundary segment.
outfolder (Union[str, Path]) – Path to folder where the model inputs will be stored.
varnames (Dict[str, str]) – Mapping between the variable/dimension names and standard naming convention of this pipeline, e.g.,
{"xq": "longitude, "yh": "latitude", "salt": "salinity", ...}
. Key “tracers” points to nested dictionary of tracers to include in boundary.segment_name (str) – Name of the segment, e.g.,
'segment_001'
.orientation (str) – Cardinal direction (lowercase) of the boundary segment, i.e.,
'east'
,'west'
,'north'
, or'south'
.startdate (str) – The starting date to use in the segment calendar.
arakawa_grid (Optional[str]) – Arakawa grid staggering type of the boundary forcing. Either
'A'
(default),'B'
, or'C'
.time_units (str) – The units used by the raw forcing files, e.g.,
hours
,days
(default).repeat_year_forcing (Optional[bool]) – When
True
the experiment runs with repeat-year forcing. WhenFalse
(default) then inter-annual forcing is used.
- encode_tidal_files_and_output(ds, filename)#
This method:
Expands the dimensions (with the segment name)
Renames some dimensions to be more specific to the segment
Provides an output file encoding
Exports the files.
- Parameters:
self.outfolder (str/path) – The output folder to save the tidal files into
dataset (xarray.Dataset) – The processed tidal dataset
filename (str) – The output file name
- Returns:
Regridded [FILENAME] files in ‘self.outfolder/[filename]_[segmentname].nc’
- Return type:
netCDF files
- General Description:
This tidal data functions are sourced from the GFDL NWA25 and changed in the following ways:
Converted code for regional-mom6 segment class
Implemented horizontal Subsetting
Combined all functions of NWA25 into a four function process (in the style of regional-mom6), i.e.,
setup_boundary_tides()
,coords()
,segment.regrid_tides()
, andsegment.encode_tidal_files_and_output()
.
Code credit:
Author(s): GFDL, James Simkins, Rob Cermak, and contributors Year: 2022 Title: "NWA25: Northwest Atlantic 1/25th Degree MOM6 Simulation" Version: N/A Type: Python Functions, Source Code Web Address: https://github.com/jsimkins2/nwa25
- regrid_tides(tpxo_v, tpxo_u, tpxo_h, times, rotational_method=RotationMethod.EXPAND_GRID, regridding_method='bilinear', fill_method=<function fill_missing_data>)#
Regrids and interpolates the tidal data for MOM6. Steps include:
Read raw tidal data (all constituents)
Perform minor transformations/conversions
Regrid the tidal elevation, and tidal velocity
Encode the output
- General Description:
The tidal data functions sourced from the GFDL’s code above were changed in the following ways:
Converted code for regional-mom6 segment class
Implemented horizontal subsetting
Combined all functions of NWA25 into a four function process (in the style of regional-mom6), i.e.,
setup_boundary_tides()
,coords()
,segment.regrid_tides()
, andsegment.encode_tidal_files_and_output()
.
- Parameters:
infile_td (str) – Raw tidal file/directory.
tpxo_v (xarray.Dataset) – Specific adjusted for MOM6 tpxo datasets (Adjusted with
setup_boundary_tides()
)tpxo_u (xarray.Dataset) – Specific adjusted for MOM6 tpxo datasets (Adjusted with
setup_boundary_tides()
)tpxo_h (xarray.Dataset) – Specific adjusted for MOM6 tpxo datasets (Adjusted with
setup_boundary_tides()
)times (pd.DateRange) – The start date of our model period.
rotational_method (rot.RotationMethod) – The method to use for rotation of the velocities. The default method,
EXPAND_GRID
, works even with non-rotated grids.regridding_method (str) – regridding method to use throughout the function. Default is
'bilinear'
fill_method (Function) – Fill method to use throughout the function. Default is
rgd.fill_missing_data
- Returns:
Regridded tidal velocity and elevation files in ‘inputdir/forcing’
- Return type:
netCDF files
Code credit:
Author(s): GFDL, James Simkins, Rob Cermak, and contributors Year: 2022 Title: "NWA25: Northwest Atlantic 1/25th Degree MOM6 Simulation" Version: N/A Type: Python Functions, Source Code Web Address: https://github.com/jsimkins2/nwa25
- regrid_velocity_tracers(rotational_method=RotationMethod.EXPAND_GRID, regridding_method='bilinear', fill_method=<function fill_missing_data>)#
Cut out and interpolate the velocities and tracers.
- Parameters:
rotational_method (rot.RotationMethod) – The method to use for rotation of the velocities. Currently, the default method,
EXPAND_GRID
, works even with non-rotated grids.regridding_method (str) – regridding method to use throughout the function. Default is
'bilinear'
fill_method (Function) – Fill method to use throughout the function. Default is
rgd.fill_missing_data
CrocoDash.rm6.regional_mom6.regridding module#
Custom-built helper methods to regrid the boundary conditions and ensure proper encoding for MOM6.
Steps:
Initial Regridding -> Find the boundary of the
hgrid
, and regrid the forcing variables to that boundary. Call (initial_regridding
) and then use thexesmf.Regridder
with any datasets you need.Ensure temperatures are in Celsius.
Fill in NaNs. This step is important for MOM6 (
fill_missing_data
) -> This diverges betweenFor tides, we split the tides into an amplitude and a phase
In some cases, here is a great place to rotate the velocities to match a curved grid (tidal_velocity), velocity is also a good place to do this.
We then add the time coordinate
- Add several depth-related coordinates to fields that are not related to the ocean’s surface (like, e.g., surface wind stress).
Add a
dz
variable in layer thicknessSome metadata issues later on
Now we do up the metadata
Rename variables to var_segment_num
(For fields with vertical dimension) Rename the vertical coordinate of the variable to
nz_segment_num_var
.(For fields with vertical dimension) Declare this new vertical coordinate as a increasing series of integers.
Re-add the “perpendicular” dimension.
….Add layer thickness
dz
to the forcing fields with vertical dimension.Add to encoding_dict a fill value(_FillValue) and zlib, dtype, for time, lat lon, … and each variable (no type needed though).
- CrocoDash.rm6.regional_mom6.regridding.add_or_update_time_dim(ds: Dataset, times) Dataset #
Add the time dimension to the dataset, in tides case can be one time step.
- Parameters:
ds (xr.Dataset) – The dataset to add the time dimension to
times (list, np.Array, xr.DataArray) – The list of times
- Returns:
The dataset with the time dimension added
- Return type:
(xr.Dataset)
- CrocoDash.rm6.regional_mom6.regridding.add_secondary_dimension(ds: Dataset, var: str, coords, segment_name: str, to_beginning=False) Dataset #
Add the perpendiciular dimension to the dataset, even if it is only one value since it is required.
- Parameters:
ds (xr.Dataset) – The dataset to add the perpendicular dimension to
var (str) – The variable to add the perpendicular dimension to
coords (xr.Dataset) – The output xarray Dataset from the coords function. Contains information required to add the perpendicular dimension.
segment_name (str) – The segment name
to_beginning (bool, optional) – Whether to add the perpendicular dimension to the beginning or to the selected position, by default False
Returns
(xr.Dataset): The dataset with the vertical dimension added
- CrocoDash.rm6.regional_mom6.regridding.coords(hgrid: Dataset, orientation: str, segment_name: str, coords_at_t_points=False, angle_variable_name='angle_dx') Dataset #
Allows us to call the coords for use in the
xesmf.Regridder
in theregrid_tides()
function.self.coords
gives us the subset of thehgrid
based on the orientation.- Parameters:
hgrid (xr.Dataset) – The horizontal grid dataset.
orientation (str) – The orientation of the boundary.
segment_name (str) – The name of the segment.
coords_at_t_points (bool, optional) – Whether to return the boundary t-points instead of the q/u/v of a general boundary for rotation. Default:
False
.
- Returns:
The correct coordinate space for the orientation
- Return type:
xr.Dataset
Code credit:
Author(s): GFDL, James Simkins, Rob Cermak, and contributors Year: 2022 Title: "NWA25: Northwest Atlantic 1/25th Degree MOM6 Simulation" Version: N/A Type: Python Functions, Source Code Web Address: https://github.com/jsimkins2/nwa25
- CrocoDash.rm6.regional_mom6.regridding.create_regridder(forcing_variables: Dataset, output_grid: Dataset, outfile: Path = None, method: str = 'bilinear', locstream_out: bool = True, periodic: bool = False) Regridder #
Basic regridder for any forcing variables. This is essentially a wrapper for the xesmf regridder with some default parameter choices.
- Parameters:
forcing_variables (xr.Dataset) – The dataset of the forcing variables.
output_grid (xr.Dataset) – The dataset of the output grid. This is the boundary of the
hgrid
outfile (Path, optional) – The path to the output file for weights; default: None
method (str, optional) – The regridding method; default:
"bilinear"
locstream_out (bool, optional) – Whether to output the locstream; default:
True
periodic (bool, optional) – Whether the grid is periodic; default:
False
- Returns:
The regridding object
- Return type:
xe.Regridder
- CrocoDash.rm6.regional_mom6.regridding.fill_missing_data(ds: Dataset, xdim: str = 'locations', zdim: str = 'z', fill: str = 'b') Dataset #
Fill in missing values.
- Parameters:
ds (xr.Dataset) – The dataset to be filled in
z_dim_name (str) – The name of the
z
dimension
- Returns:
The filled dataset
- Return type:
xr.Dataset
Code credit:
Author(s): GFDL, James Simkins, Rob Cermak, and contributors Year: 2022 Title: "NWA25: Northwest Atlantic 1/25th Degree MOM6 Simulation" Version: N/A Type: Python Functions, Source Code Web Address: https://github.com/jsimkins2/nwa25
- CrocoDash.rm6.regional_mom6.regridding.generate_dz(ds: Dataset, z_dim_name: str) Dataset #
Generate the vertical coordinate spacing.
- Parameters:
ds (xr.Dataset) – The dataset from which we extract the vertical coordinate.
z_dim_name (str) – The name of the vertical coordinate.
- Returns
(xr.Dataset): The vertical spacing variable.
- CrocoDash.rm6.regional_mom6.regridding.generate_encoding(ds: Dataset, encoding_dict, default_fill_value=9.969209968386869e+36) dict #
Generate the encoding dictionary for the dataset.
- Parameters:
ds (xr.Dataset) – The dataset to generate the encoding for
encoding_dict (dict) – The starting encoding dict with some specifications needed for time and other vars, this will be updated with encodings in this function
default_fill_value (float, optional) – The default fill value; default: 1.0e20
- Returns:
The encoding dictionary
- Return type:
(dict)
- CrocoDash.rm6.regional_mom6.regridding.generate_layer_thickness(ds: Dataset, var: str, segment_name: str, old_vert_coord_name: str) Dataset #
Generate Layer Thickness Variable, needed for vars with vertical dimensions :param ds: The dataset to generate the layer thickness for :type ds: xr.Dataset :param var: The variable to generate the layer thickness for :type var: str :param segment_name: The segment name :type segment_name: str :param old_vert_coord_name: The old vertical coordinate name :type old_vert_coord_name: str
- Returns:
The dataset with the layer thickness variable added
- Return type:
xr.Dataset
- CrocoDash.rm6.regional_mom6.regridding.get_boundary_mask(bathy: Dataset, side: str, minimum_depth=0, x_dim_name='nx', y_dim_name='ny') ndarray #
Mask out the boundary conditions based on the bathymetry. We don’t want to have boundary conditions on land. :param bathy: The bathymetry dataset :type bathy: xr.Dataset :param side: The side of the boundary, “north”, “south”, “east”, or “west” :type side: str :param minimum_depth: The minimum depth to consider land, by default 0 :type minimum_depth: float, optional
- Returns:
The boundary mask
- Return type:
np.ndarray
- CrocoDash.rm6.regional_mom6.regridding.get_hgrid_arakawa_c_points(hgrid: Dataset, point_type='t') Dataset #
Get the Arakawa C points from the hgrid.
Credit: Method originally by Fred Castruccio.
- Parameters:
hgrid (xr.Dataset) – The hgrid dataset
- Returns:
The specific points x, y, & point indexes
- Return type:
xr.Dataset
- CrocoDash.rm6.regional_mom6.regridding.mask_dataset(ds: Dataset, bathymetry: Dataset, orientation, y_dim_name='ny', x_dim_name='nx', fill_value=-1e+20) Dataset #
This function masks the dataset to the provided bathymetry. If bathymetry is not provided, it fills all NaNs with 0. :param ds: The dataset to mask :type ds: xr.Dataset :param bathymetry: The bathymetry dataset :type bathymetry: xr.Dataset :param orientation: The orientation of the boundary :type orientation: str :param fill_value: The value land points should be filled with :type fill_value: float
- CrocoDash.rm6.regional_mom6.regridding.vertical_coordinate_encoding(ds: Dataset, var: str, segment_name: str, old_vert_coord_name: str) Dataset #
Rename vertical coordinate to nz[additional-text] then change it to regular increments
- Parameters:
ds (xr.Dataset) – The dataset to rename the vertical coordinate in
var (str) – The variable to rename the vertical coordinate in
segment_name (str) – The segment name
old_vert_coord_name (str) – The old vertical coordinate name
CrocoDash.rm6.regional_mom6.rotation module#
- class CrocoDash.rm6.regional_mom6.rotation.RotationMethod(value, names=None, *, module=None, qualname=None, type=None, start=1, boundary=None)#
Bases:
Enum
Prescribes the rotational method to be used in boundary conditions when the grid does not have coordinates along lines of constant longitude-latitude. regional-mom6 main class passes this
Enum
toregrid_tides()
and toregrid_velocity_tracers()
.- EXPAND_GRID#
This method is used with the basis that we can find the angles at the q-u-v points by pretending we have another row/column of the
hgrid
with the same distances as the t-point to u/v points in the actual grid then use the four points to calculate the angle. This method replicates exactly what MOM6 does.- Type:
int
- GIVEN_ANGLE#
Expects a pre-given angle called
angle_dx
.- Type:
int
- NO_ROTATION#
Grid is along lines of constant latitude-longitude and therefore no rotation is required.
- Type:
int
- EXPAND_GRID = 1#
- GIVEN_ANGLE = 2#
- NO_ROTATION = 3#
- CrocoDash.rm6.regional_mom6.rotation.create_expanded_hgrid(hgrid: Dataset, expansion_width=1) Dataset #
Adds an additional boundary to the hgrid to allow for the calculation of the
angle_dx
for the boundary points usingmom6_angle_calculation_method()
.
- CrocoDash.rm6.regional_mom6.rotation.get_rotation_angle(rotational_method: RotationMethod, hgrid: Dataset, orientation=None)#
Returns the rotation angle - with the assumption of degrees - based on the rotational method and provided hgrid, if orientation & coords are provided, it will assume the boundary is requested.
- Parameters:
rotational_method (RotationMethod) – The rotational method to use
hgrid (xr.Dataset) – The hgrid dataset
orientation (xr.Dataset) – The orientation, which also lets us now that we are on a boundary
- Returns:
angle in degrees
- Return type:
xr.DataArray
- CrocoDash.rm6.regional_mom6.rotation.initialize_grid_rotation_angle(hgrid: Dataset) DataArray #
Calculate the
angle_dx
in degrees from the truex
direction (parallel to latitude) counter-clockwise and return as a dataarray. (Mimics MOM6 angle calculation functionmom6_angle_calculation_method()
)- Parameters:
hgrid (xr.Dataset) – The hgrid dataset
- Returns:
The t-point angles
- Return type:
xr.DataArray
- CrocoDash.rm6.regional_mom6.rotation.initialize_grid_rotation_angles_using_expanded_hgrid(hgrid: Dataset) Dataset #
Calculate the
angle_dx
(in degrees) from the truex
direction (parallel to latitude) counter-clockwise and return as a dataarray.- Parameters:
hgrid (xr.Dataset) – The hgrid dataset
- Returns:
The t-point angles
- Return type:
xr.DataArray
- CrocoDash.rm6.regional_mom6.rotation.modulo_around_point(x, x0, L)#
Returns the modulo-\(L\) value of \(x\) within the interval \([x_0 - L/2, x_0 + L/2]\). If \(L ≤ 0\), then method returns \(x\).
(Adapted from MOM6 code; mom-ocean/MOM6)
- Parameters:
x (xr.DataArray) – Value(s) to which to apply modulo arithmetic
x0 (xr.DataArray) – Center(s) of modulo range
L (float) – Modulo range width
- Returns:
x
shifted by an integer multiple ofL
to be closer tox0
, i.e., within the interval[x0 - L/2, x0 + L/2]
- Return type:
float
- CrocoDash.rm6.regional_mom6.rotation.mom6_angle_calculation_method(len_lon, top_left: DataArray, top_right: DataArray, bottom_left: DataArray, bottom_right: DataArray, point: DataArray) DataArray #
Calculate the angle of the grid point’s local x-direction compared to East-West direction using the MOM6 method adapted from: mom-ocean/MOM6
Note: this is exactly the same as the angle of the grid point’s local y-direction compared to North-South direction.
This method can handle vectorized computations.
- Parameters:
len_lon (float) – The extent of the longitude of the regional domain (in degrees).
top_left (xr.DataArray) – The four points around the point to calculate the angle from the
hgrid
; requires both anx`
andy
component (both in degrees).top_right (xr.DataArray) – The four points around the point to calculate the angle from the
hgrid
; requires both anx`
andy
component (both in degrees).bottom_left (xr.DataArray) – The four points around the point to calculate the angle from the
hgrid
; requires both anx`
andy
component (both in degrees).bottom_right (xr.DataArray) – The four points around the point to calculate the angle from the
hgrid
; requires both anx`
andy
component (both in degrees).point (xr.DataArray) – The point to calculate the angle from the
hgrid
- Returns:
The angle of the grid point’s local
x
-direction compared to East-West direction.- Return type:
xr.DataArray
CrocoDash.rm6.regional_mom6.utils module#
- CrocoDash.rm6.regional_mom6.utils.angle_between(v1, v2, v3)#
Return the angle
v2
-v1
-v3
(in radians), wherev1
,v2
,v3
are 3-vectors. That is, the angle that is formed between vectorsv2 - v1
and vectorv3 - v1
.Example
>>> from regional_mom6.utils import angle_between >>> v1 = (0, 0, 1) >>> v2 = (1, 0, 0) >>> v3 = (0, 1, 0) >>> angle_between(v1, v2, v3) 1.5707963267948966 >>> from numpy import rad2deg >>> rad2deg(angle_between(v1, v2, v3)) 90.0
- CrocoDash.rm6.regional_mom6.utils.ap2ep(uc, vc)#
Convert complex tidal u and v to tidal ellipse.
Adapted from ap2ep.m for Matlab. Copyright notice:
Authorship: The author retains the copyright of this program, while you are welcome to use and distribute it as long as you credit the author properly and respect the program name itself. Particularly, you are expected to retain the original author's name in this original version or any of its modified version that you might make. You are also expected not to essentially change the name of the programs except for adding possible extension for your own version you might create, e.g. ap2ep_xx is acceptable. Any suggestions are welcome and enjoy my program(s)! Author Info: Zhigang Xu, Ph.D. (pronounced as Tsi Gahng Hsu) Research Scientist Coastal Circulation Bedford Institute of Oceanography 1 Challenge Dr. P.O. Box 1006 Phone (902) 426-2307 (o) Dartmouth, Nova Scotia Fax (902) 426-7827 CANADA B2Y 4A2 email xuz@dfo-mpo.gc.ca Release Date: Nov. 2000, Revised on May. 2002 to adopt Foreman's northern semi major axis convention.
- Parameters:
uc – complex tidal u velocity
vc – complex tidal v velocity
- Returns:
(semi-major axis, eccentricity, inclination [radians], phase [radians])
- CrocoDash.rm6.regional_mom6.utils.ep2ap(SEMA, ECC, INC, PHA)#
Convert tidal ellipse to real u and v amplitude and phase.
Adapted from ep2ap.m for Matlab. Copyright notice:
Authorship: The author of this program retains the copyright of this program, while you are welcome to use and distribute this program as long as you credit the author properly and respect the program name itself. Particularly, you are expected to retain the original author's name in this original version of the program or any of its modified version that you might make. You are also expected not to essentially change the name of the programs except for adding possible extension for your own version you might create, e.g. app2ep_xx is acceptable. Any suggestions are welcome and enjoy my program(s)! Author Info: Zhigang Xu, Ph.D. (pronounced as Tsi Gahng Hsu) Research Scientist Coastal Circulation Bedford Institute of Oceanography 1 Challenge Dr. P.O. Box 1006 Phone (902) 426-2307 (o) Dartmouth, Nova Scotia Fax (902) 426-7827 CANADA B2Y 4A2 email xuz@dfo-mpo.gc.ca Release Date: Nov. 2000
- Parameters:
SEMA – semi-major axis
ECC – eccentricity
INC – inclination [radians]
PHA – phase [radians]
- Returns:
(u amplitude, u phase [radians], v amplitude, v phase [radians])
- CrocoDash.rm6.regional_mom6.utils.find_files_by_pattern(paths: list, patterns: list, error_message=None) list #
Function searchs paths for patterns and returns the list of the file paths with that pattern
- CrocoDash.rm6.regional_mom6.utils.get_edge(ds, edge, x_name=None, y_name=None)#
- CrocoDash.rm6.regional_mom6.utils.is_rectilinear_hgrid(hgrid: Dataset, rtol: float = 0.001) bool #
Check if the
hgrid
is a rectilinear grid by comparing the first and last rows and columns of the tlon and tlat arrays.From
mom6_bathy.grid.is_rectangular
by Alper (Altuntas).- Parameters:
hgrid (xarray.Dataset) – The horizontal grid dataset.
rtol (float) – Relative tolerance. Default is 1e-3.
- CrocoDash.rm6.regional_mom6.utils.latlon_to_cartesian(lat, lon, R=1)#
Convert latitude and longitude (in degrees) to Cartesian coordinates on a sphere of radius
R
. By defaultR = 1
.- Parameters:
lat (float) – Latitude (in degrees).
lon (float) – Longitude (in degrees).
R (float) – The radius of the sphere; default: 1.
- Returns:
Tuple with the Cartesian coordinates
x, y, z
- Return type:
tuple
Examples
Find the Cartesian coordinates that correspond to point with
(lat, lon) = (0, 0)
on a sphere with unit radius.>>> from regional_mom6.utils import latlon_to_cartesian >>> latlon_to_cartesian(0, 0) (1.0, 0.0, 0.0)
Now let’s do the same on a sphere with Earth’s radius
>>> from regional_mom6.utils import latlon_to_cartesian >>> R = 6371e3 >>> latlon_to_cartesian(0, 0, R) (6371000.0, 0.0, 0.0)
- CrocoDash.rm6.regional_mom6.utils.quadrilateral_area(v1, v2, v3, v4)#
Return the area of a spherical quadrilateral on the unit sphere that has vertices on the 3-vectors
v1
,v2
,v3
,v4
(counter-clockwise orientation is implied). The area is computed via the excess of the sum of the spherical angles of the quadrilateral from 2π.Example
Calculate the area that corresponds to half the Northern hemisphere of a sphere of radius R. This should be 1/4 of the sphere’s total area, that is π R2.
>>> from regional_mom6.utils import quadrilateral_area, latlon_to_cartesian >>> R = 434.3 >>> v1 = latlon_to_cartesian(0, 0, R) >>> v2 = latlon_to_cartesian(0, 90, R) >>> v3 = latlon_to_cartesian(90, 0, R) >>> v4 = latlon_to_cartesian(0, -90, R) >>> quadrilateral_area(v1, v2, v3, v4) 592556.1793298927 >>> from numpy import pi >>> quadrilateral_area(v1, v2, v3, v4) == pi * R**2 True
- CrocoDash.rm6.regional_mom6.utils.quadrilateral_areas(lat, lon, R=1)#
Return the area of spherical quadrilaterals on a sphere of radius
R
. By default,R = 1
. The quadrilaterals are formed by constant latitude and longitude lines on thelat
-lon
grid provided.- Parameters:
lat (numpy.array) – Array of latitude points (in degrees).
lon (numpy.array) – Array of longitude points (in degrees).
R (float) – The radius of the sphere; default: 1.
- Returns:
Array with the areas of the quadrilaterals defined by the
lat
-lon
grid provided. If the providedlat
,lon
arrays are of dimension m \(\times\) n then returned areas array is of dimension (m-1) \(\times\) (n-1).- Return type:
numpy.array
Example
Let’s construct a lat-lon grid on the sphere with 60 degree spacing. Then we compute the areas of each grid cell and confirm that the sum of the areas gives us the total area of the sphere.
>>> from regional_mom6.utils import quadrilateral_areas >>> import numpy as np >>> λ = np.linspace(0, 360, 7) >>> φ = np.linspace(-90, 90, 4) >>> lon, lat = np.meshgrid(λ, φ) >>> lon array([[ 0., 60., 120., 180., 240., 300., 360.], [ 0., 60., 120., 180., 240., 300., 360.], [ 0., 60., 120., 180., 240., 300., 360.], [ 0., 60., 120., 180., 240., 300., 360.]]) >>> lat array([[-90., -90., -90., -90., -90., -90., -90.], [-30., -30., -30., -30., -30., -30., -30.], [ 30., 30., 30., 30., 30., 30., 30.], [ 90., 90., 90., 90., 90., 90., 90.]]) >>> R = 6371e3 >>> areas = quadrilateral_areas(lat, lon, R) >>> areas array([[1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13], [4.56284230e+13, 4.56284230e+13, 4.56284230e+13, 4.56284230e+13, 4.56284230e+13, 4.56284230e+13], [1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13, 1.96911611e+13]]) >>> np.isclose(areas.sum(), 4 * np.pi * R**2, atol=np.finfo(areas.dtype).eps) True
- CrocoDash.rm6.regional_mom6.utils.rotate(u, v, radian_angle)#
Rotate the velocities counter-clockwise by an angle
radian_angle
(in radians). (Same asrotate_complex()
.)- Parameters:
u (xarray.DataArray) – The \(u\)-component of the velocity.
v (xarray.DataArray) – The \(v\)-component of the velocity.
radian_angle (xarray.DataArray) – The angle of rotation (in radians).
- Returns:
The rotated \(u\) and \(v\) components of the velocity.
- Return type:
Tuple[xarray.DataArray, xarray.DataArray]
- CrocoDash.rm6.regional_mom6.utils.rotate_complex(u, v, radian_angle)#
Rotate velocities counter-clockwise by angle
radian_angle
(in radians) using complex number math (Same asrotate()
.)- Parameters:
u (xarray.DataArray) – The \(u\)-component of the velocity.
v (xarray.DataArray) – The \(v\)-component of the velocity.
radian_angle (xarray.DataArray) – The angle of rotation (in radians).
- Returns:
The rotated \(u\) and \(v\) components of the velocity.
- Return type:
Tuple[xarray.DataArray, xarray.DataArray]
- CrocoDash.rm6.regional_mom6.utils.setup_logger(name: str, set_handler=False, log_level=20) Logger #
Setup general configuration for a logger.
- CrocoDash.rm6.regional_mom6.utils.vecdot(v1, v2)#
Return the dot product of vectors
v1
andv2
.v1
andv2
can be either numpy vectors or numpy.ndarrays in which case the last dimension is considered the dimension over which the dot product is taken.