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

#include <G4Sphere.hh>

Inheritance diagram for G4Sphere:
G4CSGSolid 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=nullptr, G4ThreeVector *n=nullptr) const
 
void DumpInfo () const
 
G4double EstimateCubicVolume (G4int nStat, G4double epsilon) const
 
G4double EstimateSurfaceArea (G4int nStat, G4double ell) const
 
 G4Sphere (__void__ &)
 
 G4Sphere (const G4Sphere &rhs)
 
 G4Sphere (const G4String &pName, G4double pRmin, G4double pRmax, G4double pSPhi, G4double pDPhi, G4double pSTheta, G4double pDTheta)
 
virtual G4VSolidGetConstituentSolid (G4int no)
 
virtual const G4VSolidGetConstituentSolid (G4int no) const
 
G4double GetCosEndPhi () const
 
G4double GetCosEndTheta () const
 
G4double GetCosStartPhi () const
 
G4double GetCosStartTheta () const
 
G4double GetCubicVolume ()
 
G4double GetDeltaPhiAngle () const
 
G4double GetDeltaThetaAngle () const
 
virtual G4DisplacedSolidGetDisplacedSolidPtr ()
 
virtual const G4DisplacedSolidGetDisplacedSolidPtr () const
 
G4GeometryType GetEntityType () const
 
G4VisExtent GetExtent () const
 
G4double GetInnerRadius () const
 
G4String GetName () const
 
G4double GetOuterRadius () const
 
G4ThreeVector GetPointOnSurface () const
 
virtual G4PolyhedronGetPolyhedron () const
 
G4double GetSinEndPhi () const
 
G4double GetSinEndTheta () const
 
G4double GetSinStartPhi () const
 
G4double GetSinStartTheta () const
 
G4double GetStartPhiAngle () const
 
G4double GetStartThetaAngle () const
 
G4double GetSurfaceArea ()
 
G4double GetTolerance () const
 
EInside Inside (const G4ThreeVector &p) const
 
G4Sphereoperator= (const G4Sphere &rhs)
 
G4bool operator== (const G4VSolid &s) const
 
void SetDeltaPhiAngle (G4double newDphi)
 
void SetDeltaThetaAngle (G4double newDTheta)
 
void SetInnerRadius (G4double newRMin)
 
void SetName (const G4String &name)
 
void SetOuterRadius (G4double newRmax)
 
void SetStartPhiAngle (G4double newSphi, G4bool trig=true)
 
void SetStartThetaAngle (G4double newSTheta)
 
std::ostream & StreamInfo (std::ostream &os) const
 
G4ThreeVector SurfaceNormal (const G4ThreeVector &p) const
 
 ~G4Sphere ()
 

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
 
G4double GetRadiusInRing (G4double rmin, G4double rmax) const
 

Protected Attributes

G4double fCubicVolume = 0.0
 
G4PolyhedronfpPolyhedron = nullptr
 
G4bool fRebuildPolyhedron = false
 
G4double fSurfaceArea = 0.0
 
G4double kCarTolerance
 

Private Types

enum  ENorm {
  kNRMin , kNRMax , kNSPhi , kNEPhi ,
  kNSTheta , kNETheta
}
 
enum  ESide {
  kNull , kRMin , kRMax , kSPhi ,
  kEPhi , kSTheta , kETheta
}
 

Private Member Functions

G4ThreeVector ApproxSurfaceNormal (const G4ThreeVector &p) const
 
void CheckDPhiAngle (G4double dPhi)
 
void CheckPhiAngles (G4double sPhi, G4double dPhi)
 
void CheckSPhiAngle (G4double sPhi)
 
void CheckThetaAngles (G4double sTheta, G4double dTheta)
 
void ClipPolygonToSimpleLimits (G4ThreeVectorList &pPolygon, G4ThreeVectorList &outputPolygon, const G4VoxelLimits &pVoxelLimit) const
 
void Initialize ()
 
void InitializePhiTrigonometry ()
 
void InitializeThetaTrigonometry ()
 

Private Attributes

G4double cosCPhi
 
G4double cosEPhi
 
G4double cosETheta
 
G4double cosHDPhi
 
G4double cosHDPhiIT
 
G4double cosHDPhiOT
 
G4double cosSPhi
 
G4double cosSTheta
 
G4double cPhi
 
G4double ePhi
 
G4double eTheta
 
G4double fDPhi
 
G4double fDTheta
 
G4double fEpsilon = 2.e-11
 
G4bool fFullPhiSphere =false
 
G4bool fFullSphere =true
 
G4bool fFullThetaSphere =false
 
G4double fRmax
 
G4double fRmaxTolerance
 
G4double fRmin
 
G4double fRminTolerance
 
G4String fshapeName
 
G4double fSPhi
 
G4double fSTheta
 
G4double halfAngTolerance
 
G4double halfCarTolerance
 
G4double hDPhi
 
G4double kAngTolerance
 
G4double kRadTolerance
 
G4double sinCPhi
 
G4double sinEPhi
 
G4double sinETheta
 
G4double sinSPhi
 
G4double sinSTheta
 
G4double tanETheta
 
G4double tanETheta2
 
G4double tanSTheta
 
G4double tanSTheta2
 

Detailed Description

Definition at line 80 of file G4Sphere.hh.

Member Enumeration Documentation

◆ ENorm

enum G4Sphere::ENorm
private
Enumerator
kNRMin 
kNRMax 
kNSPhi 
kNEPhi 
kNSTheta 
kNETheta 

Definition at line 213 of file G4Sphere.hh.

◆ ESide

enum G4Sphere::ESide
private
Enumerator
kNull 
kRMin 
kRMax 
kSPhi 
kEPhi 
kSTheta 
kETheta 

Definition at line 209 of file G4Sphere.hh.

Constructor & Destructor Documentation

◆ G4Sphere() [1/3]

G4Sphere::G4Sphere ( const G4String pName,
G4double  pRmin,
G4double  pRmax,
G4double  pSPhi,
G4double  pDPhi,
G4double  pSTheta,
G4double  pDTheta 
)

Definition at line 72 of file G4Sphere.cc.

76 : G4CSGSolid(pName), fSPhi(0.0), fFullPhiSphere(true), fFullThetaSphere(true)
77{
80
83
84 // Check radii and set radial tolerances
85
86 if ( (pRmin >= pRmax) || (pRmax < 1.1*kRadTolerance) || (pRmin < 0) )
87 {
88 std::ostringstream message;
89 message << "Invalid radii for Solid: " << GetName() << G4endl
90 << " pRmin = " << pRmin << ", pRmax = " << pRmax;
91 G4Exception("G4Sphere::G4Sphere()", "GeomSolids0002",
92 FatalException, message);
93 }
94 fRmin=pRmin; fRmax=pRmax;
97
98 // Check angles
99
100 CheckPhiAngles(pSPhi, pDPhi);
101 CheckThetaAngles(pSTheta, pDTheta);
102}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
#define G4endl
Definition: G4ios.hh:57
G4CSGSolid(const G4String &pName)
Definition: G4CSGSolid.cc:49
G4double GetRadialTolerance() const
static G4GeometryTolerance * GetInstance()
G4double GetAngularTolerance() const
G4double fRmin
Definition: G4Sphere.hh:220
G4double fEpsilon
Definition: G4Sphere.hh:216
G4bool fFullPhiSphere
Definition: G4Sphere.hh:234
void CheckPhiAngles(G4double sPhi, G4double dPhi)
G4double fSPhi
Definition: G4Sphere.hh:220
G4double halfAngTolerance
Definition: G4Sphere.hh:238
G4double halfCarTolerance
Definition: G4Sphere.hh:238
G4double fRmaxTolerance
Definition: G4Sphere.hh:215
G4double kAngTolerance
Definition: G4Sphere.hh:215
G4double fRmax
Definition: G4Sphere.hh:220
G4bool fFullThetaSphere
Definition: G4Sphere.hh:234
G4double kRadTolerance
Definition: G4Sphere.hh:216
void CheckThetaAngles(G4double sTheta, G4double dTheta)
G4double fRminTolerance
Definition: G4Sphere.hh:215
G4String GetName() const
G4double kCarTolerance
Definition: G4VSolid.hh:299
T max(const T t1, const T t2)
brief Return the largest of the two arguments

References CheckPhiAngles(), CheckThetaAngles(), FatalException, fEpsilon, fRmax, fRmaxTolerance, fRmin, fRminTolerance, G4endl, G4Exception(), G4GeometryTolerance::GetAngularTolerance(), G4GeometryTolerance::GetInstance(), G4VSolid::GetName(), G4GeometryTolerance::GetRadialTolerance(), halfAngTolerance, halfCarTolerance, kAngTolerance, G4VSolid::kCarTolerance, kRadTolerance, and G4INCL::Math::max().

Referenced by Clone().

◆ ~G4Sphere()

G4Sphere::~G4Sphere ( )

Definition at line 126 of file G4Sphere.cc.

127{
128}

◆ G4Sphere() [2/3]

G4Sphere::G4Sphere ( __void__ &  a)

Definition at line 109 of file G4Sphere.cc.

112 fRmin(0.), fRmax(0.), fSPhi(0.), fDPhi(0.), fSTheta(0.),
113 fDTheta(0.), sinCPhi(0.), cosCPhi(0.),
114 cosHDPhi(0.), cosHDPhiOT(0.), cosHDPhiIT(0.),
115 sinSPhi(0.), cosSPhi(0.), sinEPhi(0.), cosEPhi(0.), hDPhi(0.), cPhi(0.),
116 ePhi(0.), sinSTheta(0.), cosSTheta(0.), sinETheta(0.), cosETheta(0.),
117 tanSTheta(0.), tanSTheta2(0.), tanETheta(0.), tanETheta2(0.), eTheta(0.),
119{
120}
G4double cosHDPhiIT
Definition: G4Sphere.hh:224
G4double cPhi
Definition: G4Sphere.hh:225
G4double tanSTheta2
Definition: G4Sphere.hh:230
G4double tanETheta2
Definition: G4Sphere.hh:230
G4double cosSTheta
Definition: G4Sphere.hh:229
G4double cosHDPhiOT
Definition: G4Sphere.hh:224
G4double sinSTheta
Definition: G4Sphere.hh:229
G4double fDPhi
Definition: G4Sphere.hh:220
G4double fSTheta
Definition: G4Sphere.hh:220
G4double hDPhi
Definition: G4Sphere.hh:225
G4double cosHDPhi
Definition: G4Sphere.hh:224
G4double cosEPhi
Definition: G4Sphere.hh:225
G4double tanETheta
Definition: G4Sphere.hh:230
G4double sinEPhi
Definition: G4Sphere.hh:225
G4double cosSPhi
Definition: G4Sphere.hh:225
G4double sinCPhi
Definition: G4Sphere.hh:224
G4double tanSTheta
Definition: G4Sphere.hh:230
G4double ePhi
Definition: G4Sphere.hh:225
G4double cosCPhi
Definition: G4Sphere.hh:224
G4double cosETheta
Definition: G4Sphere.hh:229
G4double fDTheta
Definition: G4Sphere.hh:220
G4double sinETheta
Definition: G4Sphere.hh:229
G4double sinSPhi
Definition: G4Sphere.hh:225
G4double eTheta
Definition: G4Sphere.hh:230

◆ G4Sphere() [3/3]

G4Sphere::G4Sphere ( const G4Sphere rhs)

Member Function Documentation

◆ ApproxSurfaceNormal()

G4ThreeVector G4Sphere::ApproxSurfaceNormal ( const G4ThreeVector p) const
private

Definition at line 531 of file G4Sphere.cc.

532{
533 ENorm side;
534 G4ThreeVector norm;
535 G4double rho,rho2,radius,pPhi,pTheta;
536 G4double distRMin,distRMax,distSPhi,distEPhi,
537 distSTheta,distETheta,distMin;
538
539 rho2=p.x()*p.x()+p.y()*p.y();
540 radius=std::sqrt(rho2+p.z()*p.z());
541 rho=std::sqrt(rho2);
542
543 //
544 // Distance to r shells
545 //
546
547 distRMax=std::fabs(radius-fRmax);
548 if (fRmin)
549 {
550 distRMin=std::fabs(radius-fRmin);
551
552 if (distRMin<distRMax)
553 {
554 distMin=distRMin;
555 side=kNRMin;
556 }
557 else
558 {
559 distMin=distRMax;
560 side=kNRMax;
561 }
562 }
563 else
564 {
565 distMin=distRMax;
566 side=kNRMax;
567 }
568
569 //
570 // Distance to phi planes
571 //
572 // Protected against (0,0,z)
573
574 pPhi = std::atan2(p.y(),p.x());
575 if (pPhi<0) { pPhi += twopi; }
576
577 if (!fFullPhiSphere && rho)
578 {
579 if (fSPhi<0)
580 {
581 distSPhi=std::fabs(pPhi-(fSPhi+twopi))*rho;
582 }
583 else
584 {
585 distSPhi=std::fabs(pPhi-fSPhi)*rho;
586 }
587
588 distEPhi=std::fabs(pPhi-fSPhi-fDPhi)*rho;
589
590 // Find new minimum
591 //
592 if (distSPhi<distEPhi)
593 {
594 if (distSPhi<distMin)
595 {
596 distMin = distSPhi;
597 side = kNSPhi;
598 }
599 }
600 else
601 {
602 if (distEPhi<distMin)
603 {
604 distMin = distEPhi;
605 side = kNEPhi;
606 }
607 }
608 }
609
610 //
611 // Distance to theta planes
612 //
613
614 if (!fFullThetaSphere && radius)
615 {
616 pTheta=std::atan2(rho,p.z());
617 distSTheta=std::fabs(pTheta-fSTheta)*radius;
618 distETheta=std::fabs(pTheta-fSTheta-fDTheta)*radius;
619
620 // Find new minimum
621 //
622 if (distSTheta<distETheta)
623 {
624 if (distSTheta<distMin)
625 {
626 distMin = distSTheta ;
627 side = kNSTheta ;
628 }
629 }
630 else
631 {
632 if (distETheta<distMin)
633 {
634 distMin = distETheta ;
635 side = kNETheta ;
636 }
637 }
638 }
639
640 switch (side)
641 {
642 case kNRMin: // Inner radius
643 norm=G4ThreeVector(-p.x()/radius,-p.y()/radius,-p.z()/radius);
644 break;
645 case kNRMax: // Outer radius
646 norm=G4ThreeVector(p.x()/radius,p.y()/radius,p.z()/radius);
647 break;
648 case kNSPhi:
650 break;
651 case kNEPhi:
653 break;
654 case kNSTheta:
655 norm=G4ThreeVector(-cosSTheta*std::cos(pPhi),
656 -cosSTheta*std::sin(pPhi),
657 sinSTheta );
658 break;
659 case kNETheta:
660 norm=G4ThreeVector( cosETheta*std::cos(pPhi),
661 cosETheta*std::sin(pPhi),
662 -sinETheta );
663 break;
664 default: // Should never reach this case ...
665 DumpInfo();
666 G4Exception("G4Sphere::ApproxSurfaceNormal()",
667 "GeomSolids1002", JustWarning,
668 "Undefined side for valid surface normal to solid.");
669 break;
670 }
671
672 return norm;
673}
ENorm
Definition: G4Cons.cc:65
@ JustWarning
static constexpr double twopi
Definition: G4SIunits.hh:56
CLHEP::Hep3Vector G4ThreeVector
double G4double
Definition: G4Types.hh:83
double z() const
double x() const
double y() const
void DumpInfo() const

