Finished perf-tests of max. leaf size vs. build/query time. Added new figures.

perf-tests/README.txt

@@ -21,3 +21,13 @@ want to reproduce the experiments but want to generate the graphs yourself.
 Jose Luis Blanco
 Aug 26, 2011
 
+
+=== ADDED ON May, 2012 ===
+
+More performance tests: 
+
+$ nanoflann/test_leaf_max_size > LEAF_STATS.txt
+
+And visualize the results with the MATLAB/Octave "analyze_leafsize_stats.m". 
+Graphs online have been generated with these tools.
+

perf-tests/analyze_leafsize_stats.m

@@ -0,0 +1,54 @@
+function [] = analyze_leafsize_stats()
+	% Compute the stats from the result files of performance tests wrt 
+	% the max. leaf size
+    close all;
+    
+    %D=load('LEAF_STATS.txt');
+    %D=load('LEAF_STATS_DOUBLE.txt');
+    D=load('LEAF_STATS_DATASET.txt');
+	MAXs = unique(D(:,2));
+    
+    COLs = {'k','b','r','g'};  % Cell array of colros.
+    
+    figure(1);
+    for i=1:length(MAXs),
+        MaxLeaf = MAXs(i);
+		idxs = find(D(:,2)==MaxLeaf);
+        
+        TBs = D(idxs,3); % Time of tree Builds
+        TQs = D(idxs,4); % Time of tree Queries
+        
+        TBm   = 1e3*mean(TBs);
+        TQm   = 1e6*mean(TQs);
+        TBstd = 1e3*std(TBs);
+        TQstd = 1e6*std(TQs);
+                
+        TBms(i)=TBm;
+        TQms(i)=TQm;
+        
+        % Plot a 95% ellipse:
+        my_plot_ellipse([TBm TQm],2*[TBstd TQstd], COLs{ 1+mod(i-1,length(COLs))} );
+        
+        %plot(1e3*TBs,1e6*TQs,'.','color',COLs{ 1+mod(i-1,length(COLs))});
+        
+        text(TBm,TQm+6*TQstd,sprintf('%i',MaxLeaf));
+    end
+    
+    plot(TBms,TQms,'color',[0.4 0.4 0.4],'LineWidth',2.5);
+        
+    set(gca,'YScale','log');
+    
+    ylabel('Tree query time (us)');
+    xlabel('Tree build time (ms)');
+    title(sprintf('nanoflann performance vs. Leaf Max.Size (#points=%.01e)',D(1,1)));
+    
+end
+
+
+function [] = my_plot_ellipse(Ms,STDs, COL )
+    angs=linspace(0,2*pi,200);
+    xs= Ms(1) + cos(angs)*STDs(1);
+    ys= Ms(2) + sin(angs)*STDs(2);
+    h=plot(xs,ys,'color', COL,'LineWidth',2);
+    hold on;
+end

perf-tests/nanoflann/test_leaf_max_size.cpp

@@ -34,6 +34,9 @@
 using namespace std;
 using namespace nanoflann;
 
+// Comment-out for using random points:
+#define REAL_DATASET_FILE "scan_071_points.dat"
+
 //#define VERBOSE_OUTPUT
 
 #ifdef VERBOSE_OUTPUT
@@ -64,18 +67,18 @@ struct PointCloud
 	inline size_t kdtree_get_point_count() const { return pts.size(); }
 
 	// Returns the distance between the vector "p1[0:size-1]" and the data point with index "idx_p2" stored in the class:
-	inline float kdtree_distance(const float *p1, const size_t idx_p2,size_t size) const
+	inline T kdtree_distance(const T *p1, const size_t idx_p2,size_t size) const
 	{
-		float d0=p1[0]-pts[idx_p2].x;
-		float d1=p1[1]-pts[idx_p2].y;
-		float d2=p1[2]-pts[idx_p2].z;
+		T d0=p1[0]-pts[idx_p2].x;
+		T d1=p1[1]-pts[idx_p2].y;
+		T d2=p1[2]-pts[idx_p2].z;
 		return d0*d0+d1*d1+d2*d2;
 	}
 
 	// Returns the dim'th component of the idx'th point in the class:
 	// Since this is inlined and the "dim" argument is typically an immediate value, the
 	//  "if/else's" are actually solved at compile time.
-	inline float kdtree_get_pt(const size_t idx, int dim) const
+	inline T kdtree_get_pt(const size_t idx, int dim) const
 	{
 		if (dim==0) return pts[idx].x;
 		else if (dim==1) return pts[idx].y;
@@ -88,6 +91,32 @@ struct PointCloud
 	template <class BBOX>
 	bool kdtree_get_bbox(BBOX &bb) const { return false; }
 
+	// Default ctor.
+	PointCloud() {} 
+
+	// load dataset in Freiburg 3D scans format:
+	PointCloud(const char* sFil)
+	{
+		FILE *f = fopen(sFil,"rt");
+		if (!f) throw std::runtime_error("can't open dataset file!");
+		pts.clear();
+
+		char str[300];
+		while (fgets(str,sizeof(str),f))
+		{
+			float x,y,z;
+			if (sscanf(str,"%*f %*f %*f %f %f %f\n",&x,&y,&z)==3) 
+			{
+				pts.resize(pts.size()+1);
+				pts.rbegin()->x=x;
+				pts.rbegin()->y=y;
+				pts.rbegin()->z=z;
+			}
+		}
+		fprintf(stderr,"Read dataset: %u points\n", static_cast<unsigned int>(pts.size()));
+
+		fclose(f);
+	} 
 };
 
 double get_time()
@@ -97,6 +126,11 @@ double get_time()
 	return tv.tv_sec+tv.tv_usec/1000000.0;
 }
 
+double my_random(const double min,const double max) 
+{
+	return min+(max-min)*(rand() % 1000) / 1000.0;
+}
+
 template <typename T>
 void generateRandomPointCloud(PointCloud<T> &point, const size_t N, const T max_range = 10)
 {
@@ -115,16 +149,23 @@ void generateRandomPointCloud(PointCloud<T> &point, const size_t N, const T max_
 	VERB_COUT << "done in "<< (t1-t0)*1e3 << " ms\n";
 }
 
+// Load dataset only once:
+#ifdef REAL_DATASET_FILE
+PointCloud<float> cloud(REAL_DATASET_FILE);
+#endif
+
 template <typename num_t>
 void perf_test(const size_t N, const size_t max_leaf_elements)
 {
+#ifndef REAL_DATASET_FILE
+	// Generate random points:
 	PointCloud<num_t> cloud;
-
-	// Generate points:
 	generateRandomPointCloud(cloud, N);
-
 	num_t query_pt[3] = { 0.5, 0.5, 0.5};
-
+#else
+	// Sample dataset is [-40,40]x[-40,40]x[0,15]: Query at random:
+	num_t query_pt[3] = { my_random(-40.0,40.0), my_random(-40.0,40.0), my_random(0,10)};
+#endif
 
 	// construct an randomized kd-tree index using 4 kd-trees
 	double t0=get_time();
@@ -159,7 +200,7 @@ void perf_test(const size_t N, const size_t max_leaf_elements)
 
 	// Output for stats generation:
 	cout
-		<< N << "\t "
+		<< cloud.pts.size() << "\t "
 		<< max_leaf_elements << "\t "
 		<< At_build << "\t "
 		<< At_search << "\n";