b311480e |
1 | // Created on: 1993-03-17 |
2 | // Created by: Laurent BUCHARD |
3 | // Copyright (c) 1993-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 | |
2c26a53d |
17 | #include <IntSurf_PntOn2S.hxx> |
7fd59977 |
18 | #include <TColStd_Array1OfReal.hxx> |
19 | #include <math_FunctionSetRoot.hxx> |
20 | #include <StdFail_NotDone.hxx> |
21 | #include <GeomAbs_SurfaceType.hxx> |
22 | #include <Precision.hxx> |
23 | |
4e14c88f |
24 | //======================================================================= |
25 | //function : IsSingular |
26 | //purpose : Returns TRUE if vectors theDU || theDV or if at least one |
27 | // of them has null-magnitude. |
28 | // theSqLinTol is square of linear tolerance. |
29 | // theAngTol is angular tolerance. |
30 | //======================================================================= |
31 | static Standard_Boolean IsSingular( const gp_Vec& theDU, |
32 | const gp_Vec& theDV, |
33 | const Standard_Real theSqLinTol, |
34 | const Standard_Real theAngTol) |
35 | { |
36 | gp_Vec aDU(theDU), aDV(theDV); |
7fd59977 |
37 | |
4e14c88f |
38 | const Standard_Real aSqMagnDU = aDU.SquareMagnitude(), |
39 | aSqMagnDV = aDV.SquareMagnitude(); |
40 | |
41 | if(aSqMagnDU < theSqLinTol) |
42 | return Standard_True; |
43 | |
44 | aDU.Divide(sqrt(aSqMagnDU)); |
45 | |
46 | if(aSqMagnDV < theSqLinTol) |
47 | return Standard_True; |
48 | |
49 | aDV.Divide(sqrt(aSqMagnDV)); |
50 | |
51 | //Here aDU and aDV vectors have magnitude 1.0. |
52 | |
53 | if(aDU.Crossed(aDV).SquareMagnitude() < theAngTol*theAngTol) |
54 | return Standard_True; |
55 | |
56 | return Standard_False; |
57 | } |
58 | |
59 | //======================================================================= |
60 | //function : SingularProcessing |
61 | //purpose : Computes 2D-representation (in UV-coordinates) of |
62 | // theTg3D vector on the surface in case when |
63 | // theDU.Crossed(theDV).Magnitude() == 0.0. Stores result in |
64 | // theTg2D variable. |
65 | // theDU and theDV are vectors of 1st derivative |
66 | // (with respect to U and V variables correspondingly). |
67 | // If theIsTo3DTgCompute == TRUE then theTg3D has not been |
68 | // defined yet (it should be computed). |
69 | // theLinTol is SQUARE of the tolerance. |
70 | // |
71 | //Algorithm: |
72 | // Condition |
73 | // Tg=theDU*theTg2D.X()+theDV*theTg2D.Y() |
74 | // has to be satisfied strictly. |
75 | // More over, vector Tg has to be NORMALYZED |
76 | // (if theIsTo3DTgCompute == TRUE then new computed vector will |
77 | // always have magnitude 1.0). |
78 | //======================================================================= |
79 | static Standard_Boolean SingularProcessing( const gp_Vec& theDU, |
80 | const gp_Vec& theDV, |
81 | const Standard_Boolean theIsTo3DTgCompute, |
82 | const Standard_Real theLinTol, |
83 | const Standard_Real theAngTol, |
84 | gp_Vec& theTg3D, |
85 | gp_Vec2d& theTg2D) |
86 | { |
87 | //Attention: @ \sin theAngTol \approx theAngTol @ (for cross-product) |
88 | |
89 | //Really, vector theTg3D has to be normalyzed (if theIsTo3DTgCompute == FALSE). |
90 | const Standard_Real aSQTan = theTg3D.SquareMagnitude(); |
91 | |
92 | const Standard_Real aSqMagnDU = theDU.SquareMagnitude(), |
93 | aSqMagnDV = theDV.SquareMagnitude(); |
94 | |
95 | //There are some reasons of singularity |
96 | |
97 | //1. |
98 | if((aSqMagnDU < theLinTol) && (aSqMagnDV < theLinTol)) |
99 | { |
100 | //For future, this case can be processed as same as in case of |
101 | //osculating surfaces (expanding in Taylor series). Here, |
102 | //we return only. |
103 | |
104 | return Standard_False; |
105 | } |
106 | |
107 | //2. |
108 | if(aSqMagnDU < theLinTol) |
109 | { |
110 | //In this case, theTg3D vector will be parallel with theDV. |
111 | //Its true direction shall be precised later (the algorithm is |
112 | //based on array of Walking-points). |
113 | |
114 | if(theIsTo3DTgCompute) |
115 | { |
116 | //theTg3D will be normalyzed. Its magnitude is |
117 | const Standard_Real aTgMagn = 1.0; |
118 | |
119 | const Standard_Real aNorm = sqrt(aSqMagnDV); |
120 | theTg3D = theDV.Divided(aNorm); |
121 | theTg2D.SetCoord(0.0, aTgMagn/aNorm); |
122 | } |
123 | else |
124 | { |
125 | //theTg3D is already defined. |
126 | //Here we check only, if this tangent is parallel to theDV. |
127 | |
128 | if(theDV.Crossed(theTg3D).SquareMagnitude() < |
129 | theAngTol*theAngTol*aSqMagnDV*aSQTan) |
130 | { |
131 | //theTg3D is parallel to theDV |
132 | |
133 | //Use sign "+" if theTg3D and theDV are codirectional |
134 | //and sign "-" if opposite |
135 | const Standard_Real aDP = theTg3D.Dot(theDV); |
136 | theTg2D.SetCoord(0.0, Sign(sqrt(aSQTan/aSqMagnDV), aDP)); |
137 | } |
138 | else |
139 | { |
140 | //theTg3D is not parallel to theDV |
141 | //It is abnormal |
142 | |
143 | return Standard_False; |
144 | } |
145 | } |
146 | |
147 | return Standard_True; |
148 | } |
149 | |
150 | //3. |
151 | if(aSqMagnDV < theLinTol) |
152 | { |
153 | //In this case, theTg3D vector will be parallel with theDU. |
154 | //Its true direction shall be precised later (the algorithm is |
155 | //based on array of Walking-points). |
156 | |
157 | if(theIsTo3DTgCompute) |
158 | { |
159 | //theTg3D will be normalyzed. Its magnitude is |
160 | const Standard_Real aTgMagn = 1.0; |
161 | |
162 | const Standard_Real aNorm = sqrt(aSqMagnDU); |
163 | theTg3D = theDU.Divided(aNorm); |
164 | theTg2D.SetCoord(aTgMagn/aNorm, 0.0); |
165 | } |
166 | else |
167 | { |
168 | //theTg3D is already defined. |
169 | //Here we check only, if this tangent is parallel to theDU. |
170 | |
171 | if(theDU.Crossed(theTg3D).SquareMagnitude() < |
172 | theAngTol*theAngTol*aSqMagnDU*aSQTan) |
173 | { |
174 | //theTg3D is parallel to theDU |
175 | |
176 | //Use sign "+" if theTg3D and theDU are codirectional |
177 | //and sign "-" if opposite |
178 | const Standard_Real aDP = theTg3D.Dot(theDU); |
179 | theTg2D.SetCoord(Sign(sqrt(aSQTan/aSqMagnDU), aDP), 0.0); |
180 | } |
181 | else |
182 | { |
183 | //theTg3D is not parallel to theDU |
184 | //It is abnormal |
185 | |
186 | return Standard_False; |
187 | } |
188 | } |
189 | |
190 | return Standard_True; |
191 | } |
192 | |
193 | //4. If aSqMagnDU > 0.0 && aSqMagnDV > 0.0 but theDV || theDU. |
194 | |
195 | const Standard_Real aLenU = sqrt(aSqMagnDU), |
196 | aLenV = sqrt(aSqMagnDV); |
197 | |
198 | //aLenSum > 0.0 definitely |
199 | const Standard_Real aLenSum = aLenU + aLenV; |
200 | |
201 | if(theDV.Dot(theDU) > 0.0) |
202 | { |
203 | //Vectors theDV and theDU are codirectional. |
204 | |
205 | if(theIsTo3DTgCompute) |
206 | { |
207 | theTg2D.SetCoord(1.0/aLenSum, 1.0/aLenSum); |
208 | theTg3D = theDU*theTg2D.X() + theDV*theTg2D.Y(); |
209 | } |
210 | else |
211 | { |
212 | //theTg3D is already defined. |
213 | //Here we check only, if this tangent is parallel to theDU |
214 | //(and theDV together). |
215 | |
216 | if(theDU.Crossed(theTg3D).SquareMagnitude() < |
217 | theAngTol*theAngTol*aSqMagnDU*aSQTan) |
218 | { |
219 | //theTg3D is parallel to theDU |
220 | |
221 | const Standard_Real aDP = theTg3D.Dot(theDU); |
222 | const Standard_Real aLenTg = Sign(sqrt(aSQTan), aDP); |
223 | theTg2D.SetCoord(aLenTg/aLenSum, aLenTg/aLenSum); |
224 | } |
225 | else |
226 | { |
227 | //theTg3D is not parallel to theDU |
228 | //It is abnormal |
229 | |
230 | return Standard_False; |
231 | } |
232 | } |
233 | } |
234 | else |
235 | { |
236 | //Vectors theDV and theDU are opposite. |
237 | |
238 | if(theIsTo3DTgCompute) |
239 | { |
240 | //Here we chose theDU as direction of theTg3D. |
241 | //True direction shall be precised later (the algorithm is |
242 | //based on array of Walking-points). |
243 | |
244 | theTg2D.SetCoord(1.0/aLenSum, -1.0/aLenSum); |
245 | theTg3D = theDU*theTg2D.X() + theDV*theTg2D.Y(); |
246 | } |
247 | else |
248 | { |
249 | //theTg3D is already defined. |
250 | //Here we check only, if this tangent is parallel to theDU |
251 | //(and theDV together). |
252 | |
253 | if(theDU.Crossed(theTg3D).SquareMagnitude() < |
254 | theAngTol*theAngTol*aSqMagnDU*aSQTan) |
255 | { |
256 | //theTg3D is parallel to theDU |
257 | |
258 | const Standard_Real aDP = theTg3D.Dot(theDU); |
259 | const Standard_Real aLenTg = Sign(sqrt(aSQTan), aDP); |
260 | theTg2D.SetCoord(aLenTg/aLenSum, -aLenTg/aLenSum); |
261 | } |
262 | else |
263 | { |
264 | //theTg3D is not parallel to theDU |
265 | //It is abnormal |
266 | |
267 | return Standard_False; |
268 | } |
269 | } |
270 | } |
271 | |
272 | return Standard_True; |
273 | } |
274 | |
275 | //======================================================================= |
276 | //function : NonSingularProcessing |
277 | //purpose : Computes 2D-representation (in UV-coordinates) of |
278 | // theTg3D vector on the surface in case when |
279 | // theDU.Crossed(theDV).Magnitude() > 0.0. Stores result in |
280 | // theTg2D variable. |
281 | // theDU and theDV are vectors of 1st derivative |
282 | // (with respect to U and V variables correspondingly). |
283 | // theLinTol is SQUARE of the tolerance. |
284 | // |
285 | //Algorithm: |
286 | // Condition |
287 | // Tg=theDU*theTg2D.X()+theDV*theTg2D.Y() |
288 | // has to be satisfied strictly. |
289 | // More over, vector Tg has always to be NORMALYZED. |
290 | //======================================================================= |
291 | static Standard_Boolean NonSingularProcessing(const gp_Vec& theDU, |
292 | const gp_Vec& theDV, |
293 | const gp_Vec& theTg3D, |
294 | const Standard_Real theLinTol, |
295 | const Standard_Real theAngTol, |
296 | gp_Vec2d& theTg2D) |
297 | { |
298 | const gp_Vec aNormal = theDU.Crossed(theDV); |
299 | const Standard_Real aSQMagn = aNormal.SquareMagnitude(); |
300 | |
301 | if(IsSingular(theDU, theDV, theLinTol, theAngTol)) |
302 | { |
303 | gp_Vec aTg(theTg3D); |
304 | return |
305 | SingularProcessing(theDU, theDV, Standard_False, |
306 | theLinTol, theAngTol, aTg, theTg2D); |
307 | } |
308 | |
309 | //If @\vec{T}=\vec{A}*U+\vec{B}*V@ then |
310 | |
311 | // \left\{\begin{matrix} |
312 | // \vec{A} \times \vec{T} = (\vec{A} \times \vec{B})*V |
313 | // \vec{B} \times \vec{T} = (\vec{B} \times \vec{A})*U |
314 | // \end{matrix}\right. |
315 | |
316 | //From here, values of U and V can be found very easily |
317 | //(if @\left \| \vec{A} \times \vec{B} \right \| > 0.0 @, |
318 | //else it is singular case). |
319 | |
320 | const gp_Vec aTgU(theTg3D.Crossed(theDU)), aTgV(theTg3D.Crossed(theDV)); |
321 | const Standard_Real aDeltaU = aTgV.SquareMagnitude()/aSQMagn; |
322 | const Standard_Real aDeltaV = aTgU.SquareMagnitude()/aSQMagn; |
323 | |
324 | theTg2D.SetCoord(Sign(sqrt(aDeltaU), aTgV.Dot(aNormal)), -Sign(sqrt(aDeltaV), aTgU.Dot(aNormal))); |
325 | |
326 | return Standard_True; |
327 | } |
7fd59977 |
328 | |
329 | //-------------------------------------------------------------------------------- |
330 | ApproxInt_ImpPrmSvSurfaces::ApproxInt_ImpPrmSvSurfaces( const TheISurface& ISurf |
331 | ,const ThePSurface& PSurf): |
332 | MyHasBeenComputed(Standard_False), |
333 | MyHasBeenComputedbis(Standard_False), |
334 | MyImplicitFirst(Standard_True), |
335 | MyZerImpFunc(PSurf,ISurf) |
336 | { |
337 | } |
338 | //-------------------------------------------------------------------------------- |
339 | ApproxInt_ImpPrmSvSurfaces::ApproxInt_ImpPrmSvSurfaces( const ThePSurface& PSurf |
340 | ,const TheISurface& ISurf): |
341 | MyHasBeenComputed(Standard_False), |
342 | MyHasBeenComputedbis(Standard_False), |
343 | MyImplicitFirst(Standard_False), |
344 | MyZerImpFunc(PSurf,ISurf) |
345 | { |
346 | } |
347 | //-------------------------------------------------------------------------------- |
348 | void ApproxInt_ImpPrmSvSurfaces::Pnt(const Standard_Real u1, |
349 | const Standard_Real v1, |
350 | const Standard_Real u2, |
351 | const Standard_Real v2, |
352 | gp_Pnt& P) { |
353 | gp_Pnt aP; |
354 | gp_Vec aT; |
355 | gp_Vec2d aTS1,aTS2; |
356 | Standard_Real tu1=u1; |
357 | Standard_Real tu2=u2; |
358 | Standard_Real tv1=v1; |
359 | Standard_Real tv2=v2; |
7fd59977 |
360 | this->Compute(tu1,tv1,tu2,tv2,aP,aT,aTS1,aTS2); |
7fd59977 |
361 | P=MyPnt; |
362 | } |
363 | //-------------------------------------------------------------------------------- |
364 | Standard_Boolean ApproxInt_ImpPrmSvSurfaces::Tangency(const Standard_Real u1, |
365 | const Standard_Real v1, |
366 | const Standard_Real u2, |
367 | const Standard_Real v2, |
368 | gp_Vec& T) { |
369 | gp_Pnt aP; |
370 | gp_Vec aT; |
371 | gp_Vec2d aTS1,aTS2; |
372 | Standard_Real tu1=u1; |
373 | Standard_Real tu2=u2; |
374 | Standard_Real tv1=v1; |
375 | Standard_Real tv2=v2; |
376 | Standard_Boolean t=this->Compute(tu1,tv1,tu2,tv2,aP,aT,aTS1,aTS2); |
377 | T=MyTg; |
378 | return(t); |
379 | } |
380 | //-------------------------------------------------------------------------------- |
381 | Standard_Boolean ApproxInt_ImpPrmSvSurfaces::TangencyOnSurf1(const Standard_Real u1, |
382 | const Standard_Real v1, |
383 | const Standard_Real u2, |
384 | const Standard_Real v2, |
385 | gp_Vec2d& T) { |
386 | gp_Pnt aP; |
387 | gp_Vec aT; |
388 | gp_Vec2d aTS1,aTS2; |
389 | Standard_Real tu1=u1; |
390 | Standard_Real tu2=u2; |
391 | Standard_Real tv1=v1; |
392 | Standard_Real tv2=v2; |
393 | Standard_Boolean t=this->Compute(tu1,tv1,tu2,tv2,aP,aT,aTS1,aTS2); |
394 | T=MyTguv1; |
395 | return(t); |
396 | } |
397 | //-------------------------------------------------------------------------------- |
398 | Standard_Boolean ApproxInt_ImpPrmSvSurfaces::TangencyOnSurf2(const Standard_Real u1, |
399 | const Standard_Real v1, |
400 | const Standard_Real u2, |
401 | const Standard_Real v2, |
402 | gp_Vec2d& T) { |
403 | gp_Pnt aP; |
404 | gp_Vec aT; |
405 | gp_Vec2d aTS1,aTS2; |
406 | Standard_Real tu1=u1; |
407 | Standard_Real tu2=u2; |
408 | Standard_Real tv1=v1; |
409 | Standard_Real tv2=v2; |
410 | Standard_Boolean t=this->Compute(tu1,tv1,tu2,tv2,aP,aT,aTS1,aTS2); |
411 | T=MyTguv2; |
412 | return(t); |
413 | } |
4e14c88f |
414 | |
415 | //======================================================================= |
416 | //function : Compute |
417 | //purpose : Computes point on curve, 3D and 2D-tangents of a curve and |
418 | // parameters on the surfaces. |
419 | //======================================================================= |
420 | Standard_Boolean ApproxInt_ImpPrmSvSurfaces::Compute( Standard_Real& u1, |
421 | Standard_Real& v1, |
422 | Standard_Real& u2, |
423 | Standard_Real& v2, |
424 | gp_Pnt& P, |
425 | gp_Vec& Tg, |
426 | gp_Vec2d& Tguv1, |
427 | gp_Vec2d& Tguv2) |
428 | { |
429 | const IntSurf_Quadric& aQSurf = MyZerImpFunc.ISurface(); |
430 | const ThePSurface& aPSurf = MyZerImpFunc.PSurface(); |
431 | gp_Vec2d& aQuadTg = MyImplicitFirst ? Tguv1 : Tguv2; |
432 | gp_Vec2d& aPrmTg = MyImplicitFirst ? Tguv2 : Tguv1; |
433 | |
434 | //for square |
435 | const Standard_Real aNullValue = Precision::Approximation()* |
436 | Precision::Approximation(), |
437 | anAngTol = Precision::Angular(); |
438 | |
7fd59977 |
439 | Standard_Real tu1=u1; |
440 | Standard_Real tu2=u2; |
441 | Standard_Real tv1=v1; |
442 | Standard_Real tv2=v2; |
443 | |
444 | if(MyHasBeenComputed) { |
445 | if( (MyParOnS1.X()==u1)&&(MyParOnS1.Y()==v1) |
446 | &&(MyParOnS2.X()==u2)&&(MyParOnS2.Y()==v2)) { |
447 | return(MyIsTangent); |
448 | } |
449 | else if(MyHasBeenComputedbis == Standard_False) { |
450 | MyTgbis = MyTg; |
451 | MyTguv1bis = MyTguv1; |
452 | MyTguv2bis = MyTguv2; |
453 | MyPntbis = MyPnt; |
454 | MyParOnS1bis = MyParOnS1; |
455 | MyParOnS2bis = MyParOnS2; |
456 | MyIsTangentbis = MyIsTangent; |
457 | MyHasBeenComputedbis = MyHasBeenComputed; |
458 | } |
459 | } |
460 | |
461 | if(MyHasBeenComputedbis) { |
462 | if( (MyParOnS1bis.X()==u1)&&(MyParOnS1bis.Y()==v1) |
463 | &&(MyParOnS2bis.X()==u2)&&(MyParOnS2bis.Y()==v2)) { |
464 | |
465 | gp_Vec TV(MyTg); |
466 | gp_Vec2d TV1(MyTguv1); |
467 | gp_Vec2d TV2(MyTguv2); |
468 | gp_Pnt TP(MyPnt); |
469 | gp_Pnt2d TP1(MyParOnS1); |
470 | gp_Pnt2d TP2(MyParOnS2); |
471 | Standard_Boolean TB=MyIsTangent; |
472 | |
473 | MyTg = MyTgbis; |
474 | MyTguv1 = MyTguv1bis; |
475 | MyTguv2 = MyTguv2bis; |
476 | MyPnt = MyPntbis; |
477 | MyParOnS1 = MyParOnS1bis; |
478 | MyParOnS2 = MyParOnS2bis; |
479 | MyIsTangent = MyIsTangentbis; |
480 | |
481 | MyTgbis = TV; |
482 | MyTguv1bis = TV1; |
483 | MyTguv2bis = TV2; |
484 | MyPntbis = TP; |
485 | MyParOnS1bis = TP1; |
486 | MyParOnS2bis = TP2; |
487 | MyIsTangentbis = TB; |
488 | |
489 | return(MyIsTangent); |
490 | } |
491 | } |
492 | |
2c26a53d |
493 | math_Vector X(1,2); |
494 | math_Vector BornInf(1,2), BornSup(1,2), Tolerance(1,2); |
495 | //--- ThePSurfaceTool::GetResolution(aPSurf,Tolerance(1),Tolerance(2)); |
496 | Tolerance(1) = 1.0e-8; Tolerance(2) = 1.0e-8; |
7fd59977 |
497 | Standard_Real binfu,bsupu,binfv,bsupv; |
4e14c88f |
498 | binfu = ThePSurfaceTool::FirstUParameter(aPSurf); |
499 | binfv = ThePSurfaceTool::FirstVParameter(aPSurf); |
500 | bsupu = ThePSurfaceTool::LastUParameter(aPSurf); |
501 | bsupv = ThePSurfaceTool::LastVParameter(aPSurf); |
7fd59977 |
502 | BornInf(1) = binfu; BornSup(1) = bsupu; |
503 | BornInf(2) = binfv; BornSup(2) = bsupv; |
2c26a53d |
504 | Standard_Real TranslationU = 0., TranslationV = 0.; |
7fd59977 |
505 | |
2c26a53d |
506 | if (!FillInitialVectorOfSolution(u1, v1, u2, v2, |
507 | binfu, bsupu, binfv, bsupv, |
508 | X, |
509 | TranslationU, TranslationV)) |
510 | { |
511 | MyIsTangent=MyIsTangentbis=Standard_False; |
512 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
513 | return(Standard_False); |
514 | } |
7fd59977 |
515 | |
516 | Standard_Real PourTesterU = X(1); |
517 | Standard_Real PourTesterV = X(2); |
2c26a53d |
518 | |
519 | math_FunctionSetRoot Rsnld(MyZerImpFunc); |
520 | Rsnld.SetTolerance(Tolerance); |
7fd59977 |
521 | Rsnld.Perform(MyZerImpFunc,X,BornInf,BornSup); |
522 | if(Rsnld.IsDone()) { |
523 | MyHasBeenComputed = Standard_True; |
524 | Rsnld.Root(X); |
525 | |
526 | Standard_Real DistAvantApresU = Abs(PourTesterU-X(1)); |
527 | Standard_Real DistAvantApresV = Abs(PourTesterV-X(2)); |
7fd59977 |
528 | |
4e14c88f |
529 | MyPnt = P = ThePSurfaceTool::Value(aPSurf, X(1), X(2)); |
530 | |
531 | if( (DistAvantApresV <= 0.001 ) && |
532 | (DistAvantApresU <= 0.001 )) |
533 | { |
534 | gp_Vec aD1uPrm,aD1vPrm; |
535 | gp_Vec aD1uQuad,aD1vQuad; |
536 | |
537 | if(MyImplicitFirst) |
538 | { |
539 | u2 = X(1)-TranslationU; |
540 | v2 = X(2)-TranslationV; |
541 | |
542 | if(aQSurf.TypeQuadric() != GeomAbs_Plane) |
543 | { |
544 | while(u1-tu1>M_PI) u1-=M_PI+M_PI; |
545 | while(tu1-u1>M_PI) u1+=M_PI+M_PI; |
546 | } |
547 | |
548 | MyParOnS1.SetCoord(tu1,tv1); |
549 | MyParOnS2.SetCoord(tu2,tv2); |
550 | |
551 | gp_Pnt aP2; |
552 | |
553 | ThePSurfaceTool::D1(aPSurf, X(1), X(2), P, aD1uPrm, aD1vPrm); |
554 | aQSurf.D1(u1,v1, aP2, aD1uQuad, aD1vQuad); |
555 | |
556 | //Middle-point of P-P2 segment |
557 | P.BaryCenter(1.0, aP2, 1.0); |
7fd59977 |
558 | } |
4e14c88f |
559 | else |
560 | { |
561 | u1 = X(1)-TranslationU; |
562 | v1 = X(2)-TranslationV; |
563 | //aQSurf.Parameters(P, u2, v2); |
564 | if(aQSurf.TypeQuadric() != GeomAbs_Plane) |
565 | { |
566 | while(u2-tu2>M_PI) u2-=M_PI+M_PI; |
567 | while(tu2-u2>M_PI) u2+=M_PI+M_PI; |
568 | } |
569 | |
570 | MyParOnS1.SetCoord(tu1,tv1); |
571 | MyParOnS2.SetCoord(tu2,tu2); |
572 | |
573 | gp_Pnt aP2; |
574 | ThePSurfaceTool::D1(aPSurf, X(1), X(2), P, aD1uPrm, aD1vPrm); |
575 | |
576 | aQSurf.D1(u2, v2, aP2, aD1uQuad, aD1vQuad); |
577 | |
578 | //Middle-point of P-P2 segment |
579 | P.BaryCenter(1.0, aP2, 1.0); |
7fd59977 |
580 | } |
4e14c88f |
581 | |
582 | MyPnt = P; |
583 | |
584 | //Normals to the surfaces |
585 | gp_Vec aNormalPrm(aD1uPrm.Crossed(aD1vPrm)), |
586 | aNormalImp(aQSurf.Normale(MyPnt)); |
587 | |
588 | const Standard_Real aSQMagnPrm = aNormalPrm.SquareMagnitude(), |
589 | aSQMagnImp = aNormalImp.SquareMagnitude(); |
590 | |
591 | Standard_Boolean isPrmSingular = Standard_False, |
592 | isImpSingular = Standard_False; |
593 | |
594 | if(IsSingular(aD1uPrm, aD1vPrm, aNullValue, anAngTol)) |
595 | { |
596 | isPrmSingular = Standard_True; |
597 | if(!SingularProcessing(aD1uPrm, aD1vPrm, Standard_True, |
598 | aNullValue, anAngTol, Tg, aPrmTg)) |
599 | { |
600 | MyIsTangent = Standard_False; |
601 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
602 | return Standard_False; |
603 | } |
604 | |
605 | MyTg = Tg; |
7fd59977 |
606 | } |
4e14c88f |
607 | else |
608 | { |
609 | aNormalPrm.Divide(sqrt(aSQMagnPrm)); |
7fd59977 |
610 | } |
4e14c88f |
611 | |
612 | //Analogicaly for implicit surface |
613 | if(aSQMagnImp < aNullValue) |
614 | { |
615 | isImpSingular = Standard_True; |
616 | |
617 | if(!SingularProcessing(aD1uQuad, aD1vQuad, !isPrmSingular, |
618 | aNullValue, anAngTol, Tg, aQuadTg)) |
619 | { |
620 | MyIsTangent = Standard_False; |
621 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
622 | return Standard_False; |
623 | } |
624 | |
625 | MyTg = Tg; |
626 | } |
627 | else |
628 | { |
629 | aNormalImp.Divide(sqrt(aSQMagnImp)); |
7fd59977 |
630 | } |
631 | |
4e14c88f |
632 | if(isImpSingular && isPrmSingular) |
633 | { |
634 | //All is OK. All abnormal cases were processed above. |
7fd59977 |
635 | |
4e14c88f |
636 | MyTguv1 = Tguv1; |
637 | MyTguv2 = Tguv2; |
638 | |
639 | MyIsTangent=Standard_True; |
640 | return MyIsTangent; |
7fd59977 |
641 | } |
4e14c88f |
642 | else if(!(isImpSingular || isPrmSingular)) |
643 | { |
644 | //Processing pure non-singular case |
645 | //(3D- and 2D-tangents are still not defined) |
646 | |
647 | //Ask to pay attention to the fact that here |
648 | //aNormalImp and aNormalPrm are normalyzed. |
649 | //Therefore, @ \left \| \vec{Tg} \right \| = 0.0 @ |
650 | //if and only if (aNormalImp || aNormalPrm). |
651 | Tg = aNormalImp.Crossed(aNormalPrm); |
652 | } |
653 | |
654 | const Standard_Real aSQMagnTg = Tg.SquareMagnitude(); |
655 | |
656 | if(aSQMagnTg < aNullValue) |
657 | { |
658 | MyIsTangent = Standard_False; |
659 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
660 | return Standard_False; |
661 | } |
662 | |
663 | //Normalyze Tg vector |
664 | Tg.Divide(sqrt(aSQMagnTg)); |
665 | MyTg = Tg; |
666 | |
667 | if(!isPrmSingular) |
668 | { |
669 | //If isPrmSingular==TRUE then aPrmTg has already been computed. |
670 | |
671 | if(!NonSingularProcessing(aD1uPrm, aD1vPrm, Tg, aNullValue, anAngTol, aPrmTg)) |
672 | { |
673 | MyIsTangent = Standard_False; |
674 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
675 | return Standard_False; |
676 | } |
677 | } |
678 | |
679 | if(!isImpSingular) |
680 | { |
681 | //If isImpSingular==TRUE then aQuadTg has already been computed. |
682 | |
683 | if(!NonSingularProcessing(aD1uQuad, aD1vQuad, Tg, aNullValue, anAngTol, aQuadTg)) |
684 | { |
685 | MyIsTangent = Standard_False; |
686 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
687 | return Standard_False; |
688 | } |
7fd59977 |
689 | } |
4e14c88f |
690 | |
691 | MyTguv1 = Tguv1; |
692 | MyTguv2 = Tguv2; |
693 | |
7fd59977 |
694 | MyIsTangent=Standard_True; |
4e14c88f |
695 | |
696 | #ifdef OCCT_DEBUG |
697 | //cout << "+++++++++++++++++ ApproxInt_ImpPrmSvSurfaces::Compute(...) ++++++++++" << endl; |
698 | //printf( "P2d_1(%+10.20f, %+10.20f); P2d_2(%+10.20f, %+10.20f)\n" |
699 | // "P(%+10.20f, %+10.20f, %+10.20f);\n" |
700 | // "Tg = {%+10.20f, %+10.20f, %+10.20f};\n" |
701 | // "Tguv1 = {%+10.20f, %+10.20f};\n" |
702 | // "Tguv2 = {%+10.20f, %+10.20f}\n", |
703 | // u1, v1, u2, v2, |
704 | // P.X(), P.Y(), P.Z(), |
705 | // Tg.X(), Tg.Y(), Tg.Z(), |
706 | // Tguv1.X(), Tguv1.Y(), Tguv2.X(), Tguv2.Y()); |
707 | //cout << "-----------------------------------------------------------------------" << endl; |
708 | #endif |
709 | |
710 | return Standard_True; |
7fd59977 |
711 | } |
712 | else { |
7fd59977 |
713 | //-- cout<<" ApproxInt_ImpImpSvSurfaces.gxx : Distance apres recadrage Trop Grande "<<endl; |
714 | |
715 | MyIsTangent=Standard_False; |
716 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
4e14c88f |
717 | return Standard_False; |
7fd59977 |
718 | } |
719 | } |
720 | else { |
721 | MyIsTangent = Standard_False; |
722 | MyHasBeenComputed = MyHasBeenComputedbis = Standard_False; |
4e14c88f |
723 | return Standard_False; |
7fd59977 |
724 | } |
725 | } |
726 | |
2c26a53d |
727 | //======================================================================= |
728 | //function : SeekPoint |
729 | //purpose : Computes point on curve and |
730 | // parameters on the surfaces. |
731 | //======================================================================= |
732 | Standard_Boolean ApproxInt_ImpPrmSvSurfaces::SeekPoint(const Standard_Real u1, |
733 | const Standard_Real v1, |
734 | const Standard_Real u2, |
735 | const Standard_Real v2, |
736 | IntSurf_PntOn2S& Point) { |
737 | gp_Pnt aP; |
738 | gp_Vec aT; |
739 | gp_Vec2d aTS1,aTS2; |
740 | const IntSurf_Quadric& aQSurf = MyZerImpFunc.ISurface(); |
741 | const ThePSurface& aPSurf = MyZerImpFunc.PSurface(); |
7fd59977 |
742 | |
2c26a53d |
743 | math_Vector X(1,2); |
744 | math_Vector BornInf(1,2), BornSup(1,2), Tolerance(1,2); |
745 | //--- ThePSurfaceTool::GetResolution(aPSurf,Tolerance(1),Tolerance(2)); |
746 | Tolerance(1) = 1.0e-8; Tolerance(2) = 1.0e-8; |
747 | Standard_Real binfu,bsupu,binfv,bsupv; |
748 | binfu = ThePSurfaceTool::FirstUParameter(aPSurf); |
749 | binfv = ThePSurfaceTool::FirstVParameter(aPSurf); |
750 | bsupu = ThePSurfaceTool::LastUParameter(aPSurf); |
751 | bsupv = ThePSurfaceTool::LastVParameter(aPSurf); |
752 | BornInf(1) = binfu; BornSup(1) = bsupu; |
753 | BornInf(2) = binfv; BornSup(2) = bsupv; |
754 | Standard_Real TranslationU = 0., TranslationV = 0.; |
755 | |
756 | if (!FillInitialVectorOfSolution(u1, v1, u2, v2, |
757 | binfu, bsupu, binfv, bsupv, |
758 | X, |
759 | TranslationU, TranslationV)) |
760 | return Standard_False; |
761 | |
762 | Standard_Real NewU1, NewV1, NewU2, NewV2; |
763 | |
764 | math_FunctionSetRoot Rsnld(MyZerImpFunc); |
765 | Rsnld.SetTolerance(Tolerance); |
766 | Rsnld.Perform(MyZerImpFunc,X,BornInf,BornSup); |
767 | if(Rsnld.IsDone()) { |
768 | MyHasBeenComputed = Standard_True; |
769 | Rsnld.Root(X); |
770 | |
771 | MyPnt = ThePSurfaceTool::Value(aPSurf, X(1), X(2)); |
772 | |
773 | if(MyImplicitFirst) |
774 | { |
775 | NewU2 = X(1)-TranslationU; |
776 | NewV2 = X(2)-TranslationV; |
777 | |
778 | aQSurf.Parameters(MyPnt, NewU1, NewV1); |
779 | //adjust U |
780 | if (aQSurf.TypeQuadric() != GeomAbs_Plane) |
781 | { |
782 | Standard_Real sign = (NewU1 > u1)? -1 : 1; |
783 | while (Abs(u1 - NewU1) > M_PI) |
784 | NewU1 += sign*(M_PI+M_PI); |
785 | } |
786 | } |
787 | else |
788 | { |
789 | NewU1 = X(1)-TranslationU; |
790 | NewV1 = X(2)-TranslationV; |
791 | |
792 | aQSurf.Parameters(MyPnt, NewU2, NewV2); |
793 | //adjust U |
794 | if (aQSurf.TypeQuadric() != GeomAbs_Plane) |
795 | { |
796 | Standard_Real sign = (NewU2 > u2)? -1 : 1; |
797 | while (Abs(u2 - NewU2) > M_PI) |
798 | NewU2 += sign*(M_PI+M_PI); |
799 | } |
800 | } |
801 | } |
802 | else |
803 | return Standard_False; |
804 | |
805 | Point.SetValue(MyPnt, NewU1, NewV1, NewU2, NewV2); |
806 | return Standard_True; |
807 | } |
808 | //-------------------------------------------------------------------------------- |
7fd59977 |
809 | |
2c26a53d |
810 | Standard_Boolean |
811 | ApproxInt_ImpPrmSvSurfaces::FillInitialVectorOfSolution(const Standard_Real u1, |
812 | const Standard_Real v1, |
813 | const Standard_Real u2, |
814 | const Standard_Real v2, |
815 | const Standard_Real binfu, |
816 | const Standard_Real bsupu, |
817 | const Standard_Real binfv, |
818 | const Standard_Real bsupv, |
819 | math_Vector& X, |
820 | Standard_Real& TranslationU, |
821 | Standard_Real& TranslationV) |
822 | { |
823 | const ThePSurface& aPSurf = MyZerImpFunc.PSurface(); |
7fd59977 |
824 | |
2c26a53d |
825 | math_Vector F(1,1); |
7fd59977 |
826 | |
2c26a53d |
827 | TranslationU = 0.0; |
828 | TranslationV = 0.0; |
7fd59977 |
829 | |
2c26a53d |
830 | if(MyImplicitFirst) { |
831 | if(u2<binfu-0.0000000001) { |
832 | if(ThePSurfaceTool::IsUPeriodic(aPSurf)) { |
833 | Standard_Real d = ThePSurfaceTool::UPeriod(aPSurf); |
834 | do { TranslationU+=d; } while(u2+TranslationU < binfu); |
835 | } |
836 | else |
837 | return(Standard_False); |
838 | } |
839 | else if(u2>bsupu+0.0000000001) { |
840 | if(ThePSurfaceTool::IsUPeriodic(aPSurf)) { |
841 | Standard_Real d = ThePSurfaceTool::UPeriod(aPSurf); |
842 | do { TranslationU-=d; } while(u2+TranslationU > bsupu); |
843 | } |
844 | else |
845 | return(Standard_False); |
846 | } |
847 | if(v2<binfv-0.0000000001) { |
848 | if(ThePSurfaceTool::IsVPeriodic(aPSurf)) { |
849 | Standard_Real d = ThePSurfaceTool::VPeriod(aPSurf); |
850 | do { TranslationV+=d; } while(v2+TranslationV < binfv); |
851 | } |
852 | else |
853 | return(Standard_False); |
854 | } |
855 | else if(v2>bsupv+0.0000000001) { |
856 | if(ThePSurfaceTool::IsVPeriodic(aPSurf)) { |
857 | Standard_Real d = ThePSurfaceTool::VPeriod(aPSurf); |
858 | do { TranslationV-=d; } while(v2+TranslationV > bsupv); |
859 | } |
860 | else |
861 | return(Standard_False); |
862 | } |
863 | X(1) = u2+TranslationU; |
864 | X(2) = v2+TranslationV; |
865 | } |
866 | else { |
867 | if(u1<binfu-0.0000000001) { |
868 | if(ThePSurfaceTool::IsUPeriodic(aPSurf)) { |
869 | Standard_Real d = ThePSurfaceTool::UPeriod(aPSurf); |
870 | do { TranslationU+=d; } while(u1+TranslationU < binfu); |
871 | } |
872 | else |
873 | return(Standard_False); |
874 | } |
875 | else if(u1>bsupu+0.0000000001) { |
876 | if(ThePSurfaceTool::IsUPeriodic(aPSurf)) { |
877 | Standard_Real d = ThePSurfaceTool::UPeriod(aPSurf); |
878 | do { TranslationU-=d; } while(u1+TranslationU > bsupu); |
879 | } |
880 | else |
881 | return(Standard_False); |
882 | } |
883 | if(v1<binfv-0.0000000001) { |
884 | if(ThePSurfaceTool::IsVPeriodic(aPSurf)) { |
885 | Standard_Real d = ThePSurfaceTool::VPeriod(aPSurf); |
886 | do { TranslationV+=d; } while(v1+TranslationV < binfv); |
887 | } |
888 | else |
889 | return(Standard_False); |
890 | } |
891 | else if(v1>bsupv+0.0000000001) { |
892 | if(ThePSurfaceTool::IsVPeriodic(aPSurf)) { |
893 | Standard_Real d = ThePSurfaceTool::VPeriod(aPSurf); |
894 | do { TranslationV-=d; } while(v1+TranslationV > bsupv); |
895 | } |
896 | else |
897 | return(Standard_False); |
898 | } |
899 | X(1) = u1+TranslationU; |
900 | X(2) = v1+TranslationV; |
901 | } |
902 | |
903 | //---------------------------------------------------- |
904 | //Make a small step from boundaries in order to avoid |
905 | //finding "outboundaried" solution (Rsnld -> NotDone). |
906 | if(X(1)-0.0000000001 <= binfu) X(1)=X(1)+0.0000001; |
907 | if(X(1)+0.0000000001 >= bsupu) X(1)=X(1)-0.0000001; |
908 | if(X(2)-0.0000000001 <= binfv) X(2)=X(2)+0.0000001; |
909 | if(X(2)+0.0000000001 >= bsupv) X(2)=X(2)-0.0000001; |
910 | |
911 | return Standard_True; |
912 | } |