| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // Implementation for G4UCons wrapper class <<
27 // 26 //
28 // 30.10.13 G.Cosmo, CERN/PH << 27 // $Id:$
>> 28 //
>> 29 //
>> 30 // Implementation for G4UCons wrapper class
29 // ------------------------------------------- 31 // --------------------------------------------------------------------
30 32
31 #include "G4Cons.hh" 33 #include "G4Cons.hh"
32 #include "G4UCons.hh" 34 #include "G4UCons.hh"
33 <<
34 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G <<
35 <<
36 #include "G4GeomTools.hh" <<
37 #include "G4AffineTransform.hh" <<
38 #include "G4VPVParameterisation.hh" 35 #include "G4VPVParameterisation.hh"
39 #include "G4BoundingEnvelope.hh" << 36
40 <<
41 using namespace CLHEP; <<
42 <<
43 ////////////////////////////////////////////// 37 //////////////////////////////////////////////////////////////////////////
44 // 38 //
45 // constructor - check parameters, convert ang 39 // constructor - check parameters, convert angles so 0<sphi+dpshi<=2_PI
46 // - note if pDPhi>2PI then rese 40 // - note if pDPhi>2PI then reset to 2PI
47 41
48 G4UCons::G4UCons( const G4String& pName, 42 G4UCons::G4UCons( const G4String& pName,
49 G4double pRmin1, G4do 43 G4double pRmin1, G4double pRmax1,
50 G4double pRmin2, G4do 44 G4double pRmin2, G4double pRmax2,
51 G4double pDz, 45 G4double pDz,
52 G4double pSPhi, G4doub 46 G4double pSPhi, G4double pDPhi)
53 : Base_t(pName, pRmin1, pRmax1, pRmin2, pRma << 47 : G4USolid(pName, new UCons(pName, pRmin1, pRmax1, pRmin2, pRmax2,
>> 48 pDz, pSPhi, pDPhi))
54 { 49 {
55 } 50 }
56 51
57 ////////////////////////////////////////////// 52 ///////////////////////////////////////////////////////////////////////
58 // 53 //
59 // Fake default constructor - sets only member 54 // Fake default constructor - sets only member data and allocates memory
60 // for usage restri 55 // for usage restricted to object persistency.
61 // 56 //
62 G4UCons::G4UCons( __void__& a ) 57 G4UCons::G4UCons( __void__& a )
63 : Base_t(a) << 58 : G4USolid(a)
64 { 59 {
65 } 60 }
66 61
67 ////////////////////////////////////////////// 62 ///////////////////////////////////////////////////////////////////////
68 // 63 //
69 // Destructor 64 // Destructor
70 65
71 G4UCons::~G4UCons() = default; << 66 G4UCons::~G4UCons()
>> 67 {
>> 68 }
72 69
73 ////////////////////////////////////////////// 70 //////////////////////////////////////////////////////////////////////////
74 // 71 //
75 // Copy constructor 72 // Copy constructor
76 73
77 G4UCons::G4UCons(const G4UCons& rhs) 74 G4UCons::G4UCons(const G4UCons& rhs)
78 : Base_t(rhs) << 75 : G4USolid(rhs)
79 { 76 {
80 } 77 }
81 78
82 ////////////////////////////////////////////// 79 //////////////////////////////////////////////////////////////////////////
83 // 80 //
84 // Assignment operator 81 // Assignment operator
85 82
86 G4UCons& G4UCons::operator = (const G4UCons& r 83 G4UCons& G4UCons::operator = (const G4UCons& rhs)
87 { 84 {
88 // Check assignment to self 85 // Check assignment to self
89 // 86 //
90 if (this == &rhs) { return *this; } 87 if (this == &rhs) { return *this; }
91 88
92 // Copy base class data 89 // Copy base class data
93 // 90 //
94 Base_t::operator=(rhs); << 91 G4USolid::operator=(rhs);
95 92
96 return *this; 93 return *this;
97 } 94 }
98 95
99 ////////////////////////////////////////////// 96 /////////////////////////////////////////////////////////////////////////
100 // 97 //
101 // Accessors and modifiers <<
102 <<
103 G4double G4UCons::GetInnerRadiusMinusZ() const <<
104 { <<
105 return GetRmin1(); <<
106 } <<
107 G4double G4UCons::GetOuterRadiusMinusZ() const <<
108 { <<
109 return GetRmax1(); <<
110 } <<
111 G4double G4UCons::GetInnerRadiusPlusZ() const <<
112 { <<
113 return GetRmin2(); <<
114 } <<
115 G4double G4UCons::GetOuterRadiusPlusZ() const <<
116 { <<
117 return GetRmax2(); <<
118 } <<
119 G4double G4UCons::GetZHalfLength() const <<
120 { <<
121 return GetDz(); <<
122 } <<
123 G4double G4UCons::GetStartPhiAngle() const <<
124 { <<
125 return GetSPhi(); <<
126 } <<
127 G4double G4UCons::GetDeltaPhiAngle() const <<
128 { <<
129 return GetDPhi(); <<
130 } <<
131 G4double G4UCons::GetSinStartPhi() const <<
132 { <<
133 G4double phi = GetStartPhiAngle(); <<
134 return std::sin(phi); <<
135 } <<
136 G4double G4UCons::GetCosStartPhi() const <<
137 { <<
138 G4double phi = GetStartPhiAngle(); <<
139 return std::cos(phi); <<
140 } <<
141 G4double G4UCons::GetSinEndPhi() const <<
142 { <<
143 G4double phi = GetStartPhiAngle() + GetDelta <<
144 return std::sin(phi); <<
145 } <<
146 G4double G4UCons::GetCosEndPhi() const <<
147 { <<
148 G4double phi = GetStartPhiAngle() + GetDelta <<
149 return std::cos(phi); <<
150 } <<
151 <<
152 void G4UCons::SetInnerRadiusMinusZ(G4double Rm <<
153 { <<
154 SetRmin1(Rmin1); <<
155 fRebuildPolyhedron = true; <<
156 } <<
157 void G4UCons::SetOuterRadiusMinusZ(G4double Rm <<
158 { <<
159 SetRmax1(Rmax1); <<
160 fRebuildPolyhedron = true; <<
161 } <<
162 void G4UCons::SetInnerRadiusPlusZ(G4double Rmi <<
163 { <<
164 SetRmin2(Rmin2); <<
165 fRebuildPolyhedron = true; <<
166 } <<
167 void G4UCons::SetOuterRadiusPlusZ(G4double Rma <<
168 { <<
169 SetRmax2(Rmax2); <<
170 fRebuildPolyhedron = true; <<
171 } <<
172 void G4UCons::SetZHalfLength(G4double newDz) <<
173 { <<
174 SetDz(newDz); <<
175 fRebuildPolyhedron = true; <<
176 } <<
177 void G4UCons::SetStartPhiAngle(G4double newSPh <<
178 { <<
179 SetSPhi(newSPhi); <<
180 fRebuildPolyhedron = true; <<
181 } <<
182 void G4UCons::SetDeltaPhiAngle(G4double newDPh <<
183 { <<
184 SetDPhi(newDPhi); <<
185 fRebuildPolyhedron = true; <<
186 } <<
187 <<
188 ////////////////////////////////////////////// <<
189 // <<
190 // Dispatch to parameterisation for replicatio 98 // Dispatch to parameterisation for replication mechanism dimension
191 // computation & modification. 99 // computation & modification.
192 100
193 void G4UCons::ComputeDimensions( G4VPVPar 101 void G4UCons::ComputeDimensions( G4VPVParameterisation* p,
194 const G4int 102 const G4int n,
195 const G4VPhysi 103 const G4VPhysicalVolume* pRep )
196 { 104 {
197 p->ComputeDimensions(*(G4Cons*)this,n,pRep); 105 p->ComputeDimensions(*(G4Cons*)this,n,pRep);
198 } 106 }
199 107
200 ////////////////////////////////////////////// 108 //////////////////////////////////////////////////////////////////////////
201 // 109 //
202 // Make a clone of the object 110 // Make a clone of the object
203 111
204 G4VSolid* G4UCons::Clone() const 112 G4VSolid* G4UCons::Clone() const
205 { 113 {
206 return new G4UCons(*this); 114 return new G4UCons(*this);
207 } 115 }
208 <<
209 ////////////////////////////////////////////// <<
210 // <<
211 // Get bounding box <<
212 <<
213 void G4UCons::BoundingLimits(G4ThreeVector& pM <<
214 { <<
215 static G4bool checkBBox = true; <<
216 <<
217 G4double rmin = std::min(GetInnerRadiusMinus <<
218 G4double rmax = std::max(GetOuterRadiusMinus <<
219 G4double dz = GetZHalfLength(); <<
220 <<
221 // Find bounding box <<
222 // <<
223 if (GetDeltaPhiAngle() < twopi) <<
224 { <<
225 G4TwoVector vmin,vmax; <<
226 G4GeomTools::DiskExtent(rmin,rmax, <<
227 GetSinStartPhi(),G <<
228 GetSinEndPhi(),Get <<
229 vmin,vmax); <<
230 pMin.set(vmin.x(),vmin.y(),-dz); <<
231 pMax.set(vmax.x(),vmax.y(), dz); <<
232 } <<
233 else <<
234 { <<
235 pMin.set(-rmax,-rmax,-dz); <<
236 pMax.set( rmax, rmax, dz); <<
237 } <<
238 <<
239 // Check correctness of the bounding box <<
240 // <<
241 if (pMin.x() >= pMax.x() || pMin.y() >= pMax <<
242 { <<
243 std::ostringstream message; <<
244 message << "Bad bounding box (min >= max) <<
245 << GetName() << " !" <<
246 << "\npMin = " << pMin <<
247 << "\npMax = " << pMax; <<
248 G4Exception("G4UCons::BoundingLimits()", " <<
249 JustWarning, message); <<
250 StreamInfo(G4cout); <<
251 } <<
252 <<
253 // Check consistency of bounding boxes <<
254 // <<
255 if (checkBBox) <<
256 { <<
257 U3Vector vmin, vmax; <<
258 Extent(vmin,vmax); <<
259 if (std::abs(pMin.x()-vmin.x()) > kCarTole <<
260 std::abs(pMin.y()-vmin.y()) > kCarTole <<
261 std::abs(pMin.z()-vmin.z()) > kCarTole <<
262 std::abs(pMax.x()-vmax.x()) > kCarTole <<
263 std::abs(pMax.y()-vmax.y()) > kCarTole <<
264 std::abs(pMax.z()-vmax.z()) > kCarTole <<
265 { <<
266 std::ostringstream message; <<
267 message << "Inconsistency in bounding bo <<
268 << GetName() << " !" <<
269 << "\nBBox min: wrapper = " << p <<
270 << "\nBBox max: wrapper = " << p <<
271 G4Exception("G4UCons::BoundingLimits()", <<
272 JustWarning, message); <<
273 checkBBox = false; <<
274 } <<
275 } <<
276 } <<
277 <<
278 ////////////////////////////////////////////// <<
279 // <<
280 // Calculate extent under transform and specif <<
281 <<
282 G4bool <<
283 G4UCons::CalculateExtent(const EAxis pAxis, <<
284 const G4VoxelLimits& <<
285 const G4AffineTransfo <<
286 G4double& pMin, <<
287 { <<
288 G4ThreeVector bmin, bmax; <<
289 G4bool exist; <<
290 <<
291 // Get bounding box <<
292 BoundingLimits(bmin,bmax); <<
293 <<
294 // Check bounding box <<
295 G4BoundingEnvelope bbox(bmin,bmax); <<
296 #ifdef G4BBOX_EXTENT <<
297 if (true) return bbox.CalculateExtent(pAxis, <<
298 #endif <<
299 if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVox <<
300 { <<
301 return exist = pMin < pMax; <<
302 } <<
303 <<
304 // Get parameters of the solid <<
305 G4double rmin1 = GetInnerRadiusMinusZ(); <<
306 G4double rmax1 = GetOuterRadiusMinusZ(); <<
307 G4double rmin2 = GetInnerRadiusPlusZ(); <<
308 G4double rmax2 = GetOuterRadiusPlusZ(); <<
309 G4double dz = GetZHalfLength(); <<
310 G4double dphi = GetDeltaPhiAngle(); <<
311 <<
312 // Find bounding envelope and calculate exte <<
313 // <<
314 const G4int NSTEPS = 24; // numbe <<
315 G4double astep = twopi/NSTEPS; // max a <<
316 G4int ksteps = (dphi <= astep) ? 1 : (G4i <<
317 G4double ang = dphi/ksteps; <<
318 <<
319 G4double sinHalf = std::sin(0.5*ang); <<
320 G4double cosHalf = std::cos(0.5*ang); <<
321 G4double sinStep = 2.*sinHalf*cosHalf; <<
322 G4double cosStep = 1. - 2.*sinHalf*sinHalf; <<
323 G4double rext1 = rmax1/cosHalf; <<
324 G4double rext2 = rmax2/cosHalf; <<
325 <<
326 // bounding envelope for full cone without h <<
327 // in other cases it is a sequence of quadri <<
328 if (rmin1 == 0 && rmin2 == 0 && dphi == twop <<
329 { <<
330 G4double sinCur = sinHalf; <<
331 G4double cosCur = cosHalf; <<
332 <<
333 G4ThreeVectorList baseA(NSTEPS),baseB(NSTE <<
334 for (G4int k=0; k<NSTEPS; ++k) <<
335 { <<
336 baseA[k].set(rext1*cosCur,rext1*sinCur,- <<
337 baseB[k].set(rext2*cosCur,rext2*sinCur, <<
338 <<
339 G4double sinTmp = sinCur; <<
340 sinCur = sinCur*cosStep + cosCur*sinStep <<
341 cosCur = cosCur*cosStep - sinTmp*sinStep <<
342 } <<
343 std::vector<const G4ThreeVectorList *> pol <<
344 polygons[0] = &baseA; <<
345 polygons[1] = &baseB; <<
346 G4BoundingEnvelope benv(bmin,bmax,polygons <<
347 exist = benv.CalculateExtent(pAxis,pVoxelL <<
348 } <<
349 else <<
350 { <<
351 G4double sinStart = GetSinStartPhi(); <<
352 G4double cosStart = GetCosStartPhi(); <<
353 G4double sinEnd = GetSinEndPhi(); <<
354 G4double cosEnd = GetCosEndPhi(); <<
355 G4double sinCur = sinStart*cosHalf + cos <<
356 G4double cosCur = cosStart*cosHalf - sin <<
357 <<
358 // set quadrilaterals <<
359 G4ThreeVectorList pols[NSTEPS+2]; <<
360 for (G4int k=0; k<ksteps+2; ++k) pols[k].r <<
361 pols[0][0].set(rmin2*cosStart,rmin2*sinSta <<
362 pols[0][1].set(rmin1*cosStart,rmin1*sinSta <<
363 pols[0][2].set(rmax1*cosStart,rmax1*sinSta <<
364 pols[0][3].set(rmax2*cosStart,rmax2*sinSta <<
365 for (G4int k=1; k<ksteps+1; ++k) <<
366 { <<
367 pols[k][0].set(rmin2*cosCur,rmin2*sinCur <<
368 pols[k][1].set(rmin1*cosCur,rmin1*sinCur <<
369 pols[k][2].set(rext1*cosCur,rext1*sinCur <<
370 pols[k][3].set(rext2*cosCur,rext2*sinCur <<
371 <<
372 G4double sinTmp = sinCur; <<
373 sinCur = sinCur*cosStep + cosCur*sinStep <<
374 cosCur = cosCur*cosStep - sinTmp*sinStep <<
375 } <<
376 pols[ksteps+1][0].set(rmin2*cosEnd,rmin2*s <<
377 pols[ksteps+1][1].set(rmin1*cosEnd,rmin1*s <<
378 pols[ksteps+1][2].set(rmax1*cosEnd,rmax1*s <<
379 pols[ksteps+1][3].set(rmax2*cosEnd,rmax2*s <<
380 <<
381 // set envelope and calculate extent <<
382 std::vector<const G4ThreeVectorList *> pol <<
383 polygons.resize(ksteps+2); <<
384 for (G4int k=0; k<ksteps+2; ++k) polygons[ <<
385 G4BoundingEnvelope benv(bmin,bmax,polygons <<
386 exist = benv.CalculateExtent(pAxis,pVoxelL <<
387 } <<
388 return exist; <<
389 } <<
390 <<
391 ////////////////////////////////////////////// <<
392 // <<
393 // Create polyhedron for visualization <<
394 <<
395 G4Polyhedron* G4UCons::CreatePolyhedron() cons <<
396 { <<
397 return new G4PolyhedronCons(GetInnerRadiusMi <<
398 GetOuterRadiusMi <<
399 GetInnerRadiusPl <<
400 GetOuterRadiusPl <<
401 GetZHalfLength() <<
402 GetStartPhiAngle <<
403 GetDeltaPhiAngle <<
404 } <<
405 <<
406 #endif // G4GEOM_USE_USOLIDS <<
407 116