Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4WeightWindowProcess Class Reference

#include <G4WeightWindowProcess.hh>

+ Inheritance diagram for G4WeightWindowProcess:

Public Member Functions

 G4WeightWindowProcess (const G4VWeightWindowAlgorithm &aWeightWindowAlgorithm, const G4VWeightWindowStore &aWWStore, const G4VTrackTerminator *TrackTerminator, G4PlaceOfAction placeOfAction, const G4String &aName="WeightWindowProcess", G4bool para=false)
 
virtual ~G4WeightWindowProcess ()
 
void SetParallelWorld (G4String parallelWorldName)
 
void SetParallelWorld (G4VPhysicalVolume *parallelWorld)
 
void StartTracking (G4Track *)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual void KillTrack () const
 
virtual const G4StringGetName () const
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double, G4double, G4double &, G4GPILSelection *)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
virtual G4VParticleChangePostStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &track, const G4Step &stepData)=0
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)=0
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &track, G4ForceCondition *condition)=0
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)=0
 
G4double GetCurrentInteractionLength () const
 
void SetPILfactor (G4double value)
 
G4double GetPILfactor () const
 
G4double AlongStepGPIL (const G4Track &track, G4double previousStepSize, G4double currentMinimumStep, G4double &proposedSafety, G4GPILSelection *selection)
 
G4double AtRestGPIL (const G4Track &track, G4ForceCondition *condition)
 
