Geant4  10.00.p02
G4VCSGfaceted.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 // the GEANT4 collaboration.
27 //
28 // By copying, distributing or modifying the Program (or any work
29 // based on the Program) you indicate your acceptance of this statement,
30 // and all its terms.
31 //
32 // $Id: G4VCSGfaceted.cc 81641 2014-06-04 09:11:38Z gcosmo $
33 //
34 //
35 // --------------------------------------------------------------------
36 // GEANT 4 class source file
37 //
38 //
39 // G4VCSGfaceted.cc
40 //
41 // Implementation of the virtual class of a CSG type shape that is built
42 // entirely out of G4VCSGface faces.
43 //
44 // --------------------------------------------------------------------
45 
46 #include "G4VCSGfaceted.hh"
47 #include "G4VCSGface.hh"
48 #include "G4SolidExtentList.hh"
49 
50 #include "G4VoxelLimits.hh"
51 #include "G4AffineTransform.hh"
52 
53 #include "Randomize.hh"
54 
55 #include "G4Polyhedron.hh"
56 #include "G4VGraphicsScene.hh"
57 #include "G4VisExtent.hh"
58 
59 //
60 // Constructor
61 //
63  : G4VSolid(name),
64  numFace(0), faces(0), fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0),
65  fStatistics(1000000), fCubVolEpsilon(0.001), fAreaAccuracy(-1.)
66 {
67 }
68 
69 
70 //
71 // Fake default constructor - sets only member data and allocates memory
72 // for usage restricted to object persistency.
73 //
75  : G4VSolid(a),
76  numFace(0), faces(0), fCubicVolume(0.), fSurfaceArea(0.), fpPolyhedron(0),
77  fStatistics(1000000), fCubVolEpsilon(0.001), fAreaAccuracy(-1.)
78 {
79 }
80 
81 //
82 // Destructor
83 //
85 {
86  DeleteStuff();
87  delete fpPolyhedron;
88 }
89 
90 
91 //
92 // Copy constructor
93 //
95  : G4VSolid( source )
96 {
97  fStatistics = source.fStatistics;
100 
101  CopyStuff( source );
102 }
103 
104 
105 //
106 // Assignment operator
107 //
109 {
110  if (&source == this) { return *this; }
111 
112  // Copy base class data
113  //
114  G4VSolid::operator=(source);
115 
116  // Copy data
117  //
118  fStatistics = source.fStatistics;
120  fAreaAccuracy = source.fAreaAccuracy;
121 
122  DeleteStuff();
123  CopyStuff( source );
124 
125  return *this;
126 }
127 
128 
129 //
130 // CopyStuff (protected)
131 //
132 // Copy the contents of source
133 //
135 {
136  numFace = source.numFace;
137  if (numFace == 0) { return; } // odd, but permissable?
138 
139  faces = new G4VCSGface*[numFace];
140 
141  G4VCSGface **face = faces,
142  **sourceFace = source.faces;
143  do
144  {
145  *face = (*sourceFace)->Clone();
146  } while( ++sourceFace, ++face < faces+numFace );
147  fCubicVolume = source.fCubicVolume;
148  fSurfaceArea = source.fSurfaceArea;
149  fpPolyhedron = 0;
150 }
151 
152 
153 //
154 // DeleteStuff (protected)
155 //
156 // Delete all allocated objects
157 //
159 {
160  if (numFace)
161  {
162  G4VCSGface **face = faces;
163  do
164  {
165  delete *face;
166  } while( ++face < faces + numFace );
167 
168  delete [] faces;
169  }
170  delete fpPolyhedron; fpPolyhedron = 0;
171 }
172 
173 
174 //
175 // CalculateExtent
176 //
178  const G4VoxelLimits &voxelLimit,
179  const G4AffineTransform &transform,
180  G4double &min,
181  G4double &max ) const
182 {
183  G4SolidExtentList extentList( axis, voxelLimit );
184 
185  //
186  // Loop over all faces, checking min/max extent as we go.
187  //
188  G4VCSGface **face = faces;
189  do
190  {
191  (*face)->CalculateExtent( axis, voxelLimit, transform, extentList );
192  } while( ++face < faces + numFace );
193 
194  //
195  // Return min/max value
196  //
197  return extentList.GetExtent( min, max );
198 }
199 
200 
201 //
202 // Inside
203 //
204 // It could be a good idea to override this virtual
205 // member to add first a simple test (such as spherical
206 // test or whatnot) and to call this version only if
207 // the simplier test fails.
208 //
210 {
211  EInside answer=kOutside;
212  G4VCSGface **face = faces;
213  G4double best = kInfinity;
214  do
215  {
216  G4double distance;
217  EInside result = (*face)->Inside( p, kCarTolerance/2, &distance );
218  if (result == kSurface) { return kSurface; }
219  if (distance < best)
220  {
221  best = distance;
222  answer = result;
223  }
224  } while( ++face < faces + numFace );
225 
226  return answer;
227 }
228 
229 
230 //
231 // SurfaceNormal
232 //
234 {
235  G4ThreeVector answer;
236  G4VCSGface **face = faces;
237  G4double best = kInfinity;
238  do
239  {
240  G4double distance;
241  G4ThreeVector normal = (*face)->Normal( p, &distance );
242  if (distance < best)
243  {
244  best = distance;
245  answer = normal;
246  }
247  } while( ++face < faces + numFace );
248 
249  return answer;
250 }
251 
252 
253 //
254 // DistanceToIn(p,v)
255 //
257  const G4ThreeVector &v ) const
258 {
259  G4double distance = kInfinity;
260  G4double distFromSurface = kInfinity;
261  G4VCSGface **face = faces;
262  G4VCSGface *bestFace = *face;
263  do
264  {
265  G4double faceDistance,
266  faceDistFromSurface;
267  G4ThreeVector faceNormal;
268  G4bool faceAllBehind;
269  if ((*face)->Intersect( p, v, false, kCarTolerance/2,
270  faceDistance, faceDistFromSurface,
271  faceNormal, faceAllBehind ) )
272  {
273  //
274  // Intersecting face
275  //
276  if (faceDistance < distance)
277  {
278  distance = faceDistance;
279  distFromSurface = faceDistFromSurface;
280  bestFace = *face;
281  if (distFromSurface <= 0) { return 0; }
282  }
283  }
284  } while( ++face < faces + numFace );
285 
286  if (distance < kInfinity && distFromSurface<kCarTolerance/2)
287  {
288  if (bestFace->Distance(p,false) < kCarTolerance/2) { distance = 0; }
289  }
290 
291  return distance;
292 }
293 
294 
295 //
296 // DistanceToIn(p)
297 //
299 {
300  return DistanceTo( p, false );
301 }
302 
303 
304 //
305 // DistanceToOut(p,v)
306 //
308  const G4ThreeVector &v,
309  const G4bool calcNorm,
310  G4bool *validNorm,
311  G4ThreeVector *n ) const
312 {
313  G4bool allBehind = true;
314  G4double distance = kInfinity;
315  G4double distFromSurface = kInfinity;
317 
318  G4VCSGface **face = faces;
319  G4VCSGface *bestFace = *face;
320  do
321  {
322  G4double faceDistance,
323  faceDistFromSurface;
324  G4ThreeVector faceNormal;
325  G4bool faceAllBehind;
326  if ((*face)->Intersect( p, v, true, kCarTolerance/2,
327  faceDistance, faceDistFromSurface,
328  faceNormal, faceAllBehind ) )
329  {
330  //
331  // Intersecting face
332  //
333  if ( (distance < kInfinity) || (!faceAllBehind) ) { allBehind = false; }
334  if (faceDistance < distance)
335  {
336  distance = faceDistance;
337  distFromSurface = faceDistFromSurface;
338  normal = faceNormal;
339  bestFace = *face;
340  if (distFromSurface <= 0) { break; }
341  }
342  }
343  } while( ++face < faces + numFace );
344 
345  if (distance < kInfinity)
346  {
347  if (distFromSurface <= 0)
348  {
349  distance = 0;
350  }
351  else if (distFromSurface<kCarTolerance/2)
352  {
353  if (bestFace->Distance(p,true) < kCarTolerance/2) { distance = 0; }
354  }
355 
356  if (calcNorm)
357  {
358  *validNorm = allBehind;
359  *n = normal;
360  }
361  }
362  else
363  {
364  if (Inside(p) == kSurface) { distance = 0; }
365  if (calcNorm) { *validNorm = false; }
366  }
367 
368  return distance;
369 }
370 
371 
372 //
373 // DistanceToOut(p)
374 //
376 {
377  return DistanceTo( p, true );
378 }
379 
380 
381 //
382 // DistanceTo
383 //
384 // Protected routine called by DistanceToIn and DistanceToOut
385 //
387  const G4bool outgoing ) const
388 {
389  G4VCSGface **face = faces;
390  G4double best = kInfinity;
391  do
392  {
393  G4double distance = (*face)->Distance( p, outgoing );
394  if (distance < best) { best = distance; }
395  } while( ++face < faces + numFace );
396 
397  return (best < 0.5*kCarTolerance) ? 0 : best;
398 }
399 
400 
401 //
402 // DescribeYourselfTo
403 //
405 {
406  scene.AddSolid( *this );
407 }
408 
409 
410 //
411 // GetExtent
412 //
413 // Define the sides of the box into which our solid instance would fit.
414 //
416 {
417  static const G4ThreeVector xMax(1,0,0), xMin(-1,0,0),
418  yMax(0,1,0), yMin(0,-1,0),
419  zMax(0,0,1), zMin(0,0,-1);
420  static const G4ThreeVector *axes[6] =
421  { &xMin, &xMax, &yMin, &yMax, &zMin, &zMax };
422 
423  G4double answers[6] =
424  {-kInfinity, -kInfinity, -kInfinity, -kInfinity, -kInfinity, -kInfinity};
425 
426  G4VCSGface **face = faces;
427  do
428  {
429  const G4ThreeVector **axis = axes+5 ;
430  G4double *answer = answers+5;
431  do
432  {
433  G4double testFace = (*face)->Extent( **axis );
434  if (testFace > *answer) { *answer = testFace; }
435  }
436  while( --axis, --answer >= answers );
437 
438  } while( ++face < faces + numFace );
439 
440  return G4VisExtent( -answers[0], answers[1],
441  -answers[2], answers[3],
442  -answers[4], answers[5] );
443 }
444 
445 
446 //
447 // GetEntityType
448 //
450 {
451  return G4String("G4CSGfaceted");
452 }
453 
454 
455 //
456 // Stream object contents to an output stream
457 //
458 std::ostream& G4VCSGfaceted::StreamInfo( std::ostream& os ) const
459 {
460  os << "-----------------------------------------------------------\n"
461  << " *** Dump for solid - " << GetName() << " ***\n"
462  << " ===================================================\n"
463  << " Solid type: G4VCSGfaceted\n"
464  << " Parameters: \n"
465  << " number of faces: " << numFace << "\n"
466  << "-----------------------------------------------------------\n";
467 
468  return os;
469 }
470 
471 
472 //
473 // GetCubVolStatistics
474 //
476 {
477  return fStatistics;
478 }
479 
480 
481 //
482 // GetCubVolEpsilon
483 //
485 {
486  return fCubVolEpsilon;
487 }
488 
489 
490 //
491 // SetCubVolStatistics
492 //
494 {
495  fCubicVolume=0.;
496  fStatistics=st;
497 }
498 
499 
500 //
501 // SetCubVolEpsilon
502 //
504 {
505  fCubicVolume=0.;
506  fCubVolEpsilon=ep;
507 }
508 
509 
510 //
511 // GetAreaStatistics
512 //
514 {
515  return fStatistics;
516 }
517 
518 
519 //
520 // GetAreaAccuracy
521 //
523 {
524  return fAreaAccuracy;
525 }
526 
527 
528 //
529 // SetAreaStatistics
530 //
532 {
533  fSurfaceArea=0.;
534  fStatistics=st;
535 }
536 
537 
538 //
539 // SetAreaAccuracy
540 //
542 {
543  fSurfaceArea=0.;
544  fAreaAccuracy=ep;
545 }
546 
547 
548 //
549 // GetCubicVolume
550 //
552 {
553  if(fCubicVolume != 0.) {;}
555  return fCubicVolume;
556 }
557 
558 
559 //
560 // GetSurfaceArea
561 //
563 {
564  if(fSurfaceArea != 0.) {;}
566  return fSurfaceArea;
567 }
568 
569 
570 //
571 // GetPolyhedron
572 //
574 {
575  if (!fpPolyhedron ||
577  fpPolyhedron->GetNumberOfRotationSteps())
578  {
579  delete fpPolyhedron;
581  }
582  return fpPolyhedron;
583 }
584 
585 
586 //
587 // GetPointOnSurfaceGeneric proportional to Areas of faces
588 // in case of GenericPolycone or GenericPolyhedra
589 //
591 {
592  // Preparing variables
593  //
594  G4ThreeVector answer=G4ThreeVector(0.,0.,0.);
595  G4VCSGface **face = faces;
596  G4double area = 0;
597  G4int i;
598  std::vector<G4double> areas;
599 
600  // First step: calculate surface areas
601  //
602  do
603  {
604  G4double result = (*face)->SurfaceArea( );
605  areas.push_back(result);
606  area=area+result;
607  } while( ++face < faces + numFace );
608 
609  // Second Step: choose randomly one surface
610  //
611  G4VCSGface **face1 = faces;
612  G4double chose = area*G4UniformRand();
613  G4double Achose1, Achose2;
614  Achose1=0; Achose2=0.;
615  i=0;
616 
617  do
618  {
619  Achose2+=areas[i];
620  if(chose>=Achose1 && chose<Achose2)
621  {
622  G4ThreeVector point;
623  point= (*face1)->GetPointOnFace();
624  return point;
625  }
626  i++;
627  Achose1=Achose2;
628  } while( ++face1 < faces + numFace );
629 
630  return answer;
631 }
G4String GetName() const
virtual void DescribeYourselfTo(G4VGraphicsScene &scene) const
G4double fSurfaceArea
void SetAreaStatistics(G4int st)
G4VCSGfaceted & operator=(const G4VCSGfaceted &source)
G4Polyhedron * fpPolyhedron
static const G4double kInfinity
Definition: geomdefs.hh:42
CLHEP::Hep3Vector G4ThreeVector
virtual G4VCSGface * Clone()=0
G4String name
Definition: TRTMaterials.hh:40
G4int GetCubVolStatistics() const
G4double EstimateCubicVolume(G4int nStat, G4double epsilon) const
Definition: G4VSolid.cc:203
G4bool GetExtent(G4double &min, G4double &max) const
void SetCubVolStatistics(G4int st)
virtual G4double DistanceToIn(const G4ThreeVector &p, const G4ThreeVector &v) const
virtual ~G4VCSGfaceted()
G4VCSGfaceted(const G4String &name)
G4double a
Definition: TRTMaterials.hh:39
void CopyStuff(const G4VCSGfaceted &source)
virtual void AddSolid(const G4Box &)=0
int G4int
Definition: G4Types.hh:78
G4double fCubVolEpsilon
virtual G4ThreeVector SurfaceNormal(const G4ThreeVector &p) const
static double normal(HepRandomEngine *eptr)
Definition: RandPoisson.cc:77
virtual G4double SurfaceArea()=0
G4ThreeVector GetPointOnSurfaceGeneric() const
virtual G4double GetSurfaceArea()
G4double fCubicVolume
#define G4UniformRand()
Definition: Randomize.hh:87
virtual G4bool CalculateExtent(const EAxis pAxis, const G4VoxelLimits &pVoxelLimit, const G4AffineTransform &pTransform, G4double &pmin, G4double &pmax) const
virtual G4Polyhedron * CreatePolyhedron() const =0
G4double GetAreaAccuracy() const
bool G4bool
Definition: G4Types.hh:79
virtual G4double Extent(const G4ThreeVector axis)=0
virtual G4VisExtent GetExtent() const
const G4int n
G4double EstimateSurfaceArea(G4int nStat, G4double ell) const
Definition: G4VSolid.cc:263
virtual G4GeometryType GetEntityType() const
G4VCSGface ** faces
virtual std::ostream & StreamInfo(std::ostream &os) const
T max(const T t1, const T t2)
brief Return the largest of the two arguments
EInside
Definition: geomdefs.hh:58
EAxis
Definition: geomdefs.hh:54
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double fAreaAccuracy
virtual G4double DistanceTo(const G4ThreeVector &p, const G4bool outgoing) const
void SetAreaAccuracy(G4double ep)
G4VSolid & operator=(const G4VSolid &rhs)
Definition: G4VSolid.cc:110
G4double kCarTolerance
Definition: G4VSolid.hh:305
G4int GetNumberOfRotationStepsAtTimeOfCreation() const
G4int GetAreaStatistics() const
double G4double
Definition: G4Types.hh:76
virtual G4double DistanceToOut(const G4ThreeVector &p, const G4ThreeVector &v, const G4bool calcNorm=false, G4bool *validNorm=0, G4ThreeVector *n=0) const
virtual EInside Inside(const G4ThreeVector &p) const
virtual void CalculateExtent(const EAxis axis, const G4VoxelLimits &voxelLimit, const G4AffineTransform &tranform, G4SolidExtentList &extentList)=0
virtual G4double GetCubicVolume()
G4double GetCubVolEpsilon() const
void SetCubVolEpsilon(G4double ep)
virtual G4double Distance(const G4ThreeVector &p, G4bool outgoing)=0
virtual G4Polyhedron * GetPolyhedron() const