Skip to content
Merged
Show file tree
Hide file tree
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
8 changes: 0 additions & 8 deletions lib/iris/analysis/_interpolate_backdoor.py
Original file line number Diff line number Diff line change
Expand Up @@ -106,11 +106,3 @@ def linear(cube, sample_points, extrapolation_mode='linear'):
extrapolation_mode=extrapolation_mode)

linear.__doc__ = _interp.linear.__doc__


class Linear1dExtrapolator(six.with_metaclass(ClassWrapperSameDocstring,
_interp.Linear1dExtrapolator)):
@wraps(_interp.Linear1dExtrapolator.__init__)
def __init__(self, interpolator):
_warn_deprecated()
super(Linear1dExtrapolator, self).__init__(interpolator)
124 changes: 0 additions & 124 deletions lib/iris/analysis/_interpolate_private.py
Original file line number Diff line number Diff line change
Expand Up @@ -674,127 +674,3 @@ def linear(cube, sample_points, extrapolation_mode='linear'):

scheme = Linear(extrapolation_mode)
return cube.interpolate(sample_points, scheme)


def _interp1d_rolls_y():
"""
Determines if :class:`scipy.interpolate.interp1d` rolls its array `y` by
comparing the shape of y passed into interp1d to the shape of its internal
representation of y.

SciPy v0.13.x+ no longer rolls the axis of its internal representation
of y so we test for this occurring to prevent us subsequently
extrapolating along the wrong axis.

For further information on this change see, for example:
* https://github.com/scipy/scipy/commit/0d906d0fc54388464603c63119b9e35c9a9c4601
(the commit that introduced the change in behaviour).
* https://github.com/scipy/scipy/issues/2621
(a discussion on the change - note the issue is not resolved
at time of writing).

"""
y = np.arange(12).reshape(3, 4)
f = interp1d(np.arange(3), y, axis=0)
# If the initial shape of y and the shape internal to interp1d are *not*
# the same then scipy.interp1d rolls y.
return y.shape != f.y.shape


class Linear1dExtrapolator(object):
"""
Extension class to :class:`scipy.interpolate.interp1d` to provide linear extrapolation.

See also: :mod:`scipy.interpolate`.

.. deprecated :: 1.10

"""
roll_y = _interp1d_rolls_y()

def __init__(self, interpolator):
"""
Given an already created :class:`scipy.interpolate.interp1d` instance, return a callable object
which supports linear extrapolation.

.. deprecated :: 1.10

"""
self._interpolator = interpolator
self.x = interpolator.x
# Store the y values given to the interpolator.
self.y = interpolator.y
"""
The y values given to the interpolator object.

.. note:: These are stored with the interpolator.axis last.

"""
# Roll interpolator.axis to the end if scipy no longer does it for us.
if not self.roll_y:
self.y = np.rollaxis(self.y, self._interpolator.axis, self.y.ndim)

def all_points_in_range(self, requested_x):
"""Given the x points, do all of the points sit inside the interpolation range."""
test = (requested_x >= self.x[0]) & (requested_x <= self.x[-1])
if isinstance(test, np.ndarray):
test = test.all()
return test

def __call__(self, requested_x):
if not self.all_points_in_range(requested_x):
# cast requested_x to a numpy array if it is not already.
if not isinstance(requested_x, np.ndarray):
requested_x = np.array(requested_x)

# we need to catch the special case of providing a single value...
remember_that_i_was_0d = requested_x.ndim == 0

requested_x = requested_x.flatten()

gt = np.where(requested_x > self.x[-1])[0]
lt = np.where(requested_x < self.x[0])[0]
ok = np.where( (requested_x >= self.x[0]) & (requested_x <= self.x[-1]) )[0]

data_shape = list(self.y.shape)
data_shape[-1] = len(requested_x)
result = np.empty(data_shape, dtype=self._interpolator(self.x[0]).dtype)

# Make a variable to represent the slice into the resultant data. (This will be updated in each of gt, lt & ok)
interpolator_result_index = [slice(None, None)] * self.y.ndim

if len(ok) != 0:
interpolator_result_index[-1] = ok

r = self._interpolator(requested_x[ok])
# Reshape the properly formed array to put the interpolator.axis last i.e. dims 0, 1, 2 -> 0, 2, 1 if axis = 1
axes = list(range(r.ndim))
del axes[self._interpolator.axis]
axes.append(self._interpolator.axis)

result[interpolator_result_index] = r.transpose(axes)

if len(lt) != 0:
interpolator_result_index[-1] = lt

grad = (self.y[..., 1:2] - self.y[..., 0:1]) / (self.x[1] - self.x[0])
result[interpolator_result_index] = self.y[..., 0:1] + (requested_x[lt] - self.x[0]) * grad

if len(gt) != 0:
interpolator_result_index[-1] = gt

grad = (self.y[..., -1:] - self.y[..., -2:-1]) / (self.x[-1] - self.x[-2])
result[interpolator_result_index] = self.y[..., -1:] + (requested_x[gt] - self.x[-1]) * grad

axes = list(range(len(interpolator_result_index)))
axes.insert(self._interpolator.axis, axes.pop(axes[-1]))
result = result.transpose(axes)

if remember_that_i_was_0d:
new_shape = list(result.shape)
del new_shape[self._interpolator.axis]
result = result.reshape(new_shape)

return result
else:
return self._interpolator(requested_x)