G4Torus represents a torus or torus segment with curved sides parallel to the z-axis. The torus has a specified swept radius about which it is centered, and a given minimum and maximum radius. A minimum radius of 0 signifies a filled torus. The torus segment is specified by starting and delta angles for phi, with 0 being the +x axis, PI/2 the +y axis. A delta angle of 2PI signifies a complete, unsegmented torus/cylinder. More...

#include <G4Torus.hh>

Inheritance diagram for G4Torus:

Public Member Functions
	G4Torus (const G4String &pName, G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
	~G4Torus () override=default
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 () override
G4double	GetSurfaceArea () override
void	ComputeDimensions (G4VPVParameterisation p, const G4int n, const G4VPhysicalVolume pRep) override
void	BoundingLimits (G4ThreeVector &pMin, G4ThreeVector &pMax) const override
G4bool	CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pmin, G4double &pmax) const override
EInside	Inside (const G4ThreeVector &p) const override
G4ThreeVector	SurfaceNormal (const G4ThreeVector &p) const override
G4double	DistanceToIn (const G4ThreeVector &p, const G4ThreeVector &v) const override
G4double	DistanceToIn (const G4ThreeVector &p) const override
G4double	DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool validNorm=nullptr, G4ThreeVector n=nullptr) const override
G4double	DistanceToOut (const G4ThreeVector &p) const override
G4GeometryType	GetEntityType () const override
G4ThreeVector	GetPointOnSurface () const override
G4VSolid *	Clone () const override
std::ostream &	StreamInfo (std::ostream &os) const override
void	DescribeYourselfTo (G4VGraphicsScene &scene) const override
G4Polyhedron *	CreatePolyhedron () const override
void	SetAllParameters (G4double pRmin, G4double pRmax, G4double pRtor, G4double pSPhi, G4double pDPhi)
	G4Torus (__void__ &)
	G4Torus (const G4Torus &rhs)=default
G4Torus &	operator= (const G4Torus &rhs)
Public Member Functions inherited from G4CSGSolid
	G4CSGSolid (const G4String &pName)
	~G4CSGSolid () override
G4Polyhedron *	GetPolyhedron () const override
	G4CSGSolid (__void__ &)
	G4CSGSolid (const G4CSGSolid &rhs)
G4CSGSolid &	operator= (const G4CSGSolid &rhs)
Public Member Functions inherited from G4VSolid
	G4VSolid (const G4String &name)
virtual	~G4VSolid ()
	G4VSolid (const G4VSolid &rhs)
G4VSolid &	operator= (const G4VSolid &rhs)
G4bool	operator== (const G4VSolid &s) const
G4String	GetName () const
void	SetName (const G4String &name)
G4double	GetTolerance () const
virtual G4int	GetNumOfConstituents () const
virtual G4bool	IsFaceted () 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__ &)
G4double	EstimateCubicVolume (G4int nStat, G4double epsilon) const
G4double	EstimateSurfaceArea (G4int nStat, G4double epsilon) 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 = 0.0
G4double	fSurfaceArea = 0.0
G4bool	fRebuildPolyhedron = false
G4Polyhedron *	fpPolyhedron = nullptr
Protected Attributes inherited from G4VSolid
G4double	kCarTolerance

Detailed Description

G4Torus represents a torus or torus segment with curved sides parallel to the z-axis. The torus has a specified swept radius about which it is centered, and a given minimum and maximum radius. A minimum radius of 0 signifies a filled torus. The torus segment is specified by starting and delta angles for phi, with 0 being the +x axis, PI/2 the +y axis. A delta angle of 2PI signifies a complete, unsegmented torus/cylinder.

Definition at line 101 of file G4Torus.hh.

Constructor & Destructor Documentation

◆ G4Torus() [1/3]

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

Constructs a torus or torus segment with the given name and dimensions.

Parameters

[in]	pName	The name of the solid.
[in]	pRmin	Inner radius.
[in]	pRmax	Outer radius.
[in]	pRtor	Swept radius of torus.
[in]	pSPhi	Starting Phi angle in radians adjusted such that fSPhi+fDPhi<=2PI, fSPhi>-2PI.
[in]	pDPhi	Delta angle of the segment in radians.

Definition at line 78 of file G4Torus.cc.

  : G4CSGSolid(pName)
{
  SetAllParameters(pRmin, pRmax, pRtor, pSPhi, pDPhi);
}

Referenced by Clone(), G4Torus(), and operator=().

◆ ~G4Torus()

G4Torus::~G4Torus ( )

overridedefault

Default destructor.

◆ G4Torus() [2/3]

G4Torus::G4Torus ( __void__ & a )

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

Definition at line 179 of file G4Torus.cc.

  : G4CSGSolid(a)
{
}

◆ G4Torus() [3/3]

G4Torus::G4Torus ( const G4Torus & rhs )

default

Copy constructor and assignment operator.

Member Function Documentation

◆ BoundingLimits()

void G4Torus::BoundingLimits	(	G4ThreeVector &	pMin,
		G4ThreeVector &	pMax ) const

overridevirtual

Computes the bounding limits of the solid.

Parameters

[out]	pMin	The minimum bounding limit point.
[out]	pMax	The maximum bounding limit point.

Reimplemented from G4VSolid.

Definition at line 384 of file G4Torus.cc.

{
  G4double rmax = GetRmax();
  G4double rtor = GetRtor();
  G4double rint = rtor - rmax;
  G4double rext = rtor + rmax;
  G4double dz   = rmax;
 
  // Find bounding box
  //
  if (GetDPhi() >= twopi)
  {
    pMin.set(-rext,-rext,-dz);
    pMax.set( rext, rext, dz);
  }
  else
  {
    G4TwoVector vmin,vmax;
    G4GeomTools::DiskExtent(rint,rext,
                            GetSinStartPhi(),GetCosStartPhi(),
                            GetSinEndPhi(),GetCosEndPhi(),
                            vmin,vmax);
    pMin.set(vmin.x(),vmin.y(),-dz);
    pMax.set(vmax.x(),vmax.y(), dz);
  }
 
  // Check correctness of the bounding box
  //
  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
  {
    std::ostringstream message;
    message << "Bad bounding box (min >= max) for solid: "
            << GetName() << " !"
            << "\npMin = " << pMin
            << "\npMax = " << pMax;
    G4Exception("G4Torus::BoundingLimits()", "GeomMgt0001",
                JustWarning, message);
    DumpInfo();
  }
}

Referenced by CalculateExtent().

◆ CalculateExtent()

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

overridevirtual

Calculates the minimum and maximum extent of the solid, when under the specified transform, and within the specified limits.

Parameters

[in]	pAxis	The axis along which compute the extent.
[in]	pVoxelLimit	The limiting space dictated by voxels.
[in]	pTransform	The internal transformation applied to the solid.
[out]	pMin	The minimum extent value.
[out]	pMax	The maximum extent value.

Returns: True if the solid is intersected by the extent region.

Implements G4VSolid.

Definition at line 429 of file G4Torus.cc.

