Interpolation module reference

Interpolation routines.

References

[1]http://en.wikipedia.org/wiki/Inverse_distance_weighting
[2]Lüthi, Beat. Some Aspects of Strain, Vorticity and Material Element Dynamics as Measured with 3D Particle Tracking Velocimetry in a Turbulent Flow. PhD Thesis, ETH-Zürich (2002).
[3]http://en.wikipedia.org/wiki/Radial_basis_function

Documentation

class flowtracks.interpolation.GeneralInterpolant(method, num_neighbs=None, radius=None, param=None)[source]

Holds all parameters necessary for performing an interpolation. Use is as a callable object after initialization, see __call__().

Arguments

  • method: interpolation method. Either ‘inv’ for inverse-distance weighting, ‘rbf’ for gaussian-kernel Radial Basis Function method, or ‘corrfun’ for using a correlation function.
  • radius: of the search area for neighbours, [m]. If None, select closest neighbs.
  • neighbs: number of closest neighbours to interpolate from. If None. uses 4 neighbours for ‘inv’ method, and 7 for ‘rbf’, unless radius is not None, then neighbs is ignored.
  • param: the parameter adjusting the interpolation method. For IDW it is the inverse power (default 1), for rbf it is epsilon (default 1e5).
__call__(tracer_pos, interp_points, data, companionship=None)[source]

Sets up the necessary parameters, and performs the interpolation. Does not change the scene set by set_scene if any, so may be used for any off-scene interpolation.

Arguments

  • tracer_pos: (n,3) array, the x,y,z coordinates of one tracer per row, in [m]
  • interp_points: (m,3) array, coordinates of points where interpolation will be done.
  • data: (n,d) array, the for the d-dimensional data for tracer n. For example, in velocity interpolation this would be (n,3), each tracer having 3 components of velocity.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for analysing a simulated particle that started from a true tracer.

Returns

  • vel_interp: an (m,3) array with the interpolated value at the position of each particle, [m/s].
current_relative_positions()[source]

Returns an (m,k,3) array, the distance between interpolation point m and each of k nearest neighbours on each axis.

eulerian_jacobian(local_interp=None, eps=0.0001)[source]

A general way to calculate the velocity derivatives. It could be enhanced in the future by specific analytical derivatives of the different interpolation methods. The Jacobian is calculated for the current scene, as recorded with set_scene()

Arguments

  • local_interp: results of interpolation already performed at the position where derivatives are wanted. If not given, an interpolation of recorded scene data is automatically performed.
  • eps: the dx in each direction.

Returns

interpolate(subset=None)[source]

Performs an interpolation over the recorded scene.

Arguments

  • subset: a neighbours selection array, such as returned from which_neighbours(), to replace the recorded selection. Default value (None) uses the recorded selection. The recorded selection is not changed, so subset is forgotten after the call.

Returns

an (m,3) array with the interpolated value at the position of each of m particles.

neighb_dists(tracer_pos, interp_points, companionship=None)[source]

The distance from each interpolation point to each data point of those used for interpolation. Assumes, for now, a constant number of neighbours.

Arguments

  • tracer_pos: (n,3) array, the x,y,z coordinates of one tracer per row, in [m]
  • interp_points: (m,3) array, coordinates of points where interpolation will be done.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for analysing a simulated particle that started from a true tracer.

Returns

  • ndists: an (m,c) array, for c closest neighbours as defined during object construction.
save_config(cfg)[source]

Adds the keys necessary for recreating this interpolant into a configuration object. It is the caller’s responsibility to do a writeback to file.

Arguments

  • cfg: a ConfigParser object.
set_data_on_current_field(data)[source]

Change the data on the existing interpolation points. This enables redoing an interpolation on a scene without recalculating weights, when weights are only dependent on position, as they are for most interpolation methods.

Arguments

  • data: (n,d) array, the for the d-dimensional data for n tracers. For example, in velocity interpolation this would be (n,3), each tracer having 3 components of velocity.
set_field_positions(positions)[source]

Sets the positions of points where there is data for interpolation. This sets up the structures for efficiently finding distances etc.

Arguments

  • positions: (n,3) array, for position of n points in 3D space.
