ann.h

c++实现的KNN库:建立高维度的K-d tree,实现K邻域搜索

//		calculation) we assume that POW and # together are monotonically
//		increasing on sequences of arguments, meaning that for all
//		v0..vk and y:
//
//		POW(v0) #...# POW(vk) <= (POW(v0) #...# POW(vk)) # POW(y).
//
//	Incremental Distance Calculation:
//		The program uses an optimized method of computing distances for
//		kd-trees and bd-trees, called incremental distance calculation.
//		It is used when distances are to be updated when only a single
//		coordinate of a point has been changed.  In order to use this,
//		we assume that there is an incremental update function DIFF(x,y)
//		for #, such that if:
//
//					s = x0 # ... # xi # ... # xk 
//
//		then if s' is equal to s but with xi replaced by y, that is, 
//		
//					s' = x0 # ... # y # ... # xk
//
//		then the length of s' can be computed by:
//
//					|s'| = |s| # DIFF(xi,y).
//
//		Thus, if # is + then DIFF(xi,y) is (yi-x).  For the L_infinity
//		norm we make use of the fact that in the program this function
//		is only invoked when y > xi, and hence DIFF(xi,y)=y.
//
//		Finally, for approximate nearest neighbor queries we assume
//		that POW and ROOT are related such that
//
//					v*ROOT(x) = ROOT(POW(v)*x)
//
//		Here are the values for the various Minkowski norms:
//
//		L_p:	p even:							p odd:
//				-------------------------		------------------------
//				POW(v)			= v^p			POW(v)			= |v|^p
//				ROOT(x)			= x^(1/p)		ROOT(x)			= x^(1/p)
//				#				= +				#				= +
//				DIFF(x,y)		= y - x			DIFF(x,y)		= y - x 
//
//		L_inf:
//				POW(v)			= |v|
//				ROOT(x)			= x
//				#				= max
//				DIFF(x,y)		= y
//
//		By default the Euclidean norm is assumed.  To change the norm,
//		uncomment the appropriate set of macros below.
//----------------------------------------------------------------------

//----------------------------------------------------------------------
//	Use the following for the Euclidean norm
//----------------------------------------------------------------------
#define ANN_POW(v)			((v)*(v))
#define ANN_ROOT(x)			sqrt(x)
#define ANN_SUM(x,y)		((x) + (y))
#define ANN_DIFF(x,y)		((y) - (x))

//----------------------------------------------------------------------
//	Use the following for the L_1 (Manhattan) norm
//----------------------------------------------------------------------
// #define ANN_POW(v)		fabs(v)
// #define ANN_ROOT(x)		(x)
// #define ANN_SUM(x,y)		((x) + (y))
// #define ANN_DIFF(x,y)	((y) - (x))

//----------------------------------------------------------------------
//	Use the following for a general L_p norm
//----------------------------------------------------------------------
// #define ANN_POW(v)		pow(fabs(v),p)
// #define ANN_ROOT(x)		pow(fabs(x),1/p)
// #define ANN_SUM(x,y)		((x) + (y))
// #define ANN_DIFF(x,y)	((y) - (x))

//----------------------------------------------------------------------
//	Use the following for the L_infinity (Max) norm
//----------------------------------------------------------------------
// #define ANN_POW(v)		fabs(v)
// #define ANN_ROOT(x)		(x)
// #define ANN_SUM(x,y)		((x) > (y) ? (x) : (y))
// #define ANN_DIFF(x,y)	(y)

//----------------------------------------------------------------------
//	Array types
//		The following array types are of basic interest.  A point is
//		just a dimensionless array of coordinates, a point array is a
//		dimensionless array of points.  A distance array is a
//		dimensionless array of distances and an index array is a
//		dimensionless array of point indices.  The latter two are used
//		when returning the results of k-nearest neighbor queries.
//----------------------------------------------------------------------

typedef ANNcoord* ANNpoint;			// a point
typedef ANNpoint* ANNpointArray;	// an array of points 
typedef ANNdist*  ANNdistArray;		// an array of distances 
typedef ANNidx*   ANNidxArray;		// an array of point indices

//----------------------------------------------------------------------
//	Basic point and array utilities:
//		The following procedures are useful supplements to ANN's nearest
//		neighbor capabilities.
//
//		annDist():
//			Computes the (squared) distance between a pair of points.
//			Note that this routine is not used internally by ANN for
//			computing distance calculations.  For reasons of efficiency
//			this is done using incremental distance calculation.  Thus,
//			this routine cannot be modified as a method of changing the
//			metric.
//
//		Because points (somewhat like strings in C) are stored as
//		pointers.  Consequently, creating and destroying copies of
//		points may require storage allocation.  These procedures do
//		this.
//
//		annAllocPt() and annDeallocPt():
//				Allocate a deallocate storage for a single point, and
//				return a pointer to it.  The argument to AllocPt() is
//				used to initialize all components.
//
//		annAllocPts() and annDeallocPts():
//				Allocate and deallocate an array of points as well a
//				place to store their coordinates, and initializes the
//				points to point to their respective coordinates.  It
//				allocates point storage in a contiguous block large
//				enough to store all the points.  It performs no
//				initialization.
//
//		annCopyPt():
//				Creates a copy of a given point, allocating space for
//				the new point.  It returns a pointer to the newly
//				allocated copy.
//----------------------------------------------------------------------
   
