Skip to content

Commit

Permalink
dev
Browse files Browse the repository at this point in the history
  • Loading branch information
@davidhassell
davidhassell committed Sep 17, 2023
1 parent dd5df6a commit 0c84d90
Showing 1 changed file with 69 additions and 121 deletions.
190 changes: 69 additions & 121 deletions cf/regrid/regrid.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -55,24 +55,22 @@ class Grid:
# expected by `Data._regrid`. E.g. ['domainaxis3', 'domainaxis2']
axis_keys: list = field(default_factory=list)
# The positions of the regrid axes, in the order expected by
# `Data._regrid`. E.g. [3, 2]
# `Data._regrid`. E.g. [3, 2] or [3]
axis_indices: list = field(default_factory=list)
# The value of the *src_axes* or *dst_axes* parameter, as
# appropriate.
axes: Any = None
# The domain axis identifiers of the regridding axes. E.g. {'X':
# 'domainaxis2', 'Y': 'domainaxis1'} or ['domainaxis1',
# 'domainaxis2']
# 'domainaxis2'], or {'M': 'domainaxis3'}
axes: Any = None
# The shape of the regridding axes, in the same order as the
# 'axis_keys' attribute. E.g. (96, 73) or (1243,)
shape: tuple = None
# The regrid axis coordinates, in the order expected by
# `esmpy.Grid`. If the coordinates are 2-d (or more) then the axis
# order of each coordinate object must be as expected by
# `esmpy.Grid`. E.g. (73, 96)
# `esmpy`. If the coordinates are 2-d (or more) then the axis
# order of each coordinate object must be as expected by `esmpy`.
coords: list = field(default_factory=list)
# Only used if `mesh` is False. The regrid axis coordinate bounds,
# in the order expected by `esmpy.Grid`. If the coordinates are
# 2-d (or more) then the axis order of each bounds object must be
# as expected by `esmpy.Grid`.
# The regrid axis coordinate bounds, in the order expected by
# `esmpy`. If the coordinates are 2-d (or more) then the axis
# order of each bounds object must be as expected by `esmpy`.
bounds: list = field(default_factory=list)
# Only used if `mesh` is False. For spherical regridding, whether
# or not the longitude axis is cyclic.
Expand Down Expand Up @@ -934,15 +932,22 @@ def spherical_grid(f, name=None, method=None, cyclic=None, axes=None):
coords[dim] = coords[dim].transpose(esmpy_order)

# Set cyclicity of X axis
if cyclic is None:
if mesh:
cyclic = None
elif cyclic is None:
cyclic = f.iscyclic(x_axis)

else:
cyclic = bool(cylcic)

axis_keys = [y_axis, x_axis]
axis_sizes = (y_size, x_size)
shape = (y_size, x_size)
if mesh:
axis_keys = axis_keys[0]
axis_sizes = axis_sizes[0]

shape = axis_sizes[0]
axes = {"X": x_axis, "Y": y_axis}
else:
axes = {"M": x_axis}

if f.construct_type == "domain":
axis_indices = [0, 1]
else:
Expand All @@ -961,11 +966,11 @@ def spherical_grid(f, name=None, method=None, cyclic=None, axes=None):
return Grid(
axis_keys=axis_keys,
axis_indices=axis_indices,
axes={"X": x_axis, "Y": y_axis},
shape=axis_sizes,
axes=axes,
shape=shape,
coords=coords,
bounds=get_bounds(method, coords, mesh_location),
cyclic=bool(cyclic),
cyclic=cyclic,
coord_sys="spherical",
method=method,
name=name,
Expand All @@ -974,85 +979,8 @@ def spherical_grid(f, name=None, method=None, cyclic=None, axes=None):
mesh_location=mesh_location
)

#def spherical_mesh(f, name=None, method=None):
# """Get latitude and longitude coordinate information.
#
# Retrieve the latitude and longitude coordinates of a field, as
# well as some associated information. If 1-d lat/lon coordinates
# are found then these are returned. Otherwise if 2-d lat/lon
# coordinates found then these are returned. TODOUGRID
#
# .. versionadded:: UGRIDVER
#
# .. seealso:: `Cartesian_mesh`
#
# :Parameters:
#
# f: `Field` or `Domain`
# The construct from which to get the grid information.
#
# name: `str`
# A name to identify the mesh.
#
# method: `str`
# The regridding method.
#
# :Returns:
#
# `Grid`
# The mesh topology definition.
#
# """
# dt_key, domain_topology = f.domain_topology(item=True, default=(None, None))
# if domain_topology is None:
# raise ValueError(
# "Could not find 1-d latitude and longitude mesh coordinates: "
# "Missing domain topology construct"
# )
#
# mesh_axis = f.get_data_axes(dt_key)
#
# # Look for 1-d X and Y auxiliary coordinates
# lon = f.auxiliary_coordinate("X", filter_by_axis=mesh_axis, axis_mode='exact', default=None)
# lat = f.auxiliary_coordinate("Y", filter_by_axis=mesh_axis, axis_mode='exact', default=None)
#
# if lon is None or not (lon.Units.islongitude or lon.get_property('standard_name', None) == 'grid_longitude'):
# raise ValueError(
# "Could not find 1-d longitude mesh coordinates"
# )
#
# if lat is None or not (lat.Units.islatitude or lat.get_property('standard_name', None) == 'grid_latitude'):
# raise ValueError(
# "Could not find 1-d latitude mesh coordinates"
# )
#
# if f.construct_type == "domain":
# axis_indices = [0]
# else:
# # Make sure that the data array spans the regridding
# # axis. This might change 'f' in-place.
# data_axes = f.get_data_axes()
# if mesh_axis[0] not in f.get_data_axes():
# f.insert_dimension(mesh_axis[0], position=-1, inplace=True)
# data_axes = f.get_data_axes()
#
# # Find the index of the regridding axis
# axis_indices = [data_axes.index(mesh_axis]]
#
# return Grid(
# mesh=True,
# axis_keys=mesh_axis,
# axis_indices=axis_indices,
# shape=lon.shape,
# coords=[lon, lat], # esmpy order
# coord_sys="spherical",
# method=method,
# name=name,
# domain_topology=domain_topology,
# )


