Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | List of all members
G4Torus Class Reference

#include <G4Torus.hh>

Inheritance diagram for G4Torus:
Inheritance graph
[legend]
Collaboration diagram for G4Torus:
Collaboration graph
[legend]

Public Member Functions

 G4Torus (const G4String &pName, G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
 
 ~G4Torus ()
 
G4double GetRmin () const
 
G4double GetRmax () const
 
G4double GetRtor () const
 
G4double GetSPhi () const
 
G4double GetDPhi () const
 
G4double GetSinStartPhi () const
 
G4double GetCosStartPhi () const
 
G4double GetSinEndPhi () const
 
G4double GetCosEndPhi () const
 
G4double GetCubicVolume ()
 
G4double GetSurfaceArea ()
 
EInside Inside (const G4ThreeVector &p) const
 
void Extent (G4ThreeVector &pMin, G4ThreeVector &pMax) const
 
G4bool CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pmin, G4double &pmax) const
 
void ComputeDimensions (G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
 
G4ThreeVector SurfaceNormal (const G4ThreeVector &p) const
 
G4double DistanceToIn (const G4ThreeVector &p, const G4ThreeVector &v) const
 
G4double DistanceToIn (const G4ThreeVector &p) const
 
G4double DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=G4bool(false), G4bool *validNorm=0, G4ThreeVector *n=0) const
 
G4double DistanceToOut (const G4ThreeVector &p) const
 
G4GeometryType GetEntityType () const
 
G4ThreeVector GetPointOnSurface () const
 
G4VSolid * Clone () const
 
std::ostream & StreamInfo (std::ostream &os) const
 
void DescribeYourselfTo (G4VGraphicsScene &scene) const
 
G4Polyhedron * CreatePolyhedron () const
 
