The xmsinterp library contains classes for performing 2d (x,y) linear, natural neighbor, and idw interpolation: InterpLinear and InterpIdw. Natural neighbor interpolation is performed using the linear interpolation class. 3d (x,y,z) interpolation can be performed using the idw class. Points (x,y) must be given to create the interpolation class.
-We will learn how to use the interpolation classes by using a set of 2D points (x,y) with a measurement at each point. These points come from a csv file. Then we will generate a grid that bounds these points and use different interpolation options to visualize the data.
- -import param
-import panel as pn
-import holoviews as hv
-import numpy as np
-import pandas as pd
-import xarray as xr
-import hvplot.pandas
-import hvplot.xarray
-from holoviews.streams import Params
-
-from xms import interp
-# read in the data to do the interpolation
-# drop the id columns from the data frame so we can pass the remaining data as points to the interp class
-df = pd.read_csv("plumedat.csv").drop(['id'], axis=1)
-df.hvplot.scatter(x='x', y='y', c='c', height=600, cmap='jet', size=120, colorbar=True)
-Get the min, max (x, y) of the data and generate a grid of points.
- -# get the (x,y) range of the data and expand by 10%
-def buffer_range(range_min, range_max, buf_by=0.1):
- buf = (range_max - range_min) * buf_by
- return (range_min - buf, range_max + buf)
-
-x_range = buffer_range(df.x.min(), df.x.max())
-y_range = buffer_range(df.y.min(), df.y.max())
-
-# adjust these values for a higher/lower resolution grid
-number_of_x_points=600
-number_of_y_points=400
-x = np.linspace(x_range[0], x_range[1], number_of_x_points)
-y = np.linspace(y_range[0], y_range[1], number_of_y_points)
-grid = np.meshgrid(x, y)
-pts = np.stack([i.flatten() for i in grid], axis=1)
-Interpolate to the grid using linear interpolation and display the results with matplotlib.
- -linear = interp.interpolate.InterpLinear(df.values)
-z = linear.interpolate_to_points(pts)
-# I would prefer to be able to set up an interpolator like this:
-# linear = interp.Interpolator(kind='linear', values=df.values)
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Interpolate using natural neighbor. The first example shows a "constant" nodal function and the second example shows a "quadratic" nodals function. A nodal function can better model trends that exist in the measured values.
- -linear.set_use_natural_neighbor(on=True)
-z = linear.interpolate_to_points(pts)
-# z = linear.interpolate_to_points(pts, use='natural_neighbor')
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Now do natural neighbor with a quadratic nodal function.
- -linear.set_use_natural_neighbor(on=True, nodal_function_type="quadratic",
- nodal_function_point_search_option="nearest_points")
-z = linear.interpolate_to_points(pts)
-# z = linear.interp_to_pts(pts, use='natural_neighbor', nodal_func_type='quadratic', nd_func_pt_search_opt='nearest_pts')
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Now we will use Inverse Distance Weighted (IDW) to interpolate to the grid. Notice that unlike Linear and Natural Neighbor interpolation, IDW interpolation can extrapolate beyond the bounds of the input data points.
- -idw = interp.interpolate.InterpIdw(df.values)
-z = idw.interpolate_to_points(pts)
-# idw = interp.Interpolator(kind='idw', values=df.values)
-# z = idw.interpolate_to_points(pts)
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Now interpolate using a nodal function to try to better capture trends in the data set.
- -idw.set_nodal_function(nodal_function_type="gradient_plane")
-z = idw.interpolate_to_points(pts)
-# z = idw.interp_to_pts(pts, nodal_func_type='gradient_plane')
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Notice, in the previous interpolation that we get negative values in some locations. If we assume that there should be no negatives values (like when measuring concentration) we can use the truncation options to adjust the interpolation.
- -idw.set_truncation(maximum=150, minimum=0)
-z = idw.interpolate_to_points(pts)
-# z = idw.interpolate_to_points(pts, smax=150, smin=0)
-zz = np.reshape(z, (len(y), len(x)))
-xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
-xa.hvplot(x='x', y='y', height=500, cmap='viridis')
-Below we define a param class for the interpolation options. Then we use panel to show widgets and a dynamic plot. The user can change the widget values and then push the "Compute new plot" button to update the plot using the currently specified interpolation options.
- -class InterpOptions(param.Parameterized):
- interpolation_option = param.ObjectSelector(default='idw', objects=['linear', 'natural_neighbor', 'idw'],
- precedence=1)
- nodal_function = param.ObjectSelector(default='constant', objects=['constant', 'gradient_plane', 'quadratic'],
- precedence=1.1)
- truncation = param.Boolean(default=True, precedence=2)
- truncation_maximum = param.Number(default=df.loc[df['c'].idxmax()]['c'], precedence=2.1)
- truncation_minimum = param.Number(default=df.loc[df['c'].idxmin()]['c'], precedence=2.2)
-
- @param.depends('interpolation_option', watch=True)
- def _update_nodal_function(self):
- if self.interpolation_option == 'linear':
- self.params('nodal_function').precedence = -1
- else:
- self.params('nodal_function').precedence = 1.1
-
- @param.depends('truncation', watch=True)
- def _update_truncation_opts(self):
- if self.truncation:
- self.params('truncation_maximum').precedence = 2.1
- self.params('truncation_minimum').precedence = 2.2
- else:
- self.params('truncation_maximum').precedence = -1
- self.params('truncation_minimum').precedence = -1
-
- def generate_zz(self):
- # Get Interp Object
- interp_object = None
- if self.interpolation_option != 'idw':
- #interp_object = linear
- interp_object = interp.interpolate.InterpLinear(points=df.values)
- if self.interpolation_option == 'natural_neighbor':
- interp_object.set_use_natural_neighbor(on=True, nodal_function_type=self.nodal_function,
- nodal_function_point_search_option="nearest_points")
- else:
- interp_object = interp.interpolate.InterpIdw(points=df.values)
- interp_object.set_nodal_function(nodal_function_type=self.nodal_function)
-
- if self.truncation:
- interp_object.set_truncation(self.truncation_maximum, self.truncation_minimum)
- z = interp_object.interpolate_to_points(pts)
- zz = np.reshape(z, (len(y), len(x)))
- xa = xr.DataArray(coords={'x': x, 'y': y}, dims=['x', 'y'], data=zz.T)
- return xa
-
- def get_image(self, *args, **kwargs):
- return self.generate_zz().hvplot(x='x', y='y', height=500, cmap='viridis')
-opts = InterpOptions(name="")
-dmap = hv.util.Dynamic(opts.get_image(), operation=opts.get_image, streams=[Params(opts)])
-pn.Row(opts, dmap)
-
-
-
-