Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4ReflectedSolid.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: G4ReflectedSolid.cc 104317 2017-05-24 13:08:38Z gcosmo $
28 //
29 //
30 // Implementation for G4ReflectedSolid class
31 //
32 // Author: Vladimir Grichine, 23.07.01 (Vladimir.Grichine@cern.ch)
33 //
34 // --------------------------------------------------------------------
35 
36 #include "G4ReflectedSolid.hh"
37 
38 #include <sstream>
39 
40 #include "G4Point3D.hh"
41 #include "G4Vector3D.hh"
42 
43 #include "G4AffineTransform.hh"
44 #include "G4Transform3D.hh"
45 #include "G4VoxelLimits.hh"
46 
47 #include "G4VPVParameterisation.hh"
48 
49 #include "G4VGraphicsScene.hh"
50 #include "G4Polyhedron.hh"
51 
53 //
54 // Constructor using HepTransform3D, in fact HepReflect3D
55 
57  G4VSolid* pSolid ,
58  const G4Transform3D& transform )
59  : G4VSolid(pName), fRebuildPolyhedron(false), fpPolyhedron(0)
60 {
61  fPtrSolid = pSolid;
62  fDirectTransform3D = new G4Transform3D(transform);
63 }
64 
66 //
67 
69 {
71  delete fpPolyhedron; fpPolyhedron = 0;
72 }
73 
75 //
76 
78  : G4VSolid(rhs), fPtrSolid(rhs.fPtrSolid),
79  fRebuildPolyhedron(false), fpPolyhedron(0)
80 {
82 }
83 
85 //
86 
88 {
89  // Check assignment to self
90  //
91  if (this == &rhs) { return *this; }
92 
93  // Copy base class data
94  //
96 
97  // Copy data
98  //
99  fPtrSolid= rhs.fPtrSolid;
100  delete fDirectTransform3D;
102  fRebuildPolyhedron = false;
103  delete fpPolyhedron; fpPolyhedron= 0;
104 
105  return *this;
106 }
107 
109 //
110 
112 {
113  return G4String("G4ReflectedSolid");
114 }
115 
117 {
118  return this;
119 }
120 
122 {
123  return this;
124 }
125 
127 {
128  return fPtrSolid;
129 }
130 
132 //
133 
135 {
136  return fDirectTransform3D->inverse();
137 }
138 
140 {
141  G4Transform3D aTransform= *fDirectTransform3D;
142  return aTransform;
143 }
144 
146 {
147  fDirectTransform3D = &transform;
148  fRebuildPolyhedron = true;
149 }
150 
152 //
153 // Get bounding box
154 
156  G4ThreeVector& pMax) const
157 {
158  fPtrSolid->BoundingLimits(pMin,pMax);
159  G4double xmin = pMin.x(), ymin = pMin.y(), zmin = pMin.z();
160  G4double xmax = pMax.x(), ymax = pMax.y(), zmax = pMax.z();
164 
165  if (std::abs(xx) == 1 && std::abs(yy) == 1 && std::abs(zz) == 1)
166  {
167  // Special case of reflection in axis and pure translation
168  //
169  if (xx == -1) { G4double tmp = -xmin; xmin = -xmax; xmax = tmp; }
170  if (yy == -1) { G4double tmp = -ymin; ymin = -ymax; ymax = tmp; }
171  if (zz == -1) { G4double tmp = -zmin; zmin = -zmax; zmax = tmp; }
172  xmin += fDirectTransform3D->dx();
173  xmax += fDirectTransform3D->dx();
174  ymin += fDirectTransform3D->dy();
175  ymax += fDirectTransform3D->dy();
176  zmin += fDirectTransform3D->dz();
177  zmax += fDirectTransform3D->dz();
178  }
179  else
180  {
181  // Use additional reflection in Z to set up affine transformation
182  //
183  G4Transform3D transform3D = G4ReflectZ3D()*(*fDirectTransform3D);
184  G4AffineTransform transform(transform3D.getRotation().inverse(),
185  transform3D.getTranslation());
186 
187  // Find bounding box
188  //
189  G4VoxelLimits unLimit;
190  fPtrSolid->CalculateExtent(kXAxis,unLimit,transform,xmin,xmax);
191  fPtrSolid->CalculateExtent(kYAxis,unLimit,transform,ymin,ymax);
192  fPtrSolid->CalculateExtent(kZAxis,unLimit,transform,zmin,zmax);
193  }
194 
195  pMin.set(xmin,ymin,-zmax);
196  pMax.set(xmax,ymax,-zmin);
197 
198  // Check correctness of the bounding box
199  //
200  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
201  {
202  std::ostringstream message;
203  message << "Bad bounding box (min >= max) for solid: "
204  << GetName() << " !"
205  << "\npMin = " << pMin
206  << "\npMax = " << pMax;
207  G4Exception("G4ReflectedSolid::BoundingLimits()", "GeomMgt0001",
208  JustWarning, message);
209  DumpInfo();
210  }
211 }
212 
214 //
215 // Calculate extent under transform and specified limit
216 
217 G4bool
219  const G4VoxelLimits& pVoxelLimits,
220  const G4AffineTransform& pTransform,
221  G4double& pMin,
222  G4double& pMax ) const
223 {
224  // Separation of transformations. Calculation of the extent is done
225  // in a reflection of the global space. In such way, the voxel is
226  // reflected, but the solid is transformed just by G4AffineTransform.
227  // It allows to use CalculateExtent() of the solid.
228 
229  // Reflect voxel limits in Z
230  //
231  G4VoxelLimits limits;
232  limits.AddLimit(kXAxis, pVoxelLimits.GetMinXExtent(),
233  pVoxelLimits.GetMaxXExtent());
234  limits.AddLimit(kYAxis, pVoxelLimits.GetMinYExtent(),
235  pVoxelLimits.GetMaxYExtent());
236  limits.AddLimit(kZAxis,-pVoxelLimits.GetMaxZExtent(),
237  -pVoxelLimits.GetMinZExtent());
238 
239  // Set affine transformation
240  //
241  G4Transform3D transform3D = G4ReflectZ3D()*pTransform*(*fDirectTransform3D);
242  G4AffineTransform transform(transform3D.getRotation().inverse(),
243  transform3D.getTranslation());
244 
245  // Find extent
246  //
247  if (!fPtrSolid->CalculateExtent(pAxis, limits, transform, pMin, pMax))
248  {
249  return false;
250  }
251  if (pAxis == kZAxis)
252  {
253  G4double tmp= -pMin; pMin= -pMax; pMax= tmp;
254  }
255 
256  return true;
257 }
258 
260 //
261 //
262 
264 {
265  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
266  return fPtrSolid->Inside(newPoint);
267 }
268 
270 //
271 //
272 
275 {
276  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
278  return (*fDirectTransform3D)*normal;
279 }
280 
282 //
283 // The same algorithm as in DistanceToIn(p)
284 
285 G4double
287  const G4ThreeVector& v ) const
288 {
289  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
290  G4ThreeVector newDirection = (*fDirectTransform3D)*G4Vector3D(v);
291  return fPtrSolid->DistanceToIn(newPoint,newDirection);
292 }
293 
295 //
296 // Approximate nearest distance from the point p to the intersection of
297 // two solids
298 
299 G4double
301 {
302  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
303  return fPtrSolid->DistanceToIn(newPoint);
304 }
305 
307 //
308 // The same algorithm as DistanceToOut(p)
309 
310 G4double
312  const G4ThreeVector& v,
313  const G4bool calcNorm,
314  G4bool *validNorm,
315  G4ThreeVector *n ) const
316 {
317  G4ThreeVector solNorm;
318 
319  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
320  G4ThreeVector newDirection = (*fDirectTransform3D)*G4Vector3D(v);
321 
322  G4double dist = fPtrSolid->DistanceToOut(newPoint, newDirection,
323  calcNorm, validNorm, &solNorm);
324  if(calcNorm)
325  {
326  *n = (*fDirectTransform3D)*G4Vector3D(solNorm);
327  }
328  return dist;
329 }
330 
332 //
333 // Inverted algorithm of DistanceToIn(p)
334 
335 G4double
337 {
338  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
339  return fPtrSolid->DistanceToOut(newPoint);
340 }
341 
343 //
344 //
345 
346 void
348  const G4int,
349  const G4VPhysicalVolume* )
350 {
351  DumpInfo();
352  G4Exception("G4ReflectedSolid::ComputeDimensions()",
353  "GeomMgt0001", FatalException,
354  "Method not applicable in this context!");
355 }
356 
358 //
359 // Return a point (G4ThreeVector) randomly and uniformly selected
360 // on the solid surface
361 
363 {
365  return (*fDirectTransform3D)*G4Point3D(p);
366 }
367 
369 //
370 // Make a clone of this object
371 
373 {
374  return new G4ReflectedSolid(*this);
375 }
376 
377 
379 //
380 // Stream object contents to an output stream
381 
382 std::ostream& G4ReflectedSolid::StreamInfo(std::ostream& os) const
383 {
384  os << "-----------------------------------------------------------\n"
385  << " *** Dump for Reflected solid - " << GetName() << " ***\n"
386  << " ===================================================\n"
387  << " Solid type: " << GetEntityType() << "\n"
388  << " Parameters of constituent solid: \n"
389  << "===========================================================\n";
390  fPtrSolid->StreamInfo(os);
391  os << "===========================================================\n"
392  << " Transformations: \n"
393  << " Direct transformation - translation : \n"
394  << " " << fDirectTransform3D->getTranslation() << "\n"
395  << " - rotation : \n"
396  << " ";
398  os << "\n"
399  << "===========================================================\n";
400 
401  return os;
402 }
403 
405 //
406 //
407 
408 void
410 {
411  scene.AddSolid (*this);
412 }
413 
415 //
416 //
417 
418 G4Polyhedron*
420 {
421  G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
422  if (polyhedron)
423  {
424  polyhedron->Transform(*fDirectTransform3D);
425  return polyhedron;
426  }
427  else
428  {
429  std::ostringstream message;
430  message << "Solid - " << GetName()
431  << " - original solid has no" << G4endl
432  << "corresponding polyhedron. Returning NULL!";
433  G4Exception("G4ReflectedSolid::CreatePolyhedron()",
434  "GeomMgt1001", JustWarning, message);
435  return 0;
436  }
437 }
438 
440 //
441 //
442 
445 {
446  if (!fpPolyhedron ||
450  {
452  fRebuildPolyhedron = false;
453  }
454  return fpPolyhedron;
455 }
void set(double x, double y, double z)
CLHEP::HepRotation getRotation() const
void SetDirectTransform3D(G4Transform3D &)
G4ThreeVector GetPointOnSurface() const
virtual const G4ReflectedSolid * GetReflectedSolidPtr() const
CLHEP::Hep3Vector getTranslation() const
G4double GetMinZExtent() const
Double_t xx
virtual void AddSolid(const G4Box &)=0
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
HepGeom::Transform3D G4Transform3D
G4VSolid * Clone() const
std::ostream & print(std::ostream &os) const
Definition: RotationIO.cc:21
std::ostream & StreamInfo(std::ostream &os) const
void DescribeYourselfTo(G4VGraphicsScene &scene) const
virtual void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:625
#define G4endl
Definition: G4ios.hh:61
Transform3D inverse() const
Definition: Transform3D.cc:142
const char * p
Definition: xmltok.h:285
G4double GetMinYExtent() const
void message(RunManager *runmanager)
Definition: ts_scorers.cc:72
double dx() const
Definition: Transform3D.h:279
G4Transform3D * fDirectTransform3D
virtual G4ThreeVector GetPointOnSurface() const
Definition: G4VSolid.cc:152
double z() const
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
Double_t zz
G4VSolid * GetConstituentMovedSolid() const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
Float_t tmp
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
G4Polyhedron * CreatePolyhedron() const
virtual EInside Inside(const G4ThreeVector &p) const =0
double yy() const
Definition: Transform3D.h:264
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:110
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
double zz() const
Definition: Transform3D.h:276
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
G4double GetMaxYExtent() const
double xx() const
Definition: Transform3D.h:252
HepPolyhedron & Transform(const G4Transform3D &t)
G4double GetMaxZExtent() const
void BoundingLimits(G4ThreeVector &pMin, G4ThreeVector &pMax) const
G4double GetMaxXExtent() const
virtual std::ostream & StreamInfo(std::ostream &os) const =0
double dz() const
Definition: Transform3D.h:285
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:35
static G4int GetNumberOfRotationSteps()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.hh:65
virtual ~G4ReflectedSolid()
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
int G4int
Definition: G4Types.hh:78
EInside
Definition: geomdefs.hh:58
G4ReflectedSolid & operator=(const G4ReflectedSolid &rhs)
EInside Inside(const G4ThreeVector &p) const
EAxis
Definition: geomdefs.hh:54
G4String GetName() const
void DumpInfo() const
G4Polyhedron * GetPolyhedron() const
G4Transform3D GetTransform3D() const
virtual G4GeometryType GetEntityType() const
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
double x() const
void AddLimit(const EAxis pAxis, const G4double pMin, const G4double pMax)
double dy() const
Definition: Transform3D.h:282
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
virtual G4Polyhedron * CreatePolyhedron() const
Definition: G4VSolid.cc:660
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
Char_t n[5]
G4Polyhedron * fpPolyhedron
double y() const
G4ReflectedSolid(const G4String &pName, G4VSolid *pSolid, const G4Transform3D &transform)
G4double GetMinXExtent() const
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:35
HepRotation inverse() const
HepGeom::ReflectZ3D G4ReflectZ3D
G4Transform3D GetDirectTransform3D() const