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
9 changes: 7 additions & 2 deletions lib/iris/analysis/trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -237,12 +237,17 @@ def interpolate(cube, sample_points, method=None):
# Use a cache with _nearest_neighbour_indices_ndcoords()
cache = {}

# Cache the linear interpolator
scheme = iris.analysis.Linear()
coords, points = zip(*sample_points)
interpolator = scheme.interpolator(cube, coords)
Copy link
Member

@rhattersley rhattersley May 3, 2016

Choose a reason for hiding this comment

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

These lines need to be protected by if method in ('linear', None) - partly because it's wasted effort otherwise, but mostly because they will fail when trying to do nearest-neighbour interpolation over multi-dimensional coordinates.

NB. The docstring for iris.analysis.trajectory.interpolate is a little ambiguous about the default value of method, but it turns out that the default is not just simply "linear", but switches to "nearest" when using multi-dimensional coordinates. 👎 I definitely don't propose to change that behaviour in this PR, but it would be nice to clarify the docstring (again, possibly in a separate PR - small is beautiful).


Copy link
Member

Choose a reason for hiding this comment

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

The corresponding import statement needs to be removed.

for i in range(trajectory_size):
point = [(coord, values[i]) for coord, values in sample_points]

if method in ["linear", None]:
column = linear_regrid(cube, point)
new_cube.data[..., i] = column.data
column = interpolator([val[i] for _, val in sample_points])
new_cube.data[..., i] = column
Copy link
Member

Choose a reason for hiding this comment

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

This still needs the .data suffix, i.e.:

new_cube.data[..., i] = column.data

elif method == "nearest":
column_index = _nearest_neighbour_indices_ndcoords(cube, point, cache=cache)
column = cube[column_index]
Expand Down