| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 // class G4RegularNavigation implementation
27 //
28 // Author: Pedro Arce, May 2007
29 //
30 // -------------------------------------------
31
32 #include "G4RegularNavigation.hh"
33 #include "G4TouchableHistory.hh"
34 #include "G4PhantomParameterisation.hh"
35 #include "G4Material.hh"
36 #include "G4NormalNavigation.hh"
37 #include "G4Navigator.hh"
38 #include "G4GeometryTolerance.hh"
39 #include "G4RegularNavigationHelper.hh"
40
41 //--------------------------------------------
42 G4RegularNavigation::G4RegularNavigation()
43 {
44 kCarTolerance = G4GeometryTolerance::GetInst
45 fMinStep = 101*kCarTolerance;
46 }
47
48
49 //--------------------------------------------
50 G4RegularNavigation::~G4RegularNavigation() =
51
52
53 //--------------------------------------------
54 G4double G4RegularNavigation::
55 ComputeStep(const G4ThreeV
56 const G4ThreeV
57 const G4double
58 G4double
59 G4Naviga
60 G4bool&
61 G4ThreeV
62 G4bool&
63 G4bool&
64 G4VPhysi
65 G4int& b
66 {
67 // This method is never called because to be
68 // a regular structure. This would only happ
69 // of voxels volume. But the voxels fill com
70 // track enters the mother it automatically
71 // problems would make this method to be cal
72
73 G4ThreeVector globalPoint =
74 history.GetTopTransform().InverseTransform
75 G4ThreeVector globalDirection =
76 history.GetTopTransform().InverseTran
77
78 G4ThreeVector localPoint2 = localPoint; // t
79
80 LevelLocate( history, *pBlockedPhysical, blo
81 globalPoint, &globalDirection,
82
83
84 // Get in which voxel it is
85 //
86 G4VPhysicalVolume *motherPhysical, *daughter
87 G4LogicalVolume *motherLogical;
88 motherPhysical = history.GetTopVolume();
89 motherLogical = motherPhysical->GetLogicalVo
90 daughterPhysical = motherLogical->GetDaughte
91
92 auto daughterParam =
93 (G4PhantomParameterisation*)(daughterPh
94 G4int copyNo = daughterParam ->GetReplicaNo(
95
96 G4ThreeVector voxelTranslation = daughterPar
97 G4ThreeVector daughterPoint = localPoint - v
98
99
100 // Compute step in voxel
101 //
102 return fnormalNav->ComputeStep(daughterPoint
103 localDirectio
104 currentPropos
105 newSafety,
106 history,
107 validExitNorm
108 exitNormal,
109 exiting,
110 entering,
111 pBlockedPhysi
112 blockedReplic
113 }
114
115
116 //--------------------------------------------
117 G4double G4RegularNavigation::ComputeStepSkipp
118 G4ThreeV
119 const G4ThreeV
120 const G4double
121 G4double& newS
122 G4NavigationHi
123 G4bool& validE
124 G4ThreeVector&
125 G4bool& exitin
126 G4bool& enteri
127 G4VPhysicalVol
128 G4int& blocked
129 G4VPhysicalVol
130 {
131 G4RegularNavigationHelper::Instance()->Clear
132
133 auto param =
134 (G4PhantomParameterisation*)(pCurrentPhysi
135
136 if( !param->SkipEqualMaterials() )
137 {
138 return fnormalNav->ComputeStep(localPoint,
139 localDirect
140 currentProp
141 newSafety,
142 history,
143 validExitNo
144 exitNormal,
145 exiting,
146 entering,
147 pBlockedPhy
148 blockedRepl
149 }
150
151
152 G4double ourStep = 0.;
153
154 // To get replica No: transform local point
155 // param container volume
156 //
157 auto ide = (G4int)history.GetDepth();
158 G4ThreeVector containerPoint = history.GetTr
159 .InverseTrans
160
161 // Point in global frame
162 //
163 containerPoint = history.GetTransform(ide).I
164
165 // Point in voxel parent volume frame
166 //
167 containerPoint = history.GetTransform(ide-1)
168
169 // Store previous voxel translation to move
170 // with the new one
171 //
172 G4ThreeVector prevVoxelTranslation = contain
173
174 // Do not use the expression below: There ar
175 // fLastLocatedPointLocal does not give the
176 // (particle reaching a wall and bounced bac
177 // the wall that is deviated in an step, ...
178 // that give wrong answers in G4PhantomParam
179 //
180 // G4ThreeVector prevVoxelTranslation = para
181
182 G4int copyNo = param->GetReplicaNo(container
183
184 G4Material* currentMate = param->ComputeMate
185 G4VSolid* voxelBox = pCurrentPhysical->GetLo
186
187 G4VSolid* containerSolid = param->GetContain
188 G4Material* nextMate;
189 G4bool bFirstStep = true;
190 G4double newStep;
191 G4double totalNewStep = 0.;
192
193 // Loop while same material is found
194 //
195 //
196 fNumberZeroSteps = 0;
197 for( G4int ii = 0; ii < fNoStepsAllowed+1; +
198 {
199 if( ii == fNoStepsAllowed ) {
200 // Must kill this stuck track
201 //
202 G4ThreeVector pGlobalpoint = history.Get
203 .InverseTra
204 std::ostringstream message;
205 message << "G4RegularNavigation::Compute
206 << "Stuck Track: potential geometry or
207 << G4endl
208 << " Track stuck, moving for mo
209 << ii << " steps" << G4endl
210 << "- at point " << pGlobalpoint << G4
211 << " local direction: " << loca
212 G4Exception("G4RegularNavigation::Comput
213 "GeomRegNav1001",
214 EventMustBeAborted,
215 message);
216 }
217 newStep = voxelBox->DistanceToOut( localPo
218 fLastStepWasZero = (newStep<fMinStep);
219 if( fLastStepWasZero )
220 {
221 ++fNumberZeroSteps;
222 #ifdef G4DEBUG_NAVIGATION
223 if( fNumberZeroSteps > 1 )
224 {
225 G4ThreeVector pGlobalpoint = history.G
226 .InverseT
227 std::ostringstream message;
228 message.precision(16);
229 message << "G4RegularNavigation::Compute
230 << fNumberZeroSteps
231 << ", at " << pGlobalpoint
232 << ", nav-comp-step calls # " << i
233 << ", Step= " << newStep;
234 G4Exception("G4RegularNavigation::
235 "GeomRegNav1002", Just
236 "Potential overlap in
237 }
238 #endif
239 if( fNumberZeroSteps > fActionThreshold_
240 {
241 // Act to recover this stuck track. Pu
242 //
243 newStep = std::min(101*kCarTolerance*s
244 #ifdef G4DEBUG_NAVIGATION
245 G4ThreeVector pGlobalpoint = history.G
246 .Invers
247 std::ostringstream message;
248 message.precision(16);
249 message << "Track stuck or not moving.
250 << " Track stuc
251 << fNumberZeroSteps << "
252 << "- at point " << pGlo
253 << " (local point " << l
254 << " local direct
255 << " Potential
256 << G4endl
257 << " Trying pus
258 G4Exception("G4RegularNavigation
259 "GeomRegNav1003", Ju
260 "Potential overlap i
261 #endif
262 }
263 if( fNumberZeroSteps > fAbandonThreshold
264 {
265 // Must kill this stuck track
266 //
267 G4ThreeVector pGlobalpoint = history.G
268 .Inv
269 std::ostringstream message;
270 message << "G4RegularNavigation::Compu
271 << "Stuck Track: potential geometry
272 << G4endl
273 << " Track stuck, not moving
274 << fNumberZeroSteps << " steps" << G
275 << "- at point " << pGlobalpoint <<
276 << " local direction: " << lo
277 G4Exception("G4RegularNavigation::Comp
278 "GeomRegNav1004",
279 EventMustBeAborted,
280 message);
281 }
282 }
283 else
284 {
285 // reset the zero step counter when a no
286 fNumberZeroSteps = 0;
287 }
288 if( (bFirstStep) && (newStep < currentProp
289 {
290 exiting = true;
291 }
292 bFirstStep = false;
293
294 newStep += kCarTolerance; // Avoid preci
295 ourStep += newStep;
296 totalNewStep += newStep;
297
298 // Physical process is limiting the step,
299 //
300 if(std::fabs(totalNewStep-currentProposedS
301 {
302 return currentProposedStepLength;
303 }
304 if(totalNewStep > currentProposedStepLengt
305 {
306 G4RegularNavigationHelper::Instance()->
307 AddStepLength(copyNo, newStep-totalNew
308 return currentProposedStepLength;
309 }
310 G4RegularNavigationHelper::Instance()->Add
311
312 // Move container point until wall of voxe
313 //
314 containerPoint += newStep*localDirection;
315 if( containerSolid->Inside( containerPoint
316 {
317 break;
318 }
319
320 // Get copyNo and translation of new voxel
321 //
322 copyNo = param->GetReplicaNo(containerPoin
323 G4ThreeVector voxelTranslation = param->Ge
324
325 // Move local point until wall of voxel an
326 // local coordinates
327 //
328 localPoint += newStep*localDirection;
329 localPoint += prevVoxelTranslation - voxel
330
331 prevVoxelTranslation = voxelTranslation;
332
333 // Check if material of next voxel is the
334 nextMate = param->ComputeMaterial( copyNo,
335
336 if( currentMate != nextMate ) { break; }
337 }
338
339 return ourStep;
340 }
341
342
343 //--------------------------------------------
344 G4double
345 G4RegularNavigation::ComputeSafety(const G4Thr
346 const G4Nav
347 const G4dou
348 {
349 // This method is never called because to be
350 // regular structure. This would only happen
351 // voxels volume. But the voxels fill comple
352 // track enters the mother it automatically
353 // problems would make this method to be cal
354
355 // Compute step in voxel
356 //
357 return fnormalNav->ComputeSafety(localPoint,
358 history,
359 pMaxLength
360 }
361
362
363 //--------------------------------------------
364 G4bool
365 G4RegularNavigation::LevelLocate( G4Navigation
366 const G4VPhy
367 const G4int
368 const G4Thre
369 const G4Thre
370 const G4bool
371 G4ThreeVecto
372 {
373 G4VPhysicalVolume *motherPhysical, *pPhysica
374 G4PhantomParameterisation *pParam;
375 G4LogicalVolume *motherLogical;
376 G4ThreeVector localDir;
377 G4int replicaNo;
378
379 motherPhysical = history.GetTopVolume();
380 motherLogical = motherPhysical->GetLogicalVo
381
382 pPhysical = motherLogical->GetDaughter(0);
383 pParam = (G4PhantomParameterisation*)(pPhysi
384
385 // Save parent history in touchable history
386 // ... for use as parent t-h in ComputeMater
387 //
388 G4TouchableHistory parentTouchable( history
389
390 // Get local direction
391 //
392 if( globalDirection != nullptr )
393 {
394 localDir = history.GetTopTransform().Trans
395 }
396 else
397 {
398 localDir = G4ThreeVector(0.,0.,0.);
399 }
400
401 // Enter this daughter
402 //
403 replicaNo = pParam->GetReplicaNo( localPoint
404
405 if( replicaNo < 0 || replicaNo >= G4int(pPar
406 {
407 return false;
408 }
409
410 // Set the correct copy number in physical
411 //
412 pPhysical->SetCopyNo(replicaNo);
413 pParam->ComputeTransformation(replicaNo,pPhy
414
415 history.NewLevel(pPhysical, kParameterised,
416 localPoint = history.GetTopTransform().Trans
417
418 // Set the correct solid and material in Log
419 //
420 G4LogicalVolume *pLogical = pPhysical->GetLo
421
422 pLogical->UpdateMaterial(pParam->ComputeMate
423 pPhysical, &parentT
424 return true;
425 }
426