DLL_API ANNdist annDist(
	int				dim,		// dimension of space
	ANNpoint		p,			// points
	ANNpoint		q);

DLL_API ANNpoint annAllocPt(
	int				dim,		// dimension
	ANNcoord		c = 0);		// coordinate value (all equal)

DLL_API ANNpointArray annAllocPts(
	int				n,			// number of points
	int				dim);		// dimension

DLL_API void annDeallocPt(
	ANNpoint		&p);		// deallocate 1 point
   
DLL_API void annDeallocPts(
	ANNpointArray	&pa);		// point array

DLL_API ANNpoint annCopyPt(
	int				dim,		// dimension
	ANNpoint		source);	// point to copy

//----------------------------------------------------------------------
//Overall structure: ANN supports a number of different data structures
//for approximate and exact nearest neighbor searching.  These are:
//
//		ANNbruteForce	A simple brute-force search structure.
//		ANNkd_tree		A kd-tree tree search structure.  ANNbd_tree
//		A bd-tree tree search structure (a kd-tree with shrink
//		capabilities).
//
//		At a minimum, each of these data structures support k-nearest
//		neighbor queries.  The nearest neighbor query, annkSearch,
//		returns an integer identifier and the distance to the nearest
//		neighbor(s) and annRangeSearch returns the nearest points that
//		lie within a given query ball.
//
//		Each structure is built by invoking the appropriate constructor
//		and passing it (at a minimum) the array of points, the total
//		number of points and the dimension of the space.  Each structure
//		is also assumed to support a destructor and member functions
//		that return basic information about the point set.
//
//		Note that the array of points is not copied by the data
//		structure (for reasons of space efficiency), and it is assumed
//		to be constant throughout the lifetime of the search structure.
//
//		The search algorithm, annkSearch, is given the query point (q),
//		and the desired number of nearest neighbors to report (k), and
//		the error bound (eps) (whose default value is 0, implying exact
//		nearest neighbors).  It returns two arrays which are assumed to
//		contain at least k elements: one (nn_idx) contains the indices
//		(within the point array) of the nearest neighbors and the other
//		(dd) contains the squared distances to these nearest neighbors.
//
//		The search algorithm, annkFRSearch, is a fixed-radius kNN
//		search.  In addition to a query point, it is given a (squared)
//		radius bound.  (This is done for consistency, because the search
//		returns distances as squared quantities.) It does two things.
//		First, it computes the k nearest neighbors within the radius
//		bound, and second, it returns the total number of points lying
//		within the radius bound. It is permitted to set k = 0, in which
//		case it effectively answers a range counting query.  If the
//		error bound epsilon is positive, then the search is approximate
//		in the sense that it is free to ignore any point that lies
//		outside a ball of radius r/(1+epsilon), where r is the given
//		(unsquared) radius bound.
//
//		The generic object from which all the search structures are
//		dervied is given below.  It is a virtual object, and is useless
//		by itself.
//----------------------------------------------------------------------

class DLL_API ANNpointSet {
public:
	virtual ~ANNpointSet() {}			// virtual distructor

	virtual void annkSearch(			// approx k near neighbor search
		ANNpoint		q,				// query point
		int				k,				// number of near neighbors to return
		ANNidxArray		nn_idx,			// nearest neighbor array (modified)
		ANNdistArray	dd,				// dist to near neighbors (modified)
		double			eps=0.0			// error bound
		) = 0;							// pure virtual (defined elsewhere)

	virtual int annkFRSearch(			// approx fixed-radius kNN search
		ANNpoint		q,				// query point
		ANNdist			sqRad,			// squared radius
		int				k = 0,			// number of near neighbors to return
		ANNidxArray		nn_idx = NULL,	// nearest neighbor array (modified)
		ANNdistArray	dd = NULL,		// dist to near neighbors (modified)
		double			eps=0.0			// error bound
		) = 0;							// pure virtual (defined elsewhere)

	virtual int theDim() = 0;			// return dimension of space
	virtual int nPoints() = 0;			// return number of points
										// return pointer to points
	virtual ANNpointArray thePoints() = 0;
};

//----------------------------------------------------------------------
//	Brute-force nearest neighbor search:
//		The brute-force search structure is very simple but inefficient.
//		It has been provided primarily for the sake of comparison with
//		and validation of the more complex search structures.
//
//		Query processing is the same as described above, but the value
//		of epsilon is ignored, since all distance calculations are
//		performed exactly.
//
//		WARNING: This data structure is very slow, and should not be
//		used unless the number of points is very small.
//
//		Internal information:
//		---------------------
//		This data structure bascially consists of the array of points
//		(each a pointer to an array of coordinates).  The search is
//		performed by a simple linear scan of all the points.
//----------------------------------------------------------------------

class DLL_API ANNbruteForce: public ANNpointSet {
	int				dim;				// dimension
	int				n_pts;				// number of points
	ANNpointArray	pts;				// point array
public:
	ANNbruteForce(						// constructor from point array
		ANNpointArray	pa,				// point array
		int				n,				// number of points
		int				dd);			// dimension

	~ANNbruteForce();					// destructor

	void annkSearch(					// approx k near neighbor search
		ANNpoint		q,				// query point
		int				k,				// number of near neighbors to return
		ANNidxArray		nn_idx,			// nearest neighbor array (modified)
		ANNdistArray	dd,				// dist to near neighbors (modified)
		double			eps=0.0);		// error bound