src/GeomAPI/GeomAPI_Solid.cpp

   1 // Copyright (C) 2018-2020  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 "GeomAPI_Solid.h"
  21
  22 #include "GeomAPI_Box.h"
  23 #include "GeomAPI_Cone.h"
  24 #include "GeomAPI_Cylinder.h"
  25 #include "GeomAPI_Dir.h"
  26 #include "GeomAPI_Face.h"
  27 #include "GeomAPI_Pln.h"
  28 #include "GeomAPI_Pnt.h"
  29 #include "GeomAPI_Shell.h"
  30 #include "GeomAPI_Sphere.h"
  31 #include "GeomAPI_Torus.h"
  32 #include "GeomAPI_XYZ.h"
  33
  34 #include <BRep_Builder.hxx>
  35 #include <BRepGProp.hxx>
  36 #include <GProp_GProps.hxx>
  37 #include <Precision.hxx>
  38 #include <TopExp_Explorer.hxx>
  39 #include <TopoDS_Wire.hxx>
  40
  41 //==================================================================================================
  42 GeomAPI_Solid::GeomAPI_Solid() : GeomAPI_Shape()
  43 {
  44 }
  45
  46 //==================================================================================================
  47 GeomAPI_Solid::GeomAPI_Solid(const std::shared_ptr<GeomAPI_Shape>& theShape)
  48 {
  49   if (!theShape->isNull() && theShape->isSolid()) {
  50     setImpl(new TopoDS_Shape(theShape->impl<TopoDS_Shape>()));
  51   }
  52 }
  53
  54 //==================================================================================================
  55 std::shared_ptr<GeomAPI_Sphere> GeomAPI_Solid::getSphere() const
  56 {
  57   GeomSpherePtr aSphere;
  58   ListOfShape aShells = subShapes(SHELL);
  59   if (aShells.size() == 1)
  60     aSphere = aShells.front()->shell()->getSphere();
  61   return aSphere;
  62 }
  63
  64 //==================================================================================================
  65 std::shared_ptr<GeomAPI_Cylinder> GeomAPI_Solid::getCylinder() const
  66 {
  67   bool isCylinder = false;
  68
  69   GeomPointPtr aLocation;
  70   GeomDirPtr anAxis;
  71   double aRadius = 0.0;
  72   double aHeight = 0.0;
  73
  74   GeomPlanePtr aCaps[2];
  75
  76   for (TopExp_Explorer anExp(impl<TopoDS_Shape>(), TopAbs_FACE); anExp.More(); anExp.Next()) {
  77     GeomFacePtr aFace(new GeomAPI_Face);
  78     aFace->setImpl(new TopoDS_Shape(anExp.Current()));
  79
  80     GeomCylinderPtr aCurCyl = aFace->getCylinder();
  81     if (aCurCyl) {
  82       if (isCylinder) { // at least one cylindrical face is found
  83         // compare radii
  84         if (Abs(aRadius - aCurCyl->radius()) >= Precision::Confusion() ||
  85           // check directions are collinear
  86           !anAxis->isParallel(aCurCyl->axis()) ||
  87           // check current center is on the main axis
  88           anAxis->xyz()->cross(aLocation->xyz()->decreased(aCurCyl->location()->xyz())
  89           )->squareModulus() >= Precision::SquareConfusion()) {
  90           isCylinder = false;
  91           break;
  92         }
  93       }
  94       else { // first cylinder is found
  95         aLocation = aCurCyl->location();
  96         if (anAxis) {
  97           // the plane is found => compare directions
  98           if (!anAxis->isParallel(aCurCyl->axis()))
  99             break;
 100         }
 101         else
 102           anAxis = aCurCyl->axis();
 103         aRadius = aCurCyl->radius();
 104         aHeight = aCurCyl->height();
 105         isCylinder = true;
 106       }
 107     }
 108     else {
 109       // check the face is planar
 110       bool isPlaneApplicable = false;
 111       GeomPlanePtr aCurPln = aFace->getPlane();
 112       if (aCurPln) {
 113         // verify the plane is already exists
 114         int aLastPlanIndex = 0;
 115         while (aLastPlanIndex < 2) {
 116           if (!aCaps[aLastPlanIndex]) {
 117             // add new plane
 118             aCaps[aLastPlanIndex] = aCurPln;
 119             break;
 120           }
 121           if (aCaps[aLastPlanIndex]->isCoincident(aCurPln))
 122             break;
 123           ++aLastPlanIndex;
 124         }
 125
 126         isPlaneApplicable = aLastPlanIndex < 2;
 127       }
 128
 129       if (isPlaneApplicable) {
 130         if (!anAxis) // no cylinder is found, store the normal as further cylinder's axis
 131           anAxis = aCurPln->direction();
 132       }
 133       else {
 134         isCylinder = false;
 135         break;
 136       }
 137     }
 138   }
 139
 140   isCylinder = isCylinder && aCaps[0] && aCaps[1] &&
 141                aCaps[0]->direction()->isParallel(anAxis) &&
 142                aCaps[1]->direction()->isParallel(anAxis);
 143
 144   GeomCylinderPtr aCylinder;
 145   if (isCylinder) {
 146     // intersect planes with cylinder's axis
 147     std::shared_ptr<GeomAPI_XYZ> anAxisXYZ = anAxis->xyz();
 148     std::shared_ptr<GeomAPI_XYZ> aLocationXYZ = aLocation->xyz();
 149     double aParam0 = anAxisXYZ->dot( aCaps[0]->location()->xyz()->decreased(aLocationXYZ) );
 150     double aParam1 = anAxisXYZ->dot( aCaps[1]->location()->xyz()->decreased(aLocationXYZ) );
 151     if (aParam0 > aParam1 + Precision::Confusion()) {
 152       double tmp = aParam0;
 153       aParam0 = aParam1;
 154       aParam1 = tmp;
 155     }
 156
 157     // update location of cylinder to be coincident with one of planes
 158     aLocation->setX(aLocation->x() + aParam0 * anAxis->x());
 159     aLocation->setY(aLocation->y() + aParam0 * anAxis->y());
 160     aLocation->setZ(aLocation->z() + aParam0 * anAxis->z());
 161
 162     aHeight = aParam1 - aParam0;
 163
 164     aCylinder = GeomCylinderPtr(new GeomAPI_Cylinder(aLocation, anAxis, aRadius, aHeight));
 165   }
 166   return aCylinder;
 167 }
 168
 169 //==================================================================================================
 170 std::shared_ptr<GeomAPI_Cone> GeomAPI_Solid::getCone() const
 171 {
 172   bool isCone = false;
 173
 174   GeomPointPtr anApex;
 175   GeomDirPtr anAxis;
 176   double aSemiAngle = 0.0;
 177   double aHeight = 0.0;
 178
 179   GeomPlanePtr aCaps[2];
 180
 181   for (TopExp_Explorer anExp(impl<TopoDS_Shape>(), TopAbs_FACE); anExp.More(); anExp.Next()) {
 182     GeomFacePtr aFace(new GeomAPI_Face);
 183     aFace->setImpl(new TopoDS_Shape(anExp.Current()));
 184
 185     GeomConePtr aCurCone = aFace->getCone();
 186     if (aCurCone) {
 187       if (isCone) { // at least one conical face is found
 188         // check equal apexes
 189         if (anApex->distance(aCurCone->apex()) >= Precision::Confusion() ||
 190         // check semi-angle
 191             Abs(aSemiAngle - aCurCone->semiAngle() >= Precision::Confusion()) ||
 192         // check axes are collinear
 193             !anAxis->isParallel(aCurCone->axis())) {
 194           isCone = false;
 195           break;
 196         }
 197       }
 198       else { // first cone is found
 199         anApex = aCurCone->apex();
 200         if (anAxis) {
 201           // the plane is found => compare directions
 202           if (!anAxis->isParallel(aCurCone->axis()))
 203             break;
 204         }
 205         else
 206           anAxis = aCurCone->axis();
 207         aSemiAngle = aCurCone->semiAngle();
 208         aHeight = aCurCone->height();
 209         isCone = true;
 210       }
 211     }
 212     else {
 213       // check the face is planar
 214       bool isPlaneApplicable = false;
 215       GeomPlanePtr aCurPln = aFace->getPlane();
 216       if (aCurPln) {
 217         // verify the plane is already exists
 218         int aLastPlanIndex = 0;
 219         while (aLastPlanIndex < 2) {
 220           if (!aCaps[aLastPlanIndex]) {
 221             // add new plane
 222             aCaps[aLastPlanIndex] = aCurPln;
 223             break;
 224           }
 225           if (aCaps[aLastPlanIndex]->isCoincident(aCurPln))
 226             break;
 227           ++aLastPlanIndex;
 228         }
 229
 230         isPlaneApplicable = aLastPlanIndex < 2;
 231       }
 232
 233       if (isPlaneApplicable) {
 234         if (!anAxis) // no cone is found, store the normal as further cone's axis
 235           anAxis = aCurPln->direction();
 236       }
 237       else {
 238         isCone = false;
 239         break;
 240       }
 241     }
 242   }
 243
 244   isCone = isCone && aCaps[0] && aCaps[0]->direction()->isParallel(anAxis);
 245   if (isCone && aCaps[1]) // cone map have only one cap, if it is bounded by the apex
 246     isCone = aCaps[1]->direction()->isParallel(anAxis);
 247
 248   GeomConePtr aCone;
 249   if (isCone) {
 250     // intersect planes with cone's axis
 251     std::shared_ptr<GeomAPI_XYZ> anAxisXYZ = anAxis->xyz();
 252     std::shared_ptr<GeomAPI_XYZ> anApexXYZ = anApex->xyz();
 253     double aParam0 = anAxisXYZ->dot(aCaps[0]->location()->xyz()->decreased(anApexXYZ));
 254     double aParam1 =
 255         aCaps[1] ? anAxisXYZ->dot(aCaps[1]->location()->xyz()->decreased(anApexXYZ)) : 0.0;
 256     if (aParam0 <= 0.0 && aParam1 <= 0.0) {
 257       // reverse axis to make smaller cap be the first
 258       anAxis->reverse();
 259       aParam0 = -aParam0;
 260       aParam1 = -aParam1;
 261     }
 262     if (aParam0 > aParam1 + Precision::Confusion()) {
 263       double tmp = aParam0;
 264       aParam0 = aParam1;
 265       aParam1 = tmp;
 266     }
 267
 268     // calculate location of cone to be coincident with one of planes
 269     GeomPointPtr aLocation(new GeomAPI_Pnt(
 270         anApex->x() + aParam0 * anAxis->x(),
 271         anApex->y() + aParam0 * anAxis->y(),
 272         anApex->z() + aParam0 * anAxis->z()));
 273
 274     // calculate radii of caps
 275     aParam0 *= Tan(aSemiAngle);
 276     aParam1 *= Tan(aSemiAngle);
 277
 278     aCone = GeomConePtr(new GeomAPI_Cone(aLocation, anAxis, aSemiAngle, aParam0, aParam1));
 279   }
 280   return aCone;
 281 }
 282
 283 //==================================================================================================
 284 std::shared_ptr<GeomAPI_Torus> GeomAPI_Solid::getTorus() const
 285 {
 286   GeomTorusPtr aTorus;
 287   ListOfShape aShells = subShapes(SHELL);
 288   if (aShells.size() == 1)
 289     aTorus = aShells.front()->shell()->getTorus();
 290   return aTorus;
 291 }
 292
 293 //==================================================================================================
 294 std::shared_ptr<GeomAPI_Box> GeomAPI_Solid::getParallelepiped() const
 295 {
 296   GeomBoxPtr aBox;
 297   ListOfShape aShells = subShapes(SHELL);
 298   if (aShells.size() == 1)
 299     aBox = aShells.front()->shell()->getParallelepiped();
 300   return aBox;
 301 }
 302
 303 //==================================================================================================
 304 GeomPointPtr GeomAPI_Solid::middlePoint() const
 305 {
 306   GeomPointPtr anInnerPoint;
 307
 308   const TopoDS_Solid& aSolid = impl<TopoDS_Solid>();
 309   if (aSolid.IsNull())
 310     return anInnerPoint;
 311
 312   GProp_GProps aProps;
 313   BRepGProp::VolumeProperties(aSolid, aProps, 1.e-4);
 314
 315   gp_Pnt aPnt = aProps.CentreOfMass();
 316   anInnerPoint = GeomPointPtr(new GeomAPI_Pnt(aPnt.X(), aPnt.Y(), aPnt.Z()));
 317   return anInnerPoint;
 318 }