References cosEPhi, cosETheta, cosSPhi, cosSTheta, G4VSolid::DumpInfo(), fDPhi, fDTheta, fFullPhiSphere, fFullThetaSphere, fRmax, fRmin, fSPhi, fSTheta, G4Exception(), JustWarning, kNEPhi, kNETheta, kNRMax, kNRMin, kNSPhi, kNSTheta, sinEPhi, sinETheta, sinSPhi, sinSTheta, twopi, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by SurfaceNormal().

◆ BoundingLimits()

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

Reimplemented from G4VSolid.

Definition at line 210 of file G4Sphere.cc.

211{
212 G4double rmin = GetInnerRadius();
213 G4double rmax = GetOuterRadius();
214
215 // Find bounding box
216 //
218 {
219 pMin.set(-rmax,-rmax,-rmax);
220 pMax.set( rmax, rmax, rmax);
221 }
222 else
223 {
224 G4double sinStart = GetSinStartTheta();
225 G4double cosStart = GetCosStartTheta();
226 G4double sinEnd = GetSinEndTheta();
227 G4double cosEnd = GetCosEndTheta();
228
229 G4double stheta = GetStartThetaAngle();
230 G4double etheta = stheta + GetDeltaThetaAngle();
231 G4double rhomin = rmin*std::min(sinStart,sinEnd);
232 G4double rhomax = rmax;
233 if (stheta > halfpi) rhomax = rmax*sinStart;
234 if (etheta < halfpi) rhomax = rmax*sinEnd;
235
236 G4TwoVector xymin,xymax;
237 G4GeomTools::DiskExtent(rhomin,rhomax,
240 xymin,xymax);
241
242 G4double zmin = std::min(rmin*cosEnd,rmax*cosEnd);
243 G4double zmax = std::max(rmin*cosStart,rmax*cosStart);
244 pMin.set(xymin.x(),xymin.y(),zmin);
245 pMax.set(xymax.x(),xymax.y(),zmax);
246 }
247
248 // Check correctness of the bounding box
249 //
250 if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
251 {
252 std::ostringstream message;
253 message << "Bad bounding box (min >= max) for solid: "
254 << GetName() << " !"
255 << "\npMin = " << pMin
256 << "\npMax = " << pMax;
257 G4Exception("G4Sphere::BoundingLimits()", "GeomMgt0001",
258 JustWarning, message);
259 DumpInfo();
260 }
261}
static const G4double pMax
static const G4double pMin
static constexpr double pi
Definition: G4SIunits.hh:55
static constexpr double halfpi
Definition: G4SIunits.hh:57
double x() const
double y() const
static G4bool DiskExtent(G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
Definition: G4GeomTools.cc:390
G4double GetSinStartTheta() const
G4double GetCosStartPhi() const
G4double GetDeltaPhiAngle() const
G4double GetCosEndTheta() const
G4double GetInnerRadius() const
G4double GetOuterRadius() const
G4double GetCosEndPhi() const
G4double GetSinEndTheta() const
G4double GetDeltaThetaAngle() const
G4double GetSinEndPhi() const
G4double GetSinStartPhi() const
G4double GetStartThetaAngle() const
G4double GetCosStartTheta() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments

References G4GeomTools::DiskExtent(), G4VSolid::DumpInfo(), G4Exception(), GetCosEndPhi(), GetCosEndTheta(), GetCosStartPhi(), GetCosStartTheta(), GetDeltaPhiAngle(), GetDeltaThetaAngle(), GetInnerRadius(), G4VSolid::GetName(), GetOuterRadius(), GetSinEndPhi(), GetSinEndTheta(), GetSinStartPhi(), GetSinStartTheta(), GetStartThetaAngle(), halfpi, JustWarning, G4INCL::Math::max(), G4INCL::Math::min(), pi, pMax, pMin, twopi, CLHEP::Hep2Vector::x(), and CLHEP::Hep2Vector::y().

Referenced by CalculateExtent().

◆ 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}
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 G4Sphere::CalculateExtent ( const EAxis  pAxis,
const G4VoxelLimits pVoxelLimit,
const G4AffineTransform pTransform,
G4double pmin,
G4double pmax 
) const
virtual

Implements G4VSolid.

Definition at line 267 of file G4Sphere.cc.

271{
272 G4ThreeVector bmin, bmax;
273
274 // Get bounding box
275 BoundingLimits(bmin,bmax);
276
277 // Find extent
278 G4BoundingEnvelope bbox(bmin,bmax);
279 return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
280}
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Sphere.cc:210

References BoundingLimits(), G4BoundingEnvelope::CalculateExtent(), pMax, and pMin.

◆ CheckDPhiAngle()

void G4Sphere::CheckDPhiAngle ( G4double  dPhi)
inlineprivate

◆ CheckPhiAngles()

void G4Sphere::CheckPhiAngles ( G4double  sPhi,
G4double  dPhi 
)
inlineprivate

Referenced by G4Sphere().

◆ CheckSPhiAngle()

void G4Sphere::CheckSPhiAngle ( G4double  sPhi)
inlineprivate

◆ CheckThetaAngles()

void G4Sphere::CheckThetaAngles ( G4double  sTheta,
G4double  dTheta 
)
inlineprivate

Referenced by G4Sphere().

◆ 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 * G4Sphere::Clone ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 2742 of file G4Sphere.cc.

2743{
2744 return new G4Sphere(*this);
2745}
G4Sphere(const G4String &pName, G4double pRmin, G4double pRmax, G4double pSPhi, G4double pDPhi, G4double pSTheta, G4double pDTheta)
Definition: G4Sphere.cc:72

References G4Sphere().

◆ ComputeDimensions()

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

Reimplemented from G4VSolid.

Definition at line 199 of file G4Sphere.cc.

202{
203 p->ComputeDimensions(*this,n,pRep);
204}
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const

References G4VPVParameterisation::ComputeDimensions(), and CLHEP::detail::n.

◆ CreatePolyhedron()

G4Polyhedron * G4Sphere::CreatePolyhedron ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 2870 of file G4Sphere.cc.

References fDPhi, fDTheta, fRmax, fRmin, fSPhi, and fSTheta.

◆ DescribeYourselfTo()

void G4Sphere::DescribeYourselfTo ( G4VGraphicsScene scene) const
virtual

Implements G4VSolid.

Definition at line 2865 of file G4Sphere.cc.

2866{
2867 scene.AddSolid (*this);
2868}
virtual void AddSolid(const G4Box &)=0

References G4VGraphicsScene::AddSolid().

◆ DistanceToIn() [1/2]

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

Implements G4VSolid.

Definition at line 1653 of file G4Sphere.cc.

1654{
1655 G4double safe=0.0,safeRMin,safeRMax,safePhi,safeTheta;
1656 G4double rho2,rds,rho;
1657 G4double cosPsi;
1658 G4double pTheta,dTheta1,dTheta2;
1659 rho2=p.x()*p.x()+p.y()*p.y();
1660 rds=std::sqrt(rho2+p.z()*p.z());
1661 rho=std::sqrt(rho2);
1662
1663 //
1664 // Distance to r shells
1665 //
1666 if (fRmin)
1667 {
1668 safeRMin=fRmin-rds;
1669 safeRMax=rds-fRmax;
1670 if (safeRMin>safeRMax)
1671 {
1672 safe=safeRMin;
1673 }
1674 else
1675 {
1676 safe=safeRMax;
1677 }
1678 }
1679 else
1680 {
1681 safe=rds-fRmax;
1682 }
1683
1684 //
1685 // Distance to phi extent
1686 //
1687 if (!fFullPhiSphere && rho)
1688 {
1689 // Psi=angle from central phi to point
1690 //
1691 cosPsi=(p.x()*cosCPhi+p.y()*sinCPhi)/rho;
1692 if (cosPsi<cosHDPhi)
1693 {
1694 // Point lies outside phi range
1695 //
1696 if ((p.y()*cosCPhi-p.x()*sinCPhi)<=0)
1697 {
1698 safePhi=std::fabs(p.x()*sinSPhi-p.y()*cosSPhi);
1699 }
1700 else
1701 {
1702 safePhi=std::fabs(p.x()*sinEPhi-p.y()*cosEPhi);
1703 }
1704 if (safePhi>safe) { safe=safePhi; }
1705 }
1706 }
1707 //
1708 // Distance to Theta extent
1709 //
1710 if ((rds!=0.0) && (!fFullThetaSphere))
1711 {
1712 pTheta=std::acos(p.z()/rds);
1713 if (pTheta<0) { pTheta+=pi; }
1714 dTheta1=fSTheta-pTheta;
1715 dTheta2=pTheta-eTheta;
1716 if (dTheta1>dTheta2)
1717 {
1718 if (dTheta1>=0) // WHY ???????????
1719 {
1720 safeTheta=rds*std::sin(dTheta1);
1721 if (safe<=safeTheta)
1722 {
1723 safe=safeTheta;
1724 }
1725 }
1726 }
1727 else
1728 {
1729 if (dTheta2>=0)
1730 {
1731 safeTheta=rds*std::sin(dTheta2);
1732 if (safe<=safeTheta)
1733 {
1734 safe=safeTheta;
1735 }
1736 }
1737 }
1738 }
1739
1740 if (safe<0) { safe=0; }
1741 return safe;
1742}

References cosCPhi, cosEPhi, cosHDPhi, cosSPhi, eTheta, fFullPhiSphere, fFullThetaSphere, fRmax, fRmin, fSTheta, pi, sinCPhi, sinEPhi, sinSPhi, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

◆ DistanceToIn() [2/2]

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

Implements G4VSolid.

Definition at line 704 of file G4Sphere.cc.

