739088af9f
- aggiornamento versione.
1733 lines
54 KiB
C++
1733 lines
54 KiB
C++
/*****************************************************************************/
|
|
/* */
|
|
/* Copyright (C) 1999-2025 M. Held */
|
|
/* */
|
|
/* This code is not in the public domain. All rights reserved! Please make */
|
|
/* sure to read the full copyright statement contained in "README.txt" or in */
|
|
/* the "main" file of this code, such as "main.cc". */
|
|
/* */
|
|
/*****************************************************************************/
|
|
/* */
|
|
/* Written by: Martin Held */
|
|
/* */
|
|
/* E-Mail: held@cs.sbg.ac.at */
|
|
/* Fax Mail: (+43 662) 8044-172 */
|
|
/* Voice Mail: (+43 662) 8044-6304 */
|
|
/* Snail Mail: Martin Held */
|
|
/* FB Informatik */
|
|
/* Universitaet Salzburg */
|
|
/* A-5020 Salzburg, Austria */
|
|
/* */
|
|
/*****************************************************************************/
|
|
|
|
/* */
|
|
/* get standard libraries */
|
|
/* */
|
|
#include <math.h>
|
|
#include <assert.h>
|
|
#include <stdlib.h>
|
|
#include <stdio.h>
|
|
#include <vector>
|
|
#include <set>
|
|
#include <map>
|
|
#include <queue>
|
|
#include <algorithm>
|
|
|
|
/* */
|
|
/* get my header files */
|
|
/* */
|
|
#include "fpkernel.h"
|
|
#include "vronivector.h"
|
|
#include "vroni_object.h"
|
|
#include "numerics.h"
|
|
#include "roots.h"
|
|
#include "intersections.h"
|
|
|
|
|
|
#ifdef GRAPHICS
|
|
#define HandleSplitSeg(I, J, K) \
|
|
{\
|
|
K = SplitSeg(I, J); \
|
|
if (graphics) { \
|
|
AddToCurBuffer(I, SEG, AlertColor); \
|
|
AddToCurBuffer(K, SEG, AlertColor); \
|
|
AddToCurBuffer(J, PNT, AlertColor); \
|
|
}\
|
|
}
|
|
#else
|
|
#define HandleSplitSeg(I, J, K) \
|
|
{\
|
|
K = SplitSeg(I, J); \
|
|
}
|
|
#endif
|
|
|
|
#ifdef GRAPHICS
|
|
#define HandleSplitArc(I, J, K) \
|
|
{\
|
|
K = SplitArc(I, J); \
|
|
if (graphics) { \
|
|
AddToCurBuffer(I, ARC, AlertColor); \
|
|
AddToCurBuffer(K, ARC, AlertColor); \
|
|
AddToCurBuffer(J, PNT, AlertColor); \
|
|
}\
|
|
}
|
|
#else
|
|
#define HandleSplitArc(I, J, K) \
|
|
{\
|
|
K = SplitArc(I, J); \
|
|
}
|
|
#endif
|
|
|
|
|
|
#ifdef DBG_GRAPHICS
|
|
void SketchGrid(grid* grid);
|
|
#endif
|
|
|
|
|
|
|
|
void vroniObject::ResetIntersectionStatus(void)
|
|
{
|
|
check_completed = true;
|
|
eps_int = ZERO;
|
|
|
|
return;
|
|
}
|
|
|
|
|
|
vr_bool vroniObject::IsCheckComplete(void)
|
|
{
|
|
return check_completed;
|
|
}
|
|
|
|
|
|
vr_bool vroniObject::ThresholdReached(void)
|
|
{
|
|
return (eps_int > THRESHOLD);
|
|
}
|
|
|
|
|
|
vr_bool vroniObject::CheckIntersections(void)
|
|
{
|
|
vr_bool int_found;
|
|
|
|
int_found = ComputeAllIntersections(eps_int);
|
|
|
|
if (check_completed && (eps_int <= THRESHOLD)) {
|
|
eps_int *= 10.0;
|
|
}
|
|
|
|
return int_found;
|
|
}
|
|
|
|
|
|
double vroniObject::GetIntersectionThreshold(void)
|
|
{
|
|
return eps_int;
|
|
}
|
|
|
|
vr_bool vroniObject::AddSegIntersection(seg_entry* seg, coord intersection)
|
|
{
|
|
return seg->intersections->insert(intersection).second;
|
|
}
|
|
|
|
|
|
#ifdef GENUINE_ARCS
|
|
vr_bool vroniObject::AddArcIntersection(arc_entry* arc, coord intersection)
|
|
{
|
|
arc_intersection ai;
|
|
int i1;
|
|
coord center, start;
|
|
double alpha;
|
|
|
|
i1 = arc->num;
|
|
center = GetArcCenter(i1);
|
|
start = GetArcStartCoord(i1);
|
|
|
|
alpha = atan2(intersection.y - center.y, intersection.x - center.x) -
|
|
atan2(start.y - center.y, start.x - center.x);
|
|
if (alpha < 0)
|
|
alpha += M_2PI;
|
|
|
|
ai.c = intersection;
|
|
ai.alpha = alpha;
|
|
|
|
return arc->intersections.insert(ai).second;
|
|
}
|
|
#endif
|
|
|
|
|
|
grid* vroniObject::InitialiseGrid(int num_elements, coord bb_min, coord bb_max)
|
|
{
|
|
grid* g;
|
|
int i, j;
|
|
int n_x, n_y;
|
|
double d_x, d_y;
|
|
double min_x, min_y;
|
|
|
|
n_x = (int)(sqrt(machine_double(num_elements))) + 3;
|
|
n_y = n_x;
|
|
d_x = (2 * ZERO_MAX + bb_max.x - bb_min.x)/(n_x - 2);
|
|
d_y = (2 * ZERO_MAX + bb_max.y - bb_min.y)/(n_y - 2);
|
|
min_x = bb_min.x - d_x;
|
|
min_y = bb_min.y - d_y;
|
|
|
|
g = new grid;
|
|
|
|
g->d_x = d_x;
|
|
g->d_y = d_y;
|
|
g->n_x = n_x;
|
|
g->n_y = n_y;
|
|
|
|
g->x_min = new double[n_x];
|
|
g->x_max = new double[n_x];
|
|
g->y_min = new double[n_y];
|
|
g->y_max = new double[n_y];
|
|
|
|
for (i = 0; i < n_x; i++) {
|
|
g->x_min[i] = (i == 0 ? min_x : g->x_max[i-1]); //min_x + i*d_x
|
|
g->x_max[i] = min_x + (i+1)*d_x;
|
|
}
|
|
|
|
for (j = 0; j < n_y; j++) {
|
|
g->y_min[j] = (j == 0 ? min_y : g->y_max[j-1]); //min_y + j*d_y
|
|
g->y_max[j] = min_y + (j+1)*d_y;
|
|
}
|
|
|
|
g->cells = new grid_cell*[n_x];
|
|
for (i = 0; i < n_x; i++) {
|
|
g->cells[i] = new grid_cell[n_y];
|
|
for (j = 0; j < n_y; j++) {
|
|
g->cells[i][j].x_min = g->x_min[i];
|
|
g->cells[i][j].x_max = g->x_max[i];
|
|
g->cells[i][j].y_min = g->y_min[j];
|
|
g->cells[i][j].y_max = g->y_max[j];
|
|
g->cells[i][j].i = i;
|
|
g->cells[i][j].j = j;
|
|
}
|
|
}
|
|
|
|
return g;
|
|
}
|
|
|
|
int vroniObject::FindCellX(grid* g, double_arg x) {
|
|
return REAL_TO_INT((x - g->x_min[0])/g->d_x);
|
|
}
|
|
|
|
int vroniObject::FindCellY(grid* g, double_arg y) {
|
|
return REAL_TO_INT((y - g->y_min[0])/g->d_y);
|
|
}
|
|
|
|
|
|
vr_bool vroniObject::CheckSegIntersectsGridVerSeg(int seg, double_arg x,
|
|
double_arg y1, double_arg y2,
|
|
double_arg eps)
|
|
{
|
|
double a, b, c;
|
|
coord p, q;
|
|
double bbx_min, bbx_max, bby_min, bby_max;
|
|
|
|
assert (y1 < y2);
|
|
|
|
GetSegEqnData(seg, &a, &b, &c);
|
|
|
|
p = GetSegStartCoord(seg);
|
|
q = GetSegEndCoord(seg);
|
|
|
|
bbx_min = std::min(p.x, q.x);
|
|
bbx_max = std::max(p.x, q.x);
|
|
bby_min = std::min(p.y, q.y);
|
|
bby_max = std::max(p.y, q.y);
|
|
|
|
//check whether the bounding box of the segment contains at least part of the grid line segment
|
|
if (lt (x - bbx_min, eps) || lt(bbx_max - x, eps)) return false;
|
|
if (lt (y2 - bby_min, eps) || lt(bby_max - y1, eps)) return false;
|
|
|
|
//check whether the segment crosses the grid line at the x-coordinate x at an y-coordinate between y1 and y2
|
|
if (eq(b, eps)) {
|
|
return eq(a*x + c, eps);
|
|
} else if (gt(b, eps)) {
|
|
return le(a*x + b*y1 + c, eps) && ge(a*x + b*y2 + c, eps);
|
|
} else {
|
|
return ge(a*x + b*y1 + c, eps) && le(a*x + b*y2 + c, eps);
|
|
}
|
|
}
|
|
|
|
vr_bool vroniObject::CheckSegIntersectsGridHorSeg(int seg,
|
|
double_arg x1, double_arg x2,
|
|
double_arg y, double_arg eps)
|
|
{
|
|
double a, b, c;
|
|
coord p, q;
|
|
double bbx_min, bbx_max, bby_min, bby_max;
|
|
|
|
assert (x1 < x2);
|
|
|
|
GetSegEqnData(seg, &a, &b, &c);
|
|
|
|
p = GetSegStartCoord(seg);
|
|
q = GetSegEndCoord(seg);
|
|
|
|
bbx_min = std::min(p.x, q.x);
|
|
bbx_max = std::max(p.x, q.x);
|
|
bby_min = std::min(p.y, q.y);
|
|
bby_max = std::max(p.y, q.y);
|
|
|
|
//check whether the bounding box of the segment contains at least part of the grid line segment
|
|
if (lt (x2 - bbx_min, eps) || lt(bbx_max - x1, eps)) return false;
|
|
if (lt (y - bby_min, eps) || lt(bby_max - y, eps)) return false;
|
|
|
|
//check whether the segment crosses the grid line at the y-coordinate y at an x-coordinate between x1 and x2
|
|
if (eq(a, eps)) {
|
|
return eq(b*y + c, eps);
|
|
} else if (gt(a, eps)) {
|
|
return le((a*x1 + b*y + c), eps) && ge(a*x2 + b*y + c, eps);
|
|
} else {
|
|
return ge((a*x1 + b*y + c), eps) && le(a*x2 + b*y + c, eps);
|
|
}
|
|
}
|
|
|
|
#ifdef GENUINE_ARCS
|
|
vr_bool vroniObject::CheckArcIntersectsGridVerSeg(int arc, double_arg x,
|
|
double_arg y1, double_arg y2,
|
|
double_arg eps)
|
|
{
|
|
double r;
|
|
coord c;
|
|
double xa, y1a, y2a, yy;
|
|
double yymin, yymax;
|
|
|
|
assert (y1 < y2);
|
|
|
|
c = GetArcCenter(arc);
|
|
|
|
if (y1 < c.y && y2 > c.y) { //simplify handling of this special case
|
|
return CheckArcIntersectsGridVerSeg(arc, x, y1, c.y, eps)
|
|
|| CheckArcIntersectsGridVerSeg(arc, x, c.y, y2, eps);
|
|
}
|
|
|
|
xa = x - c.x;
|
|
y1a = y1 - c.y;
|
|
y2a = y2 - c.y;
|
|
r = GetArcRadius(arc);
|
|
|
|
yy = r*r - xa*xa;
|
|
|
|
yymin = std::min(y1a*y1a, y2a*y2a);
|
|
yymax = std::max(y1a*y1a, y2a*y2a);
|
|
|
|
if (le(yymin - yy, eps) && ge(yymax - yy, eps)) {
|
|
//grid segment intersects (or is very near to) the full circle - now check only the arc
|
|
|
|
coord p, q;
|
|
double arc_xmin, arc_xmax;
|
|
|
|
p = GetArcStartCoord(arc);
|
|
q = GetArcEndCoord(arc);
|
|
|
|
//Note: arcs are always stored CCW
|
|
|
|
if (y1 >= c.y) {
|
|
//we only consider the upper half because y1, y2 >= c.y
|
|
assert (y2 >= c.y);
|
|
|
|
//Note: we don't allow arcs that span an angle of 180 or more degrees
|
|
if (p.y < c.y) {
|
|
//arc starts in the lower half of the circle
|
|
if (q.y < c.y)
|
|
return false; //arc is completely in the lower half
|
|
arc_xmax = c.x + r;
|
|
} else {
|
|
arc_xmax = p.x;
|
|
}
|
|
|
|
if (q.y < c.y) {
|
|
arc_xmin = c.x - r;
|
|
} else {
|
|
arc_xmin = q.x;
|
|
}
|
|
} else {
|
|
//we only consider the lower half because y1, y2 <= c.y
|
|
assert (y2 <= c.y);
|
|
|
|
if (p.y > c.y) {
|
|
if (q.y > c.y)
|
|
return false; //arc is completely in the upper half
|
|
arc_xmin = c.x - r;
|
|
} else {
|
|
//Note: we don't allow arcs that span an angle of 180 or more degrees
|
|
arc_xmin = p.x;
|
|
}
|
|
|
|
if (q.y > c.y) {
|
|
arc_xmax = c.x + r;
|
|
} else {
|
|
arc_xmax = q.x;
|
|
}
|
|
}
|
|
|
|
return le(arc_xmin - x, eps) && ge(arc_xmax - x, eps);
|
|
} else {
|
|
return false;
|
|
}
|
|
}
|
|
|
|
vr_bool vroniObject::CheckArcIntersectsGridHorSeg(int arc,
|
|
double_arg x1, double_arg x2,
|
|
double_arg y, double_arg eps)
|
|
{
|
|
double r;
|
|
coord c;
|
|
double ya, x1a, x2a, xx;
|
|
double xxmin, xxmax;
|
|
|
|
assert (x1 < x2);
|
|
|
|
c = GetArcCenter(arc);
|
|
|
|
if (x1 < c.x && x2 > c.x) { //simplify handling of this special case
|
|
return CheckArcIntersectsGridHorSeg(arc, x1, c.x, y, eps)
|
|
|| CheckArcIntersectsGridHorSeg(arc, c.x, x2, y, eps);
|
|
}
|
|
|
|
ya = y - c.y;
|
|
x1a = x1 - c.x;
|
|
x2a = x2 - c.x;
|
|
r = GetArcRadius(arc);
|
|
|
|
xx = r*r - ya*ya;
|
|
|
|
xxmin = std::min(x1a*x1a, x2a*x2a);
|
|
xxmax = std::max(x1a*x1a, x2a*x2a);
|
|
|
|
if (le(xxmin - xx, eps) && ge(xxmax - xx, eps)) {
|
|
//grid segment intersects (or is very near to) the full circle - now check only the arc
|
|
|
|
coord p, q;
|
|
double arc_ymin, arc_ymax;
|
|
|
|
p = GetArcStartCoord(arc);
|
|
q = GetArcEndCoord(arc);
|
|
|
|
//Note: arcs are always stored CCW
|
|
|
|
if (x1 >= c.x) {
|
|
//we only consider the right half because x1, x2 >= c.x
|
|
assert (x2 >= c.x);
|
|
|
|
//Note: we don't allow arcs that span an angle of 180 or more degrees
|
|
if (p.x < c.x) {
|
|
//arc starts in the left half of the circle
|
|
if (q.x < c.x)
|
|
return false; //arc is completely in the left half
|
|
arc_ymin = c.y - r;
|
|
} else {
|
|
arc_ymin = p.y;
|
|
}
|
|
|
|
if (q.x < c.x) {
|
|
arc_ymax = c.y + r;
|
|
} else {
|
|
arc_ymax = q.y;
|
|
}
|
|
} else {
|
|
//we only consider the left half because x1, x2 <= c.x
|
|
assert (x2 <= c.x);
|
|
|
|
if (p.x > c.x) {
|
|
//arc starts in the right half of the circle
|
|
if (q.x > c.x)
|
|
return false; //arc is completely in the right half
|
|
arc_ymax = c.y + r;
|
|
} else {
|
|
//Note: we don't allow arcs that span an angle of 180 or more degrees
|
|
arc_ymax = p.y;
|
|
}
|
|
|
|
if (q.x > c.x) {
|
|
arc_ymin = c.y - r;
|
|
} else {
|
|
arc_ymin = q.y;
|
|
}
|
|
}
|
|
|
|
return le(arc_ymin - y, eps) && ge(arc_ymax - y, eps);
|
|
} else {
|
|
return false;
|
|
}
|
|
}
|
|
#endif
|
|
|
|
|
|
void vroniObject::AddPointToGrid(grid* g, int i1, double_arg eps)
|
|
{
|
|
int i, j;
|
|
point_entry* entry;
|
|
coord p;
|
|
|
|
assert(InPntsList(i1));
|
|
|
|
p = GetPntCoords(i1);
|
|
entry = new point_entry;
|
|
entry->num = i1;
|
|
|
|
for (i = FindCellX(g, p.x - eps); i <= FindCellX(g, p.x + eps); i++) {
|
|
for (j = FindCellY(g, p.y - eps); j <= FindCellY(g, p.y + eps); j++) {
|
|
assert (!(i < 0 || j < 0 || i >= g->n_x || j >= g->n_y));
|
|
|
|
entry->cells.insert(&(g->cells[i][j]));
|
|
g->cells[i][j].points.insert(entry);
|
|
}
|
|
}
|
|
|
|
g->points[i1] = entry;
|
|
}
|
|
|
|
void vroniObject::AddSegToGrid(grid* g, int i1, double_arg eps)
|
|
{
|
|
int i, j;
|
|
seg_entry* entry;
|
|
coord p, q;
|
|
std::queue<std::pair<int,int> > queue;
|
|
double x, y;
|
|
|
|
assert(InSegsList(i1));
|
|
|
|
entry = new seg_entry(comp_eps);
|
|
entry->num = i1;
|
|
|
|
p = GetSegStartCoord(i1);
|
|
q = GetSegEndCoord(i1);
|
|
|
|
/*
|
|
* we basically do a BFS from the seg's start vertex;
|
|
* two cells are connected if and only if the seg intersects the grid line segment between them
|
|
*/
|
|
|
|
//First find the start cell(s) -> can be up to four cells (if the start vertex is at a corner of a grid cell)
|
|
for (i = FindCellX(g, p.x - eps); i <= FindCellX(g, p.x + eps); i++) {
|
|
for (j = FindCellY(g, p.y - eps); j <= FindCellY(g, p.y + eps); j++) {
|
|
assert (!(i < 0 || j < 0 || i >= g->n_x || j >= g->n_y));
|
|
queue.push(std::make_pair(i, j));
|
|
}
|
|
}
|
|
|
|
while (!queue.empty()) {
|
|
std::pair<int,int> cell;
|
|
|
|
cell = queue.front(); queue.pop();
|
|
i = cell.first;
|
|
j = cell.second;
|
|
|
|
if (!(i < 0 || j < 0 || i >= g->n_x || j >= g->n_y)) {
|
|
|
|
if (!g->cells[i][j].segs.insert(entry).second)
|
|
continue; //we have already visited this cell
|
|
|
|
x = (q.x >= p.x ? g->x_max[i] : g->x_min[i]);
|
|
y = (q.y >= p.y ? g->y_max[j] : g->y_min[j]);
|
|
|
|
if (CheckSegIntersectsGridVerSeg(i1, x, g->y_min[j], g->y_max[j], eps)) {
|
|
queue.push(std::make_pair(q.x >= p.x ? i+1 : i-1, j));
|
|
}
|
|
if (CheckSegIntersectsGridHorSeg(i1, g->x_min[i], g->x_max[i], y, eps)) {
|
|
queue.push(std::make_pair(i, q.y >= p.y ? j+1 : j-1));
|
|
}
|
|
}
|
|
}
|
|
|
|
g->segs[i1] = entry;
|
|
}
|
|
|
|
|
|
#ifdef GENUINE_ARCS
|
|
void vroniObject::AddArcToGrid(grid* g, int i1, double_arg eps)
|
|
{
|
|
int i, j;
|
|
arc_entry* entry;
|
|
coord p;
|
|
std::queue<std::pair<int,int> > queue;
|
|
|
|
assert(InArcsList(i1));
|
|
|
|
entry = new arc_entry;
|
|
entry->num = i1;
|
|
|
|
p = GetArcStartCoord(i1);
|
|
|
|
/*
|
|
* we basically do a BFS from the arc's start vertex;
|
|
* two cells are connected if and only if the arc intersects the grid line segment between them
|
|
*/
|
|
|
|
//First find the start cell(s) -> can be up to four cells (if the start vertex is at a corner of a grid cell)
|
|
for (i = FindCellX(g, p.x - eps); i <= FindCellX(g, p.x + eps); i++) {
|
|
for (j = FindCellY(g, p.y - eps); j <= FindCellY(g, p.y + eps); j++) {
|
|
assert (!(i < 0 || j < 0 || i >= g->n_x || j >= g->n_y));
|
|
|
|
queue.push(std::make_pair(i, j));
|
|
}
|
|
}
|
|
|
|
while (!queue.empty()) {
|
|
std::pair<int,int> cell;
|
|
|
|
cell = queue.front(); queue.pop();
|
|
i = cell.first;
|
|
j = cell.second;
|
|
|
|
if (!(i < 0 || j < 0 || i >= g->n_x || j >= g->n_y)) {
|
|
|
|
if (!g->cells[i][j].arcs.insert(entry).second)
|
|
continue; //we have already visited this cell
|
|
|
|
if (CheckArcIntersectsGridVerSeg(i1, g->x_min[i], g->y_min[j], g->y_max[j], eps))
|
|
queue.push(std::make_pair(i-1, j));
|
|
if (CheckArcIntersectsGridHorSeg(i1, g->x_min[i], g->x_max[i], g->y_min[j], eps))
|
|
queue.push(std::make_pair(i, j-1));
|
|
if (CheckArcIntersectsGridVerSeg(i1, g->x_max[i], g->y_min[j], g->y_max[j], eps))
|
|
queue.push(std::make_pair(i+1, j));
|
|
if (CheckArcIntersectsGridHorSeg(i1, g->x_min[i], g->x_max[i], g->y_max[j], eps))
|
|
queue.push(std::make_pair(i, j+1));
|
|
}
|
|
}
|
|
|
|
g->arcs[i1] = entry;
|
|
}
|
|
#endif
|
|
|
|
|
|
void vroniObject::FillGrid(grid* g, double_arg eps)
|
|
{
|
|
int i;
|
|
|
|
for (i = 2; i <= num_pnts - 3; i++) {
|
|
if (!IsDeletedPnt(i)) {
|
|
AddPointToGrid(g, i, eps);
|
|
}
|
|
}
|
|
|
|
for (i = 0; i <= num_segs - 1; i++)
|
|
AddSegToGrid(g, i, eps);
|
|
|
|
#ifdef GENUINE_ARCS
|
|
for (i = 0; i <= num_arcs - 1; i++)
|
|
AddArcToGrid(g, i, eps);
|
|
#endif
|
|
}
|
|
|
|
|
|
void vroniObject::DeleteGrid(grid* grid)
|
|
{
|
|
std::map<int,point_entry*>::iterator it1;
|
|
std::map<int,seg_entry*>::iterator it2;
|
|
#ifdef GENUINE_ARCS
|
|
std::map<int,arc_entry*>::iterator it3;
|
|
#endif
|
|
int i;
|
|
|
|
for (it1 = grid->points.begin(); it1 != grid->points.end(); it1++)
|
|
delete (*it1).second;
|
|
|
|
for (it2 = grid->segs.begin(); it2 != grid->segs.end(); it2++)
|
|
delete (*it2).second;
|
|
|
|
#ifdef GENUINE_ARCS
|
|
for (it3 = grid->arcs.begin(); it3 != grid->arcs.end(); it3++)
|
|
delete (*it3).second;
|
|
#endif
|
|
|
|
delete [] grid->x_min;
|
|
delete [] grid->x_max;
|
|
delete [] grid->y_min;
|
|
delete [] grid->y_max;
|
|
|
|
for (i = 0; i < grid->n_x; i++)
|
|
delete [] grid->cells[i];
|
|
|
|
delete [] grid->cells;
|
|
|
|
delete grid;
|
|
}
|
|
|
|
|
|
#ifdef DBG_GRAPHICS
|
|
|
|
static void SketchCell(grid* grid, int i, int j)
|
|
{
|
|
coord p, q;
|
|
grid_cell* c = &grid->cells[i][j];
|
|
|
|
int color = !c->segs.empty()
|
|
#ifdef GENUINE_ARCS
|
|
|| !c->arcs.empty()
|
|
#endif
|
|
|| !c->points.empty() ? GridColor : NoColor;
|
|
|
|
p.x = c->x_min;
|
|
q.x = c->x_max;
|
|
p.y = c->y_min;
|
|
q.y = c->y_min;
|
|
AddEdgeToBuffer(p.x, p.y, q.x, q.y, color);
|
|
p.y = c->y_max;
|
|
q.y = c->y_max;
|
|
AddEdgeToBuffer(p.x, p.y, q.x, q.y, color);
|
|
p.x = c->x_min;
|
|
q.x = c->x_min;
|
|
p.y = c->y_min;
|
|
q.y = c->y_max;
|
|
AddEdgeToBuffer(p.x, p.y, q.x, q.y, color);
|
|
p.x = c->x_max;
|
|
q.x = c->x_max;
|
|
AddEdgeToBuffer(p.x, p.y, q.x, q.y, color);
|
|
}
|
|
|
|
|
|
//only for debugging purposes
|
|
static void SketchGrid(grid* grid)
|
|
{
|
|
int i, j;
|
|
for (i = 0; i < grid->n_x; ++i)
|
|
for (j = 0; j < grid->n_y; ++j)
|
|
SketchCell(grid, i, j);
|
|
}
|
|
#endif
|
|
|
|
|
|
void vroniObject::SplitAtIntersections(grid* grid, double_arg eps)
|
|
{
|
|
int i1, i2, i3;
|
|
std::set<coord,comp_seg_iscoord>::iterator it;
|
|
std::map<int,seg_entry*>::iterator it2;
|
|
vr_bool reverse;
|
|
#ifdef GENUINE_ARCS
|
|
std::set<arc_intersection,comp_arc_iscoord>::iterator it1;
|
|
std::map<int,arc_entry*>::iterator it3;
|
|
#endif
|
|
|
|
comp_eps = eps;
|
|
|
|
for (it2 = grid->segs.begin(); it2 != grid->segs.end(); it2++) {
|
|
seg_entry* sege = (*it2).second;
|
|
i1 = sege->num;
|
|
coord s, e;
|
|
|
|
s = GetSegStartCoord(i1);
|
|
e = GetSegEndCoord(i1);
|
|
|
|
reverse = lt(e.x - s.x, comp_eps) || (eq(e.x - s.x, comp_eps) && e.y < s.y);
|
|
|
|
for (it = sege->intersections->begin(); it != sege->intersections->end(); it++) {
|
|
coord c = *it;
|
|
#ifdef EXT_APPL_PNTS
|
|
// MODIF: al nuovo punto aggiunto dall'intersezione posso assegnare soltanto il loop originario, visto che non deriva
|
|
// da alcun punto della curva di partenza
|
|
i3 = StorePnt(c.x, c.y, {segs[i1].ext_appl.first, -1});
|
|
#else
|
|
i3 = StorePnt(c.x, c.y);
|
|
#endif
|
|
|
|
HandleSplitSeg(i1, i3, i2);
|
|
if (!reverse)
|
|
i1 = i2;
|
|
}
|
|
}
|
|
|
|
#ifdef GENUINE_ARCS
|
|
for (it3 = grid->arcs.begin(); it3 != grid->arcs.end(); it3++) {
|
|
arc_entry* arce = (*it3).second;
|
|
i1 = arce->num;
|
|
|
|
for (it1 = arce->intersections.begin(); it1 != arce->intersections.end(); it1++) {
|
|
coord c = (*it1).c;
|
|
#ifdef EXT_APPL_PNTS
|
|
// MODIF: al nuovo punto aggiunto dall'intersezione posso assegnare soltanto il loop originario, visto che non deriva
|
|
// da alcun punto della curva di partenza
|
|
i3 = StorePnt(c.x, c.y, {arcs[i1].ext_appl.first, -1});
|
|
#else
|
|
i3 = StorePnt(c.x, c.y);
|
|
#endif
|
|
|
|
HandleSplitArc(i1, i3, i1);
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
|
|
|
|
vr_bool vroniObject::CheckPntPnt(int i1, int i2, double_arg eps)
|
|
{
|
|
int i3;
|
|
double dist;
|
|
coord p1, p2;
|
|
|
|
assert(InPntsList(i1));
|
|
assert(InPntsList(i2));
|
|
|
|
if (i1 == i2) return false;
|
|
|
|
if (i2 < i1) VroniSwap(i1, i2, i3);
|
|
|
|
p1 = GetPntCoords(i1);
|
|
p2 = GetPntCoords(i2);
|
|
dist = PntPntDist(p1, p2);
|
|
|
|
if (eq(dist, eps)) {
|
|
IntWarning("pnt", i1, "pnt", i2);
|
|
return true;
|
|
}
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|
|
vr_bool vroniObject::CheckPntOnSeg(int seg, int pnt, double_arg eps)
|
|
{
|
|
double a, b, c, lgth, dist_q;
|
|
coord p, q;
|
|
|
|
if (!IsSegStartPnt(seg, pnt) && !IsSegEndPnt(seg, pnt)) {
|
|
GetSegEqnData(seg, &a, &b, &c);
|
|
q = GetPntCoords(pnt);
|
|
p = GetSegStartCoord(seg);
|
|
lgth = GetSegLgth(seg);
|
|
dist_q = PntSegDist(p, q, a, b, c, lgth, TINY);
|
|
|
|
if (eq(dist_q, eps)) {
|
|
IntWarning("seg", seg, "pnt", pnt);
|
|
return true;
|
|
}
|
|
}
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|
|
#ifdef GENUINE_ARCS
|
|
vr_bool vroniObject::CheckPntOnArc(int arc, int pnt, double_arg eps)
|
|
{
|
|
double r, dist;
|
|
coord p3, p, cp, ns, ne;
|
|
|
|
if (!IsArcStartPnt(arc, pnt) && !IsArcEndPnt(arc, pnt)) {
|
|
p3 = GetArcCenter(arc);
|
|
r = GetArcRadius(arc);
|
|
p = GetPntCoords(pnt);
|
|
ns = GetArcStartNormal(arc);
|
|
ne = GetArcEndNormal(arc);
|
|
cp = VecSub(p, p3);
|
|
dist = PntArcDist(ns, ne, cp, r, TINY);
|
|
|
|
if (eq(dist, eps)) {
|
|
IntWarning("arc", arc, "pnt", pnt);
|
|
return true;
|
|
}
|
|
}
|
|
|
|
return false;
|
|
}
|
|
#endif
|
|
|
|
|
|
|
|
|
|
/* */
|
|
/* check (and compute) the intersection of seg i1 and seg i2, and split */
|
|
/* the sites accordingly. */
|
|
/* */
|
|
/* note: (1) emphasis is placed on reliability rather than speed! */
|
|
/* (2) an iterative solver is used if the analytical computation of */
|
|
/* the points of intersection does not yield accurate results. */
|
|
/* (3) it is assumed that a pnt-on-seg test has been carried out prior */
|
|
/* to calling this routine. */
|
|
/* */
|
|
vr_bool vroniObject::GetSegSegIntersection(int i1, int i2, coord* result,
|
|
double_arg eps, vr_bool *ident_segs)
|
|
{
|
|
int s_p, s_q, s_u, s_v, i3, i5;
|
|
double a1, b1, c1, lgth, dist_p, dist_q, dist_u, dist_v;
|
|
double a2, b2, c2, t;
|
|
coord dir, p, q, u, v, w;
|
|
int ident = 0;
|
|
#ifndef NDEBUG
|
|
vr_bool success;
|
|
#endif
|
|
|
|
if (i1 == i2) return false;
|
|
*ident_segs = false;
|
|
|
|
if (i2 < i1) VroniSwap(i1, i2, i3);
|
|
|
|
/* */
|
|
/* compute the sidedness of the end points u,v of seg 2 w.r.t. the */
|
|
/* supporting line of seg 1. note that the segs cannot overlap as we */
|
|
/* assume that a pnt-on-seg test has been carried out prior to calling */
|
|
/* this function. thus, we may stop the test once one vertex lies on the */
|
|
/* supporting line of the other seg. */
|
|
/* */
|
|
GetSegEqnData(i1, &a1, &b1, &c1);
|
|
p = GetSegStartCoord(i1);
|
|
q = GetSegEndCoord(i1);
|
|
lgth = GetSegLgth(i1);
|
|
u = GetSegStartCoord(i2);
|
|
v = GetSegEndCoord(i2);
|
|
|
|
i5 = Get1stVtx(i2, SEG);
|
|
if (IsSegStartPnt(i1, i5) || IsSegEndPnt(i1, i5)) {
|
|
s_u = 0;
|
|
ident = 1;
|
|
}
|
|
else {
|
|
dist_u = PntLineDist(a1, b1, c1, u);
|
|
s_u = SignEps(dist_u, eps);
|
|
}
|
|
|
|
i5 = Get2ndVtx(i2, SEG);
|
|
if (IsSegStartPnt(i1, i5) || IsSegEndPnt(i1, i5)) {
|
|
s_v = 0;
|
|
++ident;
|
|
}
|
|
else {
|
|
dist_v = PntLineDist(a1, b1, c1, v);
|
|
s_v = SignEps(dist_v, eps);
|
|
}
|
|
|
|
|
|
if (ident == 2) {
|
|
//we don't allow duplicate segs
|
|
*ident_segs = true;
|
|
return false;
|
|
}
|
|
|
|
if (s_u * s_v == -1) {
|
|
/* */
|
|
/* compute the sidedness of the end points p,q of seg 1 w.r.t. the */
|
|
/* supporting line of seg 2. */
|
|
/* */
|
|
GetSegEqnData(i2, &a2, &b2, &c2);
|
|
dist_p = PntLineDist(a2, b2, c2, p);
|
|
s_p = SignEps(dist_p, eps);
|
|
dist_q = PntLineDist(a2, b2, c2, q);
|
|
s_q = SignEps(dist_q, eps);
|
|
|
|
if (s_p * s_q == -1) {
|
|
/* */
|
|
/* this is a genuine intersection. let's compute it as accurately */
|
|
/* as possible. */
|
|
/* */
|
|
dir = GetSegDirection(i1);
|
|
t = lgth * (Abs(dist_p) / (Abs(dist_p) + Abs(dist_q)));
|
|
w = RayPnt(p, dir, t);
|
|
|
|
#ifndef NDEBUG
|
|
success = IterativeIntersection(SEG, SEG, a2, b2, c2, p, 0.0, p,
|
|
0.0, lgth, dir, 0.0, &t, &w);
|
|
assert(success);
|
|
#else
|
|
(void) IterativeIntersection(SEG, SEG, a2, b2, c2, p, 0.0, p,
|
|
0.0, lgth, dir, 0.0, &t, &w);
|
|
#endif
|
|
|
|
IntWarning("seg", i1, "seg", i2);
|
|
|
|
*result = w;
|
|
return true;
|
|
}
|
|
}
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
#ifdef GENUINE_ARCS
|
|
/* */
|
|
/* check (and compute) the intersection of seg i1 and arc i2, and split */
|
|
/* the sites accordingly. */
|
|
/* */
|
|
/* note: (1) emphasis is placed on reliability rather than speed! */
|
|
/* (2) an iterative solver is used if the analytical computation of */
|
|
/* the points of intersection does not yield accurate results. */
|
|
/* (3) it is assumed that a pnt-on-seg and pnt-on-arc test has been */
|
|
/* carried out. */
|
|
/* */
|
|
int vroniObject::GetArcArcIntersection(int i1, int i2,
|
|
coord* res1, coord* res2,
|
|
double_arg eps, vr_bool *ident_arcs)
|
|
{
|
|
double R1, R2, d, alpha, alpha1, alpha2, t1, t2, f;
|
|
coord q1, q2, q3, q4, C1, C2, v, Ns1, Ne1, Ns2, Ne2, shared_point;
|
|
vr_bool success, shared = false;
|
|
int i3;
|
|
int tmp, intersections = 0;
|
|
double lo, hi, a1, a2, d1, d2;
|
|
int i, MAX_ITERATIONS = 300;
|
|
|
|
|
|
if (i1 == i2) return 0;
|
|
*ident_arcs = false;
|
|
|
|
|
|
R1 = GetArcRadius(i1);
|
|
assert(gt(R1, eps));
|
|
R2 = GetArcRadius(i2);
|
|
assert(gt(R2, eps));
|
|
|
|
if ((R1 < R2) || ((R1 == R2) && (i2 < i1))) {
|
|
VroniSwap(R1, R2, d);
|
|
VroniSwap(i1, i2, tmp);
|
|
}
|
|
|
|
/* */
|
|
/* first check whether the arcs share one or two endpoints */
|
|
/* */
|
|
i3 = Get1stVtxArc(i1);
|
|
shared_point = GetPntCoords(i3);
|
|
if (IsArcStartPnt(i2, i3) || IsArcEndPnt(i2, i3)) {
|
|
shared = true;
|
|
}
|
|
|
|
i3 = Get2ndVtxArc(i1);
|
|
if (IsArcStartPnt(i2, i3) || IsArcEndPnt(i2, i3)) {
|
|
if (shared) {
|
|
/* */
|
|
/* the arcs share both end points */
|
|
/* */
|
|
C1 = GetArcCenter(i1);
|
|
C2 = GetArcCenter(i2);
|
|
d = PntPntDist(C1, C2);
|
|
if (eq(d, eps)) {
|
|
/* */
|
|
/* the arcs are either identical or complementary. since we do */
|
|
/* not allow arcs that span an angle of 180 or more degrees, the */
|
|
/* arcs have to be identical. */
|
|
/* */
|
|
*ident_arcs = true;
|
|
return 0;
|
|
}
|
|
else {
|
|
/* */
|
|
/* there can't be a non-trivial intersection */
|
|
/* */
|
|
return 0;
|
|
}
|
|
}
|
|
shared = true;
|
|
shared_point = GetPntCoords(i3);
|
|
}
|
|
|
|
C1 = GetArcCenter(i1);
|
|
C2 = GetArcCenter(i2);
|
|
v = VecSub(C2, C1);
|
|
d = VecLen(v);
|
|
|
|
Ns1 = GetArcStartNormal(i1);
|
|
Ne1 = GetArcEndNormal(i1);
|
|
Ns2 = GetArcStartNormal(i2);
|
|
Ne2 = GetArcEndNormal(i2);
|
|
|
|
if (eq(d, eps)) {
|
|
/* */
|
|
/* the two arcs have the same center. a partial overlap of the arcs */
|
|
/* cannot occur as we assume that a pnt-on-arc test has been carried */
|
|
/* out */
|
|
/* */
|
|
return 0;
|
|
}
|
|
|
|
if (lt(d - R1 + R2, eps) || gt(d - R1 - R2 , eps)) {
|
|
/* */
|
|
/* circle 2 is completely within circle 1 or */
|
|
/* the circles don't overlap */
|
|
/* */
|
|
return 0;
|
|
}
|
|
|
|
alpha = atan2(v.y, v.x);
|
|
alpha1 = alpha - M_PI;
|
|
alpha2 = alpha + M_PI;
|
|
t1 = alpha - M_PI_2;
|
|
t2 = alpha + M_PI_2;
|
|
|
|
if (eq(d - R1 - R2, eps)) {
|
|
/* */
|
|
/* the circles touch at one point (from the outside) */
|
|
/* */
|
|
if ((d - R1 - R2) >= 0.0) {
|
|
/* */
|
|
/* in this case we don't really need an iterative search, as we */
|
|
/* know the circles intersect exactly once and the analytic result */
|
|
/* is accurate enough */
|
|
/* */
|
|
f = (R1 - R2 + d)/2/d;
|
|
q1 = VecAdd(C1, VecMult(f, v));
|
|
|
|
q3 = VecSub(q1, C1);
|
|
q4 = VecSub(q1, C2);
|
|
|
|
if ((IsInArcConeStrict(Ns1, Ne1, q3) == 0) &&
|
|
(IsInArcConeStrict(Ns2, Ne2, q4) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
return 1;
|
|
}
|
|
else {
|
|
return 0;
|
|
}
|
|
}
|
|
|
|
/* */
|
|
/* calculate the point on circle 2 that should be nearest to circle 1 */
|
|
/* */
|
|
q1 = CirclePnt(C2, R2, alpha1);
|
|
d1 = PntCircleDist(C1, R1, q1);
|
|
|
|
if (d1 >= 0) {
|
|
/* */
|
|
/* The calculated nearest point is not inside circle 1 (this might */
|
|
/* be due to numerical inaccuracies or because no such point */
|
|
/* exists): thus we can't do a iterative search using bisection... */
|
|
/* Use a ternary search to improve the calculation of the farthest */
|
|
/* point */
|
|
/* */
|
|
lo = alpha1 - M_PI_2;
|
|
hi = alpha1 + M_PI_2;
|
|
|
|
for (i = 0; i < MAX_ITERATIONS; i++) {
|
|
a1 = lo + (hi - lo)/3;
|
|
a2 = lo + (hi - lo)*2/3;
|
|
if (a1 == lo || a2 == hi)
|
|
/* we have reached machine precision, so we can stop */
|
|
break;
|
|
|
|
q1 = CirclePnt(C2, R2, a1);
|
|
q2 = CirclePnt(C2, R2, a2);
|
|
|
|
d1 = PntCircleDist(C1, R1, q1);
|
|
d2 = PntCircleDist(C1, R1, q2);
|
|
|
|
if (d1 < d2)
|
|
hi = a2;
|
|
else
|
|
lo = a1;
|
|
}
|
|
|
|
if (d1 >= 0) {
|
|
/* */
|
|
/* nearest point on circle 2 is on or outside circle 1 (but very */
|
|
/* near), so we return it as an intersection */
|
|
/* */
|
|
assert (eq(d1, 2*eps));
|
|
|
|
q3 = VecSub(q1, C1);
|
|
q4 = VecSub(q1, C2);
|
|
|
|
if ((IsInArcConeStrict(Ns1, Ne1, q3) == 0) &&
|
|
(IsInArcConeStrict(Ns2, Ne2, q4) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
return 1;
|
|
}
|
|
else {
|
|
return 0;
|
|
}
|
|
}
|
|
else {
|
|
/* */
|
|
/* we found a point on circle 2 that lies inside circle 1, so we */
|
|
/* can do the usual iterative search (bisection) */
|
|
/* */
|
|
alpha = a1;
|
|
alpha1 = alpha - M_PI;
|
|
alpha2 = alpha + M_PI;
|
|
t1 = a1;
|
|
t2 = a1;
|
|
}
|
|
}
|
|
}
|
|
else if (eq(R2 - R1 + d, eps)) {
|
|
/* */
|
|
/* circle 2's outside touches circle 1 from the inside */
|
|
/* */
|
|
/* calculate the point on circle 2 that should be farthest from the */
|
|
/* center of circle 1 */
|
|
/* */
|
|
q1 = CirclePnt(C2, R2, alpha);
|
|
d1 = PntCircleDist(C1, R1, q1);
|
|
|
|
if (d1 <= 0) {
|
|
/* */
|
|
/* the calculated farthest point is not outside circle 1 (this */
|
|
/* might be due to numerical inaccuracies or because no such point */
|
|
/* exists): thus, we can't do a iterative search using bisection... */
|
|
/* Use a ternary search to improve the calculation of the farthest */
|
|
/* point */
|
|
/* */
|
|
lo = alpha - M_PI_2;
|
|
hi = alpha + M_PI_2;
|
|
|
|
for (i = 0; i < MAX_ITERATIONS; i++) {
|
|
a1 = lo + (hi - lo)/3;
|
|
a2 = lo + (hi - lo)*2/3;
|
|
if (a1 == lo || a2 == hi)
|
|
/* we have reached machine precision, so we can stop */
|
|
break;
|
|
|
|
q1 = CirclePnt(C2, R2, a1);
|
|
q2 = CirclePnt(C2, R2, a2);
|
|
|
|
d1 = PntCircleDist(C1, R1, q1);
|
|
d2 = PntCircleDist(C1, R1, q2);
|
|
|
|
if (d1 > d2)
|
|
hi = a2;
|
|
else
|
|
lo = a1;
|
|
}
|
|
|
|
if (d1 <= 0) {
|
|
/* */
|
|
/* farthest point on circle 2 is on or inside circle 1 (but */
|
|
/* very near), so we return it as an intersection */
|
|
/* */
|
|
assert (eq(d1, 2*eps));
|
|
|
|
q3 = VecSub(q1, C1);
|
|
q4 = VecSub(q1, C2);
|
|
|
|
if ((IsInArcConeStrict(Ns1, Ne1, q3) == 0) &&
|
|
(IsInArcConeStrict(Ns2, Ne2, q4) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
return 1;
|
|
}
|
|
else {
|
|
return 0;
|
|
}
|
|
}
|
|
else {
|
|
/* */
|
|
/* we found a point on circle 2 that lies outside circle 1, so */
|
|
/* we can do the usual iterative search (bisection) */
|
|
/* */
|
|
alpha = a1;
|
|
alpha1 = alpha - M_PI;
|
|
alpha2 = alpha + M_PI;
|
|
t1 = a1;
|
|
t2 = a1;
|
|
}
|
|
}
|
|
}
|
|
|
|
success = IterativeIntersection(ARC, ARC, 0.0, 0.0, 0.0, C1, R1, C2,
|
|
alpha1, alpha, v, R2, &t1, &q1);
|
|
|
|
q3 = VecSub(q1, C1);
|
|
q4 = VecSub(q1, C2);
|
|
if ((success && IsInArcConeStrict(Ns1, Ne1, q3) == 0) &&
|
|
(IsInArcConeStrict(Ns2, Ne2, q4) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
intersections++;
|
|
}
|
|
|
|
success = IterativeIntersection(ARC, ARC, 0.0, 0.0, 0.0, C1, R1, C2,
|
|
alpha, alpha2, v, R2, &t2, &q2);
|
|
|
|
q3 = VecSub(q2, C1);
|
|
q4 = VecSub(q2, C2);
|
|
if (success &&
|
|
(IsInArcConeStrict(Ns1, Ne1, q3) == 0) &&
|
|
(IsInArcConeStrict(Ns2, Ne2, q4) == 0) &&
|
|
(!shared || ne(PntPntDist(q2, shared_point), TINY))) {
|
|
if (intersections == 1) {
|
|
if (ne(PntPntDist(q1, q2), eps)) {
|
|
*res2 = q2;
|
|
intersections++;
|
|
}
|
|
}
|
|
else {
|
|
*res1 = q2;
|
|
intersections++;
|
|
}
|
|
}
|
|
|
|
if (intersections > 0)
|
|
IntWarning("arc", i1, "arc", i2);
|
|
|
|
return intersections;
|
|
}
|
|
|
|
|
|
/* */
|
|
/* check (and compute) the intersection of seg i1 and arc i2, and split */
|
|
/* the sites accordingly. */
|
|
/* */
|
|
/* note: (1) emphasis is placed on reliability rather than speed! */
|
|
/* (2) an iterative solver is used if the analytical computation of */
|
|
/* the points of intersection does not yield accurate results. */
|
|
/* (3) it is assumed that a pnt-on-seg and pnt-on-arc test has been */
|
|
/* carried out. */
|
|
/* */
|
|
int vroniObject::GetSegArcIntersection(int seg, int arc, coord* res1,
|
|
coord* res2, double_arg eps)
|
|
{
|
|
int j1, j2;
|
|
double lgth, d, R, u1, u2, t, dist, a, b, c;
|
|
coord C, S, E, v1, v2, Ns, Ne, q1, q2, q3;
|
|
vr_bool shared, success;
|
|
coord shared_point;
|
|
int intersections;
|
|
|
|
/* */
|
|
/* first check whether the seg and the arc share one or two end points */
|
|
/* */
|
|
intersections = 0;
|
|
shared = false;
|
|
j1 = Get1stVtx(seg, SEG);
|
|
shared_point = GetSegStartCoord(seg);
|
|
if (IsArcStartPnt(arc, j1) || IsArcEndPnt(arc, j1)) {
|
|
shared = true;
|
|
}
|
|
|
|
j2 = Get2ndVtx(seg, SEG);
|
|
if (IsArcStartPnt(arc, j2) || IsArcEndPnt(arc, j2)) {
|
|
if (shared) {
|
|
/* */
|
|
/* the arc and the seg share both end points. thus, only trivial */
|
|
/* intersections exist. */
|
|
/* */
|
|
return 0;
|
|
}
|
|
shared = true;
|
|
shared_point = GetSegEndCoord(seg);
|
|
}
|
|
|
|
/* */
|
|
/* compute the parameter d of the projection p1 + v1 * d of the arc's */
|
|
/* center on the infinite line. */
|
|
/* d is important to determine the interval for the iterative computation */
|
|
/* */
|
|
C = GetArcCenter(arc);
|
|
S = GetSegStartCoord(seg);
|
|
E = GetSegEndCoord(seg);
|
|
lgth = GetSegLgth(seg);
|
|
v1 = GetSegDirection(seg);
|
|
Ns = GetArcStartNormal(arc);
|
|
Ne = GetArcEndNormal(arc);
|
|
|
|
assert(eq(VecLen(v1) - 1.0, eps));
|
|
|
|
v2 = VecSub(C, S);
|
|
d = VecDotProd(v1, v2);
|
|
|
|
GetSegEqnData(seg, &a, &b, &c);
|
|
R = GetArcRadius(arc);
|
|
dist = AbsPntLineDist(a, b, c, C);
|
|
|
|
if (lt(R - dist, eps)) { //line is fully outside of circle
|
|
return 0;
|
|
}
|
|
|
|
/* */
|
|
/* check whether the line is tangential */
|
|
/* */
|
|
/* calculate the point on the segment that is nearest to the circle's */
|
|
/* center */
|
|
/* */
|
|
q1 = VecAdd(S, VecMult(d, v1));
|
|
if (PntCircleDist(C, R, q1) >= 0) {
|
|
/* */
|
|
/* the nearest point lies on the circle (or slightly outside) */
|
|
/* -> in this case the iterative solver won't work, because we don't */
|
|
/* know a point on the segment that is inside the circle */
|
|
/* -> assume the line is tangential and use the nearest point as */
|
|
/* a candidate for the intersection */
|
|
/* */
|
|
q3 = VecSub(q1, C);
|
|
if ((IsInArcConeStrict(Ns, Ne, q3) == 0) &&
|
|
(IsInSegConeStrict(S, q1, a, b, lgth) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
// IntWarning("seg", seg, "arc", arc);
|
|
return 1;
|
|
}
|
|
else {
|
|
return 0;
|
|
}
|
|
}
|
|
|
|
/* */
|
|
/* else we can use the iterative solver */
|
|
/* */
|
|
u1 = d - R - eps;
|
|
u2 = d + R + eps;
|
|
|
|
t = (u1 + d)/2;
|
|
success = IterativeIntersection(ARC, SEG, 0.0, 0.0, 0.0,
|
|
C, R, S, u1, d, v1, 0.0,
|
|
&t, &q1);
|
|
|
|
q3 = VecSub(q1, C);
|
|
|
|
if (success &&
|
|
(IsInArcConeStrict(Ns, Ne, q3) == 0) &&
|
|
(IsInSegConeStrict(S, q1, a, b, lgth) == 0) &&
|
|
(!shared || ne(PntPntDist(q1, shared_point), TINY))) {
|
|
*res1 = q1;
|
|
intersections++;
|
|
}
|
|
|
|
t = (u2 + d)/2;
|
|
success = IterativeIntersection(ARC, SEG, 0.0, 0.0, 0.0,
|
|
C, R, S, d, u2, v1, 0.0,
|
|
&t, &q2);
|
|
|
|
q3 = VecSub(q2, C);
|
|
|
|
if (success &&
|
|
(IsInArcConeStrict(Ns, Ne, q3) == 0) &&
|
|
(IsInSegConeStrict(S, q2, a, b, lgth) == 0) &&
|
|
(!shared || ne(PntPntDist(q2, shared_point), TINY))) {
|
|
if (intersections == 1) {
|
|
if (ne(PntPntDist(q1, q2), eps)) {
|
|
*res2 = q2;
|
|
intersections++;
|
|
}
|
|
}
|
|
else {
|
|
*res1 = q2;
|
|
intersections++;
|
|
}
|
|
}
|
|
|
|
if (intersections > 0)
|
|
IntWarning("seg", seg, "arc", arc);
|
|
|
|
return intersections;
|
|
}
|
|
#endif
|
|
|
|
|
|
vr_bool vroniObject::CheckIntersectionsLocally(int i1, t_site t1,
|
|
int i2, t_site t2,
|
|
int i3, t_site t3)
|
|
{
|
|
coord res1, res2;
|
|
vr_bool duplicate;
|
|
|
|
if (t1 == SEG) {
|
|
if (t2 == SEG) {
|
|
if (GetSegSegIntersection(i1, i2, &res1, eps_int, &duplicate)) {
|
|
SetInvalidSites(i1, SEG, i2, SEG, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
if (t3 == SEG) {
|
|
if (GetSegSegIntersection(i1, i3, &res1, eps_int, &duplicate)) {
|
|
SetInvalidSites(i1, SEG, i3, SEG, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
#ifdef GENUINE_ARCS
|
|
else if (t2 == ARC) {
|
|
if (GetSegArcIntersection(i1, i2, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i1, SEG, i2, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
else if (t3 == ARC) {
|
|
if (GetSegArcIntersection(i1, i3, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i1, SEG, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
#ifdef GENUINE_ARCS
|
|
else if (t1 == ARC) {
|
|
if (t2 == SEG) {
|
|
if (GetSegArcIntersection(i2, i1, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i1, SEG, i2, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
else if (t2 == ARC) {
|
|
if (GetArcArcIntersection(i1, i2, &res1, &res2, eps_int,
|
|
&duplicate)) {
|
|
SetInvalidSites(i1, ARC, i2, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
if (t3 == SEG) {
|
|
if (GetSegArcIntersection(i3, i1, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i1, SEG, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
else if (t3 == ARC) {
|
|
if (GetArcArcIntersection(i1, i3, &res1, &res2, eps_int,
|
|
&duplicate)) {
|
|
SetInvalidSites(i1, ARC, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
if (t2 == SEG) {
|
|
if (t3 == SEG) {
|
|
if (GetSegSegIntersection(i2, i3, &res1, eps_int, &duplicate)) {
|
|
SetInvalidSites(i2, SEG, i3, SEG, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
#ifdef GENUINE_ARCS
|
|
else if (t3 == ARC) {
|
|
if (GetSegArcIntersection(i2, i3, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i2, SEG, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
#ifdef GENUINE_ARCS
|
|
else if (t2 == ARC) {
|
|
if (t3 == SEG) {
|
|
if (GetSegArcIntersection(i3, i2, &res1, &res2, eps_int)) {
|
|
SetInvalidSites(i2, SEG, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
else if (t3 == ARC) {
|
|
if (GetArcArcIntersection(i2, i3, &res1, &res2, eps_int,
|
|
&duplicate)) {
|
|
SetInvalidSites(i2, ARC, i3, ARC, eps_int);
|
|
if (IsInvalidData(false)) return true;
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
/* */
|
|
/* this function finds all intersections among */
|
|
/* */
|
|
/* pnts[2,num_pnts-3], */
|
|
/* segs[0,num_segs-1], */
|
|
/* arcs[0,num_arcs-1]. */
|
|
/* */
|
|
/* */
|
|
/* a pnt and a seg or arc are reported as intersecting if the pnt lies on */
|
|
/* the seg or arc and if it is neither an end point of the seg or arc nor a */
|
|
/* supporting point of the seg or arc. */
|
|
/* */
|
|
/* note: this routine does not delete duplicates. (all duplicates will be */
|
|
/* discarded within the next call to CleanData(). */
|
|
/* */
|
|
vr_bool vroniObject::ComputeAllIntersections(double_arg eps)
|
|
{
|
|
int i, j;
|
|
coord res1, res2;
|
|
#ifdef GENUINE_ARCS
|
|
int cnt;
|
|
#endif
|
|
int intersections = 0;
|
|
grid* grid;
|
|
vr_bool ident = false, identical_sites = false;
|
|
|
|
VD_Warning("ComputeAllIntersections() - full check for intersections!");
|
|
|
|
comp_eps = eps;
|
|
|
|
grid = InitialiseGrid(num_pnts + num_segs + num_arcs, bb_min, bb_max);
|
|
FillGrid(grid, eps);
|
|
|
|
if (check_completed) {
|
|
check_completed = false;
|
|
if (verbose)
|
|
FP_printf("...checking point intersections with precision %2.1e",
|
|
FP_PRNTARG(eps));
|
|
|
|
for (i = 0; i < grid->n_x; i++) {
|
|
for (j = 0; j < grid->n_y; j++) {
|
|
grid_cell* cell;
|
|
std::set<point_entry*>::iterator it1, it2;
|
|
std::set<seg_entry*>::iterator it3;
|
|
#ifdef GENUINE_ARCS
|
|
std::set<arc_entry*>::iterator it4;
|
|
arc_entry *arc;
|
|
#endif
|
|
point_entry *p1, *p2;
|
|
seg_entry *seg;
|
|
|
|
cell = &grid->cells[i][j];
|
|
for (it1 = cell->points.begin(); it1 != cell->points.end(); it1++){
|
|
p1 = *it1;
|
|
if (IsDeletedPnt(p1->num)) continue;
|
|
|
|
//pnt - pnt
|
|
it2 = it1;
|
|
it2++;
|
|
for (; it2 != cell->points.end(); it2++) {
|
|
p2 = *it2;
|
|
if (IsDeletedPnt(p2->num)) continue;
|
|
|
|
if (CheckPntPnt(p1->num, p2->num, eps)) {
|
|
MergePoints(p1->num, p2->num);
|
|
intersections++;
|
|
}
|
|
}
|
|
|
|
//pnt - seg
|
|
for (it3 = cell->segs.begin(); it3 != cell->segs.end(); it3++) {
|
|
seg = *it3;
|
|
|
|
if (CheckPntOnSeg(seg->num, p1->num, eps)) {
|
|
if (AddSegIntersection(seg, GetPntCoords(p1->num)))
|
|
intersections++;
|
|
}
|
|
}
|
|
|
|
#ifdef GENUINE_ARCS
|
|
//pnt - arc
|
|
for (it4 = cell->arcs.begin(); it4 != cell->arcs.end(); it4++) {
|
|
arc = *it4;
|
|
|
|
if (CheckPntOnArc(arc->num, p1->num, eps)) {
|
|
if (AddArcIntersection(arc, GetPntCoords(p1->num)))
|
|
intersections++;
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
}
|
|
}
|
|
if (verbose) {
|
|
printf(" -- done!\n");
|
|
printf(" %d point intersections found!\n", intersections);
|
|
}
|
|
|
|
if (intersections > 0) {
|
|
SplitAtIntersections(grid, eps);
|
|
DeleteGrid(grid);
|
|
return true;
|
|
}
|
|
}
|
|
|
|
if (verbose)
|
|
FP_printf("...checking seg/arc intersections with precision %2.1e",
|
|
FP_PRNTARG(eps));
|
|
|
|
for (i = 0; i < grid->n_x; i++) {
|
|
for (j = 0; j < grid->n_y; j++) {
|
|
std::set<seg_entry*>::iterator it1, it2;
|
|
seg_entry *seg1, *seg2;
|
|
#ifdef GENUINE_ARCS
|
|
std::set<arc_entry*>::iterator it3, it4;
|
|
arc_entry *arc1, *arc2;
|
|
#endif
|
|
grid_cell* cell;
|
|
|
|
cell = &grid->cells[i][j];
|
|
|
|
for (it1 = cell->segs.begin(); it1 != cell->segs.end(); it1++) {
|
|
seg1 = *it1;
|
|
|
|
//seg - seg
|
|
it2 = it1;
|
|
it2++;
|
|
for (; it2 != cell->segs.end(); it2++) {
|
|
seg2 = *it2;
|
|
|
|
if (GetSegSegIntersection(seg1->num, seg2->num, &res1, eps,
|
|
&ident)) {
|
|
if (AddSegIntersection(seg1, res1) |
|
|
AddSegIntersection(seg2, res1))
|
|
intersections++;
|
|
}
|
|
else if (ident) {
|
|
ident = false;
|
|
identical_sites = true;
|
|
}
|
|
}
|
|
|
|
|
|
#ifdef GENUINE_ARCS
|
|
//seg - arc
|
|
for (it3 = cell->arcs.begin(); it3 != cell->arcs.end(); it3++) {
|
|
arc1 = *it3;
|
|
|
|
if ((cnt = GetSegArcIntersection(seg1->num, arc1->num,
|
|
&res1, &res2, eps)) > 0) {
|
|
if ( AddSegIntersection(seg1, res1) |
|
|
AddArcIntersection(arc1, res1))
|
|
intersections++;
|
|
if (cnt > 1) {
|
|
assert(cnt == 2);
|
|
if (AddSegIntersection(seg1, res2) |
|
|
AddArcIntersection(arc1, res2))
|
|
intersections++;
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
|
|
#ifdef GENUINE_ARCS
|
|
//arc - arc
|
|
for (it3 = cell->arcs.begin(); it3 != cell->arcs.end(); it3++) {
|
|
arc1 = *it3;
|
|
|
|
it4 = it3;
|
|
it4++;
|
|
for (; it4 != cell->arcs.end(); it4++) {
|
|
arc2 = *it4;
|
|
|
|
if ((cnt = GetArcArcIntersection(arc1->num, arc2->num, &res1,
|
|
&res2, eps, &ident)) > 0){
|
|
if (AddArcIntersection(arc1, res1) |
|
|
AddArcIntersection(arc2, res1))
|
|
intersections++;
|
|
|
|
if (cnt > 1) {
|
|
assert(cnt == 2);
|
|
if (AddArcIntersection(arc1, res2) |
|
|
AddArcIntersection(arc2, res2))
|
|
intersections++;
|
|
}
|
|
}
|
|
else if (ident) {
|
|
ident = false;
|
|
identical_sites = true;
|
|
}
|
|
}
|
|
}
|
|
#endif
|
|
}
|
|
}
|
|
|
|
SplitAtIntersections(grid, eps);
|
|
DeleteGrid(grid);
|
|
|
|
check_completed = true;
|
|
|
|
if (verbose) {
|
|
printf(" -- done!\n");
|
|
if ((intersections > 0) && identical_sites)
|
|
printf(" %d seg/arc intersections and at least one duplicate!\n",
|
|
intersections);
|
|
else if (intersections > 0)
|
|
printf(" %d seg/arc intersections!\n",
|
|
intersections);
|
|
else if (identical_sites)
|
|
printf(" at least one duplicate seg/arc!\n");
|
|
else
|
|
printf(" no seg/arc intersections and no duplicates!\n");
|
|
}
|
|
|
|
return ((intersections > 0) || identical_sites);
|
|
}
|