diff --git a/docs/api.rst b/docs/api.rst index 38bac0eab..b87214e15 100644 --- a/docs/api.rst +++ b/docs/api.rst @@ -464,6 +464,7 @@ UxDataArray UxDataArray.remap.nearest_neighbor UxDataArray.remap.inverse_distance_weighted UxDataArray.remap.bilinear + UxDataArray.remap.to_rectilinear UxDataset @@ -477,6 +478,7 @@ UxDataset UxDataset.remap.nearest_neighbor UxDataset.remap.inverse_distance_weighted UxDataset.remap.bilinear + UxDataset.remap.to_rectilinear Mathematical Operators diff --git a/docs/generated/uxarray.UxDataArray.remap.to_rectilinear.rst b/docs/generated/uxarray.UxDataArray.remap.to_rectilinear.rst new file mode 100644 index 000000000..c71ba609d --- /dev/null +++ b/docs/generated/uxarray.UxDataArray.remap.to_rectilinear.rst @@ -0,0 +1,6 @@ +uxarray.UxDataArray.remap.to_rectilinear +======================================== + +.. currentmodule:: uxarray + +.. autoaccessormethod:: UxDataArray.remap.to_rectilinear diff --git a/docs/generated/uxarray.UxDataset.remap.to_rectilinear.rst b/docs/generated/uxarray.UxDataset.remap.to_rectilinear.rst new file mode 100644 index 000000000..bb972149d --- /dev/null +++ b/docs/generated/uxarray.UxDataset.remap.to_rectilinear.rst @@ -0,0 +1,6 @@ +uxarray.UxDataset.remap.to_rectilinear +====================================== + +.. currentmodule:: uxarray + +.. autoaccessormethod:: UxDataset.remap.to_rectilinear diff --git a/docs/user-guide/remapping.ipynb b/docs/user-guide/remapping.ipynb index e1d41c775..8c71dc6f4 100644 --- a/docs/user-guide/remapping.ipynb +++ b/docs/user-guide/remapping.ipynb @@ -22,7 +22,7 @@ "\n", "UXarray uses its native remapping backend by default. If [YAC](https://dkrz-sw.gitlab-pages.dkrz.de/yac/) is installed with its `yac.core` Python bindings, `.remap(...)`, `.remap.nearest_neighbor(...)`, and `.remap.bilinear(...)` can use `backend=\"yac\"` to route remapping through YAC.\n", "\n", - "When `backend=\"yac\"`, the `yac_method` parameter selects the YAC interpolation method. Supported values are `nnn`, `average`, and `conservative`. `inverse_distance_weighted()` remains UXarray-only, and `bilinear(..., backend=\"yac\")` is a convenience wrapper for `yac_method=\"average\"`. For conservative face-centered remapping, use the generic `.remap(..., backend=\"yac\", yac_method=\"conservative\")` entrypoint.\n", + "When `backend=\"yac\"`, the `yac_method` parameter selects the YAC interpolation method. Supported values are `nnn`, `average`, and `conservative`. `inverse_distance_weighted()` remains UXarray-only, and `bilinear(..., backend=\"yac\")` is a convenience wrapper for `yac_method=\"average\"`. For conservative face-centered remapping, use the generic `.remap(..., backend=\"yac\", yac_method=\"conservative\")` entrypoint. To remap onto a rectilinear latitude-longitude target and get regular xarray output, use `.remap.to_rectilinear(lon=..., lat=..., backend=\"yac\")`.\n", "\n", "See the [YAC documentation](https://dkrz-sw.gitlab-pages.dkrz.de/yac/) and [YAC installation guide](https://dkrz-sw.gitlab-pages.dkrz.de/yac/d1/d9f/installing_yac.html) for build instructions, including enabling the Python bindings.\n", "```\n" diff --git a/test/test_remap.py b/test/test_remap.py index 836fd517e..8cda30ce7 100644 --- a/test/test_remap.py +++ b/test/test_remap.py @@ -1,6 +1,7 @@ import numpy as np import numpy.testing as nt import pytest +import xarray as xr import uxarray as ux from uxarray.core.dataarray import UxDataArray @@ -196,6 +197,70 @@ def test_dataset_remap_preserves_coords(gridpath, datasetpath): assert "time" in ds_out.coords +def test_to_rectilinear_native_backend(): + """Rectilinear remap returns plain xarray output on lat/lon axes.""" + grid = ux.Grid.from_structured( + lon=np.asarray([0.0, 90.0]), + lat=np.asarray([0.0, 45.0]), + ) + da = UxDataArray( + np.asarray([1.0, 2.0, 3.0, 4.0]), + dims=["n_face"], + coords={ + "n_face": [0, 1, 2, 3], + "face_lon": ( + "n_face", + grid.face_lon.values, + {"standard_name": "longitude", "units": "degrees_east"}, + ), + "face_lat": ( + "n_face", + grid.face_lat.values, + {"standard_name": "latitude", "units": "degrees_north"}, + ), + }, + uxgrid=grid, + ) + lon = xr.DataArray( + [0.0, 90.0], + dims=["lon"], + attrs={"axis": "X", "units": "degrees_east"}, + ) + lat = xr.DataArray( + [0.0, 45.0], + dims=["lat"], + attrs={"axis": "Y", "units": "degrees_north"}, + ) + + out = da.remap.to_rectilinear(lon=lon, lat=lat, backend="uxarray") + + assert isinstance(out, xr.DataArray) + assert out.dims == ("lat", "lon") + assert out.shape == (2, 2) + nt.assert_array_equal(out.values, np.asarray([[1.0, 2.0], [3.0, 4.0]])) + assert out["lon"].attrs["units"] == "degrees_east" + assert out["lat"].attrs["units"] == "degrees_north" + + +def test_reshape_to_rectilinear_accepts_xarray_dataset(): + """Rectilinear reshaping accepts an already-open xarray Dataset.""" + from uxarray.remap.structured import ( + _normalize_rectilinear_target, + _reshape_to_rectilinear, + ) + + spec = _normalize_rectilinear_target( + lon=np.asarray([0.0, 90.0]), lat=np.asarray([0.0, 45.0]) + ) + ds = xr.Dataset({"var": ("n_face", np.asarray([1.0, 2.0, 3.0, 4.0]))}) + + out = _reshape_to_rectilinear(ds, spec) + + assert isinstance(out, xr.Dataset) + assert out["var"].dims == ("lat", "lon") + nt.assert_array_equal(out["var"].values, np.asarray([[1.0, 2.0], [3.0, 4.0]])) + + # ------------------------------------------------------------ # Bilinear tests # ------------------------------------------------------------ diff --git a/test/test_remap_yac.py b/test/test_remap_yac.py index 7d71e2656..51f2b9e33 100644 --- a/test/test_remap_yac.py +++ b/test/test_remap_yac.py @@ -256,3 +256,87 @@ def test_yac_batched_remap_with_fractional_mask(): assert out.shape == da.shape np.testing.assert_array_equal(out.values, da.values) + + +def test_yac_to_rectilinear_node_remap(): + verts = np.array([(0.0, 0.0), (90.0, 0.0), (0.0, 45.0)]) + grid = ux.open_grid(verts) + da = ux.UxDataArray( + np.asarray([1.0, 2.0, 3.0]), + dims=["n_node"], + coords={"n_node": [0, 1, 2]}, + uxgrid=grid, + ) + + out = da.remap.to_rectilinear( + lon=np.asarray([0.0, 90.0]), + lat=np.asarray([0.0, 45.0]), + backend="yac", + yac_method="nnn", + yac_options={"n": 1}, + ) + + assert out.dims == ("lat", "lon") + assert out.shape == (2, 2) + np.testing.assert_array_equal(out.values, np.asarray([[1.0, 3.0], [2.0, 3.0]])) + + +def test_yac_to_rectilinear_preserves_extra_dimensions(): + verts = np.array([(0.0, 0.0), (90.0, 0.0), (0.0, 45.0)]) + grid = ux.open_grid(verts) + da = ux.UxDataArray( + np.asarray([[1.0, 2.0, 3.0], [10.0, 20.0, 30.0]]), + dims=["time", "n_node"], + coords={"time": [0, 1], "n_node": [0, 1, 2]}, + uxgrid=grid, + ) + + out = da.remap.to_rectilinear( + lon=np.asarray([0.0, 90.0]), + lat=np.asarray([0.0, 45.0]), + backend="yac", + yac_method="nnn", + yac_options={"n": 1}, + ) + + assert out.dims == ("time", "lat", "lon") + assert out.shape == (2, 2, 2) + np.testing.assert_array_equal( + out.values, + np.asarray([[[1.0, 3.0], [2.0, 3.0]], [[10.0, 30.0], [20.0, 30.0]]]), + ) + + +def test_yac_to_rectilinear_preserves_valid_nonspatial_coords(): + verts = np.array([(0.0, 0.0), (90.0, 0.0), (0.0, 45.0)]) + grid = ux.open_grid(verts) + da = ux.UxDataArray( + np.asarray([[1.0, 2.0, 3.0], [10.0, 20.0, 30.0]]), + dims=["time", "n_node"], + coords={ + "time": [0, 1], + "n_node": [0, 1, 2], + "run": ("time", np.asarray([100, 101])), + "experiment": "demo", + "node_lon": ( + "n_node", + np.asarray([0.0, 90.0, 0.0]), + {"standard_name": "longitude", "units": "degrees_east"}, + ), + }, + uxgrid=grid, + ) + + out = da.remap.to_rectilinear( + lon=np.asarray([0.0, 90.0]), + lat=np.asarray([0.0, 45.0]), + backend="yac", + yac_method="nnn", + yac_options={"n": 1}, + ) + + assert "run" in out.coords + assert "experiment" in out.coords + assert "node_lon" not in out.coords + np.testing.assert_array_equal(out["run"].values, np.asarray([100, 101])) + assert out["experiment"].item() == "demo" diff --git a/uxarray/remap/accessor.py b/uxarray/remap/accessor.py index b1332599d..6eb0c24f3 100644 --- a/uxarray/remap/accessor.py +++ b/uxarray/remap/accessor.py @@ -34,6 +34,7 @@ def __repr__(self) -> str: + "Supported methods:\n" + " • nearest_neighbor(destination_grid, remap_to='faces')\n" + " • inverse_distance_weighted(destination_grid, remap_to='faces', power=2, k=8)\n" + + " • to_rectilinear(lon, lat, backend='yac')\n" ) def __call__( @@ -175,6 +176,64 @@ def inverse_distance_weighted( self.ux_obj, destination_grid, remap_to, power, k ) + def to_rectilinear( + self, + lon, + lat, + backend: str = "yac", + yac_method: str | None = "nnn", + yac_options: dict | None = None, + **kwargs, + ): + """ + Remap onto a rectilinear longitude/latitude grid. + + This convenience method targets 1-D longitude and latitude coordinate + arrays and returns a plain xarray object with ``lat`` and ``lon`` axes, + making the output suitable for downstream structured-grid workflows. + + Parameters + ---------- + lon : array-like or xarray.DataArray + 1-D target longitude cell-center coordinate in degrees. + lat : array-like or xarray.DataArray + 1-D target latitude cell-center coordinate in degrees. + backend : {'uxarray', 'yac'}, default='yac' + Remapping backend to use. The YAC backend uses YAC's rectilinear + grid support directly. The UXarray backend builds a temporary + structured destination grid and applies native nearest-neighbor + remapping before reshaping the result to latitude/longitude axes. + yac_method : {'nnn', 'average', 'conservative'}, optional + YAC interpolation method. Defaults to ``'nnn'`` because nearest-neighbor + works for node-, edge-, and face-centered source data. ``'conservative'`` + requires face-centered source data. + yac_options : dict, optional + YAC interpolation configuration options forwarded to the selected + YAC method. + + Returns + ------- + xarray.DataArray or xarray.Dataset + Remapped data with the source spatial dimension replaced by the + provided latitude and longitude dimensions. + """ + + _validate_backend(backend) + if backend == "yac": + from uxarray.remap.yac import _yac_remap_to_rectilinear + + return _yac_remap_to_rectilinear( + self.ux_obj, + lon, + lat, + yac_method or "nnn", + yac_options or {}, + ) + + from uxarray.remap.structured import _native_remap_to_rectilinear + + return _native_remap_to_rectilinear(self.ux_obj, lon, lat) + def bilinear( self, destination_grid: Grid, diff --git a/uxarray/remap/structured.py b/uxarray/remap/structured.py new file mode 100644 index 000000000..bc036a11d --- /dev/null +++ b/uxarray/remap/structured.py @@ -0,0 +1,216 @@ +from __future__ import annotations + +from dataclasses import dataclass + +import numpy as np +import xarray as xr + + +@dataclass(frozen=True) +class RectilinearGridSpec: + """Normalized 1-D longitude/latitude target grid metadata.""" + + lon: xr.DataArray + lat: xr.DataArray + lon_edges: np.ndarray + lat_edges: np.ndarray + cyclic_lon: bool + + @property + def shape(self) -> tuple[int, int]: + return self.lat.size, self.lon.size + + @property + def size(self) -> int: + return self.lat.size * self.lon.size + + @property + def lon_dim(self) -> str: + return self.lon.dims[0] + + @property + def lat_dim(self) -> str: + return self.lat.dims[0] + + @property + def lon_name(self) -> str: + return self.lon.name or self.lon_dim + + @property + def lat_name(self) -> str: + return self.lat.name or self.lat_dim + + @property + def lon_rad(self) -> np.ndarray: + return np.deg2rad(np.asarray(self.lon.values, dtype=np.float64)) + + @property + def lat_rad(self) -> np.ndarray: + return np.deg2rad(np.asarray(self.lat.values, dtype=np.float64)) + + @property + def lon_edges_rad(self) -> np.ndarray: + return np.deg2rad(self.lon_edges) + + @property + def lat_edges_rad(self) -> np.ndarray: + return np.deg2rad(self.lat_edges) + + def flattened_centers_rad(self) -> tuple[np.ndarray, np.ndarray]: + lon_2d, lat_2d = np.meshgrid(self.lon_rad, self.lat_rad) + return lon_2d.ravel(), lat_2d.ravel() + + +def _as_1d_coord(coord, default_name: str) -> xr.DataArray: + if isinstance(coord, xr.DataArray): + values = np.asarray(coord.values, dtype=np.float64) + if coord.ndim != 1: + raise ValueError( + f"Rectilinear {default_name!r} coordinate must be 1-D, " + f"got {coord.ndim}-D." + ) + dim = coord.dims[0] + name = coord.name or default_name + attrs = coord.attrs.copy() + else: + values = np.asarray(coord, dtype=np.float64) + if values.ndim != 1: + raise ValueError( + f"Rectilinear {default_name!r} coordinate must be 1-D, " + f"got {values.ndim}-D." + ) + dim = default_name + name = default_name + attrs = {} + + if values.size < 2: + raise ValueError( + f"Rectilinear {default_name!r} coordinate must contain at least two values." + ) + return xr.DataArray(values, dims=(dim,), name=name, attrs=attrs) + + +def _centers_to_edges(values: np.ndarray, coord_name: str) -> np.ndarray: + diffs = np.diff(values) + if np.any(diffs == 0): + raise ValueError(f"Rectilinear {coord_name!r} coordinate contains duplicates.") + if not (np.all(diffs > 0) or np.all(diffs < 0)): + raise ValueError( + f"Rectilinear {coord_name!r} coordinate must be strictly monotonic." + ) + + edges = np.empty(values.size + 1, dtype=np.float64) + edges[1:-1] = values[:-1] + 0.5 * diffs + edges[0] = values[0] - 0.5 * diffs[0] + edges[-1] = values[-1] + 0.5 * diffs[-1] + return edges + + +def _looks_global_lon(edges: np.ndarray) -> bool: + return np.isclose(abs(edges[-1] - edges[0]), 360.0, rtol=0.0, atol=1.0e-8) + + +def _normalize_rectilinear_target(lon, lat) -> RectilinearGridSpec: + lon_coord = _as_1d_coord(lon, "lon") + lat_coord = _as_1d_coord(lat, "lat") + lon_edges = _centers_to_edges(np.asarray(lon_coord.values), lon_coord.name or "lon") + lat_edges = _centers_to_edges(np.asarray(lat_coord.values), lat_coord.name or "lat") + + return RectilinearGridSpec( + lon=lon_coord, + lat=lat_coord, + lon_edges=lon_edges, + lat_edges=lat_edges, + cyclic_lon=_looks_global_lon(lon_edges), + ) + + +def _preserve_valid_coords( + da: xr.DataArray, + dropped_dim: str, + output_dims: tuple[str, ...] | list[str], +) -> dict[str, xr.DataArray]: + """Keep only coords that remain valid after replacing ``dropped_dim``.""" + + output_dims = set(output_dims) + return { + name: coord + for name, coord in da.coords.items() + if dropped_dim not in coord.dims and set(coord.dims).issubset(output_dims) + } + + +def _reshape_array_to_rectilinear( + da: xr.DataArray, spec: RectilinearGridSpec +) -> xr.DataArray: + if "n_face" not in da.dims: + return xr.DataArray(da) + + axis = da.get_axis_num("n_face") + if da.sizes["n_face"] != spec.size: + raise ValueError( + "Cannot reshape remapped data to the requested rectilinear grid. " + f"Expected {spec.size} face values, got {da.sizes['n_face']}." + ) + + shape = da.shape[:axis] + spec.shape + da.shape[axis + 1 :] + dims = da.dims[:axis] + (spec.lat_dim, spec.lon_dim) + da.dims[axis + 1 :] + coords = { + name: coord + for name, coord in da.coords.items() + if "n_face" not in coord.dims + and name not in {spec.lat_name, spec.lon_name, spec.lat_dim, spec.lon_dim} + } + coords[spec.lat_name] = spec.lat + coords[spec.lon_name] = spec.lon + + return xr.DataArray( + np.asarray(da.values).reshape(shape), + dims=dims, + coords=coords, + name=da.name, + attrs=da.attrs, + ) + + +def _reshape_to_rectilinear(obj, spec: RectilinearGridSpec): + """Convert flattened ``n_face`` remap output to lat/lon-shaped xarray.""" + + if isinstance(obj, xr.DataArray): + return _reshape_array_to_rectilinear(obj, spec) + + if isinstance(obj, xr.Dataset): + xr_obj = obj + elif hasattr(obj, "to_xarray"): + xr_obj = obj.to_xarray() + else: + xr_obj = xr.Dataset(obj) + data_vars = { + name: _reshape_array_to_rectilinear(da, spec) + for name, da in xr_obj.data_vars.items() + } + coords = { + name: coord + for name, coord in xr_obj.coords.items() + if "n_face" not in coord.dims + and name not in {spec.lat_name, spec.lon_name, spec.lat_dim, spec.lon_dim} + } + coords[spec.lat_name] = spec.lat + coords[spec.lon_name] = spec.lon + return xr.Dataset(data_vars=data_vars, coords=coords, attrs=xr_obj.attrs) + + +def _native_remap_to_rectilinear(source, lon, lat): + """Remap to a rectilinear target through UXarray's native NN backend.""" + + from uxarray.grid import Grid + from uxarray.remap.nearest_neighbor import _nearest_neighbor_remap + + spec = _normalize_rectilinear_target(lon, lat) + destination_grid = Grid.from_structured(lon=spec.lon.values, lat=spec.lat.values) + remapped = _nearest_neighbor_remap( + source, + destination_grid=destination_grid, + remap_to="faces", + ) + return _reshape_to_rectilinear(remapped, spec) diff --git a/uxarray/remap/yac.py b/uxarray/remap/yac.py index 4bbe3899b..285a89a65 100644 --- a/uxarray/remap/yac.py +++ b/uxarray/remap/yac.py @@ -10,8 +10,15 @@ from uuid import uuid4 import numpy as np +import xarray as xr import uxarray.core.dataarray +from uxarray.remap.structured import ( + RectilinearGridSpec, + _normalize_rectilinear_target, + _preserve_valid_coords, + _reshape_to_rectilinear, +) from uxarray.remap.utils import ( LABEL_TO_COORD, _assert_dimension, @@ -125,6 +132,34 @@ def _get_lon_lat(grid, dim: str) -> tuple[np.ndarray, np.ndarray]: ) +def _is_rectilinear_target(grid) -> bool: + return isinstance(grid, RectilinearGridSpec) + + +def _register_yac_grid(yac_core, name: str, grid, dim: str, define_edges: bool): + if _is_rectilinear_target(grid): + if dim != "n_face": + raise ValueError( + "Rectilinear YAC targets are face-centered; use remap_to='faces'." + ) + yac_grid = yac_core.BasicGrid.reg_2d_new( + name, + grid.lon_edges_rad, + grid.lat_edges_rad, + cyclic=[grid.cyclic_lon, False], + ) + lon, lat = grid.flattened_centers_rad() + return yac_grid, lon, lat + + yac_grid = yac_core.BasicGrid.from_uxgrid( + name, + grid, + def_edges=define_edges, + ) + lon, lat = _get_lon_lat(grid, dim) + return yac_grid, lon, lat + + def _coerce_enum(enum_type, value: Any): if not isinstance(value, str): return value @@ -160,20 +195,21 @@ def __init__( self._src_location = _get_location(yac_core, src_dim) self._tgt_location = _get_location(yac_core, tgt_dim) - define_edges = "n_edge" in (src_dim, tgt_dim) unique = uuid4().hex - self._src_grid = yac_core.BasicGrid.from_uxgrid( + self._src_grid, src_lon, src_lat = _register_yac_grid( + yac_core, f"uxarray_src_{unique}", src_grid, - def_edges=define_edges, + src_dim, + define_edges=src_dim == "n_edge", ) - self._tgt_grid = yac_core.BasicGrid.from_uxgrid( + self._tgt_grid, tgt_lon, tgt_lat = _register_yac_grid( + yac_core, f"uxarray_tgt_{unique}", tgt_grid, - def_edges=define_edges, + tgt_dim, + define_edges=tgt_dim == "n_edge", ) - src_lon, src_lat = _get_lon_lat(src_grid, src_dim) - tgt_lon, tgt_lat = _get_lon_lat(tgt_grid, tgt_dim) self._src_field = yac_core.InterpField( self._src_grid.add_coordinates(self._src_location, src_lon, src_lat) @@ -405,3 +441,76 @@ def _yac_remap(source, destination_grid, remap_to: str, yac_method: str, yac_kwa source, remapped_vars, destination_grid, remap_to ) return ds_remapped[name] if is_da else ds_remapped + + +def _yac_remap_to_rectilinear(source, lon, lat, yac_method: str, yac_kwargs): + """Remap a UXarray object to a 1-D lon/lat target through YAC.""" + + destination_dim = "n_face" + target = _normalize_rectilinear_target(lon, lat) + options = _normalize_yac_method(yac_method) + options.kwargs.update(yac_kwargs or {}) + ds, is_da, name = _to_dataset(source) + dims_to_remap = _get_remap_dims(ds) + + if options.method == "conservative": + non_face_src = dims_to_remap - {"n_face"} + if non_face_src: + raise ValueError( + "YAC conservative remapping to a rectilinear grid requires all " + "source data to be face-centered (dimension 'n_face'). " + f"Found non-face source dimension(s): {non_face_src}. " + "Use yac_method='nnn' for node- or edge-centered data." + ) + + remappers: dict[str, _YacRemapper] = { + src_dim: _YacRemapper( + ds.uxgrid, + target, + src_dim, + destination_dim, + options.method, + options.kwargs, + ) + for src_dim in dims_to_remap + } + remapped_vars = {} + + for var_name, da in ds.data_vars.items(): + src_dim = next((d for d in da.dims if d in dims_to_remap), None) + if src_dim is None: + remapped_vars[var_name] = xr.DataArray(da) + continue + + other_dims = [d for d in da.dims if d != src_dim] + da_t = da.transpose(*other_dims, src_dim) + src_values = np.asarray(da_t.values) + flat_src = src_values.reshape(-1, src_values.shape[-1]) + frac_masks = yac_kwargs.get("frac_masks") + frac_mask = ( + frac_masks.get(var_name) + if isinstance(frac_masks, dict) and var_name in frac_masks + else yac_kwargs.get("frac_mask") + ) + flat_frac_mask = None + if frac_mask is not None: + flat_frac_mask = _prepare_frac_mask(frac_mask, da_t, src_values, src_dim) + + remapper = remappers[src_dim] + out_flat = remapper.apply(flat_src, frac_mask=flat_frac_mask) + + out_shape = src_values.shape[:-1] + (remapper._tgt_size,) + out_values = out_flat.reshape(out_shape) + output_dims = tuple(other_dims + [destination_dim]) + coords = _preserve_valid_coords(da, src_dim, output_dims) + remapped_vars[var_name] = xr.DataArray( + out_values, + dims=output_dims, + coords=coords, + name=da.name, + attrs=da.attrs, + ) + + ds_remapped = xr.Dataset(remapped_vars, attrs=ds.attrs) + rectilinear = _reshape_to_rectilinear(ds_remapped, target) + return rectilinear[name] if is_da else rectilinear