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

#include <G4AdjointForcedInteractionForGamma.hh>

+ Inheritance diagram for G4AdjointForcedInteractionForGamma:

Public Member Functions

 G4AdjointForcedInteractionForGamma (G4String process_name)
 
virtual ~G4AdjointForcedInteractionForGamma ()
 
void PreparePhysicsTable (const G4ParticleDefinition &)
 
void BuildPhysicsTable (const G4ParticleDefinition &)
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &track, const G4Step &step)
 
void RegisterAdjointComptonModel (G4VEmAdjointModel *aAdjointComptonModel)
 
void RegisterAdjointBremModel (G4VEmAdjointModel *aAdjointBremModel)
 
- Public Member Functions inherited from G4VContinuousDiscreteProcess
 G4VContinuousDiscreteProcess (const G4String &, G4ProcessType aType=fNotDefined)
 
 G4VContinuousDiscreteProcess (G4VContinuousDiscreteProcess &)
 
virtual ~G4VContinuousDiscreteProcess ()
 
G4VContinuousDiscreteProcessoperator= (const G4VContinuousDiscreteProcess &)=delete
 
virtual G4double PostStepGetPhysicalInteractionLength (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4VParticleChangePostStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
 
virtual G4VParticleChangeAlongStepDoIt (const G4Track &, const G4Step &)
 
virtual G4double AtRestGetPhysicalInteractionLength (const G4Track &, G4ForceCondition *)
 
virtual G4VParticleChangeAtRestDoIt (const G4Track &, const G4Step &)
 
- Public Member Functions inherited from G4VProcess
 G4VProcess (const G4String &aName="NoName", G4ProcessType aType=fNotDefined)
 
 G4VProcess (const G4VProcess &right)
 
virtual ~G4VProcess ()
 
G4VProcessoperator= (const G4VProcess &)=delete
 
G4bool operator== (const G4VProcess &right) const
 
G4bool 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
 
virtual void ProcessDescription (std::ostream &outfile) const
 
void SetVerboseLevel (G4int value)
 
G4int GetVerboseLevel () const
 
virtual void SetMasterProcess (G4VProcess *masterP)
 
const G4VProcessGetMasterProcess () const
 
virtual void BuildWorkerPhysicsTable (const G4ParticleDefinition &part)
 
virtual void PrepareWorkerPhysicsTable (const G4ParticleDefinition &)
 

Protected Member Functions

virtual G4double GetMeanFreePath (const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
 
virtual G4double GetContinuousStepLimit (const G4Track &aTrack, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)
 
- Protected Member Functions inherited from G4VContinuousDiscreteProcess
virtual G4double GetMeanFreePath (const G4Track &aTrack, G4double previousStepSize, G4ForceCondition *condition)=0
 
virtual G4double GetContinuousStepLimit (const G4Track &aTrack, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety)=0
 
void SetGPILSelection (G4GPILSelection selection)
 
G4GPILSelection GetGPILSelection () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double prevStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- Protected Attributes inherited from G4VProcess
const G4ProcessManageraProcessManager = nullptr
 
G4VParticleChangepParticleChange = nullptr
 
G4ParticleChange aParticleChange
 
G4double theNumberOfInteractionLengthLeft = -1.0
 
G4double currentInteractionLength = -1.0
 
G4double theInitialNumberOfInteractionLength = -1.0
 
G4String theProcessName
 
G4String thePhysicsTableFileName
 
G4ProcessType theProcessType = fNotDefined
 
G4int theProcessSubType = -1
 
G4double thePILfactor = 1.0
 
G4int verboseLevel = 0
 
G4bool enableAtRestDoIt = true
 
G4bool enableAlongStepDoIt = true
 
G4bool enablePostStepDoIt = true
 

Detailed Description

Definition at line 72 of file G4AdjointForcedInteractionForGamma.hh.

Constructor & Destructor Documentation

◆ G4AdjointForcedInteractionForGamma()

G4AdjointForcedInteractionForGamma::G4AdjointForcedInteractionForGamma ( G4String  process_name)

Definition at line 37 of file G4AdjointForcedInteractionForGamma.cc.

38 :
39 G4VContinuousDiscreteProcess(process_name),theAdjointComptonModel(0),theAdjointBremModel(0)
40{ theAdjointCSManager = G4AdjointCSManager::GetAdjointCSManager();
41 fParticleChange=new G4ParticleChange();
42 lastAdjCS=0.;
43 trackid = nstep = 0;
44 is_free_flight_gamma = false;
45 copy_gamma_for_forced_interaction = false;
46 last_free_flight_trackid=1000;
47
48 theAdjointComptonModel =0;
49 theAdjointBremModel=0;
50
51 acc_track_length=0.;
52 acc_nb_adj_interaction_length=0.;
53 acc_nb_fwd_interaction_length=0.;
54 total_acc_nb_adj_interaction_length=0.;
55 total_acc_nb_fwd_interaction_length=0.;
56 continue_gamma_as_new_free_flight =false;
57
58
59
60}
static G4AdjointCSManager * GetAdjointCSManager()

◆ ~G4AdjointForcedInteractionForGamma()

G4AdjointForcedInteractionForGamma::~G4AdjointForcedInteractionForGamma ( )
virtual

Definition at line 63 of file G4AdjointForcedInteractionForGamma.cc.

65{ if (fParticleChange) delete fParticleChange;
66}

Member Function Documentation

◆ AlongStepDoIt()

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 156 of file G4AdjointForcedInteractionForGamma.cc.

157{ fParticleChange->Initialize(track);
158 //Compute nb of interactions length over step length
160 G4double stepLength = track.GetStep()->GetStepLength();
161 G4double ekin = track.GetKineticEnergy();
162 G4double nb_fwd_interaction_length_over_step=0.;
163 G4double nb_adj_interaction_length_over_step=0.;
166 ekin,track.GetMaterialCutsCouple());
167 nb_fwd_interaction_length_over_step = stepLength*lastFwdCS;
168 nb_adj_interaction_length_over_step = stepLength*lastAdjCS;
169 G4double fwd_survival_probability=std::exp(-nb_fwd_interaction_length_over_step);
170 G4double mc_induced_survival_probability=1.;
171
172 if (is_free_flight_gamma) { //for free_flight survival probability stays 1
173 //Accumulate the number of interaction lengths during free flight of gamma
174 total_acc_nb_fwd_interaction_length+=nb_fwd_interaction_length_over_step;
175 total_acc_nb_adj_interaction_length+=nb_adj_interaction_length_over_step;
176 acc_track_length+=stepLength;
177 }
178 else {
179 G4double previous_acc_nb_adj_interaction_length =acc_nb_adj_interaction_length;
180 acc_nb_fwd_interaction_length+=nb_fwd_interaction_length_over_step;
181 acc_nb_adj_interaction_length+=nb_adj_interaction_length_over_step;
182 theNumberOfInteractionLengthLeft-=nb_adj_interaction_length_over_step;
183
184 //Following condition to remove very rare FPE issue
185 //if (total_acc_nb_adj_interaction_length <= 1.e-50 && theNumberOfInteractionLengthLeft<=1.e-50) { //condition added to avoid FPE issue
186 // VI 06.11.2017 - new condition
187 if (std::abs(total_acc_nb_adj_interaction_length - previous_acc_nb_adj_interaction_length) <= 1.e-15) {
188 mc_induced_survival_probability = 1.e50;
189 /*
190 G4cout << "FPE protection: " << total_acc_nb_adj_interaction_length << " "
191 << previous_acc_nb_adj_interaction_length << " "
192 << acc_nb_fwd_interaction_length << " "
193 << acc_nb_adj_interaction_length << " "
194 << theNumberOfInteractionLengthLeft
195 << G4endl;
196 */
197 }
198 else {
199 mc_induced_survival_probability= std::exp(-acc_nb_adj_interaction_length)-std::exp(-total_acc_nb_adj_interaction_length);
200 mc_induced_survival_probability=mc_induced_survival_probability/(std::exp(-previous_acc_nb_adj_interaction_length)-std::exp(-total_acc_nb_adj_interaction_length));
201 }
202 }
203 G4double weight_correction = fwd_survival_probability/mc_induced_survival_probability;
204
205 //weight_correction = 1.;
206 //Caution!!!
207 // It is important to select the weight of the post_step_point
208 // as the current weight and not the weight of the track, as t
209 // the weight of the track is changed after having applied all
210 // the along_step_do_it.
211 G4double new_weight=weight_correction*track.GetStep()->GetPostStepPoint()->GetWeight();
212/*
213 G4cout<<"New weight "<<new_weight<<std::endl;
214 G4cout<<"Weight correction "<<weight_correction<<std::endl;
215 */
216
217 fParticleChange->SetParentWeightByProcess(false);
218 fParticleChange->SetSecondaryWeightByProcess(false);
219 fParticleChange->ProposeParentWeight(new_weight);
220
221 return fParticleChange;
222}
double G4double
Definition: G4Types.hh:83
G4double GetTotalForwardCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
G4double GetTotalAdjointCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
static G4AdjointGamma * AdjointGamma()
virtual void Initialize(const G4Track &)
G4double GetWeight() const
G4double GetStepLength() const
G4StepPoint * GetPostStepPoint() const
const G4ThreeVector & GetPosition() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
const G4Step * GetStep() const
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
void ProposeParentWeight(G4double finalWeight)
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:331

◆ BuildPhysicsTable()

void G4AdjointForcedInteractionForGamma::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 74 of file G4AdjointForcedInteractionForGamma.cc.

75{
76 theAdjointCSManager->BuildCrossSectionMatrices(); //do not worry it will be done just once
77 theAdjointCSManager->BuildTotalSigmaTables();
78}

◆ GetContinuousStepLimit()

G4double G4AdjointForcedInteractionForGamma::GetContinuousStepLimit ( const G4Track aTrack,
G4double  previousStepSize,
G4double  currentMinimumStep,
G4double currentSafety 
)
protectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 272 of file G4AdjointForcedInteractionForGamma.cc.

276{return DBL_MAX;
277}
#define DBL_MAX
Definition: templates.hh:62

◆ GetMeanFreePath()

G4double G4AdjointForcedInteractionForGamma::GetMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 280 of file G4AdjointForcedInteractionForGamma.cc.

283{ return 0.;
284}

◆ PostStepDoIt()

G4VParticleChange * G4AdjointForcedInteractionForGamma::PostStepDoIt ( const G4Track track,
const G4Step  
)
virtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 97 of file G4AdjointForcedInteractionForGamma.cc.

98{ fParticleChange->Initialize(track);
99 //For the free flight gamma no interaction occur but a gamma with same property is
100 //produces for further forced interaction
101 //It is done at the very beginning of the track such that the weight can be the same
102 if (copy_gamma_for_forced_interaction) {
103 G4ThreeVector theGammaMomentum = track.GetMomentum();
104 fParticleChange->AddSecondary(new G4DynamicParticle(G4AdjointGamma::AdjointGamma(),theGammaMomentum));
105 fParticleChange->SetParentWeightByProcess(false);
106 fParticleChange->SetSecondaryWeightByProcess(false);
107 }
108 else { //Occurrence of forced interaction
109
110 //Selection of the model to be called
111 G4VEmAdjointModel* theSelectedModel =0;
112 G4bool is_scat_proj_to_proj_case=false;
113 if (!theAdjointComptonModel && !theAdjointBremModel) return fParticleChange;
114 if (!theAdjointComptonModel) {
115 theSelectedModel = theAdjointBremModel;
116 is_scat_proj_to_proj_case=false;
117 //This is needed because the results of it will be used in the post step do it weight correction inside the model
118 theAdjointBremModel->AdjointCrossSection(
119 track.GetMaterialCutsCouple(),track.GetKineticEnergy(), false);
120
121 }
122 else if (!theAdjointBremModel) {
123 theSelectedModel = theAdjointComptonModel;
124 is_scat_proj_to_proj_case=true;
125 }
126 else { //Choose the model according to cross sections
127
128 G4double bremAdjCS = theAdjointBremModel->AdjointCrossSection(
129 track.GetMaterialCutsCouple(),track.GetKineticEnergy(), false);
130 if (G4UniformRand()*lastAdjCS<bremAdjCS) {
131 theSelectedModel = theAdjointBremModel;
132 is_scat_proj_to_proj_case=false;
133 }
134 else {
135 theSelectedModel = theAdjointComptonModel;
136 is_scat_proj_to_proj_case=true;
137 }
138 }
139
140 //Compute the weight correction factor
141 G4double one_over_effectiveAdjointCS= (1.-std::exp(acc_nb_adj_interaction_length-total_acc_nb_adj_interaction_length))/lastAdjCS;
142 G4double weight_correction_factor = lastAdjCS*one_over_effectiveAdjointCS;
143 //G4cout<<"Weight correction factor start "<<weight_correction_factor<<std::endl;
144 //Call the selected model without correction of the weight in the model
145 theSelectedModel->SetCorrectWeightForPostStepInModel(false);
146 theSelectedModel->SetAdditionalWeightCorrectionFactorForPostStepOutsideModel(weight_correction_factor);
147 theSelectedModel->SampleSecondaries(track,is_scat_proj_to_proj_case,fParticleChange);
148 theSelectedModel->SetCorrectWeightForPostStepInModel(true);
149
150 continue_gamma_as_new_free_flight =true;
151 }
152 return fParticleChange;
153}
bool G4bool
Definition: G4Types.hh:86
#define G4UniformRand()
Definition: Randomize.hh:52
void AddSecondary(G4Track *aSecondary)
G4ThreeVector GetMomentum() const
virtual void SampleSecondaries(const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)=0
void SetAdditionalWeightCorrectionFactorForPostStepOutsideModel(G4double factor)
void SetCorrectWeightForPostStepInModel(G4bool aBool)
virtual G4double AdjointCrossSection(const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)

◆ PostStepGetPhysicalInteractionLength()

G4double G4AdjointForcedInteractionForGamma::PostStepGetPhysicalInteractionLength ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
virtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 225 of file G4AdjointForcedInteractionForGamma.cc.

229{ G4int step_id = track.GetCurrentStepNumber();
231 copy_gamma_for_forced_interaction = false;
232 G4int track_id=track.GetTrackID();
233 is_free_flight_gamma = (track_id != last_free_flight_trackid+1 || continue_gamma_as_new_free_flight);
234 if (is_free_flight_gamma) {
235 if (step_id == 1 || continue_gamma_as_new_free_flight) {
237 //A gamma with same conditions will be generate at next post_step do it for the forced interaction
238 copy_gamma_for_forced_interaction = true;
239 last_free_flight_trackid = track_id;
240 acc_track_length=0.;
241 total_acc_nb_adj_interaction_length=0.;
242 total_acc_nb_fwd_interaction_length=0.;
243 continue_gamma_as_new_free_flight=false;
244 return 1.e-90;
245 }
246 else {
247 //Computation of accumulated length for
248 return DBL_MAX;
249 }
250 }
251 else { //compute the interaction length for forced interaction
252 if (step_id ==1) {
253 G4double min_val= std::exp(-total_acc_nb_adj_interaction_length);
254 theNumberOfInteractionLengthLeft = -std::log( min_val+G4UniformRand()*(1.-min_val));
256 acc_nb_adj_interaction_length=0.;
257 acc_nb_fwd_interaction_length=0.;
258 }
259 G4VPhysicalVolume* thePostPhysVolume = track.GetStep()->GetPreStepPoint()->GetPhysicalVolume();
260 G4double ekin =track.GetKineticEnergy();
261 G4double postCS=0.;
262 if (thePostPhysVolume){
264 ekin,thePostPhysVolume->GetLogicalVolume()->GetMaterialCutsCouple());
265 }
266 if (postCS>0.) return theNumberOfInteractionLengthLeft/postCS;
267 else return DBL_MAX;
268 }
269}
G4double condition(const G4ErrorSymMatrix &m)
@ NotForced
@ Forced
int G4int
Definition: G4Types.hh:85
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4VPhysicalVolume * GetPhysicalVolume() const
G4StepPoint * GetPreStepPoint() const
G4int GetTrackID() const
G4int GetCurrentStepNumber() const
G4LogicalVolume * GetLogicalVolume() const
G4double theInitialNumberOfInteractionLength
Definition: G4VProcess.hh:338

◆ PreparePhysicsTable()

void G4AdjointForcedInteractionForGamma::PreparePhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 69 of file G4AdjointForcedInteractionForGamma.cc.

70{;
71}

◆ RegisterAdjointBremModel()

void G4AdjointForcedInteractionForGamma::RegisterAdjointBremModel ( G4VEmAdjointModel aAdjointBremModel)
inline

Definition at line 91 of file G4AdjointForcedInteractionForGamma.hh.

91{theAdjointBremModel = aAdjointBremModel;}

◆ RegisterAdjointComptonModel()

void G4AdjointForcedInteractionForGamma::RegisterAdjointComptonModel ( G4VEmAdjointModel aAdjointComptonModel)
inline

Definition at line 90 of file G4AdjointForcedInteractionForGamma.hh.

90{theAdjointComptonModel = aAdjointComptonModel;}

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