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