00001 //---------------------------------------------------------------------- 00002 // File: kd_fix_rad_search.cpp 00003 // Programmer: Sunil Arya and David Mount 00004 // Description: Standard kd-tree fixed-radius kNN search 00005 // Last modified: 05/03/05 (Version 1.1) 00006 //---------------------------------------------------------------------- 00007 // Copyright (c) 1997-2005 University of Maryland and Sunil Arya and 00008 // David Mount. All Rights Reserved. 00009 // 00010 // This software and related documentation is part of the Approximate 00011 // Nearest Neighbor Library (ANN). This software is provided under 00012 // the provisions of the Lesser GNU Public License (LGPL). See the 00013 // file ../ReadMe.txt for further information. 00014 // 00015 // The University of Maryland (U.M.) and the authors make no 00016 // representations about the suitability or fitness of this software for 00017 // any purpose. It is provided "as is" without express or implied 00018 // warranty. 00019 //---------------------------------------------------------------------- 00020 // History: 00021 // Revision 1.1 05/03/05 00022 // Initial release 00023 //---------------------------------------------------------------------- 00024 00025 #include "kd_fix_rad_search.h" // kd fixed-radius search decls 00026 00027 //---------------------------------------------------------------------- 00028 // Approximate fixed-radius k nearest neighbor search 00029 // The squared radius is provided, and this procedure finds the 00030 // k nearest neighbors within the radius, and returns the total 00031 // number of points lying within the radius. 00032 // 00033 // The method used for searching the kd-tree is a variation of the 00034 // nearest neighbor search used in kd_search.cpp, except that the 00035 // radius of the search ball is known. We refer the reader to that 00036 // file for the explanation of the recursive search procedure. 00037 //---------------------------------------------------------------------- 00038 00039 //---------------------------------------------------------------------- 00040 // To keep argument lists short, a number of global variables 00041 // are maintained which are common to all the recursive calls. 00042 // These are given below. 00043 //---------------------------------------------------------------------- 00044 00045 int ANNkdFRDim; // dimension of space 00046 ANNpoint ANNkdFRQ; // query point 00047 ANNdist ANNkdFRSqRad; // squared radius search bound 00048 double ANNkdFRMaxErr; // max tolerable squared error 00049 ANNpointArray ANNkdFRPts; // the points 00050 ANNmin_k* ANNkdFRPointMK; // set of k closest points 00051 int ANNkdFRPtsVisited; // total points visited 00052 int ANNkdFRPtsInRange; // number of points in the range 00053 00054 //---------------------------------------------------------------------- 00055 // annkFRSearch - fixed radius search for k nearest neighbors 00056 //---------------------------------------------------------------------- 00057 00058 int ANNkd_tree::annkFRSearch( 00059 ANNpoint q, // the query point 00060 ANNdist sqRad, // squared radius search bound 00061 int k, // number of near neighbors to return 00062 ANNidxArray nn_idx, // nearest neighbor indices (returned) 00063 ANNdistArray dd, // the approximate nearest neighbor 00064 double eps) // the error bound 00065 { 00066 ANNkdFRDim = dim; // copy arguments to static equivs 00067 ANNkdFRQ = q; 00068 ANNkdFRSqRad = sqRad; 00069 ANNkdFRPts = pts; 00070 ANNkdFRPtsVisited = 0; // initialize count of points visited 00071 ANNkdFRPtsInRange = 0; // ...and points in the range 00072 00073 ANNkdFRMaxErr = ANN_POW(1.0 + eps); 00074 ANN_FLOP(2) // increment floating op count 00075 00076 ANNkdFRPointMK = new ANNmin_k(k); // create set for closest k points 00077 // search starting at the root 00078 root->ann_FR_search(annBoxDistance(q, bnd_box_lo, bnd_box_hi, dim)); 00079 00080 for (int i = 0; i < k; i++) { // extract the k-th closest points 00081 if (dd != NULL) 00082 dd[i] = ANNkdFRPointMK->ith_smallest_key(i); 00083 if (nn_idx != NULL) 00084 nn_idx[i] = ANNkdFRPointMK->ith_smallest_info(i); 00085 } 00086 00087 delete ANNkdFRPointMK; // deallocate closest point set 00088 return ANNkdFRPtsInRange; // return final point count 00089 } 00090 00091 //---------------------------------------------------------------------- 00092 // kd_split::ann_FR_search - search a splitting node 00093 // Note: This routine is similar in structure to the standard kNN 00094 // search. It visits the subtree that is closer to the query point 00095 // first. For fixed-radius search, there is no benefit in visiting 00096 // one subtree before the other, but we maintain the same basic 00097 // code structure for the sake of uniformity. 00098 //---------------------------------------------------------------------- 00099 00100 void ANNkd_split::ann_FR_search(ANNdist box_dist) 00101 { 00102 // check dist calc term condition 00103 if (ANNmaxPtsVisited != 0 && ANNkdFRPtsVisited > ANNmaxPtsVisited) return; 00104 00105 // distance to cutting plane 00106 ANNcoord cut_diff = ANNkdFRQ[cut_dim] - cut_val; 00107 00108 if (cut_diff < 0) { // left of cutting plane 00109 child[ANN_LO]->ann_FR_search(box_dist);// visit closer child first 00110 00111 ANNcoord box_diff = cd_bnds[ANN_LO] - ANNkdFRQ[cut_dim]; 00112 if (box_diff < 0) // within bounds - ignore 00113 box_diff = 0; 00114 // distance to further box 00115 box_dist = (ANNdist) ANN_SUM(box_dist, 00116 ANN_DIFF(ANN_POW(box_diff), ANN_POW(cut_diff))); 00117 00118 // visit further child if in range 00119 if (box_dist * ANNkdFRMaxErr <= ANNkdFRSqRad) 00120 child[ANN_HI]->ann_FR_search(box_dist); 00121 00122 } 00123 else { // right of cutting plane 00124 child[ANN_HI]->ann_FR_search(box_dist);// visit closer child first 00125 00126 ANNcoord box_diff = ANNkdFRQ[cut_dim] - cd_bnds[ANN_HI]; 00127 if (box_diff < 0) // within bounds - ignore 00128 box_diff = 0; 00129 // distance to further box 00130 box_dist = (ANNdist) ANN_SUM(box_dist, 00131 ANN_DIFF(ANN_POW(box_diff), ANN_POW(cut_diff))); 00132 00133 // visit further child if close enough 00134 if (box_dist * ANNkdFRMaxErr <= ANNkdFRSqRad) 00135 child[ANN_LO]->ann_FR_search(box_dist); 00136 00137 } 00138 ANN_FLOP(13) // increment floating ops 00139 ANN_SPL(1) // one more splitting node visited 00140 } 00141 00142 //---------------------------------------------------------------------- 00143 // kd_leaf::ann_FR_search - search points in a leaf node 00144 // Note: The unreadability of this code is the result of 00145 // some fine tuning to replace indexing by pointer operations. 00146 //---------------------------------------------------------------------- 00147 00148 void ANNkd_leaf::ann_FR_search(ANNdist box_dist) 00149 { 00150 register ANNdist dist; // distance to data point 00151 register ANNcoord* pp; // data coordinate pointer 00152 register ANNcoord* qq; // query coordinate pointer 00153 register ANNcoord t; 00154 register int d; 00155 00156 for (int i = 0; i < n_pts; i++) { // check points in bucket 00157 00158 pp = ANNkdFRPts[bkt[i]]; // first coord of next data point 00159 qq = ANNkdFRQ; // first coord of query point 00160 dist = 0; 00161 00162 for(d = 0; d < ANNkdFRDim; d++) { 00163 ANN_COORD(1) // one more coordinate hit 00164 ANN_FLOP(5) // increment floating ops 00165 00166 t = *(qq++) - *(pp++); // compute length and adv coordinate 00167 // exceeds dist to k-th smallest? 00168 if( (dist = ANN_SUM(dist, ANN_POW(t))) > ANNkdFRSqRad) { 00169 break; 00170 } 00171 } 00172 00173 if (d >= ANNkdFRDim && // among the k best? 00174 (ANN_ALLOW_SELF_MATCH || dist!=0)) { // and no self-match problem 00175 // add it to the list 00176 ANNkdFRPointMK->insert(dist, bkt[i]); 00177 ANNkdFRPtsInRange++; // increment point count 00178 } 00179 } 00180 ANN_LEAF(1) // one more leaf node visited 00181 ANN_PTS(n_pts) // increment points visited 00182 ANNkdFRPtsVisited += n_pts; // increment number of points visited 00183 }