source: trunk/source/geometry/management/src/G4ReflectedSolid.cc@ 1185

Last change on this file since 1185 was 1058, checked in by garnier, 17 years ago

file release beta
File size: 17.9 KB
Line 
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,v 1.11 2006/11/08 09:56:33 gcosmo Exp $
28//
29// GEANT4 tag $Name: geant4-09-02-ref-02 $
30//
31// Implementation for G4ReflectedSolid class for boolean
32// operations between other solids
33//
34// Author: Vladimir Grichine, 23.07.01 (Vladimir.Grichine@cern.ch)
35//
36// --------------------------------------------------------------------
37
38#include "G4ReflectedSolid.hh"
39
40#include <sstream>
41
42#include "G4Point3D.hh"
43#include "G4Normal3D.hh"
44
45#include "G4VoxelLimits.hh"
46
47#include "G4VPVParameterisation.hh"
48
49#include "G4VGraphicsScene.hh"
50#include "G4Polyhedron.hh"
51#include "G4NURBS.hh"
52// #include "G4NURBSbox.hh"
53
54
55/////////////////////////////////////////////////////////////////
56//
57// Constructor using HepTransform3D, in fact HepReflect3D
58
59G4ReflectedSolid::G4ReflectedSolid( const G4String& pName,
60 G4VSolid* pSolid ,
61 const G4Transform3D& transform )
62 : G4VSolid(pName), fpPolyhedron(0)
63{
64 fPtrSolid = pSolid ;
65 G4RotationMatrix rotMatrix ;
66
67 fDirectTransform =
68 new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
69 fPtrTransform =
70 new G4AffineTransform(rotMatrix, transform.getTranslation()) ;
71 fPtrTransform->Invert() ;
72
73 fDirectTransform3D = new G4Transform3D(transform) ;
74 fPtrTransform3D = new G4Transform3D(transform.inverse()) ;
75}
76
77///////////////////////////////////////////////////////////////////
78//
79
80G4ReflectedSolid::~G4ReflectedSolid()
81{
82 if(fPtrTransform)
83 {
84 delete fPtrTransform; fPtrTransform=0;
85 delete fDirectTransform; fDirectTransform=0;
86 }
87 if(fPtrTransform3D)
88 {
89 delete fPtrTransform3D; fPtrTransform3D=0;
90 delete fDirectTransform3D; fDirectTransform3D=0;
91 }
92 delete fpPolyhedron;
93}
94
95G4GeometryType G4ReflectedSolid::GetEntityType() const
96{
97 return G4String("G4ReflectedSolid");
98}
99
100const G4ReflectedSolid* G4ReflectedSolid::GetReflectedSolidPtr() const
101{
102 return this;
103}
104
105G4ReflectedSolid* G4ReflectedSolid::GetReflectedSolidPtr()
106{
107 return this;
108}
109
110G4VSolid* G4ReflectedSolid::GetConstituentMovedSolid() const
111{
112 return fPtrSolid;
113}
114
115/////////////////////////////////////////////////////////////////////////////
116
117G4AffineTransform G4ReflectedSolid::GetTransform() const
118{
119 G4AffineTransform aTransform = *fPtrTransform;
120 return aTransform;
121}
122
123void G4ReflectedSolid::SetTransform(G4AffineTransform& transform)
124{
125 fPtrTransform = &transform ;
126 fpPolyhedron = 0;
127}
128
129//////////////////////////////////////////////////////////////////////////////
130
131G4AffineTransform G4ReflectedSolid::GetDirectTransform() const
132{
133 G4AffineTransform aTransform= *fDirectTransform;
134 return aTransform;
135}
136
137void G4ReflectedSolid::SetDirectTransform(G4AffineTransform& transform)
138{
139 fDirectTransform = &transform ;
140 fpPolyhedron = 0;
141}
142
143/////////////////////////////////////////////////////////////////////////////
144
145G4Transform3D G4ReflectedSolid::GetTransform3D() const
146{
147 G4Transform3D aTransform = *fPtrTransform3D;
148 return aTransform;
149}
150
151void G4ReflectedSolid::SetTransform3D(G4Transform3D& transform)
152{
153 fPtrTransform3D = &transform ;
154 fpPolyhedron = 0;
155}
156
157//////////////////////////////////////////////////////////////////////////////
158
159G4Transform3D G4ReflectedSolid::GetDirectTransform3D() const
160{
161 G4Transform3D aTransform= *fDirectTransform3D;
162 return aTransform;
163}
164
165void G4ReflectedSolid::SetDirectTransform3D(G4Transform3D& transform)
166{
167 fDirectTransform3D = &transform ;
168 fpPolyhedron = 0;
169}
170
171/////////////////////////////////////////////////////////////////////////////
172
173G4RotationMatrix G4ReflectedSolid::GetFrameRotation() const
174{
175 G4RotationMatrix InvRotation= fDirectTransform->NetRotation();
176 return InvRotation;
177}
178
179void G4ReflectedSolid::SetFrameRotation(const G4RotationMatrix& matrix)
180{
181 fDirectTransform->SetNetRotation(matrix);
182}
183
184/////////////////////////////////////////////////////////////////////////////
185
186G4ThreeVector G4ReflectedSolid::GetFrameTranslation() const
187{
188 return fPtrTransform->NetTranslation();
189}
190
191void G4ReflectedSolid::SetFrameTranslation(const G4ThreeVector& vector)
192{
193 fPtrTransform->SetNetTranslation(vector);
194}
195
196///////////////////////////////////////////////////////////////
197
198G4RotationMatrix G4ReflectedSolid::GetObjectRotation() const
199{
200 G4RotationMatrix Rotation= fPtrTransform->NetRotation();
201 return Rotation;
202}
203
204void G4ReflectedSolid::SetObjectRotation(const G4RotationMatrix& matrix)
205{
206 fPtrTransform->SetNetRotation(matrix);
207}
208
209///////////////////////////////////////////////////////////////////////
210
211G4ThreeVector G4ReflectedSolid::GetObjectTranslation() const
212{
213 return fDirectTransform->NetTranslation();
214}
215
216void G4ReflectedSolid::SetObjectTranslation(const G4ThreeVector& vector)
217{
218 fDirectTransform->SetNetTranslation(vector);
219}
220
221///////////////////////////////////////////////////////////////
222//
223//
224
225G4bool
226G4ReflectedSolid::CalculateExtent( const EAxis pAxis,
227 const G4VoxelLimits& pVoxelLimit,
228 const G4AffineTransform& pTransform,
229 G4double& pMin,
230 G4double& pMax ) const
231{
232
233 G4VoxelLimits unLimit;
234 G4AffineTransform unTransform;
235
236 G4double x1 = -kInfinity, x2 = kInfinity,
237 y1 = -kInfinity, y2 = kInfinity,
238 z1 = -kInfinity, z2 = kInfinity;
239
240 G4bool existsAfterClip = false ;
241 existsAfterClip =
242 fPtrSolid->CalculateExtent(kXAxis,unLimit,unTransform,x1,x2);
243 existsAfterClip =
244 fPtrSolid->CalculateExtent(kYAxis,unLimit,unTransform,y1,y2);
245 existsAfterClip =
246 fPtrSolid->CalculateExtent(kZAxis,unLimit,unTransform,z1,z2);
247
248 existsAfterClip = false;
249 pMin = +kInfinity ;
250 pMax = -kInfinity ;
251
252 G4Transform3D pTransform3D = G4Transform3D(pTransform.NetRotation().inverse(),
253 pTransform.NetTranslation());
254
255 G4Transform3D transform3D = pTransform3D*(*fDirectTransform3D);
256
257 G4Point3D tmpPoint;
258
259// Calculate rotated vertex coordinates
260
261 G4ThreeVectorList* vertices = new G4ThreeVectorList();
262 vertices->reserve(8);
263
264 if (vertices)
265 {
266 G4ThreeVector vertex0(x1,y1,z1) ;
267 tmpPoint = transform3D*G4Point3D(vertex0);
268 vertex0 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
269 vertices->push_back(vertex0);
270
271 G4ThreeVector vertex1(x2,y1,z1) ;
272 tmpPoint = transform3D*G4Point3D(vertex1);
273 vertex1 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
274 vertices->push_back(vertex1);
275
276 G4ThreeVector vertex2(x2,y2,z1) ;
277 tmpPoint = transform3D*G4Point3D(vertex2);
278 vertex2 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
279 vertices->push_back(vertex2);
280
281 G4ThreeVector vertex3(x1,y2,z1) ;
282 tmpPoint = transform3D*G4Point3D(vertex3);
283 vertex3 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
284 vertices->push_back(vertex3);
285
286 G4ThreeVector vertex4(x1,y1,z2) ;
287 tmpPoint = transform3D*G4Point3D(vertex4);
288 vertex4 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
289 vertices->push_back(vertex4);
290
291 G4ThreeVector vertex5(x2,y1,z2) ;
292 tmpPoint = transform3D*G4Point3D(vertex5);
293 vertex5 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
294 vertices->push_back(vertex5);
295
296 G4ThreeVector vertex6(x2,y2,z2) ;
297 tmpPoint = transform3D*G4Point3D(vertex6);
298 vertex6 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
299 vertices->push_back(vertex6);
300
301 G4ThreeVector vertex7(x1,y2,z2) ;
302 tmpPoint = transform3D*G4Point3D(vertex7);
303 vertex7 = G4ThreeVector(tmpPoint.x(),tmpPoint.y(),tmpPoint.z());
304 vertices->push_back(vertex7);
305 }
306 else
307 {
308 DumpInfo();
309 G4Exception("G4ReflectedSolid::CalculateExtent()",
310 "FatalError", FatalException,
311 "Error in allocation of vertices. Out of memory !");
312 }
313
314 ClipCrossSection(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
315 ClipCrossSection(vertices,4,pVoxelLimit,pAxis,pMin,pMax) ;
316 ClipBetweenSections(vertices,0,pVoxelLimit,pAxis,pMin,pMax) ;
317
318 if (pVoxelLimit.IsLimited(pAxis) == false)
319 {
320 if ( pMin != kInfinity || pMax != -kInfinity )
321 {
322 existsAfterClip = true ;
323
324 // Add 2*tolerance to avoid precision troubles
325
326 pMin -= kCarTolerance;
327 pMax += kCarTolerance;
328 }
329 }
330 else
331 {
332 G4ThreeVector clipCentre(
333 ( pVoxelLimit.GetMinXExtent()+pVoxelLimit.GetMaxXExtent())*0.5,
334 ( pVoxelLimit.GetMinYExtent()+pVoxelLimit.GetMaxYExtent())*0.5,
335 ( pVoxelLimit.GetMinZExtent()+pVoxelLimit.GetMaxZExtent())*0.5);
336
337 if ( pMin != kInfinity || pMax != -kInfinity )
338 {
339 existsAfterClip = true ;
340
341
342 // Check to see if endpoints are in the solid
343
344 clipCentre(pAxis) = pVoxelLimit.GetMinExtent(pAxis);
345
346 if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
347 {
348 pMin = pVoxelLimit.GetMinExtent(pAxis);
349 }
350 else
351 {
352 pMin -= kCarTolerance;
353 }
354 clipCentre(pAxis) = pVoxelLimit.GetMaxExtent(pAxis);
355
356 if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
357 {
358 pMax = pVoxelLimit.GetMaxExtent(pAxis);
359 }
360 else
361 {
362 pMax += kCarTolerance;
363 }
364 }
365 // Check for case where completely enveloping clipping volume
366 // If point inside then we are confident that the solid completely
367 // envelopes the clipping volume. Hence set min/max extents according
368 // to clipping volume extents along the specified axis.
369
370 else if (Inside(transform3D.inverse()*G4Point3D(clipCentre)) != kOutside)
371 {
372 existsAfterClip = true ;
373 pMin = pVoxelLimit.GetMinExtent(pAxis) ;
374 pMax = pVoxelLimit.GetMaxExtent(pAxis) ;
375 }
376 }
377 delete vertices;
378 return existsAfterClip;
379}
380
381/////////////////////////////////////////////////////
382//
383//
384
385EInside G4ReflectedSolid::Inside(const G4ThreeVector& p) const
386{
387
388 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
389 // G4Point3D newPoint = (*fPtrTransform3D)*G4Point3D(p) ;
390
391 return fPtrSolid->Inside(G4ThreeVector(newPoint.x(),
392 newPoint.y(),
393 newPoint.z())) ;
394}
395
396//////////////////////////////////////////////////////////////
397//
398//
399
400G4ThreeVector
401G4ReflectedSolid::SurfaceNormal( const G4ThreeVector& p ) const
402{
403 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
404 G4ThreeVector normal =
405 fPtrSolid->SurfaceNormal(G4ThreeVector(newPoint.x(),
406 newPoint.y(),
407 newPoint.z() ) ) ;
408 G4Point3D newN = (*fDirectTransform3D)*G4Point3D(normal) ;
409 newN.unit() ;
410
411 return G4ThreeVector(newN.x(),newN.y(),newN.z()) ;
412}
413
414/////////////////////////////////////////////////////////////
415//
416// The same algorithm as in DistanceToIn(p)
417
418G4double
419G4ReflectedSolid::DistanceToIn( const G4ThreeVector& p,
420 const G4ThreeVector& v ) const
421{
422 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
423 G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v) ;
424 newDirection.unit() ;
425 return fPtrSolid->DistanceToIn(
426 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
427 G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z())) ;
428}
429
430////////////////////////////////////////////////////////
431//
432// Approximate nearest distance from the point p to the intersection of
433// two solids
434
435G4double
436G4ReflectedSolid::DistanceToIn( const G4ThreeVector& p) const
437{
438 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
439 return fPtrSolid->DistanceToIn(
440 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z())) ;
441}
442
443//////////////////////////////////////////////////////////
444//
445// The same algorithm as DistanceToOut(p)
446
447G4double
448G4ReflectedSolid::DistanceToOut( const G4ThreeVector& p,
449 const G4ThreeVector& v,
450 const G4bool calcNorm,
451 G4bool *validNorm,
452 G4ThreeVector *n ) const
453{
454 G4ThreeVector solNorm ;
455
456 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p) ;
457 G4Point3D newDirection = (*fDirectTransform3D)*G4Point3D(v);
458 newDirection.unit() ;
459
460 G4double dist =
461 fPtrSolid->DistanceToOut(
462 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()),
463 G4ThreeVector(newDirection.x(),newDirection.y(),newDirection.z()),
464 calcNorm, validNorm, &solNorm) ;
465 if(calcNorm)
466 {
467 G4Point3D newN = (*fDirectTransform3D)*G4Point3D(solNorm);
468 newN.unit() ;
469 *n = G4ThreeVector(newN.x(),newN.y(),newN.z());
470 }
471 return dist ;
472}
473
474//////////////////////////////////////////////////////////////
475//
476// Inverted algorithm of DistanceToIn(p)
477
478G4double
479G4ReflectedSolid::DistanceToOut( const G4ThreeVector& p ) const
480{
481 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
482 return fPtrSolid->DistanceToOut(
483 G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z()));
484}
485
486//////////////////////////////////////////////////////////////
487//
488//
489
490void
491G4ReflectedSolid::ComputeDimensions( G4VPVParameterisation*,
492 const G4int,
493 const G4VPhysicalVolume* )
494{
495 DumpInfo();
496 G4Exception("G4BooleanSolid::ComputeDimensions()",
497 "NotApplicable", FatalException,
498 "Method not applicable in this context!");
499}
500
501//////////////////////////////////////////////////////////////
502//
503// Return a point (G4ThreeVector) randomly and uniformly selected
504// on the solid surface
505
506G4ThreeVector G4ReflectedSolid::GetPointOnSurface() const
507{
508 G4ThreeVector p = fPtrSolid->GetPointOnSurface();
509 G4Point3D newPoint = (*fDirectTransform3D)*G4Point3D(p);
510
511 return G4ThreeVector(newPoint.x(),newPoint.y(),newPoint.z());
512}
513
514//////////////////////////////////////////////////////////////////////////
515//
516// Stream object contents to an output stream
517
518std::ostream& G4ReflectedSolid::StreamInfo(std::ostream& os) const
519{
520 os << "-----------------------------------------------------------\n"
521 << " *** Dump for Reflected solid - " << GetName() << " ***\n"
522 << " ===================================================\n"
523 << " Solid type: " << GetEntityType() << "\n"
524 << " Parameters of constituent solid: \n"
525 << "===========================================================\n";
526 fPtrSolid->StreamInfo(os);
527 os << "===========================================================\n"
528 << " Transformations: \n"
529 << " Direct transformation - translation : \n"
530 << " " << fDirectTransform->NetTranslation() << "\n"
531 << " - rotation : \n"
532 << " ";
533 fDirectTransform->NetRotation().print(os);
534 os << "\n"
535 << "===========================================================\n";
536
537 return os;
538}
539
540/////////////////////////////////////////////////
541//
542//
543
544void
545G4ReflectedSolid::DescribeYourselfTo ( G4VGraphicsScene& scene ) const
546{
547 scene.AddSolid (*this);
548}
549
550////////////////////////////////////////////////////
551//
552//
553
554G4Polyhedron*
555G4ReflectedSolid::CreatePolyhedron () const
556{
557 G4Polyhedron* polyhedron = fPtrSolid->CreatePolyhedron();
558 if (polyhedron)
559 {
560 polyhedron->Transform(*fDirectTransform3D);
561 return polyhedron;
562 }
563 else
564 {
565 std::ostringstream oss;
566 oss << "Solid - " << GetName()
567 << " - original solid has no" << G4endl
568 << " corresponding polyhedron. Returning NULL!";
569 G4Exception("G4ReflectedSolid::CreatePolyhedron()", "InvalidSetup",
570 JustWarning, oss.str().c_str());
571 return 0;
572 }
573}
574
575/////////////////////////////////////////////////////////
576//
577//
578
579G4NURBS*
580G4ReflectedSolid::CreateNURBS () const
581{
582 // Take into account local transformation - see CreatePolyhedron.
583 // return fPtrSolid->CreateNURBS() ;
584 return 0;
585}
586
587/////////////////////////////////////////////////////////
588//
589//
590
591G4Polyhedron*
592G4ReflectedSolid::GetPolyhedron () const
593{
594 if (!fpPolyhedron ||
595 fpPolyhedron->GetNumberOfRotationStepsAtTimeOfCreation() !=
596 fpPolyhedron->GetNumberOfRotationSteps())
597 {
598 delete fpPolyhedron;
599 fpPolyhedron = CreatePolyhedron ();
600 }
601 return fpPolyhedron;
602}
Note: See TracBrowser for help on using the repository browser.