SaraP
332 lines
13 KiB
C++
332 lines
13 KiB
C++
/*****************************************************************************/
|
|
/* */
|
|
/* Copyright (C) 2007-2025 S. Huber, 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: Stefan Huber, Martin Held */
|
|
/* */
|
|
/* E-Mail: held@cs.sbg.ac.at */
|
|
/* Fax Mail: (+43 662) 8044-611 */
|
|
/* Voice Mail: (+43 662) 8044-6304 */
|
|
/* Snail Mail: Martin Held */
|
|
/* FB Informatik */
|
|
/* Universitaet Salzburg */
|
|
/* A-5020 Salzburg, Austria */
|
|
/* */
|
|
/*****************************************************************************/
|
|
|
|
/* */
|
|
/* get standard libraries */
|
|
/* */
|
|
#include <stdlib.h>
|
|
#include <stdio.h>
|
|
#include <assert.h>
|
|
#include <math.h>
|
|
#include <float.h>
|
|
|
|
/* */
|
|
/* get my header files */
|
|
/* */
|
|
#include "fpkernel.h"
|
|
#include "vronivector.h"
|
|
#include "vroni_object.h"
|
|
#include "defs.h"
|
|
#include "numerics.h"
|
|
#include "roots.h"
|
|
|
|
|
|
|
|
/* */
|
|
/* Calculate the Voronoi center of arc and two segs */
|
|
/* */
|
|
vr_bool vroniObject::ArcSegSegCntr(int i, int j, int k, int e, coord* cntr,
|
|
double *r2, vr_bool *problematic, int *site)
|
|
{
|
|
int m;
|
|
coord c1;
|
|
double rr1;
|
|
double s1a, s1b, s1c;
|
|
double s2a, s2b, s2c;
|
|
coord centers[VRONI_MAXSOL];
|
|
double radii[VRONI_MAXSOL];
|
|
int num_sol = 0, best_sol = NIL;
|
|
|
|
double eps = ZERO;
|
|
coord v, w, sp, ep, spi, epi, spj, epj;
|
|
vr_bool counter_tangential;
|
|
#ifdef VRONI_DBG_WARN
|
|
coord spk, epk;
|
|
#endif
|
|
|
|
/* */
|
|
/* Check whether the three cites are all in lists */
|
|
/* */
|
|
assert(InArcsList(i));
|
|
assert(InSegsList(j));
|
|
assert(InSegsList(k));
|
|
|
|
*problematic = false;
|
|
|
|
/* */
|
|
/* we sort the two seg indices in order to guarantee that the center of */
|
|
/* these three sites is always computed consistently. */
|
|
/* */
|
|
SortTwoNumbers(j, k, m);
|
|
|
|
/* */
|
|
/* check whether the three sites share a common end point */
|
|
/* */
|
|
m = Get1stVtx(k, SEG);
|
|
if (IsArcStartPnt(i, m) || IsArcEndPnt(i, m)) {
|
|
if (IsSegStartPnt(j, m) || IsSegEndPnt(j, m)) {
|
|
*cntr = GetPntCoords(m);
|
|
*r2 = 0.0;
|
|
*site = m;
|
|
return true; /* common end point */
|
|
}
|
|
}
|
|
m = Get2ndVtx(k, SEG);
|
|
if (IsArcStartPnt(i, m) || IsArcEndPnt(i, m)) {
|
|
if (IsSegStartPnt(j, m) || IsSegEndPnt(j, m)) {
|
|
*cntr = GetPntCoords(m);
|
|
*r2 = 0.0;
|
|
*site = m;
|
|
return true; /* common end point */
|
|
}
|
|
}
|
|
|
|
c1 = GetArcCenter(i);
|
|
rr1 = GetArcRadius(i);
|
|
|
|
GetSegEqnData(j, &s1a, &s1b, &s1c);
|
|
s1c*= -1.0;
|
|
|
|
GetSegEqnData(k, &s2a, &s2b, &s2c);
|
|
s2c*= -1.0;
|
|
|
|
if (((s1a < -ZERO) && (s2a > ZERO)) ||
|
|
((s1a > ZERO) && (s2a < ZERO))) {
|
|
s2a *= -1.0;
|
|
s2b *= -1.0;
|
|
s2c *= -1.0;
|
|
}
|
|
else if (((s1b < -ZERO) && (s2b > ZERO)) ||
|
|
((s1b > ZERO) && (s2b < ZERO))) {
|
|
s2a *= -1.0;
|
|
s2b *= -1.0;
|
|
s2c *= -1.0;
|
|
}
|
|
|
|
while (eps <= ZERO_MAX) {
|
|
num_sol = 0;
|
|
|
|
//If second segment is tangential to arc --> let it be segment j
|
|
if ((Abs(Abs(s2a*c1.x+s2b*c1.y-s2c)-rr1) <= eps) &&
|
|
((PntPntDist(GetArcStartCoord(i), GetSegStartCoord(k)) <= eps) ||
|
|
(PntPntDist(GetArcStartCoord(i), GetSegEndCoord(k)) <= eps) ||
|
|
(PntPntDist(GetArcEndCoord(i), GetSegStartCoord(k)) <= eps) ||
|
|
(PntPntDist(GetArcEndCoord(i), GetSegEndCoord(k)) <= eps)) &&
|
|
!(IsPntInArcCone(i,GetSegStartCoord(k)) &&
|
|
IsPntInArcCone(i,GetSegEndCoord(k)))) {
|
|
int tmpi;
|
|
double tmpd;
|
|
|
|
VroniSwap(j, k, tmpi);
|
|
VroniSwap(s1a, s2a, tmpd);
|
|
VroniSwap(s1b, s2b, tmpd);
|
|
VroniSwap(s1c, s2c, tmpd);
|
|
}
|
|
|
|
spi = GetArcStartCoord(i);
|
|
epi = GetArcEndCoord(i);
|
|
spj = GetSegStartCoord(j);
|
|
epj = GetSegEndCoord(j);
|
|
|
|
#ifdef TRACE
|
|
if ((i == 1305) && (((j == 1605) && (k == 1604)) ||
|
|
((j == 1604) && (k == 1605)))) {
|
|
printf("arc%d-seg%d-seg%d: eps = %20.16f\n", i, j, k, eps);
|
|
printf("start arc: (%20.16f %20.16f)\n", spi.x, spi.y);
|
|
printf("end arc: (%20.16f %20.16f)\n", epi.x, epi.y);
|
|
printf("center arc: (%20.16f %20.16f)\n", c1.x, c1.y);
|
|
printf("start seg 1: (%20.16f %20.16f)\n", spj.x, spj.y);
|
|
printf("end seg 1: (%20.16f %20.16f)\n", epj.x, epj.y);
|
|
printf("start seg 2: (%20.16f %20.16f)\n", GetSegStartCoord(k).x,
|
|
GetSegStartCoord(k).y);
|
|
printf("end seg 2: (%20.16f %20.16f)\n", GetSegEndCoord(k).x,
|
|
GetSegEndCoord(k).y);
|
|
}
|
|
#endif
|
|
|
|
//Segment tangential to arc!
|
|
//Do not change the eps to eps_MAX here! See arcs_347
|
|
if (Abs(Abs(s1a*c1.x+s1b*c1.y-s1c)-rr1) <= eps) {
|
|
counter_tangential = true;
|
|
if (PntPntDist(spi,spj) <= eps) {
|
|
ep = spj;
|
|
sp = epj;
|
|
v = VecSub(sp, ep);
|
|
w = VecSub(epi, ep);
|
|
counter_tangential = (VecDotProd(v, w) < 0.0) ? false : true;
|
|
}
|
|
else if (PntPntDist(epi,spj) <= eps) {
|
|
ep = spj;
|
|
sp = epj;
|
|
v = VecSub(sp, ep);
|
|
w = VecSub(spi, ep);
|
|
counter_tangential = (VecDotProd(v, w) < 0.0) ? false : true;
|
|
}
|
|
else if (PntPntDist(spi,epj) <= eps) {
|
|
ep = epj;
|
|
sp = spj;
|
|
v = VecSub(sp, ep);
|
|
w = VecSub(epi, ep);
|
|
counter_tangential = (VecDotProd(v, w) < 0.0) ? false : true;
|
|
}
|
|
else if (PntPntDist(epi,epj) <= eps) {
|
|
ep = epj;
|
|
sp = spj;
|
|
v = VecSub(sp, ep);
|
|
w = VecSub(spi, ep);
|
|
counter_tangential = (VecDotProd(v, w) < 0.0) ? false : true;
|
|
}
|
|
|
|
if (!counter_tangential) {
|
|
//Intersect the ray and get all centers
|
|
v.x = s1a;
|
|
v.y = s1b;
|
|
num_sol = IntersectRaySegment(ep, v, s2a, s2b, s2c, centers, radii,
|
|
eps);
|
|
}
|
|
else {
|
|
num_sol = CircSegSegCenters(c1, rr1, s1a,s1b,s1c, s2a,s2b,s2c,
|
|
centers, radii, eps);
|
|
}
|
|
}
|
|
//parallel segments with common endpoint
|
|
else if ((Abs(s1a-s2a) <= eps) && (Abs(s1b-s2b) <= eps) &&
|
|
(Abs(s1c-s2c) <= eps) &&
|
|
((PntPntDist(spj,GetSegStartCoord(k)) <= eps) ||
|
|
(PntPntDist(spj,GetSegEndCoord(k)) <= eps) ||
|
|
(PntPntDist(epj,GetSegStartCoord(k)) <= eps) ||
|
|
(PntPntDist(epj,GetSegEndCoord(k)) <= eps))) {
|
|
//Get the common endpoint
|
|
if( PntPntDist(spj, GetSegStartCoord(k))<= eps ||
|
|
PntPntDist(spj, GetSegEndCoord(k))<= eps )
|
|
ep = spj;
|
|
else
|
|
ep = epj;
|
|
v = MakeVec(s1a, s1b);
|
|
|
|
num_sol = IntersectRayCircle(ep, v, c1, rr1, centers, radii, eps);
|
|
}
|
|
else {
|
|
//Calculate all equidistant points
|
|
num_sol = CircSegSegCenters(c1, rr1, s1a,s1b,s1c, s2a,s2b,s2c,
|
|
centers, radii, eps);
|
|
}
|
|
|
|
/* */
|
|
/* classify all solutions and pick the best solution */
|
|
/* */
|
|
best_sol = ClassifySolutions(i, j, k, e, ARC, SEG, SEG,
|
|
false, true, true, centers, radii,
|
|
&num_sol, eps, problematic);
|
|
assert((0 <= best_sol) && (best_sol < num_sol));
|
|
#ifdef TRACE
|
|
if ((i == 40) && (((j == 2) && (k == 42)) ||
|
|
((j == 42) && (k == 2)))) {
|
|
printf("arc %d, seg %d, seg %d\n", i, j, k);
|
|
for (m = 0; m < num_sol; ++m) {
|
|
printf("centers[%d]: (%20.16f %20.16f)\n", m,
|
|
centers[m].x, centers[m].y);
|
|
printf("radii[%d]: %20.16f\n", m, radii[m]);
|
|
}
|
|
printf("problematic = %d, eps = %20.16f\n", *problematic, eps);
|
|
}
|
|
#endif
|
|
if (!*problematic) {
|
|
*cntr = centers[best_sol];
|
|
*r2 = radii[best_sol];
|
|
return true;
|
|
}
|
|
else {
|
|
eps *= 10.0;
|
|
}
|
|
}
|
|
|
|
if (eq(rr1, ZERO_MAX)) {
|
|
/* */
|
|
/* this is an arc with a very small radius; we replace it */
|
|
/* by a seg */
|
|
/* */
|
|
VD_Dbg_Warning("ArcSegSegCntr() - arc with small radius!");
|
|
ReplaceArc(i);
|
|
restart = true;
|
|
return false;
|
|
}
|
|
else if (eq(PntPntDist(GetArcStartCoord(i), GetArcEndCoord(i)), ZERO_MAX)) {
|
|
/* */
|
|
/* this is an arc with a very short cord; we replace it by a seg */
|
|
/* */
|
|
VD_Dbg_Warning("ArcArcArcCntr() - arc with short cord!");
|
|
ReplaceArc(i);
|
|
restart = true;
|
|
return false;
|
|
}
|
|
else if (eq(GetSegLgth(j), ZERO_MAX)) {
|
|
/* */
|
|
/* this is a seg with a very small length; we replace it */
|
|
/* */
|
|
VD_Dbg_Warning("ArcSegSegCntr() - seg with small length!");
|
|
ReplaceSeg(j);
|
|
restart = true;
|
|
return false;
|
|
}
|
|
else if (eq(GetSegLgth(k), ZERO_MAX)) {
|
|
/* */
|
|
/* this is a seg with a very small length; we replace it */
|
|
/* */
|
|
VD_Dbg_Warning("ArcSegSegCntr() - seg with small length!");
|
|
ReplaceSeg(k);
|
|
restart = true;
|
|
return false;
|
|
}
|
|
else if (CheckIntersectionsLocally(i, ARC, j, SEG, k, SEG)) {
|
|
restart = true;
|
|
return false;
|
|
}
|
|
|
|
#ifdef VRONI_DBG_WARN
|
|
spi = GetArcStartCoord(i);
|
|
epi = GetArcEndCoord(i);
|
|
c1 = GetArcCenter(i);
|
|
printf("\nCenter not computed for arc-seg-seg: %d-%d-%d\n",
|
|
i, j, k);
|
|
printf("%20.16f %20.16f %20.16f %20.16f %20.16f %20.16f\n",
|
|
spi.x, spi.y, epi.x, epi.y, c1.x, c1.y);
|
|
epj = GetSegStartCoord(j);
|
|
spj = GetSegEndCoord(j);
|
|
epk = GetSegStartCoord(k);
|
|
spk = GetSegEndCoord(k);
|
|
printf("%20.16f %20.16f %20.16f %20.16f\n",
|
|
spj.x, spj.y, epj.x, epj.y);
|
|
printf("%20.16f %20.16f %20.16f %20.16f\n",
|
|
spk.x, spk.y, epk.x, epk.y);
|
|
#endif
|
|
|
|
*cntr = centers[best_sol];
|
|
*r2 = radii[best_sol];
|
|
restart = false;
|
|
|
|
return false;
|
|
}
|
|
|
|
|
|
|