706{
707 G4double snxt = kInfinity ; // snxt = default return value
708 G4double rho2, rad2, pDotV2d, pDotV3d, pTheta ;
709 G4double tolSTheta=0., tolETheta=0. ;
710 const G4double dRmax = 100.*fRmax;
711
712 const G4double halfRmaxTolerance = fRmaxTolerance*0.5;
713 const G4double halfRminTolerance = fRminTolerance*0.5;
714 const G4double tolORMin2 = (fRmin>halfRminTolerance)
715 ? (fRmin-halfRminTolerance)*(fRmin-halfRminTolerance) : 0;
716 const G4double tolIRMin2 =
717 (fRmin+halfRminTolerance)*(fRmin+halfRminTolerance);
718 const G4double tolORMax2 =
719 (fRmax+halfRmaxTolerance)*(fRmax+halfRmaxTolerance);
720 const G4double tolIRMax2 =
721 (fRmax-halfRmaxTolerance)*(fRmax-halfRmaxTolerance);
722
723 // Intersection point
724 //
725 G4double xi, yi, zi, rhoi, rhoi2, radi2, iTheta ;
726
727 // Phi intersection
728 //
729 G4double Comp ;
730
731 // Phi precalcs
732 //
733 G4double Dist, cosPsi ;
734
735 // Theta precalcs
736 //
737 G4double dist2STheta, dist2ETheta ;
738 G4double t1, t2, b, c, d2, d, sd = kInfinity ;
739
740 // General Precalcs
741 //
742 rho2 = p.x()*p.x() + p.y()*p.y() ;
743 rad2 = rho2 + p.z()*p.z() ;
744 pTheta = std::atan2(std::sqrt(rho2),p.z()) ;
745
746 pDotV2d = p.x()*v.x() + p.y()*v.y() ;
747 pDotV3d = pDotV2d + p.z()*v.z() ;
748
749 // Theta precalcs
750 //
751 if (!fFullThetaSphere)
752 {
753 tolSTheta = fSTheta - halfAngTolerance ;
754 tolETheta = eTheta + halfAngTolerance ;
755
756 // Special case rad2 = 0 comparing with direction
757 //
758 if ((rad2!=0.0) || (fRmin!=0.0))
759 {
760 // Keep going for computation of distance...
761 }
762 else // Positioned on the sphere's origin
763 {
764 G4double vTheta = std::atan2(std::sqrt(v.x()*v.x()+v.y()*v.y()),v.z()) ;
765 if ( (vTheta < tolSTheta) || (vTheta > tolETheta) )
766 {
767 return snxt ; // kInfinity
768 }
769 return snxt = 0.0 ;
770 }
771 }
772
773 // Outer spherical shell intersection
774 // - Only if outside tolerant fRmax
775 // - Check for if inside and outer G4Sphere heading through solid (-> 0)
776 // - No intersect -> no intersection with G4Sphere
777 //
778 // Shell eqn: x^2+y^2+z^2=RSPH^2
779 //
780 // => (px+svx)^2+(py+svy)^2+(pz+svz)^2=R^2
781 //
782 // => (px^2+py^2+pz^2) +2sd(pxvx+pyvy+pzvz)+sd^2(vx^2+vy^2+vz^2)=R^2
783 // => rad2 +2sd(pDotV3d) +sd^2 =R^2
784 //
785 // => sd=-pDotV3d+-std::sqrt(pDotV3d^2-(rad2-R^2))
786
787 c = rad2 - fRmax*fRmax ;
788
789 if (c > fRmaxTolerance*fRmax)
790 {
791 // If outside tolerant boundary of outer G4Sphere
792 // [should be std::sqrt(rad2)-fRmax > halfRmaxTolerance]
793
794 d2 = pDotV3d*pDotV3d - c ;
795
796 if ( d2 >= 0 )
797 {
798 sd = -pDotV3d - std::sqrt(d2) ;
799
800 if (sd >= 0 )
801 {
802 if ( sd>dRmax ) // Avoid rounding errors due to precision issues seen on
803 { // 64 bits systems. Split long distances and recompute
804 G4double fTerm = sd-std::fmod(sd,dRmax);
805 sd = fTerm + DistanceToIn(p+fTerm*v,v);
806 }
807 xi = p.x() + sd*v.x() ;
808 yi = p.y() + sd*v.y() ;
809 rhoi = std::sqrt(xi*xi + yi*yi) ;
810
811 if (!fFullPhiSphere && rhoi) // Check phi intersection
812 {
813 cosPsi = (xi*cosCPhi + yi*sinCPhi)/rhoi ;
814
815 if (cosPsi >= cosHDPhiOT)
816 {
817 if (!fFullThetaSphere) // Check theta intersection
818 {
819 zi = p.z() + sd*v.z() ;
820
821 // rhoi & zi can never both be 0
822 // (=>intersect at origin =>fRmax=0)
823 //
824 iTheta = std::atan2(rhoi,zi) ;
825 if ( (iTheta >= tolSTheta) && (iTheta <= tolETheta) )
826 {
827 return snxt = sd ;
828 }
829 }
830 else
831 {
832 return snxt=sd;
833 }
834 }
835 }
836 else
837 {
838 if (!fFullThetaSphere) // Check theta intersection
839 {
840 zi = p.z() + sd*v.z() ;
841
842 // rhoi & zi can never both be 0
843 // (=>intersect at origin => fRmax=0 !)
844 //
845 iTheta = std::atan2(rhoi,zi) ;
846 if ( (iTheta >= tolSTheta) && (iTheta <= tolETheta) )
847 {
848 return snxt=sd;
849 }
850 }
851 else
852 {
853 return snxt = sd;
854 }
855 }
856 }
857 }
858 else // No intersection with G4Sphere
859 {
860 return snxt=kInfinity;
861 }
862 }
863 else
864 {
865 // Inside outer radius
866 // check not inside, and heading through G4Sphere (-> 0 to in)
867
868 d2 = pDotV3d*pDotV3d - c ;
869
870 if ( (rad2 > tolIRMax2)
871 && ( (d2 >= fRmaxTolerance*fRmax) && (pDotV3d < 0) ) )
872 {
873 if (!fFullPhiSphere)
874 {
875 // Use inner phi tolerant boundary -> if on tolerant
876 // phi boundaries, phi intersect code handles leaving/entering checks
877
878 cosPsi = (p.x()*cosCPhi + p.y()*sinCPhi)/std::sqrt(rho2) ;
879
880 if (cosPsi>=cosHDPhiIT)
881 {
882 // inside radii, delta r -ve, inside phi
883
884 if ( !fFullThetaSphere )
885 {
886 if ( (pTheta >= tolSTheta + kAngTolerance)
887 && (pTheta <= tolETheta - kAngTolerance) )
888 {
889 return snxt=0;
890 }
891 }
892 else // strictly inside Theta in both cases
893 {
894 return snxt=0;
895 }
896 }
897 }
898 else
899 {
900 if ( !fFullThetaSphere )
901 {
902 if ( (pTheta >= tolSTheta + kAngTolerance)
903 && (pTheta <= tolETheta - kAngTolerance) )
904 {
905 return snxt=0;
906 }
907 }
908 else // strictly inside Theta in both cases
909 {
910 return snxt=0;
911 }
912 }
913 }
914 }
915
916 // Inner spherical shell intersection
917 // - Always farthest root, because would have passed through outer
918 // surface first.
919 // - Tolerant check if travelling through solid
920
921 if (fRmin)
922 {
923 c = rad2 - fRmin*fRmin ;
924 d2 = pDotV3d*pDotV3d - c ;
925
926 // Within tolerance inner radius of inner G4Sphere
927 // Check for immediate entry/already inside and travelling outwards
928
929 if ( (c > -halfRminTolerance) && (rad2 < tolIRMin2)
930 && ( (d2 < fRmin*kCarTolerance) || (pDotV3d >= 0) ) )
931 {
932 if ( !fFullPhiSphere )
933 {
934 // Use inner phi tolerant boundary -> if on tolerant
935 // phi boundaries, phi intersect code handles leaving/entering checks
936
937 cosPsi = (p.x()*cosCPhi+p.y()*sinCPhi)/std::sqrt(rho2) ;
938 if (cosPsi >= cosHDPhiIT)
939 {
940 // inside radii, delta r -ve, inside phi
941 //
942 if ( !fFullThetaSphere )
943 {
944 if ( (pTheta >= tolSTheta + kAngTolerance)
945 && (pTheta <= tolETheta - kAngTolerance) )
946 {
947 return snxt=0;
948 }
949 }
950 else
951 {
952 return snxt = 0 ;
953 }
954 }
955 }
956 else
957 {
958 if ( !fFullThetaSphere )
959 {
960 if ( (pTheta >= tolSTheta + kAngTolerance)
961 && (pTheta <= tolETheta - kAngTolerance) )
962 {
963 return snxt = 0 ;
964 }
965 }
966 else
967 {
968 return snxt=0;
969 }
970 }
971 }
972 else // Not special tolerant case
973 {
974 if (d2 >= 0)
975 {
976 sd = -pDotV3d + std::sqrt(d2) ;
977 if ( sd >= halfRminTolerance ) // It was >= 0 ??
978 {
979 xi = p.x() + sd*v.x() ;
980 yi = p.y() + sd*v.y() ;
981 rhoi = std::sqrt(xi*xi+yi*yi) ;
982
983 if ( !fFullPhiSphere && rhoi ) // Check phi intersection
984 {
985 cosPsi = (xi*cosCPhi + yi*sinCPhi)/rhoi ;
986
987 if (cosPsi >= cosHDPhiOT)
988 {
989 if ( !fFullThetaSphere ) // Check theta intersection
990 {
991 zi = p.z() + sd*v.z() ;
992
993 // rhoi & zi can never both be 0
994 // (=>intersect at origin =>fRmax=0)
995 //
996 iTheta = std::atan2(rhoi,zi) ;
997 if ( (iTheta >= tolSTheta) && (iTheta<=tolETheta) )
998 {
999 snxt = sd;
1000 }
1001 }
1002 else
1003 {
1004 snxt=sd;
1005 }
1006 }
1007 }
1008 else
1009 {
1010 if ( !fFullThetaSphere ) // Check theta intersection
1011 {
1012 zi = p.z() + sd*v.z() ;
1013
1014 // rhoi & zi can never both be 0
1015 // (=>intersect at origin => fRmax=0 !)
1016 //
1017 iTheta = std::atan2(rhoi,zi) ;
1018 if ( (iTheta >= tolSTheta) && (iTheta <= tolETheta) )
1019 {
1020 snxt = sd;
1021 }
1022 }
1023 else
1024 {
1025 snxt = sd;
1026 }
1027 }
1028 }
1029 }
1030 }
1031 }
1032
1033 // Phi segment intersection
1034 //
1035 // o Tolerant of points inside phi planes by up to kCarTolerance*0.5
1036 //
1037 // o NOTE: Large duplication of code between sphi & ephi checks
1038 // -> only diffs: sphi -> ephi, Comp -> -Comp and half-plane
1039 // intersection check <=0 -> >=0
1040 // -> Should use some form of loop Construct
1041 //
1042 if ( !fFullPhiSphere )
1043 {
1044 // First phi surface ('S'tarting phi)
1045 // Comp = Component in outwards normal dirn
1046 //
1047 Comp = v.x()*sinSPhi - v.y()*cosSPhi ;
1048
1049 if ( Comp < 0 )
1050 {
1051 Dist = p.y()*cosSPhi - p.x()*sinSPhi ;
1052
1053 if (Dist < halfCarTolerance)
1054 {
1055 sd = Dist/Comp ;
1056
1057 if (sd < snxt)
1058 {
1059 if ( sd > 0 )
1060 {
1061 xi = p.x() + sd*v.x() ;
1062 yi = p.y() + sd*v.y() ;
1063 zi = p.z() + sd*v.z() ;
1064 rhoi2 = xi*xi + yi*yi ;
1065 radi2 = rhoi2 + zi*zi ;
1066 }
1067 else
1068 {
1069 sd = 0 ;
1070 xi = p.x() ;
1071 yi = p.y() ;
1072 zi = p.z() ;
1073 rhoi2 = rho2 ;
1074 radi2 = rad2 ;
1075 }
1076 if ( (radi2 <= tolORMax2)
1077 && (radi2 >= tolORMin2)
1078 && ((yi*cosCPhi-xi*sinCPhi) <= 0) )
1079 {
1080 // Check theta intersection
1081 // rhoi & zi can never both be 0
1082 // (=>intersect at origin =>fRmax=0)
1083 //
1084 if ( !fFullThetaSphere )
1085 {
1086 iTheta = std::atan2(std::sqrt(rhoi2),zi) ;
1087 if ( (iTheta >= tolSTheta) && (iTheta <= tolETheta) )
1088 {
1089 // r and theta intersections good
1090 // - check intersecting with correct half-plane
1091
1092 if ((yi*cosCPhi-xi*sinCPhi) <= 0)
1093 {
1094 snxt = sd;
1095 }
1096 }
1097 }
1098 else
1099 {
1100 snxt = sd;
1101 }
1102 }
1103 }
1104 }
1105 }
1106
1107 // Second phi surface ('E'nding phi)
1108 // Component in outwards normal dirn
1109
1110 Comp = -( v.x()*sinEPhi-v.y()*cosEPhi ) ;
1111
1112 if (Comp < 0)
1113 {
1114 Dist = -(p.y()*cosEPhi-p.x()*sinEPhi) ;
1115 if ( Dist < halfCarTolerance )
1116 {
1117 sd = Dist/Comp ;
1118
1119 if ( sd < snxt )
1120 {
1121 if (sd > 0)
1122 {
1123 xi = p.x() + sd*v.x() ;
1124 yi = p.y() + sd*v.y() ;
1125 zi = p.z() + sd*v.z() ;
1126 rhoi2 = xi*xi + yi*yi ;
1127 radi2 = rhoi2 + zi*zi ;
1128 }
1129 else
1130 {
1131 sd = 0 ;
1132 xi = p.x() ;
1133 yi = p.y() ;
1134 zi = p.z() ;
1135 rhoi2 = rho2 ;
1136 radi2 = rad2 ;
1137 }
1138 if ( (radi2 <= tolORMax2)
1139 && (radi2 >= tolORMin2)
1140 && ((yi*cosCPhi-xi*sinCPhi) >= 0) )
1141 {
1142 // Check theta intersection
1143 // rhoi & zi can never both be 0
1144 // (=>intersect at origin =>fRmax=0)
1145 //
1146 if ( !fFullThetaSphere )
1147 {
1148 iTheta = std::atan2(std::sqrt(rhoi2),zi) ;
1149 if ( (iTheta >= tolSTheta) && (iTheta <= tolETheta) )
1150 {
1151 // r and theta intersections good
1152 // - check intersecting with correct half-plane
1153
1154 if ((yi*cosCPhi-xi*sinCPhi) >= 0)
1155 {
1156 snxt = sd;
1157 }
1158 }
1159 }
1160 else
1161 {
1162 snxt = sd;
1163 }
1164 }
1165 }
1166 }
1167 }
1168 }
1169
1170 // Theta segment intersection
1171
1172 if ( !fFullThetaSphere )
1173 {
1174
1175 // Intersection with theta surfaces
1176 // Known failure cases:
1177 // o Inside tolerance of stheta surface, skim
1178 // ~parallel to cone and Hit & enter etheta surface [& visa versa]
1179 //
1180 // To solve: Check 2nd root of etheta surface in addition to stheta
1181 //
1182 // o start/end theta is exactly pi/2
1183 // Intersections with cones
1184 //
1185 // Cone equation: x^2+y^2=z^2tan^2(t)
1186 //
1187 // => (px+svx)^2+(py+svy)^2=(pz+svz)^2tan^2(t)
1188 //
1189 // => (px^2+py^2-pz^2tan^2(t))+2sd(pxvx+pyvy-pzvztan^2(t))
1190 // + sd^2(vx^2+vy^2-vz^2tan^2(t)) = 0
1191 //
1192 // => sd^2(1-vz^2(1+tan^2(t))+2sd(pdotv2d-pzvztan^2(t))
1193 // + (rho2-pz^2tan^2(t)) = 0
1194
1195 if (fSTheta)
1196 {
1197 dist2STheta = rho2 - p.z()*p.z()*tanSTheta2 ;
1198 }
1199 else
1200 {
1201 dist2STheta = kInfinity ;
1202 }
1203 if ( eTheta < pi )
1204 {
1205 dist2ETheta=rho2-p.z()*p.z()*tanETheta2;
1206 }
1207 else
1208 {
1209 dist2ETheta=kInfinity;
1210 }
1211 if ( pTheta < tolSTheta )
1212 {
1213 // Inside (theta<stheta-tol) stheta cone
1214 // First root of stheta cone, second if first root -ve
1215
1216 t1 = 1 - v.z()*v.z()*(1 + tanSTheta2) ;
1217 t2 = pDotV2d - p.z()*v.z()*tanSTheta2 ;
1218 if (t1)
1219 {
1220 b = t2/t1 ;
1221 c = dist2STheta/t1 ;
1222 d2 = b*b - c ;
1223
1224 if ( d2 >= 0 )
1225 {
1226 d = std::sqrt(d2) ;
1227 sd = -b - d ; // First root
1228 zi = p.z() + sd*v.z();
1229
1230 if ( (sd < 0) || (zi*(fSTheta - halfpi) > 0) )
1231 {
1232 sd = -b+d; // Second root
1233 }
1234 if ((sd >= 0) && (sd < snxt))
1235 {
1236 xi = p.x() + sd*v.x();
1237 yi = p.y() + sd*v.y();
1238 zi = p.z() + sd*v.z();
1239 rhoi2 = xi*xi + yi*yi;
1240 radi2 = rhoi2 + zi*zi;
1241 if ( (radi2 <= tolORMax2)
1242 && (radi2 >= tolORMin2)
1243 && (zi*(fSTheta - halfpi) <= 0) )
1244 {
1245 if ( !fFullPhiSphere && rhoi2 ) // Check phi intersection
1246 {
1247 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1248 if (cosPsi >= cosHDPhiOT)
1249 {
1250 snxt = sd;
1251 }
1252 }
1253 else
1254 {
1255 snxt = sd;
1256 }
1257 }
1258 }
1259 }
1260 }
1261
1262 // Possible intersection with ETheta cone.
1263 // Second >= 0 root should be considered
1264
1265 if ( eTheta < pi )
1266 {
1267 t1 = 1 - v.z()*v.z()*(1 + tanETheta2) ;
1268 t2 = pDotV2d - p.z()*v.z()*tanETheta2 ;
1269 if (t1)
1270 {
1271 b = t2/t1 ;
1272 c = dist2ETheta/t1 ;
1273 d2 = b*b - c ;
1274
1275 if (d2 >= 0)
1276 {
1277 d = std::sqrt(d2) ;
1278 sd = -b + d ; // Second root
1279
1280 if ( (sd >= 0) && (sd < snxt) )
1281 {
1282 xi = p.x() + sd*v.x() ;
1283 yi = p.y() + sd*v.y() ;
1284 zi = p.z() + sd*v.z() ;
1285 rhoi2 = xi*xi + yi*yi ;
1286 radi2 = rhoi2 + zi*zi ;
1287
1288 if ( (radi2 <= tolORMax2)
1289 && (radi2 >= tolORMin2)
1290 && (zi*(eTheta - halfpi) <= 0) )
1291 {
1292 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1293 {
1294 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1295 if (cosPsi >= cosHDPhiOT)
1296 {
1297 snxt = sd;
1298 }
1299 }
1300 else
1301 {
1302 snxt = sd;
1303 }
1304 }
1305 }
1306 }
1307 }
1308 }
1309 }
1310 else if ( pTheta > tolETheta )
1311 {
1312 // dist2ETheta<-kRadTolerance*0.5 && dist2STheta>0)
1313 // Inside (theta > etheta+tol) e-theta cone
1314 // First root of etheta cone, second if first root 'imaginary'
1315
1316 t1 = 1 - v.z()*v.z()*(1 + tanETheta2) ;
1317 t2 = pDotV2d - p.z()*v.z()*tanETheta2 ;
1318 if (t1)
1319 {
1320 b = t2/t1 ;
1321 c = dist2ETheta/t1 ;
1322 d2 = b*b - c ;
1323
1324 if (d2 >= 0)
1325 {
1326 d = std::sqrt(d2) ;
1327 sd = -b - d ; // First root
1328 zi = p.z() + sd*v.z();
1329
1330 if ( (sd < 0) || (zi*(eTheta - halfpi) > 0) )
1331 {
1332 sd = -b + d ; // second root
1333 }
1334 if ( (sd >= 0) && (sd < snxt) )
1335 {
1336 xi = p.x() + sd*v.x() ;
1337 yi = p.y() + sd*v.y() ;
1338 zi = p.z() + sd*v.z() ;
1339 rhoi2 = xi*xi + yi*yi ;
1340 radi2 = rhoi2 + zi*zi ;
1341
1342 if ( (radi2 <= tolORMax2)
1343 && (radi2 >= tolORMin2)
1344 && (zi*(eTheta - halfpi) <= 0) )
1345 {
1346 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1347 {
1348 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1349 if (cosPsi >= cosHDPhiOT)
1350 {
1351 snxt = sd;
1352 }
1353 }
1354 else
1355 {
1356 snxt = sd;
1357 }
1358 }
1359 }
1360 }
1361 }
1362
1363 // Possible intersection with STheta cone.
1364 // Second >= 0 root should be considered
1365
1366 if ( fSTheta )
1367 {
1368 t1 = 1 - v.z()*v.z()*(1 + tanSTheta2) ;
1369 t2 = pDotV2d - p.z()*v.z()*tanSTheta2 ;
1370 if (t1)
1371 {
1372 b = t2/t1 ;
1373 c = dist2STheta/t1 ;
1374 d2 = b*b - c ;
1375
1376 if (d2 >= 0)
1377 {
1378 d = std::sqrt(d2) ;
1379 sd = -b + d ; // Second root
1380
1381 if ( (sd >= 0) && (sd < snxt) )
1382 {
1383 xi = p.x() + sd*v.x() ;
1384 yi = p.y() + sd*v.y() ;
1385 zi = p.z() + sd*v.z() ;
1386 rhoi2 = xi*xi + yi*yi ;
1387 radi2 = rhoi2 + zi*zi ;
1388
1389 if ( (radi2 <= tolORMax2)
1390 && (radi2 >= tolORMin2)
1391 && (zi*(fSTheta - halfpi) <= 0) )
1392 {
1393 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1394 {
1395 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1396 if (cosPsi >= cosHDPhiOT)
1397 {
1398 snxt = sd;
1399 }
1400 }
1401 else
1402 {
1403 snxt = sd;
1404 }
1405 }
1406 }
1407 }
1408 }
1409 }
1410 }
1411 else if ( (pTheta < tolSTheta + kAngTolerance)
1412 && (fSTheta > halfAngTolerance) )
1413 {
1414 // In tolerance of stheta
1415 // If entering through solid [r,phi] => 0 to in
1416 // else try 2nd root
1417
1418 t2 = pDotV2d - p.z()*v.z()*tanSTheta2 ;
1419 if ( (t2>=0 && tolIRMin2<rad2 && rad2<tolIRMax2 && fSTheta<halfpi)
1420 || (t2<0 && tolIRMin2<rad2 && rad2<tolIRMax2 && fSTheta>halfpi)
1421 || (v.z()<0 && tolIRMin2<rad2 && rad2<tolIRMax2 && fSTheta==halfpi) )
1422 {
1423 if (!fFullPhiSphere && rho2) // Check phi intersection
1424 {
1425 cosPsi = (p.x()*cosCPhi + p.y()*sinCPhi)/std::sqrt(rho2) ;
1426 if (cosPsi >= cosHDPhiIT)
1427 {
1428 return 0 ;
1429 }
1430 }
1431 else
1432 {
1433 return 0 ;
1434 }
1435 }
1436
1437 // Not entering immediately/travelling through
1438
1439 t1 = 1 - v.z()*v.z()*(1 + tanSTheta2) ;
1440 if (t1)
1441 {
1442 b = t2/t1 ;
1443 c = dist2STheta/t1 ;
1444 d2 = b*b - c ;
1445
1446 if (d2 >= 0)
1447 {
1448 d = std::sqrt(d2) ;
1449 sd = -b + d ;
1450 if ( (sd >= halfCarTolerance) && (sd < snxt) && (fSTheta < halfpi) )
1451 { // ^^^^^^^^^^^^^^^^^^^^^ shouldn't it be >=0 instead ?
1452 xi = p.x() + sd*v.x() ;
1453 yi = p.y() + sd*v.y() ;
1454 zi = p.z() + sd*v.z() ;
1455 rhoi2 = xi*xi + yi*yi ;
1456 radi2 = rhoi2 + zi*zi ;
1457
1458 if ( (radi2 <= tolORMax2)
1459 && (radi2 >= tolORMin2)
1460 && (zi*(fSTheta - halfpi) <= 0) )
1461 {
1462 if ( !fFullPhiSphere && rhoi2 ) // Check phi intersection
1463 {
1464 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1465 if ( cosPsi >= cosHDPhiOT )
1466 {
1467 snxt = sd;
1468 }
1469 }
1470 else
1471 {
1472 snxt = sd;
1473 }
1474 }
1475 }
1476 }
1477 }
1478 }
1479 else if ((pTheta > tolETheta-kAngTolerance) && (eTheta < pi-kAngTolerance))
1480 {
1481
1482 // In tolerance of etheta
1483 // If entering through solid [r,phi] => 0 to in
1484 // else try 2nd root
1485
1486 t2 = pDotV2d - p.z()*v.z()*tanETheta2 ;
1487
1488 if ( ((t2<0) && (eTheta < halfpi)
1489 && (tolIRMin2 < rad2) && (rad2 < tolIRMax2))
1490 || ((t2>=0) && (eTheta > halfpi)
1491 && (tolIRMin2 < rad2) && (rad2 < tolIRMax2))
1492 || ((v.z()>0) && (eTheta == halfpi)
1493 && (tolIRMin2 < rad2) && (rad2 < tolIRMax2)) )
1494 {
1495 if (!fFullPhiSphere && rho2) // Check phi intersection
1496 {
1497 cosPsi = (p.x()*cosCPhi + p.y()*sinCPhi)/std::sqrt(rho2) ;
1498 if (cosPsi >= cosHDPhiIT)
1499 {
1500 return 0 ;
1501 }
1502 }
1503 else
1504 {
1505 return 0 ;
1506 }
1507 }
1508
1509 // Not entering immediately/travelling through
1510
1511 t1 = 1 - v.z()*v.z()*(1 + tanETheta2) ;
1512 if (t1)
1513 {
1514 b = t2/t1 ;
1515 c = dist2ETheta/t1 ;
1516 d2 = b*b - c ;
1517
1518 if (d2 >= 0)
1519 {
1520 d = std::sqrt(d2) ;
1521 sd = -b + d ;
1522
1523 if ( (sd >= halfCarTolerance)
1524 && (sd < snxt) && (eTheta > halfpi) )
1525 {
1526 xi = p.x() + sd*v.x() ;
1527 yi = p.y() + sd*v.y() ;
1528 zi = p.z() + sd*v.z() ;
1529 rhoi2 = xi*xi + yi*yi ;
1530 radi2 = rhoi2 + zi*zi ;
1531
1532 if ( (radi2 <= tolORMax2)
1533 && (radi2 >= tolORMin2)
1534 && (zi*(eTheta - halfpi) <= 0) )
1535 {
1536 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1537 {
1538 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1539 if (cosPsi >= cosHDPhiOT)
1540 {
1541 snxt = sd;
1542 }
1543 }
1544 else
1545 {
1546 snxt = sd;
1547 }
1548 }
1549 }
1550 }
1551 }
1552 }
1553 else
1554 {
1555 // stheta+tol<theta<etheta-tol
1556 // For BOTH stheta & etheta check 2nd root for validity [r,phi]
1557
1558 t1 = 1 - v.z()*v.z()*(1 + tanSTheta2) ;
1559 t2 = pDotV2d - p.z()*v.z()*tanSTheta2 ;
1560 if (t1)
1561 {
1562 b = t2/t1;
1563 c = dist2STheta/t1 ;
1564 d2 = b*b - c ;
1565
1566 if (d2 >= 0)
1567 {
1568 d = std::sqrt(d2) ;
1569 sd = -b + d ; // second root
1570
1571 if ((sd >= 0) && (sd < snxt))
1572 {
1573 xi = p.x() + sd*v.x() ;
1574 yi = p.y() + sd*v.y() ;
1575 zi = p.z() + sd*v.z() ;
1576 rhoi2 = xi*xi + yi*yi ;
1577 radi2 = rhoi2 + zi*zi ;
1578
1579 if ( (radi2 <= tolORMax2)
1580 && (radi2 >= tolORMin2)
1581 && (zi*(fSTheta - halfpi) <= 0) )
1582 {
1583 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1584 {
1585 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1586 if (cosPsi >= cosHDPhiOT)
1587 {
1588 snxt = sd;
1589 }
1590 }
1591 else
1592 {
1593 snxt = sd;
1594 }
1595 }
1596 }
1597 }
1598 }
1599 t1 = 1 - v.z()*v.z()*(1 + tanETheta2) ;
1600 t2 = pDotV2d - p.z()*v.z()*tanETheta2 ;
1601 if (t1)
1602 {
1603 b = t2/t1 ;
1604 c = dist2ETheta/t1 ;
1605 d2 = b*b - c ;
1606
1607 if (d2 >= 0)
1608 {
1609 d = std::sqrt(d2) ;
1610 sd = -b + d; // second root
1611
1612 if ((sd >= 0) && (sd < snxt))
1613 {
1614 xi = p.x() + sd*v.x() ;
1615 yi = p.y() + sd*v.y() ;
1616 zi = p.z() + sd*v.z() ;
1617 rhoi2 = xi*xi + yi*yi ;
1618 radi2 = rhoi2 + zi*zi ;
1619
1620 if ( (radi2 <= tolORMax2)
1621 && (radi2 >= tolORMin2)
1622 && (zi*(eTheta - halfpi) <= 0) )
1623 {
1624 if (!fFullPhiSphere && rhoi2) // Check phi intersection
1625 {
1626 cosPsi = (xi*cosCPhi + yi*sinCPhi)/std::sqrt(rhoi2) ;
1627 if ( cosPsi >= cosHDPhiOT )
1628 {
1629 snxt = sd;
1630 }
1631 }
1632 else
1633 {
1634 snxt = sd;
1635 }
1636 }
1637 }
1638 }
1639 }
1640 }
1641 }
1642 return snxt;
1643}
static const G4double d2
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4Sphere.cc:704

