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