src/GeomAPI/GeomAPI_Solid.cpp

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