Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4Polycone Class Reference

#include <G4Polycone.hh>

Inheritance diagram for G4Polycone:
Collaboration diagram for G4Polycone:

Public Member Functions

 G4Polycone (const G4String &name, G4double phiStart, G4double phiTotal, G4int numZPlanes, const G4double zPlane[], const G4double rInner[], const G4double rOuter[])
 
 G4Polycone (const G4String &name, G4double phiStart, G4double phiTotal, G4int numRZ, const G4double r[], const G4double z[])
 
virtual ~G4Polycone ()
 
EInside Inside (const G4ThreeVector &p) const
 
G4double DistanceToIn (const G4ThreeVector &p, const G4ThreeVector &v) const
 
G4double DistanceToIn (const G4ThreeVector &p) const
 
void Extent (G4ThreeVector &pMin, G4ThreeVector &pMax) const
 
G4bool CalculateExtent (const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pmin, G4double &pmax) const
 
G4ThreeVector GetPointOnSurface () const
 
void ComputeDimensions (G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
 
G4GeometryType GetEntityType () const
 
G4VSolidClone () const
 
std::ostream & StreamInfo (std::ostream &os) const
 
G4PolyhedronCreatePolyhedron () const
 
G4bool Reset ()
 
G4double GetStartPhi () const
 
G4double GetEndPhi () const
 
G4double GetSinStartPhi () const
 
G4double GetCosStartPhi () const
 
G4double GetSinEndPhi () const
 
G4double GetCosEndPhi () const
 
G4bool IsOpen () const
 
G4int GetNumRZCorner () const
 
G4PolyconeSideRZ GetCorner (G4int index) const
 
G4PolyconeHistoricalGetOriginalParameters () const
 
void SetOriginalParameters (G4PolyconeHistorical *pars)
 
 G4Polycone (__void__ &)
 
 G4Polycone (const G4Polycone &source)
 
G4Polyconeoperator= (const G4Polycone &source)
 
- Public Member Functions inherited from G4VCSGfaceted
 G4VCSGfaceted (const G4String &name)
 
virtual ~G4VCSGfaceted ()
 
 G4VCSGfaceted (const G4VCSGfaceted &source)
 
G4VCSGfacetedoperator= (const G4VCSGfaceted &source)
 
virtual G4ThreeVector SurfaceNormal (const G4ThreeVector &p) const
 
virtual G4double DistanceToOut (const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
 
virtual G4double DistanceToOut (const G4ThreeVector &p) const
 
virtual void DescribeYourselfTo (G4VGraphicsScene &scene) const
 
virtual G4VisExtent GetExtent () const
 
virtual G4PolyhedronGetPolyhedron () const
 
G4int GetCubVolStatistics () const
 
G4double GetCubVolEpsilon () const
 
void SetCubVolStatistics (G4int st)
 
void SetCubVolEpsilon (G4double ep)
 
G4int GetAreaStatistics () const
 
G4double GetAreaAccuracy () const
 
void SetAreaStatistics (G4int st)
 
void SetAreaAccuracy (G4double ep)
 
virtual G4double GetCubicVolume ()
 
virtual G4double GetSurfaceArea ()
 
 G4VCSGfaceted (__void__ &)
 
- 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 Member Functions

G4bool SetOriginalParameters (G4ReduciblePolygon *rz)
 
void Create (G4double phiStart, G4double phiTotal, G4ReduciblePolygon *rz)
 
void CopyStuff (const G4Polycone &source)
 
G4ThreeVector GetPointOnCone (G4double fRmin1, G4double fRmax1, G4double fRmin2, G4double fRmax2, G4double zOne, G4double zTwo, G4double &totArea) const
 
G4ThreeVector GetPointOnTubs (G4double fRMin, G4double fRMax, G4double zOne, G4double zTwo, G4double &totArea) const
 
G4ThreeVector GetPointOnCut (G4double fRMin1, G4double fRMax1, G4double fRMin2, G4double fRMax2, G4double zOne, G4double zTwo, G4double &totArea) const
 
G4ThreeVector GetPointOnRing (G4double fRMin, G4double fRMax, G4double fRMin2, G4double fRMax2, G4double zOne) const
 
- Protected Member Functions inherited from G4VCSGfaceted
virtual G4double DistanceTo (const G4ThreeVector &p, const G4bool outgoing) const
 
G4ThreeVector GetPointOnSurfaceGeneric () const
 
void CopyStuff (const G4VCSGfaceted &source)
 
void DeleteStuff ()
 
- 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

G4double startPhi
 
G4double endPhi
 
G4bool phiIsOpen
 
G4int numCorner
 
G4PolyconeSideRZcorners
 
G4PolyconeHistoricaloriginal_parameters
 
G4EnclosingCylinderenclosingCylinder
 
- Protected Attributes inherited from G4VCSGfaceted
G4int numFace
 
G4VCSGface ** faces
 
G4double fCubicVolume
 
G4double fSurfaceArea
 
G4bool fRebuildPolyhedron
 
G4PolyhedronfpPolyhedron
 
- Protected Attributes inherited from G4VSolid
G4double kCarTolerance
 

Detailed Description

Definition at line 82 of file G4Polycone.hh.

Constructor & Destructor Documentation

G4Polycone::G4Polycone ( const G4String name,
G4double  phiStart,
G4double  phiTotal,
G4int  numZPlanes,
const G4double  zPlane[],
const G4double  rInner[],
const G4double  rOuter[] 
)

Definition at line 63 of file G4Polycone.cc.

70  : G4VCSGfaceted( name )
71 {
72  //
73  // Some historical ugliness
74  //
76 
77  original_parameters->Start_angle = phiStart;
79  original_parameters->Num_z_planes = numZPlanes;
80  original_parameters->Z_values = new G4double[numZPlanes];
81  original_parameters->Rmin = new G4double[numZPlanes];
82  original_parameters->Rmax = new G4double[numZPlanes];
83 
84  G4int i;
85  for (i=0; i<numZPlanes; i++)
86  {
87  if(rInner[i]>rOuter[i])
88  {
89  DumpInfo();
90  std::ostringstream message;
91  message << "Cannot create a Polycone with rInner > rOuter for the same Z"
92  << G4endl
93  << " rInner > rOuter for the same Z !" << G4endl
94  << " rMin[" << i << "] = " << rInner[i]
95  << " -- rMax[" << i << "] = " << rOuter[i];
96  G4Exception("G4Polycone::G4Polycone()", "GeomSolids0002",
97  FatalErrorInArgument, message);
98  }
99  if (( i < numZPlanes-1) && ( zPlane[i] == zPlane[i+1] ))
100  {
101  if( (rInner[i] > rOuter[i+1])
102  ||(rInner[i+1] > rOuter[i]) )
103  {
104  DumpInfo();
105  std::ostringstream message;
106  message << "Cannot create a Polycone with no contiguous segments."
107  << G4endl
108  << " Segments are not contiguous !" << G4endl
109  << " rMin[" << i << "] = " << rInner[i]
110  << " -- rMax[" << i+1 << "] = " << rOuter[i+1] << G4endl
111  << " rMin[" << i+1 << "] = " << rInner[i+1]
112  << " -- rMax[" << i << "] = " << rOuter[i];
113  G4Exception("G4Polycone::G4Polycone()", "GeomSolids0002",
114  FatalErrorInArgument, message);
115  }
116  }
117  original_parameters->Z_values[i] = zPlane[i];
118  original_parameters->Rmin[i] = rInner[i];
119  original_parameters->Rmax[i] = rOuter[i];
120  }
121 
122  //
123  // Build RZ polygon using special PCON/PGON GEANT3 constructor
124  //
125  G4ReduciblePolygon *rz =
126  new G4ReduciblePolygon( rInner, rOuter, zPlane, numZPlanes );
127 
128  //
129  // Do the real work
130  //
131  Create( phiStart, phiTotal, rz );
132 
133  delete rz;
134 }
G4VCSGfaceted(const G4String &name)
int G4int
Definition: G4Types.hh:78
void DumpInfo() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
void Create(G4double phiStart, G4double phiTotal, G4ReduciblePolygon *rz)
Definition: G4Polycone.cc:181
double G4double
Definition: G4Types.hh:76
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

Here is the caller graph for this function:

G4Polycone::G4Polycone ( const G4String name,
G4double  phiStart,
G4double  phiTotal,
G4int  numRZ,
const G4double  r[],
const G4double  z[] 
)

Definition at line 140 of file G4Polycone.cc.

146  : G4VCSGfaceted( name )
147 {
148 
149  G4ReduciblePolygon *rz = new G4ReduciblePolygon( r, z, numRZ );
150 
151  Create( phiStart, phiTotal, rz );
152 
153  // Set original_parameters struct for consistency
154  //
155 
156  G4bool convertible=SetOriginalParameters(rz);
157 
158  if(!convertible)
159  {
160  std::ostringstream message;
161  message << "Polycone " << GetName() << "cannot be converted" << G4endl
162  << "to Polycone with (Rmin,Rmaz,Z) parameters!";
163  G4Exception("G4Polycone::G4Polycone()", "GeomSolids0002",
164  FatalException, message, "Use G4GenericPolycone instead!");
165  }
166  else
167  {
168  G4cout << "INFO: Converting polycone " << GetName() << G4endl
169  << "to optimized polycone with (Rmin,Rmaz,Z) parameters !"
170  << G4endl;
171  }
172  delete rz;
173 }
G4String GetName() const
G4VCSGfaceted(const G4String &name)
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
void SetOriginalParameters(G4PolyconeHistorical *pars)
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
void Create(G4double phiStart, G4double phiTotal, G4ReduciblePolygon *rz)
Definition: G4Polycone.cc:181

Here is the call graph for this function:

G4Polycone::~G4Polycone ( )
virtual

Definition at line 362 of file G4Polycone.cc.

363 {
364  delete [] corners;
365  delete original_parameters;
366  delete enclosingCylinder;
367 }
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198
G4Polycone::G4Polycone ( __void__ &  a)

Definition at line 351 of file G4Polycone.cc.

352  : G4VCSGfaceted(a), startPhi(0.), endPhi(0.), phiIsOpen(false),
355 {
356 }
G4double endPhi
Definition: G4Polycone.hh:194
G4VCSGfaceted(const G4String &name)
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
G4double startPhi
Definition: G4Polycone.hh:193
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4bool phiIsOpen
Definition: G4Polycone.hh:195
G4int numCorner
Definition: G4Polycone.hh:196
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198
G4Polycone::G4Polycone ( const G4Polycone source)

Definition at line 373 of file G4Polycone.cc.

374  : G4VCSGfaceted( source )
375 {
376  CopyStuff( source );
377 }
G4VCSGfaceted(const G4String &name)
void CopyStuff(const G4Polycone &source)
Definition: G4Polycone.cc:403

Here is the call graph for this function:

Member Function Documentation

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

Reimplemented from G4VCSGfaceted.

Definition at line 574 of file G4Polycone.cc.

578 {
579  G4ThreeVector bmin, bmax;
580  G4bool exist;
581 
582  // Check bounding box (bbox)
583  //
584  Extent(bmin,bmax);
585  G4BoundingEnvelope bbox(bmin,bmax);
586 #ifdef G4BBOX_EXTENT
587  if (true) return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
588 #endif
589  if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax))
590  {
591  return exist = (pMin < pMax) ? true : false;
592  }
593 
594  // To find the extent, RZ contour of the polycone is subdivided
595  // in triangles. The extent is calculated as cumulative extent of
596  // all sub-polycones formed by rotation of triangles around Z
597  //
598  G4TwoVectorList contourRZ;
599  G4TwoVectorList triangles;
600  std::vector<G4int> iout;
601  G4double eminlim = pVoxelLimit.GetMinExtent(pAxis);
602  G4double emaxlim = pVoxelLimit.GetMaxExtent(pAxis);
603 
604  // get RZ contour, ensure anticlockwise order of corners
605  for (G4int i=0; i<GetNumRZCorner(); ++i)
606  {
607  G4PolyconeSideRZ corner = GetCorner(i);
608  contourRZ.push_back(G4TwoVector(corner.r,corner.z));
609  }
611  G4double area = G4GeomTools::PolygonArea(contourRZ);
612  if (area < 0.) std::reverse(contourRZ.begin(),contourRZ.end());
613 
614  // triangulate RZ countour
615  if (!G4GeomTools::TriangulatePolygon(contourRZ,triangles))
616  {
617  std::ostringstream message;
618  message << "Triangulation of RZ contour has failed for solid: "
619  << GetName() << " !"
620  << "\nExtent has been calculated using boundary box";
621  G4Exception("G4Polycone::CalculateExtent()",
622  "GeomMgt1002", JustWarning, message);
623  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
624  }
625 
626  // set trigonometric values
627  const G4int NSTEPS = 24; // number of steps for whole circle
628  G4double astep = (360/NSTEPS)*deg; // max angle for one step
629 
630  G4double sphi = GetStartPhi();
631  G4double ephi = GetEndPhi();
632  G4double dphi = IsOpen() ? ephi-sphi : twopi;
633  G4int ksteps = (dphi <= astep) ? 1 : (G4int)((dphi-deg)/astep) + 1;
634  G4double ang = dphi/ksteps;
635 
636  G4double sinHalf = std::sin(0.5*ang);
637  G4double cosHalf = std::cos(0.5*ang);
638  G4double sinStep = 2.*sinHalf*cosHalf;
639  G4double cosStep = 1. - 2.*sinHalf*sinHalf;
640 
641  G4double sinStart = GetSinStartPhi();
642  G4double cosStart = GetCosStartPhi();
643  G4double sinEnd = GetSinEndPhi();
644  G4double cosEnd = GetCosEndPhi();
645 
646  // define vectors and arrays
647  std::vector<const G4ThreeVectorList *> polygons;
648  polygons.resize(ksteps+2);
649  G4ThreeVectorList pols[NSTEPS+2];
650  for (G4int k=0; k<ksteps+2; ++k) pols[k].resize(6);
651  for (G4int k=0; k<ksteps+2; ++k) polygons[k] = &pols[k];
652  G4double r0[6],z0[6]; // contour with original edges of triangle
653  G4double r1[6]; // shifted radii of external edges of triangle
654 
655  // main loop along triangles
656  pMin = kInfinity;
657  pMax =-kInfinity;
658  G4int ntria = triangles.size()/3;
659  for (G4int i=0; i<ntria; ++i)
660  {
661  G4int i3 = i*3;
662  for (G4int k=0; k<3; ++k)
663  {
664  G4int e0 = i3+k, e1 = (k<2) ? e0+1 : i3;
665  G4int k2 = k*2;
666  // set contour with original edges of triangle
667  r0[k2+0] = triangles[e0].x(); z0[k2+0] = triangles[e0].y();
668  r0[k2+1] = triangles[e1].x(); z0[k2+1] = triangles[e1].y();
669  // set shifted radii
670  r1[k2+0] = r0[k2+0];
671  r1[k2+1] = r0[k2+1];
672  if (z0[k2+1] - z0[k2+0] <= 0) continue;
673  r1[k2+0] /= cosHalf;
674  r1[k2+1] /= cosHalf;
675  }
676 
677  // rotate countour, set sequence of 6-sided polygons
678  G4double sinCur = sinStart*cosHalf + cosStart*sinHalf;
679  G4double cosCur = cosStart*cosHalf - sinStart*sinHalf;
680  for (G4int j=0; j<6; ++j) pols[0][j].set(r0[j]*cosStart,r0[j]*sinStart,z0[j]);
681  for (G4int k=1; k<ksteps+1; ++k)
682  {
683  for (G4int j=0; j<6; ++j) pols[k][j].set(r1[j]*cosCur,r1[j]*sinCur,z0[j]);
684  G4double sinTmp = sinCur;
685  sinCur = sinCur*cosStep + cosCur*sinStep;
686  cosCur = cosCur*cosStep - sinTmp*sinStep;
687  }
688  for (G4int j=0; j<6; ++j) pols[ksteps+1][j].set(r0[j]*cosEnd,r0[j]*sinEnd,z0[j]);
689 
690  // set sub-envelope and adjust extent
691  G4double emin,emax;
692  G4BoundingEnvelope benv(polygons);
693  if (!benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,emin,emax)) continue;
694  if (emin < pMin) pMin = emin;
695  if (emax > pMax) pMax = emax;
696  if (eminlim > pMin && emaxlim < pMax) return true; // max possible extent
697  }
698  return (pMin < pMax);
699 }
G4String GetName() const
static const G4double kInfinity
Definition: geomdefs.hh:42
G4bool IsOpen() const
std::vector< G4TwoVector > G4TwoVectorList
Definition: G4GeomTools.hh:50
G4int GetNumRZCorner() const
G4double GetSinStartPhi() const
static G4double PolygonArea(const G4TwoVectorList &polygon)
Definition: G4GeomTools.cc:82
G4double GetEndPhi() const
int G4int
Definition: G4Types.hh:78
static constexpr double twopi
Definition: G4SIunits.hh:76
bool G4bool
Definition: G4Types.hh:79
void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Polycone.cc:526
static G4bool TriangulatePolygon(const G4TwoVectorList &polygon, G4TwoVectorList &result)
Definition: G4GeomTools.cc:178
G4double GetCosStartPhi() const
G4double GetStartPhi() const
std::vector< G4ThreeVector > G4ThreeVectorList
G4double GetCosEndPhi() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static const G4double emax
G4double GetSinEndPhi() const
CLHEP::Hep2Vector G4TwoVector
Definition: G4TwoVector.hh:42
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
static constexpr double deg
Definition: G4SIunits.hh:152
G4double GetMaxExtent(const EAxis pAxis) const
G4PolyconeSideRZ GetCorner(G4int index) const
G4double GetMinExtent(const EAxis pAxis) const
static void RemoveRedundantVertices(G4TwoVectorList &polygon, std::vector< G4int > &iout, G4double tolerance=0)
Definition: G4GeomTools.cc:293

