Geant4  10.03
G4Box.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: G4Box.cc 99469 2016-09-22 15:04:36Z 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 // 20.09.16 - E.Tcherniaev: added Extent(pmin,pmax)
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 
135 {
136  if(dx > 2*kCarTolerance) // limit to thickness of surfaces
137  {
138  fDx = dx;
139  }
140  else
141  {
142  std::ostringstream message;
143  message << "Dimension X too small for solid: " << GetName() << "!"
144  << G4endl
145  << " hX = " << dx;
146  G4Exception("G4Box::SetXHalfLength()", "GeomSolids0002",
147  FatalException, message);
148  }
149  fCubicVolume= 0.;
150  fSurfaceArea= 0.;
151  fRebuildPolyhedron = true;
152 }
153 
155 {
156  if(dy > 2*kCarTolerance) // limit to thickness of surfaces
157  {
158  fDy = dy;
159  }
160  else
161  {
162  std::ostringstream message;
163  message << "Dimension Y too small for solid: " << GetName() << "!"
164  << G4endl
165  << " hY = " << dy;
166  G4Exception("G4Box::SetYHalfLength()", "GeomSolids0002",
167  FatalException, message);
168  }
169  fCubicVolume= 0.;
170  fSurfaceArea= 0.;
171  fRebuildPolyhedron = true;
172 }
173 
175 {
176  if(dz > 2*kCarTolerance) // limit to thickness of surfaces
177  {
178  fDz = dz;
179  }
180  else
181  {
182  std::ostringstream message;
183  message << "Dimension Z too small for solid: " << GetName() << "!"
184  << G4endl
185  << " hZ = " << dz;
186  G4Exception("G4Box::SetZHalfLength()", "GeomSolids0002",
187  FatalException, message);
188  }
189  fCubicVolume= 0.;
190  fSurfaceArea= 0.;
191  fRebuildPolyhedron = true;
192 }
193 
195 //
196 // Dispatch to parameterisation for replication mechanism dimension
197 // computation & modification.
198 
200  const G4int n,
201  const G4VPhysicalVolume* pRep)
202 {
203  p->ComputeDimensions(*this,n,pRep);
204 }
205 
207 //
208 // Get bounding box
209 
210 void G4Box::Extent(G4ThreeVector& pMin, G4ThreeVector& pMax) const
211 {
212  pMin.set(-fDx,-fDy,-fDz);
213  pMax.set( fDx, fDy, fDz);
214 
215  // Check correctness of the bounding box
216  //
217  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
218  {
219  std::ostringstream message;
220  message << "Bad bounding box (min >= max) for solid: "
221  << GetName() << " !"
222  << "\npMin = " << pMin
223  << "\npMax = " << pMax;
224  G4Exception("G4Box::Extent()", "GeomMgt0001", JustWarning, message);
225  DumpInfo();
226  }
227 }
228 
230 //
231 // Calculate extent under transform and specified limit
232 
234  const G4VoxelLimits& pVoxelLimit,
235  const G4AffineTransform& pTransform,
236  G4double& pMin, G4double& pMax) const
237 {
238  G4ThreeVector bmin, bmax;
239 
240  // Get bounding box
241  Extent(bmin,bmax);
242 
243  // Find extent
244  G4BoundingEnvelope bbox(bmin,bmax);
245  return bbox.CalculateExtent(pAxis,pVoxelLimit,pTransform,pMin,pMax);
246 }
247 
249 //
250 // Return whether point inside/outside/on surface, using tolerance
251 
253 {
254  EInside in = kOutside ;
255  G4ThreeVector q(std::fabs(p.x()), std::fabs(p.y()), std::fabs(p.z()));
256 
257  if ( q.x() <= (fDx - delta) )
258  {
259  if (q.y() <= (fDy - delta) )
260  {
261  if ( q.z() <= (fDz - delta) ) { in = kInside ; }
262  else if ( q.z() <= (fDz + delta) ) { in = kSurface ; }
263  }
264  else if ( q.y() <= (fDy + delta) )
265  {
266  if ( q.z() <= (fDz + delta) ) { in = kSurface ; }
267  }
268  }
269  else if ( q.x() <= (fDx + delta) )
270  {
271  if ( q.y() <= (fDy + delta) )
272  {
273  if ( q.z() <= (fDz + delta) ) { in = kSurface ; }
274  }
275  }
276  return in ;
277 }
278 
280 //
281 // Calculate side nearest to p, and return normal
282 // If two sides are equidistant, normal of first side (x/y/z)
283 // encountered returned
284 
286 {
287  G4double distx, disty, distz ;
288  G4ThreeVector norm(0.,0.,0.);
289 
290  // Calculate distances as if in 1st octant
291 
292  distx = std::fabs(std::fabs(p.x()) - fDx) ;
293  disty = std::fabs(std::fabs(p.y()) - fDy) ;
294  distz = std::fabs(std::fabs(p.z()) - fDz) ;
295 
296  // New code for particle on surface including edges and corners with specific
297  // normals
298 
299  const G4ThreeVector nX = G4ThreeVector( 1.0, 0,0 );
300  const G4ThreeVector nmX = G4ThreeVector(-1.0, 0,0 );
301  const G4ThreeVector nY = G4ThreeVector( 0, 1.0,0 );
302  const G4ThreeVector nmY = G4ThreeVector( 0,-1.0,0 );
303  const G4ThreeVector nZ = G4ThreeVector( 0, 0, 1.0);
304  const G4ThreeVector nmZ = G4ThreeVector( 0, 0,- 1.0);
305 
306  G4ThreeVector normX(0.,0.,0.), normY(0.,0.,0.), normZ(0.,0.,0.);
307  G4ThreeVector sumnorm(0., 0., 0.);
308  G4int noSurfaces=0;
309 
310  if (distx <= delta) // on X/mX surface and around
311  {
312  noSurfaces ++;
313  if ( p.x() >= 0. ) { normX= nX ; } // on +X surface : (1,0,0)
314  else { normX= nmX; } // (-1,0,0)
315  sumnorm= normX;
316  }
317 
318  if (disty <= delta) // on one of the +Y or -Y surfaces
319  {
320  noSurfaces ++;
321  if ( p.y() >= 0. ) { normY= nY; } // on +Y surface
322  else { normY= nmY; }
323  sumnorm += normY;
324  }
325 
326  if (distz <= delta) // on one of the +Z or -Z surfaces
327  {
328  noSurfaces ++;
329  if ( p.z() >= 0. ) { normZ= nZ; } // on +Z surface
330  else { normZ= nmZ; }
331  sumnorm += normZ;
332  }
333 
334  static const G4double invSqrt2 = 1.0 / std::sqrt(2.0);
335  static const G4double invSqrt3 = 1.0 / std::sqrt(3.0);
336 
337  if( noSurfaces > 0 )
338  {
339  if( noSurfaces == 1 )
340  {
341  norm= sumnorm;
342  }
343  else
344  {
345  // norm = sumnorm . unit();
346  if( noSurfaces == 2 )
347  {
348  // 2 surfaces -> on edge
349  norm = invSqrt2 * sumnorm;
350  }
351  else
352  {
353  // 3 surfaces (on corner)
354  norm = invSqrt3 * sumnorm;
355  }
356  }
357  }
358  else
359  {
360 #ifdef G4CSGDEBUG
361  G4Exception("G4Box::SurfaceNormal(p)", "Notification", JustWarning,
362  "Point p is not on surface !?" );
363 #endif
364  norm = ApproxSurfaceNormal(p);
365  }
366 
367  return norm;
368 }
369 
371 //
372 // Algorithm for SurfaceNormal() following the original specification
373 // for points not on the surface
374 
376 {
377  G4double distx, disty, distz ;
378  G4ThreeVector norm(0.,0.,0.);
379 
380  // Calculate distances as if in 1st octant
381 
382  distx = std::fabs(std::fabs(p.x()) - fDx) ;
383  disty = std::fabs(std::fabs(p.y()) - fDy) ;
384  distz = std::fabs(std::fabs(p.z()) - fDz) ;
385 
386  if ( distx <= disty )
387  {
388  if ( distx <= distz ) // Closest to X
389  {
390  if ( p.x() < 0 ) { norm = G4ThreeVector(-1.0,0,0) ; }
391  else { norm = G4ThreeVector( 1.0,0,0) ; }
392  }
393  else // Closest to Z
394  {
395  if ( p.z() < 0 ) { norm = G4ThreeVector(0,0,-1.0) ; }
396  else { norm = G4ThreeVector(0,0, 1.0) ; }
397  }
398  }
399  else
400  {
401  if ( disty <= distz ) // Closest to Y
402  {
403  if ( p.y() < 0 ) { norm = G4ThreeVector(0,-1.0,0) ; }
404  else { norm = G4ThreeVector(0, 1.0,0) ; }
405  }
406  else // Closest to Z
407  {
408  if ( p.z() < 0 ) { norm = G4ThreeVector(0,0,-1.0) ; }
409  else { norm = G4ThreeVector(0,0, 1.0) ; }
410  }
411  }
412  return norm;
413 }
414 
416 //
417 // Calculate distance to box from an outside point
418 // - return kInfinity if no intersection.
419 //
420 // ALGORITHM:
421 //
422 // Check that if point lies outside x/y/z extent of box, travel is towards
423 // the box (ie. there is a possibility of an intersection)
424 //
425 // Calculate pairs of minimum and maximum distances for x/y/z travel for
426 // intersection with the box's x/y/z extent.
427 // If there is a valid intersection, it is given by the maximum min distance
428 // (ie. distance to satisfy x/y/z intersections) *if* <= minimum max distance
429 // (ie. distance after which 1+ of x/y/z intersections not satisfied)
430 //
431 // NOTE:
432 //
433 // `Inside' safe - meaningful answers given if point is inside the exact
434 // shape.
435 
437  const G4ThreeVector& v) const
438 {
439  G4double safx, safy, safz ;
440  G4double smin=0.0, sminy, sminz ; // , sminx ;
441  G4double smax=kInfinity, smaxy, smaxz ; // , smaxx ; // they always > 0
442  G4double stmp ;
443  G4double sOut=kInfinity, sOuty=kInfinity, sOutz=kInfinity ;
444 
445  safx = std::fabs(p.x()) - fDx ; // minimum distance to x surface of shape
446  safy = std::fabs(p.y()) - fDy ;
447  safz = std::fabs(p.z()) - fDz ;
448 
449  // Will we intersect?
450  // If safx/y/z is >-tol/2 the point is outside/on the box's x/y/z extent.
451  // If both p.x/y/z and v.x/y/z repectively are both positive/negative,
452  // travel is in a direction away from the shape.
453 
454  if ( ((p.x()*v.x() >= 0.0) && (safx > -delta))
455  || ((p.y()*v.y() >= 0.0) && (safy > -delta))
456  || ((p.z()*v.z() >= 0.0) && (safz > -delta)) )
457  {
458  return kInfinity ; // travel away or parallel within tolerance
459  }
460 
461  // Compute min / max distances for x/y/z travel:
462  // X Planes
463 
464  if ( v.x() ) // != 0
465  {
466  stmp = 1.0/std::fabs(v.x()) ;
467 
468  if (safx >= 0.0)
469  {
470  smin = safx*stmp ;
471  smax = (fDx+std::fabs(p.x()))*stmp ;
472  }
473  else
474  {
475  if (v.x() < 0) { sOut = (fDx + p.x())*stmp ; }
476  else { sOut = (fDx - p.x())*stmp ; }
477  }
478  }
479 
480  // Y Planes
481 
482  if ( v.y() ) // != 0
483  {
484  stmp = 1.0/std::fabs(v.y()) ;
485 
486  if (safy >= 0.0)
487  {
488  sminy = safy*stmp ;
489  smaxy = (fDy+std::fabs(p.y()))*stmp ;
490 
491  if (sminy > smin) { smin=sminy ; }
492  if (smaxy < smax) { smax=smaxy ; }
493 
494  if (smin >= (smax-delta))
495  {
496  return kInfinity ; // touch XY corner
497  }
498  }
499  else
500  {
501  if (v.y() < 0) { sOuty = (fDy + p.y())*stmp ; }
502  else { sOuty = (fDy - p.y())*stmp ; }
503  if( sOuty < sOut ) { sOut = sOuty ; }
504  }
505  }
506 
507  // Z planes
508 
509  if ( v.z() ) // != 0
510  {
511  stmp = 1.0/std::fabs(v.z()) ;
512 
513  if ( safz >= 0.0 )
514  {
515  sminz = safz*stmp ;
516  smaxz = (fDz+std::fabs(p.z()))*stmp ;
517 
518  if (sminz > smin) { smin = sminz ; }
519  if (smaxz < smax) { smax = smaxz ; }
520 
521  if (smin >= (smax-delta))
522  {
523  return kInfinity ; // touch ZX or ZY corners
524  }
525  }
526  else
527  {
528  if (v.z() < 0) { sOutz = (fDz + p.z())*stmp ; }
529  else { sOutz = (fDz - p.z())*stmp ; }
530  if( sOutz < sOut ) { sOut = sOutz ; }
531  }
532  }
533 
534  if (sOut <= (smin + delta)) // travel over edge
535  {
536  return kInfinity ;
537  }
538  if (smin < delta) { smin = 0.0 ; }
539 
540  return smin ;
541 }
542 
544 //
545 // Appoximate distance to box.
546 // Returns largest perpendicular distance to the closest x/y/z sides of
547 // the box, which is the most fast estimation of the shortest distance to box
548 // - If inside return 0
549 
551 {
552  G4double safex, safey, safez, safe = 0.0 ;
553 
554  safex = std::fabs(p.x()) - fDx ;
555  safey = std::fabs(p.y()) - fDy ;
556  safez = std::fabs(p.z()) - fDz ;
557 
558  if (safex > safe) { safe = safex ; }
559  if (safey > safe) { safe = safey ; }
560  if (safez > safe) { safe = safez ; }
561 
562  return safe ;
563 }
564 
566 //
567 // Calculate distance to surface of box from inside
568 // by calculating distances to box's x/y/z planes.
569 // Smallest distance is exact distance to exiting.
570 // - Eliminate one side of each pair by considering direction of v
571 // - when leaving a surface & v.close, return 0
572 
574  const G4bool calcNorm,
575  G4bool *validNorm,G4ThreeVector *n) const
576 {
577  ESide side = kUndefined ;
578  G4double pdist,stmp,snxt=kInfinity;
579 
580  if (calcNorm) { *validNorm = true ; } // All normals are valid
581 
582  if (v.x() > 0) // X planes
583  {
584  pdist = fDx - p.x() ;
585 
586  if (pdist > delta)
587  {
588  snxt = pdist/v.x() ;
589  side = kPX ;
590  }
591  else
592  {
593  if (calcNorm) { *n = G4ThreeVector(1,0,0) ; }
594  return snxt = 0 ;
595  }
596  }
597  else if (v.x() < 0)
598  {
599  pdist = fDx + p.x() ;
600 
601  if (pdist > delta)
602  {
603  snxt = -pdist/v.x() ;
604  side = kMX ;
605  }
606  else
607  {
608  if (calcNorm) { *n = G4ThreeVector(-1,0,0) ; }
609  return snxt = 0 ;
610  }
611  }
612 
613  if (v.y() > 0) // Y planes
614  {
615  pdist = fDy-p.y();
616 
617  if (pdist > delta)
618  {
619  stmp = pdist/v.y();
620 
621  if (stmp < snxt)
622  {
623  snxt = stmp;
624  side = kPY;
625  }
626  }
627  else
628  {
629  if (calcNorm) { *n = G4ThreeVector(0,1,0) ; }
630  return snxt = 0 ;
631  }
632  }
633  else if (v.y() < 0)
634  {
635  pdist = fDy + p.y() ;
636 
637  if (pdist > delta)
638  {
639  stmp = -pdist/v.y();
640 
641  if ( stmp < snxt )
642  {
643  snxt = stmp;
644  side = kMY;
645  }
646  }
647  else
648  {
649  if (calcNorm) { *n = G4ThreeVector(0,-1,0) ; }
650  return snxt = 0 ;
651  }
652  }
653 
654  if (v.z() > 0) // Z planes
655  {
656  pdist = fDz-p.z();
657 
658  if ( pdist > delta )
659  {
660  stmp = pdist/v.z();
661 
662  if ( stmp < snxt )
663  {
664  snxt = stmp;
665  side = kPZ;
666  }
667  }
668  else
669  {
670  if (calcNorm) { *n = G4ThreeVector(0,0,1) ; }
671  return snxt = 0 ;
672  }
673  }
674  else if (v.z() < 0)
675  {
676  pdist = fDz + p.z();
677 
678  if ( pdist > delta )
679  {
680  stmp = -pdist/v.z();
681 
682  if ( stmp < snxt )
683  {
684  snxt = stmp;
685  side = kMZ;
686  }
687  }
688  else
689  {
690  if (calcNorm) { *n = G4ThreeVector(0,0,-1) ; }
691  return snxt = 0 ;
692  }
693  }
694 
695  if (calcNorm)
696  {
697  switch (side)
698  {
699  case kPX:
700  *n=G4ThreeVector(1,0,0);
701  break;
702  case kMX:
703  *n=G4ThreeVector(-1,0,0);
704  break;
705  case kPY:
706  *n=G4ThreeVector(0,1,0);
707  break;
708  case kMY:
709  *n=G4ThreeVector(0,-1,0);
710  break;
711  case kPZ:
712  *n=G4ThreeVector(0,0,1);
713  break;
714  case kMZ:
715  *n=G4ThreeVector(0,0,-1);
716  break;
717  default:
718  G4cout << G4endl;
719  DumpInfo();
720  std::ostringstream message;
721  G4int oldprc = message.precision(16);
722  message << "Undefined side for valid surface normal to solid."
723  << G4endl
724  << "Position:" << G4endl << G4endl
725  << "p.x() = " << p.x()/mm << " mm" << G4endl
726  << "p.y() = " << p.y()/mm << " mm" << G4endl
727  << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl
728  << "Direction:" << G4endl << G4endl
729  << "v.x() = " << v.x() << G4endl
730  << "v.y() = " << v.y() << G4endl
731  << "v.z() = " << v.z() << G4endl << G4endl
732  << "Proposed distance :" << G4endl << G4endl
733  << "snxt = " << snxt/mm << " mm" << G4endl;
734  message.precision(oldprc);
735  G4Exception("G4Box::DistanceToOut(p,v,..)", "GeomSolids1002",
736  JustWarning, message);
737  break;
738  }
739  }
740  return snxt;
741 }
742 
744 //
745 // Calculate exact shortest distance to any boundary from inside
746 // - If outside return 0
747 
749 {
750  G4double safx1,safx2,safy1,safy2,safz1,safz2,safe=0.0;
751 
752 #ifdef G4CSGDEBUG
753  if( Inside(p) == kOutside )
754  {
755  G4int oldprc = G4cout.precision(16) ;
756  G4cout << G4endl ;
757  DumpInfo();
758  G4cout << "Position:" << G4endl << G4endl ;
759  G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl ;
760  G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl ;
761  G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl ;
762  G4cout.precision(oldprc) ;
763  G4Exception("G4Box::DistanceToOut(p)", "GeomSolids1002",
764  JustWarning, "Point p is outside !?" );
765  }
766 #endif
767 
768  safx1 = fDx - p.x() ;
769  safx2 = fDx + p.x() ;
770  safy1 = fDy - p.y() ;
771  safy2 = fDy + p.y() ;
772  safz1 = fDz - p.z() ;
773  safz2 = fDz + p.z() ;
774 
775  // shortest Dist to any boundary now MIN(safx1,safx2,safy1..)
776 
777  if (safx2 < safx1) { safe = safx2; }
778  else { safe = safx1; }
779  if (safy1 < safe) { safe = safy1; }
780  if (safy2 < safe) { safe = safy2; }
781  if (safz1 < safe) { safe = safz1; }
782  if (safz2 < safe) { safe = safz2; }
783 
784  if (safe < 0) { safe = 0 ; }
785  return safe ;
786 }
787 
789 //
790 // GetEntityType
791 
793 {
794  return G4String("G4Box");
795 }
796 
798 //
799 // Stream object contents to an output stream
800 
801 std::ostream& G4Box::StreamInfo(std::ostream& os) const
802 {
803  G4int oldprc = os.precision(16);
804  os << "-----------------------------------------------------------\n"
805  << " *** Dump for solid - " << GetName() << " ***\n"
806  << " ===================================================\n"
807  << " Solid type: G4Box\n"
808  << " Parameters: \n"
809  << " half length X: " << fDx/mm << " mm \n"
810  << " half length Y: " << fDy/mm << " mm \n"
811  << " half length Z: " << fDz/mm << " mm \n"
812  << "-----------------------------------------------------------\n";
813  os.precision(oldprc);
814 
815  return os;
816 }
817 
819 //
820 // GetPointOnSurface
821 //
822 // Return a point (G4ThreeVector) randomly and uniformly selected
823 // on the solid surface
824 
826 {
827  G4double px, py, pz, select, sumS;
828  G4double Sxy = fDx*fDy, Sxz = fDx*fDz, Syz = fDy*fDz;
829 
830  sumS = Sxy + Sxz + Syz;
831  select = sumS*G4UniformRand();
832 
833  if( select < Sxy )
834  {
835  px = -fDx +2*fDx*G4UniformRand();
836  py = -fDy +2*fDy*G4UniformRand();
837 
838  if(G4UniformRand() > 0.5) { pz = fDz; }
839  else { pz = -fDz; }
840  }
841  else if ( ( select - Sxy ) < Sxz )
842  {
843  px = -fDx +2*fDx*G4UniformRand();
844  pz = -fDz +2*fDz*G4UniformRand();
845 
846  if(G4UniformRand() > 0.5) { py = fDy; }
847  else { py = -fDy; }
848  }
849  else
850  {
851  py = -fDy +2*fDy*G4UniformRand();
852  pz = -fDz +2*fDz*G4UniformRand();
853 
854  if(G4UniformRand() > 0.5) { px = fDx; }
855  else { px = -fDx; }
856  }
857  return G4ThreeVector(px,py,pz);
858 }
859 
861 //
862 // Make a clone of the object
863 //
865 {
866  return new G4Box(*this);
867 }
868 
870 //
871 // Methods for visualisation
872 
874 {
875  scene.AddSolid (*this);
876 }
877 
879 {
880  return G4VisExtent (-fDx, fDx, -fDy, fDy, -fDz, fDz);
881 }
882 
884 {
885  return new G4PolyhedronBox (fDx, fDy, fDz);
886 }
887 #endif
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
Definition: G4Box.cc:199
G4ThreeVector GetPointOnSurface() const
Definition: G4Box.cc:825
G4String GetName() const
void SetZHalfLength(G4double dz)
Definition: G4Box.cc:174
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
Definition: G4Box.cc:573
G4double fDy
Definition: G4Box.hh:146
static constexpr double mm
Definition: G4SIunits.hh:115
static const G4double kInfinity
Definition: geomdefs.hh:42
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
Definition: G4Box.cc:233
CLHEP::Hep3Vector G4ThreeVector
G4Box & operator=(const G4Box &rhs)
Definition: G4Box.cc:112
G4bool fRebuildPolyhedron
Definition: G4CSGSolid.hh:80
G4double delta
Definition: G4Box.hh:147
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
Definition: G4Box.hh:64
G4GeometryType GetEntityType() const
Definition: G4Box.cc:792
void DescribeYourselfTo(G4VGraphicsScene &scene) const
Definition: G4Box.cc:873
G4double fCubicVolume
Definition: G4CSGSolid.hh:78
G4Box(const G4String &pName, G4double pX, G4double pY, G4double pZ)
Definition: G4Box.cc:63
G4VisExtent GetExtent() const
Definition: G4Box.cc:878
virtual void AddSolid(const G4Box &)=0
G4VSolid * Clone() const
Definition: G4Box.cc:864
int G4int
Definition: G4Types.hh:78
G4Polyhedron * CreatePolyhedron() const
Definition: G4Box.cc:883
std::ostream & StreamInfo(std::ostream &os) const
Definition: G4Box.cc:801
const G4int smax
void DumpInfo() const
G4double fDz
Definition: G4Box.hh:146
#define G4UniformRand()
Definition: Randomize.hh:97
G4double fSurfaceArea
Definition: G4CSGSolid.hh:79
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
G4double fDx
Definition: G4Box.hh:146
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimits, const G4Transform3D &pTransform3D, G4double &pMin, G4double &pMax) const
const G4int n
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
ESide
Definition: G4Box.hh:135
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
Definition: G4Box.cc:285
virtual void ComputeDimensions(G4Box &, const G4int, const G4VPhysicalVolume *) const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Definition: G4Box.cc:436
EInside Inside(const G4ThreeVector &p) const
Definition: G4Box.cc:252
EInside
Definition: geomdefs.hh:58
virtual ~G4Box()
Definition: G4Box.cc:95
EAxis
Definition: geomdefs.hh:54
void SetYHalfLength(G4double dy)
Definition: G4Box.cc:154
void SetXHalfLength(G4double dx)
Definition: G4Box.cc:134
#define G4endl
Definition: G4ios.hh:61
G4double kCarTolerance
Definition: G4VSolid.hh:307
double G4double
Definition: G4Types.hh:76
G4CSGSolid & operator=(const G4CSGSolid &rhs)
Definition: G4CSGSolid.cc:91
G4ThreeVector ApproxSurfaceNormal(const G4ThreeVector &p) const
Definition: G4Box.cc:375
void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4Box.cc:210