{
  G4ThreeVector bmin, bmax;
  G4bool exist;
 
  // Get bounding box
  BoundingLimits(bmin,bmax);
 
  // Check bounding box
  G4BoundingEnvelope bbox(bmin,bmax);
#ifdef G4BBOX_EXTENT
  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
#endif
  if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax))
  {
    return exist = pMin < pMax;
  }
 
  // Get parameters of the solid
  G4double rmin = GetRmin();
  G4double rmax = GetRmax();
  G4double rtor = GetRtor();
  G4double dphi = GetDPhi();
  G4double sinStart = GetSinStartPhi();
  G4double cosStart = GetCosStartPhi();
  G4double sinEnd   = GetSinEndPhi();
  G4double cosEnd   = GetCosEndPhi();
  G4double rint = rtor - rmax;
  G4double rext = rtor + rmax;
 
  // Find bounding envelope and calculate extent
  //
  static const G4int NPHI  = 24; // number of steps for whole torus
  static const G4int NDISK = 16; // number of steps for disk
  static const G4double sinHalfDisk = std::sin(pi/NDISK);
  static const G4double cosHalfDisk = std::cos(pi/NDISK);
  static const G4double sinStepDisk = 2.*sinHalfDisk*cosHalfDisk;
  static const G4double cosStepDisk = 1. - 2.*sinHalfDisk*sinHalfDisk;
 
  G4double astep = (360/NPHI)*deg; // max angle for one slice in phi
  G4int    kphi  = (dphi <= astep) ? 1 : (G4int)((dphi-deg)/astep) + 1;
  G4double ang   = dphi/kphi;
 
  G4double sinHalf = std::sin(0.5*ang);
  G4double cosHalf = std::cos(0.5*ang);
  G4double sinStep = 2.*sinHalf*cosHalf;
  G4double cosStep = 1. - 2.*sinHalf*sinHalf;
 
  // define vectors for bounding envelope
  G4ThreeVectorList pols[NDISK+1];
  for (auto & pol : pols) { pol.resize(4); }
 
  std::vector<const G4ThreeVectorList *> polygons;
  polygons.resize(NDISK+1);
  for (G4int k=0; k<NDISK+1; ++k) { polygons[k] = &pols[k]; }
 
  // set internal and external reference circles
  G4TwoVector rzmin[NDISK];
  G4TwoVector rzmax[NDISK];
 
  if ((rtor-rmin*sinHalfDisk)/cosHalf > (rtor+rmin*sinHalfDisk)) { rmin = 0; }
  rmax /= cosHalfDisk;
  G4double sinCurDisk = sinHalfDisk;
  G4double cosCurDisk = cosHalfDisk;
  for (G4int k=0; k<NDISK; ++k)
  {
    G4double rmincur = rtor + rmin*cosCurDisk;
    if (cosCurDisk < 0 && rmin > 0) { rmincur /= cosHalf; }
    rzmin[k].set(rmincur,rmin*sinCurDisk);
 
    G4double rmaxcur = rtor + rmax*cosCurDisk;
    if (cosCurDisk > 0) { rmaxcur /= cosHalf; }
    rzmax[k].set(rmaxcur,rmax*sinCurDisk);
 
    G4double sinTmpDisk = sinCurDisk;
    sinCurDisk = sinCurDisk*cosStepDisk + cosCurDisk*sinStepDisk;
    cosCurDisk = cosCurDisk*cosStepDisk - sinTmpDisk*sinStepDisk;
  }
 
  // Loop along slices in Phi. The extent is calculated as cumulative
  // extent of the slices
  pMin =  kInfinity;
  pMax = -kInfinity;
  G4double eminlim = pVoxelLimit.GetMinExtent(pAxis);
  G4double emaxlim = pVoxelLimit.GetMaxExtent(pAxis);
  G4double sinCur1 = 0, cosCur1 = 0, sinCur2 = 0, cosCur2 = 0;
  for (G4int i=0; i<kphi+1; ++i)
  {
    if (i == 0)
    {
      sinCur1 = sinStart;
      cosCur1 = cosStart;
      sinCur2 = sinCur1*cosHalf + cosCur1*sinHalf;
      cosCur2 = cosCur1*cosHalf - sinCur1*sinHalf;
    }
    else
    {
      sinCur1 = sinCur2;
      cosCur1 = cosCur2;
      sinCur2 = (i == kphi) ? sinEnd : sinCur1*cosStep + cosCur1*sinStep;
      cosCur2 = (i == kphi) ? cosEnd : cosCur1*cosStep - sinCur1*sinStep;
    }
    for (G4int k=0; k<NDISK; ++k)
    {
      G4double r1 = rzmin[k].x(), r2 = rzmax[k].x();
      G4double z1 = rzmin[k].y(), z2 = rzmax[k].y();
      pols[k][0].set(r1*cosCur1,r1*sinCur1,z1);
      pols[k][1].set(r2*cosCur1,r2*sinCur1,z2);
      pols[k][2].set(r2*cosCur2,r2*sinCur2,z2);
      pols[k][3].set(r1*cosCur2,r1*sinCur2,z1);
    }
    pols[NDISK] = pols[0];
 
    // get bounding box of current slice
    G4TwoVector vmin,vmax;
    G4GeomTools::
      DiskExtent(rint,rext,sinCur1,cosCur1,sinCur2,cosCur2,vmin,vmax);
    bmin.setX(vmin.x()); bmin.setY(vmin.y());
    bmax.setX(vmax.x()); bmax.setY(vmax.y());
 
    // set bounding envelope for current slice and adjust extent
    G4double emin,emax;
    G4BoundingEnvelope benv(bmin,bmax,polygons);
    if (!benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,emin,emax))
    {
      continue;
    }
    if (emin < pMin) { pMin = emin; }
    if (emax > pMax) { pMax = emax; }
    if (eminlim > pMin && emaxlim < pMax) { break; } // max possible extent
  }
  return (pMin < pMax);
}

◆ Clone()

G4VSolid * G4Torus::Clone ( ) const

overridevirtual

Makes a clone of the object for use in multi-treading.

Returns: A pointer to the new cloned allocated solid.

Reimplemented from G4VSolid.

Definition at line 1556 of file G4Torus.cc.

{
  return new G4Torus(*this);
}

◆ ComputeDimensions()

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

overridevirtual

Dispatch method for parameterisation replication mechanism and dimension computation.

Reimplemented from G4VSolid.

Definition at line 215 of file G4Torus.cc.

{
  p->ComputeDimensions(*this,n,pRep);
}

◆ CreatePolyhedron()

G4Polyhedron * G4Torus::CreatePolyhedron ( ) const

overridevirtual

Creates a Polyhedron used for Visualisation. It is the caller's responsibility to delete it. A null pointer means "not created".

Reimplemented from G4VSolid.

Definition at line 1666 of file G4Torus.cc.

{
  return new G4PolyhedronTorus (fRmin, fRmax, fRtor, fSPhi, fDPhi);
}

◆ DescribeYourselfTo()