Here is the call graph for this function:

G4VSolid * G4Polycone::Clone ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 722 of file G4Polycone.cc.

723 {
724  return new G4Polycone(*this);
725 }
G4Polycone(const G4String &name, G4double phiStart, G4double phiTotal, G4int numZPlanes, const G4double zPlane[], const G4double rInner[], const G4double rOuter[])
Definition: G4Polycone.cc:63

Here is the call graph for this function:

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

Reimplemented from G4VSolid.

Definition at line 704 of file G4Polycone.cc.

707 {
708  p->ComputeDimensions(*this,n,pRep);
709 }
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const

Here is the call graph for this function:

void G4Polycone::CopyStuff ( const G4Polycone source)
protected

Definition at line 403 of file G4Polycone.cc.

404 {
405  //
406  // Simple stuff
407  //
408  startPhi = source.startPhi;
409  endPhi = source.endPhi;
410  phiIsOpen = source.phiIsOpen;
411  numCorner = source.numCorner;
412 
413  //
414  // The corner array
415  //
417 
418  G4PolyconeSideRZ *corn = corners,
419  *sourceCorn = source.corners;
420  do // Loop checking, 13.08.2015, G.Cosmo
421  {
422  *corn = *sourceCorn;
423  } while( ++sourceCorn, ++corn < corners+numCorner );
424 
425  //
426  // Original parameters
427  //
428  if (source.original_parameters)
429  {
432  }
433 
434  //
435  // Enclosing cylinder
436  //
438 
439  fRebuildPolyhedron = false;
440  fpPolyhedron = 0;
441 }
G4bool fRebuildPolyhedron
G4Polyhedron * fpPolyhedron
G4double endPhi
Definition: G4Polycone.hh:194
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
G4double startPhi
Definition: G4Polycone.hh:193
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4bool phiIsOpen
Definition: G4Polycone.hh:195
G4int numCorner
Definition: G4Polycone.hh:196
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the caller graph for this function:

