Geant4_10
G4TwistedTubs.cc
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: G4TwistedTubs.cc 72937 2013-08-14 13:20:38Z gcosmo $
28 //
29 //
30 // --------------------------------------------------------------------
31 // GEANT 4 class source file
32 //
33 //
34 // G4TwistTubsSide.cc
35 //
36 // Author:
37 // 01-Aug-2002 - Kotoyo Hoshina (hoshina@hepburn.s.chiba-u.ac.jp)
38 //
39 // History:
40 // 13-Nov-2003 - O.Link (Oliver.Link@cern.ch), Integration in Geant4
41 // from original version in Jupiter-2.5.02 application.
42 // --------------------------------------------------------------------
43 
44 #include "G4TwistedTubs.hh"
45 
46 #include "G4PhysicalConstants.hh"
47 #include "G4SystemOfUnits.hh"
48 #include "G4GeometryTolerance.hh"
49 #include "G4VoxelLimits.hh"
50 #include "G4AffineTransform.hh"
51 #include "G4SolidExtentList.hh"
52 #include "G4ClippablePolygon.hh"
53 #include "G4VPVParameterisation.hh"
54 #include "meshdefs.hh"
55 
56 #include "G4VGraphicsScene.hh"
57 #include "G4Polyhedron.hh"
58 #include "G4VisExtent.hh"
59 
60 #include "Randomize.hh"
61 
62 //=====================================================================
63 //* constructors ------------------------------------------------------
64 
65 G4TwistedTubs::G4TwistedTubs(const G4String &pname,
66  G4double twistedangle,
67  G4double endinnerrad,
68  G4double endouterrad,
69  G4double halfzlen,
70  G4double dphi)
71  : G4VSolid(pname), fDPhi(dphi),
72  fLowerEndcap(0), fUpperEndcap(0), fLatterTwisted(0),
73  fFormerTwisted(0), fInnerHype(0), fOuterHype(0),
74  fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
75 {
76  if (endinnerrad < DBL_MIN)
77  {
78  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
79  FatalErrorInArgument, "Invalid end-inner-radius!");
80  }
81 
82  G4double sinhalftwist = std::sin(0.5 * twistedangle);
83 
84  G4double endinnerradX = endinnerrad * sinhalftwist;
85  G4double innerrad = std::sqrt( endinnerrad * endinnerrad
86  - endinnerradX * endinnerradX );
87 
88  G4double endouterradX = endouterrad * sinhalftwist;
89  G4double outerrad = std::sqrt( endouterrad * endouterrad
90  - endouterradX * endouterradX );
91 
92  // temporary treatment!!
93  SetFields(twistedangle, innerrad, outerrad, -halfzlen, halfzlen);
94  CreateSurfaces();
95 }
96 
97 G4TwistedTubs::G4TwistedTubs(const G4String &pname,
98  G4double twistedangle,
99  G4double endinnerrad,
100  G4double endouterrad,
101  G4double halfzlen,
102  G4int nseg,
103  G4double totphi)
104  : G4VSolid(pname),
105  fLowerEndcap(0), fUpperEndcap(0), fLatterTwisted(0),
106  fFormerTwisted(0), fInnerHype(0), fOuterHype(0),
107  fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
108 {
109 
110  if (!nseg)
111  {
112  std::ostringstream message;
113  message << "Invalid number of segments." << G4endl
114  << " nseg = " << nseg;
115  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
116  FatalErrorInArgument, message);
117  }
118  if (totphi == DBL_MIN || endinnerrad < DBL_MIN)
119  {
120  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
121  FatalErrorInArgument, "Invalid total-phi or end-inner-radius!");
122  }
123 
124  G4double sinhalftwist = std::sin(0.5 * twistedangle);
125 
126  G4double endinnerradX = endinnerrad * sinhalftwist;
127  G4double innerrad = std::sqrt( endinnerrad * endinnerrad
128  - endinnerradX * endinnerradX );
129 
130  G4double endouterradX = endouterrad * sinhalftwist;
131  G4double outerrad = std::sqrt( endouterrad * endouterrad
132  - endouterradX * endouterradX );
133 
134  // temporary treatment!!
135  fDPhi = totphi / nseg;
136  SetFields(twistedangle, innerrad, outerrad, -halfzlen, halfzlen);
137  CreateSurfaces();
138 }
139 
140 G4TwistedTubs::G4TwistedTubs(const G4String &pname,
141  G4double twistedangle,
142  G4double innerrad,
143  G4double outerrad,
144  G4double negativeEndz,
145  G4double positiveEndz,
146  G4double dphi)
147  : G4VSolid(pname), fDPhi(dphi),
148  fLowerEndcap(0), fUpperEndcap(0), fLatterTwisted(0),
149  fFormerTwisted(0), fInnerHype(0), fOuterHype(0),
150  fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
151 {
152  if (innerrad < DBL_MIN)
153  {
154  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
155  FatalErrorInArgument, "Invalid end-inner-radius!");
156  }
157 
158  SetFields(twistedangle, innerrad, outerrad, negativeEndz, positiveEndz);
159  CreateSurfaces();
160 }
161 
162 G4TwistedTubs::G4TwistedTubs(const G4String &pname,
163  G4double twistedangle,
164  G4double innerrad,
165  G4double outerrad,
166  G4double negativeEndz,
167  G4double positiveEndz,
168  G4int nseg,
169  G4double totphi)
170  : G4VSolid(pname),
171  fLowerEndcap(0), fUpperEndcap(0), fLatterTwisted(0),
172  fFormerTwisted(0), fInnerHype(0), fOuterHype(0),
173  fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
174 {
175  if (!nseg)
176  {
177  std::ostringstream message;
178  message << "Invalid number of segments." << G4endl
179  << " nseg = " << nseg;
180  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
181  FatalErrorInArgument, message);
182  }
183  if (totphi == DBL_MIN || innerrad < DBL_MIN)
184  {
185  G4Exception("G4TwistedTubs::G4TwistedTubs()", "GeomSolids0002",
186  FatalErrorInArgument, "Invalid total-phi or end-inner-radius!");
187  }
188 
189  fDPhi = totphi / nseg;
190  SetFields(twistedangle, innerrad, outerrad, negativeEndz, positiveEndz);
191  CreateSurfaces();
192 }
193 
194 //=====================================================================
195 //* Fake default constructor ------------------------------------------
196 
197 G4TwistedTubs::G4TwistedTubs( __void__& a )
198  : G4VSolid(a), fPhiTwist(0.), fInnerRadius(0.), fOuterRadius(0.), fDPhi(0.),
199  fZHalfLength(0.), fInnerStereo(0.), fOuterStereo(0.), fTanInnerStereo(0.),
200  fTanOuterStereo(0.), fKappa(0.), fInnerRadius2(0.), fOuterRadius2(0.),
201  fTanInnerStereo2(0.), fTanOuterStereo2(0.), fLowerEndcap(0), fUpperEndcap(0),
202  fLatterTwisted(0), fFormerTwisted(0), fInnerHype(0), fOuterHype(0),
203  fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0)
204 {
205  fEndZ[0] = 0.; fEndZ[1] = 0.;
206  fEndInnerRadius[0] = 0.; fEndInnerRadius[1] = 0.;
207  fEndOuterRadius[0] = 0.; fEndOuterRadius[1] = 0.;
208  fEndPhi[0] = 0.; fEndPhi[1] = 0.;
209  fEndZ2[0] = 0.; fEndZ2[1] = 0.;
210 }
211 
212 //=====================================================================
213 //* destructor --------------------------------------------------------
214 
215 G4TwistedTubs::~G4TwistedTubs()
216 {
217  if (fLowerEndcap) { delete fLowerEndcap; }
218  if (fUpperEndcap) { delete fUpperEndcap; }
219  if (fLatterTwisted) { delete fLatterTwisted; }
220  if (fFormerTwisted) { delete fFormerTwisted; }
221  if (fInnerHype) { delete fInnerHype; }
222  if (fOuterHype) { delete fOuterHype; }
223  if (fpPolyhedron) { delete fpPolyhedron; }
224 }
225 
226 //=====================================================================
227 //* Copy constructor --------------------------------------------------
228 
229 G4TwistedTubs::G4TwistedTubs(const G4TwistedTubs& rhs)
230  : G4VSolid(rhs), fPhiTwist(rhs.fPhiTwist),
231  fInnerRadius(rhs.fInnerRadius), fOuterRadius(rhs.fOuterRadius),
232  fDPhi(rhs.fDPhi), fZHalfLength(rhs.fZHalfLength),
233  fInnerStereo(rhs.fInnerStereo), fOuterStereo(rhs.fOuterStereo),
234  fTanInnerStereo(rhs.fTanInnerStereo), fTanOuterStereo(rhs.fTanOuterStereo),
235  fKappa(rhs.fKappa), fInnerRadius2(rhs.fInnerRadius2),
236  fOuterRadius2(rhs.fOuterRadius2), fTanInnerStereo2(rhs.fTanInnerStereo2),
237  fTanOuterStereo2(rhs.fTanOuterStereo2),
238  fLowerEndcap(0), fUpperEndcap(0), fLatterTwisted(0), fFormerTwisted(0),
239  fInnerHype(0), fOuterHype(0),
240  fCubicVolume(rhs.fCubicVolume), fSurfaceArea(rhs.fSurfaceArea),
241  fpPolyhedron(0), fLastInside(rhs.fLastInside), fLastNormal(rhs.fLastNormal),
242  fLastDistanceToIn(rhs.fLastDistanceToIn),
243  fLastDistanceToOut(rhs.fLastDistanceToOut),
244  fLastDistanceToInWithV(rhs.fLastDistanceToInWithV),
245  fLastDistanceToOutWithV(rhs.fLastDistanceToOutWithV)
246 {
247  for (size_t i=0; i<2; ++i)
248  {
249  fEndZ[i] = rhs.fEndZ[i];
250  fEndInnerRadius[i] = rhs.fEndInnerRadius[i];
251  fEndOuterRadius[i] = rhs.fEndOuterRadius[i];
252  fEndPhi[i] = rhs.fEndPhi[i];
253  fEndZ2[i] = rhs.fEndZ2[i];
254  }
255  CreateSurfaces();
256 }
257 
258 
259 //=====================================================================
260 //* Assignment operator -----------------------------------------------
261 
262 G4TwistedTubs& G4TwistedTubs::operator = (const G4TwistedTubs& rhs)
263 {
264  // Check assignment to self
265  //
266  if (this == &rhs) { return *this; }
267 
268  // Copy base class data
269  //
270  G4VSolid::operator=(rhs);
271 
272  // Copy data
273  //
274  fPhiTwist= rhs.fPhiTwist;
275  fInnerRadius= rhs.fInnerRadius; fOuterRadius= rhs.fOuterRadius;
276  fDPhi= rhs.fDPhi; fZHalfLength= rhs.fZHalfLength;
277  fInnerStereo= rhs.fInnerStereo; fOuterStereo= rhs.fOuterStereo;
278  fTanInnerStereo= rhs.fTanInnerStereo; fTanOuterStereo= rhs.fTanOuterStereo;
279  fKappa= rhs.fKappa; fInnerRadius2= rhs.fInnerRadius2;
280  fOuterRadius2= rhs.fOuterRadius2; fTanInnerStereo2= rhs.fTanInnerStereo2;
281  fTanOuterStereo2= rhs.fTanOuterStereo2;
282  fLowerEndcap= fUpperEndcap= fLatterTwisted= fFormerTwisted= 0;
283  fInnerHype= fOuterHype= 0;
284  fCubicVolume= rhs.fCubicVolume; fSurfaceArea= rhs.fSurfaceArea;
285  fpPolyhedron= 0;
286  fLastInside= rhs.fLastInside; fLastNormal= rhs.fLastNormal;
287  fLastDistanceToIn= rhs.fLastDistanceToIn;
288  fLastDistanceToOut= rhs.fLastDistanceToOut;
289  fLastDistanceToInWithV= rhs.fLastDistanceToInWithV;
290  fLastDistanceToOutWithV= rhs.fLastDistanceToOutWithV;
291 
292  for (size_t i=0; i<2; ++i)
293  {
294  fEndZ[i] = rhs.fEndZ[i];
295  fEndInnerRadius[i] = rhs.fEndInnerRadius[i];
296  fEndOuterRadius[i] = rhs.fEndOuterRadius[i];
297  fEndPhi[i] = rhs.fEndPhi[i];
298  fEndZ2[i] = rhs.fEndZ2[i];
299  }
300 
301  CreateSurfaces();
302 
303  return *this;
304 }
305 
306 //=====================================================================
307 //* ComputeDimensions -------------------------------------------------
308 
309 void G4TwistedTubs::ComputeDimensions(G4VPVParameterisation* /* p */ ,
310  const G4int /* n */ ,
311  const G4VPhysicalVolume* /* pRep */ )
312 {
313  G4Exception("G4TwistedTubs::ComputeDimensions()",
314  "GeomSolids0001", FatalException,
315  "G4TwistedTubs does not support Parameterisation.");
316 }
317 
318 
319 //=====================================================================
320 //* CalculateExtent ---------------------------------------------------
321 
322 G4bool G4TwistedTubs::CalculateExtent( const EAxis axis,
323  const G4VoxelLimits &voxelLimit,
324  const G4AffineTransform &transform,
325  G4double &min,
326  G4double &max ) const
327 {
328 
329  G4SolidExtentList extentList( axis, voxelLimit );
330  G4double maxEndOuterRad = (fEndOuterRadius[0] > fEndOuterRadius[1] ?
331  fEndOuterRadius[0] : fEndOuterRadius[1]);
332  G4double maxEndInnerRad = (fEndInnerRadius[0] > fEndInnerRadius[1] ?
333  fEndInnerRadius[0] : fEndInnerRadius[1]);
334  G4double maxphi = (std::fabs(fEndPhi[0]) > std::fabs(fEndPhi[1]) ?
335  std::fabs(fEndPhi[0]) : std::fabs(fEndPhi[1]));
336 
337  //
338  // Choose phi size of our segment(s) based on constants as
339  // defined in meshdefs.hh
340  //
341  // G4int numPhi = kMaxMeshSections;
342  G4double sigPhi = 2*maxphi + fDPhi;
343  G4double rFudge = 1.0/std::cos(0.5*sigPhi);
344  G4double fudgeEndOuterRad = rFudge * maxEndOuterRad;
345 
346  //
347  // We work around in phi building polygons along the way.
348  // As a reasonable compromise between accuracy and
349  // complexity (=cpu time), the following facets are chosen:
350  //
351  // 1. If fOuterRadius/maxEndOuterRad > 0.95, approximate
352  // the outer surface as a cylinder, and use one
353  // rectangular polygon (0-1) to build its mesh.
354  //
355  // Otherwise, use two trapazoidal polygons that
356  // meet at z = 0 (0-4-1)
357  //
358  // 2. If there is no inner surface, then use one
359  // polygon for each entire endcap. (0) and (1)
360  //
361  // Otherwise, use a trapazoidal polygon for each
362  // phi segment of each endcap. (0-2) and (1-3)
363  //
364  // 3. For the inner surface, if fInnerRadius/maxEndInnerRad > 0.95,
365  // approximate the inner surface as a cylinder of
366  // radius fInnerRadius and use one rectangular polygon
367  // to build each phi segment of its mesh. (2-3)
368  //
369  // Otherwise, use one rectangular polygon centered
370  // at z = 0 (5-6) and two connecting trapazoidal polygons
371  // for each phi segment (2-5) and (3-6).
372  //
373 
374  G4bool splitOuter = (fOuterRadius/maxEndOuterRad < 0.95);
375  G4bool splitInner = (fInnerRadius/maxEndInnerRad < 0.95);
376 
377  //
378  // Vertex assignments (v and w arrays)
379  // [0] and [1] are mandatory
380  // the rest are optional
381  //
382  // + -
383  // [0]------[4]------[1] <--- outer radius
384  // | |
385  // | |
386  // [2]---[5]---[6]---[3] <--- inner radius
387  //
388 
389  G4ClippablePolygon endPoly1, endPoly2;
390 
391  G4double phimax = maxphi + 0.5*fDPhi;
392  if ( phimax > pi/2) { phimax = pi-phimax; }
393  G4double phimin = - phimax;
394 
395  G4ThreeVector v0, v1, v2, v3, v4, v5, v6; // -ve phi verticies for polygon
396  G4ThreeVector w0, w1, w2, w3, w4, w5, w6; // +ve phi verticies for polygon
397 
398  //
399  // decide verticies of -ve phi boundary
400  //
401 
402  G4double cosPhi = std::cos(phimin);
403  G4double sinPhi = std::sin(phimin);
404 
405  // Outer hyperbolic surface
406 
407  v0 = transform.TransformPoint(
408  G4ThreeVector(fudgeEndOuterRad*cosPhi, fudgeEndOuterRad*sinPhi,
409  + fZHalfLength));
410  v1 = transform.TransformPoint(
411  G4ThreeVector(fudgeEndOuterRad*cosPhi, fudgeEndOuterRad*sinPhi,
412  - fZHalfLength));
413  if (splitOuter)
414  {
415  v4 = transform.TransformPoint(
416  G4ThreeVector(fudgeEndOuterRad*cosPhi, fudgeEndOuterRad*sinPhi, 0));
417  }
418 
419  // Inner hyperbolic surface
420 
421  G4double zInnerSplit = 0.;
422  if (splitInner)
423  {
424  v2 = transform.TransformPoint(
425  G4ThreeVector(maxEndInnerRad*cosPhi, maxEndInnerRad*sinPhi,
426  + fZHalfLength));
427  v3 = transform.TransformPoint(
428  G4ThreeVector(maxEndInnerRad*cosPhi, maxEndInnerRad*sinPhi,
429  - fZHalfLength));
430 
431  // Find intersection of tangential line of inner
432  // surface at z = fZHalfLength and line r=fInnerRadius.
433  G4double dr = fZHalfLength * fTanInnerStereo2;
434  G4double dz = maxEndInnerRad;
435  zInnerSplit = fZHalfLength + (fInnerRadius - maxEndInnerRad) * dz / dr;
436 
437  // Build associated vertices
438  v5 = transform.TransformPoint(
439  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
440  + zInnerSplit));
441  v6 = transform.TransformPoint(
442  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
443  - zInnerSplit));
444  }
445  else
446  {
447  v2 = transform.TransformPoint(
448  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
449  + fZHalfLength));
450  v3 = transform.TransformPoint(
451  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
452  - fZHalfLength));
453  }
454 
455  //
456  // decide vertices of +ve phi boundary
457  //
458 
459  cosPhi = std::cos(phimax);
460  sinPhi = std::sin(phimax);
461 
462  // Outer hyperbolic surface
463 
464  w0 = transform.TransformPoint(
465  G4ThreeVector(fudgeEndOuterRad*cosPhi, fudgeEndOuterRad*sinPhi,
466  + fZHalfLength));
467  w1 = transform.TransformPoint(
468  G4ThreeVector(fudgeEndOuterRad*cosPhi, fudgeEndOuterRad*sinPhi,
469  - fZHalfLength));
470  if (splitOuter)
471  {
472  G4double r = rFudge*fOuterRadius;
473 
474  w4 = transform.TransformPoint(G4ThreeVector( r*cosPhi, r*sinPhi, 0 ));
475 
476  AddPolyToExtent( v0, v4, w4, w0, voxelLimit, axis, extentList );
477  AddPolyToExtent( v4, v1, w1, w4, voxelLimit, axis, extentList );
478  }
479  else
480  {
481  AddPolyToExtent( v0, v1, w1, w0, voxelLimit, axis, extentList );
482  }
483 
484  // Inner hyperbolic surface
485 
486  if (splitInner)
487  {
488  w2 = transform.TransformPoint(
489  G4ThreeVector(maxEndInnerRad*cosPhi, maxEndInnerRad*sinPhi,
490  + fZHalfLength));
491  w3 = transform.TransformPoint(
492  G4ThreeVector(maxEndInnerRad*cosPhi, maxEndInnerRad*sinPhi,
493  - fZHalfLength));
494 
495  w5 = transform.TransformPoint(
496  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
497  + zInnerSplit));
498  w6 = transform.TransformPoint(
499  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
500  - zInnerSplit));
501 
502  AddPolyToExtent( v3, v6, w6, w3, voxelLimit, axis, extentList );
503  AddPolyToExtent( v6, v5, w5, w6, voxelLimit, axis, extentList );
504  AddPolyToExtent( v5, v2, w2, w5, voxelLimit, axis, extentList );
505 
506  }
507  else
508  {
509  w2 = transform.TransformPoint(
510  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
511  + fZHalfLength));
512  w3 = transform.TransformPoint(
513  G4ThreeVector(fInnerRadius*cosPhi, fInnerRadius*sinPhi,
514  - fZHalfLength));
515 
516  AddPolyToExtent( v3, v2, w2, w3, voxelLimit, axis, extentList );
517  }
518 
519  //
520  // Endplate segments
521  //
522  AddPolyToExtent( v1, v3, w3, w1, voxelLimit, axis, extentList );
523  AddPolyToExtent( v2, v0, w0, w2, voxelLimit, axis, extentList );
524 
525  //
526  // Return min/max value
527  //
528  return extentList.GetExtent( min, max );
529 }
530 
531 
532 //=====================================================================
533 //* AddPolyToExtent ---------------------------------------------------
534 
535 void G4TwistedTubs::AddPolyToExtent( const G4ThreeVector &v0,
536  const G4ThreeVector &v1,
537  const G4ThreeVector &w1,
538  const G4ThreeVector &w0,
539  const G4VoxelLimits &voxelLimit,
540  const EAxis axis,
541  G4SolidExtentList &extentList )
542 {
543  // Utility function for CalculateExtent
544  //
545  G4ClippablePolygon phiPoly;
546 
547  phiPoly.AddVertexInOrder( v0 );
548  phiPoly.AddVertexInOrder( v1 );
549  phiPoly.AddVertexInOrder( w1 );
550  phiPoly.AddVertexInOrder( w0 );
551 
552  if (phiPoly.PartialClip( voxelLimit, axis ))
553  {
554  phiPoly.SetNormal( (v1-v0).cross(w0-v0).unit() );
555  extentList.AddSurface( phiPoly );
556  }
557 }
558 
559 
560 //=====================================================================
561 //* Inside ------------------------------------------------------------
562 
563 EInside G4TwistedTubs::Inside(const G4ThreeVector& p) const
564 {
565 
566  const G4double halftol
567  = 0.5 * G4GeometryTolerance::GetInstance()->GetRadialTolerance();
568  // static G4int timerid = -1;
569  // G4Timer timer(timerid, "G4TwistedTubs", "Inside");
570  // timer.Start();
571 
572  G4ThreeVector *tmpp;
573  EInside *tmpinside;
574  if (fLastInside.p == p)
575  {
576  return fLastInside.inside;
577  }
578  else
579  {
580  tmpp = const_cast<G4ThreeVector*>(&(fLastInside.p));
581  tmpinside = const_cast<EInside*>(&(fLastInside.inside));
582  tmpp->set(p.x(), p.y(), p.z());
583  }
584 
585  EInside outerhypearea = ((G4TwistTubsHypeSide *)fOuterHype)->Inside(p);
586  G4double innerhyperho = ((G4TwistTubsHypeSide *)fInnerHype)->GetRhoAtPZ(p);
587  G4double distanceToOut = p.getRho() - innerhyperho; // +ve: inside
588 
589  if ((outerhypearea == kOutside) || (distanceToOut < -halftol))
590  {
591  *tmpinside = kOutside;
592  }
593  else if (outerhypearea == kSurface)
594  {
595  *tmpinside = kSurface;
596  }
597  else
598  {
599  if (distanceToOut <= halftol)
600  {
601  *tmpinside = kSurface;
602  }
603  else
604  {
605  *tmpinside = kInside;
606  }
607  }
608 
609  return fLastInside.inside;
610 }
611 
612 //=====================================================================
613 //* SurfaceNormal -----------------------------------------------------
614 
615 G4ThreeVector G4TwistedTubs::SurfaceNormal(const G4ThreeVector& p) const
616 {
617  //
618  // return the normal unit vector to the Hyperbolical Surface at a point
619  // p on (or nearly on) the surface
620  //
621  // Which of the three or four surfaces are we closest to?
622  //
623 
624  if (fLastNormal.p == p)
625  {
626  return fLastNormal.vec;
627  }
628  G4ThreeVector *tmpp =
629  const_cast<G4ThreeVector*>(&(fLastNormal.p));
630  G4ThreeVector *tmpnormal =
631  const_cast<G4ThreeVector*>(&(fLastNormal.vec));
632  G4VTwistSurface **tmpsurface =
633  const_cast<G4VTwistSurface**>(fLastNormal.surface);
634  tmpp->set(p.x(), p.y(), p.z());
635 
636  G4double distance = kInfinity;
637 
638  G4VTwistSurface *surfaces[6];
639  surfaces[0] = fLatterTwisted;
640  surfaces[1] = fFormerTwisted;
641  surfaces[2] = fInnerHype;
642  surfaces[3] = fOuterHype;
643  surfaces[4] = fLowerEndcap;
644  surfaces[5] = fUpperEndcap;
645 
646  G4ThreeVector xx;
647  G4ThreeVector bestxx;
648  G4int i;
649  G4int besti = -1;
650  for (i=0; i< 6; i++)
651  {
652  G4double tmpdistance = surfaces[i]->DistanceTo(p, xx);
653  if (tmpdistance < distance)
654  {
655  distance = tmpdistance;
656  bestxx = xx;
657  besti = i;
658  }
659  }
660 
661  tmpsurface[0] = surfaces[besti];
662  *tmpnormal = tmpsurface[0]->GetNormal(bestxx, true);
663 
664  return fLastNormal.vec;
665 }
666 
667 //=====================================================================
668 //* DistanceToIn (p, v) -----------------------------------------------
669 
670 G4double G4TwistedTubs::DistanceToIn (const G4ThreeVector& p,
671  const G4ThreeVector& v ) const
672 {
673 
674  // DistanceToIn (p, v):
675  // Calculate distance to surface of shape from `outside'
676  // along with the v, allowing for tolerance.
677  // The function returns kInfinity if no intersection or
678  // just grazing within tolerance.
679 
680  //
681  // checking last value
682  //
683 
684  G4ThreeVector *tmpp;
685  G4ThreeVector *tmpv;
686  G4double *tmpdist;
687  if ((fLastDistanceToInWithV.p == p) && (fLastDistanceToInWithV.vec == v))
688  {
689  return fLastDistanceToIn.value;
690  }
691  else
692  {
693  tmpp = const_cast<G4ThreeVector*>(&(fLastDistanceToInWithV.p));
694  tmpv = const_cast<G4ThreeVector*>(&(fLastDistanceToInWithV.vec));
695  tmpdist = const_cast<G4double*>(&(fLastDistanceToInWithV.value));
696  tmpp->set(p.x(), p.y(), p.z());
697  tmpv->set(v.x(), v.y(), v.z());
698  }
699 
700  //
701  // Calculate DistanceToIn(p,v)
702  //
703 
704  EInside currentside = Inside(p);
705 
706  if (currentside == kInside)
707  {
708  }
709  else
710  {
711  if (currentside == kSurface)
712  {
713  // particle is just on a boundary.
714  // If the particle is entering to the volume, return 0.
715  //
716  G4ThreeVector normal = SurfaceNormal(p);
717  if (normal*v < 0)
718  {
719  *tmpdist = 0;
720  return fLastDistanceToInWithV.value;
721  }
722  }
723  }
724 
725  // now, we can take smallest positive distance.
726 
727  // Initialize
728  //
729  G4double distance = kInfinity;
730 
731  // find intersections and choose nearest one.
732  //
733  G4VTwistSurface *surfaces[6];
734  surfaces[0] = fLowerEndcap;
735  surfaces[1] = fUpperEndcap;
736  surfaces[2] = fLatterTwisted;
737  surfaces[3] = fFormerTwisted;
738  surfaces[4] = fInnerHype;
739  surfaces[5] = fOuterHype;
740 
741  G4ThreeVector xx;
742  G4ThreeVector bestxx;
743  G4int i;
744  for (i=0; i< 6; i++)
745  {
746  G4double tmpdistance = surfaces[i]->DistanceToIn(p, v, xx);
747  if (tmpdistance < distance)
748  {
749  distance = tmpdistance;
750  bestxx = xx;
751  }
752  }
753  *tmpdist = distance;
754 
755  return fLastDistanceToInWithV.value;
756 }
757 
758 //=====================================================================
759 //* DistanceToIn (p) --------------------------------------------------
760 
761 G4double G4TwistedTubs::DistanceToIn (const G4ThreeVector& p) const
762 {
763  // DistanceToIn(p):
764  // Calculate distance to surface of shape from `outside',
765  // allowing for tolerance
766 
767  //
768  // checking last value
769  //
770 
771  G4ThreeVector *tmpp;
772  G4double *tmpdist;
773  if (fLastDistanceToIn.p == p)
774  {
775  return fLastDistanceToIn.value;
776  }
777  else
778  {
779  tmpp = const_cast<G4ThreeVector*>(&(fLastDistanceToIn.p));
780  tmpdist = const_cast<G4double*>(&(fLastDistanceToIn.value));
781  tmpp->set(p.x(), p.y(), p.z());
782  }
783 
784  //
785  // Calculate DistanceToIn(p)
786  //
787 
788  EInside currentside = Inside(p);
789 
790  switch (currentside)
791  {
792  case (kInside) :
793  {}
794  case (kSurface) :
795  {
796  *tmpdist = 0.;
797  return fLastDistanceToIn.value;
798  }
799  case (kOutside) :
800  {
801  // Initialize
802  G4double distance = kInfinity;
803 
804  // find intersections and choose nearest one.
805  G4VTwistSurface *surfaces[6];
806  surfaces[0] = fLowerEndcap;
807  surfaces[1] = fUpperEndcap;
808  surfaces[2] = fLatterTwisted;
809  surfaces[3] = fFormerTwisted;
810  surfaces[4] = fInnerHype;
811  surfaces[5] = fOuterHype;
812 
813  G4int i;
814  G4ThreeVector xx;
815  G4ThreeVector bestxx;
816  for (i=0; i< 6; i++)
817  {
818  G4double tmpdistance = surfaces[i]->DistanceTo(p, xx);
819  if (tmpdistance < distance)
820  {
821  distance = tmpdistance;
822  bestxx = xx;
823  }
824  }
825  *tmpdist = distance;
826  return fLastDistanceToIn.value;
827  }
828  default :
829  {
830  G4Exception("G4TwistedTubs::DistanceToIn(p)", "GeomSolids0003",
831  FatalException, "Unknown point location!");
832  }
833  } // switch end
834 
835  return kInfinity;
836 }
837 
838 //=====================================================================
839 //* DistanceToOut (p, v) ----------------------------------------------
840 
841 G4double G4TwistedTubs::DistanceToOut( const G4ThreeVector& p,
842  const G4ThreeVector& v,
843  const G4bool calcNorm,
844  G4bool *validNorm,
845  G4ThreeVector *norm ) const
846 {
847  // DistanceToOut (p, v):
848  // Calculate distance to surface of shape from `inside'
849  // along with the v, allowing for tolerance.
850  // The function returns kInfinity if no intersection or
851  // just grazing within tolerance.
852 
853  //
854  // checking last value
855  //
856 
857  G4ThreeVector *tmpp;
858  G4ThreeVector *tmpv;
859  G4double *tmpdist;
860  if ((fLastDistanceToOutWithV.p == p) && (fLastDistanceToOutWithV.vec == v) )
861  {
862  return fLastDistanceToOutWithV.value;
863  }
864  else
865  {
866  tmpp = const_cast<G4ThreeVector*>(&(fLastDistanceToOutWithV.p));
867  tmpv = const_cast<G4ThreeVector*>(&(fLastDistanceToOutWithV.vec));
868  tmpdist = const_cast<G4double*>(&(fLastDistanceToOutWithV.value));
869  tmpp->set(p.x(), p.y(), p.z());
870  tmpv->set(v.x(), v.y(), v.z());
871  }
872 
873  //
874  // Calculate DistanceToOut(p,v)
875  //
876 
877  EInside currentside = Inside(p);
878 
879  if (currentside == kOutside)
880  {
881  }
882  else
883  {
884  if (currentside == kSurface)
885  {
886  // particle is just on a boundary.
887  // If the particle is exiting from the volume, return 0.
888  //
889  G4ThreeVector normal = SurfaceNormal(p);
890  G4VTwistSurface *blockedsurface = fLastNormal.surface[0];
891  if (normal*v > 0)
892  {
893  if (calcNorm)
894  {
895  *norm = (blockedsurface->GetNormal(p, true));
896  *validNorm = blockedsurface->IsValidNorm();
897  }
898  *tmpdist = 0.;
899  return fLastDistanceToOutWithV.value;
900  }
901  }
902  }
903 
904  // now, we can take smallest positive distance.
905 
906  // Initialize
907  //
908  G4double distance = kInfinity;
909 
910  // find intersections and choose nearest one.
911  //
912  G4VTwistSurface *surfaces[6];
913  surfaces[0] = fLatterTwisted;
914  surfaces[1] = fFormerTwisted;
915  surfaces[2] = fInnerHype;
916  surfaces[3] = fOuterHype;
917  surfaces[4] = fLowerEndcap;
918  surfaces[5] = fUpperEndcap;
919 
920  G4int i;
921  G4int besti = -1;
922  G4ThreeVector xx;
923  G4ThreeVector bestxx;
924  for (i=0; i< 6; i++)
925  {
926  G4double tmpdistance = surfaces[i]->DistanceToOut(p, v, xx);
927  if (tmpdistance < distance)
928  {
929  distance = tmpdistance;
930  bestxx = xx;
931  besti = i;
932  }
933  }
934 
935  if (calcNorm)
936  {
937  if (besti != -1)
938  {
939  *norm = (surfaces[besti]->GetNormal(p, true));
940  *validNorm = surfaces[besti]->IsValidNorm();
941  }
942  }
943 
944  *tmpdist = distance;
945 
946  return fLastDistanceToOutWithV.value;
947 }
948 
949 
950 //=====================================================================
951 //* DistanceToOut (p) ----------------------------------------------
952 
953 G4double G4TwistedTubs::DistanceToOut( const G4ThreeVector& p ) const
954 {
955  // DistanceToOut(p):
956  // Calculate distance to surface of shape from `inside',
957  // allowing for tolerance
958 
959  //
960  // checking last value
961  //
962 
963  G4ThreeVector *tmpp;
964  G4double *tmpdist;
965  if (fLastDistanceToOut.p == p)
966  {
967  return fLastDistanceToOut.value;
968  }
969  else
970  {
971  tmpp = const_cast<G4ThreeVector*>(&(fLastDistanceToOut.p));
972  tmpdist = const_cast<G4double*>(&(fLastDistanceToOut.value));
973  tmpp->set(p.x(), p.y(), p.z());
974  }
975 
976  //
977  // Calculate DistanceToOut(p)
978  //
979 
980  EInside currentside = Inside(p);
981 
982  switch (currentside)
983  {
984  case (kOutside) :
985  {
986  }
987  case (kSurface) :
988  {
989  *tmpdist = 0.;
990  return fLastDistanceToOut.value;
991  }
992  case (kInside) :
993  {
994  // Initialize
995  G4double distance = kInfinity;
996 
997  // find intersections and choose nearest one.
998  G4VTwistSurface *surfaces[6];
999  surfaces[0] = fLatterTwisted;
1000  surfaces[1] = fFormerTwisted;
1001  surfaces[2] = fInnerHype;
1002  surfaces[3] = fOuterHype;
1003  surfaces[4] = fLowerEndcap;
1004  surfaces[5] = fUpperEndcap;
1005 
1006  G4int i;
1007  G4ThreeVector xx;
1008  G4ThreeVector bestxx;
1009  for (i=0; i< 6; i++)
1010  {
1011  G4double tmpdistance = surfaces[i]->DistanceTo(p, xx);
1012  if (tmpdistance < distance)
1013  {
1014  distance = tmpdistance;
1015  bestxx = xx;
1016  }
1017  }
1018  *tmpdist = distance;
1019 
1020  return fLastDistanceToOut.value;
1021  }
1022  default :
1023  {
1024  G4Exception("G4TwistedTubs::DistanceToOut(p)", "GeomSolids0003",
1025  FatalException, "Unknown point location!");
1026  }
1027  } // switch end
1028 
1029  return 0;
1030 }
1031 
1032 //=====================================================================
1033 //* StreamInfo --------------------------------------------------------
1034 
1035 std::ostream& G4TwistedTubs::StreamInfo(std::ostream& os) const
1036 {
1037  //
1038  // Stream object contents to an output stream
1039  //
1040  G4int oldprc = os.precision(16);
1041  os << "-----------------------------------------------------------\n"
1042  << " *** Dump for solid - " << GetName() << " ***\n"
1043  << " ===================================================\n"
1044  << " Solid type: G4TwistedTubs\n"
1045  << " Parameters: \n"
1046  << " -ve end Z : " << fEndZ[0]/mm << " mm \n"
1047  << " +ve end Z : " << fEndZ[1]/mm << " mm \n"
1048  << " inner end radius(-ve z): " << fEndInnerRadius[0]/mm << " mm \n"
1049  << " inner end radius(+ve z): " << fEndInnerRadius[1]/mm << " mm \n"
1050  << " outer end radius(-ve z): " << fEndOuterRadius[0]/mm << " mm \n"
1051  << " outer end radius(+ve z): " << fEndOuterRadius[1]/mm << " mm \n"
1052  << " inner radius (z=0) : " << fInnerRadius/mm << " mm \n"
1053  << " outer radius (z=0) : " << fOuterRadius/mm << " mm \n"
1054  << " twisted angle : " << fPhiTwist/degree << " degrees \n"
1055  << " inner stereo angle : " << fInnerStereo/degree << " degrees \n"
1056  << " outer stereo angle : " << fOuterStereo/degree << " degrees \n"
1057  << " phi-width of a piece : " << fDPhi/degree << " degrees \n"
1058  << "-----------------------------------------------------------\n";
1059  os.precision(oldprc);
1060 
1061  return os;
1062 }
1063 
1064 
1065 //=====================================================================
1066 //* DiscribeYourselfTo ------------------------------------------------
1067 
1068 void G4TwistedTubs::DescribeYourselfTo (G4VGraphicsScene& scene) const
1069 {
1070  scene.AddSolid (*this);
1071 }
1072 
1073 //=====================================================================
1074 //* GetExtent ---------------------------------------------------------
1075 
1076 G4VisExtent G4TwistedTubs::GetExtent() const
1077 {
1078  // Define the sides of the box into which the G4Tubs instance would fit.
1079 
1080  G4double maxEndOuterRad = (fEndOuterRadius[0] > fEndOuterRadius[1] ? 0 : 1);
1081  return G4VisExtent( -maxEndOuterRad, maxEndOuterRad,
1082  -maxEndOuterRad, maxEndOuterRad,
1083  -fZHalfLength, fZHalfLength );
1084 }
1085 
1086 //=====================================================================
1087 //* CreatePolyhedron --------------------------------------------------
1088 
1089 G4Polyhedron* G4TwistedTubs::CreatePolyhedron () const
1090 {
1091  // number of meshes
1092  //
1093  G4double dA = std::max(fDPhi,fPhiTwist);
1094  const G4int k =
1095  G4int(G4Polyhedron::GetNumberOfRotationSteps() * dA / twopi) + 2;
1096  const G4int n =
1097  G4int(G4Polyhedron::GetNumberOfRotationSteps() * fPhiTwist / twopi) + 2;
1098 
1099  const G4int nnodes = 4*(k-1)*(n-2) + 2*k*k ;
1100  const G4int nfaces = 4*(k-1)*(n-1) + 2*(k-1)*(k-1) ;
1101 
1102  G4Polyhedron *ph=new G4Polyhedron;
1103  typedef G4double G4double3[3];
1104  typedef G4int G4int4[4];
1105  G4double3* xyz = new G4double3[nnodes]; // number of nodes
1106  G4int4* faces = new G4int4[nfaces] ; // number of faces
1107  fLowerEndcap->GetFacets(k,k,xyz,faces,0) ;
1108  fUpperEndcap->GetFacets(k,k,xyz,faces,1) ;
1109  fInnerHype->GetFacets(k,n,xyz,faces,2) ;
1110  fFormerTwisted->GetFacets(k,n,xyz,faces,3) ;
1111  fOuterHype->GetFacets(k,n,xyz,faces,4) ;
1112  fLatterTwisted->GetFacets(k,n,xyz,faces,5) ;
1113 
1114  ph->createPolyhedron(nnodes,nfaces,xyz,faces);
1115 
1116  delete[] xyz;
1117  delete[] faces;
1118 
1119  return ph;
1120 }
1121 
1122 //=====================================================================
1123 //* GetPolyhedron -----------------------------------------------------
1124 
1125 G4Polyhedron* G4TwistedTubs::GetPolyhedron () const
1126 {
1127  if ((!fpPolyhedron) ||
1128  (fpPolyhedron->GetNumberOfRotationStepsAtTimeOfCreation() !=
1129  fpPolyhedron->GetNumberOfRotationSteps()))
1130  {
1131  delete fpPolyhedron;
1132  fpPolyhedron = CreatePolyhedron();
1133  }
1134  return fpPolyhedron;
1135 }
1136 
1137 //=====================================================================
1138 //* CreateSurfaces ----------------------------------------------------
1139 
1140 void G4TwistedTubs::CreateSurfaces()
1141 {
1142  // create 6 surfaces of TwistedTub
1143 
1144  G4ThreeVector x0(0, 0, fEndZ[0]);
1145  G4ThreeVector n (0, 0, -1);
1146 
1147  fLowerEndcap = new G4TwistTubsFlatSide("LowerEndcap",
1148  fEndInnerRadius, fEndOuterRadius,
1149  fDPhi, fEndPhi, fEndZ, -1) ;
1150 
1151  fUpperEndcap = new G4TwistTubsFlatSide("UpperEndcap",
1152  fEndInnerRadius, fEndOuterRadius,
1153  fDPhi, fEndPhi, fEndZ, 1) ;
1154 
1155  G4RotationMatrix rotHalfDPhi;
1156  rotHalfDPhi.rotateZ(0.5*fDPhi);
1157 
1158  fLatterTwisted = new G4TwistTubsSide("LatterTwisted",
1159  fEndInnerRadius, fEndOuterRadius,
1160  fDPhi, fEndPhi, fEndZ,
1161  fInnerRadius, fOuterRadius, fKappa,
1162  1 ) ;
1163  fFormerTwisted = new G4TwistTubsSide("FormerTwisted",
1164  fEndInnerRadius, fEndOuterRadius,
1165  fDPhi, fEndPhi, fEndZ,
1166  fInnerRadius, fOuterRadius, fKappa,
1167  -1 ) ;
1168 
1169  fInnerHype = new G4TwistTubsHypeSide("InnerHype",
1170  fEndInnerRadius, fEndOuterRadius,
1171  fDPhi, fEndPhi, fEndZ,
1172  fInnerRadius, fOuterRadius,fKappa,
1173  fTanInnerStereo, fTanOuterStereo, -1) ;
1174  fOuterHype = new G4TwistTubsHypeSide("OuterHype",
1175  fEndInnerRadius, fEndOuterRadius,
1176  fDPhi, fEndPhi, fEndZ,
1177  fInnerRadius, fOuterRadius,fKappa,
1178  fTanInnerStereo, fTanOuterStereo, 1) ;
1179 
1180 
1181  // set neighbour surfaces
1182  //
1183  fLowerEndcap->SetNeighbours(fInnerHype, fLatterTwisted,
1184  fOuterHype, fFormerTwisted);
1185  fUpperEndcap->SetNeighbours(fInnerHype, fLatterTwisted,
1186  fOuterHype, fFormerTwisted);
1187  fLatterTwisted->SetNeighbours(fInnerHype, fLowerEndcap,
1188  fOuterHype, fUpperEndcap);
1189  fFormerTwisted->SetNeighbours(fInnerHype, fLowerEndcap,
1190  fOuterHype, fUpperEndcap);
1191  fInnerHype->SetNeighbours(fLatterTwisted, fLowerEndcap,
1192  fFormerTwisted, fUpperEndcap);
1193  fOuterHype->SetNeighbours(fLatterTwisted, fLowerEndcap,
1194  fFormerTwisted, fUpperEndcap);
1195 }
1196 
1197 
1198 //=====================================================================
1199 //* GetEntityType -----------------------------------------------------
1200 
1201 G4GeometryType G4TwistedTubs::GetEntityType() const
1202 {
1203  return G4String("G4TwistedTubs");
1204 }
1205 
1206 //=====================================================================
1207 //* Clone -------------------------------------------------------------
1208 
1209 G4VSolid* G4TwistedTubs::Clone() const
1210 {
1211  return new G4TwistedTubs(*this);
1212 }
1213 
1214 //=====================================================================
1215 //* GetCubicVolume ----------------------------------------------------
1216 
1217 G4double G4TwistedTubs::GetCubicVolume()
1218 {
1219  if(fCubicVolume != 0.) {;}
1220  else { fCubicVolume = fDPhi*fZHalfLength*(fOuterRadius*fOuterRadius
1221  -fInnerRadius*fInnerRadius); }
1222  return fCubicVolume;
1223 }
1224 
1225 //=====================================================================
1226 //* GetSurfaceArea ----------------------------------------------------
1227 
1228 G4double G4TwistedTubs::GetSurfaceArea()
1229 {
1230  if(fSurfaceArea != 0.) {;}
1231  else { fSurfaceArea = G4VSolid::GetSurfaceArea(); }
1232  return fSurfaceArea;
1233 }
1234 
1235 //=====================================================================
1236 //* GetPointOnSurface -------------------------------------------------
1237 
1238 G4ThreeVector G4TwistedTubs::GetPointOnSurface() const
1239 {
1240 
1241  G4double z = G4RandFlat::shoot(fEndZ[0],fEndZ[1]);
1242  G4double phi , phimin, phimax ;
1243  G4double x , xmin, xmax ;
1244  G4double r , rmin, rmax ;
1245 
1246  G4double a1 = fOuterHype->GetSurfaceArea() ;
1247  G4double a2 = fInnerHype->GetSurfaceArea() ;
1248  G4double a3 = fLatterTwisted->GetSurfaceArea() ;
1249  G4double a4 = fFormerTwisted->GetSurfaceArea() ;
1250  G4double a5 = fLowerEndcap->GetSurfaceArea() ;
1251  G4double a6 = fUpperEndcap->GetSurfaceArea() ;
1252 
1253  G4double chose = G4RandFlat::shoot(0.,a1 + a2 + a3 + a4 + a5 + a6) ;
1254 
1255  if(chose < a1)
1256  {
1257 
1258  phimin = fOuterHype->GetBoundaryMin(z) ;
1259  phimax = fOuterHype->GetBoundaryMax(z) ;
1260  phi = G4RandFlat::shoot(phimin,phimax) ;
1261 
1262  return fOuterHype->SurfacePoint(phi,z,true) ;
1263 
1264  }
1265  else if ( (chose >= a1) && (chose < a1 + a2 ) )
1266  {
1267 
1268  phimin = fInnerHype->GetBoundaryMin(z) ;
1269  phimax = fInnerHype->GetBoundaryMax(z) ;
1270  phi = G4RandFlat::shoot(phimin,phimax) ;
1271 
1272  return fInnerHype->SurfacePoint(phi,z,true) ;
1273 
1274  }
1275  else if ( (chose >= a1 + a2 ) && (chose < a1 + a2 + a3 ) )
1276  {
1277 
1278  xmin = fLatterTwisted->GetBoundaryMin(z) ;
1279  xmax = fLatterTwisted->GetBoundaryMax(z) ;
1280  x = G4RandFlat::shoot(xmin,xmax) ;
1281 
1282  return fLatterTwisted->SurfacePoint(x,z,true) ;
1283 
1284  }
1285  else if ( (chose >= a1 + a2 + a3 ) && (chose < a1 + a2 + a3 + a4 ) )
1286  {
1287 
1288  xmin = fFormerTwisted->GetBoundaryMin(z) ;
1289  xmax = fFormerTwisted->GetBoundaryMax(z) ;
1290  x = G4RandFlat::shoot(xmin,xmax) ;
1291 
1292  return fFormerTwisted->SurfacePoint(x,z,true) ;
1293  }
1294  else if( (chose >= a1 + a2 + a3 + a4 )&&(chose < a1 + a2 + a3 + a4 + a5 ) )
1295  {
1296  rmin = GetEndInnerRadius(0) ;
1297  rmax = GetEndOuterRadius(0) ;
1298  r = std::sqrt(G4RandFlat::shoot()*(sqr(rmax)-sqr(rmin))+sqr(rmin));
1299 
1300  phimin = fLowerEndcap->GetBoundaryMin(r) ;
1301  phimax = fLowerEndcap->GetBoundaryMax(r) ;
1302  phi = G4RandFlat::shoot(phimin,phimax) ;
1303 
1304  return fLowerEndcap->SurfacePoint(phi,r,true) ;
1305  }
1306  else
1307  {
1308  rmin = GetEndInnerRadius(1) ;
1309  rmax = GetEndOuterRadius(1) ;
1310  r = rmin + (rmax-rmin)*std::sqrt(G4RandFlat::shoot());
1311 
1312  phimin = fUpperEndcap->GetBoundaryMin(r) ;
1313  phimax = fUpperEndcap->GetBoundaryMax(r) ;
1314  phi = G4RandFlat::shoot(phimin,phimax) ;
1315 
1316  return fUpperEndcap->SurfacePoint(phi,r,true) ;
1317  }
1318 }
CLHEP::Hep3Vector::set
void set(double x, double y, double z)
G4VSolid::GetName
G4String GetName() const
G4INCL::DeJongSpin::shoot
ThreeVector shoot(const G4int Ap, const G4int Af)
Definition: G4INCLDeJongSpin.cc:67
G4Polyhedron.hh
G4ClippablePolygon.hh
G4TwistedTubs.hh
G4TwistedTubs::GetEndInnerRadius
G4double GetEndInnerRadius() const
Definition: G4TwistedTubs.hh:158
test.a
tuple a
Definition: test.py:11
G4VGraphicsScene.hh
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
G4TwistedTubs::GetPointOnSurface
G4ThreeVector GetPointOnSurface() const
Definition: G4TwistedTubs.cc:1238
CLHEP::Hep3Vector::x
double x() const
G4VTwistSurface::GetNormal
virtual G4ThreeVector GetNormal(const G4ThreeVector &xx, G4bool isGlobal)=0
G4VoxelLimits.hh
HepPolyhedron::createPolyhedron
G4int createPolyhedron(G4int Nnodes, G4int Nfaces, const G4double xyz[][3], const G4int faces[][4])
Definition: HepPolyhedron.cc:1358
p
const char * p
Definition: xmltok.h:285
G4TwistedTubs::StreamInfo
std::ostream & StreamInfo(std::ostream &os) const
Definition: G4TwistedTubs.cc:1035
G4SolidExtentList::GetExtent
G4bool GetExtent(G4double &min, G4double &max) const
Definition: G4SolidExtentList.cc:131
norm
Float_t norm
Definition: comparison_ascii.C:71
G4TwistedTubs::CreatePolyhedron
G4Polyhedron * CreatePolyhedron() const
Definition: G4TwistedTubs.cc:1089
xx
Double_t xx
Definition: macro.C:10
G4ClippablePolygon::SetNormal
void SetNormal(const G4ThreeVector &newNormal)
G4ClippablePolygon::PartialClip
virtual G4bool PartialClip(const G4VoxelLimits &voxelLimit, const EAxis IgnoreMe)
Definition: G4ClippablePolygon.cc:101
G4ClippablePolygon::AddVertexInOrder
virtual void AddVertexInOrder(const G4ThreeVector vertex)
Definition: G4ClippablePolygon.cc:66
CLHEP::Hep3Vector::getRho
double getRho() const
test.x
tuple x
Definition: test.py:50
G4TwistedTubs::Inside
EInside Inside(const G4ThreeVector &p) const
Definition: G4TwistedTubs.cc:563
G4VGraphicsScene::AddSolid
virtual void AddSolid(const G4Box &)=0
G4int
int G4int
Definition: G4Types.hh:78
G4TwistedTubs::ComputeDimensions
void ComputeDimensions(G4VPVParameterisation *, const G4int, const G4VPhysicalVolume *)
Definition: G4TwistedTubs.cc:309
G4VTwistSurface::IsValidNorm
G4bool IsValidNorm() const
Definition: G4VTwistSurface.hh:156
CLHEP::Hep3Vector::z
double z() const
G4TwistedTubs::GetEntityType
G4GeometryType GetEntityType() const
Definition: G4TwistedTubs.cc:1201
G4VTwistSurface::DistanceToIn
virtual G4double DistanceToIn(const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector &gxxbest)
Definition: G4VTwistSurface.cc:289
G4TwistedTubs::GetEndOuterRadius
G4double GetEndOuterRadius() const
Definition: G4TwistedTubs.hh:161
G4TwistedTubs::SurfaceNormal
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
Definition: G4TwistedTubs.cc:615
G4VTwistSurface::GetSurfaceArea
virtual G4double GetSurfaceArea()=0
G4VisExtent.hh
n
Char_t n[5]
Definition: comparison_ascii.C:55
G4VTwistSurface::SetNeighbours
void SetNeighbours(G4VTwistSurface *axis0min, G4VTwistSurface *axis1min, G4VTwistSurface *axis0max, G4VTwistSurface *axis1max)
G4VTwistSurface::DistanceTo
virtual G4double DistanceTo(const G4ThreeVector &gp, G4ThreeVector &gxx)
Definition: G4VTwistSurface.cc:553
G4TwistedTubs::DescribeYourselfTo
void DescribeYourselfTo(G4VGraphicsScene &scene) const
Definition: G4TwistedTubs.cc:1068
python.hepunit.degree
tuple degree
Definition: hepunit.py:69
G4VTwistSurface::DistanceToOut
virtual G4double DistanceToOut(const G4ThreeVector &gp, const G4ThreeVector &gv, G4ThreeVector &gxxbest)
Definition: G4VTwistSurface.cc:478
G4VTwistSurface::SurfacePoint
virtual G4ThreeVector SurfacePoint(G4double, G4double, G4bool isGlobal=false)=0
G4GeometryTolerance::GetRadialTolerance
G4double GetRadialTolerance() const
Definition: G4GeometryTolerance.cc:97
G4bool
bool G4bool
Definition: G4Types.hh:79
G4VPVParameterisation.hh
G4SolidExtentList::AddSurface
void AddSurface(const G4ClippablePolygon &surface)
Definition: G4SolidExtentList.cc:91
G4SolidExtentList.hh
G4TwistedTubs::DistanceToOut
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcnorm=G4bool(false), G4bool *validnorm=0, G4ThreeVector *n=0) const
Definition: G4TwistedTubs.cc:841
Randomize.hh
G4TwistedTubs::Clone
G4VSolid * Clone() const
Definition: G4TwistedTubs.cc:1209
G4VTwistSurface::GetBoundaryMax
virtual G4double GetBoundaryMax(G4double)=0
test.v
tuple v
Definition: test.py:18
eplot.pname
string pname
Definition: eplot.py:33
G4TwistedTubs::GetSurfaceArea
G4double GetSurfaceArea()
Definition: G4TwistedTubs.cc:1228
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4PhysicalConstants.hh
G4VTwistSurface::GetFacets
virtual void GetFacets(G4int m, G4int n, G4double xyz[][3], G4int faces[][4], G4int iside)=0
G4VTwistSurface::GetBoundaryMin
virtual G4double GetBoundaryMin(G4double)=0
G4AffineTransform::TransformPoint
G4ThreeVector TransformPoint(const G4ThreeVector &vec) const
r
jump r
Definition: plot.C:36
FatalErrorInArgument
Definition: G4ExceptionSeverity.hh:61
G4AffineTransform.hh
G4TwistedTubs::DistanceToIn
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4TwistedTubs.cc:670
G4TwistedTubs::~G4TwistedTubs
virtual ~G4TwistedTubs()
Definition: G4TwistedTubs.cc:215
HepPolyhedron::GetNumberOfRotationSteps
static G4int GetNumberOfRotationSteps()
Definition: HepPolyhedron.cc:196
G4INCL::Math::max
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Definition: G4INCLGlobals.hh:116
EInside
EInside
Definition: geomdefs.hh:58
EAxis
EAxis
Definition: geomdefs.hh:54
CLHEP::Hep3Vector::y
double y() const
FatalException
Definition: G4ExceptionSeverity.hh:60
DBL_MIN
#define DBL_MIN
Definition: templates.hh:75
G4INCL::Math::min
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
Definition: G4INCLGlobals.hh:121
test.z
tuple z
Definition: test.py:28
CLHEP::HepRotation::rotateZ
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
G4TwistedTubs::G4TwistedTubs
G4TwistedTubs(const G4String &pname, G4double twistedangle, G4double endinnerrad, G4double endouterrad, G4double halfzlen, G4double dphi)
Definition: G4TwistedTubs.cc:65
G4endl
#define G4endl
Definition: G4ios.hh:61
G4VSolid::operator=
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:110
G4Polyhedron::GetNumberOfRotationStepsAtTimeOfCreation
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
Definition: G4Polyhedron.hh:123
sqr
T sqr(const T &x)
Definition: templates.hh:145
G4double
double G4double
Definition: G4Types.hh:76
G4SystemOfUnits.hh
G4TwistedTubs::GetCubicVolume
G4double GetCubicVolume()
Definition: G4TwistedTubs.cc:1217
G4GeometryTolerance.hh
G4TwistedTubs::GetExtent
G4VisExtent GetExtent() const
Definition: G4TwistedTubs.cc:1076
G4VSolid::GetSurfaceArea
virtual G4double GetSurfaceArea()
Definition: G4VSolid.cc:250
G4GeometryTolerance::GetInstance
static G4GeometryTolerance * GetInstance()
Definition: G4GeometryTolerance.cc:74
G4TwistedTubs::GetPolyhedron
G4Polyhedron * GetPolyhedron() const
Definition: G4TwistedTubs.cc:1125
G4TwistedTubs::CalculateExtent
G4bool CalculateExtent(const EAxis paxis, const G4VoxelLimits &pvoxellimit, const G4AffineTransform &ptransform, G4double &pmin, G4double &pmax) const
Definition: G4TwistedTubs.cc:322
G4TwistedTubs::operator=
G4TwistedTubs & operator=(const G4TwistedTubs &rhs)
Definition: G4TwistedTubs.cc:262
Generated on Sat Dec 14 2013 14:34:12 for Geant4_10 by   doxygen 1.8.5