Geant4  10.00.p03
G4USolid.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 // GEANT4 tag $Name:$
29 //
30 //
31 // G4USolid implementation
32 //
33 // --------------------------------------------------------------------
34 
35 #include "G4USolid.hh"
36 #include "G4AffineTransform.hh"
37 #include "G4VoxelLimits.hh"
38 #include "G4VGraphicsScene.hh"
39 #include "G4Polyhedron.hh"
40 #include "G4PolyhedronArbitrary.hh"
41 #include "G4VisExtent.hh"
42 #include "G4PhysicalConstants.hh"
43 
44 #include "G4AutoLock.hh"
45 
46 namespace
47 {
48  G4Mutex polyhedronMutex = G4MUTEX_INITIALIZER;
49 }
50 
52  G4VSolid(name), fShape(s), fRebuildPolyhedron(false), fPolyhedron(0)
53 {
54 }
55 
57  : G4VSolid(a), fShape(0), fRebuildPolyhedron(false), fPolyhedron(0)
58 {
59 }
60 
62 {
63  delete fPolyhedron; fPolyhedron = 0;
64 }
65 
67 {
68  return (this == &s) ? true : false;
69 }
70 
72 {
73  UVector3 pt;
74  VUSolid::EnumInside in_temp;
75  EInside in = kOutside;
76  pt.x = p.x();
77  pt.y = p.y();
78  pt.z = p.z(); // better assign at construction
79 
80  in_temp = fShape->Inside(pt);
81 
82  if (in_temp == VUSolid::eSurface)return kSurface;
83  if (in_temp == VUSolid::eInside)return kInside;
84 
85  return in;
86 }
87 
89 {
90  UVector3 p;
91  p.x = pt.x();
92  p.y = pt.y();
93  p.z = pt.z();
94  UVector3 n;
95  fShape->Normal(p, n);
96  return G4ThreeVector(n.x, n.y, n.z);
97 }
98 
100  const G4ThreeVector& d)const
101 {
102  UVector3 p;
103  p.x = pt.x();
104  p.y = pt.y();
105  p.z = pt.z(); // better assign at construction
106  UVector3 v;
107  v.x = d.x();
108  v.y = d.y();
109  v.z = d.z(); // better assign at construction
110  G4double dist = fShape->DistanceToIn(p, v);
111  if (dist > kInfinity) dist = kInfinity;
112  return dist;
113 }
114 
116 {
117  UVector3 p;
118  p.x = pt.x();
119  p.y = pt.y();
120  p.z = pt.z(); // better assign at construction
121  G4double dist = fShape->SafetyFromOutside(p); // true?
122  if (dist > kInfinity) dist = kInfinity;
123  return dist;
124 }
125 
127  const G4ThreeVector& d,
128  const G4bool calcNorm,
129  G4bool* validNorm,
130  G4ThreeVector* norm) const
131 {
132  UVector3 p;
133  p.x = pt.x();
134  p.y = pt.y();
135  p.z = pt.z(); // better assign at construction
136  UVector3 v;
137  v.x = d.x();
138  v.y = d.y();
139  v.z = d.z(); // better assign at construction
140  UVector3 n;
141  bool valid;
142  G4double dist = fShape->DistanceToOut(p, v, n,valid); // should use local variable
143  if(calcNorm)
144  {
145  if(valid){ *validNorm = true;}
146  else {* validNorm =false;}
147  if(*validNorm)
148  { norm->setX(n.x);
149  norm->setY(n.y);
150  norm->setZ(n.z);
151  } // *norm = n, but only after calcNorm check
152  }
153  if (dist > kInfinity) dist = kInfinity;
154  return dist;
155 }
156 
158 {
159  UVector3 p;
160  p.x = pt.x();
161  p.y = pt.y();
162  p.z = pt.z(); // better assign at construction
163  return fShape->SafetyFromInside(p); // true?
164 }
165 
167 {
168  return fShape->Capacity();
169 }
170 
172 {
173  return fShape->SurfaceArea();
174 }
175 
177 {
178  UVector3 p;
179  p = fShape->GetPointOnSurface();
180  return G4ThreeVector(p.x, p.y, p.z);
181 }
182 
184  const G4VoxelLimits& pVoxelLimit,
185  const G4AffineTransform& pTransform,
186  G4double& pMin, G4double& pMax) const
187 {
188  if (!pTransform.IsRotated())
189  {
191  G4double offset = pTransform.NetTranslation().x();
192  if (pAxis == kYAxis)
193  {
194  eAxis = VUSolid::eYaxis;
195  offset = pTransform.NetTranslation().y();
196  }
197  if (pAxis == kZAxis)
198  {
199  eAxis = VUSolid::eZaxis;
200  offset = pTransform.NetTranslation().z();
201  }
202  fShape->ExtentAxis(eAxis, pMin, pMax);
203 
204  pMin += offset;
205  pMax += offset;
206 
207  if (pVoxelLimit.IsLimited())
208  {
209  switch (pAxis)
210  {
211  case kXAxis:
212  if ((pMin > pVoxelLimit.GetMaxXExtent() + kCarTolerance) ||
213  (pMax < pVoxelLimit.GetMinXExtent() - kCarTolerance))
214  {
215  return false;
216  }
217  else
218  {
219  pMin = std::max(pMin, pVoxelLimit.GetMinXExtent());
220  pMax = std::min(pMax, pVoxelLimit.GetMaxXExtent());
221  }
222  break;
223  case kYAxis:
224  if ((pMin > pVoxelLimit.GetMaxYExtent() + kCarTolerance) ||
225  (pMax < pVoxelLimit.GetMinYExtent() - kCarTolerance))
226  {
227  return false;
228  }
229  else
230  {
231  pMin = std::max(pMin, pVoxelLimit.GetMinYExtent());
232  pMax = std::min(pMax, pVoxelLimit.GetMaxYExtent());
233  }
234  break;
235  case kZAxis:
236  if ((pMin > pVoxelLimit.GetMaxZExtent() + kCarTolerance) ||
237  (pMax < pVoxelLimit.GetMinZExtent() - kCarTolerance))
238  {
239  return false;
240  }
241  else
242  {
243  pMin = std::max(pMin, pVoxelLimit.GetMinZExtent());
244  pMax = std::min(pMax, pVoxelLimit.GetMaxZExtent());
245  }
246  break;
247  default:
248  break;
249  }
250  pMin -= kCarTolerance ;
251  pMax += kCarTolerance ;
252  }
253  return true;
254  }
255  else // General rotated case - create and clip mesh to boundaries
256  {
257  // Rotate BoundingBox and Calculate Extent as for BREPS
258 
259  G4bool existsAfterClip = false ;
260  G4ThreeVectorList* vertices ;
261 
262  pMin = +kInfinity ;
263  pMax = -kInfinity ;
264 
265  // Calculate rotated vertex coordinates
266 
267  vertices = CreateRotatedVertices(pTransform) ;
268  ClipCrossSection(vertices, 0, pVoxelLimit, pAxis, pMin, pMax) ;
269  ClipCrossSection(vertices, 4, pVoxelLimit, pAxis, pMin, pMax) ;
270  ClipBetweenSections(vertices, 0, pVoxelLimit, pAxis, pMin, pMax) ;
271 
272  if (pVoxelLimit.IsLimited(pAxis) == false)
273  {
274  if ((pMin != kInfinity) || (pMax != -kInfinity))
275  {
276  existsAfterClip = true ;
277 
278  // Add 2*tolerance to avoid precision troubles
279 
280  pMin -= kCarTolerance;
281  pMax += kCarTolerance;
282  }
283  }
284  else
285  {
286  G4ThreeVector clipCentre(
287  (pVoxelLimit.GetMinXExtent() + pVoxelLimit.GetMaxXExtent()) * 0.5,
288  (pVoxelLimit.GetMinYExtent() + pVoxelLimit.GetMaxYExtent()) * 0.5,
289  (pVoxelLimit.GetMinZExtent() + pVoxelLimit.GetMaxZExtent()) * 0.5);
290 
291  if ((pMin != kInfinity) || (pMax != -kInfinity))
292  {
293  existsAfterClip = true ;
294 
295 
296  // Check to see if endpoints are in the solid
297 
298  clipCentre(pAxis) = pVoxelLimit.GetMinExtent(pAxis);
299 
300  if (Inside(pTransform.Inverse().TransformPoint(clipCentre)) != kOutside)
301  {
302  pMin = pVoxelLimit.GetMinExtent(pAxis);
303  }
304  else
305  {
306  pMin -= kCarTolerance;
307  }
308  clipCentre(pAxis) = pVoxelLimit.GetMaxExtent(pAxis);
309 
310  if (Inside(pTransform.Inverse().TransformPoint(clipCentre)) != kOutside)
311  {
312  pMax = pVoxelLimit.GetMaxExtent(pAxis);
313  }
314  else
315  {
316  pMax += kCarTolerance;
317  }
318  }
319 
320  // Check for case where completely enveloping clipping volume
321  // If point inside then we are confident that the solid completely
322  // envelopes the clipping volume. Hence set min/max extents according
323  // to clipping volume extents along the specified axis.
324 
325  else if (Inside(pTransform.Inverse().TransformPoint(clipCentre))
326  != kOutside)
327  {
328  existsAfterClip = true ;
329  pMin = pVoxelLimit.GetMinExtent(pAxis) ;
330  pMax = pVoxelLimit.GetMaxExtent(pAxis) ;
331  }
332  }
333  delete vertices;
334  return existsAfterClip;
335  }
336 }
337 
339  const G4int,
340  const G4VPhysicalVolume*)
341 {
342  std::ostringstream message;
343  message << "Illegal call to G4USolid::ComputeDimensions()" << G4endl
344  << "Method not overloaded by derived class !";
345  G4Exception("G4USolid::ComputeDimensions()", "GeomSolids0003",
346  FatalException, message);
347 }
348 
350 {
351  scene.AddSolid(*this);
352 }
353 
355 {
356 
357  G4String string = fShape->GetEntityType();
358  return "G4" + string;
359 }
360 
361 std::ostream& G4USolid::StreamInfo(std::ostream& os) const
362 {
363  os << "-----------------------------------------------------------\n"
364  << " *** Dump for solid - " << fShape->GetName() << " ***\n";
365  // << " ===================================================\n"
366  // << " Solid type: " << fShape->GetEntityType() << "\n"
367  // << "-----------------------------------------------------------\n";
368 
369  return os;//fShape->StreamInfo(os);
370 }
371 
373  : G4VSolid(rhs), fRebuildPolyhedron(false), fPolyhedron(0)
374 {
375  fShape = rhs.fShape->Clone();
376 }
377 
379 {
380  // Check assignment to self
381  //
382  if (this == &rhs)
383  {
384  return *this;
385  }
386 
387  // Copy base class data
388  //
389  G4VSolid::operator=(rhs);
390 
391  // Copy data
392  //
393  fShape = rhs.fShape->Clone();
394  fRebuildPolyhedron = false;
395  delete fPolyhedron; fPolyhedron = 0;
396 
397  return *this;
398 }
399 
401 {
402  std::ostringstream message;
403  message << "Clone() method not implemented for type: "
404  << GetEntityType() << "!" << G4endl
405  << "Returning NULL pointer!";
406  G4Exception("G4USolid::Clone()", "GeomSolids1001", JustWarning, message);
407  return 0;
408 }
409 
412 {
413  G4double xMin, xMax, yMin, yMax, zMin, zMax;
414 
415  fShape->ExtentAxis(VUSolid::eXaxis, xMin, xMax);
416  fShape->ExtentAxis(VUSolid::eYaxis, yMin, yMax);
417  fShape->ExtentAxis(VUSolid::eZaxis, zMin, zMax);
418 
419  G4ThreeVectorList* vertices;
420  vertices = new G4ThreeVectorList();
421 
422  if (vertices)
423  {
424  vertices->reserve(8);
425  G4ThreeVector vertex0(xMin, yMin, zMin);
426  G4ThreeVector vertex1(xMax, yMin, zMin);
427  G4ThreeVector vertex2(xMax, yMax, zMin);
428  G4ThreeVector vertex3(xMin, yMax, zMin);
429  G4ThreeVector vertex4(xMin, yMin, zMax);
430  G4ThreeVector vertex5(xMax, yMin, zMax);
431  G4ThreeVector vertex6(xMax, yMax, zMax);
432  G4ThreeVector vertex7(xMin, yMax, zMax);
433 
434  vertices->push_back(pTransform.TransformPoint(vertex0));
435  vertices->push_back(pTransform.TransformPoint(vertex1));
436  vertices->push_back(pTransform.TransformPoint(vertex2));
437  vertices->push_back(pTransform.TransformPoint(vertex3));
438  vertices->push_back(pTransform.TransformPoint(vertex4));
439  vertices->push_back(pTransform.TransformPoint(vertex5));
440  vertices->push_back(pTransform.TransformPoint(vertex6));
441  vertices->push_back(pTransform.TransformPoint(vertex7));
442  }
443  else
444  {
445  G4Exception("G4VUSolid::CreateRotatedVertices()", "FatalError",
446  FatalException, "Out of memory - Cannot allocate vertices!");
447  }
448  return vertices;
449 }
450 
452 {
453  G4int index = 0;
454  if (fShape->GetEntityType() == "Box")
455  {
456  double array[3];
457  fShape->GetParametersList(index, array);
458  return new G4PolyhedronBox(array[0], array[1], array[2]);
459  }
460  if (fShape->GetEntityType() == "Tubs")
461  {
462  double array[5];
463  fShape->GetParametersList(index, array);
464  return new G4PolyhedronTubs(array[0], array[1], array[2], array[3], array[4]);
465  }
466  if (fShape->GetEntityType() == "Cons")
467  {
468  double array[7];
469  fShape->GetParametersList(index, array);
470  return new G4PolyhedronCons(array[0], array[1], array[2], array[3], array[4], array[5], array[6]);
471  }
472  if (fShape->GetEntityType() == "Orb")
473  {
474  double array[1];
475  fShape->GetParametersList(index, array);
476  return new G4PolyhedronSphere(0., array[0], 0., 2 * pi, 0., pi);
477  }
478  if (fShape->GetEntityType() == "Sphere")
479  {
480  double array[6];
481  fShape->GetParametersList(index, array);
482  return new G4PolyhedronSphere(array[0], array[1], array[2], array[3], array[4], array[5]);
483  }
484  if (fShape->GetEntityType() == "Tet")
485  {
486  double array[12];
487  fShape->GetParametersList(index, array);
488  G4Polyhedron* ph = new G4Polyhedron;
489  double xyz[4][3];
490  static int faces[4][4] = {{1, 3, 2, 0}, {1, 4, 3, 0}, {1, 2, 4, 0}, {2, 3, 4, 0}};
491  xyz[0][0] = array[0];
492  xyz[0][1] = array[1];
493  xyz[0][2] = array[2];
494  xyz[1][0] = array[3];
495  xyz[1][1] = array[4];
496  xyz[1][2] = array[5];
497  xyz[2][0] = array[6];
498  xyz[2][1] = array[7];
499  xyz[2][2] = array[8];
500  xyz[3][0] = array[9];
501  xyz[3][1] = array[10];
502  xyz[3][2] = array[11];
503 
504  ph->createPolyhedron(4, 4, xyz, faces);
505  return ph;
506  }
507  if (fShape->GetEntityType() == "Trd")
508  {
509  double array[5];
510  fShape->GetParametersList(index, array);
511  return new G4PolyhedronTrd2(array[0], array[1], array[2], array[3], array[4]);
512  }
513  if (fShape->GetEntityType() == "Trap")
514  {
515  double array[12];
516  fShape->GetParametersList(index, array);
517  double phi = (array[11] != 1.0) ? (std::atan(array[10] / array[9])) : (0.0);
518  double alpha1 = std::atan(array[4]);
519  double alpha2 = std::atan(array[8]);
520  double theta = std::acos(array[11]);
521 
522  return new G4PolyhedronTrap(array[0], theta, phi,
523  array[1], array[2], array[3], alpha1,
524  array[5], array[6], array[7], alpha2);
525  }
526 
527  /*
528  if(fShape->GetEntityType()=="TessellatedSolid"){
529 
530  G4Polyhedron *uPolyhedron=fShape->GetPolyhedron();
531  std::size_t nVertices = (*uPolyhedron).vertices.size();
532  std::size_t nFacets = (*uPolyhedron).facets.size();
533 
534  G4PolyhedronArbitrary *polyhedron =
535  new G4PolyhedronArbitrary (nVertices, nFacets);
536 
537  for (std::vector<UVector3>::const_iterator v = (*uPolyhedron).vertices.begin();
538  v!=(*uPolyhedron).vertices.end(); v++)
539  {
540  UVector3 p=(*v);
541  G4ThreeVector pt(p.x,p.y,p.z);
542 
543  polyhedron->AddVertex(pt);
544  }
545  for (std::vector<UFacet>::const_iterator f=(*uPolyhedron).facets.begin();
546  f != (*uPolyhedron).facets.end(); f++)
547  {
548  polyhedron->AddFacet((*f).f1,(*f).f2,(*f).f3,(*f).f4);
549  }
550 
551  return (G4Polyhedron*) polyhedron;
552  }
553  */
554 
555  return 0;
556 }
557 
559 {
560  if (!fPolyhedron ||
563  fPolyhedron->GetNumberOfRotationSteps())
564  {
565  G4AutoLock l(&polyhedronMutex);
566  delete fPolyhedron;
568  fRebuildPolyhedron = false;
569  l.unlock();
570  }
571  return fPolyhedron;
572 }
573 
575 {
576  G4VisExtent extent;
577  G4VoxelLimits voxelLimits; // Defaults to "infinite" limits.
578  G4AffineTransform affineTransform;
579  G4double vmin, vmax;
580  CalculateExtent(kXAxis, voxelLimits, affineTransform, vmin, vmax);
581  extent.SetXmin(vmin);
582  extent.SetXmax(vmax);
583  CalculateExtent(kYAxis, voxelLimits, affineTransform, vmin, vmax);
584  extent.SetYmin(vmin);
585  extent.SetYmax(vmax);
586  CalculateExtent(kZAxis, voxelLimits, affineTransform, vmin, vmax);
587  extent.SetZmin(vmin);
588  extent.SetZmax(vmax);
589  return extent;
590 }
virtual UVector3 GetPointOnSurface() const =0
void ClipCrossSection(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
Definition: G4VSolid.cc:347
G4ThreeVectorList * CreateRotatedVertices(const G4AffineTransform &pT) const
Definition: G4USolid.cc:411
virtual UGeometryType GetEntityType() const =0
virtual double SurfaceArea()=0
Definition: VUSolid.cc:255
virtual void DescribeYourselfTo(G4VGraphicsScene &scene) const
Definition: G4USolid.cc:349
static const G4double kInfinity
Definition: geomdefs.hh:42
G4double GetMinYExtent() const
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
Definition: G4USolid.cc:126
CLHEP::Hep3Vector G4ThreeVector
G4AffineTransform Inverse() const
G4Polyhedron * CreatePolyhedron() const
Definition: G4USolid.cc:451
const std::string & GetName() const
Definition: VUSolid.hh:103
G4String name
Definition: TRTMaterials.hh:40
const G4double pi
G4bool IsRotated() const
virtual std::ostream & StreamInfo(std::ostream &os) const
Definition: G4USolid.cc:361
void SetXmax(G4double xmax)
Definition: G4VisExtent.hh:102
G4ThreeVector NetTranslation() const
virtual G4double GetCubicVolume()
Definition: G4USolid.cc:166
virtual G4double GetSurfaceArea()
Definition: G4USolid.cc:171
G4double a
Definition: TRTMaterials.hh:39
virtual bool Normal(const UVector3 &aPoint, UVector3 &aNormal) const =0
virtual EnumInside Inside(const UVector3 &aPoint) const =0
virtual EInside Inside(const G4ThreeVector &p) const
Definition: G4USolid.cc:71
virtual void AddSolid(const G4Box &)=0
int G4int
Definition: G4Types.hh:78
double x
Definition: UVector3.hh:136
#define G4MUTEX_INITIALIZER
Definition: G4Threading.hh:158
G4bool operator==(const G4USolid &s) const
Definition: G4USolid.cc:66
void SetYmax(G4double ymax)
Definition: G4VisExtent.hh:106
static const double s
Definition: G4SIunits.hh:150
G4double GetMaxXExtent() const
virtual G4Polyhedron * GetPolyhedron() const
Definition: G4USolid.cc:558
G4double GetMinZExtent() const
virtual double DistanceToIn(const UVector3 &aPoint, const UVector3 &aDirection, double aPstep=UUtils::kInfinity) const =0
G4bool IsLimited() const
VUSolid * fShape
Definition: G4USolid.hh:182
bool G4bool
Definition: G4Types.hh:79
std::vector< G4ThreeVector > G4ThreeVectorList
Definition: G4VSolid.hh:79
G4bool fRebuildPolyhedron
Definition: G4USolid.hh:183
EAxisType
Definition: VUSolid.hh:27
EnumInside
Definition: VUSolid.hh:23
const G4int n
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4USolid.cc:99
virtual double SafetyFromInside(const UVector3 &aPoint, bool aAccurate=false) const =0
G4double GetMinXExtent() const
G4int G4Mutex
Definition: G4Threading.hh:156
G4ThreeVector TransformPoint(const G4ThreeVector &vec) const
G4double GetMaxZExtent() const
void SetZmin(G4double zmin)
Definition: G4VisExtent.hh:108
G4USolid(const G4String &pName, VUSolid *shape)
Definition: G4USolid.cc:51
G4Polyhedron * fPolyhedron
Definition: G4USolid.hh:184
T max(const T t1, const T t2)
brief Return the largest of the two arguments
EInside
Definition: geomdefs.hh:58
virtual void ExtentAxis(EAxisType aAxis, double &aMin, double &aMax) const
Definition: VUSolid.cc:318
EAxis
Definition: geomdefs.hh:54
virtual double SafetyFromOutside(const UVector3 &aPoint, bool aAccurate=false) const =0
virtual G4VisExtent GetExtent() const
Definition: G4USolid.cc:574
G4USolid & operator=(const G4USolid &rhs)
Definition: G4USolid.cc:378
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
Definition: G4USolid.cc:88
#define G4endl
Definition: G4ios.hh:61
G4double GetMaxYExtent() const
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:110
G4double kCarTolerance
Definition: G4VSolid.hh:305
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
virtual double DistanceToOut(const UVector3 &aPoint, const UVector3 &aDirection, UVector3 &aNormalVector, bool &aConvex, double aPstep=UUtils::kInfinity) const =0
void SetZmax(G4double zmax)
Definition: G4VisExtent.hh:110
double z
Definition: UVector3.hh:138
virtual VUSolid * Clone() const =0
virtual ~G4USolid()
Definition: G4USolid.cc:61
void SetXmin(G4double xmin)
Definition: G4VisExtent.hh:100
double G4double
Definition: G4Types.hh:76
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
Definition: G4USolid.cc:183
G4double GetMaxExtent(const EAxis pAxis) const
void SetYmin(G4double ymin)
Definition: G4VisExtent.hh:104
virtual G4VSolid * Clone() const
Definition: G4USolid.cc:400
double y
Definition: UVector3.hh:137
void ClipBetweenSections(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
Definition: G4VSolid.cc:378
G4double GetMinExtent(const EAxis pAxis) const
virtual double Capacity()=0
Definition: VUSolid.cc:199
virtual G4GeometryType GetEntityType() const
Definition: G4USolid.cc:354
virtual void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
Definition: G4USolid.cc:338
virtual G4ThreeVector GetPointOnSurface() const
Definition: G4USolid.cc:176
virtual void GetParametersList(int aNumber, double *aArray) const =0