References cosCPhi, cosEPhi, cosHDPhiIT, cosHDPhiOT, cosSPhi, d2, DistanceToIn(), eTheta, fFullPhiSphere, fFullThetaSphere, fRmax, fRmaxTolerance, fRmin, fRminTolerance, fSTheta, halfAngTolerance, halfCarTolerance, halfpi, kAngTolerance, G4VSolid::kCarTolerance, kInfinity, pi, sinCPhi, sinEPhi, sinSPhi, tanETheta2, tanSTheta2, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by DistanceToIn().

◆ DistanceToOut() [1/2]

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

Implements G4VSolid.

Definition at line 2640 of file G4Sphere.cc.

2641{
2642 G4double safe=0.0,safeRMin,safeRMax,safePhi,safeTheta;
2643 G4double rho2,rds,rho;
2644 G4double pTheta,dTheta1 = kInfinity,dTheta2 = kInfinity;
2645 rho2=p.x()*p.x()+p.y()*p.y();
2646 rds=std::sqrt(rho2+p.z()*p.z());
2647 rho=std::sqrt(rho2);
2648
2649#ifdef G4CSGDEBUG
2650 if( Inside(p) == kOutside )
2651 {
2652 G4int old_prc = G4cout.precision(16);
2653 G4cout << G4endl;
2654 DumpInfo();
2655 G4cout << "Position:" << G4endl << G4endl ;
2656 G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl ;
2657 G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl ;
2658 G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl ;
2659 G4cout.precision(old_prc) ;
2660 G4Exception("G4Sphere::DistanceToOut(p)",
2661 "GeomSolids1002", JustWarning, "Point p is outside !?" );
2662 }
2663#endif
2664
2665 // Distance to r shells
2666 //
2667 safeRMax = fRmax-rds;
2668 safe = safeRMax;
2669 if (fRmin)
2670 {
2671 safeRMin = rds-fRmin;
2672 safe = std::min( safeRMin, safeRMax );
2673 }
2674
2675 // Distance to phi extent
2676 //
2677 if ( !fFullPhiSphere )
2678 {
2679 if (rho>0.0)
2680 {
2681 if ((p.y()*cosCPhi-p.x()*sinCPhi)<=0)
2682 {
2683 safePhi=-(p.x()*sinSPhi-p.y()*cosSPhi);
2684 }
2685 else
2686 {
2687 safePhi=(p.x()*sinEPhi-p.y()*cosEPhi);
2688 }
2689 }
2690 else
2691 {
2692 safePhi = 0.0; // Distance to both Phi surfaces (extended)
2693 }
2694 // Both cases above can be improved - in case fRMin > 0.0
2695 // although it may be costlier (good for precise, not fast version)
2696
2697 safe= std::min(safe, safePhi);
2698 }
2699
2700 // Distance to Theta extent
2701 //
2702 if ( !fFullThetaSphere )
2703 {
2704 if( rds > 0.0 )
2705 {
2706 pTheta=std::acos(p.z()/rds);
2707 if (pTheta<0) { pTheta+=pi; }
2708 if(fSTheta>0.)
2709 { dTheta1=pTheta-fSTheta;}
2710 if(eTheta<pi)
2711 { dTheta2=eTheta-pTheta;}
2712
2713 safeTheta=rds*std::sin(std::min(dTheta1, dTheta2) );
2714 }
2715 else
2716 {
2717 safeTheta= 0.0;
2718 // An improvement will be to return negative answer if outside (TODO)
2719 }
2720 safe = std::min( safe, safeTheta );
2721 }
2722
2723 if (safe<0.0) { safe=0; }
2724 // An improvement to return negative answer if outside (TODO)
2725
2726 return safe;
2727}
static constexpr double mm
Definition: G4SIunits.hh:95
G4GLOB_DLL std::ostream G4cout
EInside Inside(const G4ThreeVector &p) const
Definition: G4Sphere.cc:288
@ kOutside
Definition: geomdefs.hh:68

References cosCPhi, cosEPhi, cosSPhi, G4VSolid::DumpInfo(), eTheta, fFullPhiSphere, fFullThetaSphere, fRmax, fRmin, fSTheta, G4cout, G4endl, G4Exception(), Inside(), JustWarning, kInfinity, kOutside, G4INCL::Math::min(), mm, pi, sinCPhi, sinEPhi, sinSPhi, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

◆ DistanceToOut() [2/2]

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

Implements G4VSolid.

Definition at line 1749 of file G4Sphere.cc.