void SetAllParameters (G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
 
 G4Torus (__void__ &)
 
 G4Torus (const G4Torus &rhs)
 
G4Torus & operator= (const G4Torus &rhs)
 
- Public Member Functions inherited from G4CSGSolid
 G4CSGSolid (const G4String &pName)
 
virtual ~G4CSGSolid ()
 
virtual G4Polyhedron * GetPolyhedron () const
 
 G4CSGSolid (__void__ &)
 
 G4CSGSolid (const G4CSGSolid &rhs)
 
G4CSGSolid & operator= (const G4CSGSolid &rhs)
 
- Public Member Functions inherited from G4VSolid
 G4VSolid (const G4String &name)
 
virtual ~G4VSolid ()
 
G4bool operator== (const G4VSolid &s) const
 
G4String GetName () const
 
void SetName (const G4String &name)
 
G4double GetTolerance () const
 
void DumpInfo () const
 
virtual G4VisExtent GetExtent () const
 
virtual const G4VSolid * GetConstituentSolid (G4int no) const
 
virtual G4VSolid * GetConstituentSolid (G4int no)
 
virtual const G4DisplacedSolid * GetDisplacedSolidPtr () const
 
virtual G4DisplacedSolid * GetDisplacedSolidPtr ()
 
 G4VSolid (__void__ &)
 
 G4VSolid (const G4VSolid &rhs)
 
G4VSolid & operator= (const G4VSolid &rhs)
 
G4double EstimateCubicVolume (G4int nStat, G4double epsilon) const
 
G4double EstimateSurfaceArea (G4int nStat, G4double ell) const
 

Additional Inherited Members

- Protected Member Functions inherited from G4CSGSolid
G4double GetRadiusInRing (G4double rmin, G4double rmax) const
 
- Protected Member Functions inherited from G4VSolid
void CalculateClippedPolygonExtent (G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
 
void ClipCrossSection (G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
 
void ClipBetweenSections (G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
 
void ClipPolygon (G4ThreeVectorList &pPolygon, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis) const
 
- Protected Attributes inherited from G4CSGSolid
G4double fCubicVolume
 
G4double fSurfaceArea
 
G4bool fRebuildPolyhedron
 
G4Polyhedron * fpPolyhedron
 
- Protected Attributes inherited from G4VSolid
G4double kCarTolerance
 

Detailed Description

Definition at line 102 of file G4Torus.hh.

Constructor & Destructor Documentation

G4Torus::G4Torus ( const G4String &  pName,
G4double  pRmin,
G4double  pRmax,
G4double  pRtor,
G4double  pSPhi,
G4double  pDPhi 
)

Definition at line 87 of file G4Torus.cc.

93  : G4CSGSolid(pName)
94 {
95  SetAllParameters(pRmin, pRmax, pRtor, pSPhi, pDPhi);
96 }
G4CSGSolid::G4CSGSolid
G4CSGSolid(const G4String &pName)
Definition: G4CSGSolid.cc:49
G4Torus::SetAllParameters
void SetAllParameters(G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
Definition: G4Torus.cc:103

Here is the call graph for this function:

Here is the caller graph for this function:

G4Torus::~G4Torus ( )

Definition at line 200 of file G4Torus.cc.

201 {}
G4Torus::G4Torus ( __void__ &  a)

Definition at line 188 of file G4Torus.cc.

189  : G4CSGSolid(a), fRmin(0.), fRmax(0.), fRtor(0.), fSPhi(0.),
190  fDPhi(0.), fRminTolerance(0.), fRmaxTolerance(0. ),
191  kRadTolerance(0.), kAngTolerance(0.),
192  halfCarTolerance(0.), halfAngTolerance(0.)
193 {
194 }
G4CSGSolid::G4CSGSolid
G4CSGSolid(const G4String &pName)
Definition: G4CSGSolid.cc:49
G4Torus::G4Torus ( const G4Torus &  rhs)

Definition at line 207 of file G4Torus.cc.

208  : G4CSGSolid(rhs), fRmin(rhs.fRmin),fRmax(rhs.fRmax),
209  fRtor(rhs.fRtor),fSPhi(rhs.fSPhi),fDPhi(rhs.fDPhi),
210  fRminTolerance(rhs.fRminTolerance), fRmaxTolerance(rhs.fRmaxTolerance),
211  kRadTolerance(rhs.kRadTolerance), kAngTolerance(rhs.kAngTolerance),
212  halfCarTolerance(rhs.halfCarTolerance),
213  halfAngTolerance(rhs.halfAngTolerance)
214 {
215 }
G4CSGSolid::G4CSGSolid
G4CSGSolid(const G4String &pName)
Definition: G4CSGSolid.cc:49

Member Function Documentation

G4bool G4Torus::CalculateExtent ( const EAxis  pAxis,
const G4VoxelLimits &  pVoxelLimit,
const G4AffineTransform &  pTransform,
G4double &  pmin,
G4double &  pmax 
) const
virtual

Implements G4VSolid.

Definition at line 461 of file G4Torus.cc.

465 {
466  G4ThreeVector bmin, bmax;
467  G4bool exist;
468 
469  // Get bounding box
470  Extent(bmin,bmax);
471 
472  // Check bounding box
473  G4BoundingEnvelope bbox(bmin,bmax);
474 #ifdef G4BBOX_EXTENT
475  if (true) return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
476 #endif
477  if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax))
478  {
479  return exist = (pMin < pMax) ? true : false;
480  }
481 
482  // Get parameters of the solid
483  G4double rmin = GetRmin();
484  G4double rmax = GetRmax();
485  G4double rtor = GetRtor();
486  G4double dphi = GetDPhi();
487  G4double sinStart = GetSinStartPhi();
488  G4double cosStart = GetCosStartPhi();
489  G4double sinEnd = GetSinEndPhi();
490  G4double cosEnd = GetCosEndPhi();
491  G4double rint = rtor - rmax;
492  G4double rext = rtor + rmax;
493 
494  // Find bounding envelope and calculate extent
495  //
496  static const G4int NPHI = 24; // number of steps for whole torus
497  static const G4int NDISK = 16; // number of steps for disk
498  static const G4double sinHalfDisk = std::sin(pi/NDISK);
499  static const G4double cosHalfDisk = std::cos(pi/NDISK);
500  static const G4double sinStepDisk = 2.*sinHalfDisk*cosHalfDisk;
501  static const G4double cosStepDisk = 1. - 2.*sinHalfDisk*sinHalfDisk;
502 
503  G4double astep = (360/NPHI)*deg; // max angle for one slice in phi
504  G4int kphi = (dphi <= astep) ? 1 : (G4int)((dphi-deg)/astep) + 1;
505  G4double ang = dphi/kphi;
506 
507  G4double sinHalf = std::sin(0.5*ang);
508  G4double cosHalf = std::cos(0.5*ang);
509  G4double sinStep = 2.*sinHalf*cosHalf;
510  G4double cosStep = 1. - 2.*sinHalf*sinHalf;
511 
512  // define vectors for bounding envelope
513  G4ThreeVectorList pols[NDISK+1];
514  for (G4int k=0; k<NDISK+1; ++k) pols[k].resize(4);
515 
516  std::vector<const G4ThreeVectorList *> polygons;
517  polygons.resize(NDISK+1);
518  for (G4int k=0; k<NDISK+1; ++k) polygons[k] = &pols[k];
519 
520  // set internal and external reference circles
521  G4TwoVector rzmin[NDISK];
522  G4TwoVector rzmax[NDISK];
523 
524  if ((rtor-rmin*sinHalfDisk)/cosHalf > (rtor+rmin*sinHalfDisk)) rmin = 0;
525  rmax /= cosHalfDisk;
526  G4double sinCurDisk = sinHalfDisk;
527  G4double cosCurDisk = cosHalfDisk;
528  for (G4int k=0; k<NDISK; ++k)
529  {
530  G4double rmincur = rtor + rmin*cosCurDisk;
531  if (cosCurDisk < 0 && rmin > 0) rmincur /= cosHalf;
532  rzmin[k].set(rmincur,rmin*sinCurDisk);
533 
534  G4double rmaxcur = rtor + rmax*cosCurDisk;
535  if (cosCurDisk > 0) rmaxcur /= cosHalf;
536  rzmax[k].set(rmaxcur,rmax*sinCurDisk);
537 
538  G4double sinTmpDisk = sinCurDisk;
539  sinCurDisk = sinCurDisk*cosStepDisk + cosCurDisk*sinStepDisk;
540  cosCurDisk = cosCurDisk*cosStepDisk - sinTmpDisk*sinStepDisk;
541  }
542 
543  // Loop along slices in Phi. The extent is calculated as cumulative
544  // extent of the slices
545  pMin = kInfinity;
546  pMax = -kInfinity;
547  G4double eminlim = pVoxelLimit.GetMinExtent(pAxis);
548  G4double emaxlim = pVoxelLimit.GetMaxExtent(pAxis);
549  G4double sinCur1 = 0, cosCur1 = 0, sinCur2 = 0, cosCur2 = 0;
550  for (G4int i=0; i<kphi+1; ++i)
551  {
552  if (i == 0)
553  {
554  sinCur1 = sinStart;
555  cosCur1 = cosStart;
556  sinCur2 = sinCur1*cosHalf + cosCur1*sinHalf;
557  cosCur2 = cosCur1*cosHalf - sinCur1*sinHalf;
558  }
559  else
560  {
561  sinCur1 = sinCur2;
562  cosCur1 = cosCur2;
563  sinCur2 = (i == kphi) ? sinEnd : sinCur1*cosStep + cosCur1*sinStep;
564  cosCur2 = (i == kphi) ? cosEnd : cosCur1*cosStep - sinCur1*sinStep;
565  }
566  for (G4int k=0; k<NDISK; ++k)
567  {
568  G4double r1 = rzmin[k].x(), r2 = rzmax[k].x();
569  G4double z1 = rzmin[k].y(), z2 = rzmax[k].y();
570  pols[k][0].set(r1*cosCur1,r1*sinCur1,z1);
571  pols[k][1].set(r2*cosCur1,r2*sinCur1,z2);
572  pols[k][2].set(r2*cosCur2,r2*sinCur2,z2);
573  pols[k][3].set(r1*cosCur2,r1*sinCur2,z1);
574  }
575  pols[NDISK] = pols[0];
576 
577  // get bounding box of current slice
578  G4TwoVector vmin,vmax;
579  G4GeomTools::
580  DiskExtent(rint,rext,sinCur1,cosCur1,sinCur2,cosCur2,vmin,vmax);
581  bmin.setX(vmin.x()); bmin.setY(vmin.y());
582  bmax.setX(vmax.x()); bmax.setY(vmax.y());
583 
584  // set bounding envelope for current slice and adjust extent
585  G4double emin,emax;
586  G4BoundingEnvelope benv(bmin,bmax,polygons);
587  if (!benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,emin,emax)) continue;
588  if (emin < pMin) pMin = emin;
589  if (emax > pMax) pMax = emax;
590  if (eminlim > pMin && emaxlim < pMax) break; // max possible extent
591  }
592  return (pMin < pMax);
593 }
CLHEP::Hep2Vector::y
double y() const
CLHEP::Hep2Vector::x
double x() const
kInfinity
static const G4double kInfinity
Definition: geomdefs.hh:42
G4Torus::GetRmax
G4double GetRmax() const
G4Torus::GetSinStartPhi
G4double GetSinStartPhi() const
G4int
int G4int
Definition: G4Types.hh:78
CLHEP::Hep3Vector::setY
void setY(double)
CLHEP::Hep3Vector::setX
void setX(double)
G4Torus::GetRtor
G4double GetRtor() const
G4Torus::GetSinEndPhi
G4double GetSinEndPhi() const
G4Torus::Extent
void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Torus.cc:417
G4Torus::GetRmin
G4double GetRmin() const
CLHEP::Hep2Vector::set
void set(double x, double y)
G4Torus::GetCosEndPhi
G4double GetCosEndPhi() const
G4bool
bool G4bool
Definition: G4Types.hh:79
G4Torus::GetDPhi
G4double GetDPhi() const
G4ThreeVectorList
std::vector< G4ThreeVector > G4ThreeVectorList
Definition: G4BoundingEnvelope.hh:61
emax
static const G4double emax
Definition: G4ChatterjeeCrossSection.hh:39
G4Torus::GetCosStartPhi
G4double GetCosStartPhi() const
pi
static constexpr double pi
Definition: G4SIunits.hh:75
G4double
double G4double
Definition: G4Types.hh:76
deg
static constexpr double deg
Definition: G4SIunits.hh:152
G4VoxelLimits::GetMaxExtent
G4double GetMaxExtent(const EAxis pAxis) const
G4VoxelLimits::GetMinExtent
G4double GetMinExtent(const EAxis pAxis) const
G4GeomTools::DiskExtent
static G4bool DiskExtent(G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
Definition: G4GeomTools.cc:378

Here is the call graph for this function:

G4VSolid * G4Torus::Clone ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1556 of file G4Torus.cc.

1557 {
1558  return new G4Torus(*this);
1559 }
G4Torus::G4Torus
G4Torus(const G4String &pName, G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
Definition: G4Torus.cc:87

Here is the call graph for this function:

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

Reimplemented from G4VSolid.

Definition at line 248 of file G4Torus.cc.

251 {
252  p->ComputeDimensions(*this,n,pRep);
253 }
G4VPVParameterisation::ComputeDimensions
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
Definition: G4VPVParameterisation.hh:95

Here is the call graph for this function:

G4Polyhedron * G4Torus::CreatePolyhedron ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1641 of file G4Torus.cc.

1642 {
1643  return new G4PolyhedronTorus (fRmin, fRmax, fRtor, fSPhi, fDPhi);
1644 }
void G4Torus::DescribeYourselfTo ( G4VGraphicsScene &  scene) const
virtual

Implements G4VSolid.

Definition at line 1636 of file G4Torus.cc.

1637 {
1638  scene.AddSolid (*this);
1639 }
G4VGraphicsScene::AddSolid
virtual void AddSolid(const G4Box &)=0

Here is the call graph for this function:

G4double G4Torus::DistanceToIn ( const G4ThreeVector &  p,
const G4ThreeVector &  v 
) const
virtual

Implements G4VSolid.

Definition at line 947 of file G4Torus.cc.

949 {
950 
951  G4double snxt=kInfinity, sphi=kInfinity; // snxt = default return value
952 
953  G4double sd[4] ;
954 
955  // Precalculated trig for phi intersections - used by r,z intersections to
956  // check validity
957 
958  G4bool seg; // true if segmented
959  G4double hDPhi; // half dphi
960  G4double cPhi,sinCPhi=0.,cosCPhi=0.; // central phi
961 
962  G4double tolORMin2; // `generous' radii squared
963  G4double tolORMax2;
964 
965  G4double Dist,xi,yi,zi,rhoi,it2; // Intersection point variables
966 
967  G4double Comp;
968  G4double cosSPhi,sinSPhi; // Trig for phi start intersect
969  G4double ePhi,cosEPhi,sinEPhi; // for phi end intersect
970 
971  // Set phi divided flag and precalcs
972  //
973  if ( fDPhi < twopi )
974  {
975  seg = true ;
976  hDPhi = 0.5*fDPhi ; // half delta phi
977  cPhi = fSPhi + hDPhi ;
978  sinCPhi = std::sin(cPhi) ;
979  cosCPhi = std::cos(cPhi) ;
980  }
981  else
982  {
983  seg = false ;
984  }
985 
986  if (fRmin > fRminTolerance) // Calculate tolerant rmin and rmax
987  {
988  tolORMin2 = (fRmin - fRminTolerance)*(fRmin - fRminTolerance) ;
989  }
990  else
991  {
992  tolORMin2 = 0 ;
993  }
994  tolORMax2 = (fRmax + fRmaxTolerance)*(fRmax + fRmaxTolerance) ;
995 
996  // Intersection with Rmax (possible return) and Rmin (must also check phi)
997 
998  snxt = SolveNumericJT(p,v,fRmax,true);
999 
1000  if (fRmin) // Possible Rmin intersection
1001  {
1002  sd[0] = SolveNumericJT(p,v,fRmin,true);
1003  if ( sd[0] < snxt ) { snxt = sd[0] ; }
1004  }
1005 
1006  //
1007  // Phi segment intersection
1008  //
1009  // o Tolerant of points inside phi planes by up to kCarTolerance*0.5
1010  //
1011  // o NOTE: Large duplication of code between sphi & ephi checks
1012  // -> only diffs: sphi -> ephi, Comp -> -Comp and half-plane
1013  // intersection check <=0 -> >=0
1014  // -> use some form of loop Construct ?
1015 
1016  if (seg)
1017  {
1018  sinSPhi = std::sin(fSPhi) ; // First phi surface ('S'tarting phi)
1019  cosSPhi = std::cos(fSPhi) ;
1020  Comp = v.x()*sinSPhi - v.y()*cosSPhi ; // Component in outwards
1021  // normal direction
1022  if (Comp < 0 )
1023  {
1024  Dist = (p.y()*cosSPhi - p.x()*sinSPhi) ;
1025 
1026  if (Dist < halfCarTolerance)
1027  {
1028  sphi = Dist/Comp ;
1029  if (sphi < snxt)
1030  {
1031  if ( sphi < 0 ) { sphi = 0 ; }
1032 
1033  xi = p.x() + sphi*v.x() ;
1034  yi = p.y() + sphi*v.y() ;
1035  zi = p.z() + sphi*v.z() ;
1036  rhoi = std::hypot(xi,yi);
1037  it2 = zi*zi + (rhoi-fRtor)*(rhoi-fRtor);
1038 
1039  if ( it2 >= tolORMin2 && it2 <= tolORMax2 )
1040  {
1041  // r intersection is good - check intersecting
1042  // with correct half-plane
1043  //
1044  if ((yi*cosCPhi-xi*sinCPhi)<=0) { snxt=sphi; }
1045  }
1046  }
1047  }
1048  }
1049  ePhi=fSPhi+fDPhi; // Second phi surface ('E'nding phi)
1050  sinEPhi=std::sin(ePhi);
1051  cosEPhi=std::cos(ePhi);
1052  Comp=-(v.x()*sinEPhi-v.y()*cosEPhi);
1053 
1054  if ( Comp < 0 ) // Component in outwards normal dirn
1055  {
1056  Dist = -(p.y()*cosEPhi - p.x()*sinEPhi) ;
1057 
1058  if (Dist < halfCarTolerance )
1059  {
1060  sphi = Dist/Comp ;
1061 
1062  if (sphi < snxt )
1063  {
1064  if (sphi < 0 ) { sphi = 0 ; }
1065 
1066  xi = p.x() + sphi*v.x() ;
1067  yi = p.y() + sphi*v.y() ;
1068  zi = p.z() + sphi*v.z() ;
1069  rhoi = std::hypot(xi,yi);
1070  it2 = zi*zi + (rhoi-fRtor)*(rhoi-fRtor);
1071 
1072  if (it2 >= tolORMin2 && it2 <= tolORMax2)
1073  {
1074  // z and r intersections good - check intersecting
1075  // with correct half-plane
1076  //
1077  if ((yi*cosCPhi-xi*sinCPhi)>=0) { snxt=sphi; }
1078  }
1079  }
1080  }
1081  }
1082  }
1083  if(snxt < halfCarTolerance) { snxt = 0.0 ; }
1084 
1085  return snxt ;
1086 }
kInfinity
static const G4double kInfinity
Definition: geomdefs.hh:42
CLHEP::Hep3Vector::x
double x() const
CLHEP::Hep3Vector::z
double z() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4bool
bool G4bool
Definition: G4Types.hh:79
CLHEP::Hep3Vector::y
double y() const
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4Torus::DistanceToIn ( const G4ThreeVector &  p) const
virtual

Implements G4VSolid.

Definition at line 1095 of file G4Torus.cc.

1096 {
1097  G4double safe=0.0, safe1, safe2 ;
1098  G4double phiC, cosPhiC, sinPhiC, safePhi, ePhi, cosPsi ;
1099  G4double rho, pt ;
1100 
1101  rho = std::hypot(p.x(),p.y());
1102  pt = std::hypot(p.z(),rho-fRtor);
1103  safe1 = fRmin - pt ;
1104  safe2 = pt - fRmax ;
1105 
1106  if (safe1 > safe2) { safe = safe1; }
1107  else { safe = safe2; }
1108 
1109  if ( fDPhi < twopi && rho )
1110  {
1111  phiC = fSPhi + fDPhi*0.5 ;
1112  cosPhiC = std::cos(phiC) ;
1113  sinPhiC = std::sin(phiC) ;
1114  cosPsi = (p.x()*cosPhiC + p.y()*sinPhiC)/rho ;
1115 
1116  if (cosPsi < std::cos(fDPhi*0.5) ) // Psi=angle from central phi to point
1117  { // Point lies outside phi range
1118  if ((p.y()*cosPhiC - p.x()*sinPhiC) <= 0 )
1119  {
1120  safePhi = std::fabs(p.x()*std::sin(fSPhi) - p.y()*std::cos(fSPhi)) ;
1121  }
1122  else
1123  {
1124  ePhi = fSPhi + fDPhi ;
1125  safePhi = std::fabs(p.x()*std::sin(ePhi) - p.y()*std::cos(ePhi)) ;
1126  }
1127  if (safePhi > safe) { safe = safePhi ; }
1128  }
1129  }
1130  if (safe < 0 ) { safe = 0 ; }
1131  return safe;
1132 }
CLHEP::Hep3Vector::x
double x() const
CLHEP::Hep3Vector::z
double z() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
CLHEP::Hep3Vector::y
double y() const
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

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

Implements G4VSolid.

Definition at line 1140 of file G4Torus.cc.

1145 {
1146  ESide side = kNull, sidephi = kNull ;
1147  G4double snxt = kInfinity, sphi, sd[4] ;
1148 
1149  // Vars for phi intersection
1150  //
1151  G4double sinSPhi, cosSPhi, ePhi, sinEPhi, cosEPhi;
1152  G4double cPhi, sinCPhi, cosCPhi ;
1153  G4double pDistS, compS, pDistE, compE, sphi2, xi, yi, zi, vphi ;
1154 
1155  // Radial Intersections Defenitions & General Precals
1156 
1158 
1159 #if 1
1160 
1161  // This is the version with the calculation of CalcNorm = true
1162  // To be done: Check the precision of this calculation.
1163  // If you want return always validNorm = false, then take the version below
1164 
1165 
1166  G4double rho = std::hypot(p.x(),p.y());
1167  G4double pt = hypot(p.z(),rho-fRtor);
1168 
1169  G4double pDotV = p.x()*v.x() + p.y()*v.y() + p.z()*v.z() ;
1170 
1171  G4double tolRMax = fRmax - fRmaxTolerance ;
1172 
1173  G4double vDotNmax = pDotV - fRtor*(v.x()*p.x() + v.y()*p.y())/rho ;
1174  G4double pDotxyNmax = (1 - fRtor/rho) ;
1175 
1176  if( (pt*pt > tolRMax*tolRMax) && (vDotNmax >= 0) )
1177  {
1178  // On tolerant boundary & heading outwards (or perpendicular to) outer
1179  // radial surface -> leaving immediately with *n for really convex part
1180  // only
1181 
1182  if ( calcNorm && (pDotxyNmax >= -2.*fRmaxTolerance) )
1183  {
1184  *n = G4ThreeVector( p.x()*(1 - fRtor/rho)/pt,
1185  p.y()*(1 - fRtor/rho)/pt,
1186  p.z()/pt ) ;
1187  *validNorm = true ;
1188  }
1189 
1190  return snxt = 0 ; // Leaving by Rmax immediately
1191  }
1192 
1193  snxt = SolveNumericJT(p,v,fRmax,false);
1194  side = kRMax ;
1195 
1196  // rmin
1197 
1198  if ( fRmin )
1199  {
1200  G4double tolRMin = fRmin + fRminTolerance ;
1201 
1202  if ( (pt*pt < tolRMin*tolRMin) && (vDotNmax < 0) )
1203  {
1204  if (calcNorm) { *validNorm = false ; } // Concave surface of the torus
1205  return snxt = 0 ; // Leaving by Rmin immediately
1206  }
1207 
1208  sd[0] = SolveNumericJT(p,v,fRmin,false);
1209  if ( sd[0] < snxt )
1210  {
1211  snxt = sd[0] ;
1212  side = kRMin ;
1213  }
1214  }
1215 
1216 #else
1217 
1218  // this is the "conservative" version which return always validnorm = false
1219  // NOTE: using this version the unit test testG4Torus will break
1220 
1221  snxt = SolveNumericJT(p,v,fRmax,false);
1222  side = kRMax ;
1223 
1224  if ( fRmin )
1225  {
1226  sd[0] = SolveNumericJT(p,v,fRmin,false);
1227  if ( sd[0] < snxt )
1228  {
1229  snxt = sd[0] ;
1230  side = kRMin ;
1231  }
1232  }
1233 
1234  if ( calcNorm && (snxt == 0.0) )
1235  {
1236  *validNorm = false ; // Leaving solid, but possible re-intersection
1237  return snxt ;
1238  }
1239 
1240 #endif
1241 
1242  if (fDPhi < twopi) // Phi Intersections
1243  {
1244  sinSPhi = std::sin(fSPhi) ;
1245  cosSPhi = std::cos(fSPhi) ;
1246  ePhi = fSPhi + fDPhi ;
1247  sinEPhi = std::sin(ePhi) ;
1248  cosEPhi = std::cos(ePhi) ;
1249  cPhi = fSPhi + fDPhi*0.5 ;
1250  sinCPhi = std::sin(cPhi) ;
1251  cosCPhi = std::cos(cPhi) ;
1252 
1253  // angle calculation with correction
1254  // of difference in domain of atan2 and Sphi
1255  //
1256  vphi = std::atan2(v.y(),v.x()) ;
1257 
1258  if ( vphi < fSPhi - halfAngTolerance ) { vphi += twopi; }
1259  else if ( vphi > ePhi + halfAngTolerance ) { vphi -= twopi; }
1260 
1261  if ( p.x() || p.y() ) // Check if on z axis (rho not needed later)
1262  {
1263  pDistS = p.x()*sinSPhi - p.y()*cosSPhi ; // pDist -ve when inside
1264  pDistE = -p.x()*sinEPhi + p.y()*cosEPhi ;
1265 
1266  // Comp -ve when in direction of outwards normal
1267  //
1268  compS = -sinSPhi*v.x() + cosSPhi*v.y() ;
1269  compE = sinEPhi*v.x() - cosEPhi*v.y() ;
1270  sidephi = kNull ;
1271 
1272  if( ( (fDPhi <= pi) && ( (pDistS <= halfCarTolerance)
1273  && (pDistE <= halfCarTolerance) ) )
1274  || ( (fDPhi > pi) && !((pDistS > halfCarTolerance)
1275  && (pDistE > halfCarTolerance) ) ) )
1276  {
1277  // Inside both phi *full* planes
1278 
1279  if ( compS < 0 )
1280  {
1281  sphi = pDistS/compS ;
1282 
1283  if (sphi >= -halfCarTolerance)
1284  {
1285  xi = p.x() + sphi*v.x() ;
1286  yi = p.y() + sphi*v.y() ;
1287 
1288  // Check intersecting with correct half-plane
1289  // (if not -> no intersect)
1290  //
1291  if ( (std::fabs(xi)<=kCarTolerance)
1292  && (std::fabs(yi)<=kCarTolerance) )
1293  {
1294  sidephi = kSPhi;
1295  if ( ((fSPhi-halfAngTolerance)<=vphi)
1296  && ((ePhi+halfAngTolerance)>=vphi) )
1297  {
1298  sphi = kInfinity;
1299  }
1300  }
1301  else if ( yi*cosCPhi-xi*sinCPhi >=0 )
1302  {
1303  sphi = kInfinity ;
1304  }
1305  else
1306  {
1307  sidephi = kSPhi ;
1308  }
1309  }
1310  else
1311  {
1312  sphi = kInfinity ;
1313  }
1314  }
1315  else
1316  {
1317  sphi = kInfinity ;
1318  }
1319 
1320  if ( compE < 0 )
1321  {
1322  sphi2 = pDistE/compE ;
1323 
1324  // Only check further if < starting phi intersection
1325  //
1326  if ( (sphi2 > -kCarTolerance) && (sphi2 < sphi) )
1327  {
1328  xi = p.x() + sphi2*v.x() ;
1329  yi = p.y() + sphi2*v.y() ;
1330 
1331  if ( (std::fabs(xi)<=kCarTolerance)
1332  && (std::fabs(yi)<=kCarTolerance) )
1333  {
1334  // Leaving via ending phi
1335  //
1336  if( !( (fSPhi-halfAngTolerance <= vphi)
1337  && (ePhi+halfAngTolerance >= vphi) ) )
1338  {
1339  sidephi = kEPhi ;
1340  sphi = sphi2;
1341  }
1342  }
1343  else // Check intersecting with correct half-plane
1344  {
1345  if ( (yi*cosCPhi-xi*sinCPhi) >= 0)
1346  {
1347  // Leaving via ending phi
1348  //
1349  sidephi = kEPhi ;
1350  sphi = sphi2;
1351 
1352  }
1353  }
1354  }
1355  }
1356  }
1357  else
1358  {
1359  sphi = kInfinity ;
1360  }
1361  }
1362  else
1363  {
1364  // On z axis + travel not || to z axis -> if phi of vector direction
1365  // within phi of shape, Step limited by rmax, else Step =0
1366 
1367  vphi = std::atan2(v.y(),v.x());
1368 
1369  if ( ( fSPhi-halfAngTolerance <= vphi ) &&
1370  ( vphi <= ( ePhi+halfAngTolerance ) ) )
1371  {
1372  sphi = kInfinity;
1373  }
1374  else
1375  {
1376  sidephi = kSPhi ; // arbitrary
1377  sphi=0;
1378  }
1379  }
1380 
1381  // Order intersections
1382 
1383  if (sphi<snxt)
1384  {
1385  snxt=sphi;
1386  side=sidephi;
1387  }
1388  }
1389 
1390  G4double rhoi,it,iDotxyNmax ;
1391  // Note: by numerical computation we know where the ray hits the torus
1392  // So I propose to return the side where the ray hits
1393 
1394  if (calcNorm)
1395  {
1396  switch(side)
1397  {
1398  case kRMax: // n is unit vector
1399  xi = p.x() + snxt*v.x() ;
1400  yi = p.y() + snxt*v.y() ;
1401  zi = p.z() + snxt*v.z() ;
1402  rhoi = std::hypot(xi,yi);
1403  it = hypot(zi,rhoi-fRtor);
1404 
1405  iDotxyNmax = (1-fRtor/rhoi) ;
1406  if(iDotxyNmax >= -2.*fRmaxTolerance) // really convex part of Rmax
1407  {
1408  *n = G4ThreeVector( xi*(1-fRtor/rhoi)/it,
1409  yi*(1-fRtor/rhoi)/it,
1410  zi/it ) ;
1411  *validNorm = true ;
1412  }
1413  else
1414  {
1415  *validNorm = false ; // concave-convex part of Rmax
1416  }
1417  break ;
1418 
1419  case kRMin:
1420  *validNorm = false ; // Rmin is concave or concave-convex
1421  break;
1422 
1423  case kSPhi:
1424  if (fDPhi <= pi )
1425  {
1426  *n=G4ThreeVector(std::sin(fSPhi),-std::cos(fSPhi),0);
1427  *validNorm=true;
1428  }
1429  else
1430  {
1431  *validNorm = false ;
1432  }
1433  break ;
1434 
1435  case kEPhi:
1436  if (fDPhi <= pi)
1437  {
1438  *n=G4ThreeVector(-std::sin(fSPhi+fDPhi),std::cos(fSPhi+fDPhi),0);
1439  *validNorm=true;
1440  }
1441  else
1442  {
1443  *validNorm = false ;
1444  }
1445  break;
1446 
1447  default:
1448 
1449  // It seems we go here from time to time ...
1450 
1451  G4cout << G4endl;
1452  DumpInfo();
1453  std::ostringstream message;
1454  G4int oldprc = message.precision(16);
1455  message << "Undefined side for valid surface normal to solid."
1456  << G4endl
1457  << "Position:" << G4endl << G4endl
1458  << "p.x() = " << p.x()/mm << " mm" << G4endl
1459  << "p.y() = " << p.y()/mm << " mm" << G4endl
1460  << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl
1461  << "Direction:" << G4endl << G4endl
1462  << "v.x() = " << v.x() << G4endl
1463  << "v.y() = " << v.y() << G4endl
1464  << "v.z() = " << v.z() << G4endl << G4endl
1465  << "Proposed distance :" << G4endl << G4endl
1466  << "snxt = " << snxt/mm << " mm" << G4endl;
1467  message.precision(oldprc);
1468  G4Exception("G4Torus::DistanceToOut(p,v,..)",
1469  "GeomSolids1002",JustWarning, message);
1470  break;
1471  }
1472  }
1473  if ( snxt<halfCarTolerance ) { snxt=0 ; }
1474 
1475  return snxt;
1476 }
mm
static constexpr double mm
Definition: G4SIunits.hh:115
kInfinity
static const G4double kInfinity
Definition: geomdefs.hh:42
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
CLHEP::Hep3Vector::x
double x() const
JustWarning
Definition: G4ExceptionSeverity.hh:64
G4int
int G4int
Definition: G4Types.hh:78
CLHEP::Hep3Vector::z
double z() const
G4VSolid::DumpInfo
void DumpInfo() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4cout
G4GLOB_DLL std::ostream G4cout
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
CLHEP::Hep3Vector::y
double y() const
G4endl
#define G4endl
Definition: G4ios.hh:61
G4VSolid::kCarTolerance
G4double kCarTolerance
Definition: G4VSolid.hh:307
pi
static constexpr double pi
Definition: G4SIunits.hh:75
ESide
ESide
Definition: G4Cons.cc:76
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4Torus::DistanceToOut ( const G4ThreeVector &  p) const
virtual

Implements G4VSolid.

Definition at line 1482 of file G4Torus.cc.

1483 {
1484  G4double safe=0.0,safeR1,safeR2;
1485  G4double rho,pt ;
1486  G4double safePhi,phiC,cosPhiC,sinPhiC,ePhi;
1487 
1488  rho = std::hypot(p.x(),p.y());
1489  pt = std::hypot(p.z(),rho-fRtor);
1490 
1491 #ifdef G4CSGDEBUG
1492  if( Inside(p) == kOutside )
1493  {
1494  G4int oldprc = G4cout.precision(16) ;
1495  G4cout << G4endl ;
1496  DumpInfo();
1497  G4cout << "Position:" << G4endl << G4endl ;
1498  G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl ;
1499  G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl ;
1500  G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl ;
1501  G4cout.precision(oldprc);
1502  G4Exception("G4Torus::DistanceToOut(p)", "GeomSolids1002",
1503  JustWarning, "Point p is outside !?" );
1504  }
1505 #endif
1506 
1507  if (fRmin)
1508  {
1509  safeR1 = pt - fRmin ;
1510  safeR2 = fRmax - pt ;
1511 
1512  if (safeR1 < safeR2) { safe = safeR1 ; }
1513  else { safe = safeR2 ; }
1514  }
1515  else
1516  {
1517  safe = fRmax - pt ;
1518  }
1519 
1520  // Check if phi divided, Calc distances closest phi plane
1521  //
1522  if (fDPhi<twopi) // Above/below central phi of Torus?
1523  {
1524  phiC = fSPhi + fDPhi*0.5 ;
1525  cosPhiC = std::cos(phiC) ;
1526  sinPhiC = std::sin(phiC) ;
1527 
1528  if ((p.y()*cosPhiC-p.x()*sinPhiC)<=0)
1529  {
1530  safePhi = -(p.x()*std::sin(fSPhi) - p.y()*std::cos(fSPhi)) ;
1531  }
1532  else
1533  {
1534  ePhi = fSPhi + fDPhi ;
1535  safePhi = (p.x()*std::sin(ePhi) - p.y()*std::cos(ePhi)) ;
1536  }
1537  if (safePhi < safe) { safe = safePhi ; }
1538  }
1539  if (safe < 0) { safe = 0 ; }
1540  return safe ;
1541 }
G4Torus::Inside
EInside Inside(const G4ThreeVector &p) const
Definition: G4Torus.cc:599
mm
static constexpr double mm
Definition: G4SIunits.hh:115
CLHEP::Hep3Vector::x
double x() const
JustWarning
Definition: G4ExceptionSeverity.hh:64
G4int
int G4int
Definition: G4Types.hh:78
CLHEP::Hep3Vector::z
double z() const
G4VSolid::DumpInfo
void DumpInfo() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4cout
G4GLOB_DLL std::ostream G4cout
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
CLHEP::Hep3Vector::y
double y() const
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

void G4Torus::Extent ( G4ThreeVector &  pMin,
G4ThreeVector &  pMax 
) const
virtual

Reimplemented from G4VSolid.

Definition at line 417 of file G4Torus.cc.

418 {
419  G4double rmax = GetRmax();
420  G4double rtor = GetRtor();
421  G4double rint = rtor - rmax;
422  G4double rext = rtor + rmax;
423  G4double dz = rmax;
424 
425  // Find bounding box
426  //
427  if (GetDPhi() >= twopi)
428  {
429  pMin.set(-rext,-rext,-dz);
430  pMax.set( rext, rext, dz);
431  }
432  else
433  {
434  G4TwoVector vmin,vmax;
435  G4GeomTools::DiskExtent(rint,rext,
436  GetSinStartPhi(),GetCosStartPhi(),
437  GetSinEndPhi(),GetCosEndPhi(),
438  vmin,vmax);
439  pMin.set(vmin.x(),vmin.y(),-dz);
440  pMax.set(vmax.x(),vmax.y(), dz);
441  }
442 
443  // Check correctness of the bounding box
444  //
445  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
446  {
447  std::ostringstream message;
448  message << "Bad bounding box (min >= max) for solid: "
449  << GetName() << " !"
450  << "\npMin = " << pMin
451  << "\npMax = " << pMax;
452  G4Exception("G4Torus::Extent()", "GeomMgt0001", JustWarning, message);
453  DumpInfo();
454  }
455 }
CLHEP::Hep3Vector::set
void set(double x, double y, double z)
G4VSolid::GetName
G4String GetName() const
CLHEP::Hep2Vector::y
double y() const
CLHEP::Hep2Vector::x
double x() const
CLHEP::Hep3Vector::x
double x() const
G4Torus::GetRmax
G4double GetRmax() const
G4Torus::GetSinStartPhi
G4double GetSinStartPhi() const
JustWarning
Definition: G4ExceptionSeverity.hh:64
CLHEP::Hep3Vector::z
double z() const
G4VSolid::DumpInfo
void DumpInfo() const
G4Torus::GetRtor
G4double GetRtor() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4Torus::GetSinEndPhi
G4double GetSinEndPhi() const
G4Torus::GetCosEndPhi
G4double GetCosEndPhi() const
G4Torus::GetDPhi
G4double GetDPhi() const
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4Torus::GetCosStartPhi
G4double GetCosStartPhi() const
CLHEP::Hep3Vector::y
double y() const
G4double
double G4double
Definition: G4Types.hh:76
G4GeomTools::DiskExtent
static G4bool DiskExtent(G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
Definition: G4GeomTools.cc:378

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4Torus::GetCosEndPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetCosStartPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetCubicVolume ( )
inlinevirtual

Reimplemented from G4VSolid.

G4double G4Torus::GetDPhi ( ) const
inline

Here is the caller graph for this function:

G4GeometryType G4Torus::GetEntityType ( ) const
virtual

Implements G4VSolid.

Definition at line 1547 of file G4Torus.cc.

1548 {
1549  return G4String("G4Torus");
1550 }
G4ThreeVector G4Torus::GetPointOnSurface ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1588 of file G4Torus.cc.

1589 {
1590  G4double cosu, sinu,cosv, sinv, aOut, aIn, aSide, chose, phi, theta, rRand;
1591 
1592  phi = G4RandFlat::shoot(fSPhi,fSPhi+fDPhi);
1593  theta = G4RandFlat::shoot(0.,twopi);
1594 
1595  cosu = std::cos(phi); sinu = std::sin(phi);
1596  cosv = std::cos(theta); sinv = std::sin(theta);
1597 
1598  // compute the areas
1599 
1600  aOut = (fDPhi)*twopi*fRtor*fRmax;
1601  aIn = (fDPhi)*twopi*fRtor*fRmin;
1602  aSide = pi*(fRmax*fRmax-fRmin*fRmin);
1603 
1604  if ((fSPhi == 0) && (fDPhi == twopi)){ aSide = 0; }
1605  chose = G4RandFlat::shoot(0.,aOut + aIn + 2.*aSide);
1606 
1607  if(chose < aOut)
1608  {
1609  return G4ThreeVector ((fRtor+fRmax*cosv)*cosu,
1610  (fRtor+fRmax*cosv)*sinu, fRmax*sinv);
1611  }
1612  else if( (chose >= aOut) && (chose < aOut + aIn) )
1613  {
1614  return G4ThreeVector ((fRtor+fRmin*cosv)*cosu,
1615  (fRtor+fRmin*cosv)*sinu, fRmin*sinv);
1616  }
1617  else if( (chose >= aOut + aIn) && (chose < aOut + aIn + aSide) )
1618  {
1619  rRand = GetRadiusInRing(fRmin,fRmax);
1620  return G4ThreeVector ((fRtor+rRand*cosv)*std::cos(fSPhi),
1621  (fRtor+rRand*cosv)*std::sin(fSPhi), rRand*sinv);
1622  }
1623  else
1624  {
1625  rRand = GetRadiusInRing(fRmin,fRmax);
1626  return G4ThreeVector ((fRtor+rRand*cosv)*std::cos(fSPhi+fDPhi),
1627  (fRtor+rRand*cosv)*std::sin(fSPhi+fDPhi),
1628  rRand*sinv);
1629  }
1630 }
G4INCL::DeJongSpin::shoot
ThreeVector shoot(const G4int Ap, const G4int Af)
Definition: G4INCLDeJongSpin.cc:68
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
G4CSGSolid::GetRadiusInRing
G4double GetRadiusInRing(G4double rmin, G4double rmax) const
Definition: G4CSGSolid.cc:111
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
pi
static constexpr double pi
Definition: G4SIunits.hh:75
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4Torus::GetRmax ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetRmin ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetRtor ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetSinEndPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetSinStartPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetSPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Torus::GetSurfaceArea ( )
inlinevirtual

Reimplemented from G4VSolid.

EInside G4Torus::Inside ( const G4ThreeVector &  p) const
virtual

Implements G4VSolid.

Definition at line 599 of file G4Torus.cc.

600 {
601  G4double r, pt2, pPhi, tolRMin, tolRMax ;
602 
603  EInside in = kOutside ;
604 
605  // General precals
606  //
607  r = std::hypot(p.x(),p.y());
608  pt2 = p.z()*p.z() + (r-fRtor)*(r-fRtor);
609 
610  if (fRmin) tolRMin = fRmin + fRminTolerance ;
611  else tolRMin = 0 ;
612 
613  tolRMax = fRmax - fRmaxTolerance;
614 
615  if (pt2 >= tolRMin*tolRMin && pt2 <= tolRMax*tolRMax )
616  {
617  if ( fDPhi == twopi || pt2 == 0 ) // on torus swept axis
618  {
619  in = kInside ;
620  }
621  else
622  {
623  // Try inner tolerant phi boundaries (=>inside)
624  // if not inside, try outer tolerant phi boundaries
625 
626  pPhi = std::atan2(p.y(),p.x()) ;
627 
628  if ( pPhi < -halfAngTolerance ) { pPhi += twopi ; } // 0<=pPhi<2pi
629  if ( fSPhi >= 0 )
630  {
631  if ( (std::fabs(pPhi) < halfAngTolerance)
632  && (std::fabs(fSPhi + fDPhi - twopi) < halfAngTolerance) )
633  {
634  pPhi += twopi ; // 0 <= pPhi < 2pi
635  }
636  if ( (pPhi >= fSPhi + halfAngTolerance)
637  && (pPhi <= fSPhi + fDPhi - halfAngTolerance) )
638  {
639  in = kInside ;
640  }
641  else if ( (pPhi >= fSPhi - halfAngTolerance)
642  && (pPhi <= fSPhi + fDPhi + halfAngTolerance) )
643  {
644  in = kSurface ;
645  }
646  }
647  else // fSPhi < 0
648  {
649  if ( (pPhi <= fSPhi + twopi - halfAngTolerance)
650  && (pPhi >= fSPhi + fDPhi + halfAngTolerance) ) {;}
651  else
652  {
653  in = kSurface ;
654  }
655  }
656  }
657  }
658  else // Try generous boundaries
659  {
660  tolRMin = fRmin - fRminTolerance ;
661  tolRMax = fRmax + fRmaxTolerance ;
662 
663  if (tolRMin < 0 ) { tolRMin = 0 ; }
664 
665  if ( (pt2 >= tolRMin*tolRMin) && (pt2 <= tolRMax*tolRMax) )
666  {
667  if ( (fDPhi == twopi) || (pt2 == 0) ) // Continuous in phi or on z-axis
668  {
669  in = kSurface ;
670  }
671  else // Try outer tolerant phi boundaries only
672  {
673  pPhi = std::atan2(p.y(),p.x()) ;
674 
675  if ( pPhi < -halfAngTolerance ) { pPhi += twopi ; } // 0<=pPhi<2pi
676  if ( fSPhi >= 0 )
677  {
678  if ( (std::fabs(pPhi) < halfAngTolerance)
679  && (std::fabs(fSPhi + fDPhi - twopi) < halfAngTolerance) )
680  {
681  pPhi += twopi ; // 0 <= pPhi < 2pi
682  }
683  if ( (pPhi >= fSPhi - halfAngTolerance)
684  && (pPhi <= fSPhi + fDPhi + halfAngTolerance) )
685  {
686  in = kSurface;
687  }
688  }
689  else // fSPhi < 0
690  {
691  if ( (pPhi <= fSPhi + twopi - halfAngTolerance)
692  && (pPhi >= fSPhi + fDPhi + halfAngTolerance) ) {;}
693  else
694  {
695  in = kSurface ;
696  }
697  }
698  }
699  }
700  }
701  return in ;
702 }
CLHEP::Hep3Vector::x
double x() const
CLHEP::Hep3Vector::z
double z() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
EInside
EInside
Definition: geomdefs.hh:58
CLHEP::Hep3Vector::y
double y() const
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 221 of file G4Torus.cc.

222 {
223  // Check assignment to self
224  //
225  if (this == &rhs) { return *this; }
226 
227  // Copy base class data
228  //
229  G4CSGSolid::operator=(rhs);
230 
231  // Copy data
232  //
233  fRmin = rhs.fRmin; fRmax = rhs.fRmax;
234  fRtor = rhs.fRtor; fSPhi = rhs.fSPhi; fDPhi = rhs.fDPhi;
235  fRminTolerance = rhs.fRminTolerance; fRmaxTolerance = rhs.fRmaxTolerance;
236  kRadTolerance = rhs.kRadTolerance; kAngTolerance = rhs.kAngTolerance;
237  halfCarTolerance = rhs.halfCarTolerance;
238  halfAngTolerance = rhs.halfAngTolerance;
239 
240  return *this;
241 }
G4CSGSolid::operator=
G4CSGSolid & operator=(const G4CSGSolid &rhs)
Definition: G4CSGSolid.cc:91

Here is the call graph for this function:

void G4Torus::SetAllParameters ( G4double  pRmin,
G4double  pRmax,
G4double  pRtor,
G4double  pSPhi,
G4double  pDPhi 
)

Definition at line 103 of file G4Torus.cc.

108 {
109  const G4double fEpsilon = 4.e-11; // relative tolerance of radii
110 
111  fCubicVolume = 0.;
112  fSurfaceArea = 0.;
113  fRebuildPolyhedron = true;
114 
115  kRadTolerance = G4GeometryTolerance::GetInstance()->GetRadialTolerance();
116  kAngTolerance = G4GeometryTolerance::GetInstance()->GetAngularTolerance();
117 
118  halfCarTolerance = 0.5*kCarTolerance;
119  halfAngTolerance = 0.5*kAngTolerance;
120 
121  if ( pRtor >= pRmax+1.e3*kCarTolerance ) // Check swept radius, as in G4Cons
122  {
123  fRtor = pRtor ;
124  }
125  else
126  {
127  std::ostringstream message;
128  message << "Invalid swept radius for Solid: " << GetName() << G4endl
129  << " pRtor = " << pRtor << ", pRmax = " << pRmax;
130  G4Exception("G4Torus::SetAllParameters()",
131  "GeomSolids0002", FatalException, message);
132  }
133 
134  // Check radii, as in G4Cons
135  //
136  if ( pRmin < pRmax - 1.e2*kCarTolerance && pRmin >= 0 )
137  {
138  if (pRmin >= 1.e2*kCarTolerance) { fRmin = pRmin ; }
139  else { fRmin = 0.0 ; }
140  fRmax = pRmax ;
141  }
142  else
143  {
144  std::ostringstream message;
145  message << "Invalid values of radii for Solid: " << GetName() << G4endl
146  << " pRmin = " << pRmin << ", pRmax = " << pRmax;
147  G4Exception("G4Torus::SetAllParameters()",
148  "GeomSolids0002", FatalException, message);
149  }
150 
151  // Relative tolerances
152  //
153  fRminTolerance = (fRmin)
154  ? 0.5*std::max( kRadTolerance, fEpsilon*(fRtor-fRmin )) : 0;
155  fRmaxTolerance = 0.5*std::max( kRadTolerance, fEpsilon*(fRtor+fRmax) );
156 
157  // Check angles
158  //
159  if ( pDPhi >= twopi ) { fDPhi = twopi ; }
160  else
161  {
162  if (pDPhi > 0) { fDPhi = pDPhi ; }
163  else
164  {
165  std::ostringstream message;
166  message << "Invalid Z delta-Phi for Solid: " << GetName() << G4endl
167  << " pDPhi = " << pDPhi;
168  G4Exception("G4Torus::SetAllParameters()",
169  "GeomSolids0002", FatalException, message);
170  }
171  }
172 
173  // Ensure psphi in 0-2PI or -2PI-0 range if shape crosses 0
174  //
175  fSPhi = pSPhi;
176 
177  if (fSPhi < 0) { fSPhi = twopi-std::fmod(std::fabs(fSPhi),twopi) ; }
178  else { fSPhi = std::fmod(fSPhi,twopi) ; }
179 
180  if (fSPhi+fDPhi > twopi) { fSPhi-=twopi ; }
181 }
G4VSolid::GetName
G4String GetName() const
G4CSGSolid::fRebuildPolyhedron
G4bool fRebuildPolyhedron
Definition: G4CSGSolid.hh:80
G4CSGSolid::fCubicVolume
G4double fCubicVolume
Definition: G4CSGSolid.hh:78
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4CSGSolid::fSurfaceArea
G4double fSurfaceArea
Definition: G4CSGSolid.hh:79
G4GeometryTolerance::GetRadialTolerance
G4double GetRadialTolerance() const
Definition: G4GeometryTolerance.cc:97
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4INCL::Math::max
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Definition: G4INCLGlobals.hh:112
FatalException
Definition: G4ExceptionSeverity.hh:60
G4endl
#define G4endl
Definition: G4ios.hh:61
G4VSolid::kCarTolerance
G4double kCarTolerance
Definition: G4VSolid.hh:307
G4double
double G4double
Definition: G4Types.hh:76
G4GeometryTolerance::GetAngularTolerance
G4double GetAngularTolerance() const
Definition: G4GeometryTolerance.cc:92
G4GeometryTolerance::GetInstance
static G4GeometryTolerance * GetInstance()
Definition: G4GeometryTolerance.cc:73

Here is the call graph for this function:

Here is the caller graph for this function:

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

Reimplemented from G4CSGSolid.

Definition at line 1565 of file G4Torus.cc.

1566 {
1567  G4int oldprc = os.precision(16);
1568  os << "-----------------------------------------------------------\n"
1569  << " *** Dump for solid - " << GetName() << " ***\n"
1570  << " ===================================================\n"
1571  << " Solid type: G4Torus\n"
1572  << " Parameters: \n"
1573  << " inner radius: " << fRmin/mm << " mm \n"
1574  << " outer radius: " << fRmax/mm << " mm \n"
1575  << " swept radius: " << fRtor/mm << " mm \n"
1576  << " starting phi: " << fSPhi/degree << " degrees \n"
1577  << " delta phi : " << fDPhi/degree << " degrees \n"
1578  << "-----------------------------------------------------------\n";
1579  os.precision(oldprc);
1580 
1581  return os;
1582 }
G4VSolid::GetName
G4String GetName() const
mm
static constexpr double mm
Definition: G4SIunits.hh:115
G4int
int G4int
Definition: G4Types.hh:78
degree
static constexpr double degree
Definition: G4SIunits.hh:144

Here is the call graph for this function:

G4ThreeVector G4Torus::SurfaceNormal ( const G4ThreeVector &  p) const
virtual

Implements G4VSolid.

Definition at line 710 of file G4Torus.cc.

711 {
712  G4int noSurfaces = 0;
713  G4double rho, pt, pPhi;
714  G4double distRMin = kInfinity;
715  G4double distSPhi = kInfinity, distEPhi = kInfinity;
716 
717  // To cope with precision loss
718  //
719  const G4double delta = std::max(10.0*kCarTolerance,
720  1.0e-8*(fRtor+fRmax));
721  const G4double dAngle = 10.0*kAngTolerance;
722 
723  G4ThreeVector nR, nPs, nPe;
724  G4ThreeVector norm, sumnorm(0.,0.,0.);
725 
726  rho = std::hypot(p.x(),p.y());
727  pt = std::hypot(p.z(),rho-fRtor);
728 
729  G4double distRMax = std::fabs(pt - fRmax);
730  if(fRmin) distRMin = std::fabs(pt - fRmin);
731 
732  if( rho > delta && pt != 0.0 )
733  {
734  G4double redFactor= (rho-fRtor)/rho;
735  nR = G4ThreeVector( p.x()*redFactor, // p.x()*(1.-fRtor/rho),
736  p.y()*redFactor, // p.y()*(1.-fRtor/rho),
737  p.z() );
738  nR *= 1.0/pt;
739  }
740 
741  if ( fDPhi < twopi ) // && rho ) // old limitation against (0,0,z)
742  {
743  if ( rho )
744  {
745  pPhi = std::atan2(p.y(),p.x());
746 
747  if(pPhi < fSPhi-delta) { pPhi += twopi; }
748  else if(pPhi > fSPhi+fDPhi+delta) { pPhi -= twopi; }
749 
750  distSPhi = std::fabs( pPhi - fSPhi );
751  distEPhi = std::fabs(pPhi-fSPhi-fDPhi);
752  }
753  nPs = G4ThreeVector(std::sin(fSPhi),-std::cos(fSPhi),0);
754  nPe = G4ThreeVector(-std::sin(fSPhi+fDPhi),std::cos(fSPhi+fDPhi),0);
755  }
756  if( distRMax <= delta )
757  {
758  noSurfaces ++;
759  sumnorm += nR;
760  }
761  else if( fRmin && (distRMin <= delta) ) // Must not be on both Outer and Inner
762  {
763  noSurfaces ++;
764  sumnorm -= nR;
765  }
766 
767  // To be on one of the 'phi' surfaces,
768  // it must be within the 'tube' - with tolerance
769 
770  if( (fDPhi < twopi) && (fRmin-delta <= pt) && (pt <= (fRmax+delta)) )
771  {
772  if (distSPhi <= dAngle)
773  {
774  noSurfaces ++;
775  sumnorm += nPs;
776  }
777  if (distEPhi <= dAngle)
778  {
779  noSurfaces ++;
780  sumnorm += nPe;
781  }
782  }
783  if ( noSurfaces == 0 )
784  {
785 #ifdef G4CSGDEBUG
786  G4ExceptionDescription ed;
787  ed.precision(16);
788 
789  EInside inIt= Inside( p );
790 
791  if( inIt != kSurface )
792  {
793  ed << " ERROR> Surface Normal was called for Torus,"
794  << " with point not on surface." << G4endl;
795  }
796  else
797  {
798  ed << " ERROR> Surface Normal has not found a surface, "
799  << " despite the point being on the surface. " <<G4endl;
800  }
801 
802  if( inIt != kInside)
803  {
804  ed << " Safety (Dist To In) = " << DistanceToIn(p) << G4endl;
805  }
806  if( inIt != kOutside)
807  {
808  ed << " Safety (Dist to Out) = " << DistanceToOut(p) << G4endl;
809  }
810  ed << " Coordinates of point : " << p << G4endl;
811  ed << " Parameters of solid : " << G4endl << *this << G4endl;
812 
813  if( inIt == kSurface )
814  {
815  G4Exception("G4Torus::SurfaceNormal(p)", "GeomSolids1002",
816  JustWarning, ed,
817  "Failing to find normal, even though point is on surface!");
818  }
819  else
820  {
821  static const char* NameInside[3]= { "Inside", "Surface", "Outside" };
822  ed << " The point is " << NameInside[inIt] << " the solid. "<< G4endl;
823  G4Exception("G4Torus::SurfaceNormal(p)", "GeomSolids1002",
824  JustWarning, ed, "Point p is not on surface !?" );
825  }
826 #endif
827  norm = ApproxSurfaceNormal(p);
828  }
829  else if ( noSurfaces == 1 ) { norm = sumnorm; }
830  else { norm = sumnorm.unit(); }
831 
832  return norm ;
833 }
G4Torus::Inside
EInside Inside(const G4ThreeVector &p) const
Definition: G4Torus.cc:599
kInfinity
static const G4double kInfinity
Definition: geomdefs.hh:42
G4ExceptionDescription
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
CLHEP::Hep3Vector::x
double x() const
G4Torus::DistanceToIn
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4Torus.cc:947
JustWarning
Definition: G4ExceptionSeverity.hh:64
G4int
int G4int
Definition: G4Types.hh:78
CLHEP::Hep3Vector::z
double z() const
twopi
static constexpr double twopi
Definition: G4SIunits.hh:76
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4Torus::DistanceToOut
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=G4bool(false), G4bool *validNorm=0, G4ThreeVector *n=0) const
Definition: G4Torus.cc:1140
G4INCL::Math::max
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Definition: G4INCLGlobals.hh:112
EInside
EInside
Definition: geomdefs.hh:58
CLHEP::Hep3Vector::unit
Hep3Vector unit() const
CLHEP::Hep3Vector::y
double y() const
G4endl
#define G4endl
Definition: G4ios.hh:61
G4VSolid::kCarTolerance
G4double kCarTolerance
Definition: G4VSolid.hh:307
G4double
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

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