Commit | Line | Data |
---|---|---|
b311480e | 1 | // Created on: 1995-06-06 |
2 | // Created by: Xavier BENVENISTE | |
3 | // Copyright (c) 1995-1999 Matra Datavision | |
973c2be1 | 4 | // Copyright (c) 1999-2014 OPEN CASCADE SAS |
b311480e | 5 | // |
973c2be1 | 6 | // This file is part of Open CASCADE Technology software library. |
b311480e | 7 | // |
d5f74e42 | 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 | |
973c2be1 | 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. | |
b311480e | 13 | // |
973c2be1 | 14 | // Alternatively, this file may be used under the terms of Open CASCADE |
15 | // commercial license or contractual agreement. | |
7fd59977 | 16 | |
17 | // Modified by skv - Wed Jun 2 11:49:59 2004 OCC5898 | |
18 | ||
42cf5bc1 | 19 | #include <Adaptor2d_HCurve2d.hxx> |
7fd59977 | 20 | #include <Adaptor3d_CurveOnSurface.hxx> |
42cf5bc1 | 21 | #include <Adaptor3d_HCurve.hxx> |
22 | #include <Adaptor3d_HSurface.hxx> | |
23 | #include <AdvApprox_ApproxAFunction.hxx> | |
24 | #include <Approx_SameParameter.hxx> | |
25 | #include <BSplCLib.hxx> | |
26 | #include <Extrema_ExtPC.hxx> | |
27 | #include <Extrema_LocateExtPC.hxx> | |
28 | #include <GCPnts_QuasiUniformDeflection.hxx> | |
29 | #include <Geom2d_BSplineCurve.hxx> | |
30 | #include <Geom2d_Curve.hxx> | |
7fd59977 | 31 | #include <Geom2dAdaptor_Curve.hxx> |
32 | #include <Geom2dAdaptor_HCurve.hxx> | |
42cf5bc1 | 33 | #include <Geom_Curve.hxx> |
34 | #include <Geom_Surface.hxx> | |
7fd59977 | 35 | #include <GeomAdaptor_Curve.hxx> |
36 | #include <GeomAdaptor_HCurve.hxx> | |
7fd59977 | 37 | #include <GeomAdaptor_HSurface.hxx> |
42cf5bc1 | 38 | #include <GeomAdaptor_Surface.hxx> |
7fd59977 | 39 | #include <GeomLib_MakeCurvefromApprox.hxx> |
40 | #include <Precision.hxx> | |
42cf5bc1 | 41 | #include <Standard_ConstructionError.hxx> |
42 | #include <Standard_OutOfRange.hxx> | |
43 | #include <TColStd_Array1OfReal.hxx> | |
7fd59977 | 44 | |
7fd59977 | 45 | //======================================================================= |
52db4751 | 46 | //class : Approx_SameParameter_Evaluator |
47 | //purpose : Used in same parameterization curve approximation. | |
7fd59977 | 48 | //======================================================================= |
7fd59977 | 49 | class Approx_SameParameter_Evaluator : public AdvApprox_EvaluatorFunction |
50 | { | |
a86d3ec0 | 51 | public: |
52db4751 | 52 | Approx_SameParameter_Evaluator (const TColStd_Array1OfReal& theFlatKnots, |
53 | const TColStd_Array1OfReal& thePoles, | |
54 | const Handle(Adaptor2d_HCurve2d)& theHCurve2d) | |
55 | : FlatKnots(theFlatKnots), | |
56 | Poles(thePoles), | |
57 | HCurve2d(theHCurve2d) {} | |
7fd59977 | 58 | |
59 | virtual void Evaluate (Standard_Integer *Dimension, | |
52db4751 | 60 | Standard_Real StartEnd[2], |
61 | Standard_Real *Parameter, | |
62 | Standard_Integer *DerivativeRequest, | |
63 | Standard_Real *Result, // [Dimension] | |
64 | Standard_Integer *ErrorCode); | |
a86d3ec0 | 65 | |
66 | private: | |
7fd59977 | 67 | const TColStd_Array1OfReal& FlatKnots; |
68 | const TColStd_Array1OfReal& Poles; | |
69 | Handle(Adaptor2d_HCurve2d) HCurve2d; | |
70 | }; | |
71 | ||
52db4751 | 72 | //======================================================================= |
73 | //function : Evaluate | |
74 | //purpose : | |
75 | //======================================================================= | |
7fd59977 | 76 | void Approx_SameParameter_Evaluator::Evaluate (Standard_Integer *,/*Dimension*/ |
52db4751 | 77 | Standard_Real /*StartEnd*/[2], |
78 | Standard_Real *Parameter, | |
79 | Standard_Integer *DerivativeRequest, | |
80 | Standard_Real *Result, | |
81 | Standard_Integer *ReturnCode) | |
82 | { | |
83 | const Standard_Integer aDegree = 3; | |
84 | Standard_Integer extrap_mode[2] = {aDegree, aDegree}; | |
85 | Standard_Real eval_result[2]; | |
86 | Standard_Real *PolesArray = (Standard_Real *) &Poles(Poles.Lower()) ; | |
87 | ||
88 | // Evaluate the 1D B-Spline that represents the change in parameterization. | |
7fd59977 | 89 | BSplCLib::Eval(*Parameter, |
52db4751 | 90 | Standard_False, |
91 | *DerivativeRequest, | |
92 | extrap_mode[0], | |
93 | aDegree, | |
94 | FlatKnots, | |
95 | 1, | |
96 | PolesArray[0], | |
97 | eval_result[0]); | |
98 | ||
99 | gp_Pnt2d aPoint; | |
100 | gp_Vec2d aVector; | |
101 | if (*DerivativeRequest == 0) | |
102 | { | |
103 | HCurve2d->D0(eval_result[0], aPoint); | |
104 | aPoint.Coord(Result[0],Result[1]); | |
7fd59977 | 105 | } |
52db4751 | 106 | else if (*DerivativeRequest == 1) |
107 | { | |
108 | HCurve2d->D1(eval_result[0], aPoint, aVector); | |
109 | aVector.Multiply(eval_result[1]); | |
110 | aVector.Coord(Result[0],Result[1]); | |
7fd59977 | 111 | } |
52db4751 | 112 | |
113 | ReturnCode[0] = 0; | |
114 | } | |
115 | ||
116 | //======================================================================= | |
117 | //function : ProjectPointOnCurve | |
118 | //purpose : | |
119 | //======================================================================= | |
120 | static void ProjectPointOnCurve(const Standard_Real InitValue, | |
121 | const gp_Pnt APoint, | |
122 | const Standard_Real Tolerance, | |
123 | const Standard_Integer NumIteration, | |
124 | const Adaptor3d_Curve& Curve, | |
125 | Standard_Boolean& Status, | |
126 | Standard_Real& Result) | |
127 | { | |
128 | Standard_Integer num_iter = 0, not_done = 1, ii; | |
129 | ||
130 | gp_Pnt a_point; | |
131 | gp_Vec vector, d1, d2; | |
132 | Standard_Real func, func_derivative, | |
133 | param = InitValue; | |
134 | Status = Standard_False; | |
135 | do | |
136 | { | |
137 | num_iter++; | |
138 | Curve.D2(param, a_point, d1, d2); | |
139 | for (ii = 1 ; ii <= 3 ; ii++) | |
140 | vector.SetCoord(ii, APoint.Coord(ii) - a_point.Coord(ii)); | |
141 | ||
142 | func = vector.Dot(d1); | |
143 | if ( Abs(func) < Tolerance * d1.Magnitude()) | |
144 | { | |
145 | not_done = 0; | |
146 | Status = Standard_True; | |
147 | } | |
148 | else | |
149 | { | |
150 | func_derivative = vector.Dot(d2) - d1.Dot(d1); | |
151 | ||
152 | // Avoid division by zero. | |
153 | const Standard_Real Toler = 1.0e-12; | |
154 | if( Abs(func_derivative) > Toler ) | |
155 | param -= func / func_derivative; | |
156 | ||
157 | param = Max(param,Curve.FirstParameter()); | |
158 | param = Min(param,Curve.LastParameter()); | |
159 | } | |
160 | } while (not_done && num_iter <= NumIteration); | |
161 | ||
162 | Result = param; | |
7fd59977 | 163 | } |
164 | ||
52db4751 | 165 | //======================================================================= |
166 | //function : ComputeTolReached | |
167 | //purpose : | |
168 | //======================================================================= | |
7fd59977 | 169 | static Standard_Real ComputeTolReached(const Handle(Adaptor3d_HCurve)& c3d, |
52db4751 | 170 | const Adaptor3d_CurveOnSurface& cons, |
171 | const Standard_Integer nbp) | |
7fd59977 | 172 | { |
52db4751 | 173 | Standard_Real d2 = 0.0; // Square max discrete deviation. |
a86d3ec0 | 174 | const Standard_Real first = c3d->FirstParameter(); |
175 | const Standard_Real last = c3d->LastParameter(); | |
52db4751 | 176 | for(Standard_Integer i = 0; i <= nbp; i++) |
177 | { | |
178 | Standard_Real t = IntToReal(i) / IntToReal(nbp); | |
179 | Standard_Real u = first * (1.0 - t) + last * t; | |
7fd59977 | 180 | gp_Pnt Pc3d = c3d->Value(u); |
181 | gp_Pnt Pcons = cons.Value(u); | |
182 | if (Precision::IsInfinite(Pcons.X()) || | |
52db4751 | 183 | Precision::IsInfinite(Pcons.Y()) || |
184 | Precision::IsInfinite(Pcons.Z())) | |
185 | { | |
a86d3ec0 | 186 | d2=Precision::Infinite(); |
187 | break; | |
7fd59977 | 188 | } |
52db4751 | 189 | d2 = Max(d2, Pc3d.SquareDistance(Pcons)); |
7fd59977 | 190 | } |
52db4751 | 191 | |
192 | const Standard_Real aMult = 1.5; // To be tolerant to discrete tolerance computing. | |
193 | Standard_Real aDeviation = aMult * sqrt(d2); | |
194 | aDeviation = Max(aDeviation, Precision::Confusion()); // Tolerance in modeling space. | |
195 | return aDeviation; | |
7fd59977 | 196 | } |
197 | ||
52db4751 | 198 | //======================================================================= |
199 | //function : Check | |
200 | //purpose : Check current interpolation for validity. | |
201 | //======================================================================= | |
7fd59977 | 202 | static Standard_Boolean Check(const TColStd_Array1OfReal& FlatKnots, |
52db4751 | 203 | const TColStd_Array1OfReal& Poles, |
204 | const Standard_Integer nbp, | |
205 | const TColStd_Array1OfReal& pc3d, | |
206 | const TColStd_Array1OfReal& , | |
207 | const Handle(Adaptor3d_HCurve)& c3d, | |
208 | const Adaptor3d_CurveOnSurface& cons, | |
209 | Standard_Real& tol, | |
210 | const Standard_Real oldtol) | |
7fd59977 | 211 | { |
52db4751 | 212 | const Standard_Integer aDegree = 3; |
213 | Standard_Integer extrap_mode[2] = {aDegree, aDegree}; | |
7fd59977 | 214 | |
a86d3ec0 | 215 | // Correction of the interval of valid values. This condition has no sensible |
216 | // grounds. But it is better then the old one (which is commented out) because | |
217 | // it fixes the bug OCC5898. To develop more or less sensible criterion it is | |
218 | // necessary to deeply investigate this problem which is not possible in frames | |
219 | // of debugging. | |
52db4751 | 220 | Standard_Real aParamFirst = 3.0 * pc3d(1) - 2.0 * pc3d(nbp); |
221 | Standard_Real aParamLast = 3.0 * pc3d(nbp) - 2.0 * pc3d(1); | |
a86d3ec0 | 222 | |
4590b551 | 223 | Standard_Real FirstPar = cons.FirstParameter(); |
224 | Standard_Real LastPar = cons.LastParameter(); | |
52db4751 | 225 | if (aParamFirst < FirstPar) |
226 | aParamFirst = FirstPar; | |
227 | if (aParamLast > LastPar) | |
228 | aParamLast = LastPar; | |
229 | ||
230 | ||
231 | Standard_Real d2 = 0.0; // Maximum square deviation on the samples. | |
232 | const Standard_Real d = tol; | |
233 | const Standard_Integer nn = 2 * nbp; | |
234 | const Standard_Real unsurnn = 1.0/nn; | |
235 | for(Standard_Integer i = 0; i <= nn; i++) | |
236 | { | |
237 | // Compute corresponding parameter on 2d curve. | |
238 | // It should be inside of 3d curve parameter space. | |
7fd59977 | 239 | Standard_Real t = unsurnn*i; |
240 | Standard_Real tc3d = pc3d(1)*(1.-t) + pc3d(nbp)*t; | |
241 | gp_Pnt Pc3d = c3d->Value(tc3d); | |
242 | Standard_Real tcons; | |
243 | BSplCLib::Eval(tc3d,Standard_False,0,extrap_mode[0], | |
52db4751 | 244 | aDegree,FlatKnots,1, (Standard_Real&)Poles(1),tcons); |
245 | if (tcons < aParamFirst || | |
246 | tcons > aParamLast) | |
247 | { | |
248 | tol = Precision::Infinite(); | |
7fd59977 | 249 | return Standard_False; |
250 | } | |
251 | gp_Pnt Pcons = cons.Value(tcons); | |
252 | Standard_Real temp = Pc3d.SquareDistance(Pcons); | |
253 | if(temp > d2) d2 = temp; | |
254 | } | |
255 | tol = sqrt(d2); | |
7fd59977 | 256 | |
52db4751 | 257 | // Check poles parameters to be ordered. |
258 | for(Standard_Integer i = Poles.Lower() + 1; i <= Poles.Upper(); ++i) | |
259 | { | |
260 | const Standard_Real aPreviousParam = Poles(i - 1); | |
261 | const Standard_Real aCurrentParam = Poles(i); | |
262 | ||
263 | if (aPreviousParam > aCurrentParam) | |
264 | return Standard_False; | |
265 | } | |
266 | ||
267 | return (tol <= d || tol > 0.8 * oldtol); | |
268 | } | |
7fd59977 | 269 | |
270 | //======================================================================= | |
271 | //function : Approx_SameParameter | |
272 | //purpose : | |
273 | //======================================================================= | |
7fd59977 | 274 | Approx_SameParameter::Approx_SameParameter(const Handle(Geom_Curve)& C3D, |
52db4751 | 275 | const Handle(Geom2d_Curve)& C2D, |
276 | const Handle(Geom_Surface)& S, | |
277 | const Standard_Real Tol) | |
278 | : mySameParameter(Standard_True), | |
279 | myDone(Standard_False) | |
7fd59977 | 280 | { |
281 | myHCurve2d = new Geom2dAdaptor_HCurve(C2D); | |
282 | myC3d = new GeomAdaptor_HCurve(C3D); | |
283 | mySurf = new GeomAdaptor_HSurface(S); | |
284 | Build(Tol); | |
285 | } | |
286 | ||
7fd59977 | 287 | //======================================================================= |
288 | //function : Approx_SameParameter | |
289 | //purpose : | |
290 | //======================================================================= | |
7fd59977 | 291 | Approx_SameParameter::Approx_SameParameter(const Handle(Adaptor3d_HCurve)& C3D, |
52db4751 | 292 | const Handle(Geom2d_Curve)& C2D, |
293 | const Handle(Adaptor3d_HSurface)& S, | |
294 | const Standard_Real Tol) | |
295 | : mySameParameter(Standard_True), | |
296 | myDone(Standard_False) | |
7fd59977 | 297 | { |
298 | myC3d = C3D; | |
299 | mySurf = S; | |
300 | myHCurve2d = new Geom2dAdaptor_HCurve(C2D); | |
301 | Build(Tol); | |
302 | } | |
303 | ||
7fd59977 | 304 | //======================================================================= |
305 | //function : Approx_SameParameter | |
306 | //purpose : | |
307 | //======================================================================= | |
7fd59977 | 308 | Approx_SameParameter::Approx_SameParameter(const Handle(Adaptor3d_HCurve)& C3D, |
52db4751 | 309 | const Handle(Adaptor2d_HCurve2d)& C2D, |
310 | const Handle(Adaptor3d_HSurface)& S, | |
311 | const Standard_Real Tol) | |
312 | : mySameParameter(Standard_True), | |
313 | myDone(Standard_False) | |
7fd59977 | 314 | { |
315 | myC3d = C3D; | |
316 | mySurf = S; | |
317 | myHCurve2d = C2D; | |
318 | Build(Tol); | |
319 | } | |
320 | ||
7fd59977 | 321 | //======================================================================= |
322 | //function : Build | |
323 | //purpose : | |
324 | //======================================================================= | |
7fd59977 | 325 | void Approx_SameParameter::Build(const Standard_Real Tolerance) |
326 | { | |
a86d3ec0 | 327 | const Standard_Real anErrorMAX = 1.0e15; |
328 | const Standard_Integer aMaxArraySize = 1000; | |
329 | const Standard_Integer NCONTROL = 22; | |
330 | ||
7fd59977 | 331 | Standard_Integer ii ; |
332 | Adaptor3d_CurveOnSurface CurveOnSurface(myHCurve2d,mySurf); | |
333 | Standard_Real fcons = CurveOnSurface.FirstParameter(); | |
334 | Standard_Real lcons = CurveOnSurface.LastParameter(); | |
335 | Standard_Real fc3d = myC3d->FirstParameter(); | |
336 | Standard_Real lc3d = myC3d->LastParameter(); | |
337 | ||
0d969553 Y |
338 | //Control tangents at the extremities to know if the |
339 | //reparametring is possible and calculate the tangents | |
340 | //at the extremities of the function of change of variable. | |
d20d815b | 341 | Standard_Real tangent[2] = { 0.0, 0.0 }; |
7fd59977 | 342 | gp_Pnt Pcons,Pc3d; |
343 | gp_Vec Vcons,Vc3d; | |
344 | ||
52db4751 | 345 | const Standard_Real Tol = Tolerance; |
346 | const Standard_Real Tol2 = Tol * Tol; | |
347 | Standard_Real deltamin = Precision::PConfusion(); | |
7fd59977 | 348 | |
349 | Standard_Real besttol2 = Tol2; | |
7fd59977 | 350 | |
52db4751 | 351 | // Check tangency on curve border. |
352 | Standard_Boolean extrok = 1; | |
7fd59977 | 353 | CurveOnSurface.D1(fcons,Pcons,Vcons); |
354 | myC3d->D1(fc3d,Pc3d,Vc3d); | |
355 | Standard_Real dist2 = Pcons.SquareDistance(Pc3d); | |
356 | Standard_Real dmax2 = dist2; | |
357 | ||
358 | Standard_Real magVcons = Vcons.Magnitude(); | |
52db4751 | 359 | if (magVcons > 1.e-12) |
7fd59977 | 360 | tangent[0] = Vc3d.Magnitude() / magVcons; |
7fd59977 | 361 | else extrok = 0; |
362 | ||
363 | CurveOnSurface.D1(lcons,Pcons,Vcons); | |
364 | myC3d->D1(lc3d,Pc3d,Vc3d); | |
365 | dist2 = Pcons.SquareDistance(Pc3d); | |
366 | ||
52db4751 | 367 | dmax2 = Max(dmax2, dist2); |
7fd59977 | 368 | magVcons = Vcons.Magnitude(); |
52db4751 | 369 | if (magVcons > 1.e-12) |
7fd59977 | 370 | tangent[1] = Vc3d.Magnitude() / magVcons; |
7fd59977 | 371 | else extrok = 0; |
372 | ||
373 | ||
0d969553 Y |
374 | //Take a multiple of the sample pof CheckShape, |
375 | //at least the control points will be correct. No comment!!! | |
7fd59977 | 376 | |
7fd59977 | 377 | Standard_Boolean interpolok = 0; |
378 | Standard_Real tolsov = 1.e200; | |
0d969553 Y |
379 | //Take parameters with constant step on the curve on surface |
380 | //and on curve 3d. | |
7fd59977 | 381 | Standard_Real deltacons = lcons - fcons; |
382 | deltacons /= (NCONTROL); | |
383 | Standard_Real deltac3d = lc3d - fc3d; | |
384 | deltac3d /= (NCONTROL); | |
385 | ||
386 | Standard_Real wcons = fcons; | |
387 | Standard_Real wc3d = fc3d; | |
a86d3ec0 | 388 | |
389 | Standard_Real qpcons[aMaxArraySize], qnewpcons[aMaxArraySize], | |
52db4751 | 390 | qpc3d[aMaxArraySize], qnewpc3d[aMaxArraySize]; |
7fd59977 | 391 | Standard_Real * pcons = qpcons; Standard_Real * newpcons = qnewpcons; |
392 | Standard_Real * pc3d = qpc3d; Standard_Real * newpc3d = qnewpc3d; | |
393 | ||
394 | for ( ii = 0 ; ii < NCONTROL; ii++) { | |
395 | pcons[ii] = wcons; | |
396 | pc3d[ii] = wc3d; | |
397 | wcons += deltacons; | |
398 | wc3d += deltac3d; | |
399 | } | |
400 | pcons[NCONTROL] = lcons; | |
401 | pc3d[NCONTROL] = lc3d; | |
402 | ||
52db4751 | 403 | // Change number of points in case of C0 continuity. |
7fd59977 | 404 | Standard_Integer New_NCONTROL = NCONTROL; |
52db4751 | 405 | GeomAbs_Shape Continuity = myHCurve2d->Continuity(); |
406 | if(Continuity > GeomAbs_C1) Continuity = GeomAbs_C1; | |
407 | if(Continuity < GeomAbs_C1) | |
408 | { | |
a86d3ec0 | 409 | Standard_Integer NbInt = myHCurve2d->NbIntervals(GeomAbs_C1) + 1; |
410 | TColStd_Array1OfReal Param_de_decoupeC1 (1, NbInt); | |
411 | myHCurve2d->Intervals(Param_de_decoupeC1, GeomAbs_C1); | |
412 | TColStd_SequenceOfReal new_par; | |
413 | Standard_Integer inter = 1; | |
414 | ii =1; | |
415 | new_par.Append(fcons); | |
416 | ||
3f5bebe8 | 417 | while(inter <= NbInt && Param_de_decoupeC1(inter) <= fcons + deltamin) inter++; |
418 | while(NbInt > 0 && Param_de_decoupeC1(NbInt) >= lcons - deltamin) NbInt--; | |
a86d3ec0 | 419 | |
3f5bebe8 | 420 | while(inter <= NbInt || (ii < NCONTROL && inter <= Param_de_decoupeC1.Length()) ) { |
a86d3ec0 | 421 | if(Param_de_decoupeC1(inter) < pcons[ii]) { |
422 | new_par.Append(Param_de_decoupeC1(inter)); | |
423 | if((pcons[ii] - Param_de_decoupeC1(inter)) <= deltamin) { | |
424 | ii++; | |
425 | if(ii > NCONTROL) {ii = NCONTROL;} | |
426 | } | |
427 | inter++; | |
428 | } | |
429 | else { | |
430 | if((Param_de_decoupeC1(inter) - pcons[ii]) > deltamin) { | |
431 | new_par.Append(pcons[ii]); | |
432 | } | |
433 | ii++; | |
434 | } | |
435 | } | |
436 | ||
437 | new_par.Append(lcons); | |
438 | New_NCONTROL = new_par.Length() - 1; | |
31e0b8e8 | 439 | // Simple protection if New_NCONTROL > allocated elements in array but one |
440 | // aMaxArraySize - 1 index may be filled after projection. | |
441 | if (New_NCONTROL > aMaxArraySize - 1) { | |
a86d3ec0 | 442 | mySameParameter = Standard_False; |
443 | return; | |
444 | } | |
445 | for(ii = 1; ii <= New_NCONTROL; ii++){ | |
446 | pcons[ii] = pc3d[ii] = new_par.Value(ii + 1); | |
447 | } | |
448 | pc3d[New_NCONTROL] = lc3d; | |
449 | } | |
450 | ||
52db4751 | 451 | // Check existing same parameter state. |
7fd59977 | 452 | Extrema_LocateExtPC Projector; |
453 | Projector.Initialize(myC3d->Curve(),fc3d,lc3d,Tol); | |
a86d3ec0 | 454 | |
7fd59977 | 455 | Standard_Integer count = 1; |
52db4751 | 456 | Standard_Real previousp = fc3d, initp=0, curp; |
7fd59977 | 457 | Standard_Real bornesup = lc3d - deltamin; |
458 | Standard_Boolean projok = 0, | |
459 | use_parameter ; | |
460 | for (ii = 1; ii < New_NCONTROL; ii++){ | |
461 | CurveOnSurface.D0(pcons[ii],Pcons); | |
462 | myC3d->D0(pc3d[ii],Pc3d); | |
463 | dist2 = Pcons.SquareDistance(Pc3d); | |
464 | use_parameter = (dist2 <= Tol2 && (pc3d[ii] > pc3d[count-1] + deltamin)) ; | |
2a739b6d | 465 | Standard_Real aDistMin = RealLast();; |
7fd59977 | 466 | if(use_parameter) { |
a86d3ec0 | 467 | |
7fd59977 | 468 | if(dist2 > dmax2) dmax2 = dist2; |
469 | initp = previousp = pc3d[count] = pc3d[ii]; | |
470 | pcons[count] = pcons[ii]; | |
471 | count++; | |
2a739b6d | 472 | |
7fd59977 | 473 | } |
474 | else { | |
475 | if(!projok) initp = pc3d[ii]; | |
476 | projok = mySameParameter = Standard_False; | |
477 | Projector.Perform(Pcons, initp); | |
478 | if (Projector.IsDone()) { | |
a86d3ec0 | 479 | curp = Projector.Point().Parameter(); |
480 | Standard_Real dist_2 = Projector.SquareDistance(); | |
2a739b6d | 481 | projok = Standard_True; |
482 | aDistMin = dist_2; | |
7fd59977 | 483 | } |
a86d3ec0 | 484 | else |
485 | { | |
486 | ProjectPointOnCurve(initp,Pcons,Tol,30,myC3d->Curve(),projok,curp); | |
2a739b6d | 487 | if(projok) |
488 | { | |
489 | const gp_Pnt& ap1 =myC3d->Value(curp); | |
490 | aDistMin = Pcons.SquareDistance(ap1); | |
491 | } | |
7fd59977 | 492 | } |
2a739b6d | 493 | projok = (projok && (curp > previousp + deltamin && curp < bornesup)); |
a86d3ec0 | 494 | if(projok) |
495 | { | |
2a739b6d | 496 | initp = previousp = pc3d[count] = curp; |
497 | pcons[count] = pcons[ii]; | |
498 | count++; | |
499 | ||
7fd59977 | 500 | } |
a86d3ec0 | 501 | else |
502 | { | |
503 | Extrema_ExtPC PR(Pcons,myC3d->Curve(),fc3d,lc3d,Tol); | |
504 | if(PR.IsDone()) | |
505 | { | |
506 | const Standard_Integer aNbExt = PR.NbExt(); | |
507 | if(aNbExt > 0) | |
508 | { | |
509 | Standard_Integer anIndMin = 0; | |
2a739b6d | 510 | Standard_Real aCurDistMin = RealLast(); |
a86d3ec0 | 511 | for(Standard_Integer i = 1; i <= aNbExt; i++) |
512 | { | |
513 | const gp_Pnt &aP = PR.Point(i).Value(); | |
514 | Standard_Real aDist2 = aP.SquareDistance(Pcons); | |
2a739b6d | 515 | if(aDist2 < aCurDistMin) |
a86d3ec0 | 516 | { |
2a739b6d | 517 | aCurDistMin = aDist2; |
a86d3ec0 | 518 | anIndMin = i; |
519 | } | |
520 | } | |
2a739b6d | 521 | if(anIndMin) |
a86d3ec0 | 522 | { |
2a739b6d | 523 | curp = PR.Point(anIndMin).Parameter(); |
524 | if( curp > previousp + deltamin && curp < bornesup) | |
525 | { | |
526 | aDistMin = aCurDistMin; | |
527 | initp = previousp = pc3d[count] = curp; | |
528 | pcons[count] = pcons[ii]; | |
529 | count++; | |
530 | projok = Standard_True; | |
531 | ||
532 | } | |
a86d3ec0 | 533 | } |
2a739b6d | 534 | |
a86d3ec0 | 535 | } |
536 | } | |
537 | } | |
2a739b6d | 538 | if(projok && besttol2 < aDistMin) |
539 | besttol2 = aDistMin; | |
540 | ||
541 | else if(!projok) | |
a86d3ec0 | 542 | { |
543 | //Projector | |
0797d9d3 | 544 | #ifdef OCCT_DEBUG |
04232180 | 545 | std::cout << "Projection not done" << std::endl; |
7fd59977 | 546 | #endif |
547 | } | |
548 | } | |
549 | } | |
a86d3ec0 | 550 | |
7fd59977 | 551 | if(mySameParameter){ |
552 | myTolReached = 1.5*sqrt(dmax2); | |
553 | return; | |
554 | } | |
a86d3ec0 | 555 | |
52db4751 | 556 | if(!extrok) |
557 | { | |
558 | // If not already SameP and tangent to mill, abandon. | |
7fd59977 | 559 | mySameParameter = Standard_False; |
0797d9d3 | 560 | #ifdef OCCT_DEBUG |
04232180 | 561 | std::cout<<"SameParameter problem : zero tangent to extremities"<<std::endl; |
7fd59977 | 562 | #endif |
563 | return; | |
564 | } | |
565 | ||
566 | pcons[count] = lcons; | |
567 | pc3d[count] = lc3d; | |
568 | ||
52db4751 | 569 | // There is at least one point where same parameter is broken. |
570 | // Try to build B-spline interpolation curve with degree 3. | |
571 | // The loop is organized over number of poles. | |
a86d3ec0 | 572 | Standard_Boolean hasCountChanged = Standard_False; |
52db4751 | 573 | do |
a86d3ec0 | 574 | { |
0d969553 | 575 | // The tables and their limits for the interpolation. |
7fd59977 | 576 | Standard_Integer num_knots = count + 7; |
577 | Standard_Integer num_poles = count + 3; | |
578 | TColStd_Array1OfReal Paramc3d(*pc3d,1,count+1); | |
579 | TColStd_Array1OfReal Paramcons(*pcons,1,count+1); | |
580 | TColStd_Array1OfInteger ContactOrder(1,num_poles) ; | |
581 | TColStd_Array1OfReal Poles(1,num_poles) ; | |
582 | TColStd_Array1OfReal InterpolationParameters(1,num_poles) ; | |
583 | TColStd_Array1OfReal FlatKnots(1,num_knots) ; | |
a86d3ec0 | 584 | |
0d969553 | 585 | // Fill tables taking attention to end values. |
7fd59977 | 586 | ContactOrder.Init(0); |
587 | ContactOrder(2) = ContactOrder(num_poles - 1) = 1; | |
a86d3ec0 | 588 | |
7fd59977 | 589 | FlatKnots(1) = FlatKnots(2) = FlatKnots(3) = FlatKnots(4) = fc3d; |
590 | FlatKnots(num_poles + 1) = FlatKnots(num_poles + 2) = | |
591 | FlatKnots(num_poles + 3) = FlatKnots(num_poles + 4) = lc3d; | |
a86d3ec0 | 592 | |
7fd59977 | 593 | Poles(1) = fcons; Poles(num_poles) = lcons; |
594 | Poles(2) = tangent[0]; Poles(num_poles - 1) = tangent[1]; | |
a86d3ec0 | 595 | |
7fd59977 | 596 | InterpolationParameters(1) = InterpolationParameters(2) = fc3d; |
597 | InterpolationParameters(num_poles - 1) = InterpolationParameters(num_poles) = lc3d; | |
a86d3ec0 | 598 | |
7fd59977 | 599 | for (ii = 3; ii <= num_poles - 2; ii++) { |
600 | Poles(ii) = Paramcons(ii - 1); | |
601 | InterpolationParameters(ii) = FlatKnots(ii+2) = Paramc3d(ii - 1); | |
602 | } | |
603 | Standard_Integer inversion_problem; | |
604 | BSplCLib::Interpolate(3,FlatKnots,InterpolationParameters,ContactOrder, | |
a86d3ec0 | 605 | 1,Poles(1),inversion_problem); |
7fd59977 | 606 | if(inversion_problem) { |
9775fa61 | 607 | throw Standard_ConstructionError(); |
7fd59977 | 608 | } |
609 | ||
7fd59977 | 610 | // Test if par2d(par3d) is monotonous function or not ----- IFV, Jan 2000 |
611 | // and try to insert new point to improve BSpline interpolation | |
612 | ||
613 | Standard_Integer extrap_mode[2] ; | |
614 | extrap_mode[0] = extrap_mode[1] = 3; | |
615 | Standard_Real eval_result[2] ; | |
616 | Standard_Integer DerivativeRequest = 0; | |
617 | Standard_Real *PolesArray = | |
618 | (Standard_Real *) &Poles(Poles.Lower()) ; | |
619 | ||
620 | Standard_Integer newcount = 0; | |
621 | for (ii = 0; ii < count; ii++) { | |
a86d3ec0 | 622 | |
7fd59977 | 623 | newpcons[newcount] = pcons[ii]; |
624 | newpc3d[newcount] = pc3d[ii]; | |
625 | newcount++; | |
626 | ||
a86d3ec0 | 627 | if(count - ii + newcount == aMaxArraySize) continue; |
7fd59977 | 628 | |
629 | BSplCLib::Eval(.5*(pc3d[ii]+pc3d[ii+1]), Standard_False, DerivativeRequest, | |
a86d3ec0 | 630 | extrap_mode[0], 3, FlatKnots, 1, PolesArray[0], eval_result[0]); |
631 | ||
7fd59977 | 632 | if(eval_result[0] < pcons[ii] || eval_result[0] > pcons[ii+1]) { |
a86d3ec0 | 633 | Standard_Real ucons = 0.5*(pcons[ii]+pcons[ii+1]); |
634 | Standard_Real uc3d = 0.5*(pc3d[ii]+pc3d[ii+1]); | |
635 | ||
636 | CurveOnSurface.D0(ucons,Pcons); | |
637 | Projector.Perform(Pcons, uc3d); | |
638 | if (Projector.IsDone()) { | |
639 | curp = Projector.Point().Parameter(); | |
640 | Standard_Real dist_2 = Projector.SquareDistance(); | |
641 | if(dist_2 > besttol2) besttol2 = dist_2; | |
642 | projok = 1; | |
643 | } | |
644 | else { | |
645 | ProjectPointOnCurve(uc3d,Pcons,Tol,30,myC3d->Curve(),projok,curp); | |
646 | } | |
647 | if(projok){ | |
648 | if(curp > pc3d[ii] + deltamin && curp < pc3d[ii+1] - deltamin){ | |
649 | newpc3d[newcount] = curp; | |
650 | newpcons[newcount] = ucons; | |
651 | newcount ++; | |
652 | } | |
653 | } | |
654 | else { | |
0797d9d3 | 655 | #ifdef OCCT_DEBUG |
04232180 | 656 | std::cout << "Projection not done" << std::endl; |
7fd59977 | 657 | #endif |
a86d3ec0 | 658 | } |
7fd59977 | 659 | } |
a86d3ec0 | 660 | |
7fd59977 | 661 | } |
662 | ||
663 | newpc3d[newcount] = pc3d[count]; | |
664 | newpcons[newcount] = pcons[count]; | |
665 | Standard_Real * temp; | |
666 | temp = pc3d; | |
667 | pc3d = newpc3d; | |
668 | newpc3d = temp; | |
669 | temp = pcons; | |
670 | pcons = newpcons; | |
671 | newpcons = temp; | |
672 | ||
52db4751 | 673 | if((count != newcount) && newcount < aMaxArraySize) |
674 | { | |
675 | hasCountChanged = Standard_True; | |
676 | count = newcount; | |
677 | continue; | |
678 | } | |
7fd59977 | 679 | |
680 | count = newcount; | |
681 | ||
682 | Standard_Real algtol = sqrt(besttol2); | |
683 | ||
684 | interpolok = Check (FlatKnots, Poles, count+1, Paramc3d, Paramcons, | |
52db4751 | 685 | myC3d, CurveOnSurface, algtol, tolsov); |
7fd59977 | 686 | |
687 | if (Precision::IsInfinite(algtol)) { | |
688 | mySameParameter = Standard_False; | |
0797d9d3 | 689 | #ifdef OCCT_DEBUG |
04232180 | 690 | std::cout<<"SameParameter problem : function of interpolation of parametration at mills !!"<<std::endl; |
7fd59977 | 691 | #endif |
692 | return; | |
693 | } | |
694 | ||
695 | tolsov = algtol; | |
696 | ||
52db4751 | 697 | interpolok = (interpolok || // Good result. |
698 | count >= aMaxArraySize - 1 ); // Number of points. | |
7fd59977 | 699 | |
700 | if(interpolok) { | |
a86d3ec0 | 701 | Standard_Real besttol = sqrt(besttol2); |
52db4751 | 702 | |
7fd59977 | 703 | Handle(TColStd_HArray1OfReal) tol1d,tol2d,tol3d; |
704 | tol1d = new TColStd_HArray1OfReal(1,2) ; | |
705 | tol1d->SetValue(1, mySurf->UResolution(besttol)); | |
706 | tol1d->SetValue(2, mySurf->VResolution(besttol)); | |
707 | ||
708 | Approx_SameParameter_Evaluator ev (FlatKnots, Poles, myHCurve2d); | |
709 | AdvApprox_ApproxAFunction anApproximator(2,0,0,tol1d,tol2d,tol3d,fc3d,lc3d, | |
a86d3ec0 | 710 | Continuity,11,40,ev); |
7fd59977 | 711 | |
712 | if (anApproximator.IsDone() || anApproximator.HasResult()) { | |
a86d3ec0 | 713 | Adaptor3d_CurveOnSurface ACS = CurveOnSurface; |
714 | GeomLib_MakeCurvefromApprox aCurveBuilder(anApproximator) ; | |
715 | Handle(Geom2d_BSplineCurve) aC2d = aCurveBuilder.Curve2dFromTwo1d(1,2) ; | |
716 | Handle(Adaptor2d_HCurve2d) aHCurve2d = new Geom2dAdaptor_HCurve(aC2d); | |
717 | CurveOnSurface.Load(aHCurve2d); | |
718 | ||
719 | myTolReached = ComputeTolReached(myC3d,CurveOnSurface,NCONTROL); | |
720 | ||
721 | if(myTolReached > anErrorMAX) | |
722 | { | |
52db4751 | 723 | //This tolerance is too big. Probably, we will not be able to get |
a86d3ec0 | 724 | //edge with sameparameter in this case. |
725 | ||
726 | myDone = Standard_False; | |
727 | return; | |
728 | } | |
729 | ||
730 | if( (myTolReached < 250.0*besttol) || | |
731 | (count >= aMaxArraySize-2) || | |
732 | !hasCountChanged) //if count does not change after adding new point | |
733 | //(else we can have circularity) | |
734 | { | |
735 | myCurve2d = aC2d; | |
52db4751 | 736 | myHCurve2d = aHCurve2d; |
a86d3ec0 | 737 | myDone = Standard_True; |
738 | } | |
739 | else | |
740 | { | |
741 | interpolok = Standard_False; | |
742 | CurveOnSurface = ACS; | |
743 | } | |
7fd59977 | 744 | } |
745 | } | |
a86d3ec0 | 746 | |
747 | if(!interpolok) | |
748 | { | |
7fd59977 | 749 | |
a86d3ec0 | 750 | newcount = 0; |
7fd59977 | 751 | for(Standard_Integer n = 0; n < count; n++){ |
a86d3ec0 | 752 | newpc3d[newcount] = pc3d[n]; |
753 | newpcons[newcount] = pcons[n]; | |
754 | newcount ++; | |
755 | ||
756 | if(count - n + newcount == aMaxArraySize) continue; | |
757 | ||
758 | Standard_Real ucons = 0.5*(pcons[n]+pcons[n+1]); | |
759 | Standard_Real uc3d = 0.5*(pc3d[n]+pc3d[n+1]); | |
760 | ||
761 | CurveOnSurface.D0(ucons,Pcons); | |
762 | Projector.Perform(Pcons, uc3d); | |
763 | if (Projector.IsDone()) { | |
764 | curp = Projector.Point().Parameter(); | |
765 | Standard_Real dist_2 = Projector.SquareDistance(); | |
766 | if(dist_2 > besttol2) besttol2 = dist_2; | |
767 | projok = 1; | |
768 | } | |
769 | else { | |
770 | ProjectPointOnCurve(uc3d,Pcons,Tol,30,myC3d->Curve(),projok,curp); | |
771 | } | |
772 | if(projok){ | |
773 | if(curp > pc3d[n] + deltamin && curp < pc3d[n+1] - deltamin){ | |
774 | newpc3d[newcount] = curp; | |
775 | newpcons[newcount] = ucons; | |
776 | newcount ++; | |
777 | } | |
778 | } | |
779 | else { | |
0797d9d3 | 780 | #ifdef OCCT_DEBUG |
04232180 | 781 | std::cout << "Projection not done" << std::endl; |
7fd59977 | 782 | #endif |
a86d3ec0 | 783 | } |
7fd59977 | 784 | } |
785 | newpc3d[newcount] = pc3d[count]; | |
786 | newpcons[newcount] = pcons[count]; | |
787 | Standard_Real * tempx; | |
788 | tempx = pc3d; | |
789 | pc3d = newpc3d; | |
790 | newpc3d = tempx; | |
791 | tempx = pcons; | |
792 | pcons = newpcons; | |
793 | newpcons = tempx; | |
a86d3ec0 | 794 | |
795 | if(count != newcount) | |
796 | { | |
797 | count = newcount; | |
798 | hasCountChanged = Standard_True; | |
799 | } | |
800 | else | |
801 | { | |
802 | hasCountChanged = Standard_False; | |
803 | } | |
7fd59977 | 804 | } |
52db4751 | 805 | } while(!interpolok && hasCountChanged); |
806 | ||
807 | if (!myDone) | |
808 | { | |
809 | // Loop is finished unsuccessfully. Fix tolerance by maximal deviation, | |
810 | // using data from the last loop iteration. | |
811 | Standard_Integer num_knots = count + 7; | |
812 | Standard_Integer num_poles = count + 3; | |
813 | TColStd_Array1OfReal Paramc3d(*pc3d,1,count + 1); | |
814 | TColStd_Array1OfReal Paramcons(*pcons,1,count + 1); | |
815 | TColStd_Array1OfInteger ContactOrder(1,num_poles) ; | |
816 | TColStd_Array1OfReal Poles(1,num_poles) ; | |
817 | TColStd_Array1OfReal InterpolationParameters(1,num_poles) ; | |
818 | TColStd_Array1OfReal FlatKnots(1,num_knots) ; | |
819 | ||
820 | // Fill tables taking attention to end values. | |
821 | ContactOrder.Init(0); | |
822 | ContactOrder(2) = ContactOrder(num_poles - 1) = 1; | |
823 | ||
824 | FlatKnots(1) = FlatKnots(2) = FlatKnots(3) = FlatKnots(4) = fc3d; | |
825 | FlatKnots(num_poles + 1) = FlatKnots(num_poles + 2) = | |
826 | FlatKnots(num_poles + 3) = FlatKnots(num_poles + 4) = lc3d; | |
827 | ||
828 | Poles(1) = fcons; Poles(num_poles) = lcons; | |
829 | Poles(2) = tangent[0]; Poles(num_poles - 1) = tangent[1]; | |
830 | ||
831 | InterpolationParameters(1) = InterpolationParameters(2) = fc3d; | |
832 | InterpolationParameters(num_poles - 1) = InterpolationParameters(num_poles) = lc3d; | |
833 | ||
834 | for (ii = 3; ii <= num_poles - 2; ii++) | |
835 | { | |
836 | Poles(ii) = Paramcons(ii - 1); | |
837 | InterpolationParameters(ii) = FlatKnots(ii+2) = Paramc3d(ii - 1); | |
838 | } | |
839 | Standard_Integer inversion_problem; | |
840 | BSplCLib::Interpolate(3,FlatKnots,InterpolationParameters,ContactOrder, | |
841 | 1,Poles(1),inversion_problem); | |
842 | if(inversion_problem) | |
843 | { | |
9775fa61 | 844 | throw Standard_ConstructionError(); |
52db4751 | 845 | } |
846 | ||
847 | Standard_Real besttol = sqrt(besttol2); | |
848 | Handle(TColStd_HArray1OfReal) tol1d,tol2d,tol3d; | |
849 | tol1d = new TColStd_HArray1OfReal(1,2) ; | |
850 | tol1d->SetValue(1, mySurf->UResolution(besttol)); | |
851 | tol1d->SetValue(2, mySurf->VResolution(besttol)); | |
852 | ||
853 | Approx_SameParameter_Evaluator ev (FlatKnots, Poles, myHCurve2d); | |
854 | AdvApprox_ApproxAFunction anApproximator(2,0,0,tol1d,tol2d,tol3d,fc3d,lc3d, | |
855 | Continuity,11,40,ev); | |
856 | ||
857 | if (!anApproximator.IsDone() && | |
858 | !anApproximator.HasResult() ) | |
859 | { | |
860 | myDone = Standard_False; | |
861 | return; | |
862 | } | |
863 | ||
864 | GeomLib_MakeCurvefromApprox aCurveBuilder(anApproximator) ; | |
865 | Handle(Geom2d_BSplineCurve) aC2d = aCurveBuilder.Curve2dFromTwo1d(1,2) ; | |
866 | Handle(Adaptor2d_HCurve2d) aHCurve2d = new Geom2dAdaptor_HCurve(aC2d); | |
867 | CurveOnSurface.Load(aHCurve2d); | |
868 | ||
869 | myTolReached = ComputeTolReached(myC3d,CurveOnSurface,NCONTROL); | |
870 | ||
871 | if(myTolReached > anErrorMAX) | |
872 | { | |
873 | //This tolerance is too big. Probably, we will not be able to get | |
874 | //edge with sameparameter in this case. | |
875 | myDone = Standard_False; | |
876 | return; | |
877 | } | |
878 | ||
879 | myCurve2d = aC2d; | |
880 | myHCurve2d = aHCurve2d; | |
881 | myDone = Standard_True; | |
7fd59977 | 882 | } |
883 | } |