0023830: BRepExtrema_DistShapeShape does not find intersection of face with edge
[occt.git] / src / Extrema / Extrema_GenExtCS.cxx
CommitLineData
b311480e 1// Created on: 1995-07-18
2// Created by: Modelistation
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
23#include <Extrema_GenExtCS.ixx>
24
25#include <math_Vector.hxx>
26#include <math_FunctionSetRoot.hxx>
27#include <Precision.hxx>
28
29#include <Geom_Line.hxx>
30#include <Adaptor3d_HCurve.hxx>
31#include <GeomAdaptor_Curve.hxx>
32
33#include <Extrema_ExtCC.hxx>
34#include <Extrema_POnCurv.hxx>
35#include <Extrema_ExtPS.hxx>
36
37//=======================================================================
38//function : Extrema_GenExtCS
39//purpose :
40//=======================================================================
41
b311480e 42Extrema_GenExtCS::Extrema_GenExtCS()
7fd59977 43{
44 myDone = Standard_False;
45 myInit = Standard_False;
46}
47
48//=======================================================================
49//function : Extrema_GenExtCS
50//purpose :
51//=======================================================================
52
53 Extrema_GenExtCS::Extrema_GenExtCS(const Adaptor3d_Curve& C,
54 const Adaptor3d_Surface& S,
55 const Standard_Integer NbT,
56 const Standard_Integer NbU,
57 const Standard_Integer NbV,
58 const Standard_Real Tol1,
59 const Standard_Real Tol2)
60{
61 Initialize(S, NbU, NbV, Tol2);
62 Perform(C, NbT, Tol1);
63}
64
65//=======================================================================
66//function : Extrema_GenExtCS
67//purpose :
68//=======================================================================
69
70 Extrema_GenExtCS::Extrema_GenExtCS(const Adaptor3d_Curve& C,
71 const Adaptor3d_Surface& S,
72 const Standard_Integer NbT,
73 const Standard_Integer NbU,
74 const Standard_Integer NbV,
75 const Standard_Real tmin,
76 const Standard_Real tsup,
77 const Standard_Real Umin,
78 const Standard_Real Usup,
79 const Standard_Real Vmin,
80 const Standard_Real Vsup,
81 const Standard_Real Tol1,
82 const Standard_Real Tol2)
83{
84 Initialize(S, NbU, NbV, Umin,Usup,Vmin,Vsup,Tol2);
85 Perform(C, NbT, tmin, tsup, Tol1);
86}
87
88//=======================================================================
89//function : Initialize
90//purpose :
91//=======================================================================
92
93void Extrema_GenExtCS::Initialize(const Adaptor3d_Surface& S,
94 const Standard_Integer NbU,
95 const Standard_Integer NbV,
96 const Standard_Real Tol2)
97{
98 myumin = S.FirstUParameter();
99 myusup = S.LastUParameter();
100 myvmin = S.FirstVParameter();
101 myvsup = S.LastVParameter();
102 Initialize(S,NbU,NbV,myumin,myusup,myvmin,myvsup,Tol2);
103}
104
105//=======================================================================
106//function : Initialize
107//purpose :
108//=======================================================================
109
110void Extrema_GenExtCS::Initialize(const Adaptor3d_Surface& S,
111 const Standard_Integer NbU,
112 const Standard_Integer NbV,
113 const Standard_Real Umin,
114 const Standard_Real Usup,
115 const Standard_Real Vmin,
116 const Standard_Real Vsup,
117 const Standard_Real Tol2)
118{
119 myS = (Adaptor3d_SurfacePtr)&S;
7fd59977 120 myusample = NbU;
121 myvsample = NbV;
122 myumin = Umin;
123 myusup = Usup;
124 myvmin = Vmin;
125 myvsup = Vsup;
126 mytol2 = Tol2;
7fd59977 127}
128
129//=======================================================================
130//function : Perform
131//purpose :
132//=======================================================================
133
134void Extrema_GenExtCS::Perform(const Adaptor3d_Curve& C,
135 const Standard_Integer NbT,
136 const Standard_Real Tol1)
137{
138 mytmin = C.FirstParameter();
139 mytsup = C.LastParameter();
140 Perform(C, NbT, mytmin, mytsup,Tol1);
141}
142
143//=======================================================================
144//function : Perform
145//purpose :
146//=======================================================================
147
148void Extrema_GenExtCS::Perform(const Adaptor3d_Curve& C,
149 const Standard_Integer NbT,
150 const Standard_Real tmin,
151 const Standard_Real tsup,
152 const Standard_Real Tol1)
153{
154 myDone = Standard_False;
155 myF.Initialize(C,*myS);
156 mytmin = tmin;
157 mytsup = tsup;
158 mytol1 = Tol1;
159 mytsample = NbT;
b659a6dc 160 // Modif de lvt pour trimer la surface non pas aux infinis mais a +/- 10000
7fd59977 161
162 Standard_Real trimusup = myusup, trimumin = myumin,trimvsup = myvsup,trimvmin = myvmin;
163 if (Precision::IsInfinite(trimusup)){
164 trimusup = 1.0e+10;
165 }
166 if (Precision::IsInfinite(trimvsup)){
167 trimvsup = 1.0e+10;
168 }
169 if (Precision::IsInfinite(trimumin)){
170 trimumin = -1.0e+10;
171 }
172 if (Precision::IsInfinite(trimvmin)){
173 trimvmin = -1.0e+10;
174 }
b659a6dc 175 //
7fd59977 176 math_Vector Tol(1, 3);
177 Tol(1) = mytol1;
178 Tol(2) = mytol2;
179 Tol(3) = mytol2;
180 math_Vector UV(1,3), UVinf(1,3), UVsup(1,3);
181 UVinf(1) = mytmin;
182 UVinf(2) = trimumin;
183 UVinf(3) = trimvmin;
184
185 UVsup(1) = mytsup;
186 UVsup(2) = trimusup;
187 UVsup(3) = trimvsup;
7fd59977 188 // 18/02/02 akm vvv : (OCC163) bad extremas - extrusion surface versus the line.
189 // Try to preset the initial solution as extrema between
190 // extrusion direction and the curve.
191 if (myS->GetType() == GeomAbs_SurfaceOfExtrusion)
192 {
193 gp_Dir aDir = myS->Direction();
194 Handle(Adaptor3d_HCurve) aCurve = myS->BasisCurve();
195 Standard_Real dfUFirst = aCurve->FirstParameter();
196 // Create iso line of U=U0
197 GeomAdaptor_Curve anAx(new Geom_Line(aCurve->Value(dfUFirst), aDir),
198 trimvmin, trimvsup);
199 Extrema_ExtCC aLocator(C, anAx);
200 if (aLocator.IsDone() && aLocator.NbExt()>0)
201 {
202 Standard_Integer iExt;
203 // Try to find all extremas
204 Extrema_POnCurv aP1, aP2;
205 for (iExt=1; iExt<=aLocator.NbExt(); iExt++)
206 {
b659a6dc 207 aLocator.Points (iExt, aP1, aP2);
208 // Parameter on curve
209 UV(1) = aP1.Parameter();
210 // To find parameters on surf, try ExtPS
211 Extrema_ExtPS aPreciser (aP1.Value(), *myS, mytol2, mytol2);
212 if (aPreciser.IsDone())
213 {
214 // Managed to find extremas between point and surface
215 Standard_Integer iPExt;
216 for (iPExt=1; iPExt<=aPreciser.NbExt(); iPExt++)
217 {
218 aPreciser.Point(iPExt).Parameter(UV(2),UV(3));
219 math_FunctionSetRoot S1 (myF,UV,Tol,UVinf,UVsup);
220 }
221 }
222 else
223 {
224 // Failed... try the point on iso line
225 UV(2) = dfUFirst;
226 UV(3) = aP2.Parameter();
227 math_FunctionSetRoot S1 (myF,UV,Tol,UVinf,UVsup);
228 }
7fd59977 229 } // for (iExt=1; iExt<=aLocator.NbExt(); iExt++)
230 } // if (aLocator.IsDone() && aLocator.NbExt()>0)
231 } // if (myS.Type() == GeomAbs_ExtrusionSurface)
b659a6dc 232 else
7fd59977 233 {
b659a6dc 234 Standard_Real aCUAdd = (mytsup - mytmin) / mytsample;
235 Standard_Real aSUAdd = (myusup - myumin) / myusample;
236 Standard_Real aSVAdd = (myvsup - myvmin) / myvsample;
237 TColgp_HArray1OfPnt aCPs(1, mytsample);
238 TColgp_HArray2OfPnt aSPs(1, myusample, 1, myvsample);
239 Standard_Integer aRestIterCount = 3;
240 // The value is calculated by the bug CR23830.
241 Standard_Integer aCUDen = 2, aSUDen = 2, aSVDen = 2;
242 Standard_Boolean anAreAvSqsInited = Standard_False;
243 Standard_Real aCUSq = 0, aSUSq = 0, aSVSq = 0;
244 while (aRestIterCount--)
245 {
246 Standard_Real aMinCU, aMinSU, aMinSV, aMaxCU, aMaxSU, aMaxSV;
247 Standard_Real aMinSqDist = DBL_MAX, aMaxSqDist = -DBL_MAX;
248 for (Standard_Integer aSUNom = 1; aSUNom < aSUDen; aSUNom += 2)
249 {
250 Standard_Real aSU0 = myumin + (aSUNom * aSUAdd) / aSUDen;
251 for (Standard_Integer aSVNom = 1; aSVNom < aSVDen; aSVNom += 2)
252 {
253 Standard_Real aSV0 = myvmin + (aSVNom * aSVAdd) / aSVDen;
254 for (Standard_Integer aCUNom = 1; aCUNom < aCUDen; aCUNom += 2)
255 {
256 Standard_Real aCU0 = mytmin + (aCUNom * aCUAdd) / aCUDen;
257 Standard_Real aCU = aCU0;
258 for (Standard_Integer aCUI = 1; aCUI <= mytsample;
259 aCUI++, aCU += aCUAdd)
260 {
261 aCPs.SetValue(aCUI, C.Value(aCU));
262 }
263 //
264 aCU = aCU0;
265 Standard_Real aSU = aSU0;
266 for (Standard_Integer aSUI = 1; aSUI <= myusample;
267 aSUI++, aSU += aSUAdd)
268 {
269 Standard_Real aSV = aSV0;
270 for (Standard_Integer aSVI = 1; aSVI <= myvsample;
271 aSVI++, aSV += aSVAdd)
272 {
273 gp_Pnt aSP = myS->Value(aSU, aSV);
274 aSPs.ChangeValue(aSUI, aSVI) = aSP;
275 Standard_Real aCU = aCU0;
276 for (Standard_Integer aCUI = 1; aCUI <= mytsample;
277 aCUI++, aCU += aCUAdd)
278 {
279 Standard_Real aSqDist = aSP.SquareDistance(aCPs.Value(aCUI));
280 if (aSqDist < aMinSqDist)
281 {
282 aMinSqDist = aSqDist;
283 aMinCU = aCU;
284 aMinSU = aSU;
285 aMinSV = aSV;
286 }
287 if (aMaxSqDist < aSqDist)
288 {
289 aMaxSqDist = aSqDist;
290 aMaxCU = aCU;
291 aMaxSU = aSU;
292 aMaxSV = aSV;
293 }
294 }
295 }
296 }
297 }
298 }
299 }
300 // Find min approximation.
301 UV(1) = aMinCU;
302 UV(2) = aMinSU;
303 UV(3) = aMinSV;
304 math_FunctionSetRoot anA(myF, UV, Tol, UVinf, UVsup, 100, aRestIterCount);
305 // Find max approximation.
306 if (!anA.IsDivergent() || !aRestIterCount)
307 {
308 UV(1) = aMaxCU;
309 UV(2) = aMaxSU;
310 UV(3) = aMaxSV;
311 math_FunctionSetRoot(myF, UV, Tol, UVinf, UVsup);
312 break;
313 }
314 //
315 if (!anAreAvSqsInited)
316 {
317 anAreAvSqsInited = Standard_True;
318 //
319 const gp_Pnt * aCP1 = &aCPs.Value(1);
320 for (Standard_Integer aCUI = 2; aCUI <= mytsample; aCUI++)
321 {
322 const gp_Pnt & aCP2 = aCPs.Value(aCUI);
323 aCUSq += aCP1->SquareDistance(aCP2);
324 aCP1 = &aCP2;
325 }
326 aCUSq /= mytsample - 1;
327 //
328 for (Standard_Integer aSVI = 1; aSVI <= myvsample; aSVI++)
329 {
330 const gp_Pnt * aSP1 = &aSPs.Value(1, aSVI);
331 for (Standard_Integer aSUI = 2; aSUI <= myusample; aSUI++)
332 {
333 const gp_Pnt & aSP2 = aSPs.Value(aSUI, aSVI);
334 aSUSq += aSP1->SquareDistance(aSP2);
335 aSP1 = &aSP2;
336 }
337 }
338 aSUSq /= myvsample * (myusample - 1);
339 //
340 for (Standard_Integer aSUI = 1; aSUI <= myusample; aSUI++)
341 {
342 const gp_Pnt * aSP1 = &aSPs.Value(aSUI, 1);
343 for (Standard_Integer aSVI = 2; aSVI <= myvsample; aSVI++)
344 {
345 const gp_Pnt & aSP2 = aSPs.Value(aSUI, aSVI);
346 aSVSq += aSP1->SquareDistance(aSP2);
347 aSP1 = &aSP2;
348 }
349 }
350 aSVSq /= myusample * (myvsample - 1);
351 }
352 //
353 if ((aSUSq <= aCUSq) && (aSVSq <= aCUSq))
354 {
355 aCUDen += aCUDen;
356 aCUSq /= 4;
357 }
358 else if ((aCUSq <= aSUSq) && (aSVSq <= aSUSq))
359 {
360 aSUDen += aSUDen;
361 aSUSq /= 4;
362 }
363 else
364 {
365 aSVDen += aSVDen;
366 aSVSq /= 4;
7fd59977 367 }
368 }
7fd59977 369 }
370
371 myDone = Standard_True;
372}
373
374//=======================================================================
375//function : IsDone
376//purpose :
377//=======================================================================
378
379Standard_Boolean Extrema_GenExtCS::IsDone() const
380{
381 return myDone;
382}
383
384//=======================================================================
385//function : NbExt
386//purpose :
387//=======================================================================
388
389Standard_Integer Extrema_GenExtCS::NbExt() const
390{
391 if (!IsDone()) { StdFail_NotDone::Raise(); }
392 return myF.NbExt();
393}
394
395//=======================================================================
396//function : Value
397//purpose :
398//=======================================================================
399
400Standard_Real Extrema_GenExtCS::SquareDistance(const Standard_Integer N) const
401{
402 if (!IsDone()) { StdFail_NotDone::Raise(); }
403 return myF.SquareDistance(N);
404}
405
406//=======================================================================
407//function : PointOnCurve
408//purpose :
409//=======================================================================
410
411const Extrema_POnCurv& Extrema_GenExtCS::PointOnCurve(const Standard_Integer N) const
412{
413 if (!IsDone()) { StdFail_NotDone::Raise(); }
414 return myF.PointOnCurve(N);
415}
416
417//=======================================================================
418//function : PointOnSurface
419//purpose :
420//=======================================================================
421
422const Extrema_POnSurf& Extrema_GenExtCS::PointOnSurface(const Standard_Integer N) const
423{
424 if (!IsDone()) { StdFail_NotDone::Raise(); }
425 return myF.PointOnSurface(N);
426}
427
428//=======================================================================
429//function : BidonSurface
430//purpose :
431//=======================================================================
432
433Adaptor3d_SurfacePtr Extrema_GenExtCS::BidonSurface() const
434{
435 return (Adaptor3d_SurfacePtr)0L;
436}
437
438//=======================================================================
439//function : BidonCurve
440//purpose :
441//=======================================================================
442
443Adaptor3d_CurvePtr Extrema_GenExtCS::BidonCurve() const
444{
445 return (Adaptor3d_CurvePtr)0L;
446}
447