void G4Polycone::Create ( G4double  phiStart,
G4double  phiTotal,
G4ReduciblePolygon rz 
)
protected

Definition at line 181 of file G4Polycone.cc.

184 {
185  //
186  // Perform checks of rz values
187  //
188  if (rz->Amin() < 0.0)
189  {
190  std::ostringstream message;
191  message << "Illegal input parameters - " << GetName() << G4endl
192  << " All R values must be >= 0 !";
193  G4Exception("G4Polycone::Create()", "GeomSolids0002",
194  FatalErrorInArgument, message);
195  }
196 
197  G4double rzArea = rz->Area();
198  if (rzArea < -kCarTolerance)
199  {
200  rz->ReverseOrder();
201  }
202  else if (rzArea < kCarTolerance)
203  {
204  std::ostringstream message;
205  message << "Illegal input parameters - " << GetName() << G4endl
206  << " R/Z cross section is zero or near zero: " << rzArea;
207  G4Exception("G4Polycone::Create()", "GeomSolids0002",
208  FatalErrorInArgument, message);
209  }
210 
212  || (!rz->RemoveRedundantVertices( kCarTolerance )) )
213  {
214  std::ostringstream message;
215  message << "Illegal input parameters - " << GetName() << G4endl
216  << " Too few unique R/Z values !";
217  G4Exception("G4Polycone::Create()", "GeomSolids0002",
218  FatalErrorInArgument, message);
219  }
220 
221  if (rz->CrossesItself(1/kInfinity))
222  {
223  std::ostringstream message;
224  message << "Illegal input parameters - " << GetName() << G4endl
225  << " R/Z segments cross !";
226  G4Exception("G4Polycone::Create()", "GeomSolids0002",
227  FatalErrorInArgument, message);
228  }
229 
230  numCorner = rz->NumVertices();
231 
232  //
233  // Phi opening? Account for some possible roundoff, and interpret
234  // nonsense value as representing no phi opening
235  //
236  if (phiTotal <= 0 || phiTotal > twopi-1E-10)
237  {
238  phiIsOpen = false;
239  startPhi = 0;
240  endPhi = twopi;
241  }
242  else
243  {
244  phiIsOpen = true;
245 
246  //
247  // Convert phi into our convention
248  //
249  startPhi = phiStart;
250  while( startPhi < 0 ) // Loop checking, 13.08.2015, G.Cosmo
251  startPhi += twopi;
252 
253  endPhi = phiStart+phiTotal;
254  while( endPhi < startPhi ) // Loop checking, 13.08.2015, G.Cosmo
255  endPhi += twopi;
256  }
257 
258  //
259  // Allocate corner array.
260  //
262 
263  //
264  // Copy corners
265  //
266  G4ReduciblePolygonIterator iterRZ(rz);
267 
268  G4PolyconeSideRZ *next = corners;
269  iterRZ.Begin();
270  do // Loop checking, 13.08.2015, G.Cosmo
271  {
272  next->r = iterRZ.GetA();
273  next->z = iterRZ.GetB();
274  } while( ++next, iterRZ.Next() );
275 
276  //
277  // Allocate face pointer array
278  //
280  faces = new G4VCSGface*[numFace];
281 
282  //
283  // Construct conical faces
284  //
285  // But! Don't construct a face if both points are at zero radius!
286  //
287  G4PolyconeSideRZ *corner = corners,
288  *prev = corners + numCorner-1,
289  *nextNext;
290  G4VCSGface **face = faces;
291  do // Loop checking, 13.08.2015, G.Cosmo
292  {
293  next = corner+1;
294  if (next >= corners+numCorner) next = corners;
295  nextNext = next+1;
296  if (nextNext >= corners+numCorner) nextNext = corners;
297 
298  if (corner->r < 1/kInfinity && next->r < 1/kInfinity) continue;
299 
300  //
301  // We must decide here if we can dare declare one of our faces
302  // as having a "valid" normal (i.e. allBehind = true). This
303  // is never possible if the face faces "inward" in r.
304  //
305  G4bool allBehind;
306  if (corner->z > next->z)
307  {
308  allBehind = false;
309  }
310  else
311  {
312  //
313  // Otherwise, it is only true if the line passing
314  // through the two points of the segment do not
315  // split the r/z cross section
316  //
317  allBehind = !rz->BisectedBy( corner->r, corner->z,
318  next->r, next->z, kCarTolerance );
319  }
320 
321  *face++ = new G4PolyconeSide( prev, corner, next, nextNext,
322  startPhi, endPhi-startPhi, phiIsOpen, allBehind );
323  } while( prev=corner, corner=next, corner > corners );
324 
325  if (phiIsOpen)
326  {
327  //
328  // Construct phi open edges
329  //
330  *face++ = new G4PolyPhiFace( rz, startPhi, 0, endPhi );
331  *face++ = new G4PolyPhiFace( rz, endPhi, 0, startPhi );
332  }
333 
334  //
335  // We might have dropped a face or two: recalculate numFace
336  //
337  numFace = face-faces;
338 
339  //
340  // Make enclosingCylinder
341  //
343  new G4EnclosingCylinder( rz, phiIsOpen, phiStart, phiTotal );
344 }
G4String GetName() const
G4bool CrossesItself(G4double tolerance)
static const G4double kInfinity
Definition: geomdefs.hh:42
G4double Amin() const
G4double endPhi
Definition: G4Polycone.hh:194
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
static constexpr double twopi
Definition: G4SIunits.hh:76
G4bool RemoveDuplicateVertices(G4double tolerance)
G4bool RemoveRedundantVertices(G4double tolerance)
G4int NumVertices() const
G4double startPhi
Definition: G4Polycone.hh:193
bool G4bool
Definition: G4Types.hh:79
G4VCSGface ** faces
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4bool BisectedBy(G4double a1, G4double b1, G4double a2, G4double b2, G4double tolerance)
#define G4endl
Definition: G4ios.hh:61
G4double kCarTolerance
Definition: G4VSolid.hh:307
G4bool phiIsOpen
Definition: G4Polycone.hh:195
G4int numCorner
Definition: G4Polycone.hh:196
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4Polyhedron * G4Polycone::CreatePolyhedron ( ) const
virtual

