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-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
  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
  31 Extrema_GenExtCC::Extrema_GenExtCC () : myDone (Standard_False)
  32 {
  33 }
  34
  35 Extrema_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
  48 Extrema_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
  64 void 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
  73 void Extrema_GenExtCC::SetTolerance (const Standard_Real Tol)
  74 {
  75   myF.SetTolerance (Tol);
  76 }
  77
  78
  79 //=============================================================================
  80 void Extrema_GenExtCC::Perform ()
  81 /*-----------------------------------------------------------------------------
  82 Fonction:
  83    Recherche de toutes les distances extremales entre les courbes C1 et C2.
  84   a partir de 2 echantillonnages (NbU,NbV).
  85
  86 Methode:
  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
 103 Traitement:
 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 /*
 138 a- 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   const Curve1& aCurve1 = *((Curve1*)(myF.CurvePtr (1)));
 150   const Curve2& aCurve2 = *((Curve1*)(myF.CurvePtr (2)));
 151   const Handle(ArrayOfPnt)& aPntArray1 = aCache1->Points();
 152   const Handle(ArrayOfPnt)& aPntArray2 = aCache2->Points();
 153   Standard_Integer NoU, NoV;
 154   TColStd_Array2OfReal TbDist2(0, aNbU + 1, 0, aNbV + 1);
 155   for (NoU = 1; NoU <= aNbU; NoU++) {
 156     const Pnt& P1 = aPntArray1->Value (NoU);
 157     for (NoV = 1; NoV <= aNbV; NoV++) {
 158       const Pnt& P2 = aPntArray2->Value (NoV);
 159       TbDist2(NoU,NoV) = P1.SquareDistance(P2);
 160     }
 161   }
 162
 163 /*
 164 b- Calcul des minima:
 165    -----------------
 166    b.a) Initialisations:
 167 */
 168 //     - generales
 169   math_Vector Tol(1, 2);
 170   Tol(1) = myF.Tolerance();
 171   Tol(2) = myF.Tolerance();
 172   math_Vector UV(1,2), UVinf(1,2), UVsup(1,2);
 173   UVinf(1) = aCache1->TrimFirstParameter();
 174   UVinf(2) = aCache2->TrimFirstParameter();
 175   UVsup(1) = aCache1->TrimLastParameter();
 176   UVsup(2) = aCache2->TrimLastParameter();
 177
 178   myF.SubIntervalInitialize(UVinf,UVsup);
 179
 180 //     - des 'bords' du tableau TbDist2
 181   for (NoV = 0; NoV <= aNbV+1; NoV++) {
 182     TbDist2(0,NoV) = RealLast();
 183     TbDist2(aNbU+1,NoV) = RealLast();
 184   }
 185   for (NoU = 1; NoU <= aNbU; NoU++) {
 186     TbDist2(NoU,0) = RealLast();
 187     TbDist2(NoU,aNbV+1) = RealLast();
 188   }
 189
 190 //     - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
 191   TColStd_Array2OfInteger TbSel(0,aNbU+1,0,aNbV+1);
 192   TbSel.Init(0);
 193
 194 /*
 195    b.b) Calcul des minima:
 196 */
 197 //     - recherche d un minimum sur la grille
 198   Standard_Integer Nu, Nv;
 199   Standard_Real Dist2;
 200   const Handle(TColStd_HArray1OfReal)& aParamArray1 = aCache1->Parameters();
 201   const Handle(TColStd_HArray1OfReal)& aParamArray2 = aCache2->Parameters();
 202   for (NoU = 1; NoU <= aNbU; NoU++) {
 203     for (NoV = 1; NoV <= aNbV; NoV++) {
 204       if (TbSel(NoU,NoV) == 0) {
 205         Dist2 = TbDist2(NoU,NoV);
 206         if ((TbDist2(NoU-1,NoV-1) >= Dist2) &&
 207             (TbDist2(NoU-1,NoV  ) >= Dist2) &&
 208             (TbDist2(NoU-1,NoV+1) >= Dist2) &&
 209             (TbDist2(NoU  ,NoV-1) >= Dist2) &&
 210             (TbDist2(NoU  ,NoV+1) >= Dist2) &&
 211             (TbDist2(NoU+1,NoV-1) >= Dist2) &&
 212             (TbDist2(NoU+1,NoV  ) >= Dist2) &&
 213             (TbDist2(NoU+1,NoV+1) >= Dist2)) {
 214
 215 //     - calcul de l extremum sur la surface:
 216
 217           UV(1) = aParamArray1->Value(NoU);
 218           UV(2) = aParamArray2->Value(NoV);
 219
 220           math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
 221
 222
 223 //     - mise a jour du tableau TbSel
 224           for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
 225             for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
 226               TbSel(Nu,Nv) = 1;
 227             }
 228           }
 229         }
 230       } // if (TbSel(NoU,NoV)
 231     } // for (NoV = 1; ...
 232   } // for (NoU = 1; ...
 233 /*
 234 c- Calcul des maxima:
 235    -----------------
 236    c.a) Initialisations:
 237 */
 238 //     - des 'bords' du tableau TbDist2
 239   for (NoV = 0; NoV <= aNbV+1; NoV++) {
 240     TbDist2(0,NoV) = RealFirst();
 241     TbDist2(aNbU+1,NoV) = RealFirst();
 242   }
 243   for (NoU = 1; NoU <= aNbU; NoU++) {
 244     TbDist2(NoU,0) = RealFirst();
 245     TbDist2(NoU,aNbV+1) = RealFirst();
 246   }
 247
 248 //     - du tableau TbSel(0,aNbU+1,0,aNbV+1) de selection des points
 249   TbSel.Init(0);
 250   /*for (NoU = 0; NoU <= aNbU+1; NoU++) {
 251     for (NoV = 0; NoV <= aNbV+1; NoV++) {
 252       TbSel(NoU,NoV) = 0;
 253     }
 254   }*/
 255 /*
 256    c.b) Calcul des maxima:
 257 */
 258 //     - recherche d un maximum sur la grille
 259   for (NoU = 1; NoU <= aNbU; NoU++) {
 260     for (NoV = 1; NoV <= aNbV; NoV++) {
 261       if (TbSel(NoU,NoV) == 0) {
 262         Dist2 = TbDist2(NoU,NoV);
 263         if ((TbDist2(NoU-1,NoV-1) <= Dist2) &&
 264             (TbDist2(NoU-1,NoV  ) <= Dist2) &&
 265             (TbDist2(NoU-1,NoV+1) <= Dist2) &&
 266             (TbDist2(NoU  ,NoV-1) <= Dist2) &&
 267             (TbDist2(NoU  ,NoV+1) <= Dist2) &&
 268             (TbDist2(NoU+1,NoV-1) <= Dist2) &&
 269             (TbDist2(NoU+1,NoV  ) <= Dist2) &&
 270             (TbDist2(NoU+1,NoV+1) <= Dist2)) {
 271
 272 //     - calcul de l extremum sur la surface:
 273
 274           UV(1) = aParamArray1->Value(NoU);
 275           UV(2) = aParamArray2->Value(NoV);
 276
 277           math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup);
 278
 279 //     - mise a jour du tableau TbSel
 280           for (Nu = NoU-1; Nu <= NoU+1; Nu++) {
 281             for (Nv = NoV-1; Nv <= NoV+1; Nv++) {
 282               TbSel(Nu,Nv) = 1;
 283             }
 284           }
 285         }
 286       } // if (TbSel(NoU,NoV))
 287     } // for (NoV = 1; ...)
 288   } // for (NoU = 1; ...)
 289   myDone = Standard_True;
 290 }
 291 //=============================================================================
 292
 293 Standard_Boolean Extrema_GenExtCC::IsDone () const { return myDone; }
 294 //=============================================================================
 295
 296 Standard_Integer Extrema_GenExtCC::NbExt () const
 297 {
 298   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::NbExt()")
 299   return myF.NbExt();
 300 }
 301 //=============================================================================
 302
 303 Standard_Real Extrema_GenExtCC::SquareDistance (const Standard_Integer N) const
 304 {
 305   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
 306   Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
 307   return myF.SquareDistance(N);
 308 }
 309 //=============================================================================
 310
 311 void Extrema_GenExtCC::Points (const Standard_Integer N,
 312                             POnC& P1, POnC& P2) const
 313 {
 314   StdFail_NotDone_Raise_if (!myDone, "Extrema_GenExtCC::SquareDistance()")
 315   Standard_OutOfRange_Raise_if ((N < 1 || N > NbExt()), "Extrema_GenExtCC::SquareDistance()")
 316   myF.Points(N,P1,P2);
 317 }