projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0024239: Eliminate GCC compiler warning (multi-line comment)
[occt.git] / src / Approx / Approx_SameParameter.cxx
CommitLineData
b311480e 1// Created on: 1995-06-06
2// Created by: Xavier BENVENISTE
3// Copyright (c) 1995-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
7fd59977 21
22// Modified by skv - Wed Jun 2 11:49:59 2004 OCC5898
23
24#include <Approx_SameParameter.ixx>
25#include <TColStd_Array1OfReal.hxx>
26#include <BSplCLib.hxx>
27#include <Adaptor3d_CurveOnSurface.hxx>
28#include <Geom2dAdaptor_Curve.hxx>
29#include <Geom2dAdaptor_HCurve.hxx>
30#include <GeomAdaptor_Curve.hxx>
31#include <GeomAdaptor_HCurve.hxx>
32#include <GeomAdaptor_Surface.hxx>
33#include <GeomAdaptor_HSurface.hxx>
34//#include <GCPnts_UniformDeflection.hxx>
35#include <GCPnts_QuasiUniformDeflection.hxx>
36#include <Extrema_LocateExtPC.hxx>
37#include <AdvApprox_ApproxAFunction.hxx>
38#include <GeomLib_MakeCurvefromApprox.hxx>
39#include <Precision.hxx>
40
41#define MAX_ARRAY_SIZE 1000 // IFV, Jan 2000
42
43#ifdef DEB
44#ifdef DRAW
45#include <DrawTrSurf.hxx>
46#endif
47#include <Geom2d_BSplineCurve.hxx>
48#include <stdio.h>
49static Standard_Boolean Voir = Standard_False;
50static Standard_Boolean AffichFw = Standard_False;
51static Standard_Integer NbCurve = 0;
52#endif
53//
0d969553 54// allows testing if Extrema produces correct results/
7fd59977 55
56
57static void ProjectPointOnCurve(const Standard_Real InitValue,
58 const gp_Pnt APoint,
59 const Standard_Real Tolerance,
60 const Standard_Integer NumIteration,
61 const Adaptor3d_Curve& Curve,
62 Standard_Boolean& Status,
63 Standard_Real& Result)
64{
65 Standard_Integer num_iter = 0,
66 not_done = 1,
67 ii ;
68
69 gp_Pnt a_point ;
70 gp_Vec vector,
71 d1,
72 d2 ;
73 Standard_Real func,
74 func_derivative,
75 param = InitValue ;
76 Status = Standard_False ;
77 Standard_Real Toler = 1.0e-12;
78 do {
79 num_iter += 1 ;
80 Curve.D2(param,
81 a_point,
82 d1,
83 d2) ;
84 for (ii = 1 ; ii <= 3 ; ii++) {
85 vector.SetCoord(ii, APoint.Coord(ii) - a_point.Coord(ii)) ;
86 }
87 func = vector.Dot(d1) ;
88 func_derivative = vector.Dot(d2) ;
89 func_derivative -= d1.Dot(d1) ;
90 if ( Abs(func) < Tolerance * d1.Magnitude()) {
91 not_done = 0 ;
92 Status = Standard_True ;
93 }
94 else
95 { // fixing a bug PRO18577 : avoid divizion by zero
96 if( Abs(func_derivative) > Toler ) {
97 param -= func / func_derivative ;
98 }
99 param = Max(param,Curve.FirstParameter()) ;
100 param = Min(param,Curve.LastParameter()) ;
101 Status = Standard_True ;
102 }
103 }
104 while (not_done && num_iter <= NumIteration) ;
105 Result = param ;
106}
107
108
109
110//=======================================================================
111//class : Approx_SameParameter_Evaluator
112//purpose :
113//=======================================================================
114
115class Approx_SameParameter_Evaluator : public AdvApprox_EvaluatorFunction
116{
117 public:
118 Approx_SameParameter_Evaluator (const TColStd_Array1OfReal& theFlatKnots,
119 const TColStd_Array1OfReal& thePoles,
120 const Handle(Adaptor2d_HCurve2d)& theHCurve2d)
121 : FlatKnots(theFlatKnots), Poles(thePoles), HCurve2d(theHCurve2d) {}
122
123 virtual void Evaluate (Standard_Integer *Dimension,
124 Standard_Real StartEnd[2],
125 Standard_Real *Parameter,
126 Standard_Integer *DerivativeRequest,
127 Standard_Real *Result, // [Dimension]
128 Standard_Integer *ErrorCode);
129
130 private:
131 const TColStd_Array1OfReal& FlatKnots;
132 const TColStd_Array1OfReal& Poles;
133 Handle(Adaptor2d_HCurve2d) HCurve2d;
134};
135
136void Approx_SameParameter_Evaluator::Evaluate (Standard_Integer *,/*Dimension*/
137 Standard_Real /*StartEnd*/[2],
138 Standard_Real *Parameter,
139 Standard_Integer *DerivativeRequest,
140 Standard_Real *Result,
141 Standard_Integer *ReturnCode)
142{
143 gp_Pnt2d Point ;
144 gp_Vec2d Vector ;
145 Standard_Integer extrap_mode[2] ;
146 extrap_mode[0] = extrap_mode[1] = 3;
147 Standard_Real eval_result[2] ;
148 Standard_Real *PolesArray =
149 (Standard_Real *) &Poles(Poles.Lower()) ;
150 //
151 // evaluate the 1D bspline that represents the change in parameterization
152 //
153 BSplCLib::Eval(*Parameter,
154 Standard_False,
155 *DerivativeRequest,
156 extrap_mode[0],
157 3,
158 FlatKnots,
159 1,
160 PolesArray[0],
161 eval_result[0]) ;
162
163
164 if (*DerivativeRequest == 0){
165 HCurve2d->D0(eval_result[0],Point);
166 Point.Coord(Result[0],Result[1]);
167 }
168 else if (*DerivativeRequest == 1){
169 HCurve2d->D1(eval_result[0], Point, Vector);
170 Vector.Multiply(eval_result[1]);
171 Vector.Coord(Result[0],Result[1]);
172 }
173 ReturnCode[0] = 0 ;
174}
175
176static Standard_Real ComputeTolReached(const Handle(Adaptor3d_HCurve)& c3d,
177 const Adaptor3d_CurveOnSurface& cons,
178 const Standard_Integer nbp)
179{
180 Standard_Real d2 = 0.;
181 Standard_Integer nn = nbp;
182 Standard_Real unsurnn = 1./nn;
183 Standard_Real first = c3d->FirstParameter();
184 Standard_Real last = c3d->LastParameter();
185 for(Standard_Integer i = 0; i <= nn; i++){
186 Standard_Real t = unsurnn*i;
187 Standard_Real u = first*(1.-t) + last*t;
188 gp_Pnt Pc3d = c3d->Value(u);
189 gp_Pnt Pcons = cons.Value(u);
190 if (Precision::IsInfinite(Pcons.X()) ||
191 Precision::IsInfinite(Pcons.Y()) ||
192 Precision::IsInfinite(Pcons.Z())) {
193 d2=Precision::Infinite();
194 break;
195 }
196 Standard_Real temp = Pc3d.SquareDistance(Pcons);
197 if(temp > d2) d2 = temp;
198 }
199 d2 = 1.5*sqrt(d2);
200 if(d2<1.e-7) d2 = 1.e-7;
201 return d2;
202}
203
204static Standard_Boolean Check(const TColStd_Array1OfReal& FlatKnots,
205 const TColStd_Array1OfReal& Poles,
206 const Standard_Integer nbp,
207 const TColStd_Array1OfReal& pc3d,
208// const TColStd_Array1OfReal& pcons,
209 const TColStd_Array1OfReal& ,
210 const Handle(Adaptor3d_HCurve)& c3d,
211 const Adaptor3d_CurveOnSurface& cons,
212 Standard_Real& tol,
213 const Standard_Real oldtol)
214{
215 Standard_Real d = tol;
216 Standard_Integer extrap_mode[2] ;
217 extrap_mode[0] = extrap_mode[1] = 3;
218 Standard_Integer i;
219#ifdef DEB
220 if (Voir) {
221 cout<<endl;
0d969553
Y		 222 cout<<"Control the change of variable : "<<endl;
223 cout<<"yawn mesured by projection : "<<d<<endl;
224 cout<<"Number of points : "<<nbp<<endl;
7fd59977 225 }
226#endif
227#if 0
228 Standard_Real glis = 0., dglis = 0.;
229 for(i = 1; i <= nbp; i++){
230 Standard_Real tc3d = pc3d(i);
231 gp_Pnt Pc3d = c3d->Value(tc3d);
232 Standard_Real tcons;
233 BSplCLib::Eval(tc3d,Standard_False,0,extrap_mode[0],
234 3,FlatKnots,1, (Standard_Real&)Poles(1),tcons);
235 gp_Pnt Pcons = cons.Value(tcons);
236 Standard_Real temp = Pc3d.SquareDistance(Pcons);
237 if(temp >= dglis) dglis = temp;
238 temp = Abs(tcons-pcons(i));
239 if(temp >= glis) glis = temp;
240 }
241 dglis = sqrt(dglis);
242#ifdef DEB
243 if ( Voir) {
0d969553
Y		 244 cout<<"shift of parameter to the imposed points : "<<glis<<endl;
245 cout<<"shift distance at the imposed points : "<<dglis<<endl;
7fd59977 246 }
247#endif
248 dglis = 0.;
249 for(i = 1; i < nbp; i++){
250 Standard_Real tc3d = 0.5*(pc3d(i)+pc3d(i+1));
251 gp_Pnt Pc3d = c3d->Value(tc3d);
252 Standard_Real tcons;
253 BSplCLib::Eval(tc3d,Standard_False,0,extrap_mode[0],
254 3,FlatKnots,1, (Standard_Real&)Poles(1),tcons);
255 gp_Pnt Pcons = cons.Value(tcons);
256 Standard_Real temp = Pc3d.SquareDistance(Pcons);
257 if(temp >= dglis) dglis = temp;
258 }
259 dglis = sqrt(dglis);
260#ifdef DEB
261 if (Voir)
262 cout<<"distance de glissement en milieu d intervals : "<<dglis<<endl;
263#endif
264#endif
265
266 Standard_Real d2 = 0.;
267 Standard_Integer nn = 2*nbp;
268 Standard_Real unsurnn = 1./nn;
269// Modified by skv - Wed Jun 2 11:49:59 2004 OCC5898 Begin
270// Correction of the interval of valid values. This condition has no sensible
271// grounds. But it is better then the old one (which is commented out) because
272// it fixes the bug OCC5898. To develop more or less sensible criterion it is
273// necessary to deeply investigate this problem which is not possible in frames
274// of debugging.
275
276// Standard_Real firstborne= 2*pc3d(1)-pc3d(nbp);
277// Standard_Real lastborne= 2*pc3d(nbp)-pc3d(1);
278 Standard_Real firstborne= 3.*pc3d(1) - 2.*pc3d(nbp);
279 Standard_Real lastborne = 3.*pc3d(nbp) - 2.*pc3d(1);
280// Modified by skv - Wed Jun 2 11:50:03 2004 OCC5898 End
281 for(i = 0; i <= nn; i++){
282 Standard_Real t = unsurnn*i;
283 Standard_Real tc3d = pc3d(1)*(1.-t) + pc3d(nbp)*t;
284 gp_Pnt Pc3d = c3d->Value(tc3d);
285 Standard_Real tcons;
286 BSplCLib::Eval(tc3d,Standard_False,0,extrap_mode[0],
287 3,FlatKnots,1, (Standard_Real&)Poles(1),tcons);
288 if (tcons < firstborne || tcons > lastborne) {
289 tol=Precision::Infinite();
290 return Standard_False;
291 }
292 gp_Pnt Pcons = cons.Value(tcons);
293 Standard_Real temp = Pc3d.SquareDistance(Pcons);
294 if(temp > d2) d2 = temp;
295 }
296 tol = sqrt(d2);
297#ifdef DEB
298 if (Voir)
0d969553 299 cout<<"distance max on "<<nn<<" points : "<<tol<<endl<<endl;
7fd59977 300#endif
301 return ((tol <= d) || (tol > 0.8 * oldtol));
302}
303
304
305//=======================================================================
306//function : Approx_SameParameter
307//purpose :
308//=======================================================================
309
310Approx_SameParameter::Approx_SameParameter(const Handle(Geom_Curve)& C3D,
311 const Handle(Geom2d_Curve)& C2D,
312 const Handle(Geom_Surface)& S,
313 const Standard_Real Tol):
314 mySameParameter(Standard_True), myDone(Standard_False)
315{
316 myHCurve2d = new Geom2dAdaptor_HCurve(C2D);
317 myC3d = new GeomAdaptor_HCurve(C3D);
318 mySurf = new GeomAdaptor_HSurface(S);
319 Build(Tol);
320}
321
322
323//=======================================================================
324//function : Approx_SameParameter
325//purpose :
326//=======================================================================
327
328Approx_SameParameter::Approx_SameParameter(const Handle(Adaptor3d_HCurve)& C3D,
329 const Handle(Geom2d_Curve)& C2D,
330 const Handle(Adaptor3d_HSurface)& S,
331 const Standard_Real Tol):
332 mySameParameter(Standard_True), myDone(Standard_False)
333{
334 myC3d = C3D;
335 mySurf = S;
336 myHCurve2d = new Geom2dAdaptor_HCurve(C2D);
337 Build(Tol);
338}
339
340
341//=======================================================================
342//function : Approx_SameParameter
343//purpose :
344//=======================================================================
345
346Approx_SameParameter::Approx_SameParameter(const Handle(Adaptor3d_HCurve)& C3D,
347 const Handle(Adaptor2d_HCurve2d)& C2D,
348 const Handle(Adaptor3d_HSurface)& S,
349 const Standard_Real Tol):
350 mySameParameter(Standard_True), myDone(Standard_False)
351{
352 myC3d = C3D;
353 mySurf = S;
354 myHCurve2d = C2D;
355 Build(Tol);
356}
357
358
359//=======================================================================
360//function : Build
361//purpose :
362//=======================================================================
363
364void Approx_SameParameter::Build(const Standard_Real Tolerance)
365{
366 Standard_Integer ii ;
367 Adaptor3d_CurveOnSurface CurveOnSurface(myHCurve2d,mySurf);
368 Standard_Real fcons = CurveOnSurface.FirstParameter();
369 Standard_Real lcons = CurveOnSurface.LastParameter();
370 Standard_Real fc3d = myC3d->FirstParameter();
371 Standard_Real lc3d = myC3d->LastParameter();
372
373 GeomAbs_Shape Continuity = myHCurve2d->Continuity();
374
375 if(Continuity > GeomAbs_C1) Continuity = GeomAbs_C1;
376
0d969553
Y		 377 //Control tangents at the extremities to know if the
378 //reparametring is possible and calculate the tangents
379 //at the extremities of the function of change of variable.
7fd59977 380 Standard_Real tangent[2];
381 gp_Pnt Pcons,Pc3d;
382 gp_Vec Vcons,Vc3d;
383
384 Standard_Real Tol = Tolerance;
385 Standard_Real Tol2 = Tol * Tol;
386 Standard_Real Tolp = myC3d->Resolution(Tol), deltamin = 50*Tolp;
387
388 Standard_Real besttol2 = Tol2;
389 Standard_Boolean extrok = 0;
390
391 extrok = 1;
392 CurveOnSurface.D1(fcons,Pcons,Vcons);
393 myC3d->D1(fc3d,Pc3d,Vc3d);
394 Standard_Real dist2 = Pcons.SquareDistance(Pc3d);
395 Standard_Real dmax2 = dist2;
396
397 Standard_Real magVcons = Vcons.Magnitude();
398 if (magVcons > 1.e-12){
399 tangent[0] = Vc3d.Magnitude() / magVcons;
400 }
401 else extrok = 0;
402
403 CurveOnSurface.D1(lcons,Pcons,Vcons);
404 myC3d->D1(lc3d,Pc3d,Vc3d);
405 dist2 = Pcons.SquareDistance(Pc3d);
406
407 if(dist2 > dmax2) dmax2 = dist2;
408 magVcons = Vcons.Magnitude();
409 if (magVcons > 1.e-12){
410 tangent[1] = Vc3d.Magnitude() / magVcons;
411 }
412 else extrok = 0;
413
414
415 if(dmax2 > besttol2) besttol2 = dmax2;
416
0d969553
Y		 417 //Take a multiple of the sample pof CheckShape,
418 //at least the control points will be correct. No comment!!!
7fd59977 419
420 Standard_Integer NCONTROL = 22;
421#ifdef DEB
422 Standard_Integer nbcoups = 0;
423#endif
424
425 Standard_Boolean interpolok = 0;
426 Standard_Real tolsov = 1.e200;
0d969553
Y		 427 //Take parameters with constant step on the curve on surface
428 //and on curve 3d.
7fd59977 429 Standard_Real deltacons = lcons - fcons;
430 deltacons /= (NCONTROL);
431 Standard_Real deltac3d = lc3d - fc3d;
432 deltac3d /= (NCONTROL);
433
434 Standard_Real wcons = fcons;
435 Standard_Real wc3d = fc3d;
436
437 Standard_Real qpcons[MAX_ARRAY_SIZE], qnewpcons[MAX_ARRAY_SIZE],
438 qpc3d[MAX_ARRAY_SIZE], qnewpc3d[MAX_ARRAY_SIZE];
439 Standard_Real * pcons = qpcons; Standard_Real * newpcons = qnewpcons;
440 Standard_Real * pc3d = qpc3d; Standard_Real * newpc3d = qnewpc3d;
441
442 for ( ii = 0 ; ii < NCONTROL; ii++) {
443 pcons[ii] = wcons;
444 pc3d[ii] = wc3d;
445 wcons += deltacons;
446 wc3d += deltac3d;
447 }
448 pcons[NCONTROL] = lcons;
449 pc3d[NCONTROL] = lc3d;
450
451 Standard_Integer New_NCONTROL = NCONTROL;
452 if(Continuity < GeomAbs_C1) {
453 Standard_Integer NbInt = myHCurve2d->NbIntervals(GeomAbs_C1) + 1;
454 TColStd_Array1OfReal Param_de_decoupeC1 (1, NbInt);
455 myHCurve2d->Intervals(Param_de_decoupeC1, GeomAbs_C1);
456 TColStd_SequenceOfReal new_par;
457 Standard_Integer inter = 1;
458 ii =1;
459 new_par.Append(fcons);
460
461 while(Param_de_decoupeC1(inter) <= fcons + deltamin) inter++;
462 while(Param_de_decoupeC1(NbInt) >= lcons - deltamin) NbInt--;
463
464 while(inter <= NbInt || ii < NCONTROL) {
465 if(Param_de_decoupeC1(inter) < pcons[ii]) {
466 new_par.Append(Param_de_decoupeC1(inter));
467 if((pcons[ii] - Param_de_decoupeC1(inter)) <= deltamin) {
468 ii++;
469 if(ii > NCONTROL) {ii = NCONTROL;}
470 }
471 inter++;
472 }
473 else {
474 if((Param_de_decoupeC1(inter) - pcons[ii]) > deltamin) {
475 new_par.Append(pcons[ii]);
476 }
477 ii++;
478 }
479 }
480
481 new_par.Append(lcons);
482 New_NCONTROL = new_par.Length() - 1;
483 //simple protection if New_NCONTROL > allocated elements in array
484 if (New_NCONTROL > MAX_ARRAY_SIZE) {
485 mySameParameter = Standard_False;
486 return;
487 }
488 for(ii = 1; ii <= New_NCONTROL; ii++){
489 pcons[ii] = pc3d[ii] = new_par.Value(ii + 1);
490 }
491 pc3d[New_NCONTROL] = lc3d;
492 }
493
494
495 Extrema_LocateExtPC Projector;
496 Projector.Initialize(myC3d->Curve(),fc3d,lc3d,Tol);
497
498 Standard_Integer count = 1;
499 Standard_Real previousp = fc3d, initp=0, curp;//, deltamin = 50*Tolp;
500 Standard_Real bornesup = lc3d - deltamin;
501 Standard_Boolean projok = 0,
502 use_parameter ;
503 for (ii = 1; ii < New_NCONTROL; ii++){
504 CurveOnSurface.D0(pcons[ii],Pcons);
505 myC3d->D0(pc3d[ii],Pc3d);
506 dist2 = Pcons.SquareDistance(Pc3d);
507 use_parameter = (dist2 <= Tol2 && (pc3d[ii] > pc3d[count-1] + deltamin)) ;
508 if(use_parameter) {
509
510 if(dist2 > dmax2) dmax2 = dist2;
511 initp = previousp = pc3d[count] = pc3d[ii];
512 pcons[count] = pcons[ii];
513 count++;
514 }
515 else {
516 if(!projok) initp = pc3d[ii];
517 projok = mySameParameter = Standard_False;
518 Projector.Perform(Pcons, initp);
519 if (Projector.IsDone()) {
520 curp = Projector.Point().Parameter();
521 Standard_Real dist_2 = Projector.SquareDistance();
522 if(dist_2 > besttol2) besttol2 = dist_2;
523 projok = 1;
524 }
525 else {
526 ProjectPointOnCurve(initp,Pcons,Tol,30,myC3d->Curve(),projok,curp);
527 }
528 if(projok){
529 if(curp > previousp + deltamin && curp < bornesup){
530 initp = previousp = pc3d[count] = curp;
531 pcons[count] = pcons[ii];
532 count++;
533 }
534 }
535 else {
536#ifdef DEB
537 // JAG
538 cout << "Projection not done" << endl;
539#endif
540 }
541 }
542 }
543 if(mySameParameter){
544 myTolReached = 1.5*sqrt(dmax2);
545 return;
546 }
547
0d969553 548 if(!extrok) { // If not already SameP and tangent to mill, abandon.
7fd59977 549 mySameParameter = Standard_False;
550#ifdef DEB
0d969553 551 cout<<"SameParameter problem : zero tangent to extremities"<<endl;
7fd59977 552#endif
553 return;
554 }
555
556 pcons[count] = lcons;
557 pc3d[count] = lc3d;
558
559#ifdef DEB
560 if (AffichFw) {
561 char Name[17];
562 Name[0]='\0';
563 TColgp_Array1OfPnt2d DEBP2d (0,count);
564 TColStd_Array1OfInteger DEBMults(0,count);
565 DEBMults.Init(1); DEBMults(0) = 2; DEBMults(count) = 2;
566 TColStd_Array1OfReal DEBKnots(0,count);
567 for (Standard_Integer DEBi = 0; DEBi <= count; DEBi++) {
568 DEBKnots(DEBi) = DEBi;
569 DEBP2d (DEBi) = gp_Pnt2d(pc3d[DEBi],pcons[DEBi]);
570 }
571 Handle(Geom2d_BSplineCurve) DEBBS =
572 new Geom2d_BSplineCurve(DEBP2d,DEBKnots,DEBMults,1);
91322f44 573 Sprintf(Name,"DEBC2d_%d",++NbCurve);
7fd59977 574#ifdef DRAW
575 DrawTrSurf::Set(Name,DEBBS);
576#endif
577 }
578#endif
579
580 while(!interpolok){
0d969553 581 // The tables and their limits for the interpolation.
7fd59977 582 Standard_Integer num_knots = count + 7;
583 Standard_Integer num_poles = count + 3;
584 TColStd_Array1OfReal Paramc3d(*pc3d,1,count+1);
585 TColStd_Array1OfReal Paramcons(*pcons,1,count+1);
586 TColStd_Array1OfInteger ContactOrder(1,num_poles) ;
587 TColStd_Array1OfReal Poles(1,num_poles) ;
588 TColStd_Array1OfReal InterpolationParameters(1,num_poles) ;
589 TColStd_Array1OfReal FlatKnots(1,num_knots) ;
590
0d969553 591 // Fill tables taking attention to end values.
7fd59977 592 ContactOrder.Init(0);
593 ContactOrder(2) = ContactOrder(num_poles - 1) = 1;
594
595 FlatKnots(1) = FlatKnots(2) = FlatKnots(3) = FlatKnots(4) = fc3d;
596 FlatKnots(num_poles + 1) = FlatKnots(num_poles + 2) =
597 FlatKnots(num_poles + 3) = FlatKnots(num_poles + 4) = lc3d;
598
599 Poles(1) = fcons; Poles(num_poles) = lcons;
600 Poles(2) = tangent[0]; Poles(num_poles - 1) = tangent[1];
601
602 InterpolationParameters(1) = InterpolationParameters(2) = fc3d;
603 InterpolationParameters(num_poles - 1) = InterpolationParameters(num_poles) = lc3d;
604
605 for (ii = 3; ii <= num_poles - 2; ii++) {
606 Poles(ii) = Paramcons(ii - 1);
607 InterpolationParameters(ii) = FlatKnots(ii+2) = Paramc3d(ii - 1);
608 }
609 Standard_Integer inversion_problem;
610 BSplCLib::Interpolate(3,FlatKnots,InterpolationParameters,ContactOrder,
611 1,Poles(1),inversion_problem);
612 if(inversion_problem) {
613 Standard_ConstructionError::Raise();
614 }
615
616 //-------------------------------------------
617#ifdef DEB
618 if (AffichFw) {
619 nbcoups ++;
620 char Name[17];
621 Name[0] = '\0';
622 Standard_Integer nnn = 100;
623 TColgp_Array1OfPnt2d DEBP2d (0,nnn);
624 TColStd_Array1OfInteger DEBMults(0,nnn);
625 DEBMults.Init(1); DEBMults(0) = 2; DEBMults(nnn) = 2;
626 TColStd_Array1OfReal DEBKnots(0,nnn);
627 Standard_Real du = (lc3d - fc3d) / nnn;
628 Standard_Real u3d = fc3d;
629 Standard_Integer extrap_mode[2] ;
630 extrap_mode[0] = extrap_mode[1] = 3;
631 Standard_Real eval_result[2] ;
632 Standard_Integer DerivativeRequest = 0;
633 Standard_Real *PolesArray =
634 (Standard_Real *) &Poles(Poles.Lower()) ;
635
636 for (Standard_Integer DEBi = 0; DEBi <= nnn; DEBi++) {
637 DEBKnots(DEBi) = DEBi;
638 BSplCLib::Eval(u3d,
639 Standard_False,
640 DerivativeRequest,
641 extrap_mode[0],
642 3,
643 FlatKnots,
644 1,
645 PolesArray[0],
646 eval_result[0]) ;
647
648 DEBP2d (DEBi) = gp_Pnt2d(u3d,eval_result[0]);
649 u3d += du;
650 }
651
652 Handle(Geom2d_BSplineCurve) DEBBS =
653 new Geom2d_BSplineCurve(DEBP2d,DEBKnots,DEBMults,1);
91322f44 654 Sprintf(Name,"DEBC2d_%d_%d",NbCurve,nbcoups );
7fd59977 655#ifdef DRAW
656 DrawTrSurf::Set(Name,DEBBS);
657#endif
658 }
659#endif
660//-------------------------------------------
661
662//-------------------------------------------
663 // Test if par2d(par3d) is monotonous function or not ----- IFV, Jan 2000
664 // and try to insert new point to improve BSpline interpolation
665
666 Standard_Integer extrap_mode[2] ;
667 extrap_mode[0] = extrap_mode[1] = 3;
668 Standard_Real eval_result[2] ;
669 Standard_Integer DerivativeRequest = 0;
670 Standard_Real *PolesArray =
671 (Standard_Real *) &Poles(Poles.Lower()) ;
672
673 Standard_Integer newcount = 0;
674 for (ii = 0; ii < count; ii++) {
675
676 newpcons[newcount] = pcons[ii];
677 newpc3d[newcount] = pc3d[ii];
678 newcount++;
679
680 if(count - ii + newcount == MAX_ARRAY_SIZE) continue;
681
682 BSplCLib::Eval(.5*(pc3d[ii]+pc3d[ii+1]), Standard_False, DerivativeRequest,
683 extrap_mode[0], 3, FlatKnots, 1, PolesArray[0], eval_result[0]);
684
685 if(eval_result[0] < pcons[ii] || eval_result[0] > pcons[ii+1]) {
686 Standard_Real ucons = 0.5*(pcons[ii]+pcons[ii+1]);
687 Standard_Real uc3d = 0.5*(pc3d[ii]+pc3d[ii+1]);
688
689 CurveOnSurface.D0(ucons,Pcons);
690 Projector.Perform(Pcons, uc3d);
691 if (Projector.IsDone()) {
692 curp = Projector.Point().Parameter();
693 Standard_Real dist_2 = Projector.SquareDistance();
694 if(dist_2 > besttol2) besttol2 = dist_2;
695 projok = 1;
696 }
697 else {
698 ProjectPointOnCurve(uc3d,Pcons,Tol,30,myC3d->Curve(),projok,curp);
699 }
700 if(projok){
701 if(curp > pc3d[ii] + deltamin && curp < pc3d[ii+1] - deltamin){
702 newpc3d[newcount] = curp;
703 newpcons[newcount] = ucons;
704 newcount ++;
705 }
706 }
707 else {
708#ifdef DEB
709 // JAG
710 cout << "Projection not done" << endl;
711#endif
712 }
713 }
714
715 }
716
717 newpc3d[newcount] = pc3d[count];
718 newpcons[newcount] = pcons[count];
719 Standard_Real * temp;
720 temp = pc3d;
721 pc3d = newpc3d;
722 newpc3d = temp;
723 temp = pcons;
724 pcons = newpcons;
725 newpcons = temp;
726
727 if((count != newcount) && newcount < MAX_ARRAY_SIZE) { count = newcount; continue;}
728
729 count = newcount;
730
731 Standard_Real algtol = sqrt(besttol2);
732
733 interpolok = Check (FlatKnots, Poles, count+1, Paramc3d, Paramcons,
734 myC3d, CurveOnSurface, algtol, tolsov);
735
736 if (Precision::IsInfinite(algtol)) {
737 mySameParameter = Standard_False;
738#ifdef DEB
0d969553 739 cout<<"SameParameter problem : function of interpolation of parametration at mills !!"<<endl;
7fd59977 740#endif
741 return;
742 }
743
744 tolsov = algtol;
745
746 interpolok = (interpolok || count >= MAX_ARRAY_SIZE);
747
748 if(interpolok) {
749 Standard_Real besttol = sqrt(besttol2);
750#ifdef DEB
751 if (Voir) {
752 if(algtol > besttol){
0d969553 753 cout<<"SameParameter : Tol can't be reached before approx"<<endl;
7fd59977 754 }
755 }
756#endif
757 Handle(TColStd_HArray1OfReal) tol1d,tol2d,tol3d;
758 tol1d = new TColStd_HArray1OfReal(1,2) ;
759 tol1d->SetValue(1, mySurf->UResolution(besttol));
760 tol1d->SetValue(2, mySurf->VResolution(besttol));
761
762 Approx_SameParameter_Evaluator ev (FlatKnots, Poles, myHCurve2d);
763 AdvApprox_ApproxAFunction anApproximator(2,0,0,tol1d,tol2d,tol3d,fc3d,lc3d,
764 Continuity,11,40,ev);
765
766 if (anApproximator.IsDone() || anApproximator.HasResult()) {
767 GeomLib_MakeCurvefromApprox aCurveBuilder(anApproximator) ;
768 myCurve2d = aCurveBuilder.Curve2dFromTwo1d(1,2) ;
769 myHCurve2d = new Geom2dAdaptor_HCurve(myCurve2d);
770 CurveOnSurface.Load(myHCurve2d);
771 myTolReached = ComputeTolReached(myC3d,CurveOnSurface,NCONTROL);
772 myDone = Standard_True;
773 }
774 }
775 else {
776#ifdef DEB
777 if (Voir)
0d969553 778 cout<<"SameParameter : Not enough points, enrich"<<endl;
7fd59977 779#endif
780
781 Standard_Integer newcount = 0;
782 for(Standard_Integer n = 0; n < count; n++){
783 newpc3d[newcount] = pc3d[n];
784 newpcons[newcount] = pcons[n];
785 newcount ++;
786
787 if(count - n + newcount == MAX_ARRAY_SIZE) continue;
788
789 Standard_Real ucons = 0.5*(pcons[n]+pcons[n+1]);
790 Standard_Real uc3d = 0.5*(pc3d[n]+pc3d[n+1]);
791
792 CurveOnSurface.D0(ucons,Pcons);
793 Projector.Perform(Pcons, uc3d);
794 if (Projector.IsDone()) {
795 curp = Projector.Point().Parameter();
796 Standard_Real dist_2 = Projector.SquareDistance();
797 if(dist_2 > besttol2) besttol2 = dist_2;
798 projok = 1;
799 }
800 else {
801 ProjectPointOnCurve(uc3d,Pcons,Tol,30,myC3d->Curve(),projok,curp);
802 }
803 if(projok){
804 if(curp > pc3d[n] + deltamin && curp < pc3d[n+1] - deltamin){
805 newpc3d[newcount] = curp;
806 newpcons[newcount] = ucons;
807 newcount ++;
808 }
809 }
810 else {
811#ifdef DEB
812 // JAG
813 cout << "Projection not done" << endl;
814#endif
815 }
816 }
817 newpc3d[newcount] = pc3d[count];
818 newpcons[newcount] = pcons[count];
819 Standard_Real * tempx;
820 tempx = pc3d;
821 pc3d = newpc3d;
822 newpc3d = tempx;
823 tempx = pcons;
824 pcons = newpcons;
825 newpcons = tempx;
826 count = newcount;
827 }
828 }
829}
Open CASCADE Technology repository