Salome HOME
Copyright update 2021
[modules/shaper.git] / src / GeomAlgoAPI / GeomAlgoAPI_ShapeTools.cpp
index 48e3b6cf5054c72c56b4edcda66f97e989115ae0..98d3cfc4c9252e0c2a84236a0c646cea9ccd03c1 100644 (file)
@@ -1,4 +1,4 @@
-// Copyright (C) 2014-2017  CEA/DEN, EDF R&D
+// Copyright (C) 2014-2021  CEA/DEN, EDF R&D
 //
 // This library is free software; you can redistribute it and/or
 // modify it under the terms of the GNU Lesser General Public
 //
 // You should have received a copy of the GNU Lesser General Public
 // License along with this library; if not, write to the Free Software
-// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
 //
-// See http://www.salome-platform.org/ or
-// email : webmaster.salome@opencascade.com<mailto:webmaster.salome@opencascade.com>
+// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
 //
 
 #include "GeomAlgoAPI_ShapeTools.h"
 #include <GeomAPI_Wire.h>
 
 #include <Bnd_Box.hxx>
+
+#include <BRep_Tool.hxx>
 #include <BRep_Builder.hxx>
-#include <BRepAdaptor_Curve.hxx>
 #include <BRepAlgo.hxx>
 #include <BRepAlgo_FaceRestrictor.hxx>
+#include <BRepAdaptor_Curve.hxx>
 #include <BRepBndLib.hxx>
+#include <BRepBuilderAPI_Copy.hxx>
 #include <BRepBuilderAPI_FindPlane.hxx>
 #include <BRepBuilderAPI_MakeEdge.hxx>
 #include <BRepBuilderAPI_MakeFace.hxx>
+#include <BRepBuilderAPI_MakeVertex.hxx>
+#include <BRepBuilderAPI_MakeWire.hxx>
 #include <BRepCheck_Analyzer.hxx>
 #include <BRepExtrema_DistShapeShape.hxx>
 #include <BRepExtrema_ExtCF.hxx>
 #include <BRepTools_WireExplorer.hxx>
 #include <BRepTopAdaptor_FClass2d.hxx>
 #include <BRepClass_FaceClassifier.hxx>
+#include <BRepLib_CheckCurveOnSurface.hxx>
+#include <BRepLProp.hxx>
+
+#include <BOPAlgo_Builder.hxx>
+
 #include <Geom2d_Curve.hxx>
 #include <Geom2d_Curve.hxx>
-#include <BRepLib_CheckCurveOnSurface.hxx>
-#include <BRep_Tool.hxx>
-#include  <Geom_CylindricalSurface.hxx>
+
+#include <Geom_CylindricalSurface.hxx>
 #include <Geom_Line.hxx>
 #include <Geom_Plane.hxx>
+#include <Geom_RectangularTrimmedSurface.hxx>
+
 #include <GeomAPI_ProjectPointOnCurve.hxx>
 #include <GeomAPI_ShapeIterator.h>
+
 #include <GeomLib_IsPlanarSurface.hxx>
 #include <GeomLib_Tool.hxx>
 #include <GeomAPI_IntCS.hxx>
+
 #include <gp_Pln.hxx>
 #include <GProp_GProps.hxx>
+
 #include <IntAna_IntConicQuad.hxx>
 #include <IntAna_Quadric.hxx>
-#include <NCollection_Vector.hxx>
+
 #include <ShapeAnalysis.hxx>
 #include <ShapeAnalysis_Surface.hxx>
-#include <TopoDS_Builder.hxx>
+
+#include <TopoDS.hxx>
 #include <TopoDS_Edge.hxx>
 #include <TopoDS_Face.hxx>
 #include <TopoDS_Shape.hxx>
 #include <TopoDS_Shell.hxx>
 #include <TopoDS_Vertex.hxx>
-#include <TopoDS.hxx>
+#include <TopoDS_Builder.hxx>
+
 #include <TopExp.hxx>
 #include <TopExp_Explorer.hxx>
-#include <TopTools_ListIteratorOfListOfShape.hxx>
 
+#include <TopTools_ListIteratorOfListOfShape.hxx>
 
-#include <BOPAlgo_Builder.hxx>
-#include <BRepBuilderAPI_MakeVertex.hxx>
-#include <TopoDS_Edge.hxx>
+#include <NCollection_Vector.hxx>
 
 //==================================================================================================
 static GProp_GProps props(const TopoDS_Shape& theShape)
