Blender V2.61 - r43446

pbvh.c

Go to the documentation of this file.
00001 /*
00002  * ***** BEGIN GPL LICENSE BLOCK *****
00003  *
00004  * This program is free software; you can redistribute it and/or
00005  * modify it under the terms of the GNU General Public License
00006  * as published by the Free Software Foundation; either version 2
00007  * of the License, or (at your option) any later version.
00008  *
00009  * This program is distributed in the hope that it will be useful,
00010  * but WITHOUT ANY WARRANTY; without even the implied warranty of
00011  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00012  * GNU General Public License for more details.
00013  *
00014  * You should have received a copy of the GNU General Public License
00015  * along with this program; if not, write to the Free Software Foundation,
00016  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
00017  *
00018  * ***** END GPL LICENSE BLOCK *****
00019  */
00020 
00027 #include "DNA_meshdata_types.h"
00028 
00029 #include "MEM_guardedalloc.h"
00030 
00031 #include "BLI_math.h"
00032 #include "BLI_utildefines.h"
00033 #include "BLI_ghash.h"
00034 #include "BLI_pbvh.h"
00035 
00036 #include "BKE_DerivedMesh.h"
00037 #include "BKE_mesh.h" /* for mesh_calc_normals */
00038 #include "BKE_global.h" /* for mesh_calc_normals */
00039 
00040 #include "GPU_buffers.h"
00041 
00042 #define LEAF_LIMIT 10000
00043 
00044 //#define PERFCNTRS
00045 
00046 /* Bitmap */
00047 typedef char* BLI_bitmap;
00048 
00049 static BLI_bitmap BLI_bitmap_new(int tot)
00050 {
00051     return MEM_callocN((tot >> 3) + 1, "BLI bitmap");
00052 }
00053 
00054 static int BLI_bitmap_get(BLI_bitmap b, int index)
00055 {
00056     return b[index >> 3] & (1 << (index & 7));
00057 }
00058 
00059 static void BLI_bitmap_set(BLI_bitmap b, int index)
00060 {
00061     b[index >> 3] |= (1 << (index & 7));
00062 }
00063 
00064 #if 0 /* UNUSED */
00065 static void BLI_bitmap_clear(BLI_bitmap b, int index)
00066 {
00067     b[index >> 3] &= ~(1 << (index & 7));
00068 }
00069 #endif
00070 
00071 /* Axis-aligned bounding box */
00072 typedef struct {
00073     float bmin[3], bmax[3];
00074 } BB;
00075 
00076 /* Axis-aligned bounding box with centroid */
00077 typedef struct {
00078     float bmin[3], bmax[3], bcentroid[3];
00079 } BBC;
00080 
00081 struct PBVHNode {
00082     /* Opaque handle for drawing code */
00083     GPU_Buffers *draw_buffers;
00084 
00085     /* Voxel bounds */
00086     BB vb;
00087     BB orig_vb;
00088 
00089     /* For internal nodes, the offset of the children in the PBVH
00090        'nodes' array. */
00091     int children_offset;
00092 
00093     /* Pointer into the PBVH prim_indices array and the number of
00094        primitives used by this leaf node.
00095 
00096        Used for leaf nodes in both mesh- and multires-based PBVHs.
00097     */
00098     int *prim_indices;
00099     unsigned int totprim;
00100 
00101     /* Array of indices into the mesh's MVert array. Contains the
00102        indices of all vertices used by faces that are within this
00103        node's bounding box.
00104 
00105        Note that a vertex might be used by a multiple faces, and
00106        these faces might be in different leaf nodes. Such a vertex
00107        will appear in the vert_indices array of each of those leaf
00108        nodes.
00109 
00110        In order to support cases where you want access to multiple
00111        nodes' vertices without duplication, the vert_indices array
00112        is ordered such that the first part of the array, up to
00113        index 'uniq_verts', contains "unique" vertex indices. These
00114        vertices might not be truly unique to this node, but if
00115        they appear in another node's vert_indices array, they will
00116        be above that node's 'uniq_verts' value.
00117 
00118        Used for leaf nodes in a mesh-based PBVH (not multires.)
00119     */
00120     int *vert_indices;
00121     unsigned int uniq_verts, face_verts;
00122 
00123     /* An array mapping face corners into the vert_indices
00124        array. The array is sized to match 'totprim', and each of
00125        the face's corners gets an index into the vert_indices
00126        array, in the same order as the corners in the original
00127        MFace. The fourth value should not be used if the original
00128        face is a triangle.
00129 
00130        Used for leaf nodes in a mesh-based PBVH (not multires.)
00131     */
00132     int (*face_vert_indices)[4];
00133 
00134     /* Indicates whether this node is a leaf or not; also used for
00135        marking various updates that need to be applied. */
00136     PBVHNodeFlags flag : 8;
00137 
00138     /* Used for raycasting: how close bb is to the ray point. */
00139     float tmin;
00140 
00141     int proxy_count;
00142     PBVHProxyNode* proxies;
00143 };
00144 
00145 struct PBVH {
00146     PBVHNode *nodes;
00147     int node_mem_count, totnode;
00148 
00149     int *prim_indices;
00150     int totprim;
00151     int totvert;
00152 
00153     int leaf_limit;
00154 
00155     /* Mesh data */
00156     MVert *verts;
00157     MFace *faces;
00158 
00159     /* Grid Data */
00160     DMGridData **grids;
00161     DMGridAdjacency *gridadj;
00162     void **gridfaces;
00163     int totgrid;
00164     int gridsize;
00165 
00166     /* Only used during BVH build and update,
00167        don't need to remain valid after */
00168     BLI_bitmap vert_bitmap;
00169 
00170 #ifdef PERFCNTRS
00171     int perf_modified;
00172 #endif
00173 
00174     /* flag are verts/faces deformed */
00175     int deformed;
00176 };
00177 
00178 #define STACK_FIXED_DEPTH   100
00179 
00180 typedef struct PBVHStack {
00181     PBVHNode *node;
00182     int revisiting;
00183 } PBVHStack;
00184 
00185 typedef struct PBVHIter {
00186     PBVH *bvh;
00187     BLI_pbvh_SearchCallback scb;
00188     void *search_data;
00189 
00190     PBVHStack *stack;
00191     int stacksize;
00192 
00193     PBVHStack stackfixed[STACK_FIXED_DEPTH];
00194     int stackspace;
00195 } PBVHIter;
00196 
00197 static void BB_reset(BB *bb)
00198 {
00199     bb->bmin[0] = bb->bmin[1] = bb->bmin[2] = FLT_MAX;
00200     bb->bmax[0] = bb->bmax[1] = bb->bmax[2] = -FLT_MAX;
00201 }
00202 
00203 /* Expand the bounding box to include a new coordinate */
00204 static void BB_expand(BB *bb, float co[3])
00205 {
00206     int i;
00207     for(i = 0; i < 3; ++i) {
00208         bb->bmin[i] = MIN2(bb->bmin[i], co[i]);
00209         bb->bmax[i] = MAX2(bb->bmax[i], co[i]);
00210     }
00211 }
00212 
00213 /* Expand the bounding box to include another bounding box */
00214 static void BB_expand_with_bb(BB *bb, BB *bb2)
00215 {
00216     int i;
00217     for(i = 0; i < 3; ++i) {
00218         bb->bmin[i] = MIN2(bb->bmin[i], bb2->bmin[i]);
00219         bb->bmax[i] = MAX2(bb->bmax[i], bb2->bmax[i]);
00220     }
00221 }
00222 
00223 /* Return 0, 1, or 2 to indicate the widest axis of the bounding box */
00224 static int BB_widest_axis(BB *bb)
00225 {
00226     float dim[3];
00227     int i;
00228 
00229     for(i = 0; i < 3; ++i)
00230         dim[i] = bb->bmax[i] - bb->bmin[i];
00231 
00232     if(dim[0] > dim[1]) {
00233         if(dim[0] > dim[2])
00234             return 0;
00235         else
00236             return 2;
00237     }
00238     else {
00239         if(dim[1] > dim[2])
00240             return 1;
00241         else
00242             return 2;
00243     }
00244 }
00245 
00246 static void BBC_update_centroid(BBC *bbc)
00247 {
00248     int i;
00249     for(i = 0; i < 3; ++i)
00250         bbc->bcentroid[i] = (bbc->bmin[i] + bbc->bmax[i]) * 0.5f;
00251 }
00252 
00253 /* Not recursive */
00254 static void update_node_vb(PBVH *bvh, PBVHNode *node)
00255 {
00256     BB vb;
00257 
00258     BB_reset(&vb);
00259     
00260     if(node->flag & PBVH_Leaf) {
00261         PBVHVertexIter vd;
00262 
00263         BLI_pbvh_vertex_iter_begin(bvh, node, vd, PBVH_ITER_ALL) {
00264             BB_expand(&vb, vd.co);
00265         }
00266         BLI_pbvh_vertex_iter_end;
00267     }
00268     else {
00269         BB_expand_with_bb(&vb,
00270                   &bvh->nodes[node->children_offset].vb);
00271         BB_expand_with_bb(&vb,
00272                   &bvh->nodes[node->children_offset + 1].vb);
00273     }
00274 
00275     node->vb= vb;
00276 }
00277 
00278 //void BLI_pbvh_node_BB_reset(PBVHNode* node)
00279 //{
00280 //  BB_reset(&node->vb);
00281 //}
00282 //
00283 //void BLI_pbvh_node_BB_expand(PBVHNode* node, float co[3])
00284 //{
00285 //  BB_expand(&node->vb, co);
00286 //}
00287 
00288 
00289 /* Adapted from BLI_kdopbvh.c */
00290 /* Returns the index of the first element on the right of the partition */
00291 static int partition_indices(int *prim_indices, int lo, int hi, int axis,
00292                  float mid, BBC *prim_bbc)
00293 {
00294     int i=lo, j=hi;
00295     for(;;) {
00296         for(; prim_bbc[prim_indices[i]].bcentroid[axis] < mid; i++);
00297         for(; mid < prim_bbc[prim_indices[j]].bcentroid[axis]; j--);
00298         
00299         if(!(i < j))
00300             return i;
00301         
00302         SWAP(int, prim_indices[i], prim_indices[j]);
00303         i++;
00304     }
00305 }
00306 
00307 static void check_partitioning(int *prim_indices, int lo, int hi, int axis,
00308                    float mid, BBC *prim_bbc, int index_of_2nd_partition)
00309 {
00310     int i;
00311     for(i = lo; i <= hi; ++i) {
00312         const float c = prim_bbc[prim_indices[i]].bcentroid[axis];
00313 
00314         if((i < index_of_2nd_partition && c > mid) ||
00315            (i > index_of_2nd_partition && c < mid)) {
00316             printf("fail\n");
00317         }
00318     }
00319 }
00320 
00321 static void grow_nodes(PBVH *bvh, int totnode)
00322 {
00323     if(totnode > bvh->node_mem_count) {
00324         PBVHNode *prev = bvh->nodes;
00325         bvh->node_mem_count *= 1.33;
00326         if(bvh->node_mem_count < totnode)
00327             bvh->node_mem_count = totnode;
00328         bvh->nodes = MEM_callocN(sizeof(PBVHNode) * bvh->node_mem_count,
00329                      "bvh nodes");
00330         memcpy(bvh->nodes, prev, bvh->totnode * sizeof(PBVHNode));
00331         MEM_freeN(prev);
00332     }
00333 
00334     bvh->totnode = totnode;
00335 }
00336 
00337 /* Add a vertex to the map, with a positive value for unique vertices and
00338    a negative value for additional vertices */
00339 static int map_insert_vert(PBVH *bvh, GHash *map,
00340                 unsigned int *face_verts,
00341                 unsigned int *uniq_verts, int vertex)
00342 {
00343     void *value, *key = SET_INT_IN_POINTER(vertex);
00344 
00345     if(!BLI_ghash_haskey(map, key)) {
00346         if(BLI_bitmap_get(bvh->vert_bitmap, vertex)) {
00347             value = SET_INT_IN_POINTER(~(*face_verts));
00348             ++(*face_verts);
00349         }
00350         else {
00351             BLI_bitmap_set(bvh->vert_bitmap, vertex);
00352             value = SET_INT_IN_POINTER(*uniq_verts);
00353             ++(*uniq_verts);
00354         }
00355         
00356         BLI_ghash_insert(map, key, value);
00357         return GET_INT_FROM_POINTER(value);
00358     }
00359     else
00360         return GET_INT_FROM_POINTER(BLI_ghash_lookup(map, key));
00361 }
00362 
00363 /* Find vertices used by the faces in this node and update the draw buffers */
00364 static void build_mesh_leaf_node(PBVH *bvh, PBVHNode *node)
00365 {
00366     GHashIterator *iter;
00367     GHash *map;
00368     int i, j, totface;
00369 
00370     map = BLI_ghash_new(BLI_ghashutil_inthash, BLI_ghashutil_intcmp, "build_mesh_leaf_node gh");
00371     
00372     node->uniq_verts = node->face_verts = 0;
00373     totface= node->totprim;
00374 
00375     node->face_vert_indices = MEM_callocN(sizeof(int) * 4*totface,
00376                           "bvh node face vert indices");
00377 
00378     for(i = 0; i < totface; ++i) {
00379         MFace *f = bvh->faces + node->prim_indices[i];
00380         int sides = f->v4 ? 4 : 3;
00381 
00382         for(j = 0; j < sides; ++j) {
00383             node->face_vert_indices[i][j]= 
00384                 map_insert_vert(bvh, map, &node->face_verts,
00385                         &node->uniq_verts, (&f->v1)[j]);
00386         }
00387     }
00388 
00389     node->vert_indices = MEM_callocN(sizeof(int) *
00390                      (node->uniq_verts + node->face_verts),
00391                      "bvh node vert indices");
00392 
00393     /* Build the vertex list, unique verts first */
00394     for(iter = BLI_ghashIterator_new(map), i = 0;
00395         !BLI_ghashIterator_isDone(iter);
00396         BLI_ghashIterator_step(iter), ++i) {
00397         void *value = BLI_ghashIterator_getValue(iter);
00398         int ndx = GET_INT_FROM_POINTER(value);
00399 
00400         if(ndx < 0)
00401             ndx = -ndx + node->uniq_verts - 1;
00402 
00403         node->vert_indices[ndx] =
00404             GET_INT_FROM_POINTER(BLI_ghashIterator_getKey(iter));
00405     }
00406 
00407     BLI_ghashIterator_free(iter);
00408 
00409     for(i = 0; i < totface; ++i) {
00410         MFace *f = bvh->faces + node->prim_indices[i];
00411         int sides = f->v4 ? 4 : 3;
00412         
00413         for(j = 0; j < sides; ++j) {
00414             if(node->face_vert_indices[i][j] < 0)
00415                 node->face_vert_indices[i][j]=
00416                     -node->face_vert_indices[i][j] +
00417                     node->uniq_verts - 1;
00418         }
00419     }
00420 
00421     if(!G.background) {
00422         node->draw_buffers =
00423             GPU_build_mesh_buffers(map, bvh->verts, bvh->faces,
00424                     node->prim_indices,
00425                     node->totprim, node->vert_indices,
00426                     node->uniq_verts,
00427                     node->uniq_verts + node->face_verts);
00428     }
00429 
00430     node->flag |= PBVH_UpdateDrawBuffers;
00431 
00432     BLI_ghash_free(map, NULL, NULL);
00433 }
00434 
00435 static void build_grids_leaf_node(PBVH *bvh, PBVHNode *node)
00436 {
00437     if(!G.background) {
00438         node->draw_buffers =
00439             GPU_build_grid_buffers(bvh->grids, node->prim_indices,
00440                 node->totprim, bvh->gridsize);
00441     }
00442     node->flag |= PBVH_UpdateDrawBuffers;
00443 }
00444 
00445 /* Recursively build a node in the tree
00446 
00447    vb is the voxel box around all of the primitives contained in
00448    this node.
00449 
00450    cb is the bounding box around all the centroids of the primitives
00451    contained in this node
00452 
00453    offset and start indicate a range in the array of primitive indices
00454 */
00455 
00456 static void build_sub(PBVH *bvh, int node_index, BB *cb, BBC *prim_bbc,
00457            int offset, int count)
00458 {
00459     int i, axis, end;
00460     BB cb_backing;
00461 
00462     /* Decide whether this is a leaf or not */
00463     // XXX adapt leaf limit for grids
00464     if(count <= bvh->leaf_limit) {
00465         bvh->nodes[node_index].flag |= PBVH_Leaf;
00466 
00467         bvh->nodes[node_index].prim_indices = bvh->prim_indices + offset;
00468         bvh->nodes[node_index].totprim = count;
00469 
00470         /* Still need vb for searches */
00471         BB_reset(&bvh->nodes[node_index].vb);
00472         for(i = offset + count - 1; i >= offset; --i) {
00473             BB_expand_with_bb(&bvh->nodes[node_index].vb,
00474                       (BB*)(prim_bbc +
00475                         bvh->prim_indices[i]));
00476         }
00477         
00478         if(bvh->faces)
00479             build_mesh_leaf_node(bvh, bvh->nodes + node_index);
00480         else
00481             build_grids_leaf_node(bvh, bvh->nodes + node_index);
00482         bvh->nodes[node_index].orig_vb= bvh->nodes[node_index].vb;
00483 
00484         /* Done with this subtree */
00485         return;
00486     }
00487     else {
00488         BB_reset(&bvh->nodes[node_index].vb);
00489         bvh->nodes[node_index].children_offset = bvh->totnode;
00490         grow_nodes(bvh, bvh->totnode + 2);
00491 
00492         if(!cb) {
00493             cb = &cb_backing;
00494             BB_reset(cb);
00495             for(i = offset + count - 1; i >= offset; --i)
00496                 BB_expand(cb, prim_bbc[bvh->prim_indices[i]].bcentroid);
00497         }
00498     }
00499 
00500     axis = BB_widest_axis(cb);
00501 
00502     for(i = offset + count - 1; i >= offset; --i) {
00503         BB_expand_with_bb(&bvh->nodes[node_index].vb,
00504                   (BB*)(prim_bbc + bvh->prim_indices[i]));
00505     }
00506 
00507     bvh->nodes[node_index].orig_vb= bvh->nodes[node_index].vb;
00508 
00509     end = partition_indices(bvh->prim_indices, offset, offset + count - 1,
00510                             axis,
00511                             (cb->bmax[axis] + cb->bmin[axis]) * 0.5f,
00512                             prim_bbc);
00513     check_partitioning(bvh->prim_indices, offset, offset + count - 1,
00514                        axis,
00515                        (cb->bmax[axis] + cb->bmin[axis]) * 0.5f,
00516                        prim_bbc, end);
00517 
00518     build_sub(bvh, bvh->nodes[node_index].children_offset, NULL,
00519           prim_bbc, offset, end - offset);
00520     build_sub(bvh, bvh->nodes[node_index].children_offset + 1, NULL,
00521           prim_bbc, end, offset + count - end);
00522 }
00523 
00524 static void pbvh_build(PBVH *bvh, BB *cb, BBC *prim_bbc, int totprim)
00525 {
00526     int i;
00527 
00528     if(totprim != bvh->totprim) {
00529         bvh->totprim = totprim;
00530         if(bvh->nodes) MEM_freeN(bvh->nodes);
00531         if(bvh->prim_indices) MEM_freeN(bvh->prim_indices);
00532         bvh->prim_indices = MEM_callocN(sizeof(int) * totprim,
00533                         "bvh prim indices");
00534         for(i = 0; i < totprim; ++i)
00535             bvh->prim_indices[i] = i;
00536         bvh->totnode = 0;
00537         if(bvh->node_mem_count < 100) {
00538             bvh->node_mem_count = 100;
00539             bvh->nodes = MEM_callocN(sizeof(PBVHNode) *
00540                          bvh->node_mem_count,
00541                          "bvh initial nodes");
00542         }
00543     }
00544 
00545     bvh->totnode = 1;
00546     build_sub(bvh, 0, cb, prim_bbc, 0, totprim);
00547 }
00548 
00549 /* Do a full rebuild with on Mesh data structure */
00550 void BLI_pbvh_build_mesh(PBVH *bvh, MFace *faces, MVert *verts, int totface, int totvert)
00551 {
00552     BBC *prim_bbc = NULL;
00553     BB cb;
00554     int i, j;
00555 
00556     bvh->faces = faces;
00557     bvh->verts = verts;
00558     bvh->vert_bitmap = BLI_bitmap_new(totvert);
00559     bvh->totvert = totvert;
00560     bvh->leaf_limit = LEAF_LIMIT;
00561 
00562     BB_reset(&cb);
00563 
00564     /* For each face, store the AABB and the AABB centroid */
00565     prim_bbc = MEM_mallocN(sizeof(BBC) * totface, "prim_bbc");
00566 
00567     for(i = 0; i < totface; ++i) {
00568         MFace *f = faces + i;
00569         const int sides = f->v4 ? 4 : 3;
00570         BBC *bbc = prim_bbc + i;
00571 
00572         BB_reset((BB*)bbc);
00573 
00574         for(j = 0; j < sides; ++j)
00575             BB_expand((BB*)bbc, verts[(&f->v1)[j]].co);
00576 
00577         BBC_update_centroid(bbc);
00578 
00579         BB_expand(&cb, bbc->bcentroid);
00580     }
00581 
00582     if(totface)
00583         pbvh_build(bvh, &cb, prim_bbc, totface);
00584 
00585     MEM_freeN(prim_bbc);
00586     MEM_freeN(bvh->vert_bitmap);
00587 }
00588 
00589 /* Do a full rebuild with on Grids data structure */
00590 void BLI_pbvh_build_grids(PBVH *bvh, DMGridData **grids, DMGridAdjacency *gridadj,
00591     int totgrid, int gridsize, void **gridfaces)
00592 {
00593     BBC *prim_bbc = NULL;
00594     BB cb;
00595     int i, j;
00596 
00597     bvh->grids= grids;
00598     bvh->gridadj= gridadj;
00599     bvh->gridfaces= gridfaces;
00600     bvh->totgrid= totgrid;
00601     bvh->gridsize= gridsize;
00602     bvh->leaf_limit = MAX2(LEAF_LIMIT/((gridsize-1)*(gridsize-1)), 1);
00603 
00604     BB_reset(&cb);
00605 
00606     /* For each grid, store the AABB and the AABB centroid */
00607     prim_bbc = MEM_mallocN(sizeof(BBC) * totgrid, "prim_bbc");
00608 
00609     for(i = 0; i < totgrid; ++i) {
00610         DMGridData *grid= grids[i];
00611         BBC *bbc = prim_bbc + i;
00612 
00613         BB_reset((BB*)bbc);
00614 
00615         for(j = 0; j < gridsize*gridsize; ++j)
00616             BB_expand((BB*)bbc, grid[j].co);
00617 
00618         BBC_update_centroid(bbc);
00619 
00620         BB_expand(&cb, bbc->bcentroid);
00621     }
00622 
00623     if(totgrid)
00624         pbvh_build(bvh, &cb, prim_bbc, totgrid);
00625 
00626     MEM_freeN(prim_bbc);
00627 }
00628 
00629 PBVH *BLI_pbvh_new(void)
00630 {
00631     PBVH *bvh = MEM_callocN(sizeof(PBVH), "pbvh");
00632 
00633     return bvh;
00634 }
00635 
00636 void BLI_pbvh_free(PBVH *bvh)
00637 {
00638     PBVHNode *node;
00639     int i;
00640 
00641     for(i = 0; i < bvh->totnode; ++i) {
00642         node= &bvh->nodes[i];
00643 
00644         if(node->flag & PBVH_Leaf) {
00645             if(node->draw_buffers)
00646                 GPU_free_buffers(node->draw_buffers);
00647             if(node->vert_indices)
00648                 MEM_freeN(node->vert_indices);
00649             if(node->face_vert_indices)
00650                 MEM_freeN(node->face_vert_indices);
00651         }
00652     }
00653 
00654     if (bvh->deformed) {
00655         if (bvh->verts) {
00656             /* if pbvh was deformed, new memory was allocated for verts/faces -- free it */
00657 
00658             MEM_freeN(bvh->verts);
00659             if(bvh->faces)
00660                 MEM_freeN(bvh->faces);
00661         }
00662     }
00663 
00664     if(bvh->nodes)
00665         MEM_freeN(bvh->nodes);
00666 
00667     if(bvh->prim_indices)
00668         MEM_freeN(bvh->prim_indices);
00669 
00670     MEM_freeN(bvh);
00671 }
00672 
00673 static void pbvh_iter_begin(PBVHIter *iter, PBVH *bvh, BLI_pbvh_SearchCallback scb, void *search_data)
00674 {
00675     iter->bvh= bvh;
00676     iter->scb= scb;
00677     iter->search_data= search_data;
00678 
00679     iter->stack= iter->stackfixed;
00680     iter->stackspace= STACK_FIXED_DEPTH;
00681 
00682     iter->stack[0].node= bvh->nodes;
00683     iter->stack[0].revisiting= 0;
00684     iter->stacksize= 1;
00685 }
00686 
00687 static void pbvh_iter_end(PBVHIter *iter)
00688 {
00689     if(iter->stackspace > STACK_FIXED_DEPTH)
00690         MEM_freeN(iter->stack);
00691 }
00692 
00693 static void pbvh_stack_push(PBVHIter *iter, PBVHNode *node, int revisiting)
00694 {
00695     if(iter->stacksize == iter->stackspace) {
00696         PBVHStack *newstack;
00697 
00698         iter->stackspace *= 2;
00699         newstack= MEM_callocN(sizeof(PBVHStack)*iter->stackspace, "PBVHStack");
00700         memcpy(newstack, iter->stack, sizeof(PBVHStack)*iter->stacksize);
00701 
00702         if(iter->stackspace > STACK_FIXED_DEPTH)
00703             MEM_freeN(iter->stack);
00704         iter->stack= newstack;
00705     }
00706 
00707     iter->stack[iter->stacksize].node= node;
00708     iter->stack[iter->stacksize].revisiting= revisiting;
00709     iter->stacksize++;
00710 }
00711 
00712 static PBVHNode *pbvh_iter_next(PBVHIter *iter)
00713 {
00714     PBVHNode *node;
00715     int revisiting;
00716 
00717     /* purpose here is to traverse tree, visiting child nodes before their
00718        parents, this order is necessary for e.g. computing bounding boxes */
00719 
00720     while(iter->stacksize) {
00721         /* pop node */
00722         iter->stacksize--;
00723         node= iter->stack[iter->stacksize].node;
00724 
00725         /* on a mesh with no faces this can happen
00726          * can remove this check if we know meshes have at least 1 face */
00727         if(node==NULL)
00728             return NULL;
00729 
00730         revisiting= iter->stack[iter->stacksize].revisiting;
00731 
00732         /* revisiting node already checked */
00733         if(revisiting)
00734             return node;
00735 
00736         if(iter->scb && !iter->scb(node, iter->search_data))
00737             continue; /* don't traverse, outside of search zone */
00738 
00739         if(node->flag & PBVH_Leaf) {
00740             /* immediately hit leaf node */
00741             return node;
00742         }
00743         else {
00744             /* come back later when children are done */
00745             pbvh_stack_push(iter, node, 1);
00746 
00747             /* push two child nodes on the stack */
00748             pbvh_stack_push(iter, iter->bvh->nodes+node->children_offset+1, 0);
00749             pbvh_stack_push(iter, iter->bvh->nodes+node->children_offset, 0);
00750         }
00751     }
00752 
00753     return NULL;
00754 }
00755 
00756 static PBVHNode *pbvh_iter_next_occluded(PBVHIter *iter)
00757 {
00758     PBVHNode *node;
00759 
00760     while(iter->stacksize) {
00761         /* pop node */
00762         iter->stacksize--;
00763         node= iter->stack[iter->stacksize].node;
00764 
00765         /* on a mesh with no faces this can happen
00766          * can remove this check if we know meshes have at least 1 face */
00767         if(node==NULL) return NULL;
00768 
00769         if(iter->scb && !iter->scb(node, iter->search_data)) continue; /* don't traverse, outside of search zone */
00770 
00771         if(node->flag & PBVH_Leaf) {
00772             /* immediately hit leaf node */
00773             return node;
00774         }
00775         else {
00776             pbvh_stack_push(iter, iter->bvh->nodes+node->children_offset+1, 0);
00777             pbvh_stack_push(iter, iter->bvh->nodes+node->children_offset, 0);
00778         }
00779     }
00780 
00781     return NULL;
00782 }
00783 
00784 void BLI_pbvh_search_gather(PBVH *bvh,
00785     BLI_pbvh_SearchCallback scb, void *search_data,
00786     PBVHNode ***r_array, int *r_tot)
00787 {
00788     PBVHIter iter;
00789     PBVHNode **array= NULL, **newarray, *node;
00790     int tot= 0, space= 0;
00791 
00792     pbvh_iter_begin(&iter, bvh, scb, search_data);
00793 
00794     while((node=pbvh_iter_next(&iter))) {
00795         if(node->flag & PBVH_Leaf) {
00796             if(tot == space) {
00797                 /* resize array if needed */
00798                 space= (tot == 0)? 32: space*2;
00799                 newarray= MEM_callocN(sizeof(PBVHNode)*space, "PBVHNodeSearch");
00800 
00801                 if(array) {
00802                     memcpy(newarray, array, sizeof(PBVHNode)*tot);
00803                     MEM_freeN(array);
00804                 }
00805 
00806                 array= newarray;
00807             }
00808 
00809             array[tot]= node;
00810             tot++;
00811         }
00812     }
00813 
00814     pbvh_iter_end(&iter);
00815 
00816     if(tot == 0 && array) {
00817         MEM_freeN(array);
00818         array= NULL;
00819     }
00820 
00821     *r_array= array;
00822     *r_tot= tot;
00823 }
00824 
00825 void BLI_pbvh_search_callback(PBVH *bvh,
00826     BLI_pbvh_SearchCallback scb, void *search_data,
00827     BLI_pbvh_HitCallback hcb, void *hit_data)
00828 {
00829     PBVHIter iter;
00830     PBVHNode *node;
00831 
00832     pbvh_iter_begin(&iter, bvh, scb, search_data);
00833 
00834     while((node=pbvh_iter_next(&iter)))
00835         if (node->flag & PBVH_Leaf)
00836             hcb(node, hit_data);
00837 
00838     pbvh_iter_end(&iter);
00839 }
00840 
00841 typedef struct node_tree {
00842     PBVHNode* data;
00843 
00844     struct node_tree* left;
00845     struct node_tree* right;
00846 } node_tree;
00847 
00848 static void node_tree_insert(node_tree* tree, node_tree* new_node)
00849 {
00850     if (new_node->data->tmin < tree->data->tmin) {
00851         if (tree->left) {
00852             node_tree_insert(tree->left, new_node);
00853         }
00854         else {
00855             tree->left = new_node;
00856         }
00857     }
00858     else {
00859         if (tree->right) {
00860             node_tree_insert(tree->right, new_node);
00861         }
00862         else {
00863             tree->right = new_node;
00864         }
00865     }
00866 }
00867 
00868 static void traverse_tree(node_tree* tree, BLI_pbvh_HitOccludedCallback hcb, void* hit_data, float* tmin)
00869 {
00870     if (tree->left) traverse_tree(tree->left, hcb, hit_data, tmin);
00871 
00872     hcb(tree->data, hit_data, tmin);
00873 
00874     if (tree->right) traverse_tree(tree->right, hcb, hit_data, tmin);
00875 }
00876 
00877 static void free_tree(node_tree* tree)
00878 {
00879     if (tree->left) {
00880         free_tree(tree->left);
00881         tree->left = 0;
00882     }
00883 
00884     if (tree->right) {
00885         free_tree(tree->right);
00886         tree->right = 0;
00887     }
00888 
00889     free(tree);
00890 }
00891 
00892 float BLI_pbvh_node_get_tmin(PBVHNode* node)
00893 {
00894     return node->tmin;
00895 }
00896 
00897 static void BLI_pbvh_search_callback_occluded(PBVH *bvh,
00898     BLI_pbvh_SearchCallback scb, void *search_data,
00899     BLI_pbvh_HitOccludedCallback hcb, void *hit_data)
00900 {
00901     PBVHIter iter;
00902     PBVHNode *node;
00903     node_tree *tree = 0;
00904 
00905     pbvh_iter_begin(&iter, bvh, scb, search_data);
00906 
00907     while((node=pbvh_iter_next_occluded(&iter))) {
00908         if(node->flag & PBVH_Leaf) {
00909             node_tree* new_node = malloc(sizeof(node_tree));
00910 
00911             new_node->data = node;
00912 
00913             new_node->left  = NULL;
00914             new_node->right = NULL;
00915 
00916             if (tree) {
00917                 node_tree_insert(tree, new_node);
00918             }
00919             else {
00920                 tree = new_node;
00921             }
00922         }
00923     }
00924 
00925     pbvh_iter_end(&iter);
00926 
00927     if (tree) {
00928         float tmin = FLT_MAX;
00929         traverse_tree(tree, hcb, hit_data, &tmin);
00930         free_tree(tree);
00931     }
00932 }
00933 
00934 static int update_search_cb(PBVHNode *node, void *data_v)
00935 {
00936     int flag= GET_INT_FROM_POINTER(data_v);
00937 
00938     if(node->flag & PBVH_Leaf)
00939         return (node->flag & flag);
00940     
00941     return 1;
00942 }
00943 
00944 static void pbvh_update_normals(PBVH *bvh, PBVHNode **nodes,
00945     int totnode, float (*face_nors)[3])
00946 {
00947     float (*vnor)[3];
00948     int n;
00949 
00950     if(bvh->grids)
00951         return;
00952 
00953     /* could be per node to save some memory, but also means
00954        we have to store for each vertex which node it is in */
00955     vnor= MEM_callocN(sizeof(float)*3*bvh->totvert, "bvh temp vnors");
00956 
00957     /* subtle assumptions:
00958        - We know that for all edited vertices, the nodes with faces
00959          adjacent to these vertices have been marked with PBVH_UpdateNormals.
00960          This is true because if the vertex is inside the brush radius, the
00961          bounding box of it's adjacent faces will be as well.
00962        - However this is only true for the vertices that have actually been
00963          edited, not for all vertices in the nodes marked for update, so we
00964          can only update vertices marked with ME_VERT_PBVH_UPDATE.
00965     */
00966 
00967     #pragma omp parallel for private(n) schedule(static)
00968     for(n = 0; n < totnode; n++) {
00969         PBVHNode *node= nodes[n];
00970 
00971         if((node->flag & PBVH_UpdateNormals)) {
00972             int i, j, totface, *faces;
00973 
00974             faces= node->prim_indices;
00975             totface= node->totprim;
00976 
00977             for(i = 0; i < totface; ++i) {
00978                 MFace *f= bvh->faces + faces[i];
00979                 float fn[3];
00980                 unsigned int *fv = &f->v1;
00981                 int sides= (f->v4)? 4: 3;
00982 
00983                 if(f->v4)
00984                     normal_quad_v3(fn, bvh->verts[f->v1].co, bvh->verts[f->v2].co,
00985                                    bvh->verts[f->v3].co, bvh->verts[f->v4].co);
00986                 else
00987                     normal_tri_v3(fn, bvh->verts[f->v1].co, bvh->verts[f->v2].co,
00988                                   bvh->verts[f->v3].co);
00989 
00990                 for(j = 0; j < sides; ++j) {
00991                     int v= fv[j];
00992 
00993                     if(bvh->verts[v].flag & ME_VERT_PBVH_UPDATE) {
00994                         /* this seems like it could be very slow but profile
00995                            does not show this, so just leave it for now? */
00996                         #pragma omp atomic
00997                         vnor[v][0] += fn[0];
00998                         #pragma omp atomic
00999                         vnor[v][1] += fn[1];
01000                         #pragma omp atomic
01001                         vnor[v][2] += fn[2];
01002                     }
01003                 }
01004 
01005                 if(face_nors)
01006                     copy_v3_v3(face_nors[faces[i]], fn);
01007             }
01008         }
01009     }
01010 
01011     #pragma omp parallel for private(n) schedule(static)
01012     for(n = 0; n < totnode; n++) {
01013         PBVHNode *node= nodes[n];
01014 
01015         if(node->flag & PBVH_UpdateNormals) {
01016             int i, *verts, totvert;
01017 
01018             verts= node->vert_indices;
01019             totvert= node->uniq_verts;
01020 
01021             for(i = 0; i < totvert; ++i) {
01022                 const int v = verts[i];
01023                 MVert *mvert= &bvh->verts[v];
01024 
01025                 if(mvert->flag & ME_VERT_PBVH_UPDATE) {
01026                     float no[3];
01027 
01028                     copy_v3_v3(no, vnor[v]);
01029                     normalize_v3(no);
01030                     
01031                     mvert->no[0] = (short)(no[0]*32767.0f);
01032                     mvert->no[1] = (short)(no[1]*32767.0f);
01033                     mvert->no[2] = (short)(no[2]*32767.0f);
01034                     
01035                     mvert->flag &= ~ME_VERT_PBVH_UPDATE;
01036                 }
01037             }
01038 
01039             node->flag &= ~PBVH_UpdateNormals;
01040         }
01041     }
01042 
01043     MEM_freeN(vnor);
01044 }
01045 
01046 static void pbvh_update_BB_redraw(PBVH *bvh, PBVHNode **nodes,
01047     int totnode, int flag)
01048 {
01049     int n;
01050 
01051     /* update BB, redraw flag */
01052     #pragma omp parallel for private(n) schedule(static)
01053     for(n = 0; n < totnode; n++) {
01054         PBVHNode *node= nodes[n];
01055 
01056         if((flag & PBVH_UpdateBB) && (node->flag & PBVH_UpdateBB))
01057             /* don't clear flag yet, leave it for flushing later */
01058             update_node_vb(bvh, node);
01059 
01060         if((flag & PBVH_UpdateOriginalBB) && (node->flag & PBVH_UpdateOriginalBB))
01061             node->orig_vb= node->vb;
01062 
01063         if((flag & PBVH_UpdateRedraw) && (node->flag & PBVH_UpdateRedraw))
01064             node->flag &= ~PBVH_UpdateRedraw;
01065     }
01066 }
01067 
01068 static void pbvh_update_draw_buffers(PBVH *bvh, PBVHNode **nodes, int totnode, int smooth)
01069 {
01070     PBVHNode *node;
01071     int n;
01072 
01073     /* can't be done in parallel with OpenGL */
01074     for(n = 0; n < totnode; n++) {
01075         node= nodes[n];
01076 
01077         if(node->flag & PBVH_UpdateDrawBuffers) {
01078             if(bvh->grids) {
01079                 GPU_update_grid_buffers(node->draw_buffers,
01080                            bvh->grids,
01081                            node->prim_indices,
01082                            node->totprim,
01083                            bvh->gridsize,
01084                            smooth);
01085             }
01086             else {
01087                 GPU_update_mesh_buffers(node->draw_buffers,
01088                            bvh->verts,
01089                            node->vert_indices,
01090                            node->uniq_verts +
01091                            node->face_verts);
01092             }
01093 
01094             node->flag &= ~PBVH_UpdateDrawBuffers;
01095         }
01096     }
01097 }
01098 
01099 static int pbvh_flush_bb(PBVH *bvh, PBVHNode *node, int flag)
01100 {
01101     int update= 0;
01102 
01103     /* difficult to multithread well, we just do single threaded recursive */
01104     if(node->flag & PBVH_Leaf) {
01105         if(flag & PBVH_UpdateBB) {
01106             update |= (node->flag & PBVH_UpdateBB);
01107             node->flag &= ~PBVH_UpdateBB;
01108         }
01109 
01110         if(flag & PBVH_UpdateOriginalBB) {
01111             update |= (node->flag & PBVH_UpdateOriginalBB);
01112             node->flag &= ~PBVH_UpdateOriginalBB;
01113         }
01114 
01115         return update;
01116     }
01117     else {
01118         update |= pbvh_flush_bb(bvh, bvh->nodes + node->children_offset, flag);
01119         update |= pbvh_flush_bb(bvh, bvh->nodes + node->children_offset + 1, flag);
01120 
01121         if(update & PBVH_UpdateBB)
01122             update_node_vb(bvh, node);
01123         if(update & PBVH_UpdateOriginalBB)
01124             node->orig_vb= node->vb;
01125     }
01126 
01127     return update;
01128 }
01129 
01130 void BLI_pbvh_update(PBVH *bvh, int flag, float (*face_nors)[3])
01131 {
01132     PBVHNode **nodes;
01133     int totnode;
01134 
01135     if(!bvh->nodes)
01136         return;
01137 
01138     BLI_pbvh_search_gather(bvh, update_search_cb, SET_INT_IN_POINTER(flag),
01139         &nodes, &totnode);
01140 
01141     if(flag & PBVH_UpdateNormals)
01142         pbvh_update_normals(bvh, nodes, totnode, face_nors);
01143 
01144     if(flag & (PBVH_UpdateBB|PBVH_UpdateOriginalBB|PBVH_UpdateRedraw))
01145         pbvh_update_BB_redraw(bvh, nodes, totnode, flag);
01146 
01147     if(flag & (PBVH_UpdateBB|PBVH_UpdateOriginalBB))
01148         pbvh_flush_bb(bvh, bvh->nodes, flag);
01149 
01150     if(nodes) MEM_freeN(nodes);
01151 }
01152 
01153 void BLI_pbvh_redraw_BB(PBVH *bvh, float bb_min[3], float bb_max[3])
01154 {
01155     PBVHIter iter;
01156     PBVHNode *node;
01157     BB bb;
01158 
01159     BB_reset(&bb);
01160 
01161     pbvh_iter_begin(&iter, bvh, NULL, NULL);
01162 
01163     while((node=pbvh_iter_next(&iter)))
01164         if(node->flag & PBVH_UpdateRedraw)
01165             BB_expand_with_bb(&bb, &node->vb);
01166 
01167     pbvh_iter_end(&iter);
01168 
01169     copy_v3_v3(bb_min, bb.bmin);
01170     copy_v3_v3(bb_max, bb.bmax);
01171 }
01172 
01173 void BLI_pbvh_get_grid_updates(PBVH *bvh, int clear, void ***gridfaces, int *totface)
01174 {
01175     PBVHIter iter;
01176     PBVHNode *node;
01177     GHashIterator *hiter;
01178     GHash *map;
01179     void *face, **faces;
01180     unsigned i;
01181     int tot;
01182 
01183     map = BLI_ghash_new(BLI_ghashutil_ptrhash, BLI_ghashutil_ptrcmp, "pbvh_get_grid_updates gh");
01184 
01185     pbvh_iter_begin(&iter, bvh, NULL, NULL);
01186 
01187     while((node=pbvh_iter_next(&iter))) {
01188         if(node->flag & PBVH_UpdateNormals) {
01189             for(i = 0; i < node->totprim; ++i) {
01190                 face= bvh->gridfaces[node->prim_indices[i]];
01191                 if(!BLI_ghash_lookup(map, face))
01192                     BLI_ghash_insert(map, face, face);
01193             }
01194 
01195             if(clear)
01196                 node->flag &= ~PBVH_UpdateNormals;
01197         }
01198     }
01199 
01200     pbvh_iter_end(&iter);
01201     
01202     tot= BLI_ghash_size(map);
01203     if(tot == 0) {
01204         *totface= 0;
01205         *gridfaces= NULL;
01206         BLI_ghash_free(map, NULL, NULL);
01207         return;
01208     }
01209 
01210     faces= MEM_callocN(sizeof(void*)*tot, "PBVH Grid Faces");
01211 
01212     for(hiter = BLI_ghashIterator_new(map), i = 0;
01213         !BLI_ghashIterator_isDone(hiter);
01214         BLI_ghashIterator_step(hiter), ++i)
01215         faces[i]= BLI_ghashIterator_getKey(hiter);
01216 
01217     BLI_ghashIterator_free(hiter);
01218 
01219     BLI_ghash_free(map, NULL, NULL);
01220 
01221     *totface= tot;
01222     *gridfaces= faces;
01223 }
01224 
01225 /***************************** Node Access ***********************************/
01226 
01227 void BLI_pbvh_node_mark_update(PBVHNode *node)
01228 {
01229     node->flag |= PBVH_UpdateNormals|PBVH_UpdateBB|PBVH_UpdateOriginalBB|PBVH_UpdateDrawBuffers|PBVH_UpdateRedraw;
01230 }
01231 
01232 void BLI_pbvh_node_get_verts(PBVH *bvh, PBVHNode *node, int **vert_indices, MVert **verts)
01233 {
01234     if(vert_indices) *vert_indices= node->vert_indices;
01235     if(verts) *verts= bvh->verts;
01236 }
01237 
01238 void BLI_pbvh_node_num_verts(PBVH *bvh, PBVHNode *node, int *uniquevert, int *totvert)
01239 {
01240     if(bvh->grids) {
01241         const int tot= node->totprim*bvh->gridsize*bvh->gridsize;
01242         if(totvert) *totvert= tot;
01243         if(uniquevert) *uniquevert= tot;
01244     }
01245     else {
01246         if(totvert) *totvert= node->uniq_verts + node->face_verts;
01247         if(uniquevert) *uniquevert= node->uniq_verts;
01248     }
01249 }
01250 
01251 void BLI_pbvh_node_get_grids(PBVH *bvh, PBVHNode *node, int **grid_indices, int *totgrid, int *maxgrid, int *gridsize, DMGridData ***griddata, DMGridAdjacency **gridadj)
01252 {
01253     if(bvh->grids) {
01254         if(grid_indices) *grid_indices= node->prim_indices;
01255         if(totgrid) *totgrid= node->totprim;
01256         if(maxgrid) *maxgrid= bvh->totgrid;
01257         if(gridsize) *gridsize= bvh->gridsize;
01258         if(griddata) *griddata= bvh->grids;
01259         if(gridadj) *gridadj= bvh->gridadj;
01260     }
01261     else {
01262         if(grid_indices) *grid_indices= NULL;
01263         if(totgrid) *totgrid= 0;
01264         if(maxgrid) *maxgrid= 0;
01265         if(gridsize) *gridsize= 0;
01266         if(griddata) *griddata= NULL;
01267         if(gridadj) *gridadj= NULL;
01268     }
01269 }
01270 
01271 void BLI_pbvh_node_get_BB(PBVHNode *node, float bb_min[3], float bb_max[3])
01272 {
01273     copy_v3_v3(bb_min, node->vb.bmin);
01274     copy_v3_v3(bb_max, node->vb.bmax);
01275 }
01276 
01277 void BLI_pbvh_node_get_original_BB(PBVHNode *node, float bb_min[3], float bb_max[3])
01278 {
01279     copy_v3_v3(bb_min, node->orig_vb.bmin);
01280     copy_v3_v3(bb_max, node->orig_vb.bmax);
01281 }
01282 
01283 void BLI_pbvh_node_get_proxies(PBVHNode* node, PBVHProxyNode** proxies, int* proxy_count)
01284 {
01285     if (node->proxy_count > 0) {
01286         if (proxies) *proxies = node->proxies;
01287         if (proxy_count) *proxy_count = node->proxy_count;
01288     }
01289     else {
01290         if (proxies) *proxies = 0;
01291         if (proxy_count) *proxy_count = 0;
01292     }
01293 }
01294 
01295 /********************************* Raycast ***********************************/
01296 
01297 typedef struct {
01298     /* Ray */
01299     float start[3];
01300     int sign[3];
01301     float inv_dir[3];
01302     int original;
01303 } RaycastData;
01304 
01305 /* Adapted from here: http://www.gamedev.net/community/forums/topic.asp?topic_id=459973 */
01306 static int ray_aabb_intersect(PBVHNode *node, void *data_v)
01307 {
01308     RaycastData *ray = data_v;
01309     float bbox[2][3];
01310     float tmin, tmax, tymin, tymax, tzmin, tzmax;
01311 
01312     if(ray->original)
01313         BLI_pbvh_node_get_original_BB(node, bbox[0], bbox[1]);
01314     else
01315         BLI_pbvh_node_get_BB(node, bbox[0], bbox[1]);
01316 
01317     tmin = (bbox[ray->sign[0]][0] - ray->start[0]) * ray->inv_dir[0];
01318     tmax = (bbox[1-ray->sign[0]][0] - ray->start[0]) * ray->inv_dir[0];
01319 
01320     tymin = (bbox[ray->sign[1]][1] - ray->start[1]) * ray->inv_dir[1];
01321     tymax = (bbox[1-ray->sign[1]][1] - ray->start[1]) * ray->inv_dir[1];
01322 
01323     if((tmin > tymax) || (tymin > tmax))
01324         return 0;
01325 
01326     if(tymin > tmin)
01327         tmin = tymin;
01328 
01329     if(tymax < tmax)
01330         tmax = tymax;
01331 
01332     tzmin = (bbox[ray->sign[2]][2] - ray->start[2]) * ray->inv_dir[2];
01333     tzmax = (bbox[1-ray->sign[2]][2] - ray->start[2]) * ray->inv_dir[2];
01334 
01335     if((tmin > tzmax) || (tzmin > tmax))
01336         return 0;
01337 
01338     if(tzmin > tmin)
01339         tmin = tzmin;
01340 
01341     // XXX jwilkins: tmax does not need to be updated since we don't use it
01342     // keeping this here for future reference
01343     //if(tzmax < tmax) tmax = tzmax; 
01344 
01345     node->tmin = tmin;
01346 
01347     return 1;
01348 }
01349 
01350 void BLI_pbvh_raycast(PBVH *bvh, BLI_pbvh_HitOccludedCallback cb, void *data,
01351               float ray_start[3], float ray_normal[3], int original)
01352 {
01353     RaycastData rcd;
01354 
01355     copy_v3_v3(rcd.start, ray_start);
01356     rcd.inv_dir[0] = 1.0f / ray_normal[0];
01357     rcd.inv_dir[1] = 1.0f / ray_normal[1];
01358     rcd.inv_dir[2] = 1.0f / ray_normal[2];
01359     rcd.sign[0] = rcd.inv_dir[0] < 0;
01360     rcd.sign[1] = rcd.inv_dir[1] < 0;
01361     rcd.sign[2] = rcd.inv_dir[2] < 0;
01362     rcd.original = original;
01363 
01364     BLI_pbvh_search_callback_occluded(bvh, ray_aabb_intersect, &rcd, cb, data);
01365 }
01366 
01367 static int ray_face_intersection(float ray_start[3], float ray_normal[3],
01368                  float *t0, float *t1, float *t2, float *t3,
01369                  float *fdist)
01370 {
01371     float dist;
01372 
01373     if ((isect_ray_tri_epsilon_v3(ray_start, ray_normal, t0, t1, t2, &dist, NULL, 0.1f) && dist < *fdist) ||
01374        (t3 && isect_ray_tri_epsilon_v3(ray_start, ray_normal, t0, t2, t3, &dist, NULL, 0.1f) && dist < *fdist))
01375     {
01376         *fdist = dist;
01377         return 1;
01378     }
01379     else {
01380         return 0;
01381     }
01382 }
01383 
01384 int BLI_pbvh_node_raycast(PBVH *bvh, PBVHNode *node, float (*origco)[3],
01385     float ray_start[3], float ray_normal[3], float *dist)
01386 {
01387     int hit= 0;
01388 
01389     if(bvh->faces) {
01390         MVert *vert = bvh->verts;
01391         int *faces= node->prim_indices;
01392         int totface= node->totprim;
01393         int i;
01394 
01395         for(i = 0; i < totface; ++i) {
01396             MFace *f = bvh->faces + faces[i];
01397             int *face_verts = node->face_vert_indices[i];
01398 
01399             if(origco) {
01400                 /* intersect with backuped original coordinates */
01401                 hit |= ray_face_intersection(ray_start, ray_normal,
01402                              origco[face_verts[0]],
01403                              origco[face_verts[1]],
01404                              origco[face_verts[2]],
01405                              f->v4? origco[face_verts[3]]: NULL,
01406                              dist);
01407             }
01408             else {
01409                 /* intersect with current coordinates */
01410                 hit |= ray_face_intersection(ray_start, ray_normal,
01411                              vert[f->v1].co,
01412                              vert[f->v2].co,
01413                              vert[f->v3].co,
01414                              f->v4 ? vert[f->v4].co : NULL,
01415                              dist);
01416             }
01417         }
01418     }
01419     else {
01420         int totgrid= node->totprim;
01421         int gridsize= bvh->gridsize;
01422         int i, x, y;
01423 
01424         for(i = 0; i < totgrid; ++i) {
01425             DMGridData *grid= bvh->grids[node->prim_indices[i]];
01426             if (!grid)
01427                 continue;
01428 
01429             for(y = 0; y < gridsize-1; ++y) {
01430                 for(x = 0; x < gridsize-1; ++x) {
01431                     if(origco) {
01432                         hit |= ray_face_intersection(ray_start, ray_normal,
01433                                      origco[y*gridsize + x],
01434                                      origco[y*gridsize + x+1],
01435                                      origco[(y+1)*gridsize + x+1],
01436                                      origco[(y+1)*gridsize + x],
01437                                      dist);
01438                     }
01439                     else {
01440                         hit |= ray_face_intersection(ray_start, ray_normal,
01441                                      grid[y*gridsize + x].co,
01442                                      grid[y*gridsize + x+1].co,
01443                                      grid[(y+1)*gridsize + x+1].co,
01444                                      grid[(y+1)*gridsize + x].co,
01445                                      dist);
01446                     }
01447                 }
01448             }
01449 
01450             if(origco)
01451                 origco += gridsize*gridsize;
01452         }
01453     }
01454 
01455     return hit;
01456 }
01457 
01458 //#include <GL/glew.h>
01459 
01460 void BLI_pbvh_node_draw(PBVHNode *node, void *UNUSED(data))
01461 {
01462 #if 0
01463     /* XXX: Just some quick code to show leaf nodes in different colors */
01464     float col[3]; int i;
01465 
01466     if(0) { //is_partial) {
01467         col[0] = (rand() / (float)RAND_MAX); col[1] = col[2] = 0.6;
01468     }
01469     else {
01470         srand((long long)node);
01471         for(i = 0; i < 3; ++i)
01472             col[i] = (rand() / (float)RAND_MAX) * 0.3 + 0.7;
01473     }
01474     glMaterialfv(GL_FRONT_AND_BACK, GL_DIFFUSE, col);
01475 
01476     glColor3f(1, 0, 0);
01477 #endif
01478     GPU_draw_buffers(node->draw_buffers);
01479 }
01480 
01481 /* Adapted from:
01482    http://www.gamedev.net/community/forums/topic.asp?topic_id=512123
01483    Returns true if the AABB is at least partially within the frustum
01484    (ok, not a real frustum), false otherwise.
01485 */
01486 int BLI_pbvh_node_planes_contain_AABB(PBVHNode *node, void *data)
01487 {
01488     float (*planes)[4] = data;
01489     int i, axis;
01490     float vmin[3] /*, vmax[3]*/, bb_min[3], bb_max[3];
01491 
01492     BLI_pbvh_node_get_BB(node, bb_min, bb_max);
01493 
01494     for(i = 0; i < 4; ++i) { 
01495         for(axis = 0; axis < 3; ++axis) {
01496             if(planes[i][axis] > 0) { 
01497                 vmin[axis] = bb_min[axis];
01498                 /*vmax[axis] = bb_max[axis];*/ /*UNUSED*/
01499             }
01500             else {
01501                 vmin[axis] = bb_max[axis];
01502                 /*vmax[axis] = bb_min[axis];*/ /*UNUSED*/
01503             }
01504         }
01505         
01506         if(dot_v3v3(planes[i], vmin) + planes[i][3] > 0)
01507             return 0;
01508     } 
01509 
01510     return 1;
01511 }
01512 
01513 void BLI_pbvh_draw(PBVH *bvh, float (*planes)[4], float (*face_nors)[3], int smooth)
01514 {
01515     PBVHNode **nodes;
01516     int totnode;
01517 
01518     BLI_pbvh_search_gather(bvh, update_search_cb, SET_INT_IN_POINTER(PBVH_UpdateNormals|PBVH_UpdateDrawBuffers),
01519         &nodes, &totnode);
01520 
01521     pbvh_update_normals(bvh, nodes, totnode, face_nors);
01522     pbvh_update_draw_buffers(bvh, nodes, totnode, smooth);
01523 
01524     if(nodes) MEM_freeN(nodes);
01525 
01526     if(planes) {
01527         BLI_pbvh_search_callback(bvh, BLI_pbvh_node_planes_contain_AABB,
01528                 planes, BLI_pbvh_node_draw, NULL);
01529     }
01530     else {
01531         BLI_pbvh_search_callback(bvh, NULL, NULL, BLI_pbvh_node_draw, NULL);
01532     }
01533 }
01534 
01535 void BLI_pbvh_grids_update(PBVH *bvh, DMGridData **grids, DMGridAdjacency *gridadj, void **gridfaces)
01536 {
01537     bvh->grids= grids;
01538     bvh->gridadj= gridadj;
01539     bvh->gridfaces= gridfaces;
01540 }
01541 
01542 float (*BLI_pbvh_get_vertCos(PBVH *pbvh))[3]
01543 {
01544     int a;
01545     float (*vertCos)[3]= NULL;
01546 
01547     if (pbvh->verts) {
01548         float *co;
01549         MVert *mvert= pbvh->verts;
01550 
01551         vertCos= MEM_callocN(3*pbvh->totvert*sizeof(float), "BLI_pbvh_get_vertCoords");
01552         co= (float*)vertCos;
01553 
01554         for (a= 0; a<pbvh->totvert; a++, mvert++, co+= 3) {
01555             copy_v3_v3(co, mvert->co);
01556         }
01557     }
01558 
01559     return vertCos;
01560 }
01561 
01562 void BLI_pbvh_apply_vertCos(PBVH *pbvh, float (*vertCos)[3])
01563 {
01564     int a;
01565 
01566     if (!pbvh->deformed) {
01567         if (pbvh->verts) {
01568             /* if pbvh is not already deformed, verts/faces points to the */
01569             /* original data and applying new coords to this arrays would lead to */
01570             /* unneeded deformation -- duplicate verts/faces to avoid this */
01571 
01572             pbvh->verts= MEM_dupallocN(pbvh->verts);
01573             pbvh->faces= MEM_dupallocN(pbvh->faces);
01574 
01575             pbvh->deformed= 1;
01576         }
01577     }
01578 
01579     if (pbvh->verts) {
01580         MVert *mvert= pbvh->verts;
01581         /* copy new verts coords */
01582         for (a= 0; a < pbvh->totvert; ++a, ++mvert) {
01583             copy_v3_v3(mvert->co, vertCos[a]);
01584             mvert->flag |= ME_VERT_PBVH_UPDATE;
01585         }
01586 
01587         /* coordinates are new -- normals should also be updated */
01588         mesh_calc_normals(pbvh->verts, pbvh->totvert, pbvh->faces, pbvh->totprim, NULL);
01589 
01590         for (a= 0; a < pbvh->totnode; ++a)
01591             BLI_pbvh_node_mark_update(&pbvh->nodes[a]);
01592 
01593         BLI_pbvh_update(pbvh, PBVH_UpdateBB, NULL);
01594         BLI_pbvh_update(pbvh, PBVH_UpdateOriginalBB, NULL);
01595 
01596     }
01597 }
01598 
01599 int BLI_pbvh_isDeformed(PBVH *pbvh)
01600 {
01601     return pbvh->deformed;
01602 }
01603 /* Proxies */
01604 
01605 PBVHProxyNode* BLI_pbvh_node_add_proxy(PBVH* bvh, PBVHNode* node)
01606 {
01607     int index, totverts;
01608 
01609     #pragma omp critical
01610     {
01611 
01612         index = node->proxy_count;
01613 
01614         node->proxy_count++;
01615 
01616         if (node->proxies)
01617             node->proxies= MEM_reallocN(node->proxies, node->proxy_count*sizeof(PBVHProxyNode));
01618         else
01619             node->proxies= MEM_mallocN(sizeof(PBVHProxyNode), "PBVHNodeProxy");
01620 
01621         if (bvh->grids)
01622             totverts = node->totprim*bvh->gridsize*bvh->gridsize;
01623         else
01624             totverts = node->uniq_verts;
01625 
01626         node->proxies[index].co= MEM_callocN(sizeof(float[3])*totverts, "PBVHNodeProxy.co");
01627     }
01628 
01629     return node->proxies + index;
01630 }
01631 
01632 void BLI_pbvh_node_free_proxies(PBVHNode* node)
01633 {
01634     #pragma omp critical
01635     {
01636         int p;
01637 
01638         for (p= 0; p < node->proxy_count; p++) {
01639             MEM_freeN(node->proxies[p].co);
01640             node->proxies[p].co= 0;
01641         }
01642 
01643         MEM_freeN(node->proxies);
01644         node->proxies = 0;
01645 
01646         node->proxy_count= 0;
01647     }
01648 }
01649 
01650 void BLI_pbvh_gather_proxies(PBVH* pbvh, PBVHNode*** r_array,  int* r_tot)
01651 {
01652     PBVHNode **array= NULL, **newarray, *node;
01653     int tot= 0, space= 0;
01654     int n;
01655 
01656     for (n= 0; n < pbvh->totnode; n++) {
01657         node = pbvh->nodes + n;
01658 
01659         if(node->proxy_count > 0) {
01660             if(tot == space) {
01661                 /* resize array if needed */
01662                 space= (tot == 0)? 32: space*2;
01663                 newarray= MEM_callocN(sizeof(PBVHNode)*space, "BLI_pbvh_gather_proxies");
01664 
01665                 if (array) {
01666                     memcpy(newarray, array, sizeof(PBVHNode)*tot);
01667                     MEM_freeN(array);
01668                 }
01669 
01670                 array= newarray;
01671             }
01672 
01673             array[tot]= node;
01674             tot++;
01675         }
01676     }
01677 
01678     if(tot == 0 && array) {
01679         MEM_freeN(array);
01680         array= NULL;
01681     }
01682 
01683     *r_array= array;
01684     *r_tot= tot;
01685 }