Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4USphere.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 // Implementation for G4USphere wrapper class
31 // --------------------------------------------------------------------
32 
33 #include "G4Sphere.hh"
34 #include "G4USphere.hh"
35 
36 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G4GEOM_USE_PARTIAL_USOLIDS) )
37 
38 #include "G4GeomTools.hh"
39 #include "G4AffineTransform.hh"
40 #include "G4VPVParameterisation.hh"
41 #include "G4BoundingEnvelope.hh"
42 
43 using namespace CLHEP;
44 
46 //
47 // constructor - check parameters, convert angles so 0<sphi+dpshi<=2_PI
48 // - note if pDPhi>2PI then reset to 2PI
49 
50 G4USphere::G4USphere( const G4String& pName,
51  G4double pRmin, G4double pRmax,
52  G4double pSPhi, G4double pDPhi,
53  G4double pSTheta, G4double pDTheta )
54  : G4USolid(pName, new USphere(pName, pRmin, pRmax, pSPhi, pDPhi,
55  pSTheta, pDTheta))
56 {
57 }
58 
60 //
61 // Fake default constructor - sets only member data and allocates memory
62 // for usage restricted to object persistency.
63 //
64 G4USphere::G4USphere( __void__& a )
65  : G4USolid(a)
66 {
67 }
68 
70 //
71 // Destructor
72 
73 G4USphere::~G4USphere()
74 {
75 }
76 
78 //
79 // Copy constructor
80 
81 G4USphere::G4USphere(const G4USphere& rhs)
82  : G4USolid(rhs)
83 {
84 }
85 
87 //
88 // Assignment operator
89 
90 G4USphere& G4USphere::operator = (const G4USphere& rhs)
91 {
92  // Check assignment to self
93  //
94  if (this == &rhs) { return *this; }
95 
96  // Copy base class data
97  //
98  G4USolid::operator=(rhs);
99 
100  return *this;
101 }
102 
104 //
105 // Accessors & modifiers
106 
107 G4double G4USphere::GetInnerRadius() const
108 {
109  return GetShape()->GetInnerRadius();
110 }
111 G4double G4USphere::GetOuterRadius() const
112 {
113  return GetShape()->GetOuterRadius();
114 }
115 G4double G4USphere::GetStartPhiAngle() const
116 {
117  return GetShape()->GetStartPhiAngle();
118 }
119 G4double G4USphere::GetDeltaPhiAngle() const
120 {
121  return GetShape()->GetDeltaPhiAngle();
122 }
123 G4double G4USphere::GetStartThetaAngle() const
124 {
125  return GetShape()->GetStartThetaAngle();
126 }
127 G4double G4USphere::GetDeltaThetaAngle() const
128 {
129  return GetShape()->GetDeltaThetaAngle();
130 }
131 G4double G4USphere::GetSinStartPhi() const
132 {
133  G4double phi = GetShape()->GetStartPhiAngle();
134  return std::sin(phi);
135 }
136 G4double G4USphere::GetCosStartPhi() const
137 {
138  G4double phi = GetShape()->GetStartPhiAngle();
139  return std::cos(phi);
140 }
141 G4double G4USphere::GetSinEndPhi() const
142 {
143  G4double phi = GetShape()->GetStartPhiAngle() +
144  GetShape()->GetDeltaPhiAngle();
145  return std::sin(phi);
146 }
147 G4double G4USphere::GetCosEndPhi() const
148 {
149  G4double phi = GetShape()->GetStartPhiAngle() +
150  GetShape()->GetDeltaPhiAngle();
151  return std::cos(phi);
152 }
153 G4double G4USphere::GetSinStartTheta() const
154 {
155  G4double theta = GetShape()->GetStartThetaAngle();
156  return std::sin(theta);
157 }
158 G4double G4USphere::GetCosStartTheta() const
159 {
160  G4double theta = GetShape()->GetStartThetaAngle();
161  return std::cos(theta);
162 }
163 G4double G4USphere::GetSinEndTheta() const
164 {
165  G4double theta = GetShape()->GetStartThetaAngle() +
166  GetShape()->GetDeltaThetaAngle();
167  return std::sin(theta);
168 }
169 G4double G4USphere::GetCosEndTheta() const
170 {
171  G4double theta = GetShape()->GetStartThetaAngle() +
172  GetShape()->GetDeltaThetaAngle();
173  return std::cos(theta);
174 }
175 
176 void G4USphere::SetInnerRadius(G4double newRMin)
177 {
178  GetShape()->SetInnerRadius(newRMin);
179  fRebuildPolyhedron = true;
180 }
181 void G4USphere::SetOuterRadius(G4double newRmax)
182 {
183  GetShape()->SetOuterRadius(newRmax);
184  fRebuildPolyhedron = true;
185 }
186 void G4USphere::SetStartPhiAngle(G4double newSphi, G4bool trig)
187 {
188  GetShape()->SetStartPhiAngle(newSphi, trig);
189  fRebuildPolyhedron = true;
190 }
191 void G4USphere::SetDeltaPhiAngle(G4double newDphi)
192 {
193  GetShape()->SetDeltaPhiAngle(newDphi);
194  fRebuildPolyhedron = true;
195 }
196 void G4USphere::SetStartThetaAngle(G4double newSTheta)
197 {
198  GetShape()->SetStartThetaAngle(newSTheta);
199  fRebuildPolyhedron = true;
200 }
201 void G4USphere::SetDeltaThetaAngle(G4double newDTheta)
202 {
203  GetShape()->SetDeltaThetaAngle(newDTheta);
204  fRebuildPolyhedron = true;
205 }
206 
208 //
209 // Dispatch to parameterisation for replication mechanism dimension
210 // computation & modification.
211 
212 void G4USphere::ComputeDimensions( G4VPVParameterisation* p,
213  const G4int n,
214  const G4VPhysicalVolume* pRep)
215 {
216  p->ComputeDimensions(*(G4Sphere*)this,n,pRep);
217 }
218 
220 //
221 // Make a clone of the object
222 
223 G4VSolid* G4USphere::Clone() const
224 {
225  return new G4USphere(*this);
226 }
227 
229 //
230 // Get bounding box
231 
232 void G4USphere::Extent(G4ThreeVector& pMin, G4ThreeVector& pMax) const
233 {
234  static G4bool checkBBox = true;
235 
236  G4double rmin = GetInnerRadius();
237  G4double rmax = GetOuterRadius();
238 
239  // Find bounding box
240  //
241  if (GetDeltaThetaAngle() >= pi && GetDeltaPhiAngle() >= twopi)
242  {
243  pMin.set(-rmax,-rmax,-rmax);
244  pMax.set( rmax, rmax, rmax);
245  }
246  else
247  {
248  G4double sinStart = GetSinStartTheta();
249  G4double cosStart = GetCosStartTheta();
250  G4double sinEnd = GetSinEndTheta();
251  G4double cosEnd = GetCosEndTheta();
252 
253  G4double stheta = GetStartThetaAngle();
254  G4double etheta = stheta + GetDeltaThetaAngle();
255  G4double rhomin = rmin*std::min(sinStart,sinEnd);
256  G4double rhomax = rmax;
257  if (stheta > halfpi) rhomax = rmax*sinStart;
258  if (etheta < halfpi) rhomax = rmax*sinEnd;
259 
260  G4TwoVector xymin,xymax;
261  G4GeomTools::DiskExtent(rhomin,rhomax,
262  GetSinStartPhi(),GetCosStartPhi(),
263  GetSinEndPhi(),GetCosEndPhi(),
264  xymin,xymax);
265 
266  G4double zmin = std::min(rmin*cosEnd,rmax*cosEnd);
267  G4double zmax = std::max(rmin*cosStart,rmax*cosStart);
268  pMin.set(xymin.x(),xymin.y(),zmin);
269  pMax.set(xymax.x(),xymax.y(),zmax);
270  }
271 
272  // Check correctness of the bounding box
273  //
274  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
275  {
276  std::ostringstream message;
277  message << "Bad bounding box (min >= max) for solid: "
278  << GetName() << " !"
279  << "\npMin = " << pMin
280  << "\npMax = " << pMax;
281  G4Exception("G4USphere::Extent()", "GeomMgt0001", JustWarning, message);
282  StreamInfo(G4cout);
283  }
284 
285  // Check consistency of bounding boxes
286  //
287  if (checkBBox)
288  {
289  UVector3 vmin, vmax;
290  GetShape()->Extent(vmin,vmax);
291  if (std::abs(pMin.x()-vmin.x()) > kCarTolerance ||
292  std::abs(pMin.y()-vmin.y()) > kCarTolerance ||
293  std::abs(pMin.z()-vmin.z()) > kCarTolerance ||
294  std::abs(pMax.x()-vmax.x()) > kCarTolerance ||
295  std::abs(pMax.y()-vmax.y()) > kCarTolerance ||
296  std::abs(pMax.z()-vmax.z()) > kCarTolerance)
297  {
298  std::ostringstream message;
299  message << "Inconsistency in bounding boxes for solid: "
300  << GetName() << " !"
301  << "\nBBox min: wrapper = " << pMin << " solid = " << vmin
302  << "\nBBox max: wrapper = " << pMax << " solid = " << vmax;
303  G4Exception("G4USphere::Extent()", "GeomMgt0001", JustWarning, message);
304  checkBBox = false;
305  }
306  }
307 }
308 
310 //
311 // Calculate extent under transform and specified limit
312 
313 G4bool G4USphere::CalculateExtent(const EAxis pAxis,
314  const G4VoxelLimits& pVoxelLimit,
315  const G4AffineTransform& pTransform,
316  G4double& pMin, G4double& pMax) const
317 {
318  G4ThreeVector bmin, bmax;
319 
320  // Get bounding box
321  Extent(bmin,bmax);
322 
323  // Find extent
324  G4BoundingEnvelope bbox(bmin,bmax);
325  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
326 }
327 
329 //
330 // Create polyhedron for visualization
331 
332 G4Polyhedron* G4USphere::CreatePolyhedron() const
333 {
334  return new G4PolyhedronSphere(GetInnerRadius(),
335  GetOuterRadius(),
336  GetStartPhiAngle(),
337  GetDeltaPhiAngle(),
338  GetStartThetaAngle(),
339  GetDeltaThetaAngle());
340 }
341 
342 #endif // G4GEOM_USE_USOLIDS
void set(double x, double y, double z)
double y() const
double x() const
double x() const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
const char * p
Definition: xmltok.h:285
int G4int
Definition: G4Types.hh:78
double z() const
static constexpr double twopi
Definition: G4SIunits.hh:76
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
const G4double kCarTolerance
const G4int n
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
EAxis
Definition: geomdefs.hh:54
double y() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static constexpr double pi
Definition: G4SIunits.hh:75
static constexpr double halfpi
Definition: G4SIunits.hh:77
double G4double
Definition: G4Types.hh:76
static G4bool DiskExtent(G4double rmin, G4double rmax, G4double startPhi, G4double delPhi, G4TwoVector &pmin, G4TwoVector &pmax)
Definition: G4GeomTools.cc:378