0023824: Bad results of sweep operation when a path curve has unpredictable torsion...
[occt.git] / src / GeomFill / GeomFill_Sweep.cxx
1 // Created on: 1997-11-21
2 // Created by: Philippe MANGIN
3 // Copyright (c) 1997-1999 Matra Datavision
4 // Copyright (c) 1999-2012 OPEN CASCADE SAS
5 //
6 // The content of this file is subject to the Open CASCADE Technology Public
7 // License Version 6.5 (the "License"). You may not use the content of this file
8 // except in compliance with the License. Please obtain a copy of the License
9 // at http://www.opencascade.org and read it completely before using this file.
10 //
11 // The Initial Developer of the Original Code is Open CASCADE S.A.S., having its
12 // main offices at: 1, place des Freres Montgolfier, 78280 Guyancourt, France.
13 //
14 // The Original Code and all software distributed under the License is
15 // distributed on an "AS IS" basis, without warranty of any kind, and the
16 // Initial Developer hereby disclaims all such warranties, including without
17 // limitation, any warranties of merchantability, fitness for a particular
18 // purpose or non-infringement. Please see the License for the specific terms
19 // and conditions governing the rights and limitations under the License.
20
21
22 //  Modified by skv - Fri Feb  6 11:44:48 2004 OCC5073
23
24 #include <GeomFill_Sweep.ixx>
25 #include <GeomFill_SweepFunction.hxx>
26 #include <GeomFill_LocFunction.hxx>
27
28 #include <Standard_ErrorHandler.hxx>
29
30 #include <gp_Pnt2d.hxx>
31 #include <gp_Dir2d.hxx>
32 #include <gp_Pnt.hxx>
33 #include <gp_Dir.hxx>
34 #include <gp_Lin.hxx>
35 #include <gp_Circ.hxx>
36 #include <gp_GTrsf.hxx>
37 #include <gp_Mat.hxx>
38 #include <gp_Ax2.hxx>
39
40 #include <TColgp_Array1OfPnt.hxx>
41 #include <TColgp_Array2OfPnt.hxx>
42 #include <TColgp_HArray2OfPnt.hxx>
43 //#include <GeomLib_Array1OfMat.hxx>
44 #include <TColStd_Array1OfInteger.hxx>
45 #include <TColStd_Array1OfReal.hxx>
46 #include <TColStd_Array2OfReal.hxx>
47
48 #include <GeomAbs_CurveType.hxx>
49 #include <GeomAdaptor_Curve.hxx>
50 #include <GeomLib.hxx>
51
52 #include <Geom2d_Line.hxx>
53 #include <Geom2d_BSplineCurve.hxx>
54 #include <Geom2d_TrimmedCurve.hxx>
55
56 #include <Geom_Circle.hxx>
57 #include <Geom_Line.hxx>
58 #include <Geom_BSplineSurface.hxx>
59 #include <Geom_Plane.hxx>
60 #include <Geom_SurfaceOfLinearExtrusion.hxx>
61 #include <Geom_CylindricalSurface.hxx>
62 #include <Geom_ConicalSurface.hxx>
63 #include <Geom_ToroidalSurface.hxx>
64 #include <Geom_SphericalSurface.hxx>
65 #include <Geom_SurfaceOfRevolution.hxx>
66 #include <Geom_RectangularTrimmedSurface.hxx>
67
68 #include <Approx_SweepApproximation.hxx>
69 #include <AdvApprox_PrefAndRec.hxx>
70 #include <AdvApprox_ApproxAFunction.hxx>
71 #include <GeomConvert_ApproxSurface.hxx>
72
73 #include <Precision.hxx>
74 #include <ElCLib.hxx>
75
76 //=======================================================================
77 //class : GeomFill_Sweep_Eval
78 //purpose: The evaluator for curve approximation
79 //=======================================================================
80
81 class GeomFill_Sweep_Eval : public AdvApprox_EvaluatorFunction
82 {
83  public:
84   GeomFill_Sweep_Eval (GeomFill_LocFunction& theTool)
85     : theAncore(theTool) {}
86   
87   virtual void Evaluate (Standard_Integer *Dimension,
88                          Standard_Real     StartEnd[2],
89                          Standard_Real    *Parameter,
90                          Standard_Integer *DerivativeRequest,
91                          Standard_Real    *Result, // [Dimension]
92                          Standard_Integer *ErrorCode);
93   
94  private:
95   GeomFill_LocFunction& theAncore;
96 };
97
98 void GeomFill_Sweep_Eval::Evaluate (Standard_Integer *,/*Dimension*/
99                                     Standard_Real     StartEnd[2],
100                                     Standard_Real    *Parameter,
101                                     Standard_Integer *DerivativeRequest,
102                                     Standard_Real    *Result,// [Dimension]
103                                     Standard_Integer *ErrorCode)
104 {
105   theAncore.DN (*Parameter,
106                 StartEnd[0],
107                 StartEnd[1],
108                 *DerivativeRequest, 
109                 Result[0],
110                 ErrorCode[0]);
111 }
112
113 //===============================================================
114 // Function : Create
115 // Purpose :
116 //===============================================================
117 GeomFill_Sweep::GeomFill_Sweep(const Handle(GeomFill_LocationLaw)& Location,
118                                const Standard_Boolean WithKpart)
119 {
120   done = Standard_False;
121
122   myLoc =  Location;
123   myKPart =  WithKpart;
124   SetTolerance(1.e-4);
125   myForceApproxC1 = Standard_False;
126
127   myLoc->GetDomain(First, Last);
128   SFirst = SLast = 30.081996;
129   SError = RealLast();
130 }
131
132 //===============================================================
133 // Function : SetDomain
134 // Purpose :
135 //===============================================================
136  void GeomFill_Sweep::SetDomain(const Standard_Real LocFirst,
137                                 const Standard_Real LocLast,
138                                 const Standard_Real SectionFirst,
139                                 const Standard_Real SectionLast) 
140 {
141   First = LocFirst;
142   Last =  LocLast;
143   SFirst = SectionFirst;
144   SLast = SectionLast;
145 }
146
147 //===============================================================
148 // Function : SetTolerance
149 // Purpose :
150 //===============================================================
151  void GeomFill_Sweep::SetTolerance(const Standard_Real Tolerance3d,
152                                    const Standard_Real BoundTolerance,
153                                    const Standard_Real Tolerance2d,
154                                    const Standard_Real ToleranceAngular) 
155 {
156   Tol3d = Tolerance3d;
157   BoundTol = BoundTolerance;
158   Tol2d =Tolerance2d;
159   TolAngular = ToleranceAngular;
160 }
161
162 //=======================================================================
163 //Function : SetForceApproxC1
164 //Purpose  : Set the flag that indicates attempt to approximate
165 //           a C1-continuous surface if a swept surface proved
166 //           to be C0.
167 //=======================================================================
168  void GeomFill_Sweep::SetForceApproxC1(const Standard_Boolean ForceApproxC1)
169 {
170   myForceApproxC1 = ForceApproxC1;
171 }
172
173
174 //===============================================================
175 // Function : ExchangeUV
176 // Purpose :
177 //===============================================================
178  Standard_Boolean GeomFill_Sweep::ExchangeUV() const
179 {
180   return myExchUV;
181 }
182
183 //===============================================================
184 // Function : UReversed
185 // Purpose :
186 //===============================================================
187  Standard_Boolean GeomFill_Sweep::UReversed() const
188 {
189   return isUReversed;
190 }
191
192 //===============================================================
193 // Function : VReversed
194 // Purpose :
195 //===============================================================
196  Standard_Boolean GeomFill_Sweep::VReversed() const
197 {
198   return isVReversed;
199 }
200
201 //===============================================================
202 // Function : Build
203 // Purpose :
204 //===============================================================
205  void GeomFill_Sweep::Build(const Handle(GeomFill_SectionLaw)& Section,        
206                             const GeomFill_ApproxStyle Methode,
207                             const GeomAbs_Shape Continuity,
208                             const Standard_Integer Degmax,
209                             const Standard_Integer Segmax) 
210 {
211   // Inits
212   done = Standard_False;
213   myExchUV = Standard_False;
214   isUReversed = isVReversed = Standard_False;
215   mySec =  Section;
216
217   if ((SFirst == SLast) &&  (SLast == 30.081996)) {
218     mySec->GetDomain(SFirst, SLast);
219   }
220    
221   Standard_Boolean isKPart = Standard_False, 
222                    isProduct  = Standard_False;
223  
224  // Traitement des KPart
225  if (myKPart)  isKPart = BuildKPart();
226  
227  // Traitement des produits Formelles
228  if ((!isKPart) && (Methode == GeomFill_Location)) {
229      Handle(Geom_BSplineSurface) BS;
230      BS = mySec->BSplineSurface();
231      if (! BS.IsNull()) {
232       // Approx de la loi
233 //    isProduct = BuildProduct(Continuity, Degmax, Segmax);
234      }
235    }
236
237  if (isKPart || isProduct) {
238  // Approx du 2d
239    done =  Build2d(Continuity, Degmax, Segmax);
240  }
241   else {
242  // Approx globale
243     done =  BuildAll(Continuity, Degmax, Segmax);
244   }
245 }
246
247 //===============================================================
248 // Function ::Build2d
249 // Purpose :A venir...
250 //===============================================================
251 // Standard_Boolean GeomFill_Sweep::Build2d(const GeomAbs_Shape Continuity,
252  Standard_Boolean GeomFill_Sweep::Build2d(const GeomAbs_Shape ,
253 //                                        const Standard_Integer Degmax,
254                                           const Standard_Integer ,
255 //                                        const Standard_Integer Segmax) 
256                                           const Standard_Integer ) 
257 {
258   Standard_Boolean Ok = Standard_False;
259   if (myLoc->Nb2dCurves() == 0) {
260     Ok = Standard_True;
261   }
262   return Ok;
263 }
264
265 //===============================================================
266 // Function : BuildAll
267 // Purpose :
268 //===============================================================
269  Standard_Boolean GeomFill_Sweep::BuildAll(const GeomAbs_Shape Continuity,
270                                            const Standard_Integer Degmax,
271                                            const Standard_Integer Segmax) 
272 {
273   Standard_Boolean Ok = Standard_False;
274   Standard_Integer nbsegmax = Segmax, nbspan = myLoc->NbIntervals(GeomAbs_C1);
275   if (Segmax < nbspan)  nbsegmax = nbspan;
276
277   Handle(GeomFill_SweepFunction) Func 
278     = new (GeomFill_SweepFunction) (mySec, myLoc, First, SFirst,
279                                    (SLast-SFirst)/(Last-First) );
280   Approx_SweepApproximation Approx( Func );
281
282   Approx.Perform(First,  Last,
283                  Tol3d,  BoundTol,  Tol2d,  TolAngular,
284                  Continuity, Degmax,  Segmax);
285
286   if (Approx.IsDone()) {
287     Ok = Standard_True;
288
289 #if DEB
290     Approx.Dump(cout);
291 #endif
292     
293     // La surface
294     Standard_Integer UDegree,VDegree,NbUPoles,
295                      NbVPoles,NbUKnots,NbVKnots;
296     Approx.SurfShape(UDegree,VDegree,NbUPoles,
297                      NbVPoles,NbUKnots,NbVKnots);
298
299     TColgp_Array2OfPnt Poles(1,NbUPoles, 1,NbVPoles);
300     TColStd_Array2OfReal Weights(1,NbUPoles, 1,NbVPoles);
301     TColStd_Array1OfReal UKnots(1, NbUKnots),VKnots(1, NbVKnots); 
302     TColStd_Array1OfInteger UMults(1, NbUKnots), VMults(1, NbVKnots);
303
304     Approx.Surface(Poles, Weights,
305                    UKnots,VKnots,
306                    UMults,VMults);
307
308     mySurface = new (Geom_BSplineSurface)
309                    (Poles, Weights,
310                     UKnots,VKnots,
311                     UMults,VMults,
312                     Approx.UDegree(),  Approx.VDegree(),
313                     mySec->IsUPeriodic());
314     SError = Approx. MaxErrorOnSurf();
315
316     if (myForceApproxC1 && !mySurface->IsCNv(1))
317     {
318       Standard_Real theTol = 1.e-4;
319       GeomAbs_Shape theUCont = GeomAbs_C1, theVCont = GeomAbs_C1;
320       Standard_Integer degU = 14, degV = 14;
321       Standard_Integer nmax = 16;
322       Standard_Integer thePrec = 1;  
323       
324       GeomConvert_ApproxSurface ConvertApprox(mySurface,theTol,theUCont,theVCont,
325                                               degU,degV,nmax,thePrec);
326       if (ConvertApprox.HasResult())
327       {
328         mySurface = ConvertApprox.Surface();
329         myCurve2d = new  (TColGeom2d_HArray1OfCurve) (1, 2);
330         CError =  new (TColStd_HArray2OfReal) (1,2, 1,2);
331
332         const Handle(Geom_BSplineSurface)& BSplSurf =
333           Handle(Geom_BSplineSurface)::DownCast(mySurface);
334         
335         gp_Dir2d D(0., 1.);
336         gp_Pnt2d P(BSplSurf->UKnot(1), 0);
337         Handle(Geom2d_Line) LC1 = new (Geom2d_Line) (P, D);
338         Handle(Geom2d_TrimmedCurve) TC1 =
339           new (Geom2d_TrimmedCurve) (LC1, 0, BSplSurf->VKnot(BSplSurf->NbVKnots()));
340         
341         myCurve2d->SetValue(1, TC1);
342         CError->SetValue(1, 1, 0.);
343         CError->SetValue(2, 1, 0.);
344  
345         P.SetCoord(BSplSurf->UKnot(BSplSurf->NbUKnots()), 0);
346         Handle(Geom2d_Line) LC2 = new (Geom2d_Line) (P, D);
347         Handle(Geom2d_TrimmedCurve) TC2 = 
348           new (Geom2d_TrimmedCurve) (LC2, 0, BSplSurf->VKnot(BSplSurf->NbVKnots()));
349         
350         myCurve2d->SetValue(myCurve2d->Length(), TC2);
351         CError->SetValue(1, myCurve2d->Length(), 0.);
352         CError->SetValue(2, myCurve2d->Length(), 0.);
353         
354         SError = theTol;
355       }
356     } //if (!mySurface->IsCNv(1))
357     
358     // Les Courbes 2d
359     if (myCurve2d.IsNull())
360     {
361       myCurve2d = new  (TColGeom2d_HArray1OfCurve) (1, 2+myLoc->TraceNumber());
362       CError =  new (TColStd_HArray2OfReal) (1,2, 1, 2+myLoc->TraceNumber());
363       Standard_Integer kk,ii, ifin = 1, ideb;
364       
365       if (myLoc->HasFirstRestriction()) {
366         ideb = 1;
367       }
368       else {
369         ideb = 2;
370       }
371       ifin += myLoc->TraceNumber();
372       if (myLoc->HasLastRestriction()) ifin++;
373       
374       for (ii=ideb, kk=1; ii<=ifin; ii++, kk++) {
375         Handle(Geom2d_BSplineCurve) C 
376           = new (Geom2d_BSplineCurve) (Approx.Curve2dPoles(kk),
377                                        Approx.Curves2dKnots(),
378                                        Approx.Curves2dMults(),
379                                        Approx.Curves2dDegree());
380         myCurve2d->SetValue(ii, C);
381         CError->SetValue(1, ii,  Approx.Max2dError(kk));
382         CError->SetValue(2, ii,  Approx.Max2dError(kk));
383       }
384       
385       // Si les courbes de restriction, ne sont pas calcules, on prend
386       // les iso Bords.
387       if (! myLoc->HasFirstRestriction()) {
388         gp_Dir2d D(0., 1.);
389         gp_Pnt2d P(UKnots(UKnots.Lower()), 0);
390         Handle(Geom2d_Line) LC = new (Geom2d_Line) (P, D);
391         Handle(Geom2d_TrimmedCurve) TC = new (Geom2d_TrimmedCurve)
392           (LC, First, Last);
393         
394         myCurve2d->SetValue(1, TC);
395         CError->SetValue(1, 1, 0.);
396         CError->SetValue(2, 1, 0.);
397       }
398       
399       if (! myLoc->HasLastRestriction()) {
400         gp_Dir2d D(0., 1.);
401         gp_Pnt2d P(UKnots(UKnots.Upper()), 0);
402         Handle(Geom2d_Line) LC = new (Geom2d_Line) (P, D);
403         Handle(Geom2d_TrimmedCurve) TC = 
404           new (Geom2d_TrimmedCurve) (LC, First, Last);
405         myCurve2d->SetValue(myCurve2d->Length(), TC);
406         CError->SetValue(1, myCurve2d->Length(), 0.);
407         CError->SetValue(2, myCurve2d->Length(), 0.);
408       }
409     } //if (myCurve2d.IsNull())
410   } 
411   return Ok;
412 }
413
414 //===============================================================
415 // Function : BuildProduct
416 // Purpose : A venir...
417 //===============================================================
418  Standard_Boolean GeomFill_Sweep::BuildProduct(const GeomAbs_Shape Continuity,
419                                           const Standard_Integer Degmax,
420                                           const Standard_Integer Segmax) 
421 {
422   Standard_Boolean Ok = Standard_False;
423
424   Handle(Geom_BSplineSurface) BSurf;
425   BSurf = Handle(Geom_BSplineSurface)::DownCast(
426            mySec->BSplineSurface()->Copy());
427   if (BSurf.IsNull()) return Ok; // Ce mode de construction est impossible  
428
429
430   Standard_Integer NbIntervalC2,  NbIntervalC3;
431   GeomFill_LocFunction Func(myLoc);
432
433   NbIntervalC2 = myLoc->NbIntervals(GeomAbs_C2);
434   NbIntervalC3 = myLoc->NbIntervals(GeomAbs_C3);
435   TColStd_Array1OfReal Param_de_decoupeC2 (1, NbIntervalC2+1);
436   myLoc->Intervals(Param_de_decoupeC2, GeomAbs_C2);
437   TColStd_Array1OfReal Param_de_decoupeC3 (1, NbIntervalC3+1);
438   myLoc->Intervals(Param_de_decoupeC3, GeomAbs_C3);
439
440
441   AdvApprox_PrefAndRec Preferentiel(Param_de_decoupeC2,
442                                      Param_de_decoupeC3);
443    
444   Handle(TColStd_HArray1OfReal) ThreeDTol = new (TColStd_HArray1OfReal) (1,4);
445   ThreeDTol->Init(Tol3d); // A Affiner...
446
447   GeomFill_Sweep_Eval eval (Func);
448   AdvApprox_ApproxAFunction Approx(0, 0, 4,
449                                    ThreeDTol,
450                                    ThreeDTol,
451                                    ThreeDTol,
452                                    First,
453                                    Last, 
454                                    Continuity,
455                                    Degmax,
456                                    Segmax, 
457                                    eval,
458                                    Preferentiel);
459 #if DEB
460   Approx.Dump(cout);
461 #endif
462
463   Ok = Approx.HasResult();
464   if (Ok) {
465 /*    TColgp_Array1OfMat TM(1, nbpoles);
466     Handle(TColgp_HArray2OfPnt) ResPoles;
467     ResPoles = Approx.Poles();
468
469     // Produit Tensoriel
470     for (ii=1; ii<=nbpoles; ii++) {
471       TM(ii).SetCols(ResPoles->Value(ii,2).XYZ(), 
472                      ResPoles->Value(ii,3).XYZ(),  
473                      ResPoles->Value(ii,4).XYZ());
474       TR(ii) = ResPoles->Value(ii,1);
475     }
476     GeomLib::TensorialProduct(BSurf, TM, TR,
477                               Approx.Knots()->Array1(), 
478                               Approx.Multiplicities()->Array1());
479
480     // Somme
481     TColgp_Array1OfPnt TPoles(1, nbpoles);
482     for (ii=1; ii<=nbpoles; ii++) {
483       TPoles(ii) = ResPoles->Value(ii,1);
484     }
485     Handle(Geom_BsplineCurve) BS = 
486       new (Geom_BsplineCurve) (Poles, 
487                                Approx.Knots()->Array1(), 
488                                Approx.Multiplicities()->Array1(),
489                                Approx.Degree());
490     for (ii=1; ii<=BSurf->NbVKnots(); ii++)
491       BS->InsertKnot( BSurf->VKnot(ii), 
492                      BSurf->VMultiplicity(ii), 
493                      Precision::Confusion());
494    TColgp_Array2OfPnt SurfPoles (1, BSurf->NbUPoles());
495    for (ii=1; 
496         
497 */
498     mySurface = BSurf;
499   }
500   return Ok;
501 }
502
503 //  Modified by skv - Thu Feb  5 18:05:03 2004 OCC5073 Begin
504 //  Conditions:
505 //     * theSec should be constant
506 //     * the type of section should be a line
507 //     * theLoc should represent a translation.
508
509 static Standard_Boolean IsSweepParallelSpine (const Handle(GeomFill_LocationLaw) &theLoc,
510                              const Handle(GeomFill_SectionLaw)  &theSec,
511                              const Standard_Real                 theTol)
512 {
513   // Get the first and last transformations of the location
514   Standard_Real aFirst;
515   Standard_Real aLast;
516   gp_Vec        VBegin;
517   gp_Vec        VEnd;
518   gp_Mat        M;
519   gp_GTrsf      GTfBegin;
520   gp_Trsf       TfBegin;
521   gp_GTrsf      GTfEnd;
522   gp_Trsf       TfEnd;
523
524   theLoc->GetDomain(aFirst, aLast);
525
526 // Get the first transformation
527   theLoc->D0(aFirst, M, VBegin);
528
529   GTfBegin.SetVectorialPart(M);
530   GTfBegin.SetTranslationPart(VBegin.XYZ());
531
532   TfBegin.SetValues(GTfBegin(1,1), GTfBegin(1,2), GTfBegin(1,3), GTfBegin(1,4),
533                     GTfBegin(2,1), GTfBegin(2,2), GTfBegin(2,3), GTfBegin(2,4),
534                     GTfBegin(3,1), GTfBegin(3,2), GTfBegin(3,3), GTfBegin(3,4),
535                     1.e-12, 1.e-14);
536
537 // Get the last transformation
538   theLoc->D0(aLast, M, VEnd);
539
540   GTfEnd.SetVectorialPart(M);
541   GTfEnd.SetTranslationPart(VEnd.XYZ());
542
543   TfEnd.SetValues(GTfEnd(1,1), GTfEnd(1,2), GTfEnd(1,3), GTfEnd(1,4),
544                   GTfEnd(2,1), GTfEnd(2,2), GTfEnd(2,3), GTfEnd(2,4),
545                   GTfEnd(3,1), GTfEnd(3,2), GTfEnd(3,3), GTfEnd(3,4),
546                   1.e-12, 1.e-14);
547
548   Handle(Geom_Surface) aSurf = theSec->BSplineSurface();
549   Standard_Real Umin;
550   Standard_Real Umax;
551   Standard_Real Vmin;
552   Standard_Real Vmax;
553
554   aSurf->Bounds(Umin, Umax, Vmin, Vmax);
555
556   // Get and transform the first section
557   Handle(Geom_Curve) FirstSection = theSec->ConstantSection();
558   GeomAdaptor_Curve  ACFirst(FirstSection);
559
560   Standard_Real UFirst = ACFirst.FirstParameter();
561   gp_Lin        L      = ACFirst.Line();
562
563   L.Transform(TfBegin);
564
565   // Get and transform the last section
566   Handle(Geom_Curve) aLastSection    = aSurf->VIso(Vmax);
567   Standard_Real      aFirstParameter = aLastSection->FirstParameter();
568   gp_Pnt             aPntLastSec     = aLastSection->Value(aFirstParameter);
569
570   aPntLastSec.Transform(TfEnd);
571
572   gp_Pnt        aPntFirstSec = ElCLib::Value( UFirst, L );
573   gp_Vec        aVecSec( aPntFirstSec, aPntLastSec );
574   gp_Vec        aVecSpine = VEnd - VBegin;
575
576   Standard_Boolean isParallel = aVecSec.IsParallel(aVecSpine, theTol);
577
578   return isParallel;
579 }
580 //  Modified by skv - Thu Feb  5 18:05:01 2004 OCC5073 End
581
582 //===============================================================
583 // Function : BuildKPart
584 // Purpose :
585 //===============================================================
586  Standard_Boolean GeomFill_Sweep::BuildKPart() 
587 {
588   Standard_Boolean Ok = Standard_False;
589   Standard_Boolean isUPeriodic = Standard_False;
590   Standard_Boolean isVPeriodic = Standard_False;
591   Standard_Boolean IsTrsf = Standard_True;
592
593   isUPeriodic = mySec->IsUPeriodic();
594   Handle(Geom_Surface) S;
595   GeomAbs_CurveType SectionType;
596   gp_Vec V;
597   gp_Mat M;
598   Standard_Real levier, error = 0 ;
599   Standard_Real UFirst=0, VFirst=First, ULast=0, VLast=Last;
600   Standard_Real Tol = Min (Tol3d, BoundTol);
601
602   // (1) Trajectoire Rectilignes -------------------------
603   if (myLoc->IsTranslation(error)) {
604     // Donne de la translation
605     gp_Vec DP, DS;
606     myLoc->D0(1, M, DS);
607     myLoc->D0(0, M, V);
608     DP = DS - V;
609     DP.Normalize();
610     gp_GTrsf Tf;
611     gp_Trsf Tf2;
612     Tf.SetVectorialPart(M);
613     Tf.SetTranslationPart(V.XYZ());
614     try { // Pas joli mais il n'y as pas d'autre moyens de tester SetValues
615       OCC_CATCH_SIGNALS
616       Tf2.SetValues(Tf(1,1), Tf(1,2), Tf(1,3), Tf(1,4),
617                     Tf(2,1), Tf(2,2), Tf(2,3), Tf(2,4),
618                     Tf(3,1), Tf(3,2), Tf(3,3), Tf(3,4),
619                     1.e-12, 1.e-14);
620     }
621     catch (Standard_ConstructionError) {
622       IsTrsf = Standard_False;
623     }
624     if (!IsTrsf) {
625       return Standard_False;
626     }
627
628     // (1.1) Cas Extrusion
629     if (mySec->IsConstant(error)) {
630       Handle(Geom_Curve) Section;
631       Section = mySec->ConstantSection();
632       GeomAdaptor_Curve AC(Section);
633       SectionType = AC.GetType();
634       UFirst = AC.FirstParameter();
635       ULast  = AC.LastParameter();
636   // (1.1.a) Cas Plan
637       if ( (SectionType == GeomAbs_Line) && IsTrsf) {
638 //  Modified by skv - Thu Feb  5 11:39:06 2004 OCC5073 Begin
639         if (!IsSweepParallelSpine(myLoc, mySec, Tol))
640           return Standard_False;
641 //  Modified by skv - Thu Feb  5 11:39:08 2004 OCC5073 End
642         gp_Lin L = AC.Line();
643         L.Transform(Tf2);
644         DS.SetXYZ(L.Position().Direction().XYZ());
645         DS.Normalize();
646         levier = Abs(DS.Dot(DP));
647         SError = error + levier * Abs(Last-First);
648         if (SError <= Tol) {
649           Ok = Standard_True;
650           gp_Ax2 AxisOfPlane (L.Location(), DS^DP, DS);
651           S = new (Geom_Plane) (AxisOfPlane);
652         }
653         else SError = 0.;
654       }
655     
656   // (1.1.b) Cas Cylindrique
657       if ( (SectionType == GeomAbs_Circle) && IsTrsf) {
658         gp_Circ C = AC.Circle();
659         C.Transform(Tf2);
660         
661         DS.SetXYZ (C.Position().Direction().XYZ());
662         DS.Normalize();
663         levier = Abs(DS.CrossMagnitude(DP)) * C.Radius();
664         SError = levier * Abs(Last - First);
665         if (SError <= Tol) {
666           Ok = Standard_True;
667           gp_Ax3 axe (C.Location(), DP, C.Position().XDirection());
668           S = new (Geom_CylindricalSurface) 
669             (axe, C.Radius());
670             if (C.Position().Direction().
671                 IsOpposite(axe.Direction(), 0.1) ) {
672               Standard_Real f, l;
673               // L'orientation parametrique est inversee
674               l = 2*M_PI - UFirst;
675               f = 2*M_PI - ULast;
676               UFirst = f;
677               ULast  = l;
678               isUReversed = Standard_True;
679             }
680         }
681         else SError = 0.;
682       }
683
684   // (1.1.c) C'est bien une extrusion
685       if (!Ok) {
686         if (IsTrsf) {
687           Section->Transform(Tf2);
688           S = new (Geom_SurfaceOfLinearExtrusion)
689             (Section, DP);
690           SError = 0.;
691           Ok = Standard_True;
692         }
693         else { // extrusion sur BSpline
694           
695         }
696       }
697     }
698   
699   // (1.2) Cas conique
700    else if (mySec->IsConicalLaw(error)) {
701
702      gp_Pnt P1, P2, Centre0, Centre1, Centre2;
703      gp_Vec dsection;
704      Handle(Geom_Curve) Section;
705      GeomAdaptor_Curve AC;
706      gp_Circ C;
707      Standard_Real R1, R2;
708
709
710      Section = mySec->CirclSection(SLast);
711      Section->Transform(Tf2);
712      Section->Translate(Last*DP);
713      AC.Load(Section);
714      C = AC.Circle();
715      Centre2 = C.Location();
716      AC.D1(0, P2, dsection);
717      R2 = C.Radius();
718     
719      Section = mySec->CirclSection(SFirst);
720      Section->Transform(Tf2);
721      Section->Translate(First*DP);
722      AC.Load(Section);
723      C =  AC.Circle();
724      Centre1 = C.Location();
725      P1 = AC.Value(0);
726      R1 = C.Radius();
727
728      Section = mySec->CirclSection(SFirst - First*(SLast-SFirst)/(Last-First));
729      Section->Transform(Tf2);
730      AC.Load(Section);
731      C =  AC.Circle();     
732      Centre0 = C.Location();
733      
734      Standard_Real Angle;
735      gp_Vec  N(Centre1, P1);
736      if (N.Magnitude() < 1.e-9) {
737        gp_Vec Bis(Centre2, P2);
738        N = Bis;
739      }
740      gp_Vec  L(P1, P2), Dir(Centre1,Centre2);
741
742      Angle = L.Angle(Dir);
743      if ((Angle > 0.01) && (Angle < M_PI/2-0.01)) {
744        if (R2<R1) Angle = -Angle;
745        SError = error;
746        gp_Ax3 Axis(Centre0, Dir, N);
747        S = new (Geom_ConicalSurface) 
748          (Axis, Angle, C.Radius());
749        // Calcul du glissement parametrique
750        VFirst = First / Cos(Angle);
751        VLast  = Last  / Cos(Angle);
752
753        // Bornes en U
754        UFirst = AC.FirstParameter();
755        ULast  = AC.LastParameter();
756        gp_Vec diso;
757        gp_Pnt pbis;
758        S->VIso(VLast)->D1(0, pbis, diso);
759        if (diso.Magnitude()>1.e-9 && dsection.Magnitude()>1.e-9)
760          isUReversed = diso.IsOpposite(dsection, 0.1);
761        if (isUReversed ) {
762          Standard_Real f, l;
763          // L'orientation parametrique est inversee
764          l = 2*M_PI - UFirst;
765          f = 2*M_PI - ULast;
766          UFirst = f;
767          ULast  = l;
768        } 
769
770        // C'est un cone
771        Ok = Standard_True;
772      }
773    }
774   }
775
776   // (2) Trajectoire Circulaire
777   if (myLoc->IsRotation(error)) {
778     if (mySec->IsConstant(error)) {
779       // La trajectoire
780       gp_Pnt Centre;
781       isVPeriodic = (Abs(Last-First -2*M_PI) < 1.e-15);
782       Standard_Real RotRadius;
783       gp_Vec DP, DS, DN;
784       myLoc->D0(0.1, M, DS);
785       myLoc->D0(0, M, V);
786       myLoc->Rotation(Centre);
787
788       DP = DS - V;
789       DS.SetXYZ(V.XYZ() - Centre.XYZ());   
790       RotRadius = DS.Magnitude();
791       if (RotRadius > 1.e-15) DS.Normalize();
792       else return Standard_False; // Pas de KPart, rotation degeneree
793       DN = DS ^ DP;
794       DN.Normalize();
795       DP = DN ^ DS;
796       DP.Normalize();
797
798       gp_GTrsf Tf;
799       gp_Trsf Tf2;
800       Tf.SetVectorialPart(M);
801       Tf.SetTranslationPart(V.XYZ());
802 //      try { // Pas joli mais il n'y as pas d'autre moyens de tester SetValues
803 //        OCC_CATCH_SIGNALS
804         Tf2.SetValues(Tf(1,1), Tf(1,2), Tf(1,3), Tf(1,4),
805                       Tf(2,1), Tf(2,2), Tf(2,3), Tf(2,4),
806                       Tf(3,1), Tf(3,2), Tf(3,3), Tf(3,4),
807                       1.e-14, 1.e-15);
808 //      }
809 //      catch (Standard_ConstructionError) {
810 //      IsTrsf = Standard_False;
811 //      }
812       // La section
813       Handle(Geom_Curve) Section;
814       Section = mySec->ConstantSection();
815       GeomAdaptor_Curve AC(Section);
816       SectionType = AC.GetType();
817       UFirst = AC.FirstParameter();
818       ULast  = AC.LastParameter();
819
820       // (2.1) Tore/Sphere ?
821       if ((SectionType == GeomAbs_Circle) && IsTrsf) {
822         gp_Circ C = AC.Circle();
823         Standard_Real Radius;
824         Standard_Boolean IsGoodSide = Standard_True;;
825         C.Transform(Tf2);
826         gp_Vec DC;
827         // On calcul le centre eventuel
828         DC.SetXYZ(C.Location().XYZ() - Centre.XYZ());
829         Centre.ChangeCoord() += (DC.Dot(DN))*DN.XYZ();
830         DC.SetXYZ(C.Location().XYZ() - Centre.XYZ());
831         Radius = DC.Magnitude(); //grand Rayon du tore
832         if ((Radius > Tol) && (DC.Dot(DS) < 0)) IsGoodSide = Standard_False;
833         if (Radius < Tol/100) DC = DS; // Pour definir le tore
834
835         // On verifie d'abord que le plan de la section est // a 
836         // l'axe de rotation
837         gp_Vec NC;
838         NC.SetXYZ (C.Position().Direction().XYZ());
839         NC.Normalize();
840         error =  Abs(NC.Dot(DN));
841         // Puis on evalue l'erreur commise sur la section, 
842         // en pivotant son plan ( pour contenir l'axe de rotation)
843         error += Abs(NC.Dot(DS));
844         error *= C.Radius();
845         if (error <= Tol) {
846           SError = error;
847           error += Radius + Abs(RotRadius - C.Radius())/2;
848           if (error <= Tol) {
849             // (2.1.a) Sphere
850             Standard_Real f = UFirst , l =  ULast;
851             SError = error;
852             Centre.BaryCenter(1.0, C.Location(), 1.0); 
853             gp_Ax3 AxisOfSphere(Centre, DN, DS);  
854             S = new (Geom_SphericalSurface) 
855               (AxisOfSphere, (RotRadius + C.Radius())/2 );
856             // Pour les spheres on ne peut pas controler le parametre
857             // V (donc U car  myExchUV = Standard_True)
858             // Il faut donc modifier UFirst, ULast...
859             if (C.Position().Direction().
860                 IsOpposite(AxisOfSphere.YDirection(), 0.1) ) {
861               // L'orientation parametrique est inversee
862               l = 2*M_PI - UFirst;
863               f = 2*M_PI - ULast;
864               isUReversed = Standard_True;
865             }
866             // On calcul le "glissement" parametrique.
867             Standard_Real rot; 
868             rot = C.Position().XDirection().AngleWithRef
869               (AxisOfSphere.XDirection(), AxisOfSphere.YDirection());
870             f -= rot;
871             l  -= rot;
872
873             if ( (f >= -M_PI/2) && (l <= M_PI/2)) {
874               Ok = Standard_True;
875               myExchUV = Standard_True;
876               UFirst = f;
877               ULast  = l;
878             }
879             else { // On restaure ce qu'il faut
880               isUReversed = Standard_False;
881             }
882           }
883           else if (IsGoodSide) {         
884             // (2.1.b) Tore
885             gp_Ax3 AxisOfTore(Centre, DN, DC);
886             S = new (Geom_ToroidalSurface) (AxisOfTore, 
887                    Radius , C.Radius());
888         
889             // Pour les tores on ne peut pas controler le parametre
890             // V (donc U car  myExchUV = Standard_True)
891             // Il faut donc modifier UFirst, ULast...
892             Handle(Geom_Circle) Iso;
893             Iso =  Handle(Geom_Circle)::DownCast(S->UIso(0.));
894             gp_Ax2 axeiso;
895             axeiso = Iso->Circ().Position();
896               
897             if (C.Position().Direction().
898                 IsOpposite(axeiso.Direction(), 0.1) ) {
899               Standard_Real f, l;
900               // L'orientation parametrique est inversee
901               l = 2*M_PI - UFirst;
902               f = 2*M_PI - ULast;
903               UFirst = f;
904               ULast  = l;
905               isUReversed = Standard_True;
906             }
907             // On calcul le "glissement" parametrique.
908             Standard_Real rot; 
909             rot = C.Position().XDirection().AngleWithRef
910               (axeiso.XDirection(), axeiso.Direction());
911             UFirst -= rot;
912             ULast  -= rot;
913
914             myExchUV = Standard_True;
915             // Attention l'arete de couture dans le cas periodique 
916             // n'est peut etre pas a la bonne place...
917             if (isUPeriodic && Abs(UFirst)>Precision::PConfusion()) 
918               isUPeriodic = Standard_False; //Pour trimmer la surface...
919             Ok = Standard_True;
920           }
921         }
922         else {
923           SError = 0.;
924         }
925       }
926       // (2.2) Cone / Cylindre
927       if ((SectionType == GeomAbs_Line) && IsTrsf)  {
928         gp_Lin L =  AC.Line();
929         L.Transform(Tf2);
930         gp_Vec DL;
931         DL.SetXYZ(L.Direction().XYZ());
932         levier = Max(Abs(AC.FirstParameter()), AC.LastParameter());
933         // si la line est ortogonale au cercle de rotation  
934         SError = error + levier * Abs(DL.Dot(DP));
935         if (SError <= Tol) {
936           Standard_Boolean reverse;
937           gp_Lin Dir(Centre, DN);
938           Standard_Real aux;
939           aux = DL.Dot(DN);
940           reverse = (aux < 0);  // On choisit ici le sens de parametrisation
941             
942           // Calcul du centre du vecteur supportant la "XDirection"
943           gp_Pnt CentreOfSurf;
944           gp_Vec O1O2(Centre, L.Location()), trans;
945           trans = DN;
946           trans *= DN.Dot(O1O2);
947           CentreOfSurf = Centre.Translated(trans);
948           DS.SetXYZ(L.Location().XYZ() - CentreOfSurf.XYZ());
949
950           error = SError;
951           error += (DL.XYZ()).CrossMagnitude(DN.XYZ())*levier;
952           if (error <= Tol) {
953             // (2.2.a) Cylindre
954             // si la line est orthogonale au plan de rotation
955             SError = error;
956             gp_Ax3 Axis(CentreOfSurf, Dir.Direction(), DS);
957             S = new (Geom_CylindricalSurface) 
958                     (Axis, L.Distance(CentreOfSurf));
959             Ok = Standard_True;
960             myExchUV = Standard_True;
961           }
962           else {
963             // On evalue l'angle du cone
964             Standard_Real Angle = Abs(Dir.Angle(L));
965             if (Angle > M_PI/2) Angle = M_PI -Angle;
966             if (reverse) Angle = -Angle;
967             aux = DS.Dot(DL);
968             if (aux < 0) {
969               Angle = - Angle;
970             }
971             if (Abs(Abs(Angle) - M_PI/2) > 0.01) {
972               // (2.2.b) Cone
973               // si les 2 droites ne sont pas orthogonales
974               Standard_Real Radius = CentreOfSurf.Distance(L.Location());
975               gp_Ax3 Axis(CentreOfSurf, Dir.Direction(), DS);
976               S = new (Geom_ConicalSurface) 
977                     (Axis, Angle, Radius);
978               myExchUV = Standard_True;
979               Ok = Standard_True;
980             }
981             else {
982               // On n'as pas conclue, on remet l'erreur a 0.
983               SError = 0.;
984             }
985           }
986           if (Ok && reverse) {
987             // On reverse le parametre
988             Standard_Real uf, ul;
989             Handle(Geom_Line) CL = new (Geom_Line)(L);
990             uf = CL->ReversedParameter(ULast);
991             ul = CL->ReversedParameter(UFirst);
992             UFirst = uf;
993             ULast = ul;
994             isUReversed = Standard_True;
995           }
996         }
997         else SError = 0.;
998       }
999   
1000       // (2.3) Revolution
1001       if (!Ok) {
1002         if (IsTrsf) { 
1003           Section->Transform(Tf2);
1004           gp_Ax1 Axis (Centre, DN);
1005           S = new (Geom_SurfaceOfRevolution)
1006             (Section, Axis);
1007           myExchUV = Standard_True;
1008           SError = 0.;
1009           Ok = Standard_True;
1010         }
1011       }
1012     }
1013   }
1014
1015
1016   if (Ok) { // On trimme la surface
1017     if (myExchUV) {
1018       Standard_Boolean b;
1019       b = isUPeriodic; isUPeriodic = isVPeriodic;  isVPeriodic = b;
1020       Standard_Real r;
1021       r = UFirst; UFirst = VFirst; VFirst = r;
1022       r = ULast; ULast = VLast; VLast = r;
1023     }
1024
1025     if (!isUPeriodic && !isVPeriodic)
1026       mySurface = new (Geom_RectangularTrimmedSurface)
1027         (S, UFirst, ULast, VFirst, VLast);
1028     else if (isUPeriodic) {
1029       if (isVPeriodic) mySurface = S;
1030       else mySurface = new (Geom_RectangularTrimmedSurface)
1031         (S, VFirst, VLast, Standard_False);
1032     }
1033     else
1034       mySurface = new (Geom_RectangularTrimmedSurface)
1035         (S,UFirst, ULast, Standard_True);
1036
1037 #if DEB
1038   if (isUPeriodic && !mySurface->IsUPeriodic()) 
1039     cout<<"Pb de periodicite en U" << endl;
1040   if (isUPeriodic && !mySurface->IsUClosed())
1041     cout<<"Pb de fermeture en U" << endl;
1042   if (isVPeriodic && !mySurface->IsVPeriodic()) 
1043     cout << "Pb de periodicite en V" << endl;
1044   if (isVPeriodic && !mySurface->IsVClosed())
1045     cout<<"Pb de fermeture en V" << endl;
1046 #endif
1047   }
1048
1049
1050   return Ok;
1051 }
1052
1053 //===============================================================
1054 // Function : IsDone
1055 // Purpose :
1056 //===============================================================
1057  Standard_Boolean GeomFill_Sweep::IsDone() const
1058 {
1059   return done;
1060 }
1061
1062 //===============================================================
1063 // Function :ErrorOnSurface
1064 // Purpose :
1065 //===============================================================
1066  Standard_Real GeomFill_Sweep::ErrorOnSurface() const
1067 {
1068   return SError;
1069 }
1070
1071 //===============================================================
1072 // Function ::ErrorOnRestriction
1073 // Purpose :
1074 //===============================================================
1075  void GeomFill_Sweep::ErrorOnRestriction(const Standard_Boolean IsFirst,
1076                                          Standard_Real& UError,
1077                                          Standard_Real& VError) const
1078 {
1079   Standard_Integer ind;
1080   if (IsFirst) ind=1;
1081   else ind = myCurve2d->Length();
1082
1083   UError =  CError->Value(1, ind);
1084   VError =  CError->Value(2, ind);
1085 }
1086
1087 //===============================================================
1088 // Function :ErrorOnTrace
1089 // Purpose :
1090 //===============================================================
1091  void GeomFill_Sweep::ErrorOnTrace(const Standard_Integer IndexOfTrace,
1092                                    Standard_Real& UError,
1093                                    Standard_Real& VError) const
1094 {
1095   Standard_Integer ind = IndexOfTrace+1;
1096   if (IndexOfTrace > myLoc->TraceNumber())
1097     Standard_OutOfRange::Raise(" GeomFill_Sweep::ErrorOnTrace");
1098
1099   UError =  CError->Value(1, ind);
1100   VError =  CError->Value(2, ind);
1101 }
1102
1103 //===============================================================
1104 // Function :Surface
1105 // Purpose :
1106 //===============================================================
1107  Handle(Geom_Surface) GeomFill_Sweep::Surface() const
1108 {
1109   return mySurface;
1110 }
1111
1112 //===============================================================
1113 // Function ::Restriction
1114 // Purpose :
1115 //===============================================================
1116  Handle(Geom2d_Curve) GeomFill_Sweep::Restriction(const Standard_Boolean IsFirst) const
1117 {
1118   if (IsFirst)
1119     return  myCurve2d->Value(1);
1120   return  myCurve2d->Value(myCurve2d->Length());
1121  
1122 }
1123
1124 //===============================================================
1125 // Function :
1126 // Purpose :
1127 //===============================================================
1128  Standard_Integer GeomFill_Sweep::NumberOfTrace() const
1129 {
1130   return myLoc->TraceNumber();
1131 }
1132
1133 //===============================================================
1134 // Function : 
1135 // Purpose :
1136 //===============================================================
1137  Handle(Geom2d_Curve) 
1138      GeomFill_Sweep::Trace(const Standard_Integer IndexOfTrace) const
1139 {
1140   Standard_Integer ind = IndexOfTrace+1;
1141   if (IndexOfTrace > myLoc->TraceNumber())
1142     Standard_OutOfRange::Raise(" GeomFill_Sweep::Trace");
1143   return  myCurve2d->Value(ind);  
1144 }