Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4Hype Class Reference

#include <G4Hype.hh>

Inheritance diagram for G4Hype:
Collaboration diagram for G4Hype:

Public Member Functions

 G4Hype (const G4String &pName, G4double newInnerRadius, G4double newOuterRadius, G4double newInnerStereo, G4double newOuterStereo, G4double newHalfLenZ)
 
virtual ~G4Hype ()
 
void ComputeDimensions (G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
 
void Extent (G4ThreeVector &pMin, G4ThreeVector &pMax) const
 
G4bool CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
 
G4double GetInnerRadius () const
 
G4double GetOuterRadius () const
 
G4double GetZHalfLength () const
 
G4double GetInnerStereo () const
 
G4double GetOuterStereo () const
 
void SetInnerRadius (G4double newIRad)
 
void SetOuterRadius (G4double newORad)
 
void SetZHalfLength (G4double newHLZ)
 
void SetInnerStereo (G4double newISte)
 
void SetOuterStereo (G4double newOSte)
 
EInside Inside (const G4ThreeVector &p) const
 
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
 
G4VSolidClone () const
 
std::ostream & StreamInfo (std::ostream &os) const
 
G4double GetCubicVolume ()
 
G4double GetSurfaceArea ()
 
G4ThreeVector GetPointOnSurface () const
 
void DescribeYourselfTo (G4VGraphicsScene &scene) const
 
G4VisExtent GetExtent () const
 
G4PolyhedronCreatePolyhedron () const
 
G4PolyhedronGetPolyhedron () const
 
 G4Hype (__void__ &)
 
 G4Hype (const G4Hype &rhs)
 
G4Hypeoperator= (const G4Hype &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 const G4VSolidGetConstituentSolid (G4int no) const
 
virtual G4VSolidGetConstituentSolid (G4int no)
 
virtual const G4DisplacedSolidGetDisplacedSolidPtr () const
 
virtual G4DisplacedSolidGetDisplacedSolidPtr ()
 
 G4VSolid (__void__ &)
 
 G4VSolid (const G4VSolid &rhs)
 
G4VSolidoperator= (const G4VSolid &rhs)
 
G4double EstimateCubicVolume (G4int nStat, G4double epsilon) const
 
G4double EstimateSurfaceArea (G4int nStat, G4double ell) const
 

Protected Types

enum  ESide { outerFace, innerFace, leftCap, rightCap }
 

Protected Member Functions

G4bool InnerSurfaceExists () const
 
G4double HypeInnerRadius2 (G4double zVal) const
 
G4double HypeOuterRadius2 (G4double zVal) 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
 

Static Protected Member Functions

static G4double ApproxDistOutside (G4double pr, G4double pz, G4double r0, G4double tanPhi)
 
static G4double ApproxDistInside (G4double pr, G4double pz, G4double r0, G4double tan2Phi)
 
static G4int IntersectHype (const G4ThreeVector &p, const G4ThreeVector &v, G4double r2, G4double tan2Phi, G4double s[2])
 

Protected Attributes

G4double innerRadius
 
G4double outerRadius
 
G4double halfLenZ
 
G4double innerStereo
 
G4double outerStereo
 
G4double tanInnerStereo
 
G4double tanOuterStereo
 
G4double tanInnerStereo2
 
G4double tanOuterStereo2
 
G4double innerRadius2
 
G4double outerRadius2
 
G4double endInnerRadius2
 
G4double endOuterRadius2
 
G4double endInnerRadius
 
G4double endOuterRadius
 
- Protected Attributes inherited from G4VSolid
G4double kCarTolerance
 

Detailed Description

Definition at line 67 of file G4Hype.hh.

Member Enumeration Documentation

enum G4Hype::ESide
protected
Enumerator
outerFace 
innerFace 
leftCap 
rightCap 

Definition at line 183 of file G4Hype.hh.

Constructor & Destructor Documentation

G4Hype::G4Hype ( const G4String pName,
G4double  newInnerRadius,
G4double  newOuterRadius,
G4double  newInnerStereo,
G4double  newOuterStereo,
G4double  newHalfLenZ 
)

Definition at line 74 of file G4Hype.cc.

80  : G4VSolid(pName), fCubicVolume(0.), fSurfaceArea(0.),
81  fRebuildPolyhedron(false), fpPolyhedron(0)
82 {
83  fHalfTol = 0.5*kCarTolerance;
84 
85  // Check z-len
86  //
87  if (newHalfLenZ<=0)
88  {
89  std::ostringstream message;
90  message << "Invalid Z half-length - " << GetName() << G4endl
91  << " Invalid Z half-length: "
92  << newHalfLenZ/mm << " mm";
93  G4Exception("G4Hype::G4Hype()", "GeomSolids0002",
94  FatalErrorInArgument, message);
95  }
96  halfLenZ=newHalfLenZ;
97 
98  // Check radii
99  //
100  if (newInnerRadius<0 || newOuterRadius<0)
101  {
102  std::ostringstream message;
103  message << "Invalid radii - " << GetName() << G4endl
104  << " Invalid radii ! Inner radius: "
105  << newInnerRadius/mm << " mm" << G4endl
106  << " Outer radius: "
107  << newOuterRadius/mm << " mm";
108  G4Exception("G4Hype::G4Hype()", "GeomSolids0002",
109  FatalErrorInArgument, message);
110  }
111  if (newInnerRadius >= newOuterRadius)
112  {
113  std::ostringstream message;
114  message << "Outer > inner radius - " << GetName() << G4endl
115  << " Invalid radii ! Inner radius: "
116  << newInnerRadius/mm << " mm" << G4endl
117  << " Outer radius: "
118  << newOuterRadius/mm << " mm";
119  G4Exception("G4Hype::G4Hype()", "GeomSolids0002",
120  FatalErrorInArgument, message);
121  }
122 
123  innerRadius=newInnerRadius;
124  outerRadius=newOuterRadius;
125 
128 
129  SetInnerStereo( newInnerStereo );
130  SetOuterStereo( newOuterStereo );
131 }
G4String GetName() const
G4double outerRadius2
Definition: G4Hype.hh:175
static constexpr double mm
Definition: G4SIunits.hh:115
G4double innerRadius2
Definition: G4Hype.hh:174
G4double outerRadius
Definition: G4Hype.hh:163
G4double halfLenZ
Definition: G4Hype.hh:164
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
G4VSolid(const G4String &name)
Definition: G4VSolid.cc:61
G4double kCarTolerance
Definition: G4VSolid.hh:307
void SetOuterStereo(G4double newOSte)
G4double innerRadius
Definition: G4Hype.hh:162
void SetInnerStereo(G4double newISte)

Here is the call graph for this function:

Here is the caller graph for this function:

G4Hype::~G4Hype ( )
virtual

Definition at line 152 of file G4Hype.cc.

153 {
154  delete fpPolyhedron; fpPolyhedron = 0;
155 }
G4Hype::G4Hype ( __void__ &  a)

Definition at line 138 of file G4Hype.cc.

139  : G4VSolid(a), innerRadius(0.), outerRadius(0.), halfLenZ(0.), innerStereo(0.),
143  fCubicVolume(0.), fSurfaceArea(0.), fHalfTol(0.),
144  fRebuildPolyhedron(false), fpPolyhedron(0)
145 {
146 }
G4double outerStereo
Definition: G4Hype.hh:166
G4double outerRadius2
Definition: G4Hype.hh:175
G4double innerRadius2
Definition: G4Hype.hh:174
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4double innerStereo
Definition: G4Hype.hh:165
G4double endOuterRadius2
Definition: G4Hype.hh:177
G4double tanOuterStereo
Definition: G4Hype.hh:171
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
G4double tanInnerStereo
Definition: G4Hype.hh:170
G4VSolid(const G4String &name)
Definition: G4VSolid.cc:61
G4double endInnerRadius2
Definition: G4Hype.hh:176
G4double innerRadius
Definition: G4Hype.hh:162
G4double endInnerRadius
Definition: G4Hype.hh:178
G4Hype::G4Hype ( const G4Hype rhs)

Definition at line 161 of file G4Hype.cc.

162  : G4VSolid(rhs), innerRadius(rhs.innerRadius),
170  fCubicVolume(rhs.fCubicVolume), fSurfaceArea(rhs.fSurfaceArea),
171  fHalfTol(rhs.fHalfTol), fRebuildPolyhedron(false), fpPolyhedron(0)
172 {
173 }
G4double outerStereo
Definition: G4Hype.hh:166
G4double outerRadius2
Definition: G4Hype.hh:175
G4double innerRadius2
Definition: G4Hype.hh:174
G4double innerStereo
Definition: G4Hype.hh:165
G4double endOuterRadius2
Definition: G4Hype.hh:177
G4double tanOuterStereo
Definition: G4Hype.hh:171
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
G4double tanInnerStereo
Definition: G4Hype.hh:170
G4VSolid(const G4String &name)
Definition: G4VSolid.cc:61
G4double endInnerRadius2
Definition: G4Hype.hh:176
G4double innerRadius
Definition: G4Hype.hh:162
G4double endInnerRadius
Definition: G4Hype.hh:178

Member Function Documentation

G4double G4Hype::ApproxDistInside ( G4double  pr,
G4double  pz,
G4double  r0,
G4double  tan2Phi 
)
staticprotected

Definition at line 1091 of file G4Hype.cc.

1093 {
1094  if (tan2Phi < DBL_MIN) return r0 - pr;
1095 
1096  //
1097  // Corresponding position and normal on hyperbolic
1098  //
1099  G4double rh = std::sqrt( r0*r0 + pz*pz*tan2Phi );
1100 
1101  G4double dr = -rh;
1102  G4double dz = pz*tan2Phi;
1103  G4double len = std::sqrt(dr*dr + dz*dz);
1104 
1105  //
1106  // Answer
1107  //
1108  return std::fabs((pr-rh)*dr)/len;
1109 }
const XML_Char int len
Definition: expat.h:262
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4Hype::ApproxDistOutside ( G4double  pr,
G4double  pz,
G4double  r0,
G4double  tanPhi 
)
staticprotected

Definition at line 1033 of file G4Hype.cc.

1035 {
1036  if (tanPhi < DBL_MIN) return pr-r0;
1037 
1038  G4double tan2Phi = tanPhi*tanPhi;
1039 
1040  //
1041  // First point
1042  //
1043  G4double z1 = pz;
1044  G4double r1 = std::sqrt( r0*r0 + z1*z1*tan2Phi );
1045 
1046  //
1047  // Second point
1048  //
1049  G4double z2 = (pr*tanPhi + pz)/(1 + tan2Phi);
1050  G4double r2 = std::sqrt( r0*r0 + z2*z2*tan2Phi );
1051 
1052  //
1053  // Line between them
1054  //
1055  G4double dr = r2-r1;
1056  G4double dz = z2-z1;
1057 
1058  G4double len = std::sqrt(dr*dr + dz*dz);
1059  if (len < DBL_MIN)
1060  {
1061  //
1062  // The two points are the same?? I guess we
1063  // must have really bracketed the normal
1064  //
1065  dr = pr-r1;
1066  dz = pz-z1;
1067  return std::sqrt( dr*dr + dz*dz );
1068  }
1069 
1070  //
1071  // Distance
1072  //
1073  return std::fabs((pr-r1)*dz - (pz-z1)*dr)/len;
1074 }
const XML_Char int len
Definition: expat.h:262
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

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

Implements G4VSolid.

Definition at line 246 of file G4Hype.cc.

250 {
251  G4ThreeVector bmin, bmax;
252 
253  // Get bounding box
254  Extent(bmin,bmax);
255 
256  // Find extent
257  G4BoundingEnvelope bbox(bmin,bmax);
258  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
259 }
void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Hype.cc:223

Here is the call graph for this function:

G4VSolid * G4Hype::Clone ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1124 of file G4Hype.cc.

1125 {
1126  return new G4Hype(*this);
1127 }
G4Hype(const G4String &pName, G4double newInnerRadius, G4double newOuterRadius, G4double newInnerStereo, G4double newOuterStereo, G4double newHalfLenZ)
Definition: G4Hype.cc:74

Here is the call graph for this function:

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

Reimplemented from G4VSolid.

Definition at line 212 of file G4Hype.cc.

215 {
216  p->ComputeDimensions(*this,n,pRep);
217 }
const G4int n
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const

Here is the call graph for this function:

G4Polyhedron * G4Hype::CreatePolyhedron ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1306 of file G4Hype.cc.

1307 {
1310 }
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double halfLenZ
Definition: G4Hype.hh:164
G4double innerRadius
Definition: G4Hype.hh:162

Here is the caller graph for this function:

void G4Hype::DescribeYourselfTo ( G4VGraphicsScene scene) const
virtual

Implements G4VSolid.

Definition at line 1284 of file G4Hype.cc.

1285 {
1286  scene.AddSolid (*this);
1287 }
virtual void AddSolid(const G4Box &)=0

Here is the call graph for this function:

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

Implements G4VSolid.

Definition at line 355 of file G4Hype.cc.

357 {
358  //
359  // Quick test. Beware! This assumes v is a unit vector!
360  //
361  if (std::fabs(p.x()*v.y() - p.y()*v.x()) > endOuterRadius+kCarTolerance)
362  return kInfinity;
363 
364  //
365  // Take advantage of z symmetry, and reflect throught the
366  // z=0 plane so that pz is always positive
367  //
368  G4double pz(p.z()), vz(v.z());
369  if (pz < 0)
370  {
371  pz = -pz;
372  vz = -vz;
373  }
374 
375  //
376  // We must be very careful if we don't want to
377  // create subtle leaks at the edges where the
378  // hyperbolic surfaces connect to the endplate.
379  // The only reliable way to do so is to make sure
380  // that the decision as to when a track passes
381  // over the edge of one surface is exactly the
382  // same decision as to when a track passes into the
383  // other surface. By "exact", we don't mean algebraicly
384  // exact, but we mean the same machine instructions
385  // should be used.
386  //
387  G4bool couldMissOuter(true),
388  couldMissInner(true),
389  cantMissInnerCylinder(false);
390 
391  //
392  // Check endplate intersection
393  //
394  G4double sigz = pz-halfLenZ;
395 
396  if (sigz > -fHalfTol) // equivalent to: if (pz > halfLenZ - fHalfTol)
397  {
398  //
399  // We start in front of the endplate (within roundoff)
400  // Correct direction to intersect endplate?
401  //
402  if (vz >= 0)
403  {
404  //
405  // Nope. As long as we are far enough away, we
406  // can't intersect anything
407  //
408  if (sigz > 0) return kInfinity;
409 
410  //
411  // Otherwise, we may still hit a hyperbolic surface
412  // if the point is on the hyperbolic surface (within tolerance)
413  //
414  G4double pr2 = p.x()*p.x() + p.y()*p.y();
416  return kInfinity;
417 
418  if (InnerSurfaceExists())
419  {
421  return kInfinity;
422  if ( (pr2 < endOuterRadius2 - kCarTolerance*endOuterRadius)
423  && (pr2 > endInnerRadius2 + kCarTolerance*endInnerRadius) )
424  return kInfinity;
425  }
426  else
427  {
428  if (pr2 < endOuterRadius2 - kCarTolerance*endOuterRadius)
429  return kInfinity;
430  }
431  }
432  else
433  {
434  //
435  // Where do we intersect at z = halfLenZ?
436  //
437  G4double q = -sigz/vz;
438  G4double xi = p.x() + q*v.x(),
439  yi = p.y() + q*v.y();
440 
441  //
442  // Is this on the endplate? If so, return s, unless
443  // we are on the tolerant surface, in which case return 0
444  //
445  G4double pr2 = xi*xi + yi*yi;
446  if (pr2 <= endOuterRadius2)
447  {
448  if (InnerSurfaceExists())
449  {
450  if (pr2 >= endInnerRadius2) return (sigz < fHalfTol) ? 0 : q;
451  //
452  // This test is sufficient to ensure that the
453  // trajectory cannot miss the inner hyperbolic surface
454  // for z > 0, if the normal is correct.
455  //
456  G4double dot1 = (xi*v.x() + yi*v.y())*endInnerRadius/std::sqrt(pr2);
457  couldMissInner = (dot1 - halfLenZ*tanInnerStereo2*vz <= 0);
458 
459  if (pr2 > endInnerRadius2*(1 - 2*DBL_EPSILON) )
460  {
461  //
462  // There is a potential leak if the inner
463  // surface is a cylinder
464  //
465  if ( (innerStereo < DBL_MIN)
466  && ((std::fabs(v.x()) > DBL_MIN) || (std::fabs(v.y()) > DBL_MIN)))
467  cantMissInnerCylinder = true;
468  }
469  }
470  else
471  {
472  return (sigz < fHalfTol) ? 0 : q;
473  }
474  }
475  else
476  {
477  G4double dotR( xi*v.x() + yi*v.y() );
478  if (dotR >= 0)
479  {
480  //
481  // Otherwise, if we are traveling outwards, we know
482  // we must miss the hyperbolic surfaces also, so
483  // we need not bother checking
484  //
485  return kInfinity;
486  }
487  else
488  {
489  //
490  // This test is sufficient to ensure that the
491  // trajectory cannot miss the outer hyperbolic surface
492  // for z > 0, if the normal is correct.
493  //
494  G4double dot1 = dotR*endOuterRadius/std::sqrt(pr2);
495  couldMissOuter = (dot1 - halfLenZ*tanOuterStereo2*vz>= 0);
496  }
497  }
498  }
499  }
500 
501  //
502  // Check intersection with outer hyperbolic surface, save
503  // distance to valid intersection into "best".
504  //
505  G4double best = kInfinity;
506 
507  G4double q[2];
509 
510  if (n > 0)
511  {
512  //
513  // Potential intersection: is p on this surface?
514  //
515  if (pz < halfLenZ+fHalfTol)
516  {
517  G4double dr2 = p.x()*p.x() + p.y()*p.y() - HypeOuterRadius2(pz);
518  if (std::fabs(dr2) < kCarTolerance*endOuterRadius)
519  {
520  //
521  // Sure, but make sure we're traveling inwards at
522  // this point
523  //
524  if (p.x()*v.x() + p.y()*v.y() - pz*tanOuterStereo2*vz < 0)
525  return 0;
526  }
527  }
528 
529  //
530  // We are now certain that p is not on the tolerant surface.
531  // Accept only position distance q
532  //
533  G4int i;
534  for( i=0; i<n; i++ )
535  {
536  if (q[i] >= 0)
537  {
538  //
539  // Check to make sure this intersection point is
540  // on the surface, but only do so if we haven't
541  // checked the endplate intersection already
542  //
543  G4double zi = pz + q[i]*vz;
544 
545  if (zi < -halfLenZ) continue;
546  if (zi > +halfLenZ && couldMissOuter) continue;
547 
548  //
549  // Check normal
550  //
551  G4double xi = p.x() + q[i]*v.x(),
552  yi = p.y() + q[i]*v.y();
553 
554  if (xi*v.x() + yi*v.y() - zi*tanOuterStereo2*vz > 0) continue;
555 
556  best = q[i];
557  break;
558  }
559  }
560  }
561 
562  if (!InnerSurfaceExists()) return best;
563 
564  //
565  // Check intersection with inner hyperbolic surface
566  //
567  n = IntersectHype( p, v, innerRadius2, tanInnerStereo2, q );
568  if (n == 0)
569  {
570  if (cantMissInnerCylinder) return (sigz < fHalfTol) ? 0 : -sigz/vz;
571 
572  return best;
573  }
574 
575  //
576  // P on this surface?
577  //
578  if (pz < halfLenZ+fHalfTol)
579  {
580  G4double dr2 = p.x()*p.x() + p.y()*p.y() - HypeInnerRadius2(pz);
581  if (std::fabs(dr2) < kCarTolerance*endInnerRadius)
582  {
583  //
584  // Sure, but make sure we're traveling outwards at
585  // this point
586  //
587  if (p.x()*v.x() + p.y()*v.y() - pz*tanInnerStereo2*vz > 0) return 0;
588  }
589  }
590 
591  //
592  // No, so only positive q is valid. Search for a valid intersection
593  // that is closer than the outer intersection (if it exists)
594  //
595  G4int i;
596  for( i=0; i<n; i++ )
597  {
598  if (q[i] > best) break;
599  if (q[i] >= 0)
600  {
601  //
602  // Check to make sure this intersection point is
603  // on the surface, but only do so if we haven't
604  // checked the endplate intersection already
605  //
606  G4double zi = pz + q[i]*vz;
607 
608  if (zi < -halfLenZ) continue;
609  if (zi > +halfLenZ && couldMissInner) continue;
610 
611  //
612  // Check normal
613  //
614  G4double xi = p.x() + q[i]*v.x(),
615  yi = p.y() + q[i]*v.y();
616 
617  if (xi*v.x() + yi*v.y() - zi*tanOuterStereo2*vz < 0) continue;
618 
619  best = q[i];
620  break;
621  }
622  }
623 
624  //
625  // Done
626  //
627  return best;
628 }
G4double outerRadius2
Definition: G4Hype.hh:175
static const G4double kInfinity
Definition: geomdefs.hh:42
double x() const
G4double innerRadius2
Definition: G4Hype.hh:174
G4double innerStereo
Definition: G4Hype.hh:165
G4double endOuterRadius2
Definition: G4Hype.hh:177
int G4int
Definition: G4Types.hh:78
G4bool InnerSurfaceExists() const
G4double HypeOuterRadius2(G4double zVal) const
double z() const
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double HypeInnerRadius2(G4double zVal) const
G4double endOuterRadius
Definition: G4Hype.hh:179
bool G4bool
Definition: G4Types.hh:79
#define DBL_EPSILON
Definition: templates.hh:87
const G4int n
G4double halfLenZ
Definition: G4Hype.hh:164
static G4int IntersectHype(const G4ThreeVector &p, const G4ThreeVector &v, G4double r2, G4double tan2Phi, G4double s[2])
Definition: G4Hype.cc:965
double y() const
#define DBL_MIN
Definition: templates.hh:75
G4double endInnerRadius2
Definition: G4Hype.hh:176
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4double endInnerRadius
Definition: G4Hype.hh:178

Here is the call graph for this function:

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

Implements G4VSolid.

Definition at line 650 of file G4Hype.cc.

651 {
652  G4double absZ(std::fabs(p.z()));
653 
654  //
655  // Check region
656  //
657  G4double r2 = p.x()*p.x() + p.y()*p.y();
658  G4double r = std::sqrt(r2);
659 
660  G4double sigz = absZ - halfLenZ;
661 
662  if (r < endOuterRadius)
663  {
664  if (sigz > -fHalfTol)
665  {
666  if (InnerSurfaceExists())
667  {
668  if (r > endInnerRadius)
669  return sigz < fHalfTol ? 0 : sigz; // Region 1
670 
671  G4double dr = endInnerRadius - r;
672  if (sigz > dr*tanInnerStereo2)
673  {
674  //
675  // In region 5
676  //
677  G4double answer = std::sqrt( dr*dr + sigz*sigz );
678  return answer < fHalfTol ? 0 : answer;
679  }
680  }
681  else
682  {
683  //
684  // In region 1 (no inner surface)
685  //
686  return sigz < fHalfTol ? 0 : sigz;
687  }
688  }
689  }
690  else
691  {
692  G4double dr = r - endOuterRadius;
693  if (sigz > -dr*tanOuterStereo2)
694  {
695  //
696  // In region 2
697  //
698  G4double answer = std::sqrt( dr*dr + sigz*sigz );
699  return answer < fHalfTol ? 0 : answer;
700  }
701  }
702 
703  if (InnerSurfaceExists())
704  {
706  {
707  //
708  // In region 4
709  //
711  return answer < fHalfTol ? 0 : answer;
712  }
713  }
714 
715  //
716  // We are left by elimination with region 3
717  //
719  return answer < fHalfTol ? 0 : answer;
720 }
double x() const
static G4double ApproxDistInside(G4double pr, G4double pz, G4double r0, G4double tan2Phi)
Definition: G4Hype.cc:1091
G4double tanOuterStereo
Definition: G4Hype.hh:171
G4bool InnerSurfaceExists() const
double z() const
G4double tanOuterStereo2
Definition: G4Hype.hh:173
static G4double ApproxDistOutside(G4double pr, G4double pz, G4double r0, G4double tanPhi)
Definition: G4Hype.cc:1033
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double HypeInnerRadius2(G4double zVal) const
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
double y() const
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4double innerRadius
Definition: G4Hype.hh:162
G4double endInnerRadius
Definition: G4Hype.hh:178

Here is the call graph for this function:

G4double G4Hype::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 731 of file G4Hype.cc.

734 {
735  static const G4ThreeVector normEnd1(0.0,0.0,+1.0);
736  static const G4ThreeVector normEnd2(0.0,0.0,-1.0);
737 
738  //
739  // Keep track of closest surface
740  //
741  G4double sBest; // distance to
742  const G4ThreeVector *nBest; // normal vector
743  G4bool vBest; // whether "valid"
744 
745  //
746  // Check endplate, taking advantage of symmetry.
747  // Note that the endcap is the only surface which
748  // has a "valid" normal, i.e. is a surface of which
749  // the entire solid is behind.
750  //
751  G4double pz(p.z()), vz(v.z());
752  if (vz < 0)
753  {
754  pz = -pz;
755  vz = -vz;
756  nBest = &normEnd2;
757  }
758  else
759  nBest = &normEnd1;
760 
761  //
762  // Possible intercept. Are we on the surface?
763  //
764  if (pz > halfLenZ-fHalfTol)
765  {
766  if (calcNorm) { *norm = *nBest; *validNorm = true; }
767  return 0;
768  }
769 
770  //
771  // Nope. Get distance. Beware of zero vz.
772  //
773  sBest = (vz > DBL_MIN) ? (halfLenZ - pz)/vz : kInfinity;
774  vBest = true;
775 
776  //
777  // Check outer surface
778  //
779  G4double r2 = p.x()*p.x() + p.y()*p.y();
780 
781  G4double q[2];
783 
784  G4ThreeVector norm1, norm2;
785 
786  if (n > 0)
787  {
788  //
789  // We hit somewhere. Are we on the surface?
790  //
791  G4double dr2 = r2 - HypeOuterRadius2(pz);
792  if (std::fabs(dr2) < endOuterRadius*kCarTolerance)
793  {
794  G4ThreeVector normHere( p.x(), p.y(), -p.z()*tanOuterStereo2 );
795  //
796  // Sure. But are we going the right way?
797  //
798  if (normHere.dot(v) > 0)
799  {
800  if (calcNorm) { *norm = normHere.unit(); *validNorm = false; }
801  return 0;
802  }
803  }
804 
805  //
806  // Nope. Check closest positive intercept.
807  //
808  G4int i;
809  for( i=0; i<n; i++ )
810  {
811  if (q[i] > sBest) break;
812  if (q[i] > 0)
813  {
814  //
815  // Make sure normal is correct (that this
816  // solution is an outgoing solution)
817  //
818  G4ThreeVector pk(p+q[i]*v);
819  norm1 = G4ThreeVector( pk.x(), pk.y(), -pk.z()*tanOuterStereo2 );
820  if (norm1.dot(v) > 0)
821  {
822  sBest = q[i];
823  nBest = &norm1;
824  vBest = false;
825  break;
826  }
827  }
828  }
829  }
830 
831  if (InnerSurfaceExists())
832  {
833  //
834  // Check inner surface
835  //
836  n = IntersectHype( p, v, innerRadius2, tanInnerStereo2, q );
837  if (n > 0)
838  {
839  //
840  // On surface?
841  //
842  G4double dr2 = r2 - HypeInnerRadius2(pz);
843  if (std::fabs(dr2) < endInnerRadius*kCarTolerance)
844  {
845  G4ThreeVector normHere( -p.x(), -p.y(), p.z()*tanInnerStereo2 );
846  if (normHere.dot(v) > 0)
847  {
848  if (calcNorm)
849  {
850  *norm = normHere.unit();
851  *validNorm = false;
852  }
853  return 0;
854  }
855  }
856 
857  //
858  // Check closest positive
859  //
860  G4int i;
861  for( i=0; i<n; i++ )
862  {
863  if (q[i] > sBest) break;
864  if (q[i] > 0)
865  {
866  G4ThreeVector pk(p+q[i]*v);
867  norm2 = G4ThreeVector( -pk.x(), -pk.y(), pk.z()*tanInnerStereo2 );
868  if (norm2.dot(v) > 0)
869  {
870  sBest = q[i];
871  nBest = &norm2;
872  vBest = false;
873  break;
874  }
875  }
876  }
877  }
878  }
879 
880  //
881  // Done!
882  //
883  if (calcNorm)
884  {
885  *validNorm = vBest;
886 
887  if (nBest == &norm1 || nBest == &norm2)
888  *norm = nBest->unit();
889  else
890  *norm = *nBest;
891  }
892 
893  return sBest;
894 }
G4double outerRadius2
Definition: G4Hype.hh:175
static const G4double kInfinity
Definition: geomdefs.hh:42
CLHEP::Hep3Vector G4ThreeVector
double x() const
double dot(const Hep3Vector &) const
G4double innerRadius2
Definition: G4Hype.hh:174
int G4int
Definition: G4Types.hh:78
G4bool InnerSurfaceExists() const
G4double HypeOuterRadius2(G4double zVal) const
double z() const
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double HypeInnerRadius2(G4double zVal) const
G4double endOuterRadius
Definition: G4Hype.hh:179
bool G4bool
Definition: G4Types.hh:79
const G4int n
G4double halfLenZ
Definition: G4Hype.hh:164
static G4int IntersectHype(const G4ThreeVector &p, const G4ThreeVector &v, G4double r2, G4double tan2Phi, G4double s[2])
Definition: G4Hype.cc:965
Hep3Vector unit() const
double y() const
#define DBL_MIN
Definition: templates.hh:75
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4double endInnerRadius
Definition: G4Hype.hh:178

Here is the call graph for this function:

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

Implements G4VSolid.

Definition at line 902 of file G4Hype.cc.

903 {
904  //
905  // Try each surface and remember the closest
906  //
907  G4double absZ(std::fabs(p.z()));
908  G4double r(p.perp());
909 
910  G4double sBest = halfLenZ - absZ;
911 
912  G4double tryOuter = ApproxDistInside( r, absZ, outerRadius, tanOuterStereo2 );
913  if (tryOuter < sBest)
914  sBest = tryOuter;
915 
916  if (InnerSurfaceExists())
917  {
919  if (tryInner < sBest) sBest = tryInner;
920  }
921 
922  return sBest < 0.5*kCarTolerance ? 0 : sBest;
923 }
static G4double ApproxDistInside(G4double pr, G4double pz, G4double r0, G4double tan2Phi)
Definition: G4Hype.cc:1091
G4bool InnerSurfaceExists() const
double z() const
G4double tanOuterStereo2
Definition: G4Hype.hh:173
static G4double ApproxDistOutside(G4double pr, G4double pz, G4double r0, G4double tanPhi)
Definition: G4Hype.cc:1033
G4double outerRadius
Definition: G4Hype.hh:163
G4double halfLenZ
Definition: G4Hype.hh:164
G4double tanInnerStereo
Definition: G4Hype.hh:170
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4double innerRadius
Definition: G4Hype.hh:162
double perp() const

Here is the call graph for this function:

void G4Hype::Extent ( G4ThreeVector pMin,
G4ThreeVector pMax 
) const
virtual

Reimplemented from G4VSolid.

Definition at line 223 of file G4Hype.cc.

224 {
227 
228  // Check correctness of the bounding box
229  //
230  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
231  {
232  std::ostringstream message;
233  message << "Bad bounding box (min >= max) for solid: "
234  << GetName() << " !"
235  << "\npMin = " << pMin
236  << "\npMax = " << pMax;
237  G4Exception("G4Hype::Extent()", "GeomMgt0001", JustWarning, message);
238  DumpInfo();
239  }
240 }
void set(double x, double y, double z)
G4String GetName() const
double x() const
double z() const
void DumpInfo() const
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double y() const

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4Hype::GetCubicVolume ( )
virtual

Reimplemented from G4VSolid.

Definition at line 1133 of file G4Hype.cc.

1134 {
1135  if(fCubicVolume != 0.) {;}
1136  else { fCubicVolume = G4VSolid::GetCubicVolume(); }
1137  return fCubicVolume;
1138 }
virtual G4double GetCubicVolume()
Definition: G4VSolid.cc:189

Here is the call graph for this function:

G4GeometryType G4Hype::GetEntityType ( ) const
virtual

Implements G4VSolid.

Definition at line 1115 of file G4Hype.cc.

1116 {
1117  return G4String("G4Hype");
1118 }
G4VisExtent G4Hype::GetExtent ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1293 of file G4Hype.cc.

1294 {
1295  // Define the sides of the box into which the G4Tubs instance would fit.
1296  //
1299  -halfLenZ, halfLenZ );
1300 }
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
G4double G4Hype::GetInnerRadius ( ) const
inline

Here is the caller graph for this function:

G4double G4Hype::GetInnerStereo ( ) const
inline

Here is the caller graph for this function:

G4double G4Hype::GetOuterRadius ( ) const
inline

Here is the caller graph for this function:

G4double G4Hype::GetOuterStereo ( ) const
inline

Here is the caller graph for this function:

G4ThreeVector G4Hype::GetPointOnSurface ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1179 of file G4Hype.cc.

1180 {
1181  G4double xRand, yRand, zRand, r2 , aOne, aTwo, aThree, chose, sinhu;
1182  G4double phi, cosphi, sinphi, rBar2Out, rBar2In, alpha, t, rOut, rIn2, rOut2;
1183 
1184  // we use the formula of the area of a surface of revolution to compute
1185  // the areas, using the equation of the hyperbola:
1186  // x^2 + y^2 = (z*tanphi)^2 + r^2
1187 
1188  rBar2Out = outerRadius2;
1189  alpha = 2.*pi*rBar2Out*std::cos(outerStereo)/tanOuterStereo;
1191  t = std::log(t+std::sqrt(sqr(t)+1));
1192  aOne = std::fabs(2.*alpha*(std::sinh(2.*t)/4.+t/2.));
1193 
1194 
1195  rBar2In = innerRadius2;
1196  alpha = 2.*pi*rBar2In*std::cos(innerStereo)/tanInnerStereo;
1198  t = std::log(t+std::sqrt(sqr(t)+1));
1199  aTwo = std::fabs(2.*alpha*(std::sinh(2.*t)/4.+t/2.));
1200 
1203 
1204  if(outerStereo == 0.) {aOne = std::fabs(2.*pi*outerRadius*2.*halfLenZ);}
1205  if(innerStereo == 0.) {aTwo = std::fabs(2.*pi*innerRadius*2.*halfLenZ);}
1206 
1207  phi = G4RandFlat::shoot(0.,2.*pi);
1208  cosphi = std::cos(phi);
1209  sinphi = std::sin(phi);
1212 
1213  chose = G4RandFlat::shoot(0.,aOne+aTwo+2.*aThree);
1214  if(chose>=0. && chose < aOne)
1215  {
1216  if(outerStereo != 0.)
1217  {
1218  zRand = outerRadius*sinhu/tanOuterStereo;
1219  xRand = std::sqrt(sqr(sinhu)+1)*outerRadius*cosphi;
1220  yRand = std::sqrt(sqr(sinhu)+1)*outerRadius*sinphi;
1221  return G4ThreeVector (xRand, yRand, zRand);
1222  }
1223  else
1224  {
1225  return G4ThreeVector(outerRadius*cosphi,outerRadius*sinphi,
1227  }
1228  }
1229  else if(chose>=aOne && chose<aOne+aTwo)
1230  {
1231  if(innerStereo != 0.)
1232  {
1235  zRand = innerRadius*sinhu/tanInnerStereo;
1236  xRand = std::sqrt(sqr(sinhu)+1)*innerRadius*cosphi;
1237  yRand = std::sqrt(sqr(sinhu)+1)*innerRadius*sinphi;
1238  return G4ThreeVector (xRand, yRand, zRand);
1239  }
1240  else
1241  {
1242  return G4ThreeVector(innerRadius*cosphi,innerRadius*sinphi,
1244  }
1245  }
1246  else if(chose>=aOne+aTwo && chose<aOne+aTwo+aThree)
1247  {
1249  rOut2 = outerRadius2+tanOuterStereo2*halfLenZ*halfLenZ;
1250  rOut = std::sqrt(rOut2) ;
1251 
1252  do // Loop checking, 13.08.2015, G.Cosmo
1253  {
1254  xRand = G4RandFlat::shoot(-rOut,rOut) ;
1255  yRand = G4RandFlat::shoot(-rOut,rOut) ;
1256  r2 = xRand*xRand + yRand*yRand ;
1257  } while ( ! ( r2 >= rIn2 && r2 <= rOut2 ) ) ;
1258 
1259  zRand = halfLenZ;
1260  return G4ThreeVector (xRand, yRand, zRand);
1261  }
1262  else
1263  {
1264  rIn2 = innerRadius2+tanInnerStereo2*halfLenZ*halfLenZ;
1265  rOut2 = outerRadius2+tanOuterStereo2*halfLenZ*halfLenZ;
1266  rOut = std::sqrt(rOut2) ;
1267 
1268  do // Loop checking, 13.08.2015, G.Cosmo
1269  {
1270  xRand = G4RandFlat::shoot(-rOut,rOut) ;
1271  yRand = G4RandFlat::shoot(-rOut,rOut) ;
1272  r2 = xRand*xRand + yRand*yRand ;
1273  } while ( ! ( r2 >= rIn2 && r2 <= rOut2 ) ) ;
1274 
1275  zRand = -1.*halfLenZ;
1276  return G4ThreeVector (xRand, yRand, zRand);
1277  }
1278 }
G4double outerStereo
Definition: G4Hype.hh:166
G4double outerRadius2
Definition: G4Hype.hh:175
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4double innerRadius2
Definition: G4Hype.hh:174
G4double innerStereo
Definition: G4Hype.hh:165
G4double tanOuterStereo
Definition: G4Hype.hh:171
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double halfLenZ
Definition: G4Hype.hh:164
G4double tanInnerStereo
Definition: G4Hype.hh:170
static constexpr double pi
Definition: G4SIunits.hh:75
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
G4double innerRadius
Definition: G4Hype.hh:162
static const G4double alpha

Here is the call graph for this function:

G4Polyhedron * G4Hype::GetPolyhedron ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 1316 of file G4Hype.cc.

1317 {
1318  if (!fpPolyhedron ||
1319  fRebuildPolyhedron ||
1320  fpPolyhedron->GetNumberOfRotationStepsAtTimeOfCreation() !=
1321  fpPolyhedron->GetNumberOfRotationSteps())
1322  {
1323  G4AutoLock l(&polyhedronMutex);
1324  delete fpPolyhedron;
1325  fpPolyhedron = CreatePolyhedron();
1326  fRebuildPolyhedron = false;
1327  l.unlock();
1328  }
1329  return fpPolyhedron;
1330 }
G4Polyhedron * CreatePolyhedron() const
Definition: G4Hype.cc:1306
static G4int GetNumberOfRotationSteps()
G4int GetNumberOfRotationStepsAtTimeOfCreation() const

Here is the call graph for this function:

G4double G4Hype::GetSurfaceArea ( )
virtual

Reimplemented from G4VSolid.

Definition at line 1144 of file G4Hype.cc.

1145 {
1146  if(fSurfaceArea != 0.) {;}
1147  else { fSurfaceArea = G4VSolid::GetSurfaceArea(); }
1148  return fSurfaceArea;
1149 }
virtual G4double GetSurfaceArea()
Definition: G4VSolid.cc:251

Here is the call graph for this function:

G4double G4Hype::GetZHalfLength ( ) const
inline

Here is the caller graph for this function:

G4double G4Hype::HypeInnerRadius2 ( G4double  zVal) const
inlineprotected

Here is the caller graph for this function:

G4double G4Hype::HypeOuterRadius2 ( G4double  zVal) const
inlineprotected

Here is the caller graph for this function:

G4bool G4Hype::InnerSurfaceExists ( ) const
inlineprotected

Here is the caller graph for this function:

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

Implements G4VSolid.

Definition at line 265 of file G4Hype.cc.

266 {
267  //
268  // Check z extents: are we outside?
269  //
270  const G4double absZ(std::fabs(p.z()));
271  if (absZ > halfLenZ + fHalfTol) return kOutside;
272 
273  //
274  // Check outer radius
275  //
276  const G4double oRad2(HypeOuterRadius2(absZ));
277  const G4double xR2( p.x()*p.x()+p.y()*p.y() );
278 
279  if (xR2 > oRad2 + kCarTolerance*endOuterRadius) return kOutside;
280 
281  if (xR2 > oRad2 - kCarTolerance*endOuterRadius) return kSurface;
282 
283  if (InnerSurfaceExists())
284  {
285  //
286  // Check inner radius
287  //
288  const G4double iRad2(HypeInnerRadius2(absZ));
289 
290  if (xR2 < iRad2 - kCarTolerance*endInnerRadius) return kOutside;
291 
292  if (xR2 < iRad2 + kCarTolerance*endInnerRadius) return kSurface;
293  }
294 
295  //
296  // We are inside in radius, now check endplate surface
297  //
298  if (absZ > halfLenZ - fHalfTol) return kSurface;
299 
300  return kInside;
301 }
double x() const
G4bool InnerSurfaceExists() const
G4double HypeOuterRadius2(G4double zVal) const
double z() const
G4double HypeInnerRadius2(G4double zVal) const
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
double y() const
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4double endInnerRadius
Definition: G4Hype.hh:178

Here is the call graph for this function:

G4int G4Hype::IntersectHype ( const G4ThreeVector p,
const G4ThreeVector v,
G4double  r2,
G4double  tan2Phi,
G4double  s[2] 
)
staticprotected

Definition at line 965 of file G4Hype.cc.

967 {
968  G4double x0 = p.x(), y0 = p.y(), z0 = p.z();
969  G4double tx = v.x(), ty = v.y(), tz = v.z();
970 
971  G4double a = tx*tx + ty*ty - tz*tz*tan2Phi;
972  G4double b = 2*( x0*tx + y0*ty - z0*tz*tan2Phi );
973  G4double c = x0*x0 + y0*y0 - r2 - z0*z0*tan2Phi;
974 
975  if (std::fabs(a) < DBL_MIN)
976  {
977  //
978  // The trajectory is parallel to the asympotic limit of
979  // the surface: single solution
980  //
981  if (std::fabs(b) < DBL_MIN) return 0;
982  // Unless we travel through exact center
983 
984  ss[0] = c/b;
985  return 1;
986  }
987 
988 
989  G4double radical = b*b - 4*a*c;
990 
991  if (radical < -DBL_MIN) return 0; // No solution
992 
993  if (radical < DBL_MIN)
994  {
995  //
996  // Grazes surface
997  //
998  ss[0] = -b/a/2.0;
999  return 1;
1000  }
1001 
1002  radical = std::sqrt(radical);
1003 
1004  G4double q = -0.5*( b + (b < 0 ? -radical : +radical) );
1005  G4double sa = q/a;
1006  G4double sb = c/q;
1007  if (sa < sb) { ss[0] = sa; ss[1] = sb; } else { ss[0] = sb; ss[1] = sa; }
1008  return 2;
1009 }
double x() const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
double z() const
tuple b
Definition: test.py:12
double y() const
#define DBL_MIN
Definition: templates.hh:75
double G4double
Definition: G4Types.hh:76
tuple c
Definition: test.py:13

Here is the call graph for this function:

Here is the caller graph for this function:

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

Definition at line 179 of file G4Hype.cc.

180 {
181  // Check assignment to self
182  //
183  if (this == &rhs) { return *this; }
184 
185  // Copy base class data
186  //
187  G4VSolid::operator=(rhs);
188 
189  // Copy data
190  //
192  halfLenZ = rhs.halfLenZ;
199  fCubicVolume = rhs.fCubicVolume; fSurfaceArea = rhs.fSurfaceArea;
200  fHalfTol = rhs.fHalfTol;
201  fRebuildPolyhedron = false;
202  delete fpPolyhedron; fpPolyhedron = 0;
203 
204  return *this;
205 }
G4double outerStereo
Definition: G4Hype.hh:166
G4double outerRadius2
Definition: G4Hype.hh:175
G4double innerRadius2
Definition: G4Hype.hh:174
G4double innerStereo
Definition: G4Hype.hh:165
G4double endOuterRadius2
Definition: G4Hype.hh:177
G4double tanOuterStereo
Definition: G4Hype.hh:171
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double outerRadius
Definition: G4Hype.hh:163
G4double endOuterRadius
Definition: G4Hype.hh:179
G4double halfLenZ
Definition: G4Hype.hh:164
G4double tanInnerStereo
Definition: G4Hype.hh:170
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:111
G4double endInnerRadius2
Definition: G4Hype.hh:176
G4double innerRadius
Definition: G4Hype.hh:162
G4double endInnerRadius
Definition: G4Hype.hh:178

Here is the call graph for this function:

void G4Hype::SetInnerRadius ( G4double  newIRad)
inline

Here is the caller graph for this function:

void G4Hype::SetInnerStereo ( G4double  newISte)
inline

Here is the caller graph for this function:

void G4Hype::SetOuterRadius ( G4double  newORad)
inline

Here is the caller graph for this function:

void G4Hype::SetOuterStereo ( G4double  newOSte)
inline

Here is the caller graph for this function:

void G4Hype::SetZHalfLength ( G4double  newHLZ)
inline

Here is the caller graph for this function:

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

Implements G4VSolid.

Definition at line 1155 of file G4Hype.cc.

1156 {
1157  G4int oldprc = os.precision(16);
1158  os << "-----------------------------------------------------------\n"
1159  << " *** Dump for solid - " << GetName() << " ***\n"
1160  << " ===================================================\n"
1161  << " Solid type: G4Hype\n"
1162  << " Parameters: \n"
1163  << " half length Z: " << halfLenZ/mm << " mm \n"
1164  << " inner radius : " << innerRadius/mm << " mm \n"
1165  << " outer radius : " << outerRadius/mm << " mm \n"
1166  << " inner stereo angle : " << innerStereo/degree << " degrees \n"
1167  << " outer stereo angle : " << outerStereo/degree << " degrees \n"
1168  << "-----------------------------------------------------------\n";
1169  os.precision(oldprc);
1170 
1171  return os;
1172 }
G4double outerStereo
Definition: G4Hype.hh:166
G4String GetName() const
static constexpr double mm
Definition: G4SIunits.hh:115
G4double innerStereo
Definition: G4Hype.hh:165
int G4int
Definition: G4Types.hh:78
G4double outerRadius
Definition: G4Hype.hh:163
static constexpr double degree
Definition: G4SIunits.hh:144
G4double halfLenZ
Definition: G4Hype.hh:164
G4double innerRadius
Definition: G4Hype.hh:162

Here is the call graph for this function:

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

Implements G4VSolid.

Definition at line 309 of file G4Hype.cc.

310 {
311  //
312  // Which of the three or four surfaces are we closest to?
313  //
314  const G4double absZ(std::fabs(p.z()));
315  const G4double distZ(absZ - halfLenZ);
316  const G4double dist2Z(distZ*distZ);
317 
318  const G4double xR2( p.x()*p.x()+p.y()*p.y() );
319  const G4double dist2Outer( std::fabs(xR2 - HypeOuterRadius2(absZ)) );
320 
321  if (InnerSurfaceExists())
322  {
323  //
324  // Has inner surface: is this closest?
325  //
326  const G4double dist2Inner( std::fabs(xR2 - HypeInnerRadius2(absZ)) );
327  if (dist2Inner < dist2Z && dist2Inner < dist2Outer)
328  return G4ThreeVector( -p.x(), -p.y(), p.z()*tanInnerStereo2 ).unit();
329  }
330 
331  //
332  // Do the "endcaps" win?
333  //
334  if (dist2Z < dist2Outer)
335  return G4ThreeVector( 0.0, 0.0, p.z() < 0 ? -1.0 : 1.0 );
336 
337 
338  //
339  // Outer surface wins
340  //
341  return G4ThreeVector( p.x(), p.y(), -p.z()*tanOuterStereo2 ).unit();
342 }
CLHEP::Hep3Vector G4ThreeVector
double x() const
G4bool InnerSurfaceExists() const
G4double HypeOuterRadius2(G4double zVal) const
double z() const
G4double tanOuterStereo2
Definition: G4Hype.hh:173
G4double tanInnerStereo2
Definition: G4Hype.hh:172
G4double HypeInnerRadius2(G4double zVal) const
G4double halfLenZ
Definition: G4Hype.hh:164
double y() const
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Member Data Documentation

G4double G4Hype::endInnerRadius
protected

Definition at line 178 of file G4Hype.hh.

G4double G4Hype::endInnerRadius2
protected

Definition at line 176 of file G4Hype.hh.

G4double G4Hype::endOuterRadius
protected

Definition at line 179 of file G4Hype.hh.

G4double G4Hype::endOuterRadius2
protected

Definition at line 177 of file G4Hype.hh.

G4double G4Hype::halfLenZ
protected

Definition at line 164 of file G4Hype.hh.

G4double G4Hype::innerRadius
protected

Definition at line 162 of file G4Hype.hh.

G4double G4Hype::innerRadius2
protected

Definition at line 174 of file G4Hype.hh.

G4double G4Hype::innerStereo
protected

Definition at line 165 of file G4Hype.hh.

G4double G4Hype::outerRadius
protected

Definition at line 163 of file G4Hype.hh.

G4double G4Hype::outerRadius2
protected

Definition at line 175 of file G4Hype.hh.

G4double G4Hype::outerStereo
protected

Definition at line 166 of file G4Hype.hh.

G4double G4Hype::tanInnerStereo
protected

Definition at line 170 of file G4Hype.hh.

G4double G4Hype::tanInnerStereo2
protected

Definition at line 172 of file G4Hype.hh.

G4double G4Hype::tanOuterStereo
protected

Definition at line 171 of file G4Hype.hh.

G4double G4Hype::tanOuterStereo2
protected

Definition at line 173 of file G4Hype.hh.


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