basegfx/source/curve/b2dcubicbezier.cxx

   1 /*************************************************************************
   2  *
   3  * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
   4  *
   5  * Copyright 2008 by Sun Microsystems, Inc.
   6  *
   7  * OpenOffice.org - a multi-platform office productivity suite
   8  *
   9  * $RCSfile: b2dcubicbezier.cxx,v $
  10  * $Revision: 1.16 $
  11  *
  12  * This file is part of OpenOffice.org.
  13  *
  14  * OpenOffice.org is free software: you can redistribute it and/or modify
  15  * it under the terms of the GNU Lesser General Public License version 3
  16  * only, as published by the Free Software Foundation.
  17  *
  18  * OpenOffice.org is distributed in the hope that it will be useful,
  19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  21  * GNU Lesser General Public License version 3 for more details
  22  * (a copy is included in the LICENSE file that accompanied this code).
  23  *
  24  * You should have received a copy of the GNU Lesser General Public License
  25  * version 3 along with OpenOffice.org.  If not, see
  26  * <http://www.openoffice.org/license.html>
  27  * for a copy of the LGPLv3 License.
  28  *
  29  ************************************************************************/
  30
  31 // MARKER(update_precomp.py): autogen include statement, do not remove
  32 #include "precompiled_basegfx.hxx"
  33 #include <basegfx/curve/b2dcubicbezier.hxx>
  34 #include <basegfx/vector/b2dvector.hxx>
  35 #include <basegfx/polygon/b2dpolygon.hxx>
  36 #include <basegfx/numeric/ftools.hxx>
  37
  38 #include <limits>
  39
  40 // #i37443#
  41 #define FACTOR_FOR_UNSHARPEN    (1.6)
  42 #ifdef DBG_UTIL
  43 static double fMultFactUnsharpen = FACTOR_FOR_UNSHARPEN;
  44 #endif
  45
  46 //////////////////////////////////////////////////////////////////////////////
  47
  48 namespace basegfx
  49 {
  50     namespace
  51     {
  52         void ImpSubDivAngle(
  53             const B2DPoint& rfPA,                       // start point
  54             const B2DPoint& rfEA,                       // edge on A
  55             const B2DPoint& rfEB,                       // edge on B
  56             const B2DPoint& rfPB,                       // end point
  57             B2DPolygon& rTarget,                        // target polygon
  58             double fAngleBound,                         // angle bound in [0.0 .. 2PI]
  59             bool bAllowUnsharpen,                       // #i37443# allow the criteria to get unsharp in recursions
  60             sal_uInt16 nMaxRecursionDepth)      // endless loop protection
  61         {
  62             if(nMaxRecursionDepth)
  63             {
  64                 // do angle test
  65                 B2DVector aLeft(rfEA - rfPA);
  66                 B2DVector aRight(rfEB - rfPB);
  67
  68                 // #i72104#
  69                 if(aLeft.equalZero())
  70                 {
  71                     aLeft = rfEB - rfPA;
  72                 }
  73
  74                 if(aRight.equalZero())
  75                 {
  76                     aRight = rfEA - rfPB;
  77                 }
  78
  79                 const double fCurrentAngle(aLeft.angle(aRight));
  80
  81                 if(fabs(fCurrentAngle) > (F_PI - fAngleBound))
  82                 {
  83                     // end recursion
  84                     nMaxRecursionDepth = 0;
  85                 }
  86                 else
  87                 {
  88                     if(bAllowUnsharpen)
  89                     {
  90                         // #i37443# unsharpen criteria
  91 #ifdef DBG_UTIL
  92                         fAngleBound *= fMultFactUnsharpen;
  93 #else
  94                         fAngleBound *= FACTOR_FOR_UNSHARPEN;
  95 #endif
  96                     }
  97                 }
  98             }
  99
 100             if(nMaxRecursionDepth)
 101             {
 102                 // divide at 0.5
 103                 const B2DPoint aS1L(average(rfPA, rfEA));
 104                 const B2DPoint aS1C(average(rfEA, rfEB));
 105                 const B2DPoint aS1R(average(rfEB, rfPB));
 106                 const B2DPoint aS2L(average(aS1L, aS1C));
 107                 const B2DPoint aS2R(average(aS1C, aS1R));
 108                 const B2DPoint aS3C(average(aS2L, aS2R));
 109
 110                 // left recursion
 111                 ImpSubDivAngle(rfPA, aS1L, aS2L, aS3C, rTarget, fAngleBound, bAllowUnsharpen, nMaxRecursionDepth - 1);
 112
 113                 // right recursion
 114                 ImpSubDivAngle(aS3C, aS2R, aS1R, rfPB, rTarget, fAngleBound, bAllowUnsharpen, nMaxRecursionDepth - 1);
 115             }
 116             else
 117             {
 118                 rTarget.append(rfPB);
 119             }
 120         }
 121
 122         void ImpSubDivAngleStart(
 123             const B2DPoint& rfPA,                       // start point
 124             const B2DPoint& rfEA,                       // edge on A
 125             const B2DPoint& rfEB,                       // edge on B
 126             const B2DPoint& rfPB,                       // end point
 127             B2DPolygon& rTarget,                        // target polygon
 128             const double& rfAngleBound,         // angle bound in [0.0 .. 2PI]
 129             bool bAllowUnsharpen)                       // #i37443# allow the criteria to get unsharp in recursions
 130         {
 131             sal_uInt16 nMaxRecursionDepth(8);
 132             const B2DVector aLeft(rfEA - rfPA);
 133             const B2DVector aRight(rfEB - rfPB);
 134             bool bLeftEqualZero(aLeft.equalZero());
 135             bool bRightEqualZero(aRight.equalZero());
 136             bool bAllParallel(false);
 137
 138             if(bLeftEqualZero && bRightEqualZero)
 139             {
 140                 nMaxRecursionDepth = 0;
 141             }
 142             else
 143             {
 144                 const B2DVector aBase(rfPB - rfPA);
 145                 const bool bBaseEqualZero(aBase.equalZero()); // #i72104#
 146
 147                 if(!bBaseEqualZero)
 148                 {
 149                     const bool bLeftParallel(bLeftEqualZero ? true : areParallel(aLeft, aBase));
 150                     const bool bRightParallel(bRightEqualZero ? true : areParallel(aRight, aBase));
 151
 152                     if(bLeftParallel && bRightParallel)
 153                     {
 154                         bAllParallel = true;
 155
 156                         if(!bLeftEqualZero)
 157                         {
 158                             double fFactor;
 159
 160                             if(fabs(aBase.getX()) > fabs(aBase.getY()))
 161                             {
 162                                 fFactor = aLeft.getX() / aBase.getX();
 163                             }
 164                             else
 165                             {
 166                                 fFactor = aLeft.getY() / aBase.getY();
 167                             }
 168
 169                             if(fFactor >= 0.0 && fFactor <= 1.0)
 170                             {
 171                                 bLeftEqualZero = true;
 172                             }
 173                         }
 174
 175                         if(!bRightEqualZero)
 176                         {
 177                             double fFactor;
 178
 179                             if(fabs(aBase.getX()) > fabs(aBase.getY()))
 180                             {
 181                                 fFactor = aRight.getX() / -aBase.getX();
 182                             }
 183                             else
 184                             {
 185                                 fFactor = aRight.getY() / -aBase.getY();
 186                             }
 187
 188                             if(fFactor >= 0.0 && fFactor <= 1.0)
 189                             {
 190                                 bRightEqualZero = true;
 191                             }
 192                         }
 193
 194                         if(bLeftEqualZero && bRightEqualZero)
 195                         {
 196                             nMaxRecursionDepth = 0;
 197                         }
 198                     }
 199                 }
 200             }
 201
 202             if(nMaxRecursionDepth)
 203             {
 204                 // divide at 0.5 ad test both edges for angle criteria
 205                 const B2DPoint aS1L(average(rfPA, rfEA));
 206                 const B2DPoint aS1C(average(rfEA, rfEB));
 207                 const B2DPoint aS1R(average(rfEB, rfPB));
 208                 const B2DPoint aS2L(average(aS1L, aS1C));
 209                 const B2DPoint aS2R(average(aS1C, aS1R));
 210                 const B2DPoint aS3C(average(aS2L, aS2R));
 211
 212                 // test left
 213                 bool bAngleIsSmallerLeft(bAllParallel && bLeftEqualZero);
 214                 if(!bAngleIsSmallerLeft)
 215                 {
 216                     const B2DVector aLeftLeft(bLeftEqualZero ? aS2L - aS1L : aS1L - rfPA); // #i72104#
 217                     const B2DVector aRightLeft(aS2L - aS3C);
 218                     const double fCurrentAngleLeft(aLeftLeft.angle(aRightLeft));
 219                     bAngleIsSmallerLeft = (fabs(fCurrentAngleLeft) > (F_PI - rfAngleBound));
 220                 }
 221
 222                 // test right
 223                 bool bAngleIsSmallerRight(bAllParallel && bRightEqualZero);
 224                 if(!bAngleIsSmallerRight)
 225                 {
 226                     const B2DVector aLeftRight(aS2R - aS3C);
 227                     const B2DVector aRightRight(bRightEqualZero ? aS2R - aS1R : aS1R - rfPB); // #i72104#
 228                     const double fCurrentAngleRight(aLeftRight.angle(aRightRight));
 229                     bAngleIsSmallerRight = (fabs(fCurrentAngleRight) > (F_PI - rfAngleBound));
 230                 }
 231
 232                 if(bAngleIsSmallerLeft && bAngleIsSmallerRight)
 233                 {
 234                     // no recursion necessary at all
 235                     nMaxRecursionDepth = 0;
 236                 }
 237                 else
 238                 {
 239                     // left
 240                     if(bAngleIsSmallerLeft)
 241                     {
 242                         rTarget.append(aS3C);
 243                     }
 244                     else
 245                     {
 246                         ImpSubDivAngle(rfPA, aS1L, aS2L, aS3C, rTarget, rfAngleBound, bAllowUnsharpen, nMaxRecursionDepth);
 247                     }
 248
 249                     // right
 250                     if(bAngleIsSmallerRight)
 251                     {
 252                         rTarget.append(rfPB);
 253                     }
 254                     else
 255                     {
 256                         ImpSubDivAngle(aS3C, aS2R, aS1R, rfPB, rTarget, rfAngleBound, bAllowUnsharpen, nMaxRecursionDepth);
 257                     }
 258                 }
 259             }
 260
 261             if(!nMaxRecursionDepth)
 262             {
 263                 rTarget.append(rfPB);
 264             }
 265         }
 266
 267         void ImpSubDivDistance(
 268             const B2DPoint& rfPA,                       // start point
 269             const B2DPoint& rfEA,                       // edge on A
 270             const B2DPoint& rfEB,                       // edge on B
 271             const B2DPoint& rfPB,                       // end point
 272             B2DPolygon& rTarget,                        // target polygon
 273             double fDistanceBound2,                     // quadratic distance criteria
 274             double fLastDistanceError2,         // the last quadratic distance error
 275             sal_uInt16 nMaxRecursionDepth)      // endless loop protection
 276         {
 277             if(nMaxRecursionDepth)
 278             {
 279                 // decide if another recursion is needed. If not, set
 280                 // nMaxRecursionDepth to zero
 281
 282                 // Perform bezier flatness test (lecture notes from R. Schaback,
 283                 // Mathematics of Computer-Aided Design, Uni Goettingen, 2000)
 284                 //
 285                 // ||P(t) - L(t)|| <= max     ||b_j - b_0 - j/n(b_n - b_0)||
 286                 //                    0<=j<=n
 287                 //
 288                 // What is calculated here is an upper bound to the distance from
 289                 // a line through b_0 and b_3 (rfPA and P4 in our notation) and the
 290                 // curve. We can drop 0 and n from the running indices, since the
 291                 // argument of max becomes zero for those cases.
 292                 const double fJ1x(rfEA.getX() - rfPA.getX() - 1.0/3.0*(rfPB.getX() - rfPA.getX()));
 293                 const double fJ1y(rfEA.getY() - rfPA.getY() - 1.0/3.0*(rfPB.getY() - rfPA.getY()));
 294                 const double fJ2x(rfEB.getX() - rfPA.getX() - 2.0/3.0*(rfPB.getX() - rfPA.getX()));
 295                 const double fJ2y(rfEB.getY() - rfPA.getY() - 2.0/3.0*(rfPB.getY() - rfPA.getY()));
 296                 const double fDistanceError2(::std::max(fJ1x*fJ1x + fJ1y*fJ1y, fJ2x*fJ2x + fJ2y*fJ2y));
 297
 298                 // stop if error measure does not improve anymore. This is a
 299                 // safety guard against floating point inaccuracies.
 300                 // stop if distance from line is guaranteed to be bounded by d
 301                 const bool bFurtherDivision(fLastDistanceError2 > fDistanceError2 && fDistanceError2 >= fDistanceBound2);
 302
 303                 if(bFurtherDivision)
 304                 {
 305                     // remember last error value
 306                     fLastDistanceError2 = fDistanceError2;
 307                 }
 308                 else
 309                 {
 310                     // stop recustion
 311                     nMaxRecursionDepth = 0;
 312                 }
 313             }
 314
 315             if(nMaxRecursionDepth)
 316             {
 317                 // divide at 0.5
 318                 const B2DPoint aS1L(average(rfPA, rfEA));
 319                 const B2DPoint aS1C(average(rfEA, rfEB));
 320                 const B2DPoint aS1R(average(rfEB, rfPB));
 321                 const B2DPoint aS2L(average(aS1L, aS1C));
 322                 const B2DPoint aS2R(average(aS1C, aS1R));
 323                 const B2DPoint aS3C(average(aS2L, aS2R));
 324
 325                 // left recursion
 326                 ImpSubDivDistance(rfPA, aS1L, aS2L, aS3C, rTarget, fDistanceBound2, fLastDistanceError2, nMaxRecursionDepth - 1);
 327
 328                 // right recursion
 329                 ImpSubDivDistance(aS3C, aS2R, aS1R, rfPB, rTarget, fDistanceBound2, fLastDistanceError2, nMaxRecursionDepth - 1);
 330             }
 331             else
 332             {
 333                 rTarget.append(rfPB);
 334             }
 335         }
 336     } // end of anonymous namespace
 337 } // end of namespace basegfx
 338
 339 //////////////////////////////////////////////////////////////////////////////
 340
 341 namespace basegfx
 342 {
 343     B2DCubicBezier::B2DCubicBezier(const B2DCubicBezier& rBezier)
 344     :   maStartPoint(rBezier.maStartPoint),
 345         maEndPoint(rBezier.maEndPoint),
 346         maControlPointA(rBezier.maControlPointA),
 347         maControlPointB(rBezier.maControlPointB)
 348     {
 349     }
 350
 351     B2DCubicBezier::B2DCubicBezier()
 352     {
 353     }
 354
 355     B2DCubicBezier::B2DCubicBezier(const B2DPoint& rStart, const B2DPoint& rEnd)
 356     :   maStartPoint(rStart),
 357         maEndPoint(rEnd),
 358         maControlPointA(rStart),
 359         maControlPointB(rEnd)
 360     {
 361     }
 362
 363     B2DCubicBezier::B2DCubicBezier(const B2DPoint& rStart, const B2DPoint& rControlPointA, const B2DPoint& rControlPointB, const B2DPoint& rEnd)
 364     :   maStartPoint(rStart),
 365         maEndPoint(rEnd),
 366         maControlPointA(rControlPointA),
 367         maControlPointB(rControlPointB)
 368     {
 369     }
 370
 371     B2DCubicBezier::~B2DCubicBezier()
 372     {
 373     }
 374
 375     // assignment operator
 376     B2DCubicBezier& B2DCubicBezier::operator=(const B2DCubicBezier& rBezier)
 377     {
 378         maStartPoint = rBezier.maStartPoint;
 379         maEndPoint = rBezier.maEndPoint;
 380         maControlPointA = rBezier.maControlPointA;
 381         maControlPointB = rBezier.maControlPointB;
 382
 383         return *this;
 384     }
 385
 386     // compare operators
 387     bool B2DCubicBezier::operator==(const B2DCubicBezier& rBezier) const
 388     {
 389         return (
 390             maStartPoint == rBezier.maStartPoint
 391             && maEndPoint == rBezier.maEndPoint
 392             && maControlPointA == rBezier.maControlPointA
 393             && maControlPointB == rBezier.maControlPointB
 394         );
 395     }
 396
 397     bool B2DCubicBezier::operator!=(const B2DCubicBezier& rBezier) const
 398     {
 399         return (
 400             maStartPoint != rBezier.maStartPoint
 401             || maEndPoint != rBezier.maEndPoint
 402             || maControlPointA != rBezier.maControlPointA
 403             || maControlPointB != rBezier.maControlPointB
 404         );
 405     }
 406
 407     bool B2DCubicBezier::equal(const B2DCubicBezier& rBezier) const
 408     {
 409         return (
 410             maStartPoint.equal(rBezier.maStartPoint)
 411             && maEndPoint.equal(rBezier.maEndPoint)
 412             && maControlPointA.equal(rBezier.maControlPointA)
 413             && maControlPointB.equal(rBezier.maControlPointB)
 414         );
 415     }
 416
 417     // test if vectors are used
 418     bool B2DCubicBezier::isBezier() const
 419     {
 420         if(maControlPointA != maStartPoint || maControlPointB != maEndPoint)
 421         {
 422             return true;
 423         }
 424
 425         return false;
 426     }
 427
 428     void B2DCubicBezier::testAndSolveTrivialBezier()
 429     {
 430         if(maControlPointA != maStartPoint || maControlPointB != maEndPoint)
 431         {
 432             const B2DVector aEdge(maEndPoint - maStartPoint);
 433
 434             // controls parallel to edge can be trivial. No edge -> not parallel -> control can
 435             // still not be trivial (e.g. ballon loop)
 436             if(!aEdge.equalZero())
 437             {
 438                 // get control vectors
 439                 const B2DVector aVecA(maControlPointA - maStartPoint);
 440                 const B2DVector aVecB(maControlPointB - maEndPoint);
 441
 442                 // check if trivial per se
 443                 bool bAIsTrivial(aVecA.equalZero());
 444                 bool bBIsTrivial(aVecB.equalZero());
 445
 446                 // #i102241# prepare inverse edge length to normalize cross values;
 447                 // else the small compare value used in fTools::equalZero
 448                 // will be length dependent and this detection will work as less
 449                 // precise as longer the edge is. In principle, the length of the control
 450                 // vector would need to be used too, but to be trivial it is assumed to
 451                 // be of roughly equal length to the edge, so edge length can be used
 452                 // for both. Only needed when one of both is not trivial per se.
 453                 const double fInverseEdgeLength(bAIsTrivial && bBIsTrivial
 454                     ? 1.0
 455                     : 1.0 / aEdge.getLength());
 456
 457                 // if A is not zero, check if it could be
 458                 if(!bAIsTrivial)
 459                 {
 460                     // #i102241# parallel to edge? Check aVecA, aEdge. Use cross() which does what
 461                     // we need here with the precision we need
 462                     const double fCross(aVecA.cross(aEdge) * fInverseEdgeLength);
 463
 464                     if(fTools::equalZero(fCross))
 465                     {
 466                         // get scale to edge. Use bigger distance for numeric quality
 467                         const double fScale(fabs(aEdge.getX()) > fabs(aEdge.getY())
 468                             ? aVecA.getX() / aEdge.getX()
 469                             : aVecA.getY() / aEdge.getY());
 470
 471                         // relative end point of vector in edge range?
 472                         if(fTools::moreOrEqual(fScale, 0.0) && fTools::lessOrEqual(fScale, 1.0))
 473                         {
 474                             bAIsTrivial = true;
 475                         }
 476                     }
 477                 }
 478
 479                 // if B is not zero, check if it could be, but only if A is already trivial;
 480                 // else solve to trivial will not be possible for whole edge
 481                 if(bAIsTrivial && !bBIsTrivial)
 482                 {
 483                     // parallel to edge? Check aVecB, aEdge
 484                     const double fCross(aVecB.cross(aEdge) * fInverseEdgeLength);
 485
 486                     if(fTools::equalZero(fCross))
 487                     {
 488                         // get scale to edge. Use bigger distance for numeric quality
 489                         const double fScale(fabs(aEdge.getX()) > fabs(aEdge.getY())
 490                             ? aVecB.getX() / aEdge.getX()
 491                             : aVecB.getY() / aEdge.getY());
 492
 493                         // end point of vector in edge range? Caution: controlB is directed AGAINST edge
 494                         if(fTools::lessOrEqual(fScale, 0.0) && fTools::moreOrEqual(fScale, -1.0))
 495                         {
 496                             bBIsTrivial = true;
 497                         }
 498                     }
 499                 }
 500
 501                 // if both are/can be reduced, do it.
 502                 // Not possible if only one is/can be reduced (!)
 503                 if(bAIsTrivial && bBIsTrivial)
 504                 {
 505                     maControlPointA = maStartPoint;
 506                     maControlPointB = maEndPoint;
 507                 }
 508             }
 509         }
 510     }
 511
 512     namespace {
 513         double impGetLength(const B2DCubicBezier& rEdge, double fDeviation, sal_uInt32 nRecursionWatch)
 514         {
 515             const double fEdgeLength(rEdge.getEdgeLength());
 516             const double fControlPolygonLength(rEdge.getControlPolygonLength());
 517             const double fCurrentDeviation(fTools::equalZero(fControlPolygonLength) ? 0.0 : 1.0 - (fEdgeLength / fControlPolygonLength));
 518
 519             if(!nRecursionWatch || fTools:: lessOrEqual(fCurrentDeviation, fDeviation))
 520             {
 521                 return (fEdgeLength + fControlPolygonLength) * 0.5;
 522             }
 523             else
 524             {
 525                 B2DCubicBezier aLeft, aRight;
 526                 const double fNewDeviation(fDeviation * 0.5);
 527                 const sal_uInt32 nNewRecursionWatch(nRecursionWatch - 1);
 528
 529                 rEdge.split(0.5, &aLeft, &aRight);
 530
 531                 return impGetLength(aLeft, fNewDeviation, nNewRecursionWatch)
 532                     + impGetLength(aRight, fNewDeviation, nNewRecursionWatch);
 533             }
 534         }
 535     }
 536
 537     double B2DCubicBezier::getLength(double fDeviation) const
 538     {
 539         if(isBezier())
 540         {
 541             if(fDeviation < 0.00000001)
 542             {
 543                 fDeviation = 0.00000001;
 544             }
 545
 546             return impGetLength(*this, fDeviation, 6);
 547         }
 548         else
 549         {
 550             return B2DVector(getEndPoint() - getStartPoint()).getLength();
 551         }
 552     }
 553
 554     double B2DCubicBezier::getEdgeLength() const
 555     {
 556         const B2DVector aEdge(maEndPoint - maStartPoint);
 557         return aEdge.getLength();
 558     }
 559
 560     double B2DCubicBezier::getControlPolygonLength() const
 561     {
 562         const B2DVector aVectorA(maControlPointA - maStartPoint);
 563         const B2DVector aVectorB(maEndPoint - maControlPointB);
 564
 565         if(!aVectorA.equalZero() || !aVectorB.equalZero())
 566         {
 567             const B2DVector aTop(maControlPointB - maControlPointA);
 568             return (aVectorA.getLength() + aVectorB.getLength() + aTop.getLength());
 569         }
 570         else
 571         {
 572             return getEdgeLength();
 573         }
 574     }
 575
 576     void B2DCubicBezier::adaptiveSubdivideByAngle(B2DPolygon& rTarget, double fAngleBound, bool bAllowUnsharpen) const
 577     {
 578         if(isBezier())
 579         {
 580             // use support method #i37443# and allow unsharpen the criteria
 581             ImpSubDivAngleStart(maStartPoint, maControlPointA, maControlPointB, maEndPoint, rTarget, fAngleBound * F_PI180, bAllowUnsharpen);
 582         }
 583         else
 584         {
 585             rTarget.append(getEndPoint());
 586         }
 587     }
 588
 589     B2DVector B2DCubicBezier::getTangent(double t) const
 590     {
 591         if(fTools::lessOrEqual(t, 0.0))
 592         {
 593             // tangent in start point
 594             B2DVector aTangent(getControlPointA() - getStartPoint());
 595
 596             if(!aTangent.equalZero())
 597             {
 598                 return aTangent;
 599             }
 600
 601             // start point and control vector are the same, fallback
 602             // to implicit start vector to control point B
 603             aTangent = (getControlPointB() - getStartPoint()) * 0.3;
 604
 605             if(!aTangent.equalZero())
 606             {
 607                 return aTangent;
 608             }
 609
 610             // not a bezier at all, return edge vector
 611             return (getEndPoint() - getStartPoint()) * 0.3;
 612         }
 613         else if(fTools::moreOrEqual(t, 1.0))
 614         {
 615             // tangent in end point
 616             B2DVector aTangent(getEndPoint() - getControlPointB());
 617
 618             if(!aTangent.equalZero())
 619             {
 620                 return aTangent;
 621             }
 622
 623             // end point and control vector are the same, fallback
 624             // to implicit start vector from control point A
 625             aTangent = (getEndPoint() - getControlPointA()) * 0.3;
 626
 627             if(!aTangent.equalZero())
 628             {
 629                 return aTangent;
 630             }
 631
 632             // not a bezier at all, return edge vector
 633             return (getEndPoint() - getStartPoint()) * 0.3;
 634         }
 635         else
 636         {
 637             // t is in ]0.0 .. 1.0[. Split and extract
 638             B2DCubicBezier aRight;
 639             split(t, 0, &aRight);
 640
 641             return aRight.getControlPointA() - aRight.getStartPoint();
 642         }
 643     }
 644
 645     // #i37443# adaptive subdivide by nCount subdivisions
 646     void B2DCubicBezier::adaptiveSubdivideByCount(B2DPolygon& rTarget, sal_uInt32 nCount) const
 647     {
 648         const double fLenFact(1.0 / static_cast< double >(nCount + 1));
 649
 650         for(sal_uInt32 a(1); a <= nCount; a++)
 651         {
 652             const double fPos(static_cast< double >(a) * fLenFact);
 653             rTarget.append(interpolatePoint(fPos));
 654         }
 655
 656         rTarget.append(getEndPoint());
 657     }
 658
 659     // adaptive subdivide by distance
 660     void B2DCubicBezier::adaptiveSubdivideByDistance(B2DPolygon& rTarget, double fDistanceBound) const
 661     {
 662         if(isBezier())
 663         {
 664             ImpSubDivDistance(maStartPoint, maControlPointA, maControlPointB, maEndPoint, rTarget,
 665                 fDistanceBound * fDistanceBound, ::std::numeric_limits<double>::max(), 30);
 666         }
 667         else
 668         {
 669             rTarget.append(getEndPoint());
 670         }
 671     }
 672
 673     B2DPoint B2DCubicBezier::interpolatePoint(double t) const
 674     {
 675         OSL_ENSURE(t >= 0.0 && t <= 1.0, "B2DCubicBezier::interpolatePoint: Access out of range (!)");
 676
 677         if(isBezier())
 678         {
 679             const B2DPoint aS1L(interpolate(maStartPoint, maControlPointA, t));
 680             const B2DPoint aS1C(interpolate(maControlPointA, maControlPointB, t));
 681             const B2DPoint aS1R(interpolate(maControlPointB, maEndPoint, t));
 682             const B2DPoint aS2L(interpolate(aS1L, aS1C, t));
 683             const B2DPoint aS2R(interpolate(aS1C, aS1R, t));
 684
 685             return interpolate(aS2L, aS2R, t);
 686         }
 687         else
 688         {
 689             return interpolate(maStartPoint, maEndPoint, t);
 690         }
 691     }
 692
 693     double B2DCubicBezier::getSmallestDistancePointToBezierSegment(const B2DPoint& rTestPoint, double& rCut) const
 694     {
 695         const sal_uInt32 nInitialDivisions(3L);
 696         B2DPolygon aInitialPolygon;
 697
 698         // as start make a fix division, creates nInitialDivisions + 2L points
 699         aInitialPolygon.append(getStartPoint());
 700         adaptiveSubdivideByCount(aInitialPolygon, nInitialDivisions);
 701
 702         // now look for the closest point
 703         const sal_uInt32 nPointCount(aInitialPolygon.count());
 704         B2DVector aVector(rTestPoint - aInitialPolygon.getB2DPoint(0L));
 705         double fQuadDist(aVector.getX() * aVector.getX() + aVector.getY() * aVector.getY());
 706         double fNewQuadDist;
 707         sal_uInt32 nSmallestIndex(0L);
 708
 709         for(sal_uInt32 a(1L); a < nPointCount; a++)
 710         {
 711             aVector = B2DVector(rTestPoint - aInitialPolygon.getB2DPoint(a));
 712             fNewQuadDist = aVector.getX() * aVector.getX() + aVector.getY() * aVector.getY();
 713
 714             if(fNewQuadDist < fQuadDist)
 715             {
 716                 fQuadDist = fNewQuadDist;
 717                 nSmallestIndex = a;
 718             }
 719         }
 720
 721         // look right and left for even smaller distances
 722         double fStepValue(1.0 / (double)((nPointCount - 1L) * 2L)); // half the edge step width
 723         double fPosition((double)nSmallestIndex / (double)(nPointCount - 1L));
 724         bool bDone(false);
 725
 726         while(!bDone)
 727         {
 728             if(!bDone)
 729             {
 730                 // test left
 731                 double fPosLeft(fPosition - fStepValue);
 732
 733                 if(fPosLeft < 0.0)
 734                 {
 735                     fPosLeft = 0.0;
 736                     aVector = B2DVector(rTestPoint - maStartPoint);
 737                 }
 738                 else
 739                 {
 740                     aVector = B2DVector(rTestPoint - interpolatePoint(fPosLeft));
 741                 }
 742
 743                 fNewQuadDist = aVector.getX() * aVector.getX() + aVector.getY() * aVector.getY();
 744
 745                 if(fTools::less(fNewQuadDist, fQuadDist))
 746                 {
 747                     fQuadDist = fNewQuadDist;
 748                     fPosition = fPosLeft;
 749                 }
 750                 else
 751                 {
 752                     // test right
 753                     double fPosRight(fPosition + fStepValue);
 754
 755                     if(fPosRight > 1.0)
 756                     {
 757                         fPosRight = 1.0;
 758                         aVector = B2DVector(rTestPoint - maEndPoint);
 759                     }
 760                     else
 761                     {
 762                         aVector = B2DVector(rTestPoint - interpolatePoint(fPosRight));
 763                     }
 764
 765                     fNewQuadDist = aVector.getX() * aVector.getX() + aVector.getY() * aVector.getY();
 766
 767                     if(fTools::less(fNewQuadDist, fQuadDist))
 768                     {
 769                         fQuadDist = fNewQuadDist;
 770                         fPosition = fPosRight;
 771                     }
 772                     else
 773                     {
 774                         // not less left or right, done
 775                         bDone = true;
 776                     }
 777                 }
 778             }
 779
 780             if(0.0 == fPosition || 1.0 == fPosition)
 781             {
 782                 // if we are completely left or right, we are done
 783                 bDone = true;
 784             }
 785
 786             if(!bDone)
 787             {
 788                 // prepare next step value
 789                 fStepValue /= 2.0;
 790             }
 791         }
 792
 793         rCut = fPosition;
 794         return sqrt(fQuadDist);
 795     }
 796
 797     void B2DCubicBezier::split(double t, B2DCubicBezier* pBezierA, B2DCubicBezier* pBezierB) const
 798     {
 799         OSL_ENSURE(t >= 0.0 && t <= 1.0, "B2DCubicBezier::split: Access out of range (!)");
 800
 801         if(!pBezierA && !pBezierB)
 802         {
 803             return;
 804         }
 805
 806         if(isBezier())
 807         {
 808             const B2DPoint aS1L(interpolate(maStartPoint, maControlPointA, t));
 809             const B2DPoint aS1C(interpolate(maControlPointA, maControlPointB, t));
 810             const B2DPoint aS1R(interpolate(maControlPointB, maEndPoint, t));
 811             const B2DPoint aS2L(interpolate(aS1L, aS1C, t));
 812             const B2DPoint aS2R(interpolate(aS1C, aS1R, t));
 813             const B2DPoint aS3C(interpolate(aS2L, aS2R, t));
 814
 815             if(pBezierA)
 816             {
 817                 pBezierA->setStartPoint(maStartPoint);
 818                 pBezierA->setEndPoint(aS3C);
 819                 pBezierA->setControlPointA(aS1L);
 820                 pBezierA->setControlPointB(aS2L);
 821             }
 822
 823             if(pBezierB)
 824             {
 825                 pBezierB->setStartPoint(aS3C);
 826                 pBezierB->setEndPoint(maEndPoint);
 827                 pBezierB->setControlPointA(aS2R);
 828                 pBezierB->setControlPointB(aS1R);
 829             }
 830         }
 831         else
 832         {
 833             const B2DPoint aSplit(interpolate(maStartPoint, maEndPoint, t));
 834
 835             if(pBezierA)
 836             {
 837                 pBezierA->setStartPoint(maStartPoint);
 838                 pBezierA->setEndPoint(aSplit);
 839                 pBezierA->setControlPointA(maStartPoint);
 840                 pBezierA->setControlPointB(aSplit);
 841             }
 842
 843             if(pBezierB)
 844             {
 845                 pBezierB->setStartPoint(aSplit);
 846                 pBezierB->setEndPoint(maEndPoint);
 847                 pBezierB->setControlPointA(aSplit);
 848                 pBezierB->setControlPointB(maEndPoint);
 849             }
 850         }
 851     }
 852
 853     B2DCubicBezier B2DCubicBezier::snippet(double fStart, double fEnd) const
 854     {
 855         B2DCubicBezier aRetval;
 856
 857         if(fTools::more(fStart, 1.0))
 858         {
 859             fStart = 1.0;
 860         }
 861         else if(fTools::less(fStart, 0.0))
 862         {
 863             fStart = 0.0;
 864         }
 865
 866         if(fTools::more(fEnd, 1.0))
 867         {
 868             fEnd = 1.0;
 869         }
 870         else if(fTools::less(fEnd, 0.0))
 871         {
 872             fEnd = 0.0;
 873         }
 874
 875         if(fEnd <= fStart)
 876         {
 877             // empty or NULL, create single point at center
 878             const double fSplit((fEnd + fStart) * 0.5);
 879             const B2DPoint aPoint(interpolate(getStartPoint(), getEndPoint(), fSplit));
 880             aRetval.setStartPoint(aPoint);
 881             aRetval.setEndPoint(aPoint);
 882             aRetval.setControlPointA(aPoint);
 883             aRetval.setControlPointB(aPoint);
 884         }
 885         else
 886         {
 887             if(isBezier())
 888             {
 889                 // copy bezier; cut off right, then cut off left. Do not forget to
 890                 // adapt cut value when both cuts happen
 891                 const bool bEndIsOne(fTools::equal(fEnd, 1.0));
 892                 const bool bStartIsZero(fTools::equalZero(fStart));
 893                 aRetval = *this;
 894
 895                 if(!bEndIsOne)
 896                 {
 897                     aRetval.split(fEnd, &aRetval, 0);
 898
 899                     if(!bStartIsZero)
 900                     {
 901                         fStart /= fEnd;
 902                     }
 903                 }
 904
 905                 if(!bStartIsZero)
 906                 {
 907                     aRetval.split(fStart, 0, &aRetval);
 908                 }
 909             }
 910             else
 911             {
 912                 // no bezier, create simple edge
 913                 const B2DPoint aPointA(interpolate(getStartPoint(), getEndPoint(), fStart));
 914                 const B2DPoint aPointB(interpolate(getStartPoint(), getEndPoint(), fEnd));
 915                 aRetval.setStartPoint(aPointA);
 916                 aRetval.setEndPoint(aPointB);
 917                 aRetval.setControlPointA(aPointA);
 918                 aRetval.setControlPointB(aPointB);
 919             }
 920         }
 921
 922         return aRetval;
 923     }
 924
 925     B2DRange B2DCubicBezier::getRange() const
 926     {
 927         B2DRange aRetval(maStartPoint, maEndPoint);
 928
 929         aRetval.expand(maControlPointA);
 930         aRetval.expand(maControlPointB);
 931
 932         return aRetval;
 933     }
 934
 935     bool B2DCubicBezier::getMinimumExtremumPosition(double& rfResult) const
 936     {
 937         ::std::vector< double > aAllResults;
 938
 939         aAllResults.reserve(4);
 940         getAllExtremumPositions(aAllResults);
 941
 942         const sal_uInt32 nCount(aAllResults.size());
 943
 944         if(!nCount)
 945         {
 946             return false;
 947         }
 948         else if(1 == nCount)
 949         {
 950             rfResult = aAllResults[0];
 951             return true;
 952         }
 953         else
 954         {
 955             rfResult = *(::std::min_element(aAllResults.begin(), aAllResults.end()));
 956             return true;
 957         }
 958     }
 959
 960     namespace
 961     {
 962         inline void impCheckExtremumResult(double fCandidate, ::std::vector< double >& rResult)
 963         {
 964             // check for range ]0.0 .. 1.0[ with excluding 1.0 and 0.0 clearly
 965             // by using the equalZero test, NOT ::more or ::less which will use the
 966             // ApproxEqual() which is too exact here
 967             if(fCandidate > 0.0 && !fTools::equalZero(fCandidate))
 968             {
 969                 if(fCandidate < 1.0 && !fTools::equalZero(fCandidate - 1.0))
 970                 {
 971                     rResult.push_back(fCandidate);
 972                 }
 973             }
 974         }
 975     }
 976
 977     void B2DCubicBezier::getAllExtremumPositions(::std::vector< double >& rResults) const
 978     {
 979         rResults.clear();
 980
 981         // calculate the x-extrema parameters by zeroing first x-derivative
 982         // of the cubic bezier's parametric formula, which results in a
 983         // quadratic equation: dBezier/dt = t*t*fAX - 2*t*fBX + fCX
 984         const B2DPoint aRelativeEndPoint(maEndPoint-maStartPoint);
 985         const double fAX = 3 * (maControlPointA.getX() - maControlPointB.getX()) + aRelativeEndPoint.getX();
 986         const double fBX = 2 * maControlPointA.getX() - maControlPointB.getX() - maStartPoint.getX();
 987         double fCX(maControlPointA.getX() - maStartPoint.getX());
 988
 989         if(fTools::equalZero(fCX))
 990         {
 991             // detect fCX equal zero and truncate to real zero value in that case
 992             fCX = 0.0;
 993         }
 994
 995         if( !fTools::equalZero(fAX) )
 996         {
 997             // derivative is polynomial of order 2 => use binomial formula
 998             const double fD = fBX*fBX - fAX*fCX;
 999             if( fD >= 0.0 )
1000             {
1001                 const double fS = sqrt(fD);
1002                 // same as above but for very small fAX and/or fCX
1003                 // this has much better numerical stability
1004                 // see NRC chapter 5-6 (thanks THB!)
1005                 const double fQ = fBX + ((fBX >= 0) ? +fS : -fS);
1006                 impCheckExtremumResult(fQ / fAX, rResults);
1007                 impCheckExtremumResult(fCX / fQ, rResults);
1008             }
1009         }
1010         else if( !fTools::equalZero(fBX) )
1011         {
1012             // derivative is polynomial of order 1 => one extrema
1013             impCheckExtremumResult(fCX / (2 * fBX), rResults);
1014         }
1015
1016         // calculate the y-extrema parameters by zeroing first y-derivative
1017         const double fAY = 3 * (maControlPointA.getY() - maControlPointB.getY()) + aRelativeEndPoint.getY();
1018         const double fBY = 2 * maControlPointA.getY() - maControlPointB.getY() - maStartPoint.getY();
1019         double fCY(maControlPointA.getY() - maStartPoint.getY());
1020
1021         if(fTools::equalZero(fCY))
1022         {
1023             // detect fCY equal zero and truncate to real zero value in that case
1024             fCY = 0.0;
1025         }
1026
1027         if( !fTools::equalZero(fAY) )
1028         {
1029             // derivative is polynomial of order 2 => use binomial formula
1030             const double fD = fBY*fBY - fAY*fCY;
1031             if( fD >= 0 )
1032             {
1033                 const double fS = sqrt(fD);
1034                 // same as above but for very small fAX and/or fCX
1035                 // this has much better numerical stability
1036                 // see NRC chapter 5-6 (thanks THB!)
1037                 const double fQ = fBY + ((fBY >= 0) ? +fS : -fS);
1038                 impCheckExtremumResult(fQ / fAY, rResults);
1039                 impCheckExtremumResult(fCY / fQ, rResults);
1040             }
1041         }
1042         else if( !fTools::equalZero(fBY) )
1043         {
1044             // derivative is polynomial of order 1 => one extrema
1045             impCheckExtremumResult(fCY / (2 * fBY), rResults);
1046         }
1047     }
1048 } // end of namespace basegfx
1049
1050 // eof