Geant4_10
G4VTwistedFaceted.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: G4VTwistedFaceted.hh 66356 2012-12-18 09:02:32Z gcosmo $
28 //
29 // --------------------------------------------------------------------
30 // GEANT 4 class header file
31 //
32 //
33 // G4VTwistedFaceted
34 //
35 // Class description:
36 //
37 // G4VTwistedFaceted is an abstract base class for twisted boxoids:
38 // G4TwistedTrd, G4TwistedTrap and G4TwistedBox
39 
40 // Author:
41 //
42 // 27-Oct-2004 - O.Link (Oliver.Link@cern.ch)
43 //
44 // --------------------------------------------------------------------
45 
46 #ifndef __G4VTWISTEDFACETED__
47 #define __G4VTWISTEDFACETED__
48 
49 #include "G4VSolid.hh"
50 #include "G4TwistTrapAlphaSide.hh"
52 #include "G4TwistBoxSide.hh"
53 #include "G4TwistTrapFlatSide.hh"
54 
55 class G4SolidExtentList;
56 class G4ClippablePolygon;
57 
59 {
60  public: // with description
61 
62  G4VTwistedFaceted(const G4String &pname, // Name of instance
63  G4double PhiTwist, // twist angle
64  G4double pDz, // half z lenght
65  G4double pTheta, // direction between end planes
66  G4double pPhi, // defined by polar & azim. angles
67  G4double pDy1, // half y length at -pDz
68  G4double pDx1, // half x length at -pDz,-pDy
69  G4double pDx2, // half x length at -pDz,+pDy
70  G4double pDy2, // half y length at +pDz
71  G4double pDx3, // half x length at +pDz,-pDy
72  G4double pDx4, // half x length at +pDz,+pDy
73  G4double pAlph // tilt angle at +pDz
74  );
75 
76  virtual ~G4VTwistedFaceted();
77 
79  const G4int,
80  const G4VPhysicalVolume* );
81 
82  virtual G4bool CalculateExtent(const EAxis pAxis,
83  const G4VoxelLimits &pVoxelLimit,
84  const G4AffineTransform &pTransform,
85  G4double &pMin,
86  G4double &pMax ) const;
87 
88  virtual G4double DistanceToIn (const G4ThreeVector &p,
89  const G4ThreeVector &v ) const;
90 
91  virtual G4double DistanceToIn (const G4ThreeVector &p ) const;
92 
93  virtual G4double DistanceToOut(const G4ThreeVector &p,
94  const G4ThreeVector &v,
95  const G4bool calcnorm = false,
96  G4bool *validnorm = 0,
97  G4ThreeVector *n=0 ) const;
98 
99  virtual G4double DistanceToOut(const G4ThreeVector &p) const;
100 
101  virtual EInside Inside (const G4ThreeVector &p) const;
102 
103  virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const;
104 
107 
108  virtual inline G4double GetCubicVolume() ;
109  virtual inline G4double GetSurfaceArea() ;
110 
111  virtual void DescribeYourselfTo (G4VGraphicsScene &scene) const;
112  virtual G4Polyhedron *CreatePolyhedron () const ;
113  virtual G4Polyhedron *GetPolyhedron () const;
114 
115  virtual std::ostream &StreamInfo(std::ostream& os) const;
116 
117  // accessors
118 
119  inline G4double GetTwistAngle () const { return fPhiTwist; }
120 
121  inline G4double GetDx1 () const { return fDx1 ; }
122  inline G4double GetDx2 () const { return fDx2 ; }
123  inline G4double GetDx3 () const { return fDx3 ; }
124  inline G4double GetDx4 () const { return fDx4 ; }
125  inline G4double GetDy1 () const { return fDy1 ; }
126  inline G4double GetDy2 () const { return fDy2 ; }
127  inline G4double GetDz () const { return fDz ; }
128  inline G4double GetPhi () const { return fPhi ; }
129  inline G4double GetTheta () const { return fTheta ; }
130  inline G4double GetAlpha () const { return fAlph ; }
131 
132  inline G4double Xcoef(G4double u,G4double phi, G4double ftg) const ;
133  // For calculating the w(u) function
134 
135  inline G4double GetValueA(G4double phi) const;
136  inline G4double GetValueB(G4double phi) const;
137  inline G4double GetValueD(G4double phi) const;
138 
139  virtual G4VisExtent GetExtent () const;
140  virtual G4GeometryType GetEntityType() const;
141 
142  public: // without description
143 
144  G4VTwistedFaceted(__void__&);
145  // Fake default constructor for usage restricted to direct object
146  // persistency for clients requiring preallocation of memory for
147  // persistifiable objects.
148 
151  // Copy constructor and assignment operator.
152 
153  protected: // with description
154 
156  CreateRotatedVertices(const G4AffineTransform& pTransform) const;
157  // Create the List of transformed vertices in the format required
158  // for G4VSolid:: ClipCrossSection and ClipBetweenSections.
159 
160  private:
161 
162  void CreateSurfaces();
163 
164  private:
165 
166  G4double fTheta;
167  G4double fPhi ;
168 
169  G4double fDy1;
170  G4double fDx1;
171  G4double fDx2;
172 
173  G4double fDy2;
174  G4double fDx3;
175  G4double fDx4;
176 
177  G4double fDz; // Half-length along the z axis
178 
179  G4double fDx ; // maximum side in x
180  G4double fDy ; // maximum side in y
181 
182  G4double fAlph ;
183  G4double fTAlph ; // std::tan(fAlph)
184 
185  G4double fdeltaX ;
186  G4double fdeltaY ;
187 
188  G4double fPhiTwist; // twist angle ( dphi in surface equation)
189 
190  G4VTwistSurface *fLowerEndcap ; // surface of -ve z
191  G4VTwistSurface *fUpperEndcap ; // surface of +ve z
192 
193  G4VTwistSurface *fSide0 ; // Twisted Side at phi = 0 deg
194  G4VTwistSurface *fSide90 ; // Twisted Side at phi = 90 deg
195  G4VTwistSurface *fSide180 ; // Twisted Side at phi = 180 deg
196  G4VTwistSurface *fSide270 ; // Twisted Side at phi = 270 deg
197 
198  G4double fCubicVolume ; // volume of the solid
199  G4double fSurfaceArea ; // area of the solid
200 
201  mutable G4Polyhedron* fpPolyhedron; // pointer to polyhedron for vis
202 
203  class LastState // last Inside result
204  {
205  public:
206  LastState()
207  {
208  p.set(kInfinity,kInfinity,kInfinity); inside = kOutside;
209  }
210  ~LastState(){}
211  LastState(const LastState& r) : p(r.p), inside(r.inside){}
212  LastState& operator=(const LastState& r)
213  {
214  if (this == &r) { return *this; }
215  p = r.p; inside = r.inside;
216  return *this;
217  }
218  public:
220  EInside inside;
221  };
222 
223  class LastVector // last SurfaceNormal result
224  {
225  public:
226  LastVector()
227  {
228  p.set(kInfinity,kInfinity,kInfinity);
229  vec.set(kInfinity,kInfinity,kInfinity);
230  surface = new G4VTwistSurface*[1];
231  }
232  ~LastVector()
233  {
234  delete [] surface;
235  }
236  LastVector(const LastVector& r) : p(r.p), vec(r.vec)
237  {
238  surface = new G4VTwistSurface*[1];
239  surface[0] = r.surface[0];
240  }
241  LastVector& operator=(const LastVector& r)
242  {
243  if (&r == this) { return *this; }
244  p = r.p; vec = r.vec;
245  delete [] surface; surface = new G4VTwistSurface*[1];
246  surface[0] = r.surface[0];
247  return *this;
248  }
249  public:
251  G4ThreeVector vec;
252  G4VTwistSurface **surface;
253  };
254 
255  class LastValue // last G4double value
256  {
257  public:
258  LastValue()
259  {
260  p.set(kInfinity,kInfinity,kInfinity);
261  value = DBL_MAX;
262  }
263  ~LastValue(){}
264  LastValue(const LastValue& r) : p(r.p), value(r.value){}
265  LastValue& operator=(const LastValue& r)
266  {
267  if (this == &r) { return *this; }
268  p = r.p; value = r.value;
269  return *this;
270  }
271  public:
273  G4double value;
274  };
275 
276  class LastValueWithDoubleVector // last G4double value
277  {
278  public:
279  LastValueWithDoubleVector()
280  {
281  p.set(kInfinity,kInfinity,kInfinity);
282  vec.set(kInfinity,kInfinity,kInfinity);
283  value = DBL_MAX;
284  }
285  ~LastValueWithDoubleVector(){}
286  LastValueWithDoubleVector(const LastValueWithDoubleVector& r)
287  : p(r.p), vec(r.vec), value(r.value){}
288  LastValueWithDoubleVector& operator=(const LastValueWithDoubleVector& r)
289  {
290  if (this == &r) { return *this; }
291  p = r.p; vec = r.vec; value = r.value;
292  return *this;
293  }
294  public:
296  G4ThreeVector vec;
297  G4double value;
298  };
299 
300  LastState fLastInside;
301  LastVector fLastNormal;
302  LastValue fLastDistanceToIn;
303  LastValue fLastDistanceToOut;
304  LastValueWithDoubleVector fLastDistanceToInWithV;
305  LastValueWithDoubleVector fLastDistanceToOutWithV;
306 
307  };
308 
309 //=====================================================================
310 
311 inline
313 {
314  if(fCubicVolume != 0.) ;
315  else fCubicVolume = 2 * fDz
316  * ( ( fDx1 + fDx2 ) * fDy1 + ( fDx3 + fDx4 ) * fDy2 );
317  return fCubicVolume;
318 }
319 
320 inline
322 {
323  if(fSurfaceArea != 0.) ;
324  else fSurfaceArea = G4VSolid::GetSurfaceArea();
325  return fSurfaceArea;
326 }
327 
328 inline
330 {
331  return ( fDx4 + fDx2 + ( fDx4 - fDx2 ) * ( 2 * phi ) / fPhiTwist ) ;
332 }
333 
334 inline
336 {
337  return ( fDx3 + fDx1 + ( fDx3 - fDx1 ) * ( 2 * phi ) / fPhiTwist ) ;
338 }
339 
340 inline
342 {
343  return ( fDy2 + fDy1 + ( fDy2 - fDy1 ) * ( 2 * phi ) / fPhiTwist ) ;
344 }
345 
346 inline
348 {
349  return GetValueA(phi)/2. + (GetValueD(phi)-GetValueA(phi))/4.
350  - u*( ( GetValueD(phi)-GetValueA(phi) ) / ( 2 * GetValueB(phi) ) - ftg );
351 }
352 
353 #endif
void set(double x, double y, double z)
virtual void ComputeDimensions(G4VPVParameterisation *, const G4int, const G4VPhysicalVolume *)
virtual G4double GetSurfaceArea()
G4ThreeVectorList * CreateRotatedVertices(const G4AffineTransform &pTransform) const
virtual G4double GetCubicVolume()
G4VTwistedFaceted(const G4String &pname, G4double PhiTwist, G4double pDz, G4double pTheta, G4double pPhi, G4double pDy1, G4double pDx1, G4double pDx2, G4double pDy2, G4double pDx3, G4double pDx4, G4double pAlph)
virtual G4GeometryType GetEntityType() const
const char * p
Definition: xmltok.h:285
G4ThreeVector GetPointOnSurface() const
G4double GetDy2() const
G4ThreeVector GetPointInSolid(G4double z) const
virtual G4Polyhedron * CreatePolyhedron() const
G4double GetDy1() const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4double GetDx3() const
G4VTwistedFaceted & operator=(const G4VTwistedFaceted &rhs)
G4double GetTwistAngle() const
int G4int
Definition: G4Types.hh:78
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
G4double Xcoef(G4double u, G4double phi, G4double ftg) const
G4double GetValueD(G4double phi) const
Char_t n[5]
G4double GetTheta() const
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
bool G4bool
Definition: G4Types.hh:79
G4double GetValueB(G4double phi) const
std::vector< G4ThreeVector > G4ThreeVectorList
Definition: G4VSolid.hh:79
tuple v
Definition: test.py:18
string pname
Definition: eplot.py:33
G4double GetDx2() const
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcnorm=false, G4bool *validnorm=0, G4ThreeVector *n=0) const
G4double GetPhi() const
jump r
Definition: plot.C:36
virtual G4Polyhedron * GetPolyhedron() const
EInside
Definition: geomdefs.hh:58
virtual G4VisExtent GetExtent() const
virtual EInside Inside(const G4ThreeVector &p) const
EAxis
Definition: geomdefs.hh:54
virtual void DescribeYourselfTo(G4VGraphicsScene &scene) const
tuple z
Definition: test.py:28
G4double GetAlpha() const
const XML_Char int const XML_Char * value
Definition: expat.h:331
G4double GetDx1() const
double G4double
Definition: G4Types.hh:76
virtual std::ostream & StreamInfo(std::ostream &os) const
G4double GetDx4() const
#define DBL_MAX
Definition: templates.hh:83
virtual G4double GetSurfaceArea()
Definition: G4VSolid.cc:250
G4double GetValueA(G4double phi) const
G4double GetDz() const