3e4c7db594605632cbe0241f8f56c95e9ef73383
[occt.git] / src / IntPatch / IntPatch_Intersection.cxx
1 // Created by: Modelization
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
3 //
4 // This file is part of Open CASCADE Technology software library.
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
11 //
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
14
15
16 #include <Adaptor3d_HSurface.hxx>
17 #include <Adaptor3d_TopolTool.hxx>
18 #include <IntPatch_ALine.hxx>
19 #include <IntPatch_ALineToWLine.hxx>
20 #include <IntPatch_GLine.hxx>
21 #include <IntPatch_ImpImpIntersection.hxx>
22 #include <IntPatch_ImpPrmIntersection.hxx>
23 #include <IntPatch_Intersection.hxx>
24 #include <IntPatch_Line.hxx>
25 #include <IntPatch_Point.hxx>
26 #include <IntPatch_PrmPrmIntersection.hxx>
27 #include <IntPatch_RLine.hxx>
28 #include <IntPatch_WLine.hxx>
29 #include <IntSurf_Quadric.hxx>
30 #include <Standard_ConstructionError.hxx>
31 #include <Standard_DomainError.hxx>
32 #include <Standard_OutOfRange.hxx>
33 #include <StdFail_NotDone.hxx>
34
35 #include <stdio.h>
36 #define DEBUG 0 
37 static const Standard_Integer aNbPointsInALine = 200;
38
39 //=======================================================================
40 //function : MinMax
41 //purpose  : Replaces theParMIN = MIN(theParMIN, theParMAX),
42 //                    theParMAX = MAX(theParMIN, theParMAX).
43 //=======================================================================
44 static inline void MinMax(Standard_Real& theParMIN, Standard_Real& theParMAX)
45 {
46   if(theParMIN > theParMAX)
47   {
48     const Standard_Real aTmp = theParMAX;
49     theParMAX = theParMIN;
50     theParMIN = aTmp;
51   }
52 }
53
54 //=======================================================================
55 //function : IsSeam
56 //purpose  : Returns:
57 //            0 - if interval [theU1, theU2] does not intersect the "seam-edge"
58 //                or if "seam-edge" do not exist;
59 //            1 - if interval (theU1, theU2) intersect the "seam-edge".
60 //            2 - if theU1 or/and theU2 lie ON the "seam-edge"
61 //
62 //ATTENTION!!!
63 //  If (theU1 == theU2) then this function will return only both 0 or 2.
64 //=======================================================================
65 static Standard_Integer IsSeam( const Standard_Real theU1,
66                                 const Standard_Real theU2,
67                                 const Standard_Real thePeriod)
68 {
69   if(IsEqual(thePeriod, 0.0))
70     return 0;
71
72   //If interval [theU1, theU2] intersect seam-edge then there exists an integer
73   //number N such as 
74   //    (theU1 <= T*N <= theU2) <=> (theU1/T <= N <= theU2/T),
75   //where T is the period.
76   //I.e. the inerval [theU1/T, theU2/T] must contain at least one
77   //integer number. In this case, Floor(theU1/T) and Floor(theU2/T)
78   //return different values or theU1/T is strictly integer number.
79   //Examples:
80   //  1. theU1/T==2.8, theU2/T==3.5 => Floor(theU1/T) == 2, Floor(theU2/T) == 3.
81   //  2. theU1/T==2.0, theU2/T==2.6 => Floor(theU1/T) == Floor(theU2/T) == 2.
82
83   const Standard_Real aVal1 = theU1/thePeriod,
84                       aVal2 = theU2/thePeriod;
85   const Standard_Integer aPar1 = static_cast<Standard_Integer>(Floor(aVal1));
86   const Standard_Integer aPar2 = static_cast<Standard_Integer>(Floor(aVal2));
87   if(aPar1 != aPar2)
88   {//Interval (theU1, theU2] intersects seam-edge
89     if(IsEqual(aVal2, static_cast<Standard_Real>(aPar2)))
90     {//aVal2 is an integer number => theU2 lies ON the "seam-edge"
91       return 2;
92     }
93
94     return 1;
95   }
96
97   //Here, aPar1 == aPar2. 
98
99   if(IsEqual(aVal1, static_cast<Standard_Real>(aPar1)))
100   {//aVal1 is an integer number => theU1 lies ON the "seam-edge"
101     return 2;
102   }
103
104   //If aVal2 is a true integer number then always (aPar1 != aPar2).
105
106   return 0;
107 }
108
109 //=======================================================================
110 //function : IsSeamOrBound
111 //purpose  : Returns TRUE if segment [thePtf, thePtl] intersects "seam-edge"
112 //            (if it exist) or surface boundaries and both thePtf and thePtl do
113 //            not match "seam-edge" or boundaries.
114 //           Point thePtmid lies in this segment. If thePtmid match
115 //            "seam-edge" or boundaries strictly (without any tolerance) then
116 //            the function will return TRUE.
117 //            See comments in function body for detail information.
118 //=======================================================================
119 static Standard_Boolean IsSeamOrBound(const IntSurf_PntOn2S& thePtf,
120                                       const IntSurf_PntOn2S& thePtl,
121                                       const IntSurf_PntOn2S& thePtmid,
122                                       const Standard_Real theU1Period,
123                                       const Standard_Real theU2Period,
124                                       const Standard_Real theV1Period,
125                                       const Standard_Real theV2Period,
126                                       const Standard_Real theUfSurf1,
127                                       const Standard_Real theUlSurf1,
128                                       const Standard_Real theVfSurf1,
129                                       const Standard_Real theVlSurf1,
130                                       const Standard_Real theUfSurf2,
131                                       const Standard_Real theUlSurf2,
132                                       const Standard_Real theVfSurf2,
133                                       const Standard_Real theVlSurf2)
134 {
135   Standard_Real aU11 = 0.0, aU12 = 0.0, aV11 = 0.0, aV12 = 0.0;
136   Standard_Real aU21 = 0.0, aU22 = 0.0, aV21 = 0.0, aV22 = 0.0;
137   thePtf.Parameters(aU11, aV11, aU12, aV12);
138   thePtl.Parameters(aU21, aV21, aU22, aV22);
139
140   MinMax(aU11, aU21);
141   MinMax(aV11, aV21);
142   MinMax(aU12, aU22);
143   MinMax(aV12, aV22);
144
145   if((aU11 - theUfSurf1)*(aU21 - theUfSurf1) < 0.0)
146   {//Interval [aU11, aU21] intersects theUfSurf1
147     return Standard_True;
148   }
149
150   if((aU11 - theUlSurf1)*(aU21 - theUlSurf1) < 0.0)
151   {//Interval [aU11, aU21] intersects theUlSurf1
152     return Standard_True;
153   }
154
155   if((aV11 - theVfSurf1)*(aV21 - theVfSurf1) < 0.0)
156   {//Interval [aV11, aV21] intersects theVfSurf1
157     return Standard_True;
158   }
159
160   if((aV11 - theVlSurf1)*(aV21 - theVlSurf1) < 0.0)
161   {//Interval [aV11, aV21] intersects theVlSurf1
162     return Standard_True;
163   }
164
165   if((aU12 - theUfSurf2)*(aU22 - theUfSurf2) < 0.0)
166   {//Interval [aU12, aU22] intersects theUfSurf2
167     return Standard_True;
168   }
169
170   if((aU12 - theUlSurf2)*(aU22 - theUlSurf2) < 0.0)
171   {//Interval [aU12, aU22] intersects theUlSurf2
172     return Standard_True;
173   }
174
175   if((aV12 - theVfSurf2)*(aV22 - theVfSurf2) < 0.0)
176   {//Interval [aV12, aV22] intersects theVfSurf2
177     return Standard_True;
178   }
179
180   if((aV12 - theVlSurf2)*(aV22 - theVlSurf2) < 0.0)
181   {//Interval [aV12, aV22] intersects theVlSurf2
182     return Standard_True;
183   }
184
185   if(IsSeam(aU11, aU21, theU1Period))
186     return Standard_True;
187
188   if(IsSeam(aV11, aV21, theV1Period))
189     return Standard_True;
190
191   if(IsSeam(aU12, aU22, theU2Period))
192     return Standard_True;
193
194   if(IsSeam(aV12, aV22, theV2Period))
195     return Standard_True;
196
197   /*
198     The segment [thePtf, thePtl] does not intersect the boundaries and
199     the seam-edge of the surfaces.
200     Nevertheless, following situation is possible:
201
202                   seam or
203                    bound
204                      |
205         thePtf  *    |
206                      |
207                      * thePtmid
208           thePtl  *  |
209                      |
210
211     This case must be processed, too.
212   */
213
214   Standard_Real aU1 = 0.0, aU2 = 0.0, aV1 = 0.0, aV2 = 0.0;
215   thePtmid.Parameters(aU1, aV1, aU2, aV2);
216
217   if(IsEqual(aU1, theUfSurf1) || IsEqual(aU1, theUlSurf1))
218     return Standard_True;
219
220   if(IsEqual(aU2, theUfSurf2) || IsEqual(aU2, theUlSurf2))
221     return Standard_True;
222
223   if(IsEqual(aV1, theVfSurf1) || IsEqual(aV1, theVlSurf1))
224     return Standard_True;
225
226   if(IsEqual(aV2, theVfSurf2) || IsEqual(aV2, theVlSurf2))
227     return Standard_True;
228
229   if(IsSeam(aU1, aU1, theU1Period))
230     return Standard_True;
231
232   if(IsSeam(aU2, aU2, theU2Period))
233     return Standard_True;
234
235   if(IsSeam(aV1, aV1, theV1Period))
236     return Standard_True;
237
238   if(IsSeam(aV2, aV2, theV2Period))
239     return Standard_True;
240
241   return Standard_False;
242 }
243
244 //=======================================================================
245 //function : JoinWLines
246 //purpose  : joins all WLines from theSlin to one if it is possible and
247 //            records the result into theSlin again.
248 //            Lines will be kept to be splitted if:
249 //              a) they are separated (has no common points);
250 //              b) resulted line (after joining) go through
251 //                 seam-edges or surface boundaries.
252 //
253 //          In addition, if points in theSPnt lies at least in one of 
254 //          the line in theSlin, this point will be deleted.
255 //=======================================================================
256 static void JoinWLines(IntPatch_SequenceOfLine& theSlin,
257                 IntPatch_SequenceOfPoint& theSPnt,
258                 const Standard_Real theTol3D,
259                 const Standard_Real theU1Period,
260                 const Standard_Real theU2Period,
261                 const Standard_Real theV1Period,
262                 const Standard_Real theV2Period,
263                 const Standard_Real theUfSurf1,
264                 const Standard_Real theUlSurf1,
265                 const Standard_Real theVfSurf1,
266                 const Standard_Real theVlSurf1,
267                 const Standard_Real theUfSurf2,
268                 const Standard_Real theUlSurf2,
269                 const Standard_Real theVfSurf2,
270                 const Standard_Real theVlSurf2)
271 {
272   if(theSlin.Length() == 0)
273     return;
274
275   for(Standard_Integer aNumOfLine1 = 1; aNumOfLine1 <= theSlin.Length(); aNumOfLine1++)
276   {
277     Handle(IntPatch_WLine) aWLine1 (Handle(IntPatch_WLine)::DownCast(theSlin.Value(aNumOfLine1)));
278
279     if(aWLine1.IsNull())
280     {//We must have failed to join not-point-lines
281       return;
282     }
283
284     const Standard_Integer aNbPntsWL1 = aWLine1->NbPnts();
285     const IntSurf_PntOn2S& aPntFW1 = aWLine1->Point(1);
286     const IntSurf_PntOn2S& aPntLW1 = aWLine1->Point(aNbPntsWL1);
287
288     for(Standard_Integer aNPt = 1; aNPt <= theSPnt.Length(); aNPt++)
289     {
290       const IntSurf_PntOn2S aPntCur = theSPnt.Value(aNPt).PntOn2S();
291
292       if( aPntCur.IsSame(aPntFW1, Precision::Confusion()) ||
293         aPntCur.IsSame(aPntLW1, Precision::Confusion()))
294       {
295         theSPnt.Remove(aNPt);
296         aNPt--;
297       }
298     }
299
300     Standard_Boolean hasBeenRemoved = Standard_False;
301     for(Standard_Integer aNumOfLine2 = aNumOfLine1 + 1; aNumOfLine2 <= theSlin.Length(); aNumOfLine2++)
302     {
303       Handle(IntPatch_WLine) aWLine2 (Handle(IntPatch_WLine)::DownCast(theSlin.Value(aNumOfLine2)));
304
305       const Standard_Integer aNbPntsWL2 = aWLine2->NbPnts();
306
307       const IntSurf_PntOn2S& aPntFWL1 = aWLine1->Point(1);
308       const IntSurf_PntOn2S& aPntLWL1 = aWLine1->Point(aNbPntsWL1);
309
310       const IntSurf_PntOn2S& aPntFWL2 = aWLine2->Point(1);
311       const IntSurf_PntOn2S& aPntLWL2 = aWLine2->Point(aNbPntsWL2);
312
313       if(aPntFWL1.IsSame(aPntFWL2, Precision::Confusion()))
314       {
315         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(2);
316         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(2);
317         if(!IsSeamOrBound(aPt1, aPt2, aPntFWL1, theU1Period, theU2Period,
318                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
319                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
320                           theVfSurf2, theVlSurf2))
321         {
322           aWLine1->ClearVertexes();
323           for(Standard_Integer aNPt = 1; aNPt <= aNbPntsWL2; aNPt++)
324           {
325             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
326             aWLine1->Curve()->InsertBefore(1, aPt);
327           }
328
329           aWLine1->ComputeVertexParameters(theTol3D);
330
331           theSlin.Remove(aNumOfLine2);
332           aNumOfLine2--;
333           hasBeenRemoved = Standard_True;
334
335           continue;
336         }
337       }
338
339       if(aPntFWL1.IsSame(aPntLWL2, Precision::Confusion()))
340       {
341         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(2);
342         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(aNbPntsWL2-1);
343         if(!IsSeamOrBound(aPt1, aPt2, aPntFWL1, theU1Period, theU2Period,
344                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
345                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
346                           theVfSurf2, theVlSurf2))
347         {
348           aWLine1->ClearVertexes();
349           for(Standard_Integer aNPt = aNbPntsWL2; aNPt >= 1; aNPt--)
350           {
351             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
352             aWLine1->Curve()->InsertBefore(1, aPt);
353           }
354
355           aWLine1->ComputeVertexParameters(theTol3D);
356
357           theSlin.Remove(aNumOfLine2);
358           aNumOfLine2--;
359           hasBeenRemoved = Standard_True;
360
361           continue;
362         }
363       }
364
365       if(aPntLWL1.IsSame(aPntFWL2, Precision::Confusion()))
366       {
367         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(aNbPntsWL1-1);
368         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(2);
369         if(!IsSeamOrBound(aPt1, aPt2, aPntLWL1, theU1Period, theU2Period,
370                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
371                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
372                           theVfSurf2, theVlSurf2))
373         {
374           aWLine1->ClearVertexes();
375           for(Standard_Integer aNPt = 1; aNPt <= aNbPntsWL2; aNPt++)
376           {
377             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
378             aWLine1->Curve()->Add(aPt);
379           }
380
381           aWLine1->ComputeVertexParameters(theTol3D);
382
383           theSlin.Remove(aNumOfLine2);
384           aNumOfLine2--;
385           hasBeenRemoved = Standard_True;
386
387           continue;
388         }
389       }
390
391       if(aPntLWL1.IsSame(aPntLWL2, Precision::Confusion()))
392       {
393         const IntSurf_PntOn2S& aPt1 = aWLine1->Point(aNbPntsWL1-1);
394         const IntSurf_PntOn2S& aPt2 = aWLine2->Point(aNbPntsWL2-1);
395         if(!IsSeamOrBound(aPt1, aPt2, aPntLWL1, theU1Period, theU2Period,
396                           theV1Period, theV2Period, theUfSurf1, theUlSurf1,
397                           theVfSurf1, theVlSurf1, theUfSurf2, theUlSurf2,
398                           theVfSurf2, theVlSurf2))
399         {
400           aWLine1->ClearVertexes();
401           for(Standard_Integer aNPt = aNbPntsWL2; aNPt >= 1; aNPt--)
402           {
403             const IntSurf_PntOn2S& aPt = aWLine2->Point(aNPt);
404             aWLine1->Curve()->Add(aPt);
405           }
406
407           aWLine1->ComputeVertexParameters(theTol3D);
408
409           theSlin.Remove(aNumOfLine2);
410           aNumOfLine2--;
411           hasBeenRemoved = Standard_True;
412
413           continue;
414         }
415       }
416     }
417
418     if(hasBeenRemoved)
419       aNumOfLine1--;
420   }
421 }
422
423 //======================================================================
424 // function: SequenceOfLine
425 //======================================================================
426 const IntPatch_SequenceOfLine& IntPatch_Intersection::SequenceOfLine() const { return(slin); }
427
428 //======================================================================
429 // function: IntPatch_Intersection
430 //======================================================================
431 IntPatch_Intersection::IntPatch_Intersection ()
432  : done(Standard_False),
433    //empt, tgte, oppo,
434    myTolArc(0.0), myTolTang(0.0),
435    myUVMaxStep(0.0), myFleche(0.0),
436    myIsStartPnt(Standard_False)
437    //myU1Start, myV1Start, myU2Start, myV2Start
438 {
439 }
440
441 //======================================================================
442 // function: IntPatch_Intersection
443 //======================================================================
444 IntPatch_Intersection::IntPatch_Intersection(const Handle(Adaptor3d_HSurface)&  S1,
445                                              const Handle(Adaptor3d_TopolTool)& D1,
446                                              const Handle(Adaptor3d_HSurface)&  S2,
447                                              const Handle(Adaptor3d_TopolTool)& D2,
448                                              const Standard_Real TolArc,
449                                              const Standard_Real TolTang)
450  : done(Standard_False),
451    //empt, tgte, oppo,
452    myTolArc(TolArc), myTolTang(TolTang),
453    myUVMaxStep(0.0), myFleche(0.0),
454    myIsStartPnt(Standard_False)
455    //myU1Start, myV1Start, myU2Start, myV2Start
456 {
457   if(myTolArc<1e-8) myTolArc=1e-8;
458   if(myTolTang<1e-8) myTolTang=1e-8;
459   if(myTolArc>0.5) myTolArc=0.5;
460   if(myTolTang>0.5) myTolTang=0.5;
461   Perform(S1,D1,S2,D2,TolArc,TolTang);
462 }
463
464 //======================================================================
465 // function: IntPatch_Intersection
466 //======================================================================
467 IntPatch_Intersection::IntPatch_Intersection(const Handle(Adaptor3d_HSurface)&  S1,
468                                              const Handle(Adaptor3d_TopolTool)& D1,
469                                              const Standard_Real TolArc,
470                                              const Standard_Real TolTang)
471  : done(Standard_False),
472    //empt, tgte, oppo,
473    myTolArc(TolArc), myTolTang(TolTang),
474    myUVMaxStep(0.0), myFleche(0.0),
475    myIsStartPnt(Standard_False)
476    //myU1Start, myV1Start, myU2Start, myV2Start
477 {
478   Perform(S1,D1,TolArc,TolTang);
479 }
480
481 //======================================================================
482 // function: SetTolerances
483 //======================================================================
484 void IntPatch_Intersection::SetTolerances(const Standard_Real TolArc,
485                                           const Standard_Real TolTang,
486                                           const Standard_Real UVMaxStep,
487                                           const Standard_Real Fleche)
488
489   myTolArc     = TolArc;
490   myTolTang    = TolTang;
491   myUVMaxStep  = UVMaxStep;
492   myFleche     = Fleche;
493   if(myTolArc<1e-8) myTolArc=1e-8;
494   if(myTolTang<1e-8) myTolTang=1e-8;
495   if(myTolArc>0.5) myTolArc=0.5;
496   if(myTolTang>0.5) myTolTang=0.5;  
497   if(myFleche<1.0e-3) myFleche=1e-3;
498   if(myUVMaxStep<1.0e-3) myUVMaxStep=1e-3;
499   if(myFleche>10) myFleche=10;
500   if(myUVMaxStep>0.5) myUVMaxStep=0.5;
501 }
502
503 //======================================================================
504 // function: Perform
505 //======================================================================
506 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  S1,
507                                     const Handle(Adaptor3d_TopolTool)& D1,
508                                     const Standard_Real TolArc,
509                                     const Standard_Real TolTang)
510 {
511   myTolArc = TolArc;
512   myTolTang = TolTang;
513   if(myFleche == 0.0)  myFleche = 0.01;
514   if(myUVMaxStep==0.0) myUVMaxStep = 0.01;
515
516   done = Standard_True;
517   spnt.Clear();
518   slin.Clear();
519   
520   empt = Standard_True;
521   tgte = Standard_False;
522   oppo = Standard_False;
523
524   switch (S1->GetType())
525   { 
526     case GeomAbs_Plane:
527     case GeomAbs_Cylinder:
528     case GeomAbs_Sphere:
529     case GeomAbs_Cone:
530     case GeomAbs_Torus: break;
531     default:
532     {
533       IntPatch_PrmPrmIntersection interpp;
534       interpp.Perform(S1,D1,TolArc,TolTang,myFleche,myUVMaxStep);
535       if (interpp.IsDone())
536       {
537         done = Standard_True;
538         tgte = Standard_False;
539         empt = interpp.IsEmpty();
540         const Standard_Integer nblm = interpp.NbLines();
541         for (Standard_Integer i=1; i<=nblm; i++) slin.Append(interpp.Line(i));
542       }
543     }
544     break;
545   }
546 }
547
548 /////////////////////////////////////////////////////////////////////////////
549 //  These several support functions provide methods which can help basic   //
550 //  algorithm to intersect infinite surfaces of the following types:       //
551 //                                                                         //
552 //  a.) SurfaceOfExtrusion;                                                //
553 //  b.) SurfaceOfRevolution;                                               //
554 //  c.) OffsetSurface.                                                     //
555 //                                                                         //
556 /////////////////////////////////////////////////////////////////////////////
557 #include <TColgp_Array1OfXYZ.hxx>
558 #include <TColgp_Array1OfPnt.hxx>
559 #include <TColgp_SequenceOfPnt.hxx>
560 #include <Extrema_ExtPS.hxx>
561 #include <Extrema_POnSurf.hxx>
562 #include <Geom2d_Curve.hxx>
563 #include <Geom2dAPI_InterCurveCurve.hxx>
564 #include <GeomAdaptor.hxx>
565 #include <GeomAdaptor_HCurve.hxx>
566 #include <GeomAdaptor_Curve.hxx>
567 #include <GeomAdaptor_Surface.hxx>
568 #include <GeomAdaptor_HSurface.hxx>
569 #include <Geom_Plane.hxx>
570 #include <ProjLib_ProjectOnPlane.hxx>
571 #include <GeomProjLib.hxx>
572 #include <ElCLib.hxx>
573 #include <Geom_TrimmedCurve.hxx>
574 #include <Geom_Surface.hxx>
575 #include <Geom_SurfaceOfLinearExtrusion.hxx>
576 #include <Geom_OffsetSurface.hxx>
577 #include <Geom_SurfaceOfRevolution.hxx>
578 #include <Geom_RectangularTrimmedSurface.hxx>
579
580 //===============================================================
581 //function: FUN_GetMinMaxXYZPnt
582 //===============================================================
583 static void FUN_GetMinMaxXYZPnt( const Handle(Adaptor3d_HSurface)& S,
584                                  gp_Pnt& pMin, gp_Pnt& pMax )
585 {
586   const Standard_Real DU = 0.25 * Abs(S->LastUParameter() - S->FirstUParameter());
587   const Standard_Real DV = 0.25 * Abs(S->LastVParameter() - S->FirstVParameter());
588   Standard_Real tMinXYZ = RealLast();
589   Standard_Real tMaxXYZ = -tMinXYZ;
590   gp_Pnt PUV, ptMax, ptMin;
591   for(Standard_Real U = S->FirstUParameter(); U <= S->LastUParameter(); U += DU)
592   {
593     for(Standard_Real V = S->FirstVParameter(); V <= S->LastVParameter(); V += DV)
594     {
595       S->D0(U,V,PUV);
596       const Standard_Real cXYZ = PUV.XYZ().Modulus();
597       if(cXYZ > tMaxXYZ) { tMaxXYZ = cXYZ; ptMax = PUV; }
598       if(cXYZ < tMinXYZ) { tMinXYZ = cXYZ; ptMin = PUV; }
599     }
600   }
601   pMin = ptMin;
602   pMax = ptMax;
603 }
604 //==========================================================================
605 //function: FUN_TrimInfSurf
606 //==========================================================================
607 static void FUN_TrimInfSurf(const gp_Pnt& Pmin,
608                             const gp_Pnt& Pmax,
609                             const Handle(Adaptor3d_HSurface)& InfSurf,
610                             const Standard_Real& AlternativeTrimPrm,
611                             Handle(Adaptor3d_HSurface)& TrimS)
612 {
613   Standard_Real TP = AlternativeTrimPrm;
614   Extrema_ExtPS ext1(Pmin, InfSurf->Surface(), 1.e-7, 1.e-7);
615   Extrema_ExtPS ext2(Pmax, InfSurf->Surface(), 1.e-7, 1.e-7);
616   if(ext1.IsDone() || ext2.IsDone())
617   {
618     Standard_Real Umax = -1.e+100, Umin = 1.e+100, Vmax = -1.e+100, Vmin = 1.e+100, cU, cV;
619     if(ext1.IsDone())
620     {
621       for(Standard_Integer i = 1; i <= ext1.NbExt(); i++)
622       {
623         const Extrema_POnSurf & pons = ext1.Point(i);
624         pons.Parameter(cU,cV);
625         if(cU > Umax) Umax = cU;
626         if(cU < Umin) Umin = cU;
627         if(cV > Vmax) Vmax = cV;
628         if(cV < Vmin) Vmin = cV;
629       }
630     }
631     if(ext2.IsDone())
632     {
633       for(Standard_Integer i = 1; i <= ext2.NbExt(); i++)
634       {
635         const Extrema_POnSurf & pons = ext2.Point(i);
636         pons.Parameter(cU,cV);
637         if(cU > Umax) Umax = cU;
638         if(cU < Umin) Umin = cU;
639         if(cV > Vmax) Vmax = cV;
640         if(cV < Vmin) Vmin = cV;
641       }
642     }
643     TP = Max(Abs(Umin),Max(Abs(Umax),Max(Abs(Vmin),Abs(Vmax))));
644   }
645   if(TP == 0.) { TrimS = InfSurf; return; }
646   else
647   {
648     const Standard_Boolean Uinf = Precision::IsNegativeInfinite(InfSurf->FirstUParameter()); 
649     const Standard_Boolean Usup = Precision::IsPositiveInfinite(InfSurf->LastUParameter());
650     const Standard_Boolean Vinf = Precision::IsNegativeInfinite(InfSurf->FirstVParameter()); 
651     const Standard_Boolean Vsup = Precision::IsPositiveInfinite(InfSurf->LastVParameter());
652     Handle(Adaptor3d_HSurface) TmpSS;
653     Standard_Integer IsTrimed = 0;
654     const Standard_Real tp = 1000.0 * TP;
655     if(Vinf && Vsup) { TrimS = InfSurf->VTrim(-tp, tp, 1.0e-7); IsTrimed = 1; }
656     if(Vinf && !Vsup){ TrimS = InfSurf->VTrim(-tp, InfSurf->LastVParameter(), 1.0e-7); IsTrimed = 1; }
657     if(!Vinf && Vsup){ TrimS = InfSurf->VTrim(InfSurf->FirstVParameter(), tp, 1.0e-7); IsTrimed = 1; }
658     if(IsTrimed)
659     {
660       TmpSS = TrimS;
661       if(Uinf && Usup)  TrimS = TmpSS->UTrim(-tp, tp, 1.0e-7);
662       if(Uinf && !Usup) TrimS = TmpSS->UTrim(-tp, InfSurf->LastUParameter(), 1.0e-7);
663       if(!Uinf && Usup) TrimS = TmpSS->UTrim(InfSurf->FirstUParameter(), tp, 1.0e-7);
664     }
665     else
666     {
667       if(Uinf && Usup)  TrimS = InfSurf->UTrim(-tp, tp, 1.0e-7);
668       if(Uinf && !Usup) TrimS = InfSurf->UTrim(-tp, InfSurf->LastUParameter(), 1.0e-7);
669       if(!Uinf && Usup) TrimS = InfSurf->UTrim(InfSurf->FirstUParameter(), tp, 1.0e-7);
670     }
671   }
672 }
673 //================================================================================
674 //function: FUN_GetUiso
675 //================================================================================
676 static void FUN_GetUiso(const Handle(Geom_Surface)& GS,
677                         const GeomAbs_SurfaceType&  T,
678                         const Standard_Real&        FirstV,
679                         const Standard_Real&        LastV,
680                         const Standard_Boolean&     IsVC,
681                         const Standard_Boolean&     IsVP,
682                         const Standard_Real&        U,
683                         Handle(Geom_Curve)&         I)
684 {
685   if(T !=  GeomAbs_OffsetSurface)
686   {
687     Handle(Geom_Curve) gc = GS->UIso(U);
688     if(IsVP && (FirstV == 0.0 && LastV == (2.*M_PI))) I = gc;
689     else
690     {
691       Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstV,LastV);
692       //szv:I = Handle(Geom_Curve)::DownCast(gtc);
693       I = gtc;
694     }
695   }
696   else//OffsetSurface
697   {
698     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (GS);
699     const Handle(Geom_Surface) bs = gos->BasisSurface();
700     Handle(Geom_Curve) gcbs = bs->UIso(U);
701     GeomAdaptor_Curve gac(gcbs);
702     const GeomAbs_CurveType GACT = gac.GetType();
703     if(IsVP || IsVC || GACT == GeomAbs_BSplineCurve || GACT == GeomAbs_BezierCurve || Abs(LastV - FirstV) < 1.e+5)
704     {
705       Handle(Geom_Curve) gc = gos->UIso(U);
706       if(IsVP && (FirstV == 0.0 && LastV == (2*M_PI))) I = gc;
707       else
708       {
709         Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstV,LastV);
710         //szv:I = Handle(Geom_Curve)::DownCast(gtc);
711         I = gtc;
712       }
713     }
714     else//Offset Line, Parab, Hyperb
715     {
716       Standard_Real VmTr, VMTr;
717       if(GACT != GeomAbs_Hyperbola)
718       {
719         if(FirstV >= 0. && LastV >= 0.){ VmTr = FirstV; VMTr = ((LastV - FirstV) > 1.e+4) ? (FirstV + 1.e+4) : LastV; }
720         else if(FirstV < 0. && LastV < 0.){ VMTr = LastV; VmTr = ((FirstV - LastV) < -1.e+4) ? (LastV - 1.e+4) : FirstV; }
721         else { VmTr = (FirstV < -1.e+4) ? -1.e+4 : FirstV; VMTr = (LastV > 1.e+4) ? 1.e+4 : LastV; }
722       }
723       else//Hyperbola
724       {
725         if(FirstV >= 0. && LastV >= 0.)
726         {
727           if(FirstV > 4.) return;
728           VmTr = FirstV; VMTr = (LastV > 4.) ? 4. : LastV;
729         }
730         else if(FirstV < 0. && LastV < 0.)
731         {
732           if(LastV < -4.) return;
733           VMTr = LastV; VmTr = (FirstV < -4.) ? -4. : FirstV;
734         }
735         else { VmTr = (FirstV < -4.) ? -4. : FirstV; VMTr = (LastV > 4.) ? 4. : LastV; }
736       }
737       //Make trimmed surface
738       Handle(Geom_RectangularTrimmedSurface) rts = new Geom_RectangularTrimmedSurface(gos,VmTr,VMTr,Standard_True);
739       I = rts->UIso(U);
740     }
741   }
742 }
743 //================================================================================
744 //function: FUN_GetViso
745 //================================================================================
746 static void FUN_GetViso(const Handle(Geom_Surface)& GS,
747                         const GeomAbs_SurfaceType&  T,
748                         const Standard_Real&        FirstU,
749                         const Standard_Real&        LastU,
750                         const Standard_Boolean&     IsUC,
751                         const Standard_Boolean&     IsUP,
752                         const Standard_Real&        V,
753                         Handle(Geom_Curve)&         I)
754 {
755   if(T !=  GeomAbs_OffsetSurface)
756   {
757     Handle(Geom_Curve) gc = GS->VIso(V);
758     if(IsUP && (FirstU == 0.0 && LastU == (2*M_PI))) I = gc;
759     else
760     {
761       Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstU,LastU);
762       //szv:I = Handle(Geom_Curve)::DownCast(gtc);
763       I = gtc;
764     }
765   }
766   else//OffsetSurface
767   {
768     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (GS);
769     const Handle(Geom_Surface) bs = gos->BasisSurface();
770     Handle(Geom_Curve) gcbs = bs->VIso(V);
771     GeomAdaptor_Curve gac(gcbs);
772     const GeomAbs_CurveType GACT = gac.GetType();
773     if(IsUP || IsUC || GACT == GeomAbs_BSplineCurve || GACT == GeomAbs_BezierCurve || Abs(LastU - FirstU) < 1.e+5)
774     {
775       Handle(Geom_Curve) gc = gos->VIso(V);
776       if(IsUP && (FirstU == 0.0 && LastU == (2*M_PI))) I = gc;
777       else
778       {
779         Handle(Geom_TrimmedCurve) gtc = new Geom_TrimmedCurve(gc,FirstU,LastU);
780         //szv:I = Handle(Geom_Curve)::DownCast(gtc);
781         I = gtc;
782       }
783     }
784     else//Offset Line, Parab, Hyperb
785     {
786       Standard_Real UmTr, UMTr;
787       if(GACT != GeomAbs_Hyperbola)
788       {
789         if(FirstU >= 0. && LastU >= 0.){ UmTr = FirstU; UMTr = ((LastU - FirstU) > 1.e+4) ? (FirstU + 1.e+4) : LastU; }
790         else if(FirstU < 0. && LastU < 0.){ UMTr = LastU; UmTr = ((FirstU - LastU) < -1.e+4) ? (LastU - 1.e+4) : FirstU; }
791         else { UmTr = (FirstU < -1.e+4) ? -1.e+4 : FirstU; UMTr = (LastU > 1.e+4) ? 1.e+4 : LastU; }
792       }
793       else//Hyperbola
794       {
795         if(FirstU >= 0. && LastU >= 0.)
796         {
797           if(FirstU > 4.) return;
798           UmTr = FirstU; UMTr = (LastU > 4.) ? 4. : LastU;
799         }
800         else if(FirstU < 0. && LastU < 0.)
801         {
802           if(LastU < -4.) return;
803           UMTr = LastU; UmTr = (FirstU < -4.) ? -4. : FirstU;
804         }
805         else { UmTr = (FirstU < -4.) ? -4. : FirstU; UMTr = (LastU > 4.) ? 4. : LastU; }
806       }
807       //Make trimmed surface
808       Handle(Geom_RectangularTrimmedSurface) rts = new Geom_RectangularTrimmedSurface(gos,UmTr,UMTr,Standard_True);
809       I = rts->VIso(V);
810     }
811   }
812 }
813 //================================================================================
814 //function: FUN_PL_Intersection
815 //================================================================================
816 static void FUN_PL_Intersection(const Handle(Adaptor3d_HSurface)& S1,
817                                 const GeomAbs_SurfaceType&        T1,
818                                 const Handle(Adaptor3d_HSurface)& S2,
819                                 const GeomAbs_SurfaceType&        T2,
820                                 Standard_Boolean&                 IsOk,
821                                 TColgp_SequenceOfPnt&             SP,
822                                 gp_Vec&                           DV)
823 {
824   IsOk = Standard_False;
825   // 1. Check: both surfaces have U(V)isos - lines.
826   DV = gp_Vec(0.,0.,1.);
827   Standard_Boolean isoS1isLine[2] = {0, 0};
828   Standard_Boolean isoS2isLine[2] = {0, 0};
829   Handle(Geom_Curve) C1, C2;
830   const GeomAdaptor_Surface & gas1 = *(GeomAdaptor_Surface*)(&(S1->Surface()));
831   const GeomAdaptor_Surface & gas2 = *(GeomAdaptor_Surface*)(&(S2->Surface()));
832   const Handle(Geom_Surface) gs1 = gas1.Surface();
833   const Handle(Geom_Surface) gs2 = gas2.Surface();
834   Standard_Real MS1[2], MS2[2];
835   MS1[0] = 0.5 * (S1->LastUParameter() + S1->FirstUParameter());
836   MS1[1] = 0.5 * (S1->LastVParameter() + S1->FirstVParameter());
837   MS2[0] = 0.5 * (S2->LastUParameter() + S2->FirstUParameter());
838   MS2[1] = 0.5 * (S2->LastVParameter() + S2->FirstVParameter());
839   if(T1 == GeomAbs_SurfaceOfExtrusion) isoS1isLine[0] = Standard_True;
840   else if(!S1->IsVPeriodic() && !S1->IsVClosed()) {
841     if(T1 != GeomAbs_OffsetSurface) C1 = gs1->UIso(MS1[0]);
842     else {
843       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
844       const Handle(Geom_Surface) bs = gos->BasisSurface();
845       C1 = bs->UIso(MS1[0]);
846     }
847     GeomAdaptor_Curve gac(C1);
848     if(gac.GetType() == GeomAbs_Line) isoS1isLine[0] = Standard_True;
849   }
850   if(!S1->IsUPeriodic() && !S1->IsUClosed()) {
851     if(T1 != GeomAbs_OffsetSurface) C1 = gs1->VIso(MS1[1]);
852     else {
853       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
854       const Handle(Geom_Surface) bs = gos->BasisSurface();
855       C1 = bs->VIso(MS1[1]);
856     }
857     GeomAdaptor_Curve gac(C1);
858     if(gac.GetType() == GeomAbs_Line) isoS1isLine[1] = Standard_True;
859   }
860   if(T2 == GeomAbs_SurfaceOfExtrusion) isoS2isLine[0] = Standard_True;
861   else if(!S2->IsVPeriodic() && !S2->IsVClosed()) {
862     if(T2 != GeomAbs_OffsetSurface) C2 = gs2->UIso(MS2[0]);
863     else {
864       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
865       const Handle(Geom_Surface) bs = gos->BasisSurface();
866       C2 = bs->UIso(MS2[0]);
867     }
868     GeomAdaptor_Curve gac(C2);
869     if(gac.GetType() == GeomAbs_Line) isoS2isLine[0] = Standard_True;
870   }
871   if(!S2->IsUPeriodic() && !S2->IsUClosed()) {
872     if(T2 != GeomAbs_OffsetSurface) C2 = gs2->VIso(MS2[1]);
873     else {
874       const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
875       const Handle(Geom_Surface) bs = gos->BasisSurface();
876       C2 = bs->VIso(MS2[1]);
877     }
878     GeomAdaptor_Curve gac(C2);
879     if(gac.GetType() == GeomAbs_Line) isoS2isLine[1] = Standard_True;
880   }
881   Standard_Boolean IsBothLines = ((isoS1isLine[0] || isoS1isLine[1]) &&
882                                   (isoS2isLine[0] || isoS2isLine[1]));
883   if(!IsBothLines){
884     return;
885   }
886   // 2. Check: Uiso lines of both surfaces are collinear.
887   gp_Pnt puvS1, puvS2;
888   gp_Vec derS1[2], derS2[2];
889   S1->D1(MS1[0], MS1[1], puvS1, derS1[0], derS1[1]);
890   S2->D1(MS2[0], MS2[1], puvS2, derS2[0], derS2[1]);
891   C1.Nullify(); C2.Nullify();
892   Standard_Integer iso = 0;
893   if(isoS1isLine[0] && isoS2isLine[0] &&
894      derS1[1].IsParallel(derS2[1],Precision::Angular())) {
895     iso = 1;
896     FUN_GetViso(gs1,T1,S1->FirstUParameter(),S1->LastUParameter(),
897                 S1->IsUClosed(),S1->IsUPeriodic(),MS1[1],C1);
898     FUN_GetViso(gs2,T2,S2->FirstUParameter(),S2->LastUParameter(),
899                 S2->IsUClosed(),S2->IsUPeriodic(),MS2[1],C2);
900   }
901   else if(isoS1isLine[0] && isoS2isLine[1] &&
902           derS1[1].IsParallel(derS2[0],Precision::Angular())) {
903     iso = 1;
904     FUN_GetViso(gs1,T1,S1->FirstUParameter(),S1->LastUParameter(),
905                 S1->IsUClosed(),S1->IsUPeriodic(),MS1[1],C1);
906     FUN_GetUiso(gs2,T2,S2->FirstVParameter(),S2->LastVParameter(),
907                 S2->IsVClosed(),S2->IsVPeriodic(),MS2[0],C2);
908   }
909   else if(isoS1isLine[1] && isoS2isLine[0] &&
910           derS1[0].IsParallel(derS2[1],Precision::Angular())) {
911     iso = 0;
912     FUN_GetUiso(gs1,T1,S1->FirstVParameter(),S1->LastVParameter(),
913                 S1->IsVClosed(),S1->IsVPeriodic(),MS1[0],C1);
914     FUN_GetViso(gs2,T2,S2->FirstUParameter(),S2->LastUParameter(),
915                 S2->IsUClosed(),S2->IsUPeriodic(),MS2[1],C2);
916   }
917   else if(isoS1isLine[1] && isoS2isLine[1] &&
918           derS1[0].IsParallel(derS2[0],Precision::Angular())) {
919     iso = 0;
920     FUN_GetUiso(gs1,T1,S1->FirstVParameter(),S1->LastVParameter(),
921                 S1->IsVClosed(),S1->IsVPeriodic(),MS1[0],C1);
922     FUN_GetUiso(gs2,T2,S2->FirstVParameter(),S2->LastVParameter(),
923                 S2->IsVClosed(),S2->IsVPeriodic(),MS2[0],C2);
924   }
925   else {
926     IsOk = Standard_False;
927     return;
928   }
929   IsOk = Standard_True;
930   // 3. Make intersections of V(U)isos
931   if(C1.IsNull() || C2.IsNull()) return;
932   DV = derS1[iso];
933   Handle(Geom_Plane) GPln = new Geom_Plane(gp_Pln(puvS1,gp_Dir(DV)));
934   Handle(Geom_Curve) C1Prj =
935     GeomProjLib::ProjectOnPlane(C1,GPln,gp_Dir(DV),Standard_True);
936   Handle(Geom_Curve) C2Prj =
937     GeomProjLib::ProjectOnPlane(C2,GPln,gp_Dir(DV),Standard_True);
938   if(C1Prj.IsNull() || C2Prj.IsNull()) return;
939   Handle(Geom2d_Curve) C1Prj2d =
940     GeomProjLib::Curve2d(C1Prj,Handle(Geom_Surface)::DownCast (GPln));
941   Handle(Geom2d_Curve) C2Prj2d =
942     GeomProjLib::Curve2d(C2Prj,Handle(Geom_Surface)::DownCast (GPln));
943   Geom2dAPI_InterCurveCurve ICC(C1Prj2d,C2Prj2d,1.0e-7);
944   if(ICC.NbPoints() > 0 )
945   {
946     for(Standard_Integer ip = 1; ip <= ICC.NbPoints(); ip++)
947     {
948       gp_Pnt2d P = ICC.Point(ip);
949       gp_Pnt P3d = ElCLib::To3d(gp_Ax2(puvS1,gp_Dir(DV)),P);
950       SP.Append(P3d);
951     }
952   }
953 }
954 //================================================================================
955 //function: FUN_NewFirstLast
956 //================================================================================
957 static void FUN_NewFirstLast(const GeomAbs_CurveType& ga_ct,
958                              const Standard_Real&     Fst,
959                              const Standard_Real&     Lst,
960                              const Standard_Real&     TrVal,
961                              Standard_Real&           NewFst,
962                              Standard_Real&           NewLst,
963                              Standard_Boolean&        NeedTr)
964 {
965   NewFst = Fst; NewLst = Lst; NeedTr = Standard_False;
966   switch (ga_ct)
967   {
968     case GeomAbs_Line:
969     case GeomAbs_Parabola:
970     {
971       if(Abs(Lst - Fst) > TrVal)
972       {
973         if(Fst >= 0. && Lst >= 0.)
974         {
975           NewFst = Fst;
976           NewLst = ((Fst + TrVal) < Lst) ? (Fst + TrVal) : Lst;
977         }
978         if(Fst < 0. && Lst < 0.)
979         {
980           NewLst = Lst;
981           NewFst = ((Lst - TrVal) > Fst) ? (Lst - TrVal) : Fst;
982         }
983         else
984         {
985           NewFst = (Fst < -TrVal) ? -TrVal : Fst;
986           NewLst = (Lst > TrVal) ? TrVal : Lst;
987         }
988         NeedTr = Standard_True;
989       }
990       break;
991     }
992     case GeomAbs_Hyperbola:
993     {
994       if(Abs(Lst - Fst) > 10.)
995       { 
996         if(Fst >= 0. && Lst >= 0.)
997         {
998           if(Fst > 4.) return;
999           NewFst = Fst;
1000           NewLst = (Lst > 4.) ? 4. : Lst;
1001         }
1002         if(Fst < 0. && Lst < 0.)
1003         {
1004           if(Lst < -4.) return;
1005           NewLst = Lst;
1006           NewFst = (Fst < -4.) ? -4. : Fst;
1007         }
1008         else
1009         {
1010           NewFst = (Fst < -4.) ? -4. : Fst;
1011           NewLst = (Lst > 4.) ? 4. : Lst;
1012         }
1013         NeedTr = Standard_True;
1014       }
1015       break;
1016     }
1017     default:
1018     break;
1019   }
1020 }
1021 //================================================================================
1022 //function: FUN_TrimBothSurf
1023 //================================================================================
1024 static void FUN_TrimBothSurf(const Handle(Adaptor3d_HSurface)& S1,
1025                              const GeomAbs_SurfaceType&        T1,
1026                              const Handle(Adaptor3d_HSurface)& S2,
1027                              const GeomAbs_SurfaceType&        T2,
1028                              const Standard_Real&              TV,
1029                              Handle(Adaptor3d_HSurface)&       NS1,
1030                              Handle(Adaptor3d_HSurface)&       NS2)
1031 {
1032   const GeomAdaptor_Surface & gas1 = *(GeomAdaptor_Surface*)(&(S1->Surface()));
1033   const GeomAdaptor_Surface & gas2 = *(GeomAdaptor_Surface*)(&(S2->Surface()));
1034   const Handle(Geom_Surface) gs1 = gas1.Surface();
1035   const Handle(Geom_Surface) gs2 = gas2.Surface();
1036   const Standard_Real UM1 = 0.5 * (S1->LastUParameter() + S1->FirstUParameter());
1037   const Standard_Real UM2 = 0.5 * (S2->LastUParameter() + S2->FirstUParameter());
1038   const Standard_Real VM1 = 0.5 * (S1->LastVParameter() + S1->FirstVParameter());
1039   const Standard_Real VM2 = 0.5 * (S2->LastVParameter() + S2->FirstVParameter());
1040   Handle(Geom_Curve) visoS1, visoS2, uisoS1, uisoS2;
1041   if(T1 != GeomAbs_OffsetSurface){ visoS1 = gs1->VIso(VM1); uisoS1 = gs1->UIso(UM1); }
1042   else
1043   {
1044     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs1);
1045     const Handle(Geom_Surface) bs = gos->BasisSurface();
1046     visoS1 = bs->VIso(VM1); uisoS1 = bs->UIso(UM1);
1047   }
1048   if(T2 != GeomAbs_OffsetSurface){ visoS2 = gs2->VIso(VM2); uisoS2 = gs2->UIso(UM2); }
1049   else
1050   {
1051     const Handle(Geom_OffsetSurface) gos = Handle(Geom_OffsetSurface)::DownCast (gs2);
1052     const Handle(Geom_Surface) bs = gos->BasisSurface();
1053     visoS2 = bs->VIso(VM2); uisoS2 = bs->UIso(UM2);
1054   }
1055   if(uisoS1.IsNull() || uisoS2.IsNull() || visoS1.IsNull() || visoS2.IsNull()){ NS1 = S1; NS2 = S2; return; }
1056   GeomAdaptor_Curve gau1(uisoS1);
1057   GeomAdaptor_Curve gav1(visoS1);
1058   GeomAdaptor_Curve gau2(uisoS2);
1059   GeomAdaptor_Curve gav2(visoS2);
1060   GeomAbs_CurveType GA_U1 = gau1.GetType();
1061   GeomAbs_CurveType GA_V1 = gav1.GetType();
1062   GeomAbs_CurveType GA_U2 = gau2.GetType();
1063   GeomAbs_CurveType GA_V2 = gav2.GetType();
1064   Standard_Boolean TrmU1 = Standard_False;
1065   Standard_Boolean TrmV1 = Standard_False;
1066   Standard_Boolean TrmU2 = Standard_False;
1067   Standard_Boolean TrmV2 = Standard_False;
1068   Standard_Real V1S1,V2S1,U1S1,U2S1, V1S2,V2S2,U1S2,U2S2;
1069   FUN_NewFirstLast(GA_U1,S1->FirstVParameter(),S1->LastVParameter(),TV,V1S1,V2S1,TrmV1);
1070   FUN_NewFirstLast(GA_V1,S1->FirstUParameter(),S1->LastUParameter(),TV,U1S1,U2S1,TrmU1);
1071   FUN_NewFirstLast(GA_U2,S2->FirstVParameter(),S2->LastVParameter(),TV,V1S2,V2S2,TrmV2);
1072   FUN_NewFirstLast(GA_V2,S2->FirstUParameter(),S2->LastUParameter(),TV,U1S2,U2S2,TrmU2);
1073   if(TrmV1) NS1 = S1->VTrim(V1S1, V2S1, 1.0e-7);
1074   if(TrmV2) NS2 = S2->VTrim(V1S2, V2S2, 1.0e-7);
1075   if(TrmU1)
1076   {
1077     if(TrmV1)
1078     {
1079       Handle(Adaptor3d_HSurface) TS = NS1;
1080       NS1 = TS->UTrim(U1S1, U2S1, 1.0e-7);
1081     }
1082     else NS1 = S1->UTrim(U1S1, U2S1, 1.0e-7);
1083   }
1084   if(TrmU2)
1085   {
1086     if(TrmV2)
1087     {
1088       Handle(Adaptor3d_HSurface) TS = NS2;
1089       NS2 = TS->UTrim(U1S2, U2S2, 1.0e-7);
1090     }
1091     else NS2 = S2->UTrim(U1S2, U2S2, 1.0e-7);
1092   }
1093 }
1094
1095 //=======================================================================
1096 //function : Perform
1097 //purpose  : 
1098 //=======================================================================
1099 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  theS1,
1100                                     const Handle(Adaptor3d_TopolTool)& theD1,
1101                                     const Handle(Adaptor3d_HSurface)&  theS2,
1102                                     const Handle(Adaptor3d_TopolTool)& theD2,
1103                                     const Standard_Real TolArc,
1104                                     const Standard_Real TolTang,
1105                                     const Standard_Boolean isGeomInt,
1106                                     const Standard_Boolean theIsReqToKeepRLine)
1107 {
1108   myTolArc = TolArc;
1109   myTolTang = TolTang;
1110   if(myFleche <= Precision::PConfusion())
1111     myFleche = 0.01;
1112   if(myUVMaxStep <= Precision::PConfusion())
1113     myUVMaxStep = 0.01;
1114
1115   done = Standard_False;
1116   spnt.Clear();
1117   slin.Clear();
1118   empt = Standard_True;
1119   tgte = Standard_False;
1120   oppo = Standard_False;
1121
1122   GeomAbs_SurfaceType typs1 = theS1->GetType();
1123   GeomAbs_SurfaceType typs2 = theS2->GetType();
1124   
1125   //treatment of the cases with cone or torus
1126   Standard_Boolean TreatAsBiParametric = Standard_False;
1127   Standard_Integer bGeomGeom = 0;
1128   //
1129   if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone ||
1130       typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1131     gp_Ax1 aCTAx, aGeomAx;
1132     GeomAbs_SurfaceType aCTType;
1133     Standard_Boolean bToCheck;
1134     //
1135     const Handle(Adaptor3d_HSurface)& aCTSurf = 
1136       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS1 : theS2;
1137     const Handle(Adaptor3d_HSurface)& aGeomSurf = 
1138       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS2 : theS1;
1139     //
1140     aCTType = aCTSurf->GetType();
1141     bToCheck = Standard_False;
1142     //
1143     if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone) {
1144       const gp_Cone aCon1 = (aCTType == GeomAbs_Cone) ? 
1145         aCTSurf->Cone() : aGeomSurf->Cone();
1146       Standard_Real a1 = Abs(aCon1.SemiAngle());
1147       bToCheck = (a1 < 0.02) || (a1 > 1.55);
1148       //
1149       if (typs1 == typs2) {
1150         const gp_Cone aCon2 = aGeomSurf->Cone();
1151         Standard_Real a2 = Abs(aCon2.SemiAngle());
1152         bToCheck = bToCheck || (a2 < 0.02) || (a2 > 1.55);
1153         //
1154         if (a1 > 1.55 && a2 > 1.55) {//quasi-planes: if same domain, treat as canonic
1155           const gp_Ax1 A1 = aCon1.Axis(), A2 = aCon2.Axis();
1156           if (A1.IsParallel(A2,Precision::Angular())) {
1157             const gp_Pnt Apex1 = aCon1.Apex(), Apex2 = aCon2.Apex();
1158             const gp_Pln Plan1( Apex1, A1.Direction() );
1159             if (Plan1.Distance( Apex2 ) <= Precision::Confusion()) {
1160               bToCheck = Standard_False;
1161             }
1162           }
1163         }
1164       }
1165       //
1166       TreatAsBiParametric = bToCheck;
1167       if (aCTType == GeomAbs_Cone) {
1168         aCTAx = aCon1.Axis();
1169       }
1170     }
1171     //
1172     if (typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1173       const gp_Torus aTor1 = (aCTType == GeomAbs_Torus) ? 
1174         aCTSurf->Torus() : aGeomSurf->Torus();
1175       bToCheck = aTor1.MajorRadius() > aTor1.MinorRadius();
1176       if (typs1 == typs2) {
1177         const gp_Torus aTor2 = aGeomSurf->Torus();
1178         bToCheck = aTor2.MajorRadius() > aTor2.MinorRadius();
1179       }
1180       //
1181       if (aCTType == GeomAbs_Torus) {
1182         aCTAx = aTor1.Axis();
1183       }
1184     }
1185     //
1186     if (bToCheck) {
1187       const gp_Lin aL1(aCTAx);
1188       //
1189       switch (aGeomSurf->GetType()) {
1190       case GeomAbs_Plane: {
1191         aGeomAx = aGeomSurf->Plane().Axis();
1192         if (aCTType == GeomAbs_Cone) {
1193           bGeomGeom = 1;
1194           if (Abs(aCTSurf->Cone().SemiAngle()) < 0.02) {
1195             Standard_Real ps = Abs(aCTAx.Direction().Dot(aGeomAx.Direction()));
1196             if(ps < 0.015) {
1197               bGeomGeom = 0;
1198             }
1199           }
1200         }
1201         else {
1202           if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) ||
1203               (aCTAx.IsNormal(aGeomAx, Precision::Angular()) && 
1204                (aGeomSurf->Plane().Distance(aCTAx.Location()) < Precision::Confusion()))) {
1205             bGeomGeom = 1;
1206           }
1207         }
1208         bToCheck = Standard_False;
1209         break;
1210       }
1211       case GeomAbs_Sphere: {
1212         if (aL1.Distance(aGeomSurf->Sphere().Location()) < Precision::Confusion()) {
1213           bGeomGeom = 1;
1214         }
1215         bToCheck = Standard_False;
1216         break;
1217       }
1218       case GeomAbs_Cylinder:
1219         aGeomAx = aGeomSurf->Cylinder().Axis();
1220         break;
1221       case GeomAbs_Cone: 
1222         aGeomAx = aGeomSurf->Cone().Axis();
1223         break;
1224       case GeomAbs_Torus: 
1225         aGeomAx = aGeomSurf->Torus().Axis();
1226         break;
1227       default: 
1228         bToCheck = Standard_False;
1229         break;
1230       }
1231       //
1232       if (bToCheck) {
1233         if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) &&
1234             (aL1.Distance(aGeomAx.Location()) <= Precision::Confusion())) {
1235           bGeomGeom = 1;
1236         }
1237       }
1238       //
1239       if (bGeomGeom == 1) {
1240         TreatAsBiParametric = Standard_False;
1241       }
1242     }
1243   }
1244   //
1245
1246   if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite()) {
1247     TreatAsBiParametric= Standard_False;
1248   }
1249
1250 //  Modified by skv - Mon Sep 26 14:58:30 2005 Begin
1251 //   if(TreatAsBiParametric) { typs1 = typs2 = GeomAbs_BezierSurface; }
1252   if(TreatAsBiParametric)
1253   {
1254     if (typs1 == GeomAbs_Cone && typs2 == GeomAbs_Plane)
1255       typs1 = GeomAbs_BezierSurface; // Using Imp-Prm Intersector
1256     else if (typs1 == GeomAbs_Plane && typs2 == GeomAbs_Cone)
1257       typs2 = GeomAbs_BezierSurface; // Using Imp-Prm Intersector
1258     else {
1259       // Using Prm-Prm Intersector
1260       typs1 = GeomAbs_BezierSurface;
1261       typs2 = GeomAbs_BezierSurface;
1262     }
1263   }
1264 //  Modified by skv - Mon Sep 26 14:58:30 2005 End
1265
1266   // Surface type definition
1267   Standard_Integer ts1 = 0;
1268   switch (typs1)
1269   {
1270     case GeomAbs_Plane:
1271     case GeomAbs_Cylinder:
1272     case GeomAbs_Sphere:
1273     case GeomAbs_Cone: ts1 = 1; break;
1274     case GeomAbs_Torus: ts1 = bGeomGeom; break;
1275     default: break;
1276   }
1277
1278   Standard_Integer ts2 = 0;
1279   switch (typs2)
1280   {
1281     case GeomAbs_Plane:
1282     case GeomAbs_Cylinder:
1283     case GeomAbs_Sphere:
1284     case GeomAbs_Cone: ts2 = 1; break;
1285     case GeomAbs_Torus: ts2 = bGeomGeom; break;
1286     default: break;
1287   }
1288   //
1289   // treatment of the cases with torus and any other geom surface
1290   //
1291   // Possible intersection types: 1. ts1 == ts2 == 1 <Geom-Geom>
1292   //                              2. ts1 != ts2      <Geom-Param>
1293   //                              3. ts1 == ts2 == 0 <Param-Param>
1294
1295   // Geom - Geom
1296   if(ts1 == ts2 && ts1 == 1)
1297   {
1298     const Standard_Boolean RestrictLine = Standard_True;
1299     IntSurf_ListOfPntOn2S ListOfPnts;
1300     ListOfPnts.Clear();
1301     if(isGeomInt)
1302     {
1303       if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite())
1304       {
1305         GeomGeomPerfom( theS1, theD1, theS2, theD2, TolArc, 
1306                         TolTang, ListOfPnts, RestrictLine,
1307                         typs1, typs2, theIsReqToKeepRLine);
1308       }
1309       else
1310       {
1311         GeomGeomPerfomTrimSurf( theS1, theD1, theS2, theD2,
1312                                 TolArc, TolTang, ListOfPnts, RestrictLine,
1313                                 typs1, typs2, theIsReqToKeepRLine);
1314       }
1315     }
1316     else
1317     {
1318       ParamParamPerfom(theS1, theD1, theS2, theD2, 
1319               TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1320     }
1321   }
1322
1323   // Geom - Param
1324   if(ts1 != ts2)
1325   {
1326     GeomParamPerfom(theS1, theD1, theS2, theD2, ts1 == 0, typs1, typs2);
1327   }
1328
1329   // Param - Param 
1330   if(ts1 == ts2 && ts1 == 0)
1331   {
1332     const Standard_Boolean RestrictLine = Standard_True;
1333     IntSurf_ListOfPntOn2S ListOfPnts;
1334     ListOfPnts.Clear();
1335
1336     ParamParamPerfom(theS1, theD1, theS2, theD2, TolArc,
1337                         TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1338   }
1339 }
1340
1341 //=======================================================================
1342 //function : Perform
1343 //purpose  : 
1344 //=======================================================================
1345 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  theS1,
1346                                     const Handle(Adaptor3d_TopolTool)& theD1,
1347                                     const Handle(Adaptor3d_HSurface)&  theS2,
1348                                     const Handle(Adaptor3d_TopolTool)& theD2,
1349                                     const Standard_Real TolArc,
1350                                     const Standard_Real TolTang,
1351                                     IntSurf_ListOfPntOn2S& ListOfPnts,
1352                                     const Standard_Boolean RestrictLine,
1353                                     const Standard_Boolean isGeomInt)
1354 {
1355   myTolArc = TolArc;
1356   myTolTang = TolTang;
1357   if(myFleche <= Precision::PConfusion())
1358     myFleche = 0.01;
1359   if(myUVMaxStep <= Precision::PConfusion())
1360     myUVMaxStep = 0.01;
1361     
1362   done = Standard_False;
1363   spnt.Clear();
1364   slin.Clear();
1365   empt = Standard_True;
1366   tgte = Standard_False;
1367   oppo = Standard_False;
1368
1369   GeomAbs_SurfaceType typs1 = theS1->GetType();
1370   GeomAbs_SurfaceType typs2 = theS2->GetType();
1371   //
1372   //treatment of the cases with cone or torus
1373   Standard_Boolean TreatAsBiParametric = Standard_False;
1374   Standard_Integer bGeomGeom = 0;
1375   //
1376   if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone ||
1377       typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1378     gp_Ax1 aCTAx, aGeomAx;
1379     GeomAbs_SurfaceType aCTType;
1380     Standard_Boolean bToCheck;
1381     //
1382     const Handle(Adaptor3d_HSurface)& aCTSurf = 
1383       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS1 : theS2;
1384     const Handle(Adaptor3d_HSurface)& aGeomSurf = 
1385       (typs1 == GeomAbs_Cone || typs1 == GeomAbs_Torus) ? theS2 : theS1;
1386     //
1387     aCTType = aCTSurf->GetType();
1388     bToCheck = Standard_False;
1389     //
1390     if (typs1 == GeomAbs_Cone  || typs2 == GeomAbs_Cone) {
1391       const gp_Cone aCon1 = (aCTType == GeomAbs_Cone) ? 
1392         aCTSurf->Cone() : aGeomSurf->Cone();
1393       Standard_Real a1 = Abs(aCon1.SemiAngle());
1394       bToCheck = (a1 < 0.02) || (a1 > 1.55);
1395       //
1396       if (typs1 == typs2) {
1397         const gp_Cone aCon2 = aGeomSurf->Cone();
1398         Standard_Real a2 = Abs(aCon2.SemiAngle());
1399         bToCheck = bToCheck || (a2 < 0.02) || (a2 > 1.55);
1400         //
1401         if (a1 > 1.55 && a2 > 1.55) {//quasi-planes: if same domain, treat as canonic
1402           const gp_Ax1 A1 = aCon1.Axis(), A2 = aCon2.Axis();
1403           if (A1.IsParallel(A2,Precision::Angular())) {
1404             const gp_Pnt Apex1 = aCon1.Apex(), Apex2 = aCon2.Apex();
1405             const gp_Pln Plan1( Apex1, A1.Direction() );
1406             if (Plan1.Distance( Apex2 ) <= Precision::Confusion()) {
1407               bToCheck = Standard_False;
1408             }
1409           }
1410         }
1411       }
1412       //
1413       TreatAsBiParametric = bToCheck;
1414       if (aCTType == GeomAbs_Cone) {
1415         aCTAx = aCon1.Axis();
1416       }
1417     }
1418     //
1419     if (typs1 == GeomAbs_Torus || typs2 == GeomAbs_Torus) {
1420       const gp_Torus aTor1 = (aCTType == GeomAbs_Torus) ? 
1421         aCTSurf->Torus() : aGeomSurf->Torus();
1422       bToCheck = aTor1.MajorRadius() > aTor1.MinorRadius();
1423       if (typs1 == typs2) {
1424         const gp_Torus aTor2 = aGeomSurf->Torus();
1425         bToCheck = aTor2.MajorRadius() > aTor2.MinorRadius();
1426       }
1427       //
1428       if (aCTType == GeomAbs_Torus) {
1429         aCTAx = aTor1.Axis();
1430       }
1431     }
1432     //
1433     if (bToCheck) {
1434       const gp_Lin aL1(aCTAx);
1435       //
1436       switch (aGeomSurf->GetType()) {
1437       case GeomAbs_Plane: {
1438         aGeomAx = aGeomSurf->Plane().Axis();
1439         if (aCTType == GeomAbs_Cone) {
1440           bGeomGeom = 1;
1441           if (Abs(aCTSurf->Cone().SemiAngle()) < 0.02) {
1442             Standard_Real ps = Abs(aCTAx.Direction().Dot(aGeomAx.Direction()));
1443             if(ps < 0.015) {
1444               bGeomGeom = 0;
1445             }
1446           }
1447         }
1448         else {
1449           if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) ||
1450               (aCTAx.IsNormal(aGeomAx, Precision::Angular()) && 
1451                (aGeomSurf->Plane().Distance(aCTAx.Location()) < Precision::Confusion()))) {
1452             bGeomGeom = 1;
1453           }
1454         }
1455         bToCheck = Standard_False;
1456         break;
1457       }
1458       case GeomAbs_Sphere: {
1459         if (aL1.Distance(aGeomSurf->Sphere().Location()) < Precision::Confusion()) {
1460           bGeomGeom = 1;
1461         }
1462         bToCheck = Standard_False;
1463         break;
1464       }
1465       case GeomAbs_Cylinder:
1466         aGeomAx = aGeomSurf->Cylinder().Axis();
1467         break;
1468       case GeomAbs_Cone: 
1469         aGeomAx = aGeomSurf->Cone().Axis();
1470         break;
1471       case GeomAbs_Torus: 
1472         aGeomAx = aGeomSurf->Torus().Axis();
1473         break;
1474       default: 
1475         bToCheck = Standard_False;
1476         break;
1477       }
1478       //
1479       if (bToCheck) {
1480         if (aCTAx.IsParallel(aGeomAx, Precision::Angular()) &&
1481             (aL1.Distance(aGeomAx.Location()) <= Precision::Confusion())) {
1482           bGeomGeom = 1;
1483         }
1484       }
1485       //
1486       if (bGeomGeom == 1) {
1487         TreatAsBiParametric = Standard_False;
1488       }
1489     }
1490   }
1491   //
1492
1493   if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite()) {
1494     TreatAsBiParametric= Standard_False;
1495   }
1496
1497   if(TreatAsBiParametric)
1498   {
1499     // Using Prm-Prm Intersector
1500     typs1 = GeomAbs_BezierSurface;
1501     typs2 = GeomAbs_BezierSurface;
1502   }
1503
1504   // Surface type definition
1505   Standard_Integer ts1 = 0;
1506   switch (typs1)
1507   {
1508     case GeomAbs_Plane:
1509     case GeomAbs_Cylinder:
1510     case GeomAbs_Sphere:
1511     case GeomAbs_Cone: ts1 = 1; break;
1512     case GeomAbs_Torus: ts1 = bGeomGeom; break;
1513     default: break;
1514   }
1515
1516   Standard_Integer ts2 = 0;
1517   switch (typs2)
1518   {
1519     case GeomAbs_Plane:
1520     case GeomAbs_Cylinder:
1521     case GeomAbs_Sphere:
1522     case GeomAbs_Cone: ts2 = 1; break;
1523     case GeomAbs_Torus: ts2 = bGeomGeom; break;
1524     default: break;
1525   }
1526   //
1527   // Possible intersection types: 1. ts1 == ts2 == 1 <Geom-Geom>
1528   //                              2. ts1 != ts2      <Geom-Param>
1529   //                              3. ts1 == ts2 == 0 <Param-Param>
1530
1531   if(!isGeomInt)
1532   {
1533     ParamParamPerfom(theS1, theD1, theS2, theD2, 
1534                 TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1535   }
1536   else if(ts1 != ts2)
1537   {
1538     GeomParamPerfom(theS1, theD1, theS2, theD2, ts1 == 0, typs1, typs2);
1539   }
1540   else if (ts1 == 0)
1541   {
1542     ParamParamPerfom(theS1, theD1, theS2, theD2,
1543                 TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1544   }
1545   else if(ts1 == 1)
1546   {
1547     if(theD1->DomainIsInfinite() || theD2->DomainIsInfinite())
1548     {
1549       GeomGeomPerfom(theS1, theD1, theS2, theD2, TolArc, 
1550                       TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1551     }
1552     else
1553     {
1554       GeomGeomPerfomTrimSurf(theS1, theD1, theS2, theD2,
1555               TolArc, TolTang, ListOfPnts, RestrictLine, typs1, typs2);
1556     }
1557   }
1558 }
1559
1560 //=======================================================================
1561 //function : ParamParamPerfom
1562 //purpose  : 
1563 //=======================================================================
1564 void IntPatch_Intersection::ParamParamPerfom(const Handle(Adaptor3d_HSurface)&  theS1,
1565                                              const Handle(Adaptor3d_TopolTool)& theD1,
1566                                              const Handle(Adaptor3d_HSurface)&  theS2,
1567                                              const Handle(Adaptor3d_TopolTool)& theD2,
1568                                              const Standard_Real TolArc,
1569                                              const Standard_Real TolTang,
1570                                              IntSurf_ListOfPntOn2S& ListOfPnts,
1571                                              const Standard_Boolean RestrictLine,
1572                                              const GeomAbs_SurfaceType typs1,
1573                                              const GeomAbs_SurfaceType typs2)
1574 {
1575   IntPatch_PrmPrmIntersection interpp;
1576   if(!theD1->DomainIsInfinite() && !theD2->DomainIsInfinite())
1577   {
1578     Standard_Boolean ClearFlag = Standard_True;
1579     if(!ListOfPnts.IsEmpty())
1580     {
1581       interpp.Perform(theS1,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep, ListOfPnts, RestrictLine);
1582       ClearFlag = Standard_False;
1583     }
1584     interpp.Perform(theS1,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep,ClearFlag);   //double call!!!!!!!
1585   }
1586   else if((theD1->DomainIsInfinite()) ^ (theD2->DomainIsInfinite()))
1587   {
1588     gp_Pnt pMaxXYZ, pMinXYZ;
1589     if(theD1->DomainIsInfinite())
1590     {
1591       FUN_GetMinMaxXYZPnt( theS2, pMinXYZ, pMaxXYZ );
1592       const Standard_Real MU = Max(Abs(theS2->FirstUParameter()),Abs(theS2->LastUParameter()));
1593       const Standard_Real MV = Max(Abs(theS2->FirstVParameter()),Abs(theS2->LastVParameter()));
1594       const Standard_Real AP = Max(MU, MV);
1595       Handle(Adaptor3d_HSurface) SS;
1596       FUN_TrimInfSurf(pMinXYZ, pMaxXYZ, theS1, AP, SS);
1597       interpp.Perform(SS,theD1,theS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1598     }
1599     else
1600     {
1601       FUN_GetMinMaxXYZPnt( theS1, pMinXYZ, pMaxXYZ );
1602       const Standard_Real MU = Max(Abs(theS1->FirstUParameter()),Abs(theS1->LastUParameter()));
1603       const Standard_Real MV = Max(Abs(theS1->FirstVParameter()),Abs(theS1->LastVParameter()));
1604       const Standard_Real AP = Max(MU, MV);
1605       Handle(Adaptor3d_HSurface) SS;
1606       FUN_TrimInfSurf(pMinXYZ, pMaxXYZ, theS2, AP, SS);
1607       interpp.Perform(theS1, theD1, SS, theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1608     }
1609   }//(theD1->DomainIsInfinite()) ^ (theD2->DomainIsInfinite())
1610   else
1611   {
1612     if(typs1 == GeomAbs_OtherSurface || typs2 == GeomAbs_OtherSurface)
1613     {
1614       done = Standard_False;
1615       return;
1616     }
1617
1618     Standard_Boolean IsPLInt = Standard_False;
1619     TColgp_SequenceOfPnt sop;
1620     gp_Vec v;
1621     FUN_PL_Intersection(theS1,typs1,theS2,typs2,IsPLInt,sop,v);
1622
1623     if(IsPLInt)
1624     {
1625       if(sop.Length() > 0)
1626       {
1627         for(Standard_Integer ip = 1; ip <= sop.Length(); ip++)
1628         {
1629           gp_Lin lin(sop.Value(ip),gp_Dir(v));
1630           Handle(IntPatch_GLine) gl = new IntPatch_GLine(lin,Standard_False);
1631           slin.Append(Handle(IntPatch_Line)::DownCast (gl));
1632         }
1633
1634         done = Standard_True;
1635       }
1636       else
1637         done = Standard_False;
1638
1639       return;
1640     }// 'COLLINEAR LINES'
1641     else
1642     {
1643       Handle(Adaptor3d_HSurface) nS1 = theS1;
1644       Handle(Adaptor3d_HSurface) nS2 = theS2;
1645       FUN_TrimBothSurf(theS1,typs1,theS2,typs2,1.e+8,nS1,nS2);
1646       interpp.Perform(nS1,theD1,nS2,theD2,TolArc,TolTang,myFleche,myUVMaxStep);
1647     }// 'NON - COLLINEAR LINES'
1648   }// both domains are infinite
1649
1650   if (interpp.IsDone())
1651   {
1652     done = Standard_True;
1653     tgte = Standard_False;
1654     empt = interpp.IsEmpty();
1655
1656     for(Standard_Integer i = 1; i <= interpp.NbLines(); i++)
1657     {
1658       if(interpp.Line(i)->ArcType() != IntPatch_Walking)
1659         slin.Append(interpp.Line(i));
1660     }
1661
1662     for (Standard_Integer i = 1; i <= interpp.NbLines(); i++)
1663     {
1664       if(interpp.Line(i)->ArcType() == IntPatch_Walking)
1665         slin.Append(interpp.Line(i));
1666     }
1667   }
1668 }
1669
1670 //=======================================================================
1671 ////function : GeomGeomPerfom
1672 //purpose  : 
1673 //=======================================================================
1674 void IntPatch_Intersection::GeomGeomPerfom(const Handle(Adaptor3d_HSurface)& theS1,
1675                                            const Handle(Adaptor3d_TopolTool)& theD1,
1676                                            const Handle(Adaptor3d_HSurface)& theS2,
1677                                            const Handle(Adaptor3d_TopolTool)& theD2,
1678                                            const Standard_Real TolArc,
1679                                            const Standard_Real TolTang,
1680                                            IntSurf_ListOfPntOn2S& ListOfPnts,
1681                                            const Standard_Boolean RestrictLine,
1682                                            const GeomAbs_SurfaceType theTyps1,
1683                                            const GeomAbs_SurfaceType theTyps2,
1684                                            const Standard_Boolean theIsReqToKeepRLine)
1685 {
1686   IntPatch_ImpImpIntersection interii(theS1,theD1,theS2,theD2,
1687                                       myTolArc,myTolTang, theIsReqToKeepRLine);
1688   const Standard_Boolean anIS = interii.IsDone();
1689   if (anIS)
1690   {
1691     done = anIS;
1692     empt = interii.IsEmpty();
1693     if (!empt)
1694     {
1695       tgte = interii.TangentFaces();
1696       if (tgte)
1697         oppo = interii.OppositeFaces();
1698
1699       for (Standard_Integer i = 1; i <= interii.NbLines(); i++)
1700       {
1701         const Handle(IntPatch_Line)& line = interii.Line(i);
1702         if (line->ArcType() == IntPatch_Analytic)
1703         {
1704           const GeomAbs_SurfaceType aTyps1 = theS1->GetType();
1705           const GeomAbs_SurfaceType aTyps2 = theS2->GetType();
1706           IntSurf_Quadric Quad1,Quad2;
1707           
1708           switch(aTyps1)
1709           {
1710           case GeomAbs_Plane:
1711             Quad1.SetValue(theS1->Plane());
1712             break;
1713
1714           case GeomAbs_Cylinder:
1715             Quad1.SetValue(theS1->Cylinder());
1716             break;
1717
1718           case GeomAbs_Sphere:
1719             Quad1.SetValue(theS1->Sphere());
1720             break;
1721
1722           case GeomAbs_Cone:
1723             Quad1.SetValue(theS1->Cone());
1724             break;
1725
1726           case GeomAbs_Torus:
1727             Quad1.SetValue(theS1->Torus());
1728             break;
1729
1730           default:
1731             break;
1732           }
1733
1734           switch(aTyps2)
1735           {
1736           case GeomAbs_Plane:
1737             Quad2.SetValue(theS2->Plane());
1738             break;
1739           case GeomAbs_Cylinder:
1740             Quad2.SetValue(theS2->Cylinder());
1741             break;
1742
1743           case GeomAbs_Sphere:
1744             Quad2.SetValue(theS2->Sphere());
1745             break;
1746
1747           case GeomAbs_Cone:
1748             Quad2.SetValue(theS2->Cone());
1749             break;
1750
1751           case GeomAbs_Torus:
1752             Quad2.SetValue(theS2->Torus());
1753             break;
1754
1755           default:
1756             break;
1757           }
1758
1759           IntPatch_ALineToWLine AToW(Quad1,Quad2,0.01,0.05,aNbPointsInALine);
1760           Handle(IntPatch_Line) wlin=AToW.MakeWLine((*((Handle(IntPatch_ALine) *)(&line))));
1761           slin.Append(wlin);
1762         }
1763         else
1764           slin.Append(interii.Line(i));
1765       }
1766
1767       for (Standard_Integer i = 1; i <= interii.NbPnts(); i++)
1768       {
1769         spnt.Append(interii.Point(i));
1770       }
1771     }
1772   }
1773   else
1774     ParamParamPerfom(theS1, theD1, theS2, theD2, 
1775                 TolArc, TolTang, ListOfPnts, RestrictLine, theTyps1, theTyps2);
1776 }
1777
1778 //=======================================================================
1779 //function : GeomParamPerfom
1780 //purpose  : 
1781 //=======================================================================
1782 void IntPatch_Intersection::
1783   GeomParamPerfom(const Handle(Adaptor3d_HSurface)&  theS1,
1784                   const Handle(Adaptor3d_TopolTool)& theD1,
1785                   const Handle(Adaptor3d_HSurface)&  theS2,
1786                   const Handle(Adaptor3d_TopolTool)& theD2,
1787                   const Standard_Boolean isNotAnalitical,
1788                   const GeomAbs_SurfaceType typs1,
1789                   const GeomAbs_SurfaceType typs2)
1790 {
1791   IntPatch_ImpPrmIntersection interip;
1792   if (myIsStartPnt)
1793   {
1794     if (isNotAnalitical/*ts1 == 0*/)
1795       interip.SetStartPoint(myU1Start,myV1Start);
1796     else
1797       interip.SetStartPoint(myU2Start,myV2Start);
1798   }
1799
1800   if(theD1->DomainIsInfinite() && theD2->DomainIsInfinite())
1801   {
1802     Standard_Boolean IsPLInt = Standard_False;
1803     TColgp_SequenceOfPnt sop;
1804     gp_Vec v;
1805     FUN_PL_Intersection(theS1,typs1,theS2,typs2,IsPLInt,sop,v);
1806     
1807     if(IsPLInt)
1808     {
1809       if(sop.Length() > 0)
1810       {
1811         for(Standard_Integer ip = 1; ip <= sop.Length(); ip++)
1812         {
1813           gp_Lin lin(sop.Value(ip),gp_Dir(v));
1814           Handle(IntPatch_GLine) gl = new IntPatch_GLine(lin,Standard_False);
1815           slin.Append(Handle(IntPatch_Line)::DownCast (gl));
1816         }
1817
1818         done = Standard_True;
1819       }
1820       else
1821         done = Standard_False;
1822
1823       return;
1824     }
1825     else
1826     {
1827       Handle(Adaptor3d_HSurface) nS1 = theS1;
1828       Handle(Adaptor3d_HSurface) nS2 = theS2;
1829       FUN_TrimBothSurf(theS1,typs1,theS2,typs2,1.e+5,nS1,nS2);
1830       interip.Perform(nS1,theD1,nS2,theD2,myTolArc,myTolTang,myFleche,myUVMaxStep);
1831     }
1832   }
1833   else
1834     interip.Perform(theS1,theD1,theS2,theD2,myTolArc,myTolTang,myFleche,myUVMaxStep);
1835
1836   if (interip.IsDone()) 
1837   {
1838     done = Standard_True;
1839     empt = interip.IsEmpty();
1840
1841     if (!empt)
1842     {
1843       const Standard_Integer aNbLines = interip.NbLines();
1844       for(Standard_Integer i = 1; i <= aNbLines; i++)
1845       {
1846         if(interip.Line(i)->ArcType() != IntPatch_Walking)
1847           slin.Append(interip.Line(i));
1848       }
1849
1850       for(Standard_Integer i = 1; i <= aNbLines; i++)
1851       {
1852         if(interip.Line(i)->ArcType() == IntPatch_Walking)
1853           slin.Append(interip.Line(i));
1854       }
1855
1856       for (Standard_Integer i = 1; i <= interip.NbPnts(); i++)
1857         spnt.Append(interip.Point(i));
1858     }
1859   }
1860 }
1861
1862 //=======================================================================
1863 //function : GeomGeomPerfomTrimSurf
1864 //purpose  : This function returns ready walking-line (which is not need
1865 //            in convertation) as an intersection line between two
1866 //            trimmed surfaces.
1867 //=======================================================================
1868 void IntPatch_Intersection::
1869   GeomGeomPerfomTrimSurf( const Handle(Adaptor3d_HSurface)& theS1,
1870                           const Handle(Adaptor3d_TopolTool)& theD1,
1871                           const Handle(Adaptor3d_HSurface)& theS2,
1872                           const Handle(Adaptor3d_TopolTool)& theD2,
1873                           const Standard_Real theTolArc,
1874                           const Standard_Real theTolTang,
1875                           IntSurf_ListOfPntOn2S& theListOfPnts,
1876                           const Standard_Boolean RestrictLine,
1877                           const GeomAbs_SurfaceType theTyps1,
1878                           const GeomAbs_SurfaceType theTyps2,
1879                           const Standard_Boolean theIsReqToKeepRLine)
1880 {
1881   IntSurf_Quadric Quad1,Quad2;
1882
1883   if((theTyps1 == GeomAbs_Cylinder) && (theTyps2 == GeomAbs_Cylinder))
1884   {
1885     IntPatch_ImpImpIntersection anInt;
1886     anInt.Perform(theS1, theD1, theS2, theD2, myTolArc,
1887                   myTolTang, Standard_True, theIsReqToKeepRLine);
1888
1889     done = anInt.IsDone();
1890
1891     if(done)
1892     {
1893       empt = anInt.IsEmpty();
1894       if (!empt)
1895       {
1896         tgte = anInt.TangentFaces();
1897         if (tgte)
1898           oppo = anInt.OppositeFaces();
1899
1900         const Standard_Integer aNbLin = anInt.NbLines();
1901         const Standard_Integer aNbPts = anInt.NbPnts();
1902
1903         for(Standard_Integer aLID = 1; aLID <= aNbLin; aLID++)
1904         {
1905           const Handle(IntPatch_Line)& aLine = anInt.Line(aLID);
1906           slin.Append(aLine);
1907         }
1908
1909         for(Standard_Integer aPID = 1; aPID <= aNbPts; aPID++)
1910         {
1911           const IntPatch_Point& aPoint = anInt.Point(aPID);
1912           spnt.Append(aPoint);
1913         }
1914
1915         JoinWLines( slin, spnt, theTolTang,
1916                     theS1->IsUPeriodic()? theS1->UPeriod() : 0.0,
1917                     theS2->IsUPeriodic()? theS2->UPeriod() : 0.0,
1918                     theS1->IsVPeriodic()? theS1->VPeriod() : 0.0,
1919                     theS2->IsVPeriodic()? theS2->VPeriod() : 0.0,
1920                     theS1->FirstUParameter(), theS1->LastUParameter(),
1921                     theS1->FirstVParameter(), theS1->LastVParameter(),
1922                     theS2->FirstUParameter(), theS2->LastUParameter(),
1923                     theS2->FirstVParameter(), theS2->LastVParameter());
1924       }
1925     }
1926   }
1927   else
1928   {
1929     GeomGeomPerfom(theS1, theD1, theS2, theD2,
1930             theTolArc, theTolTang, theListOfPnts,
1931             RestrictLine, theTyps1, theTyps2, theIsReqToKeepRLine);
1932   }
1933 }
1934
1935
1936 void IntPatch_Intersection::Perform(const Handle(Adaptor3d_HSurface)&  S1,
1937                                     const Handle(Adaptor3d_TopolTool)& D1,
1938                                     const Handle(Adaptor3d_HSurface)&  S2,
1939                                     const Handle(Adaptor3d_TopolTool)& D2,
1940                                     const Standard_Real U1,
1941                                     const Standard_Real V1,
1942                                     const Standard_Real U2,
1943                                     const Standard_Real V2,
1944                                     const Standard_Real TolArc,
1945                                     const Standard_Real TolTang)
1946 {
1947   myTolArc = TolArc;
1948   myTolTang = TolTang;
1949   if(myFleche == 0.0) {
1950 #if DEBUG
1951     //cout<<" -- IntPatch_Intersection::myFleche fixe par defaut a 0.01 --"<<endl;
1952     //cout<<" -- Utiliser la Methode SetTolerances( ... ) "<<endl;
1953 #endif
1954     myFleche = 0.01;
1955   }
1956   if(myUVMaxStep==0.0) {
1957 #if DEBUG
1958     //cout<<" -- IntPatch_Intersection::myUVMaxStep fixe par defaut a 0.01 --"<<endl;
1959     //cout<<" -- Utiliser la Methode SetTolerances( ... ) "<<endl;
1960 #endif
1961     myUVMaxStep = 0.01;
1962   }
1963
1964   done = Standard_False;
1965   spnt.Clear();
1966   slin.Clear();
1967
1968   empt = Standard_True;
1969   tgte = Standard_False;
1970   oppo = Standard_False;
1971
1972   const GeomAbs_SurfaceType typs1 = S1->GetType();
1973   const GeomAbs_SurfaceType typs2 = S2->GetType();
1974   
1975   if(   typs1==GeomAbs_Plane 
1976      || typs1==GeomAbs_Cylinder
1977      || typs1==GeomAbs_Sphere
1978      || typs1==GeomAbs_Cone
1979      || typs2==GeomAbs_Plane 
1980      || typs2==GeomAbs_Cylinder
1981      || typs2==GeomAbs_Sphere
1982      || typs2==GeomAbs_Cone)
1983   {
1984     myIsStartPnt = Standard_True;
1985     myU1Start = U1; myV1Start = V1; myU2Start = U2; myV2Start = V2;
1986     Perform(S1,D1,S2,D2,TolArc,TolTang);
1987     myIsStartPnt = Standard_False;
1988   }
1989   else
1990   {
1991     IntPatch_PrmPrmIntersection interpp;
1992     interpp.Perform(S1,D1,S2,D2,U1,V1,U2,V2,TolArc,TolTang,myFleche,myUVMaxStep);
1993     if (interpp.IsDone())
1994     {
1995       done = Standard_True;
1996       tgte = Standard_False;
1997       empt = interpp.IsEmpty();
1998       const Standard_Integer nblm = interpp.NbLines();
1999       Standard_Integer i = 1;
2000       for (; i<=nblm; i++) slin.Append(interpp.Line(i));
2001     }
2002   }
2003 }
2004 //======================================================================
2005 #include <IntPatch_IType.hxx>
2006 #include <IntPatch_LineConstructor.hxx>
2007 #include <Adaptor2d_HCurve2d.hxx>
2008 #include <Geom_Curve.hxx>
2009 #define MAXR 200
2010
2011
2012 //void IntPatch_Intersection__MAJ_R(Handle(Adaptor2d_HCurve2d) *R1,
2013 //                                     Handle(Adaptor2d_HCurve2d) *R2,
2014 //                                     int *NR1,
2015 //                                     int *NR2,
2016 //                                     Standard_Integer nbR1,
2017 //                                     Standard_Integer nbR2,
2018 //                                     const IntPatch_Point& VTX)
2019 void IntPatch_Intersection__MAJ_R(Handle(Adaptor2d_HCurve2d) *,
2020                                   Handle(Adaptor2d_HCurve2d) *,
2021                                   int *,
2022                                   int *,
2023                                   Standard_Integer ,
2024                                   Standard_Integer ,
2025                                   const IntPatch_Point& )
2026
2027   /*
2028   if(VTX.IsOnDomS1()) { 
2029     
2030     //-- long unsigned ptr= *((long unsigned *)(((Handle(Standard_Transient) *)(&(VTX.ArcOnS1())))));
2031     for(Standard_Integer i=0; i<nbR1;i++) { 
2032       if(VTX.ArcOnS1()==R1[i]) { 
2033         NR1[i]++;
2034         printf("\n ******************************");
2035         return;
2036       }
2037     }
2038     printf("\n R Pas trouvee  (IntPatch)\n");
2039     
2040   }
2041   */
2042 }
2043
2044
2045 //void IntPatch_Intersection::Dump(const Standard_Integer Mode,
2046 void IntPatch_Intersection::Dump(const Standard_Integer ,
2047                                  const Handle(Adaptor3d_HSurface)&  S1,
2048                                  const Handle(Adaptor3d_TopolTool)& D1,
2049                                  const Handle(Adaptor3d_HSurface)&  S2,
2050                                  const Handle(Adaptor3d_TopolTool)& D2) const 
2051
2052   
2053   //-- ----------------------------------------------------------------------
2054   //--  construction de la liste des restrictions & vertex 
2055   //--
2056   int NR1[MAXR],NR2[MAXR];
2057   Handle(Adaptor2d_HCurve2d) R1[MAXR],R2[MAXR];
2058   Standard_Integer nbR1=0,nbR2=0;
2059   for(D1->Init();D1->More() && nbR1<MAXR; D1->Next()) { 
2060     R1[nbR1]=D1->Value(); 
2061     NR1[nbR1]=0;
2062     nbR1++;
2063   }
2064   for(D2->Init();D2->More() && nbR2<MAXR; D2->Next()) { 
2065     R2[nbR2]=D2->Value();
2066     NR2[nbR2]=0;
2067     nbR2++;
2068   }
2069   
2070   printf("\nDUMP_INT:  ----empt:%2ud  tgte:%2ud  oppo:%2ud ---------------------------------",empt,tgte,empt);
2071   Standard_Integer i,j,nbr1,nbr2,nbgl,nbgc,nbge,nbgp,nbgh,nbl,nbr,nbg,nbw,nba;
2072   nbl=nbr=nbg=nbw=nba=nbgl=nbge=nbr1=nbr2=nbgc=nbgp=nbgh=0;
2073   nbl=NbLines();
2074   for(i=1;i<=nbl;i++) { 
2075     const Handle(IntPatch_Line)& line=Line(i);
2076     const IntPatch_IType IType=line->ArcType();
2077     if(IType == IntPatch_Walking) nbw++;
2078     else     if(IType == IntPatch_Restriction) { 
2079       nbr++;
2080       Handle(IntPatch_RLine) rlin (Handle(IntPatch_RLine)::DownCast (line));
2081       if(rlin->IsArcOnS1()) nbr1++;
2082       if(rlin->IsArcOnS2()) nbr2++;
2083     }
2084     else     if(IType == IntPatch_Analytic) nba++;
2085     else     { 
2086       nbg++; 
2087       if(IType == IntPatch_Lin) nbgl++;
2088       else if(IType == IntPatch_Circle) nbgc++;
2089       else if(IType == IntPatch_Parabola) nbgp++;
2090       else if(IType == IntPatch_Hyperbola) nbgh++;
2091       else if(IType == IntPatch_Ellipse) nbge++;
2092     }
2093   }
2094   
2095   
2096   printf("\nDUMP_INT:Lines:%2d Wlin:%2d Restr:%2d(On1:%2d On2:%2d) Ana:%2d Geom:%2d(L:%2d C:%2d E:%2d H:%2d P:%2d)",
2097          nbl,nbw,nbr,nbr1,nbr2,nba,nbg,nbgl,nbgc,nbge,nbgh,nbgp);
2098   
2099   IntPatch_LineConstructor LineConstructor(2);
2100   
2101   Standard_Integer nbllc=0;
2102   nbw=nbr=nbg=nba=0;
2103   Standard_Integer nbva,nbvw,nbvr,nbvg;
2104   nbva=nbvr=nbvw=nbvg=0;
2105   for (j=1; j<=nbl; j++) {
2106     Standard_Integer v,nbvtx;
2107     const Handle(IntPatch_Line)& intersLinej = Line(j);
2108     Standard_Integer NbLines;
2109     LineConstructor.Perform(SequenceOfLine(),intersLinej,S1,D1,S2,D2,1e-7);
2110     NbLines = LineConstructor.NbLines();
2111     
2112     for(Standard_Integer k=1;k<=NbLines;k++) { 
2113       nbllc++;
2114       const Handle(IntPatch_Line)& LineK = LineConstructor.Line(k);
2115       if (LineK->ArcType() == IntPatch_Analytic) { 
2116         Handle(IntPatch_ALine) alin (Handle(IntPatch_ALine)::DownCast (LineK));
2117         nbvtx=alin->NbVertex();
2118         nbva+=nbvtx;        nba++;
2119         for(v=1;v<=nbvtx;v++) { 
2120           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,alin->Vertex(v));
2121         }
2122       }
2123       else if (LineK->ArcType() == IntPatch_Restriction) {
2124         Handle(IntPatch_RLine) rlin (Handle(IntPatch_RLine)::DownCast (LineK));
2125         nbvtx=rlin->NbVertex();
2126         nbvr+=nbvtx;        nbr++;
2127         for(v=1;v<=nbvtx;v++) { 
2128           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,rlin->Vertex(v));
2129         }
2130       }
2131       else if (LineK->ArcType() == IntPatch_Walking) {
2132         Handle(IntPatch_WLine) wlin (Handle(IntPatch_WLine)::DownCast (LineK));
2133         nbvtx=wlin->NbVertex();
2134         nbvw+=nbvtx;        nbw++;
2135         for(v=1;v<=nbvtx;v++) { 
2136           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,wlin->Vertex(v));
2137         }
2138       }
2139       else { 
2140         Handle(IntPatch_GLine) glin (Handle(IntPatch_GLine)::DownCast (LineK));
2141         nbvtx=glin->NbVertex();
2142         nbvg+=nbvtx;        nbg++;
2143         for(v=1;v<=nbvtx;v++) { 
2144           IntPatch_Intersection__MAJ_R(R1,R2,NR1,NR2,nbR1,nbR2,glin->Vertex(v));
2145         }
2146       }
2147     }
2148   }
2149   printf("\nDUMP_LC :Lines:%2d WLin:%2d Restr:%2d Ana:%2d Geom:%2d",
2150          nbllc,nbw,nbr,nba,nbg);
2151   printf("\nDUMP_LC :vtx          :%2d     r:%2d    :%2d     :%2d",
2152          nbvw,nbvr,nbva,nbvg);
2153
2154
2155
2156    printf("\n");
2157 }