forked from bartvdbraak/blender
Fix crash in delaunay triangulation due to epsilon issues.
This commit is contained in:
Howard Trickey
parent
51d8d790d7
commit
a0892bb690
@ -53,6 +53,7 @@ typedef struct CDTVert {
|
||||
SymEdge *symedge; /* Some edge attached to it. */
|
||||
LinkNode *input_ids; /* List of corresponding vertex input ids. */
|
||||
int index; /* Index into array that cdt keeps. */
|
||||
int visit_index; /* Which visit epoch has this been seen. */
|
||||
} CDTVert;
|
||||
|
||||
typedef struct CDTEdge {
|
||||
@ -110,7 +111,9 @@ static void dump_v(const CDTVert *v, const char *lab);
|
||||
static void dump_se_cycle(const SymEdge *se, const char *lab, const int limit);
|
||||
static void dump_id_list(const LinkNode *id_list, const char *lab);
|
||||
static void dump_cdt(const CDT_state *cdt, const char *lab);
|
||||
static void dump_cdt_vert_neighborhood(CDT_state *cdt, int v, int maxdist, const char *lab);
|
||||
static void cdt_draw(CDT_state *cdt, const char *lab);
|
||||
static void write_cdt_input_to_file(const CDT_input *inp);
|
||||
static void validate_face_centroid(SymEdge *se);
|
||||
static void validate_cdt(CDT_state *cdt, bool check_all_tris);
|
||||
#endif
|
||||
@ -151,14 +154,11 @@ static int CCW_test(const double a[2], const double b[2], const double c[2], con
|
||||
return 0;
|
||||
}
|
||||
|
||||
/** return true if a -- b -- c are in that order, assuming they are on a straight line. */
|
||||
static bool in_line(const double a[2], const double b[2], const double c[2])
|
||||
/** return true if a -- b -- c are in that order, assuming they are on a straight line according to
|
||||
* CCW_test. */
|
||||
static bool in_line(const double a[2], const double b[2], const double c[2], double eps)
|
||||
{
|
||||
double dir_ab[2], dir_ac[2];
|
||||
|
||||
sub_v2_v2v2_db(dir_ab, a, b);
|
||||
sub_v2_v2v2_db(dir_ac, a, c);
|
||||
return dot_v2v2_db(dir_ab, dir_ac) >= 0.0;
|
||||
return fabs(len_v2v2_db(a, c) - (len_v2v2_db(a, b) + len_v2v2_db(b, c))) <= eps;
|
||||
}
|
||||
|
||||
#ifndef NDEBUG
|
||||
@ -208,6 +208,7 @@ static CDTVert *add_cdtvert(CDT_state *cdt, double x, double y)
|
||||
}
|
||||
BLI_assert(cdt->vert_array_len < cdt->vert_array_len_alloc);
|
||||
v->index = cdt->vert_array_len;
|
||||
v->visit_index = 0;
|
||||
cdt->vert_array[cdt->vert_array_len++] = v;
|
||||
return v;
|
||||
}
|
||||
@ -365,8 +366,8 @@ static CDTEdge *add_diagonal(CDT_state *cdt, SymEdge *s1, SymEdge *s2)
|
||||
CDTFace *fold, *fnew;
|
||||
SymEdge *sdiag, *sdiagsym;
|
||||
SymEdge *s1prev, *s1prevsym, *s2prev, *s2prevsym, *se;
|
||||
BLI_assert(reachable(s1, s2, 20));
|
||||
BLI_assert(reachable(s2, s1, 20));
|
||||
BLI_assert(reachable(s1, s2, 2000));
|
||||
BLI_assert(reachable(s2, s1, 2000));
|
||||
fold = s1->face;
|
||||
fnew = add_cdtface(cdt);
|
||||
s1prev = prev(s1);
|
||||
@ -385,7 +386,7 @@ static CDTEdge *add_diagonal(CDT_state *cdt, SymEdge *s1, SymEdge *s2)
|
||||
s2->rot = sdiagsym;
|
||||
sdiagsym->rot = s2prevsym;
|
||||
#ifdef DEBUG_CDT
|
||||
BLI_assert(reachable(s2, sdiag, 20));
|
||||
BLI_assert(reachable(s2, sdiag, 2000));
|
||||
#endif
|
||||
for (se = s2; se != sdiag; se = se->next) {
|
||||
se->face = fnew;
|
||||
@ -1505,7 +1506,8 @@ static void add_edge_constraint(
|
||||
SymEdge *se;
|
||||
CDTEdge *edge;
|
||||
int ccw1, ccw2, isect;
|
||||
int i, search_count;
|
||||
int i, search_count, visit;
|
||||
double curco[2];
|
||||
double lambda;
|
||||
const double epsilon = cdt->epsilon;
|
||||
bool done, state_through_vert;
|
||||
@ -1553,7 +1555,7 @@ static void add_edge_constraint(
|
||||
* in = a
|
||||
* out = b
|
||||
*
|
||||
* crossings[0] will have in = NULL, and crossings[last] will have out = NULL
|
||||
* crossings[0] will have in = NULL, and crossings[last] will have out = NULL.
|
||||
*/
|
||||
if (r_edges) {
|
||||
*r_edges = NULL;
|
||||
@ -1569,6 +1571,9 @@ static void add_edge_constraint(
|
||||
dump_se(vse1, " se1");
|
||||
dump_se(vse2, " se2");
|
||||
}
|
||||
if (dbg_level > 2) {
|
||||
dump_cdt(cdt, "before insert_segment");
|
||||
}
|
||||
}
|
||||
#endif
|
||||
if (v1 == v2) {
|
||||
@ -1579,6 +1584,37 @@ static void add_edge_constraint(
|
||||
#endif
|
||||
return;
|
||||
}
|
||||
/* Quick check for special case where segment is one of edges from v1. */
|
||||
tstart = t = vse1;
|
||||
done = false;
|
||||
do {
|
||||
if (t->next->vert == v2) {
|
||||
cdata.in = NULL;
|
||||
cdata.out = t->next;
|
||||
cdata.lambda = 0.0;
|
||||
cdata.vert = t->vert;
|
||||
BLI_array_append(crossings, cdata);
|
||||
cdata.in = t->next->rot;
|
||||
cdata.out = NULL;
|
||||
cdata.lambda = 0.0;
|
||||
cdata.vert = v2;
|
||||
BLI_array_append(crossings, cdata);
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "special one segment case\n");
|
||||
dump_se(t, " ");
|
||||
}
|
||||
#endif
|
||||
done = true;
|
||||
break;
|
||||
}
|
||||
t = t->rot;
|
||||
} while (t != tstart);
|
||||
if (!done) {
|
||||
/* To prevent infinite loop in the face of epsilon tests that might lead us back to
|
||||
* an already-visited (vertex, face) pair, use visit indices.
|
||||
*/
|
||||
visit = ++cdt->visit_count;
|
||||
state_through_vert = true;
|
||||
done = false;
|
||||
t = vse1;
|
||||
@ -1588,16 +1624,23 @@ static void add_edge_constraint(
|
||||
* but not including the crossing of edge t, which will either be through a vert
|
||||
* (if state_through_vert is true) or through edge t not at either end.
|
||||
* In the latter case, t->face is the face that ray v1--v2 goes through after path-so-far.
|
||||
* Rather than continually looking for intersection of v1--v2, we keep track of
|
||||
* last vertex or intersection point in curco,because that may be slightly off the ray
|
||||
* v1--v2.
|
||||
*/
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(
|
||||
stderr, "top of insert_segment main loop, state_through_vert=%d\n", state_through_vert);
|
||||
fprintf(stderr,
|
||||
"top of insert_segment main loop, state_through_vert=%d\n",
|
||||
state_through_vert);
|
||||
dump_se_cycle(t, "current t ", 4);
|
||||
}
|
||||
#endif
|
||||
BLI_assert(t->vert->visit_index != visit || t->face->visit_index != visit);
|
||||
t->vert->visit_index = visit;
|
||||
t->face->visit_index = visit;
|
||||
if (state_through_vert) {
|
||||
/* Invariant: ray v1--v2 contains t->vert. */
|
||||
/* Invariant: ray vcur--v2 contains t->vert. */
|
||||
cdata.in = (BLI_array_len(crossings) == 0) ? NULL : t;
|
||||
cdata.out = NULL; /* To be filled in if this isn't final. */
|
||||
cdata.lambda = 0.0;
|
||||
@ -1612,7 +1655,7 @@ static void add_edge_constraint(
|
||||
done = true;
|
||||
}
|
||||
else {
|
||||
/* Do ccw scan of triangles around t->vert to find exit triangle for ray v1--v2. */
|
||||
/* Do ccw scan of triangles around t->vert to find exit triangle for ray vcur--v2. */
|
||||
tstart = t;
|
||||
tout = NULL;
|
||||
do {
|
||||
@ -1628,7 +1671,8 @@ static void add_edge_constraint(
|
||||
fprintf(stderr, "ccw1=%d, ccw2=%d\n", ccw1, ccw2);
|
||||
}
|
||||
#endif
|
||||
if (ccw1 == 0 && in_line(t->vert->co, va->co, v2->co)) {
|
||||
if (ccw1 == 0 && in_line(t->vert->co, va->co, v2->co, epsilon) &&
|
||||
va->visit_index != visit) {
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "ray goes through va\n");
|
||||
@ -1639,7 +1683,8 @@ static void add_edge_constraint(
|
||||
t = t->next;
|
||||
break;
|
||||
}
|
||||
else if (ccw2 == 0 && in_line(t->vert->co, vb->co, v2->co)) {
|
||||
else if (ccw2 == 0 && in_line(t->vert->co, vb->co, v2->co, epsilon) &&
|
||||
vb->visit_index != visit) {
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "ray goes through vb\n");
|
||||
@ -1650,7 +1695,9 @@ static void add_edge_constraint(
|
||||
tout = sym(t);
|
||||
break;
|
||||
}
|
||||
else if (ccw1 > 0 && ccw2 < 0) {
|
||||
else if (ccw1 > 0 && ccw2 < 0 &&
|
||||
(t->next->vert->visit_index != visit ||
|
||||
t->next->face->visit_index != visit)) {
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "segment intersects\n");
|
||||
@ -1668,26 +1715,84 @@ static void add_edge_constraint(
|
||||
}
|
||||
#endif
|
||||
} while (t != tstart);
|
||||
BLI_assert(tout != NULL); /* TODO: something sensible for "this can't happen" */
|
||||
if (tout == NULL) {
|
||||
/* With exact arithmetic this shouldn't happen, but maybe the epsilon tests made it so
|
||||
* that we want to go back to a previous vertex.
|
||||
* As desperation measure, pick unvisited vertex that is closest in line with
|
||||
* destination.
|
||||
*/
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "add_edge_constraint desperation search\n");
|
||||
}
|
||||
#endif
|
||||
SymEdge *bestt = NULL;
|
||||
double dot, bestdot = -2.0;
|
||||
double dir_tv_v2[2], dir_tvnext_v2[2];
|
||||
sub_v2_v2v2_db(dir_tv_v2, v2->co, t->vert->co);
|
||||
do {
|
||||
if (t->next->vert->visit_index != visit) {
|
||||
sub_v2_v2v2_db(dir_tvnext_v2, v2->co, t->next->vert->co);
|
||||
dot = dot_v2v2_db(dir_tv_v2, dir_tvnext_v2);
|
||||
if (dot > bestdot) {
|
||||
bestdot = dot;
|
||||
bestt = t->next;
|
||||
}
|
||||
}
|
||||
t = t->rot;
|
||||
} while (t != tstart);
|
||||
if (bestt == NULL) {
|
||||
/* No unvisited place to go! Give up on adding this edge constraint. */
|
||||
#ifdef DEBUG_CDT
|
||||
fprintf(stderr, "could not add edge constraint\n");
|
||||
#endif
|
||||
return;
|
||||
}
|
||||
else {
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "add_edge_constraint desperation search chose to go through\n");
|
||||
dump_v(bestt->vert, "desperation vert");
|
||||
}
|
||||
#endif
|
||||
tout = bestt;
|
||||
t = t->next;
|
||||
}
|
||||
}
|
||||
crossings[BLI_array_len(crossings) - 1].out = tout;
|
||||
}
|
||||
}
|
||||
else { /* State is "through edge", not "through vert" */
|
||||
/* Invariant: ray v1--v2 intersects segment t->edge, not at either end.
|
||||
* and t->face is the face we have just passed through. */
|
||||
* and t->face is the face we have just passed through.
|
||||
* Whatever we make t next should not have both vert and face visited.
|
||||
*/
|
||||
va = t->vert;
|
||||
vb = t->next->vert;
|
||||
/* Get curco; cdata should have data from last time through the loop still. */
|
||||
if (cdata.lambda == 0.0) {
|
||||
copy_v3_v3_db(curco, cdata.vert->co);
|
||||
}
|
||||
else {
|
||||
interp_v2_v2v2_db(curco, cdata.in->vert->co, cdata.in->next->vert->co, lambda);
|
||||
}
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "through edge case\n");
|
||||
fprintf(stderr, "through edge case, curco=(%f,%f)\n", F2(curco));
|
||||
dump_v(va, " va");
|
||||
dump_v(vb, " vb");
|
||||
}
|
||||
#endif
|
||||
isect = isect_seg_seg_v2_lambda_mu_db(va->co, vb->co, v1->co, v2->co, &lambda, NULL);
|
||||
/* TODO: something sensible for "this can't happen" */
|
||||
BLI_assert(isect == ISECT_LINE_LINE_CROSS);
|
||||
UNUSED_VARS_NDEBUG(isect);
|
||||
isect = isect_seg_seg_v2_lambda_mu_db(va->co, vb->co, curco, v2->co, &lambda, NULL);
|
||||
if (isect != ISECT_LINE_LINE_CROSS) {
|
||||
/* Shouldn't happen. Just pick something. */
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "add_edge_constraint no intersect found, using lambda = 0.5\n");
|
||||
}
|
||||
#endif
|
||||
lambda = 0.5;
|
||||
}
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "intersect point at %f along va--vb\n", lambda);
|
||||
@ -1706,41 +1811,41 @@ static void add_edge_constraint(
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
dump_se_cycle(tout, "next search tri", 4);
|
||||
dump_se(tout, "tout");
|
||||
}
|
||||
#endif
|
||||
/* 'tout' is 'symedge' from 'vb' to third vertex, 'vc'. */
|
||||
/* 'tout' is 'symedge' from 'va' to third vertex, 'vc'. */
|
||||
BLI_assert(tout->vert == va);
|
||||
vc = tout->next->vert;
|
||||
ccw1 = CCW_test(v1->co, v2->co, vc->co, epsilon);
|
||||
ccw1 = CCW_test(curco, v2->co, vc->co, epsilon);
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "now searching with third vertex ");
|
||||
dump_v(vc, "vc");
|
||||
fprintf(stderr, "ccw(v1, v2, vc) = %d\n", ccw1);
|
||||
fprintf(stderr, "ccw(vcur, v2, vc) = %d\n", ccw1);
|
||||
}
|
||||
#endif
|
||||
if (ccw1 == -1) {
|
||||
/* v1--v2 should intersect vb--vc. */
|
||||
if (ccw1 == -1 && (vc->visit_index != visit || tout->next->face->visit_index != visit)) {
|
||||
/* vcur--v2 should intersect vb--vc. */
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "v1--v2 intersects vb--vc\n");
|
||||
fprintf(stderr, "curco--v2 intersects vb--vc\n");
|
||||
}
|
||||
#endif
|
||||
t = tout->next;
|
||||
state_through_vert = false;
|
||||
}
|
||||
else if (ccw1 == 1) {
|
||||
/* v1--v2 should intersect va--vc. */
|
||||
else if (ccw1 == 1 && tout->face->visit_index != visit) {
|
||||
/* vcur--v2 should intersect va--vc. */
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "v1--v2 intersects va--vc\n");
|
||||
fprintf(stderr, "curco--v2 intersects va--vc\n");
|
||||
}
|
||||
#endif
|
||||
t = tout;
|
||||
state_through_vert = false;
|
||||
}
|
||||
else {
|
||||
/* ccw1 == 0. */
|
||||
else if (ccw1 == 0 && vc->visit_index != visit) {
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 1) {
|
||||
fprintf(stderr, "ccw==0 case, so going through or to vc\n");
|
||||
@ -1749,6 +1854,16 @@ static void add_edge_constraint(
|
||||
t = tout->next;
|
||||
state_through_vert = true;
|
||||
}
|
||||
else {
|
||||
#ifdef DEBUG_CDT
|
||||
fprintf(stderr, "add_edge_constraint desperation search 2\n");
|
||||
#endif
|
||||
if (tout->face->visit_index != visit) {
|
||||
/* Treat as if an intersection of va-vc. */
|
||||
t = tout;
|
||||
state_through_vert = false;
|
||||
}
|
||||
}
|
||||
}
|
||||
if (++search_count > 10000) {
|
||||
fprintf(stderr, "infinite loop? bailing out\n");
|
||||
@ -1756,6 +1871,7 @@ static void add_edge_constraint(
|
||||
break;
|
||||
}
|
||||
}
|
||||
}
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 0) {
|
||||
fprintf(stderr, "Crossing info gathered:\n");
|
||||
@ -1968,7 +2084,7 @@ static void dissolve_symedge(CDT_state *cdt, SymEdge *se)
|
||||
if (se->face->symedge == se) {
|
||||
se->face->symedge = se->next;
|
||||
}
|
||||
if (symse->face->symedge == se) {
|
||||
if (symse->face->symedge == symse) {
|
||||
symse->face->symedge = symse->next;
|
||||
}
|
||||
}
|
||||
@ -2380,7 +2496,13 @@ CDT_result *BLI_delaunay_2d_cdt_calc(const CDT_input *input, const CDT_output_ty
|
||||
CDTEdge *face_edge;
|
||||
SymEdge *face_symedge;
|
||||
#ifdef DEBUG_CDT
|
||||
int dbg_level = 1;
|
||||
int dbg_level = 0;
|
||||
#endif
|
||||
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level == -1) {
|
||||
write_cdt_input_to_file(input);
|
||||
}
|
||||
#endif
|
||||
|
||||
if ((nv > 0 && input->vert_coords == NULL) || (ne > 0 && input->edges == NULL) ||
|
||||
@ -2446,6 +2568,15 @@ CDT_result *BLI_delaunay_2d_cdt_calc(const CDT_input *input, const CDT_output_ty
|
||||
continue;
|
||||
}
|
||||
add_edge_constraint(cdt, verts[v1], verts[v2], i, NULL);
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 3) {
|
||||
char namebuf[60];
|
||||
sprintf(namebuf, "after edge constraint %d = (%d,%d)\n", i, v1, v2);
|
||||
cdt_draw(cdt, namebuf);
|
||||
dump_cdt(cdt, namebuf);
|
||||
validate_cdt(cdt, true);
|
||||
}
|
||||
#endif
|
||||
}
|
||||
#ifdef DEBUG_CDT
|
||||
if (dbg_level > 2) {
|
||||
@ -2579,20 +2710,37 @@ void BLI_delaunay_2d_cdt_free(CDT_result *result)
|
||||
|
||||
#ifdef DEBUG_CDT
|
||||
|
||||
static void dump_se(const SymEdge *se, const char *lab)
|
||||
static const char *vertname(const CDTVert *v)
|
||||
{
|
||||
if (se->next) {
|
||||
fprintf(
|
||||
stderr, "%s((%.3f,%.3f)->(%.3f,%.3f))\n", lab, F2(se->vert->co), F2(se->next->vert->co));
|
||||
static char vertnamebuf[20];
|
||||
|
||||
if (v->index < 4) {
|
||||
sprintf(vertnamebuf, "[%c]", "ABCD"[v->index]);
|
||||
}
|
||||
else {
|
||||
fprintf(stderr, "%s((%.3f,%.3f)->NULL)\n", lab, F2(se->vert->co));
|
||||
sprintf(vertnamebuf, "[%d]", v->index - 4);
|
||||
}
|
||||
return vertnamebuf;
|
||||
}
|
||||
|
||||
static void dump_v(const CDTVert *v, const char *lab)
|
||||
{
|
||||
fprintf(stderr, "%s(%.3f,%.3f)\n", lab, F2(v->co));
|
||||
fprintf(stderr, "%s%s(%.3f,%.3f)\n", lab, vertname(v), F2(v->co));
|
||||
}
|
||||
|
||||
static void dump_se(const SymEdge *se, const char *lab)
|
||||
{
|
||||
if (se->next) {
|
||||
fprintf(stderr,
|
||||
"%s%s((%.3f,%.3f)->(%.3f,%.3f))\n",
|
||||
lab,
|
||||
vertname(se->vert),
|
||||
F2(se->vert->co),
|
||||
F2(se->next->vert->co));
|
||||
}
|
||||
else {
|
||||
fprintf(stderr, "%s%s((%.3f,%.3f)->NULL)\n", lab, vertname(se->vert), F2(se->vert->co));
|
||||
}
|
||||
}
|
||||
|
||||
static void dump_se_cycle(const SymEdge *se, const char *lab, const int limit)
|
||||
@ -2622,21 +2770,12 @@ static void dump_id_list(const LinkNode *id_list, const char *lab)
|
||||
}
|
||||
}
|
||||
|
||||
static const char *vertname(CDTVert *v)
|
||||
{
|
||||
static char vertnamebuf[20];
|
||||
|
||||
if (v->index < 4) {
|
||||
sprintf(vertnamebuf, "[%c]", "ABCD"[v->index]);
|
||||
}
|
||||
else {
|
||||
sprintf(vertnamebuf, "[%d]", v->index - 4);
|
||||
}
|
||||
return vertnamebuf;
|
||||
}
|
||||
|
||||
/* If filter_fn != NULL, only dump vert v its edges when filter_fn(cdt, v, filter_data) is true. */
|
||||
# define PL(p) (POINTER_AS_UINT(p) & 0xFFFF)
|
||||
static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
static void dump_cdt_filtered(const CDT_state *cdt,
|
||||
bool (*filter_fn)(const CDT_state *, int, void *),
|
||||
void *filter_data,
|
||||
const char *lab)
|
||||
{
|
||||
LinkNode *ln;
|
||||
CDTVert *v, *vother;
|
||||
@ -2648,6 +2787,9 @@ static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
fprintf(stderr, "\nCDT %s\n", lab);
|
||||
fprintf(stderr, "\nVERTS\n");
|
||||
for (i = 0; i < cdt->vert_array_len; i++) {
|
||||
if (filter_fn && !filter_fn(cdt, i, filter_data)) {
|
||||
continue;
|
||||
}
|
||||
v = cdt->vert_array[i];
|
||||
fprintf(stderr, "%s %x: (%f,%f) symedge=%x\n", vertname(v), PL(v), F2(v->co), PL(v->symedge));
|
||||
dump_id_list(v->input_ids, " ");
|
||||
@ -2664,6 +2806,9 @@ static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
fprintf(stderr, "\n");
|
||||
}
|
||||
}
|
||||
if (filter_fn) {
|
||||
return;
|
||||
}
|
||||
fprintf(stderr, "\nEDGES\n");
|
||||
for (ln = cdt->edges; ln; ln = ln->next) {
|
||||
e = (CDTEdge *)ln->link;
|
||||
@ -2709,6 +2854,55 @@ static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
}
|
||||
# undef PL
|
||||
|
||||
static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
{
|
||||
dump_cdt_filtered(cdt, NULL, NULL, lab);
|
||||
}
|
||||
|
||||
typedef struct ReachableFilterData {
|
||||
int vstart_index;
|
||||
int maxdist;
|
||||
} ReachableFilterData;
|
||||
|
||||
/* Stupid algorithm will repeat itself. Don't use for large n. */
|
||||
static bool reachable_filter(const CDT_state *cdt, int v_index, void *filter_data)
|
||||
{
|
||||
CDTVert *v;
|
||||
SymEdge *se;
|
||||
ReachableFilterData *rfd_in = (ReachableFilterData *)filter_data;
|
||||
ReachableFilterData rfd_next;
|
||||
|
||||
if (v_index == rfd_in->vstart_index) {
|
||||
return true;
|
||||
}
|
||||
if (rfd_in->maxdist <= 0 || v_index < 0 || v_index >= cdt->vert_array_len) {
|
||||
return false;
|
||||
}
|
||||
else {
|
||||
v = cdt->vert_array[v_index];
|
||||
se = v->symedge;
|
||||
if (se != NULL) {
|
||||
rfd_next.vstart_index = rfd_in->vstart_index;
|
||||
rfd_next.maxdist = rfd_in->maxdist - 1;
|
||||
do {
|
||||
if (reachable_filter(cdt, se->next->vert->index, &rfd_next)) {
|
||||
return true;
|
||||
}
|
||||
se = se->rot;
|
||||
} while (se != v->symedge);
|
||||
}
|
||||
}
|
||||
return false;
|
||||
}
|
||||
|
||||
static void dump_cdt_vert_neighborhood(CDT_state *cdt, int v, int maxdist, const char *lab)
|
||||
{
|
||||
ReachableFilterData rfd;
|
||||
rfd.vstart_index = v;
|
||||
rfd.maxdist = maxdist;
|
||||
dump_cdt_filtered(cdt, reachable_filter, &rfd, lab);
|
||||
}
|
||||
|
||||
/**
|
||||
* Make an html file with svg in it to display the argument cdt.
|
||||
* Mouse-overs will reveal the coordinates of vertices and edges.
|
||||
@ -2716,8 +2910,12 @@ static void dump_cdt(const CDT_state *cdt, const char *lab)
|
||||
* The first call creates DRAWFILE; subsequent calls append to it.
|
||||
*/
|
||||
# define DRAWFILE "/tmp/debug_draw.html"
|
||||
# define MAX_DRAW_WIDTH 1000
|
||||
# define MAX_DRAW_HEIGHT 700
|
||||
# define MAX_DRAW_WIDTH 2000
|
||||
# define MAX_DRAW_HEIGHT 1400
|
||||
# define THIN_LINE 1
|
||||
# define THICK_LINE 3
|
||||
# define VERT_RADIUS 3
|
||||
# define DRAW_VERT_LABELS 1
|
||||
static void cdt_draw(CDT_state *cdt, const char *lab)
|
||||
{
|
||||
static bool append = false;
|
||||
@ -2727,9 +2925,7 @@ static void cdt_draw(CDT_state *cdt, const char *lab)
|
||||
double maxx = cdt->maxx + draw_margin;
|
||||
double miny = cdt->miny - draw_margin;
|
||||
double maxy = cdt->maxy + draw_margin;
|
||||
double width = maxx - minx;
|
||||
double height = maxy - miny;
|
||||
double aspect = height / width;
|
||||
double width, height, aspect;
|
||||
int view_width, view_height;
|
||||
double scale;
|
||||
LinkNode *ln;
|
||||
@ -2737,6 +2933,10 @@ static void cdt_draw(CDT_state *cdt, const char *lab)
|
||||
CDTEdge *e;
|
||||
int i, strokew;
|
||||
|
||||
/* Note: to debug a small area: assign custom min's/max's here. */
|
||||
width = maxx - minx;
|
||||
height = maxy - miny;
|
||||
aspect = height / width;
|
||||
view_width = MAX_DRAW_WIDTH;
|
||||
view_height = (int)(view_width * aspect);
|
||||
if (view_height > MAX_DRAW_HEIGHT) {
|
||||
@ -2767,28 +2967,36 @@ static void cdt_draw(CDT_state *cdt, const char *lab)
|
||||
}
|
||||
u = e->symedges[0].vert;
|
||||
v = e->symedges[1].vert;
|
||||
strokew = is_constrained_edge(e) ? 5 : 2;
|
||||
strokew = is_constrained_edge(e) ? THICK_LINE : THIN_LINE;
|
||||
fprintf(f,
|
||||
"<line fill=\"none\" stroke=\"black\" stroke-width=\"%d\" "
|
||||
"x1=\"%.3f\" y1=\"%.3f\" x2=\"%.3f\" y2=\"%.3f\">\n",
|
||||
"x1=\"%f\" y1=\"%f\" x2=\"%f\" y2=\"%f\">\n",
|
||||
strokew,
|
||||
SX(u->co[0]),
|
||||
SY(u->co[1]),
|
||||
SX(v->co[0]),
|
||||
SY(v->co[1]));
|
||||
fprintf(
|
||||
f, " <title>(%.3f,%.3f)(%.3f,%.3f)</title>\n", u->co[0], u->co[1], v->co[0], v->co[1]);
|
||||
fprintf(f, " <title>%s", vertname(u));
|
||||
fprintf(f, "%s</title>\n", vertname(v));
|
||||
fprintf(f, "</line>\n");
|
||||
}
|
||||
i = cdt->output_prepared ? NUM_BOUND_VERTS : 0;
|
||||
for (; i < cdt->vert_array_len; i++) {
|
||||
v = cdt->vert_array[i];
|
||||
fprintf(f,
|
||||
"<circle fill=\"black\" cx=\"%.3f\" cy=\"%.3f\" r=\"5\">\n",
|
||||
"<circle fill=\"black\" cx=\"%f\" cy=\"%f\" r=\"%d\">\n",
|
||||
SX(v->co[0]),
|
||||
SY(v->co[1]));
|
||||
fprintf(f, " <title>(%.3f,%.3f)</title>\n", v->co[0], v->co[1]);
|
||||
SY(v->co[1]),
|
||||
VERT_RADIUS);
|
||||
fprintf(f, " <title>%s(%.3f,%.3f)</title>\n", vertname(v), v->co[0], v->co[1]);
|
||||
fprintf(f, "</circle>\n");
|
||||
# if DRAW_VERT_LABELS
|
||||
fprintf(f,
|
||||
"<text x=\"%f\" y=\"%f\" font-size=\"small\">%s</text>\n",
|
||||
SX(v->co[0]) + VERT_RADIUS,
|
||||
SY(v->co[1]) - VERT_RADIUS,
|
||||
vertname(v));
|
||||
# endif
|
||||
}
|
||||
|
||||
fprintf(f, "</svg>\n</div>\n");
|
||||
@ -2798,6 +3006,28 @@ static void cdt_draw(CDT_state *cdt, const char *lab)
|
||||
# undef SY
|
||||
}
|
||||
|
||||
# define CDTFILE "/tmp/cdtinput.txt"
|
||||
static void write_cdt_input_to_file(const CDT_input *inp)
|
||||
{
|
||||
int i, j;
|
||||
FILE *f = fopen(CDTFILE, "w");
|
||||
|
||||
fprintf(f, "%d %d %d\n", inp->verts_len, inp->edges_len, inp->faces_len);
|
||||
for (i = 0; i < inp->verts_len; i++) {
|
||||
fprintf(f, "%.17f %.17f\n", inp->vert_coords[i][0], inp->vert_coords[i][1]);
|
||||
}
|
||||
for (i = 0; i < inp->edges_len; i++) {
|
||||
fprintf(f, "%d %d\n", inp->edges[i][0], inp->edges[i][1]);
|
||||
}
|
||||
for (i = 0; i < inp->faces_len; i++) {
|
||||
for (j = 0; j < inp->faces_len_table[i]; j++) {
|
||||
fprintf(f, "%d ", inp->faces[j + inp->faces_start_table[i]]);
|
||||
}
|
||||
fprintf(f, "\n");
|
||||
}
|
||||
fclose(f);
|
||||
}
|
||||
|
||||
# ifndef NDEBUG /* Only used in assert. */
|
||||
/**
|
||||
* Is a visible from b: i.e., ab crosses no edge of cdt?
|
||||
|
||||
@ -15,6 +15,10 @@ extern "C" {
|
||||
#include <fstream>
|
||||
#include <sstream>
|
||||
|
||||
#define DO_REGULAR_TESTS 1
|
||||
#define DO_RANDOM_TESTS 0
|
||||
#define DO_FILE_TESTS 0
|
||||
|
||||
static void fill_input_verts(CDT_input *r_input, float (*vcos)[2], int nverts)
|
||||
{
|
||||
r_input->verts_len = nverts;
|
||||
@ -43,6 +47,117 @@ static void add_input_faces(
|
||||
r_input->faces_len_table = faces_len_table;
|
||||
}
|
||||
|
||||
/* The spec should have the form:
|
||||
* #verts #edges #faces
|
||||
* <float> <float> [#verts lines)
|
||||
* <int> <int> [#edges lines]
|
||||
* <int> <int> ... <int> [#faces lines]
|
||||
*/
|
||||
static void fill_input_from_string(CDT_input *r_input, const char *spec)
|
||||
{
|
||||
std::string line;
|
||||
std::vector<std::vector<int>> faces;
|
||||
int i, j;
|
||||
int nverts, nedges, nfaces;
|
||||
float(*p)[2];
|
||||
int(*e)[2];
|
||||
int *farr;
|
||||
int *flen;
|
||||
int *fstart;
|
||||
|
||||
std::istringstream ss(spec);
|
||||
getline(ss, line);
|
||||
std::istringstream hdrss(line);
|
||||
hdrss >> nverts >> nedges >> nfaces;
|
||||
if (nverts == 0) {
|
||||
return;
|
||||
}
|
||||
p = (float(*)[2])MEM_malloc_arrayN(nverts, 2 * sizeof(float), __func__);
|
||||
if (nedges > 0) {
|
||||
e = (int(*)[2])MEM_malloc_arrayN(nedges, 2 * sizeof(int), __func__);
|
||||
}
|
||||
if (nfaces > 0) {
|
||||
flen = (int *)MEM_malloc_arrayN(nfaces, sizeof(int), __func__);
|
||||
fstart = (int *)MEM_malloc_arrayN(nfaces, sizeof(int), __func__);
|
||||
}
|
||||
i = 0;
|
||||
while (i < nverts && getline(ss, line)) {
|
||||
std::istringstream iss(line);
|
||||
iss >> p[i][0] >> p[i][1];
|
||||
i++;
|
||||
}
|
||||
i = 0;
|
||||
while (i < nedges && getline(ss, line)) {
|
||||
std::istringstream ess(line);
|
||||
ess >> e[i][0] >> e[i][1];
|
||||
i++;
|
||||
}
|
||||
i = 0;
|
||||
while (i < nfaces && getline(ss, line)) {
|
||||
std::istringstream fss(line);
|
||||
int v;
|
||||
faces.push_back(std::vector<int>());
|
||||
while (fss >> v) {
|
||||
faces[i].push_back(v);
|
||||
}
|
||||
i++;
|
||||
}
|
||||
fill_input_verts(r_input, p, nverts);
|
||||
if (nedges > 0) {
|
||||
add_input_edges(r_input, e, nedges);
|
||||
}
|
||||
if (nfaces > 0) {
|
||||
for (i = 0; i < nfaces; i++) {
|
||||
flen[i] = (int)faces[i].size();
|
||||
if (i == 0) {
|
||||
fstart[i] = 0;
|
||||
}
|
||||
else {
|
||||
fstart[i] = fstart[i - 1] + flen[i - 1];
|
||||
}
|
||||
}
|
||||
farr = (int *)MEM_malloc_arrayN(fstart[nfaces - 1] + flen[nfaces - 1], sizeof(int), __func__);
|
||||
for (i = 0; i < nfaces; i++) {
|
||||
for (j = 0; j < (int)faces[i].size(); j++) {
|
||||
farr[fstart[i] + j] = faces[i][j];
|
||||
}
|
||||
}
|
||||
add_input_faces(r_input, farr, fstart, flen, nfaces);
|
||||
}
|
||||
}
|
||||
|
||||
static void fill_input_from_file(CDT_input *in, const char *filename)
|
||||
{
|
||||
std::FILE *fp = std::fopen(filename, "rb");
|
||||
if (fp) {
|
||||
std::string contents;
|
||||
std::fseek(fp, 0, SEEK_END);
|
||||
contents.resize(std::ftell(fp));
|
||||
std::rewind(fp);
|
||||
std::fread(&contents[0], 1, contents.size(), fp);
|
||||
std::fclose(fp);
|
||||
fill_input_from_string(in, contents.c_str());
|
||||
}
|
||||
else {
|
||||
printf("couldn't open file %s\n", filename);
|
||||
}
|
||||
}
|
||||
|
||||
static void free_spec_arrays(CDT_input *in)
|
||||
{
|
||||
if (in->vert_coords) {
|
||||
MEM_freeN(in->vert_coords);
|
||||
}
|
||||
if (in->edges) {
|
||||
MEM_freeN(in->edges);
|
||||
}
|
||||
if (in->faces_len_table) {
|
||||
MEM_freeN(in->faces_len_table);
|
||||
MEM_freeN(in->faces_start_table);
|
||||
MEM_freeN(in->faces);
|
||||
}
|
||||
}
|
||||
|
||||
/* which output vert index goes with given input vertex? -1 if not found */
|
||||
static int get_output_vert_index(const CDT_result *r, int in_index)
|
||||
{
|
||||
@ -187,6 +302,7 @@ static void dump_result(CDT_result *r)
|
||||
}
|
||||
}
|
||||
|
||||
#if DO_REGULAR_TESTS
|
||||
TEST(delaunay, Empty)
|
||||
{
|
||||
CDT_input in;
|
||||
@ -692,7 +808,108 @@ TEST(delaunay, TriInTri)
|
||||
BLI_delaunay_2d_cdt_free(out);
|
||||
}
|
||||
|
||||
#if 0
|
||||
TEST(delaunay, DiamondInSquare)
|
||||
{
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
const char *spec = R"(8 0 2
|
||||
0.0 0.0
|
||||
1.0 0.0
|
||||
1.0 1.0
|
||||
0.0 1.0
|
||||
0.14644660940672627 0.5
|
||||
0.5 0.14644660940672627
|
||||
0.8535533905932737 0.5
|
||||
0.5 0.8535533905932737
|
||||
0 1 2 3
|
||||
4 5 6 7
|
||||
)";
|
||||
fill_input_from_string(&in, spec);
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_CONSTRAINTS_VALID_BMESH);
|
||||
EXPECT_EQ(out->verts_len, 8);
|
||||
EXPECT_EQ(out->edges_len, 10);
|
||||
EXPECT_EQ(out->faces_len, 3);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
|
||||
TEST(delaunay, DiamondInSquareWire)
|
||||
{
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
const char *spec = R"(8 8 0
|
||||
0.0 0.0
|
||||
1.0 0.0
|
||||
1.0 1.0
|
||||
0.0 1.0
|
||||
0.14644660940672627 0.5
|
||||
0.5 0.14644660940672627
|
||||
0.8535533905932737 0.5
|
||||
0.5 0.8535533905932737
|
||||
0 1
|
||||
1 2
|
||||
2 3
|
||||
3 0
|
||||
4 5
|
||||
5 6
|
||||
6 7
|
||||
7 4
|
||||
)";
|
||||
fill_input_from_string(&in, spec);
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_CONSTRAINTS);
|
||||
EXPECT_EQ(out->verts_len, 8);
|
||||
EXPECT_EQ(out->edges_len, 8);
|
||||
EXPECT_EQ(out->faces_len, 2);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
|
||||
TEST(delaunay, ClosePts)
|
||||
{
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
const char *spec = R"(7 2 1
|
||||
0.46876350045204163 0.06087132915854454
|
||||
0.46865847706794739 0.03632887825369835
|
||||
0.49176687002182007 0.03632888197898865
|
||||
0.49166208505630493 0.06087132543325424
|
||||
0.49171400070190430 0.04841339960694313
|
||||
0.49171534180641174 0.04839951172471046
|
||||
0.49045535922050476 0.06087132915854454
|
||||
4 5
|
||||
6 4
|
||||
0 1 2 3
|
||||
)";
|
||||
fill_input_from_string(&in, spec);
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_FULL);
|
||||
EXPECT_EQ(out->verts_len, 7);
|
||||
EXPECT_EQ(out->edges_len, 12);
|
||||
EXPECT_EQ(out->faces_len, 6);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
|
||||
TEST(delaunay, ClosePts2)
|
||||
{
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
const char *spec = R"(6 1 1
|
||||
-0.17878936231136322 -0.44374340772628784
|
||||
-0.17871695756912231 -0.45601493120193481
|
||||
-0.17544384300708771 -0.45601493120193481
|
||||
-0.17537136375904083 -0.44374340772628784
|
||||
-0.17544738948345184 -0.45602506399154663
|
||||
-0.17872454226016998 -0.45472940802574158
|
||||
4 5
|
||||
0 1 2 3
|
||||
)";
|
||||
fill_input_from_string(&in, spec);
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_FULL);
|
||||
EXPECT_EQ(out->verts_len, 6);
|
||||
EXPECT_EQ(out->edges_len, 10);
|
||||
EXPECT_EQ(out->faces_len, 5);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
#endif
|
||||
|
||||
#if DO_RANDOM_TESTS
|
||||
enum {
|
||||
RANDOM_PTS,
|
||||
RANDOM_SEGS,
|
||||
@ -780,17 +997,17 @@ static void rand_delaunay_test(int test_kind,
|
||||
|
||||
TEST(delaunay, randompts)
|
||||
{
|
||||
rand_delaunay_test(RANDOM_PTS, 7, 1, CDT_FULL);
|
||||
rand_delaunay_test(RANDOM_PTS, 7, 100, CDT_FULL);
|
||||
}
|
||||
|
||||
TEST(delaunay, randomsegs)
|
||||
{
|
||||
rand_delaunay_test(RANDOM_SEGS, 7, 1, CDT_FULL);
|
||||
rand_delaunay_test(RANDOM_SEGS, 7, 100, CDT_FULL);
|
||||
}
|
||||
|
||||
TEST(delaunay, randompoly)
|
||||
{
|
||||
rand_delaunay_test(RANDOM_POLY, 7, 1, CDT_FULL);
|
||||
rand_delaunay_test(RANDOM_POLY, 7, 100, CDT_FULL);
|
||||
}
|
||||
|
||||
TEST(delaunay, randompoly_inside)
|
||||
@ -809,47 +1026,35 @@ TEST(delaunay, randompoly_validbmesh)
|
||||
}
|
||||
#endif
|
||||
|
||||
#if 0
|
||||
/* For debugging or timing large examples.
|
||||
* The given file should have one point per line, as a space-separated pair of floats
|
||||
#if DO_FILE_TESTS
|
||||
/* For timing large examples of points only.
|
||||
* The given file should have one point per line, as a space-separated pair of floats.
|
||||
*/
|
||||
static void points_from_file_test(const char *filename)
|
||||
{
|
||||
std::ifstream f;
|
||||
std::string line;
|
||||
struct XY {
|
||||
float x;
|
||||
float y;
|
||||
} xy;
|
||||
std::vector<XY> pts;
|
||||
int i, n;
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
double tstart;
|
||||
float (*p)[2];
|
||||
|
||||
f.open(filename);
|
||||
while (getline(f, line)) {
|
||||
std::istringstream iss(line);
|
||||
iss >> xy.x >> xy.y;
|
||||
pts.push_back(xy);
|
||||
}
|
||||
n = (int)pts.size();
|
||||
fprintf(stderr, "read %d points\n", (int)pts.size());
|
||||
p = (float (*)[2])MEM_malloc_arrayN(n, 2 * sizeof(float), "delaunay");
|
||||
for (i = 0; i < n; i++) {
|
||||
xy = pts[i];
|
||||
p[i][0] = xy.x;
|
||||
p[i][1] = xy.y;
|
||||
}
|
||||
fill_input_from_file(&in, filename);
|
||||
tstart = PIL_check_seconds_timer();
|
||||
fill_input_verts(&in, p, n);
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_FULL);
|
||||
fprintf(stderr, "time to triangulate=%f seconds\n", PIL_check_seconds_timer() - tstart);
|
||||
BLI_delaunay_2d_cdt_free(out);
|
||||
MEM_freeN(p);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
|
||||
TEST(delaunay, debug)
|
||||
{
|
||||
CDT_input in;
|
||||
CDT_result *out;
|
||||
fill_input_from_file(&in, "/tmp/cdtinput.txt");
|
||||
out = BLI_delaunay_2d_cdt_calc(&in, CDT_FULL);
|
||||
BLI_delaunay_2d_cdt_free(out);
|
||||
free_spec_arrays(&in);
|
||||
}
|
||||
|
||||
# if 0
|
||||
# define POINTFILEROOT "/tmp/"
|
||||
|
||||
TEST(delaunay, terrain1)
|
||||
@ -867,3 +1072,4 @@ TEST(delaunay, terrain3)
|
||||
points_from_file_test(POINTFILEROOT "points3.txt");
|
||||
}
|
||||
# endif
|
||||
#endif
|
||||
|
||||
Loading…
Reference in New Issue
Block a user