projects / occt.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
0022610: The algorithm GeomAPI_ProjectPointOnSurf produces wrong results
[occt.git] / src / Extrema / Extrema_GenExtCC.gxx
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#include <StdFail_NotDone.hxx>
23#include <math_FunctionSetRoot.hxx>
24#include <math_NewtonFunctionSetRoot.hxx>
25#include <TColStd_Array2OfInteger.hxx>
26#include <TColStd_Array2OfReal.hxx>
27#include <Standard_NullObject.hxx>
28#include <Standard_OutOfRange.hxx>
29#include <Precision.hxx>
30
31Extrema_GenExtCC::Extrema_GenExtCC () : myDone (Standard_False)
32{
33}
34
35Extrema_GenExtCC::Extrema_GenExtCC (const Curve1& C1,
36 const Curve2& C2,
37 const Standard_Integer NbU,
38 const Standard_Integer NbV,
39 const Standard_Real TolU,
40 const Standard_Real TolV) : myF (C1,C2, Min(TolU, TolV)), myDone (Standard_False)
41{
42 SetCurveCache (1, new Cache (C1, C1.FirstParameter(), C1.LastParameter(), NbU, Standard_True));
43 SetCurveCache (2, new Cache (C2, C2.FirstParameter(), C2.LastParameter(), NbV, Standard_True));
44 Perform();
45}
46
47
48Extrema_GenExtCC::Extrema_GenExtCC (const Curve1& C1,
49 const Curve2& C2,
50 const Standard_Real Uinf,
51 const Standard_Real Usup,
52 const Standard_Real Vinf,
53 const Standard_Real Vsup,
54 const Standard_Integer NbU,
55 const Standard_Integer NbV,
56 const Standard_Real TolU,
57 const Standard_Real TolV) : myF (C1,C2), myDone (Standard_False)
58{
59 SetCurveCache (1, new Cache (C1, Uinf, Usup, NbU, Standard_True));
60 SetCurveCache (2, new Cache (C2, Vinf, Vsup, NbV, Standard_True));
61 Perform();
62}
63
64void Extrema_GenExtCC::SetCurveCache (const Standard_Integer theRank,
65 const Handle(Cache)& theCache)
66{
67 Standard_OutOfRange_Raise_if (theRank < 1 || theRank > 2, "Extrema_GenExtCC::SetCurveCache()")
68 myF.SetCurve (theRank, *(Curve1*)theCache->CurvePtr());
69 Standard_Integer anInd = theRank - 1;
70 myCache[anInd] = theCache;
71}
72
73void Extrema_GenExtCC::SetTolerance (const Standard_Real Tol)
74{
75 myF.SetTolerance (Tol);
76}
77
78
79//=============================================================================
80void Extrema_GenExtCC::Perform ()
81/*-----------------------------------------------------------------------------
82Fonction:
83 Recherche de toutes les distances extremales entre les courbes C1 et C2.
84 a partir de 2 echantillonnages (NbU,NbV).
85
86Methode:
87 L'algorithme part de l'hypothese que les echantillonnages sont suffisamment
88 fins pour que, s'il existe N distances extremales entre les 2 courbes,
89 alors il existe aussi N extrema entre les 2 ensembles de points.
90 Ainsi, l'algorithme consiste a partir des extrema des echantillons
91 pour trouver les extrema des courbes.
92 Les extrema sont calcules par l'algorithme math_FunctionSetRoot avec les
93 arguments suivants:
94 - myF: Extrema_FuncExtCC cree a partir de C1 et C2,
95 - UV: math_Vector dont les composantes sont les parametres des points de
96 l'extremum sur les ensembles de points,
97 - Tol: Min(TolU,TolV), (Prov.:math_FunctionSetRoot n'autorise pas un vecteur)
98 - UVinf: math_Vector dont les composantes sont les bornes inferieures de u et
99 v,
100 - UVsup: math_Vector dont les composantes sont les bornes superieures de u et
101 v.
102
103Traitement:
104 a- Constitution du tableau des square distances (TbDist2(0,NbU+1,0,NbV+1)):
105 Le tableau est volontairement etendu; les lignes 0 et NbU+1 et les
106 colonnes 0 et NbV+1 seront initialisees a RealFirst() ou RealLast()
107 pour simplifier les tests effectues dans l'etape b
108 (on n'a pas besoin de tester si le point est sur une extremite).
109 b- Calcul des extrema:
110 On recherche d'abord les minima et ensuite les maxima. Ces 2 traitements
111 se passent de facon similaire:
112 b.a- Initialisations:
113 - des 'bords' du tableau TbDist2 (a RealLast() dans le cas des minima
114 et a RealLast() dans le cas des maxima),
115 - du tableau TbSel(0,NbU+1,0,NbV+1) de selection des points pour un
116 calcul d'extremum local (a 0). Lorsqu'un couple de points sera
117 selectionne, il ne sera plus selectionnable, ainsi que les couples
118 adjacents (8 au maximum).
119 Les adresses correspondantes seront mises a 1.
120 b.b- Calcul des minima (ou maxima):
121 On boucle sur toutes les square distances du tableau TbDist2:
122 - recherche d'un minimum (ou maximum) sur les echantillons,
123 - calcul de l'extremum sur les courbes,
124 - mise a jour du tableau TbSel.
125-----------------------------------------------------------------------------*/
126{
127 myDone = Standard_False;
128
129 const Handle(Cache)& aCache1 = myCache[0];
130 const Handle(Cache)& aCache2 = myCache[1];
131 Standard_NullObject_Raise_if ((aCache1.IsNull() || aCache2.IsNull()),
132 "Extrema_GenExtCC::Perform()")
133
134 Standard_Integer aNbU = aCache1->NbSamples(), aNbV = aCache2->NbSamples();
135 Standard_OutOfRange_Raise_if ((aNbU < 2 ||aNbV < 2), "Extrema_GenExtCC::Perform()")
136
137/*
138a- Constitution du tableau des distances (TbDist2(0,NbU+1,0,NbV+1)):
139 ---------------------------------------------------------------
140*/
141
142 //ensure that caches have been calculated
143 if (!aCache1->IsValid())
144 aCache1->CalculatePoints();
145 if (!aCache2->IsValid())
146 aCache2->CalculatePoints();
147
148// Calcul des distances
149 //initialization of variables in the same way as in CalculatePoints()
150 Standard_Real PasU = aCache1->TrimLastParameter() - aCache1->TrimFirstParameter();
151 Standard_Real PasV = aCache2->TrimLastParameter() - aCache2->TrimFirstParameter();
152 Standard_Real U0 = PasU / aNbU / 100.;
153 Standard_Real V0 = PasV / aNbV / 100.;
154 PasU = (PasU - U0) / (aNbU - 1);
155 PasV = (PasV - V0) / (aNbV - 1);
156 U0 = aCache1->TrimFirstParameter() + (U0/2.);
157 V0 = aCache2->TrimFirstParameter() + (V0/2.);
158
159 const Curve1& aCurve1 = *((Curve1*)(myF.CurvePtr (1)));
160 const Curve2& aCurve2 = *((Curve1*)(myF.CurvePtr (2)));
161 const Handle(ArrayOfPnt)& aPntArray1 = aCache1->Points();
162 const Handle(ArrayOfPnt)& aPntArray2 = aCache2->Points();
163 Standard_Integer NoU, NoV;
164 TColStd_Array2OfReal TbDist2(0, aNbU + 1, 0, aNbV + 1);
165 for (NoU = 1; NoU <= aNbU; NoU++) {
166 const Pnt& P1 = aPntArray1->Value (NoU);
167 for (NoV = 1; NoV <= aNbV; NoV++) {
168 const Pnt& P2 = aPntArray2->Value (NoV);
169 TbDist2(NoU,NoV) = P1.SquareDistance(P2);
170 }
171 }
172
173/*
174b- Calcul des minima:
175 -----------------
176 b.a) Initialisations:
177*/
178// - generales
179 math_Vector Tol(1, 2);
180 Tol(1) = myF.Tolerance();
181 Tol(2) = myF.Tolerance();
182 math_Vector UV(1,2), UVinf(1,2), UVsup(1,2);
183 UVinf(1) = aCache1->TrimFirstParameter();
184 UVinf(2) = aCache2->TrimFirstParameter();
185 UVsup(1) = aCache1->TrimLastParameter();
186 UVsup(2) = aCache2->TrimLastParameter();
187
188// - des 'bords' du tableau TbDist2
189 for (NoV = 0; NoV <= aNbV+1; NoV++) {
190 TbDist2(0,NoV) = RealLast();
191 TbDist2(aNbU+1,NoV) = RealLast();
192 }
193 for (NoU = 1; NoU <= aNbU; NoU++) {
194 TbDist2(NoU,0) = RealLast();
195 TbDist2(NoU,aNbV+1) = RealLast();
196 }
197
198// - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
199 TColStd_Array2OfInteger TbSel(0,aNbU+1,0,aNbV+1);
200 TbSel.Init(0);
201
202/*
203 b.b) Calcul des minima:
204*/
205// - recherche d un minimum sur la grille
206 Standard_Integer Nu, Nv;
207 Standard_Real Dist2;
208 for (NoU = 1; NoU <= aNbU; NoU++) {
209 for (NoV = 1; NoV <= aNbV; NoV++) {
210 if (TbSel(NoU,NoV) == 0) {
211 Dist2 = TbDist2(NoU,NoV);
212 if ((TbDist2(NoU-1,NoV-1) >= Dist2) &&
213 (TbDist2(NoU-1,NoV ) >= Dist2) &&
214 (TbDist2(NoU-1,NoV+1) >= Dist2) &&
215 (TbDist2(NoU ,NoV-1) >= Dist2) &&
216 (TbDist2(NoU ,NoV+1) >= Dist2) &&
217 (TbDist2(NoU+1,NoV-1) >= Dist2) &&
218 (TbDist2(NoU+1,NoV ) >= Dist2) &&
219 (TbDist2(NoU+1,NoV+1) >= Dist2)) {
220
221// - calcul de l extremum sur la surface:
222 UV(1) = U0 + (NoU - 1) * PasU;
223 UV(2) = V0 + (NoV - 1) * PasV;
224
225
226 math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
227
228
229// - mise a jour du tableau TbSel
230 for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
231 for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
232 TbSel(Nu,Nv) = 1;
233 }
234 }
235 }
236 } // if (TbSel(NoU,NoV)
237 } // for (NoV = 1; ...
238 } // for (NoU = 1; ...
239/*
240c- Calcul des maxima:
241 -----------------
242 c.a) Initialisations:
243*/
244// - des 'bords' du tableau TbDist2
245 for (NoV = 0; NoV <= aNbV+1; NoV++) {
246 TbDist2(0,NoV) = RealFirst();
247 TbDist2(aNbU+1,NoV) = RealFirst();
248 }
249 for (NoU = 1; NoU <= aNbU; NoU++) {
250 TbDist2(NoU,0) = RealFirst();
251 TbDist2(NoU,aNbV+1) = RealFirst();
252 }
253
254// - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
255 TbSel.Init(0);
256 /*for (NoU = 0; NoU <= aNbU+1; NoU++) {
257 for (NoV = 0; NoV <= aNbV+1; NoV++) {
258 TbSel(NoU,NoV) = 0;
259 }
260 }*/
261/*
262 c.b) Calcul des maxima:
263*/
264// - recherche d un maximum sur la grille
265 for (NoU = 1; NoU <= aNbU; NoU++) {
266 for (NoV = 1; NoV <= aNbV; NoV++) {
267 if (TbSel(NoU,NoV) == 0) {
268 Dist2 = TbDist2(NoU,NoV);
269 if ((TbDist2(NoU-1,NoV-1) <= Dist2) &&
270 (TbDist2(NoU-1,NoV ) <= Dist2) &&
271 (TbDist2(NoU-1,NoV+1) <= Dist2) &&
272 (TbDist2(NoU ,NoV-1) <= Dist2) &&
273 (TbDist2(NoU ,NoV+1) <= Dist2) &&
274 (TbDist2(NoU+1,NoV-1) <= Dist2) &&
275 (TbDist2(NoU+1,NoV ) <= Dist2) &&
276 (TbDist2(NoU+1,NoV+1) <= Dist2)) {
277
278// - calcul de l extremum sur la surface:
279 UV(1) = U0 + (NoU - 1) * PasU;
280 UV(2) = V0 + (NoV - 1) * PasV;
281
282 math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
283
284// - mise a jour du tableau TbSel
285 for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
286 for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
287 TbSel(Nu,Nv) = 1;
288 }
289 }
290 }
291 } // if (TbSel(NoU,NoV))
292 } // for (NoV = 1; ...)
293 } // for (NoU = 1; ...)
294 myDone = Standard_True;
295}
296//=============================================================================
297
298Standard_Boolean Extrema_GenExtCC::IsDone () const { return myDone; }
299//=============================================================================
300
301Standard_Integer Extrema_GenExtCC::NbExt () const
302{
303 StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::NbExt()")
304 return myF.NbExt();
305}
306//=============================================================================
307
308Standard_Real Extrema_GenExtCC::SquareDistance (const Standard_Integer N) const
309{
310 StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
311 Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
312 return myF.SquareDistance(N);
313}
314//=============================================================================
315
316void Extrema_GenExtCC::Points (const Standard_Integer N,
317 POnC& P1, POnC& P2) const
318{
319 StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
320 Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
321 myF.Points(N,P1,P2);
322}
Open CASCADE Technology repository