Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4UnionSolid.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: G4UnionSolid.cc 109479 2018-04-24 14:33:49Z gcosmo $
28 //
29 // Implementation of methods for the class G4UnionSolid
30 //
31 // History:
32 //
33 // 12.09.98 V.Grichine: first implementation
34 // 28.11.98 V.Grichine: fix while loops in DistToIn/Out
35 // 27.07.99 V.Grichine: modifications in DistToOut(p,v,...), while -> do-while
36 // 16.03.01 V.Grichine: modifications in CalculateExtent()
37 // 17.03.17 E.Tcherniaev: revision of SurfaceNormal()
38 // 23.04.18 E.Tcherniaev: added extended BBox, yearly return in Inside()
39 //
40 // --------------------------------------------------------------------
41 
42 #include <sstream>
43 
44 #include "G4UnionSolid.hh"
45 
46 #include "G4SystemOfUnits.hh"
47 #include "G4VoxelLimits.hh"
48 #include "G4VPVParameterisation.hh"
49 #include "G4GeometryTolerance.hh"
50 
51 #include "G4VGraphicsScene.hh"
52 #include "G4Polyhedron.hh"
53 #include "HepPolyhedronProcessor.h"
54 
56 //
57 // Transfer all data members to G4BooleanSolid which is responsible
58 // for them. pName will be in turn sent to G4VSolid
59 
61  G4VSolid* pSolidA ,
62  G4VSolid* pSolidB )
63  : G4BooleanSolid(pName,pSolidA,pSolidB)
64 {
66  G4ThreeVector pmin, pmax;
67  BoundingLimits(pmin, pmax);
68  fPMin = pmin - pdelta;
69  fPMax = pmax + pdelta;
70 }
71 
73 //
74 // Constructor
75 
77  G4VSolid* pSolidA ,
78  G4VSolid* pSolidB ,
79  G4RotationMatrix* rotMatrix,
80  const G4ThreeVector& transVector )
81  : G4BooleanSolid(pName,pSolidA,pSolidB,rotMatrix,transVector)
82 
83 {
85  G4ThreeVector pmin, pmax;
86  BoundingLimits(pmin, pmax);
87  fPMin = pmin - pdelta;
88  fPMax = pmax + pdelta;
89 }
90 
92 //
93 // Constructor
94 
96  G4VSolid* pSolidA ,
97  G4VSolid* pSolidB ,
98  const G4Transform3D& transform )
99  : G4BooleanSolid(pName,pSolidA,pSolidB,transform)
100 {
102  G4ThreeVector pmin, pmax;
103  BoundingLimits(pmin, pmax);
104  fPMin = pmin - pdelta;
105  fPMax = pmax + pdelta;
106 }
107 
109 //
110 // Fake default constructor - sets only member data and allocates memory
111 // for usage restricted to object persistency.
112 
114  : G4BooleanSolid(a)
115 {
116 }
117 
119 //
120 // Destructor
121 
123 {
124 }
125 
127 //
128 // Copy constructor
129 
131  : G4BooleanSolid (rhs)
132 {
133  fPMin = rhs.fPMin;
134  fPMax = rhs.fPMax;
135 }
136 
138 //
139 // Assignment operator
140 
142 {
143  // Check assignment to self
144  //
145  if (this == &rhs) { return *this; }
146 
147  // Copy base class data
148  //
150 
151  fPMin = rhs.fPMin;
152  fPMax = rhs.fPMax;
153  return *this;
154 }
155 
157 //
158 // Get bounding box
159 
161  G4ThreeVector& pMax) const
162 {
163  G4ThreeVector minA,maxA, minB,maxB;
164  fPtrSolidA->BoundingLimits(minA,maxA);
165  fPtrSolidB->BoundingLimits(minB,maxB);
166 
167  pMin.set(std::min(minA.x(),minB.x()),
168  std::min(minA.y(),minB.y()),
169  std::min(minA.z(),minB.z()));
170 
171  pMax.set(std::max(maxA.x(),maxB.x()),
172  std::max(maxA.y(),maxB.y()),
173  std::max(maxA.z(),maxB.z()));
174 
175  // Check correctness of the bounding box
176  //
177  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
178  {
179  std::ostringstream message;
180  message << "Bad bounding box (min >= max) for solid: "
181  << GetName() << " !"
182  << "\npMin = " << pMin
183  << "\npMax = " << pMax;
184  G4Exception("G4UnionSolid::BoundingLimits()", "GeomMgt0001",
185  JustWarning, message);
186  DumpInfo();
187  }
188 }
189 
191 //
192 // Calculate extent under transform and specified limit
193 
194 G4bool
196  const G4VoxelLimits& pVoxelLimit,
197  const G4AffineTransform& pTransform,
198  G4double& pMin,
199  G4double& pMax ) const
200 {
201  G4bool touchesA, touchesB, out ;
202  G4double minA = kInfinity, minB = kInfinity,
203  maxA = -kInfinity, maxB = -kInfinity;
204 
205  touchesA = fPtrSolidA->CalculateExtent( pAxis, pVoxelLimit,
206  pTransform, minA, maxA);
207  touchesB = fPtrSolidB->CalculateExtent( pAxis, pVoxelLimit,
208  pTransform, minB, maxB);
209  if( touchesA || touchesB )
210  {
211  pMin = std::min( minA, minB );
212  pMax = std::max( maxA, maxB );
213  out = true ;
214  }
215  else
216  {
217  out = false ;
218  }
219 
220  return out ; // It exists in this slice if either one does.
221 }
222 
224 //
225 // Important comment: When solids A and B touch together along flat
226 // surface the surface points will be considered as kSurface, while points
227 // located around will correspond to kInside
228 
230 {
231  if (std::max(p.z()-fPMax.z(),fPMin.z()-p.z()) > 0) return kOutside;
232 
233  EInside positionA = fPtrSolidA->Inside(p);
234  if (positionA == kInside) { return positionA; } // inside A
235  EInside positionB = fPtrSolidB->Inside(p);
236  if (positionA == kOutside) { return positionB; }
237 
238  if (positionB == kInside) { return positionB; } // inside B
239  if (positionB == kOutside) { return positionA; } // surface A
240 
241  // Both points are on surface
242  //
243  static const G4double rtol
245 
246  return ((fPtrSolidA->SurfaceNormal(p) +
247  fPtrSolidB->SurfaceNormal(p)).mag2() < rtol) ? kInside : kSurface;
248 }
249 
251 //
252 // Get surface normal
253 
256 {
257  EInside positionA = fPtrSolidA->Inside(p);
258  EInside positionB = fPtrSolidB->Inside(p);
259 
260  if (positionA == kSurface &&
261  positionB == kOutside) return fPtrSolidA->SurfaceNormal(p);
262 
263  if (positionA == kOutside &&
264  positionB == kSurface) return fPtrSolidB->SurfaceNormal(p);
265 
266  if (positionA == kSurface &&
267  positionB == kSurface)
268  {
269  if (Inside(p) == kSurface)
270  {
271  G4ThreeVector normalA = fPtrSolidA->SurfaceNormal(p);
272  G4ThreeVector normalB = fPtrSolidB->SurfaceNormal(p);
273  return (normalA + normalB).unit();
274  }
275  }
276 #ifdef G4BOOLDEBUG
277  G4String surf[3] = { "OUTSIDE", "SURFACE", "INSIDE" };
278  std::ostringstream message;
279  G4int oldprc = message.precision(16);
280  message << "Invalid call of SurfaceNormal(p) for union solid: "
281  << GetName() << " !"
282  << "\nPoint p" << p << " is " << surf[Inside(p)] << " !!!";
283  message.precision(oldprc);
284  G4Exception("G4UnionSolid::SurfaceNormal()", "GeomMgt0001",
285  JustWarning, message);
286 #endif
287  return fPtrSolidA->SurfaceNormal(p);
288 }
289 
291 //
292 // The same algorithm as in DistanceToIn(p)
293 
294 G4double
296  const G4ThreeVector& v ) const
297 {
298 #ifdef G4BOOLDEBUG
299  if( Inside(p) == kInside )
300  {
301  G4cout << "WARNING - Invalid call in "
302  << "G4UnionSolid::DistanceToIn(p,v)" << G4endl
303  << " Point p is inside !" << G4endl;
304  G4cout << " p = " << p << G4endl;
305  G4cout << " v = " << v << G4endl;
306  G4cerr << "WARNING - Invalid call in "
307  << "G4UnionSolid::DistanceToIn(p,v)" << G4endl
308  << " Point p is inside !" << G4endl;
309  G4cerr << " p = " << p << G4endl;
310  G4cerr << " v = " << v << G4endl;
311  }
312 #endif
313 
314  return std::min(fPtrSolidA->DistanceToIn(p,v),
315  fPtrSolidB->DistanceToIn(p,v) ) ;
316 }
317 
319 //
320 // Approximate nearest distance from the point p to the union of
321 // two solids
322 
323 G4double
325 {
326 #ifdef G4BOOLDEBUG
327  if( Inside(p) == kInside )
328  {
329  G4cout << "WARNING - Invalid call in "
330  << "G4UnionSolid::DistanceToIn(p)" << G4endl
331  << " Point p is inside !" << G4endl;
332  G4cout << " p = " << p << G4endl;
333  G4cerr << "WARNING - Invalid call in "
334  << "G4UnionSolid::DistanceToIn(p)" << G4endl
335  << " Point p is inside !" << G4endl;
336  G4cerr << " p = " << p << G4endl;
337  }
338 #endif
341  G4double safety = std::min(distA,distB) ;
342  if(safety < 0.0) safety = 0.0 ;
343  return safety ;
344 }
345 
347 //
348 // The same algorithm as DistanceToOut(p)
349 
350 G4double
352  const G4ThreeVector& v,
353  const G4bool calcNorm,
354  G4bool *validNorm,
355  G4ThreeVector *n ) const
356 {
357  G4double dist = 0.0, disTmp = 0.0 ;
358  G4ThreeVector normTmp;
359  G4ThreeVector* nTmp= &normTmp;
360 
361  if( Inside(p) == kOutside )
362  {
363 #ifdef G4BOOLDEBUG
364  G4cout << "Position:" << G4endl << G4endl;
365  G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl;
366  G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl;
367  G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl;
368  G4cout << "Direction:" << G4endl << G4endl;
369  G4cout << "v.x() = " << v.x() << G4endl;
370  G4cout << "v.y() = " << v.y() << G4endl;
371  G4cout << "v.z() = " << v.z() << G4endl << G4endl;
372  G4cout << "WARNING - Invalid call in "
373  << "G4UnionSolid::DistanceToOut(p,v)" << G4endl
374  << " Point p is outside !" << G4endl;
375  G4cout << " p = " << p << G4endl;
376  G4cout << " v = " << v << G4endl;
377  G4cerr << "WARNING - Invalid call in "
378  << "G4UnionSolid::DistanceToOut(p,v)" << G4endl
379  << " Point p is outside !" << G4endl;
380  G4cerr << " p = " << p << G4endl;
381  G4cerr << " v = " << v << G4endl;
382 #endif
383  }
384  else
385  {
386  EInside positionA = fPtrSolidA->Inside(p) ;
387  // EInside positionB = fPtrSolidB->Inside(p) ;
388 
389  if( positionA != kOutside )
390  {
391  do // Loop checking, 13.08.2015, G.Cosmo
392  {
393  disTmp = fPtrSolidA->DistanceToOut(p+dist*v,v,calcNorm,
394  validNorm,nTmp);
395  dist += disTmp ;
396 
397  if(fPtrSolidB->Inside(p+dist*v) != kOutside)
398  {
399  disTmp = fPtrSolidB->DistanceToOut(p+dist*v,v,calcNorm,
400  validNorm,nTmp);
401  dist += disTmp ;
402  }
403  }
404  while( (fPtrSolidA->Inside(p+dist*v) != kOutside)
405  && (disTmp > 0.5*kCarTolerance) );
406  }
407  else // if( positionB != kOutside )
408  {
409  do // Loop checking, 13.08.2015, G.Cosmo
410  {
411  disTmp = fPtrSolidB->DistanceToOut(p+dist*v,v,calcNorm,
412  validNorm,nTmp);
413  dist += disTmp ;
414 
415  if(fPtrSolidA->Inside(p+dist*v) != kOutside)
416  {
417  disTmp = fPtrSolidA->DistanceToOut(p+dist*v,v,calcNorm,
418  validNorm,nTmp);
419  dist += disTmp ;
420  }
421  }
422  while( (fPtrSolidB->Inside(p+dist*v) != kOutside)
423  && (disTmp > 0.5*kCarTolerance) );
424  }
425  }
426  if( calcNorm )
427  {
428  *validNorm = false ;
429  *n = *nTmp ;
430  }
431  return dist ;
432 }
433 
435 //
436 // Inverted algorithm of DistanceToIn(p)
437 
438 G4double
440 {
441  G4double distout = 0.0;
442  if( Inside(p) == kOutside )
443  {
444 #ifdef G4BOOLDEBUG
445  G4cout << "WARNING - Invalid call in "
446  << "G4UnionSolid::DistanceToOut(p)" << G4endl
447  << " Point p is outside !" << G4endl;
448  G4cout << " p = " << p << G4endl;
449  G4cerr << "WARNING - Invalid call in "
450  << "G4UnionSolid::DistanceToOut(p)" << G4endl
451  << " Point p is outside !" << G4endl;
452  G4cerr << " p = " << p << G4endl;
453 #endif
454  }
455  else
456  {
457  EInside positionA = fPtrSolidA->Inside(p) ;
458  EInside positionB = fPtrSolidB->Inside(p) ;
459 
460  // Is this equivalent ??
461  // if( ! ( (positionA == kOutside)) &&
462  // (positionB == kOutside)) )
463  if((positionA == kInside && positionB == kInside ) ||
464  (positionA == kInside && positionB == kSurface ) ||
465  (positionA == kSurface && positionB == kInside ) )
466  {
467  distout= std::max(fPtrSolidA->DistanceToOut(p),
468  fPtrSolidB->DistanceToOut(p) ) ;
469  }
470  else
471  {
472  if(positionA == kOutside)
473  {
474  distout= fPtrSolidB->DistanceToOut(p) ;
475  }
476  else
477  {
478  distout= fPtrSolidA->DistanceToOut(p) ;
479  }
480  }
481  }
482  return distout;
483 }
484 
486 //
487 //
488 
490 {
491  return G4String("G4UnionSolid");
492 }
493 
495 //
496 // Make a clone of the object
497 
499 {
500  return new G4UnionSolid(*this);
501 }
502 
504 //
505 //
506 
507 void
509  const G4int,
510  const G4VPhysicalVolume* )
511 {
512 }
513 
515 //
516 //
517 
518 void
520 {
521  scene.AddSolid (*this);
522 }
523 
525 //
526 //
527 
528 G4Polyhedron*
530 {
532  // Stack components and components of components recursively
533  // See G4BooleanSolid::StackPolyhedron
534  G4Polyhedron* top = StackPolyhedron(processor, this);
535  G4Polyhedron* result = new G4Polyhedron(*top);
536  if (processor.execute(*result)) { return result; }
537  else { return 0; }
538 }
void set(double x, double y, double z)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4UnionSolid & operator=(const G4UnionSolid &rhs)
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
virtual void AddSolid(const G4Box &)=0
G4Polyhedron * CreatePolyhedron() const
static const G4double kInfinity
Definition: geomdefs.hh:42
G4ThreeVector fPMax
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
static constexpr double mm
Definition: G4SIunits.hh:115
virtual void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:625
#define G4endl
Definition: G4ios.hh:61
const char * p
Definition: xmltok.h:285
G4UnionSolid(const G4String &pName, G4VSolid *pSolidA, G4VSolid *pSolidB)
Definition: G4UnionSolid.cc:60
void message(RunManager *runmanager)
Definition: ts_scorers.cc:72
G4BooleanSolid & operator=(const G4BooleanSolid &rhs)
double z() const
G4double GetRadialTolerance() const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
virtual EInside Inside(const G4ThreeVector &p) const =0
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
G4double kCarTolerance
Definition: G4VSolid.hh:307
bool execute(HepPolyhedron &)
static const G4int maxA
Float_t distB
Definition: plotDend.C:22
G4GLOB_DLL std::ostream G4cerr
Hep3Vector unit() const
G4double G4ParticleHPJENDLHEData::G4double result
#define processor
Definition: xmlparse.cc:617
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
G4ThreeVector fPMin
G4Polyhedron * StackPolyhedron(HepPolyhedronProcessor &, const G4VSolid *) const
int G4int
Definition: G4Types.hh:78
EInside
Definition: geomdefs.hh:58
G4VSolid * fPtrSolidA
G4VSolid * fPtrSolidB
EAxis
Definition: geomdefs.hh:54
G4String GetName() const
void DumpInfo() const
G4GeometryType GetEntityType() const
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4GLOB_DLL std::ostream G4cout
double x() const
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
static G4GeometryTolerance * GetInstance()
Char_t n[5]
double y() const
Float_t distA
Definition: plotDend.C:22
virtual ~G4UnionSolid()
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
EInside Inside(const G4ThreeVector &p) const
G4VSolid * Clone() const