| 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 G4UOrb wrapper class << 26 // $Id:$
27 // 27 //
28 // 30.10.13 G.Cosmo, CERN/PH << 28 //
>> 29 // Implementation for G4UOrb wrapper class
29 // ------------------------------------------- 30 // --------------------------------------------------------------------
30 31
31 #include "G4Orb.hh" 32 #include "G4Orb.hh"
32 #include "G4UOrb.hh" 33 #include "G4UOrb.hh"
33 34
34 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G 35 #if ( defined(G4GEOM_USE_USOLIDS) || defined(G4GEOM_USE_PARTIAL_USOLIDS) )
35 36
36 #include "G4TwoVector.hh" 37 #include "G4TwoVector.hh"
37 #include "G4AffineTransform.hh" 38 #include "G4AffineTransform.hh"
38 #include "G4GeometryTolerance.hh" 39 #include "G4GeometryTolerance.hh"
39 #include "G4BoundingEnvelope.hh" 40 #include "G4BoundingEnvelope.hh"
40 41
41 #include "G4VPVParameterisation.hh" 42 #include "G4VPVParameterisation.hh"
42 #include "G4PhysicalConstants.hh" 43 #include "G4PhysicalConstants.hh"
43 44
44 using namespace CLHEP; 45 using namespace CLHEP;
45 46
46 ////////////////////////////////////////////// 47 ////////////////////////////////////////////////////////////////////////
47 // 48 //
48 // constructor - check positive radius 49 // constructor - check positive radius
49 // 50 //
50 51
51 G4UOrb::G4UOrb( const G4String& pName, G4doubl 52 G4UOrb::G4UOrb( const G4String& pName, G4double pRmax )
52 : Base_t(pName, pRmax) << 53 : G4USolid(pName, new UOrb(pName, pRmax))
53 { 54 {
54 } 55 }
55 56
56 ////////////////////////////////////////////// 57 ///////////////////////////////////////////////////////////////////////
57 // 58 //
58 // Fake default constructor - sets only member 59 // Fake default constructor - sets only member data and allocates memory
59 // for usage restri 60 // for usage restricted to object persistency.
60 // 61 //
61 G4UOrb::G4UOrb( __void__& a ) 62 G4UOrb::G4UOrb( __void__& a )
62 : Base_t(a) << 63 : G4USolid(a)
63 { 64 {
64 } 65 }
65 66
66 ////////////////////////////////////////////// 67 /////////////////////////////////////////////////////////////////////
67 // 68 //
68 // Destructor 69 // Destructor
69 70
70 G4UOrb::~G4UOrb() = default; << 71 G4UOrb::~G4UOrb()
>> 72 {
>> 73 }
71 74
72 ////////////////////////////////////////////// 75 //////////////////////////////////////////////////////////////////////////
73 // 76 //
74 // Copy constructor 77 // Copy constructor
75 78
76 G4UOrb::G4UOrb(const G4UOrb& rhs) 79 G4UOrb::G4UOrb(const G4UOrb& rhs)
77 : Base_t(rhs) << 80 : G4USolid(rhs)
78 { 81 {
79 } 82 }
80 83
81 ////////////////////////////////////////////// 84 //////////////////////////////////////////////////////////////////////////
82 // 85 //
83 // Assignment operator 86 // Assignment operator
84 87
85 G4UOrb& G4UOrb::operator = (const G4UOrb& rhs) 88 G4UOrb& G4UOrb::operator = (const G4UOrb& rhs)
86 { 89 {
87 // Check assignment to self 90 // Check assignment to self
88 // 91 //
89 if (this == &rhs) { return *this; } 92 if (this == &rhs) { return *this; }
90 93
91 // Copy base class data 94 // Copy base class data
92 // 95 //
93 Base_t::operator=(rhs); << 96 G4USolid::operator=(rhs);
94 97
95 return *this; 98 return *this;
96 } 99 }
97 100
98 ////////////////////////////////////////////// 101 //////////////////////////////////////////////////////////////////////////
99 // 102 //
100 // Accessors & modifiers 103 // Accessors & modifiers
101 104
102 G4double G4UOrb::GetRadius() const 105 G4double G4UOrb::GetRadius() const
103 { 106 {
104 return Base_t::GetRadius(); << 107 return GetShape()->GetRadius();
105 } 108 }
106 109
107 void G4UOrb::SetRadius(G4double newRmax) 110 void G4UOrb::SetRadius(G4double newRmax)
108 { 111 {
109 Base_t::SetRadius(newRmax); << 112 GetShape()->SetRadius(newRmax);
110 fRebuildPolyhedron = true; 113 fRebuildPolyhedron = true;
111 } 114 }
112 115
113 G4double G4UOrb::GetRadialTolerance() const <<
114 { <<
115 return Base_t::GetRadialTolerance(); <<
116 } <<
117 <<
118 ////////////////////////////////////////////// 116 //////////////////////////////////////////////////////////////////////////
119 // 117 //
120 // Dispatch to parameterisation for replicatio 118 // Dispatch to parameterisation for replication mechanism dimension
121 // computation & modification. 119 // computation & modification.
122 120
123 void G4UOrb::ComputeDimensions( G4VPVPar 121 void G4UOrb::ComputeDimensions( G4VPVParameterisation* p,
124 const G4int n, 122 const G4int n,
125 const G4VPhysic 123 const G4VPhysicalVolume* pRep )
126 { 124 {
127 p->ComputeDimensions(*(G4Orb*)this,n,pRep); 125 p->ComputeDimensions(*(G4Orb*)this,n,pRep);
128 } 126 }
129 127
130 ////////////////////////////////////////////// 128 //////////////////////////////////////////////////////////////////////////
131 // 129 //
132 // Make a clone of the object 130 // Make a clone of the object
133 131
134 G4VSolid* G4UOrb::Clone() const 132 G4VSolid* G4UOrb::Clone() const
135 { 133 {
136 return new G4UOrb(*this); 134 return new G4UOrb(*this);
137 } 135 }
138 136
139 ////////////////////////////////////////////// 137 //////////////////////////////////////////////////////////////////////////
140 // 138 //
141 // Get bounding box 139 // Get bounding box
142 140
143 void G4UOrb::BoundingLimits(G4ThreeVector& pMi << 141 void G4UOrb::Extent(G4ThreeVector& pMin, G4ThreeVector& pMax) const
144 { 142 {
145 G4double radius = GetRadius(); 143 G4double radius = GetRadius();
146 pMin.set(-radius,-radius,-radius); 144 pMin.set(-radius,-radius,-radius);
147 pMax.set( radius, radius, radius); 145 pMax.set( radius, radius, radius);
148 146
149 // Check correctness of the bounding box 147 // Check correctness of the bounding box
150 // 148 //
151 if (pMin.x() >= pMax.x() || pMin.y() >= pMax 149 if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
152 { 150 {
153 std::ostringstream message; 151 std::ostringstream message;
154 message << "Bad bounding box (min >= max) 152 message << "Bad bounding box (min >= max) for solid: "
155 << GetName() << " !" 153 << GetName() << " !"
156 << "\npMin = " << pMin 154 << "\npMin = " << pMin
157 << "\npMax = " << pMax; 155 << "\npMax = " << pMax;
158 G4Exception("G4UOrb::BoundingLimits()", "G << 156 G4Exception("G4UOrb::Extent()", "GeomMgt0001", JustWarning, message);
159 JustWarning, message); <<
160 StreamInfo(G4cout); 157 StreamInfo(G4cout);
161 } 158 }
162 } 159 }
163 160
164 ////////////////////////////////////////////// 161 //////////////////////////////////////////////////////////////////////////
165 // 162 //
166 // Calculate extent under transform and specif 163 // Calculate extent under transform and specified limit
167 164
168 G4bool 165 G4bool
169 G4UOrb::CalculateExtent(const EAxis pAxis, 166 G4UOrb::CalculateExtent(const EAxis pAxis,
170 const G4VoxelLimits& p 167 const G4VoxelLimits& pVoxelLimit,
171 const G4AffineTransfor 168 const G4AffineTransform& pTransform,
172 G4double& pMin, G4doub 169 G4double& pMin, G4double& pMax) const
173 { 170 {
174 G4ThreeVector bmin, bmax; 171 G4ThreeVector bmin, bmax;
175 G4bool exist; 172 G4bool exist;
176 173
177 // Get bounding box 174 // Get bounding box
178 BoundingLimits(bmin,bmax); << 175 Extent(bmin,bmax);
179 176
180 // Check bounding box 177 // Check bounding box
181 G4BoundingEnvelope bbox(bmin,bmax); 178 G4BoundingEnvelope bbox(bmin,bmax);
182 #ifdef G4BBOX_EXTENT 179 #ifdef G4BBOX_EXTENT
183 if (true) return bbox.CalculateExtent(pAxis, 180 if (true) return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
184 #endif 181 #endif
185 if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVox 182 if (bbox.BoundingBoxVsVoxelLimits(pAxis,pVoxelLimit,pTransform,pMin,pMax))
186 { 183 {
187 return exist = pMin < pMax; << 184 return exist = (pMin < pMax) ? true : false;
188 } 185 }
189 186
190 // Find bounding envelope and calculate exte 187 // Find bounding envelope and calculate extent
191 // 188 //
192 static const G4int NTHETA = 8; // number of 189 static const G4int NTHETA = 8; // number of steps along Theta
193 static const G4int NPHI = 16; // number of 190 static const G4int NPHI = 16; // number of steps along Phi
194 static const G4double sinHalfTheta = std::si 191 static const G4double sinHalfTheta = std::sin(halfpi/NTHETA);
195 static const G4double cosHalfTheta = std::co 192 static const G4double cosHalfTheta = std::cos(halfpi/NTHETA);
196 static const G4double sinHalfPhi = std::si 193 static const G4double sinHalfPhi = std::sin(pi/NPHI);
197 static const G4double cosHalfPhi = std::co 194 static const G4double cosHalfPhi = std::cos(pi/NPHI);
198 static const G4double sinStepTheta = 2.*sinH 195 static const G4double sinStepTheta = 2.*sinHalfTheta*cosHalfTheta;
199 static const G4double cosStepTheta = 1. - 2. 196 static const G4double cosStepTheta = 1. - 2.*sinHalfTheta*sinHalfTheta;
200 static const G4double sinStepPhi = 2.*sinH 197 static const G4double sinStepPhi = 2.*sinHalfPhi*cosHalfPhi;
201 static const G4double cosStepPhi = 1. - 2. 198 static const G4double cosStepPhi = 1. - 2.*sinHalfPhi*sinHalfPhi;
202 199
203 G4double radius = GetRadius(); 200 G4double radius = GetRadius();
204 G4double rtheta = radius/cosHalfTheta; 201 G4double rtheta = radius/cosHalfTheta;
205 G4double rphi = rtheta/cosHalfPhi; 202 G4double rphi = rtheta/cosHalfPhi;
206 203
207 // set reference circle 204 // set reference circle
208 G4TwoVector xy[NPHI]; 205 G4TwoVector xy[NPHI];
209 G4double sinCurPhi = sinHalfPhi; 206 G4double sinCurPhi = sinHalfPhi;
210 G4double cosCurPhi = cosHalfPhi; 207 G4double cosCurPhi = cosHalfPhi;
211 for (auto & k : xy) << 208 for (G4int k=0; k<NPHI; ++k)
212 { 209 {
213 k.set(cosCurPhi,sinCurPhi); << 210 xy[k].set(cosCurPhi,sinCurPhi);
214 G4double sinTmpPhi = sinCurPhi; 211 G4double sinTmpPhi = sinCurPhi;
215 sinCurPhi = sinCurPhi*cosStepPhi + cosCurP 212 sinCurPhi = sinCurPhi*cosStepPhi + cosCurPhi*sinStepPhi;
216 cosCurPhi = cosCurPhi*cosStepPhi - sinTmpP 213 cosCurPhi = cosCurPhi*cosStepPhi - sinTmpPhi*sinStepPhi;
217 } 214 }
218 215
219 // set bounding circles 216 // set bounding circles
220 G4ThreeVectorList circles[NTHETA]; 217 G4ThreeVectorList circles[NTHETA];
221 for (auto & circle : circles) circle.resize( << 218 for (G4int i=0; i<NTHETA; ++i) circles[i].resize(NPHI);
222 219
223 G4double sinCurTheta = sinHalfTheta; 220 G4double sinCurTheta = sinHalfTheta;
224 G4double cosCurTheta = cosHalfTheta; 221 G4double cosCurTheta = cosHalfTheta;
225 for (auto & circle : circles) << 222 for (G4int i=0; i<NTHETA; ++i)
226 { 223 {
227 G4double z = rtheta*cosCurTheta; 224 G4double z = rtheta*cosCurTheta;
228 G4double rho = rphi*sinCurTheta; 225 G4double rho = rphi*sinCurTheta;
229 for (G4int k=0; k<NPHI; ++k) 226 for (G4int k=0; k<NPHI; ++k)
230 { 227 {
231 circle[k].set(rho*xy[k].x(),rho*xy[k].y( << 228 circles[i][k].set(rho*xy[k].x(),rho*xy[k].y(),z);
232 } 229 }
233 G4double sinTmpTheta = sinCurTheta; 230 G4double sinTmpTheta = sinCurTheta;
234 sinCurTheta = sinCurTheta*cosStepTheta + c 231 sinCurTheta = sinCurTheta*cosStepTheta + cosCurTheta*sinStepTheta;
235 cosCurTheta = cosCurTheta*cosStepTheta - s 232 cosCurTheta = cosCurTheta*cosStepTheta - sinTmpTheta*sinStepTheta;
236 } 233 }
237 234
238 // set envelope and calculate extent 235 // set envelope and calculate extent
239 std::vector<const G4ThreeVectorList *> polyg 236 std::vector<const G4ThreeVectorList *> polygons;
240 polygons.resize(NTHETA); 237 polygons.resize(NTHETA);
241 for (G4int i=0; i<NTHETA; ++i) polygons[i] = 238 for (G4int i=0; i<NTHETA; ++i) polygons[i] = &circles[i];
242 239
243 G4BoundingEnvelope benv(bmin,bmax,polygons); 240 G4BoundingEnvelope benv(bmin,bmax,polygons);
244 exist = benv.CalculateExtent(pAxis,pVoxelLim 241 exist = benv.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
245 return exist; 242 return exist;
246 } 243 }
247 244
248 ////////////////////////////////////////////// 245 //////////////////////////////////////////////////////////////////////////
249 // 246 //
250 // Create polyhedron for visualization 247 // Create polyhedron for visualization
251 248
252 G4Polyhedron* G4UOrb::CreatePolyhedron() const 249 G4Polyhedron* G4UOrb::CreatePolyhedron() const
253 { 250 {
254 return new G4PolyhedronSphere(0., GetRadius( 251 return new G4PolyhedronSphere(0., GetRadius(), 0., twopi, 0., pi);
255 } 252 }
256 253
257 #endif // G4GEOM_USE_USOLIDS 254 #endif // G4GEOM_USE_USOLIDS
258 255