diff options
Diffstat (limited to 'basegfx/source/curve')
| -rw-r--r-- | basegfx/source/curve/b2dbeziertools.cxx | 533 |
1 files changed, 5 insertions, 528 deletions
diff --git a/basegfx/source/curve/b2dbeziertools.cxx b/basegfx/source/curve/b2dbeziertools.cxx index 92008a6158a3..3cf7484169e9 100644 --- a/basegfx/source/curve/b2dbeziertools.cxx +++ b/basegfx/source/curve/b2dbeziertools.cxx @@ -4,9 +4,9 @@ * * $RCSfile: b2dbeziertools.cxx,v $ * - * $Revision: 1.11 $ + * $Revision: 1.12 $ * - * last change: $Author: obo $ $Date: 2006-09-17 07:58:15 $ + * last change: $Author: obo $ $Date: 2007-07-18 11:05:04 $ * * The Contents of this file are made available subject to * the terms of GNU Lesser General Public License Version 2.1. @@ -36,534 +36,11 @@ // MARKER(update_precomp.py): autogen include statement, do not remove #include "precompiled_basegfx.hxx" -#include <limits> -#include <algorithm> -#ifndef _USE_MATH_DEFINES -#define _USE_MATH_DEFINES // needed by Visual C++ for math constants -#endif -#include <math.h> // M_PI definition - #ifndef _BGFX_CURVE_B2DCUBICBEZIER_HXX #include <basegfx/curve/b2dcubicbezier.hxx> #endif -#ifndef _BGFX_CURVE_B2DQUADRATICBEZIER_HXX -#include <basegfx/curve/b2dquadraticbezier.hxx> -#endif - -#ifndef _BGFX_POLYGON_B2DPOLYGON_HXX -#include <basegfx/polygon/b2dpolygon.hxx> -#endif - -#ifndef _BGFX_VECTOR_B2DVECTOR_HXX -#include <basegfx/vector/b2dvector.hxx> -#endif - -#ifndef _BGFX_POINT_B2DPOINT_HXX -#include <basegfx/point/b2dpoint.hxx> -#endif - -#ifndef _BGFX_NUMERIC_FTOOLS_HXX -#include <basegfx/numeric/ftools.hxx> -#endif - -#include <basegfx/curve/b2dbeziertools.hxx> - - -namespace -{ - class DistanceErrorFunctor - { - public: - DistanceErrorFunctor( const double& distance ) : - mfDistance2( distance*distance ), - mfLastDistanceError2( ::std::numeric_limits<double>::max() ) - { - } - - bool needsFurtherSubdivision( const double& P1x, const double& P1y, - const double& P2x, const double& P2y, - const double& P3x, const double& P3y, - const double& P4x, const double& P4y, - const double&, const double& ) // last two values not used here - { - // Perform bezier flatness test (lecture notes from R. Schaback, - // Mathematics of Computer-Aided Design, Uni Goettingen, 2000) - // - // ||P(t) - L(t)|| <= max ||b_j - b_0 - j/n(b_n - b_0)|| - // 0<=j<=n - // - // What is calculated here is an upper bound to the distance from - // a line through b_0 and b_3 (P1 and P4 in our notation) and the - // curve. We can drop 0 and n from the running indices, since the - // argument of max becomes zero for those cases. - const double fJ1x( P2x - P1x - 1.0/3.0*(P4x - P1x) ); - const double fJ1y( P2y - P1y - 1.0/3.0*(P4y - P1y) ); - const double fJ2x( P3x - P1x - 2.0/3.0*(P4x - P1x) ); - const double fJ2y( P3y - P1y - 2.0/3.0*(P4y - P1y) ); - const double distanceError2( ::std::max( fJ1x*fJ1x + fJ1y*fJ1y, - fJ2x*fJ2x + fJ2y*fJ2y) ); - - // stop if error measure does not improve anymore. This is a - // safety guard against floating point inaccuracies. - // stop if distance from line is guaranteed to be bounded by d - bool bRet( mfLastDistanceError2 > distanceError2 && - distanceError2 >= mfDistance2 ); - - mfLastDistanceError2 = distanceError2; - - return bRet; - } - - private: - double mfDistance2; - double mfLastDistanceError2; - }; - - - class AngleErrorFunctor - { - public: - AngleErrorFunctor( const double& angleBounds ) : - // take positive branch, convert from degree to radiant - mfTanAngle( tan( fabs(angleBounds)/180.0*M_PI ) ), - mfLastTanAngle( ::std::numeric_limits<double>::max() ) - { - } - - bool needsFurtherSubdivision( const double P1x, const double P1y, - const double P2x, const double P2y, - const double P3x, const double P3y, - const double P4x, const double P4y, - const double Pdx, const double Pdy ) - { - // Test angle differences between two lines (ad - // and bd), meeting in the t=0.5 division point - // (d), and the angle from the other ends of those - // lines (b and a, resp.) to the tangents to the - // curve at this points: - // - // *__________ - // ......*b - // ... - // .. - // . - // * *d - // | . - // | . - // | . - // | . - // |. - // |. - // * - // a - // - // The idea of this test is based on the fact that - // small angle deviations between straight lines - // go by mostly unnoticed by the human eye. When - // we stop subdivision according to this test, the - // resulting line polygon is guaranteed to have - // angle differences between adjacent lines lower - // than the given threshold. - // - // When using half of the angle bound for the - // difference to the tangents at a or b, resp., - // this procedure guarantees that no angle in the - // resulting line polygon is larger than the - // specified angle bound. This is because during - // subdivision, adjacent curve segments will have - // collinear tangent vectors, thus, when each - // side's line segments differs by at most angle/2 - // from that tangent, the summed difference will - // be at most angle (this was modeled after an - // idea from Armin Weiss). - - // To stay within the notation above, a equals P1, - // the other end point of the tangent starting at - // a is P2, d is Pd, and so forth. - const ::basegfx::B2DVector vecAD( Pdx - P1x, Pdy - P1y ); - const ::basegfx::B2DVector vecDB( P4x - Pdx, P4y - Pdy ); - - const double scalarVecADDB( vecAD.scalar( vecDB ) ); - const double crossVecADDB( vecAD.cross( vecDB ) ); - - const ::basegfx::B2DVector vecStartTangent( P2x - P1x, P2y - P1y ); - const ::basegfx::B2DVector vecEndTangent( P4x - P3x, P4y - P3y ); - - const double scalarVecStartTangentAD( vecStartTangent.scalar( vecAD ) ); - const double crossVecStartTangentAD( vecStartTangent.cross( vecAD ) ); - - const double scalarVecDBEndTangent( vecDB.scalar( vecEndTangent ) ); - const double crossVecDBEndTangent( vecDB.cross( vecEndTangent ) ); - - double fCurrAngle( ::std::numeric_limits<double>::max() ); - - // are vecAD, vecDB, start and end tangent collinear? If - // yes, we're already done, and no further subdivision can - // bring us any closer to a straight line... - if( ::basegfx::fTools::equalZero( crossVecADDB ) && - ::basegfx::fTools::equalZero( crossVecStartTangentAD ) && - ::basegfx::fTools::equalZero( crossVecDBEndTangent ) ) - { - mfLastTanAngle = 0.0; - return false; - } - - // anyone has zero denominator? then we're at - // +infinity, anyway, and should better keep on - // subdividing - if( ::basegfx::fTools::equalZero( scalarVecADDB ) || - ::basegfx::fTools::equalZero( scalarVecStartTangentAD ) || - ::basegfx::fTools::equalZero( scalarVecDBEndTangent ) ) - { - mfLastTanAngle = fCurrAngle; - return true; - } - - // now, sieve out quadrants which have - // deviating angles of more than 90 - // degrees. By inspection, this is everything - // with a negative scalar product. If we - // encounter such a negative scalar product, - // we can simply keep on subdividing, since at - // least one angle is then >90 degrees. - if( ::basegfx::fTools::less( scalarVecADDB, 0.0 ) || - ::basegfx::fTools::less( scalarVecStartTangentAD, 0.0 ) || - ::basegfx::fTools::less( scalarVecDBEndTangent, 0.0 ) ) - { - mfLastTanAngle = fCurrAngle; // TODO: Do we need the - // correct value here? - return true; - } - - // take the maximum tangens of angle - // deviation, to compare against the threshold - // below - fCurrAngle = ::std::max( fabs( crossVecADDB / scalarVecADDB ), - ::std::max( fabs( crossVecStartTangentAD / scalarVecStartTangentAD ), - fabs( crossVecDBEndTangent / scalarVecDBEndTangent ) ) ); - - // stop if error measure does not improve anymore. This is a - // safety guard against floating point inaccuracies. - // stop if angle difference is guaranteed to be bounded by mfTanAngle - bool bRet( mfLastTanAngle > fCurrAngle && - fCurrAngle >= mfTanAngle ); - - mfLastTanAngle = fCurrAngle; - - return bRet; - } - - private: - double mfTanAngle; - double mfLastTanAngle; - }; - - /* Recursively subdivide cubic bezier curve via deCasteljau. - - @param rPoly - Polygon to append generated points to - - @param d2 - Maximal squared difference of curve to a straight line - - @param P* - Exactly four points, interpreted as support and control points of - a cubic bezier curve. - - @param old_distance2 - Last squared distance to line for this recursion - path. Used as an end condition, if it is no longer - improving. - - @param recursionDepth - Depth of recursion. Used as a termination criterion, to - prevent endless looping. - */ - template < class ErrorFunctor > int ImplAdaptiveSubdivide( ::basegfx::B2DPolygon& rPoly, - ErrorFunctor rErrorFunctor, - const double P1x, const double P1y, - const double P2x, const double P2y, - const double P3x, const double P3y, - const double P4x, const double P4y, - int recursionDepth ) - { - // Hard limit on recursion depth, empiric number. Note that if - // for some obscure reason, we're subdividing always, this - // would lead to 2^maxRecursionDepth points generated - enum {maxRecursionDepth=30}; - - // deCasteljau bezier arc, split at t=0.5 - // Foley/vanDam, p. 508 - - // Note that for the pure distance error method, this - // subdivision could be moved into the if-branch. But - // since this accounts for saved work only for the - // very last subdivision step, and we need the - // subdivided curve for the angle criterium, I think - // it's justified here. - const double L1x( P1x ), L1y( P1y ); - const double L2x( (P1x + P2x)*0.5 ), L2y( (P1y + P2y)*0.5 ); - const double Hx ( (P2x + P3x)*0.5 ), Hy ( (P2y + P3y)*0.5 ); - const double L3x( (L2x + Hx)*0.5 ), L3y( (L2y + Hy)*0.5 ); - const double R4x( P4x ), R4y( P4y ); - const double R3x( (P3x + P4x)*0.5 ), R3y( (P3y + P4y)*0.5 ); - const double R2x( (Hx + R3x)*0.5 ), R2y( (Hy + R3y)*0.5 ); - const double R1x( (L3x + R2x)*0.5 ), R1y( (L3y + R2y)*0.5 ); - const double L4x( R1x ), L4y( R1y ); - - // stop at recursion level 128. This is a safety guard against - // floating point inaccuracies. - if( recursionDepth < maxRecursionDepth && - rErrorFunctor.needsFurtherSubdivision( P1x, P1y, - P2x, P2y, - P3x, P3y, - P4x, P4y, - R1x, R1y ) ) - { - // subdivide further - ++recursionDepth; - - int nGeneratedPoints(0); - - nGeneratedPoints += ImplAdaptiveSubdivide(rPoly, rErrorFunctor, L1x, L1y, L2x, L2y, L3x, L3y, L4x, L4y, recursionDepth); - nGeneratedPoints += ImplAdaptiveSubdivide(rPoly, rErrorFunctor, R1x, R1y, R2x, R2y, R3x, R3y, R4x, R4y, recursionDepth); - - // return number of points generated in this - // recursion branch - return nGeneratedPoints; - } - else - { - // requested resolution reached. Add end points to - // output iterator. order is preserved, since - // this is so to say depth first traversal. - rPoly.append( ::basegfx::B2DPoint( P1x, P1y ) ); - - // return number of points generated in this - // recursion branch - return 1; - } - } -} - -namespace basegfx -{ - -// TODO: For bezier length calculations, use Jens Gravesens algorithm (e.g. Graphic Gems V, p. 199) -// LATER -#if 0 - /* Approximate given cubic bezier curve by quadratic bezier segments */ - void ImplQuadBezierApprox( B2DPolygon& rPoly, - BitStream& rBits, - Point& rLastPoint, - const double d2, - const double P1x, const double P1y, - const double P2x, const double P2y, - const double P3x, const double P3y, - const double P4x, const double P4y ) - { - // Check for degenerate case, where the given cubic bezier curve - // is already quadratic: P4 == 3P3 - 3P2 + P1 - if( P4x == 3.0*P3x - 3.0*P2x + P1x && - P4y == 3.0*P3y - 3.0*P2y + P1y ) - { - Impl_addQuadBezier( rBits, rLastPoint, - 3.0/2.0*P2x - 1.0/2.0*P1x, 3.0/2.0*P2y - 1.0/2.0*P1y, - P4x, P4y); - } - else - { - // Create quadratic segment for given cubic: - // Start and end point must coincide, determine quadratic control - // point in such a way that it lies on the intersection of the - // tangents at start and end point, resp. Thus, both cubic and - // quadratic curve segments will match in 0th and 1st derivative - // at the start and end points - - // Intersection of P2P1 and P4P3 - // (P2y-P4y)(P3x-P4x)-(P2x-P4x)(P3y-P4y) - // lambda = ------------------------------------- - // (P1x-P2x)(P3y-P4y)-(P1y-P2y)(P3x-P4x) - // - // Intersection point IP is now - // IP = P2 + lambda(P1-P2) - // - const double nominator( (P2y-P4y)*(P3x-P4x) - (P2x-P4x)*(P3y-P4y) ); - const double denominator( (P1x-P2x)*(P3y-P4y) - (P1y-P2y)*(P3x-P4x) ); - const double lambda( nominator / denominator ); - - const double IPx( P2x + lambda*( P1x - P2x) ); - const double IPy( P2y + lambda*( P1y - P2y) ); - - // Introduce some alias names: quadratic start point is P1, end - // point is P4, control point is IP - const double QP1x( P1x ); - const double QP1y( P1y ); - const double QP2x( IPx ); - const double QP2y( IPy ); - const double QP3x( P4x ); - const double QP3y( P4y ); - - // Adapted bezier flatness test (lecture notes from R. Schaback, - // Mathematics of Computer-Aided Design, Uni Goettingen, 2000) - // - // ||C(t) - Q(t)|| <= max ||c_j - q_j|| - // 0<=j<=n - // - // In this case, we don't need the distance from the cubic bezier - // to a straight line, but to a quadratic bezier. The c_j's are - // the cubic bezier's bernstein coefficients, the q_j's the - // quadratic bezier's. We have the c_j's given, the q_j's can be - // calculated from QPi like this (sorry, mixed index notation, we - // use [1,n], formulas use [0,n-1]): - // - // q_0 = QP1 = P1 - // q_1 = 1/3 QP1 + 2/3 QP2 - // q_2 = 2/3 QP2 + 1/3 QP3 - // q_3 = QP3 = P4 - // - // We can drop case 0 and 3, since there the curves coincide - // (distance is zero) - - // calculate argument of max for j=1 and j=2 - const double fJ1x( P2x - 1.0/3.0*QP1x - 2.0/3.0*QP2x ); - const double fJ1y( P2y - 1.0/3.0*QP1y - 2.0/3.0*QP2y ); - const double fJ2x( P3x - 2.0/3.0*QP2x - 1.0/3.0*QP3x ); - const double fJ2y( P3y - 2.0/3.0*QP2y - 1.0/3.0*QP3y ); - - // stop if distance from cubic curve is guaranteed to be bounded by d - // Should denominator be 0: then P1P2 and P3P4 are parallel (P1P2^T R[90,P3P4] = 0.0), - // meaning that either we have a straight line or an inflexion point (see else block below) - if( 0.0 != denominator && - ::std::max( fJ1x*fJ1x + fJ1y*fJ1y, - fJ2x*fJ2x + fJ2y*fJ2y) < d2 ) - { - // requested resolution reached. - // Add end points to output file. - // order is preserved, since this is so to say depth first traversal. - Impl_addQuadBezier( rBits, rLastPoint, - QP2x, QP2y, - QP3x, QP3y); - } - else - { - // Maybe subdivide further - - // This is for robustness reasons, since the line intersection - // method below gets instable if the curve gets closer to a - // straight line. If the given cubic bezier does not deviate by - // more than d/4 from a straight line, either: - // - take the line (that's what we do here) - // - express the line by a quadratic bezier - - // Perform bezier flatness test (lecture notes from R. Schaback, - // Mathematics of Computer-Aided Design, Uni Goettingen, 2000) - // - // ||P(t) - L(t)|| <= max ||b_j - b_0 - j/n(b_n - b_0)|| - // 0<=j<=n - // - // What is calculated here is an upper bound to the distance from - // a line through b_0 and b_3 (P1 and P4 in our notation) and the - // curve. We can drop 0 and n from the running indices, since the - // argument of max becomes zero for those cases. - const double fJ1x( P2x - P1x - 1.0/3.0*(P4x - P1x) ); - const double fJ1y( P2y - P1y - 1.0/3.0*(P4y - P1y) ); - const double fJ2x( P3x - P1x - 2.0/3.0*(P4x - P1x) ); - const double fJ2y( P3y - P1y - 2.0/3.0*(P4y - P1y) ); - - // stop if distance from line is guaranteed to be bounded by d/4 - if( ::std::max( fJ1x*fJ1x + fJ1y*fJ1y, - fJ2x*fJ2x + fJ2y*fJ2y) < d2/16.0 ) - { - // do not subdivide further, add straight line instead - Impl_addStraightLine( rBits, rLastPoint, P4x, P4y); - } - else - { - // deCasteljau bezier arc, split at t=0.5 - // Foley/vanDam, p. 508 - const double L1x( P1x ), L1y( P1y ); - const double L2x( (P1x + P2x)*0.5 ), L2y( (P1y + P2y)*0.5 ); - const double Hx ( (P2x + P3x)*0.5 ), Hy ( (P2y + P3y)*0.5 ); - const double L3x( (L2x + Hx)*0.5 ), L3y( (L2y + Hy)*0.5 ); - const double R4x( P4x ), R4y( P4y ); - const double R3x( (P3x + P4x)*0.5 ), R3y( (P3y + P4y)*0.5 ); - const double R2x( (Hx + R3x)*0.5 ), R2y( (Hy + R3y)*0.5 ); - const double R1x( (L3x + R2x)*0.5 ), R1y( (L3y + R2y)*0.5 ); - const double L4x( R1x ), L4y( R1y ); - - // subdivide further - Impl_quadBezierApprox(rBits, rLastPoint, d2, L1x, L1y, L2x, L2y, L3x, L3y, L4x, L4y); - Impl_quadBezierApprox(rBits, rLastPoint, d2, R1x, R1y, R2x, R2y, R3x, R3y, R4x, R4y); - } - } - } - } - } -#endif - - sal_Int32 adaptiveSubdivideByDistance( B2DPolygon& rPoly, - const B2DCubicBezier& rCurve, - double distanceBounds, - bool bAddEndPoint ) - { - const B2DPoint start( rCurve.getStartPoint() ); - const B2DPoint control1( rCurve.getControlPointA() ); - const B2DPoint control2( rCurve.getControlPointB() ); - const B2DPoint end( rCurve.getEndPoint() ); - - sal_Int32 nPoints( ImplAdaptiveSubdivide( rPoly, - DistanceErrorFunctor( distanceBounds ), - start.getX(), start.getY(), - control1.getX(), control1.getY(), - control2.getX(), control2.getY(), - end.getX(), end.getY(), - 0 ) ); - // finish polygon - if ( bAddEndPoint ) - { - rPoly.append( end ); - nPoints++; - } - - return nPoints; - } - - sal_Int32 adaptiveSubdivideByAngle( B2DPolygon& rPoly, - const B2DCubicBezier& rCurve, - double angleBounds, - bool bAddEndPoint ) - { - const B2DPoint start( rCurve.getStartPoint() ); - const B2DPoint control1( rCurve.getControlPointA() ); - const B2DPoint control2( rCurve.getControlPointB() ); - const B2DPoint end( rCurve.getEndPoint() ); - - sal_Int32 nPoints( ImplAdaptiveSubdivide( rPoly, - AngleErrorFunctor( angleBounds ), - start.getX(), start.getY(), - control1.getX(), control1.getY(), - control2.getX(), control2.getY(), - end.getX(), end.getY(), - 0 ) ); - // finish polygon - if ( bAddEndPoint ) - { - rPoly.append( end ); - nPoints++; - } - - return nPoints; - } - - sal_Int32 adaptiveSubdivideByDistance( B2DPolygon& /*rPoly*/, - const B2DQuadraticBezier& /*rCurve*/, - double /*distanceBounds*/, - bool /*bAddEndPoint*/ ) - { - // TODO - return 0; - } -} +////////////////////////////////////////////////////////////////////////////// +// necessary stuff moved to B2DCubicBezier +// eof |