Geant4 11.4.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4BooleanSolid Class Reference

G4BooleanSolid is the base class for solids created by Boolean operations between other solids. More...

#include <G4BooleanSolid.hh>

Inheritance diagram for G4BooleanSolid:

Public Member Functions

 G4BooleanSolid (const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB)
 G4BooleanSolid (const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB, G4RotationMatrix *rotMatrix, const G4ThreeVector &transVector)
 G4BooleanSolid (const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB, const G4Transform3D &transform)
 ~G4BooleanSolid () override
 G4BooleanSolid (__void__ &)
 G4BooleanSolid (const G4BooleanSolid &rhs)
G4BooleanSolidoperator= (const G4BooleanSolid &rhs)
const G4VSolidGetConstituentSolid (G4int no) const override
G4VSolidGetConstituentSolid (G4int no) override
G4double GetCubicVolume () override
G4double GetSurfaceArea () override
G4GeometryType GetEntityType () const override
G4PolyhedronGetPolyhedron () const override
std::ostream & StreamInfo (std::ostream &os) const override
G4int GetCubVolStatistics () const
void SetCubVolStatistics (G4int st)
G4double GetCubVolEpsilon () const
void SetCubVolEpsilon (G4double ep)
G4int GetAreaStatistics () const
void SetAreaStatistics (G4int st)
G4double GetAreaAccuracy () const
void SetAreaAccuracy (G4double ep)
G4ThreeVector GetPointOnSurface () const override
G4int GetNumOfConstituents () const override
G4bool IsFaceted () const override
Public Member Functions inherited from G4VSolid
 G4VSolid (const G4String &name)
virtual ~G4VSolid ()
 G4VSolid (const G4VSolid &rhs)
