Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
#include "tcurveutil.h"
Toshihiro Shimizu 890ddd
#include "tcurves.h"
Toshihiro Shimizu 890ddd
#include "tmathutil.h"
Toshihiro Shimizu 890ddd
#include "tbezier.h"
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
/*
Toshihiro Shimizu 890ddd
Questa funzione ritorna un vettore di
Toshihiro Shimizu 890ddd
coppie di double (DoublePair) che individua i parametri 
Toshihiro Shimizu 890ddd
dei punti d'intersezione.
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  L'intero restituito indica il numero d'intersezioni che 
Toshihiro Shimizu 890ddd
  sono state individuate (per due segmenti una).
Toshihiro Shimizu 890ddd
  
Toshihiro Shimizu 890ddd
    Se i segmenti sono paralleli il parametro viene posto a -1.
Toshihiro Shimizu 890ddd
*/
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
int intersect(const TSegment &first,
Toshihiro Shimizu 890ddd
			  const TSegment &second,
Toshihiro Shimizu 890ddd
			  std::vector<doublepair> &intersections)</doublepair>
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	return intersect(first.getP0(), first.getP1(), second.getP0(), second.getP1(), intersections);
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
int intersect(const TPointD &p1, const TPointD &p2, const TPointD &p3, const TPointD &p4,
Toshihiro Shimizu 890ddd
			  std::vector<doublepair> &intersections)</doublepair>
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	// This algorithm is presented in Graphics Geems III pag 199
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	static double Ax, Bx, Ay, By, Cx, Cy, d, f, e;
Toshihiro Shimizu 890ddd
	static double x1lo, x1hi, y1lo, y1hi;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	Ax = p2.x - p1.x;
Toshihiro Shimizu 890ddd
	Bx = p3.x - p4.x;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	//test delle BBox
