src/FrontTrack/FrontTrack_Projector.cxx

   1 // File      : FrontTrack_Projector.cxx
   2 // Created   : Wed Apr 26 20:33:55 2017
   3 // Author    : Edward AGAPOV (eap)
   4
   5 #include "FrontTrack_Projector.hxx"
   6
   7 #include <BRepAdaptor_Curve.hxx>
   8 #include <BRepBndLib.hxx>
   9 #include <BRepTopAdaptor_FClass2d.hxx>
  10 #include <BRep_Tool.hxx>
  11 #include <ElCLib.hxx>
  12 #include <ElSLib.hxx>
  13 #include <GCPnts_UniformDeflection.hxx>
  14 #include <GeomAdaptor_Curve.hxx>
  15 #include <GeomLib_IsPlanarSurface.hxx>
  16 #include <ShapeAnalysis_Curve.hxx>
  17 #include <ShapeAnalysis_Surface.hxx>
  18 #include <TopExp.hxx>
  19 #include <TopoDS.hxx>
  20 #include <TopoDS_Edge.hxx>
  21 #include <TopoDS_Face.hxx>
  22 #include <TopoDS_Vertex.hxx>
  23 #include <gp_Circ.hxx>
  24 #include <gp_Cylinder.hxx>
  25 #include <gp_Dir.hxx>
  26 #include <gp_Pln.hxx>
  27 #include <gp_Pnt.hxx>
  28 #include <gp_Sphere.hxx>
  29 #include <gp_Vec.hxx>
  30
  31 #include <limits>
  32
  33 //-----------------------------------------------------------------------------
  34 /*!
  35  * \brief Root class of a projector of a point to a boundary shape
  36  */
  37 struct FT_RealProjector
  38 {
  39   virtual ~FT_RealProjector() {}
  40
  41   /*!
  42    * \brief Project a point to a boundary shape
  43    *  \param [in] point - the point to project
  44    *  \param [out] newSolution - position on the shape (U or UV) found during the projection
  45    *  \param [in] prevSolution - position already found during the projection of a neighbor point
  46    *  \return gp_Pnt - the projection point
  47    */
  48   virtual gp_Pnt project( const gp_Pnt& point,
  49                           double*       newSolution,
  50                           const double* prevSolution = 0) = 0;
  51
  52   /*!
  53    * \brief Project a point to a boundary shape and check if the projection is within
  54    *        the shape boundary
  55    *  \param [in] point - the point to project
  56    *  \param [in] maxDist2 - the maximal allowed square distance between point and projection
  57    *  \param [out] projection - the projection point
  58    *  \param [out] newSolution - position on the shape (U or UV) found during the projection
  59    *  \param [in] prevSolution - position already found during the projection of a neighbor point
  60    *  \return bool - false if the projection point lies out of the shape boundary or
  61    the distance the point and the projection is more than sqrt(maxDist2)
  62   */
  63   virtual bool projectAndClassify( const gp_Pnt& point,
  64                                    const double  maxDist2,
  65                                    gp_Pnt&       projection,
  66                                    double*       newSolution,
  67                                    const double* prevSolution = 0) = 0;
  68
  69   // return true if a previously found solution can be used to speed up the projection
  70
  71   virtual bool canUsePrevSolution() const { return false; }
  72
  73
  74   double _dist; // distance between the point being projected and its projection
  75 };
  76
  77 namespace // actual projection algorithms
  78 {
  79   const double theEPS = 1e-12;
  80
  81   //================================================================================
  82   /*!
  83    * \brief Projector to any curve
  84    */
  85   //================================================================================
  86
  87   struct CurveProjector : public FT_RealProjector
  88   {
  89     BRepAdaptor_Curve   _curve;
  90     double              _tol;
  91     ShapeAnalysis_Curve _projector;
  92     double              _uRange[2];
  93
  94     //-----------------------------------------------------------------------------
  95     CurveProjector( const TopoDS_Edge& e, const double tol ):
  96       _curve( e ), _tol( tol )
  97     {
  98       BRep_Tool::Range( e, _uRange[0], _uRange[1] );
  99     }
 100
 101     //-----------------------------------------------------------------------------
 102     // project a point to the curve
 103     virtual gp_Pnt project( const gp_Pnt& P,
 104                             double*       newSolution,
 105                             const double* prevSolution = 0)
 106     {
 107 #ifdef _DEBUG_
 108     std::cout << ".. project a point to the curve prevSolution = " << prevSolution << std::endl;
 109 #endif
 110       gp_Pnt         proj;
 111       Standard_Real param;
 112
 113       if ( prevSolution )
 114       {
 115         _dist = _projector.NextProject( prevSolution[0], _curve, P, _tol, proj, param );
 116       }
 117       else
 118       {
 119         _dist = _projector.Project( _curve, P, _tol, proj, param, false );
 120       }
 121 #ifdef _DEBUG_
 122     std::cout << "..    _dist : " << _dist << std::endl;
 123 #endif
 124       proj = _curve.Value( param );
 125
 126       newSolution[0] = param;
 127
 128       return proj;
 129     }
 130
 131     //-----------------------------------------------------------------------------
 132     // project a point to a curve and check if the projection is within the curve boundary
 133     virtual bool projectAndClassify( const gp_Pnt& point,
 134                                      const double  maxDist2,
 135                                      gp_Pnt&       projection,
 136                                      double*       newSolution,
 137                                      const double* prevSolution = 0)
 138     {
 139 #ifdef _DEBUG_
 140     std::cout << ".. project a point to a curve and check " << std::endl;
 141 #endif
 142       projection = project( point, newSolution, prevSolution );
 143       return ( _uRange[0] < newSolution[0] && newSolution[0] < _uRange[1] &&
 144                _dist * _dist < maxDist2 );
 145     }
 146
 147     //-----------------------------------------------------------------------------
 148     // return true if a previously found solution can be used to speed up the projection
 149     virtual bool canUsePrevSolution() const { return true; }
 150   };
 151
 152   //================================================================================
 153   /*!
 154    * \brief Projector to a straight curve. Don't project, classify only
 155    */
 156   //================================================================================
 157
 158   struct LineProjector : public FT_RealProjector
 159   {
 160     gp_Pnt _p0, _p1;
 161
 162     //-----------------------------------------------------------------------------
 163     LineProjector( TopoDS_Edge e )
 164     {
 165       e.Orientation( TopAbs_FORWARD );
 166       _p0 = BRep_Tool::Pnt( TopExp::FirstVertex( e ));
 167       _p1 = BRep_Tool::Pnt( TopExp::LastVertex ( e ));
 168     }
 169
 170     //-----------------------------------------------------------------------------
 171     // does nothing
 172     virtual gp_Pnt project( const gp_Pnt& P,
 173                             double*       newSolution,
 174                             const double* prevSolution = 0)
 175     {
 176       return P;
 177     }
 178     //-----------------------------------------------------------------------------
 179     // check if a point lies within the line segment
 180     virtual bool projectAndClassify( const gp_Pnt& point,
 181                                      const double  maxDist2,
 182                                      gp_Pnt&       projection,
 183                                      double*       newSolution,
 184                                      const double* prevSolution = 0)
 185     {
 186       gp_Vec edge( _p0, _p1 );
 187       gp_Vec p0p ( _p0, point  );
 188       double u = ( edge * p0p ) / edge.SquareMagnitude();  // param [0,1] on the edge
 189       projection = ( 1. - u ) * _p0.XYZ() + u * _p1.XYZ(); // projection of the point on the edge
 190       if ( u < 0 || 1 < u )
 191         return false;
 192
 193       // check distance
 194       return point.SquareDistance( projection ) < theEPS * theEPS;
 195     }
 196   };
 197
 198   //================================================================================
 199   /*!
 200    * \brief Projector to a circular edge
 201    */
 202   //================================================================================
 203
 204   struct CircleProjector : public FT_RealProjector
 205   {
 206     gp_Circ _circle;
 207     double _uRange[2];
 208
 209     //-----------------------------------------------------------------------------
 210     CircleProjector( const gp_Circ& c, const double f, const double l ):
 211       _circle( c )
 212     {
 213       _uRange[0] = f;
 214       _uRange[1] = l;
 215     }
 216
 217     //-----------------------------------------------------------------------------
 218     // project a point to the circle
 219     virtual gp_Pnt project( const gp_Pnt& P,
 220                             double*       newSolution,
 221                             const double* prevSolution = 0)
 222     {
 223       // assume that P is already on the the plane of circle, since
 224       // it is in the middle of two points lying on the circle
 225
 226       // move P to the circle
 227       const gp_Pnt& O = _circle.Location();
 228       gp_Vec radiusVec( O, P );
 229       double radius = radiusVec.Magnitude();
 230       if ( radius < std::numeric_limits<double>::min() )
 231         return P; // P in on the axe
 232
 233       gp_Pnt proj = O.Translated( radiusVec.Multiplied( _circle.Radius() / radius ));
 234
 235       _dist = _circle.Radius() - radius;
 236
 237       return proj;
 238     }
 239
 240     //-----------------------------------------------------------------------------
 241     // project and check if a projection lies within the circular edge
 242     virtual bool projectAndClassify( const gp_Pnt& point,
 243                                      const double  maxDist2,
 244                                      gp_Pnt&       projection,
 245                                      double*       newSolution,
 246                                      const double* prevSolution = 0)
 247     {
 248       _dist = -1;
 249       projection = project( point, newSolution );
 250       if ( _dist < 0 || // ?
 251            _dist * _dist > maxDist2 )
 252         return false;
 253
 254       newSolution[0] = ElCLib::Parameter( _circle, projection );
 255       return ( _uRange[0] < newSolution[0] && newSolution[0] < _uRange[1] );
 256     }
 257   };
 258
 259   //================================================================================
 260   /*!
 261    * \brief Projector to any surface
 262    */
 263   //================================================================================
 264
 265   struct SurfaceProjector : public FT_RealProjector
 266   {
 267     ShapeAnalysis_Surface    _projector;
 268     double                   _tol;
 269     BRepTopAdaptor_FClass2d* _classifier;
 270
 271     //-----------------------------------------------------------------------------
 272     SurfaceProjector( const TopoDS_Face& face, const double tol, BRepTopAdaptor_FClass2d* cls ):
 273       _projector( BRep_Tool::Surface( face )),
 274       _tol( tol ),
 275       _classifier( cls )
 276     {
 277     }
 278     //-----------------------------------------------------------------------------
 279     // delete _classifier
 280     ~SurfaceProjector()
 281     {
 282       delete _classifier;
 283     }
 284
 285     //-----------------------------------------------------------------------------
 286     // project a point to a surface
 287     virtual gp_Pnt project( const gp_Pnt& P,
 288                             double*       newSolution,
 289                             const double* prevSolution = 0)
 290     {
 291       gp_Pnt2d uv;
 292
 293       if ( prevSolution )
 294       {
 295         gp_Pnt2d prevUV( prevSolution[0], prevSolution[1] );
 296         uv = _projector.NextValueOfUV( prevUV, P, _tol );
 297       }
 298       else
 299       {
 300         uv = _projector.ValueOfUV( P, _tol );
 301       }
 302
 303       uv.Coord( newSolution[0], newSolution[1] );
 304
 305       gp_Pnt proj = _projector.Value( uv );
 306
 307       _dist = _projector.Gap();
 308
 309       return proj;
 310     }
 311
 312     //-----------------------------------------------------------------------------
 313     // project a point to a surface and check if the projection is within the surface boundary
 314     virtual bool projectAndClassify( const gp_Pnt& point,
 315                                      const double  maxDist2,
 316                                      gp_Pnt&       projection,
 317                                      double*       newSolution,
 318                                      const double* prevSolution = 0)
 319     {
 320       projection = project( point, newSolution, prevSolution );
 321       return ( _dist * _dist < maxDist2 )  &&  classify( newSolution );
 322     }
 323
 324     //-----------------------------------------------------------------------------
 325     // check if the projection is within the shape boundary
 326     bool classify( const double* newSolution )
 327     {
 328       TopAbs_State state = _classifier->Perform( gp_Pnt2d( newSolution[0], newSolution[1]) );
 329       return ( state != TopAbs_OUT );
 330     }
 331
 332     //-----------------------------------------------------------------------------
 333     // return true if a previously found solution can be used to speed up the projection
 334     virtual bool canUsePrevSolution() const { return true; }
 335   };
 336
 337   //================================================================================
 338   /*!
 339    * \brief Projector to a plane. Don't project, classify only
 340    */
 341   //================================================================================
 342
 343   struct PlaneProjector : public SurfaceProjector
 344   {
 345     gp_Pln _plane;
 346     bool   _isRealPlane; // false means that a surface is planar but parametrization is different
 347
 348     //-----------------------------------------------------------------------------
 349     PlaneProjector( const gp_Pln&            pln,
 350                     const TopoDS_Face&       face,
 351                     BRepTopAdaptor_FClass2d* cls,
 352                     bool                     isRealPlane=true):
 353       SurfaceProjector( face, 0, cls ),
 354       _plane( pln ),
 355       _isRealPlane( isRealPlane )
 356     {}
 357
 358     //-----------------------------------------------------------------------------
 359     // does nothing
 360     virtual gp_Pnt project( const gp_Pnt& P,
 361                             double*       newSolution,
 362                             const double* prevSolution = 0)
 363     {
 364       return P;
 365     }
 366     //-----------------------------------------------------------------------------
 367     // check if a point lies within the boundry of the planar face
 368     virtual bool projectAndClassify( const gp_Pnt& point,
 369                                      const double  maxDist2,
 370                                      gp_Pnt&       projection,
 371                                      double*       newSolution,
 372                                      const double* prevSolution = 0)
 373     {
 374       if ( _isRealPlane )
 375       {
 376         ElSLib::PlaneParameters( _plane.Position(), point, newSolution[0], newSolution[1]);
 377         projection = ElSLib::PlaneValue ( newSolution[0], newSolution[1], _plane.Position() );
 378         if ( projection.SquareDistance( point ) > theEPS * theEPS )
 379           return false;
 380
 381         return SurfaceProjector::classify( newSolution );
 382       }
 383       else
 384       {
 385         return SurfaceProjector::projectAndClassify( point, maxDist2, projection,
 386                                                      newSolution, prevSolution );
 387       }
 388     }
 389     //-----------------------------------------------------------------------------
 390     // return true if a previously found solution can be used to speed up the projection
 391     virtual bool canUsePrevSolution() const { return false; }
 392   };
 393
 394   //================================================================================
 395   /*!
 396    * \brief Projector to a cylinder
 397    */
 398   //================================================================================
 399
 400   struct CylinderProjector : public SurfaceProjector
 401   {
 402     gp_Cylinder _cylinder;
 403
 404     //-----------------------------------------------------------------------------
 405     CylinderProjector( const gp_Cylinder&       c,
 406                        const TopoDS_Face&       face,
 407                        BRepTopAdaptor_FClass2d* cls ):
 408       SurfaceProjector( face, 0, cls ),
 409       _cylinder( c )
 410     {}
 411
 412     //-----------------------------------------------------------------------------
 413     // project a point to the cylinder
 414     virtual gp_Pnt project( const gp_Pnt& P,
 415                             double*       newSolution,
 416                             const double* prevSolution = 0)
 417     {
 418       // project the point P to the cylinder axis -> Pp
 419       const gp_Pnt& O   = _cylinder.Position().Location();
 420       const gp_Dir& axe = _cylinder.Position().Direction();
 421       gp_Vec       trsl = gp_Vec( axe ).Multiplied( gp_Vec( O, P ).Dot( axe ));
 422       gp_Pnt       Pp   = O.Translated( trsl );
 423
 424       // move Pp to the cylinder
 425       gp_Vec radiusVec( Pp, P );
 426       double radius = radiusVec.Magnitude();
 427       if ( radius < std::numeric_limits<double>::min() )
 428         return P; // P in on the axe
 429
 430       gp_Pnt proj = Pp.Translated( radiusVec.Multiplied( _cylinder.Radius() / radius ));
 431
 432       _dist = _cylinder.Radius() - radius;
 433
 434       return proj;
 435     }
 436     //-----------------------------------------------------------------------------
 437     // project a point to the cylinder and check if the projection is within the surface boundary
 438     virtual bool projectAndClassify( const gp_Pnt& point,
 439                                      const double  maxDist2,
 440                                      gp_Pnt&       projection,
 441                                      double*       newSolution,
 442                                      const double* prevSolution = 0)
 443     {
 444       ElSLib::CylinderParameters( _cylinder.Position(), _cylinder.Radius(), point,
 445                                   newSolution[0], newSolution[1]);
 446       projection = ElSLib::CylinderValue( newSolution[0], newSolution[1],
 447                                           _cylinder.Position(), _cylinder.Radius() );
 448
 449       return ( _dist * _dist < maxDist2 )  &&  SurfaceProjector::classify( newSolution );
 450     }
 451     //-----------------------------------------------------------------------------
 452     // return true if a previously found solution can be used to speed up the projection
 453     virtual bool canUsePrevSolution() const { return false; }
 454   };
 455
 456   //================================================================================
 457   /*!
 458    * \brief Projector to a cone
 459    */
 460   //================================================================================
 461
 462   struct ConeProjector : public SurfaceProjector
 463   {
 464     gp_Cone _cone;
 465
 466     //-----------------------------------------------------------------------------
 467     ConeProjector( const gp_Cone&           c,
 468                    const TopoDS_Face&       face,
 469                    BRepTopAdaptor_FClass2d* cls ):
 470       SurfaceProjector( face, 0, cls ),
 471       _cone( c )
 472     {}
 473
 474     //-----------------------------------------------------------------------------
 475     // project a point to the cone
 476     virtual gp_Pnt project( const gp_Pnt& point,
 477                             double*       newSolution,
 478                             const double* prevSolution = 0)
 479     {
 480       ElSLib::ConeParameters( _cone.Position(), _cone.RefRadius(), _cone.SemiAngle(),
 481                               point, newSolution[0], newSolution[1]);
 482       gp_Pnt proj = ElSLib::ConeValue( newSolution[0], newSolution[1],
 483                                        _cone.Position(), _cone.RefRadius(), _cone.SemiAngle() );
 484       _dist = point.Distance( proj );
 485
 486       return proj;
 487     }
 488
 489     //-----------------------------------------------------------------------------
 490     // project a point to the cone and check if the projection is within the surface boundary
 491     virtual bool projectAndClassify( const gp_Pnt& point,
 492                                      const double  maxDist2,
 493                                      gp_Pnt&       projection,
 494                                      double*       newSolution,
 495                                      const double* prevSolution = 0)
 496     {
 497       projection = project( point, newSolution, prevSolution );
 498
 499       return ( _dist * _dist < maxDist2 )  &&  SurfaceProjector::classify( newSolution );
 500     }
 501     //-----------------------------------------------------------------------------
 502     // return true if a previously found solution can be used to speed up the projection
 503     virtual bool canUsePrevSolution() const { return false; }
 504   };
 505
 506   //================================================================================
 507   /*!
 508    * \brief Projector to a sphere
 509    */
 510   //================================================================================
 511
 512   struct SphereProjector : public SurfaceProjector
 513   {
 514     gp_Sphere _sphere;
 515
 516     //-----------------------------------------------------------------------------
 517     SphereProjector( const gp_Sphere&         s,
 518                      const TopoDS_Face&       face,
 519                      BRepTopAdaptor_FClass2d* cls ):
 520       SurfaceProjector( face, 0, cls ),
 521       _sphere( s )
 522     {}
 523
 524     //-----------------------------------------------------------------------------
 525     // project a point to the sphere
 526     virtual gp_Pnt project( const gp_Pnt& P,
 527                             double*       newSolution,
 528                             const double* prevSolution = 0)
 529     {
 530       // move Pp to the Sphere
 531       const gp_Pnt& O = _sphere.Location();
 532       gp_Vec radiusVec( O, P );
 533       double radius = radiusVec.Magnitude();
 534       if ( radius < std::numeric_limits<double>::min() )
 535         return P; // P is on O
 536
 537       gp_Pnt proj = O.Translated( radiusVec.Multiplied( _sphere.Radius() / radius ));
 538
 539       _dist = _sphere.Radius() - radius;
 540
 541       return proj;
 542     }
 543
 544     //-----------------------------------------------------------------------------
 545     // project a point to the sphere and check if the projection is within the surface boundary
 546     virtual bool projectAndClassify( const gp_Pnt& point,
 547                                      const double  maxDist2,
 548                                      gp_Pnt&       projection,
 549                                      double*       newSolution,
 550                                      const double* prevSolution = 0)
 551     {
 552       ElSLib::SphereParameters( _sphere.Position(), _sphere.Radius(), point,
 553                                   newSolution[0], newSolution[1]);
 554       projection = ElSLib::SphereValue( newSolution[0], newSolution[1],
 555                                         _sphere.Position(), _sphere.Radius() );
 556
 557       return ( _dist * _dist < maxDist2 )  &&  SurfaceProjector::classify( newSolution );
 558     }
 559     //-----------------------------------------------------------------------------
 560     // return true if a previously found solution can be used to speed up the projection
 561     virtual bool canUsePrevSolution() const { return false; }
 562   };
 563
 564   //================================================================================
 565   /*!
 566    * \brief Projector to a torus
 567    */
 568   //================================================================================
 569
 570   struct TorusProjector : public SurfaceProjector
 571   {
 572     gp_Torus _torus;
 573
 574     //-----------------------------------------------------------------------------
 575     TorusProjector( const gp_Torus&          t,
 576                     const TopoDS_Face&       face,
 577                     BRepTopAdaptor_FClass2d* cls ):
 578       SurfaceProjector( face, 0, cls ),
 579       _torus( t )
 580     {}
 581
 582     //-----------------------------------------------------------------------------
 583     // project a point to the torus
 584     virtual gp_Pnt project( const gp_Pnt& point,
 585                             double*       newSolution,
 586                             const double* prevSolution = 0)
 587     {
 588       ElSLib::TorusParameters( _torus.Position(), _torus.MajorRadius(), _torus.MinorRadius(),
 589                                point, newSolution[0], newSolution[1]);
 590       gp_Pnt proj = ElSLib::TorusValue( newSolution[0], newSolution[1],
 591                                         _torus.Position(), _torus.MajorRadius(), _torus.MinorRadius() );
 592       _dist = point.Distance( proj );
 593
 594       return proj;
 595     }
 596
 597     //-----------------------------------------------------------------------------
 598     // project a point to the torus and check if the projection is within the surface boundary
 599     virtual bool projectAndClassify( const gp_Pnt& point,
 600                                      const double  maxDist2,
 601                                      gp_Pnt&       projection,
 602                                      double*       newSolution,
 603                                      const double* prevSolution = 0)
 604     {
 605       projection = project( point, newSolution, prevSolution );
 606
 607       return ( _dist * _dist < maxDist2 )  &&  SurfaceProjector::classify( newSolution );
 608     }
 609     //-----------------------------------------------------------------------------
 610     // return true if a previously found solution can be used to speed up the projection
 611     virtual bool canUsePrevSolution() const { return false; }
 612   };
 613
 614   //================================================================================
 615   /*!
 616    * \brief Check if a curve can be considered straight
 617    */
 618   //================================================================================
 619
 620   bool isStraight( const GeomAdaptor_Curve& curve, const double tol )
 621   {
 622     // rough check: evaluate how far from a straight line connecting the curve ends
 623     // stand several internal points of the curve
 624
 625     const double  f = curve.FirstParameter();
 626     const double  l = curve.LastParameter();
 627     const gp_Pnt pf = curve.Value( f );
 628     const gp_Pnt pl = curve.Value( l );
 629     const gp_Vec lineVec( pf, pl );
 630     const double lineLen2 = lineVec.SquareMagnitude();
 631     if ( lineLen2 < std::numeric_limits< double >::min() )
 632       return false; // E seems closed
 633
 634     const double nbSamples = 7;
 635     for ( int i = 0; i < nbSamples; ++i )
 636     {
 637       const double  r = ( i + 1 ) / nbSamples;
 638       const gp_Pnt pi = curve.Value( f * r + l * ( 1 - r ));
 639       const gp_Vec vi( pf, pi );
 640       const double h2 = lineVec.Crossed( vi ).SquareMagnitude() / lineLen2;
 641       if ( h2 > tol * tol )
 642         return false;
 643     }
 644
 645     // thorough check
 646     GCPnts_UniformDeflection divider( curve, tol );
 647     return ( divider.IsDone() && divider.NbPoints() < 3 );
 648   }
 649 }
 650
 651 //================================================================================
 652 /*!
 653  * \brief Initialize with a boundary shape
 654  */
 655 //================================================================================
 656
 657 FT_Projector::FT_Projector(const TopoDS_Shape& shape)
 658 {
 659   _realProjector = 0;
 660   setBoundaryShape( shape );
 661   _tryWOPrevSolution = false;
 662 }
 663
 664 //================================================================================
 665 /*!
 666  * \brief Copy another projector
 667  */
 668 //================================================================================
 669
 670 FT_Projector::FT_Projector(const FT_Projector& other)
 671 {
 672   _realProjector = 0;
 673   _shape = other._shape;
 674   _bndBox = other._bndBox;
 675   _tryWOPrevSolution = false;
 676 }
 677
 678 //================================================================================
 679 /*!
 680  * \brief Destructor. Delete _realProjector
 681  */
 682 //================================================================================
 683
 684 FT_Projector::~FT_Projector()
 685 {
 686   delete _realProjector;
 687 }
 688
 689 //================================================================================
 690 /*!
 691  * \brief Initialize with a boundary shape. Compute the bounding box
 692  */
 693 //================================================================================
 694
 695 void FT_Projector::setBoundaryShape(const TopoDS_Shape& shape)
 696 {
 697   delete _realProjector; _realProjector = 0;
 698   _shape = shape;
 699   if ( shape.IsNull() )
 700     return;
 701
 702   BRepBndLib::Add( shape, _bndBox );
 703   _bndBox.Enlarge( 1e-5 * sqrt( _bndBox.SquareExtent() ));
 704 }
 705
 706 //================================================================================
 707 /*!
 708  * \brief Create a real projector
 709  */
 710 //================================================================================
 711
 712 void FT_Projector::prepareForProjection()
 713 {
 714   if ( _shape.IsNull() || _realProjector )
 715     return;
 716
 717   if ( _shape.ShapeType() == TopAbs_EDGE )
 718   {
 719     const TopoDS_Edge& edge = TopoDS::Edge( _shape );
 720
 721     double tol = 1e-6 * sqrt( _bndBox.SquareExtent() );
 722
 723     double f,l;
 724     Handle(Geom_Curve) curve = BRep_Tool::Curve( edge, f,l );
 725     if ( curve.IsNull() )
 726       return; // degenerated edge
 727
 728     GeomAdaptor_Curve acurve( curve, f, l );
 729     switch ( acurve.GetType() )
 730     {
 731     case GeomAbs_Line:
 732       _realProjector = new LineProjector( edge );
 733       break;
 734     case GeomAbs_Circle:
 735       _realProjector = new CircleProjector( acurve.Circle(), f, l );
 736       break;
 737     case GeomAbs_BezierCurve:
 738     case GeomAbs_BSplineCurve:
 739     case GeomAbs_OffsetCurve:
 740     case GeomAbs_OtherCurve:
 741       if ( isStraight( acurve, tol ))
 742       {
 743         _realProjector = new LineProjector( edge );
 744         break;
 745       }
 746     case GeomAbs_Ellipse:
 747     case GeomAbs_Hyperbola:
 748     case GeomAbs_Parabola:
 749       _realProjector = new CurveProjector( edge, tol );
 750     }
 751   }
 752   else if ( _shape.ShapeType() == TopAbs_FACE )
 753   {
 754     TopoDS_Face face = TopoDS::Face( _shape );
 755
 756     Handle(Geom_Surface) surface = BRep_Tool::Surface( face );
 757     if ( surface.IsNull() )
 758       return;
 759
 760     GeomAdaptor_Surface asurface( surface );
 761     Standard_Real tol   = BRep_Tool::Tolerance( face );
 762     Standard_Real toluv = Min( asurface.UResolution( tol ), asurface.VResolution( tol ));
 763     BRepTopAdaptor_FClass2d* classifier = new BRepTopAdaptor_FClass2d( face, toluv );
 764
 765     switch ( asurface.GetType() )
 766     {
 767     case GeomAbs_Plane:
 768       _realProjector = new PlaneProjector( asurface.Plane(), face, classifier );
 769       break;
 770     case GeomAbs_Cylinder:
 771       _realProjector = new CylinderProjector( asurface.Cylinder(), face, classifier );
 772       break;
 773     case GeomAbs_Sphere:
 774       _realProjector = new SphereProjector( asurface.Sphere(), face, classifier );
 775       break;
 776     case GeomAbs_Cone:
 777       _realProjector = new ConeProjector( asurface.Cone(), face, classifier );
 778       break;
 779     case GeomAbs_Torus:
 780       _realProjector = new TorusProjector( asurface.Torus(), face, classifier );
 781       break;
 782     case GeomAbs_BezierSurface:
 783     case GeomAbs_BSplineSurface:
 784     case GeomAbs_SurfaceOfRevolution:
 785     case GeomAbs_SurfaceOfExtrusion:
 786     case GeomAbs_OffsetSurface:
 787     case GeomAbs_OtherSurface:
 788       GeomLib_IsPlanarSurface isPlaneCheck( surface, tol );
 789       if ( isPlaneCheck.IsPlanar() )
 790       {
 791         _realProjector = new PlaneProjector( isPlaneCheck.Plan(), face, classifier,
 792                                              /*isRealPlane=*/false);
 793       }
 794       else
 795       {
 796         _realProjector = new SurfaceProjector( face, tol, classifier );
 797       }
 798       break;
 799     }
 800
 801     if ( !_realProjector )
 802       delete classifier;
 803   }
 804 }
 805
 806 //================================================================================
 807 /*!
 808  * \brief Return true if projection is not needed
 809  */
 810 //================================================================================
 811
 812 bool FT_Projector::isPlanarBoundary() const
 813 {
 814   return ( dynamic_cast< LineProjector*  >( _realProjector ) ||
 815            dynamic_cast< PlaneProjector* >( _realProjector ) );
 816 }
 817
 818 //================================================================================
 819 /*!
 820  * \brief Check if a point lies on the boundary shape
 821  *  \param [in] point - the point to check
 822  *  \param [in] tol2 - a square tolerance allowing to decide whether a point is on the shape
 823  *  \param [in] newSolution - position on the shape (U or UV) of the point found
 824  *         during projecting
 825  *  \param [in] prevSolution - position on the shape (U or UV) of a neighbor point
 826  *  \return bool - \c true if the point lies on the boundary shape
 827  *
 828  * This method is used to select a shape by checking if all neighbor nodes of a node to move
 829  * lie on a shape.
 830  */
 831 //================================================================================
 832
 833 bool FT_Projector::isOnShape( const gp_Pnt& point,
 834                               const double  tol2,
 835                               double*       newSolution,
 836                               const double* prevSolution)
 837 {
 838   if ( _bndBox.IsOut( point ) || !_realProjector )
 839     return false;
 840
 841   gp_Pnt proj;
 842   if ( isPlanarBoundary() )
 843     return projectAndClassify( point, tol2, proj, newSolution, prevSolution );
 844
 845   return project( point, tol2, proj, newSolution, prevSolution );
 846 }
 847
 848 //================================================================================
 849 /*!
 850  * \brief Project a point to the boundary shape
 851  *  \param [in] point - the point to project
 852  *  \param [in] maxDist2 - the maximal square distance between the point and the projection
 853  *  \param [out] projection - the projection
 854  *  \param [out] newSolution - position on the shape (U or UV) of the point found
 855  *         during projecting
 856  *  \param [in] prevSolution - already found position on the shape (U or UV) of a neighbor point
 857  *  \return bool - false if the distance between the point and the projection
 858  *         is more than sqrt(maxDist2)
 859  *
 860  * This method is used to project a node in the case where only one shape is found by name
 861  */
 862 //================================================================================
 863
 864 bool FT_Projector::project( const gp_Pnt& point,
 865                             const double  maxDist2,
 866                             gp_Pnt&       projection,
 867                             double*       newSolution,
 868                             const double* prevSolution)
 869 {
 870   if ( !_realProjector )
 871     return false;
 872
 873   _realProjector->_dist = 1e100;
 874   projection = _realProjector->project( point, newSolution, prevSolution );
 875
 876   bool ok = ( _realProjector->_dist * _realProjector->_dist < maxDist2 );
 877   if ( !ok && _tryWOPrevSolution && prevSolution )
 878   {
 879     projection = _realProjector->project( point, newSolution );
 880     ok = ( _realProjector->_dist * _realProjector->_dist < maxDist2 );
 881   }
 882   return ok;
 883 }
 884
 885 //================================================================================
 886 /*!
 887  * \brief Project a point to the boundary shape and check if the projection lies within
 888  *        the shape boundary
 889  *  \param [in] point - the point to project
 890  *  \param [in] maxDist2 - the maximal square distance between the point and the projection
 891  *  \param [out] projection - the projection
 892  *  \param [out] newSolution - position on the shape (U or UV) of the point found
 893  *         during projecting
 894  *  \param [in] prevSolution - already found position on the shape (U or UV) of a neighbor point
 895  *  \return bool - false if the projection point lies out of the shape boundary or
 896  *          the distance between the point and the projection is more than sqrt(maxDist2)
 897  *
 898  * This method is used to project a node in the case where several shapes are selected for
 899  * projection of a node group
 900  */
 901 //================================================================================
 902
 903 bool FT_Projector::projectAndClassify( const gp_Pnt& point,
 904                                        const double  maxDist2,
 905                                        gp_Pnt&       projection,
 906                                        double*       newSolution,
 907                                        const double* prevSolution)
 908 {
 909   if ( _bndBox.IsOut( point ) || !_realProjector )
 910     return false;
 911
 912   bool ok = _realProjector->projectAndClassify( point, maxDist2, projection,
 913                                                 newSolution, prevSolution );
 914   if ( !ok && _tryWOPrevSolution && prevSolution )
 915     ok = _realProjector->projectAndClassify( point, maxDist2, projection, newSolution );
 916
 917   return ok;
 918 }
 919
 920 //================================================================================
 921 /*!
 922  * \brief Return true if a previously found solution can be used to speed up the projection
 923  */
 924 //================================================================================
 925
 926 bool FT_Projector::canUsePrevSolution() const
 927 {
 928   return ( _realProjector && _realProjector->canUsePrevSolution() );
 929 }