Skip to content

Commit 5669ae3

Browse files
davidhasselldavidhassell
committed
dev
1 parent 750875e commit 5669ae3

File tree

4 files changed

+200
-398
lines changed

4 files changed

+200
-398
lines changed

‎cf/data/dask_regrid.py‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -488,7 +488,8 @@ def _regrid(
488488
continue
489489

490490
w = data[i0:i1]
491-
D_j = 1 - w[mask].sum()
491+
D_j = w[~mask].sum()
492+
# D_j = 1 - w[mask].sum()
492493
w = w / D_j
493494
w[mask] = 0
494495
data[i0:i1] = w

‎cf/regrid/regrid.py‎

Lines changed: 42 additions & 21 deletions
Original file line numberDiff line numberDiff line change
@@ -92,6 +92,10 @@ class Grid:
9292
# for a 2-d grid for which one of the dimensions is a size 2 dummy
9393
# dimension.
9494
dummy_size_2_dimension: bool = False
95+
# Whether or not the grid is a structured grid.
96+
is_grid: bool = False
97+
# Whether or not the grid is a UGRID mesh.
98+
is_mesh: bool = False
9599
# The location on a UGRID mesh topology of the grid cells. An
96100
# empty string means that the grid is not a UGRID mesh
97101
# topology. E.g. '' or 'face'.
@@ -402,12 +406,6 @@ def regrid(
402406
raise ValueError(
403407
f"{method!r} regridding is only available for 2-d regridding"
404408
)
405-
# elif method in ("linear", "bilinear"):
406-
# if src_grid.mesh_location or dst_grid.mesh_location:
407-
# raise ValueError(
408-
# f"{method!r} regridding is only available when neither the "
409-
# "source nor destination grid is a UGRID mesh"
410-
# )
411409
elif method in ("nearest_dtos", "nearest_stod"):
412410
if not has_coordinate_arrays(src_grid) and not has_coordinate_arrays(
413411
dst_grid
@@ -416,8 +414,23 @@ def regrid(
416414
f"{method!r} regridding is only available when both the "
417415
"source and destination grids have coordinate arrays"
418416
)
419-
420-
417+
418+
if method == "nearest_dtos" and not (
419+
(src_grid.is_mesh and dst_grid.is_mesh)
420+
or (src_grid.is_grid and dst_grid.is_grid)
421+
):
422+
raise ValueError(
423+
f"{method!r} regridding is (at the moment) only available "
424+
"when neither or both the source and destination grids are "
425+
"a UGRID mesh"
426+
)
427+
428+
elif cartesian and (src_grid.is_mesh or dst_grid.is_mesh):
429+
raise ValueError(
430+
"Cartesian regridding is (at the moment) only available when "
431+
"neither the source nor destination grid is a UGRID mesh"
432+
)
433+
421434
if create_regrid_operator:
422435
# ------------------------------------------------------------
423436
# Create a new regrid operator
@@ -474,7 +487,7 @@ def regrid(
474487
if is_log_level_debug(logger):
475488
logger.debug(
476489
f"Source ESMF Grid:\n{src_esmpy_grid}\nDestination ESMF Grid:\n{dst_esmpy_grid}\n"
477-
) # pragma: no cover
490+
) # pragma: no cover
478491

479492
esmpy_regrid_operator = [] if return_esmpy_regrid_operator else None
480493

@@ -590,7 +603,7 @@ def regrid(
590603
# ----------------------------------------------------------------
591604
update_data(src, regridded_data, src_grid)
592605

593-
if coord_sys == "spherical" and not dst_grid.mesh_location:
606+
if coord_sys == "spherical" and dst_grid.is_grid:
594607
# Set the cyclicity of the longitude axis of the new field
595608
key, x = src.dimension_coordinate("X", default=(None, None), item=True)
596609
if x is not None and x.Units.equivalent(Units("degrees")):
@@ -1018,6 +1031,8 @@ def spherical_grid(f, name=None, method=None, cyclic=None, axes=None):
10181031
data_axes = f.get_data_axes()
10191032
axis_indices = [data_axes.index(key) for key in axis_keys]
10201033

1034+
is_mesh = bool(mesh_location)
1035+
10211036
return Grid(
10221037
axis_keys=axis_keys,
10231038
axis_indices=axis_indices,
@@ -1031,6 +1046,8 @@ def spherical_grid(f, name=None, method=None, cyclic=None, axes=None):
10311046
coord_sys="spherical",
10321047
method=method,
10331048
name=name,
1049+
is_mesh=is_mesh,
1050+
is_grid=not is_mesh,
10341051
mesh_location=mesh_location,
10351052
domain_topology=domain_topology,
10361053
)
@@ -1193,6 +1210,8 @@ def Cartesian_grid(f, name=None, method=None, axes=None):
11931210
if mesh_location:
11941211
regridding_dimensionality += 1
11951212

1213+
is_mesh = bool(mesh_location)
1214+
11961215
return Grid(
11971216
axis_keys=axis_keys,
11981217
axis_indices=axis_indices,
@@ -1207,6 +1226,8 @@ def Cartesian_grid(f, name=None, method=None, axes=None):
12071226
method=method,
12081227
name=name,
12091228
dummy_size_2_dimension=dummy_size_2_dimension,
1229+
is_mesh=is_mesh,
1230+
is_grid=not is_mesh,
12101231
mesh_location=mesh_location,
12111232
domain_topology=domain_topology,
12121233
)
@@ -1422,7 +1443,7 @@ def create_esmpy_grid(grid, mask=None):
14221443
The `esmpy.Grid` or `esmpy.Mesh` derived from *grid*.
14231444
14241445
"""
1425-
if grid.mesh_location:
1446+
if grid.is_mesh:
14261447
# Create an `esmpy.Mesh`
14271448
return create_esmpy_mesh(grid, mask)
14281449

@@ -1690,13 +1711,13 @@ def create_esmpy_mesh(grid, mask=None):
16901711
node_count = node_ids.size
16911712
node_owners = np.zeros(node_count)
16921713

1693-
# print (
1694-
# 'index=',index,
1695-
# 'index.shape=',index.shape,
1696-
# 'node_ids=',node_ids,
1697-
# 'node_ids.shape=',node_ids.shape,
1698-
# 'node_coords=',node_coords,
1699-
# 'node_coords.shape=',node_coords.shape)
1714+
# print (
1715+
# 'index=',index,
1716+
# 'index.shape=',index.shape,
1717+
# 'node_ids=',node_ids,
1718+
# 'node_ids.shape=',node_ids.shape,
1719+
# 'node_coords=',node_coords,
1720+
# 'node_coords.shape=',node_coords.shape)
17001721
# Add nodes. This must be done before `add_elements`.
17011722
esmpy_mesh.add_nodes(
17021723
node_count=node_count,
@@ -1862,10 +1883,10 @@ def create_esmpy_weights(
18621883
dst_meshloc = esmpy.api.constants.MeshLoc.NODE
18631884

18641885
src_esmpy_field = esmpy.Field(
1865-
src_esmpy_grid, "src", meshloc=src_meshloc
1886+
src_esmpy_grid, name="src", meshloc=src_meshloc
18661887
)
18671888
dst_esmpy_field = esmpy.Field(
1868-
dst_esmpy_grid, "dst", meshloc=dst_meshloc
1889+
dst_esmpy_grid, name="dst", meshloc=dst_meshloc
18691890
)
18701891

18711892
mask_values = np.array([0], dtype="int32")
@@ -2234,7 +2255,7 @@ def update_coordinates(src, dst, src_grid, dst_grid):
22342255

22352256
# Copy domain topology and cell connectivity constructs from the
22362257
# destination grid
2237-
if dst_grid.mesh_location:
2258+
if dst_grid.is_mesh:
22382259
for key, topology in dst.constructs(
22392260
filter_by_type=("domain_topology", "cell_connectivity"),
22402261
filter_by_axis=dst_axis_keys,

‎cf/regrid/regridoperator.py‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -623,13 +623,13 @@ def tosparse(self):
623623
# Note: It is much more efficient to access 'weights.indptr'
624624
# and 'weights.data' directly, rather than iterating
625625
# over rows of 'weights' and using 'weights.getrow'.
626-
# count_nonzero = np.count_nonzero
626+
# count_nonzero = np.count_nonzero
627627
indptr = weights.indptr.tolist()
628628
data = weights.data
629629
for j, (i0, i1) in enumerate(zip(indptr[:-1], indptr[1:])):
630-
# print (data[i0:i1])
631-
if not data[i0:i1].size: #count_nonzero(data[i0:i1]):
632-
# print (j, repr(data[i0:i1]))
630+
# print (data[i0:i1])
631+
if not data[i0:i1].size: # count_nonzero(data[i0:i1]):
632+
# print (j, repr(data[i0:i1]))
633633
dst_mask[j] = True
634634

635635
if not dst_mask.any():

0 commit comments

Comments
 (0)