src/Extrema/Extrema_GenExtPS.cxx

   1 // File:        Extrema_GenExtPS.cxx
   2 // Created:     Tue Jul 18 08:21:34 1995
   3 // Author:      Modelistation
   4 //              <model@metrox>
   5
   6 //  Modified by skv - Thu Sep 30 15:21:07 2004 OCC593
   7
   8
   9 #include <Extrema_GenExtPS.ixx>
  10 #include <StdFail_NotDone.hxx>
  11 #include <Standard_OutOfRange.hxx>
  12 #include <TColStd_Array2OfInteger.hxx>
  13 #include <TColStd_Array2OfReal.hxx>
  14 #include <TColgp_Array2OfPnt.hxx>
  15 #include <math_FunctionSetRoot.hxx>
  16 #include <math_Vector.hxx>
  17 #include <math_NewtonFunctionSetRoot.hxx>
  18 #include <GeomAbs_IsoType.hxx>
  19 #include <Bnd_Sphere.hxx>
  20 #include <Extrema_HUBTreeOfSphere.hxx>
  21 #include <Extrema_ExtFlag.hxx>
  22 #include <Bnd_Array1OfSphere.hxx>
  23 #include <Bnd_HArray1OfSphere.hxx>
  24
  25 //IMPLEMENT_HARRAY1(Extrema_HArray1OfSphere)
  26
  27
  28 class Bnd_SphereUBTreeSelector : public Extrema_UBTreeOfSphere::Selector
  29 {
  30  public:
  31
  32   Bnd_SphereUBTreeSelector (const Handle(Bnd_HArray1OfSphere)& theSphereArray,
  33 Bnd_Sphere& theSol)
  34     : myXYZ(0,0,0),
  35       mySphereArray(theSphereArray),
  36       mySol(theSol)
  37   {
  38     //myXYZ = gp_Pnt(0, 0, 0);
  39   }
  40
  41   void DefineCheckPoint( const gp_Pnt& theXYZ )
  42   { myXYZ = theXYZ; }
  43
  44   Bnd_Sphere& Sphere() const
  45   { return mySol; }
  46
  47   virtual Standard_Boolean Reject( const Bnd_Sphere &theBnd ) const = 0;
  48
  49   virtual Standard_Boolean Accept(const Standard_Integer& theObj) = 0;
  50
  51  protected:
  52   gp_Pnt                      myXYZ;
  53   const Handle(Bnd_HArray1OfSphere)& mySphereArray;
  54   Bnd_Sphere&                                            mySol;
  55
  56 };
  57
  58 class Bnd_SphereUBTreeSelectorMin : public Bnd_SphereUBTreeSelector
  59 {
  60 public:
  61   Bnd_SphereUBTreeSelectorMin (const Handle(Bnd_HArray1OfSphere)& theSphereArray,
  62                 Bnd_Sphere& theSol)
  63                 : Bnd_SphereUBTreeSelector(theSphereArray, theSol),
  64                   myMinDist(RealLast())
  65   {}
  66
  67   void SetMinDist( const Standard_Real theMinDist )
  68   { myMinDist = theMinDist; }
  69
  70   Standard_Real MinDist() const
  71   { return myMinDist; }
  72
  73   Standard_Boolean Reject( const Bnd_Sphere &theBnd ) const
  74   {
  75     Bnd_SphereUBTreeSelectorMin* me =
  76       const_cast<Bnd_SphereUBTreeSelectorMin*>(this);
  77     // myMinDist is decreased each time a nearer object is found
  78     return theBnd.IsOut( myXYZ.XYZ(), me->myMinDist );
  79   }
  80
  81   Standard_Boolean Accept(const Standard_Integer&);
  82
  83 private:
  84         Standard_Real   myMinDist;
  85 };
  86
  87 Standard_Boolean Bnd_SphereUBTreeSelectorMin::Accept(const Standard_Integer& theInd)
  88 {
  89   const Bnd_Sphere& aSph = mySphereArray->Value(theInd);
  90   Standard_Real aCurDist;
  91
  92     if ( (aCurDist = aSph.SquareDistance(myXYZ.XYZ())) < mySol.SquareDistance(myXYZ.XYZ()) )
  93     {
  94       mySol = aSph;
  95       if ( aCurDist < myMinDist )
  96         myMinDist = aCurDist;
  97
  98       return Standard_True;
  99     }
 100
 101   return Standard_False;
 102 }
 103
 104 class Bnd_SphereUBTreeSelectorMax : public Bnd_SphereUBTreeSelector
 105 {
 106 public:
 107   Bnd_SphereUBTreeSelectorMax (const Handle(Bnd_HArray1OfSphere)& theSphereArray,
 108                 Bnd_Sphere& theSol)
 109                 : Bnd_SphereUBTreeSelector(theSphereArray, theSol),
 110                   myMaxDist(0)
 111   {}
 112
 113   void SetMaxDist( const Standard_Real theMaxDist )
 114   { myMaxDist = theMaxDist; }
 115
 116   Standard_Real MaxDist() const
 117   { return myMaxDist; }
 118
 119   Standard_Boolean Reject( const Bnd_Sphere &theBnd ) const
 120   {
 121     Bnd_SphereUBTreeSelectorMax* me =
 122       const_cast<Bnd_SphereUBTreeSelectorMax*>(this);
 123     // myMaxDist is decreased each time a nearer object is found
 124     return theBnd.IsOut( myXYZ.XYZ(), me->myMaxDist );
 125   }
 126
 127   Standard_Boolean Accept(const Standard_Integer&);
 128
 129 private:
 130         Standard_Real   myMaxDist;
 131 };
 132
 133 Standard_Boolean Bnd_SphereUBTreeSelectorMax::Accept(const Standard_Integer& theInd)
 134 {
 135   const Bnd_Sphere& aSph = mySphereArray->Value(theInd);
 136   Standard_Real aCurDist;
 137
 138     if ( (aCurDist = aSph.SquareDistance(myXYZ.XYZ())) > mySol.SquareDistance(myXYZ.XYZ()) )
 139     {
 140       mySol = aSph;
 141       if ( aCurDist > myMaxDist )
 142         myMaxDist = aCurDist;
 143
 144       return Standard_True;
 145     }
 146
 147   return Standard_False;
 148 }
 149
 150
 151
 152 //=============================================================================
 153
 154 /*-----------------------------------------------------------------------------
 155 Function:
 156    Find all extremum distances between point P and surface
 157   S using sampling (NbU,NbV).
 158
 159 Method:
 160    The algorithm bases on the hypothesis that sampling is precise enough
 161   pour que, s'il existe N distances extremales entre le point et la surface,
 162   alors il existe aussi N extrema entre le point et la grille.
 163   So, the algorithm consists in starting from extrema of the grid to find the
 164   extrema of the surface.
 165   The extrema are calculated by the algorithm math_FunctionSetRoot with the
 166   following arguments:
 167   - F: Extrema_FuncExtPS created from P and S,
 168   - UV: math_Vector the components which of are parameters of the extremum on the
 169     grid,
 170   - Tol: Min(TolU,TolV), (Prov.:math_FunctionSetRoot does not autorize a vector)
 171   - UVinf: math_Vector the components which of are lower limits of u and v,
 172   - UVsup: math_Vector the components which of are upper limits of u and v.
 173
 174 Processing:
 175   a- Creation of the table of distances (TbDist(0,NbU+1,0,NbV+1)):
 176      The table is expanded at will; lines 0 and NbU+1 and
 177      columns 0 and NbV+1 are initialized at RealFirst() or RealLast()
 178      to simplify the tests carried out at stage b
 179      (there is no need to test if the point is on border of the grid).
 180   b- Calculation of extrema:
 181      First the minimums and then the maximums are found. These 2 procedured
 182      pass in a similar way:
 183   b.a- Initialization:
 184       - 'borders' of table  TbDist (RealLast() in case of minimums
 185         and  RealLast() in case of maximums),
 186       - table TbSel(0,NbU+1,0,NbV+1) of selection of points for
 187         calculation of local extremum (0). When a point will selected,
 188         it will not be selectable, as well as the ajacent points
 189         (8 at least). The corresponding addresses will be set to 1.
 190   b.b- Calculation of minimums (or maximums):
 191        All distances from table TbDist are parsed in a loop:
 192       - search minimum (or maximum) in the grid,
 193       - calculate extremum on the surface,
 194       - update table TbSel.
 195 -----------------------------------------------------------------------------*/
 196
 197 static Standard_Boolean IsoIsDeg  (const Adaptor3d_Surface& S,
 198                                    const Standard_Real      Param,
 199                                    const GeomAbs_IsoType    IT,
 200                                    const Standard_Real      TolMin,
 201                                    const Standard_Real      TolMax)
 202 {
 203     Standard_Real U1=0.,U2=0.,V1=0.,V2=0.,T;
 204     Standard_Boolean Along = Standard_True;
 205     U1 = S.FirstUParameter();
 206     U2 = S.LastUParameter();
 207     V1 = S.FirstVParameter();
 208     V2 = S.LastVParameter();
 209     gp_Vec D1U,D1V;
 210     gp_Pnt P;
 211     Standard_Real Step,D1NormMax;
 212     if (IT == GeomAbs_IsoV)
 213     {
 214       Step = (U2 - U1)/10;
 215       D1NormMax=0.;
 216       for (T=U1;T<=U2;T=T+Step)
 217       {
 218         S.D1(T,Param,P,D1U,D1V);
 219         D1NormMax=Max(D1NormMax,D1U.Magnitude());
 220       }
 221
 222       if (D1NormMax >TolMax || D1NormMax < TolMin )
 223            Along = Standard_False;
 224     }
 225     else
 226     {
 227       Step = (V2 - V1)/10;
 228       D1NormMax=0.;
 229       for (T=V1;T<=V2;T=T+Step)
 230       {
 231         S.D1(Param,T,P,D1U,D1V);
 232         D1NormMax=Max(D1NormMax,D1V.Magnitude());
 233       }
 234
 235       if (D1NormMax >TolMax || D1NormMax < TolMin )
 236            Along = Standard_False;
 237
 238
 239     }
 240     return Along;
 241 }
 242 //----------------------------------------------------------
 243 Extrema_GenExtPS::Extrema_GenExtPS()
 244 {
 245   myDone = Standard_False;
 246   myInit = Standard_False;
 247   myFlag = Extrema_ExtFlag_MINMAX;
 248   myAlgo = Extrema_ExtAlgo_Grad;
 249 }
 250
 251
 252 Extrema_GenExtPS::Extrema_GenExtPS (const gp_Pnt&          P,
 253                               const Adaptor3d_Surface& S,
 254                               const Standard_Integer NbU,
 255                               const Standard_Integer NbV,
 256                               const Standard_Real    TolU,
 257                               const Standard_Real    TolV,
 258                               const Extrema_ExtFlag F,
 259                               const Extrema_ExtAlgo A)
 260   : myF (P,S), myFlag(F), myAlgo(A)
 261 {
 262   Initialize(S, NbU, NbV, TolU, TolV);
 263   Perform(P);
 264 }
 265
 266 Extrema_GenExtPS::Extrema_GenExtPS (const gp_Pnt&          P,
 267                               const Adaptor3d_Surface& S,
 268                               const Standard_Integer NbU,
 269                               const Standard_Integer NbV,
 270                               const Standard_Real    Umin,
 271                               const Standard_Real    Usup,
 272                               const Standard_Real    Vmin,
 273                               const Standard_Real    Vsup,
 274                               const Standard_Real    TolU,
 275                               const Standard_Real    TolV,
 276                               const Extrema_ExtFlag F,
 277                               const Extrema_ExtAlgo A)
 278   : myF (P,S), myFlag(F), myAlgo(A)
 279 {
 280   Initialize(S, NbU, NbV, Umin, Usup, Vmin, Vsup, TolU, TolV);
 281   Perform(P);
 282 }
 283
 284
 285 void Extrema_GenExtPS::Initialize(const Adaptor3d_Surface& S,
 286                                   const Standard_Integer NbU,
 287                                   const Standard_Integer NbV,
 288                                   const Standard_Real    TolU,
 289                                   const Standard_Real    TolV)
 290 {
 291   myumin = S.FirstUParameter();
 292   myusup = S.LastUParameter();
 293   myvmin = S.FirstVParameter();
 294   myvsup = S.LastVParameter();
 295   Initialize(S,NbU,NbV,myumin,myusup,myvmin,myvsup,TolU,TolV);
 296 }
 297
 298
 299 void Extrema_GenExtPS::Initialize(const Adaptor3d_Surface& S,
 300                                   const Standard_Integer NbU,
 301                                   const Standard_Integer NbV,
 302                                   const Standard_Real    Umin,
 303                                   const Standard_Real    Usup,
 304                                   const Standard_Real    Vmin,
 305                                   const Standard_Real    Vsup,
 306                                   const Standard_Real    TolU,
 307                                   const Standard_Real    TolV)
 308 {
 309   myInit = Standard_True;
 310   myS = (Adaptor3d_SurfacePtr)&S;
 311   myusample = NbU;
 312   myvsample = NbV;
 313   mytolu = TolU;
 314   mytolv = TolV;
 315   myumin = Umin;
 316   myusup = Usup;
 317   myvmin = Vmin;
 318   myvsup = Vsup;
 319
 320   if ((myusample < 2) ||
 321       (myvsample < 2)) { Standard_OutOfRange::Raise(); }
 322
 323   myF.Initialize(S);
 324
 325   mySphereUBTree.Clear();
 326
 327   if(myAlgo == Extrema_ExtAlgo_Grad)
 328   {
 329           //If flag was changed and extrema not reinitialized Extrema would fail
 330     mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
 331   Standard_Real PasU = myusup - myumin;
 332   Standard_Real PasV = myvsup - myvmin;
 333   Standard_Real U0 = PasU / myusample / 100.;
 334   Standard_Real V0 = PasV / myvsample / 100.;
 335   gp_Pnt P1;
 336   PasU = (PasU - U0) / (myusample - 1);
 337   PasV = (PasV - V0) / (myvsample - 1);
 338   U0 = U0/2. + myumin;
 339   V0 = V0/2. + myvmin;
 340
 341 // Calcul des distances
 342
 343   Standard_Integer NoU, NoV;
 344   Standard_Real U, V;
 345   for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
 346     for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
 347       P1 = myS->Value(U, V);
 348       mypoints->SetValue(NoU,NoV,P1);
 349     }
 350   }
 351   }
 352
 353   //mypoints = new TColgp_HArray2OfPnt(0,myusample+1,0,myvsample+1);
 354
 355 /*
 356 a- Constitution du tableau des distances (TbDist(0,myusample+1,0,myvsample+1)):
 357    ---------------------------------------------------------------
 358 */
 359
 360 // Parametrage de l echantillon
 361
 362
 363 }
 364
 365 void Extrema_GenExtPS::BuildTree()
 366 {
 367         if(mySphereUBTree.IsEmpty())
 368         {
 369           Standard_Real PasU = myusup - myumin;
 370           Standard_Real PasV = myvsup - myvmin;
 371           Standard_Real U0 = PasU / myusample / 100.;
 372           Standard_Real V0 = PasV / myvsample / 100.;
 373           gp_Pnt P1;
 374           PasU = (PasU - U0) / (myusample - 1);
 375           PasV = (PasV - V0) / (myvsample - 1);
 376           U0 = U0/2. + myumin;
 377           V0 = V0/2. + myvmin;
 378
 379         // Calcul des distances
 380
 381           Standard_Integer NoU, NoV;
 382           //mySphereUBTree.Clear();
 383           Extrema_UBTreeFillerOfSphere aFiller(mySphereUBTree.ChangeTree());
 384           Standard_Integer i = myusample * myvsample;
 385           Standard_Real U, V;
 386           /*for ( NoU = 1; NoU <= myusample; NoU++) {
 387                 for ( NoV = 1; NoV <= myvsample; NoV++) {
 388                         i++;
 389                 }
 390           }*/
 391           mySphereArray = new Bnd_HArray1OfSphere(0, i);
 392           i = 0;
 393           for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
 394                 for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
 395                   P1 = myS->Value(U, V);
 396                   //mypoints->SetValue(NoU, NoV, P1);
 397                   Bnd_Sphere aSph(P1.XYZ(), 0/*mytolu < mytolv ? mytolu : mytolv*/, NoU, NoV);
 398                   aFiller.Add(i, aSph);
 399                   mySphereArray->SetValue( i, aSph );
 400                   i++;
 401                 }
 402           }
 403           aFiller.Fill();
 404         }
 405 }
 406
 407 void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, const Standard_Real PasU, const Standard_Real PasV, const Extrema_ExtFlag f)
 408 {
 409   math_Vector Tol(1,2);
 410   Tol(1) = mytolu;
 411   Tol(2) = mytolv;
 412
 413   math_Vector UVinf(1,2), UVsup(1,2);
 414   UVinf(1) = myumin;
 415   UVinf(2) = myvmin;
 416   UVsup(1) = myusup;
 417   UVsup(2) = myvsup;
 418
 419   math_Vector errors(1,2);
 420   math_Vector root(1, 2);
 421   Standard_Real eps = 1.e-9;
 422   Standard_Integer nbsubsample = 11;
 423
 424   Standard_Integer aNbMaxIter = 100;
 425
 426   if (myF.HasDegIso())
 427     aNbMaxIter = 150;
 428
 429   math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
 430   if(S.IsDone()) {
 431     root = S.Root();
 432     myF.Value(root, errors);
 433     if(f == Extrema_ExtFlag_MIN)
 434     {
 435       if(Abs(errors(1)) > eps || Abs(errors(2)) > eps) {
 436         //try to improve solution on subgrid of sample points
 437         gp_Pnt PSol = myS->Value(root(1), root(2));
 438         Standard_Real DistSol = P.SquareDistance(PSol);
 439
 440         Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup);
 441         Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup);
 442
 443         if(u2 - u1 < 2.*PasU) {
 444           if(Abs(u1 - myumin) < 1.e-9) {
 445             u2 = u1 + 2.*PasU;
 446             u2 = Min(u2, myusup);
 447           }
 448           if(Abs(u2 - myusup) < 1.e-9) {
 449             u1 = u2 - 2.*PasU;
 450             u1 = Max(u1, myumin);
 451           }
 452         }
 453
 454         if(v2 - v1 < 2.*PasV) {
 455           if(Abs(v1 - myvmin) < 1.e-9) {
 456             v2 = v1 + 2.*PasV;
 457             v2 = Min(v2, myvsup);
 458           }
 459           if(Abs(v2 - myvsup) < 1.e-9) {
 460             v1 = v2 - 2.*PasV;
 461             v1 = Max(v1, myvmin);
 462           }
 463         }
 464
 465         Standard_Real du = (u2 - u1)/(nbsubsample-1);
 466         Standard_Real dv = (v2 - v1)/(nbsubsample-1);
 467         Standard_Real u, v;
 468         Standard_Real dist;
 469
 470         Standard_Boolean NewSolution = Standard_False;
 471         Standard_Integer Nu, Nv;
 472         for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
 473           for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
 474             gp_Pnt Puv = myS->Value(u, v);
 475             dist = P.SquareDistance(Puv);
 476
 477             if(dist < DistSol) {
 478               UV(1) = u;
 479               UV(2) = v;
 480               NewSolution = Standard_True;
 481               DistSol = dist;
 482             }
 483           }
 484         }
 485
 486         if(NewSolution) {
 487           //try to precise
 488           math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
 489         }
 490
 491       }
 492     }
 493   }
 494   myDone = Standard_True;
 495 }
 496
 497 void Extrema_GenExtPS::SetFlag(const Extrema_ExtFlag F)
 498 {
 499   myFlag = F;
 500 }
 501
 502 void Extrema_GenExtPS::SetAlgo(const Extrema_ExtAlgo A)
 503 {
 504   myAlgo = A;
 505 }
 506
 507 void Extrema_GenExtPS::Perform(const gp_Pnt& P)
 508 {
 509   //BuildTree();
 510   myDone = Standard_False;
 511   myF.SetPoint(P);
 512
 513   Standard_Real PasU = myusup - myumin;
 514   Standard_Real PasV = myvsup - myvmin;
 515   Standard_Real U, U0 = PasU / myusample / 100.;
 516   Standard_Real V, V0 = PasV / myvsample / 100.;
 517   PasU = (PasU - U0) / (myusample - 1);
 518   PasV = (PasV - V0) / (myvsample - 1);
 519   U0 = U0/2.+myumin;
 520   V0 = V0/2.+myvmin;
 521
 522   //math_Vector Tol(1, 2);
 523   math_Vector UV(1,2);
 524
 525   if(myAlgo == Extrema_ExtAlgo_Grad)
 526   {
 527     Standard_Integer NoU,NoV;
 528
 529     TColStd_Array2OfReal TheDist(0, myusample+1, 0, myvsample+1);
 530     for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
 531       for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
 532         TheDist(NoU, NoV) = P.SquareDistance(mypoints->Value(NoU, NoV));
 533       }
 534     }
 535
 536     TColStd_Array2OfInteger TbSel(0,myusample+1,0,myvsample+1);
 537     TbSel.Init(0);
 538
 539     Standard_Real Dist;
 540
 541     if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
 542     {
 543       for (NoV = 0; NoV <= myvsample+1; NoV++) {
 544         TheDist(0,NoV) = RealLast();
 545         TheDist(myusample+1,NoV) = RealLast();
 546       }
 547       for (NoU = 1; NoU <= myusample; NoU++) {
 548         TheDist(NoU,0) = RealLast();
 549         TheDist(NoU,myvsample+1) = RealLast();
 550       }
 551       for (NoU = 1; NoU <= myusample; NoU++) {
 552         for (NoV = 1; NoV <= myvsample; NoV++) {
 553           if (TbSel(NoU,NoV) == 0) {
 554             Dist = TheDist(NoU,NoV);
 555             if ((TheDist(NoU-1,NoV-1) >= Dist) &&
 556               (TheDist(NoU-1,NoV  ) >= Dist) &&
 557               (TheDist(NoU-1,NoV+1) >= Dist) &&
 558               (TheDist(NoU  ,NoV-1) >= Dist) &&
 559               (TheDist(NoU  ,NoV+1) >= Dist) &&
 560               (TheDist(NoU+1,NoV-1) >= Dist) &&
 561               (TheDist(NoU+1,NoV  ) >= Dist) &&
 562               (TheDist(NoU+1,NoV+1) >= Dist)) {
 563                 //Create array of UV vectors to calculate min
 564                 UV(1) = U0 + (NoU - 1) * PasU;
 565                 UV(2) = V0 + (NoV - 1) * PasV;
 566                 FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
 567             }
 568           }
 569         }
 570       }
 571     }
 572
 573     if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
 574     {
 575       for (NoV = 0; NoV <= myvsample+1; NoV++) {
 576         TheDist(0,NoV) = RealFirst();
 577         TheDist(myusample+1,NoV) = RealFirst();
 578       }
 579       for (NoU = 1; NoU <= myusample; NoU++) {
 580         TheDist(NoU,0) = RealFirst();
 581         TheDist(NoU,myvsample+1) = RealFirst();
 582       }
 583       for (NoU = 1; NoU <= myusample; NoU++) {
 584         for (NoV = 1; NoV <= myvsample; NoV++) {
 585           if (TbSel(NoU,NoV) == 0) {
 586             Dist = TheDist(NoU,NoV);
 587             if ((TheDist(NoU-1,NoV-1) <= Dist) &&
 588               (TheDist(NoU-1,NoV  ) <= Dist) &&
 589               (TheDist(NoU-1,NoV+1) <= Dist) &&
 590               (TheDist(NoU  ,NoV-1) <= Dist) &&
 591               (TheDist(NoU  ,NoV+1) <= Dist) &&
 592               (TheDist(NoU+1,NoV-1) <= Dist) &&
 593               (TheDist(NoU+1,NoV  ) <= Dist) &&
 594               (TheDist(NoU+1,NoV+1) <= Dist)) {
 595                 //Create array of UV vectors to calculate max
 596                 UV(1) = U0 + (NoU - 1) * PasU;
 597                 UV(2) = V0 + (NoV - 1) * PasV;
 598                 FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
 599             }
 600           }
 601         }
 602       }
 603     }
 604   }
 605   else
 606   {
 607     BuildTree();
 608     if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
 609     {
 610       Bnd_Sphere aSol = mySphereArray->Value(0);
 611       Bnd_SphereUBTreeSelectorMin aSelector(mySphereArray, aSol);
 612       //aSelector.SetMaxDist( RealLast() );
 613       aSelector.DefineCheckPoint( P );
 614       Standard_Integer aNbSel = mySphereUBTree.Select( aSelector );
 615       //TODO: check if no solution in binary tree
 616       Bnd_Sphere& aSph = aSelector.Sphere();
 617
 618       UV(1) = U0 + (aSph.U() - 1) * PasU;
 619       UV(2) = V0 + (aSph.V() - 1) * PasV;
 620
 621       FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
 622     }
 623     if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
 624     {
 625       Bnd_Sphere aSol = mySphereArray->Value(0);
 626       Bnd_SphereUBTreeSelectorMax aSelector(mySphereArray, aSol);
 627       //aSelector.SetMaxDist( RealLast() );
 628       aSelector.DefineCheckPoint( P );
 629       Standard_Integer aNbSel = mySphereUBTree.Select( aSelector );
 630       //TODO: check if no solution in binary tree
 631       Bnd_Sphere& aSph = aSelector.Sphere();
 632
 633       UV(1) = U0 + (aSph.U() - 1) * PasU;
 634       UV(2) = V0 + (aSph.V() - 1) * PasV;
 635
 636       FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
 637     }
 638   }
 639 }
 640 //=============================================================================
 641
 642 Standard_Boolean Extrema_GenExtPS::IsDone () const { return myDone; }
 643 //=============================================================================
 644
 645 Standard_Integer Extrema_GenExtPS::NbExt () const
 646 {
 647   if (!IsDone()) { StdFail_NotDone::Raise(); }
 648   return myF.NbExt();
 649 }
 650 //=============================================================================
 651
 652 Standard_Real Extrema_GenExtPS::SquareDistance (const Standard_Integer N) const
 653 {
 654   if (!IsDone()) { StdFail_NotDone::Raise(); }
 655   return myF.SquareDistance(N);
 656 }
 657 //=============================================================================
 658
 659 Extrema_POnSurf Extrema_GenExtPS::Point (const Standard_Integer N) const
 660 {
 661   if (!IsDone()) { StdFail_NotDone::Raise(); }
 662   return myF.Point(N);
 663 }
 664 //=============================================================================