@@ -112,7 +125,7 @@ static GProp_GProps props(const TopoDS_Shape& theShape)
 }
 
 //==================================================================================================
-double GeomAlgoAPI_ShapeTools::volume(const std::shared_ptr<GeomAPI_Shape> theShape)
+double GeomAlgoAPI_ShapeTools::length(const std::shared_ptr<GeomAPI_Shape> theShape)
 {
   GProp_GProps aGProps;
   if(!theShape.get()) {
@@ -122,12 +135,29 @@ double GeomAlgoAPI_ShapeTools::volume(const std::shared_ptr<GeomAPI_Shape> theSh
   if(aShape.IsNull()) {
     return 0.0;
   }
+
+  BRepGProp::LinearProperties(aShape, aGProps, Standard_True);
+  return  aGProps.Mass();
+}
+
+//==================================================================================================
+double GeomAlgoAPI_ShapeTools::volume(const std::shared_ptr<GeomAPI_Shape> theShape)
+{
+  if(!theShape.get()) {
+    return 0.0;
+  }
+  const TopoDS_Shape& aShape = theShape->impl<TopoDS_Shape>();
+  if(aShape.IsNull()) {
+    return 0.0;
+  }
   const Standard_Real anEps = 1.e-6;
-  if (aShape.ShapeType() <= TopAbs_SOLID)
-    BRepGProp::VolumeProperties(aShape, aGProps, anEps);
-  else
-    BRepGProp::SurfaceProperties(aShape, aGProps, anEps);
-  return aGProps.Mass();
+  double aVolume = 0.0;
+  for (TopExp_Explorer anExp(aShape, TopAbs_SOLID); anExp.More(); anExp.Next()) {
+    GProp_GProps aGProps;
+    BRepGProp::VolumeProperties(anExp.Current(), aGProps, anEps);
+    aVolume += aGProps.Mass();
+  }
+  return aVolume;
 }
 
 //==================================================================================================
@@ -200,9 +230,12 @@ double GeomAlgoAPI_ShapeTools::minimalDistance(const GeomShapePtr& theShape1,
 std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
   const std::shared_ptr<GeomAPI_Shape> theCompound,
   const GeomAPI_Shape::ShapeType theType,
-  ListOfShape& theCombinedShapes,
-  ListOfShape& theFreeShapes)
+  ListOfShape& theResuts)
 {
+
+  ListOfShape aResCombinedShapes;
+  ListOfShape aResFreeShapes;
+
   GeomShapePtr aResult = theCompound;
 
   if(!theCompound.get()) {
@@ -220,17 +253,25 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
     aTA = TopAbs_SOLID;
   }
 
-  theCombinedShapes.clear();
-  theFreeShapes.clear();
+  // map from the resulting shapes to minimal index of the used shape from theCompound list
+  std::map<GeomShapePtr, int> anInputOrder;
+  // map from ancestors-shapes to the index of shapes in theCompound
+  NCollection_DataMap<TopoDS_Shape, int> anAncestorsOrder;
 
   // Get free shapes.
+  int anOrder = 0;
   const TopoDS_Shape& aShapesComp = theCompound->impl<TopoDS_Shape>();
-  for(TopoDS_Iterator anIter(aShapesComp); anIter.More(); anIter.Next() ) {
+  for(TopoDS_Iterator anIter(aShapesComp); anIter.More(); anIter.Next(), anOrder++) {
     const TopoDS_Shape& aShape = anIter.Value();
     if(aShape.ShapeType() > aTA) {
       std::shared_ptr<GeomAPI_Shape> aGeomShape(new GeomAPI_Shape);
       aGeomShape->setImpl<TopoDS_Shape>(new TopoDS_Shape(aShape));
-      theFreeShapes.push_back(aGeomShape);
+      aResFreeShapes.push_back(aGeomShape);
+      anInputOrder[aGeomShape] = anOrder;
+    } else {
+      for(TopExp_Explorer anExp(aShape, aTA); anExp.More(); anExp.Next()) {
+        anAncestorsOrder.Bind(anExp.Current(), anOrder);
+      }
     }
   }
 
@@ -240,13 +281,13 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
   if(aMapSA.IsEmpty()) {
     return aResult;
   }
+  theResuts.clear();
 
   // Get all shapes with common sub-shapes and free shapes.
   NCollection_Map<TopoDS_Shape> aFreeShapes;
   NCollection_Vector<NCollection_Map<TopoDS_Shape>> aShapesWithCommonSubshapes;
   for(TopTools_IndexedDataMapOfShapeListOfShape::Iterator
       anIter(aMapSA); anIter.More(); anIter.Next()) {
-    const TopoDS_Shape& aShape = anIter.Key();
     TopTools_ListOfShape& aListOfShape = anIter.ChangeValue();
     if(aListOfShape.IsEmpty()) {
       continue;
@@ -268,8 +309,8 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
           aTempIter(aTempList); aTempIter.More(); aTempIter.Next()) {
         const TopoDS_Shape& aTempShape = aTempIter.Value();
         for(TopTools_IndexedDataMapOfShapeListOfShape::Iterator
-            anIter(aMapSA); anIter.More(); anIter.Next()) {
-          TopTools_ListOfShape& aTempListOfShape = anIter.ChangeValue();
+            anIter2(aMapSA); anIter2.More(); anIter2.Next()) {
+          TopTools_ListOfShape& aTempListOfShape = anIter2.ChangeValue();
           if(aTempListOfShape.IsEmpty()) {
             continue;
           } else if(aTempListOfShape.Size() == 1 && aTempListOfShape.First() == aTempShape) {
@@ -305,6 +346,7 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
     TopoDS_Shell aShell;
     TopoDS_CompSolid aCSolid;
     TopoDS_Builder aBuilder;
+    anOrder = -1;
     theType ==
       GeomAPI_Shape::COMPSOLID ? aBuilder.MakeCompSolid(aCSolid) : aBuilder.MakeShell(aShell);
     NCollection_Map<TopoDS_Shape>& aShapesMap = anIter.ChangeValue();
@@ -314,13 +356,17 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
         theType ==
           GeomAPI_Shape::COMPSOLID ? aBuilder.Add(aCSolid, aShape) : aBuilder.Add(aShell, aShape);
         aShapesMap.Remove(aShape);
+        int aThisOrder = anAncestorsOrder.Find(aShape);
+        if (anOrder == -1 || aThisOrder < anOrder)
+          anOrder = aThisOrder; // take the minimum order position
       }
     }
     std::shared_ptr<GeomAPI_Shape> aGeomShape(new GeomAPI_Shape);
     TopoDS_Shape* aSh = theType == GeomAPI_Shape::COMPSOLID ? new TopoDS_Shape(aCSolid) :
                                                               new TopoDS_Shape(aShell);
     aGeomShape->setImpl<TopoDS_Shape>(aSh);
-    theCombinedShapes.push_back(aGeomShape);
+    aResCombinedShapes.push_back(aGeomShape);
+    anInputOrder[aGeomShape] = anOrder;
   }
 
   // Adding free shapes.
@@ -329,25 +375,30 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::combineShapes(
     if(aFreeShapes.Contains(aShape)) {
       std::shared_ptr<GeomAPI_Shape> aGeomShape(new GeomAPI_Shape);
       aGeomShape->setImpl<TopoDS_Shape>(new TopoDS_Shape(aShape));
-      theFreeShapes.push_back(aGeomShape);
+      aResFreeShapes.push_back(aGeomShape);
+      anInputOrder[aGeomShape] = anAncestorsOrder.Find(aShape);
     }
   }
 
-  if(theCombinedShapes.size() == 1 && theFreeShapes.size() == 0) {
-    aResult = theCombinedShapes.front();
-  } else if(theCombinedShapes.size() == 0 && theFreeShapes.size() == 1) {
-    aResult = theFreeShapes.front();
+  if(aResCombinedShapes.size() == 1 && aResFreeShapes.size() == 0) {
+    aResult = aResCombinedShapes.front();
+    theResuts.push_back(aResult);
+  } else if(aResCombinedShapes.size() == 0 && aResFreeShapes.size() == 1) {
+    aResult = aResFreeShapes.front();
+    theResuts.push_back(aResult);
   } else {
     TopoDS_Compound aResultComp;
     TopoDS_Builder aBuilder;
     aBuilder.MakeCompound(aResultComp);
-    for(ListOfShape::const_iterator anIter = theCombinedShapes.cbegin();
-        anIter != theCombinedShapes.cend(); anIter++) {
-      aBuilder.Add(aResultComp, (*anIter)->impl<TopoDS_Shape>());
-    }
-    for(ListOfShape::const_iterator anIter = theFreeShapes.cbegin();
-        anIter != theFreeShapes.cend(); anIter++) {
-      aBuilder.Add(aResultComp, (*anIter)->impl<TopoDS_Shape>());
+    // put to result compound and result list in accordance to the order numbers
+    std::map<GeomShapePtr, int>::iterator anInputIter = anInputOrder.begin();
+    std::map<int, GeomShapePtr> aNums;
+    for(; anInputIter != anInputOrder.end(); anInputIter++)
+      aNums[anInputIter->second] = anInputIter->first;
+    std::map<int, GeomShapePtr>::iterator aNumsIter = aNums.begin();
+    for(; aNumsIter != aNums.end(); aNumsIter++) {
+      aBuilder.Add(aResultComp, (aNumsIter->second)->impl<TopoDS_Shape>());
+      theResuts.push_back(aNumsIter->second);
     }
     aResult->setImpl(new TopoDS_Shape(aResultComp));
   }
@@ -487,7 +538,7 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::groupSharedTopology(
   TopoDS_Compound aCompound;
   BRep_Builder aBuilder;
   aBuilder.MakeCompound(aCompound);
-  ListOfShape aCompSolids, aFreeSolids;
+  ListOfShape aSolids;
   for (NCollection_Vector<TopTools_MapOfShape>::Iterator anIt(aGroups); anIt.More(); anIt.Next()) {
     const TopTools_MapOfShape& aGroup = anIt.ChangeValue();
     GeomShapePtr aGeomShape(new GeomAPI_Shape());
@@ -505,8 +556,7 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::groupSharedTopology(
       aGeomShape->setImpl(new TopoDS_Shape(makeCompound(anOrderedGoup)));
       aGeomShape = GeomAlgoAPI_ShapeTools::combineShapes(aGeomShape,
                                                          GeomAPI_Shape::COMPSOLID,
-                                                         aCompSolids,
-                                                         aFreeSolids);
+                                                         aSolids);
     }
     aBuilder.Add(aCompound, aGeomShape->impl<TopoDS_Shape>());
   }
@@ -518,6 +568,24 @@ std::shared_ptr<GeomAPI_Shape> GeomAlgoAPI_ShapeTools::groupSharedTopology(
   return aResult;
 }
 
+//==================================================================================================
+bool GeomAlgoAPI_ShapeTools::hasSharedTopology(const ListOfShape& theShapes,
+                                               const GeomAPI_Shape::ShapeType theShapeType)
+{
+  TopTools_IndexedMapOfShape aSubs;
+  for (ListOfShape::const_iterator anIt = theShapes.begin(); anIt != theShapes.end(); ++anIt) {
+    TopTools_IndexedMapOfShape aCurSubs;
+    TopExp::MapShapes((*anIt)->impl<TopoDS_Shape>(), (TopAbs_ShapeEnum)theShapeType, aCurSubs);
+    for (TopTools_IndexedMapOfShape::Iterator aSubIt(aCurSubs); aSubIt.More(); aSubIt.Next()) {
+      if (aSubs.Contains(aSubIt.Value()))
+        return true;
+      else
+        aSubs.Add(aSubIt.Value());
+    }
+  }
+  return false;
+}
+
 //==================================================================================================
 std::list<std::shared_ptr<GeomAPI_Pnt> >
   GeomAlgoAPI_ShapeTools::getBoundingBox(const ListOfShape& theShapes, const double theEnlarge)
@@ -541,7 +609,6 @@ std::list<std::shared_ptr<GeomAPI_Pnt> >
   Standard_Real aYArr[2] = {aBndBox.CornerMin().Y(), aBndBox.CornerMax().Y()};
   Standard_Real aZArr[2] = {aBndBox.CornerMin().Z(), aBndBox.CornerMax().Z()};
   std::list<std::shared_ptr<GeomAPI_Pnt> > aResultPoints;
-  int aNum = 0;
   for(int i = 0; i < 2; i++) {
     for(int j = 0; j < 2; j++) {
       for(int k = 0; k < 2; k++) {
@@ -800,6 +867,44 @@ std::shared_ptr<GeomAPI_Shape>
   return anOuterWire;
 }
 
+//==================================================================================================
+static bool boundaryOfEdge(const std::shared_ptr<GeomAPI_Edge> theEdge,
+                          const std::shared_ptr<GeomAPI_Vertex> theVertex,
+                          double& theParam)
+{
+  GeomPointPtr aPoint = theVertex->point();
+  GeomPointPtr aFirstPnt = theEdge->firstPoint();
+  double aFirstPntTol = theEdge->firstPointTolerance();
+  GeomPointPtr aLastPnt = theEdge->lastPoint();
+  double aLastPntTol = theEdge->lastPointTolerance();
+
+  double aFirst, aLast;
+  theEdge->getRange(aFirst, aLast);
+
+  bool isFirst = aPoint->distance(aFirstPnt) <= aFirstPntTol;
+  bool isLast = aPoint->distance(aLastPnt) <= aLastPntTol;
+  if (isFirst)
+    theParam = aFirst;
+  else if (isLast)
+    theParam = aLast;
+
+  return isFirst != isLast;
+}
+
+bool GeomAlgoAPI_ShapeTools::isTangent(const std::shared_ptr<GeomAPI_Edge> theEdge1,
+                                       const std::shared_ptr<GeomAPI_Edge> theEdge2,
+                                       const std::shared_ptr<GeomAPI_Vertex> theTgPoint)
+{
+  double aParE1 = 0, aParE2 = 0;
+  if (!boundaryOfEdge(theEdge1, theTgPoint, aParE1) ||
+      !boundaryOfEdge(theEdge2, theTgPoint, aParE2))
+    return false;
+
+  BRepAdaptor_Curve aC1(theEdge1->impl<TopoDS_Edge>());
+  BRepAdaptor_Curve aC2(theEdge2->impl<TopoDS_Edge>());
+  return BRepLProp::Continuity(aC1, aC2, aParE1, aParE2) >= GeomAbs_G1;
+}
+
 //==================================================================================================
 bool GeomAlgoAPI_ShapeTools::isParallel(const std::shared_ptr<GeomAPI_Edge> theEdge,
                                         const std::shared_ptr<GeomAPI_Face> theFace)
@@ -817,8 +922,7 @@ bool GeomAlgoAPI_ShapeTools::isParallel(const std::shared_ptr<GeomAPI_Edge> theE
 
 //==================================================================================================
 std::list<std::shared_ptr<GeomAPI_Vertex> > GeomAlgoAPI_ShapeTools::intersect(
-  const std::shared_ptr<GeomAPI_Edge> theEdge, const std::shared_ptr<GeomAPI_Face> theFace,
-  const bool thePointsOutsideFace)
+  const std::shared_ptr<GeomAPI_Edge> theEdge, const std::shared_ptr<GeomAPI_Face> theFace)
 {
   std::list<std::shared_ptr<GeomAPI_Vertex> > aResult;
   if(!theEdge.get() || !theFace.get()) {
@@ -1000,7 +1104,7 @@ static TopoDS_Wire fixParametricGaps(const TopoDS_Wire& theWire)
 {
   TopoDS_Wire aFixedWire;
   Handle(Geom_Curve) aPrevCurve;
-  double aPrevLastParam = 0.0;
+  double aPrevLastParam = -Precision::Infinite();
 
   BRep_Builder aBuilder;
   aBuilder.MakeWire(aFixedWire);
@@ -1010,7 +1114,7 @@ static TopoDS_Wire fixParametricGaps(const TopoDS_Wire& theWire)
     TopoDS_Edge anEdge = aWExp.Current();
     double aFirst, aLast;
     Handle(Geom_Curve) aCurve = BRep_Tool::Curve(anEdge, aFirst, aLast);
-    if (aCurve == aPrevCurve) {
+    if (aCurve == aPrevCurve && Abs(aFirst - aPrevLastParam) > Precision::Confusion()) {
       // if parametric gap occurs, create new edge based on the copied curve
       aCurve = Handle(Geom_Curve)::DownCast(aCurve->Copy());
       TopoDS_Vertex aV1, aV2;
@@ -1030,24 +1134,47 @@ static TopoDS_Wire fixParametricGaps(const TopoDS_Wire& theWire)
   return aFixedWire;
 }
 
+//==================================================================================================
 std::shared_ptr<GeomAPI_Edge> GeomAlgoAPI_ShapeTools::wireToEdge(
       const std::shared_ptr<GeomAPI_Wire>& theWire)
 {
   GeomEdgePtr anEdge;
   if (theWire) {
     TopoDS_Wire aWire = theWire->impl<TopoDS_Wire>();
-    // Workaround: when concatenate a wire consisting of two edges based on the same B-spline curve
-    // (non-periodic, but having equal start and end points), first of which is placed at the end
-    // on the curve and second is placed at the start, this workaround copies second curve to avoid
-    // treating these edges as a single curve by setting trim parameters.
-    aWire = fixParametricGaps(aWire);
-    TopoDS_Edge aNewEdge = BRepAlgo::ConcatenateWireC0(aWire);
+    BRepTools_WireExplorer aWExp(aWire);
+    TopoDS_Edge aNewEdge = aWExp.Current();
+    aWExp.Next();
+    if (aWExp.More()) {
+      // Workaround for the closed wire to avoid jumping of its start point:
+      // split this wire for two parts, convert them to edges, then compose together
+      if (BRep_Tool::IsClosed(aWire)) {
+        aWire = TopoDS::Wire(BRepBuilderAPI_Copy(aWire).Shape());
+        aWExp.Init(aWire);
+        aNewEdge = aWExp.Current();
+
+        BRep_Builder().Remove(aWire, aNewEdge);
+        GeomWirePtr aSplitWire(new GeomAPI_Wire);
+        aSplitWire->setImpl(new TopoDS_Wire(aWire));
+        GeomEdgePtr aMergedEdge = wireToEdge(aSplitWire);
+
+        aWire = BRepBuilderAPI_MakeWire(aNewEdge, aMergedEdge->impl<TopoDS_Edge>());
+      }
+
+      // Workaround: when concatenate a wire consisting of two edges based on the same B-spline
+      // curve (non-periodic, but having equal start and end points), first of which is placed
+      // at the end on the curve and second is placed at the start, this workaround copies
+      // second curve to avoid treating these edges as a single curve by setting trim parameters.
+      aWire = fixParametricGaps(aWire);
+      aWire = BRepAlgo::ConcatenateWire(aWire, GeomAbs_G1); // join smooth parts of wire
+      aNewEdge = BRepAlgo::ConcatenateWireC0(aWire); // join C0 parts of wire
+    }
     anEdge = GeomEdgePtr(new GeomAPI_Edge);
     anEdge->setImpl(new TopoDS_Edge(aNewEdge));
   }
   return anEdge;
 }
 
+//==================================================================================================
 ListOfShape GeomAlgoAPI_ShapeTools::getLowLevelSubShapes(const GeomShapePtr& theShape)
 {
   ListOfShape aSubShapes;
@@ -1068,4 +1195,131 @@ ListOfShape GeomAlgoAPI_ShapeTools::getLowLevelSubShapes(const GeomShapePtr& the
   }
 
   return aSubShapes;
-}
\ No newline at end of file
+}
+
+//==================================================================================================
+static void getMinMaxPointsOnLine(const std::list<std::shared_ptr<GeomAPI_Pnt> >& thePoints,
+                                  const gp_Dir theDir,
+                                  double& theMin, double& theMax)
+{
+  theMin = RealLast();
+  theMax = RealFirst();
+  // Project bounding points on theDir
+  for (std::list<std::shared_ptr<GeomAPI_Pnt> >::const_iterator
+         aPointsIt = thePoints.begin(); aPointsIt != thePoints.end(); aPointsIt++) {
+    const gp_Pnt& aPnt = (*aPointsIt)->impl<gp_Pnt>();
+    gp_Dir aPntDir (aPnt.XYZ());
+    Standard_Real proj = (theDir*aPntDir) * aPnt.XYZ().Modulus();
+    if (proj < theMin) theMin = proj;
+    if (proj > theMax) theMax = proj;
+  }
+}
+
+//==================================================================================================
+void GeomAlgoAPI_ShapeTools::computeThroughAll(const ListOfShape& theObjects,
+                                               const ListOfShape& theBaseShapes,
+                                               const std::shared_ptr<GeomAPI_Dir> theDir,
+                                               double& theToSize, double& theFromSize)
+{
+  // Bounding box of objects
+  std::list<std::shared_ptr<GeomAPI_Pnt> > aBndObjs =
+      GeomAlgoAPI_ShapeTools::getBoundingBox(theObjects);
+  if (aBndObjs.size() != 8) {
+    return;
+  }
+
+  // the value to enlarge the bounding box of each object to make the extruded shape
+  // a little bit larger than overall objects to get the correct result of Boolean CUT operation
+  double anEnlargement = 0.1 * aBndObjs.front()->distance(aBndObjs.back());
+
+  // Prism direction
+  if (theDir.get()) {
+    // One direction for all prisms
+    gp_Dir aDir = theDir->impl<gp_Dir>();
+
+    // Bounding box of the base
+    std::list<std::shared_ptr<GeomAPI_Pnt> > aBndBases =
+        GeomAlgoAPI_ShapeTools::getBoundingBox(theBaseShapes);
+    if (aBndBases.size() != 8) {
+      return;
+    }
+
+    // Objects bounds
+    Standard_Real lowBnd, upperBnd;
+    getMinMaxPointsOnLine(aBndObjs, aDir, lowBnd, upperBnd);
+
+    // Base bounds
+    Standard_Real lowBase, upperBase;
+    getMinMaxPointsOnLine(aBndBases, aDir, lowBase, upperBase);
+
+    // ----------.-----.---------.--------------.-----------> theDir
+    //       lowBnd   lowBase   upperBase    upperBnd
+
+    theToSize = upperBnd - lowBase;
+    theFromSize = upperBase - lowBnd;
+  } else {
+    // Direction is a normal to each base shape (different normals to bases)
+    // So we calculate own sizes for each base shape
+    theToSize = 0.0;
+    theFromSize = 0.0;
+
+    for (ListOfShape::const_iterator anIt = theBaseShapes.begin();
+         anIt != theBaseShapes.end(); ++anIt) {
+      const GeomShapePtr& aBaseShape_i = (*anIt);
+      ListOfShape aBaseShapes_i;
+      aBaseShapes_i.push_back(aBaseShape_i);
+
+      // Bounding box of the base
+      std::list<std::shared_ptr<GeomAPI_Pnt> > aBndBases =
+          GeomAlgoAPI_ShapeTools::getBoundingBox(aBaseShapes_i, anEnlargement);
+      if (aBndBases.size() != 8) {
+        return;
+      }
+
+      // Direction (normal to aBaseShapes_i)
+      // Code like in GeomAlgoAPI_Prism
+      gp_Dir aDir;
+      const TopoDS_Shape& aBaseShape = aBaseShape_i->impl<TopoDS_Shape>();
+      BRepBuilderAPI_FindPlane aFindPlane(aBaseShape);
+      if (aFindPlane.Found() == Standard_True) {
+        Handle(Geom_Plane) aPlane;
+        if (aBaseShape.ShapeType() == TopAbs_FACE || aBaseShape.ShapeType() == TopAbs_SHELL) {
+          TopExp_Explorer anExp(aBaseShape, TopAbs_FACE);
+          const TopoDS_Shape& aFace = anExp.Current();
+          Handle(Geom_Surface) aSurface = BRep_Tool::Surface(TopoDS::Face(aFace));
+          if(aSurface->DynamicType() == STANDARD_TYPE(Geom_RectangularTrimmedSurface)) {
+            Handle(Geom_RectangularTrimmedSurface) aTrimSurface =
+              Handle(Geom_RectangularTrimmedSurface)::DownCast(aSurface);
+            aSurface = aTrimSurface->BasisSurface();
+          }
+          if(aSurface->DynamicType() != STANDARD_TYPE(Geom_Plane)) {
+            return;
+          }
+          aPlane = Handle(Geom_Plane)::DownCast(aSurface);
+        } else {
+          aPlane = aFindPlane.Plane();
+        }
+        aDir = aPlane->Axis().Direction();
+      } else {
+        return;
+      }
+
+      // Objects bounds
+      Standard_Real lowBnd, upperBnd;
+      getMinMaxPointsOnLine(aBndObjs, aDir, lowBnd, upperBnd);
+
+      // Base bounds
+      Standard_Real lowBase, upperBase;
+      getMinMaxPointsOnLine(aBndBases, aDir, lowBase, upperBase);
+
+      // ----------.-----.---------.--------------.-----------> theDir
+      //       lowBnd   lowBase   upperBase    upperBnd
+
+      double aToSize_i = upperBnd - lowBase;
+      double aFromSize_i = upperBase - lowBnd;
+
+      if (aToSize_i > theToSize) theToSize = aToSize_i;
+      if (aFromSize_i > theFromSize) theFromSize = aFromSize_i;
+    }
+  }
+}