set_interp_points(points, companionship=None)[source]

Sets the points into which interpolation will be done. It is possible to set this once and then do interpolation of several datasets into these points.

Arguments

  • positions: (m,3) array, for position of m target points in 3D space.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for analysing a simulated particle that started from a true tracer.
set_scene(tracer_pos, interp_points, data=None, companionship=None)[source]

Records scene data for future interpolation using the same scene.

Arguments

  • tracer_pos: (n,3) array, the x,y,z coordinates of one tracer per row, in [m]
  • interp_points: (m,3) array, coordinates of points where interpolation will be done.
  • data: (n,d) array, the for the d-dimensional data for tracer n. For example, in velocity interpolation this would be (n,3), each tracer having 3 components of velocity.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for analysing a simulated particle that started from a true tracer.
trim_points(which)[source]

Remove interpolation points from the scene.

Arguments

  • which: a boolean array, length is number of current particle list (as given in set_scene), True to trim a point, False to keep.
flowtracks.interpolation.Interpolant(method, num_neighbs=None, radius=None, param=None)

Factory function. Returns an object of the interpolant class that matches the given method. All classes are subclassed from GeneralInterpolant.

Arguments

  • method: interpolation method. Either ‘inv’ for inverse-distance weighting, ‘rbf’ for gaussian-kernel Radial Basis Function method, or ‘corrfun’ for using a correlation function.
  • radius: of the search area for neighbours, [m]. If None, select closest neighbs.
  • neighbs: number of closest neighbours to interpolate from. If None. uses 4 neighbours for ‘inv’ method, and 7 for ‘rbf’, unless radius is not None, then neighbs is ignored.
  • param: the parameter adjusting the interpolation method. For IDW it is the inverse power (default 1), for rbf it is epsilon (default 1e5).
class flowtracks.interpolation.InverseDistanceWeighter(num_neighbs=None, radius=None, param=None)[source]

Holds all parameters necessary for performing an inverse-distance interpolation [1]. Use is either as a callable object after initialization, see __call__(), or by setting a scene for repeated interpolation, see set_scene() and interpolate()

Arguments

  • num_neighbs: number of closest neighbours to interpolate from. If None uses 4 neighbours, unless radius is not None, then neighbs is ignored.
  • radius: of the search area for neighbours, [m]. If None, select closest neighbs.
  • param: the inverse power of distance to use (default 1).
__call__(tracer_pos, interp_points, data, companionship=None)[source]

Sets up the necessary parameters, and performs the interpolation. Does not change the scene set by set_scene if any, so may be used for any off-scene interpolation.

Arguments

  • tracer_pos: (n,3) array, the x,y,z coordinates of one tracer per row, in [m]
  • interp_points: (m,3) array, coordinates of points where interpolation will be done.
  • data: (n,d) array, the for the d-dimensional data for tracer n. For example, in velocity interpolation this would be (n,3), each tracer having 3 components of velocity.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for analysing a simulated particle that started from a true tracer.

Returns

  • vel_interp: an (m,3) array with the interpolated value at the position of each particle, [m/s].
eulerian_jacobian(local_interp=None, eps=None)[source]

Velocity derivatives. The Jacobian is calculated for the current scene, as recorded with set_scene()

Arguments

  • local_interp: results of interpolation already performed at the position where derivatives are wanted. If not given, an interpolation of recorded scene data is automatically performed.
  • eps: unused, here for compatibility with base class.

Returns

(m,d,3) array, for m interpolation points and d interpolation dimentions. For each point, [i,j] = du_i/dx_j

set_scene(tracer_pos, interp_points, data=None, companionship=None)[source]

Adds to the base class only a precalculation of weights.

trim_points(which)[source]

Remove interpolation points from the scene.

Arguments

  • which: a boolean array, length is number of current particle list (as given in set_scene), True to trim a point, False to keep.
weights(dists, use_parts)[source]

Calculate the respective weight of each tracer j=1..n in the interpolation point i=1..m. The actual weight is normalized to the sum of weights in the interpolation, not here.

Arguments

  • dists: (m,n) array, the distance of interpolation_point i=1...m from tracer j=1...n, for (row,col) (i,j) [m] where n is the number of nearest neighbours.
  • use_parts: (m,n) boolean array, whether tracer j is a neighbour of particle i, same indexing as dists.

