From 037c4c3708e40a1f7569b4ea176380958fd1389a Mon Sep 17 00:00:00 2001 From: johanengelen Date: Sat, 14 Jun 2008 15:01:19 +0000 Subject: [PATCH] update 2geom --- src/2geom/convex-cover.cpp | 6 +- src/2geom/curve.h | 6 + src/2geom/elliptical-arc.h | 5 +- src/2geom/hvlinesegment.h | 10 + src/2geom/nearest-point.cpp | 556 +++++++++++++++--------------- src/2geom/nearest-point.h | 266 +++++++------- src/2geom/path.cpp | 57 ++- src/2geom/path.h | 16 + src/2geom/pathvector.cpp | 53 +++ src/2geom/pathvector.h | 2 + src/2geom/piecewise.cpp | 16 + src/2geom/piecewise.h | 1 + src/2geom/poly-laguerre-solve.cpp | 9 +- src/2geom/poly.cpp | 6 +- src/2geom/quadtree.cpp | 14 + src/2geom/quadtree.h | 24 ++ src/2geom/sbasis-geometric.cpp | 11 + src/2geom/sbasis-geometric.h | 8 + src/2geom/sbasis.h | 1 + src/2geom/solve-bezier-one-d.cpp | 6 + 20 files changed, 652 insertions(+), 421 deletions(-) diff --git a/src/2geom/convex-cover.cpp b/src/2geom/convex-cover.cpp index 2d11b7419..50c43e6d1 100644 --- a/src/2geom/convex-cover.cpp +++ b/src/2geom/convex-cover.cpp @@ -114,7 +114,7 @@ ConvexHull::angle_sort() { void ConvexHull::graham_scan() { unsigned stac = 2; - for(int i = 2; i < boundary.size(); i++) { + for(unsigned int i = 2; i < boundary.size(); i++) { double o = SignedTriangleArea(boundary[stac-2], boundary[stac-1], boundary[i]); @@ -350,6 +350,7 @@ ConvexHull sweepline_intersection(ConvexHull const &a, ConvexHull const &b) { // al and bl now point to the top of the first pair of edges that overlap in y value double sweep_y = std::min(a.boundary[al][Y], b.boundary[bl][Y]); + return ret; } /*** ConvexHull intersection(ConvexHull a, ConvexHull b); @@ -359,10 +360,11 @@ ConvexHull sweepline_intersection(ConvexHull const &a, ConvexHull const &b) { */ ConvexHull intersection(ConvexHull a, ConvexHull b) { ConvexHull ret; + /* int ai = 0, bi = 0; int aj = a.boundary.size() - 1; int bj = b.boundary.size() - 1; - + */ /*while (true) { if(a[ai] }*/ diff --git a/src/2geom/curve.h b/src/2geom/curve.h index b76c0d328..d201d5d5e 100644 --- a/src/2geom/curve.h +++ b/src/2geom/curve.h @@ -105,10 +105,16 @@ public: virtual Point pointAt(Coord t) const { return pointAndDerivatives(t, 1).front(); } virtual Coord valueAt(Coord t, Dim2 d) const { return pointAt(t)[d]; } + virtual Point operator() (double t) const { return pointAt(t); } virtual std::vector pointAndDerivatives(Coord t, unsigned n) const = 0; virtual D2 toSBasis() const = 0; }; +inline +Coord nearest_point(Point const& p, Curve const& c) +{ + return c.nearestPoint(p); +} } // end namespace Geom diff --git a/src/2geom/elliptical-arc.h b/src/2geom/elliptical-arc.h index f3d4f4931..3a3b90670 100644 --- a/src/2geom/elliptical-arc.h +++ b/src/2geom/elliptical-arc.h @@ -166,7 +166,7 @@ class EllipticalArc : public Curve return ( are_near(ray(X), 0) || are_near(ray(Y), 0) ); } - // TODO: native implementation of the following methods + Rect boundsFast() const { return boundsExact(); @@ -174,6 +174,7 @@ class EllipticalArc : public Curve Rect boundsExact() const; + // TODO: native implementation of the following methods Rect boundsLocal(Interval i, unsigned int deg) const { return SBasisCurve(toSBasis()).boundsLocal(i, deg); @@ -193,6 +194,7 @@ class EllipticalArc : public Curve return allNearestPoints(p, from, to).front(); } + // TODO: native implementation of the following methods int winding(Point p) const { return SBasisCurve(toSBasis()).winding(p); @@ -200,6 +202,7 @@ class EllipticalArc : public Curve Curve *derivative() const; + // TODO: native implementation of the following methods Curve *transformed(Matrix const &m) const { return SBasisCurve(toSBasis()).transformed(m); diff --git a/src/2geom/hvlinesegment.h b/src/2geom/hvlinesegment.h index 8c479858b..b3c200378 100644 --- a/src/2geom/hvlinesegment.h +++ b/src/2geom/hvlinesegment.h @@ -170,6 +170,11 @@ class HLineSegment : public Curve if ( from > to ) std::swap(from, to); double xfrom = pointAt(from)[X]; double xto = pointAt(to)[X]; + if ( xfrom > xto ) + { + std::swap(xfrom, xto); + std::swap(from, to); + } if ( p[X] > xfrom && p[X] < xto ) { return (p[X] - initialPoint()[X]) / (finalPoint()[X] - initialPoint()[X]); @@ -395,6 +400,11 @@ class VLineSegment : public Curve if ( from > to ) std::swap(from, to); double yfrom = pointAt(from)[Y]; double yto = pointAt(to)[Y]; + if (yfrom > yto) + { + std::swap(yfrom, yto); + std::swap(from, to); + } if ( p[Y] > yfrom && p[Y] < yto ) { return (p[Y] - initialPoint()[Y]) / (finalPoint()[Y] - initialPoint()[Y]); diff --git a/src/2geom/nearest-point.cpp b/src/2geom/nearest-point.cpp index de880713f..074de1cb3 100644 --- a/src/2geom/nearest-point.cpp +++ b/src/2geom/nearest-point.cpp @@ -1,278 +1,278 @@ -/* - * nearest point routines for D2 and Piecewise> - * - * Authors: - * - * Marco Cecchetti - * - * Copyright 2007-2008 authors - * - * This library is free software; you can redistribute it and/or - * modify it either under the terms of the GNU Lesser General Public - * License version 2.1 as published by the Free Software Foundation - * (the "LGPL") or, at your option, under the terms of the Mozilla - * Public License Version 1.1 (the "MPL"). If you do not alter this - * notice, a recipient may use your version of this file under either - * the MPL or the LGPL. - * - * You should have received a copy of the LGPL along with this library - * in the file COPYING-LGPL-2.1; if not, write to the Free Software - * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA - * You should have received a copy of the MPL along with this library - * in the file COPYING-MPL-1.1 - * - * The contents of this file are subject to the Mozilla Public License - * Version 1.1 (the "License"); you may not use this file except in - * compliance with the License. You may obtain a copy of the License at - * http://www.mozilla.org/MPL/ - * - * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY - * OF ANY KIND, either express or implied. See the LGPL or the MPL for - * the specific language governing rights and limitations. - */ - - -#include "nearest-point.h" - -namespace Geom -{ - -//////////////////////////////////////////////////////////////////////////////// -// D2 versions - -/* - * Return the parameter t of a nearest point on the portion of the curve "c", - * related to the interval [from, to], to the point "p". - * The needed curve derivative "dc" is passed as parameter. - * The function return the first nearest point to "p" that is found. - */ - -double nearest_point( Point const& p, - D2 const& c, - D2 const& dc, - double from, double to ) -{ - if ( from > to ) std::swap(from, to); - if ( from < 0 || to > 1 ) - { - THROW_RANGEERROR("[from,to] interval out of bounds"); - } - - SBasis dd = dot(c - p, dc); - std::vector zeros = Geom::roots(dd); - - double closest = from; - double min_dist_sq = L2sq(c(from) - p); - double distsq; - for ( unsigned int i = 0; i < zeros.size(); ++i ) - { - distsq = L2sq(c(zeros[i]) - p); - if ( min_dist_sq > L2sq(c(zeros[i]) - p) ) - { - closest = zeros[i]; - min_dist_sq = distsq; - } - } - if ( min_dist_sq > L2sq( c(to) - p ) ) - closest = to; - return closest; - -} - -/* - * Return the parameters t of all the nearest points on the portion of - * the curve "c", related to the interval [from, to], to the point "p". - * The needed curve derivative "dc" is passed as parameter. - */ - -std::vector -all_nearest_points( Point const& p, - D2 const& c, - D2 const& dc, - double from, double to ) -{ - std::swap(from, to); - if ( from > to ) std::swap(from, to); - if ( from < 0 || to > 1 ) - { - THROW_RANGEERROR("[from,to] interval out of bounds"); - } - - std::vector result; - SBasis dd = dot(c - p, Geom::derivative(c)); - - std::vector zeros = Geom::roots(dd); - std::vector candidates; - candidates.push_back(from); - candidates.insert(candidates.end(), zeros.begin(), zeros.end()); - candidates.push_back(to); - std::vector distsq; - distsq.reserve(candidates.size()); - for ( unsigned int i = 0; i < candidates.size(); ++i ) - { - distsq.push_back( L2sq(c(candidates[i]) - p) ); - } - unsigned int closest = 0; - double dsq = distsq[0]; - for ( unsigned int i = 1; i < candidates.size(); ++i ) - { - if ( dsq > distsq[i] ) - { - closest = i; - dsq = distsq[i]; - } - } - for ( unsigned int i = 0; i < candidates.size(); ++i ) - { - if( distsq[closest] == distsq[i] ) - { - result.push_back(candidates[i]); - } - } - return result; -} - - -//////////////////////////////////////////////////////////////////////////////// -// Piecewise< D2 > versions - - -double nearest_point( Point const& p, - Piecewise< D2 > const& c, - double from, double to ) -{ - if ( from > to ) std::swap(from, to); - if ( from < c.cuts[0] || to > c.cuts[c.size()] ) - { - THROW_RANGEERROR("[from,to] interval out of bounds"); - } - - unsigned int si = c.segN(from); - unsigned int ei = c.segN(to); - if ( si == ei ) - { - double nearest= - nearest_point(p, c[si], c.segT(from, si), c.segT(to, si)); - return c.mapToDomain(nearest, si); - } - double t; - double nearest = nearest_point(p, c[si], c.segT(from, si)); - unsigned int ni = si; - double dsq; - double mindistsq = distanceSq(p, c[si](nearest)); - Rect bb; - for ( unsigned int i = si + 1; i < ei; ++i ) - { - bb = bounds_fast(c[i]); - dsq = distanceSq(p, bb); - if ( mindistsq <= dsq ) continue; - t = nearest_point(p, c[i]); - dsq = distanceSq(p, c[i](t)); - if ( mindistsq > dsq ) - { - nearest = t; - ni = i; - mindistsq = dsq; - } - } - bb = bounds_fast(c[ei]); - dsq = distanceSq(p, bb); - if ( mindistsq > dsq ) - { - t = nearest_point(p, c[ei], 0, c.segT(to, ei)); - dsq = distanceSq(p, c[ei](t)); - if ( mindistsq > dsq ) - { - nearest = t; - ni = ei; - } - } - return c.mapToDomain(nearest, ni); -} - -std::vector -all_nearest_points( Point const& p, - Piecewise< D2 > const& c, - double from, double to ) -{ - if ( from > to ) std::swap(from, to); - if ( from < c.cuts[0] || to > c.cuts[c.size()] ) - { - THROW_RANGEERROR("[from,to] interval out of bounds"); - } - - unsigned int si = c.segN(from); - unsigned int ei = c.segN(to); - if ( si == ei ) - { - std::vector all_nearest = - all_nearest_points(p, c[si], c.segT(from, si), c.segT(to, si)); - for ( unsigned int i = 0; i < all_nearest.size(); ++i ) - { - all_nearest[i] = c.mapToDomain(all_nearest[i], si); - } - return all_nearest; - } - std::vector all_t; - std::vector< std::vector > all_np; - all_np.push_back( all_nearest_points(p, c[si], c.segT(from, si)) ); - std::vector ni; - ni.push_back(si); - double dsq; - double mindistsq = distanceSq( p, c[si](all_np.front().front()) ); - Rect bb; - for ( unsigned int i = si + 1; i < ei; ++i ) - { - bb = bounds_fast(c[i]); - dsq = distanceSq(p, bb); - if ( mindistsq < dsq ) continue; - all_t = all_nearest_points(p, c[i]); - dsq = distanceSq( p, c[i](all_t.front()) ); - if ( mindistsq > dsq ) - { - all_np.clear(); - all_np.push_back(all_t); - ni.clear(); - ni.push_back(i); - mindistsq = dsq; - } - else if ( mindistsq == dsq ) - { - all_np.push_back(all_t); - ni.push_back(i); - } - } - bb = bounds_fast(c[ei]); - dsq = distanceSq(p, bb); - if ( mindistsq >= dsq ) - { - all_t = all_nearest_points(p, c[ei], 0, c.segT(to, ei)); - dsq = distanceSq( p, c[ei](all_t.front()) ); - if ( mindistsq > dsq ) - { - for ( unsigned int i = 0; i < all_t.size(); ++i ) - { - all_t[i] = c.mapToDomain(all_t[i], ei); - } - return all_t; - } - else if ( mindistsq == dsq ) - { - all_np.push_back(all_t); - ni.push_back(ei); - } - } - std::vector all_nearest; - for ( unsigned int i = 0; i < all_np.size(); ++i ) - { - for ( unsigned int j = 0; j < all_np[i].size(); ++j ) - { - all_nearest.push_back( c.mapToDomain(all_np[i][j], ni[i]) ); - } - } - return all_nearest; -} - -} // end namespace Geom - - +/* + * nearest point routines for D2 and Piecewise> + * + * Authors: + * + * Marco Cecchetti + * + * Copyright 2007-2008 authors + * + * This library is free software; you can redistribute it and/or + * modify it either under the terms of the GNU Lesser General Public + * License version 2.1 as published by the Free Software Foundation + * (the "LGPL") or, at your option, under the terms of the Mozilla + * Public License Version 1.1 (the "MPL"). If you do not alter this + * notice, a recipient may use your version of this file under either + * the MPL or the LGPL. + * + * You should have received a copy of the LGPL along with this library + * in the file COPYING-LGPL-2.1; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + * You should have received a copy of the MPL along with this library + * in the file COPYING-MPL-1.1 + * + * The contents of this file are subject to the Mozilla Public License + * Version 1.1 (the "License"); you may not use this file except in + * compliance with the License. You may obtain a copy of the License at + * http://www.mozilla.org/MPL/ + * + * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY + * OF ANY KIND, either express or implied. See the LGPL or the MPL for + * the specific language governing rights and limitations. + */ + + +#include "nearest-point.h" + +namespace Geom +{ + +//////////////////////////////////////////////////////////////////////////////// +// D2 versions + +/* + * Return the parameter t of a nearest point on the portion of the curve "c", + * related to the interval [from, to], to the point "p". + * The needed curve derivative "dc" is passed as parameter. + * The function return the first nearest point to "p" that is found. + */ + +double nearest_point( Point const& p, + D2 const& c, + D2 const& dc, + double from, double to ) +{ + if ( from > to ) std::swap(from, to); + if ( from < 0 || to > 1 ) + { + THROW_RANGEERROR("[from,to] interval out of bounds"); + } + + SBasis dd = dot(c - p, dc); + std::vector zeros = Geom::roots(dd); + + double closest = from; + double min_dist_sq = L2sq(c(from) - p); + double distsq; + for ( unsigned int i = 0; i < zeros.size(); ++i ) + { + distsq = L2sq(c(zeros[i]) - p); + if ( min_dist_sq > L2sq(c(zeros[i]) - p) ) + { + closest = zeros[i]; + min_dist_sq = distsq; + } + } + if ( min_dist_sq > L2sq( c(to) - p ) ) + closest = to; + return closest; + +} + +/* + * Return the parameters t of all the nearest points on the portion of + * the curve "c", related to the interval [from, to], to the point "p". + * The needed curve derivative "dc" is passed as parameter. + */ + +std::vector +all_nearest_points( Point const& p, + D2 const& c, + D2 const& dc, + double from, double to ) +{ + std::swap(from, to); + if ( from > to ) std::swap(from, to); + if ( from < 0 || to > 1 ) + { + THROW_RANGEERROR("[from,to] interval out of bounds"); + } + + std::vector result; + SBasis dd = dot(c - p, Geom::derivative(c)); + + std::vector zeros = Geom::roots(dd); + std::vector candidates; + candidates.push_back(from); + candidates.insert(candidates.end(), zeros.begin(), zeros.end()); + candidates.push_back(to); + std::vector distsq; + distsq.reserve(candidates.size()); + for ( unsigned int i = 0; i < candidates.size(); ++i ) + { + distsq.push_back( L2sq(c(candidates[i]) - p) ); + } + unsigned int closest = 0; + double dsq = distsq[0]; + for ( unsigned int i = 1; i < candidates.size(); ++i ) + { + if ( dsq > distsq[i] ) + { + closest = i; + dsq = distsq[i]; + } + } + for ( unsigned int i = 0; i < candidates.size(); ++i ) + { + if( distsq[closest] == distsq[i] ) + { + result.push_back(candidates[i]); + } + } + return result; +} + + +//////////////////////////////////////////////////////////////////////////////// +// Piecewise< D2 > versions + + +double nearest_point( Point const& p, + Piecewise< D2 > const& c, + double from, double to ) +{ + if ( from > to ) std::swap(from, to); + if ( from < c.cuts[0] || to > c.cuts[c.size()] ) + { + THROW_RANGEERROR("[from,to] interval out of bounds"); + } + + unsigned int si = c.segN(from); + unsigned int ei = c.segN(to); + if ( si == ei ) + { + double nearest= + nearest_point(p, c[si], c.segT(from, si), c.segT(to, si)); + return c.mapToDomain(nearest, si); + } + double t; + double nearest = nearest_point(p, c[si], c.segT(from, si)); + unsigned int ni = si; + double dsq; + double mindistsq = distanceSq(p, c[si](nearest)); + Rect bb; + for ( unsigned int i = si + 1; i < ei; ++i ) + { + bb = bounds_fast(c[i]); + dsq = distanceSq(p, bb); + if ( mindistsq <= dsq ) continue; + t = nearest_point(p, c[i]); + dsq = distanceSq(p, c[i](t)); + if ( mindistsq > dsq ) + { + nearest = t; + ni = i; + mindistsq = dsq; + } + } + bb = bounds_fast(c[ei]); + dsq = distanceSq(p, bb); + if ( mindistsq > dsq ) + { + t = nearest_point(p, c[ei], 0, c.segT(to, ei)); + dsq = distanceSq(p, c[ei](t)); + if ( mindistsq > dsq ) + { + nearest = t; + ni = ei; + } + } + return c.mapToDomain(nearest, ni); +} + +std::vector +all_nearest_points( Point const& p, + Piecewise< D2 > const& c, + double from, double to ) +{ + if ( from > to ) std::swap(from, to); + if ( from < c.cuts[0] || to > c.cuts[c.size()] ) + { + THROW_RANGEERROR("[from,to] interval out of bounds"); + } + + unsigned int si = c.segN(from); + unsigned int ei = c.segN(to); + if ( si == ei ) + { + std::vector all_nearest = + all_nearest_points(p, c[si], c.segT(from, si), c.segT(to, si)); + for ( unsigned int i = 0; i < all_nearest.size(); ++i ) + { + all_nearest[i] = c.mapToDomain(all_nearest[i], si); + } + return all_nearest; + } + std::vector all_t; + std::vector< std::vector > all_np; + all_np.push_back( all_nearest_points(p, c[si], c.segT(from, si)) ); + std::vector ni; + ni.push_back(si); + double dsq; + double mindistsq = distanceSq( p, c[si](all_np.front().front()) ); + Rect bb; + for ( unsigned int i = si + 1; i < ei; ++i ) + { + bb = bounds_fast(c[i]); + dsq = distanceSq(p, bb); + if ( mindistsq < dsq ) continue; + all_t = all_nearest_points(p, c[i]); + dsq = distanceSq( p, c[i](all_t.front()) ); + if ( mindistsq > dsq ) + { + all_np.clear(); + all_np.push_back(all_t); + ni.clear(); + ni.push_back(i); + mindistsq = dsq; + } + else if ( mindistsq == dsq ) + { + all_np.push_back(all_t); + ni.push_back(i); + } + } + bb = bounds_fast(c[ei]); + dsq = distanceSq(p, bb); + if ( mindistsq >= dsq ) + { + all_t = all_nearest_points(p, c[ei], 0, c.segT(to, ei)); + dsq = distanceSq( p, c[ei](all_t.front()) ); + if ( mindistsq > dsq ) + { + for ( unsigned int i = 0; i < all_t.size(); ++i ) + { + all_t[i] = c.mapToDomain(all_t[i], ei); + } + return all_t; + } + else if ( mindistsq == dsq ) + { + all_np.push_back(all_t); + ni.push_back(ei); + } + } + std::vector all_nearest; + for ( unsigned int i = 0; i < all_np.size(); ++i ) + { + for ( unsigned int j = 0; j < all_np[i].size(); ++j ) + { + all_nearest.push_back( c.mapToDomain(all_np[i][j], ni[i]) ); + } + } + return all_nearest; +} + +} // end namespace Geom + + diff --git a/src/2geom/nearest-point.h b/src/2geom/nearest-point.h index 8c1a6f7d7..73ac0c3ce 100644 --- a/src/2geom/nearest-point.h +++ b/src/2geom/nearest-point.h @@ -1,133 +1,133 @@ -/* - * nearest point routines for D2 and Piecewise> - * - * - * Authors: - * - * Marco Cecchetti - * - * Copyright 2007-2008 authors - * - * This library is free software; you can redistribute it and/or - * modify it either under the terms of the GNU Lesser General Public - * License version 2.1 as published by the Free Software Foundation - * (the "LGPL") or, at your option, under the terms of the Mozilla - * Public License Version 1.1 (the "MPL"). If you do not alter this - * notice, a recipient may use your version of this file under either - * the MPL or the LGPL. - * - * You should have received a copy of the LGPL along with this library - * in the file COPYING-LGPL-2.1; if not, write to the Free Software - * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA - * You should have received a copy of the MPL along with this library - * in the file COPYING-MPL-1.1 - * - * The contents of this file are subject to the Mozilla Public License - * Version 1.1 (the "License"); you may not use this file except in - * compliance with the License. You may obtain a copy of the License at - * http://www.mozilla.org/MPL/ - * - * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY - * OF ANY KIND, either express or implied. See the LGPL or the MPL for - * the specific language governing rights and limitations. - */ - - -#ifndef _NEAREST_POINT_H_ -#define _NEAREST_POINT_H_ - - -#include - -#include "d2.h" -#include "piecewise.h" -#include "exception.h" - - - -namespace Geom -{ - -/* - * Given a line L specified by a point A and direction vector v, - * return the point on L nearest to p. Note that the returned value - * is with respect to the _normalized_ direction of v! - */ -inline double nearest_point(Point const &p, Point const &A, Point const &v) -{ - Point d(p - A); - return d[0] * v[0] + d[1] * v[1]; -} - -//////////////////////////////////////////////////////////////////////////////// -// D2 versions - -/* - * Return the parameter t of a nearest point on the portion of the curve "c", - * related to the interval [from, to], to the point "p". - * The needed curve derivative "dc" is passed as parameter. - * The function return the first nearest point to "p" that is found. - */ -double nearest_point( Point const& p, - D2 const& c, D2 const& dc, - double from = 0, double to = 1 ); - -inline -double nearest_point( Point const& p, - D2 const& c, - double from = 0, double to = 1 ) -{ - return nearest_point(p, c, Geom::derivative(c), from, to); -} - -/* - * Return the parameters t of all the nearest points on the portion of - * the curve "c", related to the interval [from, to], to the point "p". - * The needed curve derivative "dc" is passed as parameter. - */ -std::vector -all_nearest_points( Point const& p, - D2 const& c, D2 const& dc, - double from = 0, double to = 1 ); - -inline -std::vector -all_nearest_points( Point const& p, - D2 const& c, - double from = 0, double to = 1 ) -{ - return all_nearest_points(p, c, Geom::derivative(c), from, to); -} - - -//////////////////////////////////////////////////////////////////////////////// -// Piecewise< D2 > versions - -double nearest_point( Point const& p, - Piecewise< D2 > const& c, - double from, double to ); - -inline -double nearest_point( Point const& p, Piecewise< D2 > const& c ) -{ - return nearest_point(p, c, c.cuts[0], c.cuts[c.size()]); -} - - -std::vector -all_nearest_points( Point const& p, - Piecewise< D2 > const& c, - double from, double to ); - -inline -std::vector -all_nearest_points( Point const& p, Piecewise< D2 > const& c ) -{ - return all_nearest_points(p, c, c.cuts[0], c.cuts[c.size()]); -} - -} // end namespace Geom - - - -#endif /*_NEAREST_POINT_H_*/ +/* + * nearest point routines for D2 and Piecewise> + * + * + * Authors: + * + * Marco Cecchetti + * + * Copyright 2007-2008 authors + * + * This library is free software; you can redistribute it and/or + * modify it either under the terms of the GNU Lesser General Public + * License version 2.1 as published by the Free Software Foundation + * (the "LGPL") or, at your option, under the terms of the Mozilla + * Public License Version 1.1 (the "MPL"). If you do not alter this + * notice, a recipient may use your version of this file under either + * the MPL or the LGPL. + * + * You should have received a copy of the LGPL along with this library + * in the file COPYING-LGPL-2.1; if not, write to the Free Software + * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA + * You should have received a copy of the MPL along with this library + * in the file COPYING-MPL-1.1 + * + * The contents of this file are subject to the Mozilla Public License + * Version 1.1 (the "License"); you may not use this file except in + * compliance with the License. You may obtain a copy of the License at + * http://www.mozilla.org/MPL/ + * + * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY + * OF ANY KIND, either express or implied. See the LGPL or the MPL for + * the specific language governing rights and limitations. + */ + + +#ifndef _NEAREST_POINT_H_ +#define _NEAREST_POINT_H_ + + +#include + +#include "d2.h" +#include "piecewise.h" +#include "exception.h" + + + +namespace Geom +{ + +/* + * Given a line L specified by a point A and direction vector v, + * return the point on L nearest to p. Note that the returned value + * is with respect to the _normalized_ direction of v! + */ +inline double nearest_point(Point const &p, Point const &A, Point const &v) +{ + Point d(p - A); + return d[0] * v[0] + d[1] * v[1]; +} + +//////////////////////////////////////////////////////////////////////////////// +// D2 versions + +/* + * Return the parameter t of a nearest point on the portion of the curve "c", + * related to the interval [from, to], to the point "p". + * The needed curve derivative "dc" is passed as parameter. + * The function return the first nearest point to "p" that is found. + */ +double nearest_point( Point const& p, + D2 const& c, D2 const& dc, + double from = 0, double to = 1 ); + +inline +double nearest_point( Point const& p, + D2 const& c, + double from = 0, double to = 1 ) +{ + return nearest_point(p, c, Geom::derivative(c), from, to); +} + +/* + * Return the parameters t of all the nearest points on the portion of + * the curve "c", related to the interval [from, to], to the point "p". + * The needed curve derivative "dc" is passed as parameter. + */ +std::vector +all_nearest_points( Point const& p, + D2 const& c, D2 const& dc, + double from = 0, double to = 1 ); + +inline +std::vector +all_nearest_points( Point const& p, + D2 const& c, + double from = 0, double to = 1 ) +{ + return all_nearest_points(p, c, Geom::derivative(c), from, to); +} + + +//////////////////////////////////////////////////////////////////////////////// +// Piecewise< D2 > versions + +double nearest_point( Point const& p, + Piecewise< D2 > const& c, + double from, double to ); + +inline +double nearest_point( Point const& p, Piecewise< D2 > const& c ) +{ + return nearest_point(p, c, c.cuts[0], c.cuts[c.size()]); +} + + +std::vector +all_nearest_points( Point const& p, + Piecewise< D2 > const& c, + double from, double to ); + +inline +std::vector +all_nearest_points( Point const& p, Piecewise< D2 > const& c ) +{ + return all_nearest_points(p, c, c.cuts[0], c.cuts[c.size()]); +} + +} // end namespace Geom + + + +#endif /*_NEAREST_POINT_H_*/ diff --git a/src/2geom/path.cpp b/src/2geom/path.cpp index 35c7b76bc..0ff1f3b0c 100644 --- a/src/2geom/path.cpp +++ b/src/2geom/path.cpp @@ -50,8 +50,11 @@ void Path::swap(Path &other) { Rect Path::boundsFast() const { Rect bounds=front().boundsFast(); - for ( const_iterator iter=++begin(); iter != end() ; ++iter ) { - bounds.unionWith(iter->boundsFast()); + const_iterator iter = begin(); + if ( iter != end() ) { + for ( ++iter; iter != end() ; ++iter ) { + bounds.unionWith(iter->boundsFast()); + } } return bounds; } @@ -236,6 +239,56 @@ double Path::nearestPoint(Point const& _point, double from, double to) const return ni + nearest; } +Rect Path::boundsFast() +{ + Rect bound; + if (empty()) return bound; + + bound = begin()->boundsFast(); + double top = bound.top(); + double bottom = bound.bottom(); + double left = bound.left(); + double right = bound.right(); + for (iterator it = ++begin(); it != end(); ++it) + { + bound = it->boundsFast(); + if ( top > bound.top() ) top = bound.top(); + if ( bottom < bound.bottom() ) bottom = bound.bottom(); + if ( left > bound.left() ) left = bound.left(); + if ( right < bound.right() ) right = bound.right(); + } + bound[0][0] = left; + bound[0][1] = right; + bound[1][0] = top; + bound[1][1] = bottom; + return bound; +} + +Rect Path::boundsExact() +{ + Rect bound; + if (empty()) return bound; + + bound = begin()->boundsExact(); + double top = bound.top(); + double bottom = bound.bottom(); + double left = bound.left(); + double right = bound.right(); + for (iterator it = ++begin(); it != end(); ++it) + { + bound = it->boundsExact(); + if ( top > bound.top() ) top = bound.top(); + if ( bottom < bound.bottom() ) bottom = bound.bottom(); + if ( left > bound.left() ) left = bound.left(); + if ( right < bound.right() ) right = bound.right(); + } + bound[0][0] = left; + bound[0][1] = right; + bound[1][0] = top; + bound[1][1] = bottom; + return bound; +} + //This assumes that you can't be perfect in your t-vals, and as such, tweaks the start void Path::appendPortionTo(Path &ret, double from, double to) const { assert(from >= 0 && to >= 0); diff --git a/src/2geom/path.h b/src/2geom/path.h index 670d4c125..8f51c46c2 100644 --- a/src/2geom/path.h +++ b/src/2geom/path.h @@ -266,6 +266,12 @@ public: return (*this)[i].valueAt(lt, d); } + + Point operator() (double t) const + { + return pointAt(t); + } + std::vector roots(double v, Dim2 d) const { std::vector res; for(unsigned i = 0; i <= size(); i++) { @@ -297,6 +303,9 @@ public: return nearestPoint(_point, 0, sz); } + Rect boundsFast(); + Rect boundsExact(); + void appendPortionTo(Path &p, double f, double t) const; Path portion(double f, double t) const { @@ -543,6 +552,13 @@ inline static Piecewise > paths_to_pw(std::vector paths) { return ret; } +inline +Coord nearest_point(Point const& p, Path const& c) +{ + return c.nearestPoint(p); +} + + /* class PathPortion : public Curve { Path *source; diff --git a/src/2geom/pathvector.cpp b/src/2geom/pathvector.cpp index 696ccad77..9c8111767 100644 --- a/src/2geom/pathvector.cpp +++ b/src/2geom/pathvector.cpp @@ -65,6 +65,59 @@ PathVector reverse_paths_and_order (PathVector const & path_in) return path_out; } +Rect bounds_fast( PathVector const& pv ) +{ + typedef PathVector::const_iterator const_iterator; + + Rect bound; + if (pv.empty()) return bound; + + bound = (pv.begin())->boundsFast(); + double top = bound.top(); + double bottom = bound.bottom(); + double left = bound.left(); + double right = bound.right(); + for (const_iterator it = ++(pv.begin()); it != pv.end(); ++it) + { + bound = it->boundsFast(); + if ( top > bound.top() ) top = bound.top(); + if ( bottom < bound.bottom() ) bottom = bound.bottom(); + if ( left > bound.left() ) left = bound.left(); + if ( right < bound.right() ) right = bound.right(); + } + bound[0][0] = left; + bound[0][1] = right; + bound[1][0] = top; + bound[1][1] = bottom; + return bound; +} + +Rect bounds_exact( PathVector const& pv ) +{ + typedef PathVector::const_iterator const_iterator; + + Rect bound; + if (pv.empty()) return bound; + + bound = (pv.begin())->boundsExact(); + double top = bound.top(); + double bottom = bound.bottom(); + double left = bound.left(); + double right = bound.right(); + for (const_iterator it = ++(pv.begin()); it != pv.end(); ++it) + { + bound = it->boundsExact(); + if ( top > bound.top() ) top = bound.top(); + if ( bottom < bound.bottom() ) bottom = bound.bottom(); + if ( left > bound.left() ) left = bound.left(); + if ( right < bound.right() ) right = bound.right(); + } + bound[0][0] = left; + bound[0][1] = right; + bound[1][0] = top; + bound[1][1] = bottom; + return bound; +} } // namespace Geom diff --git a/src/2geom/pathvector.h b/src/2geom/pathvector.h index 97a0d19f1..a9f4d88ed 100644 --- a/src/2geom/pathvector.h +++ b/src/2geom/pathvector.h @@ -46,6 +46,8 @@ PathVector operator* (PathVector const & path_in, Matrix const &m); PathVector reverse_paths_and_order (PathVector const & path_in); +Rect bounds_fast( PathVector const & pv ); +Rect bounds_exact( PathVector const & pv ); } #endif // SEEN_GEOM_PATHVECTOR_H diff --git a/src/2geom/piecewise.cpp b/src/2geom/piecewise.cpp index dc91ab4a9..222152aa3 100644 --- a/src/2geom/piecewise.cpp +++ b/src/2geom/piecewise.cpp @@ -167,6 +167,22 @@ std::vector roots(Piecewise const &f){ return result; } +Piecewise reverse(Piecewise const &f) { + Piecewise ret = Piecewise(); + ret.cuts.resize(f.cuts.size()); + ret.segs.resize(f.segs.size()); + double start = f.cuts[0]; + double end = f.cuts.back(); + for (unsigned i = 0; i < f.cuts.size(); i++) { + double x = f.cuts[f.cuts.size() - 1 - i]; + ret.cuts[i] = end - (x - start); + } + for (unsigned i = 0; i < f.segs.size(); i++) + ret.segs[i] = reverse(f[f.segs.size() - i - 1]); + return ret; +} + + } /* Local Variables: diff --git a/src/2geom/piecewise.h b/src/2geom/piecewise.h index 574d6de62..3f13d15b0 100644 --- a/src/2geom/piecewise.h +++ b/src/2geom/piecewise.h @@ -704,6 +704,7 @@ Piecewise derivative(Piecewise const &a) { } std::vector roots(Piecewise const &f); +Piecewise reverse(Piecewise const &f); } diff --git a/src/2geom/poly-laguerre-solve.cpp b/src/2geom/poly-laguerre-solve.cpp index 766f16eda..266e24c2b 100644 --- a/src/2geom/poly-laguerre-solve.cpp +++ b/src/2geom/poly-laguerre-solve.cpp @@ -12,7 +12,7 @@ cdouble laguerre_internal_complex(Poly const & p, double n = p.degree(); quad_root = false; const unsigned shuffle_rate = 10; - static double shuffle[] = {0, 0.5, 0.25, 0.75, 0.125, 0.375, 0.625, 0.875, 1.0}; + //static double shuffle[] = {0, 0.5, 0.25, 0.75, 0.125, 0.375, 0.625, 0.875, 1.0}; unsigned shuffle_counter = 0; while(std::norm(a) > (tol*tol)) { //std::cout << "xk = " << xk << std::endl; @@ -132,7 +132,12 @@ laguerre(Poly p, const double tol) { std::vector laguerre_real_interval(Poly const & ply, const double lo, const double hi, - const double tol) { + const double tol) +{ + /* not implemented*/ + assert(false); + std::vector result; + return result; } /* diff --git a/src/2geom/poly.cpp b/src/2geom/poly.cpp index 4f5ba6c55..27f98596b 100644 --- a/src/2geom/poly.cpp +++ b/src/2geom/poly.cpp @@ -45,7 +45,7 @@ std::vector > solve(Poly const & pp) { gsl_complex_packed_ptr z = new double[p.degree()*2]; double* a = new double[p.size()]; - for(int i = 0; i < p.size(); i++) + for(unsigned int i = 0; i < p.size(); i++) a[i] = p[i]; std::vector > roots; //roots.resize(p.degree()); @@ -55,7 +55,7 @@ std::vector > solve(Poly const & pp) { gsl_poly_complex_workspace_free (w); - for (int i = 0; i < p.degree(); i++) { + for (unsigned int i = 0; i < p.degree(); i++) { roots.push_back(std::complex (z[2*i] ,z[2*i+1])); //printf ("z%d = %+.18f %+.18f\n", i, z[2*i], z[2*i+1]); } @@ -67,7 +67,7 @@ std::vector solve_reals(Poly const & p) { std::vector > roots = solve(p); std::vector real_roots; - for(int i = 0; i < roots.size(); i++) { + for(unsigned int i = 0; i < roots.size(); i++) { if(roots[i].imag() == 0) // should be more lenient perhaps real_roots.push_back(roots[i].real()); } diff --git a/src/2geom/quadtree.cpp b/src/2geom/quadtree.cpp index a84a5a7d4..bb041edbe 100644 --- a/src/2geom/quadtree.cpp +++ b/src/2geom/quadtree.cpp @@ -1,5 +1,17 @@ #include "quadtree.h" +namespace Geom{ +Quad* QuadTree::search(Rect const &r) { + return search(r[0].min(), r[1].min(), + r[0].max(), r[1].max()); +} + +void QuadTree::insert(Rect const &r, int shape) { + insert(r[0].min(), r[1].min(), + r[0].max(), r[1].max(), shape); +} + + Quad* QuadTree::search(double x0, double y0, double x1, double y1) { Quad *q = root; @@ -117,6 +129,8 @@ void QuadTree::erase(Quad *q, int shape) { return; } +}; + /* Local Variables: mode:c++ diff --git a/src/2geom/quadtree.h b/src/2geom/quadtree.h index 9b5b75e62..716c56673 100644 --- a/src/2geom/quadtree.h +++ b/src/2geom/quadtree.h @@ -1,6 +1,10 @@ #include #include +#include "d2.h" + +namespace Geom{ + class Quad{ public: Quad* children[4]; @@ -10,6 +14,23 @@ public: children[i] = 0; } typedef std::vector::iterator iterator; + Rect bounds(unsigned i, double x, double y, double d) { + double dd = d/2; + switch(i % 4) { + case 0: + return Rect(Interval(x, x+dd), Interval(y, y+dd)); + case 1: + return Rect(Interval(x+dd, x+d), Interval(y, y+dd)); + case 2: + return Rect(Interval(x, x+dd), Interval(y+dd, y+d)); + case 3: + return Rect(Interval(x+dd, x+d), Interval(y+dd, y+d)); + default: + /* just to suppress warning message + * this case should be never reached */ + assert(false); + } + } }; class QuadTree{ @@ -23,9 +44,12 @@ public: Quad* search(double x0, double y0, double x1, double y1); void insert(double x0, double y0, double x1, double y1, int shape); + Quad* search(Geom::Rect const &r); + void insert(Geom::Rect const &r, int shape); void erase(Quad *q, int shape); }; +}; /* Local Variables: diff --git a/src/2geom/sbasis-geometric.cpp b/src/2geom/sbasis-geometric.cpp index b30c3a655..655830389 100644 --- a/src/2geom/sbasis-geometric.cpp +++ b/src/2geom/sbasis-geometric.cpp @@ -157,6 +157,17 @@ Geom::atan2(D2 const &vect, double tol, unsigned order){ return atan2(Piecewise >(vect),tol,order); } +//tan2 is the pseudo-inverse of atan2. It takes an angle and returns a unit_vector that points in the direction of angle. +D2 > +Geom::tan2(SBasis const &angle, double tol, unsigned order){ + return tan2(Piecewise(angle), tol, order); +} + +D2 > +Geom::tan2(Piecewise const &angle, double tol, unsigned order){ + return D2 >(cos(angle, tol, order), sin(angle, tol, order)); +} + //unitVector(x,y) is computed as (b,-a) where a and b are solutions of: // ax+by=0 (eqn1) and a^2+b^2=1 (eqn2) Piecewise > diff --git a/src/2geom/sbasis-geometric.h b/src/2geom/sbasis-geometric.h index c4f139aa6..6e21a6395 100644 --- a/src/2geom/sbasis-geometric.h +++ b/src/2geom/sbasis-geometric.h @@ -33,6 +33,14 @@ Piecewise atan2(Piecewise >const &vect, double tol=.01, unsigned order=3); +D2 > +tan2(SBasis const &angle, + double tol=.01, unsigned order=3); + +D2 > +tan2(Piecewise const &angle, + double tol=.01, unsigned order=3); + Piecewise > unitVector(D2 const &vect, double tol=.01, unsigned order=3); diff --git a/src/2geom/sbasis.h b/src/2geom/sbasis.h index 0a48ff1cf..d6f55c8d8 100644 --- a/src/2geom/sbasis.h +++ b/src/2geom/sbasis.h @@ -296,6 +296,7 @@ inline std::ostream &operator<< (std::ostream &out_file, const SBasis & p) { return out_file; } +// These are deprecated, use sbasis-math versions if possible SBasis sin(Linear bo, int k); SBasis cos(Linear bo, int k); diff --git a/src/2geom/solve-bezier-one-d.cpp b/src/2geom/solve-bezier-one-d.cpp index cd7c36cc3..3ec738df1 100644 --- a/src/2geom/solve-bezier-one-d.cpp +++ b/src/2geom/solve-bezier-one-d.cpp @@ -86,6 +86,11 @@ Bernsteins::find_bernstein_roots(double const *w, /* The control points */ /* if deep enough, return 1 solution at midpoint */ if (depth >= MAXDEPTH) { //printf("bottom out %d\n", depth); + const double Ax = right_t - left_t; + const double Ay = w[degree] - w[0]; + + solutions.push_back(left_t - Ax*w[0] / Ay); + return; solutions.push_back((left_t + right_t) / 2.0); return; } @@ -101,6 +106,7 @@ Bernsteins::find_bernstein_roots(double const *w, /* The control points */ solutions.push_back(left_t - Ax*w[0] / Ay); return; } + } /* Otherwise, solve recursively after subdividing control polygon */ -- 2.30.2