source: trunk/source/geometry/navigation/src/G4RegularNavigation.cc@ 881

Last change on this file since 881 was 831, checked in by garnier, 17 years ago

import all except CVS
  • Property svn:executable set to *
File size: 12.6 KB
RevLine 
[831]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: G4RegularNavigation.cc,v 1.7 2007/12/10 16:30:01 gunter Exp $
28// GEANT4 tag $ Name:$
29//
30// class G4RegularNavigation implementation
31//
32// Author: Pedro Arce, May 2007
33//
34// --------------------------------------------------------------------
35
36#include "G4RegularNavigation.hh"
37#include "G4TouchableHistory.hh"
38#include "G4PhantomParameterisation.hh"
39#include "G4Material.hh"
40#include "G4NormalNavigation.hh"
41#include "G4Navigator.hh"
42#include "G4GeometryTolerance.hh"
43
44//------------------------------------------------------------------
45G4RegularNavigation::G4RegularNavigation()
46 : fVerbose(1), fCheck(true)
47{
48 kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
49}
50
51
52//------------------------------------------------------------------
53G4RegularNavigation::~G4RegularNavigation()
54{
55}
56
57
58//------------------------------------------------------------------
59G4double G4RegularNavigation::
60 ComputeStep(const G4ThreeVector& localPoint,
61 const G4ThreeVector& localDirection,
62 const G4double currentProposedStepLength,
63 G4double& newSafety,
64 G4NavigationHistory& history,
65 G4bool& validExitNormal,
66 G4ThreeVector& exitNormal,
67 G4bool& exiting,
68 G4bool& entering,
69 G4VPhysicalVolume *(*pBlockedPhysical),
70 G4int& blockedReplicaNo)
71{
72 // This method is never called because to be called the daughter has to be
73 // a regular structure. This would only happen if the track is in the mother
74 // of voxels volume. But the voxels fill completely their mother, so when a
75 // track enters the mother it automatically enters a voxel. Only precision
76 // problems would make this method to be called
77
78 G4ThreeVector globalPoint =
79 history.GetTopTransform().Inverse().TransformPoint(localPoint);
80 G4ThreeVector globalDirection =
81 history.GetTopTransform().Inverse().TransformAxis(localDirection);
82
83 G4ThreeVector localPoint2 = localPoint; // take away constantness
84
85 LevelLocate( history, *pBlockedPhysical, blockedReplicaNo,
86 globalPoint, &globalDirection, true, localPoint2 );
87
88
89 // Get in which voxel it is
90 //
91 G4VPhysicalVolume *motherPhysical, *daughterPhysical;
92 G4LogicalVolume *motherLogical;
93 motherPhysical = history.GetTopVolume();
94 motherLogical = motherPhysical->GetLogicalVolume();
95 daughterPhysical = motherLogical->GetDaughter(0);
96
97 G4PhantomParameterisation * daughterParam =
98 (G4PhantomParameterisation*)(daughterPhysical->GetParameterisation());
99 G4int copyNo = daughterParam ->GetReplicaNo(localPoint,localDirection);
100
101 G4ThreeVector voxelTranslation = daughterParam->GetTranslation( copyNo );
102 G4ThreeVector daughterPoint = localPoint - voxelTranslation;
103
104
105 // Compute step in voxel
106 //
107 return fnormalNav->ComputeStep(daughterPoint,
108 localDirection,
109 currentProposedStepLength,
110 newSafety,
111 history,
112 validExitNormal,
113 exitNormal,
114 exiting,
115 entering,
116 pBlockedPhysical,
117 blockedReplicaNo);
118}
119
120
121//------------------------------------------------------------------
122G4double G4RegularNavigation::ComputeStepSkippingEqualMaterials(
123 G4ThreeVector localPoint,
124 const G4ThreeVector& localDirection,
125 const G4double currentProposedStepLength,
126 G4double& newSafety,
127 G4NavigationHistory& history,
128 G4bool& validExitNormal,
129 G4ThreeVector& exitNormal,
130 G4bool& exiting,
131 G4bool& entering,
132 G4VPhysicalVolume *(*pBlockedPhysical),
133 G4int& blockedReplicaNo,
134 G4VPhysicalVolume* pCurrentPhysical)
135{
136 G4PhantomParameterisation *param =
137 (G4PhantomParameterisation*)(pCurrentPhysical->GetParameterisation());
138
139 if( !param->SkipEqualMaterials() ) {
140 return fnormalNav->ComputeStep(localPoint,
141 localDirection,
142 currentProposedStepLength,
143 newSafety,
144 history,
145 validExitNormal,
146 exitNormal,
147 exiting,
148 entering,
149 pBlockedPhysical,
150 blockedReplicaNo);
151 }
152
153
154 G4double ourStep = 0.;
155
156 // To get replica No: transform local point to the reference system of the
157 // param container volume
158 //
159 G4int ide = history.GetDepth();
160 G4ThreeVector containerPoint = history.GetTransform(ide).Inverse().TransformPoint(localPoint);
161
162 // Point in global frame
163 //
164 containerPoint = history.GetTransform(ide).Inverse().TransformPoint(localPoint);
165
166 // Point in voxel parent volume frame
167 //
168 containerPoint = history.GetTransform(ide-1).TransformPoint(containerPoint);
169
170 // Store previous voxel translation to move localPoint by the difference
171 // with the new one
172 //
173 G4ThreeVector prevVoxelTranslation = containerPoint - localPoint;
174
175 // Do not use the expression below: There are cases where the
176 // fLastLocatedPointLocal does not give the correct answer
177 // (particle reaching a wall and bounced back, particle travelling through
178 // the wall that is deviated in an step, ...; these are pathological cases
179 // that give wrong answers in G4PhantomParameterisation::GetReplicaNo()
180 //
181 // G4ThreeVector prevVoxelTranslation = param->GetTranslation( copyNo );
182
183 G4int copyNo = param->GetReplicaNo(containerPoint,localDirection);
184
185 G4Material* currentMate = param->ComputeMaterial( copyNo, 0, 0 );
186 G4VSolid* voxelBox = pCurrentPhysical->GetLogicalVolume()->GetSolid();
187
188 G4VSolid* containerSolid = param->GetContainerSolid();
189 G4Material* nextMate;
190 G4bool bLocatedOnEdge = false;
191 G4bool bFirstStep = true;
192 G4double newStep;
193 G4double totalNewStep = 0.;
194
195 // Loop while same material is found
196 //
197 for( ;; )
198 {
199 newStep = voxelBox->DistanceToOut( localPoint, localDirection );
200 if( (bFirstStep) && (newStep < currentProposedStepLength) )
201 {
202 exiting = true;
203 }
204 bFirstStep = false;
205
206 newStep += kCarTolerance; // Avoid precision problems
207 ourStep += newStep;
208 totalNewStep += newStep;
209
210 // Physical process is limiting the step, don't continue
211 //
212 if(std::fabs(totalNewStep-currentProposedStepLength) < kCarTolerance)
213 {
214 return currentProposedStepLength;
215 }
216
217 // Move container point until wall of voxel
218 //
219 containerPoint += newStep*localDirection;
220 if( containerSolid->Inside( containerPoint ) != kInside )
221 {
222 bLocatedOnEdge = true; // NEEDED??
223 break;
224 }
225
226 // Get copyNo and translation of new voxel
227 //
228 copyNo = param->GetReplicaNo(containerPoint,localDirection);
229 G4ThreeVector voxelTranslation = param->GetTranslation( copyNo );
230
231 // Move local point until wall of voxel and then put it in the new voxel
232 // local coordinates
233 //
234 localPoint += newStep*localDirection;
235 localPoint += prevVoxelTranslation - voxelTranslation;
236
237 prevVoxelTranslation = voxelTranslation;
238
239 // Check if material of next voxel is the same as that of the current voxel
240 nextMate = param->ComputeMaterial( copyNo, 0, 0 );
241
242 if( currentMate != nextMate ) { break; }
243 }
244
245 return ourStep;
246}
247
248
249//------------------------------------------------------------------
250G4double
251G4RegularNavigation::ComputeSafety(const G4ThreeVector& localPoint,
252 const G4NavigationHistory& history,
253 const G4double pMaxLength)
254{
255 // This method is never called because to be called the daughter has to be a
256 // regular structure. This would only happen if the track is in the mother of
257 // voxels volume. But the voxels fill completely their mother, so when a
258 // track enters the mother it automatically enters a voxel. Only precision
259 // problems would make this method to be called
260
261 // Compute step in voxel
262 //
263 return fnormalNav->ComputeSafety(localPoint,
264 history,
265 pMaxLength );
266}
267
268
269//------------------------------------------------------------------
270G4bool
271G4RegularNavigation::LevelLocate( G4NavigationHistory& history,
272 const G4VPhysicalVolume* ,
273 const G4int ,
274 const G4ThreeVector& globalPoint,
275 const G4ThreeVector* globalDirection,
276 const G4bool pLocatedOnEdge,
277 G4ThreeVector& localPoint )
278{
279 G4SmartVoxelHeader *motherVoxelHeader;
280 G4VPhysicalVolume *motherPhysical, *pPhysical;
281 G4PhantomParameterisation *pParam;
282 G4LogicalVolume *motherLogical;
283 G4VSolid *pSolid;
284 G4ThreeVector localDir;
285 G4int replicaNo;
286
287 motherPhysical = history.GetTopVolume();
288 motherLogical = motherPhysical->GetLogicalVolume();
289 motherVoxelHeader = motherLogical->GetVoxelHeader();
290
291 pPhysical = motherLogical->GetDaughter(0);
292 pParam = (G4PhantomParameterisation*)(pPhysical->GetParameterisation());
293
294 pSolid = pParam->GetContainerSolid();
295
296 // Save parent history in touchable history
297 // ... for use as parent t-h in ComputeMaterial method of param
298 //
299 G4TouchableHistory parentTouchable( history );
300
301 // Get local direction
302 //
303 if( globalDirection )
304 {
305 localDir = history.GetTopTransform().TransformAxis(*globalDirection);
306 }
307 else
308 {
309 localDir = G4ThreeVector(0.,0.,0.);
310 }
311
312 // Check that track is not on the surface and check that track is not
313 // exiting the voxel parent volume
314 //
315 if ( !G4AuxiliaryNavServices::CheckPointOnSurface(pSolid, localPoint,
316 globalDirection, history.GetTopTransform(), pLocatedOnEdge)
317 || G4AuxiliaryNavServices::CheckPointExiting(pSolid, localPoint,
318 globalDirection, history.GetTopTransform() ) )
319 {
320 }
321 else
322 {
323 // Enter this daughter
324 //
325 replicaNo = pParam->GetReplicaNo( localPoint, localDir );
326
327 if( replicaNo < 0 || replicaNo >= G4int(pParam->GetNoVoxel()) )
328 {
329 return false;
330 }
331
332 // Set the correct copy number in physical
333 //
334 pPhysical->SetCopyNo(replicaNo);
335 pParam->ComputeTransformation(replicaNo,pPhysical);
336
337 history.NewLevel(pPhysical, kParameterised, replicaNo );
338 localPoint = history.GetTopTransform().TransformPoint(globalPoint);
339
340 // Set the correct solid and material in Logical Volume
341 //
342 G4LogicalVolume *pLogical = pPhysical->GetLogicalVolume();
343
344 pLogical->UpdateMaterial(pParam->ComputeMaterial(replicaNo,
345 pPhysical, &parentTouchable) );
346 return true;
347 }
348
349 return false;
350}
Note: See TracBrowser for help on using the repository browser.