G4double PostStepGPIL (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4bool IsApplicable (const G4ParticleDefinition &)
 
virtual void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual void PreparePhysicsTable (const G4ParticleDefinition &)
 
virtual G4bool StorePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
virtual G4bool RetrievePhysicsTable (const G4ParticleDefinition *, const G4String &, G4bool)
 
const G4StringGetPhysicsTableFileName (const G4ParticleDefinition *, const G4String &directory, const G4String &tableName, G4bool ascii=false)
 
const G4StringGetProcessName () const
 
G4ProcessType GetProcessType () const
 
void SetProcessType (G4ProcessType)
 
G4int GetProcessSubType () const
 
void SetProcessSubType (G4int)
 
virtual void StartTracking (G4Track *)
 
virtual void EndTracking ()
 
virtual void SetProcessManager (const G4ProcessManager *)
 
virtual const G4ProcessManagerGetProcessManager ()
 
virtual void ResetNumberOfInteractionLengthLeft ()
 
G4double GetNumberOfInteractionLengthLeft () const
 
G4double GetTotalNumberOfInteractionLengthTraversed () const
 
G4bool isAtRestDoItIsEnabled () const
 
G4bool isAlongStepDoItIsEnabled () const
 
G4bool isPostStepDoItIsEnabled () const
 
virtual void DumpInfo () const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
- Public Member Functions inherited from G4VTrackTerminator
 G4VTrackTerminator ()
 
virtual ~G4VTrackTerminator ()
 
virtual void KillTrack () const =0
 
virtual const G4StringGetName () const =0
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager
 
G4VParticleChangepParticleChange
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft
 
G4double currentInteractionLength
 
G4double theInitialNumberOfInteractionLength
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType
 
G4int theProcessSubType
 
G4double thePILfactor
 
G4bool enableAtRestDoIt
 
G4bool enableAlongStepDoIt
 
G4bool enablePostStepDoIt
 
G4int verboseLevel
 
- Protected Attributes inherited from G4VTrackTerminator
G4double kCarTolerance
 

Detailed Description

Definition at line 62 of file G4WeightWindowProcess.hh.

Constructor & Destructor Documentation

◆ G4WeightWindowProcess()

G4WeightWindowProcess::G4WeightWindowProcess ( const G4VWeightWindowAlgorithm aWeightWindowAlgorithm,
const G4VWeightWindowStore aWWStore,
const G4VTrackTerminator TrackTerminator,
G4PlaceOfAction  placeOfAction,
const G4String aName = "WeightWindowProcess",
G4bool  para = false 
)

Definition at line 55 of file G4WeightWindowProcess.cc.

61 : G4VProcess(aName),
62 fParticleChange(new G4ParticleChange),
63 fWeightWindowAlgorithm(aWeightWindowAlgorithm),
64 fWeightWindowStore(aWWStore),
65 fPostStepAction(0),
66 fPlaceOfAction(placeOfAction),
67 fGhostNavigator(0), fNavigatorID(-1), fFieldTrack('0')
68{
69 if (TrackTerminator)
70 {
71 fPostStepAction = new G4SamplingPostStepAction(*TrackTerminator);
72 }
73 else
74 {
75 fPostStepAction = new G4SamplingPostStepAction(*this);
76 }
77 if (!fParticleChange)
78 {
79 G4Exception("G4WeightWindowProcess::G4WeightWindowProcess()",
80 "FatalError", FatalException,
81 "Failed allocation of G4ParticleChange !");
82 }
83 G4VProcess::pParticleChange = fParticleChange;
84
85 fGhostStep = new G4Step();
86 fGhostPreStepPoint = fGhostStep->GetPreStepPoint();
87 fGhostPostStepPoint = fGhostStep->GetPostStepPoint();
88
90 fPathFinder = G4PathFinder::GetInstance();
91
92 if (verboseLevel>0)
93 {
94 G4cout << GetProcessName() << " is created " << G4endl;
95 }
96
97 paraflag = para;
98
99}
@ FatalException
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
static G4PathFinder * GetInstance()
Definition: G4PathFinder.cc:57
Definition: G4Step.hh:78
G4StepPoint * GetPreStepPoint() const
G4StepPoint * GetPostStepPoint() const
static G4TransportationManager * GetTransportationManager()
G4int verboseLevel
Definition: G4VProcess.hh:368
G4VParticleChange * pParticleChange
Definition: G4VProcess.hh:283
const G4String & GetProcessName() const
Definition: G4VProcess.hh:379
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41

◆ ~G4WeightWindowProcess()

G4WeightWindowProcess::~G4WeightWindowProcess ( )
virtual

Definition at line 101 of file G4WeightWindowProcess.cc.

102{
103
104 delete fPostStepAction;
105 delete fParticleChange;
106 delete fGhostStep;
107
108}

Member Function Documentation

◆ AlongStepDoIt()

G4VParticleChange * G4WeightWindowProcess::AlongStepDoIt ( const G4Track track,
const G4Step  
)
virtual

Implements G4VProcess.

Definition at line 387 of file G4WeightWindowProcess.cc.

389{
390 // Dummy ParticleChange ie: does nothing
391 // Expecting G4Transportation to move the track
393 return pParticleChange;
394
395 // return 0;
396}
virtual void Initialize(const G4Track &)

◆ AlongStepGetPhysicalInteractionLength()

G4double G4WeightWindowProcess::AlongStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double proposedSafety,
G4GPILSelection selection 
)
virtual

Implements G4VProcess.

Definition at line 299 of file G4WeightWindowProcess.cc.