Returns

  • weights: an (m,n) array.
flowtracks.interpolation.corrfun_interp(dists, use_parts, data, corrs_hist, corrs_bins)[source]

For each of n particle, generate the velocity interpolated to its position from all neighbours as selected by caller. The weighting of neighbours is by the correlation function, e.g. if the distance at neighbor i is \(r_i\), then it adds \(\rho(r_i)*v_i\) to the interpolated velocity. This is done for each component separately.

Arguemnts

  • dists: (m,n) array, the distance of interpolation_point \(i=1...m\) from tracer \(j=1...n\), for (row,col) (i,j) [m]
  • use_parts: (m,n) boolean array, whether tracer j is a neighbour of particle i, same indexing as dists.
  • data: (n,d) array, the d components of the data that is interpolated from, for each of n tracers.
  • corrs_hist: the correlation function histogram, an array of b bins.
  • corrs_bins: same size array, the bin start point for each bin.

Returns

  • vel_avg: an (m,3) array with the interpolated velocity at each interpolation point, [units of data].
flowtracks.interpolation.interpolant(method, num_neighbs=None, radius=None, param=None)[source]

Factory function. Returns an object of the interpolant class that matches the given method. All classes are subclassed from GeneralInterpolant.

Arguments

  • method: interpolation method. Either ‘inv’ for inverse-distance weighting, ‘rbf’ for gaussian-kernel Radial Basis Function method, or ‘corrfun’ for using a correlation function.
  • radius: of the search area for neighbours, [m]. If None, select closest neighbs.
  • neighbs: number of closest neighbours to interpolate from. If None. uses 4 neighbours for ‘inv’ method, and 7 for ‘rbf’, unless radius is not None, then neighbs is ignored.
  • param: the parameter adjusting the interpolation method. For IDW it is the inverse power (default 1), for rbf it is epsilon (default 1e5).
flowtracks.interpolation.rbf_interp(tracer_dists, dists, use_parts, data, epsilon=0.01)[source]

Radial-basis interpolation [3] for each particle, from all neighbours selected by caller. The difference from inv_dist_interp is that the weights are independent of interpolation point, among other differences.

Arguments

  • tracer_dists: (n,n) array, the distance of tracer \(i=1...n\) from tracer \(j=1...n\), for (row,col) (i,j) [m]
  • dists: (m,n) array, the distance from interpolation point \(i=1...m\) to tracer j. [m]
  • use_parts: (m,n) boolean array, True where tracer \(j=1...n\) is a neighbour of interpolation point \(i=1...m\).
  • data: (n,d) array, the d components of the data for each of n tracers.

Returns

  • vel_interp: an (m,3) array with the interpolated velocity at the position of each particle, [m/s].
flowtracks.interpolation.read_interpolant(conf_fname)[source]

Builds an Interpolant object based on values in an INI-formatted file.

Arguments

  • conf_fname: path to configuration file.

Returns

an Interpolant object constructed from values in the configuration file.

flowtracks.interpolation.select_neighbs(tracer_pos, interp_points, radius=None, num_neighbs=None, companionship=None)[source]

For each of m interpolation points, find its distance to all tracers. Use result to decide which tracers are the neighbours of each interpolation point, based on either a fixed radius or the closest num_neighbs.

Arguments

  • tracer_pos: (n,3) array, the x,y,z coordinates of one tracer per row, [m]
  • interp_points: (m,3) array, coordinates of points where interpolation will be done.
  • radius: of the search area for neighbours, [m]. If None, select closest num_neighbs.
  • num_neighbs: number of closest neighbours to interpolate from. If None. uses all neighbours in a given radius. radius has precedence.
  • companionship: an optional array denoting for each interpolation point the index of a tracer that should be excluded from it (“companion tracer”), useful esp. for interpolating tracers unto themselves and for analysing a simulated particle that started from a true tracer.

Returns

  • dists: (m,n) array, the distance from each interpolation point to each tracer.
  • use_parts: (m,n) boolean array, True where tracer \(j=1...n\) is a neighbour of interpolation point \(i=1...m\).