src/StdPrs/StdPrs_Isolines.cxx

   1 // Created on: 2014-10-14
   2 // Created by: Anton POLETAEV
   3 // Copyright (c) 2013-2014 OPEN CASCADE SAS
   4 //
   5 // This file is part of Open CASCADE Technology software library.
   6 //
   7 // This library is free software; you can redistribute it and/or modify it under
   8 // the terms of the GNU Lesser General Public License version 2.1 as published
   9 // by the Free Software Foundation, with special exception defined in the file
  10 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
  11 // distribution for complete text of the license and disclaimer of any warranty.
  12 //
  13 // Alternatively, this file may be used under the terms of Open CASCADE
  14 // commercial license or contractual agreement.
  15
  16 #include <StdPrs_Isolines.hxx>
  17
  18 #include <Adaptor3d_IsoCurve.hxx>
  19 #include <BRepTools.hxx>
  20 #include <BRep_Tool.hxx>
  21 #include <GCPnts_AbscissaPoint.hxx>
  22 #include <GCPnts_QuasiUniformDeflection.hxx>
  23 #include <Geom_BezierSurface.hxx>
  24 #include <Geom_BSplineSurface.hxx>
  25 #include <Geom_Plane.hxx>
  26 #include <Geom_Line.hxx>
  27 #include <Geom2d_Line.hxx>
  28 #include <Geom2dAdaptor_Curve.hxx>
  29 #include <GeomAdaptor_Curve.hxx>
  30 #include <GeomLib_Tool.hxx>
  31 #include <gp_Lin2d.hxx>
  32 #include <Hatch_Hatcher.hxx>
  33 #include <NCollection_Shared.hxx>
  34 #include <Prs3d.hxx>
  35 #include <Prs3d_IsoAspect.hxx>
  36 #include <Poly_Array1OfTriangle.hxx>
  37 #include <Poly_Triangulation.hxx>
  38 #include <StdPrs_DeflectionCurve.hxx>
  39 #include <StdPrs_ToolRFace.hxx>
  40 #include <TColgp_SequenceOfPnt2d.hxx>
  41 #include <Standard_ErrorHandler.hxx>
  42 #include <Geom_Surface.hxx>
  43 #include <Geom_OffsetSurface.hxx>
  44 #include <Geom_RectangularTrimmedSurface.hxx>
  45 #include <algorithm>
  46
  47 namespace
  48 {
  49   const gp_Lin2d isoU (const Standard_Real theU) { return gp_Lin2d (gp_Pnt2d (theU, 0.0), gp::DY2d()); }
  50   const gp_Lin2d isoV (const Standard_Real theV) { return gp_Lin2d (gp_Pnt2d (0.0, theV), gp::DX2d()); }
  51
  52   typedef NCollection_Shared< NCollection_Vector<StdPrs_Isolines::SegOnIso> > VecOfSegments;
  53   typedef NCollection_Sequence<Handle(VecOfSegments)> SeqOfVecOfSegments;
  54
  55   //! Pack isoline segments into polylines.
  56   static void sortSegments (const SeqOfVecOfSegments&   theSegments,
  57                             Prs3d_NListOfSequenceOfPnt& thePolylines)
  58   {
  59     for (SeqOfVecOfSegments::Iterator aLineIter (theSegments); aLineIter.More(); aLineIter.Next())
  60     {
  61       Handle(VecOfSegments)& anIsoSegs = aLineIter.ChangeValue();
  62       std::stable_sort (anIsoSegs->begin(), anIsoSegs->end());
  63
  64       Handle(TColgp_HSequenceOfPnt) aPolyline = new TColgp_HSequenceOfPnt();
  65       thePolylines.Append (aPolyline);
  66       Standard_Real aLast = 0.0;
  67       for (VecOfSegments::Iterator aSegIter (*anIsoSegs); aSegIter.More(); aSegIter.Next())
  68       {
  69         if (!aPolyline->IsEmpty()
  70          && Abs (aSegIter.Value()[0].Param - aLast) > Precision::PConfusion())
  71         {
  72           aPolyline = new TColgp_HSequenceOfPnt();
  73           thePolylines.Append (aPolyline);
  74         }
  75
  76         aPolyline->Append (aSegIter.Value()[0].Pnt);
  77         aPolyline->Append (aSegIter.Value()[1].Pnt);
  78         aLast = aSegIter.Value()[1].Param;
  79       }
  80     }
  81   }
  82
  83   //! Reoder and adjust to the limit a curve's parameter values.
  84   //! @param theCurve [in] the curve.
  85   //! @param theLimit [in] the parameter limit value.
  86   //! @param theFirst [in/out] the first parameter value.
  87   //! @param theLast  [in/out] the last parameter value.
  88   static void findLimits (const Adaptor3d_Curve& theCurve,
  89                           const Standard_Real    theLimit,
  90                           Standard_Real&         theFirst,
  91                           Standard_Real&         theLast)
  92   {
  93     theFirst = Max (theCurve.FirstParameter(), theFirst);
  94     theLast  = Min (theCurve.LastParameter(), theLast);
  95
  96     Standard_Boolean isFirstInf = Precision::IsNegativeInfinite (theFirst);
  97     Standard_Boolean isLastInf  = Precision::IsPositiveInfinite (theLast);
  98
  99     if (!isFirstInf && !isLastInf)
 100     {
 101       return;
 102     }
 103
 104     gp_Pnt aP1, aP2;
 105     Standard_Real aDelta = 1.0;
 106     if (isFirstInf && isLastInf)
 107     {
 108       do
 109       {
 110         aDelta *= 2.0;
 111         theFirst = -aDelta;
 112         theLast  =  aDelta;
 113         theCurve.D0 (theFirst, aP1);
 114         theCurve.D0 (theLast,  aP2);
 115       }
 116       while (aP1.Distance (aP2) < theLimit);
 117     }
 118     else if (isFirstInf)
 119     {
 120       theCurve.D0 (theLast, aP2);
 121       do
 122       {
 123         aDelta *= 2.0;
 124         theFirst = theLast - aDelta;
 125         theCurve.D0 (theFirst, aP1);
 126       }
 127       while (aP1.Distance (aP2) < theLimit);
 128     }
 129     else if (isLastInf)
 130     {
 131       theCurve.D0 (theFirst, aP1);
 132       do
 133       {
 134         aDelta *= 2.0;
 135         theLast = theFirst + aDelta;
 136         theCurve.D0 (theLast, aP2);
 137       }
 138       while (aP1.Distance (aP2) < theLimit);
 139     }
 140   }
 141
 142 }
 143
 144 //==================================================================
 145 // function : AddOnTriangulation
 146 // purpose  :
 147 //==================================================================
 148 void StdPrs_Isolines::AddOnTriangulation (const Handle(Prs3d_Presentation)& thePresentation,
 149                                           const TopoDS_Face&                theFace,
 150                                           const Handle(Prs3d_Drawer)&       theDrawer)
 151 {
 152   Prs3d_NListOfSequenceOfPnt aUPolylines, aVPolylines;
 153   AddOnTriangulation (theFace, theDrawer, aUPolylines, aVPolylines);
 154   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->UIsoAspect(), aUPolylines);
 155   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->VIsoAspect(), aVPolylines);
 156 }
 157
 158 //==================================================================
 159 // function : AddOnTriangulation
 160 // purpose  :
 161 //==================================================================
 162 void StdPrs_Isolines::AddOnTriangulation (const TopoDS_Face&          theFace,
 163                                           const Handle(Prs3d_Drawer)& theDrawer,
 164                                           Prs3d_NListOfSequenceOfPnt& theUPolylines,
 165                                           Prs3d_NListOfSequenceOfPnt& theVPolylines)
 166 {
 167   const Standard_Integer aNbIsoU = theDrawer->UIsoAspect()->Number();
 168   const Standard_Integer aNbIsoV = theDrawer->VIsoAspect()->Number();
 169   if (aNbIsoU < 1 && aNbIsoV < 1)
 170   {
 171     return;
 172   }
 173
 174   // Evaluate parameters for uv isolines.
 175   TColStd_SequenceOfReal aUIsoParams;
 176   TColStd_SequenceOfReal aVIsoParams;
 177   Standard_Real aUmin = 0., aUmax = 0., aVmin = 0., aVmax = 0.;
 178   UVIsoParameters (theFace, aNbIsoU, aNbIsoV, theDrawer->MaximalParameterValue(), aUIsoParams, aVIsoParams,
 179                    aUmin, aUmax, aVmin, aVmax);
 180
 181   // Access surface definition.
 182   TopLoc_Location aLocSurface;
 183   Handle(Geom_Surface) aSurface = BRep_Tool::Surface (theFace, aLocSurface);
 184   if (aSurface.IsNull())
 185   {
 186     return;
 187   }
 188
 189   // Access triangulation.
 190   TopLoc_Location aLocTriangulation;
 191   const Handle(Poly_Triangulation)& aTriangulation = BRep_Tool::Triangulation (theFace, aLocTriangulation);
 192   if (aTriangulation.IsNull())
 193   {
 194     return;
 195   }
 196
 197   // Setup equal location for surface and triangulation.
 198   if (!aLocTriangulation.IsEqual (aLocSurface))
 199   {
 200     aSurface = Handle (Geom_Surface)::DownCast (
 201       aSurface->Transformed ((aLocSurface / aLocTriangulation).Transformation()));
 202   }
 203
 204   const Handle(Standard_Type)& TheType = aSurface->DynamicType();
 205   if (TheType == STANDARD_TYPE(Geom_OffsetSurface))
 206   {
 207     Standard_Real u1, u2, v1, v2;
 208     aSurface->Bounds(u1, u2, v1, v2);
 209     //Isolines of Offset surfaces are calculated by approximation and
 210     //cannot be calculated for infinite limits.
 211     if (Precision::IsInfinite(u1) || Precision::IsInfinite(u2) ||
 212       Precision::IsInfinite(v1) || Precision::IsInfinite(v2))
 213     {
 214       u1 = Max(aUmin, u1);
 215       u2 = Min(aUmax, u2);
 216       v1 = Max(aVmin, v1);
 217       v2 = Min(aVmax, v2);
 218       aSurface = new Geom_RectangularTrimmedSurface(aSurface, u1, u2, v1, v2);
 219     }
 220   }
 221
 222   addOnTriangulation (aTriangulation, aSurface, aLocTriangulation, aUIsoParams, aVIsoParams, theUPolylines, theVPolylines);
 223 }
 224
 225 //==================================================================
 226 // function : AddOnTriangulation
 227 // purpose  :
 228 //==================================================================
 229 void StdPrs_Isolines::AddOnTriangulation (const Handle(Prs3d_Presentation)& thePresentation,
 230                                           const Handle(Poly_Triangulation)& theTriangulation,
 231                                           const Handle(Geom_Surface)&       theSurface,
 232                                           const TopLoc_Location&            theLocation,
 233                                           const Handle(Prs3d_Drawer)&       theDrawer,
 234                                           const TColStd_SequenceOfReal&     theUIsoParams,
 235                                           const TColStd_SequenceOfReal&     theVIsoParams)
 236 {
 237   Prs3d_NListOfSequenceOfPnt aUPolylines, aVPolylines;
 238   addOnTriangulation (theTriangulation, theSurface, theLocation, theUIsoParams, theVIsoParams, aUPolylines, aVPolylines);
 239   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->UIsoAspect(), aUPolylines);
 240   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->VIsoAspect(), aVPolylines);
 241 }
 242
 243 //==================================================================
 244 // function : addOnTriangulation
 245 // purpose  :
 246 //==================================================================
 247 void StdPrs_Isolines::addOnTriangulation (const Handle(Poly_Triangulation)& theTriangulation,
 248                                           const Handle(Geom_Surface)&       theSurface,
 249                                           const TopLoc_Location&            theLocation,
 250                                           const TColStd_SequenceOfReal&     theUIsoParams,
 251                                           const TColStd_SequenceOfReal&     theVIsoParams,
 252                                           Prs3d_NListOfSequenceOfPnt&       theUPolylines,
 253                                           Prs3d_NListOfSequenceOfPnt&       theVPolylines)
 254 {
 255   const Poly_Array1OfTriangle& aTriangles = theTriangulation->Triangles();
 256   const TColgp_Array1OfPnt&    aNodes     = theTriangulation->Nodes();
 257   const TColgp_Array1OfPnt2d&  aUVNodes   = theTriangulation->UVNodes();
 258   for (Standard_Integer anUVIter = 0; anUVIter < 2; ++anUVIter)
 259   {
 260     const Standard_Boolean        isUIso      = anUVIter == 0;
 261     const TColStd_SequenceOfReal& anIsoParams = isUIso ? theUIsoParams : theVIsoParams;
 262     const Standard_Integer aNbIsolines = anIsoParams.Length();
 263     if (aNbIsolines == 0)
 264     {
 265       continue;
 266     }
 267
 268     SeqOfVecOfSegments aPolylines;
 269     TColStd_Array1OfInteger anIsoIndexes (1, aNbIsolines);
 270     anIsoIndexes.Init (-1);
 271     for (Standard_Integer anIsoIdx = 1; anIsoIdx <= aNbIsolines; ++anIsoIdx)
 272     {
 273       const gp_Lin2d anIsolineUV = isUIso ? isoU (anIsoParams.Value (anIsoIdx)) : isoV (anIsoParams.Value (anIsoIdx));
 274       Handle(VecOfSegments) anIsoPnts;
 275       if (anIsoIndexes.Value (anIsoIdx) != -1)
 276       {
 277         anIsoPnts = aPolylines.ChangeValue (anIsoIndexes.Value (anIsoIdx));
 278       }
 279
 280       for (Standard_Integer aTriIter = aTriangles.Lower(); aTriIter <= aTriangles.Upper(); ++aTriIter)
 281       {
 282         Standard_Integer aNodeIdxs[3];
 283         aTriangles.Value (aTriIter).Get (aNodeIdxs[0], aNodeIdxs[1],aNodeIdxs[2]);
 284         const gp_Pnt aNodesXYZ[3] = { aNodes.Value (aNodeIdxs[0]),
 285                                       aNodes.Value (aNodeIdxs[1]),
 286                                       aNodes.Value (aNodeIdxs[2]) };
 287         const gp_Pnt2d aNodesUV[3] = { aUVNodes.Value (aNodeIdxs[0]),
 288                                        aUVNodes.Value (aNodeIdxs[1]),
 289                                        aUVNodes.Value (aNodeIdxs[2]) };
 290
 291         // Find intersections with triangle in uv space and its projection on triangulation.
 292         SegOnIso aSegment;
 293         if (!findSegmentOnTriangulation (theSurface, isUIso, anIsolineUV, aNodesXYZ, aNodesUV, aSegment))
 294         {
 295           continue;
 296         }
 297
 298         if (anIsoPnts.IsNull())
 299         {
 300           aPolylines.Append (new VecOfSegments());
 301           anIsoIndexes.SetValue (anIsoIdx, aPolylines.Size());
 302           anIsoPnts = aPolylines.ChangeValue (anIsoIndexes.Value (anIsoIdx));
 303         }
 304
 305         if (!theLocation.IsIdentity())
 306         {
 307           aSegment[0].Pnt.Transform (theLocation);
 308           aSegment[1].Pnt.Transform (theLocation);
 309         }
 310         anIsoPnts->Append (aSegment);
 311       }
 312     }
 313
 314     sortSegments (aPolylines, isUIso ? theUPolylines : theVPolylines);
 315   }
 316 }
 317
 318 //==================================================================
 319 // function : AddOnSurface
 320 // purpose  :
 321 //==================================================================
 322 void StdPrs_Isolines::AddOnSurface (const Handle(Prs3d_Presentation)& thePresentation,
 323                                     const TopoDS_Face&                theFace,
 324                                     const Handle(Prs3d_Drawer)&       theDrawer,
 325                                     const Standard_Real               theDeflection)
 326 {
 327   Prs3d_NListOfSequenceOfPnt aUPolylines, aVPolylines;
 328   AddOnSurface (theFace, theDrawer, theDeflection, aUPolylines, aVPolylines);
 329   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->UIsoAspect(), aUPolylines);
 330   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->VIsoAspect(), aVPolylines);
 331 }
 332
 333 //==================================================================
 334 // function : AddOnSurface
 335 // purpose  :
 336 //==================================================================
 337 void StdPrs_Isolines::AddOnSurface (const TopoDS_Face&          theFace,
 338                                     const Handle(Prs3d_Drawer)& theDrawer,
 339                                     const Standard_Real         theDeflection,
 340                                     Prs3d_NListOfSequenceOfPnt& theUPolylines,
 341                                     Prs3d_NListOfSequenceOfPnt& theVPolylines)
 342 {
 343   const Standard_Integer aNbIsoU = theDrawer->UIsoAspect()->Number();
 344   const Standard_Integer aNbIsoV = theDrawer->VIsoAspect()->Number();
 345   if (aNbIsoU < 1 && aNbIsoV < 1)
 346   {
 347     return;
 348   }
 349
 350   // Evaluate parameters for uv isolines.
 351   TColStd_SequenceOfReal aUIsoParams, aVIsoParams;
 352   Standard_Real aUmin = 0., aUmax = 0., aVmin = 0., aVmax = 0.;
 353   UVIsoParameters (theFace, aNbIsoU, aNbIsoV, theDrawer->MaximalParameterValue(), aUIsoParams, aVIsoParams,
 354                    aUmin, aUmax, aVmin, aVmax);
 355
 356   BRepAdaptor_Surface aSurface (theFace);
 357   addOnSurface (new BRepAdaptor_HSurface (aSurface),
 358                 theDrawer,
 359                 theDeflection,
 360                 aUIsoParams,
 361                 aVIsoParams,
 362                 theUPolylines,
 363                 theVPolylines);
 364 }
 365
 366 //==================================================================
 367 // function : AddOnSurface
 368 // purpose  :
 369 //==================================================================
 370 void StdPrs_Isolines::AddOnSurface (const Handle(Prs3d_Presentation)&   thePresentation,
 371                                     const Handle(BRepAdaptor_HSurface)& theSurface,
 372                                     const Handle(Prs3d_Drawer)&         theDrawer,
 373                                     const Standard_Real                 theDeflection,
 374                                     const TColStd_SequenceOfReal&       theUIsoParams,
 375                                     const TColStd_SequenceOfReal&       theVIsoParams)
 376 {
 377   Prs3d_NListOfSequenceOfPnt aUPolylines, aVPolylines;
 378   addOnSurface (theSurface, theDrawer, theDeflection, theUIsoParams, theVIsoParams, aUPolylines, aVPolylines);
 379   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->UIsoAspect(), aUPolylines);
 380   Prs3d::AddPrimitivesGroup (thePresentation, theDrawer->VIsoAspect(), aVPolylines);
 381 }
 382
 383 //==================================================================
 384 // function : addOnSurface
 385 // purpose  :
 386 //==================================================================
 387 void StdPrs_Isolines::addOnSurface (const Handle(BRepAdaptor_HSurface)& theSurface,
 388                                     const Handle(Prs3d_Drawer)&         theDrawer,
 389                                     const Standard_Real                 theDeflection,
 390                                     const TColStd_SequenceOfReal&       theUIsoParams,
 391                                     const TColStd_SequenceOfReal&       theVIsoParams,
 392                                     Prs3d_NListOfSequenceOfPnt&         theUPolylines,
 393                                     Prs3d_NListOfSequenceOfPnt&         theVPolylines)
 394 {
 395   // Choose a deflection for sampling edge uv curves.
 396   Standard_Real aUVLimit = theDrawer->MaximalParameterValue();
 397   Standard_Real aUmin  = Max (theSurface->FirstUParameter(), -aUVLimit);
 398   Standard_Real aUmax  = Min (theSurface->LastUParameter(),   aUVLimit);
 399   Standard_Real aVmin  = Max (theSurface->FirstVParameter(), -aUVLimit);
 400   Standard_Real aVmax  = Min (theSurface->LastVParameter(),   aUVLimit);
 401   Standard_Real aSamplerDeflection = Max (aUmax - aUmin, aVmax - aVmin) * theDrawer->DeviationCoefficient();
 402   Standard_Real aHatchingTolerance = RealLast();
 403
 404   try
 405   {
 406     OCC_CATCH_SIGNALS
 407     // Determine edge points for trimming uv hatch region.
 408     TColgp_SequenceOfPnt2d aTrimPoints;
 409     StdPrs_ToolRFace anEdgeTool (theSurface);
 410     for (anEdgeTool.Init(); anEdgeTool.More(); anEdgeTool.Next())
 411     {
 412       TopAbs_Orientation anOrientation = anEdgeTool.Orientation();
 413       const Adaptor2d_Curve2d* anEdgeCurve = &anEdgeTool.Value();
 414       if (anEdgeCurve->GetType() != GeomAbs_Line)
 415       {
 416         GCPnts_QuasiUniformDeflection aSampler (*anEdgeCurve, aSamplerDeflection);
 417         if (!aSampler.IsDone())
 418         {
 419 #ifdef OCCT_DEBUG
 420           std::cout << "Cannot evaluate curve on surface" << std::endl;
 421 #endif
 422           continue;
 423         }
 424
 425         Standard_Integer aNumberOfPoints = aSampler.NbPoints();
 426         if (aNumberOfPoints < 2)
 427         {
 428           continue;
 429         }
 430
 431         for (Standard_Integer anI = 1; anI < aNumberOfPoints; ++anI)
 432         {
 433           gp_Pnt2d aP1 (aSampler.Value (anI    ).X(), aSampler.Value (anI    ).Y());
 434           gp_Pnt2d aP2 (aSampler.Value (anI + 1).X(), aSampler.Value (anI + 1).Y());
 435
 436           aHatchingTolerance = Min (aP1.SquareDistance (aP2), aHatchingTolerance);
 437
 438           aTrimPoints.Append (anOrientation == TopAbs_FORWARD ? aP1 : aP2);
 439           aTrimPoints.Append (anOrientation == TopAbs_FORWARD ? aP2 : aP1);
 440         }
 441       }
 442       else
 443       {
 444         Standard_Real aU1 = anEdgeCurve->FirstParameter();
 445         Standard_Real aU2 = anEdgeCurve->LastParameter();
 446
 447         // MSV 17.08.06 OCC13144: U2 occured less than U1, to overcome it
 448         // ensure that distance U2-U1 is not greater than aLimit*2,
 449         // if greater then choose an origin and use aLimit to define
 450         // U1 and U2 anew.
 451         Standard_Real anOrigin = 0.0;
 452
 453         if (!Precision::IsNegativeInfinite (aU1) || !Precision::IsPositiveInfinite(aU2))
 454         {
 455           if (Precision::IsNegativeInfinite (aU1))
 456           {
 457             anOrigin = aU2 - aUVLimit;
 458           }
 459           else if (Precision::IsPositiveInfinite (aU2))
 460           {
 461             anOrigin = aU1 + aUVLimit;
 462           }
 463           else
 464           {
 465             anOrigin = (aU1 + aU2) * 0.5;
 466           }
 467         }
 468
 469         aU1 = Max (anOrigin - aUVLimit, aU1);
 470         aU2 = Min (anOrigin + aUVLimit, aU2);
 471
 472         gp_Pnt2d aP1 = anEdgeCurve->Value (aU1);
 473         gp_Pnt2d aP2 = anEdgeCurve->Value (aU2);
 474
 475         aHatchingTolerance = Min (aP1.SquareDistance(aP2), aHatchingTolerance);
 476
 477         aTrimPoints.Append (anOrientation == TopAbs_FORWARD ? aP1 : aP2);
 478         aTrimPoints.Append (anOrientation == TopAbs_FORWARD ? aP2 : aP1);
 479       }
 480     }
 481
 482     // Compute a hatching tolerance.
 483     aHatchingTolerance *= 0.1;
 484     aHatchingTolerance = Max (Precision::Confusion(), aHatchingTolerance);
 485     aHatchingTolerance = Min (1.0E-5, aHatchingTolerance);
 486
 487     // Load isolines into hatcher.
 488     Hatch_Hatcher aHatcher (aHatchingTolerance, anEdgeTool.IsOriented());
 489
 490     for (Standard_Integer anIso = 1; anIso <= theUIsoParams.Length(); ++anIso)
 491     {
 492       aHatcher.AddXLine (theUIsoParams.Value (anIso));
 493     }
 494     for (Standard_Integer anIso = 1; anIso <= theVIsoParams.Length(); ++anIso)
 495     {
 496       aHatcher.AddYLine (theVIsoParams.Value (anIso));
 497     }
 498
 499     // Trim hatching region.
 500     for (Standard_Integer anI = 1; anI <= aTrimPoints.Length(); anI += 2)
 501     {
 502       aHatcher.Trim (aTrimPoints (anI), aTrimPoints (anI + 1));
 503     }
 504
 505     // Use surface definition for evaluation of Bezier, B-spline surface.
 506     // Use isoline adapter for other types of surfaces.
 507     GeomAbs_SurfaceType  aSurfType = theSurface->GetType();
 508     Handle(Geom_Surface) aBSurface;
 509     GeomAdaptor_Curve    aBSurfaceCurve;
 510     Adaptor3d_IsoCurve   aCanonicalCurve;
 511     if (aSurfType == GeomAbs_BezierSurface)
 512     {
 513       aBSurface = theSurface->Bezier();
 514     }
 515     else if (aSurfType == GeomAbs_BSplineSurface)
 516     {
 517       aBSurface = theSurface->BSpline();
 518     }
 519     else
 520     {
 521       aCanonicalCurve.Load (theSurface);
 522     }
 523
 524     // For each isoline: compute its segments.
 525     for (Standard_Integer anI = 1; anI <= aHatcher.NbLines(); anI++)
 526     {
 527       Standard_Real anIsoParam = aHatcher.Coordinate (anI);
 528       Standard_Boolean isIsoU  = aHatcher.IsXLine (anI);
 529
 530       // For each isoline's segment: evaluate its points.
 531       for (Standard_Integer aJ = 1; aJ <= aHatcher.NbIntervals (anI); aJ++)
 532       {
 533         Standard_Real aSegmentP1 = aHatcher.Start (anI, aJ);
 534         Standard_Real aSegmentP2 = aHatcher.End (anI, aJ);
 535
 536         if (!aBSurface.IsNull())
 537         {
 538           aBSurfaceCurve.Load (isIsoU ? aBSurface->UIso (anIsoParam) : aBSurface->VIso (anIsoParam));
 539
 540           findLimits (aBSurfaceCurve, aUVLimit, aSegmentP1, aSegmentP2);
 541
 542           if (aSegmentP2 - aSegmentP1 <= Precision::Confusion())
 543           {
 544             continue;
 545           }
 546         }
 547         else
 548         {
 549           aCanonicalCurve.Load (isIsoU ? GeomAbs_IsoU : GeomAbs_IsoV, anIsoParam, aSegmentP1, aSegmentP2);
 550
 551           findLimits (aCanonicalCurve, aUVLimit, aSegmentP1, aSegmentP2);
 552
 553           if (aSegmentP2 - aSegmentP1 <= Precision::Confusion())
 554           {
 555             continue;
 556           }
 557         }
 558         Adaptor3d_Curve* aCurve = aBSurface.IsNull() ? (Adaptor3d_Curve*) &aCanonicalCurve
 559           : (Adaptor3d_Curve*) &aBSurfaceCurve;
 560
 561         Handle(TColgp_HSequenceOfPnt) aPoints = new TColgp_HSequenceOfPnt();
 562         StdPrs_DeflectionCurve::Add (Handle(Prs3d_Presentation)(),
 563                                      *aCurve,
 564                                      aSegmentP1,
 565                                      aSegmentP2,
 566                                      theDeflection,
 567                                      aPoints->ChangeSequence(),
 568                                      theDrawer->DeviationAngle(),
 569                                      Standard_False);
 570         if (aPoints->IsEmpty())
 571         {
 572           continue;
 573         }
 574
 575         if (isIsoU)
 576         {
 577           theUPolylines.Append (aPoints);
 578         }
 579         else
 580         {
 581           theVPolylines.Append (aPoints);
 582         }
 583       }
 584     }
 585   }
 586   catch (Standard_Failure const&)
 587   {
 588     // ...
 589   }
 590 }
 591
 592 //==================================================================
 593 // function : UVIsoParameters
 594 // purpose  :
 595 //==================================================================
 596 void StdPrs_Isolines::UVIsoParameters (const TopoDS_Face&      theFace,
 597                                        const Standard_Integer  theNbIsoU,
 598                                        const Standard_Integer  theNbIsoV,
 599                                        const Standard_Real     theUVLimit,
 600                                        TColStd_SequenceOfReal& theUIsoParams,
 601                                        TColStd_SequenceOfReal& theVIsoParams,
 602                                        Standard_Real& theUmin,
 603                                        Standard_Real& theUmax,
 604                                        Standard_Real& theVmin,
 605                                        Standard_Real& theVmax)
 606 {
 607
 608   TopLoc_Location aLocation;
 609   const Handle(Geom_Surface)& aSurface = BRep_Tool::Surface(theFace, aLocation);
 610   if (aSurface.IsNull())
 611   {
 612     return;
 613   }
 614
 615   BRepTools::UVBounds (theFace, theUmin, theUmax, theVmin, theVmax);
 616
 617   Standard_Real aUmin = theUmin;
 618   Standard_Real aUmax = theUmax;
 619   Standard_Real aVmin = theVmin;
 620   Standard_Real aVmax = theVmax;
 621
 622   if (Precision::IsInfinite (aUmin))
 623     aUmin = -theUVLimit;
 624   if (Precision::IsInfinite (aUmax))
 625     aUmax = theUVLimit;
 626   if (Precision::IsInfinite (aVmin))
 627     aVmin = -theUVLimit;
 628   if (Precision::IsInfinite (aVmax))
 629     aVmax = theUVLimit;
 630
 631
 632   const Standard_Boolean isUClosed = aSurface->IsUClosed();
 633   const Standard_Boolean isVClosed = aSurface->IsVClosed();
 634
 635   if (!isUClosed)
 636   {
 637     aUmin = aUmin + (aUmax - aUmin) / 1000.0;
 638     aUmax = aUmax - (aUmax - aUmin) / 1000.0;
 639   }
 640
 641   if (!isVClosed)
 642   {
 643     aVmin = aVmin + (aVmax - aVmin) / 1000.0;
 644     aVmax = aVmax - (aVmax - aVmin) / 1000.0;
 645   }
 646
 647   Standard_Real aUstep = (aUmax - aUmin) / (1 + theNbIsoU);
 648   Standard_Real aVstep = (aVmax - aVmin) / (1 + theNbIsoV);
 649
 650   for (Standard_Integer anIso = 1; anIso <= theNbIsoU; ++anIso)
 651   {
 652     theUIsoParams.Append (aUmin + aUstep * anIso);
 653   }
 654
 655   for (Standard_Integer anIso = 1; anIso <= theNbIsoV; ++anIso)
 656   {
 657     theVIsoParams.Append (aVmin + aVstep * anIso);
 658   }
 659 }
 660
 661 //==================================================================
 662 // function : FindSegmentOnTriangulation
 663 // purpose  :
 664 //==================================================================
 665 Standard_Boolean StdPrs_Isolines::findSegmentOnTriangulation (const Handle(Geom_Surface)& theSurface,
 666                                                               const bool                  theIsU,
 667                                                               const gp_Lin2d&             theIsoline,
 668                                                               const gp_Pnt*               theNodesXYZ,
 669                                                               const gp_Pnt2d*             theNodesUV,
 670                                                               SegOnIso&                   theSegment)
 671 {
 672   Standard_Integer aNPoints = 0;
 673
 674   for (Standard_Integer aLinkIter = 0; aLinkIter < 3 && aNPoints < 2; ++aLinkIter)
 675   {
 676     // ...
 677     // Check that uv isoline crosses the triangulation link in parametric space
 678     // ...
 679
 680     const gp_Pnt2d& aNodeUV1 = theNodesUV[aLinkIter];
 681     const gp_Pnt2d& aNodeUV2 = theNodesUV[(aLinkIter + 1) % 3];
 682     const gp_Pnt& aNode1 = theNodesXYZ[aLinkIter];
 683     const gp_Pnt& aNode2 = theNodesXYZ[(aLinkIter + 1) % 3];
 684
 685     // Compute distance of uv points to isoline taking into consideration their relative
 686     // location against the isoline (left or right). Null value for a node means that the
 687     // isoline crosses the node. Both positive or negative means that the isoline does not
 688     // cross the segment.
 689     Standard_Boolean isLeftUV1 = (theIsoline.Direction().XY() ^ gp_Vec2d (theIsoline.Location(), aNodeUV1).XY()) > 0.0;
 690     Standard_Boolean isLeftUV2 = (theIsoline.Direction().XY() ^ gp_Vec2d (theIsoline.Location(), aNodeUV2).XY()) > 0.0;
 691     Standard_Real aDistanceUV1 = isLeftUV1 ? theIsoline.Distance (aNodeUV1) : -theIsoline.Distance (aNodeUV1);
 692     Standard_Real aDistanceUV2 = isLeftUV2 ? theIsoline.Distance (aNodeUV2) : -theIsoline.Distance (aNodeUV2);
 693
 694     // Isoline crosses first point of an edge.
 695     if (Abs (aDistanceUV1) < Precision::PConfusion())
 696     {
 697       theSegment[aNPoints].Param = theIsU ? aNodeUV1.Y() : aNodeUV1.X();
 698       theSegment[aNPoints].Pnt   = aNode1;
 699       ++aNPoints;
 700       continue;
 701     }
 702
 703     // Isoline crosses second point of an edge.
 704     if (Abs (aDistanceUV2) < Precision::PConfusion())
 705     {
 706       theSegment[aNPoints].Param = theIsU ? aNodeUV2.Y() : aNodeUV2.X();
 707       theSegment[aNPoints].Pnt   = aNode2;
 708
 709       ++aNPoints;
 710       ++aLinkIter;
 711       continue;
 712     }
 713
 714     // Isoline does not cross the triangle link.
 715     if (aDistanceUV1 * aDistanceUV2 > 0.0)
 716     {
 717       continue;
 718     }
 719
 720     // Isoline crosses degenerated link.
 721     if (aNode1.SquareDistance (aNode2) < Precision::PConfusion())
 722     {
 723       theSegment[aNPoints].Param = theIsU ? aNodeUV1.Y() : aNodeUV1.X();
 724       theSegment[aNPoints].Pnt   = aNode1;
 725       ++aNPoints;
 726       continue;
 727     }
 728
 729     // ...
 730     // Derive cross-point from parametric coordinates
 731     // ...
 732
 733     Standard_Real anAlpha = Abs (aDistanceUV1) / (Abs (aDistanceUV1) + Abs (aDistanceUV2));
 734
 735     gp_Pnt aCross (0.0, 0.0, 0.0);
 736     Standard_Real aCrossU = aNodeUV1.X() + anAlpha * (aNodeUV2.X() - aNodeUV1.X());
 737     Standard_Real aCrossV = aNodeUV1.Y() + anAlpha * (aNodeUV2.Y() - aNodeUV1.Y());
 738     Standard_Real aCrossParam = theIsU ? aCrossV : aCrossU;
 739     if (theSurface.IsNull())
 740     {
 741       // Do linear interpolation of point coordinates using triangulation nodes.
 742       aCross.SetX (aNode1.X() + anAlpha * (aNode2.X() - aNode1.X()));
 743       aCross.SetY (aNode1.Y() + anAlpha * (aNode2.Y() - aNode1.Y()));
 744       aCross.SetZ (aNode1.Z() + anAlpha * (aNode2.Z() - aNode1.Z()));
 745     }
 746     else
 747     {
 748       // Do linear interpolation of point coordinates by triangulation nodes.
 749       // Get 3d point on surface.
 750       Handle(Geom_Curve) anIso1, anIso2;
 751       Standard_Real aPntOnNode1Iso = 0.0;
 752       Standard_Real aPntOnNode2Iso = 0.0;
 753       Standard_Real aPntOnNode3Iso = 0.0;
 754
 755       if (theIsoline.Direction().X() == 0.0)
 756       {
 757         aPntOnNode1Iso = aNodeUV1.X();
 758         aPntOnNode2Iso = aNodeUV2.X();
 759         aPntOnNode3Iso = aCrossU;
 760         anIso1 = theSurface->VIso (aNodeUV1.Y());
 761         anIso2 = theSurface->VIso (aNodeUV2.Y());
 762       }
 763       else if (theIsoline.Direction().Y() == 0.0)
 764       {
 765         aPntOnNode1Iso = aNodeUV1.Y();
 766         aPntOnNode2Iso = aNodeUV2.Y();
 767         aPntOnNode3Iso = aCrossV;
 768         anIso1 = theSurface->UIso (aNodeUV1.X());
 769         anIso2 = theSurface->UIso (aNodeUV2.X());
 770       }
 771
 772       GeomAdaptor_Curve aCurveAdaptor1 (anIso1);
 773       GeomAdaptor_Curve aCurveAdaptor2 (anIso2);
 774       Standard_Real aLength1 = GCPnts_AbscissaPoint::Length (aCurveAdaptor1, aPntOnNode1Iso, aPntOnNode3Iso, 1e-2);
 775       Standard_Real aLength2 = GCPnts_AbscissaPoint::Length (aCurveAdaptor2, aPntOnNode2Iso, aPntOnNode3Iso, 1e-2);
 776       if (Abs (aLength1) < Precision::Confusion() || Abs (aLength2) < Precision::Confusion())
 777       {
 778         theSegment[aNPoints].Param = aCrossParam;
 779         theSegment[aNPoints].Pnt   = (aNode2.XYZ() - aNode1.XYZ()) * anAlpha + aNode1.XYZ();
 780         ++aNPoints;
 781         continue;
 782       }
 783
 784       aCross = (aNode2.XYZ() - aNode1.XYZ()) * (aLength1 / (aLength1 + aLength2)) + aNode1.XYZ();
 785     }
 786
 787     theSegment[aNPoints].Param = aCrossParam;
 788     theSegment[aNPoints].Pnt   = aCross;
 789     ++aNPoints;
 790   }
 791
 792   if (aNPoints != 2
 793    || Abs (theSegment[1].Param - theSegment[0].Param) <= Precision::PConfusion())
 794   {
 795     return false;
 796   }
 797
 798   if (theSegment[1].Param < theSegment[0].Param)
 799   {
 800     std::swap (theSegment[0], theSegment[1]);
 801   }
 802   return true;
 803 }