[Insight-users] itkKdTree returnes wrong closest point

[Ferdinando Rodriguez y Baena]

Dear all,

I am running a few tests on the itkKdTree and seem to be getting the wrong value returned by the closest point search. If I create a 3D KD tree of floats and search for a closest point which happens to be exactly a point in the list, the value returned is wrong. I run a corresponding linear search through the same pointset (which is the one I copy into a sample and set into the KD tree )  and the right value is returned i.e. the points have been loaded correctly, but the closest point returned by the KD tree accelerated search is different from the one returned by the linear search.
I have attached  snippet of the code below.

// KD tree typedefs
typedef itk::Statistics::ListSample< PointType >             SampleType;
typedef itk::Statistics::KdTreeGenerator< SampleType > TreeGeneratorType;
typedef TreeGeneratorType::KdTreeType                        TreeType;

// Points typedefs
typedef itk::PointSet< float, 3 >       PointSetType;
typedef PointSetType::PointType     PointType;

/// Create vars
// Create a sample variable for the tree
SampleType::Pointer Sample = SampleType::New();

// Create tree generator
TreeGeneratorType::Pointer TreeGenerator = TreeGeneratorType::New();

// create tree pointer   
TreeType::Pointer PointsTree;

// load points from a mesh into a container
PointsContainer::Pointer points = m_Mesh->GetPoints();

PointsIterator point_iterator = points->Begin();

// copy points into sample
while ( point_iterator != points->End() ) 
   {
   Sample->PushBack( point_iterator->Value() );
   ++point_iterator;                       // advance to next point
   }

TreeGenerator->SetSample( Sample );
TreeGenerator->SetBucketSize( 3 );
TreeGenerator->Update();

PointsTree = m_TreeGenerator->GetOutput();

unsigned int number_of_neighbors = 1;
TreeType::InstanceIdentifierVectorType neighbors;
PointsTree->Search( query_point, number_of_neighbors, neighbors ) ; 

min_distance = query_point.SquaredEuclideanDistanceTo( PointsTree->GetMeasurementVector( neighbors[0] ) );        // WRONG VALUE

PointsContainer::Pointer points = Mesh->GetPoints();

PointsIterator point_iterator = points->Begin();
float min_err = 10000;
PointType best_pt;

while ( point_iterator != points->End() ) 
{
   if( query_point.SquaredEuclideanDistanceTo(point_iterator->Value()) < min_err)
   {
      min_err = query_point.SquaredEuclideanDistanceTo(point_iterator->Value());                    / /RIGHT VALUE
      best_pt = point_iterator->Value();
   }

   ++point_iterator;                           // advance to next point

}

if ( best_pt.SquaredEuclideanDistanceTo( PointsTree->GetMeasurementVector( neighbors[0] ) ) > 0.0001)
   // I HAVE A PROBLEM

Also, if I set the bucket size of the three to less than 3 I get an uncaught exception.

Any ideas?

Thank you all.

Ferdinando
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://public.kitware.com/pipermail/insight-users/attachments/20041025/5fce8c3c/attachment.html