Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4Box.cc
이 파일의 문서화 페이지로 가기
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: G4Box.cc 109330 2018-04-11 08:42:02Z gcosmo $
28 //
29 //
30 //
31 // Implementation for G4Box class
32 //
33 // 24.06.98 - V.Grichine: insideEdge in DistanceToIn(p,v)
34 // 20.09.98 - V.Grichine: new algorithm of DistanceToIn(p,v)
35 // 07.05.00 - V.Grichine: d= DistanceToIn(p,v), if d<e/2, d=0
36 // 09.06.00 - V.Grichine: safety in DistanceToIn(p) against Inside(p)=kOutside
37 // and information before exception in DistanceToOut(p,v,...)
38 // 15.11.00 - D.Williams, V.Grichine: bug fixed in CalculateExtent - change
39 // algorithm for rotated vertices
40 // 23.08.16 - E.Tcherniaev: use G4BoundingEnvelope for CalculateExtent()
41 // 18.04.17 - E.Tcherniaev: complete revision, speed-up
42 // --------------------------------------------------------------------
43 
44 #include "G4Box.hh"
45 
46 #if !defined(G4GEOM_USE_UBOX)
47 
48 #include "G4SystemOfUnits.hh"
49 #include "G4VoxelLimits.hh"
50 #include "G4AffineTransform.hh"
51 #include "G4BoundingEnvelope.hh"
52 #include "Randomize.hh"
53 
54 #include "G4VPVParameterisation.hh"
55 
56 #include "G4VGraphicsScene.hh"
57 #include "G4VisExtent.hh"
58 
60 //
61 // Constructor - check & set half widths
62 
63 G4Box::G4Box(const G4String& pName,
64  G4double pX,
65  G4double pY,
66  G4double pZ)
67  : G4CSGSolid(pName), fDx(pX), fDy(pY), fDz(pZ)
68 {
69  delta = 0.5*kCarTolerance;
70  if (pX < 2*kCarTolerance ||
71  pY < 2*kCarTolerance ||
72  pZ < 2*kCarTolerance) // limit to thickness of surfaces
73  {
74  std::ostringstream message;
75  message << "Dimensions too small for Solid: " << GetName() << "!" << G4endl
76  << " hX, hY, hZ = " << pX << ", " << pY << ", " << pZ;
77  G4Exception("G4Box::G4Box()", "GeomSolids0002", FatalException, message);
78  }
79 }
80 
82 //
83 // Fake default constructor - sets only member data and allocates memory
84 // for usage restricted to object persistency
85 
86 G4Box::G4Box( __void__& a )
87  : G4CSGSolid(a), fDx(0.), fDy(0.), fDz(0.), delta(0.)
88 {
89 }
90 
92 //
93 // Destructor
94 
96 {
97 }
98 
100 //
101 // Copy constructor
102 
103 G4Box::G4Box(const G4Box& rhs)
104  : G4CSGSolid(rhs), fDx(rhs.fDx), fDy(rhs.fDy), fDz(rhs.fDz), delta(rhs.delta)
105 {
106 }
107 
109 //
110 // Assignment operator
111 
113 {
114  // Check assignment to self
115  //
116  if (this == &rhs) { return *this; }
117 
118  // Copy base class data
119  //
121 
122  // Copy data
123  //
124  fDx = rhs.fDx;
125  fDy = rhs.fDy;
126  fDz = rhs.fDz;
127  delta = rhs.delta;
128 
129  return *this;
130 }
131 
133 //
134 // Set X dimension
135 
137 {
138  if(dx > 2*kCarTolerance) // limit to thickness of surfaces
139  {
140  fDx = dx;
141  }
142  else
143  {
144  std::ostringstream message;
145  message << "Dimension X too small for solid: " << GetName() << "!"
146  << G4endl
147  << " hX = " << dx;
148  G4Exception("G4Box::SetXHalfLength()", "GeomSolids0002",
149  FatalException, message);
150  }
151  fCubicVolume = 0.;
152  fSurfaceArea = 0.;
153  fRebuildPolyhedron = true;
154 }
155 
157 //
158 // Set Y dimension
159 
161 {
162  if(dy > 2*kCarTolerance) // limit to thickness of surfaces
163  {
164  fDy = dy;
165  }
166  else
167  {
168  std::ostringstream message;
169  message << "Dimension Y too small for solid: " << GetName() << "!\n"
170  << " hY = " << dy;
171  G4Exception("G4Box::SetYHalfLength()", "GeomSolids0002",
172  FatalException, message);
173  }
174  fCubicVolume = 0.;
175  fSurfaceArea = 0.;
176  fRebuildPolyhedron = true;
177 }
178 
180 //
181 // Set Z dimension
182 
184 {
185  if(dz > 2*kCarTolerance) // limit to thickness of surfaces
186  {
187  fDz = dz;
188  }
189  else
190  {
191  std::ostringstream message;
192  message << "Dimension Z too small for solid: " << GetName() << "!\n"
193  << " hZ = " << dz;
194  G4Exception("G4Box::SetZHalfLength()", "GeomSolids0002",
195  FatalException, message);
196  }
197  fCubicVolume = 0.;
198  fSurfaceArea = 0.;
199  fRebuildPolyhedron = true;
200 }
201 
203 //
204 // Dispatch to parameterisation for replication mechanism dimension
205 // computation & modification.
206 
208  const G4int n,
209  const G4VPhysicalVolume* pRep)
210 {
211  p->ComputeDimensions(*this,n,pRep);
212 }
213 
215 //
216 // Get bounding box
217 
219 {
220  pMin.set(-fDx,-fDy,-fDz);
221  pMax.set( fDx, fDy, fDz);
222 
223  // Check correctness of the bounding box
224  //
225  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
226  {
227  std::ostringstream message;
228  message << "Bad bounding box (min >= max) for solid: "
229  << GetName() << " !"
230  << "\npMin = " << pMin
231  << "\npMax = " << pMax;
232  G4Exception("G4Box::BoundingLimits()", "GeomMgt0001", JustWarning, message);
233  DumpInfo();
234  }
235 }
236 
238 //
239 // Calculate extent under transform and specified limit
240 
242  const G4VoxelLimits& pVoxelLimit,
243  const G4AffineTransform& pTransform,
244  G4double& pMin, G4double& pMax) const
245 {
246  G4ThreeVector bmin, bmax;
247 
248  // Get bounding box
249  BoundingLimits(bmin,bmax);
250 
251  // Find extent
252  G4BoundingEnvelope bbox(bmin,bmax);
253  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
254 }
255 
257 //
258 // Return whether point inside/outside/on surface, using tolerance
259 
261 {
262  G4double dist = std::max(std::max(
263  std::abs(p.x())-fDx,
264  std::abs(p.y())-fDy),
265  std::abs(p.z())-fDz);
266  if (dist > delta) return kOutside;
267  return (dist > -delta) ? kSurface : kInside;
268 }
269 
271 //
272 // Detect the side(s) and return corresponding normal
273 
275 {
276  G4ThreeVector norm(0,0,0);
277  G4double px = p.x();
278  if (std::abs(std::abs(px) - fDx) <= delta) norm.setX(px < 0 ? -1. : 1.);
279  G4double py = p.y();
280  if (std::abs(std::abs(py) - fDy) <= delta) norm.setY(py < 0 ? -1. : 1.);
281  G4double pz = p.z();
282  if (std::abs(std::abs(pz) - fDz) <= delta) norm.setZ(pz < 0 ? -1. : 1.);
283 
284  G4double nside = norm.mag2(); // number of sides = magnitude squared
285  if (nside == 1)
286  return norm;
287  else if (nside > 1)
288  return norm.unit(); // edge or corner
289  else
290  {
291  // Point is not on the surface
292  //
293 #ifdef G4CSGDEBUG
294  std::ostringstream message;
295  G4int oldprc = message.precision(16);
296  message << "Point p is not on surface (!?) of solid: "
297  << GetName() << G4endl;
298  message << "Position:\n";
299  message << " p.x() = " << p.x()/mm << " mm\n";
300  message << " p.y() = " << p.y()/mm << " mm\n";
301  message << " p.z() = " << p.z()/mm << " mm";
302  G4cout.precision(oldprc);
303  G4Exception("G4Box::SurfaceNormal(p)", "GeomSolids1002",
304  JustWarning, message );
305  DumpInfo();
306 #endif
307  return ApproxSurfaceNormal(p);
308  }
309 }
310 
312 //
313 // Algorithm for SurfaceNormal() following the original specification
314 // for points not on the surface
315 
317 {
318  G4double distx = std::abs(p.x()) - fDx;
319  G4double disty = std::abs(p.y()) - fDy;
320  G4double distz = std::abs(p.z()) - fDz;
321 
322  if (distx >= disty && distx >= distz)
323  return G4ThreeVector(std::copysign(1.,p.x()), 0., 0.);
324  if (disty >= distx && disty >= distz)
325  return G4ThreeVector(0., std::copysign(1.,p.y()), 0.);
326  else
327  return G4ThreeVector(0., 0., std::copysign(1.,p.z()));
328 }
329 
331 //
332 // Calculate distance to box from an outside point
333 // - return kInfinity if no intersection
334 //
335 
337  const G4ThreeVector& v) const
338 {
339  // Check if point is on the surface and traveling away
340  //
341  if ((std::abs(p.x()) - fDx) >= -delta && p.x()*v.x() >= 0) return kInfinity;
342  if ((std::abs(p.y()) - fDy) >= -delta && p.y()*v.y() >= 0) return kInfinity;
343  if ((std::abs(p.z()) - fDz) >= -delta && p.z()*v.z() >= 0) return kInfinity;
344 
345  // Find intersection
346  //
347  G4double invx = (v.x() == 0) ? DBL_MAX : -1./v.x();
348  G4double dx = std::copysign(fDx,invx);
349  G4double txmin = (p.x() - dx)*invx;
350  G4double txmax = (p.x() + dx)*invx;
351 
352  G4double invy = (v.y() == 0) ? DBL_MAX : -1./v.y();
353  G4double dy = std::copysign(fDy,invy);
354  G4double tymin = std::max(txmin,(p.y() - dy)*invy);
355  G4double tymax = std::min(txmax,(p.y() + dy)*invy);
356 
357  G4double invz = (v.z() == 0) ? DBL_MAX : -1./v.z();
358  G4double dz = std::copysign(fDz,invz);
359  G4double tmin = std::max(tymin,(p.z() - dz)*invz);
360  G4double tmax = std::min(tymax,(p.z() + dz)*invz);
361 
362  if (tmax <= tmin + delta) return kInfinity; // touch or no hit
363  return (tmin < delta) ? 0. : tmin;
364 }
365 
367 //
368 // Appoximate distance to box.
369 // Returns largest perpendicular distance to the closest x/y/z sides of
370 // the box, which is the most fast estimation of the shortest distance to box
371 // - If inside return 0
372 
374 {
375  G4double dist = std::max(std::max(
376  std::abs(p.x())-fDx,
377  std::abs(p.y())-fDy),
378  std::abs(p.z())-fDz);
379  return (dist > 0) ? dist : 0.;
380 }
381 
383 //
384 // Calculate distance to surface of the box from inside and
385 // find normal at exit point, if required
386 // - when leaving the surface, return 0
387 
389  const G4ThreeVector& v,
390  const G4bool calcNorm,
391  G4bool *validNorm, G4ThreeVector *n) const
392 {
393  // Check if point is on the surface and traveling away
394  //
395  if ((std::abs(p.x()) - fDx) >= -delta && p.x()*v.x() > 0)
396  {
397  if (calcNorm)
398  {
399  *validNorm = true;
400  n->set((p.x() < 0) ? -1. : 1., 0., 0.);
401  }
402  return 0.;
403  }
404  if ((std::abs(p.y()) - fDy) >= -delta && p.y()*v.y() > 0)
405  {
406  if (calcNorm)
407  {
408  *validNorm = true;
409  n->set(0., (p.y() < 0) ? -1. : 1., 0.);
410  }
411  return 0.;
412  }
413  if ((std::abs(p.z()) - fDz) >= -delta && p.z()*v.z() > 0)
414  {
415  if (calcNorm)
416  {
417  *validNorm = true;
418  n->set(0., 0., (p.z() < 0) ? -1. : 1.);
419  }
420  return 0.;
421  }
422 
423  // Find intersection
424  //
425  G4double vx = v.x();
426  G4double tx = (vx == 0) ? DBL_MAX : (std::copysign(fDx,vx) - p.x())/vx;
427 
428  G4double vy = v.y();
429  G4double ty = (vy == 0) ? tx : (std::copysign(fDy,vy) - p.y())/vy;
430  G4double txy = std::min(tx,ty);
431 
432  G4double vz = v.z();
433  G4double tz = (vz == 0) ? txy : (std::copysign(fDz,vz) - p.z())/vz;
434  G4double tmax = std::min(txy,tz);
435 
436  // Set normal, if required, and return distance
437  //
438  if (calcNorm)
439  {
440  *validNorm = true;
441  if (tmax == tx) n->set((v.x() < 0) ? -1. : 1., 0., 0.);
442  else if (tmax == ty) n->set(0., (v.y() < 0) ? -1. : 1., 0.);
443  else n->set(0., 0., (v.z() < 0) ? -1. : 1.);
444  }
445  return tmax;
446 }
447 
449 //
450 // Calculate exact shortest distance to any boundary from inside
451 // - if outside return 0
452 
454 {
455 #ifdef G4CSGDEBUG
456  if( Inside(p) == kOutside )
457  {
458  std::ostringstream message;
459  G4int oldprc = message.precision(16);
460  message << "Point p is outside (!?) of solid: " << GetName() << G4endl;
461  message << "Position:\n";
462  message << " p.x() = " << p.x()/mm << " mm\n";
463  message << " p.y() = " << p.y()/mm << " mm\n";
464  message << " p.z() = " << p.z()/mm << " mm";
465  G4cout.precision(oldprc);
466  G4Exception("G4Box::DistanceToOut(p)", "GeomSolids1002",
467  JustWarning, message );
468  DumpInfo();
469  }
470 #endif
471  G4double dist = std::min(std::min(
472  fDx-std::abs(p.x()),
473  fDy-std::abs(p.y())),
474  fDz-std::abs(p.z()));
475  return (dist > 0) ? dist : 0.;
476 }
477 
479 //
480 // GetEntityType
481 
483 {
484  return G4String("G4Box");
485 }
486 
488 //
489 // Stream object contents to an output stream
490 
491 std::ostream& G4Box::StreamInfo(std::ostream& os) const
492 {
493  G4int oldprc = os.precision(16);
494  os << "-----------------------------------------------------------\n"
495  << " *** Dump for solid - " << GetName() << " ***\n"
496  << " ===================================================\n"
497  << "Solid type: G4Box\n"
498  << "Parameters: \n"
499  << " half length X: " << fDx/mm << " mm \n"
500  << " half length Y: " << fDy/mm << " mm \n"
501  << " half length Z: " << fDz/mm << " mm \n"
502  << "-----------------------------------------------------------\n";
503  os.precision(oldprc);
504  return os;
505 }
506 
508 //
509 // GetPointOnSurface
510 //
511 // Return a point (G4ThreeVector) randomly and uniformly selected
512 // on the solid surface
513 
515 {
516  G4double sxy = fDx*fDy, sxz = fDx*fDz, syz = fDy*fDz;
517  G4double select = (sxy + sxz + syz)*G4UniformRand();
518 
519  if (select < sxy)
520  return G4ThreeVector((2.*G4UniformRand() - 1.)*fDx,
521  (2.*G4UniformRand() - 1.)*fDy,
522  (select < 0.5*sxy) ? -fDz : fDz);
523 
524  if (select < sxy + sxz)
525  return G4ThreeVector((2.*G4UniformRand() - 1.)*fDx,
526  (select < sxy + 0.5*sxz) ? -fDy : fDy,
527  (2.*G4UniformRand() - 1.)*fDz);
528  else
529  return G4ThreeVector((select < sxy + sxz + 0.5*syz) ? -fDx : fDx,
530  (2.*G4UniformRand() - 1.)*fDy,
531  (2.*G4UniformRand() - 1.)*fDz);
532 }
533 
535 //
536 // Make a clone of the object
537 //
539 {
540  return new G4Box(*this);
541 }
542 
544 //
545 // Methods for visualisation
546 
548 {
549  scene.AddSolid (*this);
550 }
551 
553 {
554  return G4VisExtent (-fDx, fDx, -fDy, fDy, -fDz, fDz);
555 }
556 
558 {
559  return new G4PolyhedronBox (fDx, fDy, fDz);
560 }
561 #endif
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimits, const G4Transform3D &pTransform3D, G4double &pMin, G4double &pMax) const
void SetXHalfLength(G4double dx)
Definition: G4Box.cc:136
void set(double x, double y, double z)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
Definition: G4Box.cc:274
EInside Inside(const G4ThreeVector &p) const
Definition: G4Box.cc:260
G4Polyhedron * CreatePolyhedron() const
Definition: G4Box.cc:557
CLHEP::Hep3Vector G4ThreeVector
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4Box.cc:336
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4Box(const G4String &pName, G4double pX, G4double pY, G4double pZ)
Definition: G4Box.cc:63
virtual void AddSolid(const G4Box &)=0
static const G4double kInfinity
Definition: geomdefs.hh:42
G4VSolid * Clone() const
Definition: G4Box.cc:538
static constexpr double mm
Definition: G4SIunits.hh:115
G4CSGSolid & operator=(const G4CSGSolid &rhs)
Definition: G4CSGSolid.cc:92
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
#define G4endl
Definition: G4ios.hh:61
std::ostream & StreamInfo(std::ostream &os) const
Definition: G4Box.cc:491
const char * p
Definition: xmltok.h:285
void SetZHalfLength(G4double dz)
Definition: G4Box.cc:183
void message(RunManager *runmanager)
Definition: ts_scorers.cc:72
double z() const
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
Definition: G4Box.cc:388
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
Definition: G4Box.cc:207
G4Box & operator=(const G4Box &rhs)
Definition: G4Box.cc:112
void setX(double)
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
Definition: G4Box.cc:241
G4VisExtent GetExtent() const
Definition: G4Box.cc:552
G4double fDx
Definition: G4Box.hh:141
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
void setZ(double)
G4ThreeVector ApproxSurfaceNormal(const G4ThreeVector &p) const
Definition: G4Box.cc:316
Definition: G4Box.hh:64
G4double kCarTolerance
Definition: G4VSolid.hh:307
#define G4UniformRand()
Definition: Randomize.hh:53
double mag2() const
void DescribeYourselfTo(G4VGraphicsScene &scene) const
Definition: G4Box.cc:547
G4GeometryType GetEntityType() const
Definition: G4Box.cc:482
Hep3Vector unit() const
void SetYHalfLength(G4double dy)
Definition: G4Box.cc:160
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
int G4int
Definition: G4Types.hh:78
EInside
Definition: geomdefs.hh:58
G4double fDz
Definition: G4Box.hh:141
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Box.cc:218
EAxis
Definition: geomdefs.hh:54
G4String GetName() const
void DumpInfo() const
Float_t norm
G4double delta
Definition: G4Box.hh:142
G4GLOB_DLL std::ostream G4cout
double x() const
G4ThreeVector GetPointOnSurface() const
Definition: G4Box.cc:514
Char_t n[5]
double y() const
G4double fSurfaceArea
Definition: G4CSGSolid.hh:79
virtual ~G4Box()
Definition: G4Box.cc:95
#define DBL_MAX
Definition: templates.hh:83
G4double fCubicVolume
Definition: G4CSGSolid.hh:78
G4bool fRebuildPolyhedron
Definition: G4CSGSolid.hh:80
G4double fDy
Definition: G4Box.hh:141
void setY(double)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments