bmesh: improve limited dissolve result

iteratively dissolve the best edge/vert, updating the heap as the dissolve runs.

1 files changed, 202 insertions, 146 deletions

diff --git a/source/blender/bmesh/tools/bmesh_decimate_dissolve.c b/source/blender/bmesh/tools/bmesh_decimate_dissolve.c
index 9f97d8f4287..310357453a3 100644
--- a/source/blender/bmesh/tools/bmesh_decimate_dissolve.c
+++ b/source/blender/bmesh/tools/bmesh_decimate_dissolve.c
@@ -30,18 +30,22 @@
 #include "MEM_guardedalloc.h"
 
 #include "BLI_math.h"
+#include "BLI_heap.h"
 
 #include "bmesh.h"
 #include "bmesh_decimate.h"  /* own include */
 
-#define UNIT_TO_ANGLE DEG2RADF(90.0f)
-#define ANGLE_TO_UNIT (1.0f / UNIT_TO_ANGLE)
+#define COST_INVALID FLT_MAX
+
 
 /* multiply vertex edge angle by face angle
  * this means we are not left with sharp corners between _almost_ planer faces
  * convert angles [0-PI/2] -> [0-1], multiply together, then convert back to radians. */
 static float bm_vert_edge_face_angle(BMVert *v)
 {
+#define UNIT_TO_ANGLE DEG2RADF(90.0f)
+#define ANGLE_TO_UNIT (1.0f / UNIT_TO_ANGLE)
+
 	const float angle = BM_vert_calc_edge_angle(v);
 	/* note: could be either edge, it doesn't matter */
 	if (v->e && BM_edge_is_manifold(v->e)) {
@@ -50,167 +54,184 @@ static float bm_vert_edge_face_angle(BMVert *v)
 	else {
 		return angle;
 	}
-}
 
 #undef UNIT_TO_ANGLE
 #undef ANGLE_TO_UNIT
+}
 
-typedef struct DissolveElemWeight {
-	BMHeader *ele;
-	float weight;
-} DissolveElemWeight;
-
-static int dissolve_elem_cmp(const void *a1, const void *a2)
+static float bm_edge_calc_dissolve_error(const BMEdge *e, const BMO_Delimit delimit)
 {
-	const struct DissolveElemWeight *d1 = a1, *d2 = a2;
+	const bool is_contig = BM_edge_is_contiguous(e);
+	float angle;
+
+	if (!BM_edge_is_manifold(e)) {
+		goto fail;
+	}
+
+	if ((delimit & BMO_DELIM_SEAM) &&
+	    (BM_elem_flag_test(e, BM_ELEM_SEAM)))
+	{
+		goto fail;
+	}
+
+	if ((delimit & BMO_DELIM_MATERIAL) &&
+	    (e->l->f->mat_nr != e->l->radial_next->f->mat_nr))
+	{
+		goto fail;
+	}
+
+	if ((delimit & BMO_DELIM_NORMAL) &&
+	    (is_contig == false))
+	{
+		goto fail;
+	}
+
+	angle = BM_edge_calc_face_angle(e);
+	if (is_contig == false) {
+		angle = (float)M_PI - angle;
+	}
 
-	if      (d1->weight > d2->weight) return  1;
-	else if (d1->weight < d2->weight) return -1;
-	return 0;
+	return angle;
+
+fail:
+	return COST_INVALID;
 }
 
