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.
Interpolant
(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 closestneighbs
.neighbs
: number of closest neighbours to interpolate from. If None. uses 4 neighbours for ‘inv’ method, and 7 for ‘rbf’, unlessradius
is not None, thenneighbs
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)[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.
Returns
vel_interp
: an (m,3) array with the interpolated value at the position of each particle, [m/s].
-
interpolate
(subset=None)[source]¶ Performs an interpolation over the recorded scene.
Arguments
subset
: a neighbours selection array, such as returned fromwhich_neighbours()
, to replace the recorded selection. Default value (None) uses the recorded selection. The recorded selection is not changed, sosubset
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)[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.
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_scene
(tracer_pos, interp_points, data)[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.
-
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.
-
which_neighbours
()[source]¶ Finds the neighbours that would be selected for use at each interpolation point, given the current scene as set by set_scene().
Returns
(m,n) boolean array, True where tracer \(j=1...n\) is a neighbour of interpolation point \(i=1...m\) under the reigning selection criteria.
-
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 asdists
.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 ofdata
].
-
flowtracks.interpolation.
inv_dist_interp
(dists, use_parts, velocity, p=1)[source]¶ For each of n particle, generate the velocity interpolated to its position from all neighbours as selected by caller. Interpolation method is inverse-distance weighting, [1]
Arguments
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 asdists
.velocity
: (n,3) array, the u,v,w velocity components for each of n tracers, [m/s]p
: the power of inverse distance weight, w = r^(-p). default 1. Use 0 for simple averaging.
Returns
vel_avg
: an (m,3) array with the interpolated velocity at each interpolation point, [m/s].
-
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)[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.
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\).