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