| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // G4ParallelGeometriesLimiterProcess << 26 //
27 // ------------------------------------------- << 27 // $Id: G4BiasingProcessLimiterForParallelGeometries.cc $
>> 28 //
>> 29 //
28 30
29 #include "G4ios.hh" 31 #include "G4ios.hh"
30 #include "G4ParallelGeometriesLimiterProcess.h 32 #include "G4ParallelGeometriesLimiterProcess.hh"
31 #include "G4BiasingProcessSharedData.hh" 33 #include "G4BiasingProcessSharedData.hh"
32 #include "G4ProcessManager.hh" 34 #include "G4ProcessManager.hh"
33 #include "G4TransportationManager.hh" 35 #include "G4TransportationManager.hh"
34 #include "G4PathFinder.hh" 36 #include "G4PathFinder.hh"
35 #include "G4FieldTrackUpdator.hh" 37 #include "G4FieldTrackUpdator.hh"
36 38
37 #include "G4SystemOfUnits.hh" 39 #include "G4SystemOfUnits.hh"
38 40
39 G4ParallelGeometriesLimiterProcess:: << 41 G4ParallelGeometriesLimiterProcess::G4ParallelGeometriesLimiterProcess(const G4String& processName) :
40 G4ParallelGeometriesLimiterProcess(const G4Str << 42 G4VProcess(processName, fParallel),
41 : G4VProcess(processName, fParallel) << 43 fParallelWorldSafety( 0.0 ),
>> 44 fIsTrackingTime ( false ),
>> 45 fFieldTrack ( '0' )
42 { 46 {
43 // -- Process Sub Type ? << 47 // -- Process Sub Type ? ยงยง
44 48
45 fPathFinder = G4PathFinder::GetInstance(); << 49 fPathFinder = G4PathFinder::GetInstance();
46 fTransportationManager = G4TransportationMan 50 fTransportationManager = G4TransportationManager::GetTransportationManager();
47 } 51 }
48 52
>> 53
49 // ---------------------------- 54 // ----------------------------
50 // -- Add/Remove world volumes: 55 // -- Add/Remove world volumes:
51 // ---------------------------- 56 // ----------------------------
52 void G4ParallelGeometriesLimiterProcess:: << 57 void G4ParallelGeometriesLimiterProcess::AddParallelWorld(const G4String& parallelWorldName)
53 AddParallelWorld(const G4String& parallelWorld <<
54 { 58 {
>> 59
55 // -- Refuse adding parallel geometry during 60 // -- Refuse adding parallel geometry during tracking time:
56 if (fIsTrackingTime) 61 if (fIsTrackingTime)
57 { <<
58 G4ExceptionDescription ed; <<
59 ed << "G4ParallelGeometriesLimiterProcess <<
60 << "': adding a parallel world volume a <<
61 << G4endl; <<
62 G4Exception("G4ParallelGeometriesLimiterPr <<
63 "BIAS.GEN.21", JustWarning, ed <<
64 return; <<
65 } <<
66 else <<
67 { <<
68 G4VPhysicalVolume* newWorld = fTransportat <<
69 <<
70 // -- Fatal exception if requested world d <<
71 if (newWorld == nullptr) <<
72 { <<
73 G4ExceptionDescription tellWhatIsWrong; <<
74 tellWhatIsWrong << "Volume `" << parall <<
75 << "' is not a parallel <<
76 << G4endl; <<
77 G4Exception("G4ParallelGeometriesLimiter <<
78 "BIAS.GEN.22", FatalExceptio <<
79 } <<
80 <<
81 // -- Protection against adding the mass g <<
82 if ( newWorld == fTransportationManager-> <<
83 { 62 {
84 G4ExceptionDescription ed; 63 G4ExceptionDescription ed;
85 ed << "G4ParallelGeometriesLimiterProces 64 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
86 << "': trying to add the world volume << 65 << "': adding a parallel world volume at tracking time is not allowed." << G4endl;
87 << G4endl; << 66 G4Exception("G4ParallelGeometriesLimiterProcess::AddParallelWorld(const G4String& parallelWorldName)",
88 G4Exception("G4ParallelGeometriesLimiter << 67 "BIAS.GEN.21",
89 "BIAS.GEN.23", JustWarning, << 68 JustWarning, ed,
>> 69 "Call ignored.");
90 return; 70 return;
91 } 71 }
92 << 72
93 // -- Add parallel world, taking care it i << 73 else
94 G4bool isNew = true; << 74
95 for ( auto knownWorld : fParallelWorlds ) <<
96 { <<
97 if ( knownWorld == newWorld ) { isNew = <<
98 } <<
99 if ( isNew ) <<
100 { <<
101 fParallelWorlds.push_back( newWorld ); <<
102 } <<
103 else <<
104 { 75 {
105 G4ExceptionDescription ed; << 76 G4VPhysicalVolume* newWorld = fTransportationManager->IsWorldExisting( parallelWorldName );
106 ed << "G4ParallelGeometriesLimiterProces << 77
107 << "': trying to re-add the parallel << 78 // -- Fatal exception if requested world does not exist:
108 << parallelWorldName << "'." << G4end << 79 if (newWorld == 0)
109 G4Exception("G4ParallelGeometriesLimiter << 80 {
110 "BIAS.GEN.24", JustWarning, << 81 G4ExceptionDescription tellWhatIsWrong;
111 return; << 82 tellWhatIsWrong << "Volume `" << parallelWorldName
>> 83 << "' is not a parallel world nor the mass world volume."
>> 84 << G4endl;
>> 85 G4Exception("G4ParallelGeometriesLimiterProcess::SetWorldVolume(const G4String)",
>> 86 "BIAS.GEN.22",
>> 87 FatalException,
>> 88 tellWhatIsWrong);
>> 89 }
>> 90
>> 91 // -- Protection against adding the mass geometry world as parallel world:
>> 92 if ( newWorld == fTransportationManager->GetNavigatorForTracking()->GetWorldVolume() )
>> 93 {
>> 94 G4ExceptionDescription ed;
>> 95 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
>> 96 << "': trying to add the world volume for tracking as a parallel world." << G4endl;
>> 97 G4Exception("G4ParallelGeometriesLimiterProcess::AddParallelWorld(const G4String& parallelWorldName)",
>> 98 "BIAS.GEN.23",
>> 99 JustWarning, ed,
>> 100 "Call ignored.");
>> 101 return;
>> 102 }
>> 103
>> 104 // -- Add parallel world, taking care it is not in the list yet:
>> 105 G4bool isNew = true;
>> 106 for ( auto knownWorld : fParallelWorlds )
>> 107 {
>> 108 if ( knownWorld == newWorld ) isNew = false;
>> 109 }
>> 110 if ( isNew ) fParallelWorlds.push_back( newWorld );
>> 111 else
>> 112 {
>> 113 G4ExceptionDescription ed;
>> 114 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
>> 115 << "': trying to re-add the parallel world volume `" << parallelWorldName << "'." << G4endl;
>> 116 G4Exception("G4ParallelGeometriesLimiterProcess::AddParallelWorld(const G4String& parallelWorldName)",
>> 117 "BIAS.GEN.24",
>> 118 JustWarning, ed,
>> 119 "Call ignored.");
>> 120 return;
>> 121 }
112 } 122 }
113 } << 123
114 } 124 }
115 125
116 void G4ParallelGeometriesLimiterProcess:: << 126
117 RemoveParallelWorld(const G4String& parallelWo << 127 void G4ParallelGeometriesLimiterProcess::RemoveParallelWorld(const G4String& parallelWorldName)
118 { 128 {
>> 129
119 // -- Refuse refuse removing parallel geomet 130 // -- Refuse refuse removing parallel geometry during tracking time:
120 if (fIsTrackingTime) 131 if (fIsTrackingTime)
121 { <<
122 G4ExceptionDescription ed; <<
123 ed << "G4ParallelGeometriesLimiterProcess <<
124 << "': removing a parallel world volume <<
125 << G4endl; <<
126 G4Exception("G4ParallelGeometriesLimiterPr <<
127 "BIAS.GEN.25", JustWarning, ed <<
128 return; <<
129 } <<
130 else <<
131 { <<
132 G4VPhysicalVolume* newWorld = fTransportat <<
133 if (newWorld == nullptr) <<
134 { 132 {
135 G4ExceptionDescription ed; 133 G4ExceptionDescription ed;
136 ed << "G4ParallelGeometriesLimiterProces 134 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
137 << "': trying to remove an inexisting << 135 << "': removing a parallel world volume at tracking time is not allowed." << G4endl;
138 << parallelWorldName << "'." << G4end << 136 G4Exception("G4ParallelGeometriesLimiterProcess::RemoveParallelWorld(const G4String& parallelWorldName)",
139 G4Exception("G4ParallelGeometriesLimiter << 137 "BIAS.GEN.25",
140 "BIAS.GEN.26", JustWarning, << 138 JustWarning, ed,
>> 139 "Call ignored.");
141 return; 140 return;
142 } 141 }
143 << 142
144 // -- get position of world volume in list << 143 else
145 std::size_t iWorld = 0; << 144
146 for ( auto knownWorld : fParallelWorlds ) <<
147 { <<
148 if ( knownWorld == newWorld ) break; <<
149 ++iWorld; <<
150 } <<
151 <<
152 if ( iWorld == fParallelWorlds.size() ) <<
153 { 145 {
154 G4ExceptionDescription ed; << 146 G4VPhysicalVolume* newWorld = fTransportationManager->IsWorldExisting( parallelWorldName );
155 ed << "G4ParallelGeometriesLimiterProces << 147
156 << "': trying to remove an non-regist << 148 if (newWorld == 0)
157 << parallelWorldName << "'." << G4end << 149 {
158 G4Exception("G4ParallelGeometriesLimiter << 150
159 "BIAS.GEN.27", JustWarning, << 151 G4ExceptionDescription ed;
160 return; << 152 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
>> 153 << "': trying to remove an inexisting parallel world '" << parallelWorldName << "'." << G4endl;
>> 154 G4Exception("G4ParallelGeometriesLimiterProcess::RemoveParallelWorld(const G4String& parallelWorldName)",
>> 155 "BIAS.GEN.26",
>> 156 JustWarning, ed,
>> 157 "Call ignored.");
>> 158 return;
>> 159 }
>> 160
>> 161 // -- get position of world volume in list:
>> 162 size_t iWorld = 0;
>> 163 for ( auto knownWorld : fParallelWorlds )
>> 164 {
>> 165 if ( knownWorld == newWorld ) break;
>> 166 iWorld++;
>> 167 }
>> 168
>> 169 if ( iWorld == fParallelWorlds.size() )
>> 170 {
>> 171 G4ExceptionDescription ed;
>> 172 ed << "G4ParallelGeometriesLimiterProcess `" << GetProcessName()
>> 173 << "': trying to remove an non-registerered parallel world '" << parallelWorldName << "'." << G4endl;
>> 174 G4Exception("G4ParallelGeometriesLimiterProcess::RemoveParallelWorld(const G4String& parallelWorldName)",
>> 175 "BIAS.GEN.27",
>> 176 JustWarning, ed,
>> 177 "Call ignored.");
>> 178 return;
>> 179 }
>> 180
>> 181 // -- remove from vector:
>> 182 fParallelWorlds.erase( fParallelWorlds.begin() + iWorld );
>> 183
161 } 184 }
162 // -- remove from vector: << 185
163 fParallelWorlds.erase( fParallelWorlds.beg << 186
164 } << 187
165 } 188 }
166 189
>> 190
167 // -------------------- 191 // --------------------
168 // Start/End tracking: 192 // Start/End tracking:
169 // -------------------- 193 // --------------------
170 void G4ParallelGeometriesLimiterProcess::Start 194 void G4ParallelGeometriesLimiterProcess::StartTracking(G4Track* track)
171 { 195 {
172 fIsTrackingTime = true; 196 fIsTrackingTime = true;
173 197
174 // -- fetch the navigators, their indeces, a 198 // -- fetch the navigators, their indeces, and activate:
175 fParallelWorldNavigators.clear(); << 199 fParallelWorldNavigators .clear();
176 fParallelWorldNavigatorIndeces.clear(); 200 fParallelWorldNavigatorIndeces.clear();
177 fParallelWorldSafeties.clear(); << 201 fParallelWorldSafeties .clear();
178 fParallelWorldIsLimiting.clear(); << 202 fParallelWorldIsLimiting .clear();
179 fParallelWorldWasLimiting.clear(); << 203 fParallelWorldWasLimiting .clear();
180 fCurrentVolumes.clear(); << 204 fCurrentVolumes .clear();
181 fPreviousVolumes.clear(); << 205 fPreviousVolumes .clear();
182 for ( auto parallelWorld : fParallelWorlds ) 206 for ( auto parallelWorld : fParallelWorlds )
183 { << 207 {
184 fParallelWorldNavigators.push_back( fTrans << 208 fParallelWorldNavigators .push_back( fTransportationManager-> GetNavigator( parallelWorld ) );
185 fParallelWorldNavigatorIndeces.push_back( << 209 fParallelWorldNavigatorIndeces.push_back( fTransportationManager->ActivateNavigator( fParallelWorldNavigators.back() ) );
186 fParallelWorldSafeties.push_back( 0.0 ); << 210 fParallelWorldSafeties .push_back( 0.0 );
187 fParallelWorldIsLimiting.push_back( false << 211 fParallelWorldIsLimiting .push_back( false );
188 fParallelWorldWasLimiting.push_back( false << 212 fParallelWorldWasLimiting .push_back( false );
189 } << 213 }
190 214
191 fPathFinder->PrepareNewTrack( track->GetPosi 215 fPathFinder->PrepareNewTrack( track->GetPosition(), track->GetMomentumDirection() );
192 // -- Does it work at this level, after "Pre << 216 // -- ยงยง does it work at this level, after "PrepareNewTrack" above ?
193 for ( auto navigatorIndex : fParallelWorldNa 217 for ( auto navigatorIndex : fParallelWorldNavigatorIndeces )
194 { << 218 {
195 fPreviousVolumes.push_back( nullptr ); << 219 fPreviousVolumes.push_back( nullptr );
196 fCurrentVolumes .push_back( fPathFinder->G << 220 fCurrentVolumes .push_back( fPathFinder->GetLocatedVolume( navigatorIndex ) );
197 } 221 }
198 222
199 // -- will force updating safety: 223 // -- will force updating safety:
200 fParallelWorldSafety = 0.0; 224 fParallelWorldSafety = 0.0;
201 for ( std::size_t i = 0 ; i < fParallelWorld << 225 for ( size_t i = 0 ; i < fParallelWorldNavigatorIndeces.size() ; i++ ) fParallelWorldSafeties[i] = 0.0;
202 { <<
203 fParallelWorldSafeties[i] = 0.0; <<
204 } <<
205 } 226 }
206 227
>> 228
207 void G4ParallelGeometriesLimiterProcess::EndTr 229 void G4ParallelGeometriesLimiterProcess::EndTracking()
208 { 230 {
209 fIsTrackingTime = false; 231 fIsTrackingTime = false;
210 for ( auto parallelWorldNavigator : fParalle 232 for ( auto parallelWorldNavigator : fParallelWorldNavigators )
211 { <<
212 fTransportationManager->DeActivateNavigato 233 fTransportationManager->DeActivateNavigator( parallelWorldNavigator );
213 } <<
214 } 234 }
215 235
216 G4double G4ParallelGeometriesLimiterProcess:: << 236
217 PostStepGetPhysicalInteractionLength(const G4T << 237 G4double G4ParallelGeometriesLimiterProcess::PostStepGetPhysicalInteractionLength(const G4Track&, G4double, G4ForceCondition* condition)
218 G4ForceCo <<
219 { 238 {
>> 239
220 // -- push previous step limitation flags an 240 // -- push previous step limitation flags and volumes:
221 // -- consider switching pointers insteads o << 241 // -- ยงยง consider switching pointers insteads of making copies of std::vector's:
222 fParallelWorldWasLimiting = fParallelWorldIs 242 fParallelWorldWasLimiting = fParallelWorldIsLimiting;
223 fPreviousVolumes = fCurrentVolumes; << 243 fPreviousVolumes = fCurrentVolumes;
224 244
225 // -- update volumes: 245 // -- update volumes:
226 std::size_t i = 0; << 246 size_t i = 0;
227 for ( auto navigatorIndex : fParallelWorldNa << 247 for ( auto navigatorIndex : fParallelWorldNavigatorIndeces ) fCurrentVolumes[i++] = fPathFinder->GetLocatedVolume( navigatorIndex );
228 { <<
229 fCurrentVolumes[i++] = fPathFinder->GetLoc <<
230 } <<
231 248
232 *condition = NotForced; 249 *condition = NotForced;
233 return DBL_MAX; 250 return DBL_MAX;
234 } 251 }
235 252
236 G4double G4ParallelGeometriesLimiterProcess:: << 253
237 AlongStepGetPhysicalInteractionLength(const G4 << 254 G4double G4ParallelGeometriesLimiterProcess::AlongStepGetPhysicalInteractionLength(const G4Track& track,
238 G4double << 255 G4double previousStepSize,
239 G4double << 256 G4double currentMinimumStep,
240 G4double << 257 G4double& proposedSafety,
241 G4GPILSe << 258 G4GPILSelection* selection)
242 { 259 {
>> 260
243 // -- Init: 261 // -- Init:
244 // -- Note that the returnedStep must be phy << 262 // -- Note that the returnedStep must be physically meaningful, even if we return NotCandidateForSelection as condition;
245 // -- even if we return NotCandidateForSelec << 263 // -- the reason is that the stepping manager always takes the smallest alongstep among the returned ones (point related
246 // -- the reason is that the stepping manage <<
247 // -- alongstep among the returned ones (poi <<
248 // -- to geometry step length wrt to true pa 264 // -- to geometry step length wrt to true path length).
249 *selection = NotCandidateForSelection; << 265 *selection = NotCandidateForSelection;
250 G4double returnedStep = DBL_MAX; 266 G4double returnedStep = DBL_MAX;
251 267
252 // -- G4FieldTrack and ELimited: 268 // -- G4FieldTrack and ELimited:
253 static G4ThreadLocal G4FieldTrack* endTrack_ << 269 static G4ThreadLocal G4FieldTrack *endTrack_G4MT_TLS_ = 0 ;
254 if (!endTrack_MT) endTrack_MT = new G4FieldT << 270 if (!endTrack_G4MT_TLS_) endTrack_G4MT_TLS_ = new G4FieldTrack ('0') ;
255 G4FieldTrack& endTrack = *endTrack_MT; << 271 G4FieldTrack &endTrack = *endTrack_G4MT_TLS_;
256 << 272
257 static G4ThreadLocal ELimited* eLimited_MT = << 273 static G4ThreadLocal ELimited *eLimited_G4MT_TLS_ = 0 ;
258 if (!eLimited_MT) eLimited_MT = new ELimited << 274 if (!eLimited_G4MT_TLS_) eLimited_G4MT_TLS_ = new ELimited ;
259 ELimited &eLimited = *eLimited_MT; << 275 ELimited &eLimited = *eLimited_G4MT_TLS_;
>> 276
260 277
261 // ------------------- 278 // -------------------
262 // -- Update safeties: 279 // -- Update safeties:
263 // ------------------- 280 // -------------------
264 if ( previousStepSize > 0.0 ) 281 if ( previousStepSize > 0.0 )
265 { <<
266 for ( auto& parallelWorldSafety : fParalle <<
267 { 282 {
268 parallelWorldSafety -= previousStepSize; << 283 for ( auto& parallelWorldSafety : fParallelWorldSafeties )
269 if ( parallelWorldSafety < 0. ) { paral << 284 {
270 fParallelWorldSafety = parallelWorldSafe << 285 parallelWorldSafety -= previousStepSize;
271 ? parallelWorldSafe << 286 if ( parallelWorldSafety < 0. ) parallelWorldSafety = 0.0;
>> 287 fParallelWorldSafety = parallelWorldSafety < fParallelWorldSafety ? parallelWorldSafety : fParallelWorldSafety ;
>> 288 }
272 } 289 }
273 } <<
274 290
>> 291
275 // ----------------------------------------- 292 // ------------------------------------------
276 // Determination of the proposed step length 293 // Determination of the proposed step length:
277 // ----------------------------------------- 294 // ------------------------------------------
278 if ( ( currentMinimumStep <= fParallelWorldS << 295 if ( ( currentMinimumStep <= fParallelWorldSafety ) && ( currentMinimumStep > 0. ) )
279 && ( currentMinimumStep > 0. ) ) << 296 {
280 { << 297 // -- No chance to limit the step, as proposed move inside safety
281 // -- No chance to limit the step, as prop << 298
282 << 299 returnedStep = currentMinimumStep;
283 returnedStep = currentMinimumStep; << 300 proposedSafety = fParallelWorldSafety - currentMinimumStep;
284 proposedSafety = fParallelWorldSafety - cu << 301 }
285 } <<
286 else 302 else
287 { << 303 {
288 // -- Proposed move exceeds common safety, << 304 // -- Proposed move exceeds common safety, need to state
289 G4double smallestReturnedStep = -1.0; << 305 G4double smallestReturnedStep = -1.0;
290 ELimited eLimitedForSmallestStep = kDoNot; << 306 ELimited eLimitedForSmallestStep = kDoNot;
291 for ( std::size_t i = 0 ; i < fParallelWor << 307 for ( size_t i = 0 ; i < fParallelWorldNavigatorIndeces.size() ; i++ )
292 { << 308 {
293 // -- Update safety of geometries having << 309 // -- Update safety of geometries having safety smaller than current minimum step
294 if ( currentMinimumStep >= fParallelWor << 310 if ( currentMinimumStep >= fParallelWorldSafeties[i] )
295 { << 311 {
296 G4FieldTrackUpdator::Update(&fFieldTra << 312 G4FieldTrackUpdator::Update(&fFieldTrack, &track);
297 G4double tmpReturnedStep = fPathFinder << 313 G4double tmpReturnedStep = fPathFinder->ComputeStep(fFieldTrack,
298 currentMi << 314 currentMinimumStep,
299 fParallel << 315 fParallelWorldNavigatorIndeces[i],
300 track.Get << 316 track.GetCurrentStepNumber(),
301 fParallel << 317 fParallelWorldSafeties[i],
302 eLimited, << 318 eLimited,
303 << 319 endTrack,
304 if ( ( smallestReturnedStep < 0.0 ) || << 320 track.GetVolume());
305 { << 321
306 smallestReturnedStep = tmpReturnedSt << 322 if ( ( smallestReturnedStep < 0.0 ) || ( tmpReturnedStep <= smallestReturnedStep ) )
307 eLimitedForSmallestStep = eLimited; << 323 {
308 } << 324 smallestReturnedStep = tmpReturnedStep;
309 << 325 eLimitedForSmallestStep = eLimited;
310 if (eLimited == kDoNot) << 326 }
311 { << 327
312 // -- Step not limited by this geome << 328 if (eLimited == kDoNot)
313 fParallelWorldSafeties[i] = fParalle << 329 {
314 fParallelWorldIsLimiting[i] = false; << 330 // -- Step not limited by this geometry
315 } << 331 fParallelWorldSafeties[i] = fParallelWorldNavigators[i]->ComputeSafety(endTrack.GetPosition());
316 else << 332 fParallelWorldIsLimiting[i] = false;
317 { << 333 }
318 fParallelWorldIsLimiting[i] = true; << 334 else
319 } << 335 {
320 } << 336 fParallelWorldIsLimiting[i] = true;
321 << 337 }
322 // -- update with smallest safety: << 338 }
323 fParallelWorldSafety = fParallelWorldSaf << 339
324 ? fParallelWorldSa << 340 // -- update with smallest safety:
325 } << 341 fParallelWorldSafety = fParallelWorldSafeties[i] < fParallelWorldSafety ? fParallelWorldSafeties[i] : fParallelWorldSafety ;
326 << 342 }
327 // -- no geometry limitation among all geo << 343
328 // -- Beware : the returnedStep must be ph << 344 // -- no geometry limitation among all geometries, can return currentMinimumStep (or DBL_MAX):
329 if ( eLimitedForSmallestStep == kDoNot ) << 345 // -- Beware : the returnedStep must be physically meaningful, even if we say "NotCandidateForSelection" !
330 { << 346 if ( eLimitedForSmallestStep == kDoNot )
331 returnedStep = currentMinimumStep; << 347 {
332 } << 348 returnedStep = currentMinimumStep;
333 // -- proposed step length of limiting geo << 349 }
334 if ( eLimitedForSmallestStep == kUnique | << 350 // -- proposed step length of limiting geometry:
335 eLimitedForSmallestStep == kSharedOth << 351 if ( eLimitedForSmallestStep == kUnique ||
336 { << 352 eLimitedForSmallestStep == kSharedOther )
337 *selection = CandidateForSelection; << 353 {
338 returnedStep = smallestReturnedStep; << 354 *selection = CandidateForSelection;
339 } << 355 returnedStep = smallestReturnedStep;
340 else if ( eLimitedForSmallestStep == kShar << 356 }
341 { << 357 else if ( eLimitedForSmallestStep == kSharedTransport)
342 // -- Expand to disable its selection in << 358 {
343 returnedStep = smallestReturnedStep* (1. << 359 returnedStep = smallestReturnedStep* (1.0 + 1.0e-9); // -- Expand to disable its selection in Step Manager comparison
344 } << 360 }
345 << 361
346 // -- and smallest safety among geometries << 362 // -- and smallest safety among geometries:
347 proposedSafety = fParallelWorldSafety ; << 363 proposedSafety = fParallelWorldSafety ;
348 } << 364 }
349 365
350 // -- returns step length, and proposedSafet 366 // -- returns step length, and proposedSafety
351 return returnedStep; 367 return returnedStep;
352 } 368 }
353 369
354 G4VParticleChange* G4ParallelGeometriesLimiter << 370
355 AlongStepDoIt( const G4Track& track, const G4S << 371 G4VParticleChange* G4ParallelGeometriesLimiterProcess::AlongStepDoIt( const G4Track& track,
>> 372 const G4Step& )
356 { 373 {
>> 374
357 fDummyParticleChange.Initialize(track); 375 fDummyParticleChange.Initialize(track);
358 return &fDummyParticleChange; 376 return &fDummyParticleChange;
359 } 377 }
360 378
361 void G4ParallelGeometriesLimiterProcess:: << 379
362 SetProcessManager(const G4ProcessManager* mgr) << 380 void G4ParallelGeometriesLimiterProcess::SetProcessManager(const G4ProcessManager* mgr)
363 { 381 {
364 G4BiasingProcessSharedData *sharedData(nullp 382 G4BiasingProcessSharedData *sharedData(nullptr);
365 383
366 // -- initialize sharedData pointer: 384 // -- initialize sharedData pointer:
367 if ( G4BiasingProcessSharedData::fSharedData << 385 if ( G4BiasingProcessSharedData::fSharedDataMap.Find(mgr) == G4BiasingProcessSharedData::fSharedDataMap.End() )
368 { << 386 {
369 sharedData = new G4BiasingProcessSharedDat << 387 sharedData = new G4BiasingProcessSharedData( mgr );
370 G4BiasingProcessSharedData::fSharedDataMap << 388 G4BiasingProcessSharedData::fSharedDataMap[mgr] = sharedData;
371 } << 389 }
372 else << 390 else sharedData = G4BiasingProcessSharedData::fSharedDataMap[mgr] ;
373 { << 391
374 sharedData = G4BiasingProcessSharedData:: <<
375 } <<
376 <<
377 // -- add itself to the shared data: 392 // -- add itself to the shared data:
378 if ( sharedData->fParallelGeometriesLimiterP << 393 if ( sharedData->fParallelGeometriesLimiterProcess == nullptr ) sharedData->fParallelGeometriesLimiterProcess = this;
379 { <<
380 sharedData->fParallelGeometriesLimiterProc <<
381 } <<
382 else 394 else
383 { << 395 {
384 G4ExceptionDescription ed; << 396 G4ExceptionDescription ed;
385 ed << " Trying to add more than one G4Para << 397 ed << " Trying to add more than one G4ParallelGeometriesLimiterProcess process to the process manager " << mgr
386 << mgr << 398 << " (process manager for `" << mgr->GetParticleType()->GetParticleName() << "'). Only one is needed. Call ignored." << G4endl;
387 << " (process manager for `" << mgr->Ge << 399 G4Exception(" G4ParallelGeometriesLimiterProcess::SetProcessManager(...)",
388 << "'). Only one is needed. Call ignore << 400 "BIAS.GEN.29",
389 G4Exception(" G4ParallelGeometriesLimiterP << 401 JustWarning,
390 "BIAS.GEN.29", JustWarning, ed << 402 ed);
391 } << 403 }
392 } 404 }
393 405
394 G4int G4ParallelGeometriesLimiterProcess:: << 406
395 GetParallelWorldIndex( const G4VPhysicalVolume << 407 G4int G4ParallelGeometriesLimiterProcess::GetParallelWorldIndex( const G4VPhysicalVolume* parallelWorld ) const
396 { 408 {
397 G4int toReturn = -1; 409 G4int toReturn = -1;
398 G4int iWorld = 0; 410 G4int iWorld = 0;
399 for ( auto world : fParallelWorlds ) 411 for ( auto world : fParallelWorlds )
400 { <<
401 if ( world == parallelWorld ) <<
402 { 412 {
403 toReturn = iWorld; << 413 if ( world == parallelWorld )
404 break; << 414 {
>> 415 toReturn = iWorld;
>> 416 break;
>> 417 }
>> 418 iWorld++;
405 } 419 }
406 ++iWorld; <<
407 } <<
408 return toReturn; 420 return toReturn;
409 } 421 }
410 422
411 G4int G4ParallelGeometriesLimiterProcess:: << 423
412 GetParallelWorldIndex( const G4String& paralle << 424 G4int G4ParallelGeometriesLimiterProcess::GetParallelWorldIndex( G4String parallelWorldName ) const
413 { 425 {
414 G4VPhysicalVolume* aWorld = fTransportationM << 426 G4VPhysicalVolume* aWorld = fTransportationManager->IsWorldExisting( parallelWorldName ); // note aWorld might be nullptr
415 // note aWorld might be nullptr <<
416 return GetParallelWorldIndex( aWorld ); 427 return GetParallelWorldIndex( aWorld );
417 } 428 }
>> 429
418 430