kd_fix_rad_search.cpp

Go to the documentation of this file.
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 }

Generated on Wed Nov 23 19:00:20 2011 for FreeCAD by  doxygen 1.6.1