Skip to content
forked from pydata/xarray

Commit

Permalink
Fix overflow/underflow warnings in interpolate_na
Browse files Browse the repository at this point in the history
These were being triggered by casting datetime64[ns] to float32.
We now rescale the co-ordinate before interpolating, except for
nearest-neighbour interpolation. The rescaling can change the
nearest neighbour, and so is avoided in this case to maintain
pandas compatibility.
  • Loading branch information
dcherian committed Aug 20, 2018
1 parent 7186644 commit 1f1ec52
Show file tree
Hide file tree
Showing 2 changed files with 19 additions and 9 deletions.
14 changes: 12 additions & 2 deletions xarray/core/missing.py
Original file line number Diff line number Diff line change
Expand Up @@ -153,7 +153,7 @@ def _apply_over_vars_with_dim(func, self, dim=None, **kwargs):
return ds


def get_clean_interp_index(arr, dim, use_coordinate=True, **kwargs):
def get_clean_interp_index(arr, dim, method, use_coordinate=True, **kwargs):
'''get index to use for x values in interpolation.
If use_coordinate is True, the coordinate that shares the name of the
Expand All @@ -176,6 +176,15 @@ def get_clean_interp_index(arr, dim, use_coordinate=True, **kwargs):
# raise if index cannot be cast to a float (e.g. MultiIndex)
try:
index = index.values.astype(np.float64)
if method != 'nearest':
# rescale index to avoid overflow/underflow
# The division can change the nearest-neighbour
# when compared to pandas (which does not divide).
# Let's keep that compatitibility
index = (index - index.min())
if len(index) > 1:
index /= index.std()

except (TypeError, ValueError):
# pandas raises a TypeError
# xarray/nuppy raise a ValueError
Expand All @@ -202,7 +211,8 @@ def interp_na(self, dim=None, use_coordinate=True, method='linear', limit=None,
valids = _get_valid_fill_mask(self, dim, limit)

# method
index = get_clean_interp_index(self, dim, use_coordinate=use_coordinate,
index = get_clean_interp_index(self, dim, method=method,
use_coordinate=use_coordinate,
**kwargs)
interp_class, kwargs = _get_interpolator(method, **kwargs)
interpolator = partial(func_interpolate_na, interp_class, **kwargs)
Expand Down
14 changes: 7 additions & 7 deletions xarray/tests/test_missing.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,8 +11,8 @@
NumpyInterpolator, ScipyInterpolator, SplineInterpolator)
from xarray.core.pycompat import dask_array_type
from xarray.tests import (
assert_array_equal, assert_equal, raises_regex, requires_bottleneck,
requires_dask, requires_scipy)
assert_array_equal, assert_equal, assert_allclose, raises_regex,
requires_bottleneck, requires_dask, requires_scipy)


@pytest.fixture
Expand Down Expand Up @@ -89,7 +89,7 @@ def test_interpolate_pd_compat():
# only checks that interpolated values are the same (not nans)
expected.values[pd.isnull(actual.values)] = np.nan

np.testing.assert_allclose(actual.values, expected.values)
np.testing.assert_almost_equal(actual.values, expected.values)


@requires_scipy
Expand Down Expand Up @@ -244,7 +244,7 @@ def test_interpolate_limits():
expected = xr.DataArray(np.array([1, 2, 3, 4, np.nan, 6],
dtype=np.float64), dims='x')

assert_equal(actual, expected)
assert_allclose(actual, expected)


@requires_scipy
Expand Down Expand Up @@ -284,17 +284,17 @@ def test_interpolate_use_coordinate():
# use_coordinate == False is same as using the default index
actual = da.interpolate_na(dim='x', use_coordinate=False)
expected = da.interpolate_na(dim='x')
assert_equal(actual, expected)
assert_allclose(actual, expected)

# possible to specify non index coordinate
actual = da.interpolate_na(dim='x', use_coordinate='xc')
expected = da.interpolate_na(dim='x')
assert_equal(actual, expected)
assert_allclose(actual, expected)

# possible to specify index coordinate by name
actual = da.interpolate_na(dim='x', use_coordinate='x')
expected = da.interpolate_na(dim='x')
assert_equal(actual, expected)
assert_allclose(actual, expected)


@requires_dask
Expand Down

0 comments on commit 1f1ec52

Please sign in to comment.