src/GeomAlgoAPI/GeomAlgoAPI_PointCloudOnFace.cpp

   1 // Copyright (C) 2014-2022  CEA/DEN, EDF R&D
   2 //
   3 // This library is free software; you can redistribute it and/or
   4 // modify it under the terms of the GNU Lesser General Public
   5 // License as published by the Free Software Foundation; either
   6 // version 2.1 of the License, or (at your option) any later version.
   7 //
   8 // This library is distributed in the hope that it will be useful,
   9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  11 // Lesser General Public License for more details.
  12 //
  13 // You should have received a copy of the GNU Lesser General Public
  14 // License along with this library; if not, write to the Free Software
  15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  16 //
  17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  18 //
  19
  20 #include "GeomAlgoAPI_PointCloudOnFace.h"
  21
  22 #include <TopoDS_Shape.hxx>
  23 #include <TopoDS_Face.hxx>
  24 #include <TopoDS.hxx>
  25 #include <BRep_Tool.hxx>
  26 #include <BRep_Builder.hxx>
  27 #include <BRepTools.hxx>
  28 #include <BRepTools_WireExplorer.hxx>
  29 #include <BRepTopAdaptor_FClass2d.hxx>
  30 #include <BRepLib_MakeVertex.hxx>
  31 #include <GeomLib.hxx>
  32 #include <gp_Ax3.hxx>
  33 #include <Geom_Plane.hxx>
  34 #include <gp_Pln.hxx>
  35 #include <TopoDS_Iterator.hxx>
  36 #include <TopExp_Explorer.hxx>
  37 #include <TopExp.hxx>
  38 #include <GProp_GProps.hxx>
  39 #include <BRepGProp.hxx>
  40 #include <ShapeAnalysis.hxx>
  41 #include <ShapeAnalysis_Surface.hxx>
  42 #include <ShapeUpgrade_UnifySameDomain.hxx>
  43 #include <ShapeUpgrade_ShapeDivideArea.hxx>
  44 #include <BRepAdaptor_Surface.hxx>
  45 #include <BRepAdaptor_Curve2d.hxx>
  46 #include <BRepBndLib.hxx>
  47 #include <GeomLProp_SLProps.hxx>
  48 #include <BRepBuilderAPI_MakeEdge.hxx>
  49 #include <TopTools_DataMapOfShapeReal.hxx>
  50
  51 #include <Standard_Version.hxx>
  52 // code from KERNEL_SRC/src/Basics/Basics_OCCTVersion.hxx
  53 #ifdef OCC_VERSION_SERVICEPACK
  54 #  define OCC_VERSION_LARGE (OCC_VERSION_MAJOR << 24 | OCC_VERSION_MINOR << 16 | OCC_VERSION_MAINTENANCE << 8 | OCC_VERSION_SERVICEPACK)
  55 #else
  56 #  ifdef OCC_VERSION_DEVELOPMENT
  57 #    define OCC_VERSION_LARGE ((OCC_VERSION_MAJOR << 24 | OCC_VERSION_MINOR << 16 | OCC_VERSION_MAINTENANCE << 8)-1)
  58 #  else
  59 #    define OCC_VERSION_LARGE (OCC_VERSION_MAJOR << 24 | OCC_VERSION_MINOR << 16 | OCC_VERSION_MAINTENANCE << 8)
  60 #  endif
  61 #endif
  62
  63 #include <algorithm>
  64
  65 Standard_Boolean comp(const std::pair<TopoDS_Shape, Standard_Real>& theA,
  66                       const std::pair<TopoDS_Shape, Standard_Real>& theB)
  67 {
  68   return (theA.second < theB.second);
  69 }
  70
  71 Standard_Boolean IsUiso (const TopoDS_Edge& theEdge,
  72                          const TopoDS_Face& theFace)
  73 {
  74   BRepAdaptor_Curve2d aBAcurve2d (theEdge, theFace);
  75   gp_Pnt2d aP2d;
  76   gp_Vec2d aVec;
  77   aBAcurve2d.D1 (aBAcurve2d.FirstParameter(), aP2d, aVec);
  78   return (Abs(aVec.Y()) > Abs(aVec.X()));
  79 }
  80
  81 void CorrectShell (const TopoDS_Shape& theShell,
  82                    const TopoDS_Face&  theFace)
  83 {
  84   BRepAdaptor_Surface aBAsurf (theFace, Standard_False);
  85   GeomAbs_SurfaceType aType = aBAsurf.GetType();
  86   if (aType <= GeomAbs_Torus) //elementary surfaces
  87     return;
  88
  89   TopLoc_Location anInputLoc;
  90   const Handle(Geom_Surface)& anInputSurf = BRep_Tool::Surface (theFace, anInputLoc);
  91
  92   BRep_Builder aBB;
  93
  94   TopoDS_Iterator anIter (theShell);
  95   for (; anIter.More(); anIter.Next())
  96   {
  97     const TopoDS_Face& aFace = TopoDS::Face (anIter.Value());
  98     TopLoc_Location aLoc;
  99     const Handle(Geom_Surface)& aSurf = BRep_Tool::Surface (aFace, aLoc);
 100     if (aSurf == anInputSurf)
 101       continue;
 102
 103     TopExp_Explorer anExplo (aFace, TopAbs_EDGE);
 104     for (; anExplo.More(); anExplo.Next())
 105     {
 106       const TopoDS_Edge& anEdge = TopoDS::Edge (anExplo.Current());
 107       Standard_Real aFirst, aLast;
 108       Handle(Geom2d_Curve) aPCurve = BRep_Tool::CurveOnSurface (anEdge, aFace, aFirst, aLast);
 109       aBB.UpdateEdge (anEdge, aPCurve, anInputSurf, anInputLoc, 0.);
 110     }
 111     Standard_Real aTol = BRep_Tool::Tolerance (aFace);
 112     aBB.UpdateFace (aFace, anInputSurf, anInputLoc, aTol);
 113   }
 114 }
 115
 116 gp_Pnt GetMidPnt2d(const TopoDS_Face&     theFace,
 117                    const Standard_Boolean theIsNaturalRestrictions)
 118 {
 119   gp_Pnt aResPnt;
 120
 121   if (theIsNaturalRestrictions)
 122   {
 123     BRepAdaptor_Surface aBAsurf (theFace);
 124     Standard_Real aUmin, aUmax, aVmin, aVmax;
 125     aUmin = aBAsurf.FirstUParameter();
 126     aUmax = aBAsurf.LastUParameter();
 127     aVmin = aBAsurf.FirstVParameter();
 128     aVmax = aBAsurf.LastVParameter();
 129     aResPnt = aBAsurf.Value ((aUmin + aUmax)/2, (aVmin + aVmax)/2);
 130   }
 131   else
 132   {
 133     const Standard_Integer aNbSamples = 4;
 134     TopoDS_Wire aWire = BRepTools::OuterWire (theFace);
 135     TopTools_IndexedMapOfShape aEmap;
 136     TopExp::MapShapes (aWire, TopAbs_EDGE, aEmap);
 137     Standard_Integer aNbPointsOnContour = aNbSamples * aEmap.Extent();
 138     TColgp_Array1OfPnt anArray (1, aNbPointsOnContour);
 139
 140     BRepTools_WireExplorer aWexp (aWire, theFace);
 141     Standard_Integer anInd = 0;
 142     for (; aWexp.More(); aWexp.Next())
 143     {
 144       const TopoDS_Edge& anEdge = aWexp.Current();
 145       BRepAdaptor_Curve2d aBAcurve2d (anEdge, theFace);
 146       Standard_Real aDelta = (aBAcurve2d.LastParameter() - aBAcurve2d.FirstParameter())/aNbSamples;
 147       for (Standard_Integer ii = 0; ii < aNbSamples; ii++)
 148       {
 149         Standard_Real aParam = aBAcurve2d.FirstParameter() + ii * aDelta;
 150         gp_Pnt2d aP2d = aBAcurve2d.Value (aParam);
 151         gp_Pnt aPnt (aP2d.X(), aP2d.Y(), 0.);
 152         anArray (++anInd) = aPnt;
 153       }
 154     }
 155
 156     gp_Ax2 anAxis;
 157     Standard_Boolean anIsSingular;
 158     GeomLib::AxeOfInertia (anArray, anAxis, anIsSingular);
 159     gp_Pnt aBaryCentre = anAxis.Location();
 160     gp_Pnt2d aCentre2d (aBaryCentre.X(), aBaryCentre.Y());
 161     BRepTopAdaptor_FClass2d aClassifier (theFace, Precision::Confusion());
 162     BRepAdaptor_Surface aBAsurf (theFace, Standard_False);
 163
 164     TopAbs_State aStatus = aClassifier.Perform (aCentre2d);
 165     gp_Pnt2d aP2d = aCentre2d;
 166     Standard_Integer anIndVertex = 0;
 167     const Standard_Integer aNbIter = 10;
 168     while (aStatus != TopAbs_IN && anIndVertex < aNbPointsOnContour)
 169     {
 170       gp_Pnt aVertexPnt = anArray (anIndVertex+1);
 171       gp_Pnt2d aVertexP2d (aVertexPnt.X(), aVertexPnt.Y());
 172       TColgp_SequenceOfPnt2d aPseq;
 173       aPseq.Append (aCentre2d);
 174       aPseq.Append (aVertexP2d);
 175       for (Standard_Integer ii = 1; ii <= aNbIter; ii++)
 176       {
 177         for (Standard_Integer jj = 1; jj < aPseq.Length(); jj++)
 178         {
 179           aP2d.SetXY ((aPseq(jj).XY() + aPseq(jj+1).XY())/2);
 180           aStatus = aClassifier.Perform (aP2d);
 181           if (aStatus == TopAbs_IN)
 182             break;
 183           else
 184           {
 185             aPseq.InsertAfter (jj, aP2d);
 186             jj++;
 187           }
 188         }
 189         if (aStatus == TopAbs_IN)
 190           break;
 191       }
 192       anIndVertex += aNbSamples;
 193     }
 194     aResPnt = aBAsurf.Value (aP2d.X(), aP2d.Y());
 195   } //case of complex boundaries
 196
 197   return aResPnt;
 198 }
 199
 200 void ModifyFacesForGlobalResult(const TopoDS_Face&     theInputFace,
 201                                 const Standard_Real    theAverageArea,
 202                                 const Standard_Boolean theIsToAddFaces,
 203                                 Standard_Integer&      theNbExtremalFaces,
 204                                 TopTools_MapOfShape&   theExtremalFaces,
 205                                 const std::vector<std::pair<TopoDS_Shape, Standard_Real>> theFacesAndAreas,
 206                                 const TopTools_DataMapOfShapeReal& theFaceAreaMap,
 207                                 const TopTools_IndexedDataMapOfShapeListOfShape& theEFmap,
 208                                 TopoDS_Shape&          theRes,
 209                                 TopoDS_Shape&          theGlobalRes,
 210                                 TopTools_MapOfShape&   theRemovedFaces)
 211 {
 212   BRep_Builder aBB;
 213   const Standard_Integer aNbFaces = (Standard_Integer) theFacesAndAreas.size();
 214
 215   const Standard_Integer aDiff = theNbExtremalFaces - theRemovedFaces.Extent();
 216
 217   Standard_Integer aSum = 0;
 218   while (aSum < aDiff) //global loop
 219   {
 220     Standard_Integer aNbFacesDone = 0, aNbFacesInTape = 0;
 221     TopoDS_Face aStartFace;
 222
 223     Standard_Integer aStartIndex = (theIsToAddFaces)? aNbFaces-1 : 0;
 224     Standard_Integer anEndIndex  = (theIsToAddFaces)? 0 : aNbFaces-1;
 225     Standard_Integer aStep = (theIsToAddFaces)? -1 : 1;
 226
 227     for (Standard_Integer ii = aStartIndex; ii != anEndIndex; ii += aStep)
 228     {
 229       const TopoDS_Face& aFace = TopoDS::Face (theFacesAndAreas[ii].first);
 230       if (!theRemovedFaces.Contains(aFace))
 231       {
 232         aStartFace = aFace;
 233         break;
 234       }
 235     }
 236     if (aStartFace.IsNull())
 237       break;
 238
 239     theRemovedFaces.Add (aStartFace);
 240
 241     TopoDS_Edge aCommonEdge;
 242     TopoDS_Face aNextFace;
 243     Standard_Real anExtremalArea = (theIsToAddFaces)? 0. : Precision::Infinite();
 244     for (TopExp_Explorer anExplo(aStartFace, TopAbs_EDGE); anExplo.More(); anExplo.Next())
 245     {
 246       const TopoDS_Edge& anEdge = TopoDS::Edge (anExplo.Current());
 247       const TopTools_ListOfShape& aFaceList = theEFmap.FindFromKey (anEdge);
 248       TopTools_ListIteratorOfListOfShape anItl (aFaceList);
 249       for (; anItl.More(); anItl.Next())
 250       {
 251         const TopoDS_Face& aFace = TopoDS::Face (anItl.Value());
 252         if (aFace.IsSame (aStartFace) ||
 253             theRemovedFaces.Contains(aFace))
 254           continue;
 255         Standard_Real anArea = theFaceAreaMap(aFace);
 256         Standard_Boolean anIsToExchange = (theIsToAddFaces)? (anArea > anExtremalArea) : (anArea < anExtremalArea);
 257         if (anIsToExchange)
 258         {
 259           anExtremalArea = anArea;
 260           aCommonEdge = anEdge;
 261           aNextFace = aFace;
 262         }
 263       }
 264     }
 265     if (aCommonEdge.IsNull()) //all adjacent faces are already removed
 266     {
 267       theExtremalFaces.Add (theFacesAndAreas[theNbExtremalFaces].first);
 268       theNbExtremalFaces++;
 269       continue;
 270     }
 271
 272     //Start filling the shell
 273     aBB.Remove (theRes, aStartFace);
 274     aNbFacesDone++;
 275     TopoDS_Shell aShell;
 276     aBB.MakeShell (aShell);
 277     Standard_Real anAreaOfTape = 0.;
 278     aBB.Add (aShell, aStartFace);
 279     aNbFacesInTape++;
 280     anAreaOfTape += theFaceAreaMap (aStartFace);
 281
 282     Standard_Boolean anIsUiso = IsUiso (aCommonEdge, aStartFace);
 283     //Find another faces on this level
 284     TopoDS_Face aCurrentFace = aNextFace;
 285     TopoDS_Edge aCurrentEdge = aCommonEdge;
 286     Standard_Boolean anIsFirstDirection = Standard_True;
 287     aBB.Remove (theRes, aCurrentFace);
 288     theRemovedFaces.Add (aCurrentFace);
 289     if (theExtremalFaces.Contains (aCurrentFace))
 290     {
 291       aNbFacesDone++;
 292     }
 293     aBB.Add (aShell, aCurrentFace);
 294     aNbFacesInTape++;
 295     anAreaOfTape += theFaceAreaMap (aCurrentFace);
 296     Standard_Boolean anIsRound = Standard_False;
 297     for (;;) //local loop
 298     {
 299       TopoDS_Edge aNextEdge;
 300       for (TopExp_Explorer anExplo(aCurrentFace, TopAbs_EDGE); anExplo.More(); anExplo.Next())
 301       {
 302         const TopoDS_Edge& anEdge = TopoDS::Edge (anExplo.Current());
 303         if (anEdge.IsSame (aCurrentEdge))
 304           continue;
 305         const TopTools_ListOfShape& aFaceList = theEFmap.FindFromKey (anEdge);
 306         TopTools_ListIteratorOfListOfShape anItl (aFaceList);
 307         for (; anItl.More(); anItl.Next())
 308         {
 309           const TopoDS_Face& aFace = TopoDS::Face (anItl.Value());
 310           if (aFace.IsSame (aCurrentFace))
 311             continue;
 312           if (aFace.IsSame (aStartFace))
 313           {
 314             anIsRound = Standard_True;
 315             break;
 316           }
 317           if (theRemovedFaces.Contains(aFace))
 318             continue;
 319           if (anIsUiso == IsUiso (anEdge, aFace))
 320           {
 321             aNextEdge = anEdge;
 322             aNextFace = aFace;
 323             break;
 324           }
 325         }
 326         if (anIsRound || !aNextEdge.IsNull())
 327           break;
 328       }
 329       if (anIsRound) //round tape: returned to start face
 330         break;
 331       if (aNextEdge.IsNull())
 332       {
 333         if (anIsFirstDirection)
 334         {
 335           aCurrentFace = aStartFace;
 336           aCurrentEdge = aCommonEdge;
 337           anIsFirstDirection = Standard_False;
 338           continue;
 339         }
 340         else
 341           break;
 342       }
 343
 344       aBB.Add (aShell, aNextFace);
 345       aNbFacesInTape++;
 346       anAreaOfTape += theFaceAreaMap (aNextFace);
 347       aBB.Remove (theRes, aNextFace);
 348       theRemovedFaces.Add (aNextFace);
 349       if (theExtremalFaces.Contains (aNextFace))
 350       {
 351         aNbFacesDone++;
 352       }
 353       aCurrentEdge = aNextEdge;
 354       aNextEdge.Nullify();
 355       aCurrentFace = aNextFace;
 356     } //end of local loop
 357
 358     //Tape is formed
 359     Standard_Integer aNumberToSplit = (theIsToAddFaces)? aNbFacesInTape + aNbFacesDone : aNbFacesInTape - aNbFacesDone;
 360     if (!theIsToAddFaces && aNbFacesDone > 1)
 361     {
 362       Standard_Integer aRealNumberToSplit = (aNumberToSplit > 0)? aNumberToSplit : 1;
 363       Standard_Real anAverageAreaInTape = anAreaOfTape / aRealNumberToSplit;
 364       if (anAverageAreaInTape > theAverageArea)
 365       {
 366         Standard_Integer aNewNumberToSplit = RealToInt(round(anAreaOfTape / theAverageArea));
 367         if (aNewNumberToSplit < aNbFacesInTape)
 368         {
 369           Standard_Integer aNumberToIncrease = aNewNumberToSplit - aNumberToSplit;
 370           for (Standard_Integer jj = theNbExtremalFaces; jj < theNbExtremalFaces + aNumberToIncrease; jj++)
 371             theExtremalFaces.Add (theFacesAndAreas[jj].first);
 372           theNbExtremalFaces += aNumberToIncrease;
 373           aNumberToSplit = aNewNumberToSplit;
 374         }
 375       }
 376     }
 377     if (anIsRound && aNumberToSplit <= 1)
 378     {
 379       Standard_Integer aNumberToIncrease = 3 - aNumberToSplit;
 380       for (Standard_Integer jj = theNbExtremalFaces; jj < theNbExtremalFaces + aNumberToIncrease; jj++)
 381         theExtremalFaces.Add (theFacesAndAreas[jj].first);
 382       theNbExtremalFaces += aNumberToIncrease;
 383       aNumberToSplit = 3;
 384     }
 385     CorrectShell (aShell, theInputFace);
 386     ShapeUpgrade_UnifySameDomain aUnifier;
 387     aUnifier.Initialize (aShell, Standard_True, Standard_True);
 388     aUnifier.Build();
 389     TopoDS_Shape aUnifiedShape = aUnifier.Shape();
 390     //Splitting
 391     TopoDS_Shape aLocalResult = aUnifiedShape;
 392     Standard_Integer aNbFacesInLocalResult = 1;
 393     if (aNumberToSplit > 1)
 394     {
 395 #if OCC_VERSION_LARGE < 0x07050304
 396       aNbFacesInLocalResult = 1;
 397 #else
 398       ShapeUpgrade_ShapeDivideArea aLocalTool (aUnifiedShape);
 399       aLocalTool.SetSplittingByNumber (Standard_True);
 400       aLocalTool.MaxArea() = -1;
 401       if (anIsUiso)
 402         aLocalTool.SetNumbersUVSplits (aNumberToSplit, 1);
 403       else
 404         aLocalTool.SetNumbersUVSplits (1, aNumberToSplit);
 405       aLocalTool.Perform();
 406       aLocalResult = aLocalTool.Result();
 407       aNbFacesInLocalResult = aNumberToSplit;
 408 #endif
 409     }
 410     else
 411     {
 412       if (aNumberToSplit == 0)
 413       {
 414         theExtremalFaces.Add (theFacesAndAreas[theNbExtremalFaces].first);
 415         theNbExtremalFaces++;
 416       }
 417     }
 418     aBB.Add (theGlobalRes, aLocalResult);
 419
 420     aSum += Abs(aNbFacesInTape - aNbFacesInLocalResult);
 421   } //end of global loop
 422
 423   //Second global loop
 424   TopoDS_Compound aSecondComp;
 425   aBB.MakeCompound (aSecondComp);
 426   while (aSum < aDiff)
 427   {
 428     TopoDS_Shape aMaxShell;
 429     Standard_Integer aMaxNbFaces = 0;
 430     TopoDS_Iterator anIter (theGlobalRes);
 431     for (; anIter.More(); anIter.Next())
 432     {
 433       const TopoDS_Shape& aShell = anIter.Value();
 434       TopTools_IndexedMapOfShape aFaceMap;
 435       TopExp::MapShapes (aShell, TopAbs_FACE, aFaceMap);
 436       if (aFaceMap.Extent() > aMaxNbFaces)
 437       {
 438         aMaxNbFaces = aFaceMap.Extent();
 439         aMaxShell = aShell;
 440       }
 441     }
 442
 443     if (aMaxNbFaces == 1)
 444       break;
 445
 446     aBB.Remove (theGlobalRes, aMaxShell);
 447     //Find iso
 448     Standard_Boolean anIsUiso = Standard_True;
 449     TopTools_IndexedDataMapOfShapeListOfShape aLocalEFmap;
 450     TopExp::MapShapesAndAncestors (aMaxShell, TopAbs_EDGE, TopAbs_FACE, aLocalEFmap);
 451     for (Standard_Integer jj = 1; jj <= aLocalEFmap.Extent(); jj++)
 452     {
 453       const TopoDS_Edge& anEdge = TopoDS::Edge (aLocalEFmap.FindKey(jj));
 454       const TopTools_ListOfShape& aFaceList = aLocalEFmap(jj);
 455       if (aFaceList.Extent() == 2)
 456       {
 457         const TopoDS_Face& aFace = TopoDS::Face (aFaceList.First());
 458         anIsUiso = IsUiso (anEdge, aFace);
 459         break;
 460       }
 461     }
 462     CorrectShell (aMaxShell, theInputFace);
 463     ShapeUpgrade_UnifySameDomain aUnifier;
 464     aUnifier.Initialize (aMaxShell, Standard_True, Standard_True);
 465     aUnifier.Build();
 466     TopoDS_Shape aUnifiedShape = aUnifier.Shape();
 467     TopoDS_Shape aLocalResult = aUnifiedShape;
 468
 469     Standard_Integer aNumberToSplit = (theIsToAddFaces)? aMaxNbFaces + (aDiff-aSum) : aMaxNbFaces - (aDiff-aSum);
 470     if (aNumberToSplit > 1)
 471     {
 472 #if OCC_VERSION_LARGE < 0x07050304
 473       aNumberToSplit = 1;
 474 #else
 475       ShapeUpgrade_ShapeDivideArea aLocalTool (aUnifiedShape);
 476       aLocalTool.SetSplittingByNumber (Standard_True);
 477       aLocalTool.MaxArea() = -1;
 478       if (anIsUiso)
 479         aLocalTool.SetNumbersUVSplits (aNumberToSplit, 1);
 480       else
 481         aLocalTool.SetNumbersUVSplits (1, aNumberToSplit);
 482       aLocalTool.Perform();
 483       aLocalResult = aLocalTool.Result();
 484 #endif
 485     }
 486     else
 487       aNumberToSplit = 1;
 488
 489     aBB.Add (aSecondComp, aLocalResult);
 490
 491     if (theIsToAddFaces)
 492       break;
 493     aSum += aMaxNbFaces - aNumberToSplit;
 494   }
 495   aBB.Add (theGlobalRes, aSecondComp);
 496 }
 497
 498 //=================================================================================================
 499 bool GeomAlgoAPI_PointCloudOnFace::PointCloud(GeomShapePtr theFace,
 500                                               const int    theNumberOfPoints,
 501                                               GeomShapePtr thePoints,
 502                                               std::string& theError)
 503 {
 504 #ifdef _DEBUG
 505   std::cout << "GeomAlgoAPI_PointCloudOnFace::PointCloud" << std::endl;
 506 #endif
 507
 508 #if OCC_VERSION_LARGE < 0x07050304
 509   theError = "Improper OCCT version: please, use OCCT 7.5.3p4 or newer.";
 510   return false;
 511 #else
 512
 513   if (!theFace.get()) {
 514     theError = "Face for point cloud calculation is null";
 515     return false;
 516   }
 517
 518   TopoDS_Shape anInputShape = theFace->impl<TopoDS_Shape>();
 519
 520   if (anInputShape.ShapeType() != TopAbs_FACE) {
 521     theError = "Shape for normale calculation is not a face";
 522     return false;
 523   }
 524
 525   TopoDS_Face anInputFace = TopoDS::Face (anInputShape);
 526
 527   ShapeUpgrade_ShapeDivideArea tool (anInputFace);
 528   tool.SetSplittingByNumber (Standard_True);
 529   tool.NbParts() = theNumberOfPoints;
 530   tool.Perform();
 531   TopoDS_Shape res = tool.Result();
 532
 533   BRep_Builder aBB;
 534   TopoDS_Compound aGlobalRes;
 535   aBB.MakeCompound (aGlobalRes);
 536
 537   TopTools_IndexedMapOfShape aFaceMap;
 538   TopExp::MapShapes (res, TopAbs_FACE, aFaceMap);
 539   Standard_Integer aNbFaces = aFaceMap.Extent();
 540
 541   TopTools_IndexedDataMapOfShapeListOfShape aEFmap;
 542   TopExp::MapShapesAndAncestors (res, TopAbs_EDGE, TopAbs_FACE, aEFmap);
 543
 544   TopTools_MapOfShape aBiggestFaces, aSmallestFaces;
 545   Standard_Boolean aUseTriangulation = Standard_True;
 546   Standard_Boolean aSkipShared = Standard_False;
 547   if (aNbFaces != theNumberOfPoints)
 548   {
 549     Standard_Real aTotalArea = 0.;
 550     std::vector<std::pair<TopoDS_Shape, Standard_Real>> aFacesAndAreas (aNbFaces);
 551     for (Standard_Integer ii = 1; ii <= aNbFaces; ii++)
 552     {
 553       GProp_GProps aProps;
 554       BRepGProp::SurfaceProperties (aFaceMap(ii), aProps, aSkipShared, aUseTriangulation);
 555       Standard_Real anArea = aProps.Mass();
 556       aTotalArea += anArea;
 557       std::pair<TopoDS_Shape, Standard_Real> aFaceWithArea (aFaceMap(ii), anArea);
 558       aFacesAndAreas[ii-1] = aFaceWithArea;
 559     }
 560     std::sort (aFacesAndAreas.begin(), aFacesAndAreas.end(), comp);
 561
 562     Standard_Real anAverageArea = aTotalArea / theNumberOfPoints;
 563
 564     TopTools_DataMapOfShapeReal aFaceAreaMap;
 565     for (Standard_Integer ii = 0; ii < aNbFaces; ii++)
 566       aFaceAreaMap.Bind (aFacesAndAreas[ii].first, aFacesAndAreas[ii].second);
 567
 568     TopTools_MapOfShape aRemovedFaces;
 569
 570     if (aNbFaces < theNumberOfPoints)
 571     {
 572       Standard_Integer aNbMissingFaces = theNumberOfPoints - aNbFaces;
 573       for (Standard_Integer ii = aNbFaces-1; ii > aNbFaces - aNbMissingFaces - 1; ii--)
 574         aBiggestFaces.Add (aFacesAndAreas[ii].first);
 575
 576       ModifyFacesForGlobalResult (anInputFace, anAverageArea,
 577                                   Standard_True, //to add faces
 578                                   aNbMissingFaces, aBiggestFaces,
 579                                   aFacesAndAreas, aFaceAreaMap, aEFmap,
 580                                   res, aGlobalRes,
 581                                   aRemovedFaces);
 582     }
 583     else //aNbFaces > theNumberOfPoints
 584     {
 585       Standard_Integer aNbExcessFaces = aNbFaces - theNumberOfPoints;
 586       for (Standard_Integer ii = 0; ii < aNbExcessFaces; ii++)
 587         aSmallestFaces.Add (aFacesAndAreas[ii].first);
 588
 589       TopTools_IndexedDataMapOfShapeListOfShape aVFmap;
 590       TopExp::MapShapesAndAncestors (res, TopAbs_VERTEX, TopAbs_FACE, aVFmap);
 591
 592       //Remove smallest faces with free boundaries
 593       for (Standard_Integer ii = 0; ii < aNbExcessFaces; ii++)
 594       {
 595         const TopoDS_Face& aFace = TopoDS::Face (aFacesAndAreas[ii].first);
 596         Standard_Boolean anIsFreeBoundFound = Standard_False;
 597         TopExp_Explorer anExplo (aFace, TopAbs_EDGE);
 598         for (; anExplo.More(); anExplo.Next())
 599         {
 600           const TopoDS_Edge& anEdge = TopoDS::Edge (anExplo.Current());
 601           if (!BRep_Tool::Degenerated (anEdge) &&
 602               aEFmap.FindFromKey(anEdge).Extent() < 2)
 603           {
 604             anIsFreeBoundFound = Standard_True;
 605             break;
 606           }
 607         }
 608         if (anIsFreeBoundFound)
 609         {
 610           Standard_Real aMaxArea = 0.;
 611           for (anExplo.Init(aFace, TopAbs_VERTEX); anExplo.More(); anExplo.Next())
 612           {
 613             const TopoDS_Shape& aVertex = anExplo.Current();
 614             const TopTools_ListOfShape& aFaceList = aVFmap.FindFromKey (aVertex);
 615             TopTools_ListIteratorOfListOfShape anItl (aFaceList);
 616             for (; anItl.More(); anItl.Next())
 617             {
 618               Standard_Real anArea = aFaceAreaMap (anItl.Value());
 619               if (anArea > aMaxArea)
 620                 aMaxArea = anArea;
 621             }
 622           }
 623           Standard_Real anArreaOfSmallestFace = aFaceAreaMap (aFace);
 624           if (anArreaOfSmallestFace < aMaxArea / 16)
 625           {
 626             aBB.Remove (res, aFace);
 627             aRemovedFaces.Add (aFace);
 628           }
 629         }
 630       }
 631
 632       ModifyFacesForGlobalResult (anInputFace, anAverageArea,
 633                                   Standard_False, //to decrease number of faces
 634                                   aNbExcessFaces, aSmallestFaces,
 635                                   aFacesAndAreas, aFaceAreaMap, aEFmap,
 636                                   res, aGlobalRes,
 637                                   aRemovedFaces);
 638     }
 639   }
 640
 641   aBB.Add (aGlobalRes, res);
 642
 643   TopoDS_Compound aCompound;
 644   aBB.MakeCompound (aCompound);
 645   for (TopExp_Explorer aGlobalExplo(aGlobalRes, TopAbs_FACE); aGlobalExplo.More(); aGlobalExplo.Next())
 646   {
 647     const TopoDS_Face& aFace = TopoDS::Face (aGlobalExplo.Current());
 648     Standard_Boolean anIsNaturalRestrictions = Standard_True;
 649     TopExp_Explorer anExplo (aFace, TopAbs_EDGE);
 650     for (; anExplo.More(); anExplo.Next())
 651     {
 652       const TopoDS_Edge& anEdge = TopoDS::Edge (anExplo.Current());
 653       if (BRep_Tool::Degenerated (anEdge))
 654         continue;
 655       if (!aEFmap.Contains(anEdge) ||
 656           aEFmap.FindFromKey(anEdge).Extent() < 2)
 657       {
 658         anIsNaturalRestrictions = Standard_False;
 659         break;
 660       }
 661     }
 662
 663     gp_Pnt aPnt = GetMidPnt2d (aFace, anIsNaturalRestrictions);
 664     TopoDS_Vertex aVertex = BRepLib_MakeVertex (aPnt);
 665     aBB.Add (aCompound, aVertex);
 666   }
 667
 668   thePoints->setImpl(new TopoDS_Shape(aCompound));
 669
 670   return true;
 671 #endif
 672 }