Geant4  10.03
G4ReflectedSolid.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: G4ReflectedSolid.cc 100906 2016-11-03 09:59:32Z 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 {
157  fPtrSolid->Extent(pMin,pMax);
158  G4double xmin = pMin.x(), ymin = pMin.y(), zmin = pMin.z();
159  G4double xmax = pMax.x(), ymax = pMax.y(), zmax = pMax.z();
160  G4double xx = fDirectTransform3D->xx();
161  G4double yy = fDirectTransform3D->yy();
162  G4double zz = fDirectTransform3D->zz();
163 
164  if (std::abs(xx) == 1 && std::abs(yy) == 1 && std::abs(zz) == 1)
165  {
166  // Special case of reflection in axis and pure translation
167  //
168  if (xx == -1) { G4double tmp = -xmin; xmin = -xmax; xmax = tmp; }
169  if (yy == -1) { G4double tmp = -ymin; ymin = -ymax; ymax = tmp; }
170  if (zz == -1) { G4double tmp = -zmin; zmin = -zmax; zmax = tmp; }
171  xmin += fDirectTransform3D->dx();
172  xmax += fDirectTransform3D->dx();
173  ymin += fDirectTransform3D->dy();
174  ymax += fDirectTransform3D->dy();
175  zmin += fDirectTransform3D->dz();
176  zmax += fDirectTransform3D->dz();
177  }
178  else
179  {
180  // Use additional reflection in Z to set up affine transformation
181  //
182  G4Transform3D transform3D = G4ReflectZ3D()*(*fDirectTransform3D);
183  G4AffineTransform transform(transform3D.getRotation().inverse(),
184  transform3D.getTranslation());
185 
186  // Find bounding box
187  //
188  G4VoxelLimits unLimit;
189  fPtrSolid->CalculateExtent(kXAxis,unLimit,transform,xmin,xmax);
190  fPtrSolid->CalculateExtent(kYAxis,unLimit,transform,ymin,ymax);
191  fPtrSolid->CalculateExtent(kZAxis,unLimit,transform,zmin,zmax);
192  }
193 
194  pMin.set(xmin,ymin,-zmax);
195  pMax.set(xmax,ymax,-zmin);
196 
197  // Check correctness of the bounding box
198  //
199  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
200  {
201  std::ostringstream message;
202  message << "Bad bounding box (min >= max) for solid: "
203  << GetName() << " !"
204  << "\npMin = " << pMin
205  << "\npMax = " << pMax;
206  G4Exception("G4ReflectedSolid::Extent()", "GeomMgt0001",
207  JustWarning, message);
208  DumpInfo();
209  }
210 }
211 
213 //
214 // Calculate extent under transform and specified limit
215 
216 G4bool
218  const G4VoxelLimits& pVoxelLimits,
219  const G4AffineTransform& pTransform,
220  G4double& pMin,
221  G4double& pMax ) const
222 {
223  // Separation of transformations. Calculation of the extent is done
224  // in a reflection of the global space. In such way, the voxel is
225  // reflected, but the solid is transformed just by G4AffineTransform.
226  // It allows to use CalculateExtent() of the solid.
227 
228  // Reflect voxel limits in Z
229  //
230  G4VoxelLimits limits;
231  limits.AddLimit(kXAxis, pVoxelLimits.GetMinXExtent(),
232  pVoxelLimits.GetMaxXExtent());
233  limits.AddLimit(kYAxis, pVoxelLimits.GetMinYExtent(),
234  pVoxelLimits.GetMaxYExtent());
235  limits.AddLimit(kZAxis,-pVoxelLimits.GetMaxZExtent(),
236  -pVoxelLimits.GetMinZExtent());
237 
238  // Set affine transformation
239  //
240  G4Transform3D transform3D = G4ReflectZ3D()*pTransform*(*fDirectTransform3D);
241  G4AffineTransform transform(transform3D.getRotation().inverse(),
242  transform3D.getTranslation());
243 
244  // Find extent
245  //
246  if (!fPtrSolid->CalculateExtent(pAxis, limits, transform, pMin, pMax))
247  {
248  return false;
249  }
250  if (pAxis == kZAxis)
251  {
252  G4double tmp= -pMin; pMin= -pMax; pMax= tmp;
253  }
254 
255  return true;
256 }
257 
259 //
260 //
261 
263 {
264  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
265  return fPtrSolid->Inside(newPoint);
266 }
267 
269 //
270 //
271 
274 {
275  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
277  return (*fDirectTransform3D)*normal;
278 }
279 
281 //
282 // The same algorithm as in DistanceToIn(p)
283 
284 G4double
286  const G4ThreeVector& v ) const
287 {
288  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
289  G4ThreeVector newDirection = (*fDirectTransform3D)*G4Vector3D(v);
290  return fPtrSolid->DistanceToIn(newPoint,newDirection);
291 }
292 
294 //
295 // Approximate nearest distance from the point p to the intersection of
296 // two solids
297 
298 G4double
300 {
301  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
302  return fPtrSolid->DistanceToIn(newPoint);
303 }
304 
306 //
307 // The same algorithm as DistanceToOut(p)
308 
309 G4double
311  const G4ThreeVector& v,
312  const G4bool calcNorm,
313  G4bool *validNorm,
314  G4ThreeVector *n ) const
315 {
316  G4ThreeVector solNorm;
317 
318  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
319  G4ThreeVector newDirection = (*fDirectTransform3D)*G4Vector3D(v);
320 
321  G4double dist = fPtrSolid->DistanceToOut(newPoint, newDirection,
322  calcNorm, validNorm, &solNorm);
323  if(calcNorm)
324  {
325  *n = (*fDirectTransform3D)*G4Vector3D(solNorm);
326  }
327  return dist;
328 }
329 
331 //
332 // Inverted algorithm of DistanceToIn(p)
333 
334 G4double
336 {
337  G4ThreeVector newPoint = (*fDirectTransform3D)*G4Point3D(p);
338  return fPtrSolid->DistanceToOut(newPoint);
339 }
340 
342 //
343 //
344 
345 void
347  const G4int,
348  const G4VPhysicalVolume* )
349 {
350  DumpInfo();
351  G4Exception("G4ReflectedSolid::ComputeDimensions()",
352  "GeomMgt0001", FatalException,
353  "Method not applicable in this context!");
354 }
355 
357 //
358 // Return a point (G4ThreeVector) randomly and uniformly selected
359 // on the solid surface
360 
362 {
364  return (*fDirectTransform3D)*G4Point3D(p);
365 }
366 
368 //
369 // Make a clone of this object
370 
372 {
373  return new G4ReflectedSolid(*this);
374 }
375 
376 
378 //
379 // Stream object contents to an output stream
380 
381 std::ostream& G4ReflectedSolid::StreamInfo(std::ostream& os) const
382 {
383  os << "-----------------------------------------------------------\n"
384  << " *** Dump for Reflected solid - " << GetName() << " ***\n"
385  << " ===================================================\n"
386  << " Solid type: " << GetEntityType() << "\n"
387  << " Parameters of constituent solid: \n"
388  << "===========================================================\n";
389  fPtrSolid->StreamInfo(os);
390  os << "===========================================================\n"
391  << " Transformations: \n"
392  << " Direct transformation - translation : \n"
393  << " " << fDirectTransform3D->getTranslation() << "\n"
394  << " - rotation : \n"
395  << " ";
396  fDirectTransform3D->getRotation().print(os);
397  os << "\n"
398  << "===========================================================\n";
399 
400  return os;
401 }
402 
404 //
405 //
406 
407 void
409 {
410  scene.AddSolid (*this);
411 }
412 
414 //
415 //
416 
417 G4Polyhedron*
419 {
420  G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
421  if (polyhedron)
422  {
423  polyhedron->Transform(*fDirectTransform3D);
424  return polyhedron;
425  }
426  else
427  {
428  std::ostringstream message;
429  message << "Solid - " << GetName()
430  << " - original solid has no" << G4endl
431  << "corresponding polyhedron. Returning NULL!";
432  G4Exception("G4ReflectedSolid::CreatePolyhedron()",
433  "GeomMgt1001", JustWarning, message);
434  return 0;
435  }
436 }
437 
439 //
440 //
441 
444 {
445  if (!fpPolyhedron ||
448  fpPolyhedron->GetNumberOfRotationSteps())
449  {
451  fRebuildPolyhedron = false;
452  }
453  return fpPolyhedron;
454 }
G4String GetName() const
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
G4double GetMinYExtent() const
void SetDirectTransform3D(G4Transform3D &)
CLHEP::Hep3Vector G4ThreeVector
std::ostream & StreamInfo(std::ostream &os) const
G4Polyhedron * GetPolyhedron() const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
EInside Inside(const G4ThreeVector &p) const
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:35
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:35
G4ReflectedSolid & operator=(const G4ReflectedSolid &rhs)
virtual void AddSolid(const G4Box &)=0
int G4int
Definition: G4Types.hh:78
G4Transform3D * fDirectTransform3D
void DumpInfo() const
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4Transform3D GetTransform3D() const
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
G4double GetMaxXExtent() const
void AddLimit(const EAxis pAxis, const G4double pMin, const G4double pMax)
G4double GetMinZExtent() const
virtual std::ostream & StreamInfo(std::ostream &os) const =0
virtual ~G4ReflectedSolid()
virtual EInside Inside(const G4ThreeVector &p) const =0
bool G4bool
Definition: G4Types.hh:79
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const =0
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
HepGeom::Transform3D G4Transform3D
virtual G4Polyhedron * CreatePolyhedron() const
Definition: G4VSolid.cc:660
const G4int n
virtual G4GeometryType GetEntityType() const
G4VSolid * GetConstituentMovedSolid() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
HepGeom::ReflectZ3D G4ReflectZ3D
G4double GetMinXExtent() const
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
G4double GetMaxZExtent() const
EInside
Definition: geomdefs.hh:58
EAxis
Definition: geomdefs.hh:54
virtual G4ThreeVector GetPointOnSurface() const
Definition: G4VSolid.cc:153
G4Polyhedron * fpPolyhedron
G4ReflectedSolid(const G4String &pName, G4VSolid *pSolid, const G4Transform3D &transform)
#define G4endl
Definition: G4ios.hh:61
G4double GetMaxYExtent() const
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:111
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
double G4double
Definition: G4Types.hh:76
G4Polyhedron * CreatePolyhedron() const
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
G4ThreeVector GetPointOnSurface() const
virtual void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:626
G4Transform3D GetDirectTransform3D() const
virtual const G4ReflectedSolid * GetReflectedSolidPtr() const
G4VSolid * Clone() const