Fix delaunay triangulation, bad indices for output faces.

If there were merged vertices, sometimes the output faces
had wrong vertex indices. Added a test for this, and fixed.

source/blender/blenlib/intern/delaunay_2d.c

@@ -2965,6 +2965,13 @@ static CDT_result *cdt_get_output(CDT_state *cdt,
   SymEdge *se, *se_start;
   CDTEdge *e;
   CDTFace *f;
+#ifdef DEBUG_CDT
+  int dbg_level = 0;
+
+  if (dbg_level > 0) {
+    fprintf(stderr, "\nCDT_GET_OUTPUT\n\n");
+  }
+#endif
 
   prepare_cdt_for_output(cdt, output_type);
 
@@ -2973,6 +2980,12 @@ static CDT_result *cdt_get_output(CDT_state *cdt,
     return result;
   }
 
+#ifdef DEBUG_CDT
+  if (dbg_level > 1) {
+    dump_cdt(cdt, "cdt to output");
+  }
+#endif
+
   /* All verts without a merge_to_index will be output.
    * vert_to_output_map[i] will hold the output vertex index
    * corresponding to the vert in position i in cdt->vert_array.
@@ -3112,7 +3125,7 @@ static CDT_result *cdt_get_output(CDT_state *cdt,
         result->faces_start_table[i] = j;
         se = se_start = f->symedge;
         do {
-          result->faces[j++] = se->vert->index;
+          result->faces[j++] = vert_to_output_map[se->vert->index];
           se = se->next;
         } while (se != se_start);
         result->faces_len_table[i] = j - result->faces_start_table[i];
@@ -3216,7 +3229,7 @@ CDT_result *BLI_delaunay_2d_cdt_calc(const CDT_input *input, const CDT_output_ty
     ne = input->edges_len;
     nf = input->faces_len;
 #ifdef DEBUG_CDT
-    if (dbg_level > 4) {
+    if (dbg_level > 0) {
       fprintf(stderr, "input modified for near ends; now ne=%d\n", ne);
     }
 #endif

tests/gtests/blenlib/BLI_delaunay_2d_test.cc

@@ -954,6 +954,70 @@ TEST(delaunay, TwoSquaresOverlap)
   BLI_delaunay_2d_cdt_free(out);
 }
 
+TEST(delaunay, TwoFaceEdgeOverlap)
+{
+  CDT_input in;
+  CDT_result *out;
+  int i, v_out[6], v_int;
+  int e01, e1i, ei2, e20, e24, e4i, ei0;
+  int f02i, f24i, f10i;
+  const char *spec = R"(6 0 2
+  5.657 0.0
+  -1.414 -5.831
+  0.0 0.0
+  5.657 0.0
+  -2.121 -2.915
+  0.0 0.0
+  2 1 0
+  5 4 3
+  )";
+
+  fill_input_from_string(&in, spec);
+  out = BLI_delaunay_2d_cdt_calc(&in, CDT_CONSTRAINTS);
+  EXPECT_EQ(out->verts_len, 5);
+  EXPECT_EQ(out->edges_len, 7);
+  EXPECT_EQ(out->faces_len, 3);
+  if (out->verts_len == 5 && out->edges_len == 7 && out->faces_len == 3) {
+    v_int = 4;
+    for (i = 0; i < 6; i++) {
+      v_out[i] = get_output_vert_index(out, i);
+      EXPECT_NE(v_out[i], -1);
+      EXPECT_NE(v_out[i], v_int);
+    }
+    EXPECT_EQ(v_out[0], v_out[3]);
+    EXPECT_EQ(v_out[2], v_out[5]);
+    e01 = get_edge(out, v_out[0], v_out[1]);
+    EXPECT_TRUE(out_edge_has_input_id(out, e01, 1));
+    e1i = get_edge(out, v_out[1], v_int);
+    EXPECT_TRUE(out_edge_has_input_id(out, e1i, 0));
+    ei2 = get_edge(out, v_int, v_out[2]);
+    EXPECT_TRUE(out_edge_has_input_id(out, ei2, 0));
+    e20 = get_edge(out, v_out[2], v_out[0]);
+    EXPECT_TRUE(out_edge_has_input_id(out, e20, 2));
+    EXPECT_TRUE(out_edge_has_input_id(out, e20, 5));
+    e24 = get_edge(out, v_out[2], v_out[4]);
+    EXPECT_TRUE(out_edge_has_input_id(out, e24, 3));
+    e4i = get_edge(out, v_out[4], v_int);
+    EXPECT_TRUE(out_edge_has_input_id(out, e4i, 4));
+    ei0 = get_edge(out, v_int, v_out[0]);
+    EXPECT_TRUE(out_edge_has_input_id(out, ei0, 4));
+    f02i = get_face_tri(out, v_out[0], v_out[2], v_int);
+    EXPECT_NE(f02i, -1);
+    EXPECT_TRUE(out_face_has_input_id(out, f02i, 0));
+    EXPECT_TRUE(out_face_has_input_id(out, f02i, 1));
+    f24i = get_face_tri(out, v_out[2], v_out[4], v_int);
+    EXPECT_NE(f24i, -1);
+    EXPECT_TRUE(out_face_has_input_id(out, f24i, 1));
+    EXPECT_FALSE(out_face_has_input_id(out, f24i, 0));
+    f10i = get_face_tri(out, v_out[1], v_out[0], v_int);
+    EXPECT_NE(f10i, -1);
+    EXPECT_TRUE(out_face_has_input_id(out, f10i, 0));
+    EXPECT_FALSE(out_face_has_input_id(out, f10i, 1));
+  }
+  free_spec_arrays(&in);
+  BLI_delaunay_2d_cdt_free(out);
+}
+
 TEST(delaunay, TriInTri)
 {
   CDT_input in;