void G4Torus::DescribeYourselfTo ( G4VGraphicsScene & scene ) const

overridevirtual

Methods for creating graphical representations (i.e. for visualisation).

Implements G4VSolid.

Definition at line 1661 of file G4Torus.cc.

{
  scene.AddSolid (*this);
}

◆ DistanceToIn() [1/2]

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

overridevirtual

Calculates the distance to the nearest surface of a shape from an outside point. The distance can be an underestimate.

Parameters

[in] p The point at offset p.

Returns: The safety distance to enter the shape.

Implements G4VSolid.

Definition at line 1095 of file G4Torus.cc.

{
  G4double safe=0.0, safe1, safe2 ;
  G4double phiC, cosPhiC, sinPhiC, safePhi, ePhi, cosPsi ;
  G4double rho, pt ;
 
  rho = std::hypot(p.x(),p.y());
  pt  = std::hypot(p.z(),rho-fRtor);
  safe1 = fRmin - pt ;
  safe2 = pt - fRmax ;
 
  if (safe1 > safe2)  { safe = safe1; }
  else                { safe = safe2; }
 
  if ( fDPhi < twopi && (rho != 0.0) )
  {
    phiC    = fSPhi + fDPhi*0.5 ;
    cosPhiC = std::cos(phiC) ;
    sinPhiC = std::sin(phiC) ;
    cosPsi  = (p.x()*cosPhiC + p.y()*sinPhiC)/rho ;
 
    if (cosPsi < std::cos(fDPhi*0.5) ) // Psi=angle from central phi to point
    {                                  // Point lies outside phi range
      if ((p.y()*cosPhiC - p.x()*sinPhiC) <= 0 )
      {
        safePhi = std::fabs(p.x()*std::sin(fSPhi) - p.y()*std::cos(fSPhi)) ;
      }
      else
      {
        ePhi    = fSPhi + fDPhi ;
        safePhi = std::fabs(p.x()*std::sin(ePhi) - p.y()*std::cos(ePhi)) ;
      }
      if (safePhi > safe)  { safe = safePhi ; }
    }
  }
  if (safe < 0 )  { safe = 0 ; }
  return safe;
}

◆ DistanceToIn() [2/2]

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

overridevirtual

Returns the distance along the normalised vector 'v' to the shape, from the point at offset 'p'. If there is no intersection, returns kInfinity. The first intersection resulting from 'leaving' a surface/volume is discarded. Hence, it is tolerant of points on the surface of the shape.

Parameters

[in]	p	The point at offset p.
[in]	v	The normalised direction vector.

Returns: The distance to enter the shape.

Implements G4VSolid.

Definition at line 917 of file G4Torus.cc.

{
  // Get bounding box of full torus
  //
  G4double boxDx  = fRtor + fRmax;
  G4double boxDy  = boxDx;
  G4double boxDz  = fRmax;
  G4double boxMax = boxDx;
  G4double boxMin = boxDz;
 
  // Check if point is traveling away
  //
  G4double distX = std::abs(p.x()) - boxDx;
  G4double distY = std::abs(p.y()) - boxDy;
  G4double distZ = std::abs(p.z()) - boxDz;
  if (distX >= -halfCarTolerance && p.x()*v.x() >= 0) { return kInfinity; }
  if (distY >= -halfCarTolerance && p.y()*v.y() >= 0) { return kInfinity; }
  if (distZ >= -halfCarTolerance && p.z()*v.z() >= 0) { return kInfinity; }
 
  // Calculate safety distance to bounding box
  // If point is too far, move it closer and calculate distance
  //
  G4double Dmax = 32*boxMax;
  G4double safe = std::max(std::max(distX,distY),distZ);
  if (safe > Dmax)
  {
    G4double dist = safe - 1.e-8*safe - boxMin; // stay outside after the move
    dist += DistanceToIn(p + dist*v, v);
    return (dist >= kInfinity) ? kInfinity : dist;
  }
 
  // Find intersection with torus
  //
  G4double snxt=kInfinity, sphi=kInfinity; // snxt = default return value
 
  G4double  sd[4] ;
 
  // Precalculated trig for phi intersections - used by r,z intersections to
  //                                            check validity
 
  G4bool seg;        // true if segmented
  G4double hDPhi;    // half dphi
  G4double cPhi,sinCPhi=0.,cosCPhi=0.;  // central phi
 
  G4double tolORMin2;  // `generous' radii squared
  G4double tolORMax2;
 
  G4double Dist,xi,yi,zi,rhoi,it2; // Intersection point variables
 
  G4double Comp;
  G4double cosSPhi,sinSPhi;       // Trig for phi start intersect
  G4double ePhi,cosEPhi,sinEPhi;  // for phi end intersect
 
  // Set phi divided flag and precalcs
  //
  if ( fDPhi < twopi )
  {
    seg        = true ;
    hDPhi      = 0.5*fDPhi ;    // half delta phi
    cPhi       = fSPhi + hDPhi ;
    sinCPhi    = std::sin(cPhi) ;
    cosCPhi    = std::cos(cPhi) ;
  }
  else
  {
    seg = false ;
  }
 
  if (fRmin > fRminTolerance) // Calculate tolerant rmin and rmax
  {
    tolORMin2 = (fRmin - fRminTolerance)*(fRmin - fRminTolerance) ;
  }
  else
  {
    tolORMin2 = 0 ;
  }
  tolORMax2 = (fRmax + fRmaxTolerance)*(fRmax + fRmaxTolerance) ;
 
  // Intersection with Rmax (possible return) and Rmin (must also check phi)
 
  snxt = SolveNumericJT(p,v,fRmax,true);
 
  if (fRmin != 0.0)  // Possible Rmin intersection
  {
    sd[0] = SolveNumericJT(p,v,fRmin,true);
    if ( sd[0] < snxt )  { snxt = sd[0] ; }
  }
 
  //
  // Phi segment intersection
  //
  // o Tolerant of points inside phi planes by up to kCarTolerance*0.5
  //
  // o NOTE: Large duplication of code between sphi & ephi checks
  //         -> only diffs: sphi -> ephi, Comp -> -Comp and half-plane
  //            intersection check <=0 -> >=0
  //         -> use some form of loop Construct ?
 
  if (seg)
  {
    sinSPhi = std::sin(fSPhi) ; // First phi surface ('S'tarting phi)
    cosSPhi = std::cos(fSPhi) ;
    Comp    = v.x()*sinSPhi - v.y()*cosSPhi ;  // Component in outwards
                                               // normal direction
    if (Comp < 0 )
    {
      Dist = (p.y()*cosSPhi - p.x()*sinSPhi) ;
 
      if (Dist < halfCarTolerance)
      {
        sphi = Dist/Comp ;
        if (sphi < snxt)
        {
          if ( sphi < 0 )  { sphi = 0 ; }
 
          xi    = p.x() + sphi*v.x() ;
          yi    = p.y() + sphi*v.y() ;
          zi    = p.z() + sphi*v.z() ;
          rhoi = std::hypot(xi,yi);
          it2 = zi*zi + (rhoi-fRtor)*(rhoi-fRtor);
 
          if ( it2 >= tolORMin2 && it2 <= tolORMax2 )
          {
            // r intersection is good - check intersecting
            // with correct half-plane
            //
            if ((yi*cosCPhi-xi*sinCPhi)<=0)  { snxt=sphi; }
          }
        }
      }
    }
    ePhi=fSPhi+fDPhi;    // Second phi surface ('E'nding phi)
    sinEPhi=std::sin(ePhi);
    cosEPhi=std::cos(ePhi);
    Comp=-(v.x()*sinEPhi-v.y()*cosEPhi);
 
    if ( Comp < 0 )   // Component in outwards normal dirn
    {
      Dist = -(p.y()*cosEPhi - p.x()*sinEPhi) ;
 
      if (Dist < halfCarTolerance )
      {
        sphi = Dist/Comp ;
 
        if (sphi < snxt )
        {
          if (sphi < 0 )  { sphi = 0 ; }
 
          xi    = p.x() + sphi*v.x() ;
          yi    = p.y() + sphi*v.y() ;
          zi    = p.z() + sphi*v.z() ;
          rhoi = std::hypot(xi,yi);
          it2 = zi*zi + (rhoi-fRtor)*(rhoi-fRtor);
 
          if (it2 >= tolORMin2 && it2 <= tolORMax2)
          {
            // z and r intersections good - check intersecting
            // with correct half-plane
            //
            if ((yi*cosCPhi-xi*sinCPhi)>=0)  { snxt=sphi; }
          }
        }
      }
    }
  }
  if(snxt < halfCarTolerance)  { snxt = 0.0 ; }
 
  return snxt ;
}