1754{
1755 G4double snxt = kInfinity; // snxt is default return value
1756 G4double sphi= kInfinity,stheta= kInfinity;
1757 ESide side=kNull,sidephi=kNull,sidetheta=kNull;
1758
1759 const G4double halfRmaxTolerance = fRmaxTolerance*0.5;
1760 const G4double halfRminTolerance = fRminTolerance*0.5;
1761 const G4double Rmax_plus = fRmax + halfRmaxTolerance;
1762 const G4double Rmin_minus = (fRmin) ? fRmin-halfRminTolerance : 0;
1763 G4double t1,t2;
1764 G4double b,c,d;
1765
1766 // Variables for phi intersection:
1767
1768 G4double pDistS,compS,pDistE,compE,sphi2,vphi;
1769
1770 G4double rho2,rad2,pDotV2d,pDotV3d;
1771
1772 G4double xi,yi,zi; // Intersection point
1773
1774 // Theta precals
1775 //
1776 G4double rhoSecTheta;
1777 G4double dist2STheta, dist2ETheta, distTheta;
1778 G4double d2,sd;
1779
1780 // General Precalcs
1781 //
1782 rho2 = p.x()*p.x()+p.y()*p.y();
1783 rad2 = rho2+p.z()*p.z();
1784
1785 pDotV2d = p.x()*v.x()+p.y()*v.y();
1786 pDotV3d = pDotV2d+p.z()*v.z();
1787
1788 // Radial Intersections from G4Sphere::DistanceToIn
1789 //
1790 // Outer spherical shell intersection
1791 // - Only if outside tolerant fRmax
1792 // - Check for if inside and outer G4Sphere heading through solid (-> 0)
1793 // - No intersect -> no intersection with G4Sphere
1794 //
1795 // Shell eqn: x^2+y^2+z^2=RSPH^2
1796 //
1797 // => (px+svx)^2+(py+svy)^2+(pz+svz)^2=R^2
1798 //
1799 // => (px^2+py^2+pz^2) +2sd(pxvx+pyvy+pzvz)+sd^2(vx^2+vy^2+vz^2)=R^2
1800 // => rad2 +2sd(pDotV3d) +sd^2 =R^2
1801 //
1802 // => sd=-pDotV3d+-std::sqrt(pDotV3d^2-(rad2-R^2))
1803
1804 if( (rad2 <= Rmax_plus*Rmax_plus) && (rad2 >= Rmin_minus*Rmin_minus) )
1805 {
1806 c = rad2 - fRmax*fRmax;
1807
1808 if (c < fRmaxTolerance*fRmax)
1809 {
1810 // Within tolerant Outer radius
1811 //
1812 // The test is
1813 // rad - fRmax < 0.5*kRadTolerance
1814 // => rad < fRmax + 0.5*kRadTol
1815 // => rad2 < (fRmax + 0.5*kRadTol)^2
1816 // => rad2 < fRmax^2 + 2.*0.5*fRmax*kRadTol + 0.25*kRadTol*kRadTol
1817 // => rad2 - fRmax^2 <~ fRmax*kRadTol
1818
1819 d2 = pDotV3d*pDotV3d - c;
1820
1821 if( (c >- fRmaxTolerance*fRmax) // on tolerant surface
1822 && ((pDotV3d >=0) || (d2 < 0)) ) // leaving outside from Rmax
1823 // not re-entering
1824 {
1825 if(calcNorm)
1826 {
1827 *validNorm = true ;
1828 *n = G4ThreeVector(p.x()/fRmax,p.y()/fRmax,p.z()/fRmax) ;
1829 }
1830 return snxt = 0;
1831 }
1832 else
1833 {
1834 snxt = -pDotV3d+std::sqrt(d2); // second root since inside Rmax
1835 side = kRMax ;
1836 }
1837 }
1838
1839 // Inner spherical shell intersection:
1840 // Always first >=0 root, because would have passed
1841 // from outside of Rmin surface .
1842
1843 if (fRmin)
1844 {
1845 c = rad2 - fRmin*fRmin;
1846 d2 = pDotV3d*pDotV3d - c;
1847
1848 if (c >- fRminTolerance*fRmin) // 2.0 * (0.5*kRadTolerance) * fRmin
1849 {
1850 if ( (c < fRminTolerance*fRmin) // leaving from Rmin
1851 && (d2 >= fRminTolerance*fRmin) && (pDotV3d < 0) )
1852 {
1853 if(calcNorm) { *validNorm = false; } // Rmin surface is concave
1854 return snxt = 0 ;
1855 }
1856 else
1857 {
1858 if ( d2 >= 0. )
1859 {
1860 sd = -pDotV3d-std::sqrt(d2);
1861
1862 if ( sd >= 0. ) // Always intersect Rmin first
1863 {
1864 snxt = sd ;
1865 side = kRMin ;
1866 }
1867 }
1868 }
1869 }
1870 }
1871 }
1872
1873 // Theta segment intersection
1874
1875 if ( !fFullThetaSphere )
1876 {
1877 // Intersection with theta surfaces
1878 //
1879 // Known failure cases:
1880 // o Inside tolerance of stheta surface, skim
1881 // ~parallel to cone and Hit & enter etheta surface [& visa versa]
1882 //
1883 // To solve: Check 2nd root of etheta surface in addition to stheta
1884 //
1885 // o start/end theta is exactly pi/2
1886 //
1887 // Intersections with cones
1888 //
1889 // Cone equation: x^2+y^2=z^2tan^2(t)
1890 //
1891 // => (px+svx)^2+(py+svy)^2=(pz+svz)^2tan^2(t)
1892 //
1893 // => (px^2+py^2-pz^2tan^2(t))+2sd(pxvx+pyvy-pzvztan^2(t))
1894 // + sd^2(vx^2+vy^2-vz^2tan^2(t)) = 0
1895 //
1896 // => sd^2(1-vz^2(1+tan^2(t))+2sd(pdotv2d-pzvztan^2(t))
1897 // + (rho2-pz^2tan^2(t)) = 0
1898 //
1899
1900 if(fSTheta) // intersection with first cons
1901 {
1902 if( std::fabs(tanSTheta) > 5./kAngTolerance ) // kons is plane z=0
1903 {
1904 if( v.z() > 0. )
1905 {
1906 if ( std::fabs( p.z() ) <= halfRmaxTolerance )
1907 {
1908 if(calcNorm)
1909 {
1910 *validNorm = true;
1911 *n = G4ThreeVector(0.,0.,1.);
1912 }
1913 return snxt = 0 ;
1914 }
1915 stheta = -p.z()/v.z();
1916 sidetheta = kSTheta;
1917 }
1918 }
1919 else // kons is not plane
1920 {
1921 t1 = 1-v.z()*v.z()*(1+tanSTheta2);
1922 t2 = pDotV2d-p.z()*v.z()*tanSTheta2; // ~vDotN if p on cons
1923 dist2STheta = rho2-p.z()*p.z()*tanSTheta2; // t3
1924
1925 distTheta = std::sqrt(rho2)-p.z()*tanSTheta;
1926
1927 if( std::fabs(t1) < halfAngTolerance ) // 1st order equation,
1928 { // v parallel to kons
1929 if( v.z() > 0. )
1930 {
1931 if(std::fabs(distTheta) < halfRmaxTolerance) // p on surface
1932 {
1933 if( (fSTheta < halfpi) && (p.z() > 0.) )
1934 {
1935 if( calcNorm ) { *validNorm = false; }
1936 return snxt = 0.;
1937 }
1938 else if( (fSTheta > halfpi) && (p.z() <= 0) )
1939 {
1940 if( calcNorm )
1941 {
1942 *validNorm = true;
1943 if (rho2)
1944 {
1945 rhoSecTheta = std::sqrt(rho2*(1+tanSTheta2));
1946
1947 *n = G4ThreeVector( p.x()/rhoSecTheta,
1948 p.y()/rhoSecTheta,
1949 std::sin(fSTheta) );
1950 }
1951 else *n = G4ThreeVector(0.,0.,1.);
1952 }
1953 return snxt = 0.;
1954 }
1955 }
1956 stheta = -0.5*dist2STheta/t2;
1957 sidetheta = kSTheta;
1958 }
1959 } // 2nd order equation, 1st root of fSTheta cone,
1960 else // 2nd if 1st root -ve
1961 {
1962 if( std::fabs(distTheta) < halfRmaxTolerance )
1963 {
1964 if( (fSTheta > halfpi) && (t2 >= 0.) ) // leave
1965 {
1966 if( calcNorm )
1967 {
1968 *validNorm = true;
1969 if (rho2)
1970 {
1971 rhoSecTheta = std::sqrt(rho2*(1+tanSTheta2));
1972
1973 *n = G4ThreeVector( p.x()/rhoSecTheta,
1974 p.y()/rhoSecTheta,
1975 std::sin(fSTheta) );
1976 }
1977 else { *n = G4ThreeVector(0.,0.,1.); }
1978 }
1979 return snxt = 0.;
1980 }
1981 else if( (fSTheta < halfpi) && (t2 < 0.) && (p.z() >=0.) ) // leave
1982 {
1983 if( calcNorm ) { *validNorm = false; }
1984 return snxt = 0.;
1985 }
1986 }
1987 b = t2/t1;
1988 c = dist2STheta/t1;
1989 d2 = b*b - c ;
1990
1991 if ( d2 >= 0. )
1992 {
1993 d = std::sqrt(d2);
1994
1995 if( fSTheta > halfpi )
1996 {
1997 sd = -b - d; // First root
1998
1999 if ( ((std::fabs(s) < halfRmaxTolerance) && (t2 < 0.))
2000 || (sd < 0.) || ( (sd > 0.) && (p.z() + sd*v.z() > 0.) ) )
2001 {
2002 sd = -b + d ; // 2nd root
2003 }
2004 if( (sd > halfRmaxTolerance) && (p.z() + sd*v.z() <= 0.) )
2005 {
2006 stheta = sd;
2007 sidetheta = kSTheta;
2008 }
2009 }
2010 else // sTheta < pi/2, concave surface, no normal
2011 {
2012 sd = -b - d; // First root
2013
2014 if ( ( (std::fabs(sd) < halfRmaxTolerance) && (t2 >= 0.) )
2015 || (sd < 0.) || ( (sd > 0.) && (p.z() + sd*v.z() < 0.) ) )
2016 {
2017 sd = -b + d ; // 2nd root
2018 }
2019 if( (sd > halfRmaxTolerance) && (p.z() + sd*v.z() >= 0.) )
2020 {
2021 stheta = sd;
2022 sidetheta = kSTheta;
2023 }
2024 }
2025 }
2026 }
2027 }
2028 }
2029 if (eTheta < pi) // intersection with second cons
2030 {
2031 if( std::fabs(tanETheta) > 5./kAngTolerance ) // kons is plane z=0
2032 {
2033 if( v.z() < 0. )
2034 {
2035 if ( std::fabs( p.z() ) <= halfRmaxTolerance )
2036 {
2037 if(calcNorm)
2038 {
2039 *validNorm = true;
2040 *n = G4ThreeVector(0.,0.,-1.);
2041 }
2042 return snxt = 0 ;
2043 }
2044 sd = -p.z()/v.z();
2045
2046 if( sd < stheta )
2047 {
2048 stheta = sd;
2049 sidetheta = kETheta;
2050 }
2051 }
2052 }
2053 else // kons is not plane
2054 {
2055 t1 = 1-v.z()*v.z()*(1+tanETheta2);
2056 t2 = pDotV2d-p.z()*v.z()*tanETheta2; // ~vDotN if p on cons
2057 dist2ETheta = rho2-p.z()*p.z()*tanETheta2; // t3
2058
2059 distTheta = std::sqrt(rho2)-p.z()*tanETheta;
2060
2061 if( std::fabs(t1) < halfAngTolerance ) // 1st order equation,
2062 { // v parallel to kons
2063 if( v.z() < 0. )
2064 {
2065 if(std::fabs(distTheta) < halfRmaxTolerance) // p on surface
2066 {
2067 if( (eTheta > halfpi) && (p.z() < 0.) )
2068 {
2069 if( calcNorm ) { *validNorm = false; }
2070 return snxt = 0.;
2071 }
2072 else if ( (eTheta < halfpi) && (p.z() >= 0) )
2073 {
2074 if( calcNorm )
2075 {
2076 *validNorm = true;
2077 if (rho2)
2078 {
2079 rhoSecTheta = std::sqrt(rho2*(1+tanETheta2));
2080 *n = G4ThreeVector( p.x()/rhoSecTheta,
2081 p.y()/rhoSecTheta,
2082 -sinETheta );
2083 }
2084 else { *n = G4ThreeVector(0.,0.,-1.); }
2085 }
2086 return snxt = 0.;
2087 }
2088 }
2089 sd = -0.5*dist2ETheta/t2;
2090
2091 if( sd < stheta )
2092 {
2093 stheta = sd;
2094 sidetheta = kETheta;
2095 }
2096 }
2097 } // 2nd order equation, 1st root of fSTheta cone
2098 else // 2nd if 1st root -ve
2099 {
2100 if ( std::fabs(distTheta) < halfRmaxTolerance )
2101 {
2102 if( (eTheta < halfpi) && (t2 >= 0.) ) // leave
2103 {
2104 if( calcNorm )
2105 {
2106 *validNorm = true;
2107 if (rho2)
2108 {
2109 rhoSecTheta = std::sqrt(rho2*(1+tanETheta2));
2110 *n = G4ThreeVector( p.x()/rhoSecTheta,
2111 p.y()/rhoSecTheta,
2112 -sinETheta );
2113 }
2114 else *n = G4ThreeVector(0.,0.,-1.);
2115 }
2116 return snxt = 0.;
2117 }
2118 else if ( (eTheta > halfpi)
2119 && (t2 < 0.) && (p.z() <=0.) ) // leave
2120 {
2121 if( calcNorm ) { *validNorm = false; }
2122 return snxt = 0.;
2123 }
2124 }
2125 b = t2/t1;
2126 c = dist2ETheta/t1;
2127 d2 = b*b - c ;
2128 if ( (d2 <halfRmaxTolerance) && (d2 > -halfRmaxTolerance) )
2129 {
2130 d2 = 0.;
2131 }
2132 if ( d2 >= 0. )
2133 {
2134 d = std::sqrt(d2);
2135
2136 if( eTheta < halfpi )
2137 {
2138 sd = -b - d; // First root
2139
2140 if( ((std::fabs(sd) < halfRmaxTolerance) && (t2 < 0.))
2141 || (sd < 0.) )
2142 {
2143 sd = -b + d ; // 2nd root
2144 }
2145 if( sd > halfRmaxTolerance )
2146 {
2147 if( sd < stheta )
2148 {
2149 stheta = sd;
2150 sidetheta = kETheta;
2151 }
2152 }
2153 }
2154 else // sTheta+fDTheta > pi/2, concave surface, no normal
2155 {
2156 sd = -b - d; // First root
2157
2158 if ( ((std::fabs(sd) < halfRmaxTolerance) && (t2 >= 0.))
2159 || (sd < 0.)
2160 || ( (sd > 0.) && (p.z() + sd*v.z() > halfRmaxTolerance) ) )
2161 {
2162 sd = -b + d ; // 2nd root
2163 }
2164 if ( ( sd>halfRmaxTolerance )
2165 && ( p.z()+sd*v.z() <= halfRmaxTolerance ) )
2166 {
2167 if( sd < stheta )
2168 {
2169 stheta = sd;
2170 sidetheta = kETheta;
2171 }
2172 }
2173 }
2174 }
2175 }
2176 }
2177 }
2178
2179 } // end theta intersections
2180
2181 // Phi Intersection
2182
2183 if ( !fFullPhiSphere )
2184 {
2185 if ( p.x() || p.y() ) // Check if on z axis (rho not needed later)
2186 {
2187 // pDist -ve when inside
2188
2189 pDistS=p.x()*sinSPhi-p.y()*cosSPhi;
2190 pDistE=-p.x()*sinEPhi+p.y()*cosEPhi;
2191
2192 // Comp -ve when in direction of outwards normal
2193
2194 compS = -sinSPhi*v.x()+cosSPhi*v.y() ;
2195 compE = sinEPhi*v.x()-cosEPhi*v.y() ;
2196 sidephi = kNull ;
2197
2198 if ( (pDistS <= 0) && (pDistE <= 0) )
2199 {
2200 // Inside both phi *full* planes
2201
2202 if ( compS < 0 )
2203 {
2204 sphi = pDistS/compS ;
2205 xi = p.x()+sphi*v.x() ;
2206 yi = p.y()+sphi*v.y() ;
2207
2208 // Check intersection with correct half-plane (if not -> no intersect)
2209 //
2210 if( (std::fabs(xi)<=kCarTolerance) && (std::fabs(yi)<=kCarTolerance) )
2211 {
2212 vphi = std::atan2(v.y(),v.x());
2213 sidephi = kSPhi;
2214 if ( ( (fSPhi-halfAngTolerance) <= vphi)
2215 && ( (ePhi+halfAngTolerance) >= vphi) )
2216 {
2217 sphi = kInfinity;
2218 }
2219 }
2220 else if ( ( yi*cosCPhi - xi*sinCPhi ) >= 0 )
2221 {
2222 sphi=kInfinity;
2223 }
2224 else
2225 {
2226 sidephi = kSPhi ;
2227 if ( pDistS > -halfCarTolerance) { sphi = 0; } // Leave by sphi
2228 }
2229 }
2230 else { sphi = kInfinity; }
2231
2232 if ( compE < 0 )
2233 {
2234 sphi2=pDistE/compE ;
2235 if (sphi2 < sphi) // Only check further if < starting phi intersection
2236 {
2237 xi = p.x()+sphi2*v.x() ;
2238 yi = p.y()+sphi2*v.y() ;
2239
2240 // Check intersection with correct half-plane
2241 //
2242 if ( (std::fabs(xi)<=kCarTolerance)
2243 && (std::fabs(yi)<=kCarTolerance))
2244 {
2245 // Leaving via ending phi
2246 //
2247 vphi = std::atan2(v.y(),v.x()) ;
2248
2249 if( !((fSPhi-halfAngTolerance <= vphi)
2250 &&(fSPhi+fDPhi+halfAngTolerance >= vphi)) )
2251 {
2252 sidephi = kEPhi;
2253 if ( pDistE <= -halfCarTolerance ) { sphi = sphi2; }
2254 else { sphi = 0.0; }
2255 }
2256 }
2257 else if ((yi*cosCPhi-xi*sinCPhi)>=0) // Leaving via ending phi
2258 {
2259 sidephi = kEPhi ;
2260 if ( pDistE <= -halfCarTolerance )
2261 {
2262 sphi=sphi2;
2263 }
2264 else
2265 {
2266 sphi = 0 ;
2267 }
2268 }
2269 }
2270 }
2271 }
2272 else if ((pDistS >= 0) && (pDistE >= 0)) // Outside both *full* phi planes
2273 {
2274 if ( pDistS <= pDistE )
2275 {
2276 sidephi = kSPhi ;
2277 }
2278 else
2279 {
2280 sidephi = kEPhi ;
2281 }
2282 if ( fDPhi > pi )
2283 {
2284 if ( (compS < 0) && (compE < 0) ) { sphi = 0; }
2285 else { sphi = kInfinity; }
2286 }
2287 else
2288 {
2289 // if towards both >=0 then once inside (after error)
2290 // will remain inside
2291
2292 if ( (compS >= 0) && (compE >= 0) ) { sphi = kInfinity; }
2293 else { sphi = 0; }
2294 }
2295 }
2296 else if ( (pDistS > 0) && (pDistE < 0) )
2297 {
2298 // Outside full starting plane, inside full ending plane
2299
2300 if ( fDPhi > pi )
2301 {
2302 if ( compE < 0 )
2303 {
2304 sphi = pDistE/compE ;
2305 xi = p.x() + sphi*v.x() ;
2306 yi = p.y() + sphi*v.y() ;
2307
2308 // Check intersection in correct half-plane
2309 // (if not -> not leaving phi extent)
2310 //
2311 if( (std::fabs(xi)<=kCarTolerance)&&(std::fabs(yi)<=kCarTolerance) )
2312 {
2313 vphi = std::atan2(v.y(),v.x());
2314 sidephi = kSPhi;
2315 if ( ( (fSPhi-halfAngTolerance) <= vphi)
2316 && ( (ePhi+halfAngTolerance) >= vphi) )
2317 {
2318 sphi = kInfinity;
2319 }
2320 }
2321 else if ( ( yi*cosCPhi - xi*sinCPhi ) <= 0 )
2322 {
2323 sphi = kInfinity ;
2324 }
2325 else // Leaving via Ending phi
2326 {
2327 sidephi = kEPhi ;
2328 if ( pDistE > -halfCarTolerance ) { sphi = 0.; }
2329 }
2330 }
2331 else
2332 {
2333 sphi = kInfinity ;
2334 }
2335 }
2336 else
2337 {
2338 if ( compS >= 0 )
2339 {
2340 if ( compE < 0 )
2341 {
2342 sphi = pDistE/compE ;
2343 xi = p.x() + sphi*v.x() ;
2344 yi = p.y() + sphi*v.y() ;
2345
2346 // Check intersection in correct half-plane
2347 // (if not -> remain in extent)
2348 //
2349 if( (std::fabs(xi)<=kCarTolerance)
2350 && (std::fabs(yi)<=kCarTolerance) )
2351 {
2352 vphi = std::atan2(v.y(),v.x());
2353 sidephi = kSPhi;
2354 if ( ( (fSPhi-halfAngTolerance) <= vphi)
2355 && ( (ePhi+halfAngTolerance) >= vphi) )
2356 {
2357 sphi = kInfinity;
2358 }
2359 }
2360 else if ( ( yi*cosCPhi - xi*sinCPhi) <= 0 )
2361 {
2362 sphi=kInfinity;
2363 }
2364 else // otherwise leaving via Ending phi
2365 {
2366 sidephi = kEPhi ;
2367 }
2368 }
2369 else sphi=kInfinity;
2370 }
2371 else // leaving immediately by starting phi
2372 {
2373 sidephi = kSPhi ;
2374 sphi = 0 ;
2375 }
2376 }
2377 }
2378 else
2379 {
2380 // Must be pDistS < 0 && pDistE > 0
2381 // Inside full starting plane, outside full ending plane
2382
2383 if ( fDPhi > pi )
2384 {
2385 if ( compS < 0 )
2386 {
2387 sphi=pDistS/compS;
2388 xi=p.x()+sphi*v.x();
2389 yi=p.y()+sphi*v.y();
2390
2391 // Check intersection in correct half-plane
2392 // (if not -> not leaving phi extent)
2393 //
2394 if( (std::fabs(xi)<=kCarTolerance)&&(std::fabs(yi)<=kCarTolerance) )
2395 {
2396 vphi = std::atan2(v.y(),v.x()) ;
2397 sidephi = kSPhi;
2398 if ( ( (fSPhi-halfAngTolerance) <= vphi)
2399 && ( (ePhi+halfAngTolerance) >= vphi) )
2400 {
2401 sphi = kInfinity;
2402 }
2403 }
2404 else if ( ( yi*cosCPhi - xi*sinCPhi ) >= 0 )
2405 {
2406 sphi = kInfinity ;
2407 }
2408 else // Leaving via Starting phi
2409 {
2410 sidephi = kSPhi ;
2411 if ( pDistS > -halfCarTolerance ) { sphi = 0; }
2412 }
2413 }
2414 else
2415 {
2416 sphi = kInfinity ;
2417 }
2418 }
2419 else
2420 {
2421 if ( compE >= 0 )
2422 {
2423 if ( compS < 0 )
2424 {
2425 sphi = pDistS/compS ;
2426 xi = p.x()+sphi*v.x() ;
2427 yi = p.y()+sphi*v.y() ;
2428
2429 // Check intersection in correct half-plane
2430 // (if not -> remain in extent)
2431 //
2432 if( (std::fabs(xi)<=kCarTolerance)
2433 && (std::fabs(yi)<=kCarTolerance))
2434 {
2435 vphi = std::atan2(v.y(),v.x()) ;
2436 sidephi = kSPhi;
2437 if ( ( (fSPhi-halfAngTolerance) <= vphi)
2438 && ( (ePhi+halfAngTolerance) >= vphi) )
2439 {
2440 sphi = kInfinity;
2441 }
2442 }
2443 else if ( ( yi*cosCPhi - xi*sinCPhi ) >= 0 )
2444 {
2445 sphi = kInfinity ;
2446 }
2447 else // otherwise leaving via Starting phi
2448 {
2449 sidephi = kSPhi ;
2450 }
2451 }
2452 else
2453 {
2454 sphi = kInfinity ;
2455 }
2456 }
2457 else // leaving immediately by ending
2458 {
2459 sidephi = kEPhi ;
2460 sphi = 0 ;
2461 }
2462 }
2463 }
2464 }
2465 else
2466 {
2467 // On z axis + travel not || to z axis -> if phi of vector direction
2468 // within phi of shape, Step limited by rmax, else Step =0
2469
2470 if ( v.x() || v.y() )
2471 {
2472 vphi = std::atan2(v.y(),v.x()) ;
2473 if ((fSPhi-halfAngTolerance < vphi) && (vphi < ePhi+halfAngTolerance))
2474 {
2475 sphi = kInfinity;
2476 }
2477 else
2478 {
2479 sidephi = kSPhi ; // arbitrary
2480 sphi = 0 ;
2481 }
2482 }
2483 else // travel along z - no phi intersection
2484 {
2485 sphi = kInfinity ;
2486 }
2487 }
2488 if ( sphi < snxt ) // Order intersecttions
2489 {
2490 snxt = sphi ;
2491 side = sidephi ;
2492 }
2493 }
2494 if (stheta < snxt ) // Order intersections
2495 {
2496 snxt = stheta ;
2497 side = sidetheta ;
2498 }
2499
2500 if (calcNorm) // Output switch operator
2501 {
2502 switch( side )
2503 {
2504 case kRMax:
2505 xi=p.x()+snxt*v.x();
2506 yi=p.y()+snxt*v.y();
2507 zi=p.z()+snxt*v.z();
2508 *n=G4ThreeVector(xi/fRmax,yi/fRmax,zi/fRmax);
2509 *validNorm=true;
2510 break;
2511
2512 case kRMin:
2513 *validNorm=false; // Rmin is concave
2514 break;
2515
2516 case kSPhi:
2517 if ( fDPhi <= pi ) // Normal to Phi-
2518 {
2520 *validNorm=true;
2521 }
2522 else { *validNorm=false; }
2523 break ;
2524
2525 case kEPhi:
2526 if ( fDPhi <= pi ) // Normal to Phi+
2527 {
2529 *validNorm=true;
2530 }
2531 else { *validNorm=false; }
2532 break;
2533
2534 case kSTheta:
2535 if( fSTheta == halfpi )
2536 {
2537 *n=G4ThreeVector(0.,0.,1.);
2538 *validNorm=true;
2539 }
2540 else if ( fSTheta > halfpi )
2541 {
2542 xi = p.x() + snxt*v.x();
2543 yi = p.y() + snxt*v.y();
2544 rho2=xi*xi+yi*yi;
2545 if (rho2)
2546 {
2547 rhoSecTheta = std::sqrt(rho2*(1+tanSTheta2));
2548 *n = G4ThreeVector( xi/rhoSecTheta, yi/rhoSecTheta,
2549 -tanSTheta/std::sqrt(1+tanSTheta2));
2550 }
2551 else
2552 {
2553 *n = G4ThreeVector(0.,0.,1.);
2554 }
2555 *validNorm=true;
2556 }
2557 else { *validNorm=false; } // Concave STheta cone
2558 break;
2559
2560 case kETheta:
2561 if( eTheta == halfpi )
2562 {
2563 *n = G4ThreeVector(0.,0.,-1.);
2564 *validNorm = true;
2565 }
2566 else if ( eTheta < halfpi )
2567 {
2568 xi=p.x()+snxt*v.x();
2569 yi=p.y()+snxt*v.y();
2570 rho2=xi*xi+yi*yi;
2571 if (rho2)
2572 {
2573 rhoSecTheta = std::sqrt(rho2*(1+tanETheta2));
2574 *n = G4ThreeVector( xi/rhoSecTheta, yi/rhoSecTheta,
2575 -tanETheta/std::sqrt(1+tanETheta2) );
2576 }
2577 else
2578 {
2579 *n = G4ThreeVector(0.,0.,-1.);
2580 }
2581 *validNorm=true;
2582 }
2583 else { *validNorm=false; } // Concave ETheta cone
2584 break;
2585
2586 default:
2587 G4cout << G4endl;
2588 DumpInfo();
2589 std::ostringstream message;
2590 G4int oldprc = message.precision(16);
2591 message << "Undefined side for valid surface normal to solid."
2592 << G4endl
2593 << "Position:" << G4endl << G4endl
2594 << "p.x() = " << p.x()/mm << " mm" << G4endl
2595 << "p.y() = " << p.y()/mm << " mm" << G4endl
2596 << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl
2597 << "Direction:" << G4endl << G4endl
2598 << "v.x() = " << v.x() << G4endl
2599 << "v.y() = " << v.y() << G4endl
2600 << "v.z() = " << v.z() << G4endl << G4endl
2601 << "Proposed distance :" << G4endl << G4endl
2602 << "snxt = " << snxt/mm << " mm" << G4endl;
2603 message.precision(oldprc);
2604 G4Exception("G4Sphere::DistanceToOut(p,v,..)",
2605 "GeomSolids1002", JustWarning, message);
2606 break;
2607 }
2608 }
2609 if (snxt == kInfinity)
2610 {
2611 G4cout << G4endl;
2612 DumpInfo();
2613 std::ostringstream message;
2614 G4int oldprc = message.precision(16);
2615 message << "Logic error: snxt = kInfinity ???" << G4endl
2616 << "Position:" << G4endl << G4endl
2617 << "p.x() = " << p.x()/mm << " mm" << G4endl
2618 << "p.y() = " << p.y()/mm << " mm" << G4endl
2619 << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl
2620 << "Rp = "<< std::sqrt( p.x()*p.x()+p.y()*p.y()+p.z()*p.z() )/mm
2621 << " mm" << G4endl << G4endl
2622 << "Direction:" << G4endl << G4endl
2623 << "v.x() = " << v.x() << G4endl
2624 << "v.y() = " << v.y() << G4endl
2625 << "v.z() = " << v.z() << G4endl << G4endl
2626 << "Proposed distance :" << G4endl << G4endl
2627 << "snxt = " << snxt/mm << " mm" << G4endl;
2628 message.precision(oldprc);
2629 G4Exception("G4Sphere::DistanceToOut(p,v,..)",
2630 "GeomSolids1002", JustWarning, message);
2631 }
2632
2633 return snxt;
2634}
ESide
Definition: G4Cons.cc:61
static constexpr double s
Definition: G4SIunits.hh:154

