65 using namespace CLHEP;
74 :
G4VSolid(pName), fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
82 std::ostringstream message;
84 <<
" Invalid Z half-length: "
85 << newHalfLenZ/
mm <<
" mm";
93 if (newInnerRadius<0 || newOuterRadius<0)
95 std::ostringstream message;
97 <<
" Invalid radii ! Inner radius: "
98 << newInnerRadius/
mm <<
" mm" <<
G4endl
100 << newOuterRadius/
mm <<
" mm";
104 if (newInnerRadius >= newOuterRadius)
106 std::ostringstream message;
108 <<
" Invalid radii ! Inner radius: "
109 << newInnerRadius/
mm <<
" mm" <<
G4endl
111 << newOuterRadius/
mm <<
" mm";
132 :
G4VSolid(a), innerRadius(0.), outerRadius(0.), halfLenZ(0.), innerStereo(0.),
133 outerStereo(0.), tanInnerStereo(0.), tanOuterStereo(0.), tanInnerStereo2(0.),
134 tanOuterStereo2(0.), innerRadius2(0.), outerRadius2(0.), endInnerRadius2(0.),
135 endOuterRadius2(0.), endInnerRadius(0.), endOuterRadius(0.),
136 fCubicVolume(0.), fSurfaceArea(0.), fHalfTol(0.), fpPolyhedron(0)
154 :
G4VSolid(rhs), innerRadius(rhs.innerRadius),
155 outerRadius(rhs.outerRadius), halfLenZ(rhs.halfLenZ),
156 innerStereo(rhs.innerStereo), outerStereo(rhs.outerStereo),
157 tanInnerStereo(rhs.tanInnerStereo), tanOuterStereo(rhs.tanOuterStereo),
158 tanInnerStereo2(rhs.tanInnerStereo2), tanOuterStereo2(rhs.tanOuterStereo2),
159 innerRadius2(rhs.innerRadius2), outerRadius2(rhs.outerRadius2),
160 endInnerRadius2(rhs.endInnerRadius2), endOuterRadius2(rhs.endOuterRadius2),
161 endInnerRadius(rhs.endInnerRadius), endOuterRadius(rhs.endOuterRadius),
162 fCubicVolume(rhs.fCubicVolume), fSurfaceArea(rhs.fSurfaceArea),
163 fHalfTol(rhs.fHalfTol), fpPolyhedron(0)
175 if (
this == &rhs) {
return *
this; }
191 fCubicVolume = rhs.fCubicVolume; fSurfaceArea = rhs.fSurfaceArea;
192 fHalfTol = rhs.fHalfTol; fpPolyhedron = 0;
226 G4double rFudge = 1.0/std::cos(0.5*sigPhi);
279 cosPhi = std::cos(phi),
280 sinPhi = std::sin(phi);
288 w0, w1, w2, w3, w4, w5, w6;
346 if (numPhi == 1) phi = 0;
347 cosPhi = std::cos(phi),
348 sinPhi = std::sin(phi);
439 if (splitOuter) v4 = w4;
441 }
while( --numPhi > 0 );
494 phiPoly.
SetNormal( (v1-v0).cross(w0-v0).unit() );
565 if (dist2Inner < dist2Z && dist2Inner < dist2Outer)
572 if (dist2Z < dist2Outer)
625 G4bool couldMissOuter(
true),
626 couldMissInner(
true),
627 cantMissInnerCylinder(
false);
634 if (sigz > -fHalfTol)
646 if (sigz > 0)
return kInfinity;
677 yi = p.
y() + q*v.
y();
705 cantMissInnerCylinder =
true;
710 return (sigz < fHalfTol) ? 0 : q;
753 if (pz < halfLenZ+fHalfTol)
783 if (zi < -halfLenZ)
continue;
784 if (zi > +halfLenZ && couldMissOuter)
continue;
790 yi = p.
y() + q[i]*v.
y();
808 if (cantMissInnerCylinder)
return (sigz < fHalfTol) ? 0 : -sigz/vz;
816 if (pz < halfLenZ+fHalfTol)
836 if (q[i] > best)
break;
846 if (zi < -halfLenZ)
continue;
847 if (zi > +halfLenZ && couldMissInner)
continue;
853 yi = p.
y() + q[i]*v.
y();
902 if (sigz > -fHalfTol)
907 return sigz < fHalfTol ? 0 : sigz;
915 G4double answer = std::sqrt( dr*dr + sigz*sigz );
916 return answer < fHalfTol ? 0 : answer;
924 return sigz < fHalfTol ? 0 : sigz;
936 G4double answer = std::sqrt( dr*dr + sigz*sigz );
937 return answer < fHalfTol ? 0 : answer;
949 return answer < fHalfTol ? 0 : answer;
957 return answer < fHalfTol ? 0 : answer;
1004 if (calcNorm) { *norm = *nBest; *validNorm =
true; }
1036 if (normHere.dot(v) > 0)
1038 if (calcNorm) { *norm = normHere.
unit(); *validNorm =
false; }
1047 for( i=0; i<
n; i++ )
1049 if (q[i] > sBest)
break;
1058 if (norm1.
dot(v) > 0)
1084 if (normHere.dot(v) > 0)
1088 *norm = normHere.
unit();
1099 for( i=0; i<
n; i++ )
1101 if (q[i] > sBest)
break;
1106 if (norm2.
dot(v) > 0)
1125 if (nBest == &norm1 || nBest == &norm2)
1126 *norm = nBest->
unit();
1151 if (tryOuter < sBest)
1157 if (tryInner < sBest) sBest = tryInner;
1209 G4double a = tx*tx + ty*ty - tz*tz*tan2Phi;
1219 if (std::fabs(b) <
DBL_MIN)
return 0;
1228 if (radical < -
DBL_MIN)
return 0;
1239 radical = std::sqrt(radical);
1241 G4double q = -0.5*( b + (b < 0 ? -radical : +radical) );
1244 if (sa < sb) { ss[0] = sa; ss[1] = sb; }
else { ss[0] = sb; ss[1] = sa; }
1273 if (tanPhi <
DBL_MIN)
return pr-r0;
1281 G4double r1 = std::sqrt( r0*r0 + z1*z1*tan2Phi );
1286 G4double z2 = (pr*tanPhi + pz)/(1 + tan2Phi);
1287 G4double r2 = std::sqrt( r0*r0 + z2*z2*tan2Phi );
1304 return std::sqrt( dr*dr + dz*dz );
1310 return std::fabs((pr-r1)*dz - (pz-z1)*dr)/
len;
1336 G4double rh = std::sqrt( r0*r0 + pz*pz*tan2Phi );
1345 return std::fabs((pr-rh)*dr)/
len;
1363 return new G4Hype(*
this);
1372 if(fCubicVolume != 0.) {;}
1374 return fCubicVolume;
1383 if(fSurfaceArea != 0.) {;}
1385 return fSurfaceArea;
1394 G4int oldprc = os.precision(16);
1395 os <<
"-----------------------------------------------------------\n"
1396 <<
" *** Dump for solid - " <<
GetName() <<
" ***\n"
1397 <<
" ===================================================\n"
1398 <<
" Solid type: G4Hype\n"
1399 <<
" Parameters: \n"
1400 <<
" half length Z: " <<
halfLenZ/
mm <<
" mm \n"
1405 <<
"-----------------------------------------------------------\n";
1406 os.precision(oldprc);
1418 G4double xRand, yRand, zRand, r2 , aOne, aTwo, aThree, chose, sinhu;
1419 G4double phi, cosphi, sinphi, rBar2Out, rBar2In, alpha, t, rOut, rIn2, rOut2;
1428 t = std::log(t+std::sqrt(
sqr(t)+1));
1429 aOne = std::fabs(2.*alpha*(std::sinh(2.*t)/4.+t/2.));
1435 t = std::log(t+std::sqrt(
sqr(t)+1));
1436 aTwo = std::fabs(2.*alpha*(std::sinh(2.*t)/4.+t/2.));
1445 cosphi = std::cos(phi);
1446 sinphi = std::sin(phi);
1451 if(chose>=0. && chose < aOne)
1466 else if(chose>=aOne && chose<aOne+aTwo)
1483 else if(chose>=aOne+aTwo && chose<aOne+aTwo+aThree)
1487 rOut = std::sqrt(rOut2) ;
1492 r2 = xRand*xRand + yRand*yRand ;
1493 }
while ( ! ( r2 >= rIn2 && r2 <= rOut2 ) ) ;
1502 rOut = std::sqrt(rOut2) ;
1507 r2 = xRand*xRand + yRand*yRand ;
1508 }
while ( ! ( r2 >= rIn2 && r2 <= rOut2 ) ) ;
1553 if (!fpPolyhedron ||
1557 delete fpPolyhedron;
1560 return fpPolyhedron;
1569 return std::log(arg+std::sqrt(
sqr(arg)+1));
G4double GetCubicVolume()
ThreeVector shoot(const G4int Ap, const G4int Af)
G4Hype(const G4String &pName, G4double newInnerRadius, G4double newOuterRadius, G4double newInnerStereo, G4double newOuterStereo, G4double newHalfLenZ)
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
CLHEP::Hep3Vector G4ThreeVector
G4double GetSurfaceArea()
double dot(const Hep3Vector &) const
G4GeometryType GetEntityType() const
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
static G4double ApproxDistInside(G4double pr, G4double pz, G4double r0, G4double tan2Phi)
G4bool GetExtent(G4double &min, G4double &max) const
virtual G4double GetCubicVolume()
void SetNormal(const G4ThreeVector &newNormal)
virtual G4bool PartialClip(const G4VoxelLimits &voxelLimit, const EAxis IgnoreMe)
G4Polyhedron * CreatePolyhedron() const
virtual void AddVertexInOrder(const G4ThreeVector vertex)
void DescribeYourselfTo(G4VGraphicsScene &scene) const
virtual void AddSolid(const G4Box &)=0
G4bool InnerSurfaceExists() const
G4double HypeOuterRadius2(G4double zVal) const
static G4double ApproxDistOutside(G4double pr, G4double pz, G4double r0, G4double tanPhi)
G4double HypeInnerRadius2(G4double zVal) const
void AddSurface(const G4ClippablePolygon &surface)
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=G4bool(false), G4bool *validNorm=0, G4ThreeVector *n=0) const
G4Hype & operator=(const G4Hype &rhs)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
G4ThreeVector GetPointOnSurface() const
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
static void AddPolyToExtent(const G4ThreeVector &v0, const G4ThreeVector &v1, const G4ThreeVector &w1, const G4ThreeVector &w0, const G4VoxelLimits &voxelLimit, const EAxis axis, G4SolidExtentList &extentList)
G4VisExtent GetExtent() const
static G4int IntersectHype(const G4ThreeVector &p, const G4ThreeVector &v, G4double r2, G4double tan2Phi, G4double s[2])
static G4int GetNumberOfRotationSteps()
T max(const T t1, const T t2)
brief Return the largest of the two arguments
EInside Inside(const G4ThreeVector &p) const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4VSolid & operator=(const G4VSolid &rhs)
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
void SetOuterStereo(G4double newOSte)
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pmin, G4double &pmax) const
virtual G4double GetSurfaceArea()
std::ostream & StreamInfo(std::ostream &os) const
G4Polyhedron * GetPolyhedron() const
void SetInnerStereo(G4double newISte)
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
const G4int kMaxMeshSections