303{
304 if(paraflag) {
305 static G4FieldTrack endTrack('0');
306 static ELimited eLimited;
307
308 *selection = NotCandidateForSelection;
309 G4double returnedStep = DBL_MAX;
310
311 if (previousStepSize > 0.)
312 { fGhostSafety -= previousStepSize; }
313 // else
314 // { fGhostSafety = -1.; }
315 if (fGhostSafety < 0.) fGhostSafety = 0.0;
316
317 // ------------------------------------------
318 // Determination of the proposed STEP LENGTH:
319 // ------------------------------------------
320 if (currentMinimumStep <= fGhostSafety && currentMinimumStep > 0.)
321 {
322 // I have no chance to limit
323 returnedStep = currentMinimumStep;
324 fOnBoundary = false;
325 proposedSafety = fGhostSafety - currentMinimumStep;
326 }
327 else // (currentMinimumStep > fGhostSafety: I may limit the Step)
328 {
329 G4FieldTrackUpdator::Update(&fFieldTrack,&track);
330 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
331 // ComputeStep
332 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
333 returnedStep
334 = fPathFinder->ComputeStep(fFieldTrack,currentMinimumStep,fNavigatorID,
335 track.GetCurrentStepNumber(),fGhostSafety,eLimited,
336 endTrack,track.GetVolume());
337 //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
338 if(eLimited == kDoNot)
339 {
340 // Track is not on the boundary
341 fOnBoundary = false;
342 fGhostSafety = fGhostNavigator->ComputeSafety(endTrack.GetPosition());
343 }
344 else
345 {
346 // Track is on the boundary
347 fOnBoundary = true;
348 proposedSafety = fGhostSafety;
349 }
350 //xbug? proposedSafety = fGhostSafety;
351 if(eLimited == kUnique || eLimited == kSharedOther) {
352 *selection = CandidateForSelection;
353 }else if (eLimited == kSharedTransport) {
354 returnedStep *= (1.0 + 1.0e-9);
355 // Expand to disable its selection in Step Manager comparison
356 }
357
358 }
359
360 // ----------------------------------------------
361 // Returns the fGhostSafety as the proposedSafety
362 // The SteppingManager will take care of keeping
363 // the smallest one.
364 // ----------------------------------------------
365 return returnedStep;
366
367 } else {
368 return DBL_MAX;
369 // not sensible - time goes backwards! return -1.0;
370 }
371
372}
@ CandidateForSelection
@ NotCandidateForSelection
ELimited
@ kDoNot
@ kUnique
@ kSharedOther
@ kSharedTransport
double G4double
Definition: G4Types.hh:64
static void Update(G4FieldTrack *, const G4Track *)
virtual G4double ComputeSafety(const G4ThreeVector &globalpoint, const G4double pProposedMaxLength=DBL_MAX, const G4bool keepState=true)
G4double ComputeStep(const G4FieldTrack &pFieldTrack, G4double pCurrentProposedStepLength, G4int navigatorId, G4int stepNo, G4double &pNewSafety, ELimited &limitedStep, G4FieldTrack &EndState, G4VPhysicalVolume *currentVolume)
G4VPhysicalVolume * GetVolume() const
G4int GetCurrentStepNumber() const
#define DBL_MAX
Definition: templates.hh:83

◆ AtRestDoIt()

G4VParticleChange * G4WeightWindowProcess::AtRestDoIt ( const G4Track ,
const G4Step  
)
virtual

Implements G4VProcess.

Definition at line 381 of file G4WeightWindowProcess.cc.

383{
384 return 0;
385}

◆ AtRestGetPhysicalInteractionLength()

G4double G4WeightWindowProcess::AtRestGetPhysicalInteractionLength ( const G4Track ,
G4ForceCondition  
)
virtual

Implements G4VProcess.

Definition at line 374 of file G4WeightWindowProcess.cc.

377{
378 return -1.0;
379}

◆ GetName()

const G4String & G4WeightWindowProcess::GetName ( ) const
virtual

Implements G4VTrackTerminator.

Definition at line 294 of file G4WeightWindowProcess.cc.

295{
297}

◆ KillTrack()

void G4WeightWindowProcess::KillTrack ( ) const
virtual

Implements G4VTrackTerminator.

Definition at line 289 of file G4WeightWindowProcess.cc.

290{
291 fParticleChange->ProposeTrackStatus(fStopAndKill);
292}
@ fStopAndKill
void ProposeTrackStatus(G4TrackStatus status)

◆ PostStepDoIt()

G4VParticleChange * G4WeightWindowProcess::PostStepDoIt ( const G4Track aTrack,
const G4Step aStep 
)
virtual

Implements G4VProcess.

Definition at line 199 of file G4WeightWindowProcess.cc.