Referenced by DistanceToIn(), and SurfaceNormal().

◆ DistanceToOut() [1/2]

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

overridevirtual

Calculates the distance to the nearest surface of a shape from an inside point 'p'. The distance can be an underestimate.

Parameters

[in] p The point at offset p.

Returns: The safety distance to exit the shape.

Implements G4VSolid.

Definition at line 1482 of file G4Torus.cc.

{
  G4double safe=0.0,safeR1,safeR2;
  G4double rho,pt ;
  G4double safePhi,phiC,cosPhiC,sinPhiC,ePhi;
 
  rho = std::hypot(p.x(),p.y());
  pt  = std::hypot(p.z(),rho-fRtor);
 
#ifdef G4CSGDEBUG
  if( Inside(p) == kOutside )
  {
     G4long oldprc = G4cout.precision(16) ;
     G4cout << G4endl ;
     DumpInfo();
     G4cout << "Position:"  << G4endl << G4endl ;
     G4cout << "p.x() = "   << p.x()/mm << " mm" << G4endl ;
     G4cout << "p.y() = "   << p.y()/mm << " mm" << G4endl ;
     G4cout << "p.z() = "   << p.z()/mm << " mm" << G4endl << G4endl ;
     G4cout.precision(oldprc);
     G4Exception("G4Torus::DistanceToOut(p)", "GeomSolids1002",
                 JustWarning, "Point p is outside !?" );
  }
#endif
 
  if (fRmin != 0.0)
  {
    safeR1 = pt - fRmin ;
    safeR2 = fRmax - pt ;
 
    if (safeR1 < safeR2)  { safe = safeR1 ; }
    else                  { safe = safeR2 ; }
  }
  else
  {
    safe = fRmax - pt ;
  }
 
  // Check if phi divided, Calc distances closest phi plane
  //
  if (fDPhi < twopi) // Above/below central phi of Torus?
  {
    phiC    = fSPhi + fDPhi*0.5 ;
    cosPhiC = std::cos(phiC) ;
    sinPhiC = std::sin(phiC) ;
 
    if ((p.y()*cosPhiC-p.x()*sinPhiC)<=0)
    {
      safePhi = -(p.x()*std::sin(fSPhi) - p.y()*std::cos(fSPhi)) ;
    }
    else
    {
      ePhi    = fSPhi + fDPhi ;
      safePhi = (p.x()*std::sin(ePhi) - p.y()*std::cos(ePhi)) ;
    }
    if (safePhi < safe)  { safe = safePhi ; }
  }
  if (safe < 0)  { safe = 0 ; }
  return safe ;
}

◆ DistanceToOut() [2/2]

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

overridevirtual

Returns the distance along the normalised vector 'v' to the shape, from a point at an offset 'p' inside or on the surface of the shape. Intersections with surfaces, when the point is less than Tolerance/2 from a surface must be ignored.

Parameters

[in]	p	The point at offset p.
[in]	v	The normalised direction vector.
[in]	calcNorm	Flag to indicate if to calculate the normal or not.
[out]	validNorm	Flag set to true if the solid lies entirely behind or on the exiting surface. It is set false if the solid does not lie entirely behind or on the exiting surface. 'calcNorm' must be true, otherwise it is unused.
[out]	n	The exiting outwards normal vector (undefined Magnitude). 'calcNorm' must be true, otherwise it is unused.

Returns: The distance to exit the shape.

Implements G4VSolid.

Definition at line 1140 of file G4Torus.cc.

{
  ESide    side = kNull, sidephi = kNull ;
  G4double snxt = kInfinity, sphi, sd[4] ;
 
  // Vars for phi intersection
  //
  G4double sinSPhi, cosSPhi, ePhi, sinEPhi, cosEPhi;
  G4double cPhi, sinCPhi, cosCPhi ;
  G4double pDistS, compS, pDistE, compE, sphi2, xi, yi, zi, vphi ;
 
  // Radial Intersections Defenitions & General Precals
 
  //////////////////////// new calculation //////////////////////
 
#if 1
 
  // This is the version with the calculation of CalcNorm = true
  // To be done: Check the precision of this calculation.
  // If you want return always validNorm = false, then take the version below
 
 
  G4double rho = std::hypot(p.x(),p.y());
  G4double pt = hypot(p.z(),rho-fRtor);
 
  G4double pDotV = p.x()*v.x() + p.y()*v.y() + p.z()*v.z() ;
 
  G4double tolRMax = fRmax - fRmaxTolerance ;
 
  G4double vDotNmax   = pDotV - fRtor*(v.x()*p.x() + v.y()*p.y())/rho ;
  G4double pDotxyNmax = (1 - fRtor/rho) ;
 
  if( (pt*pt > tolRMax*tolRMax) && (vDotNmax >= 0) )
  {
    // On tolerant boundary & heading outwards (or perpendicular to) outer
    // radial surface -> leaving immediately with *n for really convex part
    // only
 
    if ( calcNorm && (pDotxyNmax >= -2.*fRmaxTolerance) )
    {
      *n = G4ThreeVector( p.x()*(1 - fRtor/rho)/pt,
                          p.y()*(1 - fRtor/rho)/pt,
                          p.z()/pt                  ) ;
      *validNorm = true ;
    }
 
    return snxt = 0 ; // Leaving by Rmax immediately
  }
 
  snxt = SolveNumericJT(p,v,fRmax,false);
  side = kRMax ;
 
  // rmin
 
  if ( fRmin != 0.0 )
  {
    G4double tolRMin = fRmin + fRminTolerance ;
 
    if ( (pt*pt < tolRMin*tolRMin) && (vDotNmax < 0) )
    {
      if (calcNorm)  { *validNorm = false ; } // Concave surface of the torus
      return  snxt = 0 ;                      // Leaving by Rmin immediately
    }
 
    sd[0] = SolveNumericJT(p,v,fRmin,false);
    if ( sd[0] < snxt )
    {
      snxt = sd[0] ;
      side = kRMin ;
    }
  }
 
#else
 
  // this is the "conservative" version which return always validnorm = false
  // NOTE: using this version the unit test testG4Torus will break
 
  snxt = SolveNumericJT(p,v,fRmax,false);
  side = kRMax ;
 
  if ( fRmin )
  {
    sd[0] = SolveNumericJT(p,v,fRmin,false);
    if ( sd[0] < snxt )
    {
      snxt = sd[0] ;
      side = kRMin ;
    }
  }
 
  if ( calcNorm && (snxt == 0.0) )
  {
    *validNorm = false ;    // Leaving solid, but possible re-intersection
    return snxt  ;
  }
 
#endif
 
  if (fDPhi < twopi)  // Phi Intersections
  {
    sinSPhi = std::sin(fSPhi) ;
    cosSPhi = std::cos(fSPhi) ;
    ePhi    = fSPhi + fDPhi ;
    sinEPhi = std::sin(ePhi) ;
    cosEPhi = std::cos(ePhi) ;
    cPhi    = fSPhi + fDPhi*0.5 ;
    sinCPhi = std::sin(cPhi) ;
    cosCPhi = std::cos(cPhi) ;
 
    // angle calculation with correction
    // of difference in domain of atan2 and Sphi
    //
    vphi = std::atan2(v.y(),v.x()) ;
 
    if ( vphi < fSPhi - halfAngTolerance  )    { vphi += twopi; }
    else if ( vphi > ePhi + halfAngTolerance ) { vphi -= twopi; }
 
    if ( (p.x() != 0.0) || (p.y() != 0.0) ) // Check if on z axis (rho not needed later)
    {
      pDistS = p.x()*sinSPhi - p.y()*cosSPhi ; // pDist -ve when inside
      pDistE = -p.x()*sinEPhi + p.y()*cosEPhi ;
 
      // Comp -ve when in direction of outwards normal
      //
      compS   = -sinSPhi*v.x() + cosSPhi*v.y() ;
      compE   = sinEPhi*v.x() - cosEPhi*v.y() ;
      sidephi = kNull ;
 
      if( ( (fDPhi <= pi) && ( (pDistS <= halfCarTolerance)
                            && (pDistE <= halfCarTolerance) ) )
       || ( (fDPhi >  pi) && ((pDistS <=  halfCarTolerance)
                            || (pDistE <=  halfCarTolerance) ) )  )
      {
        // Inside both phi *full* planes
 
        if ( compS < 0 )
        {
          sphi = pDistS/compS ;
 
          if (sphi >= -halfCarTolerance)
          {
            xi = p.x() + sphi*v.x() ;
            yi = p.y() + sphi*v.y() ;
 
            // Check intersecting with correct half-plane
            // (if not -> no intersect)
            //
            if ( (std::fabs(xi)<=kCarTolerance)
              && (std::fabs(yi)<=kCarTolerance) )
            {
              sidephi = kSPhi;
              if ( ((fSPhi-halfAngTolerance)<=vphi)
                && ((ePhi+halfAngTolerance)>=vphi) )
              {
                sphi = kInfinity;
              }
            }
            else if ( yi*cosCPhi-xi*sinCPhi >=0 )
            {
              sphi = kInfinity ;
            }
            else
            {
              sidephi = kSPhi ;
            }
          }
          else
          {
            sphi = kInfinity ;
          }
        }
        else
        {
          sphi = kInfinity ;
        }
 
        if ( compE < 0 )
        {
          sphi2 = pDistE/compE ;
 
          // Only check further if < starting phi intersection
          //
          if ( (sphi2 > -kCarTolerance) && (sphi2 < sphi) )
          {
            xi = p.x() + sphi2*v.x() ;
            yi = p.y() + sphi2*v.y() ;
 
            if ( (std::fabs(xi)<=kCarTolerance)
              && (std::fabs(yi)<=kCarTolerance) )
            {
              // Leaving via ending phi
              //
              if( (fSPhi-halfAngTolerance > vphi)
                  || (ePhi+halfAngTolerance < vphi) )
              {
                sidephi = kEPhi ;
                sphi = sphi2;
              }
            }
            else    // Check intersecting with correct half-plane
            {
              if ( (yi*cosCPhi-xi*sinCPhi) >= 0)
              {
                // Leaving via ending phi
                //
                sidephi = kEPhi ;
                sphi = sphi2;
 
              }
            }
          }
        }
      }
      else
      {
        sphi = kInfinity ;
      }
    }
    else
    {
      // On z axis + travel not || to z axis -> if phi of vector direction
      // within phi of shape, Step limited by rmax, else Step =0
 
      vphi = std::atan2(v.y(),v.x());
 
      if ( ( fSPhi-halfAngTolerance <= vphi ) &&
           ( vphi <= ( ePhi+halfAngTolerance ) ) )
      {
        sphi = kInfinity;
      }
      else
      {
        sidephi = kSPhi ; // arbitrary
        sphi=0;
      }
    }
 
    // Order intersections
 
    if (sphi<snxt)
    {
      snxt=sphi;
      side=sidephi;
    }
  }
 
  G4double rhoi,it,iDotxyNmax ;
  // Note: by numerical computation we know where the ray hits the torus
  // So I propose to return the side where the ray hits
 
  if (calcNorm)
  {
    switch(side)
    {
      case kRMax:                     // n is unit vector
        xi    = p.x() + snxt*v.x() ;
        yi    = p.y() + snxt*v.y() ;
        zi    = p.z() + snxt*v.z() ;
        rhoi = std::hypot(xi,yi);
        it = hypot(zi,rhoi-fRtor);
 
        iDotxyNmax = (1-fRtor/rhoi) ;
        if(iDotxyNmax >= -2.*fRmaxTolerance) // really convex part of Rmax
        {
          *n = G4ThreeVector( xi*(1-fRtor/rhoi)/it,
                              yi*(1-fRtor/rhoi)/it,
                              zi/it                 ) ;
          *validNorm = true ;
        }
        else
        {
          *validNorm = false ; // concave-convex part of Rmax
        }
        break ;
 
      case kRMin:
        *validNorm = false ;  // Rmin is concave or concave-convex
        break;
 
      case kSPhi:
        if (fDPhi <= pi )
        {
          *n=G4ThreeVector(std::sin(fSPhi),-std::cos(fSPhi),0);
          *validNorm=true;
        }
        else
        {
          *validNorm = false ;
        }
        break ;
 
      case kEPhi:
        if (fDPhi <= pi)
        {
          *n=G4ThreeVector(-std::sin(fSPhi+fDPhi),std::cos(fSPhi+fDPhi),0);
          *validNorm=true;
        }
        else
        {
          *validNorm = false ;
        }
        break;
 
      default:
 
        // It seems we go here from time to time ...
 
        G4cout << G4endl;
        DumpInfo();
        std::ostringstream message;
        G4long oldprc = message.precision(16);
        message << "Undefined side for valid surface normal to solid."
                << G4endl
                << "Position:"  << G4endl << G4endl
                << "p.x() = "   << p.x()/mm << " mm" << G4endl
                << "p.y() = "   << p.y()/mm << " mm" << G4endl
                << "p.z() = "   << p.z()/mm << " mm" << G4endl << G4endl
                << "Direction:" << G4endl << G4endl
                << "v.x() = "   << v.x() << G4endl
                << "v.y() = "   << v.y() << G4endl
                << "v.z() = "   << v.z() << G4endl << G4endl
                << "Proposed distance :" << G4endl << G4endl
                << "snxt = " << snxt/mm << " mm" << G4endl;
        message.precision(oldprc);
        G4Exception("G4Torus::DistanceToOut(p,v,..)",
                    "GeomSolids1002",JustWarning, message);
        break;
    }
  }
  if ( snxt<halfCarTolerance )  { snxt=0 ; }
 
  return snxt;
}

Referenced by SurfaceNormal().

◆ GetCosEndPhi()

G4double G4Torus::GetCosEndPhi ( ) const

inline

Referenced by BoundingLimits(), and CalculateExtent().

◆ GetCosStartPhi()

G4double G4Torus::GetCosStartPhi ( ) const

inline

Referenced by BoundingLimits(), and CalculateExtent().

◆ GetCubicVolume()

G4double G4Torus::GetCubicVolume ( )

overridevirtual

Returning an estimation of the solid volume (capacity) and surface area, in internal units.

Reimplemented from G4VSolid.

Definition at line 1627 of file G4Torus.cc.

{
  if (fCubicVolume == 0)
  {
    G4AutoLock l(&torusMutex);
    fCubicVolume = fDPhi*CLHEP::pi*fRtor*(fRmax*fRmax-fRmin*fRmin);
    l.unlock();
  }
  return fCubicVolume;
}

◆ GetDPhi()

G4double G4Torus::GetDPhi ( ) const

inline

Referenced by BoundingLimits(), CalculateExtent(), G4tgbGeometryDumper::GetSolidParams(), G4GDMLWriteParamvol::Torus_dimensionsWrite(), and G4GDMLWriteSolids::TorusWrite().

◆ GetEntityType()

G4GeometryType G4Torus::GetEntityType ( ) const

overridevirtual

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

Implements G4VSolid.

Definition at line 1547 of file G4Torus.cc.

{
  return {"G4Torus"};
}

◆ GetPointOnSurface()

G4ThreeVector G4Torus::GetPointOnSurface ( ) const

overridevirtual

Returns a random point located and uniformly distributed on the surface of the solid.

Reimplemented from G4VSolid.

Definition at line 1588 of file G4Torus.cc.

{
  G4double rrmin = fRmin*fRmin;
  G4double rrmax = fRmax*fRmax;
  G4double smax = twopi*fRtor*fDPhi*fRmax;
  G4double smin = twopi*fRtor*fDPhi*fRmin;
  G4double sphi = (fDPhi >= twopi) ? 0. : pi*(rrmax - rrmin);
 
  G4double u = G4QuickRand();
  G4double v = twopi*G4QuickRand();;
  G4double select = (smax + smin + 2.*sphi)*G4QuickRand();
  G4double phi, r, ds;
  if (select < 2.*sphi)
  {
    // phi cut
    phi = fSPhi + fDPhi*(G4double)(select < sphi);
    r = std::sqrt(rrmin + (rrmax - rrmin)*u);
    ds = fRtor + r*std::cos(v);
  }
  else
  {
    // toroidal surface (rejection sampling)
    phi = fSPhi + fDPhi*u;
    r = (select < 2.*sphi + smax) ? fRmax : fRmin;
    ds = fRtor + r*std::cos(v);
    for (auto i = 0; i < 10; ++i)
    {
      if ((fRtor + r)*G4QuickRand() < ds) { break; }
      v = twopi*G4QuickRand();
      ds = fRtor + r*std::cos(v);
    }
  }
  return { ds*std::cos(phi), ds*std::sin(phi), r*std::sin(v) };
}

◆ GetRmax()

G4double G4Torus::GetRmax ( ) const

inline

Referenced by BoundingLimits(), CalculateExtent(), G4tgbGeometryDumper::GetSolidParams(), G4GDMLWriteParamvol::Torus_dimensionsWrite(), and G4GDMLWriteSolids::TorusWrite().

◆ GetRmin()

G4double G4Torus::GetRmin ( ) const

inline

Accessors.

Referenced by CalculateExtent(), G4tgbGeometryDumper::GetSolidParams(), G4GDMLWriteParamvol::Torus_dimensionsWrite(), and G4GDMLWriteSolids::TorusWrite().

◆ GetRtor()

G4double G4Torus::GetRtor ( ) const

inline

Referenced by BoundingLimits(), CalculateExtent(), G4tgbGeometryDumper::GetSolidParams(), G4GDMLWriteParamvol::Torus_dimensionsWrite(), and G4GDMLWriteSolids::TorusWrite().

◆ GetSinEndPhi()

G4double G4Torus::GetSinEndPhi ( ) const

inline

Referenced by BoundingLimits(), and CalculateExtent().

◆ GetSinStartPhi()

G4double G4Torus::GetSinStartPhi ( ) const

inline

Referenced by BoundingLimits(), and CalculateExtent().

◆ GetSPhi()

G4double G4Torus::GetSPhi ( ) const

inline

Referenced by G4tgbGeometryDumper::GetSolidParams(), G4GDMLWriteParamvol::Torus_dimensionsWrite(), and G4GDMLWriteSolids::TorusWrite().

◆ GetSurfaceArea()

G4double G4Torus::GetSurfaceArea ( )

overridevirtual

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

Reimplemented from G4VSolid.

Definition at line 1642 of file G4Torus.cc.

{
  if (fSurfaceArea == 0)
  { 
    G4AutoLock l(&torusMutex);
    fSurfaceArea = fDPhi*CLHEP::twopi*fRtor*(fRmax+fRmin);
    if(fDPhi < CLHEP::twopi)
    {
      fSurfaceArea = fSurfaceArea + CLHEP::twopi*(fRmax*fRmax-fRmin*fRmin);
    } 
    l.unlock();
  }
  return fSurfaceArea;
}

◆ Inside()

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

overridevirtual

Concrete implementations of the expected query interfaces for solids, as defined in the base class G4VSolid.

Implements G4VSolid.

Definition at line 570 of file G4Torus.cc.

