Commit 4b9b94d5 authored by Mathieu Westphal's avatar Mathieu Westphal

Adding a vtkStaticCellLocator::FindClosestPointWithinRadius

parent ee0c9688
......@@ -62,6 +62,8 @@ int TestStaticCellLocator(int,char *[])
int subId;
vtkIdType static_cellId, ref_cellId;
double static_dist2, ref_dist2, static_closest[3], ref_closest[3];
int static_inside, ref_inside;
vtkIdType static_ret, ref_ret;
int num_failed = 0;
for (int i=0; i<10; ++i)
......@@ -84,6 +86,34 @@ int TestStaticCellLocator(int,char *[])
num_failed++;
}
static_ret = static_loc->FindClosestPointWithinRadius(
p, 5.0, static_closest, cell, static_cellId, subId, static_dist2, static_inside);
ref_ret = ref_loc->FindClosestPointWithinRadius(
p, 5.0, ref_closest, cell, ref_cellId, subId, ref_dist2, ref_inside);
if (static_ret != ref_ret)
{
std::cerr << "different closest point within radius result\n";
num_failed++;
}
else if (static_ret)
{
// note that the cell id, inside and even closest point are not always identical
// but, the distance should be nearly identical.
ok = (std::abs(static_dist2 - ref_dist2) < 1e-12);
if (!ok)
{
std::cerr << "different closest point within radius:\n";
std::cerr << "\t" << static_cellId << " - " << ref_cellId << "\n";
std::cerr << "\t" << static_dist2 << " - " << ref_dist2 << "\n";
std::cerr << "\t(" << static_closest[0] << ", " << static_closest[1] << ", "
<< static_closest[2] << ") - (" << ref_closest[0] << ", " << ref_closest[1]
<< ", " << ref_closest[2] << ")\n";
num_failed++;
}
}
}
return (num_failed==0) ? 0 : 1;
......
......@@ -285,9 +285,9 @@ struct vtkCellProcessor
double& t, double x[3], double pcoords[3],
int &subId, vtkIdType &cellId,
vtkGenericCell *cell) = 0;
virtual void FindClosestPoint(const double x[3], double closestPoint[3],
vtkGenericCell *cell, vtkIdType &cellId,
int &subId, double& dist2) = 0;
virtual vtkIdType FindClosestPointWithinRadius(const double x[3], double radius,
double closestPoint[3], vtkGenericCell* cell, vtkIdType& cellId, int& subId, double& dist2,
int& inside) = 0;
// Convenience for computing
virtual int IsEmpty(vtkIdType binId) = 0;
......@@ -370,9 +370,8 @@ struct CellProcessor : public vtkCellProcessor
double& t, double x[3], double pcoords[3],
int &subId, vtkIdType &cellId,
vtkGenericCell *cell) override;
void FindClosestPoint(const double x[3], double closestPoint[3],
vtkGenericCell *cell, vtkIdType &cellId,
int &subId, double& dist2) override;
vtkIdType FindClosestPointWithinRadius(const double x[3], double radius, double closestPoint[3],
vtkGenericCell* cell, vtkIdType& cellId, int& subId, double& dist2, int& inside) override;
int IsEmpty(vtkIdType binId) override
{
return ( this->GetNumberOfIds(static_cast<T>(binId)) > 0 ? 0 : 1 );
......@@ -976,10 +975,10 @@ double Distance2ToBounds(const double x[3], double bounds[6])
//-----------------------------------------------------------------------------
// Return closest point (if any) AND the cell on which this closest point lies
template <typename T> void CellProcessor<T>::
FindClosestPoint(const double x[3], double closestPoint[3],
vtkGenericCell *cell, vtkIdType &closestCellId,
int &closestSubId, double& minDist2)
template<typename T>
vtkIdType CellProcessor<T>::FindClosestPointWithinRadius(const double x[3], double radius,
double closestPoint[3], vtkGenericCell* cell, vtkIdType& closestCellId, int& closestSubId,
double& minDist2, int& inside)
{
std::vector<bool> binHasBeenQueued(this->NumBins, false);
std::vector<bool> cellHasBeenVisited(this->NumCells, false);
......@@ -988,6 +987,7 @@ FindClosestPoint(const double x[3], double closestPoint[3],
double distance2ToCellBounds, dist2;
int subId;
int ijk[3];
vtkIdType retVal = 0;
using node = std::pair<double, vtkIdType>;
std::priority_queue< node, std::vector<node>, std::greater<node> > queue;
......@@ -998,7 +998,7 @@ FindClosestPoint(const double x[3], double closestPoint[3],
binHasBeenQueued[binId] = true;
// distance to closest point
minDist2 = VTK_DOUBLE_MAX;
minDist2 = radius * radius;
while (!queue.empty())
{
......@@ -1056,6 +1056,8 @@ FindClosestPoint(const double x[3], double closestPoint[3],
if (stat != -1 && dist2 < minDist2)
{
retVal = 1;
inside = stat;
minDist2 = dist2;
closestCellId = cellId;
closestSubId = subId;
......@@ -1107,8 +1109,7 @@ FindClosestPoint(const double x[3], double closestPoint[3],
}
}
}
// any post-processing? don't think so.
return retVal;
}
//-----------------------------------------------------------------------------
......@@ -1356,13 +1357,26 @@ void vtkStaticCellLocator::
FindClosestPoint(const double x[3], double closestPoint[3],
vtkGenericCell *cell, vtkIdType &cellId,
int &subId, double& dist2)
{
int inside;
double radius = vtkMath::Inf();
double point[3] = { x[0], x[1], x[2] };
this->FindClosestPointWithinRadius(
point, radius, closestPoint, cell, cellId, subId, dist2, inside);
}
//----------------------------------------------------------------------------
vtkIdType vtkStaticCellLocator::FindClosestPointWithinRadius(double x[3], double radius,
double closestPoint[3], vtkGenericCell* cell, vtkIdType& cellId, int& subId, double& dist2,
int& inside)
{
this->BuildLocator();
if ( ! this->Processor )
{
return;
return 0;
}
this->Processor->FindClosestPoint(x,closestPoint,cell,cellId,subId,dist2);
return this->Processor->FindClosestPointWithinRadius(
x, radius, closestPoint, cell, cellId, subId, dist2, inside);
}
//-----------------------------------------------------------------------------
......
......@@ -75,6 +75,7 @@ public:
//@}
using vtkAbstractCellLocator::FindClosestPoint;
using vtkAbstractCellLocator::FindClosestPointWithinRadius;
/**
* Test a point to find if it is inside a cell. Returns the cellId if inside
......@@ -134,6 +135,25 @@ public:
vtkGenericCell *cell, vtkIdType &cellId,
int &subId, double& dist2) override;
/**
* Return the closest point within a specified radius and the cell which is
* closest to the point x. The closest point is somewhere on a cell, it
* need not be one of the vertices of the cell. This method returns 1 if a
* point is found within the specified radius. If there are no cells within
* the specified radius, the method returns 0 and the values of
* closestPoint, cellId, subId, and dist2 are undefined. This version takes
* in a vtkGenericCell to avoid allocating and deallocating the cell. This
* is much faster than the version which does not take a *cell, especially
* when this function is called many times in a row such as by a for loop,
* where the allocation and dealloction can be done only once outside the
* for loop. If a closest point is found, "cell" contains the points and
* ptIds for the cell "cellId" upon exit. If a closest point is found,
* inside returns the return value of the EvaluatePosition call to the
* closest cell; inside(=1) or outside(=0).
*/
vtkIdType FindClosestPointWithinRadius(double x[3], double radius, double closestPoint[3],
vtkGenericCell* cell, vtkIdType& cellId, int& subId, double& dist2, int& inside) override;
/**
* Return intersection point (if any) AND the cell which was intersected by
* the finite line. The cell is returned as a cell id and as a generic cell.
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment