intern/bsp/intern/BOP_CarveInterface.cpp

   1 /*
   2  * ***** BEGIN GPL LICENSE BLOCK *****
   3  *
   4  * This program is free software; you can redistribute it and/or
   5  * modify it under the terms of the GNU General Public License
   6  * as published by the Free Software Foundation; either version 2
   7  * of the License, or (at your option) any later version.
   8  *
   9  * This program is distributed in the hope that it will be useful,
  10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12  * GNU General Public License for more details.
  13  *
  14  * You should have received a copy of the GNU General Public License
  15  * along with this program; if not, write to the Free Software Foundation,
  16  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  17  *
  18  * The Original Code is Copyright (C) 2010 by the Blender Foundation.
  19  * All rights reserved.
  20  *
  21  * The Original Code is: all of this file.
  22  *
  23  * Contributor(s): Ken Hughes,
  24  *                 Sergey Sharybin.
  25  *
  26  * ***** END GPL LICENSE BLOCK *****
  27  */
  28
  29 /** \file bsp/intern/BOP_CarveInterface.cpp
  30  *  \ingroup bsp
  31  */
  32
  33 #include "BOP_Interface.h"
  34 #include "BSP_CSGMesh_CFIterator.h"
  35
  36 #include <carve/csg_triangulator.hpp>
  37 #include <carve/interpolator.hpp>
  38 #include <carve/rescale.hpp>
  39
  40 #include <iostream>
  41
  42 using namespace carve::mesh;
  43 using namespace carve::geom;
  44 typedef unsigned int uint;
  45
  46 #define MAX(x,y) ((x)>(y)?(x):(y))
  47 #define MIN(x,y) ((x)<(y)?(x):(y))
  48
  49 static bool isQuadPlanar(carve::geom3d::Vector &v1, carve::geom3d::Vector &v2,
  50                          carve::geom3d::Vector &v3, carve::geom3d::Vector &v4)
  51 {
  52         carve::geom3d::Vector vec1, vec2, vec3, cross;
  53
  54         vec1 = v2 - v1;
  55         vec2 = v4 - v1;
  56         vec3 = v3 - v1;
  57
  58         cross = carve::geom::cross(vec1, vec2);
  59
  60         float production = carve::geom::dot(cross, vec3);
  61         float magnitude = 1e-6 * cross.length();
  62
  63         return fabs(production) < magnitude;
  64 }
  65
  66 static bool isFacePlanar(CSG_IFace &face, std::vector<carve::geom3d::Vector> &vertices)
  67 {
  68         if (face.vertex_number == 4) {
  69                 return isQuadPlanar(vertices[face.vertex_index[0]], vertices[face.vertex_index[1]],
  70                                     vertices[face.vertex_index[2]], vertices[face.vertex_index[3]]);
  71         }
  72
  73         return true;
  74 }
  75
  76 static void Carve_copyMeshes(std::vector<MeshSet<3>::mesh_t*> &meshes, std::vector<MeshSet<3>::mesh_t*> &new_meshes)
  77 {
  78         std::vector<MeshSet<3>::mesh_t*>::iterator it = meshes.begin();
  79
  80         for(; it!=meshes.end(); it++) {
  81                 MeshSet<3>::mesh_t *mesh = *it;
  82                 MeshSet<3>::mesh_t *new_mesh = new MeshSet<3>::mesh_t(mesh->faces);
  83
  84                 new_meshes.push_back(new_mesh);
  85         }
  86 }
  87
  88 static MeshSet<3> *Carve_meshSetFromMeshes(std::vector<MeshSet<3>::mesh_t*> &meshes)
  89 {
  90         std::vector<MeshSet<3>::mesh_t*> new_meshes;
  91
  92         Carve_copyMeshes(meshes, new_meshes);
  93
  94         return new MeshSet<3>(new_meshes);
  95 }
  96
  97 static MeshSet<3> *Carve_meshSetFromTwoMeshes(std::vector<MeshSet<3>::mesh_t*> &left_meshes,
  98                                               std::vector<MeshSet<3>::mesh_t*> &right_meshes)
  99 {
 100         std::vector<MeshSet<3>::mesh_t*> new_meshes;
 101
 102         Carve_copyMeshes(left_meshes, new_meshes);
 103         Carve_copyMeshes(right_meshes, new_meshes);
 104
 105         return new MeshSet<3>(new_meshes);
 106 }
 107
 108 static bool Carve_checkEdgeFaceIntersections_do(carve::csg::Intersections &intersections,
 109                                                 MeshSet<3>::face_t *face_a, MeshSet<3>::edge_t *edge_b)
 110 {
 111         if(intersections.intersects(edge_b, face_a))
 112                 return true;
 113
 114         carve::mesh::MeshSet<3>::vertex_t::vector_t _p;
 115         if(face_a->simpleLineSegmentIntersection(carve::geom3d::LineSegment(edge_b->v1()->v, edge_b->v2()->v), _p))
 116                 return true;
 117
 118         return false;
 119 }
 120
 121 static bool Carve_checkEdgeFaceIntersections(carve::csg::Intersections &intersections,
 122                                              MeshSet<3>::face_t *face_a, MeshSet<3>::face_t *face_b)
 123 {
 124         MeshSet<3>::edge_t *edge_b;
 125
 126         edge_b = face_b->edge;
 127         do {
 128                 if(Carve_checkEdgeFaceIntersections_do(intersections, face_a, edge_b))
 129                         return true;
 130                 edge_b = edge_b->next;
 131         } while (edge_b != face_b->edge);
 132
 133         return false;
 134 }
 135
 136 static inline bool Carve_facesAreCoplanar(const MeshSet<3>::face_t *a, const MeshSet<3>::face_t *b)
 137 {
 138   carve::geom3d::Ray temp;
 139   // XXX: Find a better definition. This may be a source of problems
 140   // if floating point inaccuracies cause an incorrect answer.
 141   return !carve::geom3d::planeIntersection(a->plane, b->plane, temp);
 142 }
 143
 144 static bool Carve_checkMeshSetInterseciton_do(carve::csg::Intersections &intersections,
 145                                               const RTreeNode<3, Face<3> *> *a_node,
 146                                               const RTreeNode<3, Face<3> *> *b_node,
 147                                               bool descend_a = true)
 148 {
 149         if(!a_node->bbox.intersects(b_node->bbox))
 150                 return false;
 151
 152         if(a_node->child && (descend_a || !b_node->child)) {
 153                 for(RTreeNode<3, Face<3> *> *node = a_node->child; node; node = node->sibling) {
 154                         if(Carve_checkMeshSetInterseciton_do(intersections, node, b_node, false))
 155                                 return true;
 156                 }
 157         }
 158         else if(b_node->child) {
 159                 for(RTreeNode<3, Face<3> *> *node = b_node->child; node; node = node->sibling) {
 160                         if(Carve_checkMeshSetInterseciton_do(intersections, a_node, node, true))
 161                                 return true;
 162                 }
 163         }
 164         else {
 165                 for(size_t i = 0; i < a_node->data.size(); ++i) {
 166                         MeshSet<3>::face_t *fa = a_node->data[i];
 167                         aabb<3> aabb_a = fa->getAABB();
 168                         if(aabb_a.maxAxisSeparation(b_node->bbox) > carve::EPSILON) continue;
 169
 170                         for(size_t j = 0; j < b_node->data.size(); ++j) {
 171                                 MeshSet<3>::face_t *fb = b_node->data[j];
 172                                 aabb<3> aabb_b = fb->getAABB();
 173                                 if(aabb_b.maxAxisSeparation(aabb_a) > carve::EPSILON) continue;
 174
 175                                 std::pair<double, double> a_ra = fa->rangeInDirection(fa->plane.N, fa->edge->vert->v);
 176                                 std::pair<double, double> b_ra = fb->rangeInDirection(fa->plane.N, fa->edge->vert->v);
 177                                 if(carve::rangeSeparation(a_ra, b_ra) > carve::EPSILON) continue;
 178
 179                                 std::pair<double, double> a_rb = fa->rangeInDirection(fb->plane.N, fb->edge->vert->v);
 180                                 std::pair<double, double> b_rb = fb->rangeInDirection(fb->plane.N, fb->edge->vert->v);
 181                                 if(carve::rangeSeparation(a_rb, b_rb) > carve::EPSILON) continue;
 182
 183                                 if(!Carve_facesAreCoplanar(fa, fb)) {
 184                                         if(Carve_checkEdgeFaceIntersections(intersections, fa, fb)) {
 185                                                 return true;
 186                                         }
 187                                 }
 188                         }
 189                 }
 190         }
 191
 192         return false;
 193 }
 194
 195 static bool Carve_checkMeshSetInterseciton(RTreeNode<3, Face<3> *> *rtree_a, RTreeNode<3, Face<3> *> *rtree_b)
 196 {
 197         carve::csg::Intersections intersections;
 198
 199         return Carve_checkMeshSetInterseciton_do(intersections, rtree_a, rtree_b);
 200 }
 201
 202 static void Carve_getIntersectedOperandMeshes(std::vector<MeshSet<3>::mesh_t*> &meshes, MeshSet<3>::aabb_t &otherAABB,
 203                                               std::vector<MeshSet<3>::mesh_t*> &operandMeshes)
 204 {
 205         std::vector<MeshSet<3>::mesh_t*>::iterator it = meshes.begin();
 206         std::vector< RTreeNode<3, Face<3> *> *> meshRTree;
 207
 208         while(it != meshes.end()) {
 209                 MeshSet<3>::mesh_t *mesh = *it;
 210                 bool isAdded = false;
 211
 212                 RTreeNode<3, Face<3> *> *rtree = RTreeNode<3, Face<3> *>::construct_STR(mesh->faces.begin(), mesh->faces.end(), 4, 4);
 213
 214                 if (rtree->bbox.intersects(otherAABB)) {
 215                         bool isIntersect = false;
 216
 217                         std::vector<MeshSet<3>::mesh_t*>::iterator operand_it = operandMeshes.begin();
 218                         std::vector<RTreeNode<3, Face<3> *> *>::iterator tree_it = meshRTree.begin();
 219                         for(; operand_it!=operandMeshes.end(); operand_it++, tree_it++) {
 220                                 RTreeNode<3, Face<3> *> *operandRTree = *tree_it;
 221
 222                                 if(Carve_checkMeshSetInterseciton(rtree, operandRTree)) {
 223                                         isIntersect = true;
 224                                         break;
 225                                 }
 226                         }
 227
 228                         if(!isIntersect) {
 229                                 operandMeshes.push_back(mesh);
 230                                 meshRTree.push_back(rtree);
 231
 232                                 it = meshes.erase(it);
 233                                 isAdded = true;
 234                         }
 235                 }
 236
 237                 if (!isAdded) {
 238                         delete rtree;
 239                         it++;
 240                 }
 241         }
 242
 243         std::vector<RTreeNode<3, Face<3> *> *>::iterator tree_it = meshRTree.begin();
 244         for(; tree_it != meshRTree.end(); tree_it++) {
 245                 delete *tree_it;
 246         }
 247 }
 248
 249 static MeshSet<3> *Carve_getIntersectedOperand(std::vector<MeshSet<3>::mesh_t*> &meshes, MeshSet<3>::aabb_t &otherAABB)
 250 {
 251         std::vector<MeshSet<3>::mesh_t*> operandMeshes;
 252         Carve_getIntersectedOperandMeshes(meshes, otherAABB, operandMeshes);
 253
 254         if (operandMeshes.size() == 0)
 255                 return NULL;
 256
 257         return Carve_meshSetFromMeshes(operandMeshes);
 258 }
 259
 260 static MeshSet<3> *Carve_unionIntersectingMeshes(MeshSet<3> *poly,
 261                                                  MeshSet<3>::aabb_t &otherAABB,
 262                                                  carve::interpolate::FaceAttr<uint> &oface_num)
 263 {
 264         if(poly->meshes.size()<=1)
 265                 return poly;
 266
 267         carve::csg::CSG csg;
 268
 269         oface_num.installHooks(csg);
 270         csg.hooks.registerHook(new carve::csg::CarveTriangulator, carve::csg::CSG::Hooks::PROCESS_OUTPUT_FACE_BIT);
 271
 272         std::vector<MeshSet<3>::mesh_t*> orig_meshes =
 273                         std::vector<MeshSet<3>::mesh_t*>(poly->meshes.begin(), poly->meshes.end());
 274
 275         MeshSet<3> *left = Carve_getIntersectedOperand(orig_meshes, otherAABB);
 276
 277         if (!left) {
 278                 /* no maniforlds which intersects another object at all */
 279                 return poly;
 280         }
 281
 282         while(orig_meshes.size()) {
 283                 MeshSet<3> *right = Carve_getIntersectedOperand(orig_meshes, otherAABB);
 284
 285                 if (!right) {
 286                         /* no more intersecting manifolds which intersects other object */
 287                         break;
 288                 }
 289
 290                 try {
 291                         if(left->meshes.size()==0) {
 292                                 delete left;
 293
 294                                 left = right;
 295                         }
 296                         else {
 297                                 MeshSet<3> *result = csg.compute(left, right, carve::csg::CSG::UNION, NULL, carve::csg::CSG::CLASSIFY_EDGE);
 298
 299                                 delete left;
 300                                 delete right;
 301
 302                                 left = result;
 303                         }
 304                 }
 305                 catch(carve::exception e) {
 306                         std::cerr << "CSG failed, exception " << e.str() << std::endl;
 307
 308                         MeshSet<3> *result = Carve_meshSetFromTwoMeshes(left->meshes, right->meshes);
 309
 310                         delete left;
 311                         delete right;
 312
 313                         left = result;
 314                 }
 315                 catch(...) {
 316                         delete left;
 317                         delete right;
 318
 319                         throw "Unknown error in Carve library";
 320                 }
 321         }
 322
 323         /* append all meshes which doesn't have intersection with another operand as-is */
 324         if (orig_meshes.size()) {
 325                 MeshSet<3> *result = Carve_meshSetFromTwoMeshes(left->meshes, orig_meshes);
 326
 327                 delete left;
 328
 329                 return result;
 330         }
 331
 332         return left;
 333 }
 334
 335 static void Carve_unionIntersections(MeshSet<3> **left_r, MeshSet<3> **right_r,
 336                                      carve::interpolate::FaceAttr<uint> &oface_num)
 337 {
 338         MeshSet<3> *left, *right;
 339
 340         MeshSet<3>::aabb_t leftAABB = (*left_r)->getAABB();
 341         MeshSet<3>::aabb_t rightAABB = (*right_r)->getAABB();
 342
 343         left = Carve_unionIntersectingMeshes(*left_r, rightAABB, oface_num);
 344         right = Carve_unionIntersectingMeshes(*right_r, leftAABB, oface_num);
 345
 346         if(left != *left_r)
 347                 delete *left_r;
 348
 349         if(right != *right_r)
 350                 delete *right_r;
 351
 352         *left_r = left;
 353         *right_r = right;
 354 }
 355
 356 static MeshSet<3> *Carve_addMesh(CSG_FaceIteratorDescriptor &face_it,
 357                                  CSG_VertexIteratorDescriptor &vertex_it,
 358                                  carve::interpolate::FaceAttr<uint> &oface_num,
 359                                  uint &num_origfaces)
 360 {
 361         CSG_IVertex vertex;
 362         std::vector<carve::geom3d::Vector> vertices;
 363
 364         while (!vertex_it.Done(vertex_it.it)) {
 365                 vertex_it.Fill(vertex_it.it,&vertex);
 366                 vertices.push_back(VECTOR(vertex.position[0],
 367                                           vertex.position[1],
 368                                           vertex.position[2]));
 369                 vertex_it.Step(vertex_it.it);
 370         }
 371
 372         CSG_IFace face;
 373         std::vector<int> f;
 374         int numfaces = 0;
 375
 376         // now for the polygons.
 377         // we may need to decalare some memory for user defined face properties.
 378
 379         std::vector<int> forig;
 380         while (!face_it.Done(face_it.it)) {
 381                 face_it.Fill(face_it.it,&face);
 382
 383                 if (isFacePlanar(face, vertices)) {
 384                         f.push_back(face.vertex_number);
 385                         f.push_back(face.vertex_index[0]);
 386                         f.push_back(face.vertex_index[1]);
 387                         f.push_back(face.vertex_index[2]);
 388
 389                         if (face.vertex_number == 4)
 390                                 f.push_back(face.vertex_index[3]);
 391
 392                         forig.push_back(face.orig_face);
 393                         ++numfaces;
 394                         face_it.Step(face_it.it);
 395                         ++num_origfaces;
 396                 }
 397                 else {
 398                         f.push_back(3);
 399                         f.push_back(face.vertex_index[0]);
 400                         f.push_back(face.vertex_index[1]);
 401                         f.push_back(face.vertex_index[2]);
 402
 403                         forig.push_back(face.orig_face);
 404                         ++numfaces;
 405
 406                         if (face.vertex_number == 4) {
 407                                 f.push_back(3);
 408                                 f.push_back(face.vertex_index[0]);
 409                                 f.push_back(face.vertex_index[2]);
 410                                 f.push_back(face.vertex_index[3]);
 411
 412                                 forig.push_back(face.orig_face);
 413                                 ++numfaces;
 414                         }
 415
 416                         face_it.Step(face_it.it);
 417                         ++num_origfaces;
 418                 }
 419         }
 420
 421         MeshSet<3> *poly = new MeshSet<3> (vertices, numfaces, f);
 422
 423         uint i;
 424         MeshSet<3>::face_iter face_iter = poly->faceBegin();
 425         for (i = 0; face_iter != poly->faceEnd(); ++face_iter, ++i) {
 426                 MeshSet<3>::face_t *face = *face_iter;
 427                 oface_num.setAttribute(face, forig[i]);
 428         }
 429
 430         return poly;
 431 }
 432
 433 static double triangleArea(carve::geom3d::Vector &v1, carve::geom3d::Vector &v2, carve::geom3d::Vector &v3)
 434 {
 435         carve::geom3d::Vector a = v2 - v1;
 436         carve::geom3d::Vector b = v3 - v1;
 437
 438         return carve::geom::cross(a, b).length();
 439 }
 440
 441 static bool checkValidQuad(std::vector<MeshSet<3>::vertex_t> &vertex_storage, uint quad[4])
 442 {
 443         carve::geom3d::Vector &v1 = vertex_storage[quad[0]].v;
 444         carve::geom3d::Vector &v2 = vertex_storage[quad[1]].v;
 445         carve::geom3d::Vector &v3 = vertex_storage[quad[2]].v;
 446         carve::geom3d::Vector &v4 = vertex_storage[quad[3]].v;
 447
 448 #if 0
 449         /* disabled for now to prevent initially non-planar be triangulated
 450          * in theory this might cause some artifacts if intersections happens by non-planar
 451          * non-concave quad, but in practice it's acceptable */
 452         if (!isQuadPlanar(v1, v2, v3, v4)) {
 453                 /* non-planar faces better not be merged because of possible differences in triangulation
 454                  * of non-planar faces in opengl and renderer */
 455                 return false;
 456         }
 457 #endif
 458
 459         carve::geom3d::Vector edges[4];
 460         carve::geom3d::Vector normal;
 461         bool normal_set = false;
 462
 463         edges[0] = v2 - v1;
 464         edges[1] = v3 - v2;
 465         edges[2] = v4 - v3;
 466         edges[3] = v1 - v4;
 467
 468         for (int i = 0; i < 4; i++) {
 469                 int n = i + 1;
 470
 471                 if (n == 4)
 472                         n = 0;
 473
 474                 carve::geom3d::Vector current_normal = carve::geom::cross(edges[i], edges[n]);
 475
 476                 if (current_normal.length() > DBL_EPSILON) {
 477                         if (!normal_set) {
 478                                 normal = current_normal;
 479                                 normal_set = true;
 480                         }
 481                         else if (carve::geom::dot(normal, current_normal) < 0) {
 482                                 return false;
 483                         }
 484                 }
 485         }
 486
 487         if (!normal_set) {
 488                 /* normal wasn't set means face is degraded and better merge it in such way */
 489                 return false;
 490         }
 491
 492         double area = triangleArea(v1, v2, v3) + triangleArea(v1, v3, v4);
 493         if (area <= DBL_EPSILON)
 494                 return false;
 495
 496         return true;
 497 }
 498
 499 // check whether two faces share an edge, and if so merge them
 500 static uint quadMerge(std::map<MeshSet<3>::vertex_t*, uint> *vertexToIndex_map,
 501                                           std::vector<MeshSet<3>::vertex_t> &vertex_storage,
 502                       MeshSet<3>::face_t *f1, MeshSet<3>::face_t *f2,
 503                       uint v, uint quad[4])
 504 {
 505         uint current, n1, p1, n2, p2;
 506         uint v1[3];
 507         uint v2[3];
 508
 509         // get the vertex indices for each face
 510         v1[0] = vertexToIndex_map->find(f1->edge->vert)->second;
 511         v1[1] = vertexToIndex_map->find(f1->edge->next->vert)->second;
 512         v1[2] = vertexToIndex_map->find(f1->edge->next->next->vert)->second;
 513
 514         v2[0] = vertexToIndex_map->find(f2->edge->vert)->second;
 515         v2[1] = vertexToIndex_map->find(f2->edge->next->vert)->second;
 516         v2[2] = vertexToIndex_map->find(f2->edge->next->next->vert)->second;
 517
 518         // locate the current vertex we're examining, and find the next and
 519         // previous vertices based on the face windings
 520         if (v1[0] == v) {current = 0; p1 = 2; n1 = 1;}
 521         else if (v1[1] == v) {current = 1; p1 = 0; n1 = 2;}
 522         else {current = 2; p1 = 1; n1 = 0;}
 523
 524         if (v2[0] == v) {p2 = 2; n2 = 1;}
 525         else if (v2[1] == v) {p2 = 0; n2 = 2;}
 526         else {p2 = 1; n2 = 0;}
 527
 528         // if we find a match, place indices into quad in proper order and return
 529         // success code
 530         if (v1[p1] == v2[n2]) {
 531                 quad[0] = v1[current];
 532                 quad[1] = v1[n1];
 533                 quad[2] = v1[p1];
 534                 quad[3] = v2[p2];
 535
 536                 return checkValidQuad(vertex_storage, quad);
 537         }
 538         else if (v1[n1] == v2[p2]) {
 539                 quad[0] = v1[current];
 540                 quad[1] = v2[n2];
 541                 quad[2] = v1[n1];
 542                 quad[3] = v1[p1];
 543
 544                 return checkValidQuad(vertex_storage, quad);
 545         }
 546
 547         return 0;
 548 }
 549
 550 static bool Carve_checkDegeneratedFace(std::map<MeshSet<3>::vertex_t*, uint> *vertexToIndex_map, MeshSet<3>::face_t *face)
 551 {
 552         /* only tris and quads for now */
 553         if (face->n_edges == 3) {
 554                 uint v1, v2, v3;
 555
 556                 v1 = vertexToIndex_map->find(face->edge->prev->vert)->second;
 557                 v2 = vertexToIndex_map->find(face->edge->vert)->second;
 558                 v3 = vertexToIndex_map->find(face->edge->next->vert)->second;
 559
 560                 if (v1 == v2 || v2 == v3 || v1 == v3)
 561                         return true;
 562
 563                 return triangleArea(face->edge->prev->vert->v, face->edge->vert->v, face->edge->next->vert->v) < DBL_EPSILON;
 564         }
 565         else if (face->n_edges == 4) {
 566                 uint v1, v2, v3, v4;
 567
 568                 v1 = vertexToIndex_map->find(face->edge->prev->vert)->second;
 569                 v2 = vertexToIndex_map->find(face->edge->vert)->second;
 570                 v3 = vertexToIndex_map->find(face->edge->next->vert)->second;
 571                 v4 = vertexToIndex_map->find(face->edge->next->next->vert)->second;
 572
 573                 if (v1 == v2 || v1 == v3 || v1 == v4 || v2 == v3 || v2 == v4 || v3 == v4)
 574                         return true;
 575
 576                 return triangleArea(face->edge->vert->v, face->edge->next->vert->v, face->edge->next->next->vert->v) +
 577                        triangleArea(face->edge->prev->vert->v, face->edge->vert->v, face->edge->next->next->vert->v) < DBL_EPSILON;
 578         }
 579
 580         return false;
 581 }
 582
 583 static BSP_CSGMesh *Carve_exportMesh(MeshSet<3>* &poly, carve::interpolate::FaceAttr<uint> &oface_num,
 584                                      uint num_origfaces)
 585 {
 586         uint i;
 587         BSP_CSGMesh *outputMesh = BSP_CSGMesh::New();
 588
 589         if (outputMesh == NULL)
 590                 return NULL;
 591
 592         std::vector<BSP_MVertex> *vertices = new std::vector<BSP_MVertex>;
 593
 594         outputMesh->SetVertices(vertices);
 595
 596         std::map<MeshSet<3>::vertex_t*, uint> vertexToIndex_map;
 597         std::vector<MeshSet<3>::vertex_t>::iterator it = poly->vertex_storage.begin();
 598         for (i = 0; it != poly->vertex_storage.end(); ++i, ++it) {
 599                 MeshSet<3>::vertex_t *vertex = &(*it);
 600                 vertexToIndex_map[vertex] = i;
 601         }
 602
 603         for (i = 0; i < poly->vertex_storage.size(); ++i ) {
 604                 BSP_MVertex outVtx(MT_Point3 (poly->vertex_storage[i].v[0],
 605                                               poly->vertex_storage[i].v[1],
 606                                               poly->vertex_storage[i].v[2]));
 607                 outVtx.m_edges.clear();
 608                 outputMesh->VertexSet().push_back(outVtx);
 609         }
 610
 611         // build vectors of faces for each original face and each vertex
 612         std::vector<std::vector<uint> > vi(poly->vertex_storage.size());
 613         std::vector<std::vector<uint> > ofaces(num_origfaces);
 614         MeshSet<3>::face_iter face_iter = poly->faceBegin();
 615         for (i = 0; face_iter != poly->faceEnd(); ++face_iter, ++i) {
 616                 MeshSet<3>::face_t *f = *face_iter;
 617
 618                 if (Carve_checkDegeneratedFace(&vertexToIndex_map, f))
 619                         continue;
 620
 621                 ofaces[oface_num.getAttribute(f)].push_back(i);
 622
 623                 MeshSet<3>::face_t::edge_iter_t edge_iter = f->begin();
 624
 625                 for (; edge_iter != f->end(); ++edge_iter) {
 626                         int index = vertexToIndex_map[edge_iter->vert];
 627                         vi[index].push_back(i);
 628                 }
 629         }
 630
 631         uint quadverts[4] = {0, 0, 0, 0};
 632         // go over each set of faces which belong to an original face
 633         std::vector< std::vector<uint> >::const_iterator fii;
 634         uint orig = 0;
 635         for (fii=ofaces.begin(); fii!=ofaces.end(); ++fii, ++orig) {
 636                 std::vector<uint> fl = *fii;
 637                 // go over a single set from an original face
 638                 while (fl.size() > 0) {
 639                         // remove one new face
 640                         uint findex = fl.back();
 641                         fl.pop_back();
 642
 643                         MeshSet<3>::face_t *f = *(poly->faceBegin() + findex);
 644
 645                         // for each vertex of this face, check other faces containing
 646                         // that vertex to see if there is a neighbor also belonging to
 647                         // the original face
 648                         uint result = 0;
 649
 650                         MeshSet<3>::face_t::edge_iter_t edge_iter = f->begin();
 651                         for (; edge_iter != f->end(); ++edge_iter) {
 652                                 int v = vertexToIndex_map[edge_iter->vert];
 653                                 for (uint pos2=0; !result && pos2 < vi[v].size();pos2++) {
 654
 655                                         // if we find the current face, ignore it
 656                                         uint otherf = vi[v][pos2];
 657                                         if (findex == otherf)
 658                                                 continue;
 659
 660                                         MeshSet<3>::face_t *f2 = *(poly->faceBegin() + otherf);
 661
 662                                         // if other face doesn't have the same original face,
 663                                         // ignore it also
 664                                         uint other_orig = oface_num.getAttribute(f2);
 665                                         if (orig != other_orig)
 666                                                 continue;
 667
 668                                         // if, for some reason, we don't find the other face in
 669                                         // the current set of faces, ignore it
 670                                         uint other_index = 0;
 671                                         while (other_index < fl.size() && fl[other_index] != otherf) ++other_index;
 672                                         if (other_index == fl.size()) continue;
 673
 674                                         // see if the faces share an edge
 675                                         result = quadMerge(&vertexToIndex_map, poly->vertex_storage, f, f2, v, quadverts);
 676                                         // if faces can be merged, then remove the other face
 677                                         // from the current set
 678                                         if (result) {
 679                                                 uint replace = fl.back();
 680                                                 fl.pop_back();
 681                                                 if(otherf != replace)
 682                                                         fl[other_index] = replace;
 683                                         }
 684                                 }
 685                         }
 686
 687                         // add all information except vertices to the output mesh
 688                         outputMesh->FaceSet().push_back(BSP_MFace());
 689                         BSP_MFace& outFace = outputMesh->FaceSet().back();
 690                         outFace.m_verts.clear();
 691                         outFace.m_plane.setValue(f->plane.N.v);
 692                         outFace.m_orig_face = orig;
 693
 694                         // if we merged faces, use the list of common vertices; otherwise
 695                         // use the faces's vertices
 696                         if (result) {
 697                                 // make quat using verts stored in result
 698                                 outFace.m_verts.push_back(quadverts[0]);
 699                                 outFace.m_verts.push_back(quadverts[1]);
 700                                 outFace.m_verts.push_back(quadverts[2]);
 701                                 outFace.m_verts.push_back(quadverts[3]);
 702                         } else {
 703                                 MeshSet<3>::face_t::edge_iter_t edge_iter = f->begin();
 704                                 for (; edge_iter != f->end(); ++edge_iter) {
 705                                         //int index = ofacevert_num.getAttribute(f, edge_iter.idx());
 706                                         int index = vertexToIndex_map[edge_iter->vert];
 707                                         outFace.m_verts.push_back( index );
 708                                 }
 709                         }
 710                 }
 711         }
 712
 713         // Build the mesh edges using topological informtion
 714         outputMesh->BuildEdges();
 715
 716         return outputMesh;
 717 }
 718
 719 /**
 720  * Performs a generic booleam operation, the entry point for external modules.
 721  * @param opType Boolean operation type BOP_INTERSECTION, BOP_UNION, BOP_DIFFERENCE
 722  * @param outputMesh Output mesh, the final result (the object C)
 723  * @param obAFaces Object A faces list
 724  * @param obAVertices Object A vertices list
 725  * @param obBFaces Object B faces list
 726  * @param obBVertices Object B vertices list
 727  * @param interpFunc Interpolating function
 728  * @return operation state: BOP_OK, BOP_NO_SOLID, BOP_ERROR
 729  */
 730 BoolOpState BOP_performBooleanOperation(BoolOpType                    opType,
 731                                         BSP_CSGMesh**                 outputMesh,
 732                                         CSG_FaceIteratorDescriptor    obAFaces,
 733                                         CSG_VertexIteratorDescriptor  obAVertices,
 734                                         CSG_FaceIteratorDescriptor    obBFaces,
 735                                         CSG_VertexIteratorDescriptor  obBVertices)
 736 {
 737         carve::csg::CSG::OP op;
 738         MeshSet<3> *left, *right, *output = NULL;
 739         carve::csg::CSG csg;
 740         carve::geom3d::Vector min, max;
 741         carve::interpolate::FaceAttr<uint> oface_num;
 742         uint num_origfaces = 0;
 743
 744         switch (opType) {
 745                 case BOP_UNION:
 746                         op = carve::csg::CSG::UNION;
 747                         break;
 748                 case BOP_INTERSECTION:
 749                         op = carve::csg::CSG::INTERSECTION;
 750                         break;
 751                 case BOP_DIFFERENCE:
 752                         op = carve::csg::CSG::A_MINUS_B;
 753                         break;
 754                 default:
 755                         return BOP_ERROR;
 756         }
 757
 758         left = Carve_addMesh(obAFaces, obAVertices, oface_num, num_origfaces );
 759         right = Carve_addMesh(obBFaces, obBVertices, oface_num, num_origfaces );
 760
 761         min.x = max.x = left->vertex_storage[0].v.x;
 762         min.y = max.y = left->vertex_storage[0].v.y;
 763         min.z = max.z = left->vertex_storage[0].v.z;
 764         for (uint i = 1; i < left->vertex_storage.size(); ++i) {
 765                 min.x = MIN(min.x,left->vertex_storage[i].v.x);
 766                 min.y = MIN(min.y,left->vertex_storage[i].v.y);
 767                 min.z = MIN(min.z,left->vertex_storage[i].v.z);
 768                 max.x = MAX(max.x,left->vertex_storage[i].v.x);
 769                 max.y = MAX(max.y,left->vertex_storage[i].v.y);
 770                 max.z = MAX(max.z,left->vertex_storage[i].v.z);
 771         }
 772         for (uint i = 0; i < right->vertex_storage.size(); ++i) {
 773                 min.x = MIN(min.x,right->vertex_storage[i].v.x);
 774                 min.y = MIN(min.y,right->vertex_storage[i].v.y);
 775                 min.z = MIN(min.z,right->vertex_storage[i].v.z);
 776                 max.x = MAX(max.x,right->vertex_storage[i].v.x);
 777                 max.y = MAX(max.y,right->vertex_storage[i].v.y);
 778                 max.z = MAX(max.z,right->vertex_storage[i].v.z);
 779         }
 780
 781         carve::rescale::rescale scaler(min.x, min.y, min.z, max.x, max.y, max.z);
 782         carve::rescale::fwd fwd_r(scaler);
 783         carve::rescale::rev rev_r(scaler);
 784
 785         left->transform(fwd_r);
 786         right->transform(fwd_r);
 787
 788         // prepare operands for actual boolean operation. it's needed because operands might consist of
 789         // several intersecting meshes and in case if another operands intersect an edge loop of intersecting that
 790         // meshes tessellation of operation result can't be done properly. the only way to make such situations
 791         // working is to union intersecting meshes of the same operand
 792         Carve_unionIntersections(&left, &right, oface_num);
 793
 794         if(left->meshes.size() == 0 || right->meshes.size()==0) {
 795                 // normally sohuldn't happen (zero-faces objects are handled by modifier itself), but
 796                 // unioning intersecting meshes which doesn't have consistent normals might lead to
 797                 // empty result which wouldn't work here
 798
 799                 delete left;
 800                 delete right;
 801
 802                 return BOP_ERROR;
 803         }
 804
 805         csg.hooks.registerHook(new carve::csg::CarveTriangulator, carve::csg::CSG::Hooks::PROCESS_OUTPUT_FACE_BIT);
 806
 807         oface_num.installHooks(csg);
 808
 809         try {
 810                 output = csg.compute(left, right, op, NULL, carve::csg::CSG::CLASSIFY_EDGE);
 811         }
 812         catch(carve::exception e) {
 813                 std::cerr << "CSG failed, exception " << e.str() << std::endl;
 814         }
 815         catch(...) {
 816                 delete left;
 817                 delete right;
 818
 819                 throw "Unknown error in Carve library";
 820         }
 821
 822         delete left;
 823         delete right;
 824
 825         if(!output)
 826                 return BOP_ERROR;
 827
 828         output->transform(rev_r);
 829
 830         *outputMesh = Carve_exportMesh(output, oface_num, num_origfaces);
 831         delete output;
 832
 833         return BOP_OK;
 834 }