Geant4_10
G4PolyhedraSide.hh
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 // $Id: G4PolyhedraSide.hh 67973 2013-03-13 10:16:25Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------------
31 // GEANT 4 class header file
32 //
33 //
34 // G4PolyhedraSide
35 //
36 // Class description:
37 //
38 // Class implementing a face that represents one segmented side
39 // of a polyhedra:
40 //
41 // G4PolyhedraSide( const G4PolyhedraSideRZ *prevRZ,
42 // const G4PolyhedraSideRZ *tail,
43 // const G4PolyhedraSideRZ *head,
44 // const G4PolyhedraSideRZ *nextRZ,
45 // G4int numSide,
46 // G4double phiStart, G4double phiTotal,
47 // G4bool phiIsOpen, G4bool isAllBehind=false )
48 //
49 // Values for r1,z1 and r2,z2 should be specified in clockwise
50 // order in (r,z).
51 
52 // Author:
53 // David C. Williams (davidw@scipp.ucsc.edu)
54 // --------------------------------------------------------------------
55 
56 #ifndef G4PolyhedraSide_hh
57 #define G4PolyhedraSide_hh
58 
59 #include "G4VCSGface.hh"
60 
61 class G4IntersectingCone;
62 
64 {
65  G4double r, z; // start of vector
66 };
67 
68 // ----------------------------------------------------------------------------
69 // MT-specific utility code
70 
71 #include "G4GeomSplitter.hh"
72 
73 // The class G4PhSideData is introduced to encapsulate the
74 // fields of the class G4PolyhedraSide that may not be read-only.
75 //
77 {
78  public:
79  void initialize()
80  {
81  fPhi.first = G4ThreeVector(0,0,0);
82  fPhi.second= 0.0;
83  }
84 
85  std::pair<G4ThreeVector, G4double> fPhi; // Cached value for phi
86 
87 };
88 
89 // The type G4PhSideManager is introduced to encapsulate the methods used
90 // by both the master thread and worker threads to allocate memory space
91 // for the fields encapsulated by the class G4PhSideData.
92 //
94 
95 // This macro changes the references to fields that are now encapsulated
96 // in the class G4PhSideData.
97 //
98 #define G4MT_phphi ((subInstanceManager.offset[instanceID]).fPhi)
99 //
100 // ----------------------------------------------------------------------------
101 
103 {
104 
105  public: // with description
106 
107  G4PolyhedraSide( const G4PolyhedraSideRZ *prevRZ,
108  const G4PolyhedraSideRZ *tail,
109  const G4PolyhedraSideRZ *head,
110  const G4PolyhedraSideRZ *nextRZ,
111  G4int numSide,
112  G4double phiStart, G4double phiTotal,
113  G4bool phiIsOpen, G4bool isAllBehind=false );
114  virtual ~G4PolyhedraSide();
115 
116  G4PolyhedraSide( const G4PolyhedraSide &source );
117  G4PolyhedraSide& operator=( const G4PolyhedraSide &source );
118 
119  G4bool Intersect( const G4ThreeVector &p, const G4ThreeVector &v,
120  G4bool outgoing, G4double surfTolerance,
121  G4double &distance, G4double &distFromSurface,
122  G4ThreeVector &normal, G4bool &allBehind );
123 
124  G4double Distance( const G4ThreeVector &p, G4bool outgoing );
125 
126  EInside Inside( const G4ThreeVector &p, G4double tolerance,
127  G4double *bestDistance );
128 
129  G4ThreeVector Normal( const G4ThreeVector &p, G4double *bestDistance );
130 
131  G4double Extent( const G4ThreeVector axis );
132 
133  void CalculateExtent( const EAxis axis,
134  const G4VoxelLimits &voxelLimit,
135  const G4AffineTransform &tranform,
136  G4SolidExtentList &extentList );
137 
138  G4VCSGface *Clone() { return new G4PolyhedraSide( *this ); }
139 
140  public: // without description
141 
142  // Methods used for GetPointOnSurface()
143 
145  G4ThreeVector p2,
146  G4ThreeVector p3,
147  G4ThreeVector *p4 );
150  G4double *Area );
153 
154  public: // without description
155 
156  G4PolyhedraSide(__void__&);
157  // Fake default constructor for usage restricted to direct object
158  // persistency for clients requiring preallocation of memory for
159  // persistifiable objects.
160 
161  inline G4int GetInstanceID() const { return instanceID; }
162  // Returns the instance ID.
163 
164  static const G4PhSideManager& GetSubInstanceManager();
165  // Returns the private data instance manager.
166 
167  protected:
168 
169  //
170  // A couple internal data structures
171  //
172  struct sG4PolyhedraSideVec; // Secret recipe for allowing
173  friend struct sG4PolyhedraSideVec; // protected nested structures
174 
175  typedef struct sG4PolyhedraSideEdge
176  {
177  G4ThreeVector normal; // Unit normal to this edge
178  G4ThreeVector corner[2]; // The two corners of this phi edge
179  G4ThreeVector cornNorm[2]; // The normals of these corners
181 
182  typedef struct sG4PolyhedraSideVec
183  {
184  G4ThreeVector normal, // Normal (point out of the shape)
185  center, // Point in center of side
186  surfPhi, // Unit vector on surface pointing along phi
187  surfRZ; // Unit vector on surface pointing along R/Z
188  G4PolyhedraSideEdge *edges[2]; // The phi boundary edges to this side
189  // [0]=low phi [1]=high phi
190  G4ThreeVector edgeNorm[2]; // RZ edge normals [i] at {r[i],z[i]}
192 
194  const G4PolyhedraSideVec& vec,
195  G4double normSign,
196  G4double surfTolerance,
197  G4double &distance,
198  G4double &distFromSurface );
199 
201  const G4ThreeVector &v,
202  G4int *i1, G4int *i2 );
203 
205 
206  G4int PhiSegment( G4double phi );
207 
208  G4double GetPhi( const G4ThreeVector& p );
209 
211  const G4PolyhedraSideVec &vec,
212  G4double *normDist );
213 
215  const G4PolyhedraSideVec &vec,
216  G4double *normDist );
217 
218  void CopyStuff( const G4PolyhedraSide &source );
219 
220  protected:
221 
222  G4int numSide; // Number sides
223  G4double r[2], z[2]; // r, z parameters, in specified order
224  G4double startPhi, // Start phi (0 to 2pi), if phiIsOpen
225  deltaPhi, // Delta phi (0 to 2pi), if phiIsOpen
226  endPhi; // End phi (>startPhi), if phiIsOpen
227  G4bool phiIsOpen; // True if there is a phi slice
228  G4bool allBehind; // True if the entire solid is "behind" this face
229 
230  G4IntersectingCone *cone; // Our intersecting cone
231 
232  G4PolyhedraSideVec *vecs; // Vector set for each facet of our face
233  G4PolyhedraSideEdge *edges; // The edges belong to vecs
234  G4double lenRZ, // RZ length of each side
235  lenPhi[2]; // Phi dimensions of each side
236  G4double edgeNorm; // Normal in RZ/Phi space to each side
237 
238  private:
239 
240  G4double kCarTolerance; // Geometrical surface thickness
241  G4double fSurfaceArea; // Surface Area
242 
243  G4int instanceID;
244  // This field is used as instance ID.
245  G4GEOM_DLL static G4PhSideManager subInstanceManager;
246  // This field helps to use the class G4PhSideManager introduced above.
247 };
248 
249 #endif
G4ThreeVector GetPointOnPlane(G4ThreeVector p0, G4ThreeVector p1, G4ThreeVector p2, G4ThreeVector p3, G4double *Area)
G4double GetPhi(const G4ThreeVector &p)
G4double lenPhi[2]
G4double Extent(const G4ThreeVector axis)
G4PolyhedraSide(const G4PolyhedraSideRZ *prevRZ, const G4PolyhedraSideRZ *tail, const G4PolyhedraSideRZ *head, const G4PolyhedraSideRZ *nextRZ, G4int numSide, G4double phiStart, G4double phiTotal, G4bool phiIsOpen, G4bool isAllBehind=false)
CLHEP::Hep3Vector G4ThreeVector
EInside Inside(const G4ThreeVector &p, G4double tolerance, G4double *bestDistance)
virtual ~G4PolyhedraSide()
G4int ClosestPhiSegment(G4double phi)
std::pair< G4ThreeVector, G4double > fPhi
G4double DistanceToOneSide(const G4ThreeVector &p, const G4PolyhedraSideVec &vec, G4double *normDist)
const char * p
Definition: xmltok.h:285
G4VCSGface * Clone()
static const G4PhSideManager & GetSubInstanceManager()
G4ThreeVector Normal(const G4ThreeVector &p, G4double *bestDistance)
G4int LineHitsSegments(const G4ThreeVector &p, const G4ThreeVector &v, G4int *i1, G4int *i2)
struct G4PolyhedraSide::sG4PolyhedraSideEdge G4PolyhedraSideEdge
G4PolyhedraSideVec * vecs
int G4int
Definition: G4Types.hh:78
G4double DistanceAway(const G4ThreeVector &p, const G4PolyhedraSideVec &vec, G4double *normDist)
G4double Distance(const G4ThreeVector &p, G4bool outgoing)
G4double SurfaceTriangle(G4ThreeVector p1, G4ThreeVector p2, G4ThreeVector p3, G4ThreeVector *p4)
G4bool IntersectSidePlane(const G4ThreeVector &p, const G4ThreeVector &v, const G4PolyhedraSideVec &vec, G4double normSign, G4double surfTolerance, G4double &distance, G4double &distFromSurface)
struct G4PolyhedraSide::sG4PolyhedraSideVec G4PolyhedraSideVec
bool G4bool
Definition: G4Types.hh:79
G4GeomSplitter< G4PhSideData > G4PhSideManager
G4int PhiSegment(G4double phi)
tuple v
Definition: test.py:18
G4bool Intersect(const G4ThreeVector &p, const G4ThreeVector &v, G4bool outgoing, G4double surfTolerance, G4double &distance, G4double &distFromSurface, G4ThreeVector &normal, G4bool &allBehind)
G4ThreeVector GetPointOnFace()
EInside
Definition: geomdefs.hh:58
EAxis
Definition: geomdefs.hh:54
G4double SurfaceArea()
void CopyStuff(const G4PolyhedraSide &source)
G4IntersectingCone * cone
G4PolyhedraSideEdge * edges
G4int GetInstanceID() const
double G4double
Definition: G4Types.hh:76
void CalculateExtent(const EAxis axis, const G4VoxelLimits &voxelLimit, const G4AffineTransform &tranform, G4SolidExtentList &extentList)
G4PolyhedraSide & operator=(const G4PolyhedraSide &source)
#define G4GEOM_DLL
Definition: geomwdefs.hh:48