References cosCPhi, cosEPhi, cosSPhi, d2, G4VSolid::DumpInfo(), ePhi, eTheta, fDPhi, fFullPhiSphere, fFullThetaSphere, fRmax, fRmaxTolerance, fRmin, fRminTolerance, fSPhi, fSTheta, G4cout, G4endl, G4Exception(), halfAngTolerance, halfCarTolerance, halfpi, JustWarning, kAngTolerance, G4VSolid::kCarTolerance, kEPhi, kETheta, kInfinity, kNull, kRMax, kRMin, kSPhi, kSTheta, mm, CLHEP::detail::n, pi, s, sinCPhi, sinEPhi, sinETheta, sinSPhi, tanETheta, tanETheta2, tanSTheta, tanSTheta2, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

◆ DumpInfo()

void G4VSolid::DumpInfo ( ) const
inlineinherited

Referenced by G4Cons::ApproxSurfaceNormal(), G4CutTubs::ApproxSurfaceNormal(), ApproxSurfaceNormal(), G4Torus::ApproxSurfaceNormal(), G4Tubs::ApproxSurfaceNormal(), G4ReflectedSolid::BoundingLimits(), G4DisplacedSolid::BoundingLimits(), G4IntersectionSolid::BoundingLimits(), G4ScaledSolid::BoundingLimits(), G4SubtractionSolid::BoundingLimits(), G4UnionSolid::BoundingLimits(), G4Box::BoundingLimits(), G4Cons::BoundingLimits(), G4CutTubs::BoundingLimits(), G4Orb::BoundingLimits(), G4Para::BoundingLimits(), 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(), 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
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
EInside
Definition: geomdefs.hh:67

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 G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=nullptr, G4ThreeVector *n=nullptr) const =0
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
virtual void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:665
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
@ kInside
Definition: geomdefs.hh:70

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().