G4VSolidoperator= (const G4VSolid &rhs)
G4bool operator== (const G4VSolid &s) const
G4String GetName () const
void SetName (const G4String &name)
G4double GetTolerance () const
virtual void BoundingLimits (G4ThreeVector &pMin, G4ThreeVector &pMax) const
virtual G4bool CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
virtual EInside Inside (const G4ThreeVector &p) const =0
virtual G4ThreeVector SurfaceNormal (const G4ThreeVector &p) const =0
virtual G4double DistanceToIn (const G4ThreeVector &p, const G4ThreeVector &v) const =0
virtual G4double DistanceToIn (const G4ThreeVector &p) const =0
virtual G4double DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=nullptr, G4ThreeVector *n=nullptr) const =0
virtual G4double DistanceToOut (const G4ThreeVector &p) const =0
virtual void ComputeDimensions (G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
virtual G4VSolidClone () const
void DumpInfo () const
virtual void DescribeYourselfTo (G4VGraphicsScene &scene) const =0
virtual G4VisExtent GetExtent () const
virtual G4PolyhedronCreatePolyhedron () const
virtual const G4DisplacedSolidGetDisplacedSolidPtr () const
virtual G4DisplacedSolidGetDisplacedSolidPtr ()
 G4VSolid (__void__ &)
G4double EstimateCubicVolume (G4int nStat, G4double epsilon) const
G4double EstimateSurfaceArea (G4int nStat, G4double epsilon) const

Static Public Member Functions

static G4VBooleanProcessorGetExternalBooleanProcessor ()
static void SetExternalBooleanProcessor (G4VBooleanProcessor *extProcessor)

Protected Member Functions

void GetListOfPrimitives (std::vector< std::pair< G4VSolid *, G4Transform3D > > &, const G4Transform3D &) const
G4PolyhedronStackPolyhedron (HepPolyhedronProcessor &, const G4VSolid *) const
Protected Member Functions inherited from G4VSolid
void CalculateClippedPolygonExtent (G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
void ClipCrossSection (G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
void ClipBetweenSections (G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
void ClipPolygon (G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis) const

Protected Attributes

G4VSolidfPtrSolidA = nullptr
G4VSolidfPtrSolidB = nullptr
G4double fCubicVolume = -1.0
G4double fSurfaceArea = -1.0
Protected Attributes inherited from G4VSolid
G4double kCarTolerance

Static Protected Attributes

static G4VBooleanProcessorfExternalBoolProcessor = nullptr

Detailed Description

G4BooleanSolid is the base class for solids created by Boolean operations between other solids.

Definition at line 52 of file G4BooleanSolid.hh.

Constructor & Destructor Documentation

◆ G4BooleanSolid() [1/5]

G4BooleanSolid::G4BooleanSolid ( const G4String & pName,
G4VSolid * pSolidA,
G4VSolid * pSolidB )

◆ G4BooleanSolid() [2/5]

G4BooleanSolid::G4BooleanSolid ( const G4String & pName,
G4VSolid * pSolidA,
G4VSolid * pSolidB,
G4RotationMatrix * rotMatrix,
const G4ThreeVector & transVector )

Constructor of a Boolean composition between two solids with rotation and translation, used to transform the coordinate system of the second solid to the coordinate system of the first solid.

Parameters
[in]pNameThe name of the Boolean composition.
[in]pSolidAPointer to the first reference solid.
[in]pSolidBPointer to the second solid to form the composition.
[in]rotMatrixPointer to the rotation vector.
[in]transVectorThe translation vector.

Definition at line 66 of file G4BooleanSolid.cc.

71 : G4VSolid(pName), createdDisplacedSolid(true)
72{
73 fPtrSolidA = pSolidA ;
74 fPtrSolidB = new G4DisplacedSolid("placedB",pSolidB,rotMatrix,transVector) ;
75}

◆ G4BooleanSolid() [3/5]

G4BooleanSolid::G4BooleanSolid ( const G4String & pName,
G4VSolid * pSolidA,
G4VSolid * pSolidB,
const G4Transform3D & transform )

Constructor of a Boolean composition between two solids with a transformation that moves the second solid from its desired position to its standard position.

Parameters
[in]pNameThe name of the Boolean composition.
[in]pSolidAPointer to the first reference solid.
[in]pSolidBPointer to the second solid to form the composition.
[in]transformThe composed 3D transformation.

Definition at line 81 of file G4BooleanSolid.cc.

85 : G4VSolid(pName), createdDisplacedSolid(true)
86{
87 fPtrSolidA = pSolidA ;
88 fPtrSolidB = new G4DisplacedSolid("placedB",pSolidB,transform) ;
89}

◆ ~G4BooleanSolid()

G4BooleanSolid::~G4BooleanSolid ( )
override

Destructor. If using a displaced solid, deletes all cached transformations.

Definition at line 105 of file G4BooleanSolid.cc.

106{
107 if(createdDisplacedSolid)
108 {
109 ((G4DisplacedSolid*)fPtrSolidB)->CleanTransformations();
110 }
111 delete fpPolyhedron; fpPolyhedron = nullptr;
112}

◆ G4BooleanSolid() [4/5]

G4BooleanSolid::G4BooleanSolid ( __void__ & a)

Fake default constructor for usage restricted to direct object persistency for clients requiring preallocation of memory for persistifiable objects.

Definition at line 96 of file G4BooleanSolid.cc.

97 : G4VSolid(a)
98{
99}

◆ G4BooleanSolid() [5/5]

G4BooleanSolid::G4BooleanSolid ( const G4BooleanSolid & rhs)

Copy constructor and assignment operator.

Definition at line 118 of file G4BooleanSolid.cc.

121 fCubVolStatistics(rhs.fCubVolStatistics),
122 fAreaStatistics(rhs.fAreaStatistics),
123 fCubVolEpsilon(rhs.fCubVolEpsilon),
124 fAreaAccuracy(rhs.fAreaAccuracy),
125 createdDisplacedSolid(rhs.createdDisplacedSolid)
126{
127 fPrimitives.resize(0); fPrimitivesSurfaceArea = 0.;
128}

Member Function Documentation

◆ GetAreaAccuracy()

G4double G4BooleanSolid::GetAreaAccuracy ( ) const
inline

Accessor and setter for controlling/tuning the level of accuracy used for computing the surface area.

◆ GetAreaStatistics()

G4int G4BooleanSolid::GetAreaStatistics ( ) const
inline

Accessor and setter for controlling/tuning the number of random points to be used for computing the surface area.

◆ GetConstituentSolid() [1/2]

const G4VSolid * G4BooleanSolid::GetConstituentSolid ( G4int no) const
overridevirtual

Methods returning the component solids of the Boolean composition. If the solid is made up from a Boolean operation of two solids, return the corresponding solid (for no=0 and 1). A fatal exception is thrown if the index provided is different from 0 or 1.

Parameters
[in]noIndex 0/1 of the components.
Returns
The pointer to (const or not const) of the component solid. If the solid is not a "Boolean", returns nullptr.

Reimplemented from G4VSolid.

Definition at line 165 of file G4BooleanSolid.cc.

166{
167 const G4VSolid* subSolid = nullptr;
168 if( no == 0 )
169 {
170 subSolid = fPtrSolidA;
171 }
172 else if( no == 1 )
173 {
174 subSolid = fPtrSolidB;
175 }
176 else
177 {
178 DumpInfo();
179 G4Exception("G4BooleanSolid::GetConstituentSolid()",
180 "GeomSolids0002", FatalException, "Invalid solid index.");
181 }
182 return subSolid;
183}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
void DumpInfo() const

Referenced by G4tgbGeometryDumper::DumpBooleanVolume().

◆ GetConstituentSolid() [2/2]

G4VSolid * G4BooleanSolid::GetConstituentSolid ( G4int no)
overridevirtual

Reimplemented from G4VSolid.

Definition at line 191 of file G4BooleanSolid.cc.

192{
193 G4VSolid* subSolid = nullptr;
194 if( no == 0 )
195 {
196 subSolid = fPtrSolidA;
197 }
198 else if( no == 1 )
199 {
200 subSolid = fPtrSolidB;
201 }
202 else
203 {
204 DumpInfo();
205 G4Exception("G4BooleanSolid::GetConstituentSolid()",
206 "GeomSolids0002", FatalException, "Invalid solid index.");
207 }
208 return subSolid;
209}

◆ GetCubicVolume()

G4double G4BooleanSolid::GetCubicVolume ( )
overridevirtual

Methods returning the computed capacity and surface area of the composition. The quantities returned are an estimate obtained by randomly sampling the Boolean composition and caching them for reuse.

Reimplemented from G4VSolid.

Reimplemented in G4SubtractionSolid, and G4UnionSolid.

Definition at line 590 of file G4BooleanSolid.cc.

591{
592 if(fCubicVolume < 0.)
593 {
594 G4AutoLock l(&boolMutex);
595 fCubicVolume = EstimateCubicVolume(fCubVolStatistics, fCubVolEpsilon);
596 l.unlock();
597 }
598 return fCubicVolume;
599}
G4TemplateAutoLock< G4Mutex > G4AutoLock
G4double EstimateCubicVolume(G4int nStat, G4double epsilon) const
Definition G4VSolid.cc:229

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ GetCubVolEpsilon()

G4double G4BooleanSolid::GetCubVolEpsilon ( ) const
inline

Accessor and setter for controlling/tuning the epsilon used for computing the cubic volume.

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ GetCubVolStatistics()

G4int G4BooleanSolid::GetCubVolStatistics ( ) const
inline

Accessor and setter for controlling/tuning the number of random points to be used for computing the cubic volume.

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ GetEntityType()

G4GeometryType G4BooleanSolid::GetEntityType ( ) const
overridevirtual

Returns the type ID, "G4BooleanSolid" of the solid.

Implements G4VSolid.

Reimplemented in G4IntersectionSolid, G4SubtractionSolid, and G4UnionSolid.

Definition at line 215 of file G4BooleanSolid.cc.

216{
217 return {"G4BooleanSolid"};
218}

Referenced by StreamInfo().

◆ GetExternalBooleanProcessor()

G4VBooleanProcessor * G4BooleanSolid::GetExternalBooleanProcessor ( )
static

Gets/sets the Boolean processor for polyhedron to replace the default processor.

Definition at line 630 of file G4BooleanSolid.cc.

631{
633}
static G4VBooleanProcessor * fExternalBoolProcessor

Referenced by G4MultiUnion::CreatePolyhedron().

◆ GetListOfPrimitives()

void G4BooleanSolid::GetListOfPrimitives ( std::vector< std::pair< G4VSolid *, G4Transform3D > > & primitives,
const G4Transform3D & curPlacement ) const
protected

Gets the list of constituent primitives of the solid and their placements.

Definition at line 383 of file G4BooleanSolid.cc.

386{
387 G4Transform3D transform;
388 G4VSolid* solid;
389 G4String type;
390
391 // Repeat two times, first time for fPtrSolidA and then for fPtrSolidB
392 //
393 for (auto i=0; i<2; ++i)
394 {
395 transform = curPlacement;
396 solid = (i == 0) ? fPtrSolidA : fPtrSolidB;
397 type = solid->GetEntityType();
398
399 // While current solid is a trasformed solid just modify transform
400 //
401 while (type == "G4DisplacedSolid" ||
402 type == "G4ReflectedSolid" ||
403 type == "G4ScaledSolid")
404 {
405 if (type == "G4DisplacedSolid")
406 {
407 transform = transform * G4Transform3D(
408 ((G4DisplacedSolid*)solid)->GetObjectRotation(),
409 ((G4DisplacedSolid*)solid)->GetObjectTranslation());
410 solid = ((G4DisplacedSolid*)solid)->GetConstituentMovedSolid();
411 }
412 else if (type == "G4ReflectedSolid")
413 {
414 transform= transform*((G4ReflectedSolid*)solid)->GetDirectTransform3D();
415 solid = ((G4ReflectedSolid*)solid)->GetConstituentMovedSolid();
416 }
417 else if (type == "G4ScaledSolid")
418 {
419 transform = transform * ((G4ScaledSolid*)solid)->GetScaleTransform();
420 solid = ((G4ScaledSolid*)solid)->GetUnscaledSolid();
421 }
422 type = solid->GetEntityType();
423 }
424
425 // If current solid is a Boolean solid then continue recursion,
426 // otherwise add it to the list of primitives
427 //
428 if (type == "G4UnionSolid" ||
429 type == "G4SubtractionSolid" ||
430 type == "G4IntersectionSolid" ||
431 type == "G4BooleanSolid")
432 {
433 ((G4BooleanSolid *)solid)->GetListOfPrimitives(primitives,transform);
434 }
435 else
436 {
437 primitives.emplace_back(solid,transform);
438 }
439 }
440}
HepGeom::Transform3D G4Transform3D
G4BooleanSolid(const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB)
G4VSolid * GetConstituentMovedSolid() const
G4GeometryType GetEntityType() const override
virtual G4GeometryType GetEntityType() const =0

Referenced by GetPointOnSurface().

◆ GetNumOfConstituents()

G4int G4BooleanSolid::GetNumOfConstituents ( ) const
overridevirtual

Returns the total number of constituent solids forming the Boolean composition.

Reimplemented from G4VSolid.

Definition at line 496 of file G4BooleanSolid.cc.

497{
498 return (fPtrSolidA->GetNumOfConstituents() + fPtrSolidB->GetNumOfConstituents());
499}

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ GetPointOnSurface()

G4ThreeVector G4BooleanSolid::GetPointOnSurface ( ) const
overridevirtual

Returns a point (G4ThreeVector) randomly and uniformly generated on the surface of the solid.

Reimplemented from G4VSolid.

Definition at line 447 of file G4BooleanSolid.cc.

448{
449 std::size_t nprims = fPrimitives.size();
450 std::pair<G4VSolid *, G4Transform3D> prim;
451
452 // Get list of primitives and find the total area of their surfaces
453 //
454 if (nprims == 0)
455 {
456 GetListOfPrimitives(fPrimitives, G4Transform3D());
457 nprims = fPrimitives.size();
458 fPrimitivesSurfaceArea = 0.;
459 for (std::size_t i=0; i<nprims; ++i)
460 {
461 fPrimitivesSurfaceArea += fPrimitives[i].first->GetSurfaceArea();
462 }
463 }
464
465 // Select random primitive, get random point on its surface and
466 // check that the point belongs to the surface of the solid
467 //
469 for (std::size_t k=0; k<100000; ++k) // try 100k times
470 {
471 G4double rand = fPrimitivesSurfaceArea * G4QuickRand();
472 G4double area = 0.;
473 for (std::size_t i=0; i<nprims; ++i)
474 {
475 prim = fPrimitives[i];
476 area += prim.first->GetSurfaceArea();
477 if (rand < area) { break; }
478 }
479 p = prim.first->GetPointOnSurface();
480 p = prim.second * G4Point3D(p);
481 if (Inside(p) == kSurface) { return p; }
482 }
483 std::ostringstream message;
484 message << "Solid - " << GetName() << "\n"
485 << "All 100k attempts to generate a point on the surface have failed!\n"
486 << "The solid created may be an invalid Boolean construct!";
487 G4Exception("G4BooleanSolid::GetPointOnSurface()",
488 "GeomSolids1001", JustWarning, message);
489 return p;
490}
@ JustWarning
HepGeom::Point3D< G4double > G4Point3D
Definition G4Point3D.hh:34
G4double G4QuickRand(uint32_t seed=0)
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition G4Types.hh:83
void GetListOfPrimitives(std::vector< std::pair< G4VSolid *, G4Transform3D > > &, const G4Transform3D &) const
G4String GetName() const
virtual EInside Inside(const G4ThreeVector &p) const =0
@ kSurface
Definition geomdefs.hh:69

◆ GetPolyhedron()

G4Polyhedron * G4BooleanSolid::GetPolyhedron ( ) const
overridevirtual

Returns a pointer to the generated polyhedron representation of the composition, for use in visualisation.

Reimplemented from G4VSolid.

Definition at line 514 of file G4BooleanSolid.cc.

515{
516 if (fpPolyhedron == nullptr ||
517 fRebuildPolyhedron ||
518 fpPolyhedron->GetNumberOfRotationStepsAtTimeOfCreation() !=
519 fpPolyhedron->GetNumberOfRotationSteps())
520 {
521 G4RecursiveAutoLock l(&polyhedronMutex);
522 delete fpPolyhedron;
523 fpPolyhedron = CreatePolyhedron();
524 fRebuildPolyhedron = false;
525 l.unlock();
526 }
527 return fpPolyhedron;
528}
G4TemplateAutoLock< G4RecursiveMutex > G4RecursiveAutoLock
virtual G4Polyhedron * CreatePolyhedron() const
Definition G4VSolid.cc:726

◆ GetSurfaceArea()

G4double G4BooleanSolid::GetSurfaceArea ( )
overridevirtual

Returns an estimation of the solid surface area in internal units. This method may be overloaded by derived classes to compute the exact geometrical quantity for solids where this is possible, or anyway to cache the computed value. Note: the computed value is NOT cached.

Reimplemented from G4VSolid.

Definition at line 605 of file G4BooleanSolid.cc.

606{
607 if(fSurfaceArea < 0.)
608 {
609 G4AutoLock l(&boolMutex);
610 fSurfaceArea = EstimateSurfaceArea(fAreaStatistics, fAreaAccuracy);
611 l.unlock();
612 }
613 return fSurfaceArea;
614}
G4double EstimateSurfaceArea(G4int nStat, G4double epsilon) const
Definition G4VSolid.cc:290

◆ IsFaceted()

G4bool G4BooleanSolid::IsFaceted ( ) const
overridevirtual

Return true if the resulting solid has only planar faces.

Reimplemented from G4VSolid.

Definition at line 505 of file G4BooleanSolid.cc.

506{
507 return (fPtrSolidA->IsFaceted() && fPtrSolidB->IsFaceted());
508}

◆ operator=()

G4BooleanSolid & G4BooleanSolid::operator= ( const G4BooleanSolid & rhs)

Definition at line 134 of file G4BooleanSolid.cc.

135{
136 // Check assignment to self
137 //
138 if (this == &rhs) { return *this; }
139
140 // Copy base class data
141 //
143
144 // Copy data
145 //
148 fCubVolStatistics = rhs.fCubVolStatistics; fCubVolEpsilon = rhs.fCubVolEpsilon;
149 fAreaStatistics = rhs.fAreaStatistics; fAreaAccuracy = rhs.fAreaAccuracy;
150 createdDisplacedSolid= rhs.createdDisplacedSolid;
151
152 fRebuildPolyhedron = false;
153 delete fpPolyhedron; fpPolyhedron = nullptr;
154 fPrimitives.resize(0); fPrimitivesSurfaceArea = 0.;
155
156 return *this;
157}
G4VSolid & operator=(const G4VSolid &rhs)
Definition G4VSolid.cc:108

Referenced by G4IntersectionSolid::operator=(), G4SubtractionSolid::operator=(), and G4UnionSolid::operator=().

◆ SetAreaAccuracy()

void G4BooleanSolid::SetAreaAccuracy ( G4double ep)
inline

◆ SetAreaStatistics()

void G4BooleanSolid::SetAreaStatistics ( G4int st)
inline

◆ SetCubVolEpsilon()

void G4BooleanSolid::SetCubVolEpsilon ( G4double ep)

Definition at line 294 of file G4BooleanSolid.cc.

295{
296 if (ep != fCubVolEpsilon) { fCubicVolume = -1.; }
297 fCubVolEpsilon = ep;
298
299 // Propagate ep to all components of the 1st solid
300 if (fPtrSolidA->GetNumOfConstituents() > 1)
301 {
302 G4VSolid* ptr = fPtrSolidA;
303 while(true)
304 {
305 G4String type = ptr->GetEntityType();
306 if (type == "G4DisplacedSolid")
307 {
308 ptr = ((G4DisplacedSolid*)ptr)->GetConstituentMovedSolid();
309 continue;
310 }
311 if (type == "G4ReflectedSolid")
312 {
313 ptr = ((G4ReflectedSolid*)ptr)->GetConstituentMovedSolid();
314 continue;
315 }
316 if (type == "G4ScaledSolid")
317 {
318 ptr = ((G4ScaledSolid*)ptr)->GetUnscaledSolid();
319 continue;
320 }
321 if (type != "G4MultiUnion") // G4MultiUnion doesn't have SetCubVolEpsilon()
322 {
323 ((G4BooleanSolid*)ptr)->SetCubVolEpsilon(ep);
324 }
325 break;
326 }
327 }
328
329 // Propagate ep to all components of the 2nd solid
330 if (fPtrSolidB->GetNumOfConstituents() > 1)
331 {
332 G4VSolid* ptr = fPtrSolidB;
333 while(true)
334 {
335 G4String type = ptr->GetEntityType();
336 if (type == "G4DisplacedSolid")
337 {
338 ptr = ((G4DisplacedSolid*)ptr)->GetConstituentMovedSolid();
339 continue;
340 }
341 if (type == "G4ReflectedSolid")
342 {
343 ptr = ((G4ReflectedSolid*)ptr)->GetConstituentMovedSolid();
344 continue;
345 }
346 if (type == "G4ScaledSolid")
347 {
348 ptr = ((G4ScaledSolid*)ptr)->GetUnscaledSolid();
349 continue;
350 }
351 if (type != "G4MultiUnion") // G4MultiUnion doesn't have SetCubVolEpsilon()
352 {
353 ((G4BooleanSolid*)ptr)->SetCubVolEpsilon(ep);
354 }
355 break;
356 }
357 }
358}

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ SetCubVolStatistics()

void G4BooleanSolid::SetCubVolStatistics ( G4int st)

Definition at line 224 of file G4BooleanSolid.cc.

225{
226 if (st != fCubVolStatistics) { fCubicVolume = -1.; }
227 fCubVolStatistics = st;
228
229 // Propagate st to all components of the 1st solid
230 if (fPtrSolidA->GetNumOfConstituents() > 1)
231 {
232 G4VSolid* ptr = fPtrSolidA;
233 while(true)
234 {
235 G4String type = ptr->GetEntityType();
236 if (type == "G4DisplacedSolid")
237 {
238 ptr = ((G4DisplacedSolid*)ptr)->GetConstituentMovedSolid();
239 continue;
240 }
241 if (type == "G4ReflectedSolid")
242 {
243 ptr = ((G4ReflectedSolid*)ptr)->GetConstituentMovedSolid();
244 continue;
245 }
246 if (type == "G4ScaledSolid")
247 {
248 ptr = ((G4ScaledSolid*)ptr)->GetUnscaledSolid();
249 continue;
250 }
251 if (type != "G4MultiUnion") // G4MultiUnion doesn't have SetCubVolStatistics()
252 {
253 ((G4BooleanSolid*)ptr)->SetCubVolStatistics(st);
254 }
255 break;
256 }
257 }
258
259 // Propagate st to all components of the 2nd solid
260 if (fPtrSolidB->GetNumOfConstituents() > 1)
261 {
262 G4VSolid* ptr = fPtrSolidB;
263 while(true)
264 {
265 G4String type = ptr->GetEntityType();
266 if (type == "G4DisplacedSolid")
267 {
268 ptr = ((G4DisplacedSolid*)ptr)->GetConstituentMovedSolid();
269 continue;
270 }
271 if (type == "G4ReflectedSolid")
272 {
273 ptr = ((G4ReflectedSolid*)ptr)->GetConstituentMovedSolid();
274 continue;
275 }
276 if (type == "G4ScaledSolid")
277 {
278 ptr = ((G4ScaledSolid*)ptr)->GetUnscaledSolid();
279 continue;
280 }
281 if (type != "G4MultiUnion") // G4MultiUnion doesn't have SetCubVolStatistics()
282 {
283 ((G4BooleanSolid*)ptr)->SetCubVolStatistics(st);
284 }
285 break;
286 }
287 }
288}

Referenced by G4SubtractionSolid::GetCubicVolume(), and G4UnionSolid::GetCubicVolume().

◆ SetExternalBooleanProcessor()

void G4BooleanSolid::SetExternalBooleanProcessor ( G4VBooleanProcessor * extProcessor)
static

Definition at line 621 of file G4BooleanSolid.cc.

622{
623 fExternalBoolProcessor = extProcessor;
624}

◆ StackPolyhedron()

G4Polyhedron * G4BooleanSolid::StackPolyhedron ( HepPolyhedronProcessor & processor,
const G4VSolid * solid ) const
protected

Stacks the polyhedra for processing.

Returns
A pointer to the top polyhedron.

Definition at line 535 of file G4BooleanSolid.cc.

537{
539 const G4String& type = solid->GetEntityType();
540 if (type == "G4UnionSolid")
541 { operation = HepPolyhedronProcessor::UNION; }
542 else if (type == "G4IntersectionSolid")
544 else if (type == "G4SubtractionSolid")
546 else
547 {
548 std::ostringstream message;
549 message << "Solid - " << solid->GetName()
550 << " - Unrecognised composite solid" << G4endl
551 << " Returning NULL !";
552 G4Exception("StackPolyhedron()", "GeomSolids1001", JustWarning, message);
553 return nullptr;
554 }
555
556 G4Polyhedron* top = nullptr;
557 const G4VSolid* solidA = solid->GetConstituentSolid(0);
558 const G4VSolid* solidB = solid->GetConstituentSolid(1);
559
560 if (solidA->GetConstituentSolid(0) != nullptr)
561 {
562 top = StackPolyhedron(processor, solidA);
563 }
564 else
565 {
566 top = solidA->GetPolyhedron();
567 }
568 G4Polyhedron* operand = solidB->GetPolyhedron();
569 if (operand != nullptr)
570 {
571 processor.push_back (operation, *operand);
572 }
573 else
574 {
575 std::ostringstream message;
576 message << "Solid - " << solid->GetName()
577 << " - No G4Polyhedron for Boolean component";
578 G4Exception("G4BooleanSolid::StackPolyhedron()",
579 "GeomSolids2001", JustWarning, message);
580 }
581
582 return top;
583}
#define G4endl
Definition G4ios.hh:67
G4Polyhedron * StackPolyhedron(HepPolyhedronProcessor &, const G4VSolid *) const
virtual const G4VSolid * GetConstituentSolid(G4int no) const
Definition G4VSolid.cc:185
virtual G4Polyhedron * GetPolyhedron() const
Definition G4VSolid.cc:731
void push_back(Operation, const HepPolyhedron &)

Referenced by G4IntersectionSolid::CreatePolyhedron(), G4SubtractionSolid::CreatePolyhedron(), G4UnionSolid::CreatePolyhedron(), and StackPolyhedron().

◆ StreamInfo()

std::ostream & G4BooleanSolid::StreamInfo ( std::ostream & os) const
overridevirtual

Streams the object contents to an output stream.

Implements G4VSolid.

Definition at line 364 of file G4BooleanSolid.cc.

365{
366 os << "-----------------------------------------------------------\n"
367 << " *** Dump for Boolean solid - " << GetName() << " ***\n"
368 << " ===================================================\n"
369 << " Solid type: " << GetEntityType() << "\n"
370 << " Parameters of constituent solids: \n"
371 << "===========================================================\n";
372 fPtrSolidA->StreamInfo(os);
373 fPtrSolidB->StreamInfo(os);
374 os << "===========================================================\n";
375
376 return os;
377}
G4GeometryType GetEntityType() const override

Member Data Documentation

◆ fCubicVolume

G4double G4BooleanSolid::fCubicVolume = -1.0
protected

◆ fExternalBoolProcessor

G4VBooleanProcessor * G4BooleanSolid::fExternalBoolProcessor = nullptr
staticprotected

◆ fPtrSolidA

◆ fPtrSolidB

◆ fSurfaceArea

G4double G4BooleanSolid::fSurfaceArea = -1.0
protected

Cached value of the surface area.

Definition at line 230 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid(), GetSurfaceArea(), and operator=().


The documentation for this class was generated from the following files: