Geant4  10.01.p01
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 66356 2012-12-18 09:02:32Z gcosmo $
28 //
29 //
30 // Implementation for G4ReflectedSolid class for boolean
31 // operations between other solids
32 //
33 // Author: Vladimir Grichine, 23.07.01 (Vladimir.Grichine@cern.ch)
34 //
35 // --------------------------------------------------------------------
36 
37 #include "G4ReflectedSolid.hh"
38 
39 #include <sstream>
40 
41 #include "G4Point3D.hh"
42 #include "G4Normal3D.hh"
43 
44 #include "G4VoxelLimits.hh"
45 
46 #include "G4VPVParameterisation.hh"
47 
48 #include "G4VGraphicsScene.hh"
49 #include "G4Polyhedron.hh"
50 
51 
53 //
54 // Constructor using HepTransform3D, in fact HepReflect3D
55 
57  G4VSolid* pSolid ,
58  const G4Transform3D& transform )
59  : G4VSolid(pName), fpPolyhedron(0)
60 {
61  fPtrSolid = pSolid ;
62  G4RotationMatrix rotMatrix ;
63 
65  new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
67  new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
68  fPtrTransform->Invert() ;
69 
70  fDirectTransform3D = new G4Transform3D(transform) ;
71  fPtrTransform3D = new G4Transform3D(transform.inverse()) ;
72 }
73 
75 //
76 
78 {
79  if(fPtrTransform)
80  {
81  delete fPtrTransform; fPtrTransform=0;
83  }
84  if(fPtrTransform3D)
85  {
88  }
89  delete fpPolyhedron;
90 }
91 
93 //
94 
96  : G4VSolid(rhs), fPtrSolid(rhs.fPtrSolid), fpPolyhedron(0)
97 {
102 }
103 
105 //
106 
108 {
109  // Check assignment to self
110  //
111  if (this == &rhs) { return *this; }
112 
113  // Copy base class data
114  //
115  G4VSolid::operator=(rhs);
116 
117  // Copy data
118  //
120  delete fPtrTransform;
122  delete fDirectTransform;
124  delete fPtrTransform3D;
126  delete fDirectTransform3D;
128 
129  return *this;
130 }
131 
133 //
134 
136 {
137  return G4String("G4ReflectedSolid");
138 }
139 
141 {
142  return this;
143 }
144 
146 {
147  return this;
148 }
149 
151 {
152  return fPtrSolid;
153 }
154 
156 
158 {
159  G4AffineTransform aTransform = *fPtrTransform;
160  return aTransform;
161 }
162 
164 {
165  fPtrTransform = &transform ;
166  fpPolyhedron = 0;
167 }
168 
170 
172 {
173  G4AffineTransform aTransform= *fDirectTransform;
174  return aTransform;
175 }
176 
178 {
179  fDirectTransform = &transform ;
180  fpPolyhedron = 0;
181 }
182 
184 
186 {
187  G4Transform3D aTransform = *fPtrTransform3D;
188  return aTransform;
189 }
190 
192 {
193  fPtrTransform3D = &transform ;
194  fpPolyhedron = 0;
195 }
196 
198 
200 {
201  G4Transform3D aTransform= *fDirectTransform3D;
202  return aTransform;
203 }
204 
206 {
207  fDirectTransform3D = &transform ;
208  fpPolyhedron = 0;
209 }
210 
212 
214 {
216  return InvRotation;
217 }
218 
220 {
222 }
223 
225 
227 {
228  return fPtrTransform->NetTranslation();
229 }
230 
232 {
234 }
235 
237 
239 {
241  return Rotation;
242 }
243 
245 {
246  fPtrTransform->SetNetRotation(matrix);
247 }
248 
250 
252 {
254 }
255 
257 {
259 }
260 
262 //
263 //
264 
265 G4bool
267  const G4VoxelLimits& pVoxelLimit,
268  const G4AffineTransform& pTransform,
269  G4double& pMin,
270  G4double& pMax ) const
271 {
272 
273  G4VoxelLimits unLimit;
274  G4AffineTransform unTransform;
275 
276  G4double x1 = -kInfinity, x2 = kInfinity,
277  y1 = -kInfinity, y2 = kInfinity,
278  z1 = -kInfinity, z2 = kInfinity;
279 
280  G4bool existsAfterClip = false ;
281  existsAfterClip =
282  fPtrSolid->CalculateExtent(kXAxis,unLimit,unTransform,x1,x2);
283  existsAfterClip =
284  fPtrSolid->CalculateExtent(kYAxis,unLimit,unTransform,y1,y2);
285  existsAfterClip =
286  fPtrSolid->CalculateExtent(kZAxis,unLimit,unTransform,z1,z2);
287 
288  existsAfterClip = false;
289  pMin = +kInfinity ;
290  pMax = -kInfinity ;
291 
292  G4Transform3D pTransform3D = G4Transform3D(pTransform.NetRotation().inverse(),
293  pTransform.NetTranslation());
294 
295  G4Transform3D transform3D = pTransform3D*(*fDirectTransform3D);
296 
297  G4Point3D tmpPoint;
298 
299  // Calculate rotated vertex coordinates
300 
301  G4ThreeVectorList* vertices = new G4ThreeVectorList();
302 
303  if (vertices)
304  {
305  vertices->reserve(8);
306 
307  G4ThreeVector vertex0(x1,y1,z1) ;
308  tmpPoint = transform3D*G4Point3D(vertex0);
309  vertex0 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
310  vertices->push_back(vertex0);
311 
312  G4ThreeVector vertex1(x2,y1,z1) ;
313  tmpPoint = transform3D*G4Point3D(vertex1);
314  vertex1 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
315  vertices->push_back(vertex1);
316 
317  G4ThreeVector vertex2(x2,y2,z1) ;
318  tmpPoint = transform3D*G4Point3D(vertex2);
319  vertex2 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
320  vertices->push_back(vertex2);
321 
322  G4ThreeVector vertex3(x1,y2,z1) ;
323  tmpPoint = transform3D*G4Point3D(vertex3);
324  vertex3 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
325  vertices->push_back(vertex3);
326 
327  G4ThreeVector vertex4(x1,y1,z2) ;
328  tmpPoint = transform3D*G4Point3D(vertex4);
329  vertex4 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
330  vertices->push_back(vertex4);
331 
332  G4ThreeVector vertex5(x2,y1,z2) ;
333  tmpPoint = transform3D*G4Point3D(vertex5);
334  vertex5 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
335  vertices->push_back(vertex5);
336 
337  G4ThreeVector vertex6(x2,y2,z2) ;
338  tmpPoint = transform3D*G4Point3D(vertex6);
339  vertex6 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
340  vertices->push_back(vertex6);
341 
342  G4ThreeVector vertex7(x1,y2,z2) ;
343  tmpPoint = transform3D*G4Point3D(vertex7);
344  vertex7 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
345  vertices->push_back(vertex7);
346  }
347  else
348  {
349  DumpInfo();
350  G4Exception("G4ReflectedSolid::CalculateExtent()",
351  "GeomMgt0003", FatalException,
352  "Error in allocation of vertices. Out of memory !");
353  }
354 
355  ClipCrossSection(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
356  ClipCrossSection(vertices,4,pVoxelLimit,pAxis,pMin,pMax) ;
357  ClipBetweenSections(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
358 
359  if (pVoxelLimit.IsLimited(pAxis) == false)
360  {
361  if ( pMin != kInfinity || pMax != -kInfinity )
362  {
363  existsAfterClip = true ;
364 
365  // Add 2*tolerance to avoid precision troubles
366 
367  pMin -= kCarTolerance;
368  pMax += kCarTolerance;
369  }
370  }
371  else
372  {
373  G4ThreeVector clipCentre(
374  ( pVoxelLimit.GetMinXExtent()+pVoxelLimit.GetMaxXExtent())*0.5,
375  ( pVoxelLimit.GetMinYExtent()+pVoxelLimit.GetMaxYExtent())*0.5,
376  ( pVoxelLimit.GetMinZExtent()+pVoxelLimit.GetMaxZExtent())*0.5);
377 
378  if ( pMin != kInfinity || pMax != -kInfinity )
379  {
380  existsAfterClip = true ;
381 
382 
383  // Check to see if endpoints are in the solid
384 
385  clipCentre(pAxis) = pVoxelLimit.GetMinExtent(pAxis);
386 
387  if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
388  {
389  pMin = pVoxelLimit.GetMinExtent(pAxis);
390  }
391  else
392  {
393  pMin -= kCarTolerance;
394  }
395  clipCentre(pAxis) = pVoxelLimit.GetMaxExtent(pAxis);
396 
397  if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
398  {
399  pMax = pVoxelLimit.GetMaxExtent(pAxis);
400  }
401  else
402  {
403  pMax += kCarTolerance;
404  }
405  }
406  // Check for case where completely enveloping clipping volume
407  // If point inside then we are confident that the solid completely
408  // envelopes the clipping volume. Hence set min/max extents according
409  // to clipping volume extents along the specified axis.
410 
411  else if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
412  {
413  existsAfterClip = true ;
414  pMin = pVoxelLimit.GetMinExtent(pAxis) ;
415  pMax = pVoxelLimit.GetMaxExtent(pAxis) ;
416  }
417  }
418  delete vertices;
419  return existsAfterClip;
420 }
421 
423 //
424 //
425 
427 {
428 
429  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
430  // G4Point3D newPoint = (*fPtrTransform3D)*G4Point3D(p) ;
431 
432  return fPtrSolid->Inside(G4ThreeVector(newPoint.x(),
433  newPoint.y(),
434  newPoint.z())) ;
435 }
436 
438 //
439 //
440 
443 {
444  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
446  fPtrSolid->SurfaceNormal(G4ThreeVector(newPoint.x(),
447  newPoint.y(),
448  newPoint.z() ) ) ;
449  G4Point3D newN = (*fDirectTransform3D)*G4Point3D(normal) ;
450  newN.unit() ;
451 
452  return G4ThreeVector(newN.x(),newN.y(),newN.z()) ;
453 }
454 
456 //
457 // The same algorithm as in DistanceToIn(p)
458 
459 G4double
461  const G4ThreeVector& v ) const
462 {
463  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
464  G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v) ;
465  newDirection.unit() ;
466  return fPtrSolid->DistanceToIn(
467  G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
468  G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z())) ;
469 }
470 
472 //
473 // Approximate nearest distance from the point p to the intersection of
474 // two solids
475 
476 G4double
478 {
479  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
480  return fPtrSolid->DistanceToIn(
481  G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z())) ;
482 }
483 
485 //
486 // The same algorithm as DistanceToOut(p)
487 
488 G4double
490  const G4ThreeVector& v,
491  const G4bool calcNorm,
492  G4bool *validNorm,
493  G4ThreeVector *n ) const
494 {
495  G4ThreeVector solNorm ;
496 
497  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
498  G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v);
499  newDirection.unit() ;
500 
501  G4double dist =
503  G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
504  G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z()),
505  calcNorm, validNorm, &solNorm) ;
506  if(calcNorm)
507  {
508  G4Point3D newN = (*fDirectTransform3D)*G4Point3D(solNorm);
509  newN.unit() ;
510  *n = G4ThreeVector(newN.x(),newN.y(),newN.z());
511  }
512  return dist ;
513 }
514 
516 //
517 // Inverted algorithm of DistanceToIn(p)
518 
519 G4double
521 {
522  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
523  return fPtrSolid->DistanceToOut(
524  G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()));
525 }
526 
528 //
529 //
530 
531 void
533  const G4int,
534  const G4VPhysicalVolume* )
535 {
536  DumpInfo();
537  G4Exception("G4ReflectedSolid::ComputeDimensions()",
538  "GeomMgt0001", FatalException,
539  "Method not applicable in this context!");
540 }
541 
543 //
544 // Return a point (G4ThreeVector) randomly and uniformly selected
545 // on the solid surface
546 
548 {
550  G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
551 
552  return G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z());
553 }
554 
556 //
557 // Make a clone of this object
558 
560 {
561  return new G4ReflectedSolid(*this);
562 }
563 
564 
566 //
567 // Stream object contents to an output stream
568 
569 std::ostream& G4ReflectedSolid::StreamInfo(std::ostream& os) const
570 {
571  os << "-----------------------------------------------------------\n"
572  << " *** Dump for Reflected solid - " << GetName() << " ***\n"
573  << " ===================================================\n"
574  << " Solid type: " << GetEntityType() << "\n"
575  << " Parameters of constituent solid: \n"
576  << "===========================================================\n";
577  fPtrSolid->StreamInfo(os);
578  os << "===========================================================\n"
579  << " Transformations: \n"
580  << " Direct transformation - translation : \n"
581  << " " << fDirectTransform->NetTranslation() << "\n"
582  << " - rotation : \n"
583  << " ";
584  fDirectTransform->NetRotation().print(os);
585  os << "\n"
586  << "===========================================================\n";
587 
588  return os;
589 }
590 
592 //
593 //
594 
595 void
597 {
598  scene.AddSolid (*this);
599 }
600 
602 //
603 //
604 
605 G4Polyhedron*
607 {
608  G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
609  if (polyhedron)
610  {
611  polyhedron->Transform(*fDirectTransform3D);
612  return polyhedron;
613  }
614  else
615  {
616  std::ostringstream message;
617  message << "Solid - " << GetName()
618  << " - original solid has no" << G4endl
619  << "corresponding polyhedron. Returning NULL!";
620  G4Exception("G4ReflectedSolid::CreatePolyhedron()",
621  "GeomMgt1001", JustWarning, message);
622  return 0;
623  }
624 }
625 
627 //
628 //
629 
632 {
633  if (!fpPolyhedron ||
635  fpPolyhedron->GetNumberOfRotationSteps())
636  {
637  delete fpPolyhedron;
639  }
640  return fpPolyhedron;
641 }
G4String GetName() const
void ClipCrossSection(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
Definition: G4VSolid.cc:347
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pMin, G4double &pMax) const =0
G4ThreeVector GetObjectTranslation() const
G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
void SetTransform(G4AffineTransform &)
static const G4double kInfinity
Definition: geomdefs.hh:42
G4double GetMinYExtent() const
void SetDirectTransform3D(G4Transform3D &)
G4RotationMatrix GetObjectRotation() const
CLHEP::Hep3Vector G4ThreeVector
std::ostream & StreamInfo(std::ostream &os) const
CLHEP::HepRotation G4RotationMatrix
G4RotationMatrix GetFrameRotation() const
G4Polyhedron * GetPolyhedron() const
G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
EInside Inside(const G4ThreeVector &p) const
G4ThreeVector NetTranslation() const
G4AffineTransform GetDirectTransform() const
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:35
void SetFrameTranslation(const G4ThreeVector &)
G4ReflectedSolid & operator=(const G4ReflectedSolid &rhs)
virtual void AddSolid(const G4Box &)=0
int G4int
Definition: G4Types.hh:78
G4Transform3D * fDirectTransform3D
void DumpInfo() const
void SetObjectTranslation(const G4ThreeVector &)
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
G4AffineTransform * fDirectTransform
void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4Transform3D GetTransform3D() const
void ComputeDimensions(G4VPVParameterisation *p, const G4int n, const G4VPhysicalVolume *pRep)
G4double GetMaxXExtent() const
G4AffineTransform & Invert()
G4double GetMinZExtent() const
virtual std::ostream & StreamInfo(std::ostream &os) const =0
void SetDirectTransform(G4AffineTransform &)
G4bool IsLimited() const
void SetNetTranslation(const G4ThreeVector &tlate)
virtual ~G4ReflectedSolid()
void SetNetRotation(const G4RotationMatrix &rot)
virtual EInside Inside(const G4ThreeVector &p) const =0
void SetTransform3D(G4Transform3D &)
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
G4Transform3D * fPtrTransform3D
void SetFrameRotation(const G4RotationMatrix &)
std::vector< G4ThreeVector > G4ThreeVectorList
Definition: G4VSolid.hh:79
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const =0
HepGeom::Transform3D G4Transform3D
virtual G4Polyhedron * CreatePolyhedron() const
Definition: G4VSolid.cc:639
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
G4RotationMatrix NetRotation() const
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:152
G4AffineTransform * fPtrTransform
void SetObjectRotation(const G4RotationMatrix &)
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:110
G4double kCarTolerance
Definition: G4VSolid.hh:305
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
G4AffineTransform GetTransform() 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
G4double GetMaxExtent(const EAxis pAxis) const
G4ThreeVector GetPointOnSurface() const
G4ThreeVector GetFrameTranslation() const
G4Transform3D GetDirectTransform3D() const
void ClipBetweenSections(G4ThreeVectorList *pVertices, const G4int pSectionIndex, const G4VoxelLimits &pVoxelLimit, const EAxis pAxis, G4double &pMin, G4double &pMax) const
Definition: G4VSolid.cc:378
G4double GetMinExtent(const EAxis pAxis) const
virtual const G4ReflectedSolid * GetReflectedSolidPtr() const
G4VSolid * Clone() const