{
  G4double r, pt2, pPhi, tolRMin, tolRMax ;
 
  EInside in = kOutside ;
 
  // General precals
  //
  r   = std::hypot(p.x(),p.y());
  pt2 = p.z()*p.z() + (r-fRtor)*(r-fRtor);
 
  (fRmin != 0.0) ? (tolRMin = fRmin + fRminTolerance) : (tolRMin = 0);
 
  tolRMax = fRmax - fRmaxTolerance;
 
  if (pt2 >= tolRMin*tolRMin && pt2 <= tolRMax*tolRMax )
  {
    if ( fDPhi == twopi || pt2 == 0 )  // on torus swept axis
    {
      in = kInside ;
    }
    else
    {
      // Try inner tolerant phi boundaries (=>inside)
      // if not inside, try outer tolerant phi boundaries
 
      pPhi = std::atan2(p.y(),p.x()) ;
 
      if ( pPhi < -halfAngTolerance )  { pPhi += twopi ; }  // 0<=pPhi<2pi
      if ( fSPhi >= 0 )
      {
        if ( (std::fabs(pPhi) < halfAngTolerance)
            && (std::fabs(fSPhi + fDPhi - twopi) < halfAngTolerance) )
        {
            pPhi += twopi ; // 0 <= pPhi < 2pi
        }
        if ( (pPhi >= fSPhi + halfAngTolerance)
            && (pPhi <= fSPhi + fDPhi - halfAngTolerance) )
        {
          in = kInside ;
        }
          else if ( (pPhi >= fSPhi - halfAngTolerance)
                 && (pPhi <= fSPhi + fDPhi + halfAngTolerance) )
        {
          in = kSurface ;
        }
      }
      else  // fSPhi < 0
      {
          if ( (pPhi <= fSPhi + twopi - halfAngTolerance)
            && (pPhi >= fSPhi + fDPhi  + halfAngTolerance) )  {;}
          else
          {
            in = kSurface ;
          }
      }
    }
  }
  else   // Try generous boundaries
  {
    tolRMin = fRmin - fRminTolerance ;
    tolRMax = fRmax + fRmaxTolerance ;
 
    if (tolRMin < 0 )  { tolRMin = 0 ; }
 
    if ( (pt2 >= tolRMin*tolRMin) && (pt2 <= tolRMax*tolRMax) )
    {
      if ( (fDPhi == twopi) || (pt2 == 0) ) // Continuous in phi or on z-axis
      {
        in = kSurface ;
      }
      else // Try outer tolerant phi boundaries only
      {
        pPhi = std::atan2(p.y(),p.x()) ;
 
        if ( pPhi < -halfAngTolerance )  { pPhi += twopi ; }  // 0<=pPhi<2pi
        if ( fSPhi >= 0 )
        {
          if ( (std::fabs(pPhi) < halfAngTolerance)
            && (std::fabs(fSPhi + fDPhi - twopi) < halfAngTolerance) )
          {
            pPhi += twopi ; // 0 <= pPhi < 2pi
          }
          if ( (pPhi >= fSPhi - halfAngTolerance)
            && (pPhi <= fSPhi + fDPhi + halfAngTolerance) )
          {
            in = kSurface;
          }
        }
        else  // fSPhi < 0
        {
          if ( (pPhi <= fSPhi + twopi - halfAngTolerance)
            && (pPhi >= fSPhi + fDPhi  + halfAngTolerance) )  {;}
          else
          {
            in = kSurface ;
          }
        }
      }
    }
  }
  return in ;
}

Referenced by DistanceToOut(), and SurfaceNormal().

◆ operator=()

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

Definition at line 188 of file G4Torus.cc.

{
   // Check assignment to self
   //
   if (this == &rhs)  { return *this; }
 
   // Copy base class data
   //
   G4CSGSolid::operator=(rhs);
 
   // Copy data
   //
   fRmin = rhs.fRmin; fRmax = rhs.fRmax;
   fRtor = rhs.fRtor; fSPhi = rhs.fSPhi; fDPhi = rhs.fDPhi;
   fRminTolerance = rhs.fRminTolerance; fRmaxTolerance = rhs.fRmaxTolerance;
   kRadTolerance = rhs.kRadTolerance; kAngTolerance = rhs.kAngTolerance;
   halfCarTolerance = rhs.halfCarTolerance;
   halfAngTolerance = rhs.halfAngTolerance;
 
   return *this;
}

◆ SetAllParameters()

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

Checks and sets all the parameters given in input. Used in constructor.

Definition at line 94 of file G4Torus.cc.

{
  const G4double fEpsilon = 4.e-11;  // relative tolerance of radii
 
  fCubicVolume = 0.;
  fSurfaceArea = 0.;
  fRebuildPolyhedron = true;
 
  kRadTolerance = G4GeometryTolerance::GetInstance()->GetRadialTolerance();
  kAngTolerance = G4GeometryTolerance::GetInstance()->GetAngularTolerance();
 
  halfCarTolerance = 0.5*kCarTolerance;
  halfAngTolerance = 0.5*kAngTolerance;
 
  if ( pRtor >= pRmax+1.e3*kCarTolerance )  // Check swept radius, as in G4Cons
  {
    fRtor = pRtor ;
  }
  else
  {
    std::ostringstream message;
    message << "Invalid swept radius for Solid: " << GetName() << G4endl
            << "        pRtor = " << pRtor << ", pRmax = " << pRmax;
    G4Exception("G4Torus::SetAllParameters()",
                "GeomSolids0002", FatalException, message);
  }
 
  // Check radii, as in G4Cons
  //
  if ( pRmin < pRmax - 1.e2*kCarTolerance && pRmin >= 0 )
  {
    if (pRmin >= 1.e2*kCarTolerance) { fRmin = pRmin ; }
    else                             { fRmin = 0.0   ; }
    fRmax = pRmax ;
  }
  else
  {
    std::ostringstream message;
    message << "Invalid values of radii for Solid: " << GetName() << G4endl
            << "        pRmin = " << pRmin << ", pRmax = " << pRmax;
    G4Exception("G4Torus::SetAllParameters()",
                "GeomSolids0002", FatalException, message);
  }
 
  // Relative tolerances
  //
  fRminTolerance = (fRmin) != 0.0
                 ? 0.5*std::max( kRadTolerance, fEpsilon*(fRtor-fRmin )) : 0;
  fRmaxTolerance = 0.5*std::max( kRadTolerance, fEpsilon*(fRtor+fRmax) );
 
  // Check angles
  //
  if ( pDPhi >= twopi )  { fDPhi = twopi ; }
  else
  {
    if (pDPhi > 0)       { fDPhi = pDPhi ; }
    else
    {
      std::ostringstream message;
      message << "Invalid Z delta-Phi for Solid: " << GetName() << G4endl
              << "        pDPhi = " << pDPhi;
      G4Exception("G4Torus::SetAllParameters()",
                  "GeomSolids0002", FatalException, message);
    }
  }
 
  // Ensure psphi in 0-2PI or -2PI-0 range if shape crosses 0
  //
  fSPhi = pSPhi;
 
  if (fSPhi < 0)  { fSPhi = twopi-std::fmod(std::fabs(fSPhi),twopi) ; }
  else            { fSPhi = std::fmod(fSPhi,twopi) ; }
 
  if (fSPhi+fDPhi > twopi)  { fSPhi-=twopi ; }
}

