| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // Implementation for G4Box class
27 //
28 // 30.06.95 - P.Kent: First version
29 // 20.09.98 - V.Grichine: new algorithm of Di
30 // 18.04.17 - E.Tcherniaev: complete revision
31 // -------------------------------------------
32
33 #include "G4Box.hh"
34
35 #if !defined(G4GEOM_USE_UBOX)
36
37 #include "G4SystemOfUnits.hh"
38 #include "G4VoxelLimits.hh"
39 #include "G4AffineTransform.hh"
40 #include "G4BoundingEnvelope.hh"
41 #include "G4QuickRand.hh"
42
43 #include "G4VPVParameterisation.hh"
44
45 #include "G4VGraphicsScene.hh"
46 #include "G4VisExtent.hh"
47
48 //////////////////////////////////////////////
49 //
50 // Constructor - check & set half widths
51
52 G4Box::G4Box(const G4String& pName,
53 G4double pX,
54 G4double pY,
55 G4double pZ)
56 : G4CSGSolid(pName), fDx(pX), fDy(pY), fDz(p
57 {
58 delta = 0.5*kCarTolerance;
59 if (pX < 2*kCarTolerance ||
60 pY < 2*kCarTolerance ||
61 pZ < 2*kCarTolerance) // limit to thick
62 {
63 std::ostringstream message;
64 message << "Dimensions too small for Solid
65 << " hX, hY, hZ = " << pX << "
66 G4Exception("G4Box::G4Box()", "GeomSolids0
67 }
68 }
69
70 //////////////////////////////////////////////
71 //
72 // Fake default constructor - sets only member
73 // for usage restri
74
75 G4Box::G4Box( __void__& a )
76 : G4CSGSolid(a)
77 {
78 }
79
80 //////////////////////////////////////////////
81 //
82 // Destructor
83
84 G4Box::~G4Box() = default;
85
86 //////////////////////////////////////////////
87 //
88 // Copy constructor
89
90 G4Box::G4Box(const G4Box&) = default;
91
92 //////////////////////////////////////////////
93 //
94 // Assignment operator
95
96 G4Box& G4Box::operator = (const G4Box& rhs)
97 {
98 // Check assignment to self
99 //
100 if (this == &rhs) { return *this; }
101
102 // Copy base class data
103 //
104 G4CSGSolid::operator=(rhs);
105
106 // Copy data
107 //
108 fDx = rhs.fDx;
109 fDy = rhs.fDy;
110 fDz = rhs.fDz;
111 delta = rhs.delta;
112
113 return *this;
114 }
115
116 //////////////////////////////////////////////
117 //
118 // Set X dimension
119
120 void G4Box::SetXHalfLength(G4double dx)
121 {
122 if(dx > 2*kCarTolerance) // limit to thickn
123 {
124 fDx = dx;
125 }
126 else
127 {
128 std::ostringstream message;
129 message << "Dimension X too small for soli
130 << G4endl
131 << " hX = " << dx;
132 G4Exception("G4Box::SetXHalfLength()", "Ge
133 FatalException, message);
134 }
135 fCubicVolume = 0.;
136 fSurfaceArea = 0.;
137 fRebuildPolyhedron = true;
138 }
139
140 //////////////////////////////////////////////
141 //
142 // Set Y dimension
143
144 void G4Box::SetYHalfLength(G4double dy)
145 {
146 if(dy > 2*kCarTolerance) // limit to thickn
147 {
148 fDy = dy;
149 }
150 else
151 {
152 std::ostringstream message;
153 message << "Dimension Y too small for soli
154 << " hY = " << dy;
155 G4Exception("G4Box::SetYHalfLength()", "Ge
156 FatalException, message);
157 }
158 fCubicVolume = 0.;
159 fSurfaceArea = 0.;
160 fRebuildPolyhedron = true;
161 }
162
163 //////////////////////////////////////////////
164 //
165 // Set Z dimension
166
167 void G4Box::SetZHalfLength(G4double dz)
168 {
169 if(dz > 2*kCarTolerance) // limit to thickn
170 {
171 fDz = dz;
172 }
173 else
174 {
175 std::ostringstream message;
176 message << "Dimension Z too small for soli
177 << " hZ = " << dz;
178 G4Exception("G4Box::SetZHalfLength()", "Ge
179 FatalException, message);
180 }
181 fCubicVolume = 0.;
182 fSurfaceArea = 0.;
183 fRebuildPolyhedron = true;
184 }
185
186 //////////////////////////////////////////////
187 //
188 // Dispatch to parameterisation for replicatio
189 // computation & modification.
190
191 void G4Box::ComputeDimensions(G4VPVParameteris
192 const G4int n,
193 const G4VPhysica
194 {
195 p->ComputeDimensions(*this,n,pRep);
196 }
197
198 //////////////////////////////////////////////
199 //
200 // Get bounding box
201
202 void G4Box::BoundingLimits(G4ThreeVector& pMin
203 {
204 pMin.set(-fDx,-fDy,-fDz);
205 pMax.set( fDx, fDy, fDz);
206
207 // Check correctness of the bounding box
208 //
209 if (pMin.x() >= pMax.x() || pMin.y() >= pMax
210 {
211 std::ostringstream message;
212 message << "Bad bounding box (min >= max)
213 << GetName() << " !"
214 << "\npMin = " << pMin
215 << "\npMax = " << pMax;
216 G4Exception("G4Box::BoundingLimits()", "Ge
217 DumpInfo();
218 }
219 }
220
221 //////////////////////////////////////////////
222 //
223 // Calculate extent under transform and specif
224
225 G4bool G4Box::CalculateExtent(const EAxis pAxi
226 const G4VoxelLim
227 const G4AffineTr
228 G4double&
229 {
230 G4ThreeVector bmin, bmax;
231
232 // Get bounding box
233 BoundingLimits(bmin,bmax);
234
235 // Find extent
236 G4BoundingEnvelope bbox(bmin,bmax);
237 return bbox.CalculateExtent(pAxis,pVoxelLimi
238 }
239
240 //////////////////////////////////////////////
241 //
242 // Return whether point inside/outside/on surf
243
244 EInside G4Box::Inside(const G4ThreeVector& p)
245 {
246 G4double dist = std::max(std::max(
247 std::abs(p.x())-fDx,
248 std::abs(p.y())-fDy),
249 std::abs(p.z())-fDz);
250 return (dist > delta) ? kOutside :
251 ((dist > -delta) ? kSurface : kInside);
252 }
253
254 //////////////////////////////////////////////
255 //
256 // Detect the side(s) and return corresponding
257
258 G4ThreeVector G4Box::SurfaceNormal( const G4Th
259 {
260 G4ThreeVector norm(0,0,0);
261 G4double px = p.x();
262 if (std::abs(std::abs(px) - fDx) <= delta) n
263 G4double py = p.y();
264 if (std::abs(std::abs(py) - fDy) <= delta) n
265 G4double pz = p.z();
266 if (std::abs(std::abs(pz) - fDz) <= delta) n
267
268 G4double nside = norm.mag2(); // number of s
269 if (nside == 1)
270 return norm;
271 else if (nside > 1)
272 return norm.unit(); // edge or corner
273 else
274 {
275 // Point is not on the surface
276 //
277 #ifdef G4CSGDEBUG
278 std::ostringstream message;
279 G4int oldprc = message.precision(16);
280 message << "Point p is not on surface (!?)
281 << GetName() << G4endl;
282 message << "Position:\n";
283 message << " p.x() = " << p.x()/mm << "
284 message << " p.y() = " << p.y()/mm << "
285 message << " p.z() = " << p.z()/mm << "
286 G4cout.precision(oldprc);
287 G4Exception("G4Box::SurfaceNormal(p)", "Ge
288 JustWarning, message );
289 DumpInfo();
290 #endif
291 return ApproxSurfaceNormal(p);
292 }
293 }
294
295 //////////////////////////////////////////////
296 //
297 // Algorithm for SurfaceNormal() following the
298 // for points not on the surface
299
300 G4ThreeVector G4Box::ApproxSurfaceNormal(const
301 {
302 G4double distx = std::abs(p.x()) - fDx;
303 G4double disty = std::abs(p.y()) - fDy;
304 G4double distz = std::abs(p.z()) - fDz;
305
306 if (distx >= disty && distx >= distz)
307 return {std::copysign(1.,p.x()), 0., 0.};
308 if (disty >= distx && disty >= distz)
309 return {0., std::copysign(1.,p.y()), 0.};
310 else
311 return {0., 0., std::copysign(1.,p.z())};
312 }
313
314 //////////////////////////////////////////////
315 //
316 // Calculate distance to box from an outside p
317 // - return kInfinity if no intersection
318 //
319
320 G4double G4Box::DistanceToIn(const G4ThreeVect
321 const G4ThreeVect
322 {
323 // Check if point is on the surface and trav
324 //
325 if ((std::abs(p.x()) - fDx) >= -delta && p.x
326 if ((std::abs(p.y()) - fDy) >= -delta && p.y
327 if ((std::abs(p.z()) - fDz) >= -delta && p.z
328
329 // Find intersection
330 //
331 G4double invx = (v.x() == 0) ? DBL_MAX : -1.
332 G4double dx = std::copysign(fDx,invx);
333 G4double txmin = (p.x() - dx)*invx;
334 G4double txmax = (p.x() + dx)*invx;
335
336 G4double invy = (v.y() == 0) ? DBL_MAX : -1.
337 G4double dy = std::copysign(fDy,invy);
338 G4double tymin = std::max(txmin,(p.y() - dy)
339 G4double tymax = std::min(txmax,(p.y() + dy)
340
341 G4double invz = (v.z() == 0) ? DBL_MAX : -1.
342 G4double dz = std::copysign(fDz,invz);
343 G4double tmin = std::max(tymin,(p.z() - dz)*
344 G4double tmax = std::min(tymax,(p.z() + dz)*
345
346 if (tmax <= tmin + delta) return kInfinity;
347 return (tmin < delta) ? 0. : tmin;
348 }
349
350 //////////////////////////////////////////////
351 //
352 // Appoximate distance to box.
353 // Returns largest perpendicular distance to t
354 // the box, which is the most fast estimation
355 // - If inside return 0
356
357 G4double G4Box::DistanceToIn(const G4ThreeVect
358 {
359 G4double dist = std::max(std::max(
360 std::abs(p.x())-fDx,
361 std::abs(p.y())-fDy),
362 std::abs(p.z())-fDz);
363 return (dist > 0) ? dist : 0.;
364 }
365
366 //////////////////////////////////////////////
367 //
368 // Calculate distance to surface of the box fr
369 // find normal at exit point, if required
370 // - when leaving the surface, return 0
371
372 G4double G4Box::DistanceToOut( const G4ThreeVe
373 const G4ThreeVe
374 const G4bool ca
375 G4bool* validNo
376 {
377 // Check if point is on the surface and trav
378 //
379 if ((std::abs(p.x()) - fDx) >= -delta && p.x
380 {
381 if (calcNorm)
382 {
383 *validNorm = true;
384 n->set((p.x() < 0) ? -1. : 1., 0., 0.);
385 }
386 return 0.;
387 }
388 if ((std::abs(p.y()) - fDy) >= -delta && p.y
389 {
390 if (calcNorm)
391 {
392 *validNorm = true;
393 n->set(0., (p.y() < 0) ? -1. : 1., 0.);
394 }
395 return 0.;
396 }
397 if ((std::abs(p.z()) - fDz) >= -delta && p.z
398 {
399 if (calcNorm)
400 {
401 *validNorm = true;
402 n->set(0., 0., (p.z() < 0) ? -1. : 1.);
403 }
404 return 0.;
405 }
406
407 // Find intersection
408 //
409 G4double vx = v.x();
410 G4double tx = (vx == 0) ? DBL_MAX : (std::co
411
412 G4double vy = v.y();
413 G4double ty = (vy == 0) ? tx : (std::copysig
414 G4double txy = std::min(tx,ty);
415
416 G4double vz = v.z();
417 G4double tz = (vz == 0) ? txy : (std::copysi
418 G4double tmax = std::min(txy,tz);
419
420 // Set normal, if required, and return dista
421 //
422 if (calcNorm)
423 {
424 *validNorm = true;
425 if (tmax == tx) n->set((v.x() < 0) ?
426 else if (tmax == ty) n->set(0., (v.y() < 0
427 else n->set(0., 0., (v.z()
428 }
429 return tmax;
430 }
431
432 //////////////////////////////////////////////
433 //
434 // Calculate exact shortest distance to any bo
435 // - if outside return 0
436
437 G4double G4Box::DistanceToOut(const G4ThreeVec
438 {
439 #ifdef G4CSGDEBUG
440 if( Inside(p) == kOutside )
441 {
442 std::ostringstream message;
443 G4int oldprc = message.precision(16);
444 message << "Point p is outside (!?) of sol
445 message << "Position:\n";
446 message << " p.x() = " << p.x()/mm << "
447 message << " p.y() = " << p.y()/mm << "
448 message << " p.z() = " << p.z()/mm << "
449 G4cout.precision(oldprc);
450 G4Exception("G4Box::DistanceToOut(p)", "Ge
451 JustWarning, message );
452 DumpInfo();
453 }
454 #endif
455 G4double dist = std::min(std::min(
456 fDx-std::abs(p.x()),
457 fDy-std::abs(p.y())),
458 fDz-std::abs(p.z()));
459 return (dist > 0) ? dist : 0.;
460 }
461
462 //////////////////////////////////////////////
463 //
464 // GetEntityType
465
466 G4GeometryType G4Box::GetEntityType() const
467 {
468 return {"G4Box"};
469 }
470
471 //////////////////////////////////////////////
472 //
473 // IsFaceted
474
475 G4bool G4Box::IsFaceted() const
476 {
477 return true;
478 }
479
480 //////////////////////////////////////////////
481 //
482 // Stream object contents to an output stream
483
484 std::ostream& G4Box::StreamInfo(std::ostream&
485 {
486 G4long oldprc = os.precision(16);
487 os << "-------------------------------------
488 << " *** Dump for solid - " << GetName
489 << " =================================
490 << "Solid type: G4Box\n"
491 << "Parameters: \n"
492 << " half length X: " << fDx/mm << " mm
493 << " half length Y: " << fDy/mm << " mm
494 << " half length Z: " << fDz/mm << " mm
495 << "-------------------------------------
496 os.precision(oldprc);
497 return os;
498 }
499
500 //////////////////////////////////////////////
501 //
502 // Return a point randomly and uniformly selec
503
504 G4ThreeVector G4Box::GetPointOnSurface() const
505 {
506 G4double sxy = fDx*fDy, sxz = fDx*fDz, syz =
507 G4double select = (sxy + sxz + syz)*G4QuickR
508 G4double u = 2.*G4QuickRand() - 1.;
509 G4double v = 2.*G4QuickRand() - 1.;
510
511 if (select < sxy)
512 return { u*fDx, v*fDy, ((select < 0.5*sxy)
513 else if (select < sxy + sxz)
514 return { u*fDx, ((select < sxy + 0.5*sxz)
515 else
516 return { ((select < sxy + sxz + 0.5*syz) ?
517 }
518
519 //////////////////////////////////////////////
520 //
521 // Make a clone of the object
522 //
523 G4VSolid* G4Box::Clone() const
524 {
525 return new G4Box(*this);
526 }
527
528 //////////////////////////////////////////////
529 //
530 // Methods for visualisation
531
532 void G4Box::DescribeYourselfTo (G4VGraphicsSce
533 {
534 scene.AddSolid (*this);
535 }
536
537 G4VisExtent G4Box::GetExtent() const
538 {
539 return { -fDx, fDx, -fDy, fDy, -fDz, fDz };
540 }
541
542 G4Polyhedron* G4Box::CreatePolyhedron () const
543 {
544 return new G4PolyhedronBox (fDx, fDy, fDz);
545 }
546 #endif
547