def Cartesian_grid(f, name=None, method=None, axes=None, mesh=False):

def Cartesian_grid(f, name=None, method=None, axes=None):
"""Retrieve the specified Cartesian dimension coordinates of the
field and their corresponding keys.
Expand All @@ -1073,17 +1001,12 @@ def Cartesian_grid(f, name=None, method=None, axes=None, mesh=False):
Specifiers for the dimension coordinates to be
retrieved. See `cf.Field.domain_axes` for details.
mesh: TODOUGRID
:Returns:
`Grid`
The grid definition.
"""
if mesh:
return Cartesian_mesh(f, name=name, method=method, axes=axes)

if not axes:
if name == "source":
raise ValueError(
Expand Down Expand Up @@ -1116,6 +1039,16 @@ def Cartesian_grid(f, name=None, method=None, axes=None, mesh=False):
axis_keys.append(key)
axis_sizes.append(domain_axis.size)

# mesh = False
# domain_topology, mesh_axis, mesh_location = ffff(f, name)
# if mesh_axis is not None:
# if [mesh_axis] == axis_keys:
# mesh = True
# elif mesh_axis in axis_keys:
# raise ValueError(
# "TODOUGRID Can't regrid a combo of mesh and non-mesh axes"
# )

if f.construct_type == "domain":
axis_indices = list(range(len(axis_keys)))
else:
Expand All @@ -1131,16 +1064,25 @@ def Cartesian_grid(f, name=None, method=None, axes=None, mesh=False):
data_axes = f.get_data_axes()
axis_indices = [data_axes.index(key) for key in axis_keys]

coords = []
for key in axis_keys[::-1]:
coord = f.dimension_coordinate(filter_by_axis=(key,), default=None)
if coord is None:
raise ValueError(
f"No unique {name} dimension coordinate for domain axis "
f"{key!r}."
)
mesh = False
domain_topology, mesh_axis, mesh_location = ffff(f, name)

coords.append(coord)
cyclic = False
coords = []
if domain_topology is not None:
# TODOUGRID: check for the mesh spanning the two axes
# ... somehow. Also maybe set mesh=True and cyclic=None.
#
else:
for key in axis_keys[::-1]:
coord = f.dimension_coordinate(filter_by_axis=(key,), default=None)
if coord is None:
raise ValueError(
f"No unique {name} dimension coordinate for domain axis "
f"{key!r}."
)

coords.append(coord)

bounds = get_bounds(method, coords, mesh_location)

Expand All @@ -1166,11 +1108,13 @@ def Cartesian_grid(f, name=None, method=None, axes=None, mesh=False):
shape=tuple(axis_sizes),
coords=coords,
bounds=bounds,
cyclic=False,
cyclic=cyclic,
coord_sys="Cartesian",
method=method,
name=name,
dummy_size_2_dimension=dummy_size_2_dimension,
mesh=mesh,
domain_topology = domain_topology,
)

def Cartesian_mesh(f, name=None, method=None):
Expand Down Expand Up @@ -1369,14 +1313,15 @@ def esmpy_initialise():

return esmpy.Manager(debug=bool(regrid_logging()))

def create_esmpy_grid(grid=None, mask=None):
"""Create an `esmpy` Grid or Mesh
def create_esmpy_grid(grid, mask=None):
"""Create an `esmpy.Grid` or `esmpy.Mesh`.
.. versionadded:: 3.14.0
:Parameters:
grid: `Grid`
TODOUGRID
mask: array_like, optional
The grid mask. If `None` (the default) then there are no
Expand All @@ -1387,6 +1332,7 @@ def create_esmpy_grid(grid=None, mask=None):
:Returns:
`esmpy.Grid` or `esmpy.Mesh`
TODOUGRID
"""
if g.mesh:
Expand Down Expand Up @@ -1594,8 +1540,8 @@ def create_esmpy_grid(grid=None, mask=None):
return esmpy_grid


def create_esmpy_mesh(grid=None, mask=None):
"""Create an `esmpy` Mesh.
def create_esmpy_mesh(grid, mask=None):
"""Create an `esmpy.Mesh`.
.. versionadded:: UGRIDVER
Expand All @@ -1604,6 +1550,7 @@ def create_esmpy_mesh(grid=None, mask=None):
:Parameters:
grid: `Grid`
TODOUGRID
mask: array_like, optional
The mesh mask. If `None` (the default) then there are no
Expand All @@ -1613,6 +1560,7 @@ def create_esmpy_mesh(grid=None, mask=None):
:Returns:
`esmpy.Mesh`
TODOUGRID
"""
if grid.coord_sys == "spherical":
Expand Down Expand Up @@ -2275,8 +2223,8 @@ def update_non_coordinates(
def ffff(f, name):
"""TODOUGRID"""
key, domain_topology = f.domain_topology(item=True, default=(None, None))
if domain_topology is None:
raise ValueError("TODOUGRID")
if domain_topology is None:
return (None, None, None)

mesh_location = domain_topology.get_cell(None)
if mesh_location != "face":
Expand All @@ -2285,4 +2233,4 @@ def ffff(f, name):
f"of {mesh_location!r} cells"
)

return domain_topology, f.get_data_axes(key)[0], mesh_location
return (domain_topology, f.get_data_axes(key)[0], mesh_location )

0 comments on commit 0c84d90

Please sign in to comment.