0026593: Coding rules - revert compatibility of NCollection_CellFilter constructor...
[occt.git] / src / Extrema / Extrema_GenExtCC.gxx
1 // Created on: 1995-07-18
2 // Created by: Modelistation
3 // Copyright (c) 1995-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
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
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.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 #include <algorithm>
18
19 #include <Extrema_GlobOptFuncCC.hxx>
20 #include <math_GlobOptMin.hxx>
21 #include <Standard_NullObject.hxx>
22 #include <Standard_OutOfRange.hxx>
23 #include <StdFail_NotDone.hxx>
24 #include <TColStd_Array1OfReal.hxx>
25 #include <Precision.hxx>
26 #include <NCollection_Vector.hxx>
27 #include <NCollection_CellFilter.hxx>
28
29 // Comparator, used in std::sort.
30 static Standard_Boolean comp(const gp_XY& theA,
31                              const gp_XY& theB)
32 {
33   if (theA.X() < theB.X())
34   {
35     return Standard_True;
36   }
37   else
38   {
39     if (theA.X() == theB.X())
40     {
41       if (theA.Y() <= theB.Y())
42         return Standard_True;
43     }
44   }
45
46   return Standard_False;
47 }
48
49 class Extrema_CCPointsInspector : public NCollection_CellFilter_InspectorXY
50 {
51 public:
52   typedef gp_XY Target;
53   //! Constructor; remembers the tolerance
54   Extrema_CCPointsInspector (const Standard_Real theTol)
55   {
56     myTol = theTol * theTol;
57     myIsFind = Standard_False;
58   }
59
60   void ClearFind()
61   {
62     myIsFind = Standard_False;
63   }
64
65   Standard_Boolean isFind()
66   {
67     return myIsFind;
68   }
69
70   //! Set current point to search for coincidence
71   void SetCurrent (const gp_XY& theCurPnt) 
72   { 
73     myCurrent = theCurPnt;
74   }
75
76   //! Implementation of inspection method
77   NCollection_CellFilter_Action Inspect (const Target& theObject)
78   {
79     gp_XY aPt = myCurrent.Subtracted(theObject);
80     const Standard_Real aSQDist = aPt.SquareModulus();
81     if(aSQDist < myTol)
82     {
83       myIsFind = Standard_True;
84     }
85
86     return CellFilter_Keep;
87   }
88
89 private:
90   Standard_Real myTol;
91   gp_XY myCurrent;
92   Standard_Boolean myIsFind;
93 };
94
95 //=======================================================================
96 //function : Extrema_GenExtCC
97 //purpose  : 
98 //=======================================================================
99 Extrema_GenExtCC::Extrema_GenExtCC()
100 : myParallel(Standard_False),
101   myCurveMinTol(Precision::PConfusion()),
102   myLowBorder(1,2),
103   myUppBorder(1,2),
104   myDone(Standard_False)
105 {
106   myC[0] = myC[1] = 0;
107 }
108
109 //=======================================================================
110 //function : Extrema_GenExtCC
111 //purpose  : 
112 //=======================================================================
113 Extrema_GenExtCC::Extrema_GenExtCC(const Curve1& C1,
114                                    const Curve2& C2)
115 : myParallel(Standard_False),
116   myCurveMinTol(Precision::PConfusion()),
117   myLowBorder(1,2),
118   myUppBorder(1,2),
119   myDone(Standard_False)
120 {
121   myC[0] = (Standard_Address)&C1;
122   myC[1] = (Standard_Address)&C2;
123   myLowBorder(1) = C1.FirstParameter();
124   myLowBorder(2) = C2.FirstParameter();
125   myUppBorder(1) = C1.LastParameter();
126   myUppBorder(2) = C2.LastParameter();
127 }
128
129 //=======================================================================
130 //function : Extrema_GenExtCC
131 //purpose  : 
132 //=======================================================================
133 Extrema_GenExtCC::Extrema_GenExtCC(const Curve1& C1,
134                                    const Curve2& C2,
135                                    const Standard_Real Uinf,
136                                    const Standard_Real Usup,
137                                    const Standard_Real Vinf,
138                                    const Standard_Real Vsup)
139 : myParallel(Standard_False),
140   myCurveMinTol(Precision::PConfusion()),
141   myLowBorder(1,2),
142   myUppBorder(1,2),
143   myDone(Standard_False)
144 {
145   myC[0] = (Standard_Address)&C1;
146   myC[1] = (Standard_Address)&C2;
147   myLowBorder(1) = Uinf;
148   myLowBorder(2) = Vinf;
149   myUppBorder(1) = Usup;
150   myUppBorder(2) = Vsup;
151 }
152
153 //=======================================================================
154 //function : SetParams
155 //purpose  : 
156 //=======================================================================
157 void Extrema_GenExtCC::SetParams(const Curve1& C1,
158                                  const Curve2& C2,
159                                  const Standard_Real Uinf,
160                                  const Standard_Real Usup,
161                                  const Standard_Real Vinf,
162                                  const Standard_Real Vsup)
163 {
164   myC[0] = (Standard_Address)&C1;
165   myC[1] = (Standard_Address)&C2;
166   myLowBorder(1) = Uinf;
167   myLowBorder(2) = Vinf;
168   myUppBorder(1) = Usup;
169   myUppBorder(2) = Vsup;
170 }
171
172 //=======================================================================
173 //function : SetTolerance
174 //purpose  : 
175 //=======================================================================
176 void Extrema_GenExtCC::SetTolerance(Standard_Real theTol)
177 {
178   myCurveMinTol = theTol;
179 }
180
181 //=======================================================================
182 //function : Perform
183 //purpose  : 
184 //=======================================================================
185 void Extrema_GenExtCC::Perform()
186 {
187   myDone = Standard_False;
188   myParallel = Standard_False;
189
190   Curve1 &C1 = *(Curve1*)myC[0];
191   Curve2 &C2 = *(Curve2*)myC[1];
192
193   Standard_Integer aNbInter[2];
194   aNbInter[0] = C1.NbIntervals(GeomAbs_C2);
195   aNbInter[1] = C2.NbIntervals(GeomAbs_C2);
196   TColStd_Array1OfReal anIntervals1(1, aNbInter[0] + 1);
197   TColStd_Array1OfReal anIntervals2(1, aNbInter[1] + 1);
198   C1.Intervals(anIntervals1, GeomAbs_C2);
199   C2.Intervals(anIntervals2, GeomAbs_C2);
200
201   Extrema_GlobOptFuncCCC2 aFunc (C1, C2);
202   math_GlobOptMin aFinder(&aFunc, myLowBorder, myUppBorder);
203   Standard_Real aDiscTol = 1.0e-2;
204   Standard_Real aValueTol = 1.0e-2;
205   Standard_Real aSameTol = myCurveMinTol / (aDiscTol);
206   aFinder.SetTol(aDiscTol, aSameTol);
207
208   // Size computed to have cell index inside of int32 value.
209   const Standard_Real aCellSize = Max(anIntervals1.Upper() - anIntervals1.Lower(),
210                                       anIntervals2.Upper() - anIntervals2.Lower())
211                                   * Precision::PConfusion() / (2.0 * Sqrt(2.0));
212   Extrema_CCPointsInspector anInspector(Precision::PConfusion());
213   NCollection_CellFilter<Extrema_CCPointsInspector> aFilter(aCellSize);
214   NCollection_Vector<gp_XY> aPnts;
215
216   Standard_Integer i,j,k;
217   math_Vector aFirstBorderInterval(1,2);
218   math_Vector aSecondBorderInterval(1,2);
219   Standard_Real aF = RealLast(); // Best functional value.
220   Standard_Real aCurrF = RealLast(); // Current functional value computed on current interval.
221   for(i = 1; i <= aNbInter[0]; i++)
222   {
223     for(j = 1; j <= aNbInter[1]; j++)
224     {
225       aFirstBorderInterval(1) = anIntervals1(i);
226       aFirstBorderInterval(2) = anIntervals2(j); 
227       aSecondBorderInterval(1) = anIntervals1(i + 1);
228       aSecondBorderInterval(2) = anIntervals2(j + 1);
229
230       aFinder.SetLocalParams(aFirstBorderInterval, aSecondBorderInterval);
231       aFinder.Perform();
232
233       // Check that solution found on current interval is not worse than previous.
234       aCurrF = aFinder.GetF();
235       if (aCurrF >= aF + aSameTol * aValueTol)
236       {
237         continue;
238       }
239
240       // Clean previously computed solution if current one is better.
241       if (aCurrF > aF - aSameTol * aValueTol)
242       {
243         if (aCurrF < aF)
244           aF = aCurrF;
245       }
246       else
247       {
248         aF = aCurrF;
249         aFilter.Reset(aCellSize);
250         aPnts.Clear();
251       }
252
253       // Save found solutions avoiding repetitions.
254       math_Vector sol(1,2);
255       for(k = 1; k <= aFinder.NbExtrema(); k++)
256       {
257         aFinder.Points(k, sol);
258         gp_XY aPnt2d(sol(1), sol(2));
259
260         gp_XY aXYmin = anInspector.Shift(aPnt2d, -aCellSize);
261         gp_XY aXYmax = anInspector.Shift(aPnt2d,  aCellSize);
262
263         anInspector.ClearFind();
264         anInspector.SetCurrent(aPnt2d);
265         aFilter.Inspect(aXYmin, aXYmax, anInspector);
266         if (!anInspector.isFind())
267         {
268           // Point is out of close cells, add new one.
269           aFilter.Add(aPnt2d, aPnt2d);
270           aPnts.Append(gp_XY(sol(1), sol(2)));
271         }
272       }
273     }
274   }
275
276   if (aPnts.Size() == 0)
277   {
278     // No solutions.
279     myDone = Standard_False;
280     return;
281   }
282
283   // Check for infinity solutions case, for this:
284   // Sort points lexicographically and check midpoint between each two neighboring points.
285   // If all midpoints functional value is acceptable
286   // then set myParallel flag to true and return one soulution.
287   std::sort(aPnts.begin(), aPnts.end(), comp);
288   Standard_Boolean isParallel = Standard_False;
289   Standard_Real aVal = 0.0;
290   math_Vector aVec(1,2, 0.0);
291
292   // Avoid mark parallel case when have duplicates out of tolerance.
293   // Bad conditioned task: bug25635_1, bug23706_10, bug23706_13.
294   const Standard_Integer aMinNbInfSol = 100;
295   if (aPnts.Size() >= aMinNbInfSol)
296   {
297     isParallel = Standard_True;
298     for(Standard_Integer anIdx = aPnts.Lower(); anIdx <= aPnts.Upper() - 1; anIdx++)
299     {
300       const gp_XY& aCurrent = aPnts(anIdx);
301       const gp_XY& aNext    = aPnts(anIdx + 1);
302
303       aVec(1) = (aCurrent.X() + aNext.X()) * 0.5;
304       aVec(2) = (aCurrent.Y() + aNext.Y()) * 0.5;
305
306       aFunc.Value(aVec, aVal);
307
308       if (Abs(aVal - aF) > Precision::Confusion())
309       {
310         isParallel = Standard_False;
311         break;
312       }
313     }
314   }
315
316   if (isParallel)
317   {
318     const gp_XY& aCurrent = aPnts.First();
319     myPoints1.Append(aCurrent.X());
320     myPoints2.Append(aCurrent.Y());
321     myParallel = Standard_True;
322   }
323   else
324   {
325     for(Standard_Integer anIdx = aPnts.Lower(); anIdx <= aPnts.Upper(); anIdx++)
326     {
327       const gp_XY& aCurrent = aPnts(anIdx);
328       myPoints1.Append(aCurrent.X());
329       myPoints2.Append(aCurrent.Y());
330     }
331   }
332
333   myDone = Standard_True;
334 }
335
336 //=======================================================================
337 //function : IsDone
338 //purpose  : 
339 //=======================================================================
340 Standard_Boolean Extrema_GenExtCC::IsDone() const 
341 {
342   return myDone; 
343 }
344
345 //=======================================================================
346 //function : IsParallel
347 //purpose  : 
348 //=======================================================================
349 Standard_Boolean Extrema_GenExtCC::IsParallel() const 
350 {
351   return myParallel; 
352 }
353
354 //=======================================================================
355 //function : NbExt
356 //purpose  : 
357 //=======================================================================
358 Standard_Integer Extrema_GenExtCC::NbExt() const
359 {
360   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::NbExt()")
361
362   return myPoints1.Length();
363 }
364
365 //=======================================================================
366 //function : SquareDistance
367 //purpose  : 
368 //=======================================================================
369 Standard_Real Extrema_GenExtCC::SquareDistance(const Standard_Integer N) const
370 {
371   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
372   Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
373
374   return Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)).SquareDistance(Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));
375 }
376
377 //=======================================================================
378 //function : Points
379 //purpose  : 
380 //=======================================================================
381 void Extrema_GenExtCC::Points(const Standard_Integer N,
382                               POnC& P1,
383                               POnC& P2) const
384 {
385   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::Points()")
386   Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::Points()")
387
388   P1.SetValues(myPoints1(N), Tool1::Value(*((Curve1*)myC[0]), myPoints1(N)));
389   P2.SetValues(myPoints2(N), Tool2::Value(*((Curve2*)myC[1]), myPoints2(N)));
390 }