Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4UnionSolid.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 of methods for the class G4IntersectionSolid
30 //
31 // History:
32 //
33 // 12.09.98 V.Grichine: first implementation
34 // 28.11.98 V.Grichine: fix while loops in DistToIn/Out
35 // 27.07.99 V.Grichine: modifications in DistToOut(p,v,...), while -> do-while
36 // 16.03.01 V.Grichine: modifications in CalculateExtent()
37 //
38 // --------------------------------------------------------------------
39 
40 #include <sstream>
41 
42 #include "G4UnionSolid.hh"
43 
44 #include "G4SystemOfUnits.hh"
45 #include "G4VoxelLimits.hh"
46 #include "G4VPVParameterisation.hh"
47 #include "G4GeometryTolerance.hh"
48 
49 #include "G4VGraphicsScene.hh"
50 #include "G4Polyhedron.hh"
51 #include "HepPolyhedronProcessor.h"
52 #include "G4NURBS.hh"
53 // #include "G4NURBSbox.hh"
54 
56 //
57 // Transfer all data members to G4BooleanSolid which is responsible
58 // for them. pName will be in turn sent to G4VSolid
59 
60 G4UnionSolid:: G4UnionSolid( const G4String& pName,
61  G4VSolid* pSolidA ,
62  G4VSolid* pSolidB )
63  : G4BooleanSolid(pName,pSolidA,pSolidB)
64 {
65 }
66 
68 //
69 // Constructor
70 
71 G4UnionSolid::G4UnionSolid( const G4String& pName,
72  G4VSolid* pSolidA ,
73  G4VSolid* pSolidB ,
74  G4RotationMatrix* rotMatrix,
75  const G4ThreeVector& transVector )
76  : G4BooleanSolid(pName,pSolidA,pSolidB,rotMatrix,transVector)
77 
78 {
79 }
80 
82 //
83 // Constructor
84 
85 G4UnionSolid::G4UnionSolid( const G4String& pName,
86  G4VSolid* pSolidA ,
87  G4VSolid* pSolidB ,
88  const G4Transform3D& transform )
89  : G4BooleanSolid(pName,pSolidA,pSolidB,transform)
90 {
91 }
92 
94 //
95 // Fake default constructor - sets only member data and allocates memory
96 // for usage restricted to object persistency.
97 
98 G4UnionSolid::G4UnionSolid( __void__& a )
99  : G4BooleanSolid(a)
100 {
101 }
102 
104 //
105 // Destructor
106 
107 G4UnionSolid::~G4UnionSolid()
108 {
109 }
110 
112 //
113 // Copy constructor
114 
115 G4UnionSolid::G4UnionSolid(const G4UnionSolid& rhs)
116  : G4BooleanSolid (rhs)
117 {
118 }
119 
121 //
122 // Assignment operator
123 
124 G4UnionSolid& G4UnionSolid::operator = (const G4UnionSolid& rhs)
125 {
126  // Check assignment to self
127  //
128  if (this == &rhs) { return *this; }
129 
130  // Copy base class data
131  //
132  G4BooleanSolid::operator=(rhs);
133 
134  return *this;
135 }
136 
138 //
139 //
140 
141 G4bool
142 G4UnionSolid::CalculateExtent( const EAxis pAxis,
143  const G4VoxelLimits& pVoxelLimit,
144  const G4AffineTransform& pTransform,
145  G4double& pMin,
146  G4double& pMax ) const
147 {
148  G4bool touchesA, touchesB, out ;
149  G4double minA = kInfinity, minB = kInfinity,
150  maxA = -kInfinity, maxB = -kInfinity;
151 
152  touchesA = fPtrSolidA->CalculateExtent( pAxis, pVoxelLimit,
153  pTransform, minA, maxA);
154  touchesB= fPtrSolidB->CalculateExtent( pAxis, pVoxelLimit,
155  pTransform, minB, maxB);
156  if( touchesA || touchesB )
157  {
158  pMin = std::min( minA, minB );
159  pMax = std::max( maxA, maxB );
160  out = true ;
161  }
162  else out = false ;
163 
164  return out ; // It exists in this slice if either one does.
165 }
166 
168 //
169 // Important comment: When solids A and B touch together along flat
170 // surface the surface points will be considered as kSurface, while points
171 // located around will correspond to kInside
172 
173 EInside G4UnionSolid::Inside( const G4ThreeVector& p ) const
174 {
175  EInside positionA = fPtrSolidA->Inside(p);
176  if (positionA == kInside) { return kInside; }
177 
178  EInside positionB = fPtrSolidB->Inside(p);
179 
180  if( positionB == kInside ||
181  ( positionA == kSurface && positionB == kSurface &&
182  ( fPtrSolidA->SurfaceNormal(p) +
183  fPtrSolidB->SurfaceNormal(p) ).mag2() <
184  1000*G4GeometryTolerance::GetInstance()->GetRadialTolerance() ) )
185  {
186  return kInside;
187  }
188  else
189  {
190  if( ( positionB == kSurface ) || ( positionA == kSurface ) )
191  { return kSurface; }
192  else
193  { return kOutside; }
194  }
195 }
196 
198 //
199 //
200 
201 G4ThreeVector
202 G4UnionSolid::SurfaceNormal( const G4ThreeVector& p ) const
203 {
204  G4ThreeVector normal;
205 
206 #ifdef G4BOOLDEBUG
207  if( Inside(p) == kOutside )
208  {
209  G4cout << "WARNING - Invalid call in "
210  << "G4UnionSolid::SurfaceNormal(p)" << G4endl
211  << " Point p is outside !" << G4endl;
212  G4cout << " p = " << p << G4endl;
213  G4cerr << "WARNING - Invalid call in "
214  << "G4UnionSolid::SurfaceNormal(p)" << G4endl
215  << " Point p is outside !" << G4endl;
216  G4cerr << " p = " << p << G4endl;
217  }
218 #endif
219 
220  if(fPtrSolidA->Inside(p) == kSurface && fPtrSolidB->Inside(p) != kInside)
221  {
222  normal= fPtrSolidA->SurfaceNormal(p) ;
223  }
224  else if(fPtrSolidB->Inside(p) == kSurface &&
225  fPtrSolidA->Inside(p) != kInside)
226  {
227  normal= fPtrSolidB->SurfaceNormal(p) ;
228  }
229  else
230  {
231  normal= fPtrSolidA->SurfaceNormal(p) ;
232 #ifdef G4BOOLDEBUG
233  if(Inside(p)==kInside)
234  {
235  G4cout << "WARNING - Invalid call in "
236  << "G4UnionSolid::SurfaceNormal(p)" << G4endl
237  << " Point p is inside !" << G4endl;
238  G4cout << " p = " << p << G4endl;
239  G4cerr << "WARNING - Invalid call in "
240  << "G4UnionSolid::SurfaceNormal(p)" << G4endl
241  << " Point p is inside !" << G4endl;
242  G4cerr << " p = " << p << G4endl;
243  }
244 #endif
245  }
246  return normal;
247 }
248 
250 //
251 // The same algorithm as in DistanceToIn(p)
252 
253 G4double
254 G4UnionSolid::DistanceToIn( const G4ThreeVector& p,
255  const G4ThreeVector& v ) const
256 {
257 #ifdef G4BOOLDEBUG
258  if( Inside(p) == kInside )
259  {
260  G4cout << "WARNING - Invalid call in "
261  << "G4UnionSolid::DistanceToIn(p,v)" << G4endl
262  << " Point p is inside !" << G4endl;
263  G4cout << " p = " << p << G4endl;
264  G4cout << " v = " << v << G4endl;
265  G4cerr << "WARNING - Invalid call in "
266  << "G4UnionSolid::DistanceToIn(p,v)" << G4endl
267  << " Point p is inside !" << G4endl;
268  G4cerr << " p = " << p << G4endl;
269  G4cerr << " v = " << v << G4endl;
270  }
271 #endif
272 
273  return std::min(fPtrSolidA->DistanceToIn(p,v),
274  fPtrSolidB->DistanceToIn(p,v) ) ;
275 }
276 
278 //
279 // Approximate nearest distance from the point p to the union of
280 // two solids
281 
282 G4double
283 G4UnionSolid::DistanceToIn( const G4ThreeVector& p) const
284 {
285 #ifdef G4BOOLDEBUG
286  if( Inside(p) == kInside )
287  {
288  G4cout << "WARNING - Invalid call in "
289  << "G4UnionSolid::DistanceToIn(p)" << G4endl
290  << " Point p is inside !" << G4endl;
291  G4cout << " p = " << p << G4endl;
292  G4cerr << "WARNING - Invalid call in "
293  << "G4UnionSolid::DistanceToIn(p)" << G4endl
294  << " Point p is inside !" << G4endl;
295  G4cerr << " p = " << p << G4endl;
296  }
297 #endif
298  G4double distA = fPtrSolidA->DistanceToIn(p) ;
299  G4double distB = fPtrSolidB->DistanceToIn(p) ;
300  G4double safety = std::min(distA,distB) ;
301  if(safety < 0.0) safety = 0.0 ;
302  return safety ;
303 }
304 
306 //
307 // The same algorithm as DistanceToOut(p)
308 
309 G4double
310 G4UnionSolid::DistanceToOut( const G4ThreeVector& p,
311  const G4ThreeVector& v,
312  const G4bool calcNorm,
313  G4bool *validNorm,
314  G4ThreeVector *n ) const
315 {
316  G4double dist = 0.0, disTmp = 0.0 ;
317  G4ThreeVector normTmp;
318  G4ThreeVector* nTmp= &normTmp;
319 
320  if( Inside(p) == kOutside )
321  {
322 #ifdef G4BOOLDEBUG
323  G4cout << "Position:" << G4endl << G4endl;
324  G4cout << "p.x() = " << p.x()/mm << " mm" << G4endl;
325  G4cout << "p.y() = " << p.y()/mm << " mm" << G4endl;
326  G4cout << "p.z() = " << p.z()/mm << " mm" << G4endl << G4endl;
327  G4cout << "Direction:" << G4endl << G4endl;
328  G4cout << "v.x() = " << v.x() << G4endl;
329  G4cout << "v.y() = " << v.y() << G4endl;
330  G4cout << "v.z() = " << v.z() << G4endl << G4endl;
331  G4cout << "WARNING - Invalid call in "
332  << "G4UnionSolid::DistanceToOut(p,v)" << G4endl
333  << " Point p is outside !" << G4endl;
334  G4cout << " p = " << p << G4endl;
335  G4cout << " v = " << v << G4endl;
336  G4cerr << "WARNING - Invalid call in "
337  << "G4UnionSolid::DistanceToOut(p,v)" << G4endl
338  << " Point p is outside !" << G4endl;
339  G4cerr << " p = " << p << G4endl;
340  G4cerr << " v = " << v << G4endl;
341 #endif
342  }
343  else
344  {
345  EInside positionA = fPtrSolidA->Inside(p) ;
346  // EInside positionB = fPtrSolidB->Inside(p) ;
347 
348  if( positionA != kOutside )
349  {
350  do
351  {
352  disTmp = fPtrSolidA->DistanceToOut(p+dist*v,v,calcNorm,
353  validNorm,nTmp) ;
354  dist += disTmp ;
355 
356  if(fPtrSolidB->Inside(p+dist*v) != kOutside)
357  {
358  disTmp = fPtrSolidB->DistanceToOut(p+dist*v,v,calcNorm,
359  validNorm,nTmp) ;
360  dist += disTmp ;
361  }
362  }
363  // while( Inside(p+dist*v) == kInside ) ;
364  while( fPtrSolidA->Inside(p+dist*v) != kOutside &&
365  disTmp > 0.5*kCarTolerance ) ;
366  }
367  else // if( positionB != kOutside )
368  {
369  do
370  {
371  disTmp = fPtrSolidB->DistanceToOut(p+dist*v,v,calcNorm,
372  validNorm,nTmp) ;
373  dist += disTmp ;
374 
375  if(fPtrSolidA->Inside(p+dist*v) != kOutside)
376  {
377  disTmp = fPtrSolidA->DistanceToOut(p+dist*v,v,calcNorm,
378  validNorm,nTmp) ;
379  dist += disTmp ;
380  }
381  }
382  // while( Inside(p+dist*v) == kInside ) ;
383  while( (fPtrSolidB->Inside(p+dist*v) != kOutside)
384  && (disTmp > 0.5*kCarTolerance) ) ;
385  }
386  }
387  if( calcNorm )
388  {
389  *validNorm = false ;
390  *n = *nTmp ;
391  }
392  return dist ;
393 }
394 
396 //
397 // Inverted algorithm of DistanceToIn(p)
398 
399 G4double
400 G4UnionSolid::DistanceToOut( const G4ThreeVector& p ) const
401 {
402  G4double distout = 0.0;
403  if( Inside(p) == kOutside )
404  {
405 #ifdef G4BOOLDEBUG
406  G4cout << "WARNING - Invalid call in "
407  << "G4UnionSolid::DistanceToOut(p)" << G4endl
408  << " Point p is outside !" << G4endl;
409  G4cout << " p = " << p << G4endl;
410  G4cerr << "WARNING - Invalid call in "
411  << "G4UnionSolid::DistanceToOut(p)" << G4endl
412  << " Point p is outside !" << G4endl;
413  G4cerr << " p = " << p << G4endl;
414 #endif
415  }
416  else
417  {
418  EInside positionA = fPtrSolidA->Inside(p) ;
419  EInside positionB = fPtrSolidB->Inside(p) ;
420 
421  // Is this equivalent ??
422  // if( ! ( (positionA == kOutside)) &&
423  // (positionB == kOutside)) )
424  if((positionA == kInside && positionB == kInside ) ||
425  (positionA == kInside && positionB == kSurface ) ||
426  (positionA == kSurface && positionB == kInside ) )
427  {
428  distout= std::max(fPtrSolidA->DistanceToOut(p),
429  fPtrSolidB->DistanceToOut(p) ) ;
430  }
431  else
432  {
433  if(positionA == kOutside)
434  {
435  distout= fPtrSolidB->DistanceToOut(p) ;
436  }
437  else
438  {
439  distout= fPtrSolidA->DistanceToOut(p) ;
440  }
441  }
442  }
443  return distout;
444 }
445 
447 //
448 //
449 
450 G4GeometryType G4UnionSolid::GetEntityType() const
451 {
452  return G4String("G4UnionSolid");
453 }
454 
456 //
457 // Make a clone of the object
458 
459 G4VSolid* G4UnionSolid::Clone() const
460 {
461  return new G4UnionSolid(*this);
462 }
463 
465 //
466 //
467 
468 void
469 G4UnionSolid::ComputeDimensions( G4VPVParameterisation*,
470  const G4int,
471  const G4VPhysicalVolume* )
472 {
473 }
474 
476 //
477 //
478 
479 void
480 G4UnionSolid::DescribeYourselfTo ( G4VGraphicsScene& scene ) const
481 {
482  scene.AddSolid (*this);
483 }
484 
486 //
487 //
488 
489 G4Polyhedron*
490 G4UnionSolid::CreatePolyhedron () const
491 {
492  HepPolyhedronProcessor processor;
493  // Stack components and components of components recursively
494  // See G4BooleanSolid::StackPolyhedron
495  G4Polyhedron* top = StackPolyhedron(processor, this);
496  G4Polyhedron* result = new G4Polyhedron(*top);
497  if (processor.execute(*result)) { return result; }
498  else { return 0; }
499 }
500 
502 //
503 //
504 
505 G4NURBS*
506 G4UnionSolid::CreateNURBS () const
507 {
508  // Take into account boolean operation - see CreatePolyhedron.
509  // return new G4NURBSbox (1.0, 1.0, 1.0);
510  return 0;
511 }
Generated on Sat May 25 2013 14:33:12 for Geant4 by   doxygen 1.8.4