[Insight-users] Kd-tree implementation in ITK is buggy and unreliable

[Luis Ibanez]

Hi Kevin,


Thanks a lot for contributing this test.


The test has been added to the test suite of ITK in the directory

          Insight/Testing/Code/Numerics/Statistics

and has been named

             itkKdTreeTest1.cxx

It seems that the testing of the KdTree has been done, so far,
indirectly via other classes.  Which will certainly help to
explain how the cicades have been living in this tree for so
long.

The test provided by Tobias Heimann in his bug report
http://public.kitware.com/Bug/view.php?id=5082

has also been added to the same directory as

             itkKdTreeTest2.cxx


We take note of your point that these two tests may be
revealing two different bugs in the KdTree.

The test, of course, will be failing in the Dashboard,
as we track the problem inside the KdTree.

The second test came with a dataset that is too large
to be included in the data baseline. We can make it
a conditional test or we can try to shrink the dataset
to a minimal case.... we are looking into these two
options.



     Luis



---------------------
Kevin H. Hobbs wrote:
> Here's the actual code this time :
> 
> 
> ------------------------------------------------------------------------
> 
> #include "itkVector.h"
> #include "itkMersenneTwisterRandomVariateGenerator.h"
> #include "itkListSample.h"
> #include "itkKdTree.h"
> #include "itkKdTreeGenerator.h"
> #include "itkEuclideanDistance.h"
> 
> int main(int argc, char * argv[] )
> {
> 	// Random number generator
> 	itk::Statistics::MersenneTwisterRandomVariateGenerator::Pointer mt 
> 		= itk::Statistics::MersenneTwisterRandomVariateGenerator::New();
> 	mt->Initialize();
> 	
> 	typedef itk::Array< double > MeasurementVectorType ;
> 	typedef itk::Statistics::ListSample< MeasurementVectorType > SampleType ;
> 	const SampleType::MeasurementVectorSizeType measurementVectorSize = 2;
> 
> 	SampleType::Pointer sample = SampleType::New() ;
> 	sample->SetMeasurementVectorSize( measurementVectorSize );
> 
> 	MeasurementVectorType mv( measurementVectorSize ) ;
> 	for (unsigned int i = 0 ; i < 1000 ; ++i )
> 	{
> 		mv[0] = mt->GetNormalVariate( 0.0, 1.0 );
> 		mv[1] = mt->GetNormalVariate( 0.0, 1.0 );
> 		sample->PushBack( mv ) ;
> 	}
> 
> 	typedef itk::Statistics::KdTreeGenerator< SampleType > TreeGeneratorType ;
>       	TreeGeneratorType::Pointer treeGenerator = TreeGeneratorType::New() ;
> 
> 	treeGenerator->SetSample( sample ) ;
>       	treeGenerator->SetBucketSize( 16 ) ;
>       	treeGenerator->Update() ;
> 
> 	typedef TreeGeneratorType::KdTreeType TreeType ;
>       	typedef TreeType::NearestNeighbors NeighborsType ;
>       	typedef TreeType::KdTreeNodeType NodeType ;
> 
> 	TreeType::Pointer tree = treeGenerator->GetOutput() ;
> 
> 	MeasurementVectorType queryPoint( measurementVectorSize ) ;
> 
> 	unsigned int numberOfNeighbors = 1 ;
>       	TreeType::InstanceIdentifierVectorType neighbors ;
> 
> 	MeasurementVectorType result( measurementVectorSize ) ;
> 	MeasurementVectorType test_point( measurementVectorSize ) ;
> 	for (unsigned int i = 0 ; i < 1000 ; ++i )
> 	{
> 		double min_dist = INFINITY;
> 		queryPoint[0] = mt->GetNormalVariate( 0.0, 1.0 );
> 		queryPoint[1] = mt->GetNormalVariate( 0.0, 1.0 );
> 		tree->Search( queryPoint, numberOfNeighbors, neighbors ) ;
> 		result = tree->GetMeasurementVector( neighbors[0] );
> 		double result_dist = sqrt(
> 					(result[0] - queryPoint[0]) *
> 					(result[0] - queryPoint[0]) +
> 					(result[1] - queryPoint[1]) *
> 					(result[1] - queryPoint[1])
> 					);
> 		for (unsigned int i = 0 ; i < 1000 ; ++i )
> 		{
> 			double dist;
> 			test_point = tree->GetMeasurementVector( i );
> 			dist = sqrt(
> 					(test_point[0] - queryPoint[0]) *
> 					(test_point[0] - queryPoint[0]) +
> 					(test_point[1] - queryPoint[1]) *
> 					(test_point[1] - queryPoint[1])
> 					);
> 			if (dist < min_dist )
> 				min_dist = dist;
> 
> 		}
> 		if ( min_dist < result_dist )
> 		{
> 			std::cout << min_dist << " " << result_dist << std::endl;
> 			std::cout << "Test FAILED." << std::endl;
> 			return EXIT_FAILURE;
> 
> 		}
> 
> 	}
> 	
> 	std::cout << "Test passed." << std::endl;
> 	return EXIT_SUCCESS;
> }
> 
> 
> ------------------------------------------------------------------------
> 
> cmake_minimum_required(VERSION 2.7)
> 
> project( kdtest )
> 
> find_package ( ITK )
> if ( ITK_FOUND )
> 	include( ${USE_ITK_FILE} )
> endif( ITK_FOUND )
> 
> include_directories(
> 	${kdtest_SOURCE_DIR}
> 	)
> 
> add_executable( kdtest kdtest.cxx )
> target_link_libraries( kdtest ITKCommon)