+
 void BM_mesh_decimate_dissolve_ex(BMesh *bm, const float angle_limit, const bool do_dissolve_boundaries,
                                   const BMO_Delimit delimit,
                                   BMVert **vinput_arr, const int vinput_len,
                                   BMEdge **einput_arr, const int einput_len,
                                   const short oflag_out)
 {
-	const float angle_max = (float)M_PI / 2.0f;
-	DissolveElemWeight *weight_elems = MEM_mallocN(max_ii(einput_len, vinput_len) *
-	                                               sizeof(DissolveElemWeight), __func__);
-	int i, tot_found;
-	BMIter iter;
-	BMEdge *e_iter;
-	BMEdge **earray;
+	const int eheap_table_len = do_dissolve_boundaries ? einput_len : max_ii(einput_len, vinput_len);
+	void *_heap_table = MEM_mallocN(sizeof(HeapNode *) * eheap_table_len, __func__);
 
-	int *vert_reverse_lookup;
+	int i;
 
 	/* --- first edges --- */
-
-	/* wire -> tag */
-	BM_ITER_MESH (e_iter, &iter, bm, BM_EDGES_OF_MESH) {
-		BM_elem_flag_set(e_iter, BM_ELEM_TAG, BM_edge_is_wire(e_iter));
-	}
-
-	/* go through and split edge */
-	for (i = 0, tot_found = 0; i < einput_len; i++) {
-		BMEdge *e = einput_arr[i];
-		const bool is_contig = BM_edge_is_contiguous(e);
-		float angle;
-
-		angle = BM_edge_calc_face_angle(e);
-		if (is_contig == false) {
-			angle = (float)M_PI - angle;
+	if (1) {
+		BMEdge **earray;
+		Heap *eheap;
+		HeapNode **eheap_table = _heap_table;
+		HeapNode *enode_top;
+		int *vert_reverse_lookup;
+		BMIter iter;
+		BMEdge *e_iter;
+
+		/* --- setup heap --- */
+		eheap = BLI_heap_new_ex(einput_len);
+		eheap_table = _heap_table;
+
+		/* wire -> tag */
+		BM_ITER_MESH (e_iter, &iter, bm, BM_EDGES_OF_MESH) {
+			BM_elem_flag_set(e_iter, BM_ELEM_TAG, BM_edge_is_wire(e_iter));
+			BM_elem_index_set(e_iter, -1);  /* set dirty */
 		}
-
-		if (angle < angle_limit) {
-			tot_found++;
+		bm->elem_index_dirty |= BM_EDGE;
+
+		/* build heap */
+		for (i = 0; i < einput_len; i++) {
+			BMEdge *e = einput_arr[i];
+			const float cost = bm_edge_calc_dissolve_error(e, delimit);
+			eheap_table[i] = BLI_heap_insert(eheap, cost, e);
+			BM_elem_index_set(e, i);  /* set dirty */
 		}
-		weight_elems[i].ele = (BMHeader *)e;
-		weight_elems[i].weight = angle;
-	}
-
-	if (tot_found != 0) {
-		qsort(weight_elems, einput_len, sizeof(DissolveElemWeight), dissolve_elem_cmp);
 
-		for (i = 0; i < tot_found; i++) {
-			BMEdge *e = (BMEdge *)weight_elems[i].ele;
-			const bool is_contig = BM_edge_is_contiguous(e);
-			float angle;
+		while ((BLI_heap_is_empty(eheap) == false) &&
+		       (BLI_heap_node_value((enode_top = BLI_heap_top(eheap))) < angle_limit))
+		{
+			BMFace *f_new = NULL;
+			BMEdge *e;
 
-			/* may have become non-manifold */
-			if (!BM_edge_is_manifold(e)) {
-				continue;
-			}
+			e = BLI_heap_node_ptr(enode_top);
+			i = BM_elem_index_get(e);
 
-			if ((delimit & BMO_DELIM_SEAM) &&
-			    (BM_elem_flag_test(e, BM_ELEM_SEAM)))
-			{
-				continue;
-			}
-
-			if ((delimit & BMO_DELIM_MATERIAL) &&
-			    (e->l->f->mat_nr != e->l->radial_next->f->mat_nr))
-			{
-				continue;
-			}
-
-			if ((delimit & BMO_DELIM_NORMAL) &&
-			    (is_contig == false))
-			{
-				continue;
-			}
+			if (BM_edge_is_manifold(e)) {
+				f_new = BM_faces_join_pair(bm, e->l->f,
+				                           e->l->radial_next->f, e,
+				                           false); /* join faces */
 
-			/* check twice because cumulative effect could dissolve over angle limit */
-			angle = BM_edge_calc_face_angle(e);
-			if (is_contig == false) {
-				angle = (float)M_PI - angle;
-			}
+				if (f_new) {
+					BMLoop *l_first, *l_iter;
 
-			if (angle < angle_limit) {
-				BMFace *f_new = BM_faces_join_pair(bm, e->l->f,
-				                                   e->l->radial_next->f,
-				                                   e,
-				                                   false); /* join faces */
+					BLI_heap_remove(eheap, enode_top);
+					eheap_table[i] = NULL;
 
-				/* there may be some errors, we don't mind, just move on */
-				if (f_new) {
+					/* update normal */
 					BM_face_normal_update(f_new);
 					if (oflag_out) {
 						BMO_elem_flag_enable(bm, f_new, oflag_out);
 					}
+
+					/* re-calculate costs */
+					l_iter = l_first = BM_FACE_FIRST_LOOP(f_new);
+					do {
+						const int j = BM_elem_index_get(l_iter->e);
+						if (j != -1 && eheap_table[j]) {
+							const float cost = bm_edge_calc_dissolve_error(l_iter->e, delimit);
+							BLI_heap_remove(eheap, eheap_table[j]);
+							eheap_table[j] = BLI_heap_insert(eheap, cost, l_iter->e);
+						}
+					} while ((l_iter = l_iter->next) != l_first);
 				}
 				else {
 					BMO_error_clear(bm);
 				}
 			}
+
+			if (UNLIKELY(f_new == NULL)) {
+				BLI_heap_remove(eheap, enode_top);
+				eheap_table[i] = BLI_heap_insert(eheap, COST_INVALID, e);
+			}
 		}
-	}
 
-	/* prepare for cleanup */
-	BM_mesh_elem_index_ensure(bm, BM_VERT);
-	vert_reverse_lookup = MEM_mallocN(sizeof(int) * bm->totvert, __func__);
-	fill_vn_i(vert_reverse_lookup, bm->totvert, -1);
-	for (i = 0, tot_found = 0; i < vinput_len; i++) {
-		BMVert *v = vinput_arr[i];
-		vert_reverse_lookup[BM_elem_index_get(v)] = i;
-	}
+		/* prepare for cleanup */
+		BM_mesh_elem_index_ensure(bm, BM_VERT);
+		vert_reverse_lookup = MEM_mallocN(sizeof(int) * bm->totvert, __func__);
+		fill_vn_i(vert_reverse_lookup, bm->totvert, -1);
+		for (i = 0; i < vinput_len; i++) {
+			BMVert *v = vinput_arr[i];
+			vert_reverse_lookup[BM_elem_index_get(v)] = i;
+		}
 
-	/* --- cleanup --- */
-	earray = MEM_mallocN(sizeof(BMEdge *) * bm->totedge, __func__);
-	BM_ITER_MESH_INDEX (e_iter, &iter, bm, BM_EDGES_OF_MESH, i) {
-		earray[i] = e_iter;
-	}
-	/* remove all edges/verts left behind from dissolving, NULL'ing the vertex array so we dont re-use */
-	for (i = bm->totedge - 1; i != -1; i--) {
-		e_iter = earray[i];
-
-		if (BM_edge_is_wire(e_iter) && (BM_elem_flag_test(e_iter, BM_ELEM_TAG) == false)) {
-			/* edge has become wire */
-			int vidx_reverse;
-			BMVert *v1 = e_iter->v1;
-			BMVert *v2 = e_iter->v2;
-			BM_edge_kill(bm, e_iter);
-			if (v1->e == NULL) {
-				vidx_reverse = vert_reverse_lookup[BM_elem_index_get(v1)];
-				if (vidx_reverse != -1) vinput_arr[vidx_reverse] = NULL;
-				BM_vert_kill(bm, v1);
-			}
-			if (v2->e == NULL) {
-				vidx_reverse = vert_reverse_lookup[BM_elem_index_get(v2)];
-				if (vidx_reverse != -1) vinput_arr[vidx_reverse] = NULL;
-				BM_vert_kill(bm, v2);
+		/* --- cleanup --- */
+		earray = MEM_mallocN(sizeof(BMEdge *) * bm->totedge, __func__);
+		BM_ITER_MESH_INDEX (e_iter, &iter, bm, BM_EDGES_OF_MESH, i) {
+			earray[i] = e_iter;
+		}
+		/* remove all edges/verts left behind from dissolving, NULL'ing the vertex array so we dont re-use */
+		for (i = bm->totedge - 1; i != -1; i--) {
+			e_iter = earray[i];
+
+			if (BM_edge_is_wire(e_iter) && (BM_elem_flag_test(e_iter, BM_ELEM_TAG) == false)) {
+				/* edge has become wire */
+				int vidx_reverse;
+				BMVert *v1 = e_iter->v1;
+				BMVert *v2 = e_iter->v2;
+				BM_edge_kill(bm, e_iter);
+				if (v1->e == NULL) {
+					vidx_reverse = vert_reverse_lookup[BM_elem_index_get(v1)];
+					if (vidx_reverse != -1) vinput_arr[vidx_reverse] = NULL;
+					BM_vert_kill(bm, v1);
+				}
+				if (v2->e == NULL) {
+					vidx_reverse = vert_reverse_lookup[BM_elem_index_get(v2)];
+					if (vidx_reverse != -1) vinput_arr[vidx_reverse] = NULL;
+					BM_vert_kill(bm, v2);
+				}
 			}
 		}
-	}
-	MEM_freeN(vert_reverse_lookup);
+		MEM_freeN(vert_reverse_lookup);
+		MEM_freeN(earray);
 
-	MEM_freeN(earray);
+		BLI_heap_free(eheap, NULL);
+	}
 
 
 	/* --- second verts --- */
 	if (do_dissolve_boundaries) {
-		/* simple version of the branch below, sincve we will dissolve _all_ verts that use 2 edges */
+		/* simple version of the branch below, since we will dissolve _all_ verts that use 2 edges */
 		for (i = 0; i < vinput_len; i++) {
 			BMVert *v = vinput_arr[i];
 			if (LIKELY(v != NULL) &&
@@ -221,43 +242,78 @@ void BM_mesh_decimate_dissolve_ex(BMesh *bm, const float angle_limit, const bool
 		}
 	}
 	else {
-		for (i = 0, tot_found = 0; i < vinput_len; i++) {
-			BMVert *v = vinput_arr[i];
-			const float angle = v ? bm_vert_edge_face_angle(v) : angle_limit;
+		Heap *vheap;
+		HeapNode **vheap_table = _heap_table;
+		HeapNode *vnode_top;
 
-			if (angle < angle_limit) {
-				weight_elems[i].ele = (BMHeader *)v;
-				weight_elems[i].weight = angle;
-				tot_found++;
-			}
-			else {
-				weight_elems[i].ele = NULL;
-				weight_elems[i].weight = angle_max;
+		BMVert *v_iter;
+		BMIter iter;
+
+		BM_ITER_MESH (v_iter, &iter, bm, BM_VERTS_OF_MESH) {
+			BM_elem_index_set(v_iter, -1);  /* set dirty */
+		}
+		bm->elem_index_dirty |= BM_VERT;
+
+		vheap = BLI_heap_new_ex(vinput_len);
+
+		for (i = 0; i < vinput_len; i++) {
+			BMVert *v = vinput_arr[i];
+			if (LIKELY(v != NULL)) {
+				const float cost = bm_vert_edge_face_angle(v);
+				vheap_table[i] = BLI_heap_insert(vheap, cost, v);
+				BM_elem_index_set(v, i);  /* set dirty */
 			}
 		}
 
-		if (tot_found != 0) {
-			qsort(weight_elems, vinput_len, sizeof(DissolveElemWeight), dissolve_elem_cmp);
-
-			for (i = 0; i < tot_found; i++) {
-				BMVert *v = (BMVert *)weight_elems[i].ele;
-				if (LIKELY(v != NULL) &&
-				    /* topology changes may cause this to be un-collapsable */
-				    (BM_vert_edge_count(v) == 2) &&
-				    /* check twice because cumulative effect could dissolve over angle limit */
-				    bm_vert_edge_face_angle(v) < angle_limit)
-				{
-					BMEdge *e_new = BM_vert_collapse_edge(bm, v->e, v, true); /* join edges */
-
-					if (e_new && e_new->l) {
-						BM_edge_normals_update(e_new);
+		while ((BLI_heap_is_empty(vheap) == false) &&
+		       (BLI_heap_node_value((vnode_top = BLI_heap_top(vheap))) < angle_limit))
+		{
+			BMEdge *e_new = NULL;
+			BMVert *v;
+
+			v = BLI_heap_node_ptr(vnode_top);
+			i = BM_elem_index_get(v);
+
+			if (BM_vert_edge_count(v) == 2) {
+				e_new = BM_vert_collapse_edge(bm, v->e, v, true); /* join edges */
+
+				if (e_new) {
+
+					BLI_heap_remove(vheap, vnode_top);
+					vheap_table[i] = NULL;
+
+					/* update normal */
+					if (e_new->l) {
+						BMLoop *l_first, *l_iter;
+						l_iter = l_first = e_new->l;
+						do {
+							BM_face_normal_update(l_iter->f);
+						} while ((l_iter = l_iter->radial_next) != l_first);
+
+					}
+
+					/* re-calculate costs */
+					BM_ITER_ELEM(v_iter, &iter, e_new, BM_VERTS_OF_EDGE) {
+						const int j = BM_elem_index_get(v_iter);
+						if (j != -1 && vheap_table[j]) {
+							const float cost = bm_vert_edge_face_angle(v_iter);
+							BLI_heap_remove(vheap, vheap_table[j]);
+							vheap_table[j] = BLI_heap_insert(vheap, cost, v_iter);
+						}
 					}
 				}
 			}
+
+			if (UNLIKELY(e_new == NULL)) {
+				BLI_heap_remove(vheap, vnode_top);
+				vheap_table[i] = BLI_heap_insert(vheap, COST_INVALID, v);
+			}
 		}
+
+		BLI_heap_free(vheap, NULL);
 	}
 
-	MEM_freeN(weight_elems);
+	MEM_freeN(_heap_table);
 }
 
 void BM_mesh_decimate_dissolve(BMesh *bm, const float angle_limit, const bool do_dissolve_boundaries,