ea6742084060df52398d66e191d74c91dddcd499
[occt.git] / src / AppParCurves / AppParCurves_Variational_6.gxx
1 // Created on: 1998-12-21
2 // Created by: Igor FEOKTISTOV
3 // Copyright (c) 1998-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 //  Modified by skv - Fri Apr  8 14:58:12 2005 OCC8559
18
19 #include <PLib_Base.hxx>
20 #include <PLib_JacobiPolynomial.hxx>
21 #include <PLib_HermitJacobi.hxx>
22 #include <FEmTool_Curve.hxx>
23 #include <math_Vector.hxx>
24 #include <TColStd_Array1OfReal.hxx>
25
26 //=======================================================================
27 //function : InitSmoothCriterion
28 //purpose  : Initializes the SmoothCriterion
29 //=======================================================================
30 void AppParCurves_Variational::InitSmoothCriterion()
31 {
32   
33   const Standard_Real Eps2 = 1.e-6, Eps3 = 1.e-9;
34 //  const Standard_Real J1 = .01, J2 = .001, J3 = .001;
35
36
37
38   Standard_Real Length;
39
40   InitParameters(Length);
41
42   mySmoothCriterion->SetParameters(myParameters);
43
44   Standard_Real E1, E2, E3;
45
46   InitCriterionEstimations(Length, E1,  E2,  E3);
47 /*
48   J1 = 1.e-8; J2 = J3 = (E1 + 1.e-8) * 1.e-6;
49
50   if(E1 < J1) E1 = J1;
51   if(E2 < J2) E2 = J2;
52   if(E3 < J3) E3 = J3;
53 */
54   mySmoothCriterion->EstLength() = Length;
55   mySmoothCriterion->SetEstimation(E1, E2, E3);
56
57   Standard_Real WQuadratic, WQuality;
58
59   if(!myWithMinMax && myTolerance != 0.) 
60     WQuality = myTolerance;
61   else if (myTolerance == 0.)
62     WQuality = 1.;
63   else
64     WQuality = Max(myTolerance, Eps2 * Length);
65
66   Standard_Integer NbConstr = myNbPassPoints + myNbTangPoints + myNbCurvPoints;
67   WQuadratic = Sqrt((Standard_Real)(myNbPoints - NbConstr)) * WQuality;
68   if(WQuadratic > Eps3) WQuadratic = 1./ WQuadratic;
69
70   if(WQuadratic == 0.) WQuadratic = Max(Sqrt(E1), 1.);
71
72
73   mySmoothCriterion->SetWeight(WQuadratic, WQuality, 
74                                myPercent[0], myPercent[1], myPercent[2]);
75
76
77   Handle(PLib_Base) TheBase = new PLib_HermitJacobi(myMaxDegree, myContinuity);
78   Handle(FEmTool_Curve) TheCurve;
79   Standard_Integer NbElem;
80   Standard_Real CurvTol = Eps2 * Length / myNbPoints;
81
82   // Decoupe de l'intervalle en fonction des contraintes
83   if(myWithCutting == Standard_True && NbConstr != 0) {
84
85     InitCutting(TheBase, CurvTol, TheCurve);
86
87   }
88   else {
89
90     NbElem = 1;
91     TheCurve = new FEmTool_Curve(myDimension, NbElem, TheBase, CurvTol);
92     TheCurve->Knots().SetValue(TheCurve->Knots().Lower(), 
93                                myParameters->Value(myFirstPoint));
94     TheCurve->Knots().SetValue(TheCurve->Knots().Upper(), 
95                                myParameters->Value(myLastPoint));
96
97   }
98                              
99   mySmoothCriterion->SetCurve(TheCurve);
100
101   return;
102   
103 }
104
105 //
106 //=======================================================================
107 //function : InitParameters
108 //purpose  : Calculation of initial estimation of the multicurve's length
109 //           and parameters for multipoints.
110 //=======================================================================
111 //
112 void AppParCurves_Variational::InitParameters(Standard_Real& Length)
113 {
114
115   const Standard_Real Eps1 = Precision::Confusion() * .01;
116
117   Standard_Real aux, dist;
118   Standard_Integer i, i0, i1 = 0, ipoint;
119
120
121   Length = 0.;  
122   myParameters->SetValue(myFirstPoint, Length);
123
124   for(ipoint = myFirstPoint + 1; ipoint <= myLastPoint; ipoint++) {
125     i0 = i1; i1 += myDimension;
126     dist = 0;
127     for(i = 1; i <= myDimension; i++) {
128       aux = myTabPoints->Value(i1 + i) - myTabPoints->Value(i0 + i);
129       dist += aux * aux;
130     }
131     Length += Sqrt(dist);
132     myParameters->SetValue(ipoint, Length);
133   }
134
135
136   if(Length <= Eps1) 
137     Standard_ConstructionError::Raise("AppParCurves_Variational::InitParameters");
138
139
140   for(ipoint = myFirstPoint + 1; ipoint <= myLastPoint - 1; ipoint++) 
141     myParameters->SetValue(ipoint, myParameters->Value(ipoint) / Length);
142   
143   myParameters->SetValue(myLastPoint, 1.);
144   
145   // Avec peu de point il y a sans doute sous estimation ...
146   if(myNbPoints < 10) Length *= (1. + 0.1 / (myNbPoints - 1));
147 }
148
149 //
150 //=======================================================================
151 //function : InitCriterionEstimations
152 //purpose  : Calculation of initial estimation of the linear criteria
153 //           
154 //=======================================================================
155 //
156 void AppParCurves_Variational::InitCriterionEstimations(const Standard_Real Length, 
157                                                         Standard_Real& E1,  
158                                                         Standard_Real& E2,  
159                                                         Standard_Real& E3) const
160 {
161   E1 = Length * Length;
162
163   const Standard_Real Eps1 = Precision::Confusion() * .01;
164
165   math_Vector VTang1(1, myDimension), VTang2(1, myDimension), VTang3(1, myDimension),
166               VScnd1(1, myDimension), VScnd2(1, myDimension), VScnd3(1, myDimension);
167
168   // ========== Treatment of first point =================
169
170   Standard_Integer ipnt = myFirstPoint;
171
172   EstTangent(ipnt, VTang1);
173   ipnt++;
174   EstTangent(ipnt, VTang2);
175   ipnt++;
176   EstTangent(ipnt, VTang3);
177
178   ipnt = myFirstPoint;
179   EstSecnd(ipnt, VTang1, VTang2, Length, VScnd1);
180   ipnt++;
181   EstSecnd(ipnt, VTang1, VTang3, Length, VScnd2);
182
183 //  Modified by skv - Fri Apr  8 14:58:12 2005 OCC8559 Begin
184 //   Standard_Real Delta = .5 * (myParameters->Value(ipnt) - myParameters->Value(--ipnt));
185   Standard_Integer anInd = ipnt;
186   Standard_Real    Delta = .5 * (myParameters->Value(anInd) - myParameters->Value(--ipnt));
187 //  Modified by skv - Fri Apr  8 14:58:12 2005 OCC8559 End
188
189   if(Delta <= Eps1) Delta = 1.;
190
191   E2 = VScnd1.Norm2() * Delta;
192
193   E3 = (Delta > Eps1) ? VScnd2.Subtracted(VScnd1).Norm2() / (4. * Delta) : 0.;
194   // ========== Treatment of internal points =================
195
196   Standard_Integer CurrPoint = 2;
197
198   for(ipnt = myFirstPoint + 1; ipnt < myLastPoint; ipnt++) {
199
200     Delta = .5 * (myParameters->Value(ipnt + 1) - myParameters->Value(ipnt - 1));
201
202     if(CurrPoint == 1) {
203       if(ipnt + 1 != myLastPoint) {
204         EstTangent(ipnt + 2, VTang3);
205         EstSecnd(ipnt + 1, VTang1, VTang3, Length, VScnd2);
206       }
207       else
208         EstSecnd(ipnt + 1, VTang1, VTang2, Length, VScnd2);
209
210       E2 += VScnd1.Norm2() * Delta;
211       E3 += (Delta > Eps1) ? VScnd2.Subtracted(VScnd3).Norm2() / (4. * Delta) : 0.;
212
213     }
214     else if(CurrPoint == 2) {
215       if(ipnt + 1 != myLastPoint) {
216         EstTangent(ipnt + 2, VTang1);
217         EstSecnd(ipnt + 1, VTang2, VTang1, Length, VScnd3);
218       }
219       else
220         EstSecnd(ipnt + 1, VTang2, VTang3, Length, VScnd3);
221
222       E2 += VScnd2.Norm2() * Delta;
223       E3 += (Delta > Eps1) ? VScnd3.Subtracted(VScnd1).Norm2() / (4. * Delta) : 0.;
224
225     }
226     else {
227       if(ipnt + 1 != myLastPoint) {
228         EstTangent(ipnt + 2, VTang2);
229         EstSecnd(ipnt + 1, VTang3, VTang2, Length, VScnd1);
230       }
231       else
232         EstSecnd(ipnt + 1, VTang3, VTang1, Length, VScnd1);
233         
234       E2 += VScnd3.Norm2() * Delta;
235       E3 += (Delta > Eps1) ? VScnd1.Subtracted(VScnd2).Norm2() / (4. * Delta) : 0.;
236
237     }
238
239     CurrPoint++; if(CurrPoint == 4) CurrPoint = 1;
240   }
241
242   // ========== Treatment of last point =================
243
244   Delta = .5 * (myParameters->Value(myLastPoint) - myParameters->Value(myLastPoint - 1));
245   if(Delta <= Eps1) Delta = 1.;
246
247   Standard_Real aux;
248
249   if(CurrPoint == 1) {
250
251     E2 += VScnd1.Norm2() * Delta;
252     aux = VScnd1.Subtracted(VScnd3).Norm2();
253     E3 += (Delta > Eps1) ? aux / (4. * Delta) : aux;
254     
255   }
256   else if(CurrPoint == 2) {
257
258     E2 += VScnd2.Norm2() * Delta;
259     aux = VScnd2.Subtracted(VScnd1).Norm2();
260     E3 += (Delta > Eps1) ? aux / (4. * Delta) : aux;
261
262   }
263   else {
264     
265     E2 += VScnd3.Norm2() * Delta;
266     aux = VScnd3.Subtracted(VScnd2).Norm2();
267     E3 += (Delta > Eps1) ? aux / (4. * Delta) : aux;
268     
269   }
270
271   aux = Length * Length;
272
273   E2 *= aux; E3 *= aux;
274
275   
276 }
277
278 //
279 //=======================================================================
280 //function : EstTangent
281 //purpose  : Calculation of  estimation of the Tangent
282 //           (see fortran routine MMLIPRI)
283 //=======================================================================
284 //
285
286 void AppParCurves_Variational::EstTangent(const Standard_Integer ipnt, 
287                                           math_Vector& VTang) const
288
289 {
290   Standard_Integer i ;
291   const Standard_Real Eps1 = Precision::Confusion() * .01;
292   const Standard_Real EpsNorm = 1.e-9;
293
294   Standard_Real Wpnt = 1.;
295
296
297   if(ipnt == myFirstPoint) {
298     // Estimation at first point
299     if(myNbPoints < 3) 
300       Wpnt = 0.;
301     else {
302
303       Standard_Integer adr1 = 1, adr2 = adr1 + myDimension, 
304                        adr3 = adr2 + myDimension;
305
306       math_Vector Pnt1((Standard_Real*)&myTabPoints->Value(adr1), 1, myDimension);
307       math_Vector Pnt2((Standard_Real*)&myTabPoints->Value(adr2), 1, myDimension);
308       math_Vector Pnt3((Standard_Real*)&myTabPoints->Value(adr3), 1, myDimension);
309
310       // Parabolic interpolation
311       // if we have parabolic interpolation: F(t) = A0 + A1*t + A2*t*t,
312       // first derivative for t=0 is A1 = ((d2-1)*P1 + P2 - d2*P3)/(d*(1-d))
313       //       d= |P2-P1|/(|P2-P1|+|P3-P2|), d2 = d*d
314       Standard_Real V1 = Pnt2.Subtracted(Pnt1).Norm();
315       Standard_Real V2 = 0.;
316       if(V1 > Eps1) V2 = Pnt3.Subtracted(Pnt2).Norm();
317       if(V2 > Eps1) {
318         Standard_Real d = V1 / (V1 + V2), d1;
319         d1 = 1. / (d * (1 - d)); d *= d;
320         VTang = ((d - 1.) * Pnt1 + Pnt2 - d * Pnt3) * d1;
321       }
322       else {
323         // Simple 2-point estimation
324
325         VTang = Pnt2 - Pnt1;
326       }
327     }
328   }
329   else if (ipnt == myLastPoint) {
330     // Estimation at last point
331     if(myNbPoints < 3) 
332       Wpnt = 0.;
333     else {
334
335       Standard_Integer adr1 = (myLastPoint - 3) * myDimension + 1, 
336                        adr2 = adr1 + myDimension, 
337                        adr3 = adr2 + myDimension;
338
339       math_Vector Pnt1((Standard_Real*)&myTabPoints->Value(adr1), 1, myDimension);
340       math_Vector Pnt2((Standard_Real*)&myTabPoints->Value(adr2), 1, myDimension);
341       math_Vector Pnt3((Standard_Real*)&myTabPoints->Value(adr3), 1, myDimension);
342
343       // Parabolic interpolation
344       // if we have parabolic interpolation: F(t) = A0 + A1*t + A2*t*t,
345       // first derivative for t=1 is 2*A2 + A1 = ((d2+1)*P1 - P2 - d2*P3)/(d*(1-d))
346       //       d= |P2-P1|/(|P2-P1|+|P3-P2|), d2 = d*(d-2)
347       Standard_Real V1 = Pnt2.Subtracted(Pnt1).Norm();
348       Standard_Real V2 = 0.;
349       if(V1 > Eps1) V2 = Pnt3.Subtracted(Pnt2).Norm();
350       if(V2 > Eps1) {
351         Standard_Real d = V1 / (V1 + V2), d1;
352         d1 = 1. / (d * (1 - d)); d *= d - 2;
353         VTang = ((d + 1.) * Pnt1 - Pnt2 - d * Pnt3) * d1;
354       }
355       else {
356         // Simple 2-point estimation
357         
358         VTang = Pnt3 - Pnt2;
359       }
360     }
361   }
362   else {
363
364     Standard_Integer adr1 = (ipnt - myFirstPoint - 1) * myDimension + 1, 
365                      adr2 = adr1 + 2 * myDimension;
366     
367     math_Vector Pnt1((Standard_Real*)&myTabPoints->Value(adr1), 1, myDimension);
368     math_Vector Pnt2((Standard_Real*)&myTabPoints->Value(adr2), 1, myDimension);
369
370     VTang = Pnt2 - Pnt1;
371
372   }
373
374   Standard_Real Vnorm = VTang.Norm();
375   
376   if(Vnorm <= EpsNorm)
377     VTang.Init(0.);
378   else 
379     VTang /= Vnorm;
380
381   // Estimation with constraints
382   
383   Standard_Real Wcnt = 0.;
384   Standard_Integer IdCnt = 1;
385
386 // Warning!! Here it is suppoused that all points are in range [myFirstPoint, myLastPoint]
387   
388   Standard_Integer NbConstr = myNbPassPoints + myNbTangPoints + myNbCurvPoints;
389   math_Vector VCnt(1, myDimension, 0.);
390
391   if(NbConstr > 0) {
392
393     while(myTypConstraints->Value(2 * IdCnt - 1) < ipnt &&
394           IdCnt <= NbConstr) IdCnt++;
395     if((myTypConstraints->Value(2 * IdCnt - 1) == ipnt) &&
396        (myTypConstraints->Value(2 * IdCnt) >= 1)) {
397       Wcnt = 1.;
398       Standard_Integer i0 = 2 * myDimension * (IdCnt - 1), k = 0;
399       for( i = 1; i <= myNbP3d; i++) {
400         for(Standard_Integer j = 1; j <= 3; j++) 
401           VCnt(++k) = myTabConstraints->Value(++i0);
402         i0 += 3;
403       }
404       for(i = 1; i <= myNbP2d; i++) {
405         for(Standard_Integer j = 1; j <= 2; j++) 
406           VCnt(++k) = myTabConstraints->Value(++i0);
407         i0 += 2;
408       }
409     }
410   }
411
412   // Averaging of estimation
413
414   Standard_Real Denom = Wpnt + Wcnt;
415   if(Denom == 0.) Denom = 1.;
416   else            Denom = 1. / Denom;
417
418   VTang = (Wpnt * VTang + Wcnt * VCnt) * Denom;
419
420   Vnorm = VTang.Norm();
421
422   if(Vnorm <= EpsNorm)
423     VTang.Init(0.);
424   else 
425     VTang /= Vnorm;
426
427     
428 }
429
430 //
431 //=======================================================================
432 //function : EstSecnd
433 //purpose  : Calculation of  estimation of the second derivative
434 //           (see fortran routine MLIMSCN)
435 //=======================================================================
436 //
437 void AppParCurves_Variational::EstSecnd(const Standard_Integer ipnt, 
438                                         const math_Vector& VTang1, 
439                                         const math_Vector& VTang2, 
440                                         const Standard_Real Length, 
441                                         math_Vector& VScnd) const
442 {
443   Standard_Integer i ;
444
445   const Standard_Real Eps = 1.e-9;
446
447   Standard_Real Wpnt = 1.;
448
449   Standard_Real aux;
450
451   if(ipnt == myFirstPoint)
452     aux = myParameters->Value(ipnt + 1) - myParameters->Value(ipnt);
453   else if(ipnt == myLastPoint)
454     aux = myParameters->Value(ipnt) - myParameters->Value(ipnt - 1);
455   else
456     aux = myParameters->Value(ipnt + 1) - myParameters->Value(ipnt - 1);
457
458   if(aux <= Eps)
459     aux = 1.;
460   else
461     aux = 1. / aux;
462
463   VScnd = (VTang2 - VTang1) * aux;
464     
465   // Estimation with constraints
466   
467   Standard_Real Wcnt = 0.;
468   Standard_Integer IdCnt = 1;
469
470 // Warning!! Here it is suppoused that all points are in range [myFirstPoint, myLastPoint]
471
472   Standard_Integer NbConstr = myNbPassPoints + myNbTangPoints + myNbCurvPoints;
473   math_Vector VCnt(1, myDimension, 0.);
474
475   if(NbConstr > 0) {
476
477     while(myTypConstraints->Value(2 * IdCnt - 1) < ipnt &&
478           IdCnt <= NbConstr) IdCnt++;
479
480     if((myTypConstraints->Value(2 * IdCnt - 1) == ipnt) &&
481        (myTypConstraints->Value(2 * IdCnt) >= 2))  {
482       Wcnt = 1.;
483       Standard_Integer i0 = 2 * myDimension * (IdCnt - 1) + 3, k = 0;
484       for( i = 1; i <= myNbP3d; i++) {
485         for(Standard_Integer j = 1; j <= 3; j++) 
486           VCnt(++k) = myTabConstraints->Value(++i0);
487         i0 += 3;
488       }
489       i0--;
490       for(i = 1; i <= myNbP2d; i++) {
491         for(Standard_Integer j = 1; j <= 2; j++) 
492           VCnt(++k) = myTabConstraints->Value(++i0);
493         i0 += 2;
494       }
495     }
496   }
497
498   // Averaging of estimation
499
500   Standard_Real Denom = Wpnt + Wcnt;
501   if(Denom == 0.) Denom = 1.;
502   else            Denom = 1. / Denom;
503
504   VScnd = (Wpnt * VScnd + (Wcnt * Length) * VCnt) * Denom;
505
506 }
507
508
509 //
510 //=======================================================================
511 //function : InitCutting
512 //purpose  : Realisation of curve's cutting a priory accordingly to 
513 //           constraints (see fortran routine MLICUT)
514 //=======================================================================
515 //
516 void AppParCurves_Variational::InitCutting(const Handle(PLib_Base)& aBase,
517                                            const Standard_Real CurvTol, 
518                                            Handle(FEmTool_Curve)& aCurve) const
519 {
520
521     // Definition of number of elements
522     Standard_Integer ORCMx = -1, NCont = 0, i, kk, NbElem;
523     Standard_Integer NbConstr = myNbPassPoints + myNbTangPoints + myNbCurvPoints;
524
525     for(i = 1; i <= NbConstr; i++) {
526       kk = Abs(myTypConstraints->Value(2 * i)) + 1;
527       ORCMx = Max(ORCMx, kk);
528       NCont += kk;
529     }
530
531     if(ORCMx > myMaxDegree - myNivCont) 
532       Standard_ConstructionError::Raise("AppParCurves_Variational::InitCutting");
533
534     Standard_Integer NLibre = Max(myMaxDegree - myNivCont - (myMaxDegree + 1) / 4,
535                                   myNivCont + 1);
536
537     NbElem = (NCont % NLibre == 0) ? NCont / NLibre : NCont / NLibre + 1;
538
539     while((NbElem > myMaxSegment) && (NLibre < myMaxDegree - myNivCont)) {
540
541       NLibre++;
542       NbElem = (NCont % NLibre == 0) ? NCont / NLibre : NCont / NLibre + 1;
543
544     }
545
546
547     if(NbElem > myMaxSegment) 
548       Standard_ConstructionError::Raise("AppParCurves_Variational::InitCutting");
549
550
551     aCurve = new FEmTool_Curve(myDimension, NbElem, aBase, CurvTol);
552
553     Standard_Integer NCnt = (NCont - 1) / NbElem + 1;
554     Standard_Integer NPlus = NbElem - (NCnt * NbElem - NCont);
555
556     TColStd_Array1OfReal& Knot = aCurve->Knots();
557
558     Standard_Integer IDeb = 0, IFin = NbConstr + 1,
559                      NDeb = 0, NFin = 0,
560                      IndEl = Knot.Lower(), IUpper = Knot.Upper(),
561                      NbEl = 0;
562
563
564     Knot(IndEl) = myParameters->Value(myFirstPoint);
565     Knot(IUpper) = myParameters->Value(myLastPoint);
566
567     while(NbElem - NbEl > 1) {
568
569       IndEl++; NbEl++;
570       if(NPlus == 0) NCnt--;
571
572       while(NDeb < NCnt && IDeb < IFin) {
573         IDeb++;
574         NDeb += Abs(myTypConstraints->Value(2 * IDeb)) + 1;
575       }
576
577       if(NDeb == NCnt) {
578         NDeb = 0;
579         if(NPlus == 1 && 
580            myParameters->Value(myTypConstraints->Value(2 * IDeb - 1)) > Knot(IndEl - 1))
581           
582           Knot(IndEl) = myParameters->Value(myTypConstraints->Value(2 * IDeb - 1));
583         else
584           Knot(IndEl) = (myParameters->Value(myTypConstraints->Value(2 * IDeb - 1)) +
585                          myParameters->Value(myTypConstraints->Value(2 * IDeb + 1))) / 2;
586
587       }
588       else {
589         NDeb -= NCnt;
590         Knot(IndEl) = myParameters->Value(myTypConstraints->Value(2 * IDeb - 1));
591       }
592
593       NPlus--;
594       if(NPlus == 0) NCnt--;
595
596       if(NbElem - NbEl == 1) break;
597
598       NbEl++;
599
600       while(NFin < NCnt && IDeb < IFin) {
601         IFin--;
602         NFin += Abs(myTypConstraints->Value(2 * IFin)) + 1;
603       }
604
605       if(NFin == NCnt) {
606         NFin = 0;
607         Knot(IUpper + 1 - IndEl) = (myParameters->Value(myTypConstraints->Value(2 * IFin - 1)) +
608                                 myParameters->Value(myTypConstraints->Value(2 * IFin - 3))) / 2;
609       }
610       else {
611         NFin -= NCnt;
612         if(myParameters->Value(myTypConstraints->Value(2 * IFin - 1)) < Knot(IUpper - IndEl + 1))
613           Knot(IUpper + 1 - IndEl) = myParameters->Value(myTypConstraints->Value(2 * IFin - 1));
614         else
615           Knot(IUpper + 1 - IndEl) = (myParameters->Value(myTypConstraints->Value(2 * IFin - 1)) +
616                                   myParameters->Value(myTypConstraints->Value(2 * IFin - 3))) / 2;
617       }
618
619     }
620
621
622 }