911e90fe7513afc062379396158d7c30a47cc0ad
[occt.git] / src / ApproxInt / ApproxInt_Approx.gxx
1 // Created on: 1993-03-30
2 // Created by: Laurent BUCHARD
3 // Copyright (c) 1993-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 #include <AppParCurves_Constraint.hxx>
18 #include <GeomAbs_SurfaceType.hxx>
19 #include <IntSurf_Quadric.hxx>
20 #include <gp_Trsf.hxx>
21 #include <gp_Trsf2d.hxx>
22 #include <IntSurf_PntOn2S.hxx>
23 #include <Precision.hxx>
24 #include <ApproxInt_KnotTools.hxx>
25
26 const Standard_Integer LimRajout = 5;
27 const Standard_Integer NbPntMaxDecoupage = 30;
28 const Standard_Real RatioTol = 1.5;
29
30 //=======================================================================
31 //function : MINABS3
32 //purpose  : Compute minimal absolute distance to 0 from 3 values.
33 //=======================================================================
34 static Standard_Real MINABS3(Standard_Real a, Standard_Real b,Standard_Real c)
35 {
36   if(a<0.0) a=-a;
37   if(b<0.0) b=-b;
38   if(c<0.0) c=-c;
39   if(a>c) a=c;
40   if(a>b) a=b;
41   return(a);
42 }
43
44 //=======================================================================
45 //function : MINABS4
46 //purpose  : Compute minimal absolute distance to 0 from 4 values.
47 //=======================================================================
48 static Standard_Real MINABS4(Standard_Real a, Standard_Real b,Standard_Real c,Standard_Real d)
49 {
50   if(a<0.0) a=-a;
51   if(b<0.0) b=-b;
52   if(c<0.0) c=-c;
53   if(d<0.0) d=-d;
54   if(a>c) a=c;
55   if(a>b) a=b;
56   if(a>d) a=d;
57   return(a);
58 }
59
60 //=======================================================================
61 //function : ComputeTrsf3d
62 //purpose  : 
63 //=======================================================================
64 static void ComputeTrsf3d(const Handle(TheWLine)& theline,
65                           Standard_Real& Xo, Standard_Real& Ax,
66                           Standard_Real& Yo, Standard_Real& Ay,
67                           Standard_Real& Zo, Standard_Real& Az)
68 {
69     Standard_Integer nbp = theline->NbPnts();
70     Standard_Real z0,z1,x0,x1,y0,y1;
71     z0=y0=x0=RealLast();
72     z1=y1=x1=RealFirst();
73     for(Standard_Integer i=1;i<=nbp;i++) { 
74       const gp_Pnt& P = theline->Point(i).Value();
75       Standard_Real  X = P.X();
76       Standard_Real  Y = P.Y();
77       Standard_Real  Z = P.Z();
78       if(X<x0) x0=X;
79       if(X>x1) x1=X;
80       if(Y<y0) y0=Y;
81       if(Y>y1) y1=Y;
82       if(Z<z0) z0=Z;
83       if(Z>z1) z1=Z;
84     }
85     Standard_Real dx = x1-x0;
86     Standard_Real dy = y1-y0;
87     Standard_Real dz = z1-z0;
88     Standard_Real MaxD = dx; 
89     if(MaxD < dy) MaxD=dy;
90     if(MaxD < dz) MaxD=dz;
91     Standard_Real MaxDF = 0.01*MaxD;
92
93     if(MaxDF<1e-12) 
94       MaxDF=1.0;
95
96     if(dx > MaxDF) { Ax = 1.0 / dx;    Xo = -Ax * x0;  }
97     else {     Ax = 1.0/( MaxDF) ; Xo = -Ax*x0;   }
98     if(dy > MaxDF) { Ay = 1.0 / dy;    Yo = -Ay * y0;  }
99     else {     Ay = 1.0/( MaxDF); Yo = -Ay*y0;   }
100     if(dz > MaxDF) { Az = 1.0 / dz;    Zo = -Az * z0;   }
101     else {     Az = 1.0/(MaxDF); Zo = -Az*z0;   } 
102 }
103
104 //=======================================================================
105 //function : ComputeTrsf2d
106 //purpose  :
107 //=======================================================================
108 static void ComputeTrsf2d(const Handle(TheWLine)& theline,
109                           Standard_Real& Uo, Standard_Real& Au,
110                           Standard_Real& Vo, Standard_Real& Av,
111                           const Standard_Boolean onFirst,
112                           const Standard_Real UVResRatio = 1.0)
113 {
114     Standard_Integer nbp = theline->NbPnts();
115     Standard_Real u0,u1,v0,v1;
116     u0 = v0 = RealLast();
117     u1 = v1 = RealFirst();
118     // pointer to a member-function
119     void (IntSurf_PntOn2S::* pfunc)(Standard_Real&,Standard_Real&) const;
120     if (onFirst)
121       pfunc = &IntSurf_PntOn2S::ParametersOnS1;
122     else
123       pfunc = &IntSurf_PntOn2S::ParametersOnS2;
124     for(Standard_Integer i=1;i<=nbp;i++) { 
125       const IntSurf_PntOn2S&  POn2S = theline->Point(i);
126       Standard_Real  U,V;
127       (POn2S.*pfunc)(U,V);
128       if(U<u0) u0=U;
129       if(U>u1) u1=U;
130       if(V<v0) v0=V;
131       if(V>v1) v1=V;
132     }
133
134     Standard_Real du = (u1-u0);
135     Standard_Real dv = (v1-v0);
136
137     if (UVResRatio > 1.)
138       du *= UVResRatio;
139     else if (UVResRatio < 1.)
140       dv /= UVResRatio;
141
142     Standard_Real MaxUV=du;
143     if(MaxUV<dv) MaxUV=dv;
144
145     Standard_Real MaxUVF=0.01*MaxUV;
146
147     //-- lbr le 22 fev 99 (FPE) 
148     if(MaxUVF<1e-12) 
149       MaxUVF=1.0;
150
151     if(du > MaxUVF) { Au = 1.0 / du;    Uo = -Au * u0;  }
152     else {     Au = 1.0/(MaxUVF); Uo = -Au*u0;  }
153     if(dv > MaxUVF) { Av = 1.0 / dv;    Vo = -Av * v0;  }
154     else {     Av = 1.0/(MaxUVF); Vo = -Av*v0;  }
155 }
156
157 //=======================================================================
158 //function : Parameters
159 //purpose  :
160 //=======================================================================
161 static void Parameters(const ApproxInt_TheMultiLine& Line,
162                        const Standard_Integer firstP,
163                        const Standard_Integer lastP,
164                        const Approx_ParametrizationType Par,
165                        math_Vector& TheParameters)
166 {
167   Standard_Integer i, j, nbP2d, nbP3d;
168   Standard_Real dist;
169   gp_Pnt P1, P2;
170   gp_Pnt2d P12d, P22d;
171
172   if (Par == Approx_ChordLength || Par == Approx_Centripetal) {
173     nbP3d = ApproxInt_TheMultiLineTool::NbP3d(Line);
174     nbP2d = ApproxInt_TheMultiLineTool::NbP2d(Line);
175     Standard_Integer mynbP3d=nbP3d, mynbP2d=nbP2d;
176     if (nbP3d == 0) mynbP3d = 1;
177     if (nbP2d == 0) mynbP2d = 1;
178
179     TheParameters(firstP) = 0.0;
180     dist = 0.0;
181     TColgp_Array1OfPnt tabP(1, mynbP3d);
182     TColgp_Array1OfPnt tabPP(1, mynbP3d);
183     TColgp_Array1OfPnt2d tabP2d(1, mynbP2d);
184     TColgp_Array1OfPnt2d tabPP2d(1, mynbP2d);
185
186     for (i = firstP+1; i <= lastP; i++) {
187       if (nbP3d != 0 && nbP2d != 0) ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP, tabP2d);
188       else if (nbP2d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP2d);
189       else if (nbP3d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i-1, tabP);
190
191       if (nbP3d != 0 && nbP2d != 0) ApproxInt_TheMultiLineTool::Value(Line, i, tabPP, tabPP2d);
192       else if (nbP2d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i, tabPP2d);
193       else if (nbP3d != 0)          ApproxInt_TheMultiLineTool::Value(Line, i, tabPP);
194       dist = 0;
195       for (j = 1; j <= nbP3d; j++) {
196         P1 = tabP(j);
197         P2 = tabPP(j);
198         dist += P2.Distance(P1);
199       }
200       for (j = 1; j <= nbP2d; j++) {
201         P12d = tabP2d(j);
202         P22d = tabPP2d(j);
203         dist += P22d.Distance(P12d);
204       }
205       if(Par == Approx_ChordLength)
206         TheParameters(i) = TheParameters(i-1) + dist;
207       else {// Par == Approx_Centripetal
208         TheParameters(i) = TheParameters(i-1) + Sqrt(dist);
209       }
210     }
211     for (i = firstP; i <= lastP; i++) TheParameters(i) /= TheParameters(lastP);
212   }
213   else {
214     for (i = firstP; i <= lastP; i++) {
215       TheParameters(i) = (Standard_Real(i)-firstP)/
216         (Standard_Real(lastP)-Standard_Real(firstP));
217     }
218   }
219 }
220
221 //=======================================================================
222 //function : ApproxInt_Approx
223 //purpose  : Constructor.
224 //=======================================================================
225 ApproxInt_Approx::ApproxInt_Approx()
226 : myComputeLine(4, 8, 0.001, 0.001, 5, Standard_True),
227   myComputeLineBezier(4, 8, 0.001, 0.001, 5, Standard_True)
228 {
229   myComputeLine.SetContinuity(2);
230   myData.myBezierApprox = Standard_True;
231
232   myRelativeTol = Standard_True;
233   myNbPntMax = NbPntMaxDecoupage;
234   myData.myMinFactorXYZ = 0.0;
235   myData.myMinFactorUV  = 0.0;
236   myTolReached3d = myTolReached2d = 0.0;
237   myUVRes1 = myUVRes2 = 1.0;
238 }
239
240 //=======================================================================
241 //function : Perform
242 //purpose  : Build without surfaces information.
243 //=======================================================================
244 void ApproxInt_Approx::Perform(const Handle(TheWLine)& theline,
245                                const Standard_Boolean ApproxXYZ,
246                                const Standard_Boolean ApproxU1V1,
247                                const Standard_Boolean ApproxU2V2,
248                                const Standard_Integer indicemin,
249                                const Standard_Integer indicemax)
250 {
251   // Prepare DS.
252   prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
253
254   Standard_Integer nbpntbez = indicemax-indicemin;
255   if(nbpntbez < LimRajout) 
256     myData.myBezierApprox = Standard_False;
257   else 
258     myData.myBezierApprox = Standard_True;
259
260   // Fill data structure.
261   fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
262
263   // Build knots.
264   buildKnots(theline, NULL);
265
266   Standard_Boolean cut = Standard_True;
267   if(myRelativeTol==Standard_False)
268   {
269     myComputeLine.Init(myDegMin,
270       myDegMax,
271       myTol3d*myData.myMinFactorXYZ,
272       myTol2d*myData.myMinFactorUV,
273       myNbIterMax,
274       cut,
275       myData.parametrization);
276     myComputeLineBezier.Init(myDegMin,
277       myDegMax,
278       myTol3d*myData.myMinFactorXYZ,
279       myTol2d*myData.myMinFactorUV,
280       myNbIterMax,
281       cut,
282       myData.parametrization);
283   }
284
285   buildCurve(theline, NULL);
286 }
287
288 //=======================================================================
289 //function : Perform
290 //purpose  : Param-Param perform.
291 //=======================================================================
292 void ApproxInt_Approx::Perform(const ThePSurface& Surf1,
293                                const ThePSurface& Surf2,
294                                const Handle(TheWLine)& theline,
295                                const Standard_Boolean ApproxXYZ,
296                                const Standard_Boolean ApproxU1V1,
297                                const Standard_Boolean ApproxU2V2,
298                                const Standard_Integer indicemin,
299                                const Standard_Integer indicemax) 
300 {
301
302   GeomAbs_SurfaceType typeS1 = ThePSurfaceTool::GetType(Surf1);
303   GeomAbs_SurfaceType typeS2 = ThePSurfaceTool::GetType(Surf2);
304   if ((typeS1 != GeomAbs_Plane    &&
305        typeS1 != GeomAbs_Cylinder &&
306        typeS1 != GeomAbs_Sphere   &&
307        typeS1 != GeomAbs_Cone)    &&
308       (typeS2 != GeomAbs_Plane    &&
309        typeS2 != GeomAbs_Cylinder &&
310        typeS2 != GeomAbs_Sphere   &&
311        typeS2 != GeomAbs_Cone))
312   {
313     // Prepare DS.
314     prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
315
316     // Non-analytical case: Param-Param perform.
317     ApproxInt_ThePrmPrmSvSurfaces myPrmPrmSvSurfaces(Surf1,Surf2);
318
319     Standard_Integer nbpntbez = indicemax-indicemin;
320     if(nbpntbez < LimRajout) 
321       myData.myBezierApprox = Standard_False;
322     else 
323       myData.myBezierApprox = Standard_True;
324
325     Standard_Boolean cut = Standard_True;
326     if(nbpntbez < LimRajout)
327     {
328       cut = Standard_False;
329     }
330
331     Standard_Real aS1URes = ThePSurfaceTool::UResolution(Surf1, 1.0),
332                   aS1VRes = ThePSurfaceTool::VResolution(Surf1, 1.0),
333                   aS2URes = ThePSurfaceTool::UResolution(Surf2, 1.0),
334                   aS2VRes = ThePSurfaceTool::VResolution(Surf2, 1.0);
335     if(ApproxU1V1)
336       myUVRes1 = aS1URes / aS1VRes;
337     if(ApproxU2V2)
338       myUVRes2 = aS2URes / aS2VRes;
339
340     // Fill data structure.
341     fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
342
343     // Build knots.
344     Standard_Address ptrsvsurf = &myPrmPrmSvSurfaces;
345     buildKnots(theline, ptrsvsurf);
346
347     if(myRelativeTol==Standard_False)
348     {
349       myComputeLine.Init(myDegMin,
350         myDegMax,
351         myTol3d*myData.myMinFactorXYZ,
352         myTol2d*myData.myMinFactorUV,
353         myNbIterMax,
354         cut,
355         myData.parametrization);
356       myComputeLineBezier.Init(myDegMin,
357         myDegMax,
358         myTol3d*myData.myMinFactorXYZ,
359         myTol2d*myData.myMinFactorUV,
360         myNbIterMax,
361         cut,
362         myData.parametrization);
363     }
364     else
365     {
366       myComputeLine.Init(myDegMin,
367         myDegMax,
368         myTol3d,
369         myTol2d,
370         myNbIterMax,
371         cut,
372         myData.parametrization);
373       myComputeLineBezier.Init(myDegMin,
374         myDegMax,
375         myTol3d,
376         myTol2d,
377         myNbIterMax,
378         cut,
379         myData.parametrization);
380     }
381
382     buildCurve(theline, ptrsvsurf);
383   }
384   else
385   {
386     IntSurf_Quadric Quad;
387     Standard_Boolean SecondIsImplicit=Standard_False;
388     switch (typeS1)
389     {
390
391     case GeomAbs_Plane:
392       Quad.SetValue(ThePSurfaceTool::Plane(Surf1));
393       break;
394
395     case GeomAbs_Cylinder:
396       Quad.SetValue(ThePSurfaceTool::Cylinder(Surf1));
397       break;
398
399     case GeomAbs_Sphere:
400       Quad.SetValue(ThePSurfaceTool::Sphere(Surf1));
401       break;
402
403     case GeomAbs_Cone:
404       Quad.SetValue(ThePSurfaceTool::Cone(Surf1));
405       break;
406
407     default:
408       {
409         SecondIsImplicit = Standard_True;
410         switch (typeS2)
411         {
412         case GeomAbs_Plane:
413           Quad.SetValue(ThePSurfaceTool::Plane(Surf2));
414           break;
415
416         case GeomAbs_Cylinder:
417           Quad.SetValue(ThePSurfaceTool::Cylinder(Surf2));
418           break;
419
420         case GeomAbs_Sphere:
421           Quad.SetValue(ThePSurfaceTool::Sphere(Surf2));
422           break;
423
424         case GeomAbs_Cone:
425           Quad.SetValue(ThePSurfaceTool::Cone(Surf2));
426           break;
427
428         default:
429           break;
430         }
431       }
432       break;
433     } 
434     if(SecondIsImplicit)
435     {
436       Perform(Surf1,Quad,theline,ApproxXYZ,ApproxU1V1,ApproxU2V2,indicemin,indicemax);
437     }
438     else
439     {
440       Perform(Quad,Surf2,theline,ApproxXYZ,ApproxU1V1,ApproxU2V2,indicemin,indicemax);
441     }
442   }
443 }
444
445 //=======================================================================
446 //function : SetParameters
447 //purpose  :
448 //=======================================================================
449 void ApproxInt_Approx::SetParameters(const Standard_Real     Tol3d,
450                                      const Standard_Real     Tol2d,
451                                      const Standard_Integer  DegMin,
452                                      const Standard_Integer  DegMax,
453                                      const Standard_Integer  NbIterMax,
454                                      const Standard_Boolean  ApproxWithTangency,
455                                      const Approx_ParametrizationType Parametrization)
456 {
457     myWithTangency = ApproxWithTangency;
458     myTol3d        = Tol3d / RatioTol;
459     myTol2d        = Tol2d / RatioTol;
460     myDegMin       = DegMin;
461     myDegMax       = DegMax;
462     myNbIterMax    = NbIterMax;
463     myComputeLine.Init(myDegMin,
464       myDegMax,
465       myTol3d,
466       myTol2d,
467       myNbIterMax,
468       Standard_True,
469       Parametrization);
470
471     if(!ApproxWithTangency) { 
472       myComputeLine.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
473     }
474     myComputeLineBezier.Init(myDegMin,
475       myDegMax,
476       myTol3d,
477       myTol2d,
478       myNbIterMax,
479       Standard_True,
480       Parametrization);
481     if(!ApproxWithTangency)
482     {
483       myComputeLineBezier.SetConstraints(AppParCurves_PassPoint,AppParCurves_PassPoint);
484     }
485     myData.myBezierApprox = Standard_True;
486 }
487
488 //=======================================================================
489 //function : SetParameters
490 //purpose  : Set parameters with RelativeTol flag and NbPntMax value.
491 //=======================================================================
492 void ApproxInt_Approx::SetParameters(const Standard_Real     Tol3d,
493                                      const Standard_Real     Tol2d,
494                                      const Standard_Boolean  RelativeTol,
495                                      const Standard_Integer  DegMin,
496                                      const Standard_Integer  DegMax,
497                                      const Standard_Integer  NbIterMax,
498                                      const Standard_Integer  NbPntMax,
499                                      const Standard_Boolean  ApproxWithTangency,
500                                      const Approx_ParametrizationType Parametrization) 
501 {
502   myNbPntMax = NbPntMax ;
503   myRelativeTol = RelativeTol ;
504   SetParameters (Tol3d, Tol2d, DegMin, DegMax, NbIterMax, ApproxWithTangency, Parametrization) ;
505 }
506
507 //=======================================================================
508 //function : Perform
509 //purpose  : Param-Analytic perform.
510 //=======================================================================
511 void ApproxInt_Approx::Perform(const ThePSurface& PSurf,
512                                const TheISurface& ISurf,
513                                const Handle(TheWLine)& theline,
514                                const Standard_Boolean ApproxXYZ,
515                                const Standard_Boolean ApproxU1V1,
516                                const Standard_Boolean ApproxU2V2,
517                                const Standard_Integer indicemin,
518                                const Standard_Integer indicemax)
519 {
520   // Prepare DS.
521   prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
522
523   // Non-analytical case: Param-Analytic perform.
524   ApproxInt_TheImpPrmSvSurfaces myImpPrmSvSurfaces(PSurf,ISurf);
525
526   Standard_Integer nbpntbez = indicemax-indicemin;
527   if(nbpntbez < LimRajout) 
528     myData.myBezierApprox = Standard_False;
529   else 
530     myData.myBezierApprox = Standard_True;
531
532   Standard_Boolean cut = Standard_True;
533   if(nbpntbez < LimRajout)
534   {
535     cut = Standard_False;
536   }
537
538   Standard_Real aS1URes = ThePSurfaceTool::UResolution(PSurf, 1.0),
539                 aS1VRes = ThePSurfaceTool::VResolution(PSurf, 1.0);
540   if(ApproxU1V1)
541     myUVRes1 = aS1URes / aS1VRes;
542
543   // Fill data structure.
544   fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
545
546   // Build knots.
547   Standard_Address ptrsvsurf = &myImpPrmSvSurfaces;
548   buildKnots(theline, ptrsvsurf);
549
550   if(myRelativeTol==Standard_False)
551   {
552     myComputeLine.Init(myDegMin,
553       myDegMax,
554       myTol3d*myData.myMinFactorXYZ,
555       myTol2d*myData.myMinFactorUV,
556       myNbIterMax,
557       cut,
558       myData.parametrization);
559     myComputeLineBezier.Init(myDegMin,
560       myDegMax,
561       myTol3d*myData.myMinFactorXYZ,
562       myTol2d*myData.myMinFactorUV,
563       myNbIterMax,
564       cut,
565       myData.parametrization);
566   }
567   else
568   {
569     myComputeLine.Init(myDegMin,
570       myDegMax,
571       myTol3d,
572       myTol2d,
573       myNbIterMax,
574       cut,
575       myData.parametrization);
576     myComputeLineBezier.Init(myDegMin,
577       myDegMax,
578       myTol3d,
579       myTol2d,
580       myNbIterMax,
581       cut,
582       myData.parametrization);
583   }
584
585   buildCurve(theline, ptrsvsurf);
586 }
587
588 //=======================================================================
589 //function : Perform
590 //purpose  : Analytic-Param perform.
591 //=======================================================================
592 void ApproxInt_Approx::Perform(const TheISurface& ISurf,
593                                const ThePSurface& PSurf,
594                                const Handle(TheWLine)& theline,
595                                const Standard_Boolean ApproxXYZ,
596                                const Standard_Boolean ApproxU1V1,
597                                const Standard_Boolean ApproxU2V2,
598                                const Standard_Integer indicemin,
599                                const Standard_Integer indicemax)
600 {
601   // Prepare DS.
602   prepareDS(ApproxXYZ, ApproxU1V1, ApproxU2V2, indicemin, indicemax);
603
604   // Non-analytical case: Analytic-Param perform.
605   ApproxInt_TheImpPrmSvSurfaces myImpPrmSvSurfaces(ISurf, PSurf);
606
607   Standard_Integer nbpntbez = indicemax-indicemin;
608   if(nbpntbez < LimRajout) 
609     myData.myBezierApprox = Standard_False;
610   else 
611     myData.myBezierApprox = Standard_True;
612
613   Standard_Boolean cut = Standard_True;
614   if(nbpntbez < LimRajout)
615   {
616     cut = Standard_False;
617   }
618
619   Standard_Real aS2URes = ThePSurfaceTool::UResolution(PSurf, 1.0),
620                 aS2VRes = ThePSurfaceTool::VResolution(PSurf, 1.0);
621   if(ApproxU2V2)
622     myUVRes2 = aS2URes / aS2VRes;
623
624   // Fill data structure.
625   fillData(theline, ApproxXYZ, ApproxU1V1, ApproxU2V2);
626
627   // Build knots.
628   Standard_Address ptrsvsurf = &myImpPrmSvSurfaces;
629   buildKnots(theline, ptrsvsurf);
630
631   if(myRelativeTol==Standard_False)
632   {
633     myComputeLine.Init(myDegMin,
634       myDegMax,
635       myTol3d*myData.myMinFactorXYZ,
636       myTol2d*myData.myMinFactorUV,
637       myNbIterMax,
638       cut,
639       myData.parametrization);
640     myComputeLineBezier.Init(myDegMin,
641       myDegMax,
642       myTol3d*myData.myMinFactorXYZ,
643       myTol2d*myData.myMinFactorUV,
644       myNbIterMax,
645       cut,
646       myData.parametrization);
647   }
648   else
649   {
650     myComputeLine.Init(myDegMin,
651       myDegMax,
652       myTol3d,
653       myTol2d,
654       myNbIterMax,
655       cut,
656       myData.parametrization);
657     myComputeLineBezier.Init(myDegMin,
658       myDegMax,
659       myTol3d,
660       myTol2d,
661       myNbIterMax,
662       cut,
663       myData.parametrization);
664   }
665
666   buildCurve(theline, ptrsvsurf);
667 }
668
669 //=======================================================================
670 //function : NbMultiCurves
671 //purpose  :
672 //=======================================================================
673 Standard_Integer ApproxInt_Approx::NbMultiCurves() const
674 {
675   return 1;
676 }
677
678 //=======================================================================
679 //function : UpdateTolReached
680 //purpose  :
681 //=======================================================================
682 void ApproxInt_Approx::UpdateTolReached()
683 {
684   if (myData.myBezierApprox)
685   {
686     Standard_Integer ICur;
687     Standard_Integer NbCurves = myComputeLineBezier.NbMultiCurves();
688     for (ICur = 1 ; ICur <= NbCurves ; ICur++)
689     {
690       Standard_Real Tol3D, Tol2D ;
691       myComputeLineBezier.Error (ICur, Tol3D, Tol2D) ;
692       myTolReached3d = Max(myTolReached3d, Tol3D);
693       myTolReached2d = Max(myTolReached2d, Tol2D);
694     }
695   }
696   else
697   {
698     myComputeLine.Error (myTolReached3d, myTolReached2d);
699   }
700 }
701
702 //=======================================================================
703 //function : TolReached3d
704 //purpose  :
705 //=======================================================================
706 Standard_Real ApproxInt_Approx::TolReached3d() const
707 {
708   Standard_Real TheTol3D = RatioTol * myTolReached3d ;
709
710   if (myData.myMinFactorXYZ>1.5e-7)
711     TheTol3D = TheTol3D / myData.myMinFactorXYZ;
712
713   //cout << "Tol 3D: " << TheTol3D << endl;
714   return TheTol3D ;
715 }
716
717 //=======================================================================
718 //function : TolReached2d
719 //purpose  :
720 //=======================================================================
721 Standard_Real ApproxInt_Approx::TolReached2d() const
722 {
723   Standard_Real TheTol2D = RatioTol * myTolReached2d ;
724
725   if (myData.myMinFactorUV>1.5e-7)
726     TheTol2D = TheTol2D / myData.myMinFactorUV;
727
728   //cout << "Tol 2D: " << TheTol2D << endl;
729   return TheTol2D ;
730 }
731
732 //=======================================================================
733 //function : IsDone
734 //purpose  :
735 //=======================================================================
736 Standard_Boolean ApproxInt_Approx::IsDone() const
737 {
738   if(myData.myBezierApprox)
739   { 
740     return(myComputeLineBezier.NbMultiCurves() > 0);
741     //-- Lorsque la tolerance n est pas atteinte et qu il 
742     //-- faudrait rajouter des points sur la ligne
743     //-- les approx sortent avec la meilleure tolerance
744     //-- atteinte.  ( Pas de rajout de points ds cette version)
745     //-- return(myTolReached);
746   }
747   else
748   {
749     return(myComputeLine.IsToleranceReached());
750   }
751 }
752
753 //=======================================================================
754 //function : Value
755 //purpose  :
756 //=======================================================================
757 const AppParCurves_MultiBSpCurve& ApproxInt_Approx::Value(const Standard_Integer ) const
758 {
759   if(myData.myBezierApprox)
760   {
761     return(myBezToBSpl.Value());
762   }
763   else
764   {
765     return(myComputeLine.Value());
766   }
767 }
768
769 //=======================================================================
770 //function : fillData
771 //purpose  : Fill ApproxInt data structure.
772 //=======================================================================
773 void ApproxInt_Approx::fillData(const Handle(TheWLine)& theline,
774                                 const Standard_Boolean ApproxXYZ,
775                                 const Standard_Boolean ApproxU1V1,
776                                 const Standard_Boolean ApproxU2V2)
777 {
778   if(ApproxXYZ)
779   {
780     ComputeTrsf3d(theline, myData.Xo, myData.Ax, myData.Yo, myData.Ay, myData.Zo, myData.Az);
781   }
782   else
783   {
784     myData.Xo = myData.Yo = myData.Zo = 0.0;
785     myData.Ax = myData.Ay = myData.Az = 1.0;
786   }
787
788   if(ApproxU1V1)
789   {
790     ComputeTrsf2d(theline, myData.U1o, myData.A1u, myData.V1o, myData.A1v,Standard_True, myUVRes1);
791   }
792   else
793   {
794     myData.U1o = myData.V1o = 0.0;
795     myData.A1u = myData.A1v = 1.0;
796   }
797
798   if(ApproxU2V2)
799   {
800     ComputeTrsf2d(theline, myData.U2o, myData.A2u, myData.V2o, myData.A2v, Standard_False, myUVRes2);
801   }
802   else
803   {
804     myData.U2o = myData.V2o = 0.0;
805     myData.A2u = myData.A2v = 1.0;
806   }
807
808   Standard_Real A3d = MINABS3(myData.Ax, myData.Ay, myData.Az);
809   if((A3d < myData.myMinFactorXYZ) || (myData.myMinFactorXYZ == 0.0))
810   {
811     myData.myMinFactorXYZ = A3d;
812   }
813
814   Standard_Real A2d = MINABS4(myData.A1u, myData.A1v, myData.A2u, myData.A2v);
815   if((A2d < myData.myMinFactorUV) || (myData.myMinFactorUV == 0.0))
816   {
817     myData.myMinFactorUV = A2d;
818   }
819 }
820
821 //=======================================================================
822 //function : prepareDS
823 //purpose  :
824 //=======================================================================
825 void ApproxInt_Approx::prepareDS(const Standard_Boolean theApproxXYZ,
826                                  const Standard_Boolean theApproxU1V1,
827                                  const Standard_Boolean theApproxU2V2,
828                                  const Standard_Integer theIndicemin,
829                                  const Standard_Integer theIndicemax)
830 {
831   myData.myMinFactorXYZ = 0.0;
832   myData.myMinFactorUV  = 0.0;
833   myTolReached3d = myTolReached2d = 0.0;
834   myUVRes1 = myUVRes2 = 1.0;
835   myData.ApproxU1V1 = theApproxU1V1;
836   myData.ApproxU2V2 = theApproxU2V2;
837   myData.ApproxXYZ = theApproxXYZ;
838   myData.indicemin = theIndicemin;
839   myData.indicemax = theIndicemax;
840   myData.nbpntmax = myNbPntMax;
841
842   Approx_ParametrizationType parametrization;
843   myComputeLineBezier.Parametrization(parametrization);
844   myData.parametrization = parametrization;
845 }
846
847 //=======================================================================
848 //function : buildKnots
849 //purpose  :
850 //=======================================================================
851 void ApproxInt_Approx::buildKnots(const Handle(TheWLine)& theline,
852                                   const Standard_Address thePtrSVSurf)
853 {
854   myKnots.Clear();
855   if(myData.myBezierApprox)
856   {
857     ApproxInt_TheMultiLine aTestLine(theline,
858       thePtrSVSurf,
859       ((myData.ApproxXYZ)? 1 : 0),
860       ((myData.ApproxU1V1)? 1: 0) + ((myData.ApproxU2V2)? 1: 0),
861       myData.Xo, myData.Ax, myData.Yo, myData.Ay, myData.Zo, myData.Az,
862       myData.U1o, myData.A1u, 
863       myData.V1o, myData.A1v,
864       myData.U2o, myData.A2u,
865       myData.V2o ,myData.A2v,
866       myData.ApproxU1V1,
867       myData.indicemin,
868       myData.indicemax);
869
870     Standard_Integer nbp3d = aTestLine.NbP3d();
871     Standard_Integer nbp2d = aTestLine.NbP2d();
872     TColgp_Array1OfPnt aTabPnt3d(1, Max(1, nbp3d));
873     TColgp_Array1OfPnt2d aTabPnt2d(1, Max(1, nbp2d));
874     TColgp_Array1OfPnt aPntXYZ(myData.indicemin, myData.indicemax);
875     TColgp_Array1OfPnt2d aPntU1V1(myData.indicemin, myData.indicemax);
876     TColgp_Array1OfPnt2d aPntU2V2(myData.indicemin, myData.indicemax);
877     Standard_Integer i;
878     for(i = myData.indicemin; i <= myData.indicemax; ++i)
879     {
880       if (nbp3d != 0 && nbp2d != 0) aTestLine.Value(i, aTabPnt3d, aTabPnt2d);
881       else if (nbp2d != 0)          aTestLine.Value(i, aTabPnt2d);
882       else if (nbp3d != 0)          aTestLine.Value(i, aTabPnt3d);
883       //
884       if(nbp3d > 0)
885       {
886         aPntXYZ(i) = aTabPnt3d(1);
887       }
888       if(nbp2d > 1)
889       {
890         aPntU1V1(i) = aTabPnt2d(1);
891         aPntU2V2(i) = aTabPnt2d(2);
892       }
893       else if(nbp2d > 0)
894       {
895         if(myData.ApproxU1V1)
896         {
897           aPntU1V1(i) = aTabPnt2d(1);
898         }
899         else
900         {
901           aPntU2V2(i) = aTabPnt2d(1);
902         }
903       }
904     }
905
906     Standard_Integer aMinNbPnts = myData.nbpntmax;
907
908     // Expected parametrization.
909     math_Vector aPars(myData.indicemin, myData.indicemax);
910     Parameters(aTestLine, myData.indicemin, myData.indicemax, myData.parametrization, aPars);
911
912     ApproxInt_KnotTools::BuildKnots(aPntXYZ, aPntU1V1, aPntU2V2, aPars,
913       myData.ApproxXYZ, myData.ApproxU1V1, myData.ApproxU2V2, aMinNbPnts, myKnots);
914   }
915   else
916   {
917     myKnots.Append(myData.indicemin);
918     myKnots.Append(myData.indicemax);
919   }
920 }
921
922 //=======================================================================
923 //function : buildCurve
924 //purpose  :
925 //=======================================================================
926 void ApproxInt_Approx::buildCurve(const Handle(TheWLine)& theline,
927                                   const Standard_Address thePtrSVSurf)
928 {
929   if(myData.myBezierApprox)
930   {
931     myBezToBSpl.Reset();
932   }
933
934   Standard_Integer kind = myKnots.Lower();
935   Standard_Integer imin, imax;
936   myTolReached = Standard_True;
937   Standard_Boolean OtherInter = Standard_False;
938   do
939   {
940     // Base cycle: iterate over knots.
941     imin = myKnots(kind);
942     imax = myKnots(kind+1);
943     ApproxInt_TheMultiLine myMultiLine(theline,
944       thePtrSVSurf,
945       ((myData.ApproxXYZ)? 1 : 0),
946       ((myData.ApproxU1V1)? 1: 0) + ((myData.ApproxU2V2)? 1: 0),
947       myData.Xo, myData.Ax, myData.Yo,myData.Ay, myData.Zo,myData.Az,
948       myData.U1o, myData.A1u, myData.V1o, myData.A1v,
949       myData.U2o, myData.A2u, myData.V2o, myData.A2v,
950       myData.ApproxU1V1,
951       imin,
952       imax);
953
954     if(myData.myBezierApprox)
955     {
956       myComputeLineBezier.Perform(myMultiLine);
957       if (myComputeLineBezier.NbMultiCurves() == 0)
958         return;
959       myTolReached&=myComputeLineBezier.IsToleranceReached();
960     }
961     else
962     {
963       myComputeLine.Perform(myMultiLine);
964     }
965     UpdateTolReached();
966
967     Standard_Integer indice3d,indice2d1,indice2d2;
968     indice3d = 1; 
969     indice2d1= 2;
970     indice2d2= 3;
971     if(!myData.ApproxXYZ)  { indice2d1--; indice2d2--; } 
972     if(!myData.ApproxU1V1) { indice2d2--; } 
973     if(myData.ApproxXYZ)
974     { 
975       Standard_Real ax,bx,ay,by,az,bz;
976       ax = 1.0/myData.Ax;   bx = -myData.Xo*ax;
977       ay = 1.0/myData.Ay;   by = -myData.Yo*ay;
978       az = 1.0/myData.Az;   bz = -myData.Zo*az;
979       if(myData.myBezierApprox)
980       {
981         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
982         {
983           myComputeLineBezier.ChangeValue(nbmc).Transform(indice3d,bx,ax,by,ay,bz,az);
984         }
985       }
986       else
987       {
988         myComputeLine.ChangeValue().Transform(indice3d,bx,ax,by,ay,bz,az);
989       }
990     }
991     if(myData.ApproxU1V1)
992     {
993       Standard_Real ax,bx,ay,by;
994       ax = 1.0/myData.A1u;   bx = -myData.U1o*ax;
995       ay = 1.0/myData.A1v;   by = -myData.V1o*ay;
996       if(myData.myBezierApprox) {
997         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
998         {
999           myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d1,bx,ax,by,ay);
1000         }
1001       }
1002       else
1003       {
1004         myComputeLine.ChangeValue().Transform2d(indice2d1,bx,ax,by,ay);
1005       }
1006     }
1007     if(myData.ApproxU2V2)
1008     {
1009       Standard_Real ax,bx,ay,by;
1010       ax = 1.0/myData.A2u;   bx = -myData.U2o*ax;
1011       ay = 1.0/myData.A2v;   by = -myData.V2o*ay;
1012       if(myData.myBezierApprox)
1013       {
1014         for(Standard_Integer nbmc = myComputeLineBezier.NbMultiCurves() ; nbmc>=1; nbmc--)
1015         {
1016           myComputeLineBezier.ChangeValue(nbmc).Transform2d(indice2d2,bx,ax,by,ay);
1017         }
1018       }
1019       else
1020       {
1021         myComputeLine.ChangeValue().Transform2d(indice2d2,bx,ax,by,ay);
1022       }
1023     }
1024
1025     OtherInter = Standard_False;
1026     if(myData.myBezierApprox)
1027     {
1028       for(Standard_Integer nbmc = 1; 
1029         nbmc <= myComputeLineBezier.NbMultiCurves();
1030         nbmc++)
1031       {
1032           myBezToBSpl.Append(myComputeLineBezier.Value(nbmc));
1033       }
1034       kind++;
1035       if(kind < myKnots.Upper())
1036       {
1037         OtherInter = Standard_True;
1038       }
1039     }
1040   }
1041   while(OtherInter);
1042
1043   if(myData.myBezierApprox)
1044   {
1045     myBezToBSpl.Perform();
1046   }
1047 }