63f9aa3b571e6a7d6b1be6387eedfab26be752ed
[occt.git] / src / AdvApp2Var / AdvApp2Var_ApproxAFunc2Var.cxx
1 // Created on: 1996-07-03
2 // Created by: Joelle CHAUVET
3 // Copyright (c) 1996-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 #include <AdvApp2Var_ApproxAFunc2Var.hxx>
23 #include <AdvApp2Var_EvaluatorFunc2Var.hxx>
24 #include <AdvApp2Var_Criterion.hxx>
25 #include <AdvApp2Var_Context.hxx>
26 #include <AdvApp2Var_Patch.hxx>
27 #include <AdvApp2Var_Network.hxx>
28 #include <AdvApp2Var_Node.hxx>
29 #include <AdvApp2Var_Iso.hxx>
30 #include <AdvApp2Var_Strip.hxx>
31 #include <AdvApp2Var_Framework.hxx>
32 #include <AdvApprox_Cutting.hxx>
33
34 #include <Standard_ConstructionError.hxx>
35 #include <Standard_OutOfRange.hxx>
36 #include <TColStd_HArray1OfInteger.hxx>
37 #include <TColStd_HArray2OfInteger.hxx>
38 #include <TColStd_HArray1OfReal.hxx>
39 #include <TColStd_HArray2OfReal.hxx>
40
41 #include <gp_XY.hxx>
42 #include <gp_Pnt2d.hxx>
43 #include <gp_Pnt.hxx>
44 #include <TColgp_HArray2OfPnt.hxx>
45
46 #include <Convert_GridPolynomialToPoles.hxx>
47
48 #include <Geom_BezierSurface.hxx>
49
50
51 //=======================================================================
52 //function : AdvApp2Var_ApproxAFunc2Var
53 //purpose  : 
54 //=======================================================================
55
56 AdvApp2Var_ApproxAFunc2Var::
57 AdvApp2Var_ApproxAFunc2Var(const Standard_Integer Num1DSS,
58                            const Standard_Integer Num2DSS,
59                            const Standard_Integer Num3DSS,
60                            const Handle(TColStd_HArray1OfReal)& OneDTol,
61                            const Handle(TColStd_HArray1OfReal)& TwoDTol,
62                            const Handle(TColStd_HArray1OfReal)& ThreeDTol,
63                            const Handle(TColStd_HArray2OfReal)& OneDTolFr,
64                            const Handle(TColStd_HArray2OfReal)& TwoDTolFr,
65                            const Handle(TColStd_HArray2OfReal)& ThreeDTolFr,
66                            const Standard_Real FirstInU,
67                            const Standard_Real LastInU,
68                            const Standard_Real FirstInV,
69                            const Standard_Real LastInV,
70                            const GeomAbs_IsoType FavorIso,
71                            const GeomAbs_Shape ContInU,
72                            const GeomAbs_Shape ContInV,
73                            const Standard_Integer PrecisCode,
74                            const Standard_Integer MaxDegInU,
75                            const Standard_Integer MaxDegInV,
76                            const Standard_Integer MaxPatch,
77                            const AdvApp2Var_EvaluatorFunc2Var& Func,
78                            AdvApprox_Cutting& UChoice,
79                            AdvApprox_Cutting& VChoice) :
80 my1DTolerances(OneDTol),
81 my2DTolerances(TwoDTol),
82 my3DTolerances(ThreeDTol),
83 my1DTolOnFront(OneDTolFr),
84 my2DTolOnFront(TwoDTolFr),
85 my3DTolOnFront(ThreeDTolFr),
86 myFirstParInU(FirstInU),
87 myLastParInU(LastInU),
88 myFirstParInV(FirstInV),
89 myLastParInV(LastInV),
90 myFavoriteIso(FavorIso),
91 myContInU(ContInU),
92 myContInV(ContInV),
93 myPrecisionCode(PrecisCode),
94 myMaxDegInU(MaxDegInU),
95 myMaxDegInV(MaxDegInV),
96 myMaxPatches(MaxPatch),
97 myDone(Standard_False),
98 myHasResult(Standard_False)
99 {
100   myNumSubSpaces[0] = Num1DSS;
101   myNumSubSpaces[1] = Num2DSS;
102   myNumSubSpaces[2] = Num3DSS;
103   Init();
104   Perform(UChoice, VChoice, Func);
105   ConvertBS();
106 }
107
108 //=======================================================================
109 //function : AdvApp2Var_ApproxAFunc2Var
110 //purpose  : 
111 //=======================================================================
112
113  AdvApp2Var_ApproxAFunc2Var::
114 AdvApp2Var_ApproxAFunc2Var(const Standard_Integer Num1DSS,
115                            const Standard_Integer Num2DSS,
116                            const Standard_Integer Num3DSS,
117                            const Handle(TColStd_HArray1OfReal)& OneDTol,
118                            const Handle(TColStd_HArray1OfReal)& TwoDTol,
119                            const Handle(TColStd_HArray1OfReal)& ThreeDTol,
120                            const Handle(TColStd_HArray2OfReal)& OneDTolFr,
121                            const Handle(TColStd_HArray2OfReal)& TwoDTolFr,
122                            const Handle(TColStd_HArray2OfReal)& ThreeDTolFr,
123                            const Standard_Real FirstInU,
124                            const Standard_Real LastInU,
125                            const Standard_Real FirstInV,
126                            const Standard_Real LastInV,
127                            const GeomAbs_IsoType FavorIso,
128                            const GeomAbs_Shape ContInU,
129                            const GeomAbs_Shape ContInV,
130                            const Standard_Integer PrecisCode,
131                            const Standard_Integer MaxDegInU,
132                            const Standard_Integer MaxDegInV,
133                            const Standard_Integer MaxPatch,
134                            const AdvApp2Var_EvaluatorFunc2Var& Func,
135                            const AdvApp2Var_Criterion& Crit,
136                            AdvApprox_Cutting& UChoice,
137                            AdvApprox_Cutting& VChoice) :
138 my1DTolerances(OneDTol),
139 my2DTolerances(TwoDTol),
140 my3DTolerances(ThreeDTol),
141 my1DTolOnFront(OneDTolFr),
142 my2DTolOnFront(TwoDTolFr),
143 my3DTolOnFront(ThreeDTolFr),
144 myFirstParInU(FirstInU),
145 myLastParInU(LastInU),
146 myFirstParInV(FirstInV),
147 myLastParInV(LastInV),
148 myFavoriteIso(FavorIso),
149 myContInU(ContInU),
150 myContInV(ContInV),
151 myPrecisionCode(PrecisCode),
152 myMaxDegInU(MaxDegInU),
153 myMaxDegInV(MaxDegInV),
154 myMaxPatches(MaxPatch),
155 myDone(Standard_False),
156 myHasResult(Standard_False)
157 {
158   myNumSubSpaces[0] = Num1DSS;
159   myNumSubSpaces[1] = Num2DSS;
160   myNumSubSpaces[2] = Num3DSS;
161   Init();
162   Perform(UChoice, VChoice, Func, Crit);
163   ConvertBS();
164 }
165
166 //=======================================================================
167 //function : Init
168 //purpose  : Initialisation of the approximation
169 //=======================================================================
170
171 void AdvApp2Var_ApproxAFunc2Var::Init()
172 {
173   Standard_Integer ifav,iu=0,iv=0,ndu,ndv;
174   switch (myFavoriteIso) {
175   case GeomAbs_IsoU :
176     ifav = 1;
177     break;
178   case GeomAbs_IsoV :
179     ifav = 2;
180     break;
181   default :
182     ifav = 2;
183     break;
184   }
185   switch (myContInU) {
186   case GeomAbs_C0 :
187     iu = 0;
188     break;
189   case GeomAbs_C1 :
190     iu = 1;
191     break;
192   case GeomAbs_C2 :
193     iu = 2;
194     break;
195   default :
196     Standard_ConstructionError::Raise("AdvApp2Var_ApproxAFunc2Var : UContinuity Error");
197   }
198   switch (myContInV) {
199   case GeomAbs_C0 :
200     iv = 0;
201     break;
202   case GeomAbs_C1 :
203     iv = 1;
204     break;
205   case GeomAbs_C2 :
206     iv = 2;
207     break;
208   default :
209     Standard_ConstructionError::Raise("AdvApp2Var_ApproxAFunc2Var : VContinuity Error");
210   }
211   ndu = Max(myMaxDegInU+1,2*iu+2);
212   ndv = Max(myMaxDegInV+1,2*iv+2);
213   if (ndu<2*iu+2)
214     Standard_ConstructionError::Raise("AdvApp2Var_ApproxAFunc2Var : UMaxDegree Error");
215   if (ndv<2*iv+2)
216     Standard_ConstructionError::Raise("AdvApp2Var_ApproxAFunc2Var : VMaxDegree Error");
217   myPrecisionCode = Max(0,Min(myPrecisionCode,3));
218   AdvApp2Var_Context Conditions(ifav,iu,iv,ndu,ndv,
219                                 myPrecisionCode,
220                                 myNumSubSpaces[0],
221                                 myNumSubSpaces[1],
222                                 myNumSubSpaces[2],
223                                 my1DTolerances,
224                                 my2DTolerances,
225                                 my3DTolerances,
226                                 my1DTolOnFront,
227                                 my2DTolOnFront,
228                                 my3DTolOnFront);
229   myConditions = Conditions;
230   InitGrid(1);
231 }
232
233
234 //=======================================================================
235 //function : InitGrid
236 //purpose  : Initialisation of the approximation with regular cuttings
237 //=======================================================================
238
239 void AdvApp2Var_ApproxAFunc2Var::InitGrid(const Standard_Integer NbInt)
240 {
241   Standard_Integer iu=myConditions.UOrder(),iv=myConditions.VOrder(),iint;
242
243   AdvApp2Var_Patch M0(myFirstParInU,myLastParInU,myFirstParInV,myLastParInV,iu,iv);
244
245   AdvApp2Var_SequenceOfPatch Net;
246   Net.Append(M0);
247
248   TColStd_SequenceOfReal TheU,TheV;
249   TheU.Append(myFirstParInU);
250   TheV.Append(myFirstParInV);
251   TheU.Append(myLastParInU);
252   TheV.Append(myLastParInV);
253
254   AdvApp2Var_Network Result(Net,TheU,TheV);
255
256
257   gp_XY UV1(myFirstParInU,myFirstParInV);
258   AdvApp2Var_Node C1(UV1,iu,iv);
259   gp_XY UV2(myLastParInU,myFirstParInV);
260   AdvApp2Var_Node C2(UV2,iu,iv);
261   gp_XY UV4(myLastParInU,myLastParInV);
262   AdvApp2Var_Node C4(UV4,iu,iv);
263   gp_XY UV3(myFirstParInU,myLastParInV);
264   AdvApp2Var_Node C3(UV3,iu,iv);
265   AdvApp2Var_SequenceOfNode Bag;
266   Bag.Append(C1);
267   Bag.Append(C2);
268   Bag.Append(C3);
269   Bag.Append(C4);
270
271   AdvApp2Var_Iso V0(GeomAbs_IsoV,myFirstParInV,
272                     myFirstParInU,myLastParInU,myFirstParInV,myLastParInV,
273                     1,iu,iv);
274   AdvApp2Var_Iso V1(GeomAbs_IsoV,myLastParInV,
275                     myFirstParInU,myLastParInU,myFirstParInV,myLastParInV,
276                     2,iu,iv);
277   AdvApp2Var_Iso U0(GeomAbs_IsoU,myFirstParInU,
278                     myFirstParInU,myLastParInU,myFirstParInV,myLastParInV,
279                     3,iu,iv);
280   AdvApp2Var_Iso U1(GeomAbs_IsoU,myLastParInU,
281                     myFirstParInU,myLastParInU,myFirstParInV,myLastParInV,
282                     4,iu,iv);
283
284   AdvApp2Var_Strip BU0,BV0;
285   BU0.Append(V0);
286   BU0.Append(V1);
287   BV0.Append(U0);
288   BV0.Append(U1);
289
290   AdvApp2Var_SequenceOfStrip UStrip,VStrip;
291   UStrip.Append(BU0);
292   VStrip.Append(BV0);
293
294   AdvApp2Var_Framework Constraints(Bag,UStrip,VStrip);
295
296 // regular cutting if NbInt>1
297   Standard_Real deltu = (myLastParInU-myFirstParInU)/NbInt,
298                 deltv = (myLastParInV-myFirstParInV)/NbInt;
299   for (iint=1;iint<=NbInt-1;iint++) {
300     Result.UpdateInU(myFirstParInU+iint*deltu);
301     Constraints.UpdateInU(myFirstParInU+iint*deltu);
302     Result.UpdateInV(myFirstParInV+iint*deltv);
303     Constraints.UpdateInV(myFirstParInV+iint*deltv);
304   }
305   myResult = Result;
306   myConstraints = Constraints;
307 }
308
309 //=======================================================================
310 //function : Perform
311 //purpose  : Computation of the approximation
312 //=======================================================================
313
314 void AdvApp2Var_ApproxAFunc2Var::Perform(const AdvApprox_Cutting& UChoice,
315                                          const AdvApprox_Cutting& VChoice,
316                                          const AdvApp2Var_EvaluatorFunc2Var& Func)
317 {
318   ComputePatches(UChoice,VChoice,Func);
319   myHasResult = myDone = Standard_True;
320   Compute3DErrors();
321 }
322
323 //=======================================================================
324 //function : Perform
325 //purpose  : Computation of the approximation
326 //=======================================================================
327
328 void AdvApp2Var_ApproxAFunc2Var::Perform(const AdvApprox_Cutting& UChoice,
329                                          const AdvApprox_Cutting& VChoice,
330                                          const AdvApp2Var_EvaluatorFunc2Var& Func,
331                                          const AdvApp2Var_Criterion& Crit)
332 {
333   ComputePatches(UChoice,VChoice,Func,Crit);
334   myHasResult = myDone = Standard_True;
335   Compute3DErrors();
336   ComputeCritError();
337 }
338
339 //=======================================================================
340 //function : ComputePatches
341 //purpose  : Computation of the polynomial approximations
342 //=======================================================================
343
344 void AdvApp2Var_ApproxAFunc2Var::ComputePatches(const AdvApprox_Cutting& UChoice,
345                                                 const AdvApprox_Cutting& VChoice,
346                                          const AdvApp2Var_EvaluatorFunc2Var& Func)
347 {
348   Standard_Real Udec, Vdec;
349   Standard_Boolean Umore, Vmore;
350   Standard_Integer NbPatch, NbU, NbV, NumDec;
351   Standard_Integer FirstNA;
352
353   while (myResult.FirstNotApprox(FirstNA)) {
354
355 // complete the set of constraints 
356     ComputeConstraints(UChoice, VChoice, Func);
357
358 // discretization of constraints relative to the square
359     myResult(FirstNA).Discretise(myConditions,myConstraints,Func);
360     if ( ! myResult(FirstNA).IsDiscretised() ) {
361       myHasResult =  myDone = Standard_False;
362       Standard_ConstructionError::Raise
363               ("AdvApp2Var_ApproxAFunc2Var : Surface Discretisation Error");
364     }
365
366 // calculate the number and the type of autorized cuts
367 // depending on the max number of squares and the validity of next cuts.
368     NbU = myResult.NbPatchInU();
369     NbV = myResult.NbPatchInV();
370     NbPatch = NbU*NbV;
371     Umore = UChoice.Value(myResult(FirstNA).U0(), myResult(FirstNA).U1(),Udec);
372     Vmore = VChoice.Value(myResult(FirstNA).V0(), myResult(FirstNA).V1(),Vdec);
373
374     NumDec = 0;
375     if ( ((NbPatch+NbV)<=myMaxPatches) && ((NbPatch+NbU)>myMaxPatches)
376          && (Umore) ) NumDec = 1;
377     if ( ((NbPatch+NbV)>myMaxPatches) && ((NbPatch+NbU)<=myMaxPatches) 
378          && (Vmore) ) NumDec = 2;
379     if ( ((NbPatch+NbV)<=myMaxPatches) && ((NbPatch+NbU)<=myMaxPatches) ) {
380       if ( Umore ) NumDec = 3;
381       if ( (NbV>NbU) && Vmore ) NumDec = 4;
382     }
383     if ( (NbU+1)*(NbV+1)<=myMaxPatches ) {
384       if ( !Umore && !Vmore ) NumDec=0;
385       if ( Umore && !Vmore ) NumDec=3;
386       if ( !Umore && Vmore ) NumDec=4;
387       if ( Umore && Vmore ) NumDec=5;
388     }
389
390 // approximation of the square
391     myResult(FirstNA).MakeApprox(myConditions,myConstraints,NumDec);
392
393     if ( ! myResult(FirstNA).IsApproximated() ) {
394       switch (myResult(FirstNA).CutSense()) {
395       case 0 :
396 //      It is not possible to cut : the result is preserved
397         if ( myResult(FirstNA).HasResult()) {
398           myResult(FirstNA).OverwriteApprox();
399         }
400         else {
401           myHasResult =  myDone = Standard_False;
402           Standard_ConstructionError::Raise
403               ("AdvApp2Var_ApproxAFunc2Var : Surface Approximation Error");
404         }
405         break;
406       case 1 :
407 //      It is necessary to cut in U
408         myResult.UpdateInU(Udec);
409         myConstraints.UpdateInU(Udec);
410         break;
411       case 2 :
412 //      It is necessary to cut in V
413         myResult.UpdateInV(Vdec);
414         myConstraints.UpdateInV(Vdec);
415         break;
416       case 3 :
417 //      It is necesary to cut in U and V
418         myResult.UpdateInU(Udec);
419         myConstraints.UpdateInU(Udec);
420         myResult.UpdateInV(Vdec);
421         myConstraints.UpdateInV(Vdec);
422         break;
423       default :
424         myHasResult =  myDone = Standard_False;
425         Standard_ConstructionError::Raise
426               ("AdvApp2Var_ApproxAFunc2Var : Surface Approximation Error");
427       }
428     }
429   }
430 }
431
432 //=======================================================================
433 //function : ComputePatches
434 //purpose  : Computation of the polynomial approximations
435 //=======================================================================
436
437 void AdvApp2Var_ApproxAFunc2Var::ComputePatches(const AdvApprox_Cutting& UChoice,
438                                          const AdvApprox_Cutting& VChoice,
439                                          const AdvApp2Var_EvaluatorFunc2Var& Func,
440                                          const AdvApp2Var_Criterion& Crit)
441 {
442   Standard_Real Udec, Vdec, CritValue, m0=0., m1=0.;
443   Standard_Boolean Umore, Vmore, CritAbs = (Crit.Type() == AdvApp2Var_Absolute);
444   Standard_Integer NbPatch, NbU, NbV, NbInt, NumDec;
445   Standard_Integer FirstNA, decision=0;
446
447   while (myResult.FirstNotApprox(FirstNA)) {
448
449 // complete the set of constraints 
450     ComputeConstraints(UChoice, VChoice, Func, Crit);
451     if (decision>0) {
452       m0 = m1;
453       m1 = 0.;
454     }
455
456 // discretize the constraints relative to the square
457     myResult(FirstNA).Discretise(myConditions,myConstraints,Func);
458     if ( ! myResult(FirstNA).IsDiscretised() ) {
459       myHasResult =  myDone = Standard_False;
460       Standard_ConstructionError::Raise
461               ("AdvApp2Var_ApproxAFunc2Var : Surface Discretisation Error");
462     }
463
464 // calculate the number and type of autorized cuts
465 // depending on the max number of squares and the validity of next cuts
466     NbU = myResult.NbPatchInU();
467     NbV = myResult.NbPatchInV();
468     NbPatch = NbU*NbV;
469     NbInt = NbU;
470     Umore = UChoice.Value(myResult(FirstNA).U0(), myResult(FirstNA).U1(),Udec);
471     Vmore = VChoice.Value(myResult(FirstNA).V0(), myResult(FirstNA).V1(),Vdec);
472
473     NumDec = 0;
474     if ( ((NbPatch+NbV)<=myMaxPatches) && ((NbPatch+NbU)>myMaxPatches)
475          && (Umore) ) NumDec = 1;
476     if ( ((NbPatch+NbV)>myMaxPatches) && ((NbPatch+NbU)<=myMaxPatches) 
477          && (Vmore) ) NumDec = 2;
478     if ( ((NbPatch+NbV)<=myMaxPatches) && ((NbPatch+NbU)<=myMaxPatches) ) {
479       if ( Umore ) NumDec = 3;
480       if ( (NbV>NbU) && Vmore ) NumDec = 4;
481     }
482     if ( (NbU+1)*(NbV+1)<=myMaxPatches ) {
483       if ( !Umore && !Vmore ) NumDec=0;
484       if ( Umore && !Vmore ) NumDec=1;
485       if ( !Umore && Vmore ) NumDec=2;
486       if ( Umore && Vmore ) NumDec=5;
487     }
488
489 // approximation of the square
490     if ( CritAbs ) {
491       myResult(FirstNA).MakeApprox(myConditions,myConstraints,0);
492     }
493     else {
494       myResult(FirstNA).MakeApprox(myConditions,myConstraints,NumDec);
495     }
496     if (NumDec>=3) NumDec = NumDec - 2;
497
498 // evaluation of the criterion on the square
499     if ( myResult(FirstNA).HasResult() ) {
500       Crit.Value(myResult(FirstNA),myConditions);
501       CritValue = myResult(FirstNA).CritValue();
502       if (m1<CritValue) m1 = CritValue;
503     }
504 // is it necessary to cut ?
505     decision = myResult(FirstNA).CutSense(Crit,NumDec);
506     Standard_Boolean Regular = (Crit.Repartition() ==  AdvApp2Var_Regular);
507 //    Standard_Boolean Regular = Standard_True;
508     if (Regular && decision>0) {
509       NbInt++;
510       InitGrid(NbInt);
511     }
512     else {
513       switch (decision) {
514       case 0 :
515 //      Impossible to cut : the result is preserved
516         if ( myResult(FirstNA).HasResult() ) {
517           myResult(FirstNA).OverwriteApprox();
518         }
519         else {
520           myHasResult =  myDone = Standard_False;
521           Standard_ConstructionError::Raise
522             ("AdvApp2Var_ApproxAFunc2Var : Surface Approximation Error");
523         }
524         break;
525       case 1 :
526 //      It is necessary to cut in U
527         myResult.UpdateInU(Udec);
528         myConstraints.UpdateInU(Udec);
529         break;
530       case 2 :
531 //      It is necessary to cut in V
532         myResult.UpdateInV(Vdec);
533         myConstraints.UpdateInV(Vdec);
534         break;
535       case 3 :
536 //      It is necessary to cut in U and V
537         myResult.UpdateInU(Udec);
538         myConstraints.UpdateInU(Udec);
539         myResult.UpdateInV(Vdec);
540         myConstraints.UpdateInV(Vdec);
541         break;
542         default :
543         myHasResult =  myDone = Standard_False;
544         Standard_ConstructionError::Raise
545           ("AdvApp2Var_ApproxAFunc2Var : Surface Approximation Error");
546       }
547     }
548   }
549 }
550
551 //=======================================================================
552 //function : ComputeConstraints without Criterion
553 //purpose  : Approximation of the constraints
554 //=======================================================================
555
556 void AdvApp2Var_ApproxAFunc2Var::ComputeConstraints(const AdvApprox_Cutting& UChoice,
557                                          const AdvApprox_Cutting& VChoice,
558                                          const AdvApp2Var_EvaluatorFunc2Var& Func)
559 {
560   Standard_Real dec;
561   Standard_Boolean more;
562   Standard_Integer ind1, ind2, NbPatch, NbU, NbV;
563   AdvApp2Var_Iso Is;
564   Standard_Integer indN1, indN2;
565   Standard_Integer iu = myConditions.UOrder(), iv = myConditions.VOrder();
566   AdvApp2Var_Node N1(iu,iv), N2(iu,iv);
567
568   while ( myConstraints.FirstNotApprox(ind1, ind2, Is) ) {
569
570 // approximation of iso and calculation of constraints at extremities
571     indN1 = myConstraints.FirstNode(Is.Type(),ind1,ind2);
572     N1 = myConstraints.Node(indN1);
573     indN2 = myConstraints.LastNode(Is.Type(),ind1,ind2);
574     N2 = myConstraints.Node(indN2);
575
576     Is.MakeApprox(myConditions,
577                   myFirstParInU, myLastParInU,
578                   myFirstParInV, myLastParInV,
579                   Func, N1 , N2);
580
581     if (Is.IsApproximated()) {
582 // iso is approached at the required tolerance
583       myConstraints.ChangeIso(ind1,ind2,Is);
584       myConstraints.ChangeNode(indN1) = N1;
585       myConstraints.ChangeNode(indN2) = N2;
586     }
587
588     else {
589 // Approximation is not satisfactory
590       NbU = myResult.NbPatchInU();
591       NbV = myResult.NbPatchInV();
592       if (Is.Type()==GeomAbs_IsoV) {
593         NbPatch = (NbU+1)*NbV;
594         more = UChoice.Value(Is.T0(),Is.T1(),dec);
595       }
596       else {
597         NbPatch = (NbV+1)*NbU;
598         more = VChoice.Value(Is.T0(),Is.T1(),dec);
599       }
600
601       if (NbPatch<=myMaxPatches && more) {
602 //      It is possible to cut iso
603         if (Is.Type()==GeomAbs_IsoV) {
604           myResult.UpdateInU(dec);
605           myConstraints.UpdateInU(dec);
606         }
607         else {
608           myResult.UpdateInV(dec);
609           myConstraints.UpdateInV(dec);
610         }
611       }
612
613       else {
614 //      It is not possible to cut : the result is preserved
615         if (Is.HasResult()) {
616           Is.OverwriteApprox();
617           myConstraints.ChangeIso(ind1,ind2,Is);
618           myConstraints.ChangeNode(indN1) = N1;
619           myConstraints.ChangeNode(indN2) = N2; 
620         }
621         else {
622           myHasResult =  myDone = Standard_False;
623           Standard_ConstructionError::Raise
624             ("AdvApp2Var_ApproxAFunc2Var : Curve Approximation Error");
625         }
626       }
627     }
628   }
629 }
630
631
632 //=======================================================================
633 //function : ComputeConstraints with Criterion
634 //purpose  : Approximation of the constraints
635 //=======================================================================
636
637 void AdvApp2Var_ApproxAFunc2Var::ComputeConstraints(const AdvApprox_Cutting& UChoice,
638                                          const AdvApprox_Cutting& VChoice,
639                                          const AdvApp2Var_EvaluatorFunc2Var& Func,
640                                          const AdvApp2Var_Criterion& Crit)
641 {
642   Standard_Real dec;
643   Standard_Boolean more, CritRel = (Crit.Type() == AdvApp2Var_Relative);
644   Standard_Integer ind1, ind2, NbPatch, NbU, NbV;
645   AdvApp2Var_Iso Is;
646   Standard_Integer indN1, indN2;
647   Standard_Integer iu = myConditions.UOrder(), iv = myConditions.VOrder();
648   AdvApp2Var_Node N1(iu,iv), N2(iu,iv);
649
650     while ( myConstraints.FirstNotApprox(ind1, ind2, Is) ) {
651
652 // approximation of the iso and calculation of constraints at the extremities
653       indN1 = myConstraints.FirstNode(Is.Type(),ind1,ind2);
654       N1 = myConstraints.Node(indN1);
655       indN2 = myConstraints.LastNode(Is.Type(),ind1,ind2);
656       N2 = myConstraints.Node(indN2);
657
658       Is.MakeApprox(myConditions,
659                     myFirstParInU, myLastParInU,
660                     myFirstParInV, myLastParInV,
661                     Func, N1 , N2);
662
663       if (Is.IsApproximated()) {
664 // iso is approached at the required tolerance
665         myConstraints.ChangeIso(ind1,ind2,Is);
666         myConstraints.ChangeNode(indN1) = N1;
667         myConstraints.ChangeNode(indN2) = N2;
668       }
669
670       else {
671 // Approximation is not satisfactory
672         NbU = myResult.NbPatchInU();
673         NbV = myResult.NbPatchInV();
674         if (Is.Type()==GeomAbs_IsoV) {
675           NbPatch = (NbU+1)*NbV;
676           more = UChoice.Value(Is.T0(),Is.T1(),dec);
677         }
678         else {
679           NbPatch = (NbV+1)*NbU;
680           more = VChoice.Value(Is.T0(),Is.T1(),dec);
681         }
682
683 //      To force Overwrite if the criterion is Absolute
684         more = more && (CritRel);
685
686         if (NbPatch<=myMaxPatches && more) {
687 //      It is possible to cut iso
688           if (Is.Type()==GeomAbs_IsoV) {
689             myResult.UpdateInU(dec);
690             myConstraints.UpdateInU(dec);
691           }
692           else {
693             myResult.UpdateInV(dec);
694             myConstraints.UpdateInV(dec);
695           }
696         }
697
698         else {
699 //      It is not possible to cut: the result is preserved
700           if (Is.HasResult()) {
701             Is.OverwriteApprox();
702             myConstraints.ChangeIso(ind1,ind2,Is);
703             myConstraints.ChangeNode(indN1) = N1;
704             myConstraints.ChangeNode(indN2) = N2;       
705           }
706           else {
707             myHasResult =  myDone = Standard_False;
708             Standard_ConstructionError::Raise
709               ("AdvApp2Var_ApproxAFunc2Var : Curve Approximation Error");
710           }
711         }
712       }
713     }
714 }
715
716 //=======================================================================
717 //function : Compute3DErrors
718 //purpose  : Computation of the 3D errors
719 //=======================================================================
720
721 void AdvApp2Var_ApproxAFunc2Var::Compute3DErrors()
722 {
723
724   Standard_Integer iesp,ipat;
725   Standard_Real error_max,error_moy,error_U0,error_V0,error_U1,error_V1;
726   Standard_Real Tol,F1Tol,F2Tol,F3Tol,F4Tol;
727   if ( myNumSubSpaces[2] > 0 ) {
728     my3DMaxError = new (TColStd_HArray1OfReal) (1,myNumSubSpaces[2]);
729     my3DAverageError = new (TColStd_HArray1OfReal) (1,myNumSubSpaces[2]);
730     my3DUFrontError = new (TColStd_HArray1OfReal) (1,myNumSubSpaces[2]);
731     my3DVFrontError = new (TColStd_HArray1OfReal) (1,myNumSubSpaces[2]);
732     for (iesp=1;iesp<=myNumSubSpaces[2];iesp++) {
733       error_max = 0;
734       error_moy = 0.;
735       error_U0 = 0.;
736       error_V0 = 0.;
737       error_U1 = 0.;
738       error_V1 = 0.;
739       Tol = my3DTolerances->Value(iesp);
740       F1Tol = my3DTolOnFront->Value(iesp,1);
741       F2Tol = my3DTolOnFront->Value(iesp,2);
742       F3Tol = my3DTolOnFront->Value(iesp,3);
743       F4Tol = my3DTolOnFront->Value(iesp,4);
744       for (ipat=1;ipat<=myResult.NbPatch();ipat++) {
745         error_max = Max((myResult(ipat).MaxErrors())->Value(iesp),error_max);
746         error_U0 = Max((myResult(ipat).IsoErrors())->Value(iesp,3),error_U0);
747         error_U1 = Max((myResult(ipat).IsoErrors())->Value(iesp,4),error_U1);
748         error_V0 = Max((myResult(ipat).IsoErrors())->Value(iesp,1),error_V0);
749         error_V1 = Max((myResult(ipat).IsoErrors())->Value(iesp,2),error_V1);
750         error_moy += (myResult(ipat).AverageErrors())->Value(iesp);
751       }
752       my3DMaxError->SetValue(iesp,error_max);
753       my3DUFrontError->SetValue(iesp,Max(error_U0,error_U1));
754       my3DVFrontError->SetValue(iesp,Max(error_V0,error_V1));
755       error_moy /= (Standard_Real) myResult.NbPatch();
756       my3DAverageError->SetValue(iesp,error_moy);
757       if ( error_max>Tol 
758           || error_U0>F3Tol || error_U1>F4Tol 
759           || error_V0>F1Tol || error_V1>F2Tol) {
760         myDone = Standard_False;
761       }
762     }
763   }
764 }
765
766
767 //=======================================================================
768 //function : ComputeCritError
769 //purpose  : Computation of the max value of the Criterion
770 //=======================================================================
771
772 void AdvApp2Var_ApproxAFunc2Var::ComputeCritError()
773 {
774
775   Standard_Integer iesp,ipat;
776   Standard_Real crit_max;
777   if ( myNumSubSpaces[2] > 0 ) {
778     for (iesp=1;iesp<=myNumSubSpaces[2];iesp++) {
779       crit_max = 0.;
780       for (ipat=1;ipat<=myResult.NbPatch();ipat++) {
781         crit_max = Max((myResult(ipat).CritValue()),crit_max);
782       }
783       myCriterionError = crit_max;
784     }
785   }
786 }
787
788 //=======================================================================
789 //function : ConvertBS
790 //purpose  : Convertion of the approximation in BSpline Surface
791 //=======================================================================
792
793 void AdvApp2Var_ApproxAFunc2Var::ConvertBS()
794 {
795  // Homogeneization of degrees
796   Standard_Integer iu = myConditions.UOrder(), iv = myConditions.VOrder();
797   Standard_Integer ncfu = myConditions.ULimit(), ncfv = myConditions.VLimit();
798   myResult.SameDegree(iu,iv,ncfu,ncfv);
799   myDegreeInU = ncfu - 1;
800   myDegreeInV = ncfv - 1;
801
802  // Calculate resulting surfaces 
803   mySurfaces = new ( TColGeom_HArray1OfSurface) (1,  myNumSubSpaces[2]);
804
805   Standard_Integer j;
806   TColStd_Array1OfReal UKnots (1, myResult.NbPatchInU()+1);
807   for (j=1; j<=UKnots.Length(); j++) { UKnots.SetValue(j, myResult.UParameter(j)); } 
808  
809   TColStd_Array1OfReal VKnots (1, myResult.NbPatchInV()+1);
810   for (j=1; j<=VKnots.Length(); j++) { VKnots.SetValue(j, myResult.VParameter(j)); }
811
812  // Prepare data for conversion grid of polynoms --> poles
813   Handle(TColStd_HArray1OfReal) Uint1 = 
814     new (TColStd_HArray1OfReal) (1,2);
815   Uint1->SetValue(1, -1);
816   Uint1->SetValue(2, 1);
817   Handle(TColStd_HArray1OfReal) Vint1 = 
818     new (TColStd_HArray1OfReal) (1,2);
819   Vint1->SetValue(1, -1);
820   Vint1->SetValue(2, 1);
821
822   Handle(TColStd_HArray1OfReal) Uint2 = 
823     new (TColStd_HArray1OfReal) (1,myResult.NbPatchInU()+1);
824   for (j=1; j<=Uint2->Length(); j++) { Uint2->SetValue(j, myResult.UParameter(j)); } 
825   Handle(TColStd_HArray1OfReal) Vint2 = 
826     new (TColStd_HArray1OfReal) (1,myResult.NbPatchInV()+1);
827   for (j=1; j<=Vint2->Length(); j++) { Vint2->SetValue(j, myResult.VParameter(j)); } 
828
829   Standard_Integer nmax = myResult.NbPatchInU()*myResult.NbPatchInV(),
830                    Size_eq = myConditions.ULimit() * myConditions.VLimit() * 3;
831
832   Handle(TColStd_HArray2OfInteger) NbCoeff = 
833     new (TColStd_HArray2OfInteger) (1, nmax, 1, 2);
834   Handle(TColStd_HArray1OfReal) Poly = 
835     new (TColStd_HArray1OfReal) (1, nmax * Size_eq);
836
837   Standard_Integer SSP, i;
838   for (SSP=1; SSP <= myNumSubSpaces[2]; SSP++) { 
839
840     // Creation of the grid of polynoms
841     Standard_Integer n=0,icf=1,ieq;
842     for (j=1; j<=myResult.NbPatchInV(); j++) {
843       for (i=1; i<=myResult.NbPatchInU(); i++) {
844         n++;
845         NbCoeff->SetValue(n,1, myResult.Patch(i,j).NbCoeffInU());
846         NbCoeff->SetValue(n,2, myResult.Patch(i,j).NbCoeffInV());
847         for (ieq=1; ieq<=Size_eq;ieq++) {
848           Poly->SetValue(icf,(myResult.Patch(i,j).Coefficients(SSP,myConditions))
849                                   ->Value(ieq));
850           icf++;
851         }
852       }
853     }
854   
855     // Conversion into poles
856     Convert_GridPolynomialToPoles CvP (myResult.NbPatchInU(),myResult.NbPatchInV(),
857                                        iu,iv,myMaxDegInU,myMaxDegInV,NbCoeff,
858                                        Poly,Uint1,Vint1,Uint2,Vint2);
859     if ( !CvP.IsDone() ) { myDone = Standard_False; }
860    
861     // Conversion into BSpline
862     mySurfaces->ChangeValue(SSP) = new (Geom_BSplineSurface) 
863         ( CvP.Poles()->Array2(),   
864           CvP.UKnots()->Array1(),  CvP.VKnots()->Array1(),
865           CvP.UMultiplicities()->Array1(), CvP.VMultiplicities()->Array1(),
866           CvP.UDegree(), CvP.VDegree() );
867   }
868 }
869
870 //=======================================================================
871 //function : MaxError
872 //purpose  : 
873 //=======================================================================
874
875 Handle(TColStd_HArray1OfReal)
876  AdvApp2Var_ApproxAFunc2Var::MaxError(const Standard_Integer Dimension) const
877 {
878   Handle (TColStd_HArray1OfReal) EPtr;
879   if (Dimension <1 || Dimension >3) {
880     Standard_OutOfRange::Raise
881       ("AdvApp2Var_ApproxAFunc2Var::MaxError : Dimension must be equal to 1,2 or 3 !");
882   }
883   switch (Dimension) {
884   case 1:
885     EPtr = my1DMaxError;
886     break;
887   case 2:
888     EPtr = my2DMaxError;
889     break;
890   case 3:
891     EPtr = my3DMaxError;
892     break;
893   }
894   return EPtr;
895 }
896
897 //=======================================================================
898 //function : AverageError
899 //purpose  : 
900 //=======================================================================
901
902 Handle(TColStd_HArray1OfReal)
903  AdvApp2Var_ApproxAFunc2Var::AverageError(const Standard_Integer Dimension) const 
904 {
905   Handle (TColStd_HArray1OfReal) EPtr;
906   if (Dimension <1 || Dimension >3) {
907     Standard_OutOfRange::Raise
908       ("AdvApp2Var_ApproxAFunc2Var::AverageError : Dimension must be equal to 1,2 or 3 !");
909   }
910   switch (Dimension) {
911   case 1:
912     EPtr = my1DAverageError;
913     break;
914   case 2:
915     EPtr = my2DAverageError;
916     break;
917   case 3:
918     EPtr = my3DAverageError;
919     break;
920   }
921   return EPtr;
922 }
923
924 //=======================================================================
925 //function : UFrontError
926 //purpose  : 
927 //=======================================================================
928
929 Handle(TColStd_HArray1OfReal)
930  AdvApp2Var_ApproxAFunc2Var::UFrontError(const Standard_Integer Dimension) const 
931 {
932   Handle (TColStd_HArray1OfReal) EPtr;
933   if (Dimension <1 || Dimension >3) {
934     Standard_OutOfRange::Raise
935       ("AdvApp2Var_ApproxAFunc2Var::UFrontError : Dimension must be equal to 1,2 or 3 !");
936   }
937   switch (Dimension) {
938   case 1:
939     EPtr = my1DUFrontError;
940     break;
941   case 2:
942     EPtr = my2DUFrontError;
943     break;
944   case 3:
945     EPtr = my3DUFrontError;
946     break;
947   }
948   return EPtr;
949 }
950
951 //=======================================================================
952 //function : VFrontError
953 //purpose  : 
954 //=======================================================================
955
956 Handle(TColStd_HArray1OfReal)
957  AdvApp2Var_ApproxAFunc2Var::VFrontError(const Standard_Integer Dimension) const 
958 {
959   Handle (TColStd_HArray1OfReal) EPtr;
960   if (Dimension <=0 || Dimension >3) {
961     Standard_OutOfRange::Raise
962       ("AdvApp2Var_ApproxAFunc2Var::VFrontError : Dimension must be equal to 1,2 or 3 !");
963   }
964   switch (Dimension) {
965   case 1:
966     EPtr = my1DVFrontError;
967     break;
968   case 2:
969     EPtr = my2DVFrontError;
970     break;
971   case 3:
972     EPtr = my3DVFrontError;
973     break;
974   }
975   return EPtr;
976 }
977
978 //=======================================================================
979 //function : MaxError
980 //purpose  : 
981 //=======================================================================
982
983 Standard_Real
984  AdvApp2Var_ApproxAFunc2Var::MaxError(const Standard_Integer Dimension,
985                                       const Standard_Integer SSPIndex) const 
986 {
987   if (Dimension !=3 || SSPIndex !=1) {
988     Standard_OutOfRange::Raise
989       ("AdvApp2Var_ApproxAFunc2Var::MaxError: ONE Surface 3D only !");
990   }
991   Handle (TColStd_HArray1OfReal) EPtr = MaxError(Dimension);
992   return EPtr->Value(SSPIndex);
993 }
994
995 //=======================================================================
996 //function : AverageError
997 //purpose  : 
998 //=======================================================================
999
1000 Standard_Real
1001  AdvApp2Var_ApproxAFunc2Var::AverageError(const Standard_Integer Dimension,
1002                                           const Standard_Integer SSPIndex) const 
1003 {
1004   if (Dimension !=3 || SSPIndex !=1) {
1005     Standard_OutOfRange::Raise
1006       ("AdvApp2Var_ApproxAFunc2Var::AverageError : ONE Surface 3D only !");
1007   }
1008   Handle (TColStd_HArray1OfReal) EPtr = AverageError(Dimension);
1009   return EPtr->Value(SSPIndex);
1010 }
1011
1012 //=======================================================================
1013 //function : UFrontError
1014 //purpose  : 
1015 //=======================================================================
1016
1017 Standard_Real
1018  AdvApp2Var_ApproxAFunc2Var::UFrontError(const Standard_Integer Dimension,
1019                                          const Standard_Integer SSPIndex) const 
1020 {
1021   if (Dimension !=3 || SSPIndex !=1) {
1022     Standard_OutOfRange::Raise
1023       ("AdvApp2Var_ApproxAFunc2Var::UFrontError : ONE Surface 3D only !");
1024   }
1025   Handle (TColStd_HArray1OfReal) EPtr = UFrontError(Dimension);
1026   return EPtr->Value(SSPIndex);
1027 }
1028
1029 //=======================================================================
1030 //function : VFrontError
1031 //purpose  : 
1032 //=======================================================================
1033
1034 Standard_Real
1035  AdvApp2Var_ApproxAFunc2Var::VFrontError(const Standard_Integer Dimension,
1036                                          const Standard_Integer SSPIndex) const 
1037 {
1038   if (Dimension !=3 || SSPIndex !=1) {
1039     Standard_OutOfRange::Raise
1040       ("AdvApp2Var_ApproxAFunc2Var::VFrontError : ONE Surface 3D only !");
1041   }
1042   Handle (TColStd_HArray1OfReal) EPtr = VFrontError(Dimension);
1043   return EPtr->Value(SSPIndex);
1044 }
1045
1046
1047 //=======================================================================
1048 //function : CritError
1049 //purpose  : 
1050 //=======================================================================
1051
1052 Standard_Real
1053  AdvApp2Var_ApproxAFunc2Var::CritError(const Standard_Integer Dimension,
1054                                        const Standard_Integer SSPIndex) const 
1055 {
1056   if (Dimension !=3 || SSPIndex !=1) {
1057     Standard_OutOfRange::Raise
1058       ("AdvApp2Var_ApproxAFunc2Var::CritError: ONE Surface 3D only !");
1059   }
1060   return myCriterionError;
1061 }
1062
1063 //=======================================================================
1064 //function : Dump
1065 //purpose  : 
1066 //=======================================================================
1067
1068 void AdvApp2Var_ApproxAFunc2Var::Dump(Standard_OStream& o) const 
1069 {
1070   Standard_Integer iesp=1,NbKU,NbKV,ik;
1071   o<<endl;
1072   if (!myHasResult) { o<<"No result"<<endl; }
1073   else {
1074     o<<"There is a result";
1075     if (myDone) {
1076       o<<" within the requested tolerance "<<my3DTolerances->Value(iesp)<<endl;
1077     }
1078     else if (my3DMaxError->Value(iesp)>my3DTolerances->Value(iesp)) {
1079       o<<" WITHOUT the requested tolerance "<<my3DTolerances->Value(iesp)<<endl;
1080     }
1081     else {
1082       o<<" WITHOUT the requested continuities "<<endl;
1083     }
1084     o<<endl;
1085     o<<"Result max error :"<<my3DMaxError->Value(iesp)<<endl;
1086     o<<"Result average error :"<<my3DAverageError->Value(iesp)<<endl;
1087     o<<"Result max error on U frontiers :"<<my3DUFrontError->Value(iesp)<<endl;
1088     o<<"Result max error on V frontiers :"<<my3DVFrontError->Value(iesp)<<endl;
1089     o<<endl;
1090     o<<"Degree of Bezier patches in U : "<<myDegreeInU
1091       <<"  in V : "<<myDegreeInV<<endl;
1092     o<<endl;
1093     Handle(Geom_BSplineSurface) S
1094       = Handle(Geom_BSplineSurface)::DownCast(mySurfaces->Value(iesp));
1095     o<<"Number of poles in U : "<<S->NbUPoles()
1096       <<"  in V : "<<S->NbVPoles()<<endl;
1097     o<<endl;
1098     NbKU = S->NbUKnots();
1099     NbKV = S->NbVKnots();
1100     o<<"Number of knots in U : "<<NbKU<<endl;
1101     for (ik=1;ik<=NbKU;ik++) {
1102       o<<"   "<<ik<<" : "<<S->UKnot(ik)<<"   mult : "<<S->UMultiplicity(ik)<<endl;
1103     }
1104     o<<endl;
1105     o<<"Number of knots in V : "<<NbKV<<endl;
1106     for (ik=1;ik<=NbKV;ik++) {
1107       o<<"   "<<ik<<" : "<<S->VKnot(ik)<<"   mult : "<<S->VMultiplicity(ik)<<endl;
1108     }
1109     o<<endl;
1110   }
1111 }