Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4ScaledSolid.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:$
28 //
29 // Implementation for G4ScaledSolid class
30 //
31 // History:
32 //
33 // 27.10.15 G.Cosmo: created, based on implementation also provided in Root
34 //
35 // --------------------------------------------------------------------
36 
37 #include "G4ScaledSolid.hh"
38 #include "G4BoundingEnvelope.hh"
39 
40 #include "G4VPVParameterisation.hh"
41 
42 #include "G4ScaleTransform.hh"
43 
44 #include "G4VGraphicsScene.hh"
45 #include "G4Polyhedron.hh"
46 
48 //
49 // Constructor
50 //
52  G4VSolid* pSolid ,
53  const G4Scale3D& pScale )
54  : G4VSolid(pName), fPtrSolid(pSolid),
55  fRebuildPolyhedron(false), fpPolyhedron(0)
56 {
57  fScale = new G4ScaleTransform(pScale);
58 }
59 
61 //
62 // Fake default constructor - sets only member data and allocates memory
63 // for usage restricted to object persistency.
64 //
66  : G4VSolid(a), fPtrSolid(0), fScale(0),
67  fRebuildPolyhedron(false), fpPolyhedron(0)
68 {
69 }
70 
72 //
73 // Destructor
74 //
76 {
77  delete fpPolyhedron; fpPolyhedron= 0;
78  delete fScale; fScale= 0;
79 }
80 
82 //
83 // Copy constructor
84 //
86  : G4VSolid (rhs), fPtrSolid(rhs.fPtrSolid),
87  fRebuildPolyhedron(false), fpPolyhedron(0)
88 {
89  fScale = new G4ScaleTransform(*(rhs.fScale));
90 }
91 
93 //
94 // Assignment operator
95 //
97 {
98  // Check assignment to self
99  //
100  if (this == &rhs) { return *this; }
101 
102  // Copy base class data
103  //
104  G4VSolid::operator=(rhs);
105 
106  // Copy data
107  //
108  fPtrSolid = rhs.fPtrSolid;
109  delete fScale;
110  fScale = new G4ScaleTransform(*(rhs.fScale));
111  fRebuildPolyhedron = false;
112  delete fpPolyhedron; fpPolyhedron= 0;
113 
114  return *this;
115 }
116 
118 //
119 // Return original solid not scaled
120 //
122 {
123  return fPtrSolid;
124 }
125 
127 //
128 // Get bounding box
129 
131 {
132  G4ThreeVector bmin,bmax;
133  G4ThreeVector scale = fScale->GetScale();
134 
135  fPtrSolid->Extent(bmin,bmax);
136  pMin.set(bmin.x()*scale.x(),bmin.y()*scale.y(),bmin.z()*scale.z());
137  pMax.set(bmax.x()*scale.x(),bmax.y()*scale.y(),bmax.z()*scale.z());
138 
139  // Check correctness of the bounding box
140  //
141  if (pMin.x() >= pMax.x() || pMin.y() >= pMax.y() || pMin.z() >= pMax.z())
142  {
143  std::ostringstream message;
144  message << "Bad bounding box (min >= max) for solid: "
145  << GetName() << " !"
146  << "\npMin = " << pMin
147  << "\npMax = " << pMax;
148  G4Exception("G4ScaledSolid::Extent()", "GeomMgt0001",
149  JustWarning, message);
150  DumpInfo();
151  }
152 }
153 
155 //
156 // Calculate extent under transform and specified limit
157 //
158 G4bool
160  const G4VoxelLimits& pVoxelLimit,
161  const G4AffineTransform& pTransform,
162  G4double& pMin,
163  G4double& pMax ) const
164 {
165  // Find bounding box of unscaled solid
166  G4ThreeVector bmin,bmax;
167  fPtrSolid->Extent(bmin,bmax);
168 
169  // Set combined transformation
170  G4Transform3D transform3D =
171  G4Transform3D(pTransform.NetRotation().inverse(),
172  pTransform.NetTranslation())*GetScaleTransform();
173 
174  // Find extent
175  G4BoundingEnvelope bbox(bmin,bmax);
176  return bbox.CalculateExtent(pAxis,pVoxelLimit,transform3D,pMin,pMax);
177 }
178 
180 //
181 // Inside
182 //
184 {
185  return fPtrSolid->Inside(fScale->Transform(p));
186 }
187 
189 //
190 // SurfaceNormal
191 //
194 {
195  // Transform point to unscaled shape frame
196  G4ThreeVector newPoint;
197  fScale->Transform(p, newPoint);
198 
199  // Compute normal in unscaled frame
200  G4ThreeVector newNormal = fPtrSolid->SurfaceNormal(newPoint);
202 
203  // Convert normal to scaled frame
204  fScale->InverseTransformNormal(newNormal, normal);
205  return normal/normal.mag();
206 }
207 
209 //
210 // The same algorithm as in DistanceToIn(p)
211 //
212 G4double
214  const G4ThreeVector& v ) const
215 {
216  // Transform point and direction to unscaled shape frame
217  G4ThreeVector newPoint;
218  fScale->Transform(p, newPoint);
219 
220  // Direction is un-normalized after scale transformation
221  G4ThreeVector newDirection;
222  fScale->Transform(v, newDirection);
223  newDirection = newDirection/newDirection.mag();
224 
225  // Compute distance in unscaled system
226  G4double dist = fPtrSolid->DistanceToIn(newPoint,newDirection);
227 
228  // Return converted distance to global
229  return fScale->InverseTransformDistance(dist, newDirection);
230 }
231 
233 //
234 // Approximate nearest distance from the point p to the solid from outside
235 //
236 G4double
238 {
239  // Transform point to unscaled shape frame
240  G4ThreeVector newPoint;
241  fScale->Transform(p, newPoint);
242 
243  // Compute unscaled safety, then scale it.
244  G4double dist = fPtrSolid->DistanceToIn(newPoint);
245  return fScale->InverseTransformDistance(dist);
246 }
247 
249 //
250 // The same algorithm as DistanceToOut(p)
251 //
252 G4double
254  const G4ThreeVector& v,
255  const G4bool calcNorm,
256  G4bool *validNorm,
257  G4ThreeVector *n ) const
258 {
259  // Transform point and direction to unscaled shape frame
260  G4ThreeVector newPoint;
261  fScale->Transform(p, newPoint);
262 
263  // Direction is un-normalized after scale transformation
264  G4ThreeVector newDirection;
265  fScale->Transform(v, newDirection);
266  newDirection = newDirection/newDirection.mag();
267 
268  // Compute distance in unscaled system
269  G4ThreeVector solNorm;
270  G4double dist = fPtrSolid->DistanceToOut(newPoint,newDirection,
271  calcNorm,validNorm,&solNorm);
272  if(calcNorm)
273  {
275  fScale->TransformNormal(solNorm, normal);
276  *n = normal/normal.mag();
277  }
278 
279  // Return distance converted to global
280  return fScale->InverseTransformDistance(dist, newDirection);
281 }
282 
284 //
285 // Approximate nearest distance from the point p to the solid from inside
286 //
287 G4double
289 {
290  // Transform point to unscaled shape frame
291  G4ThreeVector newPoint;
292  fScale->Transform(p, newPoint);
293 
294  // Compute unscaled safety, then scale it.
295  G4double dist = fPtrSolid->DistanceToOut(newPoint);
296  return fScale->InverseTransformDistance(dist);
297 }
298 
300 //
301 // ComputeDimensions
302 //
303 void
305  const G4int,
306  const G4VPhysicalVolume* )
307 {
308  DumpInfo();
309  G4Exception("G4ScaledSolid::ComputeDimensions()",
310  "GeomSolids0001", FatalException,
311  "Method not applicable in this context!");
312 }
313 
315 //
316 // Returns a point (G4ThreeVector) randomly and uniformly selected
317 // on the solid surface
318 //
320 {
321  return fScale->InverseTransform(fPtrSolid->GetPointOnSurface());
322 }
323 
325 //
326 // Return object type name
327 //
329 {
330  return G4String("G4ScaledSolid");
331 }
332 
334 //
335 // Make a clone of the object
336 //
338 {
339  return new G4ScaledSolid(*this);
340 }
341 
343 //
344 // Returning the scaling transformation
345 //
347 {
348  return G4Scale3D(fScale->GetScale().x(),
349  fScale->GetScale().y(),
350  fScale->GetScale().z());
351 }
352 
354 //
355 // Setting the scaling transformation
356 //
358 {
359  if (fScale) { delete fScale; }
360  fScale = new G4ScaleTransform(scale);
361  fRebuildPolyhedron = true;
362 }
363 
365 //
366 // Stream object contents to an output stream
367 //
368 std::ostream& G4ScaledSolid::StreamInfo(std::ostream& os) const
369 {
370  os << "-----------------------------------------------------------\n"
371  << " *** Dump for Scaled solid - " << GetName() << " ***\n"
372  << " ===================================================\n"
373  << " Solid type: " << GetEntityType() << "\n"
374  << " Parameters of constituent solid: \n"
375  << "===========================================================\n";
376  fPtrSolid->StreamInfo(os);
377  os << "===========================================================\n"
378  << " Scaling: \n"
379  << " Scale transformation : \n"
380  << " " << fScale->GetScale().x() << ", "
381  << fScale->GetScale().y() << ", "
382  << fScale->GetScale().z() << "\n"
383  << "===========================================================\n";
384 
385  return os;
386 }
387 
389 //
390 // DescribeYourselfTo
391 //
392 void
394 {
395  scene.AddSolid (*this);
396 }
397 
399 //
400 // CreatePolyhedron
401 //
402 G4Polyhedron*
404 {
405  G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
406  polyhedron->Transform(GetScaleTransform());
407  return polyhedron;
408 }
409 
411 //
412 // GetPolyhedron
413 //
415 {
416  if (!fpPolyhedron ||
417  fRebuildPolyhedron ||
419  fpPolyhedron->GetNumberOfRotationSteps())
420  {
421  fpPolyhedron = CreatePolyhedron();
422  fRebuildPolyhedron = false;
423  }
424  return fpPolyhedron;
425 }
void set(double x, double y, double z)
G4String GetName() const
G4ScaledSolid & operator=(const G4ScaledSolid &rhs)
G4ThreeVector GetPointOnSurface() const
G4VSolid * Clone() const
G4double InverseTransformDistance(G4double dist, const G4ThreeVector &dir) const
void TransformNormal(const G4ThreeVector &global, G4ThreeVector &local) const
G4Scale3D GetScaleTransform() const
G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
double x() const
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
const char * p
Definition: xmltok.h:285
const G4ThreeVector & GetScale() const
void InverseTransformNormal(const G4ThreeVector &local, G4ThreeVector &global) const
G4ThreeVector NetTranslation() const
void Transform(const G4ThreeVector &global, G4ThreeVector &local) const
HepPolyhedron & Transform(const G4Transform3D &t)
virtual void AddSolid(const G4Box &)=0
int G4int
Definition: G4Types.hh:78
double z() const
HepRotation inverse() const
void DumpInfo() const
virtual ~G4ScaledSolid()
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
void InverseTransform(const G4ThreeVector &local, G4ThreeVector &global) const
G4Polyhedron * CreatePolyhedron() const
virtual std::ostream & StreamInfo(std::ostream &os) const =0
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 &pVoxelLimits, const G4Transform3D &pTransform3D, G4double &pMin, G4double &pMax) const
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
HepGeom::Transform3D G4Transform3D
HepGeom::Scale3D G4Scale3D
virtual G4Polyhedron * CreatePolyhedron() const
Definition: G4VSolid.cc:660
const G4int n
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
tuple v
Definition: test.py:18
EInside Inside(const G4ThreeVector &p) const
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4Polyhedron * GetPolyhedron() const
void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4RotationMatrix NetRotation() const
static G4int GetNumberOfRotationSteps()
EInside
Definition: geomdefs.hh:58
EAxis
Definition: geomdefs.hh:54
virtual G4ThreeVector GetPointOnSurface() const
Definition: G4VSolid.cc:153
void SetScaleTransform(const G4Scale3D &scale)
double y() const
G4VSolid * GetUnscaledSolid() const
std::ostream & StreamInfo(std::ostream &os) 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
G4ScaledSolid(const G4String &pName, G4VSolid *pSolid, const G4Scale3D &pScale)
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const =0
double mag() const
G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
virtual void Extent(G4ThreeVector &pMin, G4ThreeVector &pMax) const
Definition: G4VSolid.cc:626
G4GeometryType GetEntityType() const