Implements G4VCSGfaceted.

Definition at line 1086 of file G4Polycone.cc.

1087 {
1088  //
1089  // This has to be fixed in visualization. Fake it for the moment.
1090  //
1091 
1098 }
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198
G4double G4Polycone::DistanceToIn ( const G4ThreeVector p,
const G4ThreeVector v 
) const
virtual

Reimplemented from G4VCSGfaceted.

Definition at line 498 of file G4Polycone.cc.

500 {
501  //
502  // Quick test
503  //
504  if (enclosingCylinder->ShouldMiss(p,v))
505  return kInfinity;
506 
507  //
508  // Long answer
509  //
510  return G4VCSGfaceted::DistanceToIn( p, v );
511 }
static const G4double kInfinity
Definition: geomdefs.hh:42
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4bool ShouldMiss(const G4ThreeVector &p, const G4ThreeVector &v) const

Here is the call graph for this function:

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

Reimplemented from G4VCSGfaceted.

Definition at line 517 of file G4Polycone.cc.

518 {
519  return G4VCSGfaceted::DistanceToIn(p);
520 }
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const

Here is the call graph for this function:

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

Reimplemented from G4VSolid.

Definition at line 526 of file G4Polycone.cc.

527 {
528  G4double rmin = kInfinity, rmax = -kInfinity;
529  G4double zmin = kInfinity, zmax = -kInfinity;
530 
531  for (G4int i=0; i<GetNumRZCorner(); ++i)
532  {
533  G4PolyconeSideRZ corner = GetCorner(i);
534  if (corner.r < rmin) rmin = corner.r;
535  if (corner.r > rmax) rmax = corner.r;
536  if (corner.z < zmin) zmin = corner.z;
537  if (corner.z > zmax) zmax = corner.z;
538  }
539 
540  if (IsOpen())
541  {
542  G4TwoVector vmin,vmax;
543  G4GeomTools::DiskExtent(rmin,rmax,
546  vmin,vmax);
547  pMin.set(vmin.x(),vmin.y(),zmin);
548  pMax.set(vmax.x(),vmax.y(),zmax);
549  }
550  else
551  {
552  pMin.set(-rmax,-rmax, zmin);
553  pMax.set( rmax, rmax, zmax);
554  }
555 
556  // Check correctness of the bounding box
557  //
558  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
559  {
560  std::ostringstream message;
561  message << "Bad bounding box (min >= max) for solid: "
562  << GetName() << " !"
563  << "\npMin = " << pMin
564  << "\npMax = " << pMax;
565  G4Exception("G4Polycone::Extent()", "GeomMgt0001", JustWarning, message);
566  DumpInfo();
567  }
568 }
void set(double x, double y, double z)
G4String GetName() const
double y() const
double x() const
static const G4double kInfinity
Definition: geomdefs.hh:42
G4bool IsOpen() const
double x() const
G4int GetNumRZCorner() const
G4double GetSinStartPhi() const
int G4int
Definition: G4Types.hh:78
double z() const
void DumpInfo() const
G4double GetCosStartPhi() const
G4double GetCosEndPhi() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
double y() const
G4double GetSinEndPhi() const
double G4double
Definition: G4Types.hh:76
G4PolyconeSideRZ GetCorner(G4int index) const
static G4bool DiskExtent(G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
Definition: G4GeomTools.cc:378

Here is the call graph for this function:

Here is the caller graph for this function:

G4PolyconeSideRZ G4Polycone::GetCorner ( G4int  index) const
inline

Here is the caller graph for this function:

G4double G4Polycone::GetCosEndPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Polycone::GetCosStartPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Polycone::GetEndPhi ( ) const
inline

Here is the caller graph for this function:

G4GeometryType G4Polycone::GetEntityType ( ) const
virtual

Reimplemented from G4VCSGfaceted.

Definition at line 714 of file G4Polycone.cc.

715 {
716  return G4String("G4Polycone");
717 }
G4int G4Polycone::GetNumRZCorner ( ) const
inline

Here is the caller graph for this function:

G4PolyconeHistorical* G4Polycone::GetOriginalParameters ( ) const
inline

Here is the caller graph for this function:

G4ThreeVector G4Polycone::GetPointOnCone ( G4double  fRmin1,
G4double  fRmax1,
G4double  fRmin2,
G4double  fRmax2,
G4double  zOne,
G4double  zTwo,
G4double totArea 
) const
protected

Definition at line 782 of file G4Polycone.cc.

786 {
787  // declare working variables
788  //
789  G4double Aone, Atwo, Afive, phi, zRand, fDPhi, cosu, sinu;
790  G4double rRand1, rmin, rmax, chose, rone, rtwo, qone, qtwo;
791  G4double fDz=(zTwo-zOne)/2., afDz=std::fabs(fDz);
792  G4ThreeVector point, offset=G4ThreeVector(0.,0.,0.5*(zTwo+zOne));
793  fDPhi = endPhi - startPhi;
794  rone = (fRmax1-fRmax2)/(2.*fDz);
795  rtwo = (fRmin1-fRmin2)/(2.*fDz);
796  if(fRmax1==fRmax2){qone=0.;}
797  else
798  {
799  qone = fDz*(fRmax1+fRmax2)/(fRmax1-fRmax2);
800  }
801  if(fRmin1==fRmin2){qtwo=0.;}
802  else
803  {
804  qtwo = fDz*(fRmin1+fRmin2)/(fRmin1-fRmin2);
805  }
806  Aone = 0.5*fDPhi*(fRmax2 + fRmax1)*(sqr(fRmin1-fRmin2)+sqr(zTwo-zOne));
807  Atwo = 0.5*fDPhi*(fRmin2 + fRmin1)*(sqr(fRmax1-fRmax2)+sqr(zTwo-zOne));
808  Afive = fDz*(fRmax1-fRmin1+fRmax2-fRmin2);
809  totArea = Aone+Atwo+2.*Afive;
810 
811  phi = G4RandFlat::shoot(startPhi,endPhi);
812  cosu = std::cos(phi);
813  sinu = std::sin(phi);
814 
815 
816  if( (startPhi == 0) && (endPhi == twopi) ) { Afive = 0; }
817  chose = G4RandFlat::shoot(0.,Aone+Atwo+2.*Afive);
818  if( (chose >= 0) && (chose < Aone) )
819  {
820  if(fRmax1 != fRmax2)
821  {
822  zRand = G4RandFlat::shoot(-1.*afDz,afDz);
823  point = G4ThreeVector (rone*cosu*(qone-zRand),
824  rone*sinu*(qone-zRand), zRand);
825  }
826  else
827  {
828  point = G4ThreeVector(fRmax1*cosu, fRmax1*sinu,
829  G4RandFlat::shoot(-1.*afDz,afDz));
830 
831  }
832  }
833  else if(chose >= Aone && chose < Aone + Atwo)
834  {
835  if(fRmin1 != fRmin2)
836  {
837  zRand = G4RandFlat::shoot(-1.*afDz,afDz);
838  point = G4ThreeVector (rtwo*cosu*(qtwo-zRand),
839  rtwo*sinu*(qtwo-zRand), zRand);
840 
841  }
842  else
843  {
844  point = G4ThreeVector(fRmin1*cosu, fRmin1*sinu,
845  G4RandFlat::shoot(-1.*afDz,afDz));
846  }
847  }
848  else if( (chose >= Aone + Atwo + Afive) && (chose < Aone + Atwo + 2.*Afive) )
849  {
850  zRand = G4RandFlat::shoot(-1.*afDz,afDz);
851  rmin = fRmin2-((zRand-fDz)/(2.*fDz))*(fRmin1-fRmin2);
852  rmax = fRmax2-((zRand-fDz)/(2.*fDz))*(fRmax1-fRmax2);
853  rRand1 = std::sqrt(G4RandFlat::shoot()*(sqr(rmax)-sqr(rmin))+sqr(rmin));
854  point = G4ThreeVector (rRand1*std::cos(startPhi),
855  rRand1*std::sin(startPhi), zRand);
856  }
857  else
858  {
859  zRand = G4RandFlat::shoot(-1.*afDz,afDz);
860  rmin = fRmin2-((zRand-fDz)/(2.*fDz))*(fRmin1-fRmin2);
861  rmax = fRmax2-((zRand-fDz)/(2.*fDz))*(fRmax1-fRmax2);
862  rRand1 = std::sqrt(G4RandFlat::shoot()*(sqr(rmax)-sqr(rmin))+sqr(rmin));
863  point = G4ThreeVector (rRand1*std::cos(endPhi),
864  rRand1*std::sin(endPhi), zRand);
865 
866  }
867 
868  return point+offset;
869 }
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4double endPhi
Definition: G4Polycone.hh:194
static constexpr double twopi
Definition: G4SIunits.hh:76
G4double startPhi
Definition: G4Polycone.hh:193
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4ThreeVector G4Polycone::GetPointOnCut ( G4double  fRMin1,
G4double  fRMax1,
G4double  fRMin2,
G4double  fRMax2,
G4double  zOne,
G4double  zTwo,
G4double totArea 
) const
protected

Definition at line 979 of file G4Polycone.cc.

983 { if(zOne==zTwo)
984  {
985  return GetPointOnRing(fRMin1, fRMax1,fRMin2,fRMax2,zOne);
986  }
987  if( (fRMin1 == fRMin2) && (fRMax1 == fRMax2) )
988  {
989  return GetPointOnTubs(fRMin1, fRMax1,zOne,zTwo,totArea);
990  }
991  return GetPointOnCone(fRMin1,fRMax1,fRMin2,fRMax2,zOne,zTwo,totArea);
992 }
G4ThreeVector GetPointOnRing(G4double fRMin, G4double fRMax, G4double fRMin2, G4double fRMax2, G4double zOne) const
Definition: G4Polycone.cc:936
G4ThreeVector GetPointOnCone(G4double fRmin1, G4double fRmax1, G4double fRmin2, G4double fRmax2, G4double zOne, G4double zTwo, G4double &totArea) const
Definition: G4Polycone.cc:782
G4ThreeVector GetPointOnTubs(G4double fRMin, G4double fRMax, G4double zOne, G4double zTwo, G4double &totArea) const
Definition: G4Polycone.cc:877

Here is the call graph for this function:

Here is the caller graph for this function:

G4ThreeVector G4Polycone::GetPointOnRing ( G4double  fRMin,
G4double  fRMax,
G4double  fRMin2,
G4double  fRMax2,
G4double  zOne 
) const
protected

Definition at line 936 of file G4Polycone.cc.

939 {
940  G4double xRand,yRand,phi,cosphi,sinphi,rRand1,rRand2,A1,Atot,rCh;
942  cosphi = std::cos(phi);
943  sinphi = std::sin(phi);
944 
945  if(fRMin1==fRMin2)
946  {
947  rRand1 = fRMin1; A1=0.;
948  }
949  else
950  {
951  rRand1 = G4RandFlat::shoot(fRMin1,fRMin2);
952  A1=std::fabs(fRMin2*fRMin2-fRMin1*fRMin1);
953  }
954  if(fRMax1==fRMax2)
955  {
956  rRand2=fRMax1; Atot=A1;
957  }
958  else
959  {
960  rRand2 = G4RandFlat::shoot(fRMax1,fRMax2);
961  Atot = A1+std::fabs(fRMax2*fRMax2-fRMax1*fRMax1);
962  }
963  rCh = G4RandFlat::shoot(0.,Atot);
964 
965  if(rCh>A1) { rRand1=rRand2; }
966 
967  xRand = rRand1*cosphi;
968  yRand = rRand1*sinphi;
969 
970  return G4ThreeVector(xRand, yRand, zOne);
971 }
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4double endPhi
Definition: G4Polycone.hh:194
G4double startPhi
Definition: G4Polycone.hh:193
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4ThreeVector G4Polycone::GetPointOnSurface ( ) const
virtual

Reimplemented from G4VSolid.

Definition at line 998 of file G4Polycone.cc.

999 {
1000  G4double Area=0,totArea=0,Achose1=0,Achose2=0,phi,cosphi,sinphi,rRand;
1001  G4int i=0;
1002  G4int numPlanes = original_parameters->Num_z_planes;
1003 
1005  cosphi = std::cos(phi);
1006  sinphi = std::sin(phi);
1007 
1008  rRand = original_parameters->Rmin[0] +
1010  * std::sqrt(G4RandFlat::shoot()) );
1011 
1012  std::vector<G4double> areas; // (numPlanes+1);
1013  std::vector<G4ThreeVector> points; // (numPlanes-1);
1014 
1015  areas.push_back(pi*(sqr(original_parameters->Rmax[0])
1016  -sqr(original_parameters->Rmin[0])));
1017 
1018  for(i=0; i<numPlanes-1; i++)
1019  {
1021  * std::sqrt(sqr(original_parameters->Rmin[i]
1022  -original_parameters->Rmin[i+1])+
1025 
1026  Area += (original_parameters->Rmax[i]+original_parameters->Rmax[i+1])
1027  * std::sqrt(sqr(original_parameters->Rmax[i]
1028  -original_parameters->Rmax[i+1])+
1031 
1032  Area *= 0.5*(endPhi-startPhi);
1033 
1034  if(startPhi==0.&& endPhi == twopi)
1035  {
1036  Area += std::fabs(original_parameters->Z_values[i+1]
1039  +original_parameters->Rmax[i+1]
1041  -original_parameters->Rmin[i+1]);
1042  }
1043  areas.push_back(Area);
1044  totArea += Area;
1045  }
1046 
1047  areas.push_back(pi*(sqr(original_parameters->Rmax[numPlanes-1])-
1048  sqr(original_parameters->Rmin[numPlanes-1])));
1049 
1050  totArea += (areas[0]+areas[numPlanes]);
1051  G4double chose = G4RandFlat::shoot(0.,totArea);
1052 
1053  if( (chose>=0.) && (chose<areas[0]) )
1054  {
1055  return G4ThreeVector(rRand*cosphi, rRand*sinphi,
1057  }
1058 
1059  for (i=0; i<numPlanes-1; i++)
1060  {
1061  Achose1 += areas[i];
1062  Achose2 = (Achose1+areas[i+1]);
1063  if(chose>=Achose1 && chose<Achose2)
1064  {
1067  original_parameters->Rmin[i+1],
1068  original_parameters->Rmax[i+1],
1070  original_parameters->Z_values[i+1], Area);
1071  }
1072  }
1073 
1074  rRand = original_parameters->Rmin[numPlanes-1] +
1075  ( (original_parameters->Rmax[numPlanes-1]-original_parameters->Rmin[numPlanes-1])
1076  * std::sqrt(G4RandFlat::shoot()) );
1077 
1078  return G4ThreeVector(rRand*cosphi,rRand*sinphi,
1079  original_parameters->Z_values[numPlanes-1]);
1080 
1081 }
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4double endPhi
Definition: G4Polycone.hh:194
int G4int
Definition: G4Types.hh:78
static constexpr double twopi
Definition: G4SIunits.hh:76
G4ThreeVector GetPointOnCut(G4double fRMin1, G4double fRMax1, G4double fRMin2, G4double fRMax2, G4double zOne, G4double zTwo, G4double &totArea) const
Definition: G4Polycone.cc:979
G4double startPhi
Definition: G4Polycone.hh:193
static constexpr double pi
Definition: G4SIunits.hh:75
T sqr(const T &x)
Definition: templates.hh:145
double G4double
Definition: G4Types.hh:76
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

G4ThreeVector G4Polycone::GetPointOnTubs ( G4double  fRMin,
G4double  fRMax,
G4double  zOne,
G4double  zTwo,
G4double totArea 
) const
protected

Definition at line 877 of file G4Polycone.cc.

880 {
881  G4double xRand,yRand,zRand,phi,cosphi,sinphi,chose,
882  aOne,aTwo,aFou,rRand,fDz,fSPhi,fDPhi;
883  fDz = std::fabs(0.5*(zTwo-zOne));
884  fSPhi = startPhi;
885  fDPhi = endPhi-startPhi;
886 
887  aOne = 2.*fDz*fDPhi*fRMax;
888  aTwo = 2.*fDz*fDPhi*fRMin;
889  aFou = 2.*fDz*(fRMax-fRMin);
890  totArea = aOne+aTwo+2.*aFou;
891  phi = G4RandFlat::shoot(startPhi,endPhi);
892  cosphi = std::cos(phi);
893  sinphi = std::sin(phi);
894  rRand = fRMin + (fRMax-fRMin)*std::sqrt(G4RandFlat::shoot());
895 
896  if(startPhi == 0 && endPhi == twopi)
897  aFou = 0;
898 
899  chose = G4RandFlat::shoot(0.,aOne+aTwo+2.*aFou);
900  if( (chose >= 0) && (chose < aOne) )
901  {
902  xRand = fRMax*cosphi;
903  yRand = fRMax*sinphi;
904  zRand = G4RandFlat::shoot(-1.*fDz,fDz);
905  return G4ThreeVector(xRand, yRand, zRand+0.5*(zTwo+zOne));
906  }
907  else if( (chose >= aOne) && (chose < aOne + aTwo) )
908  {
909  xRand = fRMin*cosphi;
910  yRand = fRMin*sinphi;
911  zRand = G4RandFlat::shoot(-1.*fDz,fDz);
912  return G4ThreeVector(xRand, yRand, zRand+0.5*(zTwo+zOne));
913  }
914  else if( (chose >= aOne+aTwo) && (chose <aOne+aTwo+aFou) )
915  {
916  xRand = rRand*std::cos(fSPhi+fDPhi);
917  yRand = rRand*std::sin(fSPhi+fDPhi);
918  zRand = G4RandFlat::shoot(-1.*fDz,fDz);
919  return G4ThreeVector(xRand, yRand, zRand+0.5*(zTwo+zOne));
920  }
921 
922  // else
923 
924  xRand = rRand*std::cos(fSPhi+fDPhi);
925  yRand = rRand*std::sin(fSPhi+fDPhi);
926  zRand = G4RandFlat::shoot(-1.*fDz,fDz);
927  return G4ThreeVector(xRand, yRand, zRand+0.5*(zTwo+zOne));
928 }
ThreeVector shoot(const G4int Ap, const G4int Af)
CLHEP::Hep3Vector G4ThreeVector
G4double endPhi
Definition: G4Polycone.hh:194
static constexpr double twopi
Definition: G4SIunits.hh:76
G4double startPhi
Definition: G4Polycone.hh:193
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4Polycone::GetSinEndPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Polycone::GetSinStartPhi ( ) const
inline

Here is the caller graph for this function:

G4double G4Polycone::GetStartPhi ( ) const
inline

Here is the caller graph for this function:

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

Reimplemented from G4VCSGfaceted.

Definition at line 478 of file G4Polycone.cc.

479 {
480  //
481  // Quick test
482  //
483  if (enclosingCylinder->MustBeOutside(p)) return kOutside;
484 
485  //
486  // Long answer
487  //
488  return G4VCSGfaceted::Inside(p);
489 }
G4bool MustBeOutside(const G4ThreeVector &p) const
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
virtual EInside Inside(const G4ThreeVector &p) const

Here is the call graph for this function:

G4bool G4Polycone::IsOpen ( ) const
inline

Here is the caller graph for this function:

G4Polycone & G4Polycone::operator= ( const G4Polycone source)

Definition at line 383 of file G4Polycone.cc.

384 {
385  if (this == &source) return *this;
386 
387  G4VCSGfaceted::operator=( source );
388 
389  delete [] corners;
391 
392  delete enclosingCylinder;
393 
394  CopyStuff( source );
395 
396  return *this;
397 }
G4VCSGfaceted & operator=(const G4VCSGfaceted &source)
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
void CopyStuff(const G4Polycone &source)
Definition: G4Polycone.cc:403
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

G4bool G4Polycone::Reset ( )

Definition at line 447 of file G4Polycone.cc.

448 {
449  //
450  // Clear old setup
451  //
453  delete [] corners;
454  delete enclosingCylinder;
455 
456  //
457  // Rebuild polycone
458  //
459  G4ReduciblePolygon *rz =
466  delete rz;
467 
468  return 0;
469 }
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
G4EnclosingCylinder * enclosingCylinder
Definition: G4Polycone.hh:202
void Create(G4double phiStart, G4double phiTotal, G4ReduciblePolygon *rz)
Definition: G4Polycone.cc:181
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

Here is the caller graph for this function:

void G4Polycone::SetOriginalParameters ( G4PolyconeHistorical pars)
inline

Here is the caller graph for this function:

G4bool G4Polycone::SetOriginalParameters ( G4ReduciblePolygon rz)
protected

Definition at line 1100 of file G4Polycone.cc.

1101 {
1102  G4int numPlanes = (G4int)numCorner;
1103  G4bool isConvertible=true;
1104  G4double Zmax=rz->Bmax();
1105  rz->StartWithZMin();
1106 
1107  // Prepare vectors for storage
1108  //
1109  std::vector<G4double> Z;
1110  std::vector<G4double> Rmin;
1111  std::vector<G4double> Rmax;
1112 
1113  G4int countPlanes=1;
1114  G4int icurr=0;
1115  G4int icurl=0;
1116 
1117  // first plane Z=Z[0]
1118  //
1119  Z.push_back(corners[0].z);
1120  G4double Zprev=Z[0];
1121  if (Zprev == corners[1].z)
1122  {
1123  Rmin.push_back(corners[0].r);
1124  Rmax.push_back (corners[1].r);icurr=1;
1125  }
1126  else if (Zprev == corners[numPlanes-1].z)
1127  {
1128  Rmin.push_back(corners[numPlanes-1].r);
1129  Rmax.push_back (corners[0].r);
1130  icurl=numPlanes-1;
1131  }
1132  else
1133  {
1134  Rmin.push_back(corners[0].r);
1135  Rmax.push_back (corners[0].r);
1136  }
1137 
1138  // next planes until last
1139  //
1140  G4int inextr=0, inextl=0;
1141  for (G4int i=0; i < numPlanes-2; i++)
1142  {
1143  inextr=1+icurr;
1144  inextl=(icurl <= 0)? numPlanes-1 : icurl-1;
1145 
1146  if((corners[inextr].z >= Zmax) & (corners[inextl].z >= Zmax)) { break; }
1147 
1148  G4double Zleft = corners[inextl].z;
1149  G4double Zright = corners[inextr].z;
1150  if(Zright > Zleft) // Next plane will be Zleft
1151  {
1152  Z.push_back(Zleft);
1153  countPlanes++;
1154  G4double difZr=corners[inextr].z - corners[icurr].z;
1155  G4double difZl=corners[inextl].z - corners[icurl].z;
1156 
1157  if(std::fabs(difZl) < kCarTolerance)
1158  {
1159  if(std::fabs(difZr) < kCarTolerance)
1160  {
1161  Rmin.push_back(corners[inextl].r);
1162  Rmax.push_back(corners[icurr].r);
1163  }
1164  else
1165  {
1166  Rmin.push_back(corners[inextl].r);
1167  Rmax.push_back(corners[icurr].r + (Zleft-corners[icurr].z)/difZr
1168  *(corners[inextr].r - corners[icurr].r));
1169  }
1170  }
1171  else if (difZl >= kCarTolerance)
1172  {
1173  if(std::fabs(difZr) < kCarTolerance)
1174  {
1175  Rmin.push_back(corners[icurl].r);
1176  Rmax.push_back(corners[icurr].r);
1177  }
1178  else
1179  {
1180  Rmin.push_back(corners[icurl].r);
1181  Rmax.push_back(corners[icurr].r + (Zleft-corners[icurr].z)/difZr
1182  *(corners[inextr].r - corners[icurr].r));
1183  }
1184  }
1185  else
1186  {
1187  isConvertible=false; break;
1188  }
1189  icurl=(icurl == 0)? numPlanes-1 : icurl-1;
1190  }
1191  else if(std::fabs(Zright-Zleft)<kCarTolerance) // Zright=Zleft
1192  {
1193  Z.push_back(Zleft);
1194  countPlanes++;
1195  icurr++;
1196 
1197  icurl=(icurl == 0)? numPlanes-1 : icurl-1;
1198 
1199  Rmin.push_back(corners[inextl].r);
1200  Rmax.push_back(corners[inextr].r);
1201  }
1202  else // Zright<Zleft
1203  {
1204  Z.push_back(Zright);
1205  countPlanes++;
1206 
1207  G4double difZr=corners[inextr].z - corners[icurr].z;
1208  G4double difZl=corners[inextl].z - corners[icurl].z;
1209  if(std::fabs(difZr) < kCarTolerance)
1210  {
1211  if(std::fabs(difZl) < kCarTolerance)
1212  {
1213  Rmax.push_back(corners[inextr].r);
1214  Rmin.push_back(corners[icurr].r);
1215  }
1216  else
1217  {
1218  Rmin.push_back(corners[icurl].r + (Zright-corners[icurl].z)/difZl
1219  *(corners[inextl].r - corners[icurl].r));
1220  Rmax.push_back(corners[inextr].r);
1221  }
1222  icurr++;
1223  } // plate
1224  else if (difZr >= kCarTolerance)
1225  {
1226  if(std::fabs(difZl) < kCarTolerance)
1227  {
1228  Rmax.push_back(corners[inextr].r);
1229  Rmin.push_back (corners[icurr].r);
1230  }
1231  else
1232  {
1233  Rmax.push_back(corners[inextr].r);
1234  Rmin.push_back (corners[icurl].r+(Zright-corners[icurl].z)/difZl
1235  * (corners[inextl].r - corners[icurl].r));
1236  }
1237  icurr++;
1238  }
1239  else
1240  {
1241  isConvertible=false; break;
1242  }
1243  }
1244  } // end for loop
1245 
1246  // last plane Z=Zmax
1247  //
1248  Z.push_back(Zmax);
1249  countPlanes++;
1250  inextr=1+icurr;
1251  inextl=(icurl <= 0)? numPlanes-1 : icurl-1;
1252 
1253  Rmax.push_back(corners[inextr].r);
1254  Rmin.push_back(corners[inextl].r);
1255 
1256  // Set original parameters Rmin,Rmax,Z
1257  //
1258  if(isConvertible)
1259  {
1261  original_parameters->Z_values = new G4double[countPlanes];
1262  original_parameters->Rmin = new G4double[countPlanes];
1263  original_parameters->Rmax = new G4double[countPlanes];
1264 
1265  for(G4int j=0; j < countPlanes; j++)
1266  {
1267  original_parameters->Z_values[j] = Z[j];
1268  original_parameters->Rmax[j] = Rmax[j];
1269  original_parameters->Rmin[j] = Rmin[j];
1270  }
1273  original_parameters->Num_z_planes = countPlanes;
1274 
1275  }
1276  else // Set parameters(r,z) with Rmin==0 as convention
1277  {
1278 #ifdef G4SPECSDEBUG
1279  std::ostringstream message;
1280  message << "Polycone " << GetName() << G4endl
1281  << "cannot be converted to Polycone with (Rmin,Rmaz,Z) parameters!";
1282  G4Exception("G4Polycone::SetOriginalParameters()", "GeomSolids0002",
1283  JustWarning, message);
1284 #endif
1286  original_parameters->Z_values = new G4double[numPlanes];
1287  original_parameters->Rmin = new G4double[numPlanes];
1288  original_parameters->Rmax = new G4double[numPlanes];
1289 
1290  for(G4int j=0; j < numPlanes; j++)
1291  {
1293  original_parameters->Rmax[j] = corners[j].r;
1294  original_parameters->Rmin[j] = 0.0;
1295  }
1298  original_parameters->Num_z_planes = numPlanes;
1299  }
1300  return isConvertible;
1301 }
G4String GetName() const
G4double endPhi
Definition: G4Polycone.hh:194
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
int G4int
Definition: G4Types.hh:78
G4double Bmax() const
G4double startPhi
Definition: G4Polycone.hh:193
bool G4bool
Definition: G4Types.hh:79
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
#define G4endl
Definition: G4ios.hh:61
G4double kCarTolerance
Definition: G4VSolid.hh:307
G4int numCorner
Definition: G4Polycone.hh:196
double G4double
Definition: G4Types.hh:76
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

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

Reimplemented from G4VCSGfaceted.

Definition at line 730 of file G4Polycone.cc.

731 {
732  G4int oldprc = os.precision(16);
733  os << "-----------------------------------------------------------\n"
734  << " *** Dump for solid - " << GetName() << " ***\n"
735  << " ===================================================\n"
736  << " Solid type: G4Polycone\n"
737  << " Parameters: \n"
738  << " starting phi angle : " << startPhi/degree << " degrees \n"
739  << " ending phi angle : " << endPhi/degree << " degrees \n";
740  G4int i=0;
741 
743  os << " number of Z planes: " << numPlanes << "\n"
744  << " Z values: \n";
745  for (i=0; i<numPlanes; i++)
746  {
747  os << " Z plane " << i << ": "
748  << original_parameters->Z_values[i] << "\n";
749  }
750  os << " Tangent distances to inner surface (Rmin): \n";
751  for (i=0; i<numPlanes; i++)
752  {
753  os << " Z plane " << i << ": "
754  << original_parameters->Rmin[i] << "\n";
755  }
756  os << " Tangent distances to outer surface (Rmax): \n";
757  for (i=0; i<numPlanes; i++)
758  {
759  os << " Z plane " << i << ": "
760  << original_parameters->Rmax[i] << "\n";
761  }
762 
763  os << " number of RZ points: " << numCorner << "\n"
764  << " RZ values (corners): \n";
765  for (i=0; i<numCorner; i++)
766  {
767  os << " "
768  << corners[i].r << ", " << corners[i].z << "\n";
769  }
770  os << "-----------------------------------------------------------\n";
771  os.precision(oldprc);
772 
773  return os;
774 }
G4String GetName() const
G4double endPhi
Definition: G4Polycone.hh:194
G4PolyconeSideRZ * corners
Definition: G4Polycone.hh:197
int G4int
Definition: G4Types.hh:78
static constexpr double degree
Definition: G4SIunits.hh:144
G4double startPhi
Definition: G4Polycone.hh:193
G4int numCorner
Definition: G4Polycone.hh:196
G4PolyconeHistorical * original_parameters
Definition: G4Polycone.hh:198

Here is the call graph for this function:

Member Data Documentation

G4PolyconeSideRZ* G4Polycone::corners
protected

Definition at line 197 of file G4Polycone.hh.

G4EnclosingCylinder* G4Polycone::enclosingCylinder
protected

Definition at line 202 of file G4Polycone.hh.

G4double G4Polycone::endPhi
protected

Definition at line 194 of file G4Polycone.hh.

G4int G4Polycone::numCorner
protected

Definition at line 196 of file G4Polycone.hh.

G4PolyconeHistorical* G4Polycone::original_parameters
protected

Definition at line 198 of file G4Polycone.hh.

G4bool G4Polycone::phiIsOpen
protected

Definition at line 195 of file G4Polycone.hh.

G4double G4Polycone::startPhi
protected

Definition at line 193 of file G4Polycone.hh.


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