Toshihiro Shimizu 890ddd
	if (Ax < 0.0) {
Toshihiro Shimizu 890ddd
		x1lo = p2.x;
Toshihiro Shimizu 890ddd
		x1hi = p1.x;
Toshihiro Shimizu 890ddd
	} else {
Toshihiro Shimizu 890ddd
		x1lo = p1.x;
Toshihiro Shimizu 890ddd
		x1hi = p2.x;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (Bx > 0.0) {
Toshihiro Shimizu 890ddd
		if (x1hi < p4.x || x1lo > p3.x)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
	} else if (x1hi < p3.x || x1lo > p4.x)
Toshihiro Shimizu 890ddd
		return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	Ay = p2.y - p1.y;
Toshihiro Shimizu 890ddd
	By = p3.y - p4.y;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (Ay < 0) {
Toshihiro Shimizu 890ddd
		y1lo = p2.y;
Toshihiro Shimizu 890ddd
		y1hi = p1.y;
Toshihiro Shimizu 890ddd
	} else {
Toshihiro Shimizu 890ddd
		y1lo = p1.y;
Toshihiro Shimizu 890ddd
		y1hi = p2.y;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (By > 0) {
Toshihiro Shimizu 890ddd
		if (y1hi < p4.y || y1lo > p3.y)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
	} else if (y1hi < p3.y || y1lo > p4.y)
Toshihiro Shimizu 890ddd
		return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	Cx = p1.x - p3.x;
Toshihiro Shimizu 890ddd
	Cy = p1.y - p3.y;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	d = By * Cx - Bx * Cy;
Toshihiro Shimizu 890ddd
	f = Ay * Bx - Ax * By;
Toshihiro Shimizu 890ddd
	e = Ax * Cy - Ay * Cx;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (f > 0) {
Toshihiro Shimizu 890ddd
		if (d < 0)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (!areAlmostEqual(d, f))
Toshihiro Shimizu 890ddd
			if (d > f)
Toshihiro Shimizu 890ddd
				return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (e < 0)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
		if (!areAlmostEqual(e, f))
Toshihiro Shimizu 890ddd
			if (e > f)
Toshihiro Shimizu 890ddd
				return 0;
Toshihiro Shimizu 890ddd
	} else if (f < 0) {
Toshihiro Shimizu 890ddd
		if (d > 0)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (!areAlmostEqual(d, f))
Toshihiro Shimizu 890ddd
			if (d < f)
Toshihiro Shimizu 890ddd
				return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (e > 0)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
		if (!areAlmostEqual(e, f))
Toshihiro Shimizu 890ddd
			if (e < f)
Toshihiro Shimizu 890ddd
				return 0;
Toshihiro Shimizu 890ddd
	} else {
Toshihiro Shimizu 890ddd
		if (d < 0 || d > 1 || e < 0 || e > 1)
Toshihiro Shimizu 890ddd
			return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (p1 == p2 && p3 == p4) {
Toshihiro Shimizu 890ddd
			intersections.push_back(DoublePair(0, 0));
Toshihiro Shimizu 890ddd
			return 1;
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		// controllo che i segmenti non siano sulla stessa retta
Toshihiro Shimizu 890ddd
		if (!cross(p2 - p1, p4 - p1)) {
Toshihiro Shimizu 890ddd
			// calcolo delle combinazioni baricentriche
Toshihiro Shimizu 890ddd
			double distp2p1 = norm2(p2 - p1);
Toshihiro Shimizu 890ddd
			double distp3p4 = norm2(p3 - p4);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
			double dist2_p3p1 = norm2(p3 - p1);
Toshihiro Shimizu 890ddd
			double dist2_p4p1 = norm2(p4 - p1);
Toshihiro Shimizu 890ddd
			double dist2_p3p2 = norm2(p3 - p2);
Toshihiro Shimizu 890ddd
			double dist2_p4p2 = norm2(p4 - p2);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
			int intersection = 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
			// calcolo delle prime due soluzioni
Toshihiro Shimizu 890ddd
			double vol1;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
			if (distp3p4) {
Toshihiro Shimizu 890ddd
				distp3p4 = sqrt(distp3p4);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				vol1 = (p1 - p3) * normalize(p4 - p3);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				if (vol1 >= 0 && vol1 <= distp3p4) // combinazione baricentrica valida
Toshihiro Shimizu 890ddd
				{
Toshihiro Shimizu 890ddd
					intersections.push_back(DoublePair(0.0, vol1 / distp3p4));
Toshihiro Shimizu 890ddd
					++intersection;
Toshihiro Shimizu 890ddd
				}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				vol1 = (p2 - p3) * normalize(p4 - p3);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				if (vol1 >= 0 && vol1 <= distp3p4) {
Toshihiro Shimizu 890ddd
					intersections.push_back(DoublePair(1.0, vol1 / distp3p4));
Toshihiro Shimizu 890ddd
					++intersection;
Toshihiro Shimizu 890ddd
				}
Toshihiro Shimizu 890ddd
			}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
			if (distp2p1) {
Toshihiro Shimizu 890ddd
				distp2p1 = sqrt(distp2p1);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				vol1 = (p3 - p1) * normalize(p2 - p1);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				if (dist2_p3p2 && dist2_p3p1)
Toshihiro Shimizu 890ddd
					if (vol1 >= 0 && vol1 <= distp2p1) {
Toshihiro Shimizu 890ddd
						intersections.push_back(DoublePair(vol1 / distp2p1, 0.0));
Toshihiro Shimizu 890ddd
						++intersection;
Toshihiro Shimizu 890ddd
					}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				vol1 = (p4 - p1) * normalize(p2 - p1);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				if (dist2_p4p2 && dist2_p4p1)
Toshihiro Shimizu 890ddd
					if (vol1 >= 0 && vol1 <= distp2p1) {
Toshihiro Shimizu 890ddd
						intersections.push_back(DoublePair(vol1 / distp2p1, 1.0));
Toshihiro Shimizu 890ddd
						++intersection;
Toshihiro Shimizu 890ddd
					}
Toshihiro Shimizu 890ddd
			}
Toshihiro Shimizu 890ddd
			return intersection;
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
		return -1;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double par_s = d / f;
Toshihiro Shimizu 890ddd
	double par_t = e / f;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	intersections.push_back(DoublePair(par_s, par_t));
Toshihiro Shimizu 890ddd
	return 1;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//------------------------------------------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
int intersectCloseControlPoints(const TQuadratic &c0,
Toshihiro Shimizu 890ddd
								const TQuadratic &c1,
Toshihiro Shimizu 890ddd
								std::vector<doublepair> &intersections);</doublepair>
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
int intersect(const TQuadratic &c0,
Toshihiro Shimizu 890ddd
			  const TQuadratic &c1,
Toshihiro Shimizu 890ddd
			  std::vector<doublepair> &intersections, bool checksegments)</doublepair>
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	int ret;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	// funziona male, a volte toppa le intersezioni...
Toshihiro Shimizu 890ddd
	if (checksegments) {
Toshihiro Shimizu 890ddd
		ret = intersectCloseControlPoints(c0, c1, intersections);
Toshihiro Shimizu 890ddd
		if (ret != -2)
Toshihiro Shimizu 890ddd
			return ret;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double a = c0.getP0().x - 2 * c0.getP1().x + c0.getP2().x;
Toshihiro Shimizu 890ddd
	double b = 2 * (c0.getP1().x - c0.getP0().x);
Toshihiro Shimizu 890ddd
	double d = c0.getP0().y - 2 * c0.getP1().y + c0.getP2().y;
Toshihiro Shimizu 890ddd
	double e = 2 * (c0.getP1().y - c0.getP0().y);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double coeff = b * d - a * e;
Toshihiro Shimizu 890ddd
	int i = 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (areAlmostEqual(coeff, 0.0)) //c0 is a Segment, or a single point!!!
Toshihiro Shimizu 890ddd
	{
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		TSegment s = TSegment(c0.getP0(), c0.getP2());
Toshihiro Shimizu 890ddd
		ret = intersect(s, c1, intersections);
Toshihiro Shimizu 890ddd
		if (a == 0 && d == 0) //values of t in s coincide with values of t in c0
Toshihiro Shimizu 890ddd
			return ret;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		for (i = intersections.size() - ret; i < (int)intersections.size(); i++) {
Toshihiro Shimizu 890ddd
			intersections[i].first = c0.getT(s.getPoint(intersections[i].first));
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
		return ret;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double c = c0.getP0().x;
Toshihiro Shimizu 890ddd
	double f = c0.getP0().y;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double g = c1.getP0().x - 2 * c1.getP1().x + c1.getP2().x;
Toshihiro Shimizu 890ddd
	double h = 2 * (c1.getP1().x - c1.getP0().x);
Toshihiro Shimizu 890ddd
	double k = c1.getP0().x;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double m = c1.getP0().y - 2 * c1.getP1().y + c1.getP2().y;
Toshihiro Shimizu 890ddd
	double p = 2 * (c1.getP1().y - c1.getP0().y);
Toshihiro Shimizu 890ddd
	double q = c1.getP0().y;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (areAlmostEqual(h * m - g * p, 0.0)) //c1 is a Segment, or a single point!!!
Toshihiro Shimizu 890ddd
	{
Toshihiro Shimizu 890ddd
		TSegment s = TSegment(c1.getP0(), c1.getP2());
Toshihiro Shimizu 890ddd
		ret = intersect(c0, s, intersections);
Toshihiro Shimizu 890ddd
		if (g == 0 && m == 0) //values of t in s coincide with values of t in c0
Toshihiro Shimizu 890ddd
			return ret;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		for (i = intersections.size() - ret; i < (int)intersections.size(); i++) {
Toshihiro Shimizu 890ddd
			intersections[i].second = c1.getT(s.getPoint(intersections[i].second));
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
		return ret;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double a2 = (g * d - a * m);
Toshihiro Shimizu 890ddd
	double b2 = (h * d - a * p);
Toshihiro Shimizu 890ddd
	double c2 = ((k - c) * d + (f - q) * a);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	coeff = 1.0 / coeff;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double A = (a * a + d * d) * coeff * coeff;
Toshihiro Shimizu 890ddd
	double aux = A * c2 + (a * b + d * e) * coeff;
Toshihiro Shimizu 890ddd
Shinya Kitaoka 3bfa54
	std::vector<double> t;</double>
Shinya Kitaoka 3bfa54
	std::vector<double> solutions;</double>
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	t.push_back(aux * c2 + a * c + d * f - k * a - d * q);
Toshihiro Shimizu 890ddd
	aux += A * c2;
Toshihiro Shimizu 890ddd
	t.push_back(aux * b2 - h * a - d * p);
Toshihiro Shimizu 890ddd
	t.push_back(aux * a2 + A * b2 * b2 - g * a - d * m);
Toshihiro Shimizu 890ddd
	t.push_back(2 * A * a2 * b2);
Toshihiro Shimizu 890ddd
	t.push_back(A * a2 * a2);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	rootFinding(t, solutions);
Toshihiro Shimizu 890ddd
	//  solutions.push_back(0.0); //per convenzione; un valore vale l'altro....
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	for (i = 0; i < (int)solutions.size(); i++) {
Toshihiro Shimizu 890ddd
		if (solutions[i] < 0) {
Toshihiro Shimizu 890ddd
			if (areAlmostEqual(solutions[i], 0, 1e-6))
Toshihiro Shimizu 890ddd
				solutions[i] = 0;
Toshihiro Shimizu 890ddd
			else
Toshihiro Shimizu 890ddd
				continue;
Toshihiro Shimizu 890ddd
		} else if (solutions[i] > 1) {
Toshihiro Shimizu 890ddd
			if (areAlmostEqual(solutions[i], 1, 1e-6))
Toshihiro Shimizu 890ddd
				solutions[i] = 1;
Toshihiro Shimizu 890ddd
			else
Toshihiro Shimizu 890ddd
				continue;
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		DoublePair tt;
Toshihiro Shimizu 890ddd
		tt.second = solutions[i];
Toshihiro Shimizu 890ddd
		tt.first = coeff * (tt.second * (a2 * tt.second + b2) + c2);
Toshihiro Shimizu 890ddd
		if (tt.first < 0) {
Toshihiro Shimizu 890ddd
			if (areAlmostEqual(tt.first, 0, 1e-6))
Toshihiro Shimizu 890ddd
				tt.first = 0;
Toshihiro Shimizu 890ddd
			else
Toshihiro Shimizu 890ddd
				continue;
Toshihiro Shimizu 890ddd
		} else if (tt.first > 1) {
Toshihiro Shimizu 890ddd
			if (areAlmostEqual(tt.first, 1, 1e-6))
Toshihiro Shimizu 890ddd
				tt.first = 1;
Toshihiro Shimizu 890ddd
			else
Toshihiro Shimizu 890ddd
				continue;
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		intersections.push_back(tt);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		assert(areAlmostEqual(c0.getPoint(tt.first), c1.getPoint(tt.second), 1e-1));
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
	return intersections.size();
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
//questa funzione verifica se il punto di controllo p1 e' molto vicino a p0 o a p2:
Toshihiro Shimizu 890ddd
//in tal caso, si approssima la quadratica al segmento p0-p2.
Toshihiro Shimizu 890ddd
//se p1 e' vicino a p0, la relazione che lega il t del segmento al t della quadratica originaria e' tq = sqrt(ts),
Toshihiro Shimizu 890ddd
//se p1 e' vicino a p2, invece e' tq = 1-sqrt(1-ts).
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
int intersectCloseControlPoints(const TQuadratic &c0,
Toshihiro Shimizu 890ddd
								const TQuadratic &c1,
Toshihiro Shimizu 890ddd
								std::vector<doublepair> &intersections)</doublepair>
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	int ret = -2;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double dist1 = tdistance2(c0.getP0(), c0.getP1());
Toshihiro Shimizu 890ddd
	if (dist1 == 0)
Toshihiro Shimizu 890ddd
		dist1 = 1e-20;
Toshihiro Shimizu 890ddd
	double dist2 = tdistance2(c0.getP1(), c0.getP2());
Toshihiro Shimizu 890ddd
	if (dist2 == 0)
Toshihiro Shimizu 890ddd
		dist2 = 1e-20;
Toshihiro Shimizu 890ddd
	double val0 = tmax(dist1, dist2) / tmin(dist1, dist2);
Toshihiro Shimizu 890ddd
	double dist3 = tdistance2(c1.getP0(), c1.getP1());
Toshihiro Shimizu 890ddd
	if (dist3 == 0)
Toshihiro Shimizu 890ddd
		dist3 = 1e-20;
Toshihiro Shimizu 890ddd
	double dist4 = tdistance2(c1.getP1(), c1.getP2());
Toshihiro Shimizu 890ddd
	if (dist4 == 0)
Toshihiro Shimizu 890ddd
		dist4 = 1e-20;
Toshihiro Shimizu 890ddd
	double val1 = tmax(dist3, dist4) / tmin(dist3, dist4);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (val0 > 1000000 && val1 > 1000000) //entrambe c0 e c1  approssimate a segmenti
Toshihiro Shimizu 890ddd
	{
Toshihiro Shimizu 890ddd
		TSegment s0 = TSegment(c0.getP0(), c0.getP2());
Toshihiro Shimizu 890ddd
		TSegment s1 = TSegment(c1.getP0(), c1.getP2());
Toshihiro Shimizu 890ddd
		ret = intersect(s0, s1, intersections);
Toshihiro Shimizu 890ddd
		for (UINT i = intersections.size() - ret; i < (int)intersections.size(); i++) {
Toshihiro Shimizu 890ddd
			intersections[i].first = (dist1 < dist2) ? sqrt(intersections[i].first) : 1 - sqrt(1 - intersections[i].first);
Toshihiro Shimizu 890ddd
			intersections[i].second = (dist3 < dist4) ? sqrt(intersections[i].second) : 1 - sqrt(1 - intersections[i].second);
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
		//return ret;
Toshihiro Shimizu 890ddd
	} else if (val0 > 1000000) //solo c0 approssimata  a segmento
Toshihiro Shimizu 890ddd
	{
Toshihiro Shimizu 890ddd
		TSegment s0 = TSegment(c0.getP0(), c0.getP2());
Toshihiro Shimizu 890ddd
		ret = intersect(s0, c1, intersections);
Toshihiro Shimizu 890ddd
		for (UINT i = intersections.size() - ret; i < (int)intersections.size(); i++)
Toshihiro Shimizu 890ddd
			intersections[i].first = (dist1 < dist2) ? sqrt(intersections[i].first) : 1 - sqrt(1 - intersections[i].first);
Toshihiro Shimizu 890ddd
		//return ret;
Toshihiro Shimizu 890ddd
	} else if (val1 > 1000000) //solo c1 approssimata  a segmento
Toshihiro Shimizu 890ddd
	{
Toshihiro Shimizu 890ddd
		TSegment s1 = TSegment(c1.getP0(), c1.getP2());
Toshihiro Shimizu 890ddd
		ret = intersect(c0, s1, intersections);
Toshihiro Shimizu 890ddd
		for (UINT i = intersections.size() - ret; i < (int)intersections.size(); i++)
Toshihiro Shimizu 890ddd
			intersections[i].second = (dist3 < dist4) ? sqrt(intersections[i].second) : 1 - sqrt(1 - intersections[i].second);
Toshihiro Shimizu 890ddd
		//return ret;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	/*
Toshihiro Shimizu 890ddd
if (ret!=-2)
Toshihiro Shimizu 890ddd
  {
Toshihiro Shimizu 890ddd
  std::vector<doublepair> intersections1;</doublepair>
Toshihiro Shimizu 890ddd
  int ret1 = intersect(c0, c1, intersections1, false);
Toshihiro Shimizu 890ddd
  if (ret1>ret)
Toshihiro Shimizu 890ddd
    {
Toshihiro Shimizu 890ddd
    intersections = intersections1;
Toshihiro Shimizu 890ddd
    return ret1;
Toshihiro Shimizu 890ddd
    }
Toshihiro Shimizu 890ddd
  }
Toshihiro Shimizu 890ddd
*/
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return ret;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
int intersect(const TQuadratic &q,
Toshihiro Shimizu 890ddd
			  const TSegment &s,
Toshihiro Shimizu 890ddd
			  std::vector<doublepair> &intersections,</doublepair>
Toshihiro Shimizu 890ddd
			  bool firstIsQuad)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	int solutionNumber = 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	// nota la retta a*x+b*y+c = 0 andiamo alla ricerca delle soluzioni
Toshihiro Shimizu 890ddd
	//  di a*x(t)+b*y(t)+c=0 in [0,1]
Toshihiro Shimizu 890ddd
	double
Toshihiro Shimizu 890ddd
		a = s.getP0().y - s.getP1().y,
Toshihiro Shimizu 890ddd
		b = s.getP1().x - s.getP0().x,
Toshihiro Shimizu 890ddd
		c = -(a * s.getP0().x + b * s.getP0().y);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	// se il segmento e' un punto
Toshihiro Shimizu 890ddd
	if (0.0 == a && 0.0 == b) {
Toshihiro Shimizu 890ddd
		double outParForQuad = q.getT(s.getP0());
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
		if (areAlmostEqual(q.getPoint(outParForQuad), s.getP0())) {
Toshihiro Shimizu 890ddd
			if (firstIsQuad)
Toshihiro Shimizu 890ddd
				intersections.push_back(DoublePair(outParForQuad, 0));
Toshihiro Shimizu 890ddd
			else
Toshihiro Shimizu 890ddd
				intersections.push_back(DoublePair(0, outParForQuad));
Toshihiro Shimizu 890ddd
			return 1;
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
		return 0;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Michał Janiszewski b1cc3c
	if (q.getP2() - q.getP1() == q.getP1() - q.getP0()) {//pure il secondo e' unsegmento....
Toshihiro Shimizu 890ddd
		if (firstIsQuad)
Toshihiro Shimizu 890ddd
			return intersect(TSegment(q.getP0(), q.getP2()), s, intersections);
Toshihiro Shimizu 890ddd
		else
Toshihiro Shimizu 890ddd
			return intersect(s, TSegment(q.getP0(), q.getP2()), intersections);
Michał Janiszewski b1cc3c
	}
Toshihiro Shimizu 890ddd
Shinya Kitaoka 3bfa54
	std::vector<tpointd> bez, pol;</tpointd>
Toshihiro Shimizu 890ddd
	bez.push_back(q.getP0());
Toshihiro Shimizu 890ddd
	bez.push_back(q.getP1());
Toshihiro Shimizu 890ddd
	bez.push_back(q.getP2());
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	bezier2poly(bez, pol);
Toshihiro Shimizu 890ddd
Shinya Kitaoka 3bfa54
	std::vector<double> poly_1(3, 0), sol;</double>
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	poly_1[0] = a * pol[0].x + b * pol[0].y + c;
Toshihiro Shimizu 890ddd
	poly_1[1] = a * pol[1].x + b * pol[1].y;
Toshihiro Shimizu 890ddd
	poly_1[2] = a * pol[2].x + b * pol[2].y;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (!(rootFinding(poly_1, sol)))
Toshihiro Shimizu 890ddd
		return 0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double
Toshihiro Shimizu 890ddd
		segmentPar,
Toshihiro Shimizu 890ddd
		solution;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	TPointD v10(s.getP1() - s.getP0());
Toshihiro Shimizu 890ddd
	for (UINT i = 0; i < sol.size(); ++i) {
Toshihiro Shimizu 890ddd
		solution = sol[i];
Toshihiro Shimizu 890ddd
		if ((0.0 <= solution && solution <= 1.0) ||
Toshihiro Shimizu 890ddd
			areAlmostEqual(solution, 0.0, 1e-6) ||
Toshihiro Shimizu 890ddd
			areAlmostEqual(solution, 1.0, 1e-6)) {
Toshihiro Shimizu 890ddd
			segmentPar = (q.getPoint(solution) - s.getP0()) * v10 / (v10 * v10);
Toshihiro Shimizu 890ddd
			if ((0.0 <= segmentPar && segmentPar <= 1.0) ||
Toshihiro Shimizu 890ddd
				areAlmostEqual(segmentPar, 0.0, 1e-6) ||
Toshihiro Shimizu 890ddd
				areAlmostEqual(segmentPar, 1.0, 1e-6)) {
Toshihiro Shimizu 890ddd
				TPointD p1 = q.getPoint(solution);
Toshihiro Shimizu 890ddd
				TPointD p2 = s.getPoint(segmentPar);
Toshihiro Shimizu 890ddd
				assert(areAlmostEqual(p1, p2, 1e-1));
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
				if (firstIsQuad)
Toshihiro Shimizu 890ddd
					intersections.push_back(DoublePair(solution, segmentPar));
Toshihiro Shimizu 890ddd
				else
Toshihiro Shimizu 890ddd
					intersections.push_back(DoublePair(segmentPar, solution));
Toshihiro Shimizu 890ddd
				solutionNumber++;
Toshihiro Shimizu 890ddd
			}
Toshihiro Shimizu 890ddd
		}
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return solutionNumber;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
bool isCloseToSegment(const TPointD &point, const TSegment &segment, double distance)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	TPointD a = segment.getP0();
Toshihiro Shimizu 890ddd
	TPointD b = segment.getP1();
Toshihiro Shimizu 890ddd
	double lenght2 = tdistance2(a, b);
Toshihiro Shimizu 890ddd
	if (lenght2 < tdistance2(a, point) || lenght2 < tdistance2(point, b))
Toshihiro Shimizu 890ddd
		return false;
Toshihiro Shimizu 890ddd
	if (a.x == b.x)
Toshihiro Shimizu 890ddd
		return fabs(point.x - a.x) <= distance;
Toshihiro Shimizu 890ddd
	if (a.y == b.y)
Toshihiro Shimizu 890ddd
		return fabs(point.y - a.y) <= distance;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	// y=mx+q
Toshihiro Shimizu 890ddd
	double m = (a.y - b.y) / (a.x - b.x);
Toshihiro Shimizu 890ddd
	double q = a.y - (m * a.x);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double d2 = pow(fabs(point.y - (m * point.x) - q), 2) / (1 + (m * m));
Toshihiro Shimizu 890ddd
	return d2 <= distance * distance;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
double tdistance(const TSegment &segment, const TPointD &point)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	TPointD v1 = segment.getP1() - segment.getP0();
Toshihiro Shimizu 890ddd
	TPointD v2 = point - segment.getP0();
Toshihiro Shimizu 890ddd
	TPointD v3 = point - segment.getP1();
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (v2 * v1 <= 0)
Toshihiro Shimizu 890ddd
		return tdistance(point, segment.getP0());
Toshihiro Shimizu 890ddd
	else if (v3 * v1 >= 0)
Toshihiro Shimizu 890ddd
		return tdistance(point, segment.getP1());
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return fabs(v2 * rotate90(normalize(v1)));
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//-----------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
/*
Toshihiro Shimizu 890ddd
This formule is derived from Graphic Gems pag. 600
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  e = h^2 |a|/8
Toshihiro Shimizu 890ddd
  
Toshihiro Shimizu 890ddd
    e = pixel size
Toshihiro Shimizu 890ddd
    h = step
Toshihiro Shimizu 890ddd
    a = acceleration of curve (for a quadratic is a costant value)
Toshihiro Shimizu 890ddd
*/
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
double computeStep(const TQuadratic &quad, double pixelSize)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	double step = 2;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	TPointD A = quad.getP0() - 2.0 * quad.getP1() + quad.getP2(); // 2*A is the acceleration of the curve
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double A_len = norm(A);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	/*
Toshihiro Shimizu 890ddd
  A_len is equal to 2*norm(a)
Toshihiro Shimizu 890ddd
  pixelSize will be 0.5*pixelSize
Toshihiro Shimizu 890ddd
  now h is equal to sqrt( 8 * 0.5 * pixelSize / (2*norm(a)) ) = sqrt(2) * sqrt( pixelSize/A_len )
Toshihiro Shimizu 890ddd
  */
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (A_len > 0)
Toshihiro Shimizu 890ddd
		step = TConsts::sqrt2 * sqrt(pixelSize / A_len);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return step;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//-----------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
double computeStep(const TThickQuadratic &quad, double pixelSize)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	TThickPoint
Toshihiro Shimizu 890ddd
		cp0 = quad.getThickP0(),
Toshihiro Shimizu 890ddd
		cp1 = quad.getThickP1(),
Toshihiro Shimizu 890ddd
		cp2 = quad.getThickP2();
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	TQuadratic
Toshihiro Shimizu 890ddd
		q1(TPointD(cp0.x, cp0.y), TPointD(cp1.x, cp1.y), TPointD(cp2.x, cp2.y)),
Toshihiro Shimizu 890ddd
		q2(TPointD(cp0.y, cp0.thick), TPointD(cp1.y, cp1.thick), TPointD(cp2.y, cp2.thick)),
Toshihiro Shimizu 890ddd
		q3(TPointD(cp0.x, cp0.thick), TPointD(cp1.x, cp1.thick), TPointD(cp2.x, cp2.thick));
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return tmin(computeStep(q1, pixelSize), computeStep(q2, pixelSize), computeStep(q3, pixelSize));
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//=============================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
/*
Toshihiro Shimizu 890ddd
  Explanation: The length of a Bezier quadratic can be calculated explicitly.
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  Let Q be the quadratic. The tricks to explicitly integrate | Q'(t) | are:
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    - The integrand can be reformulated as:  | Q'(t) | = sqrt(at^2 + bt + c);
Toshihiro Shimizu 890ddd
    - Complete the square beneath the sqrt (add/subtract sq(b) / 4a)
Toshihiro Shimizu 890ddd
      and perform a linear variable change. We reduce the integrand to: sqrt(kx^2 + k),
Toshihiro Shimizu 890ddd
      where k can be taken outside => sqrt(x^2 + 1)
Toshihiro Shimizu 890ddd
    - Use x = tan y. The integrand will yield sec^3 y.
Toshihiro Shimizu 890ddd
    - Integrate by parts. In short, the resulting primitive of sqrt(x^2 + 1) is:
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
        F(x) = ( x * sqrt(x^2 + 1) + log(x + sqrt(x^2 + 1)) ) / 2;
Toshihiro Shimizu 890ddd
*/
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
void TQuadraticLengthEvaluator::setQuad(const TQuadratic &quad)
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	const TPointD &p0 = quad.getP0();
Toshihiro Shimizu 890ddd
	const TPointD &p1 = quad.getP1();
Toshihiro Shimizu 890ddd
	const TPointD &p2 = quad.getP2();
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	TPointD speed0(2.0 * (p1 - p0));
Toshihiro Shimizu 890ddd
	TPointD accel(2.0 * (p2 - p1) - speed0);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double a = accel * accel;
Toshihiro Shimizu 890ddd
	double b = 2.0 * accel * speed0;
Toshihiro Shimizu 890ddd
	m_c = speed0 * speed0;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_constantSpeed = isAlmostZero(a); // => b isAlmostZero, too
Toshihiro Shimizu 890ddd
	if (m_constantSpeed) {
Toshihiro Shimizu 890ddd
		m_c = sqrt(m_c);
Toshihiro Shimizu 890ddd
		return;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_sqrt_a_div_2 = 0.5 * sqrt(a);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_noSpeed0 = isAlmostZero(m_c); // => b isAlmostZero, too
Toshihiro Shimizu 890ddd
	if (m_noSpeed0)
Toshihiro Shimizu 890ddd
		return;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_tRef = 0.5 * b / a;
Toshihiro Shimizu 890ddd
	double d = m_c - 0.5 * b * m_tRef;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_squareIntegrand = (d < TConsts::epsilon);
Toshihiro Shimizu 890ddd
	if (m_squareIntegrand) {
Toshihiro Shimizu 890ddd
		m_f = (b > 0) ? -sq(m_tRef) : sq(m_tRef);
Toshihiro Shimizu 890ddd
		return;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_e = d / a;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double sqrt_part = sqrt(sq(m_tRef) + m_e);
Toshihiro Shimizu 890ddd
	double log_arg = m_tRef + sqrt_part;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_squareIntegrand = (log_arg < TConsts::epsilon);
Toshihiro Shimizu 890ddd
	if (m_squareIntegrand) {
Toshihiro Shimizu 890ddd
		m_f = (b > 0) ? -sq(m_tRef) : sq(m_tRef);
Toshihiro Shimizu 890ddd
		return;
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	m_primitive_0 = m_sqrt_a_div_2 * (m_tRef * sqrt_part + m_e * log(log_arg));
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//-----------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
double TQuadraticLengthEvaluator::getLengthAt(double t) const
Toshihiro Shimizu 890ddd
{
Toshihiro Shimizu 890ddd
	if (m_constantSpeed)
Toshihiro Shimizu 890ddd
		return m_c * t;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (m_noSpeed0)
Toshihiro Shimizu 890ddd
		return m_sqrt_a_div_2 * sq(t);
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	if (m_squareIntegrand) {
Toshihiro Shimizu 890ddd
		double t_plus_tRef = t + m_tRef;
Toshihiro Shimizu 890ddd
		return m_sqrt_a_div_2 * (m_f + ((t_plus_tRef > 0) ? sq(t_plus_tRef) : -sq(t_plus_tRef)));
Toshihiro Shimizu 890ddd
	}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	double y = t + m_tRef;
Toshihiro Shimizu 890ddd
	double sqrt_part = sqrt(sq(y) + m_e);
Toshihiro Shimizu 890ddd
	double log_arg = y + sqrt_part; //NOTE: log_arg >= log_arg0 >= TConsts::epsilon
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
	return m_sqrt_a_div_2 * (y * sqrt_part + m_e * log(log_arg)) - m_primitive_0;
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//-----------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
//  End Of File
Toshihiro Shimizu 890ddd
//-----------------------------------------------------------------------------