Skip to content
Closed
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
65 changes: 65 additions & 0 deletions xarray/tests/test_interp.py
Original file line number Diff line number Diff line change
Expand Up @@ -662,3 +662,68 @@ def test_datetime_interp_noerror():
coords={"time": pd.date_range("01-01-2001", periods=50, freq="H")},
)
a.interp(x=xi, time=xi.time) # should not raise an error


def interp_1d(data, dim, vals):
selectors = {dim: vals}

old_coord = data[dim]

lower = data.sel(**selectors, method='ffill')
upper = data.sel(**selectors, method='bfill')
lower_x = data[dim].sel(**selectors, method='ffill')
upper_x = data[dim].sel(**selectors, method='bfill')
weight = np.abs(vals - lower_x)/np.abs(upper_x-lower_x)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would lean towards using the built in abs() here, which is slightly more generic than the NumPy function.

# needed if upper and lower coordinates are identical
weight = weight.fillna(1.0)
ans = lower * (1-weight) + upper * weight
return ans


def interp_nd(data, coords):
for dim in coords:
data = interp_1d(data, dim, coords[dim])
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The downside of this approach (indexing separately along each dimension) is that it is potentially much more expensive than doing all indexing at once, e.g., if the result of indexing along multiple dimensions is only a single point.

I think it would be interesting to see if this could be done by indexing all dimensions at once 2**len(coords) times, with all combinations of bfill/ffill along all dimensions.

This would require modifying sel to support a dict of options for method indicating different indexing methods along different axes, along the lines of the proposal from #3223 for tolerance.

return data


def test_interp_nd():
assert False


@pytest.mark.parametrize('dim',[
'x', 'y'
])
def test_new_interp(dim):

data = xr.DataArray(np.arange(10), dims=[dim]).assign_coords({dim: np.arange(10)})
new_x = xr.DataArray(np.arange(9) + .5, dims=['xp'])
ans = interp_1d(data, dim, new_x)
np.testing.assert_allclose(new_x, ans.values)


def test_interp_1d_nd_targ():
data = xr.DataArray(np.arange(10), dims=['x'], coords={'x': np.arange(10)})
new_x_2d = xr.DataArray(np.tile(np.r_[:9] + .5, (3, 1)), dims=['y', 'xp'])
ans = interp_1d(data, 'x', new_x_2d)
np.testing.assert_allclose(new_x_2d, ans.values)


def test_interp_1d_nd_targ():
# interp with actual coords should be an isel.
data = xr.DataArray(np.arange(10), dims=['x'], coords={'x': np.arange(10)})
new_x_2d = xr.DataArray(np.tile(np.r_[:9], (3, 1)), dims=['y', 'xp'])
ans = interp_1d(data, 'x', new_x_2d)
np.testing.assert_allclose(new_x_2d, ans.values)

def test_sel_nd():
npdata = np.tile(np.arange(10), (5, 1))
data = xr.DataArray(npdata, dims=['y', 'x'],
coords={'x': np.r_[:10], 'y': np.r_[:5]})
idx = xr.DataArray(npdata, dims=['z', 'xp'])

ans = data.sel(x=idx)
assert set(ans.dims) == {'z', 'y', 'xp'}
print(ans)