◆ GetConstituentSolid() [1/2]

G4VSolid * G4VSolid::GetConstituentSolid ( G4int  no)
virtualinherited

Reimplemented in G4BooleanSolid.

Definition at line 170 of file G4VSolid.cc.

171{ return nullptr; }

◆ GetConstituentSolid() [2/2]

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

Reimplemented in G4BooleanSolid.

Definition at line 167 of file G4VSolid.cc.

168{ return nullptr; }

Referenced by G4BooleanSolid::StackPolyhedron().

◆ GetCosEndPhi()

G4double G4Sphere::GetCosEndPhi ( ) const
inline

Referenced by BoundingLimits().

◆ GetCosEndTheta()

G4double G4Sphere::GetCosEndTheta ( ) const
inline

Referenced by BoundingLimits().

◆ GetCosStartPhi()

G4double G4Sphere::GetCosStartPhi ( ) const
inline

Referenced by BoundingLimits().

◆ GetCosStartTheta()

G4double G4Sphere::GetCosStartTheta ( ) const
inline

Referenced by BoundingLimits().

◆ GetCubicVolume()

G4double G4Sphere::GetCubicVolume ( )
virtual

Reimplemented from G4VSolid.

Definition at line 2775 of file G4Sphere.cc.

2776{
2777 if (fCubicVolume == 0.)
2778 {
2779 G4double RRR = fRmax*fRmax*fRmax;
2780 G4double rrr = fRmin*fRmin*fRmin;
2781 fCubicVolume = fDPhi*(cosSTheta - cosETheta)*(RRR - rrr)/3.;
2782 }
2783 return fCubicVolume;
2784}
G4double fCubicVolume
Definition: G4CSGSolid.hh:70

References cosETheta, cosSTheta, G4CSGSolid::fCubicVolume, fDPhi, fRmax, and fRmin.

◆ GetDeltaPhiAngle()

G4double G4Sphere::GetDeltaPhiAngle ( ) const
inline

◆ GetDeltaThetaAngle()

G4double G4Sphere::GetDeltaThetaAngle ( ) const
inline

◆ 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 G4Sphere::GetEntityType ( ) const
virtual

Implements G4VSolid.

Definition at line 2733 of file G4Sphere.cc.

2734{
2735 return G4String("G4Sphere");
2736}

◆ GetExtent()

G4VisExtent G4Sphere::GetExtent ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 2859 of file G4Sphere.cc.

2860{
2861 return G4VisExtent(-fRmax, fRmax,-fRmax, fRmax,-fRmax, fRmax );
2862}

References fRmax.

◆ GetInnerRadius()

G4double G4Sphere::GetInnerRadius ( ) const
inline

◆ 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(), G4IntersectionSolid::BoundingLimits(), G4ScaledSolid::BoundingLimits(), G4SubtractionSolid::BoundingLimits(), G4UnionSolid::BoundingLimits(), G4Box::BoundingLimits(), G4Cons::BoundingLimits(), G4CutTubs::BoundingLimits(), G4Orb::BoundingLimits(), G4Para::BoundingLimits(), 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(), 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(), 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().

◆ GetOuterRadius()

G4double G4Sphere::GetOuterRadius ( ) const
inline

◆ GetPointOnSurface()

G4ThreeVector G4Sphere::GetPointOnSurface ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 2808 of file G4Sphere.cc.

2809{
2810 G4double RR = fRmax*fRmax;
2811 G4double rr = fRmin*fRmin;
2812
2813 // Find surface areas
2814 //
2815 G4double aInner = fDPhi*rr*(cosSTheta - cosETheta);
2816 G4double aOuter = fDPhi*RR*(cosSTheta - cosETheta);
2817 G4double aPhi = (!fFullPhiSphere) ? fDTheta*(RR - rr) : 0.;
2818 G4double aSTheta = (fSTheta > 0) ? 0.5*fDPhi*(RR - rr)*sinSTheta : 0.;
2819 G4double aETheta = (eTheta < pi) ? 0.5*fDPhi*(RR - rr)*sinETheta : 0.;
2820 G4double aTotal = aInner + aOuter + aPhi + aSTheta + aETheta;
2821
2822 // Select surface and generate a point
2823 //
2824 G4double select = aTotal*G4QuickRand();
2825 G4double u = G4QuickRand();
2826 G4double v = G4QuickRand();
2827 if (select < aInner + aOuter) // lateral surface
2828 {
2829 G4double r = (select < aInner) ? fRmin : fRmax;
2831 G4double rho = std::sqrt(1. - z*z);
2832 G4double phi = fDPhi*v + fSPhi;
2833 return G4ThreeVector(r*rho*std::cos(phi), r*rho*std::sin(phi), r*z);
2834 }
2835 else if (select < aInner + aOuter + aPhi) // cut in phi
2836 {
2837 G4double phi = (select < aInner + aOuter + 0.5*aPhi) ? fSPhi : fSPhi + fDPhi;
2838 G4double r = std::sqrt((RR - rr)*u + rr);
2839 G4double theta = fDTheta*v + fSTheta;
2840 G4double z = std::cos(theta);
2841 G4double rho = std::sin(theta);
2842 return G4ThreeVector(r*rho*std::cos(phi), r*rho*std::sin(phi), r*z);
2843 }
2844 else // cut in theta
2845 {
2846 G4double theta = (select < aTotal - aETheta) ? fSTheta : fSTheta + fDTheta;
2847 G4double r = std::sqrt((RR - rr)*u + rr);
2848 G4double phi = fDPhi*v + fSPhi;
2849 G4double z = std::cos(theta);
2850 G4double rho = std::sin(theta);
2851 return G4ThreeVector(r*rho*std::cos(phi), r*rho*std::sin(phi), r*z);
2852 }
2853}

References cosETheta, cosSTheta, eTheta, fDPhi, fDTheta, fFullPhiSphere, fRmax, fRmin, fSPhi, fSTheta, G4QuickRand(), pi, sinETheta, and sinSTheta.

◆ GetPolyhedron()

G4Polyhedron * G4CSGSolid::GetPolyhedron ( ) const
virtualinherited

Reimplemented from G4VSolid.

Definition at line 129 of file G4CSGSolid.cc.

