Reason of all this work: Commiting my work-in-progress on reviewed collision system...
[blender.git] / source / blender / blenlib / intern / BLI_kdopbvh.c
1 /*  kdop.c      
2
3 *
4 * ***** BEGIN GPL/BL DUAL LICENSE BLOCK *****
5 *
6 * This program is free software; you can redistribute it and/or
7 * modify it under the terms of the GNU General Public License
8 * as published by the Free Software Foundation; either version 2
9 * of the License, or (at your option) any later version. The Blender
10 * Foundation also sells licenses for use in proprietary software under
11 * the Blender License.  See http://www.blender.org/BL/ for information
12 * about this.
13 *
14 * This program is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17 * GNU General Public License for more details.
18 *
19 * You should have received a copy of the GNU General Public License
20 * along with this program; if not, write to the Free Software Foundation,
21 * Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
22 *
23 * The Original Code is Copyright (C) Blender Foundation
24 * All rights reserved.
25 *
26 * The Original Code is: all of this file.
27 *
28 * Contributor(s): Daniel Genrich
29 *
30 * ***** END GPL/BL DUAL LICENSE BLOCK *****
31 */
32
33 #include "math.h"
34 #include <stdio.h>
35 #include <stdlib.h> 
36 #include <string.h>
37
38 #include "MEM_guardedalloc.h"
39
40 #include "BKE_utildefines.h"
41
42 #include "BLI_kdopbvh.h"
43
44 #ifdef _OPENMP
45 #include <omp.h>
46 #endif
47
48 typedef struct BVHNode
49 {
50         struct BVHNode *children[8]; // max 8 children
51         struct BVHNode *parent; // needed for bottom - top update
52         float bv[26]; // Bounding volume of all nodes, max 13 axis
53         int index; /* face, edge, vertex index */
54         char totnode; // how many nodes are used, used for speedup
55         char traversed;  // how many nodes already traversed until this level?
56         char main_axis;
57 } BVHNode;
58
59 struct BVHTree
60 {
61         BVHNode **nodes;
62         BVHNode *nodearray; /* pre-alloc branch nodes */
63         int     *orig_index; /* mapping for orig_index to node_index */
64         float   epsilon; /* epslion is used for inflation of the k-dop     */
65         int     totleaf; // leafs
66         int     totbranch;
67         char    tree_type; // type of tree (4 => quadtree)
68         char    axis; // kdop type (6 => OBB, 7 => AABB, ...)
69         char    start_axis, stop_axis; // KDOP_AXES array indices according to axis
70 };
71
72 typedef struct BVHOverlapData 
73 {  
74         BVHTree *tree1, *tree2; 
75         BVHTreeOverlap *overlap; 
76         int i, max_overlap; /* i is number of overlaps */
77 } BVHOverlapData;
78 ////////////////////////////////////////
79
80
81 ////////////////////////////////////////////////////////////////////////
82 // Bounding Volume Hierarchy Definition
83 // 
84 // Notes: From OBB until 26-DOP --> all bounding volumes possible, just choose type below
85 // Notes: You have to choose the type at compile time ITM
86 // Notes: You can choose the tree type --> binary, quad, octree, choose below
87 ////////////////////////////////////////////////////////////////////////
88
89 static float KDOP_AXES[13][3] =
90 { {1.0, 0, 0}, {0, 1.0, 0}, {0, 0, 1.0}, {1.0, 1.0, 1.0}, {1.0, -1.0, 1.0}, {1.0, 1.0, -1.0},
91 {1.0, -1.0, -1.0}, {1.0, 1.0, 0}, {1.0, 0, 1.0}, {0, 1.0, 1.0}, {1.0, -1.0, 0}, {1.0, 0, -1.0},
92 {0, 1.0, -1.0}
93 };
94
95 //////////////////////////////////////////////////////////////////////////////////////////////////////
96 // Introsort 
97 // with permission deriven from the following Java code:
98 // http://ralphunden.net/content/tutorials/a-guide-to-introsort/
99 // and he derived it from the SUN STL 
100 //////////////////////////////////////////////////////////////////////////////////////////////////////
101 static int size_threshold = 16;
102 /*
103 * Common methods for all algorithms
104 */
105 static void bvh_exchange(BVHNode **a, int i, int j)
106 {
107         BVHNode *t=a[i];
108         a[i]=a[j];
109         a[j]=t;
110 }
111 static int floor_lg(int a)
112 {
113         return (int)(floor(log(a)/log(2)));
114 }
115
116 /*
117 * Insertion sort algorithm
118 */
119 static void bvh_insertionsort(BVHNode **a, int lo, int hi, int axis)
120 {
121         int i,j;
122         BVHNode *t;
123         for (i=lo; i < hi; i++)
124         {
125                 j=i;
126                 t = a[i];
127                 while((j!=lo) && (t->bv[axis] < (a[j-1])->bv[axis]))
128                 {
129                         a[j] = a[j-1];
130                         j--;
131                 }
132                 a[j] = t;
133         }
134 }
135
136 static int bvh_partition(BVHNode **a, int lo, int hi, BVHNode * x, int axis)
137 {
138         int i=lo, j=hi;
139         while (1)
140         {
141                 while ((a[i])->bv[axis] < x->bv[axis]) i++;
142                 j=j-1;
143                 while (x->bv[axis] < (a[j])->bv[axis]) j=j-1;
144                 if(!(i < j))
145                         return i;
146                 bvh_exchange(a, i,j);
147                 i++;
148         }
149 }
150
151 /*
152 * Heapsort algorithm
153 */
154 static void bvh_downheap(BVHNode **a, int i, int n, int lo, int axis)
155 {
156         BVHNode * d = a[lo+i-1];
157         int child;
158         while (i<=n/2)
159         {
160                 child = 2*i;
161                 if ((child < n) && ((a[lo+child-1])->bv[axis] < (a[lo+child])->bv[axis]))
162                 {
163                         child++;
164                 }
165                 if (!(d->bv[axis] < (a[lo+child-1])->bv[axis])) break;
166                 a[lo+i-1] = a[lo+child-1];
167                 i = child;
168         }
169         a[lo+i-1] = d;
170 }
171
172 static void bvh_heapsort(BVHNode **a, int lo, int hi, int axis)
173 {
174         int n = hi-lo, i;
175         for (i=n/2; i>=1; i=i-1)
176         {
177                 bvh_downheap(a, i,n,lo, axis);
178         }
179         for (i=n; i>1; i=i-1)
180         {
181                 bvh_exchange(a, lo,lo+i-1);
182                 bvh_downheap(a, 1,i-1,lo, axis);
183         }
184 }
185
186 static BVHNode *bvh_medianof3(BVHNode **a, int lo, int mid, int hi, int axis) // returns Sortable
187 {
188         if ((a[mid])->bv[axis] < (a[lo])->bv[axis])
189         {
190                 if ((a[hi])->bv[axis] < (a[mid])->bv[axis])
191                         return a[mid];
192                 else
193                 {
194                         if ((a[hi])->bv[axis] < (a[lo])->bv[axis])
195                                 return a[hi];
196                         else
197                                 return a[lo];
198                 }
199         }
200         else
201         {
202                 if ((a[hi])->bv[axis] < (a[mid])->bv[axis])
203                 {
204                         if ((a[hi])->bv[axis] < (a[lo])->bv[axis])
205                                 return a[lo];
206                         else
207                                 return a[hi];
208                 }
209                 else
210                         return a[mid];
211         }
212 }
213 /*
214 * Quicksort algorithm modified for Introsort
215 */
216 static void bvh_introsort_loop (BVHNode **a, int lo, int hi, int depth_limit, int axis)
217 {
218         int p;
219
220         while (hi-lo > size_threshold)
221         {
222                 if (depth_limit == 0)
223                 {
224                         bvh_heapsort(a, lo, hi, axis);
225                         return;
226                 }
227                 depth_limit=depth_limit-1;
228                 p=bvh_partition(a, lo, hi, bvh_medianof3(a, lo, lo+((hi-lo)/2)+1, hi-1, axis), axis);
229                 bvh_introsort_loop(a, p, hi, depth_limit, axis);
230                 hi=p;
231         }
232 }
233
234 static void sort(BVHNode **a0, int begin, int end, int axis)
235 {
236         if (begin < end)
237         {
238                 BVHNode **a=a0;
239                 bvh_introsort_loop(a, begin, end, 2*floor_lg(end-begin), axis);
240                 bvh_insertionsort(a, begin, end, axis);
241         }
242 }
243 void sort_along_axis(BVHTree *tree, int start, int end, int axis)
244 {
245         sort(tree->nodes, start, end, axis);
246 }
247
248 //////////////////////////////////////////////////////////////////////////////////////////////////////
249
250 void BLI_bvhtree_free(BVHTree *tree)
251 {       
252         if(tree)
253         {
254                 MEM_freeN(tree->nodes);
255                 MEM_freeN(tree->nodearray);
256                 MEM_freeN(tree->orig_index);
257                 MEM_freeN(tree);
258         }
259 }
260
261 BVHTree *BLI_bvhtree_new(int maxsize, float epsilon, char tree_type, char axis)
262 {
263         BVHTree *tree;
264         int numbranches=0, i;
265         
266         // only support up to octree
267         if(tree_type > 8)
268                 return NULL;
269
270         tree = (BVHTree *)MEM_callocN(sizeof(BVHTree), "BVHTree");
271         
272         if(tree)
273         {
274                 // calculate max number of branches, our bvh kdop is "almost perfect"
275                 for(i = 1; i <= (int)ceil((float)((float)log(maxsize)/(float)log(tree_type))); i++)
276                         numbranches += (pow(tree_type, i) / tree_type);
277                 
278                 tree->nodes = (BVHNode **)MEM_callocN(sizeof(BVHNode *)*(numbranches+maxsize + tree_type), "BVHNodes");
279                 
280                 if(!tree->nodes)
281                 {
282                         MEM_freeN(tree);
283                         return NULL;
284                 }
285                 
286                 tree->nodearray = (BVHNode *)MEM_callocN(sizeof(BVHNode)*(numbranches+maxsize + tree_type), "BVHNodeArray");
287                 
288                 if(!tree->nodearray)
289                 {
290                         MEM_freeN(tree);
291                         MEM_freeN(tree->nodes);
292                         return NULL;
293                 }
294                 
295                 tree->orig_index = (int *)MEM_callocN(sizeof(int)*(numbranches+maxsize + tree_type), "BVHIndexArray");
296                 
297                 if(!tree->orig_index)
298                 {
299                         MEM_freeN(tree);
300                         MEM_freeN(tree->nodes);
301                         MEM_freeN(tree->nodearray);
302                         return NULL;
303                 }
304                 
305                 tree->epsilon = epsilon;
306                 tree->tree_type = tree_type; 
307                 tree->axis = axis;
308                 
309                 if(axis == 26)
310                 {
311                         tree->start_axis = 0;
312                         tree->stop_axis = 13;
313                 }
314                 else if(axis == 18)
315                 {
316                         tree->start_axis = 7;
317                         tree->stop_axis = 13;
318                 }
319                 else if(axis == 14)
320                 {
321                         tree->start_axis = 0;
322                         tree->stop_axis = 7;
323                 }
324                 else if(axis == 8) // AABB
325                 {
326                         tree->start_axis = 0;
327                         tree->stop_axis = 4;
328                 }
329                 else if(axis == 6) // OBB
330                 {
331                         tree->start_axis = 0;
332                         tree->stop_axis = 3;
333                 }
334                 else
335                 {
336                         BLI_bvhtree_free(tree);
337                         return NULL;
338                 }
339         }
340
341         return tree;
342 }
343
344 static void create_kdop_hull(BVHTree *tree, BVHNode *node, float *co, int numpoints, int moving)
345 {
346         float newminmax;
347         int i, k;
348         
349         // don't init boudings for the moving case
350         if(!moving)
351         {
352                 for (i = tree->start_axis; i < tree->stop_axis; i++)
353                 {
354                         node->bv[2*i] = FLT_MAX;
355                         node->bv[2*i + 1] = -FLT_MAX;
356                 }
357         }
358         
359         for(k = 0; k < numpoints; k++)
360         {
361                 // for all Axes.
362                 for (i = tree->start_axis; i < tree->stop_axis; i++)
363                 {
364                         newminmax = INPR(&co[k * 3], KDOP_AXES[i]);
365                         if (newminmax < node->bv[2 * i])
366                                 node->bv[2 * i] = newminmax;
367                         if (newminmax > node->bv[(2 * i) + 1])
368                                 node->bv[(2 * i) + 1] = newminmax;
369                 }
370         }
371 }
372
373 // depends on the fact that the BVH's for each face is already build
374 static void refit_kdop_hull(BVHTree *tree, int start, int end, float *bv)
375 {
376         float newmin,newmax;
377         int i, j;
378         
379         for (i = tree->start_axis; i < tree->stop_axis; i++)
380         {
381                 bv[2*i] = FLT_MAX;
382                 bv[2*i + 1] = -FLT_MAX;
383         }
384
385         for (j = start; j < end; j++)
386         {
387                 // for all Axes.
388                 for (i = tree->start_axis; i < tree->stop_axis; i++)
389                 {
390                         newmin = tree->nodes[j]->bv[(2 * i)];   
391                         if ((newmin < bv[(2 * i)]))
392                                 bv[(2 * i)] = newmin;
393  
394                         newmax = tree->nodes[j]->bv[(2 * i) + 1];
395                         if ((newmax > bv[(2 * i) + 1]))
396                                 bv[(2 * i) + 1] = newmax;
397                 }
398         }
399 }
400
401 int BLI_bvhtree_insert(BVHTree *tree, int index, float *co, int numpoints)
402 {
403         BVHNode *node= NULL;
404         int i;
405         
406         // insert should only possible as long as tree->totbranch is 0
407         if(tree->totbranch > 0)
408                 return 0;
409         
410         if(tree->totleaf+1 >= MEM_allocN_len(tree->nodes))
411                 return 0;
412         
413         // TODO check if have enough nodes in array
414         
415         node = tree->nodes[tree->totleaf] = &(tree->nodearray[tree->totleaf]);
416         tree->totleaf++;
417         
418         create_kdop_hull(tree, node, co, numpoints, 0);
419         
420         // inflate the bv with some epsilon
421         for (i = tree->start_axis; i < tree->stop_axis; i++)
422         {
423                 node->bv[(2 * i)] -= tree->epsilon; // minimum 
424                 node->bv[(2 * i) + 1] += tree->epsilon; // maximum 
425         }
426
427         node->index= index;
428         
429         return 1;
430 }
431
432 // only supports x,y,z axis in the moment
433 // but we should use a plain and simple function here for speed sake
434 static char get_largest_axis(float *bv)
435 {
436         float middle_point[3];
437
438         middle_point[0] = (bv[1]) - (bv[0]); // x axis
439         middle_point[1] = (bv[3]) - (bv[2]); // y axis
440         middle_point[2] = (bv[5]) - (bv[4]); // z axis
441         if (middle_point[0] > middle_point[1]) 
442         {
443                 if (middle_point[0] > middle_point[2])
444                         return 1; // max x axis
445                 else
446                         return 5; // max z axis
447         }
448         else 
449         {
450                 if (middle_point[1] > middle_point[2])
451                         return 3; // max y axis
452                 else
453                         return 5; // max z axis
454         }
455 }
456
457 static void bvh_div_nodes(BVHTree *tree, BVHNode *node, int start, int end, char lastaxis)
458 {
459         char laxis;
460         int i, tend;
461         BVHNode *tnode;
462         int slice = (end-start+tree->tree_type-1)/tree->tree_type;      //division rounded up
463         
464         // Determine which axis to split along
465         laxis = get_largest_axis(node->bv);
466
467         // Sort along longest axis
468         if(laxis!=lastaxis)
469                 sort_along_axis(tree, start, end, laxis);
470         
471         // split nodes along longest axis
472         for (i=0; start < end; start += slice, i++) //i counts the current child
473         {       
474                 tend = start + slice;
475                 
476                 if(tend > end) tend = end;
477                 
478                 if(tend-start == 1)     // ok, we have 1 left for this node
479                 {
480                         node->children[i] = tree->nodes[start];
481                         node->children[i]->parent = node;
482                 }
483                 else
484                 {
485                         tnode = node->children[i] = tree->nodes[tree->totleaf  + tree->totbranch] = &(tree->nodearray[tree->totbranch + tree->totleaf]);
486                         tree->totbranch++;
487                         tnode->parent = node;
488
489                         refit_kdop_hull(tree, start, tend, tnode->bv);
490                         bvh_div_nodes(tree, tnode, start, tend, laxis);
491                 }
492                 node->totnode++;
493         }
494         
495         return;
496 }
497
498 static void verify_tree(BVHTree *tree)
499 {
500         int i, j, check = 0;
501         
502         // check the pointer list
503         for(i = 0; i < tree->totleaf; i++)
504         {
505                 if(tree->nodes[i]->parent == NULL)
506                         printf("Leaf has no parent: %d\n", i);
507                 else
508                 {
509                         for(j = 0; j < tree->tree_type; j++)
510                         {
511                                 if(tree->nodes[i]->parent->children[j] == tree->nodes[i])
512                                         check = 1;
513                         }
514                         if(!check)
515                         {
516                                 printf("Parent child relationship doesn't match: %d\n", i);
517                         }
518                         check = 0;
519                 }
520         }
521         
522         // check the leaf list
523         for(i = 0; i < tree->totleaf; i++)
524         {
525                 if(tree->nodearray[i].parent == NULL)
526                         printf("Leaf has no parent: %d\n", i);
527                 else
528                 {
529                         for(j = 0; j < tree->tree_type; j++)
530                         {
531                                 if(tree->nodearray[i].parent->children[j] == &tree->nodearray[i])
532                                         check = 1;
533                         }
534                         if(!check)
535                         {
536                                 printf("Parent child relationship doesn't match: %d\n", i);
537                         }
538                         check = 0;
539                 }
540         }
541         
542         printf("branches: %d, leafs: %d, total: %d\n", tree->totbranch, tree->totleaf, tree->totbranch + tree->totleaf);
543 }
544         
545 void BLI_bvhtree_balance(BVHTree *tree)
546 {
547         BVHNode *node;
548         int i;
549         
550         if(tree->totleaf == 0)
551                 return;
552         
553         // create root node
554         node = tree->nodes[tree->totleaf] = &(tree->nodearray[tree->totleaf]);
555         tree->totbranch++;
556         
557         // refit root bvh node
558         refit_kdop_hull(tree, 0, tree->totleaf, tree->nodes[tree->totleaf]->bv);
559         // create + balance tree
560         bvh_div_nodes(tree, tree->nodes[tree->totleaf], 0, tree->totleaf, 0);
561         
562         // put indices into array for O(1) access
563         for(i = 0; i < tree->totleaf; i++)
564         {
565                 tree->orig_index[tree->nodes[i]->index] = i;
566         }
567         
568         verify_tree(tree);
569 }
570
571 // overlap - is it possbile for 2 bv's to collide ?
572 static int tree_overlap(float *bv1, float *bv2, int start_axis, int stop_axis)
573 {
574         float *bv1_end = bv1 + (stop_axis<<1);
575                 
576         bv1 += start_axis<<1;
577         bv2 += start_axis<<1;
578         
579         // test all axis if min + max overlap
580         for (; bv1 != bv1_end; bv1+=2, bv2+=2)
581         {
582                 if ((*(bv1) > *(bv2 + 1)) || (*(bv2) > *(bv1 + 1))) 
583                         return 0;
584         }
585         
586         return 1;
587 }
588
589 static void traverse(BVHOverlapData *data, BVHNode *node1, BVHNode *node2)
590 {
591         int j;
592         
593         if(tree_overlap(node1->bv, node2->bv, MIN2(data->tree1->start_axis, data->tree2->start_axis), MIN2(data->tree1->stop_axis, data->tree2->stop_axis)))
594         {
595                 // check if node1 is a leaf
596                 if(node1->index)
597                 {
598                         // check if node2 is a leaf
599                         if(node2->index)
600                         {
601                                 
602                                 if(node1 == node2)
603                                 {
604                                         return;
605                                 }
606                                         
607                                 if(data->i >= data->max_overlap)
608                                 {       
609                                         // try to make alloc'ed memory bigger
610                                         data->overlap = realloc(data->overlap, sizeof(BVHTreeOverlap)*data->max_overlap*2);
611                                         
612                                         if(!data->overlap)
613                                         {
614                                                 printf("Out of Memory in traverse\n");
615                                                 return;
616                                         }
617                                         data->max_overlap *= 2;
618                                 }
619                                 
620                                 // both leafs, insert overlap!
621                                 data->overlap[data->i].indexA = node1->index;
622                                 data->overlap[data->i].indexB = node2->index;
623
624                                 data->i++;
625                         }
626                         else
627                         {
628                                 for(j = 0; j < data->tree2->tree_type; j++)
629                                 {
630                                         if(node2->children[j])
631                                                 traverse(data, node1, node2->children[j]);
632                                 }
633                         }
634                 }
635                 else
636                 {
637                         
638                         for(j = 0; j < data->tree2->tree_type; j++)
639                         {
640                                 if(node1->children[j])
641                                         traverse(data, node1->children[j], node2);
642                         }
643                 }
644         }
645         return;
646 }
647
648 BVHTreeOverlap *BLI_bvhtree_overlap(BVHTree *tree1, BVHTree *tree2, int *result)
649 {
650         int j, total = 0;
651         BVHTreeOverlap *overlap = NULL, *to = NULL;
652         BVHOverlapData *data[tree1->tree_type];
653         
654         // check for compatibility of both trees (can't compare 14-DOP with 18-DOP)
655         if((tree1->axis != tree2->axis) && ((tree1->axis == 14) || tree2->axis == 14))
656                 return 0;
657         
658         // fast check root nodes for collision before doing big splitting + traversal
659         if(!tree_overlap(tree1->nodes[tree1->totleaf]->bv, tree2->nodes[tree2->totleaf]->bv, MIN2(tree1->start_axis, tree2->start_axis), MIN2(tree1->stop_axis, tree2->stop_axis)))
660                 return 0;
661         
662         for(j = 0; j < tree1->tree_type; j++)
663         {
664                 data[j] = (BVHOverlapData *)MEM_callocN(sizeof(BVHOverlapData), "BVHOverlapData");
665                 
666                 // init BVHOverlapData
667                 data[j]->overlap = (BVHTreeOverlap *)malloc(sizeof(BVHTreeOverlap)*MAX2(tree1->totleaf, tree2->totleaf));
668                 data[j]->tree1 = tree1;
669                 data[j]->tree2 = tree2;
670                 data[j]->max_overlap = MAX2(tree1->totleaf, tree2->totleaf);
671                 data[j]->i = 0;
672         }
673
674 #pragma omp parallel for private(j) schedule(static)
675         for(j = 0; j < tree1->tree_type; j++)
676         {
677                 traverse(data[j], tree1->nodes[tree1->totleaf]->children[j], tree2->nodes[tree2->totleaf]);
678         }
679         
680         for(j = 0; j < tree1->tree_type; j++)
681                 total += data[j]->i;
682         
683         to = overlap = (BVHTreeOverlap *)MEM_callocN(sizeof(BVHTreeOverlap)*total, "BVHTreeOverlap");
684         
685         for(j = 0; j < tree1->tree_type; j++)
686         {
687                 memcpy(to, data[j]->overlap, data[j]->i*sizeof(BVHTreeOverlap));
688                 to+=data[j]->i;
689         }
690         
691         for(j = 0; j < tree1->tree_type; j++)
692         {
693                 free(data[j]->overlap);
694                 MEM_freeN(data[j]);
695         }
696         
697         (*result) = total;
698         return overlap;
699 }
700
701
702 // bottom up update of bvh tree:
703 // join the 4 children here
704 static void node_join(BVHTree *tree, BVHNode *node)
705 {
706         int i, j;
707         
708         for (i = tree->start_axis; i < tree->stop_axis; i++)
709         {
710                 node->bv[2*i] = FLT_MAX;
711                 node->bv[2*i + 1] = -FLT_MAX;
712         }
713         
714         for (i = 0; i < tree->tree_type; i++)
715         {
716                 if (node->children[i]) 
717                 {
718                         for (j = tree->start_axis; j < tree->stop_axis; j++)
719                         {
720                                 // update minimum 
721                                 if (node->children[i]->bv[(2 * j)] < node->bv[(2 * j)]) 
722                                         node->bv[(2 * j)] = node->children[i]->bv[(2 * j)];
723                                 
724                                 // update maximum 
725                                 if (node->children[i]->bv[(2 * j) + 1] > node->bv[(2 * j) + 1])
726                                         node->bv[(2 * j) + 1] = node->children[i]->bv[(2 * j) + 1];
727                         }
728                 }
729                 else
730                         break;
731         }
732 }
733
734 // call before BLI_bvhtree_update_tree()
735 int BLI_bvhtree_update_node(BVHTree *tree, int index, float *co, float *co_moving, int numpoints)
736 {
737         BVHNode *node= NULL;
738         int i = 0;
739         
740         // check if index exists
741         if(index > tree->totleaf)
742                 return 0;
743         
744         node = tree->nodes[tree->orig_index[index]];
745         
746         create_kdop_hull(tree, node, co, numpoints, 0);
747         
748         if(co_moving)
749                 create_kdop_hull(tree, node, co_moving, numpoints, 1);
750         
751         // inflate the bv with some epsilon
752         for (i = tree->start_axis; i < tree->stop_axis; i++)
753         {
754                 node->bv[(2 * i)] -= tree->epsilon; // minimum 
755                 node->bv[(2 * i) + 1] += tree->epsilon; // maximum 
756         }
757         
758         return 1;
759 }
760
761 // call BLI_bvhtree_update_node() first for every node/point/triangle
762 void BLI_bvhtree_update_tree(BVHTree *tree)
763 {
764         BVHNode *leaf, *parent;
765         
766         // reset tree traversing flag
767         for (leaf = tree->nodearray + tree->totleaf; leaf != tree->nodearray + tree->totleaf + tree->totbranch; leaf++)
768                 leaf->traversed = 0;
769         
770         for (leaf = tree->nodearray; leaf != tree->nodearray + tree->totleaf; leaf++)
771         {
772                 for (parent = leaf->parent; parent; parent = parent->parent)
773                 {
774                         parent->traversed++;    // we tried to go up in hierarchy 
775                         if (parent->traversed < parent->totnode) 
776                                 break;  // we do not need to check further 
777                         else 
778                                 node_join(tree, parent);
779                 }
780         }
781 }
782
783 float BLI_bvhtree_getepsilon(BVHTree *tree)
784 {
785         return tree->epsilon;
786 }