SHOGUN
v2.0.0
|
00001 #ifndef JLCOVERTREE_H 00002 #define JLCOVERTREE_H 00003 00004 #include <shogun/lib/JLCoverTreePoint.h> 00005 #include <shogun/mathematics/Math.h> 00006 00007 #include<math.h> 00008 #include<stdio.h> 00009 #define NDEBUG 00010 #include<assert.h> 00011 00012 /* First written by John Langford jl@hunch.net 00013 Templatization by Dinoj Surendran dinojs@gmail.com 00014 Adaptation to Shogun by Fernando José Iglesias García 00015 */ 00016 00017 // the files below may not need to be included 00018 00019 /* Whatever structure/class/type is used for P, it must have the following functions defined: 00020 00021 float distance(P v1, P v2, float upper_bound); 00022 : this returns the distance between two P objects 00023 : the distance does not have to be calculated fully if it's more than upper_bound 00024 00025 v_array<P> parse_points(char *filename); 00026 : this fills up a v_array of P objects from the input file 00027 00028 void print(point &P); 00029 : this prints out the contents of a P object. 00030 */ 00031 00032 using namespace std; 00033 using namespace shogun; 00034 00038 template<class P> 00039 struct node { 00040 00042 P p; 00043 00045 float max_dist; 00046 00048 float parent_dist; 00049 00051 node<P>* children; 00052 00054 unsigned short int num_children; 00055 00057 short int scale; 00058 00059 }; 00060 00061 //template<class P> 00062 //node<P> insert(P newpoint, node<P> *top_node); // not yet implemented 00063 // 00064 //template<class P> 00065 //void remove(P byepoint, node<P> *top_node); // not yet implemented 00066 //query 00067 00071 template<class P> 00072 struct ds_node { 00073 00075 v_array<float> dist; 00076 00078 P p; 00079 00080 }; 00081 00082 static float base = 1.3; 00083 00084 00085 static float il2 = 1. / log(base); 00086 00087 inline float dist_of_scale (int s) 00088 { 00089 return CMath::pow(base, s); 00090 } 00091 00092 inline int get_scale(float d) 00093 { 00094 return (int) CMath::ceil(il2 * log(d)); 00095 } 00096 00097 template<class P> 00098 node<P> new_node(const P &p) 00099 { 00100 node<P> new_node; 00101 new_node.p = p; 00102 return new_node; 00103 } 00104 00105 template<class P> 00106 node<P> new_leaf(const P &p) 00107 { 00108 node<P> new_leaf = {p,0.,0.,NULL,0,100}; 00109 return new_leaf; 00110 } 00111 00112 template<class P> 00113 float max_set(v_array<ds_node<P> > &v) 00114 { 00115 float max = 0.; 00116 for (int i = 0; i < v.index; i++) 00117 if ( max < v[i].dist.last()) 00118 max = v[i].dist.last(); 00119 return max; 00120 } 00121 00122 void print_space(int s) 00123 { 00124 for (int i = 0; i < s; i++) 00125 SG_SPRINT(" "); 00126 } 00127 00128 template<class P> 00129 void print(int depth, node<P> &top_node) 00130 { 00131 print_space(depth); 00132 print(top_node.p); 00133 if ( top_node.num_children > 0 ) 00134 { 00135 print_space(depth); 00136 SG_SPRINT("scale = %i\n",top_node.scale); 00137 print_space(depth); 00138 SG_SPRINT("max_dist = %f\n",top_node.max_dist); 00139 print_space(depth); 00140 SG_SPRINT("num children = %i\n",top_node.num_children); 00141 for (int i = 0; i < top_node.num_children;i++) 00142 print(depth+1, top_node.children[i]); 00143 } 00144 } 00145 00146 template<class P> 00147 void split(v_array<ds_node<P> >& point_set, v_array<ds_node<P> >& far_set, int max_scale) 00148 { 00149 unsigned int new_index = 0; 00150 float fmax = dist_of_scale(max_scale); 00151 for (int i = 0; i < point_set.index; i++) 00152 { 00153 if (point_set[i].dist.last() <= fmax) 00154 { 00155 point_set[new_index++] = point_set[i]; 00156 } 00157 else 00158 push(far_set,point_set[i]); 00159 } 00160 point_set.index=new_index; 00161 } 00162 00163 template<class P> 00164 void dist_split(v_array<ds_node<P> >& point_set, 00165 v_array<ds_node<P> >& new_point_set, 00166 P new_point, 00167 int max_scale) 00168 { 00169 unsigned int new_index = 0; 00170 float fmax = dist_of_scale(max_scale); 00171 for(int i = 0; i < point_set.index; i++) 00172 { 00173 float new_d; 00174 new_d = distance(new_point, point_set[i].p, fmax); 00175 if (new_d <= fmax ) 00176 { 00177 push(point_set[i].dist, new_d); 00178 push(new_point_set,point_set[i]); 00179 } 00180 else 00181 point_set[new_index++] = point_set[i]; 00182 } 00183 point_set.index = new_index; 00184 } 00185 00186 00187 00188 00189 /* 00190 max_scale is the maximum scale of the node we might create here. 00191 point_set contains points which are 2*max_scale or less away. 00192 */ 00193 00194 template <class P> 00195 node<P> batch_insert(const P& p, 00196 int max_scale, 00197 int top_scale, 00198 v_array<ds_node<P> >& point_set, 00199 v_array<ds_node<P> >& consumed_set, 00200 v_array<v_array<ds_node<P> > >& stack) 00201 { 00202 if (point_set.index == 0) 00203 return new_leaf(p); 00204 else { 00205 float max_dist = max_set(point_set); //O(|point_set|) 00206 int next_scale = CMath::min(max_scale - 1, get_scale(max_dist)); 00207 if (next_scale == -2147483647-1) // We have points with distance 0. 00208 { 00209 v_array<node<P> > children; 00210 push(children,new_leaf(p)); 00211 while (point_set.index > 0) 00212 { 00213 push(children,new_leaf(point_set.last().p)); 00214 push(consumed_set,point_set.last()); 00215 point_set.decr(); 00216 } 00217 node<P> n = new_node(p); 00218 n.scale = 100; // A magic number meant to be larger than all scales. 00219 n.max_dist = 0; 00220 alloc(children,children.index); 00221 n.num_children = children.index; 00222 n.children = children.elements; 00223 return n; 00224 } 00225 else 00226 { 00227 v_array<ds_node<P> > far = pop(stack); 00228 split(point_set,far,max_scale); //O(|point_set|) 00229 00230 node<P> child = batch_insert(p, next_scale, top_scale, point_set, consumed_set, stack); 00231 00232 if (point_set.index == 0) 00233 { 00234 push(stack,point_set); 00235 point_set=far; 00236 return child; 00237 } 00238 else { 00239 node<P> n = new_node(p); 00240 v_array<node<P> > children; 00241 push(children, child); 00242 v_array<ds_node<P> > new_point_set = pop(stack); 00243 v_array<ds_node<P> > new_consumed_set = pop(stack); 00244 while (point_set.index != 0) { //O(|point_set| * num_children) 00245 P new_point = point_set.last().p; 00246 float new_dist = point_set.last().dist.last(); 00247 push(consumed_set, point_set.last()); 00248 point_set.decr(); 00249 00250 dist_split(point_set, new_point_set, new_point, max_scale); //O(|point_saet|) 00251 dist_split(far,new_point_set,new_point,max_scale); //O(|far|) 00252 00253 node<P> new_child = 00254 batch_insert(new_point, next_scale, top_scale, new_point_set, new_consumed_set, stack); 00255 new_child.parent_dist = new_dist; 00256 00257 push(children, new_child); 00258 00259 float fmax = dist_of_scale(max_scale); 00260 for(int i = 0; i< new_point_set.index; i++) //O(|new_point_set|) 00261 { 00262 new_point_set[i].dist.decr(); 00263 if (new_point_set[i].dist.last() <= fmax) 00264 push(point_set, new_point_set[i]); 00265 else 00266 push(far, new_point_set[i]); 00267 } 00268 for(int i = 0; i< new_consumed_set.index; i++) //O(|new_point_set|) 00269 { 00270 new_consumed_set[i].dist.decr(); 00271 push(consumed_set, new_consumed_set[i]); 00272 } 00273 new_point_set.index = 0; 00274 new_consumed_set.index = 0; 00275 } 00276 push(stack,new_point_set); 00277 push(stack,new_consumed_set); 00278 push(stack,point_set); 00279 point_set=far; 00280 n.scale = top_scale - max_scale; 00281 n.max_dist = max_set(consumed_set); 00282 alloc(children,children.index); 00283 n.num_children = children.index; 00284 n.children = children.elements; 00285 return n; 00286 } 00287 } 00288 } 00289 } 00290 00291 template<class P> 00292 node<P> batch_create(v_array<P> points) 00293 { 00294 assert(points.index > 0); 00295 v_array<ds_node<P> > point_set; 00296 v_array<v_array<ds_node<P> > > stack; 00297 00298 for (int i = 1; i < points.index; i++) { 00299 ds_node<P> temp; 00300 push(temp.dist, distance(points[0], points[i], FLT_MAX)); 00301 temp.p = points[i]; 00302 push(point_set,temp); 00303 } 00304 00305 v_array<ds_node<P> > consumed_set; 00306 00307 float max_dist = max_set(point_set); 00308 00309 node<P> top = batch_insert (points[0], 00310 get_scale(max_dist), 00311 get_scale(max_dist), 00312 point_set, 00313 consumed_set, 00314 stack); 00315 for (int i = 0; i<consumed_set.index;i++) 00316 free(consumed_set[i].dist.elements); 00317 free(consumed_set.elements); 00318 for (int i = 0; i<stack.index;i++) 00319 free(stack[i].elements); 00320 free(stack.elements); 00321 free(point_set.elements); 00322 return top; 00323 } 00324 00325 void add_height(int d, v_array<int> &heights) 00326 { 00327 if (heights.index <= d) 00328 for(;heights.index <= d;) 00329 push(heights,0); 00330 heights[d] = heights[d] + 1; 00331 } 00332 00333 template <class P> 00334 int height_dist(const node<P> top_node,v_array<int> &heights) 00335 { 00336 if (top_node.num_children == 0) 00337 { 00338 add_height(0,heights); 00339 return 0; 00340 } 00341 else 00342 { 00343 int max_v=0; 00344 for (int i = 0; i<top_node.num_children ;i++) 00345 { 00346 int d = height_dist(top_node.children[i], heights); 00347 if (d > max_v) 00348 max_v = d; 00349 } 00350 add_height(1 + max_v, heights); 00351 return (1 + max_v); 00352 } 00353 } 00354 00355 template <class P> 00356 void depth_dist(int top_scale, const node<P> top_node,v_array<int> &depths) 00357 { 00358 if (top_node.num_children > 0) 00359 for (int i = 0; i<top_node.num_children ;i++) 00360 { 00361 add_height(top_node.scale, depths); 00362 depth_dist(top_scale, top_node.children[i], depths); 00363 } 00364 } 00365 00366 template <class P> 00367 void breadth_dist(const node<P> top_node,v_array<int> &breadths) 00368 { 00369 if (top_node.num_children == 0) 00370 add_height(0,breadths); 00371 else 00372 { 00373 for (int i = 0; i<top_node.num_children ;i++) 00374 breadth_dist(top_node.children[i], breadths); 00375 add_height(top_node.num_children, breadths); 00376 } 00377 } 00378 00382 template <class P> 00383 struct d_node { 00384 00386 float dist; 00387 00389 const node<P> *n; 00390 }; 00391 00392 template <class P> 00393 inline float compare(const d_node<P> *p1, const d_node<P>* p2) 00394 { 00395 return p1 -> dist - p2 -> dist; 00396 } 00397 00398 template <class P> 00399 void halfsort (v_array<d_node<P> > cover_set) 00400 { 00401 00402 if (cover_set.index <= 1) 00403 return; 00404 register d_node<P> *base_ptr = cover_set.elements; 00405 00406 d_node<P> *hi = &base_ptr[cover_set.index - 1]; 00407 d_node<P> *right_ptr = hi; 00408 d_node<P> *left_ptr; 00409 00410 while (right_ptr > base_ptr) 00411 { 00412 d_node<P> *mid = base_ptr + ((hi - base_ptr) >> 1); 00413 00414 if (compare ( mid, base_ptr) < 0.) 00415 CMath::swap(*mid, *base_ptr); 00416 if (compare ( hi, mid) < 0.) 00417 CMath::swap(*mid, *hi); 00418 else 00419 goto jump_over; 00420 if (compare ( mid, base_ptr) < 0.) 00421 CMath::swap(*mid, *base_ptr); 00422 jump_over:; 00423 00424 left_ptr = base_ptr + 1; 00425 right_ptr = hi - 1; 00426 00427 do 00428 { 00429 while (compare (left_ptr, mid) < 0.) 00430 left_ptr ++; 00431 00432 while (compare (mid, right_ptr) < 0.) 00433 right_ptr --; 00434 00435 if (left_ptr < right_ptr) 00436 { 00437 CMath::swap(*left_ptr, *right_ptr); 00438 if (mid == left_ptr) 00439 mid = right_ptr; 00440 else if (mid == right_ptr) 00441 mid = left_ptr; 00442 left_ptr ++; 00443 right_ptr --; 00444 } 00445 else if (left_ptr == right_ptr) 00446 { 00447 left_ptr ++; 00448 right_ptr --; 00449 break; 00450 } 00451 } 00452 while (left_ptr <= right_ptr); 00453 00454 hi = right_ptr; 00455 } 00456 } 00457 00458 00459 template <class P> 00460 v_array<v_array<d_node<P> > > get_cover_sets(v_array<v_array<v_array<d_node<P> > > > &spare_cover_sets) 00461 { 00462 v_array<v_array<d_node<P> > > ret = pop(spare_cover_sets); 00463 while (ret.index < 101) 00464 { 00465 v_array<d_node<P> > temp; 00466 push(ret, temp); 00467 } 00468 return ret; 00469 } 00470 00471 inline bool shell(float parent_query_dist, float child_parent_dist, float upper_bound) 00472 { 00473 return parent_query_dist - child_parent_dist <= upper_bound; 00474 // && child_parent_dist - parent_query_dist <= upper_bound; 00475 } 00476 00477 int internal_k =1; 00478 void update_k(float *k_upper_bound, float upper_bound) 00479 { 00480 float *end = k_upper_bound + internal_k-1; 00481 float *begin = k_upper_bound; 00482 for (;end != begin; begin++) 00483 { 00484 if (upper_bound < *(begin+1)) 00485 *begin = *(begin+1); 00486 else { 00487 *begin = upper_bound; 00488 break; 00489 } 00490 } 00491 if (end == begin) 00492 *begin = upper_bound; 00493 } 00494 float *alloc_k() 00495 { 00496 return (float *)malloc(sizeof(float) * internal_k); 00497 } 00498 void set_k(float* begin, float max) 00499 { 00500 for(float *end = begin+internal_k;end != begin; begin++) 00501 *begin = max; 00502 } 00503 00504 float internal_epsilon =0.; 00505 void update_epsilon(float *upper_bound, float new_dist) {} 00506 float *alloc_epsilon() 00507 { 00508 return (float *)malloc(sizeof(float)); 00509 } 00510 void set_epsilon(float* begin, float max) 00511 { 00512 *begin = internal_epsilon; 00513 } 00514 00515 void update_unequal(float *upper_bound, float new_dist) 00516 { 00517 if (new_dist != 0.) 00518 *upper_bound = new_dist; 00519 } 00520 float* (*alloc_unequal)() = alloc_epsilon; 00521 void set_unequal(float* begin, float max) 00522 { 00523 *begin = max; 00524 } 00525 00526 void (*update)(float *foo, float bar) = update_k; 00527 void (*setter)(float *foo, float bar) = set_k; 00528 float* (*alloc_upper)() = alloc_k; 00529 00530 template <class P> 00531 inline void copy_zero_set(node<P>* query_chi, float* new_upper_bound, 00532 v_array<d_node<P> > &zero_set, v_array<d_node<P> > &new_zero_set) 00533 { 00534 new_zero_set.index = 0; 00535 d_node<P> *end = zero_set.elements + zero_set.index; 00536 for (d_node<P> *ele = zero_set.elements; ele != end ; ele++) 00537 { 00538 float upper_dist = *new_upper_bound + query_chi->max_dist; 00539 if (shell(ele->dist, query_chi->parent_dist, upper_dist)) 00540 { 00541 float d = distance(query_chi->p, ele->n->p, upper_dist); 00542 00543 if (d <= upper_dist) 00544 { 00545 if (d < *new_upper_bound) 00546 update(new_upper_bound, d); 00547 d_node<P> temp = {d, ele->n}; 00548 push(new_zero_set,temp); 00549 } 00550 } 00551 } 00552 } 00553 00554 template <class P> 00555 inline void copy_cover_sets(node<P>* query_chi, float* new_upper_bound, 00556 v_array<v_array<d_node<P> > > &cover_sets, 00557 v_array<v_array<d_node<P> > > &new_cover_sets, 00558 int current_scale, int max_scale) 00559 { 00560 for (; current_scale <= max_scale; current_scale++) 00561 { 00562 d_node<P>* ele = cover_sets[current_scale].elements; 00563 d_node<P>* end = cover_sets[current_scale].elements + cover_sets[current_scale].index; 00564 for (; ele != end; ele++) 00565 { 00566 float upper_dist = *new_upper_bound + query_chi->max_dist + ele->n->max_dist; 00567 if (shell(ele->dist, query_chi->parent_dist, upper_dist)) 00568 { 00569 float d = distance(query_chi->p, ele->n->p, upper_dist); 00570 00571 if (d <= upper_dist) 00572 { 00573 if (d < *new_upper_bound) 00574 update(new_upper_bound,d); 00575 d_node<P> temp = {d, ele->n}; 00576 push(new_cover_sets[current_scale],temp); 00577 } 00578 } 00579 } 00580 } 00581 } 00582 00583 template <class P> 00584 void print_query(const node<P> *top_node) 00585 { 00586 SG_SPRINT ("query = \n"); 00587 print(top_node->p); 00588 if ( top_node->num_children > 0 ) { 00589 SG_SPRINT("scale = %i\n",top_node->scale); 00590 SG_SPRINT("max_dist = %f\n",top_node->max_dist); 00591 SG_SPRINT("num children = %i\n",top_node->num_children); 00592 } 00593 } 00594 00595 template <class P> 00596 void print_cover_sets(v_array<v_array<d_node<P> > > &cover_sets, 00597 v_array<d_node<P> > &zero_set, 00598 int current_scale, int max_scale) 00599 { 00600 SG_SPRINT("cover set = \n"); 00601 for (; current_scale <= max_scale; current_scale++) 00602 { 00603 d_node<P> *ele = cover_sets[current_scale].elements; 00604 d_node<P> *end = cover_sets[current_scale].elements + cover_sets[current_scale].index; 00605 SG_SPRINT ("%i\n", current_scale); 00606 for (; ele != end; ele++) 00607 { 00608 node<P> *n = (node<P> *)ele->n; 00609 print(n->p); 00610 } 00611 } 00612 d_node<P> *end = zero_set.elements + zero_set.index; 00613 SG_SPRINT ("infinity\n"); 00614 for (d_node<P> *ele = zero_set.elements; ele != end ; ele++) 00615 { 00616 node<P> *n = (node<P> *)ele->n; 00617 print(n->p); 00618 } 00619 } 00620 00621 00622 /* 00623 An optimization to consider: 00624 Make all distance evaluations occur in descend. 00625 00626 Instead of passing a cover_set, pass a stack of cover sets. The 00627 last element holds d_nodes with your distance. The next lower 00628 element holds a d_node with the distance to your query parent, 00629 next = query grand parent, etc.. 00630 00631 Compute distances in the presence of the tighter upper bound. 00632 */ 00633 00634 template <class P> 00635 inline 00636 void descend(const node<P>* query, float* upper_bound, 00637 int current_scale, 00638 int &max_scale, v_array<v_array<d_node<P> > > &cover_sets, 00639 v_array<d_node<P> > &zero_set) 00640 { 00641 d_node<P> *end = cover_sets[current_scale].elements + cover_sets[current_scale].index; 00642 for (d_node<P> *parent = cover_sets[current_scale].elements; parent != end; parent++) 00643 { 00644 const node<P> *par = parent->n; 00645 float upper_dist = *upper_bound + query->max_dist + query->max_dist; 00646 if (parent->dist <= upper_dist + par->max_dist) 00647 { 00648 node<P> *chi = par->children; 00649 if (parent->dist <= upper_dist + chi->max_dist) 00650 { 00651 if (chi->num_children > 0) 00652 { 00653 if (max_scale < chi->scale) 00654 max_scale = chi->scale; 00655 d_node<P> temp = {parent->dist, chi}; 00656 push(cover_sets[chi->scale], temp); 00657 } 00658 else if (parent->dist <= upper_dist) 00659 { 00660 d_node<P> temp = {parent->dist, chi}; 00661 push(zero_set, temp); 00662 } 00663 } 00664 node<P> *child_end = par->children + par->num_children; 00665 for (chi++; chi != child_end; chi++) 00666 { 00667 float upper_chi = *upper_bound + chi->max_dist + query->max_dist + query->max_dist; 00668 if (shell(parent->dist, chi->parent_dist, upper_chi)) 00669 { 00670 float d = distance(query->p, chi->p, upper_chi); 00671 if (d <= upper_chi) 00672 { 00673 if (d < *upper_bound) 00674 update(upper_bound, d); 00675 if (chi->num_children > 0) 00676 { 00677 if (max_scale < chi->scale) 00678 max_scale = chi->scale; 00679 d_node<P> temp = {d, chi}; 00680 push(cover_sets[chi->scale],temp); 00681 } 00682 else 00683 if (d <= upper_chi - chi->max_dist) 00684 { 00685 d_node<P> temp = {d, chi}; 00686 push(zero_set, temp); 00687 } 00688 } 00689 } 00690 } 00691 } 00692 } 00693 } 00694 00695 template <class P> 00696 void brute_nearest(const node<P>* query,v_array<d_node<P> > zero_set, 00697 float* upper_bound, 00698 v_array<v_array<P> > &results, 00699 v_array<v_array<d_node<P> > > &spare_zero_sets) 00700 { 00701 if (query->num_children > 0) 00702 { 00703 v_array<d_node<P> > new_zero_set = pop(spare_zero_sets); 00704 node<P> * query_chi = query->children; 00705 brute_nearest(query_chi, zero_set, upper_bound, results, spare_zero_sets); 00706 float* new_upper_bound = alloc_upper(); 00707 00708 node<P> *child_end = query->children + query->num_children; 00709 for (query_chi++;query_chi != child_end; query_chi++) 00710 { 00711 setter(new_upper_bound,*upper_bound + query_chi->parent_dist); 00712 copy_zero_set(query_chi, new_upper_bound, zero_set, new_zero_set); 00713 brute_nearest(query_chi, new_zero_set, new_upper_bound, results, spare_zero_sets); 00714 } 00715 free (new_upper_bound); 00716 new_zero_set.index = 0; 00717 push(spare_zero_sets, new_zero_set); 00718 } 00719 else 00720 { 00721 v_array<P> temp; 00722 push(temp, query->p); 00723 d_node<P> *end = zero_set.elements + zero_set.index; 00724 for (d_node<P> *ele = zero_set.elements; ele != end ; ele++) 00725 if (ele->dist <= *upper_bound) 00726 push(temp, ele->n->p); 00727 push(results,temp); 00728 } 00729 } 00730 00731 template <class P> 00732 void internal_batch_nearest_neighbor(const node<P> *query, 00733 v_array<v_array<d_node<P> > > &cover_sets, 00734 v_array<d_node<P> > &zero_set, 00735 int current_scale, 00736 int max_scale, 00737 float* upper_bound, 00738 v_array<v_array<P> > &results, 00739 v_array<v_array<v_array<d_node<P> > > > &spare_cover_sets, 00740 v_array<v_array<d_node<P> > > &spare_zero_sets) 00741 { 00742 if (current_scale > max_scale) // All remaining points are in the zero set. 00743 brute_nearest(query, zero_set, upper_bound, results, spare_zero_sets); 00744 else 00745 if (query->scale <= current_scale && query->scale != 100) 00746 // Our query has too much scale. Reduce. 00747 { 00748 node<P> *query_chi = query->children; 00749 v_array<d_node<P> > new_zero_set = pop(spare_zero_sets); 00750 v_array<v_array<d_node<P> > > new_cover_sets = get_cover_sets(spare_cover_sets); 00751 float* new_upper_bound = alloc_upper(); 00752 00753 node<P> *child_end = query->children + query->num_children; 00754 for (query_chi++; query_chi != child_end; query_chi++) 00755 { 00756 setter(new_upper_bound,*upper_bound + query_chi->parent_dist); 00757 copy_zero_set(query_chi, new_upper_bound, zero_set, new_zero_set); 00758 copy_cover_sets(query_chi, new_upper_bound, cover_sets, new_cover_sets, 00759 current_scale, max_scale); 00760 internal_batch_nearest_neighbor(query_chi, new_cover_sets, new_zero_set, 00761 current_scale, max_scale, new_upper_bound, 00762 results, spare_cover_sets, spare_zero_sets); 00763 } 00764 free (new_upper_bound); 00765 new_zero_set.index = 0; 00766 push(spare_zero_sets, new_zero_set); 00767 push(spare_cover_sets, new_cover_sets); 00768 internal_batch_nearest_neighbor(query->children, cover_sets, zero_set, 00769 current_scale, max_scale, upper_bound, results, 00770 spare_cover_sets, spare_zero_sets); 00771 } 00772 else // reduce cover set scale 00773 { 00774 halfsort(cover_sets[current_scale]); 00775 descend(query, upper_bound, current_scale, max_scale,cover_sets, zero_set); 00776 cover_sets[current_scale++].index = 0; 00777 internal_batch_nearest_neighbor(query, cover_sets, zero_set, 00778 current_scale, max_scale, upper_bound, results, 00779 spare_cover_sets, spare_zero_sets); 00780 } 00781 } 00782 00783 template <class P> 00784 void batch_nearest_neighbor(const node<P> &top_node, const node<P> &query, 00785 v_array<v_array<P> > &results) 00786 { 00787 v_array<v_array<v_array<d_node<P> > > > spare_cover_sets; 00788 v_array<v_array<d_node<P> > > spare_zero_sets; 00789 00790 v_array<v_array<d_node<P> > > cover_sets = get_cover_sets(spare_cover_sets); 00791 v_array<d_node<P> > zero_set = pop(spare_zero_sets); 00792 00793 float* upper_bound = alloc_upper(); 00794 setter(upper_bound,FLT_MAX); 00795 00796 float top_dist = distance(query.p, top_node.p, FLT_MAX); 00797 update(upper_bound, top_dist); 00798 00799 d_node<P> temp = {top_dist, &top_node}; 00800 push(cover_sets[0], temp); 00801 00802 internal_batch_nearest_neighbor(&query,cover_sets,zero_set,0,0,upper_bound,results, 00803 spare_cover_sets,spare_zero_sets); 00804 00805 free(upper_bound); 00806 push(spare_cover_sets, cover_sets); 00807 00808 for (int i = 0; i < spare_cover_sets.index; i++) 00809 { 00810 v_array<v_array<d_node<P> > > cover_sets2 = spare_cover_sets[i]; 00811 for (int j = 0; j < cover_sets2.index; j++) 00812 free (cover_sets2[j].elements); 00813 free(cover_sets2.elements); 00814 } 00815 free(spare_cover_sets.elements); 00816 00817 push(spare_zero_sets, zero_set); 00818 00819 for (int i = 0; i < spare_zero_sets.index; i++) 00820 free(spare_zero_sets[i].elements); 00821 free(spare_zero_sets.elements); 00822 } 00823 00824 template <class P> 00825 void k_nearest_neighbor(const node<P> &top_node, const node<P> &query, 00826 v_array<v_array<P> > &results, int k) 00827 { 00828 internal_k = k; 00829 update = update_k; 00830 setter = set_k; 00831 alloc_upper = alloc_k; 00832 00833 batch_nearest_neighbor(top_node, query,results); 00834 } 00835 00836 template <class P> 00837 void epsilon_nearest_neighbor(const node<P> &top_node, const node<P> &query, 00838 v_array<v_array<P> > &results, float epsilon) 00839 { 00840 internal_epsilon = epsilon; 00841 update = update_epsilon; 00842 setter = set_epsilon; 00843 alloc_upper = alloc_epsilon; 00844 00845 batch_nearest_neighbor(top_node, query,results); 00846 } 00847 00848 template <class P> 00849 void unequal_nearest_neighbor(const node<P> &top_node, const node<P> &query, 00850 v_array<v_array<P> > &results) 00851 { 00852 update = update_unequal; 00853 setter = set_unequal; 00854 alloc_upper = alloc_unequal; 00855 00856 batch_nearest_neighbor(top_node, query, results); 00857 } 00858 00859 #endif