ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4EllipticalTube.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4EllipticalTube.cc
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 // G4EllipticalTube implementation
27 //
28 // Author: David C. Williams (davidw@scipp.ucsc.edu)
29 // Revision: Evgueni Tcherniaev (evgueni.tcherniaev@cern.ch), 23.12.2019
30 // --------------------------------------------------------------------
31 
32 #include "G4EllipticalTube.hh"
33 
34 #if !(defined(G4GEOM_USE_UELLIPTICALTUBE) && defined(G4GEOM_USE_SYS_USOLIDS))
35 
36 #include "G4GeomTools.hh"
37 #include "G4RandomTools.hh"
38 #include "G4ClippablePolygon.hh"
39 #include "G4AffineTransform.hh"
40 #include "G4VoxelLimits.hh"
41 #include "G4BoundingEnvelope.hh"
42 
43 #include "Randomize.hh"
44 
45 #include "G4VGraphicsScene.hh"
46 #include "G4VisExtent.hh"
47 
48 #include "G4AutoLock.hh"
49 
50 namespace
51 {
52  G4Mutex polyhedronMutex = G4MUTEX_INITIALIZER;
53 }
54 
55 using namespace CLHEP;
56 
58 //
59 // Constructor
60 
62  G4double Dx,
63  G4double Dy,
64  G4double Dz )
65  : G4VSolid(name), fDx(Dx), fDy(Dy), fDz(Dz)
66 {
68 }
69 
71 //
72 // Fake default constructor - sets only member data and allocates memory
73 // for usage restricted to object persistency.
74 
76  : G4VSolid(a), halfTolerance(0.), fDx(0.), fDy(0.), fDz(0.),
77  fRsph(0.), fDDx(0.), fDDy(0.), fSx(0.), fSy(0.), fR(0.),
78  fQ1(0.), fQ2(0.), fScratch(0.)
79 {
80 }
81 
83 //
84 // Destructor
85 
87 {
88  delete fpPolyhedron; fpPolyhedron = nullptr;
89 }
90 
92 //
93 // Copy constructor
94 
96  : G4VSolid(rhs), halfTolerance(rhs.halfTolerance),
97  fDx(rhs.fDx), fDy(rhs.fDy), fDz(rhs.fDz),
98  fCubicVolume(rhs.fCubicVolume), fSurfaceArea(rhs.fSurfaceArea),
99  fRsph(rhs.fRsph), fDDx(rhs.fDDx), fDDy(rhs.fDDy),
100  fSx(rhs.fSx), fSy(rhs.fSy), fR(rhs.fR),
101  fQ1(rhs.fQ1), fQ2(rhs.fQ2), fScratch(rhs.fScratch)
102 {
103 }
104 
106 //
107 // Assignment operator
108 
110 {
111  // Check assignment to self
112  //
113  if (this == &rhs) { return *this; }
114 
115  // Copy base class data
116  //
117  G4VSolid::operator=(rhs);
118 
119  // Copy data
120  //
122  fDx = rhs.fDx;
123  fDy = rhs.fDy;
124  fDz = rhs.fDz;
127 
128  fRsph = rhs.fRsph;
129  fDDx = rhs.fDDx;
130  fDDy = rhs.fDDy;
131  fSx = rhs.fSx;
132  fSy = rhs.fSy;
133  fR = rhs.fR;
134  fQ1 = rhs.fQ1;
135  fQ2 = rhs.fQ2;
136  fScratch = rhs.fScratch;
137 
138  fRebuildPolyhedron = false;
139  delete fpPolyhedron; fpPolyhedron = nullptr;
140 
141  return *this;
142 }
143 
145 //
146 // Check dimensions
147 
149 {
150  // Check dimensions
151  //
152  halfTolerance = 0.5*kCarTolerance; // half tolerance
153  G4double dmin = 2*kCarTolerance;
154  if (fDx < dmin || fDy < dmin || fDz < dmin)
155  {
156  std::ostringstream message;
157  message << "Invalid (too small or negative) dimensions for Solid: "
158  << GetName()
159  << "\n Dx = " << fDx
160  << "\n Dy = " << fDy
161  << "\n Dz = " << fDz;
162  G4Exception("G4EllipticalTube::CheckParameters()", "GeomSolids0002",
163  FatalException, message);
164  }
165 
166  // Set pre-calculatated values
167  //
168  halfTolerance = 0.5*kCarTolerance; // half tolerance
169  fRsph = std::sqrt(fDx * fDx + fDy * fDy + fDz * fDz); // radius of surrounding sphere
170  fDDx = fDx * fDx; // X semi-axis squared
171  fDDy = fDy * fDy; // Y semi-axis squared
172 
173  fR = std::min(fDx, fDy); // resulting radius, after scaling elipse to circle
174  fSx = fR / fDx; // X scale factor
175  fSy = fR / fDy; // Y scale factor
176 
177  fQ1 = 0.5 / fR; // distance approxiamtion dist = Q1 * (x^2 + y^2) - Q2
178  fQ2 = 0.5 * (fR + halfTolerance * halfTolerance / fR);
179  fScratch = 2. * fR * fR * DBL_EPSILON; // scratch within calculation error thickness
180  // fScratch = (B * B / A) * (2. + halfTolerance / A) * halfTolerance; // alternative
181 }
182 
184 //
185 // Get bounding box
186 
188  G4ThreeVector& pMax ) const
189 {
190  pMin.set(-fDx,-fDy,-fDz);
191  pMax.set( fDx, fDy, fDz);
192 }
193 
195 //
196 // Calculate extent under transform and specified limit
197 
198 G4bool
200  const G4VoxelLimits& pVoxelLimit,
201  const G4AffineTransform& pTransform,
202  G4double& pMin, G4double& pMax ) const
203 {
204  G4ThreeVector bmin, bmax;
205  G4bool exist;
206 
207  // Check bounding box (bbox)
208  //
209  BoundingLimits(bmin,bmax);
210  G4BoundingEnvelope bbox(bmin,bmax);
211 #ifdef G4BBOX_EXTENT
212  return bbox.CalculateExtent(pAxis,pVoxelLimit, pTransform, pMin, pMax);
213 #endif
214  if (bbox.BoundingBoxVsVoxelLimits(pAxis, pVoxelLimit, pTransform, pMin, pMax))
215  {
216  return exist = (pMin < pMax) ? true : false;
217  }
218 
219  G4double dx = fDx;
220  G4double dy = fDy;
221  G4double dz = fDz;
222 
223  // Set bounding envelope (benv) and calculate extent
224  //
225  const G4int NSTEPS = 24; // number of steps for whole circle
226  G4double ang = twopi/NSTEPS;
227 
228  G4double sinHalf = std::sin(0.5*ang);
229  G4double cosHalf = std::cos(0.5*ang);
230  G4double sinStep = 2.*sinHalf*cosHalf;
231  G4double cosStep = 1. - 2.*sinHalf*sinHalf;
232  G4double sx = dx/cosHalf;
233  G4double sy = dy/cosHalf;
234 
235  G4double sinCur = sinHalf;
236  G4double cosCur = cosHalf;
237  G4ThreeVectorList baseA(NSTEPS),baseB(NSTEPS);
238  for (G4int k=0; k<NSTEPS; ++k)
239  {
240  baseA[k].set(sx*cosCur,sy*sinCur,-dz);
241  baseB[k].set(sx*cosCur,sy*sinCur, dz);
242 
243  G4double sinTmp = sinCur;
244  sinCur = sinCur*cosStep + cosCur*sinStep;
245  cosCur = cosCur*cosStep - sinTmp*sinStep;
246  }
247 
248  std::vector<const G4ThreeVectorList *> polygons(2);
249  polygons[0] = &baseA;
250  polygons[1] = &baseB;
251  G4BoundingEnvelope benv(bmin, bmax, polygons);
252  exist = benv.CalculateExtent(pAxis, pVoxelLimit, pTransform, pMin, pMax);
253  return exist;
254 }
255 
257 //
258 // Determine where is point: inside, outside or on surface
259 //
260 
262 {
263  G4double x = p.x() * fSx;
264  G4double y = p.y() * fSy;
265  G4double distR = fQ1 * (x * x + y * y) - fQ2;
266  G4double distZ = std::abs(p.z()) - fDz;
267  G4double dist = std::max(distR, distZ);
268 
269  if (dist > halfTolerance) return kOutside;
270  return (dist > -halfTolerance) ? kSurface : kInside;
271 }
272 
274 //
275 // Return unit normal at surface closest to p
276 
278 {
279  G4ThreeVector norm(0, 0, 0);
280  G4int nsurf = 0;
281 
282  // check lateral surface
283  G4double x = p.x() * fSx;
284  G4double y = p.y() * fSy;
285  G4double distR = fQ1 * (x * x + y * y) - fQ2;
286  if (std::abs(distR) <= halfTolerance)
287  {
288  norm = G4ThreeVector(p.x() * fDDy, p.y() * fDDx, 0.).unit();
289  ++nsurf;
290  }
291 
292  // check lateral bases
293  G4double distZ = std::abs(p.z()) - fDz;
294  if (std::abs(distZ) <= halfTolerance)
295  {
296  norm.setZ(p.z() < 0 ? -1. : 1.);
297  ++nsurf;
298  }
299 
300  // return normal
301  if (nsurf == 1) return norm;
302  else if (nsurf > 1) return norm.unit(); // edge
303  else
304  {
305  // Point is not on the surface
306  //
307 #ifdef G4SPECDEBUG
308  std::ostringstream message;
309  G4int oldprc = message.precision(16);
310  message << "Point p is not on surface (!?) of solid: "
311  << GetName() << G4endl;
312  message << "Position:\n";
313  message << " p.x() = " << p.x()/mm << " mm\n";
314  message << " p.y() = " << p.y()/mm << " mm\n";
315  message << " p.z() = " << p.z()/mm << " mm";
316  G4cout.precision(oldprc);
317  G4Exception("G4EllipticalTube::SurfaceNormal(p)", "GeomSolids1002",
318  JustWarning, message );
319  DumpInfo();
320 #endif
321  return ApproxSurfaceNormal(p);
322  }
323 }
324 
326 //
327 // Find surface nearest to point and return corresponding normal.
328 // The algorithm is similar to the algorithm used in Inside().
329 // This method normally should not be called.
330 
333 {
334  G4double x = p.x() * fSx;
335  G4double y = p.y() * fSy;
336  G4double distR = fQ1 * (x * x + y * y) - fQ2;
337  G4double distZ = std::abs(p.z()) - fDz;
338  if (distR > distZ && (x * x + y * y) > 0)
339  return G4ThreeVector(p.x() * fDDy, p.y() * fDDx, 0.).unit();
340  else
341  return G4ThreeVector(0, 0, (p.z() < 0 ? -1. : 1.));
342 }
343 
345 //
346 // Calculate distance to shape from outside, along normalised vector,
347 // return kInfinity if no intersection, or distance < halfTolerance
348 
350  const G4ThreeVector& v ) const
351 {
352  G4double offset = 0.;
353  G4ThreeVector pcur = p;
354 
355  // Check if point is flying away
356  //
357  G4double safex = std::abs(pcur.x()) - fDx;
358  G4double safey = std::abs(pcur.y()) - fDy;
359  G4double safez = std::abs(pcur.z()) - fDz;
360 
361  if (safez >= -halfTolerance && pcur.z() * v.z() >= 0.) return kInfinity;
362  if (safey >= -halfTolerance && pcur.y() * v.y() >= 0.) return kInfinity;
363  if (safex >= -halfTolerance && pcur.x() * v.x() >= 0.) return kInfinity;
364 
365  // Relocate point, if required
366  //
367  G4double Dmax = 32. * fRsph;
368  if (std::max(std::max(safex, safey), safez) > Dmax)
369  {
370  offset = (1. - 1.e-08) * pcur.mag() - 2. * fRsph;
371  pcur += offset * v;
372  G4double dist = DistanceToIn(pcur, v);
373  return (dist == kInfinity) ? kInfinity : dist + offset;
374  }
375 
376  // Scale elliptical tube to cylinder
377  //
378  G4double px = pcur.x() * fSx;
379  G4double py = pcur.y() * fSy;
380  G4double pz = pcur.z();
381  G4double vx = v.x() * fSx;
382  G4double vy = v.y() * fSy;
383  G4double vz = v.z();
384 
385  // Set coefficients of quadratic equation: A t^2 + 2B t + C = 0
386  //
387  G4double rr = px * px + py * py;
388  G4double A = vx * vx + vy * vy;
389  G4double B = px * vx + py * vy;
390  G4double C = rr - fR * fR;
391  G4double D = B * B - A * C;
392 
393  // Check if point is flying away relative to lateral surface
394  //
395  G4double distR = fQ1 * rr - fQ2;
396  G4bool parallelToZ = (A < DBL_EPSILON || std::abs(vz) >= 1.);
397  if (distR >= -halfTolerance && (B >= 0. || parallelToZ)) return kInfinity;
398 
399  // Find intersection with Z planes
400  //
401  G4double invz = (vz == 0) ? DBL_MAX : -1./vz;
402  G4double dz = std::copysign(fDz, invz);
403  G4double tzmin = (pz - dz) * invz;
404  G4double tzmax = (pz + dz) * invz;
405 
406  // Solve qudratic equation. There are two cases special where D <= 0:
407  // 1) trajectory parallel to Z axis (A = 0, B = 0, C - any, D = 0)
408  // 2) touch (D = 0) or no intersection (D < 0) with lateral surface
409  //
410  if (parallelToZ) return (tzmin<halfTolerance) ? offset : tzmin + offset; // 1)
411  if (D <= A * A * fScratch) return kInfinity; // 2)
412 
413  // Find roots of quadratic equation
414  G4double tmp = -B - std::copysign(std::sqrt(D), B);
415  G4double t1 = tmp / A;
416  G4double t2 = C / tmp;
417  G4double trmin = std::min(t1, t2);
418  G4double trmax = std::max(t1, t2);
419 
420  // Return distance
421  G4double tin = std::max(tzmin, trmin);
422  G4double tout = std::min(tzmax, trmax);
423 
424  if (tout <= tin + halfTolerance) return kInfinity; // touch or no hit
425  return (tin<halfTolerance) ? offset : tin + offset;
426 }
427 
429 //
430 // Estimate distance to the surface from outside,
431 // returns 0 if point is inside
432 
434 {
435  // safety distance to bounding box
436  G4double distX = std::abs(p.x()) - fDx;
437  G4double distY = std::abs(p.y()) - fDy;
438  G4double distZ = std::abs(p.z()) - fDz;
439  G4double distB = std::max(std::max(distX, distY), distZ);
440  // return (distB < 0) ? 0 : distB;
441 
442  // safety distance to lateral surface
443  G4double x = p.x() * fSx;
444  G4double y = p.y() * fSy;
445  G4double distR = std::sqrt(x * x + y * y) - fR;
446 
447  // return SafetyToIn
448  G4double dist = std::max(distB, distR);
449  return (dist < 0) ? 0 : dist;
450 }
451 
453 //
454 // Calculate distance to shape from inside and find normal
455 // at exit point, if required
456 // - when leaving the surface, return 0
457 
459  const G4ThreeVector& v,
460  const G4bool calcNorm,
461  G4bool* validNorm,
462  G4ThreeVector* n ) const
463 {
464  // Check if point flying away relative to Z planes
465  //
466  G4double pz = p.z();
467  G4double vz = v.z();
468  G4double distZ = std::abs(pz) - fDz;
469  if (distZ >= -halfTolerance && pz * vz > 0)
470  {
471  if (calcNorm)
472  {
473  *validNorm = true;
474  n->set(0, 0, (pz < 0) ? -1. : 1.);
475  }
476  return 0.;
477  }
478  G4double tzmax = (vz == 0) ? DBL_MAX : (std::copysign(fDz, vz) - pz) / vz;
479 
480  // Scale elliptical tube to cylinder
481  //
482  G4double px = p.x() * fSx;
483  G4double py = p.y() * fSy;
484  G4double vx = v.x() * fSx;
485  G4double vy = v.y() * fSy;
486 
487  // Check if point is flying away relative to lateral surface
488  //
489  G4double rr = px * px + py * py;
490  G4double B = px * vx + py * vy;
491  G4double distR = fQ1 * rr - fQ2;
492  if (distR >= -halfTolerance && B > 0.)
493  {
494  if (calcNorm)
495  {
496  *validNorm = true;
497  *n = G4ThreeVector(px * fDDy, py * fDDx, 0.).unit();
498  }
499  return 0.;
500  }
501 
502  // Just in case check if point is outside, normally it should never be
503  //
504  if (std::max(distZ, distR) > halfTolerance)
505  {
506 #ifdef G4SPECDEBUG
507  std::ostringstream message;
508  G4int oldprc = message.precision(16);
509  message << "Point p is outside (!?) of solid: "
510  << GetName() << G4endl;
511  message << "Position: " << p << G4endl;;
512  message << "Direction: " << v;
513  G4cout.precision(oldprc);
514  G4Exception("G4EllipticalTube::DistanceToOut(p,v)", "GeomSolids1002",
515  JustWarning, message );
516  DumpInfo();
517 #endif
518  if (calcNorm)
519  {
520  *validNorm = true;
521  *n = ApproxSurfaceNormal(p);
522  }
523  return 0.;
524  }
525 
526  // Set coefficients of quadratic equation: A t^2 + 2B t + C = 0
527  //
528  G4double A = vx * vx + vy * vy;
529  G4double C = rr - fR * fR;
530  G4double D = B * B - A * C;
531 
532  // Solve qudratic equation. There are two special cases where D <= 0:
533  // 1) trajectory parallel to Z axis (A = 0, B = 0, C - any, D = 0)
534  // 2) touch (D = 0) or no intersection (D < 0) with lateral surface
535  //
536  G4bool parallelToZ = (A < DBL_EPSILON || std::abs(vz) >= 1.);
537  if (parallelToZ) // 1)
538  {
539  if (calcNorm)
540  {
541  *validNorm = true;
542  n->set(0, 0, (vz < 0) ? -1. : 1.);
543  }
544  return tzmax;
545  }
546  if (D <= A * A * fScratch) // 2)
547  {
548  if (calcNorm)
549  {
550  *validNorm = true;
551  *n = G4ThreeVector(px * fDDy, py * fDDx, 0.).unit();
552  }
553  return 0.;
554  }
555 
556  // Find roots of quadratic equation
557  G4double tmp = -B - std::copysign(std::sqrt(D), B);
558  G4double t1 = tmp / A;
559  G4double t2 = C / tmp;
560  G4double trmax = std::max(t1, t2);
561 
562  // Return distance
563  G4double tmax = std::min(tzmax, trmax);
564 
565  // Set normal, if required, and return distance
566  //
567  if (calcNorm)
568  {
569  *validNorm = true;
570  G4ThreeVector pnew = p + tmax * v;
571  if (tmax == tzmax)
572  n->set(0, 0, (pnew.z() < 0) ? -1. : 1.);
573  else
574  *n = G4ThreeVector(pnew.x() * fDDy, pnew.y() * fDDx, 0.).unit();
575  }
576  return tmax;
577 }
578 
580 //
581 // Estimate distance to the surface from inside,
582 // returns 0 if point is outside
583 //
584 
586 {
587 #ifdef G4SPECDEBUG
588  if( Inside(p) == kOutside )
589  {
590  std::ostringstream message;
591  G4int oldprc = message.precision(16);
592  message << "Point p is outside (!?) of solid: " << GetName() << "\n"
593  << "Position:\n"
594  << " p.x() = " << p.x()/mm << " mm\n"
595  << " p.y() = " << p.y()/mm << " mm\n"
596  << " p.z() = " << p.z()/mm << " mm";
597  message.precision(oldprc) ;
598  G4Exception("G4ElliptocalTube::DistanceToOut(p)", "GeomSolids1002",
599  JustWarning, message);
600  DumpInfo();
601  }
602 #endif
603  // safety distance to Z-bases
604  G4double distZ = fDz - std::abs(p.z());
605 
606  // safety distance lateral surface
607  G4double x = p.x() * fSx;
608  G4double y = p.y() * fSy;
609  G4double distR = fR - std::sqrt(x * x + y * y);
610 
611  // return SafetyToOut
612  G4double dist = std::min(distZ, distR);
613  return (dist < 0) ? 0 : dist;
614 }
615 
617 //
618 // GetEntityType
619 
621 {
622  return G4String("G4EllipticalTube");
623 }
624 
626 //
627 // Make a clone of the object
628 
630 {
631  return new G4EllipticalTube(*this);
632 }
633 
635 //
636 // Return volume
637 
639 {
640  if (fCubicVolume == 0.)
641  {
642  fCubicVolume = twopi * fDx * fDy * fDz;
643  }
644  return fCubicVolume;
645 }
646 
648 //
649 // Return cached surface area
650 
652 {
653  G4ThreadLocalStatic G4double cached_Dx = 0;
654  G4ThreadLocalStatic G4double cached_Dy = 0;
655  G4ThreadLocalStatic G4double cached_Dz = 0;
656  G4ThreadLocalStatic G4double cached_area = 0;
657  if (cached_Dx != fDx || cached_Dy != fDy || cached_Dz != fDz)
658  {
659  cached_Dx = fDx;
660  cached_Dy = fDy;
661  cached_Dz = fDz;
662  cached_area = 2.*(pi*fDx*fDy + G4GeomTools::EllipsePerimeter(fDx, fDy)*fDz);
663  }
664  return cached_area;
665 }
666 
668 //
669 // Return surface area
670 
672 {
673  if(fSurfaceArea == 0.)
674  {
676  }
677  return fSurfaceArea;
678 }
679 
681 //
682 // Stream object contents to output stream
683 
684 std::ostream& G4EllipticalTube::StreamInfo(std::ostream& os) const
685 {
686  G4int oldprc = os.precision(16);
687  os << "-----------------------------------------------------------\n"
688  << " *** Dump for solid - " << GetName() << " ***\n"
689  << " ===================================================\n"
690  << " Solid type: G4EllipticalTube\n"
691  << " Parameters: \n"
692  << " length Z: " << fDz/mm << " mm \n"
693  << " lateral surface equation: \n"
694  << " (X / " << fDx << ")^2 + (Y / " << fDy << ")^2 = 1 \n"
695  << "-----------------------------------------------------------\n";
696  os.precision(oldprc);
697 
698  return os;
699 }
700 
701 
703 //
704 // Pick up a random point on the surface
705 
707 {
708  // Select surface (0 - base at -Z, 1 - base at +Z, 2 - lateral surface)
709  //
710  G4double sbase = pi * fDx * fDy;
711  G4double ssurf = GetCachedSurfaceArea();
712  G4double select = ssurf * G4UniformRand();
713 
714  G4int k = 0;
715  if (select > sbase) k = 1;
716  if (select > 2. * sbase) k = 2;
717 
718  // Pick random point on selected surface (rejection sampling)
719  //
721  switch (k) {
722  case 0: // base at -Z
723  {
725  p.set(rho.x(), rho.y(), -fDz);
726  break;
727  }
728  case 1: // base at +Z
729  {
731  p.set(rho.x(), rho.y(), fDz);
732  break;
733  }
734  case 2: // lateral surface
735  {
737  p.set(rho.x(), rho.y(), (2. * G4UniformRand() - 1.) * fDz);
738  break;
739  }
740  }
741  return p;
742 }
743 
744 
746 //
747 // CreatePolyhedron
748 
750 {
751  // create cylinder with radius=1...
752  //
753  G4Polyhedron* eTube = new G4PolyhedronTube(0., 1., fDz);
754 
755  // apply non-uniform scaling...
756  //
757  eTube->Transform(G4Scale3D(fDx, fDy, 1.));
758  return eTube;
759 }
760 
762 //
763 // GetPolyhedron
764 
766 {
767  if (fpPolyhedron == nullptr ||
771  {
772  G4AutoLock l(&polyhedronMutex);
773  delete fpPolyhedron;
775  fRebuildPolyhedron = false;
776  l.unlock();
777  }
778  return fpPolyhedron;
779 }
780 
782 //
783 // DescribeYourselfTo
784 
786 {
787  scene.AddSolid (*this);
788 }
789 
791 //
792 // GetExtent
793 
795 {
796  return G4VisExtent( -fDx, fDx, -fDy, fDy, -fDz, fDz );
797 }
798 
799 #endif // !defined(G4GEOM_USE_UELLIPTICALTUBE) || !defined(G4GEOM_USE_SYS_USOLIDS)