Geant4  10.00.p02
G4ImportanceProcess.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: G4ImportanceProcess.cc 77999 2013-12-02 08:18:01Z gcosmo $
28 //
29 // ----------------------------------------------------------------------
30 // GEANT 4 class source file
31 //
32 // G4ImportanceProcess.cc
33 //
34 // ----------------------------------------------------------------------
35 
36 #include "G4ImportanceProcess.hh"
38 #include "G4GeometryCell.hh"
40 #include "G4VTrackTerminator.hh"
41 #include "G4VIStore.hh"
42 
43 #include "G4Step.hh"
44 #include "G4Navigator.hh"
45 #include "G4VTouchable.hh"
46 #include "G4VPhysicalVolume.hh"
47 #include "G4ParticleChange.hh"
48 #include "G4PathFinder.hh"
50 #include "G4StepPoint.hh"
51 #include "G4FieldTrackUpdator.hh"
52 
53 
55 G4ImportanceProcess(const G4VImportanceAlgorithm &aImportanceAlgorithm,
56  const G4VIStore &aIstore,
57  const G4VTrackTerminator *TrackTerminator,
58  const G4String &aName, G4bool para)
59  : G4VProcess(aName),
60  fParticleChange(new G4ParticleChange),
61  fImportanceAlgorithm(aImportanceAlgorithm),
62  fIStore(aIstore),
63  fPostStepAction(0),
64  fGhostWorldName("NoParallelWorld"),fGhostWorld(0),
65  fGhostNavigator(0), fNavigatorID(-1), fFieldTrack('0'),
66  paraflag(para)
67 
68 {
69  G4cout << G4endl << G4endl << G4endl;
70  G4cout << "G4ImportanceProcess:: Creating " << G4endl;
71  if (TrackTerminator)
72  {
73  fPostStepAction = new G4SamplingPostStepAction(*TrackTerminator);
74  }
75  else
76  {
78  }
79  if (!fParticleChange)
80  {
81  G4Exception("G4ImportanceProcess::G4ImportanceProcess()",
82  "FatalError", FatalException,
83  "Failed allocation of G4ParticleChange !");
84  }
86 
87 
88  fGhostStep = new G4Step();
91 
94 
95  if (verboseLevel>0)
96  {
97  G4cout << GetProcessName() << " is created " << G4endl;
98  }
99 
100  G4cout << "G4ImportanceProcess:: importance process paraflag is: " << paraflag << G4endl;
101 
102 }
103 
105 {
106 
107  delete fPostStepAction;
108  delete fParticleChange;
109  // delete fGhostStep;
110  // delete fGhostWorld;
111  // delete fGhostNavigator;
112 
113 }
114 
115 
116 
117 //------------------------------------------------------
118 //
119 // SetParallelWorld
120 //
121 //------------------------------------------------------
123 SetParallelWorld(G4String parallelWorldName)
124 {
125  G4cout << G4endl << G4endl << G4endl;
126  G4cout << "G4ImportanceProcess:: SetParallelWorld name = " << parallelWorldName << G4endl;
127 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
128 // Get pointers of the parallel world and its navigator
129 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
130  fGhostWorldName = parallelWorldName;
133 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
134 }
135 
136 // void G4ImportanceProcess::
137 // SetParallelWorld(const G4VPhysicalVolume* parallelWorld)
138 // {
139 // //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
140 // // Get pointer of navigator
141 // //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
142 // // G4cout << " G4ImportanceProcess:: Got here 1 " << G4endl;
143 // // fGhostWorldName = parallelWorld->GetName();
144 // // G4cout << " G4ImportanceProcess:: Got here 2 ghostName:" << fGhostWorldName << G4endl;
145 // fGhostWorld = parallelWorld;
146 // G4cout << " G4ImportanceProcess:: Got here 3 " << G4endl;
147 // fGhostNavigator = fTransportationManager->GetNavigator(parallelWorld);
148 // G4cout << " G4ImportanceProcess:: Got here 4 " << G4endl;
149 // //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
150 // }
151 
152 //------------------------------------------------------
153 //
154 // StartTracking
155 //
156 //------------------------------------------------------
158 {
159 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
160 // Activate navigator and get the navigator ID
161 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
162 // G4cout << " G4ParallelWorldScoringProcess::StartTracking" << G4endl;
163 
164  if(paraflag) {
165  if(fGhostNavigator)
167  else
168  {
169  G4Exception("G4ImportanceProcess::StartTracking",
170  "ProcParaWorld000",FatalException,
171  "G4ImportanceProcess is used for tracking without having a parallel world assigned");
172  }
173 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
174 
175 // G4cout << "G4ParallelWorldScoringProcess::StartTracking <<<<<<<<<<<<<<<<<< " << G4endl;
176 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
177 // Let PathFinder initialize
178 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
180 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
181 
182 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
183 // Setup initial touchables for the first step
184 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
189 
190  // Initialize
191  fGhostSafety = -1.;
192  fOnBoundary = false;
193  }
194 
195 }
196 
197 
200  G4double ,
202 {
203 // *condition = Forced;
204 // return kInfinity;
205 
206 // *condition = StronglyForced;
207  *condition = Forced;
208  return DBL_MAX;
209 }
210 
213  const G4Step &aStep)
214 {
215  fParticleChange->Initialize(aTrack);
216 
217  if(paraflag) {
219  //xbug? fOnBoundary = false;
220  CopyStep(aStep);
221 
222  if(fOnBoundary)
223  {
224 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
225 // Locate the point and get new touchable
226 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
227  //?? fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(),
228  //?? step.GetPostStepPoint()->GetMomentumDirection());
230  //AH G4cout << " on boundary " << G4endl;
231 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
232  }
233  else
234  {
235  // Do I need this ??????????????????????????????????????????????????????????
236  // fGhostNavigator->LocateGlobalPointWithinVolume(track.GetPosition());
237  // ?????????????????????????????????????????????????????????????????????????
238 
239  // fPathFinder->ReLocate(track.GetPosition());
240 
241  // reuse the touchable
243  //AH G4cout << " NOT on boundary " << G4endl;
244  }
245 
248 
249 // if ( (aStep.GetPostStepPoint()->GetStepStatus() == fGeomBoundary)
250 // && (aStep.GetStepLength() > kCarTolerance) )
251 // {
252 // if (aTrack.GetTrackStatus()==fStopAndKill)
253 // {
254 // G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
255 // << " StopAndKill track." << G4endl;
256 // }
257 
258 // G4StepPoint *prepoint = aStep.GetPreStepPoint();
259 // G4StepPoint *postpoint = aStep.GetPostStepPoint();
260 // G4GeometryCell prekey(*(prepoint->GetPhysicalVolume()),
261 // prepoint->GetTouchable()->GetReplicaNumber());
262 // G4GeometryCell postkey(*(postpoint->GetPhysicalVolume()),
263 // postpoint->GetTouchable()->GetReplicaNumber());
264 
265 
267  && (aStep.GetStepLength() > kCarTolerance) )
268  {
269  if (aTrack.GetTrackStatus()==fStopAndKill)
270  {
271  G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
272  << " StopAndKill track. on boundary" << G4endl;
273  }
274 
279 
280  //AH
281  /*
282  G4cout << G4endl;
283  G4cout << G4endl;
284  G4cout << " inside parallel importance process " << aTrack.GetCurrentStepNumber() << G4endl;
285  G4cout << G4endl;
286  G4cout << G4endl;
287  G4cout << " prekey: " << fGhostPreStepPoint->GetPhysicalVolume()->GetName() << " replica:"
288  << fGhostPreStepPoint->GetTouchable()->GetReplicaNumber() << G4endl;
289  G4cout << " prekey ISTORE: " << fIStore.GetImportance(prekey) << G4endl;
290  G4cout << " postkey: " << G4endl;
291  G4cout << " postkey ISTORE: " << fIStore.GetImportance(postkey) << G4endl;
292  */
293  //AH
295  Calculate(fIStore.GetImportance(prekey),
296  fIStore.GetImportance(postkey),
297  aTrack.GetWeight());
298  //AH
299  /*
300  G4cout << " prekey weight: " << fIStore.GetImportance(prekey)
301  << " postkey weight: " << fIStore.GetImportance(postkey)
302  << " split weight: " << nw << G4endl;
303  G4cout << " before poststepaction " << G4endl;
304  */
305  //AH
306  fPostStepAction->DoIt(aTrack, fParticleChange, nw);
307  //AH G4cout << " after post step do it " << G4endl;
308  }
309  } else {
310  if ( (aStep.GetPostStepPoint()->GetStepStatus() == fGeomBoundary)
311  && (aStep.GetStepLength() > kCarTolerance) )
312  {
313  //AH G4cout << " inside non-parallel importance process " << G4endl;
314  if (aTrack.GetTrackStatus()==fStopAndKill)
315  {
316  G4cout << "WARNING - G4ImportanceProcess::PostStepDoIt()"
317  << " StopAndKill track. on boundary non-parallel" << G4endl;
318  }
319 
320  G4StepPoint *prepoint = aStep.GetPreStepPoint();
321  G4StepPoint *postpoint = aStep.GetPostStepPoint();
322 
323  G4GeometryCell prekey(*(prepoint->GetPhysicalVolume()),
324  prepoint->GetTouchable()->GetReplicaNumber());
325  G4GeometryCell postkey(*(postpoint->GetPhysicalVolume()),
326  postpoint->GetTouchable()->GetReplicaNumber());
327 
329  Calculate(fIStore.GetImportance(prekey),
330  fIStore.GetImportance(postkey),
331  aTrack.GetWeight());
332  //AH
333  /*
334  G4cout << " prekey weight: " << fIStore.GetImportance(prekey)
335  << " postkey weight: " << fIStore.GetImportance(postkey)
336  << " split weight: " << nw << G4endl;
337  G4cout << " before poststepaction 2 " << G4endl;
338  */
339  //AH
340  fPostStepAction->DoIt(aTrack, fParticleChange, nw);
341  //AH G4cout << " after poststepaction 2 " << G4endl;
342  }
343  }
344  return fParticleChange;
345 }
346 
348 {
350 }
351 
353 {
354  return theProcessName;
355 }
356 
359  const G4Track& track, G4double previousStepSize, G4double currentMinimumStep,
360  G4double& proposedSafety, G4GPILSelection* selection)
361 {
362 
363  //AH G4cout << " along step physical interaction length " << G4endl;
364 
365  if(paraflag) {
366  static G4ThreadLocal G4FieldTrack *endTrack_G4MT_TLS_ = 0 ; if (!endTrack_G4MT_TLS_) endTrack_G4MT_TLS_ = new G4FieldTrack ('0') ; G4FieldTrack &endTrack = *endTrack_G4MT_TLS_;
367  static G4ThreadLocal ELimited *eLimited_G4MT_TLS_ = 0 ; if (!eLimited_G4MT_TLS_) eLimited_G4MT_TLS_ = new ELimited ; ELimited &eLimited = *eLimited_G4MT_TLS_;
368 
369  *selection = NotCandidateForSelection;
370  G4double returnedStep = DBL_MAX;
371 
372  if (previousStepSize > 0.)
373  { fGhostSafety -= previousStepSize; }
374  // else
375  // { fGhostSafety = -1.; }
376  if (fGhostSafety < 0.) fGhostSafety = 0.0;
377 
378  // ------------------------------------------
379  // Determination of the proposed STEP LENGTH:
380  // ------------------------------------------
381  if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
382  {
383  // I have no chance to limit
384  returnedStep = currentMinimumStep;
385  fOnBoundary = false;
386  proposedSafety = fGhostSafety - currentMinimumStep;
387  //AH G4cout << " step not limited, why? " << G4endl;
388  }
389  else // (currentMinimumStep > fGhostSafety: I may limit the Step)
390  {
392  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
393  // ComputeStep
394  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
395  returnedStep
396  = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
397  track.GetCurrentStepNumber(),fGhostSafety,eLimited,
398  endTrack,track.GetVolume());
399  //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
400  if(eLimited == kDoNot)
401  {
402  //AH G4cout << " computestep came back with not-boundary " << G4endl;
403  // Track is not on the boundary
404  fOnBoundary = false;
405  fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition());
406  }
407  else
408  {
409  // Track is on the boundary
410  //AH G4cout << " FOUND A BOUNDARY ! " << track.GetCurrentStepNumber() << G4endl;
411  fOnBoundary = true;
412  // proposedSafety = fGhostSafety;
413  }
414  proposedSafety = fGhostSafety;
415  if(eLimited == kUnique || eLimited == kSharedOther) {
416  *selection = CandidateForSelection;
417  }else if (eLimited == kSharedTransport) {
418  returnedStep *= (1.0 + 1.0e-9);
419  // Expand to disable its selection in Step Manager comparison
420  }
421 
422  }
423 
424  // ----------------------------------------------
425  // Returns the fGhostSafety as the proposedSafety
426  // The SteppingManager will take care of keeping
427  // the smallest one.
428  // ----------------------------------------------
429  return returnedStep;
430 
431  } else {
432 
433  return DBL_MAX;
434 
435  }
436 
437 }
438 
442 {
443  return -1.0;
444 }
445 
447 AtRestDoIt(const G4Track&, const G4Step&)
448 {
449  return 0;
450 }
451 
453 AlongStepDoIt(const G4Track& aTrack, const G4Step& )
454 {
455  // Dummy ParticleChange ie: does nothing
456  // Expecting G4Transportation to move the track
457  //AH G4cout << " along step do it " << G4endl;
458  pParticleChange->Initialize(aTrack);
459  return pParticleChange;
460  // return 0;
461 }
462 
464 {
465  fGhostStep->SetTrack(step.GetTrack());
469 
470  *fGhostPreStepPoint = *(step.GetPreStepPoint());
472 
473 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
474 // Set StepStatus for ghost world
475 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
476  if(fOnBoundary)
480 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
481 }
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
G4double condition(const G4ErrorSymMatrix &m)
void SetParallelWorld(G4String parallelWorldName)
void SetStepLength(G4double value)
virtual void Initialize(const G4Track &)
virtual G4double AlongStepGetPhysicalInteractionLength(const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
void SetTrack(G4Track *value)
G4int verboseLevel
Definition: G4VProcess.hh:368
G4StepPoint * fGhostPreStepPoint
virtual G4VParticleChange * AlongStepDoIt(const G4Track &, const G4Step &)
G4double GetStepLength() const
const G4VIStore & fIStore
G4StepStatus GetStepStatus() const
ELimited
const G4ThreeVector & GetPosition() const
void CopyStep(const G4Step &step)
G4TrackStatus GetTrackStatus() const
virtual G4VParticleChange * AtRestDoIt(const G4Track &, const G4Step &)
G4TouchableHandle CreateTouchableHandle(G4int navId) const
G4SteppingControl GetControlFlag() const
const G4VTouchable * GetTouchable() const
#define G4ThreadLocal
Definition: tls.hh:52
void SetStepStatus(const G4StepStatus aValue)
virtual G4double PostStepGetPhysicalInteractionLength(const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
G4TouchableHandle fNewGhostTouchable
G4StepPoint * GetPreStepPoint() const
G4VPhysicalVolume * fGhostWorld
G4GLOB_DLL std::ostream G4cout
G4TouchableHandle fOldGhostTouchable
G4int GetCurrentStepNumber() const
virtual const G4String & GetName() const
G4VPhysicalVolume * GetPhysicalVolume() const
void SetControlFlag(G4SteppingControl StepControlFlag)
bool G4bool
Definition: G4Types.hh:79
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
Definition: G4Step.hh:76
G4TransportationManager * fTransportationManager
G4ParticleChange * fParticleChange
const G4String & GetProcessName() const
Definition: G4VProcess.hh:408
virtual G4double GetImportance(const G4GeometryCell &gCell) const =0
G4double GetTotalEnergyDeposit() const
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
static G4TransportationManager * GetTransportationManager()
virtual void KillTrack() const
void StartTracking(G4Track *)
virtual void Initialize(const G4Track &)
G4VPhysicalVolume * GetParallelWorld(const G4String &worldName)
G4int ActivateNavigator(G4Navigator *aNavigator)
void DoIt(const G4Track &aTrack, G4ParticleChange *aParticleChange, const G4Nsplit_Weight &nw)
const G4ThreeVector & GetMomentumDirection() const
G4Navigator * GetNavigator(const G4String &worldName)
G4StepPoint * GetPostStepPoint() const
G4StepPoint * fGhostPostStepPoint
virtual G4int GetReplicaNumber(G4int depth=0) const
Definition: G4VTouchable.cc:58
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
G4VPhysicalVolume * GetVolume() const
G4double GetWeight() const
G4ImportanceProcess(const G4VImportanceAlgorithm &aImportanceAlgorithm, const G4VIStore &aIstore, const G4VTrackTerminator *TrackTerminator, const G4String &aName="ImportanceProcess", G4bool para=false)
#define G4endl
Definition: G4ios.hh:61
G4String theProcessName
Definition: G4VProcess.hh:335
void SetTotalEnergyDeposit(G4double value)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
G4ForceCondition
G4Track * GetTrack() const
virtual G4double AtRestGetPhysicalInteractionLength(const G4Track &, G4ForceCondition *)
void SetTouchableHandle(const G4TouchableHandle &apValue)
#define DBL_MAX
Definition: templates.hh:83
const G4TouchableHandle & GetTouchableHandle() const
static void Update(G4FieldTrack *, const G4Track *)
G4SamplingPostStepAction * fPostStepAction
G4GPILSelection
const G4VImportanceAlgorithm & fImportanceAlgorithm
virtual G4VParticleChange * PostStepDoIt(const G4Track &, const G4Step &)