Geant4-11
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes
G4IntersectionSolid Class Reference

#include <G4IntersectionSolid.hh>

Inheritance diagram for G4IntersectionSolid:
G4BooleanSolid G4VSolid

Public Member Functions

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

Protected Member Functions

void CalculateClippedPolygonExtent (G4ThreeVectorList &pPolygon, 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 ClipCrossSection (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
 
void GetListOfPrimitives (std::vector< std::pair< G4VSolid *, G4Transform3D > > &, const G4Transform3D &) const
 
G4PolyhedronStackPolyhedron (HepPolyhedronProcessor &, const G4VSolid *) const
 

Protected Attributes

G4double fCubicVolume = -1.0
 
G4VSolidfPtrSolidA = nullptr
 
G4VSolidfPtrSolidB = nullptr
 
G4double kCarTolerance
 

Private Member Functions

void ClipPolygonToSimpleLimits (G4ThreeVectorList &pPolygon, G4ThreeVectorList &outputPolygon, const G4VoxelLimits &pVoxelLimit) const
 

Private Attributes

G4bool createdDisplacedSolid = false
 
G4double fAreaAccuracy = -1
 
G4double fCubVolEpsilon = 0.001
 
G4PolyhedronfpPolyhedron = nullptr
 
std::vector< std::pair< G4VSolid *, G4Transform3D > > fPrimitives
 
G4double fPrimitivesSurfaceArea = 0.0
 
G4bool fRebuildPolyhedron = false
 
G4String fshapeName
 
G4int fStatistics = 1000000
 
G4double fSurfaceArea = -1.0
 

Detailed Description

Definition at line 45 of file G4IntersectionSolid.hh.

Constructor & Destructor Documentation

◆ G4IntersectionSolid() [1/5]

G4IntersectionSolid::G4IntersectionSolid ( const G4String pName,
G4VSolid pSolidA,
G4VSolid pSolidB 
)

Definition at line 49 of file G4IntersectionSolid.cc.

52 : G4BooleanSolid(pName,pSolidA,pSolidB)
53{
54}
G4BooleanSolid(const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB)

Referenced by Clone().

◆ G4IntersectionSolid() [2/5]

G4IntersectionSolid::G4IntersectionSolid ( const G4String pName,
G4VSolid pSolidA,
G4VSolid pSolidB,
G4RotationMatrix rotMatrix,
const G4ThreeVector transVector 
)

Definition at line 59 of file G4IntersectionSolid.cc.

64 : G4BooleanSolid(pName,pSolidA,pSolidB,rotMatrix,transVector)
65{
66}

◆ G4IntersectionSolid() [3/5]

G4IntersectionSolid::G4IntersectionSolid ( const G4String pName,
G4VSolid pSolidA,
G4VSolid pSolidB,
const G4Transform3D transform 
)

Definition at line 72 of file G4IntersectionSolid.cc.

76 : G4BooleanSolid(pName,pSolidA,pSolidB,transform)
77{
78}
G4bool transform(G4String &input, const G4String &type)

◆ ~G4IntersectionSolid()

G4IntersectionSolid::~G4IntersectionSolid ( )
virtual

Definition at line 94 of file G4IntersectionSolid.cc.

95{
96}

◆ G4IntersectionSolid() [4/5]

G4IntersectionSolid::G4IntersectionSolid ( __void__ &  a)

Definition at line 85 of file G4IntersectionSolid.cc.

87{
88}

◆ G4IntersectionSolid() [5/5]

G4IntersectionSolid::G4IntersectionSolid ( const G4IntersectionSolid rhs)

Definition at line 102 of file G4IntersectionSolid.cc.

103 : G4BooleanSolid (rhs)
104{
105}

Member Function Documentation

◆ BoundingLimits()

void G4IntersectionSolid::BoundingLimits ( G4ThreeVector pMin,
G4ThreeVector pMax 
) const
virtual

Reimplemented from G4VSolid.

Definition at line 130 of file G4IntersectionSolid.cc.

132{
133 G4ThreeVector minA,maxA, minB,maxB;
135 fPtrSolidB->BoundingLimits(minB,maxB);
136
137 pMin.set(std::max(minA.x(),minB.x()),
138 std::max(minA.y(),minB.y()),
139 std::max(minA.z(),minB.z()));
140
141 pMax.set(std::min(maxA.x(),maxB.x()),
142 std::min(maxA.y(),maxB.y()),
143 std::min(maxA.z(),maxB.z()));
144
145 // Check correctness of the bounding box
146 //
147 if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
148 {
149 std::ostringstream message;
150 message << "Bad bounding box (min >= max) for solid: "
151 << GetName() << " !"
152 << "\npMin = " << pMin
153 << "\npMax = " << pMax;
154 G4Exception("G4IntersectionSolid::BoundingLimits()", "GeomMgt0001",
155 JustWarning, message);
156 DumpInfo();
157 }
158}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
static const G4int maxA
static const G4double pMax
static const G4double pMin
double z() const
double x() const
double y() const
G4VSolid * fPtrSolidA
G4VSolid * fPtrSolidB
G4String GetName() const
void DumpInfo() const
virtual void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:665
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

References G4VSolid::BoundingLimits(), G4VSolid::DumpInfo(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4Exception(), G4VSolid::GetName(), JustWarning, G4INCL::Math::max(), maxA, G4INCL::Math::min(), pMax, pMin, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

◆ CalculateClippedPolygonExtent()

void G4VSolid::CalculateClippedPolygonExtent ( G4ThreeVectorList pPolygon,
const G4VoxelLimits pVoxelLimit,
const EAxis  pAxis,
G4double pMin,
G4double pMax 
) const
protectedinherited

Definition at line 489 of file G4VSolid.cc.

494{
495 G4int noLeft,i;
496 G4double component;
497
498 ClipPolygon(pPolygon,pVoxelLimit,pAxis);
499 noLeft = pPolygon.size();
500
501 if ( noLeft )
502 {
503 for (i=0; i<noLeft; ++i)
504 {
505 component = pPolygon[i].operator()(pAxis);
506
507 if (component < pMin)
508 {
509 pMin = component;
510 }
511 if (component > pMax)
512 {
513 pMax = component;
514 }
515 }
516 }
517}
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
void ClipPolygon(G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis) const
Definition: G4VSolid.cc:539

References G4VSolid::ClipPolygon(), pMax, and pMin.

Referenced by G4VSolid::ClipBetweenSections(), and G4VSolid::ClipCrossSection().

◆ CalculateExtent()

G4bool G4IntersectionSolid::CalculateExtent ( const EAxis  pAxis,
const G4VoxelLimits pVoxelLimit,
const G4AffineTransform pTransform,
G4double pMin,
G4double pMax 
) const
virtual

Implements G4VSolid.

Definition at line 165 of file G4IntersectionSolid.cc.

170{
171 G4bool retA, retB, out;
172 G4double minA, minB, maxA, maxB;
173
174 retA = fPtrSolidA
175 ->CalculateExtent( pAxis, pVoxelLimit, pTransform, minA, maxA);
176 retB = fPtrSolidB
177 ->CalculateExtent( pAxis, pVoxelLimit, pTransform, minB, maxB);
178
179 if( retA && retB )
180 {
181 pMin = std::max( minA, minB );
182 pMax = std::min( maxA, maxB );
183 out = (pMax > pMin); // true;
184 }
185 else
186 {
187 out = false;
188 }
189
190 return out; // It exists in this slice only if both exist in it.
191}
bool G4bool
Definition: G4Types.hh:86
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0

References G4VSolid::CalculateExtent(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4INCL::Math::max(), maxA, G4INCL::Math::min(), pMax, and pMin.

◆ ClipBetweenSections()

void G4VSolid::ClipBetweenSections ( G4ThreeVectorList pVertices,
const G4int  pSectionIndex,
const G4VoxelLimits pVoxelLimit,
const EAxis  pAxis,
G4double pMin,
G4double pMax 
) const
protectedinherited

Definition at line 444 of file G4VSolid.cc.

449{
450 G4ThreeVectorList polygon;
451 polygon.reserve(4);
452 polygon.push_back((*pVertices)[pSectionIndex]);
453 polygon.push_back((*pVertices)[pSectionIndex+4]);
454 polygon.push_back((*pVertices)[pSectionIndex+5]);
455 polygon.push_back((*pVertices)[pSectionIndex+1]);
456 CalculateClippedPolygonExtent(polygon,pVoxelLimit,pAxis,pMin,pMax);
457 polygon.clear();
458
459 polygon.push_back((*pVertices)[pSectionIndex+1]);
460 polygon.push_back((*pVertices)[pSectionIndex+5]);
461 polygon.push_back((*pVertices)[pSectionIndex+6]);
462 polygon.push_back((*pVertices)[pSectionIndex+2]);
463 CalculateClippedPolygonExtent(polygon,pVoxelLimit,pAxis,pMin,pMax);
464 polygon.clear();
465
466 polygon.push_back((*pVertices)[pSectionIndex+2]);
467 polygon.push_back((*pVertices)[pSectionIndex+6]);
468 polygon.push_back((*pVertices)[pSectionIndex+7]);
469 polygon.push_back((*pVertices)[pSectionIndex+3]);
470 CalculateClippedPolygonExtent(polygon,pVoxelLimit,pAxis,pMin,pMax);
471 polygon.clear();
472
473 polygon.push_back((*pVertices)[pSectionIndex+3]);
474 polygon.push_back((*pVertices)[pSectionIndex+7]);
475 polygon.push_back((*pVertices)[pSectionIndex+4]);
476 polygon.push_back((*pVertices)[pSectionIndex]);
477 CalculateClippedPolygonExtent(polygon,pVoxelLimit,pAxis,pMin,pMax);
478 return;
479}
std::vector< G4ThreeVector > G4ThreeVectorList
void CalculateClippedPolygonExtent(G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
Definition: G4VSolid.cc:489

References G4VSolid::CalculateClippedPolygonExtent(), pMax, and pMin.

◆ ClipCrossSection()

void G4VSolid::ClipCrossSection ( G4ThreeVectorList pVertices,
const G4int  pSectionIndex,
const G4VoxelLimits pVoxelLimit,
const EAxis  pAxis,
G4double pMin,
G4double pMax 
) const
protectedinherited

Definition at line 414 of file G4VSolid.cc.

419{
420
421 G4ThreeVectorList polygon;
422 polygon.reserve(4);
423 polygon.push_back((*pVertices)[pSectionIndex]);
424 polygon.push_back((*pVertices)[pSectionIndex+1]);
425 polygon.push_back((*pVertices)[pSectionIndex+2]);
426 polygon.push_back((*pVertices)[pSectionIndex+3]);
427 CalculateClippedPolygonExtent(polygon,pVoxelLimit,pAxis,pMin,pMax);
428 return;
429}

References G4VSolid::CalculateClippedPolygonExtent(), pMax, and pMin.

◆ ClipPolygon()

void G4VSolid::ClipPolygon ( G4ThreeVectorList pPolygon,
const G4VoxelLimits pVoxelLimit,
const EAxis  pAxis 
) const
protectedinherited

Definition at line 539 of file G4VSolid.cc.

542{
543 G4ThreeVectorList outputPolygon;
544
545 if ( pVoxelLimit.IsLimited() )
546 {
547 if (pVoxelLimit.IsXLimited() ) // && pAxis != kXAxis)
548 {
549 G4VoxelLimits simpleLimit1;
550 simpleLimit1.AddLimit(kXAxis,pVoxelLimit.GetMinXExtent(),kInfinity);
551 ClipPolygonToSimpleLimits(pPolygon,outputPolygon,simpleLimit1);
552
553 pPolygon.clear();
554
555 if ( !outputPolygon.size() ) return;
556
557 G4VoxelLimits simpleLimit2;
558 simpleLimit2.AddLimit(kXAxis,-kInfinity,pVoxelLimit.GetMaxXExtent());
559 ClipPolygonToSimpleLimits(outputPolygon,pPolygon,simpleLimit2);
560
561 if ( !pPolygon.size() ) return;
562 else outputPolygon.clear();
563 }
564 if ( pVoxelLimit.IsYLimited() ) // && pAxis != kYAxis)
565 {
566 G4VoxelLimits simpleLimit1;
567 simpleLimit1.AddLimit(kYAxis,pVoxelLimit.GetMinYExtent(),kInfinity);
568 ClipPolygonToSimpleLimits(pPolygon,outputPolygon,simpleLimit1);
569
570 // Must always clear pPolygon - for clip to simpleLimit2 and in case of
571 // early exit
572
573 pPolygon.clear();
574
575 if ( !outputPolygon.size() ) return;
576
577 G4VoxelLimits simpleLimit2;
578 simpleLimit2.AddLimit(kYAxis,-kInfinity,pVoxelLimit.GetMaxYExtent());
579 ClipPolygonToSimpleLimits(outputPolygon,pPolygon,simpleLimit2);
580
581 if ( !pPolygon.size() ) return;
582 else outputPolygon.clear();
583 }
584 if ( pVoxelLimit.IsZLimited() ) // && pAxis != kZAxis)
585 {
586 G4VoxelLimits simpleLimit1;
587 simpleLimit1.AddLimit(kZAxis,pVoxelLimit.GetMinZExtent(),kInfinity);
588 ClipPolygonToSimpleLimits(pPolygon,outputPolygon,simpleLimit1);
589
590 // Must always clear pPolygon - for clip to simpleLimit2 and in case of
591 // early exit
592
593 pPolygon.clear();
594
595 if ( !outputPolygon.size() ) return;
596
597 G4VoxelLimits simpleLimit2;
598 simpleLimit2.AddLimit(kZAxis,-kInfinity,pVoxelLimit.GetMaxZExtent());
599 ClipPolygonToSimpleLimits(outputPolygon,pPolygon,simpleLimit2);
600
601 // Return after final clip - no cleanup
602 }
603 }
604}
void ClipPolygonToSimpleLimits(G4ThreeVectorList &pPolygon, G4ThreeVectorList &outputPolygon, const G4VoxelLimits &pVoxelLimit) const
Definition: G4VSolid.cc:612
G4bool IsYLimited() const
G4double GetMinZExtent() const
void AddLimit(const EAxis pAxis, const G4double pMin, const G4double pMax)
G4bool IsXLimited() const
G4double GetMaxYExtent() const
G4double GetMaxZExtent() const
G4double GetMinYExtent() const
G4double GetMinXExtent() const
G4bool IsZLimited() const
G4bool IsLimited() const
G4double GetMaxXExtent() const
@ kYAxis
Definition: geomdefs.hh:56
@ kXAxis
Definition: geomdefs.hh:55
@ kZAxis
Definition: geomdefs.hh:57
static const G4double kInfinity
Definition: geomdefs.hh:41

References G4VoxelLimits::AddLimit(), G4VSolid::ClipPolygonToSimpleLimits(), G4VoxelLimits::GetMaxXExtent(), G4VoxelLimits::GetMaxYExtent(), G4VoxelLimits::GetMaxZExtent(), G4VoxelLimits::GetMinXExtent(), G4VoxelLimits::GetMinYExtent(), G4VoxelLimits::GetMinZExtent(), G4VoxelLimits::IsLimited(), G4VoxelLimits::IsXLimited(), G4VoxelLimits::IsYLimited(), G4VoxelLimits::IsZLimited(), kInfinity, kXAxis, kYAxis, and kZAxis.

Referenced by G4VSolid::CalculateClippedPolygonExtent().

◆ ClipPolygonToSimpleLimits()

void G4VSolid::ClipPolygonToSimpleLimits ( G4ThreeVectorList pPolygon,
G4ThreeVectorList outputPolygon,
const G4VoxelLimits pVoxelLimit 
) const
privateinherited

Definition at line 612 of file G4VSolid.cc.

615{
616 G4int i;
617 G4int noVertices=pPolygon.size();
618 G4ThreeVector vEnd,vStart;
619
620 for (i = 0 ; i < noVertices ; ++i )
621 {
622 vStart = pPolygon[i];
623 if ( i == noVertices-1 ) vEnd = pPolygon[0];
624 else vEnd = pPolygon[i+1];
625
626 if ( pVoxelLimit.Inside(vStart) )
627 {
628 if (pVoxelLimit.Inside(vEnd))
629 {
630 // vStart and vEnd inside -> output end point
631 //
632 outputPolygon.push_back(vEnd);
633 }
634 else
635 {
636 // vStart inside, vEnd outside -> output crossing point
637 //
638 pVoxelLimit.ClipToLimits(vStart,vEnd);
639 outputPolygon.push_back(vEnd);
640 }
641 }
642 else
643 {
644 if (pVoxelLimit.Inside(vEnd))
645 {
646 // vStart outside, vEnd inside -> output inside section
647 //
648 pVoxelLimit.ClipToLimits(vStart,vEnd);
649 outputPolygon.push_back(vStart);
650 outputPolygon.push_back(vEnd);
651 }
652 else // Both point outside -> no output
653 {
654 // outputPolygon.push_back(vStart);
655 // outputPolygon.push_back(vEnd);
656 }
657 }
658 }
659}
G4bool ClipToLimits(G4ThreeVector &pStart, G4ThreeVector &pEnd) const
G4bool Inside(const G4ThreeVector &pVec) const

References G4VoxelLimits::ClipToLimits(), and G4VoxelLimits::Inside().

Referenced by G4VSolid::ClipPolygon().

◆ Clone()

G4VSolid * G4IntersectionSolid::Clone ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 529 of file G4IntersectionSolid.cc.

530{
531 return new G4IntersectionSolid(*this);
532}
G4IntersectionSolid(const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB)

References G4IntersectionSolid().

◆ ComputeDimensions()

void G4IntersectionSolid::ComputeDimensions ( G4VPVParameterisation p,
const G4int  n,
const G4VPhysicalVolume pRep 
)
virtual

Reimplemented from G4VSolid.

Definition at line 510 of file G4IntersectionSolid.cc.

513{
514}

◆ CreatePolyhedron()

G4Polyhedron * G4IntersectionSolid::CreatePolyhedron ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 549 of file G4IntersectionSolid.cc.

550{
552 // Stack components and components of components recursively
553 // See G4BooleanSolid::StackPolyhedron
555 G4Polyhedron* result = new G4Polyhedron(*top);
556 if (processor.execute(*result)) { return result; }
557 else { return nullptr; }
558}
G4Polyhedron * StackPolyhedron(HepPolyhedronProcessor &, const G4VSolid *) const
#define processor
Definition: xmlparse.cc:617

References processor, and G4BooleanSolid::StackPolyhedron().

◆ DescribeYourselfTo()

void G4IntersectionSolid::DescribeYourselfTo ( G4VGraphicsScene scene) const
virtual

Implements G4VSolid.

Definition at line 539 of file G4IntersectionSolid.cc.

540{
541 scene.AddSolid (*this);
542}
virtual void AddSolid(const G4Box &)=0

References G4VGraphicsScene::AddSolid().

◆ DistanceToIn() [1/2]

G4double G4IntersectionSolid::DistanceToIn ( const G4ThreeVector p) const
virtual

Implements G4VSolid.

Definition at line 382 of file G4IntersectionSolid.cc.

383{
384#ifdef G4BOOLDEBUG
385 if( Inside(p) == kInside )
386 {
387 G4cout << "WARNING - Invalid call in "
388 << "G4IntersectionSolid::DistanceToIn(p)" << G4endl
389 << " Point p is inside !" << G4endl;
390 G4cout << " p = " << p << G4endl;
391 G4cerr << "WARNING - Invalid call in "
392 << "G4IntersectionSolid::DistanceToIn(p)" << G4endl
393 << " Point p is inside !" << G4endl;
394 G4cerr << " p = " << p << G4endl;
395 }
396#endif
397 EInside sideA = fPtrSolidA->Inside(p) ;
398 EInside sideB = fPtrSolidB->Inside(p) ;
399 G4double dist=0.0 ;
400
401 if( sideA != kInside && sideB != kOutside )
402 {
403 dist = fPtrSolidA->DistanceToIn(p) ;
404 }
405 else
406 {
407 if( sideB != kInside && sideA != kOutside )
408 {
409 dist = fPtrSolidB->DistanceToIn(p) ;
410 }
411 else
412 {
414 fPtrSolidB->DistanceToIn(p) ) ;
415 }
416 }
417 return dist ;
418}
G4GLOB_DLL std::ostream G4cerr
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
EInside Inside(const G4ThreeVector &p) const
virtual EInside Inside(const G4ThreeVector &p) const =0
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
EInside
Definition: geomdefs.hh:67
@ kInside
Definition: geomdefs.hh:70
@ kOutside
Definition: geomdefs.hh:68

References G4VSolid::DistanceToIn(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4cerr, G4cout, G4endl, Inside(), G4VSolid::Inside(), kInside, kOutside, and G4INCL::Math::min().

◆ DistanceToIn() [2/2]

G4double G4IntersectionSolid::DistanceToIn ( const G4ThreeVector p,
const G4ThreeVector v 
) const
virtual

Implements G4VSolid.

Definition at line 275 of file G4IntersectionSolid.cc.

277{
278 G4double dist = 0.0;
279 if( Inside(p) == kInside )
280 {
281#ifdef G4BOOLDEBUG
282 G4cout << "WARNING - Invalid call in "
283 << "G4IntersectionSolid::DistanceToIn(p,v)" << G4endl
284 << " Point p is inside !" << G4endl;
285 G4cout << " p = " << p << G4endl;
286 G4cout << " v = " << v << G4endl;
287 G4cerr << "WARNING - Invalid call in "
288 << "G4IntersectionSolid::DistanceToIn(p,v)" << G4endl
289 << " Point p is inside !" << G4endl;
290 G4cerr << " p = " << p << G4endl;
291 G4cerr << " v = " << v << G4endl;
292#endif
293 }
294 else // if( Inside(p) == kSurface )
295 {
296 EInside wA = fPtrSolidA->Inside(p);
297 EInside wB = fPtrSolidB->Inside(p);
298
299 G4ThreeVector pA = p, pB = p;
300 G4double dA = 0., dA1=0., dA2=0.;
301 G4double dB = 0., dB1=0., dB2=0.;
302 G4bool doA = true, doB = true;
303
304 static const size_t max_trials=10000;
305 for (size_t trial=0; trial<max_trials; ++trial)
306 {
307 if(doA)
308 {
309 // find next valid range for A
310
311 dA1 = 0.;
312
313 if( wA != kInside )
314 {
315 dA1 = fPtrSolidA->DistanceToIn(pA, v);
316
317 if( dA1 == kInfinity ) return kInfinity;
318
319 pA += dA1*v;
320 }
321 dA2 = dA1 + fPtrSolidA->DistanceToOut(pA, v);
322 }
323 dA1 += dA;
324 dA2 += dA;
325
326 if(doB)
327 {
328 // find next valid range for B
329
330 dB1 = 0.;
331 if(wB != kInside)
332 {
333 dB1 = fPtrSolidB->DistanceToIn(pB, v);
334
335 if(dB1 == kInfinity) return kInfinity;
336
337 pB += dB1*v;
338 }
339 dB2 = dB1 + fPtrSolidB->DistanceToOut(pB, v);
340 }
341 dB1 += dB;
342 dB2 += dB;
343
344 // check if they overlap
345
346 if( dA1 < dB1 )
347 {
348 if( dB1 < dA2 ) return dB1;
349
350 dA = dA2;
351 pA = p + dA*v; // continue from here
352 wA = kSurface;
353 doA = true;
354 doB = false;
355 }
356 else
357 {
358 if( dA1 < dB2 ) return dA1;
359
360 dB = dB2;
361 pB = p + dB*v; // continue from here
362 wB = kSurface;
363 doB = true;
364 doA = false;
365 }
366 }
367 }
368#ifdef G4BOOLDEBUG
369 G4Exception("G4IntersectionSolid::DistanceToIn(p,v)",
370 "GeomSolids0001", JustWarning,
371 "Reached maximum number of iterations! Returning zero.");
372#endif
373 return dist ;
374}
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=nullptr, G4ThreeVector *n=nullptr) const =0
@ kSurface
Definition: geomdefs.hh:69

References G4VSolid::DistanceToIn(), G4VSolid::DistanceToOut(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4cerr, G4cout, G4endl, G4Exception(), Inside(), G4VSolid::Inside(), JustWarning, kInfinity, kInside, and kSurface.

◆ DistanceToOut() [1/2]

G4double G4IntersectionSolid::DistanceToOut ( const G4ThreeVector p) const
virtual

Implements G4VSolid.

Definition at line 484 of file G4IntersectionSolid.cc.

485{
486#ifdef G4BOOLDEBUG
487 if( Inside(p) == kOutside )
488 {
489 G4cout << "WARNING - Invalid call in "
490 << "G4IntersectionSolid::DistanceToOut(p)" << G4endl
491 << " Point p is outside !" << G4endl;
492 G4cout << " p = " << p << G4endl;
493 G4cerr << "WARNING - Invalid call in "
494 << "G4IntersectionSolid::DistanceToOut(p)" << G4endl
495 << " Point p is outside !" << G4endl;
496 G4cerr << " p = " << p << G4endl;
497 }
498#endif
499
502
503}

References G4VSolid::DistanceToOut(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4cerr, G4cout, G4endl, Inside(), kOutside, and G4INCL::Math::min().

◆ DistanceToOut() [2/2]

G4double G4IntersectionSolid::DistanceToOut ( const G4ThreeVector p,
const G4ThreeVector v,
const G4bool  calcNorm = false,
G4bool validNorm = 0,
G4ThreeVector n = 0 
) const
virtual

Implements G4VSolid.

Definition at line 425 of file G4IntersectionSolid.cc.

430{
431 G4bool validNormA, validNormB;
432 G4ThreeVector nA, nB;
433
434#ifdef G4BOOLDEBUG
435 if( Inside(p) == kOutside )
436 {
437 G4cout << "Position:" << G4endl << G4endl;
438 G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl;
439 G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl;
440 G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl;
441 G4cout << "Direction:" << G4endl << G4endl;
442 G4cout << "v.x() = " << v.x() << G4endl;
443 G4cout << "v.y() = " << v.y() << G4endl;
444 G4cout << "v.z() = " << v.z() << G4endl << G4endl;
445 G4cout << "WARNING - Invalid call in "
446 << "G4IntersectionSolid::DistanceToOut(p,v)" << G4endl
447 << " Point p is outside !" << G4endl;
448 G4cout << " p = " << p << G4endl;
449 G4cout << " v = " << v << G4endl;
450 G4cerr << "WARNING - Invalid call in "
451 << "G4IntersectionSolid::DistanceToOut(p,v)" << G4endl
452 << " Point p is outside !" << G4endl;
453 G4cerr << " p = " << p << G4endl;
454 G4cerr << " v = " << v << G4endl;
455 }
456#endif
457 G4double distA = fPtrSolidA->DistanceToOut(p,v,calcNorm,&validNormA,&nA) ;
458 G4double distB = fPtrSolidB->DistanceToOut(p,v,calcNorm,&validNormB,&nB) ;
459
460 G4double dist = std::min(distA,distB) ;
461
462 if( calcNorm )
463 {
464 if ( distA < distB )
465 {
466 *validNorm = validNormA;
467 *n = nA;
468 }
469 else
470 {
471 *validNorm = validNormB;
472 *n = nB;
473 }
474 }
475
476 return dist ;
477}
static constexpr double mm
Definition: G4SIunits.hh:95

References G4VSolid::DistanceToOut(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4cerr, G4cout, G4endl, Inside(), kOutside, G4INCL::Math::min(), mm, CLHEP::detail::n, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

◆ DumpInfo()

void G4VSolid::DumpInfo ( ) const
inlineinherited

Referenced by G4Cons::ApproxSurfaceNormal(), G4CutTubs::ApproxSurfaceNormal(), G4Sphere::ApproxSurfaceNormal(), G4Torus::ApproxSurfaceNormal(), G4Tubs::ApproxSurfaceNormal(), G4ReflectedSolid::BoundingLimits(), G4DisplacedSolid::BoundingLimits(), BoundingLimits(), G4ScaledSolid::BoundingLimits(), G4SubtractionSolid::BoundingLimits(), G4UnionSolid::BoundingLimits(), G4Box::BoundingLimits(), G4Cons::BoundingLimits(), G4CutTubs::BoundingLimits(), G4Orb::BoundingLimits(), G4Para::BoundingLimits(), G4Sphere::BoundingLimits(), G4Torus::BoundingLimits(), G4Trap::BoundingLimits(), G4Trd::BoundingLimits(), G4Tubs::BoundingLimits(), G4EllipticalCone::BoundingLimits(), G4ExtrudedSolid::BoundingLimits(), G4GenericPolycone::BoundingLimits(), G4GenericTrap::BoundingLimits(), G4Hype::BoundingLimits(), G4Paraboloid::BoundingLimits(), G4Polycone::BoundingLimits(), G4Polyhedra::BoundingLimits(), G4TessellatedSolid::BoundingLimits(), G4TwistedTubs::BoundingLimits(), G4ParameterisationBoxX::ComputeDimensions(), G4ParameterisationBoxY::ComputeDimensions(), G4ParameterisationBoxZ::ComputeDimensions(), G4ParameterisationConsRho::ComputeDimensions(), G4ParameterisationConsPhi::ComputeDimensions(), G4ParameterisationConsZ::ComputeDimensions(), G4ParameterisationParaX::ComputeDimensions(), G4ParameterisationParaY::ComputeDimensions(), G4ParameterisationParaZ::ComputeDimensions(), G4ParameterisationPolyconeRho::ComputeDimensions(), G4ParameterisationPolyconePhi::ComputeDimensions(), G4ParameterisationPolyconeZ::ComputeDimensions(), G4ParameterisationPolyhedraRho::ComputeDimensions(), G4ParameterisationPolyhedraPhi::ComputeDimensions(), G4ParameterisationPolyhedraZ::ComputeDimensions(), G4ParameterisationTrdX::ComputeDimensions(), G4ParameterisationTrdY::ComputeDimensions(), G4ParameterisationTrdZ::ComputeDimensions(), G4ParameterisationTubsRho::ComputeDimensions(), G4ParameterisationTubsPhi::ComputeDimensions(), G4ParameterisationTubsZ::ComputeDimensions(), G4ReflectedSolid::ComputeDimensions(), G4DisplacedSolid::ComputeDimensions(), G4ScaledSolid::ComputeDimensions(), G4ParameterisedNavigation::ComputeStep(), G4ReplicaNavigation::ComputeStep(), G4DisplacedSolid::CreatePolyhedron(), G4ScaledSolid::CreatePolyhedron(), G4SubtractionSolid::DistanceToIn(), G4Box::DistanceToOut(), G4Orb::DistanceToOut(), G4Para::DistanceToOut(), G4Trap::DistanceToOut(), G4Trd::DistanceToOut(), G4Paraboloid::DistanceToOut(), G4VTwistedFaceted::DistanceToOut(), G4Cons::DistanceToOut(), G4CutTubs::DistanceToOut(), G4Sphere::DistanceToOut(), G4Torus::DistanceToOut(), G4Tubs::DistanceToOut(), G4Ellipsoid::DistanceToOut(), G4EllipticalCone::DistanceToOut(), G4EllipticalTube::DistanceToOut(), G4GenericTrap::DistanceToOut(), export_G4VSolid(), G4Polycone::G4Polycone(), G4Polyhedra::G4Polyhedra(), G4BooleanSolid::GetConstituentSolid(), G4NavigationLogger::PostComputeStepLog(), G4Box::SurfaceNormal(), G4Para::SurfaceNormal(), G4Trap::SurfaceNormal(), G4Trd::SurfaceNormal(), G4Ellipsoid::SurfaceNormal(), G4EllipticalCone::SurfaceNormal(), G4EllipticalTube::SurfaceNormal(), G4ExtrudedSolid::SurfaceNormal(), and G4Tet::SurfaceNormal().

◆ EstimateCubicVolume()

G4double G4VSolid::EstimateCubicVolume ( G4int  nStat,
G4double  epsilon 
) const
inherited

Definition at line 203 of file G4VSolid.cc.

204{
205 G4int iInside=0;
206 G4double px,py,pz,minX,maxX,minY,maxY,minZ,maxZ,volume,halfepsilon;
208 EInside in;
209
210 // values needed for CalculateExtent signature
211
212 G4VoxelLimits limit; // Unlimited
213 G4AffineTransform origin;
214
215 // min max extents of pSolid along X,Y,Z
216
217 CalculateExtent(kXAxis,limit,origin,minX,maxX);
218 CalculateExtent(kYAxis,limit,origin,minY,maxY);
219 CalculateExtent(kZAxis,limit,origin,minZ,maxZ);
220
221 // limits
222
223 if(nStat < 100) nStat = 100;
224 if(epsilon > 0.01) epsilon = 0.01;
225 halfepsilon = 0.5*epsilon;
226
227 for(auto i = 0; i < nStat; ++i )
228 {
229 px = minX-halfepsilon+(maxX-minX+epsilon)*G4QuickRand();
230 py = minY-halfepsilon+(maxY-minY+epsilon)*G4QuickRand();
231 pz = minZ-halfepsilon+(maxZ-minZ+epsilon)*G4QuickRand();
232 p = G4ThreeVector(px,py,pz);
233 in = Inside(p);
234 if(in != kOutside) ++iInside;
235 }
236 volume = (maxX-minX+epsilon)*(maxY-minY+epsilon)
237 * (maxZ-minZ+epsilon)*iInside/nStat;
238 return volume;
239}
G4double epsilon(G4double density, G4double temperature)
static const G4int maxZ
G4double G4QuickRand()
Definition: G4QuickRand.hh:34
CLHEP::Hep3Vector G4ThreeVector

References G4VSolid::CalculateExtent(), epsilon(), G4QuickRand(), G4VSolid::Inside(), kOutside, kXAxis, kYAxis, kZAxis, and maxZ.

Referenced by G4VSolid::GetCubicVolume(), G4BooleanSolid::GetCubicVolume(), and G4VCSGfaceted::GetCubicVolume().

◆ EstimateSurfaceArea()

G4double G4VSolid::EstimateSurfaceArea ( G4int  nStat,
G4double  ell 
) const
inherited

Definition at line 265 of file G4VSolid.cc.

266{
267 static const G4double s2 = 1./std::sqrt(2.);
268 static const G4double s3 = 1./std::sqrt(3.);
269 static const G4ThreeVector directions[64] =
270 {
271 G4ThreeVector( 0, 0, 0), G4ThreeVector( -1, 0, 0), // ( , , ) ( -, , )
272 G4ThreeVector( 1, 0, 0), G4ThreeVector( -1, 0, 0), // ( +, , ) (-+, , )
273 G4ThreeVector( 0, -1, 0), G4ThreeVector(-s2,-s2, 0), // ( , -, ) ( -, -, )
274 G4ThreeVector( s2, -s2, 0), G4ThreeVector( 0, -1, 0), // ( +, -, ) (-+, -, )
275
276 G4ThreeVector( 0, 1, 0), G4ThreeVector( -s2, s2, 0), // ( , +, ) ( -, +, )
277 G4ThreeVector( s2, s2, 0), G4ThreeVector( 0, 1, 0), // ( +, +, ) (-+, +, )
278 G4ThreeVector( 0, -1, 0), G4ThreeVector( -1, 0, 0), // ( ,-+, ) ( -,-+, )
279 G4ThreeVector( 1, 0, 0), G4ThreeVector( -1, 0, 0), // ( +,-+, ) (-+,-+, )
280
281 G4ThreeVector( 0, 0, -1), G4ThreeVector(-s2, 0,-s2), // ( , , -) ( -, , -)
282 G4ThreeVector( s2, 0,-s2), G4ThreeVector( 0, 0, -1), // ( +, , -) (-+, , -)
283 G4ThreeVector( 0,-s2,-s2), G4ThreeVector(-s3,-s3,-s3), // ( , -, -) ( -, -, -)
284 G4ThreeVector( s3,-s3,-s3), G4ThreeVector( 0,-s2,-s2), // ( +, -, -) (-+, -, -)
285
286 G4ThreeVector( 0, s2,-s2), G4ThreeVector(-s3, s3,-s3), // ( , +, -) ( -, +, -)
287 G4ThreeVector( s3, s3,-s3), G4ThreeVector( 0, s2,-s2), // ( +, +, -) (-+, +, -)
288 G4ThreeVector( 0, 0, -1), G4ThreeVector(-s2, 0,-s2), // ( ,-+, -) ( -,-+, -)
289 G4ThreeVector( s2, 0,-s2), G4ThreeVector( 0, 0, -1), // ( +,-+, -) (-+,-+, -)
290
291 G4ThreeVector( 0, 0, 1), G4ThreeVector(-s2, 0, s2), // ( , , +) ( -, , +)
292 G4ThreeVector( s2, 0, s2), G4ThreeVector( 0, 0, 1), // ( +, , +) (-+, , +)
293 G4ThreeVector( 0,-s2, s2), G4ThreeVector(-s3,-s3, s3), // ( , -, +) ( -, -, +)
294 G4ThreeVector( s3,-s3, s3), G4ThreeVector( 0,-s2, s2), // ( +, -, +) (-+, -, +)
295
296 G4ThreeVector( 0, s2, s2), G4ThreeVector(-s3, s3, s3), // ( , +, +) ( -, +, +)
297 G4ThreeVector( s3, s3, s3), G4ThreeVector( 0, s2, s2), // ( +, +, +) (-+, +, +)
298 G4ThreeVector( 0, 0, 1), G4ThreeVector(-s2, 0, s2), // ( ,-+, +) ( -,-+, +)
299 G4ThreeVector( s2, 0, s2), G4ThreeVector( 0, 0, 1), // ( +,-+, +) (-+,-+, +)
300
301 G4ThreeVector( 0, 0, -1), G4ThreeVector( -1, 0, 0), // ( , ,-+) ( -, ,-+)
302 G4ThreeVector( 1, 0, 0), G4ThreeVector( -1, 0, 0), // ( +, ,-+) (-+, ,-+)
303 G4ThreeVector( 0, -1, 0), G4ThreeVector(-s2,-s2, 0), // ( , -,-+) ( -, -,-+)
304 G4ThreeVector( s2, -s2, 0), G4ThreeVector( 0, -1, 0), // ( +, -,-+) (-+, -,-+)
305
306 G4ThreeVector( 0, 1, 0), G4ThreeVector( -s2, s2, 0), // ( , +,-+) ( -, +,-+)
307 G4ThreeVector( s2, s2, 0), G4ThreeVector( 0, 1, 0), // ( +, +,-+) (-+, +,-+)
308 G4ThreeVector( 0, -1, 0), G4ThreeVector( -1, 0, 0), // ( ,-+,-+) ( -,-+,-+)
309 G4ThreeVector( 1, 0, 0), G4ThreeVector( -1, 0, 0), // ( +,-+,-+) (-+,-+,-+)
310 };
311
312 G4ThreeVector bmin, bmax;
313 BoundingLimits(bmin, bmax);
314
315 G4double dX = bmax.x() - bmin.x();
316 G4double dY = bmax.y() - bmin.y();
317 G4double dZ = bmax.z() - bmin.z();
318
319 // Define statistics and shell thickness
320 //
321 G4int npoints = (nstat < 1000) ? 1000 : nstat;
322 G4double coeff = 0.5 / std::cbrt(G4double(npoints));
323 G4double eps = (ell > 0) ? ell : coeff * std::min(std::min(dX, dY), dZ);
324 G4double del = 1.8 * eps; // shold be more than sqrt(3.)
325
326 G4double minX = bmin.x() - eps;
327 G4double minY = bmin.y() - eps;
328 G4double minZ = bmin.z() - eps;
329
330 G4double dd = 2. * eps;
331 dX += dd;
332 dY += dd;
333 dZ += dd;
334
335 // Calculate surface area
336 //
337 G4int icount = 0;
338 for(auto i = 0; i < npoints; ++i)
339 {
340 G4double px = minX + dX*G4QuickRand();
341 G4double py = minY + dY*G4QuickRand();
342 G4double pz = minZ + dZ*G4QuickRand();
343 G4ThreeVector p = G4ThreeVector(px, py, pz);
344 EInside in = Inside(p);
345 G4double dist = 0;
346 if (in == kInside)
347 {
348 if (DistanceToOut(p) >= eps) continue;
349 G4int icase = 0;
350 if (Inside(G4ThreeVector(px-del, py, pz)) != kInside) icase += 1;
351 if (Inside(G4ThreeVector(px+del, py, pz)) != kInside) icase += 2;
352 if (Inside(G4ThreeVector(px, py-del, pz)) != kInside) icase += 4;
353 if (Inside(G4ThreeVector(px, py+del, pz)) != kInside) icase += 8;
354 if (Inside(G4ThreeVector(px, py, pz-del)) != kInside) icase += 16;
355 if (Inside(G4ThreeVector(px, py, pz+del)) != kInside) icase += 32;
356 if (icase == 0) continue;
357 G4ThreeVector v = directions[icase];
358 dist = DistanceToOut(p, v);
359 G4ThreeVector n = SurfaceNormal(p + v*dist);
360 dist *= v.dot(n);
361 }
362 else if (in == kOutside)
363 {
364 if (DistanceToIn(p) >= eps) continue;
365 G4int icase = 0;
366 if (Inside(G4ThreeVector(px-del, py, pz)) != kOutside) icase += 1;
367 if (Inside(G4ThreeVector(px+del, py, pz)) != kOutside) icase += 2;
368 if (Inside(G4ThreeVector(px, py-del, pz)) != kOutside) icase += 4;
369 if (Inside(G4ThreeVector(px, py+del, pz)) != kOutside) icase += 8;
370 if (Inside(G4ThreeVector(px, py, pz-del)) != kOutside) icase += 16;
371 if (Inside(G4ThreeVector(px, py, pz+del)) != kOutside) icase += 32;
372 if (icase == 0) continue;
373 G4ThreeVector v = directions[icase];
374 dist = DistanceToIn(p, v);
375 if (dist == kInfinity) continue;
376 G4ThreeVector n = SurfaceNormal(p + v*dist);
377 dist *= -(v.dot(n));
378 }
379 if (dist < eps) ++icount;
380 }
381 return dX*dY*dZ*icount/npoints/dd;
382}
static const G4double eps
double dot(const Hep3Vector &) const
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0

References G4VSolid::BoundingLimits(), G4VSolid::DistanceToIn(), G4VSolid::DistanceToOut(), CLHEP::Hep3Vector::dot(), eps, G4QuickRand(), G4VSolid::Inside(), kInfinity, kInside, kOutside, G4INCL::Math::min(), CLHEP::detail::n, G4VSolid::SurfaceNormal(), CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by G4VSolid::GetSurfaceArea(), G4MultiUnion::GetSurfaceArea(), and G4VCSGfaceted::GetSurfaceArea().

◆ GetAreaAccuracy()

G4double G4BooleanSolid::GetAreaAccuracy ( ) const
inlineinherited

Referenced by export_G4BooleanSolid().

◆ GetAreaStatistics()

G4int G4BooleanSolid::GetAreaStatistics ( ) const
inlineinherited

Referenced by export_G4BooleanSolid().

◆ GetConstituentSolid() [1/2]

G4VSolid * G4BooleanSolid::GetConstituentSolid ( G4int  no)
virtualinherited

Reimplemented from G4VSolid.

Definition at line 181 of file G4BooleanSolid.cc.

182{
183 G4VSolid* subSolid = nullptr;
184 if( no == 0 )
185 subSolid = fPtrSolidA;
186 else if( no == 1 )
187 subSolid = fPtrSolidB;
188 else
189 {
190 DumpInfo();
191 G4Exception("G4BooleanSolid::GetConstituentSolid()",
192 "GeomSolids0002", FatalException, "Invalid solid index.");
193 }
194 return subSolid;
195}
@ FatalException

References G4VSolid::DumpInfo(), FatalException, G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, and G4Exception().

◆ GetConstituentSolid() [2/2]

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

Reimplemented from G4VSolid.

Definition at line 159 of file G4BooleanSolid.cc.

160{
161 const G4VSolid* subSolid = nullptr;
162 if( no == 0 )
163 subSolid = fPtrSolidA;
164 else if( no == 1 )
165 subSolid = fPtrSolidB;
166 else
167 {
168 DumpInfo();
169 G4Exception("G4BooleanSolid::GetConstituentSolid()",
170 "GeomSolids0002", FatalException, "Invalid solid index.");
171 }
172 return subSolid;
173}

References G4VSolid::DumpInfo(), FatalException, G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, and G4Exception().

Referenced by G4tgbGeometryDumper::DumpBooleanVolume().

◆ GetCubicVolume()

G4double G4BooleanSolid::GetCubicVolume ( )
virtualinherited

Reimplemented from G4VSolid.

Reimplemented in G4SubtractionSolid, and G4UnionSolid.

Definition at line 418 of file G4BooleanSolid.cc.

419{
420 if(fCubicVolume < 0.)
421 {
423 }
424 return fCubicVolume;
425}
G4double fCubicVolume
G4double fCubVolEpsilon
G4double EstimateCubicVolume(G4int nStat, G4double epsilon) const
Definition: G4VSolid.cc:203

References G4VSolid::EstimateCubicVolume(), G4BooleanSolid::fCubicVolume, G4BooleanSolid::fCubVolEpsilon, and G4BooleanSolid::fStatistics.

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

◆ GetCubVolEpsilon()

G4double G4BooleanSolid::GetCubVolEpsilon ( ) const
inlineinherited

Referenced by export_G4BooleanSolid().

◆ GetCubVolStatistics()

G4int G4BooleanSolid::GetCubVolStatistics ( ) const
inlineinherited

Referenced by export_G4BooleanSolid().

◆ GetDisplacedSolidPtr() [1/2]

G4DisplacedSolid * G4VSolid::GetDisplacedSolidPtr ( )
virtualinherited

Reimplemented in G4DisplacedSolid.

Definition at line 176 of file G4VSolid.cc.

177{ return nullptr; }

◆ GetDisplacedSolidPtr() [2/2]

const G4DisplacedSolid * G4VSolid::GetDisplacedSolidPtr ( ) const
virtualinherited

Reimplemented in G4DisplacedSolid.

Definition at line 173 of file G4VSolid.cc.

174{ return nullptr; }

◆ GetEntityType()

G4GeometryType G4IntersectionSolid::GetEntityType ( ) const
virtual

Reimplemented from G4BooleanSolid.

Definition at line 520 of file G4IntersectionSolid.cc.

521{
522 return G4String("G4IntersectionSolid");
523}

◆ GetExtent()

G4VisExtent G4VSolid::GetExtent ( ) const
virtualinherited

Reimplemented in G4Box, G4Orb, G4Sphere, G4Ellipsoid, G4EllipticalCone, G4EllipticalTube, G4GenericTrap, G4Hype, G4TessellatedSolid, G4Tet, G4TwistedTubs, G4VCSGfaceted, and G4VTwistedFaceted.

Definition at line 682 of file G4VSolid.cc.

683{
684 G4VisExtent extent;
685 G4VoxelLimits voxelLimits; // Defaults to "infinite" limits.
686 G4AffineTransform affineTransform;
687 G4double vmin, vmax;
688 CalculateExtent(kXAxis,voxelLimits,affineTransform,vmin,vmax);
689 extent.SetXmin (vmin);
690 extent.SetXmax (vmax);
691 CalculateExtent(kYAxis,voxelLimits,affineTransform,vmin,vmax);
692 extent.SetYmin (vmin);
693 extent.SetYmax (vmax);
694 CalculateExtent(kZAxis,voxelLimits,affineTransform,vmin,vmax);
695 extent.SetZmin (vmin);
696 extent.SetZmax (vmax);
697 return extent;
698}
void SetYmin(G4double ymin)
Definition: G4VisExtent.hh:114
void SetYmax(G4double ymax)
Definition: G4VisExtent.hh:116
void SetXmax(G4double xmax)
Definition: G4VisExtent.hh:112
void SetXmin(G4double xmin)
Definition: G4VisExtent.hh:110
void SetZmax(G4double zmax)
Definition: G4VisExtent.hh:120
void SetZmin(G4double zmin)
Definition: G4VisExtent.hh:118

References G4VSolid::CalculateExtent(), kXAxis, kYAxis, kZAxis, G4VisExtent::SetXmax(), G4VisExtent::SetXmin(), G4VisExtent::SetYmax(), G4VisExtent::SetYmin(), G4VisExtent::SetZmax(), and G4VisExtent::SetZmin().

Referenced by G4tgbVolume::BuildSolidForDivision(), G4BoundingExtentScene::ProcessVolume(), G4BoundingSphereScene::ProcessVolume(), and G4VisCommandsTouchable::SetNewValue().

◆ GetListOfPrimitives()

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

Definition at line 229 of file G4BooleanSolid.cc.

232{
234 G4VSolid* solid;
235 G4String type;
236
237 // Repeat two times, first time for fPtrSolidA and then for fPtrSolidB
238 //
239 for (auto i=0; i<2; ++i)
240 {
241 transform = curPlacement;
242 solid = (i == 0) ? fPtrSolidA : fPtrSolidB;
243 type = solid->GetEntityType();
244
245 // While current solid is a trasformed solid just modify transform
246 //
247 while (type == "G4DisplacedSolid" ||
248 type == "G4ReflectedSolid" ||
249 type == "G4ScaledSolid")
250 {
251 if (type == "G4DisplacedSolid")
252 {
254 ((G4DisplacedSolid*)solid)->GetObjectRotation(),
255 ((G4DisplacedSolid*)solid)->GetObjectTranslation());
256 solid = ((G4DisplacedSolid*)solid)->GetConstituentMovedSolid();
257 }
258 else if (type == "G4ReflectedSolid")
259 {
260 transform= transform*((G4ReflectedSolid*)solid)->GetDirectTransform3D();
261 solid = ((G4ReflectedSolid*)solid)->GetConstituentMovedSolid();
262 }
263 else if (type == "G4ScaledSolid")
264 {
265 transform = transform * ((G4ScaledSolid*)solid)->GetScaleTransform();
266 solid = ((G4ScaledSolid*)solid)->GetUnscaledSolid();
267 }
268 type = solid->GetEntityType();
269 }
270
271 // If current solid is a Boolean solid then continue recursion,
272 // otherwise add it to the list of primitives
273 //
274 if (type == "G4UnionSolid" ||
275 type == "G4SubtractionSolid" ||
276 type == "G4IntersectionSolid" ||
277 type == "G4BooleanSolid")
278 {
279 ((G4BooleanSolid *)solid)->GetListOfPrimitives(primitives,transform);
280 }
281 else
282 {
283 primitives.push_back(std::pair<G4VSolid*,G4Transform3D>(solid,transform));
284 }
285 }
286}
HepGeom::Transform3D G4Transform3D
G4GeometryType GetEntityType() const
G4VSolid * GetConstituentMovedSolid() const
virtual G4GeometryType GetEntityType() const =0

References G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4DisplacedSolid::GetConstituentMovedSolid(), G4VSolid::GetEntityType(), and G4coutFormatters::anonymous_namespace{G4coutFormatters.cc}::transform().

Referenced by G4BooleanSolid::GetPointOnSurface().

◆ GetName()

G4String G4VSolid::GetName ( ) const
inlineinherited

Referenced by G4GMocrenFileSceneHandler::AddDetector(), G4HepRepFileSceneHandler::AddHepRepInstance(), G4GMocrenFileSceneHandler::AddPrimitive(), G4HepRepFileSceneHandler::AddSolid(), G4GMocrenFileSceneHandler::AddSolid(), G4VtkSceneHandler::AddSolid(), G4GDMLWriteSolids::AddSolid(), G4NavigationLogger::AlongComputeStepLog(), G4GDMLWriteSolids::BooleanWrite(), G4ReflectedSolid::BoundingLimits(), G4DisplacedSolid::BoundingLimits(), BoundingLimits(), G4ScaledSolid::BoundingLimits(), G4SubtractionSolid::BoundingLimits(), G4UnionSolid::BoundingLimits(), G4Box::BoundingLimits(), G4Cons::BoundingLimits(), G4CutTubs::BoundingLimits(), G4Orb::BoundingLimits(), G4Para::BoundingLimits(), G4Sphere::BoundingLimits(), G4Torus::BoundingLimits(), G4Trap::BoundingLimits(), G4Trd::BoundingLimits(), G4Tubs::BoundingLimits(), G4EllipticalCone::BoundingLimits(), G4ExtrudedSolid::BoundingLimits(), G4GenericPolycone::BoundingLimits(), G4GenericTrap::BoundingLimits(), G4Hype::BoundingLimits(), G4Paraboloid::BoundingLimits(), G4Polycone::BoundingLimits(), G4Polyhedra::BoundingLimits(), G4TessellatedSolid::BoundingLimits(), G4TwistedTubs::BoundingLimits(), G4GDMLWriteSolids::BoxWrite(), G4ExtrudedSolid::CalculateExtent(), G4GenericPolycone::CalculateExtent(), G4Polycone::CalculateExtent(), G4Polyhedra::CalculateExtent(), G4NavigationLogger::CheckDaughterEntryPoint(), G4VDivisionParameterisation::CheckNDivAndWidth(), G4VDivisionParameterisation::CheckOffset(), G4GenericTrap::CheckOrder(), G4Para::CheckParameters(), G4Trap::CheckParameters(), G4Trd::CheckParameters(), G4Ellipsoid::CheckParameters(), G4EllipticalTube::CheckParameters(), G4ParameterisationPolyconeRho::CheckParametersValidity(), G4ParameterisationPolyconeZ::CheckParametersValidity(), G4ParameterisationPolyhedraRho::CheckParametersValidity(), G4ParameterisationPolyhedraPhi::CheckParametersValidity(), G4ParameterisationPolyhedraZ::CheckParametersValidity(), G4PhantomParameterisation::CheckVoxelsFillContainer(), G4GenericTrap::ComputeIsTwisted(), G4VoxelNavigation::ComputeSafety(), G4VoxelSafety::ComputeSafety(), G4NavigationLogger::ComputeSafetyLog(), G4ParameterisedNavigation::ComputeStep(), G4ReplicaNavigation::ComputeStep(), G4GDMLWriteSolids::ConeWrite(), G4Polyhedra::Create(), G4GenericPolycone::Create(), G4Polycone::Create(), G4PhysicalVolumeModel::CreateCurrentAttValues(), G4ReflectedSolid::CreatePolyhedron(), G4ReflectionFactory::CreateReflectedLV(), G4GenericTrap::CreateTessellatedSolid(), G4GDMLWriteSolids::CutTubeWrite(), G4SolidStore::DeRegister(), G4PhysicalVolumeModel::DescribeSolid(), G4SubtractionSolid::DistanceToIn(), G4Paraboloid::DistanceToIn(), G4TessellatedSolid::DistanceToIn(), G4Box::DistanceToOut(), G4Orb::DistanceToOut(), G4Para::DistanceToOut(), G4Trap::DistanceToOut(), G4Trd::DistanceToOut(), G4EllipticalCone::DistanceToOut(), G4TessellatedSolid::DistanceToOut(), G4Ellipsoid::DistanceToOut(), G4EllipticalTube::DistanceToOut(), G4tgbGeometryDumper::DumpMultiUnionVolume(), G4tgbGeometryDumper::DumpScaledVolume(), G4tgbGeometryDumper::DumpSolid(), G4GDMLWriteSolids::ElconeWrite(), G4GDMLWriteSolids::EllipsoidWrite(), G4GDMLWriteSolids::EltubeWrite(), G4PVDivision::ErrorInAxis(), G4ReplicatedSlice::ErrorInAxis(), export_G4VSolid(), G4Box::G4Box(), G4Cons::G4Cons(), G4CutTubs::G4CutTubs(), G4EllipticalCone::G4EllipticalCone(), G4Hype::G4Hype(), G4Para::G4Para(), G4Paraboloid::G4Paraboloid(), G4Polycone::G4Polycone(), G4Polyhedra::G4Polyhedra(), G4Sphere::G4Sphere(), G4Tet::G4Tet(), G4Trap::G4Trap(), G4Tubs::G4Tubs(), G4VParameterisationCons::G4VParameterisationCons(), G4VParameterisationPara::G4VParameterisationPara(), G4VParameterisationPolycone::G4VParameterisationPolycone(), G4VParameterisationPolyhedra::G4VParameterisationPolyhedra(), G4VParameterisationTrd::G4VParameterisationTrd(), G4VTwistedFaceted::G4VTwistedFaceted(), G4GDMLWriteSolids::GenericPolyconeWrite(), G4GDMLWriteSolids::GenTrapWrite(), G4Navigator::GetGlobalExitNormal(), G4Navigator::GetLocalExitNormal(), G4ITNavigator1::GetLocalExitNormal(), G4ITNavigator2::GetLocalExitNormal(), G4BooleanSolid::GetPointOnSurface(), G4PhantomParameterisation::GetReplicaNo(), G4GDMLWriteSolids::HypeWrite(), G4TessellatedSolid::InsideNoVoxels(), G4TessellatedSolid::InsideVoxels(), G4ITNavigator1::LocateGlobalPointAndSetup(), G4ITNavigator2::LocateGlobalPointAndSetup(), G4Navigator::LocateGlobalPointAndSetup(), G4GenericTrap::MakeDownFacet(), G4Trap::MakePlanes(), G4GenericTrap::MakeUpFacet(), G4GDMLWriteSolids::MultiUnionWrite(), G4GDMLWriteSolids::OrbWrite(), G4GDMLWriteSolids::ParaboloidWrite(), G4GDMLWriteParamvol::ParametersWrite(), G4GDMLWriteSolids::ParaWrite(), G4GDMLWriteSolids::PolyconeWrite(), G4GDMLWriteSolids::PolyhedraWrite(), G4NavigationLogger::PostComputeStepLog(), G4NavigationLogger::PreComputeStepLog(), G4NavigationLogger::PrintDaughterLog(), G4PseudoScene::ProcessVolume(), G4SolidStore::Register(), G4tgbVolumeMgr::RegisterMe(), G4NavigationLogger::ReportOutsideMother(), G4ASCIITreeSceneHandler::RequestPrimitives(), G4VSceneHandler::RequestPrimitives(), G4GenericPolycone::Reset(), G4Polyhedra::Reset(), G4VoxelSafety::SafetyForVoxelNode(), G4GDMLWriteSolids::ScaledWrite(), G4Torus::SetAllParameters(), G4Tet::SetBoundingLimits(), G4Polycone::SetOriginalParameters(), G4Polyhedra::SetOriginalParameters(), G4TessellatedSolid::SetSolidClosed(), G4Tet::SetVertices(), G4Box::SetXHalfLength(), G4Box::SetYHalfLength(), G4Box::SetZHalfLength(), G4GDMLWriteSolids::SphereWrite(), G4BooleanSolid::StackPolyhedron(), G4ReflectedSolid::StreamInfo(), G4BooleanSolid::StreamInfo(), G4DisplacedSolid::StreamInfo(), G4MultiUnion::StreamInfo(), G4ScaledSolid::StreamInfo(), G4Box::StreamInfo(), G4Cons::StreamInfo(), G4CSGSolid::StreamInfo(), G4CutTubs::StreamInfo(), G4Orb::StreamInfo(), G4Para::StreamInfo(), G4Sphere::StreamInfo(), G4Torus::StreamInfo(), G4Trap::StreamInfo(), G4Trd::StreamInfo(), G4Tubs::StreamInfo(), G4Ellipsoid::StreamInfo(), G4EllipticalCone::StreamInfo(), G4EllipticalTube::StreamInfo(), G4ExtrudedSolid::StreamInfo(), G4GenericPolycone::StreamInfo(), G4GenericTrap::StreamInfo(), G4Hype::StreamInfo(), G4Paraboloid::StreamInfo(), G4Polycone::StreamInfo(), G4Polyhedra::StreamInfo(), G4TessellatedSolid::StreamInfo(), G4Tet::StreamInfo(), G4TwistedBox::StreamInfo(), G4TwistedTrap::StreamInfo(), G4TwistedTrd::StreamInfo(), G4TwistedTubs::StreamInfo(), G4VCSGfaceted::StreamInfo(), G4VTwistedFaceted::StreamInfo(), G4GDMLRead::StripNames(), SubstractSolids(), G4UnionSolid::SurfaceNormal(), G4Box::SurfaceNormal(), G4Para::SurfaceNormal(), G4Trap::SurfaceNormal(), G4Trd::SurfaceNormal(), G4Ellipsoid::SurfaceNormal(), G4EllipticalCone::SurfaceNormal(), G4EllipticalTube::SurfaceNormal(), G4ExtrudedSolid::SurfaceNormal(), G4Tet::SurfaceNormal(), G4GDMLWriteSolids::TessellatedWrite(), G4GDMLWriteSolids::TetWrite(), G4GDMLWriteSolids::TorusWrite(), G4GDMLWriteSolids::TrapWrite(), G4GDMLWriteStructure::TraverseVolumeTree(), G4GDMLWriteSolids::TrdWrite(), G4GDMLWriteSolids::TubeWrite(), G4GDMLWriteSolids::TwistedboxWrite(), G4GDMLWriteSolids::TwistedtrapWrite(), G4GDMLWriteSolids::TwistedtrdWrite(), G4GDMLWriteSolids::TwistedtubsWrite(), G4PhysicalVolumeModel::VisitGeometryAndGetVisReps(), and G4GDMLWriteSolids::XtruWrite().

◆ GetPointOnSurface()

G4ThreeVector G4BooleanSolid::GetPointOnSurface ( ) const
virtualinherited

Reimplemented from G4VSolid.

Definition at line 293 of file G4BooleanSolid.cc.

294{
295 size_t nprims = fPrimitives.size();
296 std::pair<G4VSolid *, G4Transform3D> prim;
297
298 // Get list of primitives and find the total area of their surfaces
299 //
300 if (nprims == 0)
301 {
303 nprims = fPrimitives.size();
305 for (size_t i=0; i<nprims; ++i)
306 {
307 fPrimitivesSurfaceArea += fPrimitives[i].first->GetSurfaceArea();
308 }
309 }
310
311 // Select random primitive, get random point on its surface and
312 // check that the point belongs to the surface of the solid
313 //
315 for (size_t k=0; k<100000; ++k) // try 100k times
316 {
318 G4double area = 0.;
319 for (size_t i=0; i<nprims; ++i)
320 {
321 prim = fPrimitives[i];
322 area += prim.first->GetSurfaceArea();
323 if (rand < area) break;
324 }
325 p = prim.first->GetPointOnSurface();
326 p = prim.second * G4Point3D(p);
327 if (Inside(p) == kSurface) return p;
328 }
329 std::ostringstream message;
330 message << "Solid - " << GetName() << "\n"
331 << "All 100k attempts to generate a point on the surface have failed!\n"
332 << "The solid created may be an invalid Boolean construct!";
333 G4Exception("G4BooleanSolid::GetPointOnSurface()",
334 "GeomSolids1001", JustWarning, message);
335 return p;
336}
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:34
std::vector< std::pair< G4VSolid *, G4Transform3D > > fPrimitives
void GetListOfPrimitives(std::vector< std::pair< G4VSolid *, G4Transform3D > > &, const G4Transform3D &) const
G4double fPrimitivesSurfaceArea

References G4BooleanSolid::fPrimitives, G4BooleanSolid::fPrimitivesSurfaceArea, G4Exception(), G4QuickRand(), G4BooleanSolid::GetListOfPrimitives(), G4VSolid::GetName(), G4VSolid::Inside(), JustWarning, and kSurface.

Referenced by export_G4BooleanSolid().

◆ GetPolyhedron()

G4Polyhedron * G4BooleanSolid::GetPolyhedron ( ) const
virtualinherited

◆ GetSurfaceArea()

G4double G4BooleanSolid::GetSurfaceArea ( )
inlinevirtualinherited

Reimplemented from G4VSolid.

◆ GetTolerance()

G4double G4VSolid::GetTolerance ( ) const
inlineinherited

◆ Inside()

EInside G4IntersectionSolid::Inside ( const G4ThreeVector p) const
virtual

Implements G4VSolid.

Definition at line 197 of file G4IntersectionSolid.cc.

198{
199 EInside positionA = fPtrSolidA->Inside(p);
200 if(positionA == kOutside) return positionA; // outside A
201
202 EInside positionB = fPtrSolidB->Inside(p);
203 if(positionA == kInside) return positionB;
204
205 if(positionB == kOutside) return positionB; // outside B
206 return kSurface; // surface A & B
207}

References G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4VSolid::Inside(), kInside, kOutside, and kSurface.

Referenced by DistanceToIn(), and DistanceToOut().

◆ operator=()

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

Definition at line 112 of file G4IntersectionSolid.cc.

113{
114 // Check assignment to self
115 //
116 if (this == &rhs) { return *this; }
117
118 // Copy base class data
119 //
121
122 return *this;
123}
G4BooleanSolid & operator=(const G4BooleanSolid &rhs)

References G4BooleanSolid::operator=().

◆ operator==()

G4bool G4VSolid::operator== ( const G4VSolid s) const
inlineinherited

◆ SetAreaAccuracy()

void G4BooleanSolid::SetAreaAccuracy ( G4double  ep)
inlineinherited

Referenced by export_G4BooleanSolid().

◆ SetAreaStatistics()

void G4BooleanSolid::SetAreaStatistics ( G4int  st)
inlineinherited

Referenced by export_G4BooleanSolid().

◆ SetCubVolEpsilon()

void G4BooleanSolid::SetCubVolEpsilon ( G4double  ep)
inlineinherited

Referenced by export_G4BooleanSolid().

◆ SetCubVolStatistics()

void G4BooleanSolid::SetCubVolStatistics ( G4int  st)
inlineinherited

Referenced by export_G4BooleanSolid().

◆ SetName()

void G4VSolid::SetName ( const G4String name)
inherited

Definition at line 127 of file G4VSolid.cc.

128{
131}
void SetMapValid(G4bool val)
Definition: G4SolidStore.hh:76
static G4SolidStore * GetInstance()
G4String fshapeName
Definition: G4VSolid.hh:312
const char * name(G4int ptype)

References G4VSolid::fshapeName, G4SolidStore::GetInstance(), G4InuclParticleNames::name(), and G4SolidStore::SetMapValid().

Referenced by export_G4VSolid(), G4MultiUnion::G4MultiUnion(), and G4GDMLRead::StripNames().

◆ StackPolyhedron()

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

Definition at line 363 of file G4BooleanSolid.cc.

365{
367 const G4String& type = solid->GetEntityType();
368 if (type == "G4UnionSolid")
369 { operation = HepPolyhedronProcessor::UNION; }
370 else if (type == "G4IntersectionSolid")
372 else if (type == "G4SubtractionSolid")
374 else
375 {
376 std::ostringstream message;
377 message << "Solid - " << solid->GetName()
378 << " - Unrecognised composite solid" << G4endl
379 << " Returning NULL !";
380 G4Exception("StackPolyhedron()", "GeomSolids1001", JustWarning, message);
381 return nullptr;
382 }
383
384 G4Polyhedron* top = nullptr;
385 const G4VSolid* solidA = solid->GetConstituentSolid(0);
386 const G4VSolid* solidB = solid->GetConstituentSolid(1);
387
388 if (solidA->GetConstituentSolid(0) != nullptr)
389 {
390 top = StackPolyhedron(processor, solidA);
391 }
392 else
393 {
394 top = solidA->GetPolyhedron();
395 }
397 if (operand != nullptr)
398 {
399 processor.push_back (operation, *operand);
400 }
401 else
402 {
403 std::ostringstream message;
404 message << "Solid - " << solid->GetName()
405 << " - No G4Polyhedron for Boolean component";
406 G4Exception("G4BooleanSolid::StackPolyhedron()",
407 "GeomSolids2001", JustWarning, message);
408 }
409
410 return top;
411}
static int operand(pchar begin, pchar end, double &result, pchar &endp, const dic_type &dictionary)
Definition: Evaluator.cc:163
virtual const G4VSolid * GetConstituentSolid(G4int no) const
Definition: G4VSolid.cc:167
virtual G4Polyhedron * GetPolyhedron() const
Definition: G4VSolid.cc:705

References G4endl, G4Exception(), G4VSolid::GetConstituentSolid(), G4VSolid::GetEntityType(), G4VSolid::GetName(), G4VSolid::GetPolyhedron(), HepPolyhedronProcessor::INTERSECTION, JustWarning, operand(), processor, G4BooleanSolid::StackPolyhedron(), HepPolyhedronProcessor::SUBTRACTION, and HepPolyhedronProcessor::UNION.

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

◆ StreamInfo()

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

Implements G4VSolid.

Definition at line 210 of file G4BooleanSolid.cc.

211{
212 os << "-----------------------------------------------------------\n"
213 << " *** Dump for Boolean solid - " << GetName() << " ***\n"
214 << " ===================================================\n"
215 << " Solid type: " << GetEntityType() << "\n"
216 << " Parameters of constituent solids: \n"
217 << "===========================================================\n";
220 os << "===========================================================\n";
221
222 return os;
223}
virtual G4GeometryType GetEntityType() const
virtual std::ostream & StreamInfo(std::ostream &os) const =0

References G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4BooleanSolid::GetEntityType(), G4VSolid::GetName(), and G4VSolid::StreamInfo().

◆ SurfaceNormal()

G4ThreeVector G4IntersectionSolid::SurfaceNormal ( const G4ThreeVector p) const
virtual

Implements G4VSolid.

Definition at line 213 of file G4IntersectionSolid.cc.

214{
216 EInside insideA, insideB;
217
218 insideA = fPtrSolidA->Inside(p);
219 insideB = fPtrSolidB->Inside(p);
220
221#ifdef G4BOOLDEBUG
222 if( (insideA == kOutside) || (insideB == kOutside) )
223 {
224 G4cout << "WARNING - Invalid call in "
225 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
226 << " Point p is outside !" << G4endl;
227 G4cout << " p = " << p << G4endl;
228 G4cerr << "WARNING - Invalid call in "
229 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
230 << " Point p is outside !" << G4endl;
231 G4cerr << " p = " << p << G4endl;
232 }
233#endif
234
235 // On the surface of both is difficult ... treat it like on A now!
236 //
237 if( insideA == kSurface )
238 {
240 }
241 else if( insideB == kSurface )
242 {
244 }
245 else // We are on neither surface, so we should generate an exception
246 {
248 {
250 }
251 else
252 {
254 }
255#ifdef G4BOOLDEBUG
256 G4cout << "WARNING - Invalid call in "
257 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
258 << " Point p is out of surface !" << G4endl;
259 G4cout << " p = " << p << G4endl;
260 G4cerr << "WARNING - Invalid call in "
261 << "G4IntersectionSolid::SurfaceNormal(p)" << G4endl
262 << " Point p is out of surface !" << G4endl;
263 G4cerr << " p = " << p << G4endl;
264#endif
265 }
266
267 return normal;
268}
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:79

References G4VSolid::DistanceToOut(), G4BooleanSolid::fPtrSolidA, G4BooleanSolid::fPtrSolidB, G4cerr, G4cout, G4endl, G4VSolid::Inside(), kOutside, kSurface, CLHEP::normal(), and G4VSolid::SurfaceNormal().

Field Documentation

◆ createdDisplacedSolid

G4bool G4BooleanSolid::createdDisplacedSolid = false
privateinherited

◆ fAreaAccuracy

G4double G4BooleanSolid::fAreaAccuracy = -1
privateinherited

Definition at line 126 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid::operator=().

◆ fCubicVolume

G4double G4BooleanSolid::fCubicVolume = -1.0
protectedinherited

◆ fCubVolEpsilon

G4double G4BooleanSolid::fCubVolEpsilon = 0.001
privateinherited

Definition at line 125 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid::GetCubicVolume(), and G4BooleanSolid::operator=().

◆ fpPolyhedron

G4Polyhedron* G4BooleanSolid::fpPolyhedron = nullptr
mutableprivateinherited

◆ fPrimitives

std::vector<std::pair<G4VSolid *,G4Transform3D> > G4BooleanSolid::fPrimitives
mutableprivateinherited

◆ fPrimitivesSurfaceArea

G4double G4BooleanSolid::fPrimitivesSurfaceArea = 0.0
mutableprivateinherited

◆ fPtrSolidA

G4VSolid* G4BooleanSolid::fPtrSolidA = nullptr
protectedinherited

◆ fPtrSolidB

G4VSolid* G4BooleanSolid::fPtrSolidB = nullptr
protectedinherited

◆ fRebuildPolyhedron

G4bool G4BooleanSolid::fRebuildPolyhedron = false
mutableprivateinherited

Definition at line 129 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid::GetPolyhedron(), and G4BooleanSolid::operator=().

◆ fshapeName

G4String G4VSolid::fshapeName
privateinherited

Definition at line 312 of file G4VSolid.hh.

Referenced by G4VSolid::operator=(), and G4VSolid::SetName().

◆ fStatistics

G4int G4BooleanSolid::fStatistics = 1000000
privateinherited

Definition at line 124 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid::GetCubicVolume(), and G4BooleanSolid::operator=().

◆ fSurfaceArea

G4double G4BooleanSolid::fSurfaceArea = -1.0
privateinherited

Definition at line 127 of file G4BooleanSolid.hh.

Referenced by G4BooleanSolid::operator=().

◆ kCarTolerance

G4double G4VSolid::kCarTolerance
protectedinherited

Definition at line 299 of file G4VSolid.hh.

Referenced by G4TessellatedSolid::AddFacet(), G4Polycone::CalculateExtent(), G4Polyhedra::CalculateExtent(), G4Tet::CheckDegeneracy(), G4Para::CheckParameters(), G4Trd::CheckParameters(), G4Ellipsoid::CheckParameters(), G4EllipticalTube::CheckParameters(), G4GenericTrap::ComputeIsTwisted(), G4Polyhedra::Create(), G4GenericPolycone::Create(), G4Polycone::Create(), G4CutTubs::CreatePolyhedron(), G4TessellatedSolid::CreateVertexList(), G4VCSGfaceted::DistanceTo(), G4Sphere::DistanceToIn(), G4Ellipsoid::DistanceToIn(), G4Hype::DistanceToIn(), G4Paraboloid::DistanceToIn(), G4VCSGfaceted::DistanceToIn(), G4TessellatedSolid::DistanceToInCore(), G4Cons::DistanceToOut(), G4CutTubs::DistanceToOut(), G4Sphere::DistanceToOut(), G4Torus::DistanceToOut(), G4Tubs::DistanceToOut(), G4GenericTrap::DistanceToOut(), G4Hype::DistanceToOut(), G4Paraboloid::DistanceToOut(), G4VCSGfaceted::DistanceToOut(), G4TessellatedSolid::DistanceToOutCandidates(), G4TessellatedSolid::DistanceToOutCore(), G4TessellatedSolid::DistanceToOutNoVoxels(), G4GenericTrap::DistToPlane(), G4GenericTrap::DistToTriangle(), G4Box::G4Box(), G4Cons::G4Cons(), G4CutTubs::G4CutTubs(), G4EllipticalCone::G4EllipticalCone(), G4ExtrudedSolid::G4ExtrudedSolid(), G4GenericTrap::G4GenericTrap(), G4Hype::G4Hype(), G4Para::G4Para(), G4Sphere::G4Sphere(), G4Tet::G4Tet(), G4Trap::G4Trap(), G4Tubs::G4Tubs(), G4UnionSolid::G4UnionSolid(), G4VSolid::G4VSolid(), G4VTwistedFaceted::G4VTwistedFaceted(), G4GenericPolycone::GetPointOnSurface(), G4Polycone::GetPointOnSurface(), G4UnionSolid::Init(), G4Orb::Initialize(), G4TessellatedSolid::Initialize(), G4SubtractionSolid::Inside(), G4Hype::Inside(), G4Paraboloid::Inside(), G4VCSGfaceted::Inside(), G4VTwistedFaceted::Inside(), G4TessellatedSolid::InsideNoVoxels(), G4GenericTrap::InsidePolygone(), G4TessellatedSolid::InsideVoxels(), G4CutTubs::IsCrossingCutPlanes(), G4GenericTrap::IsSegCrossingZ(), G4Trap::MakePlane(), G4GenericTrap::NormalToPlane(), G4VSolid::operator=(), G4TessellatedSolid::SafetyFromInside(), G4TessellatedSolid::SafetyFromOutside(), G4Torus::SetAllParameters(), G4Polycone::SetOriginalParameters(), G4Polyhedra::SetOriginalParameters(), G4Box::SetXHalfLength(), G4Box::SetYHalfLength(), G4Box::SetZHalfLength(), G4Torus::SurfaceNormal(), G4GenericTrap::SurfaceNormal(), and G4Paraboloid::SurfaceNormal().


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