130{
131 if (fpPolyhedron == nullptr ||
135 {
137 delete fpPolyhedron;
139 fRebuildPolyhedron = false;
140 l.unlock();
141 }
142 return fpPolyhedron;
143}
G4Polyhedron * fpPolyhedron
Definition: G4CSGSolid.hh:73
G4bool fRebuildPolyhedron
Definition: G4CSGSolid.hh:72
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
virtual G4Polyhedron * CreatePolyhedron() const
Definition: G4VSolid.cc:700
static G4int GetNumberOfRotationSteps()

References G4VSolid::CreatePolyhedron(), G4CSGSolid::fpPolyhedron, G4CSGSolid::fRebuildPolyhedron, HepPolyhedron::GetNumberOfRotationSteps(), G4Polyhedron::GetNumberOfRotationStepsAtTimeOfCreation(), anonymous_namespace{G4CSGSolid.cc}::polyhedronMutex, and G4TemplateAutoLock< _Mutex_t >::unlock().

Referenced by G4ScoringBox::Draw(), G4ScoringCylinder::Draw(), G4ScoringBox::DrawColumn(), and G4ScoringCylinder::DrawColumn().

◆ GetRadiusInRing()

G4double G4CSGSolid::GetRadiusInRing ( G4double  rmin,
G4double  rmax 
) const
protectedinherited

Definition at line 109 of file G4CSGSolid.cc.

110{
111 G4double k = G4QuickRand();
112 return (rmin <= 0) ? rmax*std::sqrt(k)
113 : std::sqrt(k*rmax*rmax + (1. - k)*rmin*rmin);
114}

References G4QuickRand().

Referenced by G4Cons::GetPointOnSurface(), and G4Torus::GetPointOnSurface().

◆ GetSinEndPhi()

G4double G4Sphere::GetSinEndPhi ( ) const
inline

Referenced by BoundingLimits().

◆ GetSinEndTheta()

G4double G4Sphere::GetSinEndTheta ( ) const
inline

Referenced by BoundingLimits().

◆ GetSinStartPhi()

G4double G4Sphere::GetSinStartPhi ( ) const
inline

Referenced by BoundingLimits().

◆ GetSinStartTheta()

G4double G4Sphere::GetSinStartTheta ( ) const
inline

Referenced by BoundingLimits().

◆ GetStartPhiAngle()

G4double G4Sphere::GetStartPhiAngle ( ) const
inline

◆ GetStartThetaAngle()

G4double G4Sphere::GetStartThetaAngle ( ) const
inline

◆ GetSurfaceArea()

G4double G4Sphere::GetSurfaceArea ( )
virtual

Reimplemented from G4VSolid.

Definition at line 2790 of file G4Sphere.cc.

2791{
2792 if (fSurfaceArea == 0.)
2793 {
2794 G4double RR = fRmax*fRmax;
2795 G4double rr = fRmin*fRmin;
2796 fSurfaceArea = fDPhi*(RR + rr)*(cosSTheta - cosETheta);
2797 if (!fFullPhiSphere) fSurfaceArea += fDTheta*(RR - rr);
2798 if (fSTheta > 0) fSurfaceArea += 0.5*fDPhi*(RR - rr)*sinSTheta;
2799 if (eTheta < CLHEP::pi) fSurfaceArea += 0.5*fDPhi*(RR - rr)*sinETheta;
2800 }
2801 return fSurfaceArea;
2802}
G4double fSurfaceArea
Definition: G4CSGSolid.hh:71
static constexpr double pi
Definition: SystemOfUnits.h:55

References cosETheta, cosSTheta, eTheta, fDPhi, fDTheta, fFullPhiSphere, fRmax, fRmin, fSTheta, G4CSGSolid::fSurfaceArea, CLHEP::pi, sinETheta, and sinSTheta.

◆ GetTolerance()

G4double G4VSolid::GetTolerance ( ) const
inlineinherited

◆ Initialize()

void G4Sphere::Initialize ( )
inlineprivate

◆ InitializePhiTrigonometry()

void G4Sphere::InitializePhiTrigonometry ( )
inlineprivate

◆ InitializeThetaTrigonometry()

void G4Sphere::InitializeThetaTrigonometry ( )
inlineprivate

◆ Inside()

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

Implements G4VSolid.

Definition at line 288 of file G4Sphere.cc.

289{
290 G4double rho,rho2,rad2,tolRMin,tolRMax;
291 G4double pPhi,pTheta;
292 EInside in = kOutside;
293
294 const G4double halfRmaxTolerance = fRmaxTolerance*0.5;
295 const G4double halfRminTolerance = fRminTolerance*0.5;
296 const G4double Rmax_minus = fRmax - halfRmaxTolerance;
297 const G4double Rmin_plus = (fRmin > 0) ? fRmin+halfRminTolerance : 0;
298
299 rho2 = p.x()*p.x() + p.y()*p.y() ;
300 rad2 = rho2 + p.z()*p.z() ;
301
302 // Check radial surfaces. Sets 'in'
303
304 tolRMin = Rmin_plus;
305 tolRMax = Rmax_minus;
306
307 if(rad2 == 0.0)
308 {
309 if (fRmin > 0.0)
310 {
311 return in = kOutside;
312 }
313 if ( (!fFullPhiSphere) || (!fFullThetaSphere) )
314 {
315 return in = kSurface;
316 }
317 else
318 {
319 return in = kInside;
320 }
321 }
322
323 if ( (rad2 <= Rmax_minus*Rmax_minus) && (rad2 >= Rmin_plus*Rmin_plus) )
324 {
325 in = kInside;
326 }
327 else
328 {
329 tolRMax = fRmax + halfRmaxTolerance; // outside case
330 tolRMin = std::max(fRmin-halfRminTolerance, 0.); // outside case
331 if ( (rad2 <= tolRMax*tolRMax) && (rad2 >= tolRMin*tolRMin) )
332 {
333 in = kSurface;
334 }
335 else
336 {
337 return in = kOutside;
338 }
339 }
340
341 // Phi boundaries : Do not check if it has no phi boundary!
342
343 if ( !fFullPhiSphere && rho2 ) // [fDPhi < twopi] and [p.x or p.y]
344 {
345 pPhi = std::atan2(p.y(),p.x()) ;
346
347 if ( pPhi < fSPhi - halfAngTolerance ) { pPhi += twopi; }
348 else if ( pPhi > ePhi + halfAngTolerance ) { pPhi -= twopi; }
349
350 if ( (pPhi < fSPhi - halfAngTolerance)
351 || (pPhi > ePhi + halfAngTolerance) ) { return in = kOutside; }
352
353 else if (in == kInside) // else it's kSurface anyway already
354 {
355 if ( (pPhi < fSPhi + halfAngTolerance)
356 || (pPhi > ePhi - halfAngTolerance) ) { in = kSurface; }
357 }
358 }
359
360 // Theta bondaries
361
362 if ( (rho2 || p.z()) && (!fFullThetaSphere) )
363 {
364 rho = std::sqrt(rho2);
365 pTheta = std::atan2(rho,p.z());
366
367 if ( in == kInside )
368 {
369 if ( ((fSTheta > 0.0) && (pTheta < fSTheta + halfAngTolerance))
370 || ((eTheta < pi) && (pTheta > eTheta - halfAngTolerance)) )
371 {
372 if ( (( (fSTheta>0.0)&&(pTheta>=fSTheta-halfAngTolerance) )
373 || (fSTheta == 0.0) )
374 && ((eTheta==pi)||(pTheta <= eTheta + halfAngTolerance) ) )
375 {
376 in = kSurface;
377 }
378 else
379 {
380 in = kOutside;
381 }
382 }
383 }
384 else
385 {
386 if ( ((fSTheta > 0.0)&&(pTheta < fSTheta - halfAngTolerance))
387 ||((eTheta < pi )&&(pTheta > eTheta + halfAngTolerance)) )
388 {
389 in = kOutside;
390 }
391 }
392 }
393 return in;
394}
@ kSurface
Definition: geomdefs.hh:69

References ePhi, eTheta, fFullPhiSphere, fFullThetaSphere, fRmax, fRmaxTolerance, fRmin, fRminTolerance, fSPhi, fSTheta, halfAngTolerance, kInside, kOutside, kSurface, G4INCL::Math::max(), pi, twopi, CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Referenced by DistanceToOut().

◆ operator=()

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

Definition at line 160 of file G4Sphere.cc.

161{
162 // Check assignment to self
163 //
164 if (this == &rhs) { return *this; }
165
166 // Copy base class data
167 //
169
170 // Copy data
171 //
174 fEpsilon = rhs.fEpsilon; fRmin = rhs.fRmin; fRmax = rhs.fRmax;
175 fSPhi = rhs.fSPhi; fDPhi = rhs.fDPhi; fSTheta = rhs.fSTheta;
176 fDTheta = rhs.fDTheta; sinCPhi = rhs.sinCPhi; cosCPhi = rhs.cosCPhi;
177 cosHDPhi = rhs.cosHDPhi;
179 sinSPhi = rhs.sinSPhi; cosSPhi = rhs.cosSPhi;
180 sinEPhi = rhs.sinEPhi; cosEPhi = rhs.cosEPhi;
181 hDPhi = rhs.hDPhi; cPhi = rhs.cPhi; ePhi = rhs.ePhi;
190
191 return *this;
192}
G4CSGSolid & operator=(const G4CSGSolid &rhs)
Definition: G4CSGSolid.cc:89

References cosCPhi, cosEPhi, cosETheta, cosHDPhi, cosHDPhiIT, cosHDPhiOT, cosSPhi, cosSTheta, cPhi, ePhi, eTheta, fDPhi, fDTheta, fEpsilon, fFullPhiSphere, fFullSphere, fFullThetaSphere, fRmax, fRmaxTolerance, fRmin, fRminTolerance, fSPhi, fSTheta, halfAngTolerance, halfCarTolerance, hDPhi, kAngTolerance, kRadTolerance, G4CSGSolid::operator=(), sinCPhi, sinEPhi, sinETheta, sinSPhi, sinSTheta, tanETheta, tanETheta2, tanSTheta, and tanSTheta2.

◆ operator==()

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

◆ SetDeltaPhiAngle()

void G4Sphere::SetDeltaPhiAngle ( G4double  newDphi)
inline

◆ SetDeltaThetaAngle()

void G4Sphere::SetDeltaThetaAngle ( G4double  newDTheta)
inline

◆ SetInnerRadius()

void G4Sphere::SetInnerRadius ( G4double  newRMin)
inline

◆ 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().

◆ SetOuterRadius()

void G4Sphere::SetOuterRadius ( G4double  newRmax)
inline

◆ SetStartPhiAngle()

void G4Sphere::SetStartPhiAngle ( G4double  newSphi,
G4bool  trig = true 
)
inline

◆ SetStartThetaAngle()

void G4Sphere::SetStartThetaAngle ( G4double  newSTheta)
inline

◆ StreamInfo()

std::ostream & G4Sphere::StreamInfo ( std::ostream &  os) const
virtual

Reimplemented from G4CSGSolid.

Definition at line 2751 of file G4Sphere.cc.

2752{
2753 G4int oldprc = os.precision(16);
2754 os << "-----------------------------------------------------------\n"
2755 << " *** Dump for solid - " << GetName() << " ***\n"
2756 << " ===================================================\n"
2757 << " Solid type: G4Sphere\n"
2758 << " Parameters: \n"
2759 << " inner radius: " << fRmin/mm << " mm \n"
2760 << " outer radius: " << fRmax/mm << " mm \n"
2761 << " starting phi of segment : " << fSPhi/degree << " degrees \n"
2762 << " delta phi of segment : " << fDPhi/degree << " degrees \n"
2763 << " starting theta of segment: " << fSTheta/degree << " degrees \n"
2764 << " delta theta of segment : " << fDTheta/degree << " degrees \n"
2765 << "-----------------------------------------------------------\n";
2766 os.precision(oldprc);
2767
2768 return os;
2769}
static constexpr double degree
Definition: G4SIunits.hh:124

References degree, fDPhi, fDTheta, fRmax, fRmin, fSPhi, fSTheta, G4VSolid::GetName(), and mm.

◆ SurfaceNormal()

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

Implements G4VSolid.

Definition at line 402 of file G4Sphere.cc.

403{
404 G4int noSurfaces = 0;
405 G4double rho, rho2, radius, pTheta, pPhi=0.;
406 G4double distRMin = kInfinity;
407 G4double distSPhi = kInfinity, distEPhi = kInfinity;
408 G4double distSTheta = kInfinity, distETheta = kInfinity;
409 G4ThreeVector nR, nPs, nPe, nTs, nTe, nZ(0.,0.,1.);
410 G4ThreeVector norm, sumnorm(0.,0.,0.);
411
412 rho2 = p.x()*p.x()+p.y()*p.y();
413 radius = std::sqrt(rho2+p.z()*p.z());
414 rho = std::sqrt(rho2);
415
416 G4double distRMax = std::fabs(radius-fRmax);
417 if (fRmin) distRMin = std::fabs(radius-fRmin);
418
419 if ( rho && !fFullSphere )
420 {
421 pPhi = std::atan2(p.y(),p.x());
422
423 if (pPhi < fSPhi-halfAngTolerance) { pPhi += twopi; }
424 else if (pPhi > ePhi+halfAngTolerance) { pPhi -= twopi; }
425 }
426 if ( !fFullPhiSphere )
427 {
428 if ( rho )
429 {
430 distSPhi = std::fabs( pPhi-fSPhi );
431 distEPhi = std::fabs( pPhi-ePhi );
432 }
433 else if( !fRmin )
434 {
435 distSPhi = 0.;
436 distEPhi = 0.;
437 }
438 nPs = G4ThreeVector(sinSPhi,-cosSPhi,0);
439 nPe = G4ThreeVector(-sinEPhi,cosEPhi,0);
440 }
441 if ( !fFullThetaSphere )
442 {
443 if ( rho )
444 {
445 pTheta = std::atan2(rho,p.z());
446 distSTheta = std::fabs(pTheta-fSTheta);
447 distETheta = std::fabs(pTheta-eTheta);
448
449 nTs = G4ThreeVector(-cosSTheta*p.x()/rho,
450 -cosSTheta*p.y()/rho,
451 sinSTheta );
452
453 nTe = G4ThreeVector( cosETheta*p.x()/rho,
454 cosETheta*p.y()/rho,
455 -sinETheta );
456 }
457 else if( !fRmin )
458 {
459 if ( fSTheta )
460 {
461 distSTheta = 0.;
462 nTs = G4ThreeVector(0.,0.,-1.);
463 }
464 if ( eTheta < pi )
465 {
466 distETheta = 0.;
467 nTe = G4ThreeVector(0.,0.,1.);
468 }
469 }
470 }
471 if( radius ) { nR = G4ThreeVector(p.x()/radius,p.y()/radius,p.z()/radius); }
472
473 if( distRMax <= halfCarTolerance )
474 {
475 ++noSurfaces;
476 sumnorm += nR;
477 }
478 if( fRmin && (distRMin <= halfCarTolerance) )
479 {
480 ++noSurfaces;
481 sumnorm -= nR;
482 }
483 if( !fFullPhiSphere )
484 {
485 if (distSPhi <= halfAngTolerance)
486 {
487 ++noSurfaces;
488 sumnorm += nPs;
489 }
490 if (distEPhi <= halfAngTolerance)
491 {
492 ++noSurfaces;
493 sumnorm += nPe;
494 }
495 }
496 if ( !fFullThetaSphere )
497 {
498 if ((distSTheta <= halfAngTolerance) && (fSTheta > 0.))
499 {
500 ++noSurfaces;
501 if ((radius <= halfCarTolerance) && fFullPhiSphere) { sumnorm += nZ; }
502 else { sumnorm += nTs; }
503 }
504 if ((distETheta <= halfAngTolerance) && (eTheta < pi))
505 {
506 ++noSurfaces;
507 if ((radius <= halfCarTolerance) && fFullPhiSphere) { sumnorm -= nZ; }
508 else { sumnorm += nTe; }
509 if(sumnorm.z() == 0.) { sumnorm += nZ; }
510 }
511 }
512 if ( noSurfaces == 0 )
513 {
514#ifdef G4CSGDEBUG
515 G4Exception("G4Sphere::SurfaceNormal(p)", "GeomSolids1002",
516 JustWarning, "Point p is not on surface !?" );
517#endif
518 norm = ApproxSurfaceNormal(p);
519 }
520 else if ( noSurfaces == 1 ) { norm = sumnorm; }
521 else { norm = sumnorm.unit(); }
522 return norm;
523}
Hep3Vector unit() const
G4ThreeVector ApproxSurfaceNormal(const G4ThreeVector &p) const
Definition: G4Sphere.cc:531

References ApproxSurfaceNormal(), cosEPhi, cosETheta, cosSPhi, cosSTheta, ePhi, eTheta, fFullPhiSphere, fFullSphere, fFullThetaSphere, fRmax, fRmin, fSPhi, fSTheta, G4Exception(), halfAngTolerance, halfCarTolerance, JustWarning, kInfinity, pi, sinEPhi, sinETheta, sinSPhi, sinSTheta, twopi, CLHEP::Hep3Vector::unit(), CLHEP::Hep3Vector::x(), CLHEP::Hep3Vector::y(), and CLHEP::Hep3Vector::z().

Field Documentation

◆ cosCPhi

G4double G4Sphere::cosCPhi
private

Definition at line 224 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), and operator=().

◆ cosEPhi

G4double G4Sphere::cosEPhi
private

◆ cosETheta

G4double G4Sphere::cosETheta
private

◆ cosHDPhi

G4double G4Sphere::cosHDPhi
private

Definition at line 224 of file G4Sphere.hh.

Referenced by DistanceToIn(), and operator=().

◆ cosHDPhiIT

G4double G4Sphere::cosHDPhiIT
private

Definition at line 224 of file G4Sphere.hh.

Referenced by DistanceToIn(), and operator=().

◆ cosHDPhiOT

G4double G4Sphere::cosHDPhiOT
private

Definition at line 224 of file G4Sphere.hh.

Referenced by DistanceToIn(), and operator=().

◆ cosSPhi

G4double G4Sphere::cosSPhi
private

◆ cosSTheta

G4double G4Sphere::cosSTheta
private

◆ cPhi

G4double G4Sphere::cPhi
private

Definition at line 225 of file G4Sphere.hh.

Referenced by operator=().

◆ ePhi

G4double G4Sphere::ePhi
private

Definition at line 225 of file G4Sphere.hh.

Referenced by DistanceToOut(), Inside(), operator=(), and SurfaceNormal().

◆ eTheta

G4double G4Sphere::eTheta
private

◆ fCubicVolume

G4double G4CSGSolid::fCubicVolume = 0.0
protectedinherited

◆ fDPhi

G4double G4Sphere::fDPhi
private

◆ fDTheta

G4double G4Sphere::fDTheta
private

◆ fEpsilon

G4double G4Sphere::fEpsilon = 2.e-11
private

Definition at line 216 of file G4Sphere.hh.

Referenced by G4Sphere(), and operator=().

◆ fFullPhiSphere

G4bool G4Sphere::fFullPhiSphere =false
private

◆ fFullSphere

G4bool G4Sphere::fFullSphere =true
private

Definition at line 234 of file G4Sphere.hh.

Referenced by operator=(), and SurfaceNormal().

◆ fFullThetaSphere

G4bool G4Sphere::fFullThetaSphere =false
private

◆ fpPolyhedron

G4Polyhedron* G4CSGSolid::fpPolyhedron = nullptr
mutableprotectedinherited

◆ fRebuildPolyhedron

G4bool G4CSGSolid::fRebuildPolyhedron = false
mutableprotectedinherited

◆ fRmax

G4double G4Sphere::fRmax
private

◆ fRmaxTolerance

G4double G4Sphere::fRmaxTolerance
private

Definition at line 215 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), G4Sphere(), Inside(), and operator=().

◆ fRmin

G4double G4Sphere::fRmin
private

◆ fRminTolerance

G4double G4Sphere::fRminTolerance
private

Definition at line 215 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), G4Sphere(), Inside(), and operator=().

◆ fshapeName

G4String G4VSolid::fshapeName
privateinherited

Definition at line 312 of file G4VSolid.hh.

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

◆ fSPhi

G4double G4Sphere::fSPhi
private

◆ fSTheta

G4double G4Sphere::fSTheta
private

◆ fSurfaceArea

G4double G4CSGSolid::fSurfaceArea = 0.0
protectedinherited

◆ halfAngTolerance

G4double G4Sphere::halfAngTolerance
private

Definition at line 238 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), G4Sphere(), Inside(), operator=(), and SurfaceNormal().

◆ halfCarTolerance

G4double G4Sphere::halfCarTolerance
private

Definition at line 238 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), G4Sphere(), operator=(), and SurfaceNormal().

◆ hDPhi

G4double G4Sphere::hDPhi
private

Definition at line 225 of file G4Sphere.hh.

Referenced by operator=().

◆ kAngTolerance

G4double G4Sphere::kAngTolerance
private

Definition at line 215 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), G4Sphere(), and 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(), DistanceToIn(), G4Ellipsoid::DistanceToIn(), G4Hype::DistanceToIn(), G4Paraboloid::DistanceToIn(), G4VCSGfaceted::DistanceToIn(), G4TessellatedSolid::DistanceToInCore(), G4Cons::DistanceToOut(), G4CutTubs::DistanceToOut(), 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(), 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().

◆ kRadTolerance

G4double G4Sphere::kRadTolerance
private

Definition at line 216 of file G4Sphere.hh.

Referenced by G4Sphere(), and operator=().

◆ sinCPhi

G4double G4Sphere::sinCPhi
private

Definition at line 224 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), and operator=().

◆ sinEPhi

G4double G4Sphere::sinEPhi
private

◆ sinETheta

G4double G4Sphere::sinETheta
private

◆ sinSPhi

G4double G4Sphere::sinSPhi
private

◆ sinSTheta

G4double G4Sphere::sinSTheta
private

◆ tanETheta

G4double G4Sphere::tanETheta
private

Definition at line 230 of file G4Sphere.hh.

Referenced by DistanceToOut(), and operator=().

◆ tanETheta2

G4double G4Sphere::tanETheta2
private

Definition at line 230 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), and operator=().

◆ tanSTheta

G4double G4Sphere::tanSTheta
private

Definition at line 230 of file G4Sphere.hh.

Referenced by DistanceToOut(), and operator=().

◆ tanSTheta2

G4double G4Sphere::tanSTheta2
private

Definition at line 230 of file G4Sphere.hh.

Referenced by DistanceToIn(), DistanceToOut(), and operator=().


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