37 #ifndef FIELDUTILS_INTERPOLATOR_H
38 #define FIELDUTILS_INTERPOLATOR_H
43 #include <boost/geometry/geometries/box.hpp>
44 #include <boost/geometry/geometries/point.hpp>
45 #include <boost/geometry/index/rtree.hpp>
47 #include <boost/function.hpp>
48 #include <boost/shared_ptr.hpp>
97 short int coordId = -1,
115 const std::vector<MultiRegions::ExpListSharedPtr> expInField,
116 std::vector<MultiRegions::ExpListSharedPtr> &expOutField,
121 const std::vector<MultiRegions::ExpListSharedPtr> expInField,
128 std::vector<MultiRegions::ExpListSharedPtr> &expOutField);
155 template <
typename FuncPo
interT,
typename ObjectPo
interT>
172 : idx(idx), coords(coords), dist(dist){};
176 return (dist < comp.
dist);
182 typedef boost::geometry::model::point<NekDouble, m_dim, boost::geometry::cs::cartesian>
BPoint;
184 typedef boost::geometry::index::rtree<PtsPointPair, boost::geometry::index::rstar<16> >
PtsRtree;
214 boost::function<void(const int position, const int goal)>
219 const int maxPts = 250);
233 std::vector<PtsPoint> &neighbourPts,
235 const unsigned int maxPts = 1);
238 std::vector<PtsPoint> &neighbourPts,
239 const unsigned int numPts = 1);
Array< TwoD, unsigned int > m_neighInds
Indices of the relevant neighbours for each physical point. Structure: m_neighInds[ptIdx][neighbourId...
std::pair< BPoint, unsigned int > PtsPointPair
FIELD_UTILS_EXPORT void CalcW_Gauss(const PtsPoint &searchPt, const NekDouble sigma, const int maxPts=250)
Computes interpolation weights using gaussian interpolation.
A class that contains algorithms for interpolation between pts fields, expansions and different meshe...
static const int m_dim
dimension of this interpolator. Hardcoded to 3
void SetProgressCallback(FuncPointerT func, ObjectPointerT obj)
sets a callback funtion which gets called every time the interpolation progresses ...
FIELD_UTILS_EXPORT void FindNeighbours(const PtsPoint &searchPt, std::vector< PtsPoint > &neighbourPts, const NekDouble dist, const unsigned int maxPts=1)
Finds the neares neighbours of a point.
short int m_coordId
coordinate id along which the interpolation should be performed
PtsPoint(int idx, Array< OneD, NekDouble > coords, NekDouble dist)
boost::geometry::index::rtree< PtsPointPair, boost::geometry::index::rstar< 16 > > PtsRtree
FIELD_UTILS_EXPORT InterpMethod GetInterpMethod() const
Returns the interpolation method used by this interpolator.
FIELD_UTILS_EXPORT void Interpolate(const LibUtilities::PtsFieldSharedPtr ptsInField, LibUtilities::PtsFieldSharedPtr &ptsOutField)
Interpolate from a pts field to a pts field.
FIELD_UTILS_EXPORT void CalcW_Linear(const PtsPoint &searchPt, int coordId)
Computes interpolation weights using linear interpolation.
bool operator<(const PtsPoint &comp) const
InterpMethod m_method
Interpolation Method.
static InterpolatorSharedPtr NullInterpolator
FIELD_UTILS_EXPORT NekDouble GetFiltWidth() const
Returns the filter width.
Array< OneD, NekDouble > coords
FIELD_UTILS_EXPORT LibUtilities::PtsFieldSharedPtr GetOutField() const
Returns the output field.
LibUtilities::PtsFieldSharedPtr m_ptsInField
input field
boost::shared_ptr< PtsRtree > m_rtree
A tree structure to speed up the neighbour search. Note that we fill it with an iterator, so instead of rstar, the packing algorithm is used.
boost::geometry::model::point< NekDouble, m_dim, boost::geometry::cs::cartesian > BPoint
boost::shared_ptr< PtsField > PtsFieldSharedPtr
FIELD_UTILS_EXPORT void CalcW_NNeighbour(const PtsPoint &searchPt)
Computes interpolation weights using nearest neighbour interpolation.
FIELD_UTILS_EXPORT void CalcWeights(const LibUtilities::PtsFieldSharedPtr ptsInField, LibUtilities::PtsFieldSharedPtr &ptsOutField)
Compute interpolation weights without doing any interpolation.
boost::function< void(const int position, const int goal)> m_progressCallback
FIELD_UTILS_EXPORT void PrintStatistics()
Print statics of the interpolation weights.
FIELD_UTILS_EXPORT void FindNNeighbours(const PtsPoint &searchPt, std::vector< PtsPoint > &neighbourPts, const unsigned int numPts=1)
Finds the neares neighbours of a point.
FIELD_UTILS_EXPORT void SetupTree()
int m_maxPts
Max number of interpolation points.
Interpolator(InterpMethod method=eNoMethod, short int coordId=-1, NekDouble filtWidth=0.0, int maxPts=1000)
Constructor of the Interpolator class.
LibUtilities::PtsFieldSharedPtr m_ptsOutField
output field
FIELD_UTILS_EXPORT int GetDim() const
returns the dimension of the Interpolator. Should be higher than the dimensions of the interpolated f...
Array< TwoD, float > m_weights
Interpolation weights for each neighbour. Structure: m_weights[physPtIdx][neighbourIdx].
FIELD_UTILS_EXPORT LibUtilities::PtsFieldSharedPtr GetInField() const
Returns the input field.
std::vector< MultiRegions::ExpListSharedPtr > m_expOutField
output field
FIELD_UTILS_EXPORT void CalcW_Quadratic(const PtsPoint &searchPt, int coordId)
Computes interpolation weights using quadratic interpolation.
#define FIELD_UTILS_EXPORT
boost::shared_ptr< Interpolator > InterpolatorSharedPtr
FIELD_UTILS_EXPORT int GetCoordId() const
Returns the coordinate id along which the interpolation should be performed.
NekDouble m_filtWidth
Filter width used for some interpolation algorithms.
std::vector< MultiRegions::ExpListSharedPtr > m_expInField
input field
FIELD_UTILS_EXPORT void CalcW_Shepard(const PtsPoint &searchPt, int numPts)
Computes interpolation weights using linear interpolation.