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