SaraP
cbec90699f
- creato il progetto in Visual Studio per compilare come libreria statica - modifiche al codice originale per integrarlo nelle nostre librerie.
574 lines
21 KiB
C++
574 lines
21 KiB
C++
/*****************************************************************************/
|
|
/* */
|
|
/* F I S T : Fast, Industrial-Strength Triangulation */
|
|
/* */
|
|
/*****************************************************************************/
|
|
/* */
|
|
/* (C) Martin Held */
|
|
/* (C) Universitaet Salzburg, Salzburg, Austria */
|
|
/* */
|
|
/* This code is not in the public domain. All rights reserved! Please make */
|
|
/* sure to read the full copyright statement contained in api_functions.cpp. */
|
|
/* */
|
|
/*****************************************************************************/
|
|
|
|
/* */
|
|
/* get standard libraries */
|
|
/* */
|
|
#include <stdio.h>
|
|
#include <stdlib.h>
|
|
#include <math.h>
|
|
#include <assert.h>
|
|
|
|
/* */
|
|
/* get my header files */
|
|
/* */
|
|
#include "fpkernel.h"
|
|
#include "martin.h"
|
|
#include "defs.h"
|
|
#include "list.h"
|
|
#include "header.h"
|
|
#include "basic.h"
|
|
#include "numerics.h"
|
|
#include "bv_tree.h"
|
|
#ifdef GRAPHICS
|
|
#include "graphics.h"
|
|
#endif
|
|
|
|
/* */
|
|
/* function prototypes of functions provided in this file */
|
|
/* */
|
|
boolean Desperate(global_struct *all, heapdef *hp, list_ind ind, int i, boolean *splitted);
|
|
boolean ExistsCrossOver(global_struct *all, list_ind ind, list_ind *ind1,
|
|
int *i1, list_ind *ind2, int *i2, list_ind *ind3,
|
|
int *i3, list_ind *ind4, int *i4);
|
|
void HandleCrossOver(global_struct *all, heapdef *hp, list_ind ind1, int i1, list_ind ind2,
|
|
int i2, list_ind ind3, int i3, list_ind ind4, int i4);
|
|
boolean LetsHope(global_struct *all, heapdef *hp, list_ind ind);
|
|
void HandleSplit(global_struct *all, list_ind ind1, int i1, list_ind ind3,
|
|
int i3);
|
|
boolean FoundSplit(global_struct *all, list_ind ind5, int i5, list_ind ind,
|
|
list_ind ind1, int i1, int i3, int i4, list_ind *ind2,
|
|
int *i2);
|
|
boolean ExistsSplit(global_struct *all, list_ind ind, list_ind *ind1, int *i1, list_ind *ind2,
|
|
int *i2);
|
|
void FreeDistances(bridgedef *bridge);
|
|
int WindingNumber(global_struct *all, list_ind ind, point p);
|
|
|
|
|
|
/* */
|
|
/* the functions in this file try to ensure that we always end up with */
|
|
/* something that (topologically) is a triangulation. */
|
|
/* */
|
|
/* the more desperate we get, the more aggressive means we choose for making */
|
|
/* diagonals "valid". */
|
|
/* */
|
|
|
|
|
|
boolean Desperate(global_struct *all, heapdef *hp, list_ind ind, int i, boolean *splitted)
|
|
{
|
|
griddef *grid = &all->c_grid;
|
|
|
|
int i1, i2, i3, i4;
|
|
list_ind ind1, ind2, ind3, ind4;
|
|
|
|
*splitted = false;
|
|
|
|
/* */
|
|
/* check whether there exist consecutive vertices i1, i2, i3, i4 such */
|
|
/* that i1, i2 and i3, i4 intersect */
|
|
/* */
|
|
if (ExistsCrossOver(all, ind, &ind1, &i1, &ind2, &i2, &ind3, &i3, &ind4, &i4)) {
|
|
/* */
|
|
/* insert two new diagonals around the cross-over without checking */
|
|
/* whether they are intersection-free */
|
|
/* */
|
|
HandleCrossOver(all, hp, ind1, i1, ind2, i2, ind3, i3, ind4, i4);
|
|
return false;
|
|
}
|
|
|
|
if (grid->no_grid_yet) {
|
|
BuildGrid(all, i, i+1);
|
|
}
|
|
|
|
/* */
|
|
/* check whether there exists a valid diagonal that splits the polygon */
|
|
/* into two parts */
|
|
/* */
|
|
if (ExistsSplit(all, ind, &ind1, &i1, &ind2, &i2)) {
|
|
/* */
|
|
/* break up the polygon by inserting this diagonal (which can't be an */
|
|
/* ear -- otherwise, we would not have ended up in this part of the */
|
|
/* code). then, let's treat the two polygons separately. hopefully, */
|
|
/* this will help to handle self-overlapping polygons in the "correct" */
|
|
/* way. */
|
|
/* */
|
|
HandleSplit(all, ind1, i1, ind2, i2);
|
|
*splitted = true;
|
|
return false;
|
|
}
|
|
|
|
return true;
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
boolean ExistsCrossOver(global_struct *all, list_ind ind, list_ind *ind1, int *i1, list_ind *ind2,
|
|
int *i2, list_ind *ind3, int *i3, list_ind *ind4,
|
|
int *i4)
|
|
{
|
|
datadef *data = &all->c_data;
|
|
listdef *list = &all->c_list;
|
|
|
|
bounding_box bb1, bb2;
|
|
|
|
*ind1 = ind;
|
|
*i1 = GetIndex(list, *ind1);
|
|
*ind2 = GetNextNode(list, *ind1);
|
|
*i2 = GetIndex(list, *ind2);
|
|
*ind3 = GetNextNode(list, *ind2);
|
|
*i3 = GetIndex(list, *ind3);
|
|
*ind4 = GetNextNode(list, *ind3);
|
|
*i4 = GetIndex(list, *ind4);
|
|
do {
|
|
BBox(data, *i1, *i2, bb1);
|
|
BBox(data, *i3, *i4, bb2);
|
|
if (BBox_Overlap(bb1, bb2)) {
|
|
if (SegIntersect(all, bb1.imin, bb1.imax, bb2.imin, bb2.imax, -1))
|
|
return true;
|
|
}
|
|
*ind1 = *ind2;
|
|
*i1 = *i2;
|
|
*ind2 = *ind3;
|
|
*i2 = *i3;
|
|
*ind3 = *ind4;
|
|
*i3 = *i4;
|
|
*ind4 = GetNextNode(list, *ind3);
|
|
*i4 = GetIndex(list, *ind4);
|
|
} while (*ind1 != ind);
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|
|
|
|
void HandleCrossOver(global_struct *all, heapdef *hp, list_ind ind1, int i1,
|
|
list_ind ind2, int i2,
|
|
list_ind ind3, int i3,
|
|
list_ind ind4, int i4)
|
|
{
|
|
datadef *data = &all->c_data;
|
|
listdef *list = &all->c_list;
|
|
|
|
machine_double ratio1, ratio4;
|
|
boolean first;
|
|
int angle1, angle4;
|
|
|
|
/* */
|
|
/* which pair of triangles shall I insert?? we can use either i1, i2, i3 */
|
|
/* and i1, i3, i4, or we can use i2, i3, i4 and i1, i2, i4... */
|
|
/* */
|
|
angle1 = GetAngle(list, ind1);
|
|
angle4 = GetAngle(list, ind4);
|
|
if (angle1 < angle4) first = true;
|
|
else if (angle1 > angle4) first = false;
|
|
else if (hp->ears_sorted) {
|
|
ratio1 = GetRatio(data, i3, i4, i1);
|
|
ratio4 = GetRatio(data, i1, i2, i4);
|
|
if (ratio4 < ratio1) first = false;
|
|
else first = true;
|
|
}
|
|
else {
|
|
first = true;
|
|
}
|
|
|
|
if (first) {
|
|
/* */
|
|
/* first clip i1, i2, i3, then clip i1, i3, i4 */
|
|
/* */
|
|
DeleteLinks(list, ind2);
|
|
StoreTriangle(all, ind1, ind2, ind3,TriTriColor);
|
|
SetAngle(list, ind3, 1);
|
|
InsertIntoHeap(hp, 0.0, ind3, ind1, ind4);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
AddTriToBuffer(&all->c_redraw, i1, i2, i3, TriColor, TriFillColor);
|
|
DrawTri(data->points[i1], data->points[i2], data->points[i3],
|
|
TriColor, TriFillColor);
|
|
DrawPnt(data->points[i3], EarColor);
|
|
}
|
|
#endif
|
|
}
|
|
else {
|
|
/* */
|
|
/* first clip i2, i3, i4, then clip i1, i2, i4 */
|
|
/* */
|
|
DeleteLinks(list, ind3);
|
|
StoreTriangle(all, ind2, ind3, ind4, TriTriColor);
|
|
SetAngle(list, ind2, 1);
|
|
InsertIntoHeap(hp, 0.0, ind2, ind1, ind4);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
AddTriToBuffer(&all->c_redraw, i2, i3, i4, TriColor, TriFillColor);
|
|
DrawTri(data->points[i2], data->points[i3], data->points[i4],
|
|
TriColor, TriFillColor);
|
|
DrawPnt(data->points[i2], EarColor);
|
|
}
|
|
#endif
|
|
}
|
|
|
|
return;
|
|
}
|
|
|
|
|
|
|
|
|
|
boolean LetsHope(global_struct *all, heapdef *hp, list_ind ind)
|
|
{
|
|
listdef *list = &all->c_list;
|
|
#ifdef GRAPHICS
|
|
datadef *data = &all->c_data;
|
|
int i1;
|
|
#endif
|
|
|
|
list_ind ind0, ind1, ind2;
|
|
|
|
/* */
|
|
/* perhaps we find a 180-degree corner? such a corner is not the best of */
|
|
/* all ears. but, at least, we'd reduce the vertex count by one without */
|
|
/* causing any real harm... */
|
|
/* */
|
|
ind1 = ind;
|
|
do {
|
|
if (GetAngle(list, ind1) == 0) {
|
|
SetAngle(list, ind1, 1);
|
|
ind0 = GetPrevNode(list, ind1);
|
|
ind2 = GetNextNode(list, ind1);
|
|
InsertIntoHeap(hp, 0.0, ind1, ind0, ind2);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
i1 = GetIndex(list, ind1);
|
|
DrawPnt(data->points[i1], EarColor);
|
|
}
|
|
#endif
|
|
return true;
|
|
}
|
|
ind1 = GetNextNode(list, ind1);
|
|
} while (ind1 != ind);
|
|
|
|
/* */
|
|
/* let's clip the first convex corner. of course, we know that this is no */
|
|
/* ear in an ideal world. but this polygon isn't ideal, either! */
|
|
/* */
|
|
ind1 = ind;
|
|
do {
|
|
if (GetAngle(list, ind1) > 0) {
|
|
ind0 = GetPrevNode(list, ind1);
|
|
ind2 = GetNextNode(list, ind1);
|
|
InsertIntoHeap(hp, 0.0, ind1, ind0, ind2);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
i1 = GetIndex(list, ind1);
|
|
DrawPnt(data->points[i1], EarColor);
|
|
}
|
|
#endif
|
|
return true;
|
|
}
|
|
ind1 = GetNextNode(list, ind1);
|
|
} while (ind1 != ind);
|
|
|
|
/* */
|
|
/* no convex and tangential corners? so, let's cheat! this code ain't */
|
|
/* stop without some triangulation... ;-) g-i-g-o? right! perhaps, */
|
|
/* this is what you call a robust code?! */
|
|
/* */
|
|
SetAngle(list, ind, 1);
|
|
ind0 = GetPrevNode(list, ind);
|
|
ind2 = GetNextNode(list, ind);
|
|
InsertIntoHeap(hp, 0.0, ind, ind0, ind2);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
i1 = GetIndex(list, ind);
|
|
DrawPnt(data->points[i1], EarColor);
|
|
}
|
|
#endif
|
|
|
|
return true;
|
|
|
|
/* */
|
|
/* see, we never have to return "false"... ;-) */
|
|
/* */
|
|
/*
|
|
return false;
|
|
*/
|
|
}
|
|
|
|
|
|
|
|
|
|
boolean ExistsSplit(global_struct *all, list_ind ind, list_ind *ind1, int *i1,
|
|
list_ind *ind2, int *i2)
|
|
{
|
|
listdef *list = &all->c_list;
|
|
bridgedef *bridge = &all->c_bridge;
|
|
|
|
list_ind ind3, ind4, ind5;
|
|
int i3, i4, i5, number;
|
|
|
|
number = GetNumList(list);
|
|
if (number > bridge->max_num_dist) {
|
|
bridge->distances = CORE_ReallocateArray(bridge->memptr, bridge->distances,
|
|
bridge->max_num_dist,
|
|
number, distance,
|
|
"desperate:distances");
|
|
bridge->max_num_dist = number;
|
|
}
|
|
*ind1 = ind;
|
|
*i1 = GetIndex(list, *ind1);
|
|
ind4 = GetNextNode(list, *ind1);
|
|
i4 = GetIndex(list, ind4);
|
|
assert(*ind1 != ind4);
|
|
ind5 = GetNextNode(list, ind4);
|
|
i5 = GetIndex(list, ind5);
|
|
assert(*ind1 != *ind2);
|
|
ind3 = GetPrevNode(list, *ind1);
|
|
i3 = GetIndex(list, ind3);
|
|
assert(*ind2 != ind3);
|
|
if (FoundSplit(all, ind5, i5, ind3, *ind1, *i1, i3, i4, ind2, i2))
|
|
return true;
|
|
i3 = *i1;
|
|
*ind1 = ind4;
|
|
*i1 = i4;
|
|
ind4 = ind5;
|
|
i4 = i5;
|
|
ind5 = GetNextNode(list, ind4);
|
|
i5 = GetIndex(list, ind5);
|
|
|
|
while (ind5 != ind) {
|
|
if (FoundSplit(all, ind5, i5, ind, *ind1, *i1, i3, i4, ind2, i2))
|
|
return true;
|
|
i3 = *i1;
|
|
*ind1 = ind4;
|
|
*i1 = i4;
|
|
ind4 = ind5;
|
|
i4 = i5;
|
|
ind5 = GetNextNode(list, ind4);
|
|
i5 = GetIndex(list, ind5);
|
|
}
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|
|
/* */
|
|
/* this function computes the winding number of a polygon with respect to a */
|
|
/* point p. no care is taken to handle cases where p lies on the */
|
|
/* boundary of the polygon. (this is no issue in our application, as we will */
|
|
/* always compute the winding number with respect to the mid-point of a */
|
|
/* valid diagonal.) */
|
|
/* */
|
|
int WindingNumber(global_struct *all, list_ind ind, point p)
|
|
{
|
|
listdef *list = &all->c_list;
|
|
datadef *data = &all->c_data;
|
|
|
|
double angle;
|
|
list_ind ind2;
|
|
int i1, i2, number;
|
|
|
|
i1 = GetIndex(list, ind);
|
|
ind2 = GetNextNode(list, ind);
|
|
i2 = GetIndex(list, ind2);
|
|
angle = Angle(p, data->points[i1], data->points[i2]);
|
|
while (ind2 != ind) {
|
|
i1 = i2;
|
|
ind2 = GetNextNode(list, ind2);
|
|
i2 = GetIndex(list, ind2);
|
|
angle += Angle(p, data->points[i1], data->points[i2]);
|
|
}
|
|
|
|
angle += M_PI;
|
|
number = REAL_TO_INT(angle / M_2PI);
|
|
|
|
return number;
|
|
}
|
|
|
|
|
|
|
|
|
|
boolean FoundSplit(global_struct *all, list_ind ind5, int i5, list_ind ind, list_ind ind1,
|
|
int i1, int i3, int i4, list_ind *ind2, int *i2)
|
|
{
|
|
listdef *list = &all->c_list;
|
|
bridgedef *bridge = &all->c_bridge;
|
|
datadef *data = &all->c_data;
|
|
|
|
tmp_data_def tmp_data;
|
|
|
|
point base, center;
|
|
int num_dist = 0;
|
|
int j, i6, i7;
|
|
list_ind ind6, ind7;
|
|
bounding_box bb;
|
|
boolean convex, cone_ok;
|
|
int d_comp(const void *a, const void *);
|
|
|
|
/* */
|
|
/* sort the points according to their distance from i1 */
|
|
/* */
|
|
do {
|
|
assert(num_dist < bridge->max_num_dist);
|
|
BaseLength(data->points[i1], data->points[i5], base, bridge->distances[num_dist].dist);
|
|
bridge->distances[num_dist].ind = ind5;
|
|
++num_dist;
|
|
ind5 = GetNextNode(list, ind5);
|
|
i5 = GetIndex(list, ind5);
|
|
} while (ind5 != ind);
|
|
|
|
qsort(bridge->distances, num_dist, sizeof(distance), &d_comp);
|
|
|
|
/* */
|
|
/* find a valid diagonal. */
|
|
/* */
|
|
for (j = 0; j < num_dist; ++j) {
|
|
*ind2 = bridge->distances[j].ind;
|
|
*i2 = GetIndex(list, *ind2);
|
|
if (i1 != *i2) {
|
|
ind6 = GetPrevNode(list, *ind2);
|
|
i6 = GetIndex(list, ind6);
|
|
ind7 = GetNextNode(list, *ind2);
|
|
i7 = GetIndex(list, ind7);
|
|
convex = GetAngle(list, *ind2) > 0;
|
|
IsInCone(data, tmp_data, i6, *i2, i7, i1, convex, cone_ok);
|
|
if (cone_ok) {
|
|
convex = GetAngle(list, ind1) > 0;
|
|
IsInCone(data, tmp_data, i3, i1, i4, *i2, convex, cone_ok);
|
|
if (cone_ok) {
|
|
BBox(data, i1, *i2, bb);
|
|
if (!GridIntersectionExists(all, bb, -1, -1, ind1, -1)) {
|
|
/* */
|
|
/* check whether this is a good diagonal; we do not want a */
|
|
/* diagonal that may create figure-8's! */
|
|
/* */
|
|
VectorAdd2D(data->points[i1], data->points[*i2], center);
|
|
MultScalar2D(0.5, center);
|
|
if (WindingNumber(all, ind, center) == 1) return true;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|
|
|
|
void HandleSplit(global_struct *all, list_ind ind1, int i1, list_ind ind3, int i3)
|
|
{
|
|
listdef *list = &all->c_list;
|
|
datadef *data = &all->c_data;
|
|
|
|
list_ind ind2, ind4, prev, next;
|
|
int prv, nxt, angle;
|
|
|
|
/* */
|
|
/* duplicate nodes in order to form end points of the new diagonal */
|
|
/* */
|
|
ind2 = MakeNode(list, i1);
|
|
InsertAfter(list, ind1, ind2);
|
|
SetOriginal(list, ind2, GetOriginal(list, ind1));
|
|
ind4 = MakeNode(list, i3);
|
|
InsertAfter(list, ind3, ind4);
|
|
SetOriginal(list, ind4, GetOriginal(list, ind3));
|
|
|
|
/* */
|
|
/* insert the diagonal into the boundary loop, thus splitting the loop */
|
|
/* into two loops */
|
|
/* */
|
|
SplitSplice(list, ind1, ind2, ind3, ind4);
|
|
|
|
/* */
|
|
/* store pointers to the two new loops */
|
|
/* */
|
|
StoreChain(list, ind1);
|
|
StoreChain(list, ind3);
|
|
|
|
/* */
|
|
/* update the graphics */
|
|
/* */
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
DrawSeg(data->points[i1], data->points[i3], SplitColor);
|
|
AddEdgeToBuffer(&all->c_redraw, i1, i3, BridgeColor);
|
|
}
|
|
#endif
|
|
|
|
/* */
|
|
/* reset the angles */
|
|
/* */
|
|
next = GetNextNode(list, ind1);
|
|
nxt = GetIndex(list, next);
|
|
prev = GetPrevNode(list, ind1);
|
|
prv = GetIndex(list, prev);
|
|
angle = IsConvexAngle(list, data, prv, i1, nxt, ind1);
|
|
SetAngle(list, ind1, angle);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
if (angle == 1) DrawPnt(data->points[i1], ConvexColor);
|
|
else if (angle == 2) DrawPnt(data->points[i1], ZeroColor);
|
|
else if (angle == 0) DrawPnt(data->points[i1], TangentColor);
|
|
else DrawPnt(data->points[i1], PntsColor);
|
|
}
|
|
#endif
|
|
|
|
next = GetNextNode(list, ind2);
|
|
nxt = GetIndex(list, next);
|
|
prev = GetPrevNode(list, ind2);
|
|
prv = GetIndex(list, prev);
|
|
angle = IsConvexAngle(list, data, prv, i1, nxt, ind2);
|
|
SetAngle(list, ind2, angle);
|
|
|
|
next = GetNextNode(list, ind3);
|
|
nxt = GetIndex(list, next);
|
|
prev = GetPrevNode(list, ind3);
|
|
prv = GetIndex(list, prev);
|
|
angle = IsConvexAngle(list, data, prv, i3, nxt, ind3);
|
|
SetAngle(list, ind3, angle);
|
|
#ifdef GRAPHICS
|
|
if (all->rt_opt.graphics) {
|
|
if (angle == 1) DrawPnt(data->points[i3], ConvexColor);
|
|
else if (angle == 2) DrawPnt(data->points[i3], ZeroColor);
|
|
else if (angle == 0) DrawPnt(data->points[i3], TangentColor);
|
|
else DrawPnt(data->points[i1], PntsColor);
|
|
}
|
|
#endif
|
|
|
|
next = GetNextNode(list, ind4);
|
|
nxt = GetIndex(list, next);
|
|
prev = GetPrevNode(list, ind4);
|
|
prv = GetIndex(list, prev);
|
|
angle = IsConvexAngle(list, data, prv, i3, nxt, ind4);
|
|
SetAngle(list, ind4, angle);
|
|
|
|
return;
|
|
}
|
|
|
|
|
|
|
|
|
|
void FreeDistances(bridgedef *bridge)
|
|
{
|
|
bridge->max_num_dist = 0;
|
|
CORE_FreeMemory(bridge->memptr, &bridge->distances, "bridge:distances");
|
|
|
|
return;
|
|
}
|
|
|
|
|