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   if there exist N extreme distances between the point and the surface,
 162   so there also exist N extrema between the point and the grid.
 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.Nullify();
 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 // Calculation of 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 of the table of distances (TbDist(0,myusample+1,0,myvsample+1)):
 357    ---------------------------------------------------------------
 358 */
 359
 360 // Parameterisation of the sample
 361
 362
 363 }
 364
 365 void Extrema_GenExtPS::BuildTree()
 366 {
 367   // if tree already exists, assume it is already correctly filled
 368   if ( ! mySphereUBTree.IsNull() )
 369     return;
 370
 371   Standard_Real PasU = myusup - myumin;
 372   Standard_Real PasV = myvsup - myvmin;
 373   Standard_Real U0 = PasU / myusample / 100.;
 374   Standard_Real V0 = PasV / myvsample / 100.;
 375   gp_Pnt P1;
 376   PasU = (PasU - U0) / (myusample - 1);
 377   PasV = (PasV - V0) / (myvsample - 1);
 378   U0 = U0/2. + myumin;
 379   V0 = V0/2. + myvmin;
 380
 381   // Calculation of distances
 382   mySphereUBTree = new Extrema_UBTreeOfSphere;
 383   Extrema_UBTreeFillerOfSphere aFiller(*mySphereUBTree);
 384   Standard_Integer i = 0;
 385   Standard_Real U, V;
 386   mySphereArray = new Bnd_HArray1OfSphere(0, myusample * myvsample);
 387   Standard_Integer NoU, NoV;
 388   for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
 389     for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
 390       P1 = myS->Value(U, V);
 391       Bnd_Sphere aSph(P1.XYZ(), 0/*mytolu < mytolv ? mytolu : mytolv*/, NoU, NoV);
 392       aFiller.Add(i, aSph);
 393       mySphereArray->SetValue( i, aSph );
 394       i++;
 395     }
 396   }
 397   aFiller.Fill();
 398 }
 399
 400 void Extrema_GenExtPS::FindSolution(const gp_Pnt& P, const math_Vector& UV, const Standard_Real PasU, const Standard_Real PasV, const Extrema_ExtFlag f)
 401 {
 402   math_Vector Tol(1,2);
 403   Tol(1) = mytolu;
 404   Tol(2) = mytolv;
 405
 406   math_Vector UVinf(1,2), UVsup(1,2);
 407   UVinf(1) = myumin;
 408   UVinf(2) = myvmin;
 409   UVsup(1) = myusup;
 410   UVsup(2) = myvsup;
 411
 412   math_Vector errors(1,2);
 413   math_Vector root(1, 2);
 414   Standard_Real eps = 1.e-9;
 415   Standard_Integer nbsubsample = 11;
 416
 417   Standard_Integer aNbMaxIter = 100;
 418
 419   gp_Pnt PStart = myS->Value(UV(1), UV(2));
 420   Standard_Real DistStart = P.SquareDistance(PStart);
 421   Standard_Real DistSol = DistStart;
 422
 423   math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
 424   Standard_Boolean ToResolveOnSubgrid = Standard_False;
 425   if (f == Extrema_ExtFlag_MIN)
 426   {
 427     if(S.IsDone())
 428     {
 429       root = S.Root();
 430       myF.Value(root, errors);
 431       gp_Pnt PSol = myS->Value(root(1), root(2));
 432       DistSol = P.SquareDistance(PSol);
 433       if(Abs(errors(1)) > eps || Abs(errors(2)) > eps || DistStart < DistSol)
 434         //try to improve solution on subgrid of sample points
 435         ToResolveOnSubgrid = Standard_True;
 436     }
 437     else
 438       ToResolveOnSubgrid = Standard_True;
 439
 440     if (ToResolveOnSubgrid)
 441     {
 442       Standard_Real u1 = Max(UV(1) - PasU, myumin), u2 = Min(UV(1) + PasU, myusup);
 443       Standard_Real v1 = Max(UV(2) - PasV, myvmin), v2 = Min(UV(2) + PasV, myvsup);
 444
 445       if(u2 - u1 < 2.*PasU) {
 446         if(Abs(u1 - myumin) < 1.e-9) {
 447           u2 = u1 + 2.*PasU;
 448           u2 = Min(u2, myusup);
 449         }
 450         if(Abs(u2 - myusup) < 1.e-9) {
 451           u1 = u2 - 2.*PasU;
 452           u1 = Max(u1, myumin);
 453         }
 454       }
 455
 456       if(v2 - v1 < 2.*PasV) {
 457         if(Abs(v1 - myvmin) < 1.e-9) {
 458           v2 = v1 + 2.*PasV;
 459           v2 = Min(v2, myvsup);
 460         }
 461         if(Abs(v2 - myvsup) < 1.e-9) {
 462           v1 = v2 - 2.*PasV;
 463           v1 = Max(v1, myvmin);
 464         }
 465       }
 466
 467       Standard_Real du = (u2 - u1)/(nbsubsample-1);
 468       Standard_Real dv = (v2 - v1)/(nbsubsample-1);
 469       Standard_Real u, v;
 470       Standard_Real dist;
 471
 472       Standard_Boolean NewSolution = Standard_False;
 473       Standard_Integer Nu, Nv;
 474       for (Nu = 1, u = u1; Nu < nbsubsample; Nu++, u += du) {
 475         for (Nv = 1, v = v1; Nv < nbsubsample; Nv++, v += dv) {
 476           gp_Pnt Puv = myS->Value(u, v);
 477           dist = P.SquareDistance(Puv);
 478
 479           if(dist < DistSol) {
 480             UV(1) = u;
 481             UV(2) = v;
 482             NewSolution = Standard_True;
 483             DistSol = dist;
 484           }
 485         }
 486       }
 487
 488       if(NewSolution) {
 489         //try to precise
 490         math_FunctionSetRoot S (myF,UV,Tol,UVinf,UVsup, aNbMaxIter);
 491       }
 492     } //end of if (ToResolveOnSubgrid)
 493   } //end of if (f == Extrema_ExtFlag_MIN)
 494
 495   myDone = Standard_True;
 496 }
 497
 498 void Extrema_GenExtPS::SetFlag(const Extrema_ExtFlag F)
 499 {
 500   myFlag = F;
 501 }
 502
 503 void Extrema_GenExtPS::SetAlgo(const Extrema_ExtAlgo A)
 504 {
 505   myAlgo = A;
 506 }
 507
 508 void Extrema_GenExtPS::Perform(const gp_Pnt& P)
 509 {
 510   //BuildTree();
 511   myDone = Standard_False;
 512   myF.SetPoint(P);
 513
 514   Standard_Real PasU = myusup - myumin;
 515   Standard_Real PasV = myvsup - myvmin;
 516   Standard_Real U, U0 = PasU / myusample / 100.;
 517   Standard_Real V, V0 = PasV / myvsample / 100.;
 518   PasU = (PasU - U0) / (myusample - 1);
 519   PasV = (PasV - V0) / (myvsample - 1);
 520   U0 = U0/2.+myumin;
 521   V0 = V0/2.+myvmin;
 522
 523   //math_Vector Tol(1, 2);
 524   math_Vector UV(1,2);
 525
 526   if(myAlgo == Extrema_ExtAlgo_Grad)
 527   {
 528     Standard_Integer NoU,NoV;
 529
 530     TColStd_Array2OfReal TheDist(0, myusample+1, 0, myvsample+1);
 531     for ( NoU = 1, U = U0; NoU <= myusample; NoU++, U += PasU) {
 532       for ( NoV = 1, V = V0; NoV <= myvsample; NoV++, V += PasV) {
 533         TheDist(NoU, NoV) = P.SquareDistance(mypoints->Value(NoU, NoV));
 534       }
 535     }
 536
 537     TColStd_Array2OfInteger TbSel(0,myusample+1,0,myvsample+1);
 538     TbSel.Init(0);
 539
 540     Standard_Real Dist;
 541
 542     if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
 543     {
 544       for (NoV = 0; NoV <= myvsample+1; NoV++) {
 545         TheDist(0,NoV) = RealLast();
 546         TheDist(myusample+1,NoV) = RealLast();
 547       }
 548       for (NoU = 1; NoU <= myusample; NoU++) {
 549         TheDist(NoU,0) = RealLast();
 550         TheDist(NoU,myvsample+1) = RealLast();
 551       }
 552       for (NoU = 1; NoU <= myusample; NoU++) {
 553         for (NoV = 1; NoV <= myvsample; NoV++) {
 554           if (TbSel(NoU,NoV) == 0) {
 555             Dist = TheDist(NoU,NoV);
 556             if ((TheDist(NoU-1,NoV-1) >= Dist) &&
 557               (TheDist(NoU-1,NoV  ) >= Dist) &&
 558               (TheDist(NoU-1,NoV+1) >= Dist) &&
 559               (TheDist(NoU  ,NoV-1) >= Dist) &&
 560               (TheDist(NoU  ,NoV+1) >= Dist) &&
 561               (TheDist(NoU+1,NoV-1) >= Dist) &&
 562               (TheDist(NoU+1,NoV  ) >= Dist) &&
 563               (TheDist(NoU+1,NoV+1) >= Dist)) {
 564                 //Create array of UV vectors to calculate min
 565                 UV(1) = U0 + (NoU - 1) * PasU;
 566                 UV(2) = V0 + (NoV - 1) * PasV;
 567                 FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
 568             }
 569           }
 570         }
 571       }
 572     }
 573
 574     if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
 575     {
 576       for (NoV = 0; NoV <= myvsample+1; NoV++) {
 577         TheDist(0,NoV) = RealFirst();
 578         TheDist(myusample+1,NoV) = RealFirst();
 579       }
 580       for (NoU = 1; NoU <= myusample; NoU++) {
 581         TheDist(NoU,0) = RealFirst();
 582         TheDist(NoU,myvsample+1) = RealFirst();
 583       }
 584       for (NoU = 1; NoU <= myusample; NoU++) {
 585         for (NoV = 1; NoV <= myvsample; NoV++) {
 586           if (TbSel(NoU,NoV) == 0) {
 587             Dist = TheDist(NoU,NoV);
 588             if ((TheDist(NoU-1,NoV-1) <= Dist) &&
 589               (TheDist(NoU-1,NoV  ) <= Dist) &&
 590               (TheDist(NoU-1,NoV+1) <= Dist) &&
 591               (TheDist(NoU  ,NoV-1) <= Dist) &&
 592               (TheDist(NoU  ,NoV+1) <= Dist) &&
 593               (TheDist(NoU+1,NoV-1) <= Dist) &&
 594               (TheDist(NoU+1,NoV  ) <= Dist) &&
 595               (TheDist(NoU+1,NoV+1) <= Dist)) {
 596                 //Create array of UV vectors to calculate max
 597                 UV(1) = U0 + (NoU - 1) * PasU;
 598                 UV(2) = V0 + (NoV - 1) * PasV;
 599                 FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
 600             }
 601           }
 602         }
 603       }
 604     }
 605   }
 606   else
 607   {
 608     BuildTree();
 609     if(myFlag == Extrema_ExtFlag_MIN || myFlag == Extrema_ExtFlag_MINMAX)
 610     {
 611       Bnd_Sphere aSol = mySphereArray->Value(0);
 612       Bnd_SphereUBTreeSelectorMin aSelector(mySphereArray, aSol);
 613       //aSelector.SetMaxDist( RealLast() );
 614       aSelector.DefineCheckPoint( P );
 615       Standard_Integer aNbSel = mySphereUBTree->Select( aSelector );
 616       //TODO: check if no solution in binary tree
 617       Bnd_Sphere& aSph = aSelector.Sphere();
 618
 619       UV(1) = U0 + (aSph.U() - 1) * PasU;
 620       UV(2) = V0 + (aSph.V() - 1) * PasV;
 621
 622       FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MIN);
 623     }
 624     if(myFlag == Extrema_ExtFlag_MAX || myFlag == Extrema_ExtFlag_MINMAX)
 625     {
 626       Bnd_Sphere aSol = mySphereArray->Value(0);
 627       Bnd_SphereUBTreeSelectorMax aSelector(mySphereArray, aSol);
 628       //aSelector.SetMaxDist( RealLast() );
 629       aSelector.DefineCheckPoint( P );
 630       Standard_Integer aNbSel = mySphereUBTree->Select( aSelector );
 631       //TODO: check if no solution in binary tree
 632       Bnd_Sphere& aSph = aSelector.Sphere();
 633
 634       UV(1) = U0 + (aSph.U() - 1) * PasU;
 635       UV(2) = V0 + (aSph.V() - 1) * PasV;
 636
 637       FindSolution(P, UV, PasU, PasV, Extrema_ExtFlag_MAX);
 638     }
 639   }
 640 }
 641 //=============================================================================
 642
 643 Standard_Boolean Extrema_GenExtPS::IsDone () const { return myDone; }
 644 //=============================================================================
 645
 646 Standard_Integer Extrema_GenExtPS::NbExt () const
 647 {
 648   if (!IsDone()) { StdFail_NotDone::Raise(); }
 649   return myF.NbExt();
 650 }
 651 //=============================================================================
 652
 653 Standard_Real Extrema_GenExtPS::SquareDistance (const Standard_Integer N) const
 654 {
 655   if (!IsDone()) { StdFail_NotDone::Raise(); }
 656   return myF.SquareDistance(N);
 657 }
 658 //=============================================================================
 659
 660 Extrema_POnSurf Extrema_GenExtPS::Point (const Standard_Integer N) const
 661 {
 662   if (!IsDone()) { StdFail_NotDone::Raise(); }
 663   return myF.Point(N);
 664 }
 665 //=============================================================================