Blame c++/vector/curve.h

145120
#ifndef CURVE_H
145120
#define CURVE_H
145120
145120
145120
#include <cassert></cassert>
145120
145120
#include "vector.h"
145120
145120
145120
class Bezier;
145120
145120
class Hermite {
145120
public:
145120
	Vector p0, p1, t0, t1;
145120
145120
	explicit Hermite(
145120
		const Vector &p0 = Vector(),
145120
		const Vector &p1 = Vector(),
145120
		const Vector &t0 = Vector(),
145120
		const Vector &t1 = Vector()
145120
	):
145120
		p0(p0), p1(p1), t0(t0), t1(t1) { }
145120
145120
	bool operator< (const Hermite &a) const {
145120
		return p0 < a.p0 ? true : a.p0 < p0 ? false
145120
			 : p1 < a.p1 ? true : a.p1 < p1 ? false
145120
			 : t0 < a.t0 ? true : a.t0 < t0 ? false
145120
			 : t1 < a.t1;
145120
	}
145120
	bool operator>(const Hermite &a) const { return a < *this; }
145120
	bool operator==(const Hermite &a) const { return p0 == a.p0 && p1 == a.p1 && t0 == a.t0 && t1 == a.t1; }
145120
	bool operator!=(const Hermite &a) const { return !(*this == a); }
145120
145120
	Vector p(const Real &l) const
145120
		{ return p(l, p0, p1, t0, t1); }
145120
	Vector t(const Real &l) const
145120
		{ return t(l, p0, p1, t0, t1); }
145120
145120
	Vector operator() (const Real &l) const  { return p(l); }
145120
	
a7e4c0
	Vector d(Real l) const {
a7e4c0
		Vector t = this->t(l);
a7e4c0
		Real lensqr = t.lensqr();
a7e4c0
		if (precision*precision < lensqr) return t/sqrt(lensqr);
a7e4c0
		return (p(l + differential) - p(l - differential)).norm();
a7e4c0
	}
145120
	
145120
	Vector pp0() const { return p0 + t0/3; }
145120
	Vector pp1() const { return p1 - t1/3; }
145120
	
145120
	Hermite flip() const
145120
		{ return Hermite(p1, p0, -t1, -t0); }
145120
	
145120
	Hermite sub(const Real &a, const Real &b) const {
145120
		Real k = b - a;
145120
		return Hermite(p(a), p(b), t(a)*k, t(b)*k);
145120
	}
145120
	
145120
	void split(const Real &l, Hermite *a, Hermite *b) const {
145120
		if (!a && !b) return;
145120
		assert(a != b);
145120
		Vector point = p(l);
145120
		Vector tangent = t(l);
145120
		Hermite h = *this;
145120
		if (a) {
145120
			a->p0 = h.p0;
145120
			a->t0 = h.t0*l;
145120
			a->p1 = point;
145120
			a->t1 = tangent*l;
145120
		}
145120
		if (b && b != a) {
145120
			Real ll = 1 - l;
145120
			b->p0 = point;
145120
			b->t0 = tangent*ll;
145120
			b->p1 = h.p1;
145120
			b->t1 = h.t1*ll;
145120
		}
145120
	}
145120
145120
	Rect bounds() const
145120
		{ return Rect(p0).expand(p1).expand(pp0()).expand(pp1()); }
145120
145120
	int bends_x(Real *l) const
145120
		{ return bends(l, p0.x, p1.x, t0.x, t1.x); }
145120
	int bends_y(Real *l) const
145120
		{ return bends(l, p0.y, p1.y, t0.y, t1.y); }
145120
145120
	Range bounds_accurate_x() const
145120
		{ return bounds_accurate(p0.x, p1.x, t0.x, t1.x); }
145120
	Range bounds_accurate_y() const
145120
		{ return bounds_accurate(p0.y, p1.y, t0.y, t1.y); }
145120
	Rect bounds_accurate() const
145120
		{ return Rect(bounds_accurate_x(), bounds_accurate_y()); }
145120
a7e4c0
	int inflection_x(Real *l) const
a7e4c0
		{ return inflection(l, p0.x, p1.x, t0.x, t1.x); }
a7e4c0
	int inflection_y(Real *l) const
a7e4c0
		{ return inflection(l, p0.y, p1.y, t0.y, t1.y); }
a7e4c0
145120
	int intersections_x(Real *l, const Real &x) const
145120
		{ return intersections(l, x, p0.x, p1.x, t0.x, t1.x); }
145120
	int intersections_y(Real *l, const Real &y) const
145120
		{ return intersections(l, y, p0.y, p1.y, t0.y, t1.y); }
145120
	
145120
	inline Bezier get_bezier() const;
145120
	
145120
	static Hermite linear(const Vector &p, const Vector &t)
145120
		{ return Hermite(p, p + t, t, t);  }
145120
	static Hermite point(const Vector &p)
145120
		{ return Hermite(p, p);  }
145120
145120
	static Real p(const Real &l, const Real &p0, const Real &p1, const Real &t0, const Real &t1) {
145120
		return p0*((( 2*l - 3)*l    )*l + 1)
145120
			 + p1*(((-2*l + 3)*l    )*l    )
145120
			 + t0*(((   l - 2)*l + 1)*l    )
145120
			 + t1*(((   l - 1)*l    )*l    );
145120
	}
145120
145120
	static Real t(const Real &l, const Real &p0, const Real &p1, const Real &t0, const Real &t1) {
145120
		return p0*(( 6*l - 6)*l    )
145120
			 + p1*((-6*l + 6)*l    )
145120
			 + t0*(( 3*l - 4)*l + 1)
145120
			 + t1*(( 3*l - 2)*l    );
145120
	}
145120
145120
	static Vector p(const Real &l, const Vector &p0, const Vector &p1, const Vector &t0, const Vector &t1)
145120
		{ return Vector( p(l, p0.x, p1.x, t0.x, t1.x), p(l, p0.y, p1.y, t0.y, t1.y) ); }
145120
	static Vector t(const Real &l, const Vector &p0, const Vector &p1, const Vector &t0, const Vector &t1)
145120
		{ return Vector( t(l, p0.x, p1.x, t0.x, t1.x), t(l, p0.y, p1.y, t0.y, t1.y) ); }
145120
	
a7e4c0
	static int inflection(Real *l, const Real &p0, const Real &p1, const Real &t0, const Real &t1);
145120
	static int bends(Real *l, const Real &p0, const Real &p1, const Real &t0, const Real &t1);
145120
	static int intersections(Real *l, const Real &p, const Real &p0, const Real &p1, const Real &t0, const Real &t1);
145120
	static Range bounds_accurate(const Real &p0, const Real &p1, const Real &t0, const Real &t1);
145120
};
145120
145120
145120
class Bezier {
145120
public:
145120
	Vector p0, pp0, pp1, p1;
145120
145120
	Bezier(const Vector &p0, const Vector &p1, const Vector &pp0, const Vector &pp1):
145120
		p0(p0), pp0(pp0), pp1(pp1), p1(p1) { }
145120
	Bezier(const Vector &p0, const Vector &p1, const Vector &pp):
145120
		Bezier(p0, p1, pp, pp) { }
145120
	Bezier(const Vector &p0, const Vector &p1):
145120
		Bezier(p0, p1, (p1 - p0)/3 + p0, (p0 - p1)/3 + p1) { }
145120
	explicit Bezier(const Vector &p = Vector()):
145120
		Bezier(p, p, p, p) { }
145120
145120
	bool operator< (const Bezier &a) const {
145120
		return p0  < a.p0  ? true : a.p0  < p0  ? false
145120
			 : p1  < a.p1  ? true : a.p1  < p1  ? false
145120
			 : pp0 < a.pp0 ? true : a.pp0 < pp0 ? false
145120
			 : pp1 < a.pp1;
145120
	}
145120
	bool operator>(const Bezier &a) const { return a < *this; }
145120
	bool operator==(const Bezier &a) const { return p0 == a.p0 && p1 == a.p1 && pp0 == a.pp0 && pp1 == a.pp1; }
145120
	bool operator!=(const Bezier &a) const { return !(*this == a); }
145120
145120
	Vector p(const Real &l) const
145120
		{ return p(l, p0, p1, pp0, pp1); }
145120
	Vector t(const Real &l) const
145120
		{ return t(l, p0, p1, pp0, pp1); }
145120
145120
	Vector operator() (const Real &l) const  { return p(l); }
145120
	
a7e4c0
	Vector d(Real l) const {
a7e4c0
		Vector t = this->t(l);
a7e4c0
		Real lensqr = t.lensqr();
a7e4c0
		if (precision*precision < lensqr) return t/sqrt(lensqr);
a7e4c0
		return (p(l + differential) - p(l - differential)).norm();
a7e4c0
	}
145120
	
145120
	Vector t0() const { return 3*(pp0 - p0); }
145120
	Vector t1() const { return 3*(p1 - pp1); }
145120
	
145120
	void split(const Real &l, Bezier *a, Bezier *b) const {
145120
		if (!a && !b) return;
145120
		assert(a != b);
145120
		Real ll = 1 - l;
145120
		Vector p0 = this->p0;
145120
		Vector p1 = this->p1;
145120
		Vector a0 =  p0*ll + pp0*l;
145120
		Vector a1 = pp0*ll + pp1*l;
145120
		Vector a2 = pp1*ll +  p1*l;
145120
		Vector b0 =  a0*ll +  a1*l;
145120
		Vector b1 =  a1*ll +  a2*l;
145120
		Vector c  =  b0*ll +  b1*l;
145120
		if (a) {
a7e4c0
			a->p0  = p0;
145120
			a->pp0 = a0;
145120
			a->pp1 = b0;
145120
			a->p1  = c;
145120
		}
145120
		if (b && b != a) {
145120
			b->p0  = c;
145120
			b->pp0 = b1;
145120
			b->pp1 = a2;
a7e4c0
			b->p1  = p1;
145120
		}
145120
	}
145120
	
145120
	Bezier flip() const
145120
		{ return Bezier(p1, p0, pp1, pp0); }
145120
145120
	Bezier sub(const Real &a, const Real &b) const {
145120
		if (equal(a, b)) return Bezier( p(0.5*(a + b)) );
145120
145120
		Bezier s;
145120
		if (equal(a, 1)) {
145120
			split(b, nullptr, &s);
145120
			return s.flip();
145120
		}
145120
		
145120
		split(a, nullptr, &s);
145120
		s.split((b - a)/(1 - a), &s, nullptr);
145120
		return s;
145120
	}
145120
	
145120
	Rect bounds() const
145120
		{ return Rect(p0).expand(p1).expand(pp0).expand(pp1); }
145120
145120
	int bends_x(Real *l) const
145120
		{ return bends(l, p0.x, p1.x, pp0.x, pp1.x); }
145120
	int bends_y(Real *l) const
145120
		{ return bends(l, p0.y, p1.y, pp0.y, pp1.y); }
145120
145120
	Range bounds_accurate_x() const
145120
		{ return bounds_accurate(p0.x, p1.x, pp0.x, pp1.x); }
145120
	Range bounds_accurate_y() const
145120
		{ return bounds_accurate(p0.y, p1.y, pp0.y, pp1.y); }
145120
	Rect bounds_accurate() const
145120
		{ return Rect(bounds_accurate_x(), bounds_accurate_y()); }
145120
a7e4c0
	int inflection_x(Real *l) const
a7e4c0
		{ return inflection(l, p0.x, p1.x, pp0.x, pp1.x); }
a7e4c0
	int inflection_y(Real *l) const
a7e4c0
		{ return inflection(l, p0.y, p1.y, pp0.y, pp1.y); }
a7e4c0
145120
	int intersections_x(Real *l, const Real &x) const
145120
		{ return intersections(l, x, p0.x, p1.x, pp0.x, pp1.x); }
145120
	int intersections_y(Real *l, const Real &y) const
145120
		{ return intersections(l, y, p0.y, p1.y, pp0.y, pp1.y); }
145120
145120
	inline Hermite get_hermite() const;
145120
145120
	static Bezier linear(const Vector &p0, const Vector &p1)
145120
		{ return Bezier(p0, p1, p0 + (p1 - p0)/3, p1 + (p0 - p1)/3);  }
145120
	static Bezier point(const Vector &p)
145120
		{ return Bezier(p, p);  }
145120
145120
	static Real p(const Real &l, const Real &p0, const Real &p1, const Real &pp0, const Real &pp1) {
145120
		Real ll = 1 - l;
145120
		const Real a0 =  p0*ll + pp0*l;
145120
		const Real a1 = pp0*ll + pp1*l;
145120
		const Real a2 = pp1*ll +  p1*l;
145120
		const Real b0 =  a0*ll +  a1*l;
145120
		const Real b1 =  a1*ll +  a2*l;
145120
		return b0*ll +  b1*l;
145120
	}
145120
145120
	static Real t(const Real &l, const Real &p0, const Real &p1, const Real &pp0, const Real &pp1) {
145120
		Real ll = 1 - l;
145120
		const Real a0 =  p0*ll + pp0*l;
145120
		const Real a1 = pp0*ll + pp1*l;
145120
		const Real a2 = pp1*ll +  p1*l;
145120
		const Real b0 =  a0*ll +  a1*l;
145120
		const Real b1 =  a1*ll +  a2*l;
145120
		return 3*(b1 - b0)*l;
145120
	}
145120
145120
	static Vector p(const Real &l, const Vector &p0, const Vector &p1, const Vector &pp0, const Vector &pp1)
145120
		{ return Vector( p(l, p0.x, p1.x, pp0.x, pp1.x), p(l, p0.y, p1.y, pp0.y, pp1.y) ); }
145120
	static Vector t(const Real &l, const Vector &p0, const Vector &p1, const Vector &pp0, const Vector &pp1)
145120
		{ return Vector( t(l, p0.x, p1.x, pp0.x, pp1.x), t(l, p0.y, p1.y, pp0.y, pp1.y) ); }
145120
	
a7e4c0
	static int inflection(Real *l, const Real &p0, const Real &p1, const Real &pp0, const Real &pp1);
145120
	static int bends(Real *l, const Real &p0, const Real &p1, const Real &pp0, const Real &pp1);
145120
	static int intersections(Real *l, const Real &p, const Real &p0, const Real &p1, const Real &pp0, const Real &pp1);
145120
	static Range bounds_accurate(const Real &p0, const Real &p1, const Real &pp0, const Real &pp1);
145120
};
145120
145120
145120
inline Bezier Hermite::get_bezier() const { return Bezier(p0, p1, pp0(), pp1()); }
145120
inline Hermite Bezier::get_hermite() const { return Hermite(p0, p1, t0(), t1()); }
145120
145120
#endif