0024428: Implementation of LGPL license
[occt.git] / src / AppParCurves / AppParCurves_Variational.gxx
CommitLineData
b311480e 1// Copyright (c) 1996-1999 Matra Datavision
973c2be1 2// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 3//
973c2be1 4// This file is part of Open CASCADE Technology software library.
b311480e 5//
973c2be1 6// This library is free software; you can redistribute it and / or modify it
7// under the terms of the GNU Lesser General Public version 2.1 as published
8// by the Free Software Foundation, with special exception defined in the file
9// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10// distribution for complete text of the license and disclaimer of any warranty.
b311480e 11//
973c2be1 12// Alternatively, this file may be used under the terms of Open CASCADE
13// commercial license or contractual agreement.
b311480e 14
7fd59977 15// Jeannine PANCIATICI le 06/06/96
16// Igor FEOKTISTOV 14/12/98 - correction of Approximate() and Init().
17// Approximation d une MultiLine de points decrite par le tool MLineTool.
18// avec criteres variationnels
19
20
21#define No_Standard_RangeError
22#define No_Standard_OutOfRange
23#define No_Standard_DimensionError
24#define No_Standard_ConstructionError
25
26#include <Standard_Stream.hxx>
27
28#include <AppParCurves.hxx>
29#include <AppParCurves_Constraint.hxx>
30#include <AppParCurves_HArray1OfConstraintCouple.hxx>
31#include <AppParCurves_Array1OfMultiPoint.hxx>
32#include <AppParCurves_MultiPoint.hxx>
33#include <AppParCurves_MultiBSpCurve.hxx>
34#include <Convert_CompPolynomialToPoles.hxx>
35#include <gp_Pnt.hxx>
36#include <gp_Pnt2d.hxx>
37#include <gp_Vec.hxx>
38#include <gp_Vec2d.hxx>
39#include <TColgp_Array1OfPnt.hxx>
40#include <TColgp_Array1OfPnt2d.hxx>
41#include <TColgp_Array1OfVec.hxx>
42#include <TColgp_Array1OfVec2d.hxx>
43#include <TColStd_Array1OfInteger.hxx>
44#include <TColStd_HArray1OfInteger.hxx>
45#include <TColStd_Array1OfReal.hxx>
46#include <TColStd_HArray1OfReal.hxx>
47#include <TColStd_HArray2OfReal.hxx>
48#include <StdFail_NotDone.hxx>
49#include <Standard_SStream.hxx>
50#include <Standard_NoSuchObject.hxx>
51#include <Precision.hxx>
52//#include <Smoothing.h>
53
54#if defined(WNT)
55# include <stdio.h>
56# include <stdarg.h>
57#endif /* WNT */
58
59
60//
61//=======================================================================
62//function : AppParCurves_Variational
63//purpose : Initialization of the fields.
64//=======================================================================
65//
66 AppParCurves_Variational::AppParCurves_Variational(const MultiLine& SSP,
67 const Standard_Integer FirstPoint,
68 const Standard_Integer LastPoint,
69 const Handle(AppParCurves_HArray1OfConstraintCouple)& TheConstraints,
70 const Standard_Integer MaxDegree,
71 const Standard_Integer MaxSegment,
72 const GeomAbs_Shape Continuity,
73 const Standard_Boolean WithMinMax,
74 const Standard_Boolean WithCutting,
75 const Standard_Real Tolerance,
76 const Standard_Integer NbIterations):
77 mySSP(SSP),
78 myFirstPoint(FirstPoint),
79 myLastPoint(LastPoint),
80 myConstraints(TheConstraints),
81 myMaxDegree(MaxDegree),
82 myMaxSegment(MaxSegment),
83 myNbIterations(NbIterations),
84 myTolerance(Tolerance),
85 myContinuity(Continuity),
86 myWithMinMax(WithMinMax),
87 myWithCutting(WithCutting)
88
89
90{
91// Verifications:
92 if (myMaxDegree < 1) Standard_DomainError::Raise();
93 myMaxDegree = Min (30, myMaxDegree);
94//
95 if (myMaxSegment < 1) Standard_DomainError::Raise();
96//
97 if (myWithMinMax != 0 && myWithMinMax !=1 ) Standard_DomainError::Raise();
98 if (myWithCutting != 0 && myWithCutting !=1 ) Standard_DomainError::Raise();
99//
100 myIsOverConstr = Standard_False;
101 myIsCreated = Standard_False;
102 myIsDone = Standard_False;
103 switch (myContinuity) {
104 case GeomAbs_C0:
105 myNivCont=0;
106 break ;
107 case GeomAbs_C1:
108 myNivCont=1;
109 break ;
110 case GeomAbs_C2:
111 myNivCont=2;
112 break ;
113 default:
114 Standard_ConstructionError::Raise();
115 }
116//
117 myNbP2d = ToolLine::NbP2d(SSP);
118 myNbP3d = ToolLine::NbP3d(SSP);
119 myDimension = 2 * myNbP2d + 3* myNbP3d ;
120//
121 myPercent[0]=0.4;
122 myPercent[1]=0.2;
123 myPercent[2]=0.4;
124 myKnots= new TColStd_HArray1OfReal(1,2);
125 myKnots->SetValue(1,0.);
126 myKnots->SetValue(2,1.);
127
128// Declaration
129//
130 mySmoothCriterion = new AppParCurves_MyCriterion(mySSP, myFirstPoint, myLastPoint);
131 myParameters = new TColStd_HArray1OfReal(myFirstPoint, myLastPoint);
132 myNbPoints=myLastPoint-myFirstPoint+1;
133 if (myNbPoints <= 0) Standard_ConstructionError::Raise();
134//
135 myTabPoints= new TColStd_HArray1OfReal(1,myDimension*myNbPoints);
136//
137// Table of Points initialization
138//
139 Standard_Integer ipoint,jp2d,jp3d,index;
140 TColgp_Array1OfPnt TabP3d(1, Max(1,myNbP3d));
141 TColgp_Array1OfPnt2d TabP2d(1, Max(1,myNbP2d));
142 gp_Pnt2d P2d;
143 gp_Pnt P3d;
144 index=1;
145
146 for ( ipoint = myFirstPoint ; ipoint <= myLastPoint ; ipoint++)
147 {
148
149 if(myNbP2d !=0 && myNbP3d ==0 ) {
150 ToolLine::Value(mySSP,ipoint,TabP2d);
151
152 for ( jp2d = 1 ; jp2d <= myNbP2d ;jp2d++)
153 { P2d = TabP2d.Value(jp2d);
154
155 myTabPoints->SetValue(index++,P2d.X());
156 myTabPoints->SetValue(index++,P2d.Y());
157 }
158 }
159 if(myNbP3d !=0 && myNbP2d == 0 ) {
160 ToolLine::Value(mySSP,ipoint,TabP3d);
161
162 for ( jp3d = 1 ; jp3d <= myNbP3d ; jp3d++)
163
164 { P3d=TabP3d.Value(jp3d);
165
166 myTabPoints->SetValue(index++,P3d.X());
167 myTabPoints->SetValue(index++,P3d.Y());
168 myTabPoints->SetValue(index++,P3d.Z());
169
170 }
171 }
172 if(myNbP3d !=0 && myNbP2d != 0 ) {
173 ToolLine::Value(mySSP,ipoint,TabP3d,TabP2d);
174
175 for ( jp3d = 1 ; jp3d <= myNbP3d ; jp3d++)
176
177 { P3d=TabP3d.Value(jp3d);
178
179 myTabPoints->SetValue(index++,P3d.X());
180 myTabPoints->SetValue(index++,P3d.Y());
181 myTabPoints->SetValue(index++,P3d.Z());
182
183 }
184 for ( jp2d = 1 ; jp2d <= myNbP2d ; jp2d++)
185
186 { P2d=TabP2d.Value(jp2d);
187
188 myTabPoints->SetValue(index++,P2d.X());
189 myTabPoints->SetValue(index++,P2d.Y());
190 }
191 }
192 }
193 Init();
194}
195//
196//=======================================================================
197//function : Init
198//purpose : Initializes the tables of constraints
199// and verifies if the problem is not over-constrained.
200// This method is used in the Create and the method SetConstraint.
201//=======================================================================
202//
203void AppParCurves_Variational::Init()
204{
205
206 Standard_Integer ipoint,jp2d,jp3d,index,jndex;
207 Standard_Integer CurMultyPoint;
208 TColgp_Array1OfVec TabV3d(1, Max(1,myNbP3d));
209 TColgp_Array1OfVec2d TabV2d(1, Max(1,myNbP2d));
210 TColgp_Array1OfVec TabV3dcurv(1, Max(1,myNbP3d));
211 TColgp_Array1OfVec2d TabV2dcurv(1, Max(1,myNbP2d));
212
213 gp_Vec Vt3d, Vc3d;
214 gp_Vec2d Vt2d, Vc2d;
215
216 myNbConstraints=myConstraints->Length();
217 if (myNbConstraints < 0) Standard_ConstructionError::Raise();
218
219 myTypConstraints = new TColStd_HArray1OfInteger(1,Max(1,2*myNbConstraints));
220 myTabConstraints = new TColStd_HArray1OfReal(1,Max(1,2*myDimension*myNbConstraints));
221 myTtheta = new TColStd_HArray1OfReal(1,Max(1,(2 * myNbP2d + 6 * myNbP3d) * myNbConstraints));
222 myTfthet = new TColStd_HArray1OfReal(1,Max(1,(2 * myNbP2d + 6 * myNbP3d) * myNbConstraints));
223
224
225//
226// Table of types initialization
227 Standard_Integer iconstr;
228 index=1;
229 jndex=1;
230 CurMultyPoint = 1;
231 myNbPassPoints=0;
232 myNbTangPoints=0;
233 myNbCurvPoints=0;
234 AppParCurves_Constraint valcontr;
235
236 for ( iconstr = myConstraints->Lower() ; iconstr <= myConstraints->Upper() ; iconstr++)
237 {
238 ipoint=(myConstraints->Value(iconstr)).Index();
239
240 valcontr=(myConstraints->Value(iconstr)).Constraint();
241 switch (valcontr) {
242 case AppParCurves_NoConstraint:
243 CurMultyPoint -= myNbP3d * 6 + myNbP2d * 2;
244 break ;
245 case AppParCurves_PassPoint:
246 myTypConstraints->SetValue(index++,ipoint);
247 myTypConstraints->SetValue(index++,0);
248 myNbPassPoints++;
249 if(myNbP2d !=0 ) jndex=jndex+4*myNbP2d;
250 if(myNbP3d !=0 ) jndex=jndex+6*myNbP3d;
251 break ;
252 case AppParCurves_TangencyPoint:
253 myTypConstraints->SetValue(index++,ipoint);
254 myTypConstraints->SetValue(index++,1);
255 myNbTangPoints++;
256 if(myNbP2d !=0 && myNbP3d == 0 )
257 {
258 if (ToolLine::Tangency(mySSP,ipoint,TabV2d) == Standard_False)
259 Standard_ConstructionError::Raise();
260 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
261 {
262 Vt2d=TabV2d.Value(jp2d);
263 Vt2d.Normalize();
264 myTabConstraints->SetValue(jndex++,Vt2d.X());
265 myTabConstraints->SetValue(jndex++,Vt2d.Y());
266 jndex=jndex+2;
267 InitTthetaF(2, valcontr, CurMultyPoint + (jp2d - 1) * 2, jndex - 4);
268 }
269 }
270 if(myNbP3d !=0 && myNbP2d == 0)
271 {
272 if (ToolLine::Tangency(mySSP,ipoint,TabV3d) == Standard_False)
273 Standard_ConstructionError::Raise();
274 for (jp3d=1;jp3d<=myNbP3d;jp3d++)
275 {
276 Vt3d=TabV3d.Value(jp3d);
277 Vt3d.Normalize();
278 myTabConstraints->SetValue(jndex++,Vt3d.X());
279
280 myTabConstraints->SetValue(jndex++,Vt3d.Y());
281
282 myTabConstraints->SetValue(jndex++,Vt3d.Z());
283 jndex=jndex+3;
284 InitTthetaF(3, valcontr, CurMultyPoint + (jp3d - 1) * 6, jndex - 6);
285 }
286 }
287 if(myNbP3d !=0 && myNbP2d != 0)
288 {
289 if (ToolLine::Tangency(mySSP,ipoint,TabV3d,TabV2d) == Standard_False)
290 Standard_ConstructionError::Raise();
291 for (jp3d=1;jp3d<=myNbP3d;jp3d++)
292 {
293 Vt3d=TabV3d.Value(jp3d);
294 Vt3d.Normalize();
295 myTabConstraints->SetValue(jndex++,Vt3d.X());
296 myTabConstraints->SetValue(jndex++,Vt3d.Y());
297 myTabConstraints->SetValue(jndex++,Vt3d.Z());
298 jndex=jndex+3;
299 InitTthetaF(3, valcontr, CurMultyPoint + (jp3d - 1) * 6, jndex - 6);
300 }
301
302 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
303 {
304 Vt2d=TabV2d.Value(jp2d);
305 Vt2d.Normalize();
306 myTabConstraints->SetValue(jndex++,Vt2d.X());
307 myTabConstraints->SetValue(jndex++,Vt2d.Y());
308 jndex=jndex+2;
309 InitTthetaF(2, valcontr, CurMultyPoint + (jp2d - 1) * 2 + myNbP3d * 6, jndex - 4);
310 }
311 }
312
313
314 break ;
315
316 case AppParCurves_CurvaturePoint:
317 myTypConstraints->SetValue(index++,ipoint);
318 myTypConstraints->SetValue(index++,2);
319 myNbCurvPoints++;
320 if(myNbP2d !=0 && myNbP3d == 0)
321 {
322 if (ToolLine::Tangency(mySSP,ipoint,TabV2d) == Standard_False )
323 Standard_ConstructionError::Raise();
324 if (ToolLine::Curvature(mySSP,ipoint,TabV2dcurv) == Standard_False)
325 Standard_ConstructionError::Raise();
326 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
327 {
328 Vt2d=TabV2d.Value(jp2d);
329 Vt2d.Normalize();
330 Vc2d=TabV2dcurv.Value(jp2d);
c6541a0c 331 if (Abs(Abs(Vc2d.Angle(Vt2d)) - M_PI/2.) > Precision::Angular())
7fd59977 332 Standard_ConstructionError::Raise();
333 myTabConstraints->SetValue(jndex++,Vt2d.X());
334 myTabConstraints->SetValue(jndex++,Vt2d.Y());
335 myTabConstraints->SetValue(jndex++,Vc2d.X());
336 myTabConstraints->SetValue(jndex++,Vc2d.Y());
337 InitTthetaF(2, valcontr, CurMultyPoint + (jp2d - 1) * 2, jndex - 4);
338 }
339 }
340
341 if(myNbP3d !=0 && myNbP2d == 0 )
342 {
343 if (ToolLine::Tangency(mySSP,ipoint,TabV3d) == Standard_False )
344 Standard_ConstructionError::Raise();
345 if (ToolLine::Curvature(mySSP,ipoint,TabV3dcurv) == Standard_False)
346 Standard_ConstructionError::Raise();
347 for (jp3d=1;jp3d<=myNbP3d;jp3d++)
348 {
349 Vt3d=TabV3d.Value(jp3d);
350 Vt3d.Normalize();
351 Vc3d=TabV3dcurv.Value(jp3d);
352 if ( (Vc3d.Normalized()).IsNormal(Vt3d,Precision::Angular()) == Standard_False)
353 Standard_ConstructionError::Raise();
354 myTabConstraints->SetValue(jndex++,Vt3d.X());
355 myTabConstraints->SetValue(jndex++,Vt3d.Y());
356 myTabConstraints->SetValue(jndex++,Vt3d.Z());
357 myTabConstraints->SetValue(jndex++,Vc3d.X());
358 myTabConstraints->SetValue(jndex++,Vc3d.Y());
359 myTabConstraints->SetValue(jndex++,Vc3d.Z());
360 InitTthetaF(3, valcontr, CurMultyPoint + (jp3d - 1) * 6, jndex - 6);
361 }
362 }
363 if(myNbP3d !=0 && myNbP2d != 0 )
364 {
365 if (ToolLine::Tangency(mySSP,ipoint,TabV3d,TabV2d) == Standard_False )
366 Standard_ConstructionError::Raise();
367 if (ToolLine::Curvature(mySSP,ipoint,TabV3dcurv,TabV2dcurv) == Standard_False)
368 Standard_ConstructionError::Raise();
369 for (jp3d=1;jp3d<=myNbP3d;jp3d++)
370 {
371 Vt3d=TabV3d.Value(jp3d);
372 Vt3d.Normalize();
373 Vc3d=TabV3dcurv.Value(jp3d);
374 if ( (Vc3d.Normalized()).IsNormal(Vt3d,Precision::Angular()) == Standard_False)
375 Standard_ConstructionError::Raise();
376 myTabConstraints->SetValue(jndex++,Vt3d.X());
377 myTabConstraints->SetValue(jndex++,Vt3d.Y());
378 myTabConstraints->SetValue(jndex++,Vt3d.Z());
379 myTabConstraints->SetValue(jndex++,Vc3d.X());
380 myTabConstraints->SetValue(jndex++,Vc3d.Y());
381 myTabConstraints->SetValue(jndex++,Vc3d.Z());
382 InitTthetaF(3, valcontr, CurMultyPoint + (jp3d - 1) * 6, jndex - 6);
383 }
384 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
385 {
386 Vt2d=TabV2d.Value(jp2d);
387 Vt2d.Normalize();
388 Vc2d=TabV2dcurv.Value(jp2d);
c6541a0c 389 if (Abs(Abs(Vc2d.Angle(Vt2d)) - M_PI/2.) > Precision::Angular())
7fd59977 390 Standard_ConstructionError::Raise();
391 myTabConstraints->SetValue(jndex++,Vt2d.X());
392 myTabConstraints->SetValue(jndex++,Vt2d.Y());
393 myTabConstraints->SetValue(jndex++,Vc2d.X());
394 myTabConstraints->SetValue(jndex++,Vc2d.Y());
395 InitTthetaF(2, valcontr, CurMultyPoint + (jp2d - 1) * 2 + myNbP3d * 6, jndex - 4);
396 }
397
398 }
399 break ;
400 default:
401 Standard_ConstructionError::Raise();
402 }
403 CurMultyPoint += myNbP3d * 6 + myNbP2d * 2;
404 }
405// OverConstraint Detection
406 Standard_Integer MaxSeg;
407 if(myWithCutting == Standard_True) MaxSeg = myMaxSegment ;
408 else MaxSeg = 1;
409 if (((myMaxDegree-myNivCont)*MaxSeg-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
410 {
411 myIsOverConstr =Standard_True;
412 myIsCreated = Standard_False;
413 }
414 else {
415 InitSmoothCriterion();
416 myIsCreated = Standard_True;
417 }
418
419
420}
421//
422//=======================================================================
423//function : Approximate
424//purpose : Makes the approximation with the current fields.
425//=======================================================================
426//
427void AppParCurves_Variational::Approximate()
428
429{
430 if (myIsCreated == Standard_False ) StdFail_NotDone:: Raise();
431
432
433 Standard_Real WQuadratic, WQuality;
434
435 TColStd_Array1OfReal Ecarts(myFirstPoint, myLastPoint);
436
437 mySmoothCriterion->GetWeight(WQuadratic, WQuality);
438
439 Handle(FEmTool_Curve) TheCurve;
440
441 mySmoothCriterion->GetCurve(TheCurve);
442
443//---------------------------------------------------------------------
444
445 TheMotor(mySmoothCriterion, WQuadratic, WQuality, TheCurve, Ecarts);
446
447
448 if(myWithMinMax && myTolerance < myMaxError)
449 Adjusting(mySmoothCriterion, WQuadratic, WQuality, TheCurve, Ecarts);
450
451//---------------------------------------------------------------------
452
453 Standard_Integer jp2d,jp3d,index,ipole,
454 NbElem = TheCurve->NbElements();
455
456 TColgp_Array1OfPnt TabP3d(1, Max(1,myNbP3d));
457 TColgp_Array1OfPnt2d TabP2d(1, Max(1,myNbP2d));
458 Standard_Real debfin[2] = {-1., 1};
459
460 gp_Pnt2d P2d;
461 gp_Pnt P3d;
462
463 index=0;
464
465 {
466 Handle(TColStd_HArray2OfReal) PolynomialIntervalsPtr =
467 new TColStd_HArray2OfReal(1,NbElem,1,2) ;
468
469 Handle(TColStd_HArray1OfInteger) NbCoeffPtr =
470 new TColStd_HArray1OfInteger(1, myMaxSegment);
471
472 Standard_Integer size = myMaxSegment * (myMaxDegree + 1) * myDimension ;
473 Handle(TColStd_HArray1OfReal) CoeffPtr = new TColStd_HArray1OfReal(1,size);
474
475 CoeffPtr->Init(0.);
476
477 Handle(TColStd_HArray1OfReal) IntervallesPtr =
478 new TColStd_HArray1OfReal(1, NbElem + 1);
479
480 IntervallesPtr->ChangeArray1() = TheCurve->Knots();
481
482 TheCurve->GetPolynom(CoeffPtr->ChangeArray1());
483
484 Standard_Integer ii;
485
486 for(ii = 1; ii <= NbElem; ii++)
487 NbCoeffPtr->SetValue(ii, TheCurve->Degree(ii)+1);
488
489
490 for (ii = PolynomialIntervalsPtr->LowerRow() ;
491 ii <= PolynomialIntervalsPtr->UpperRow() ;ii++)
492 {
493 PolynomialIntervalsPtr->SetValue(ii,1,debfin[0]) ;
494 PolynomialIntervalsPtr->SetValue(ii,2,debfin[1]) ;
495 }
496/*
497 printf("\n =========== Parameters for converting\n");
498 printf("nb_courbes : %d \n", NbElem);
499 printf("tab_option[4] : %d \n", myNivCont);
500 printf("myDimension : %d \n", myDimension);
501 printf("myMaxDegree : %d \n", myMaxDegree);
502 printf("\n NbcoeffPtr :\n");
503 for(ii = 1; ii <= NbElem; ii++) printf("NbCoeffPtr(%d) = %d \n", ii, NbCoeffPtr->Value(ii));
504 size = NbElem*(myMaxDegree + 1) * myDimension;
505 printf("\n CoeffPtr :\n");
506 for(ii = 1; ii <= size; ii++) printf("CoeffPtr(%d) = %.8e \n", ii, CoeffPtr->Value(ii));
507 printf("\n PolinimialIntervalsPtr :\n");
508 for (ii = PolynomialIntervalsPtr->LowerRow() ;
509 ii <= PolynomialIntervalsPtr->UpperRow() ;ii++)
510 {
511 printf(" %d : [%f, %f] \n", ii, PolynomialIntervalsPtr->Value(ii,1),
512 PolynomialIntervalsPtr->Value(ii,2)) ;
513 }
514 printf("\n IntervallesPtr :\n");
515 for (ii = IntervallesPtr->Lower();
516 ii <= IntervallesPtr->Upper() - 1; ii++)
517 {
518 printf(" %d : [%f, %f] \n", ii, IntervallesPtr->Value(ii),
519 IntervallesPtr->Value(ii+1)) ;
520 }
521*/
522 Convert_CompPolynomialToPoles AConverter(NbElem,
523 myNivCont,
524 myDimension,
525 myMaxDegree,
526 NbCoeffPtr,
527 CoeffPtr,
528 PolynomialIntervalsPtr,
529 IntervallesPtr) ;
530 if (AConverter.IsDone())
531 {
532 Handle(TColStd_HArray2OfReal) PolesPtr ;
533 Handle(TColStd_HArray1OfInteger) Mults;
534 Standard_Integer NbPoles=AConverter.NbPoles();
535// Standard_Integer Deg=AConverter.Degree();
536 AppParCurves_Array1OfMultiPoint TabMU(1,NbPoles);
537 AConverter.Poles(PolesPtr) ;
538 AConverter.Knots(myKnots) ;
539 AConverter.Multiplicities(Mults) ;
540
541 for (ipole=PolesPtr->LowerRow();ipole<=PolesPtr->UpperRow();ipole++)
542 {
543 index=PolesPtr->LowerCol();
544/* if(myNbP2d !=0 )
545 {
546 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
547 {
548 P2d.SetX(PolesPtr->Value(ipole,index++));
549 P2d.SetY(PolesPtr->Value(ipole,index++));
550 TabP2d.SetValue(jp2d,P2d);
551 }
552 }*/
553 if(myNbP3d !=0 )
554 {
555 for (jp3d=1;jp3d<=myNbP3d;jp3d++)
556 {
557 // cout << "\n Poles(ipole,1)" << PolesPtr->Value(ipole,index);
558 P3d.SetX(PolesPtr->Value(ipole,index++));
559 // cout << "\n Poles(ipole,1)" << PolesPtr->Value(ipole,index);
560 P3d.SetY(PolesPtr->Value(ipole,index++));
561 // cout << "\n Poles(ipole,1)" << PolesPtr->Value(ipole,index);
562 P3d.SetZ(PolesPtr->Value(ipole,index++));
563 TabP3d.SetValue(jp3d,P3d);
564 }
565 }
566 if(myNbP2d !=0 )
567 {
568 for (jp2d=1;jp2d<=myNbP2d;jp2d++)
569 {
570 P2d.SetX(PolesPtr->Value(ipole,index++));
571 P2d.SetY(PolesPtr->Value(ipole,index++));
572 TabP2d.SetValue(jp2d,P2d);
573 }
574 }
575 if(myNbP2d !=0 && myNbP3d !=0)
576 {
577 AppParCurves_MultiPoint aMultiPoint(TabP3d,TabP2d);
578 TabMU.SetValue(ipole,aMultiPoint);
579 }
580 else if (myNbP2d !=0)
581 {
582 AppParCurves_MultiPoint aMultiPoint(TabP2d);
583 TabMU.SetValue(ipole,aMultiPoint);
584 }
585 else
586 {
587
588 AppParCurves_MultiPoint aMultiPoint(TabP3d);
589 TabMU.SetValue(ipole,aMultiPoint);
590 }
591
592 }
593 AppParCurves_MultiBSpCurve aCurve(TabMU,myKnots->Array1(),Mults->Array1());
594 myMBSpCurve=aCurve;
595 myIsDone = Standard_True;
596 }
597 }
598
599}
600//
601//=======================================================================
602//function : IsCreated
603//purpose : returns True if the creation is done
604//=======================================================================
605//
606Standard_Boolean AppParCurves_Variational::IsCreated() const
607{
608 return myIsCreated;
609}
610//
611//=======================================================================
612//function : IsDone
613//purpose : returns True if the approximation is ok
614//=======================================================================
615//
616Standard_Boolean AppParCurves_Variational::IsDone() const
617{
618 return myIsDone;
619}
620//
621//=======================================================================
622//function : IsOverConstrained
623//purpose : returns True if the problem is overconstrained
624// in this case, approximation cannot be done.
625//=======================================================================
626//
627Standard_Boolean AppParCurves_Variational::IsOverConstrained() const
628{
629 return myIsOverConstr;
630}
631//
632//=======================================================================
633//function : Value
634//purpose : returns all the BSpline curves approximating the
635// MultiLine SSP after minimization of the parameter.
636
637//=======================================================================
638//
639AppParCurves_MultiBSpCurve AppParCurves_Variational::Value() const
640{
641 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
642 return myMBSpCurve;
643
644}
645//
646//=======================================================================
647//function : MaxError
648//purpose : returns the maximum of the distances between
649// the points of the multiline and the approximation
650// curves.
651//=======================================================================
652//
653Standard_Real AppParCurves_Variational::MaxError() const
654{
655 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
656 return myMaxError;
657}
658//
659//=======================================================================
660//function : MaxErrorIndex
661//purpose : returns the index of the MultiPoint of ErrorMax
662//=======================================================================
663//
664Standard_Integer AppParCurves_Variational::MaxErrorIndex() const
665{
666 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
667 return myMaxErrorIndex;
668}
669//
670//=======================================================================
671//function : QuadraticError
672//purpose : returns the quadratic average of the distances between
673// the points of the multiline and the approximation
674// curves.
675//=======================================================================
676//
677Standard_Real AppParCurves_Variational::QuadraticError() const
678{
679 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
680 return myCriterium[0];
681}
682//
683//=======================================================================
684//function : Distance
685//purpose : returns the distances between the points of the
686// multiline and the approximation curves.
687//=======================================================================
688//
689void AppParCurves_Variational::Distance(math_Matrix& mat)
690
691{
692 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
693 Standard_Integer ipoint,jp2d,jp3d,index;
694 TColgp_Array1OfPnt TabP3d(1,Max(1,myNbP3d));
695 TColgp_Array1OfPnt2d TabP2d(1, Max(1,myNbP2d));
696 Standard_Integer j0 = mat.LowerCol() - myFirstPoint;
697
698 gp_Pnt2d P2d;
699 gp_Pnt P3d;
700
701
702 gp_Pnt Pt3d;
703 gp_Pnt2d Pt2d;
704
705 for ( ipoint = myFirstPoint ; ipoint <= myLastPoint ; ipoint++)
706 {
707 index=1;
708 if(myNbP3d !=0 ) {
709 ToolLine::Value(mySSP,ipoint,TabP3d);
710
711 for ( jp3d = 1 ; jp3d <= myNbP3d ; jp3d++)
712
713 { P3d=TabP3d.Value(jp3d);
714 myMBSpCurve.Value(index,myParameters->Value(ipoint),Pt3d);
715 mat(index++, j0 + ipoint)=P3d.Distance(Pt3d);
716
717 }
718 }
719 if(myNbP2d !=0 ) {
720 if(myNbP3d == 0 ) ToolLine::Value(mySSP,ipoint,TabP2d);
721 else ToolLine::Value(mySSP,ipoint,TabP3d,TabP2d);
722 for ( jp2d = 1 ; jp2d <= myNbP2d ;jp2d++)
723
724 { P2d = TabP2d.Value(jp2d);
725 myMBSpCurve.Value(index,myParameters->Value(ipoint),Pt2d);
726 mat(index++, j0 + ipoint)=P2d.Distance(Pt2d);
727 }
728 }
729 }
730
731}
732//
733//=======================================================================
734//function : AverageError
735//purpose : returns the average error between
736// the MultiLine and the approximation.
737//=======================================================================
738//
739Standard_Real AppParCurves_Variational::AverageError() const
740{
741 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
742 return myAverageError;
743}
744//
745//=======================================================================
746//function : Parameters
747//purpose : returns the parameters uses to the approximations
748//=======================================================================
749//
750const Handle(TColStd_HArray1OfReal)& AppParCurves_Variational::Parameters() const
751{
752 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
753 return myParameters;
754}
755//
756//=======================================================================
757//function : Knots
758//purpose : returns the knots uses to the approximations
759//=======================================================================
760//
761const Handle(TColStd_HArray1OfReal)& AppParCurves_Variational::Knots() const
762{
763 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
764 return myKnots;
765}
766//
767//=======================================================================
768//function : Criterium
769//purpose : returns the values of the quality criterium.
770//=======================================================================
771//
772void AppParCurves_Variational::Criterium(Standard_Real& VFirstOrder, Standard_Real& VSecondOrder, Standard_Real& VThirdOrder) const
773{
774 if (myIsDone == Standard_False) StdFail_NotDone::Raise();
775 VFirstOrder=myCriterium[1] ;
776 VSecondOrder=myCriterium[2];
777 VThirdOrder=myCriterium[3];
778}
779//
780//=======================================================================
781//function : CriteriumWeight
782//purpose : returns the Weights (as percent) associed to the criterium used in
783// the optimization.
784//=======================================================================
785//
786void AppParCurves_Variational::CriteriumWeight(Standard_Real& Percent1, Standard_Real& Percent2, Standard_Real& Percent3) const
787{
788 Percent1 = myPercent[0];
789 Percent2 = myPercent[1];
790 Percent3 = myPercent[2] ;
791}
792//
793//=======================================================================
794//function : MaxDegree
795//purpose : returns the Maximum Degree used in the approximation
796//=======================================================================
797//
798Standard_Integer AppParCurves_Variational::MaxDegree() const
799{
800 return myMaxDegree;
801}
802//
803//=======================================================================
804//function : MaxSegment
805//purpose : returns the Maximum of segment used in the approximation
806//=======================================================================
807//
808Standard_Integer AppParCurves_Variational::MaxSegment() const
809{
810 return myMaxSegment;
811 }
812//
813//=======================================================================
814//function : Continuity
815//purpose : returns the Continuity used in the approximation
816//=======================================================================
817//
818GeomAbs_Shape AppParCurves_Variational::Continuity() const
819{
820 return myContinuity;
821}
822//
823//=======================================================================
824//function : WithMinMax
825//purpose : returns if the approximation search to minimize the
826// maximum Error or not.
827//=======================================================================
828//
829Standard_Boolean AppParCurves_Variational::WithMinMax() const
830{
831 return myWithMinMax;
832}
833//
834//=======================================================================
835//function : WithCutting
836//purpose : returns if the approximation can insert new Knots or not.
837//=======================================================================
838//
839Standard_Boolean AppParCurves_Variational::WithCutting() const
840{
841 return myWithCutting;
842}
843//
844//=======================================================================
845//function : Tolerance
846//purpose : returns the tolerance used in the approximation.
847//=======================================================================
848//
849Standard_Real AppParCurves_Variational::Tolerance() const
850{
851 return myTolerance;
852}
853//
854//=======================================================================
855//function : NbIterations
856//purpose : returns the number of iterations used in the approximation.
857//=======================================================================
858//
859Standard_Integer AppParCurves_Variational::NbIterations() const
860{
861 return myNbIterations;
862}
863//
864//=======================================================================
865//function : Dump
866//purpose : Prints on the stream o information on the current state
867// of the object.
868//=======================================================================
869//
870void AppParCurves_Variational::Dump(Standard_OStream& o) const
871{
872 o << " \nVariational Smoothing " << endl;
873 o << " Number of multipoints " << myNbPoints << endl;
874 o << " Number of 2d par multipoint " << myNbP2d << endl;
875 o << " Nombre of 3d par multipoint " << myNbP3d << endl;
876 o << " Number of PassagePoint " << myNbPassPoints << endl;
877 o << " Number of TangencyPoints " << myNbTangPoints << endl;
878 o << " Number of CurvaturePoints " << myNbCurvPoints << endl;
879 o << " \nTolerance " << o.setf(ios::scientific) << setprecision(3) << setw(9) << myTolerance;
880 if ( WithMinMax()) { o << " as Max Error." << endl;}
881 else { o << " as size Error." << endl;}
882 o << "CriteriumWeights : " << myPercent[0] << " , "
883 << myPercent[1] << " , " << myPercent[2] << endl;
884
885 if (myIsDone ) {
886 o << " MaxError " << setprecision(3) << setw(9) << myMaxError << endl;
887 o << " Index of MaxError " << myMaxErrorIndex << endl;
888 o << " Average Error " << setprecision(3) << setw(9) << myAverageError << endl;
889 o << " Quadratic Error " << setprecision(3) << setw(9) << myCriterium[0] << endl;
890 o << " Tension Criterium " << setprecision(3) << setw(9) << myCriterium[1] << endl;
891 o << " Flexion Criterium " << setprecision(3) << setw(9) << myCriterium[2] << endl;
892 o << " Jerk Criterium " << setprecision(3) << setw(9) << myCriterium[3] << endl;
893 o << " NbSegments " << myKnots->Length()-1 << endl;
894 }
895 else
896 { if (myIsOverConstr) o << "The probleme is overconstraint " << endl;
897 else o << " Erreur dans l''approximation" << endl;
898 }
899}
900//
901//=======================================================================
902//function : SetConstraints
903//purpose : Define the constraints to approximate
904// If this value is incompatible with the others fields
905// this method modify nothing and returns false
906//=======================================================================
907//
908Standard_Boolean AppParCurves_Variational::SetConstraints(const Handle(AppParCurves_HArray1OfConstraintCouple)& aConstraint)
909
910{
911 myConstraints=aConstraint;
912 Init();
913 if (myIsOverConstr ) return Standard_False;
914 else return Standard_True;
915}
916//
917//=======================================================================
918//function : SetParameters
919//purpose : Defines the parameters used by the approximations.
920//=======================================================================
921//
922void AppParCurves_Variational::SetParameters(const Handle(TColStd_HArray1OfReal)& param)
923{
924 myParameters->ChangeArray1() = param->Array1();
925}
926//
927//=======================================================================
928//function : SetKnots
929//purpose : Defines the knots used by the approximations
930// -- If this value is incompatible with the others fields
931// -- this method modify nothing and returns false
932//=======================================================================
933//
934Standard_Boolean AppParCurves_Variational::SetKnots(const Handle(TColStd_HArray1OfReal)& knots)
935{
936 myKnots->ChangeArray1() = knots->Array1();
937 return Standard_True;
938}
939//
940//=======================================================================
941//function : SetMaxDegree
942//purpose : Define the Maximum Degree used in the approximation
943// If this value is incompatible with the others fields
944// this method modify nothing and returns false
945//=======================================================================
946//
947Standard_Boolean AppParCurves_Variational::SetMaxDegree(const Standard_Integer Degree)
948{
949 if (((Degree-myNivCont)*myMaxSegment-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
950 return Standard_False;
951 else
952 {
953 myMaxDegree=Degree;
954
955 InitSmoothCriterion();
956
957 return Standard_True;
958 }
959
960}
961//
962//=======================================================================
963//function : SetMaxSegment
964//purpose : Define the maximum number of segments used in the approximation
965// If this value is incompatible with the others fields
966// this method modify nothing and returns false
967//=======================================================================
968//
969Standard_Boolean AppParCurves_Variational::SetMaxSegment(const Standard_Integer NbSegment)
970{
971 if ( myWithCutting == Standard_True &&
972 ((myMaxDegree-myNivCont)*NbSegment-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
973 return Standard_False;
974 else
975 {
976 myMaxSegment=NbSegment;
977 return Standard_True;
978 }
979}
980//
981//=======================================================================
982//function : SetContinuity
983//purpose : Define the Continuity used in the approximation
984// If this value is incompatible with the others fields
985// this method modify nothing and returns false
986//=======================================================================
987//
988Standard_Boolean AppParCurves_Variational::SetContinuity(const GeomAbs_Shape C)
989{
990 Standard_Integer NivCont=0;
991 switch (C) {
992 case GeomAbs_C0:
993 NivCont=0;
994 break ;
995 case GeomAbs_C1:
996 NivCont=1;
997 break ;
998 case GeomAbs_C2:
999 NivCont=2;
1000 break ;
1001 default:
1002 Standard_ConstructionError::Raise();
1003 }
1004 if (((myMaxDegree-NivCont)*myMaxSegment-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
1005 return Standard_False;
1006 else
1007 {
1008 myContinuity= C;
1009 myNivCont=NivCont;
1010
1011 InitSmoothCriterion();
1012 return Standard_True;
1013 }
1014}
1015//
1016//=======================================================================
1017//function : SetWithMinMax
1018//purpose : Define if the approximation search to minimize the
1019// maximum Error or not.
1020//=======================================================================
1021//
1022void AppParCurves_Variational::SetWithMinMax(const Standard_Boolean MinMax)
1023{
1024 myWithMinMax=MinMax;
1025
1026 InitSmoothCriterion();
1027}
1028//
1029//=======================================================================
1030//function : SetWithCutting
1031//purpose : Define if the approximation can insert new Knots or not.
1032// If this value is incompatible with the others fields
1033// this method modify nothing and returns false
1034//=======================================================================
1035//
1036Standard_Boolean AppParCurves_Variational::SetWithCutting(const Standard_Boolean Cutting)
1037{
1038 if (Cutting == Standard_False)
1039 {
1040 if (((myMaxDegree-myNivCont)*myKnots->Length()-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
1041 return Standard_False;
1042
1043 else
1044 {
1045 myWithCutting=Cutting;
1046 InitSmoothCriterion();
1047 return Standard_True;
1048 }
1049 }
1050 else
1051 {
1052 if (((myMaxDegree-myNivCont)*myMaxSegment-myNbPassPoints-2*myNbTangPoints-3*myNbCurvPoints) < 0 )
1053 return Standard_False;
1054
1055 else
1056 {
1057 myWithCutting=Cutting;
1058 InitSmoothCriterion();
1059 return Standard_True;
1060 }
1061 }
1062}
1063//
1064//=======================================================================
1065//function : SetCriteriumWeight
1066//purpose : define the Weights (as percent) associed to the criterium used in
1067// the optimization.
1068//=======================================================================
1069//
1070void AppParCurves_Variational::SetCriteriumWeight(const Standard_Real Percent1, const Standard_Real Percent2, const Standard_Real Percent3)
1071{
1072 if (Percent1 < 0 || Percent2 < 0 || Percent3 < 0 ) Standard_DomainError::Raise();
1073 Standard_Real Total = Percent1 + Percent2 + Percent3;
1074 myPercent[0] = Percent1/Total;
1075 myPercent[1] = Percent2/Total;
1076 myPercent[2] = Percent3/Total;
1077
1078 InitSmoothCriterion();
1079}
1080//
1081//=======================================================================
1082//function : SetCriteriumWeight
1083//purpose : define the Weight (as percent) associed to the
1084// criterium Order used in the optimization : Others
1085// weights are updated.
1086//=======================================================================
1087//
1088void AppParCurves_Variational::SetCriteriumWeight(const Standard_Integer Order, const Standard_Real Percent)
1089{
1090 if ( Percent < 0 ) Standard_DomainError::Raise();
1091 if ( Order < 1 || Order > 3 ) Standard_ConstructionError::Raise();
1092 myPercent[Order-1] = Percent;
1093 Standard_Real Total = myPercent[0] + myPercent[1] + myPercent[2];
1094 myPercent[0] = myPercent[0] / Total;
1095 myPercent[1] = myPercent[1] / Total;
1096 myPercent[2] = myPercent[2] / Total;
1097
1098 InitSmoothCriterion();
1099
1100}
1101//
1102//=======================================================================
1103//function : SetTolerance
1104//purpose : define the tolerance used in the approximation.
1105//=======================================================================
1106//
1107void AppParCurves_Variational::SetTolerance(const Standard_Real Tol)
1108{
1109 myTolerance=Tol;
1110 InitSmoothCriterion();
1111}
1112//
1113//=======================================================================
1114//function : SetNbIterations
1115//purpose : define the number of iterations used in the approximation.
1116//=======================================================================
1117//
1118void AppParCurves_Variational::SetNbIterations(const Standard_Integer Iter)
1119{
1120 myNbIterations=Iter;
1121}
1122
1123
1124// Private methods
1125#include <AppParCurves_Variational_1.gxx> // TheMotor
1126#include <AppParCurves_Variational_2.gxx> // Optimization
1127#include <AppParCurves_Variational_3.gxx> // Project
1128#include <AppParCurves_Variational_4.gxx> // ACR
1129#include <AppParCurves_Variational_5.gxx> // SplitCurve
1130#include <AppParCurves_Variational_6.gxx> // InitSmoothCriterion and other Init... methods
1131#include <AppParCurves_Variational_7.gxx> // Adjusting
1132#include <AppParCurves_Variational_8.gxx> // AssemblingConstraints
1133#include <AppParCurves_Variational_9.gxx> // InitTtheta and InitTfthet methods