source/blender/bmesh/tools/bmesh_decimate_collapse.c

   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  * Contributor(s): Campbell Barton
  19  *
  20  * ***** END GPL LICENSE BLOCK *****
  21  */
  22
  23 /** \file blender/bmesh/tools/bmesh_decimate_collapse.c
  24  *  \ingroup bmesh
  25  *
  26  * BMesh decimator that uses an edge collapse method.
  27  */
  28
  29 #include <stddef.h>
  30
  31 #include "MEM_guardedalloc.h"
  32
  33 #include "DNA_scene_types.h"
  34
  35 #include "BLI_math.h"
  36 #include "BLI_quadric.h"
  37 #include "BLI_heap.h"
  38
  39 #include "BKE_customdata.h"
  40
  41 #include "bmesh.h"
  42 #include "bmesh_decimate.h"  /* own include */
  43
  44 #include "../intern/bmesh_structure.h"
  45
  46 /* defines for testing */
  47 #define USE_CUSTOMDATA
  48 #define USE_TRIANGULATE
  49 #define USE_VERT_NORMAL_INTERP  /* has the advantage that flipped faces don't mess up vertex normals */
  50
  51 /* these checks are for rare cases that we can't avoid since they are valid meshes still */
  52 #define USE_SAFETY_CHECKS
  53
  54 #define BOUNDARY_PRESERVE_WEIGHT 100.0f
  55 #define OPTIMIZE_EPS 0.01f  /* FLT_EPSILON is too small, see [#33106] */
  56 #define COST_INVALID FLT_MAX
  57
  58 typedef enum CD_UseFlag {
  59         CD_DO_VERT = (1 << 0),
  60         CD_DO_EDGE = (1 << 1),
  61         CD_DO_LOOP = (1 << 2)
  62 } CD_UseFlag;
  63
  64
  65 /* BMesh Helper Functions
  66  * ********************** */
  67
  68 /**
  69  * \param vquadrics must be calloc'd
  70  */
  71 static void bm_decim_build_quadrics(BMesh *bm, Quadric *vquadrics)
  72 {
  73         BMIter iter;
  74         BMFace *f;
  75         BMEdge *e;
  76
  77         BM_ITER_MESH (f, &iter, bm, BM_FACES_OF_MESH) {
  78                 BMLoop *l_first;
  79                 BMLoop *l_iter;
  80
  81                 const float *co = BM_FACE_FIRST_LOOP(f)->v->co;
  82                 const float *no = f->no;
  83                 const float offset = -dot_v3v3(no, co);
  84                 Quadric q;
  85
  86                 BLI_quadric_from_v3_dist(&q, no, offset);
  87
  88                 l_iter = l_first = BM_FACE_FIRST_LOOP(f);
  89                 do {
  90                         BLI_quadric_add_qu_qu(&vquadrics[BM_elem_index_get(l_iter->v)], &q);
  91                 } while ((l_iter = l_iter->next) != l_first);
  92         }
  93
  94         /* boundary edges */
  95         BM_ITER_MESH (e, &iter, bm, BM_EDGES_OF_MESH) {
  96                 if (UNLIKELY(BM_edge_is_boundary(e))) {
  97                         float edge_vector[3];
  98                         float edge_cross[3];
  99                         sub_v3_v3v3(edge_vector, e->v2->co, e->v1->co);
 100                         f = e->l->f;
 101                         cross_v3_v3v3(edge_cross, edge_vector, f->no);
 102
 103                         if (normalize_v3(edge_cross) > FLT_EPSILON) {
 104                                 Quadric q;
 105                                 BLI_quadric_from_v3_dist(&q, edge_cross, -dot_v3v3(edge_cross, e->v1->co));
 106                                 BLI_quadric_mul(&q, BOUNDARY_PRESERVE_WEIGHT);
 107
 108                                 BLI_quadric_add_qu_qu(&vquadrics[BM_elem_index_get(e->v1)], &q);
 109                                 BLI_quadric_add_qu_qu(&vquadrics[BM_elem_index_get(e->v2)], &q);
 110                         }
 111                 }
 112         }
 113 }
 114
 115
 116 static void bm_decim_calc_target_co(BMEdge *e, float optimize_co[3],
 117                                     const Quadric *vquadrics)
 118 {
 119         /* compute an edge contraction target for edge 'e'
 120          * this is computed by summing it's vertices quadrics and
 121          * optimizing the result. */
 122         Quadric q;
 123
 124         BLI_quadric_add_qu_ququ(&q,
 125                                 &vquadrics[BM_elem_index_get(e->v1)],
 126                                 &vquadrics[BM_elem_index_get(e->v2)]);
 127
 128
 129         if (BLI_quadric_optimize(&q, optimize_co, OPTIMIZE_EPS)) {
 130                 return;  /* all is good */
 131         }
 132         else {
 133                 mid_v3_v3v3(optimize_co, e->v1->co, e->v2->co);
 134         }
 135 }
 136
 137 static bool bm_edge_collapse_is_degenerate_flip(BMEdge *e, const float optimize_co[3])
 138 {
 139         BMIter liter;
 140         BMLoop *l;
 141         unsigned int i;
 142
 143         for (i = 0; i < 2; i++) {
 144                 /* loop over both verts */
 145                 BMVert *v = *((&e->v1) + i);
 146
 147                 BM_ITER_ELEM (l, &liter, v, BM_LOOPS_OF_VERT) {
 148                         if (l->e != e && l->prev->e != e) {
 149                                 float *co_prev = l->prev->v->co;
 150                                 float *co_next = l->next->v->co;
 151                                 float cross_exist[3];
 152                                 float cross_optim[3];
 153
 154 #if 1
 155                                 float vec_other[3];  /* line between the two outer verts, re-use for both cross products */
 156                                 float vec_exist[3];  /* before collapse */
 157                                 float vec_optim[3];  /* after collapse */
 158
 159                                 sub_v3_v3v3(vec_other, co_prev, co_next);
 160                                 sub_v3_v3v3(vec_exist, co_prev, v->co);
 161                                 sub_v3_v3v3(vec_optim, co_prev, optimize_co);
 162
 163                                 cross_v3_v3v3(cross_exist, vec_other, vec_exist);
 164                                 cross_v3_v3v3(cross_optim, vec_other, vec_optim);
 165
 166                                 /* normalize isn't really needed, but ensures the value at a unit we can compare against */
 167                                 normalize_v3(cross_exist);
 168                                 normalize_v3(cross_optim);
 169 #else
 170                                 normal_tri_v3(cross_exist, v->co,       co_prev, co_next);
 171                                 normal_tri_v3(cross_optim, optimize_co, co_prev, co_next);
 172 #endif
 173
 174                                 /* use a small value rather then zero so we don't flip a face in multiple steps
 175                                  * (first making it zero area, then flipping again) */
 176                                 if (dot_v3v3(cross_exist, cross_optim) <= FLT_EPSILON) {
 177                                         //printf("no flip\n");
 178                                         return true;
 179                                 }
 180                         }
 181                 }
 182         }
 183
 184         return false;
 185 }
 186
 187 static void bm_decim_build_edge_cost_single(BMEdge *e,
 188                                             const Quadric *vquadrics, const float *vweights,
 189                                             Heap *eheap, HeapNode **eheap_table)
 190 {
 191         const Quadric *q1, *q2;
 192         float optimize_co[3];
 193         float cost;
 194
 195         if (eheap_table[BM_elem_index_get(e)]) {
 196                 BLI_heap_remove(eheap, eheap_table[BM_elem_index_get(e)]);
 197         }
 198
 199         /* check we can collapse, some edges we better not touch */
 200         if (BM_edge_is_boundary(e)) {
 201                 if (e->l->f->len == 3) {
 202                         /* pass */
 203                 }
 204                 else {
 205                         /* only collapse tri's */
 206                         eheap_table[BM_elem_index_get(e)] = NULL;
 207                         return;
 208                 }
 209         }
 210         else if (BM_edge_is_manifold(e)) {
 211                 if ((e->l->f->len == 3) && (e->l->radial_next->f->len == 3)) {
 212                         /* pass */
 213                 }
 214                 else {
 215                         /* only collapse tri's */
 216                         eheap_table[BM_elem_index_get(e)] = NULL;
 217                         return;
 218                 }
 219         }
 220         else {
 221                 eheap_table[BM_elem_index_get(e)] = NULL;
 222                 return;
 223         }
 224
 225         if (vweights) {
 226                 if ((vweights[BM_elem_index_get(e->v1)] >= BM_MESH_DECIM_WEIGHT_MAX) &&
 227                     (vweights[BM_elem_index_get(e->v2)] >= BM_MESH_DECIM_WEIGHT_MAX))
 228                 {
 229                         /* skip collapsing this edge */
 230                         eheap_table[BM_elem_index_get(e)] = NULL;
 231                         return;
 232                 }
 233         }
 234         /* end sanity check */
 235
 236
 237         bm_decim_calc_target_co(e, optimize_co, vquadrics);
 238
 239         q1 = &vquadrics[BM_elem_index_get(e->v1)];
 240         q2 = &vquadrics[BM_elem_index_get(e->v2)];
 241
 242         if (vweights == NULL) {
 243                 cost = (BLI_quadric_evaluate(q1, optimize_co) +
 244                         BLI_quadric_evaluate(q2, optimize_co));
 245         }
 246         else {
 247                 /* add 1.0 so planar edges are still weighted against */
 248                 cost = (((BLI_quadric_evaluate(q1, optimize_co) + 1.0f) * vweights[BM_elem_index_get(e->v1)]) +
 249                         ((BLI_quadric_evaluate(q2, optimize_co) + 1.0f) * vweights[BM_elem_index_get(e->v2)]));
 250         }
 251         // print("COST %.12f\n");
 252
 253         /* note, 'cost' shouldn't be negative but happens sometimes with small values.
 254          * this can cause faces that make up a flat surface to over-collapse, see [#37121] */
 255         cost = fabsf(cost);
 256         eheap_table[BM_elem_index_get(e)] = BLI_heap_insert(eheap, cost, e);
 257 }
 258
 259
 260 /* use this for degenerate cases - add back to the heap with an invalid cost,
 261  * this way it may be calculated again if surrounding geometry changes */
 262 static void bm_decim_invalid_edge_cost_single(BMEdge *e,
 263                                               Heap *eheap, HeapNode **eheap_table)
 264 {
 265         BLI_assert(eheap_table[BM_elem_index_get(e)] == NULL);
 266         eheap_table[BM_elem_index_get(e)] = BLI_heap_insert(eheap, COST_INVALID, e);
 267 }
 268
 269 static void bm_decim_build_edge_cost(BMesh *bm,
 270                                      const Quadric *vquadrics, const float *vweights,
 271                                      Heap *eheap, HeapNode **eheap_table)
 272 {
 273         BMIter iter;
 274         BMEdge *e;
 275         unsigned int i;
 276
 277         BM_ITER_MESH_INDEX (e, &iter, bm, BM_EDGES_OF_MESH, i) {
 278                 eheap_table[i] = NULL;  /* keep sanity check happy */
 279                 bm_decim_build_edge_cost_single(e, vquadrics, vweights, eheap, eheap_table);
 280         }
 281 }
 282
 283 #ifdef USE_TRIANGULATE
 284 /* Temp Triangulation
 285  * ****************** */
 286
 287 /**
 288  * To keep things simple we can only collapse edges on triangulated data
 289  * (limitation with edge collapse and error calculation functions).
 290  *
 291  * But to avoid annoying users by only giving triangle results, we can
 292  * triangulate, keeping a reference between the faces, then join after
 293  * if the edges don't collapse, this will also allow more choices when
 294  * collapsing edges so even has some advantage over decimating quads
 295  * directly.
 296  *
 297  * \return true if any faces were triangulated.
 298  */
 299
 300 static bool bm_decim_triangulate_begin(BMesh *bm)
 301 {
 302         BMIter iter;
 303         BMFace *f;
 304         // bool has_quad;  // could optimize this a little
 305         bool has_cut = false;
 306
 307         BLI_assert((bm->elem_index_dirty & BM_VERT) == 0);
 308
 309         /* first clear loop index values */
 310         BM_ITER_MESH (f, &iter, bm, BM_FACES_OF_MESH) {
 311                 BMLoop *l_iter;
 312                 BMLoop *l_first;
 313
 314                 l_iter = l_first = BM_FACE_FIRST_LOOP(f);
 315                 do {
 316                         BM_elem_index_set(l_iter, -1);
 317                 } while ((l_iter = l_iter->next) != l_first);
 318
 319                 // has_quad |= (f->len == 4)
 320         }
 321
 322         /* adding new faces as we loop over faces
 323          * is normally best avoided, however in this case its not so bad because any face touched twice
 324          * will already be triangulated*/
 325         BM_ITER_MESH (f, &iter, bm, BM_FACES_OF_MESH) {
 326                 if (f->len == 4) {
 327                         BMLoop *f_l[4];
 328                         BMLoop *l_a, *l_b;
 329
 330                         {
 331                                 BMLoop *l_iter = BM_FACE_FIRST_LOOP(f);
 332
 333                                 f_l[0] = l_iter; l_iter = l_iter->next;
 334                                 f_l[1] = l_iter; l_iter = l_iter->next;
 335                                 f_l[2] = l_iter; l_iter = l_iter->next;
 336                                 f_l[3] = l_iter;
 337                         }
 338
 339                         if (len_squared_v3v3(f_l[0]->v->co, f_l[2]->v->co) <
 340                             len_squared_v3v3(f_l[1]->v->co, f_l[3]->v->co))
 341                         {
 342                                 l_a = f_l[0];
 343                                 l_b = f_l[2];
 344                         }
 345                         else {
 346                                 l_a = f_l[1];
 347                                 l_b = f_l[3];
 348                         }
 349
 350 #ifdef USE_SAFETY_CHECKS
 351                         if (BM_edge_exists(l_a->v, l_b->v) == NULL)
 352 #endif
 353                         {
 354                                 BMFace *f_new;
 355                                 BMLoop *l_new;
 356
 357                                 /* warning, NO_DOUBLE option here isn't handled as nice as it could be
 358                                  * - if there is a quad that has a free standing edge joining it along
 359                                  * where we want to split the face, there isnt a good way we can handle this.
 360                                  * currently that edge will get removed when joining the tris back into a quad. */
 361                                 f_new = BM_face_split(bm, f, l_a, l_b, &l_new, NULL, false);
 362
 363                                 if (f_new) {
 364                                         /* the value of this doesn't matter, only that the 2 loops match and have unique values */
 365                                         const int f_index = BM_elem_index_get(f);
 366
 367                                         /* since we just split theres only ever 2 loops */
 368                                         BLI_assert(BM_edge_is_manifold(l_new->e));
 369
 370                                         BM_elem_index_set(l_new, f_index);
 371                                         BM_elem_index_set(l_new->radial_next, f_index);
 372
 373                                         BM_face_normal_update(f);
 374                                         BM_face_normal_update(f_new);
 375
 376                                         has_cut = true;
 377                                 }
 378                         }
 379                 }
 380         }
 381
 382         BLI_assert((bm->elem_index_dirty & BM_VERT) == 0);
 383
 384         if (has_cut) {
 385                 /* now triangulation is done we need to correct index values */
 386                 BM_mesh_elem_index_ensure(bm, BM_EDGE | BM_FACE);
 387         }
 388
 389         return has_cut;
 390 }
 391
 392 static void bm_decim_triangulate_end(BMesh *bm)
 393 {
 394         /* decimation finished, now re-join */
 395         BMIter iter;
 396         BMEdge *e, *e_next;
 397
 398         /* boundary edges */
 399         BM_ITER_MESH_MUTABLE (e, e_next, &iter, bm, BM_EDGES_OF_MESH) {
 400                 BMLoop *l_a, *l_b;
 401                 if (BM_edge_loop_pair(e, &l_a, &l_b)) {
 402                         const int l_a_index = BM_elem_index_get(l_a);
 403                         if (l_a_index != -1) {
 404                                 const int l_b_index = BM_elem_index_get(l_b);
 405                                 if (l_a_index == l_b_index) {
 406                                         if (LIKELY(l_a->f->len == 3 && l_b->f->len == 3)) {
 407                                                 if (l_a->v != l_b->v) {  /* if this is the case, faces have become flipped */
 408                                                         /* check we are not making a degenerate quad */
 409                                                         BMVert *vquad[4] = {
 410                                                                 e->v1,
 411                                                                 BM_vert_in_edge(e, l_a->next->v) ? l_a->prev->v : l_a->next->v,
 412                                                                 e->v2,
 413                                                                 BM_vert_in_edge(e, l_b->next->v) ? l_b->prev->v : l_b->next->v,
 414                                                         };
 415
 416                                                         BLI_assert(ELEM3(vquad[0], vquad[1], vquad[2], vquad[3]) == false);
 417                                                         BLI_assert(ELEM3(vquad[1], vquad[0], vquad[2], vquad[3]) == false);
 418                                                         BLI_assert(ELEM3(vquad[2], vquad[1], vquad[0], vquad[3]) == false);
 419                                                         BLI_assert(ELEM3(vquad[3], vquad[1], vquad[2], vquad[0]) == false);
 420
 421                                                         if (is_quad_convex_v3(vquad[0]->co, vquad[1]->co, vquad[2]->co, vquad[3]->co)) {
 422                                                                 /* highly unlikely to fail, but prevents possible double-ups */
 423                                                                 BMFace *f[2] = {l_a->f, l_b->f};
 424                                                                 BM_faces_join(bm, f, 2, true);
 425                                                         }
 426                                                 }
 427                                         }
 428                                 }
 429                         }
 430                 }
 431         }
 432 }
 433
 434 #endif  /* USE_TRIANGULATE */
 435
 436 /* Edge Collapse Functions
 437  * *********************** */
 438
 439 #ifdef USE_CUSTOMDATA
 440
 441 /**
 442  * \param v is the target to merge into.
 443  */
 444 static void bm_edge_collapse_loop_customdata(BMesh *bm, BMLoop *l, BMVert *v_clear, BMVert *v_other,
 445                                              const float customdata_fac)
 446 {
 447         /* disable seam check - the seam check would have to be done per layer, its not really that important */
 448 //#define USE_SEAM
 449         /* these don't need to be updated, since they will get removed when the edge collapses */
 450         BMLoop *l_clear, *l_other;
 451         const bool is_manifold = BM_edge_is_manifold(l->e);
 452         int side;
 453
 454         /* l defines the vert to collapse into  */
 455
 456         /* first find the loop of 'v_other' thats attached to the face of 'l' */
 457         if (l->v == v_clear) {
 458                 l_clear = l;
 459                 l_other = l->next;
 460         }
 461         else {
 462                 l_clear = l->next;
 463                 l_other = l;
 464         }
 465
 466         BLI_assert(l_clear->v == v_clear);
 467         BLI_assert(l_other->v == v_other);
 468         (void)v_other;  /* quiet warnings for release */
 469
 470         /* now we have both corners of the face 'l->f' */
 471         for (side = 0; side < 2; side++) {
 472 #ifdef USE_SEAM
 473                 bool is_seam = false;
 474 #endif
 475                 void *src[2];
 476                 BMFace *f_exit = is_manifold ? l->radial_next->f : NULL;
 477                 BMEdge *e_prev = l->e;
 478                 BMLoop *l_first;
 479                 BMLoop *l_iter;
 480                 float w[2];
 481
 482                 if (side == 0) {
 483                         l_iter = l_first = l_clear;
 484                         src[0] = l_clear->head.data;
 485                         src[1] = l_other->head.data;
 486
 487                         w[0] = customdata_fac;
 488                         w[1] = 1.0f - customdata_fac;
 489                 }
 490                 else {
 491                         l_iter = l_first = l_other;
 492                         src[0] = l_other->head.data;
 493                         src[1] = l_clear->head.data;
 494
 495                         w[0] = 1.0f - customdata_fac;
 496                         w[1] = customdata_fac;
 497                 }
 498
 499                 // print_v2("weights", w);
 500
 501                 /* WATCH IT! - should NOT reference (_clear or _other) vars for this while loop */
 502
 503                 /* walk around the fan using 'e_prev' */
 504                 while (((l_iter = BM_vert_step_fan_loop(l_iter, &e_prev)) != l_first) && (l_iter != NULL)) {
 505                         int i;
 506                         /* quit once we hit the opposite face, if we have one */
 507                         if (f_exit && UNLIKELY(f_exit == l_iter->f)) {
 508                                 break;
 509                         }
 510
 511 #ifdef USE_SEAM
 512                         /* break out unless we find a match */
 513                         is_seam = true;
 514 #endif
 515
 516                         /* ok. we have a loop. now be smart with it! */
 517                         for (i = 0; i < bm->ldata.totlayer; i++) {
 518                                 if (CustomData_layer_has_math(&bm->ldata, i)) {
 519                                         const int offset = bm->ldata.layers[i].offset;
 520                                         const int type = bm->ldata.layers[i].type;
 521                                         void *cd_src[2] = {(char *)src[0] + offset,
 522                                                            (char *)src[1] + offset};
 523                                         void *cd_iter = (char *)l_iter->head.data + offset;
 524
 525                                         /* detect seams */
 526                                         if (CustomData_data_equals(type, cd_src[0], cd_iter)) {
 527                                                 CustomData_bmesh_interp_n(&bm->ldata, cd_src, w, NULL, 2, l_iter->head.data, i);
 528 #ifdef USE_SEAM
 529                                                 is_seam = false;
 530 #endif
 531                                         }
 532                                 }
 533                         }
 534
 535 #ifdef USE_SEAM
 536                         if (is_seam) {
 537                                 break;
 538                         }
 539 #endif
 540                 }
 541         }
 542
 543 //#undef USE_SEAM
 544
 545 }
 546 #endif  /* USE_CUSTOMDATA */
 547
 548 /**
 549  * Check if the collapse will result in a degenerate mesh,
 550  * that is - duplicate edges or faces.
 551  *
 552  * This situation could be checked for when calculating collapse cost
 553  * however its quite slow and a degenerate collapse could eventuate
 554  * after the cost is calculated, so instead, check just before collapsing.
 555  */
 556
 557 static void bm_edge_tag_enable(BMEdge *e)
 558 {
 559         BM_elem_flag_enable(e->v1, BM_ELEM_TAG);
 560         BM_elem_flag_enable(e->v2, BM_ELEM_TAG);
 561         if (e->l) {
 562                 BM_elem_flag_enable(e->l->f, BM_ELEM_TAG);
 563                 if (e->l != e->l->radial_next) {
 564                         BM_elem_flag_enable(e->l->radial_next->f, BM_ELEM_TAG);
 565                 }
 566         }
 567 }
 568
 569 static void bm_edge_tag_disable(BMEdge *e)
 570 {
 571         BM_elem_flag_disable(e->v1, BM_ELEM_TAG);
 572         BM_elem_flag_disable(e->v2, BM_ELEM_TAG);
 573         if (e->l) {
 574                 BM_elem_flag_disable(e->l->f, BM_ELEM_TAG);
 575                 if (e->l != e->l->radial_next) {
 576                         BM_elem_flag_disable(e->l->radial_next->f, BM_ELEM_TAG);
 577                 }
 578         }
 579 }
 580
 581 static int bm_edge_tag_test(BMEdge *e)
 582 {
 583         /* is the edge or one of its faces tagged? */
 584         return (BM_elem_flag_test(e->v1, BM_ELEM_TAG) ||
 585                 BM_elem_flag_test(e->v2, BM_ELEM_TAG) ||
 586                 (e->l && (BM_elem_flag_test(e->l->f, BM_ELEM_TAG) ||
 587                           (e->l != e->l->radial_next &&
 588                           BM_elem_flag_test(e->l->radial_next->f, BM_ELEM_TAG))))
 589                 );
 590 }
 591
 592 /* takes the edges loop */
 593 BLI_INLINE int bm_edge_is_manifold_or_boundary(BMLoop *l)
 594 {
 595 #if 0
 596         /* less optimized version of check below */
 597         return (BM_edge_is_manifold(l->e) || BM_edge_is_boundary(l->e);
 598 #else
 599         /* if the edge is a boundary it points to its self, else this must be a manifold */
 600         return LIKELY(l) && LIKELY(l->radial_next->radial_next == l);
 601 #endif
 602 }
 603
 604 static bool bm_edge_collapse_is_degenerate_topology(BMEdge *e_first)
 605 {
 606         /* simply check that there is no overlap between faces and edges of each vert,
 607          * (excluding the 2 faces attached to 'e' and 'e' its self) */
 608
 609         BMEdge *e_iter;
 610
 611         /* clear flags on both disks */
 612         e_iter = e_first;
 613         do {
 614                 if (!bm_edge_is_manifold_or_boundary(e_iter->l)) {
 615                         return true;
 616                 }
 617                 bm_edge_tag_disable(e_iter);
 618         } while ((e_iter = bmesh_disk_edge_next(e_iter, e_first->v1)) != e_first);
 619
 620         e_iter = e_first;
 621         do {
 622                 if (!bm_edge_is_manifold_or_boundary(e_iter->l)) {
 623                         return true;
 624                 }
 625                 bm_edge_tag_disable(e_iter);
 626         } while ((e_iter = bmesh_disk_edge_next(e_iter, e_first->v2)) != e_first);
 627
 628         /* now enable one side... */
 629         e_iter = e_first;
 630         do {
 631                 bm_edge_tag_enable(e_iter);
 632         } while ((e_iter = bmesh_disk_edge_next(e_iter, e_first->v1)) != e_first);
 633
 634         /* ... except for the edge we will collapse, we know thats shared,
 635          * disable this to avoid false positive. We could be smart and never enable these
 636          * face/edge tags in the first place but easier to do this */
 637         // bm_edge_tag_disable(e_first);
 638         /* do inline... */
 639         {
 640 #if 0
 641                 BMIter iter;
 642                 BMIter liter;
 643                 BMLoop *l;
 644                 BMVert *v;
 645                 BM_ITER_ELEM (l, &liter, e_first, BM_LOOPS_OF_EDGE) {
 646                         BM_elem_flag_disable(l->f, BM_ELEM_TAG);
 647                         BM_ITER_ELEM (v, &iter, l->f, BM_VERTS_OF_FACE) {
 648                                 BM_elem_flag_disable(v, BM_ELEM_TAG);
 649                         }
 650                 }
 651 #else
 652                 /* we know each face is a triangle, no looping/iterators needed here */
 653
 654                 BMLoop *l_radial;
 655                 BMLoop *l_face;
 656
 657                 l_radial = e_first->l;
 658                 l_face = l_radial;
 659                 BLI_assert(l_face->f->len == 3);
 660                 BM_elem_flag_disable(l_face->f, BM_ELEM_TAG);
 661                 BM_elem_flag_disable((l_face = l_radial)->v,     BM_ELEM_TAG);
 662                 BM_elem_flag_disable((l_face = l_face->next)->v, BM_ELEM_TAG);
 663                 BM_elem_flag_disable((         l_face->next)->v, BM_ELEM_TAG);
 664                 l_face = l_radial->radial_next;
 665                 if (l_radial != l_face) {
 666                         BLI_assert(l_face->f->len == 3);
 667                         BM_elem_flag_disable(l_face->f, BM_ELEM_TAG);
 668                         BM_elem_flag_disable((l_face = l_radial->radial_next)->v, BM_ELEM_TAG);
 669                         BM_elem_flag_disable((l_face = l_face->next)->v,          BM_ELEM_TAG);
 670                         BM_elem_flag_disable((         l_face->next)->v,          BM_ELEM_TAG);
 671                 }
 672 #endif
 673         }
 674
 675         /* and check for overlap */
 676         e_iter = e_first;
 677         do {
 678                 if (bm_edge_tag_test(e_iter)) {
 679                         return true;
 680                 }
 681         } while ((e_iter = bmesh_disk_edge_next(e_iter, e_first->v2)) != e_first);
 682
 683         return false;
 684 }
 685
 686 /**
 687  * special, highly limited edge collapse function
 688  * intended for speed over flexibility.
 689  * can only collapse edges connected to (1, 2) tris.
 690  *
 691  * Important - dont add vert/edge/face data on collapsing!
 692  *
 693  * \param e_clear_other let caller know what edges we remove besides \a e_clear
 694  * \param customdata_flag merge factor, scales from 0 - 1 ('v_clear' -> 'v_other')
 695  */
 696 static bool bm_edge_collapse(BMesh *bm, BMEdge *e_clear, BMVert *v_clear, int r_e_clear_other[2],
 697 #ifdef USE_CUSTOMDATA
 698                              const CD_UseFlag customdata_flag,
 699                              const float customdata_fac
 700 #else
 701                              const CD_UseFlag UNUSED(customdata_flag),
 702                              const float UNUSED(customdata_fac)
 703 #endif
 704                             )
 705 {
 706         BMVert *v_other;
 707
 708         v_other = BM_edge_other_vert(e_clear, v_clear);
 709         BLI_assert(v_other != NULL);
 710
 711         if (BM_edge_is_manifold(e_clear)) {
 712                 BMLoop *l_a, *l_b;
 713                 BMEdge *e_a_other[2], *e_b_other[2];
 714                 bool ok;
 715
 716                 ok = BM_edge_loop_pair(e_clear, &l_a, &l_b);
 717
 718                 BLI_assert(ok == true);
 719                 BLI_assert(l_a->f->len == 3);
 720                 BLI_assert(l_b->f->len == 3);
 721
 722                 /* keep 'v_clear' 0th */
 723                 if (BM_vert_in_edge(l_a->prev->e, v_clear)) {
 724                         e_a_other[0] = l_a->prev->e;
 725                         e_a_other[1] = l_a->next->e;
 726                 }
 727                 else {
 728                         e_a_other[1] = l_a->prev->e;
 729                         e_a_other[0] = l_a->next->e;
 730                 }
 731
 732                 if (BM_vert_in_edge(l_b->prev->e, v_clear)) {
 733                         e_b_other[0] = l_b->prev->e;
 734                         e_b_other[1] = l_b->next->e;
 735                 }
 736                 else {
 737                         e_b_other[1] = l_b->prev->e;
 738                         e_b_other[0] = l_b->next->e;
 739                 }
 740
 741                 /* we could assert this case, but better just bail out */
 742 #if 0
 743                 BLI_assert(e_a_other[0] != e_b_other[0]);
 744                 BLI_assert(e_a_other[0] != e_b_other[1]);
 745                 BLI_assert(e_b_other[0] != e_a_other[0]);
 746                 BLI_assert(e_b_other[0] != e_a_other[1]);
 747 #endif
 748                 /* not totally common but we want to avoid */
 749                 if (ELEM(e_a_other[0], e_b_other[0], e_b_other[1]) ||
 750                     ELEM(e_a_other[1], e_b_other[0], e_b_other[1]))
 751                 {
 752                         return false;
 753                 }
 754
 755                 BLI_assert(BM_edge_share_vert(e_a_other[0], e_b_other[0]));
 756                 BLI_assert(BM_edge_share_vert(e_a_other[1], e_b_other[1]));
 757
 758                 r_e_clear_other[0] = BM_elem_index_get(e_a_other[0]);
 759                 r_e_clear_other[1] = BM_elem_index_get(e_b_other[0]);
 760
 761 #ifdef USE_CUSTOMDATA
 762                 /* before killing, do customdata */
 763                 if (customdata_flag & CD_DO_VERT) {
 764                         BM_data_interp_from_verts(bm, v_other, v_clear, v_other, customdata_fac);
 765                 }
 766                 if (customdata_flag & CD_DO_EDGE) {
 767                         BM_data_interp_from_edges(bm, e_a_other[1], e_a_other[0], e_a_other[1], customdata_fac);
 768                         BM_data_interp_from_edges(bm, e_b_other[1], e_b_other[0], e_b_other[1], customdata_fac);
 769                 }
 770                 if (customdata_flag & CD_DO_LOOP) {
 771                         bm_edge_collapse_loop_customdata(bm, e_clear->l,              v_clear, v_other, customdata_fac);
 772                         bm_edge_collapse_loop_customdata(bm, e_clear->l->radial_next, v_clear, v_other, customdata_fac);
 773                 }
 774 #endif
 775
 776                 BM_edge_kill(bm, e_clear);
 777
 778                 v_other->head.hflag |= v_clear->head.hflag;
 779                 BM_vert_splice(bm, v_clear, v_other);
 780
 781                 e_a_other[1]->head.hflag |= e_a_other[0]->head.hflag;
 782                 e_b_other[1]->head.hflag |= e_b_other[0]->head.hflag;
 783                 BM_edge_splice(bm, e_a_other[0], e_a_other[1]);
 784                 BM_edge_splice(bm, e_b_other[0], e_b_other[1]);
 785
 786                 // BM_mesh_validate(bm);
 787
 788                 return true;
 789         }
 790         else if (BM_edge_is_boundary(e_clear)) {
 791                 /* same as above but only one triangle */
 792                 BMLoop *l_a;
 793                 BMEdge *e_a_other[2];
 794
 795                 l_a = e_clear->l;
 796
 797                 BLI_assert(l_a->f->len == 3);
 798
 799                 /* keep 'v_clear' 0th */
 800                 if (BM_vert_in_edge(l_a->prev->e, v_clear)) {
 801                         e_a_other[0] = l_a->prev->e;
 802                         e_a_other[1] = l_a->next->e;
 803                 }
 804                 else {
 805                         e_a_other[1] = l_a->prev->e;
 806                         e_a_other[0] = l_a->next->e;
 807                 }
 808
 809                 r_e_clear_other[0] = BM_elem_index_get(e_a_other[0]);
 810                 r_e_clear_other[1] = -1;
 811
 812 #ifdef USE_CUSTOMDATA
 813                 /* before killing, do customdata */
 814                 if (customdata_flag & CD_DO_VERT) {
 815                         BM_data_interp_from_verts(bm, v_other, v_clear, v_other, customdata_fac);
 816                 }
 817                 if (customdata_flag & CD_DO_EDGE) {
 818                         BM_data_interp_from_edges(bm, e_a_other[1], e_a_other[0], e_a_other[1], customdata_fac);
 819                 }
 820                 if (customdata_flag & CD_DO_LOOP) {
 821                         bm_edge_collapse_loop_customdata(bm, e_clear->l, v_clear, v_other, customdata_fac);
 822                 }
 823 #endif
 824
 825                 BM_edge_kill(bm, e_clear);
 826
 827                 v_other->head.hflag |= v_clear->head.hflag;
 828                 BM_vert_splice(bm, v_clear, v_other);
 829
 830                 e_a_other[1]->head.hflag |= e_a_other[0]->head.hflag;
 831                 BM_edge_splice(bm, e_a_other[0], e_a_other[1]);
 832
 833                 // BM_mesh_validate(bm);
 834
 835                 return true;
 836         }
 837         else {
 838                 return false;
 839         }
 840 }
 841
 842
 843 /* collapse e the edge, removing e->v2 */
 844 static void bm_decim_edge_collapse(BMesh *bm, BMEdge *e,
 845                                    Quadric *vquadrics, float *vweights,
 846                                    Heap *eheap, HeapNode **eheap_table,
 847                                    const CD_UseFlag customdata_flag)
 848 {
 849         int e_clear_other[2];
 850         BMVert *v_other = e->v1;
 851         int v_clear_index = BM_elem_index_get(e->v2);  /* the vert is removed so only store the index */
 852         float optimize_co[3];
 853         float customdata_fac;
 854
 855 #ifdef USE_VERT_NORMAL_INTERP
 856         float v_clear_no[3];
 857         copy_v3_v3(v_clear_no, e->v2->no);
 858 #endif
 859
 860         /* disallow collapsing which results in degenerate cases */
 861         if (UNLIKELY(bm_edge_collapse_is_degenerate_topology(e))) {
 862                 bm_decim_invalid_edge_cost_single(e, eheap, eheap_table);  /* add back with a high cost */
 863                 return;
 864         }
 865
 866         bm_decim_calc_target_co(e, optimize_co, vquadrics);
 867
 868         /* check if this would result in an overlapping face */
 869         if (UNLIKELY(bm_edge_collapse_is_degenerate_flip(e, optimize_co))) {
 870                 bm_decim_invalid_edge_cost_single(e, eheap, eheap_table);  /* add back with a high cost */
 871                 return;
 872         }
 873
 874         /* use for customdata merging */
 875         if (LIKELY(compare_v3v3(e->v1->co, e->v2->co, FLT_EPSILON) == false)) {
 876                 customdata_fac = line_point_factor_v3(optimize_co, e->v1->co, e->v2->co);
 877 #if 0
 878                 /* simple test for stupid collapse */
 879                 if (customdata_fac < 0.0 - FLT_EPSILON || customdata_fac > 1.0f + FLT_EPSILON) {
 880                         return;
 881                 }
 882 #endif
 883         }
 884         else {
 885                 /* avoid divide by zero */
 886                 customdata_fac = 0.5f;
 887         }
 888
 889         if (bm_edge_collapse(bm, e, e->v2, e_clear_other, customdata_flag, customdata_fac)) {
 890                 /* update collapse info */
 891                 int i;
 892
 893                 if (vweights) {
 894                         vweights[BM_elem_index_get(v_other)] += vweights[v_clear_index];
 895                 }
 896
 897                 e = NULL;  /* paranoid safety check */
 898
 899                 copy_v3_v3(v_other->co, optimize_co);
 900
 901                 /* remove eheap */
 902                 for (i = 0; i < 2; i++) {
 903                         /* highly unlikely 'eheap_table[ke_other[i]]' would be NULL, but do for sanity sake */
 904                         if ((e_clear_other[i] != -1) && (eheap_table[e_clear_other[i]] != NULL)) {
 905                                 BLI_heap_remove(eheap, eheap_table[e_clear_other[i]]);
 906                                 eheap_table[e_clear_other[i]] = NULL;
 907                         }
 908                 }
 909
 910                 /* update vertex quadric, add kept vertex from killed vertex */
 911                 BLI_quadric_add_qu_qu(&vquadrics[BM_elem_index_get(v_other)], &vquadrics[v_clear_index]);
 912
 913                 /* update connected normals */
 914
 915                 /* in fact face normals are not used for progressive updates, no need to update them */
 916                 // BM_vert_normal_update_all(v);
 917 #ifdef USE_VERT_NORMAL_INTERP
 918                 interp_v3_v3v3(v_other->no, v_other->no, v_clear_no, customdata_fac);
 919                 normalize_v3(v_other->no);
 920 #else
 921                 BM_vert_normal_update(v_other);
 922 #endif
 923
 924
 925                 /* update error costs and the eheap */
 926                 if (LIKELY(v_other->e)) {
 927                         BMEdge *e_iter;
 928                         BMEdge *e_first;
 929                         e_iter = e_first = v_other->e;
 930                         do {
 931                                 BLI_assert(BM_edge_find_double(e_iter) == NULL);
 932                                 bm_decim_build_edge_cost_single(e_iter, vquadrics, vweights, eheap, eheap_table);
 933                         } while ((e_iter = bmesh_disk_edge_next(e_iter, v_other)) != e_first);
 934                 }
 935
 936                 /* this block used to be disabled,
 937                  * but enable now since surrounding faces may have been
 938                  * set to COST_INVALID because of a face overlap that no longer occurs */
 939 #if 1
 940                 /* optional, update edges around the vertex face fan */
 941                 {
 942                         BMIter liter;
 943                         BMLoop *l;
 944                         BM_ITER_ELEM (l, &liter, v_other, BM_LOOPS_OF_VERT) {
 945                                 if (l->f->len == 3) {
 946                                         BMEdge *e_outer;
 947                                         if (BM_vert_in_edge(l->prev->e, l->v))
 948                                                 e_outer = l->next->e;
 949                                         else
 950                                                 e_outer = l->prev->e;
 951
 952                                         BLI_assert(BM_vert_in_edge(e_outer, l->v) == false);
 953
 954                                         bm_decim_build_edge_cost_single(e_outer, vquadrics, vweights, eheap, eheap_table);
 955                                 }
 956                         }
 957                 }
 958                 /* end optional update */
 959 #endif
 960         }
 961         else {
 962                 /* add back with a high cost */
 963                 bm_decim_invalid_edge_cost_single(e, eheap, eheap_table);
 964         }
 965 }
 966
 967
 968 /* Main Decimate Function
 969  * ********************** */
 970
 971 /**
 972  * \brief BM_mesh_decimate
 973  * \param bm The mesh
 974  * \param factor face count multiplier [0 - 1]
 975  * \param vweights Optional array of vertex  aligned weights [0 - 1],
 976  *        a vertex group is the usual source for this.
 977  */
 978 void BM_mesh_decimate_collapse(BMesh *bm, const float factor, float *vweights, const bool do_triangulate)
 979 {
 980         Heap *eheap;             /* edge heap */
 981         HeapNode **eheap_table;  /* edge index aligned table pointing to the eheap */
 982         Quadric *vquadrics;      /* vert index aligned quadrics */
 983         int tot_edge_orig;
 984         int face_tot_target;
 985         bool use_triangulate;
 986
 987         CD_UseFlag customdata_flag = 0;
 988
 989 #ifdef USE_TRIANGULATE
 990         /* temp convert quads to triangles */
 991         use_triangulate = bm_decim_triangulate_begin(bm);
 992 #endif
 993
 994
 995         /* alloc vars */
 996         vquadrics = MEM_callocN(sizeof(Quadric) * bm->totvert, __func__);
 997         /* since some edges may be degenerate, we might be over allocing a little here */
 998         eheap = BLI_heap_new_ex(bm->totedge);
 999         eheap_table = MEM_mallocN(sizeof(HeapNode *) * bm->totedge, __func__);
1000         tot_edge_orig = bm->totedge;
1001
1002
1003         /* build initial edge collapse cost data */
1004         bm_decim_build_quadrics(bm, vquadrics);
1005
1006         bm_decim_build_edge_cost(bm, vquadrics, vweights, eheap, eheap_table);
1007
1008         face_tot_target = bm->totface * factor;
1009         bm->elem_index_dirty |= BM_FACE | BM_EDGE | BM_VERT;
1010
1011
1012 #ifdef USE_CUSTOMDATA
1013         /* initialize customdata flag, we only need math for loops */
1014         if (CustomData_has_interp(&bm->vdata))  customdata_flag |= CD_DO_VERT;
1015         if (CustomData_has_interp(&bm->edata))  customdata_flag |= CD_DO_EDGE;
1016         if (CustomData_has_math(&bm->ldata))    customdata_flag |= CD_DO_LOOP;
1017 #endif
1018
1019         /* iterative edge collapse and maintain the eheap */
1020         while ((bm->totface > face_tot_target) &&
1021                (BLI_heap_is_empty(eheap) == false) &&
1022                (BLI_heap_node_value(BLI_heap_top(eheap)) != COST_INVALID))
1023         {
1024                 // const float value = BLI_heap_node_value(BLI_heap_top(eheap));
1025                 BMEdge *e = BLI_heap_popmin(eheap);
1026                 BLI_assert(BM_elem_index_get(e) < tot_edge_orig);  /* handy to detect corruptions elsewhere */
1027
1028                 // printf("COST %.10f\n", value);
1029
1030                 /* under normal conditions wont be accessed again,
1031                  * but NULL just incase so we don't use freed node */
1032                 eheap_table[BM_elem_index_get(e)] = NULL;
1033
1034                 bm_decim_edge_collapse(bm, e, vquadrics, vweights, eheap, eheap_table, customdata_flag);
1035         }
1036
1037
1038 #ifdef USE_TRIANGULATE
1039         if (do_triangulate == false) {
1040                 /* its possible we only had triangles, skip this step in that case */
1041                 if (LIKELY(use_triangulate)) {
1042                         /* temp convert quads to triangles */
1043                         bm_decim_triangulate_end(bm);
1044                 }
1045         }
1046 #endif
1047
1048         /* free vars */
1049         MEM_freeN(vquadrics);
1050         MEM_freeN(eheap_table);
1051         BLI_heap_free(eheap, NULL);
1052
1053         /* testing only */
1054         // BM_mesh_validate(bm);
1055
1056         (void)tot_edge_orig;  /* quiet release build warning */
1057 }