Referenced by G4Torus().

◆ StreamInfo()

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

overridevirtual

Streams the object contents to an output stream.

Reimplemented from G4CSGSolid.

Definition at line 1565 of file G4Torus.cc.

{
  G4long oldprc = os.precision(16);
  os << "-----------------------------------------------------------\n"
     << "    *** Dump for solid - " << GetName() << " ***\n"
     << "    ===================================================\n"
     << " Solid type: G4Torus\n"
     << " Parameters: \n"
     << "    inner radius: " << fRmin/mm << " mm \n"
     << "    outer radius: " << fRmax/mm << " mm \n"
     << "    swept radius: " << fRtor/mm << " mm \n"
     << "    starting phi: " << fSPhi/degree << " degrees \n"
     << "    delta phi   : " << fDPhi/degree << " degrees \n"
     << "-----------------------------------------------------------\n";
  os.precision(oldprc);
 
  return os;
}

◆ SurfaceNormal()

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

overridevirtual

Returns the outwards pointing unit normal of the shape for the surface closest to the point at offset 'p'.

Parameters

[in] p The point at offset p.

Returns: The outwards pointing unit normal.

Implements G4VSolid.

Definition at line 680 of file G4Torus.cc.

{
  G4int noSurfaces = 0;
  G4double rho, pt, pPhi;
  G4double distRMin = kInfinity;
  G4double distSPhi = kInfinity, distEPhi = kInfinity;
 
  // To cope with precision loss
  //
  const G4double delta = std::max(10.0*kCarTolerance,
                                  1.0e-8*(fRtor+fRmax));
  const G4double dAngle = 10.0*kAngTolerance;
 
  G4ThreeVector nR, nPs, nPe;
  G4ThreeVector norm, sumnorm(0.,0.,0.);
 
  rho = std::hypot(p.x(),p.y());
  pt  = std::hypot(p.z(),rho-fRtor);
 
  G4double  distRMax = std::fabs(pt - fRmax);
  if(fRmin != 0.0) { distRMin = std::fabs(pt - fRmin); }
 
  if( rho > delta && pt != 0.0 )
  {
    G4double redFactor= (rho-fRtor)/rho;
    nR = G4ThreeVector( p.x()*redFactor,  // p.x()*(1.-fRtor/rho),
                        p.y()*redFactor,  // p.y()*(1.-fRtor/rho),
                        p.z()          );
    nR *= 1.0/pt;
  }
 
  if ( fDPhi < twopi ) // && rho ) // old limitation against (0,0,z)
  {
    if ( rho != 0.0 )
    {
      pPhi = std::atan2(p.y(),p.x());
 
      if(pPhi < fSPhi-delta)            { pPhi += twopi; }
      else if(pPhi > fSPhi+fDPhi+delta) { pPhi -= twopi; }
 
      distSPhi = std::fabs( pPhi - fSPhi );
      distEPhi = std::fabs(pPhi-fSPhi-fDPhi);
    }
    nPs = G4ThreeVector(std::sin(fSPhi),-std::cos(fSPhi),0);
    nPe = G4ThreeVector(-std::sin(fSPhi+fDPhi),std::cos(fSPhi+fDPhi),0);
  }
  if( distRMax <= delta )
  {
    ++noSurfaces;
    sumnorm += nR;
  }
  else if( (fRmin != 0.0) && (distRMin <= delta) ) // Must not be on both Outer and Inner
  {
    ++noSurfaces;
    sumnorm -= nR;
  }
 
  //  To be on one of the 'phi' surfaces,
  //  it must be within the 'tube' - with tolerance
 
  if( (fDPhi < twopi) && (fRmin-delta <= pt) && (pt <= (fRmax+delta)) )
  {
    if (distSPhi <= dAngle)
    {
      ++noSurfaces;
      sumnorm += nPs;
    }
    if (distEPhi <= dAngle)
    {
      ++noSurfaces;
      sumnorm += nPe;
    }
  }
  if ( noSurfaces == 0 )
  {
#ifdef G4CSGDEBUG
     G4ExceptionDescription ed;
     ed.precision(16);
 
     EInside  inIt= Inside( p );
 
     if( inIt != kSurface )
     {
        ed << " ERROR>  Surface Normal was called for Torus,"
           << " with point not on surface." << G4endl;
     }
     else
     {
        ed << " ERROR>  Surface Normal has not found a surface, "
           << " despite the point being on the surface. " <<G4endl;
     }
 
     if( inIt != kInside)
     {
         ed << " Safety (Dist To In)  = " << DistanceToIn(p) << G4endl;
     }
     if( inIt != kOutside)
     {
         ed << " Safety (Dist to Out) = " << DistanceToOut(p) << G4endl;
     }
     ed << " Coordinates of point : " << p << G4endl;
     ed << " Parameters  of solid : " << G4endl << *this << G4endl;
 
     if( inIt == kSurface )
     {
        G4Exception("G4Torus::SurfaceNormal(p)", "GeomSolids1002",
                    JustWarning, ed,
                    "Failing to find normal, even though point is on surface!");
     }
     else
     {
        static const char* NameInside[3]= { "Inside", "Surface", "Outside" };
        ed << "  The point is " << NameInside[inIt] << " the solid. "<< G4endl;
        G4Exception("G4Torus::SurfaceNormal(p)", "GeomSolids1002",
                    JustWarning, ed, "Point p is not on surface !?" );
     }
#endif
     norm = ApproxSurfaceNormal(p);
  }
  else if ( noSurfaces == 1 )  { norm = sumnorm; }
  else                         { norm = sumnorm.unit(); }
 
  return norm ;
}

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

Public Member Functions

Additional Inherited Members

Detailed Description

Constructor & Destructor Documentation

◆ G4Torus() [1/3]

◆ ~G4Torus()

◆ G4Torus() [2/3]

◆ G4Torus() [3/3]

Member Function Documentation

◆ BoundingLimits()

◆ CalculateExtent()

◆ Clone()

◆ ComputeDimensions()

◆ CreatePolyhedron()

◆ DescribeYourselfTo()

◆ DistanceToIn() [1/2]

◆ DistanceToIn() [2/2]

◆ DistanceToOut() [1/2]

◆ DistanceToOut() [2/2]

◆ GetCosEndPhi()

◆ GetCosStartPhi()

◆ GetCubicVolume()

◆ GetDPhi()

◆ GetEntityType()

◆ GetPointOnSurface()

◆ GetRmax()

◆ GetRmin()

◆ GetRtor()

◆ GetSinEndPhi()

◆ GetSinStartPhi()

◆ GetSPhi()

◆ GetSurfaceArea()

◆ Inside()

◆ operator=()

◆ SetAllParameters()

◆ StreamInfo()

◆ SurfaceNormal()