201{
202
203 fParticleChange->Initialize(aTrack);
204
205 if(paraflag) {
206 fOldGhostTouchable = fGhostPostStepPoint->GetTouchableHandle();
207 //xbug? fOnBoundary = false;
208 CopyStep(aStep);
209
210 if(fOnBoundary)
211 {
212//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
213// Locate the point and get new touchable
214//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
215 //?? fPathFinder->Locate(step.GetPostStepPoint()->GetPosition(),
216 //?? step.GetPostStepPoint()->GetMomentumDirection());
217 fNewGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
218//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
219 }
220 else
221 {
222 // Do I need this ??????????????????????????????????????????????????????????
223 // fGhostNavigator->LocateGlobalPointWithinVolume(track.GetPosition());
224 // ?????????????????????????????????????????????????????????????????????????
225
226 // fPathFinder->ReLocate(track.GetPosition());
227
228 // reuse the touchable
229 fNewGhostTouchable = fOldGhostTouchable;
230 }
231
232 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
233 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
234
235 }
236
237 if (aStep.GetStepLength() > kCarTolerance)
238 {
239// if ( ( fPlaceOfAction == onBoundaryAndCollision)
240// || ( (fPlaceOfAction == onBoundary) &&
241// (aStep.GetPostStepPoint()->GetStepStatus() == fGeomBoundary) )
242// || ( (fPlaceOfAction == onCollision) &&
243// (aStep.GetPostStepPoint()->GetStepStatus() != fGeomBoundary) ) )
244 if(paraflag) {
245 if ( ( fPlaceOfAction == onBoundaryAndCollision)
246 || ( (fPlaceOfAction == onBoundary) &&
247 (fGhostPostStepPoint->GetStepStatus() == fGeomBoundary) )
248 || ( (fPlaceOfAction == onCollision) &&
249 (fGhostPostStepPoint->GetStepStatus() != fGeomBoundary) ) )
250 {
251
252// G4StepPoint *postpoint = aStep.GetPostStepPoint();
253
254// G4GeometryCell postCell(*(postpoint->GetPhysicalVolume()),
255// postpoint->GetTouchable()->GetReplicaNumber());
256
257 G4GeometryCell postCell(*(fGhostPostStepPoint->GetPhysicalVolume()),
258 fGhostPostStepPoint->GetTouchable()->GetReplicaNumber());
259 G4Nsplit_Weight nw =
260 fWeightWindowAlgorithm.Calculate(aTrack.GetWeight(),
261 fWeightWindowStore.GetLowerWeight(postCell,
262 aTrack.GetKineticEnergy()));
263 fPostStepAction->DoIt(aTrack, fParticleChange, nw);
264 }
265 } else {
266 if ( ( fPlaceOfAction == onBoundaryAndCollision)
267 || ( (fPlaceOfAction == onBoundary) &&
269 || ( (fPlaceOfAction == onCollision) &&
271 {
272
273 G4StepPoint *postpoint = aStep.GetPostStepPoint();
274
275 G4GeometryCell postCell(*(postpoint->GetPhysicalVolume()),
276 postpoint->GetTouchable()->GetReplicaNumber());
277
278 G4Nsplit_Weight nw =
279 fWeightWindowAlgorithm.Calculate(aTrack.GetWeight(),
280 fWeightWindowStore.GetLowerWeight(postCell,
281 aTrack.GetKineticEnergy()));
282 fPostStepAction->DoIt(aTrack, fParticleChange, nw);
283 }
284 }
285 }
286 return fParticleChange;
287}
@ onCollision
@ onBoundaryAndCollision
@ onBoundary
@ fGeomBoundary
Definition: G4StepStatus.hh:54
virtual void Initialize(const G4Track &)
G4TouchableHandle CreateTouchableHandle(G4int navId) const
void DoIt(const G4Track &aTrack, G4ParticleChange *aParticleChange, const G4Nsplit_Weight &nw)
G4StepStatus GetStepStatus() const
const G4VTouchable * GetTouchable() const
void SetTouchableHandle(const G4TouchableHandle &apValue)
const G4TouchableHandle & GetTouchableHandle() const
G4VPhysicalVolume * GetPhysicalVolume() const
G4double GetStepLength() const
G4double GetWeight() const
G4double GetKineticEnergy() const
virtual G4int GetReplicaNumber(G4int depth=0) const
Definition: G4VTouchable.cc:58
virtual G4Nsplit_Weight Calculate(G4double init_w, G4double lowerWeightBound) const =0
virtual G4double GetLowerWeight(const G4GeometryCell &gCell, G4double partEnergy) const =0

◆ PostStepGetPhysicalInteractionLength()

G4double G4WeightWindowProcess::PostStepGetPhysicalInteractionLength ( const G4Track aTrack,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Implements G4VProcess.

Definition at line 185 of file G4WeightWindowProcess.cc.

189{
190// *condition = Forced;
191// return kInfinity;
192
193// *condition = StronglyForced;
194 *condition = Forced;
195 return DBL_MAX;
196}
G4double condition(const G4ErrorSymMatrix &m)
@ Forced

◆ SetParallelWorld() [1/2]

void G4WeightWindowProcess::SetParallelWorld ( G4String  parallelWorldName)

Definition at line 116 of file G4WeightWindowProcess.cc.

118{
119//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
120// Get pointers of the parallel world and its navigator
121//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
122 fGhostWorldName = parallelWorldName;
123 fGhostWorld = fTransportationManager->GetParallelWorld(fGhostWorldName);
124 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
125//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
126}
G4VPhysicalVolume * GetParallelWorld(const G4String &worldName)
G4Navigator * GetNavigator(const G4String &worldName)

Referenced by G4WeightWindowConfigurator::Configure().

◆ SetParallelWorld() [2/2]

void G4WeightWindowProcess::SetParallelWorld ( G4VPhysicalVolume parallelWorld)

Definition at line 128 of file G4WeightWindowProcess.cc.

130{
131//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
132// Get pointer of navigator
133//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
134 fGhostWorldName = parallelWorld->GetName();
135 fGhostWorld = parallelWorld;
136 fGhostNavigator = fTransportationManager->GetNavigator(fGhostWorld);
137//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
138}
const G4String & GetName() const

◆ StartTracking()

void G4WeightWindowProcess::StartTracking ( G4Track trk)
virtual

Reimplemented from G4VProcess.

Definition at line 145 of file G4WeightWindowProcess.cc.

146{
147//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
148// Activate navigator and get the navigator ID
149//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
150// G4cout << " G4ParallelWorldScoringProcess::StartTracking" << G4endl;
151
152 if(paraflag) {
153 if(fGhostNavigator)
154 { fNavigatorID = fTransportationManager->ActivateNavigator(fGhostNavigator); }
155 else
156 {
157 G4Exception("G4WeightWindowProcess::StartTracking",
158 "ProcParaWorld000",FatalException,
159 "G4WeightWindowProcess is used for tracking without having a parallel world assigned");
160 }
161//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
162
163// G4cout << "G4ParallelWorldScoringProcess::StartTracking <<<<<<<<<<<<<<<<<< " << G4endl;
164//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
165// Let PathFinder initialize
166//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
167 fPathFinder->PrepareNewTrack(trk->GetPosition(),trk->GetMomentumDirection());
168//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
169
170//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
171// Setup initial touchables for the first step
172//++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
173 fOldGhostTouchable = fPathFinder->CreateTouchableHandle(fNavigatorID);
174 fGhostPreStepPoint->SetTouchableHandle(fOldGhostTouchable);
175 fNewGhostTouchable = fOldGhostTouchable;
176 fGhostPostStepPoint->SetTouchableHandle(fNewGhostTouchable);
177
178 // Initialize
179 fGhostSafety = -1.;
180 fOnBoundary = false;
181 }
182}
void PrepareNewTrack(const G4ThreeVector &position, const G4ThreeVector &direction, G4VPhysicalVolume *massStartVol=0)
const G4ThreeVector & GetPosition() const
const G4ThreeVector & GetMomentumDirection() const
G4int ActivateNavigator(G4Navigator *aNavigator)

The documentation for this class was generated from the following files: