Geant4  10.03.p01
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
G4AdjointForcedInteractionForGamma Class Reference

#include <G4AdjointForcedInteractionForGamma.hh>

Inheritance diagram for G4AdjointForcedInteractionForGamma:
Collaboration 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 ()
 
virtual G4double AlongStepGetPhysicalInteractionLength (const G4Track &, G4double previousStepSize, G4double currentMinimumStep, G4double &currentSafety, G4GPILSelection *selection)
 
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 ()
 
G4int operator== (const G4VProcess &right) const
 
G4int operator!= (const G4VProcess &right) const
 
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 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
 
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
void SetGPILSelection (G4GPILSelection selection)
 
G4GPILSelection GetGPILSelection () const
 
- Protected Member Functions inherited from G4VProcess
void SubtractNumberOfInteractionLengthLeft (G4double previousStepSize)
 
void ClearNumberOfInteractionLengthLeft ()
 

Additional Inherited Members

- Static Public Member Functions inherited from G4VProcess
static const G4StringGetProcessTypeName (G4ProcessType)
 
- 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
 

Detailed Description

Definition at line 73 of file G4AdjointForcedInteractionForGamma.hh.

Constructor & Destructor Documentation

G4AdjointForcedInteractionForGamma::G4AdjointForcedInteractionForGamma ( G4String  process_name)

Definition at line 39 of file G4AdjointForcedInteractionForGamma.cc.

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

Here is the call graph for this function:

G4AdjointForcedInteractionForGamma::~G4AdjointForcedInteractionForGamma ( )
virtual

Definition at line 65 of file G4AdjointForcedInteractionForGamma.cc.

66 { if (fParticleChange) delete fParticleChange;
67 }

Member Function Documentation

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 157 of file G4AdjointForcedInteractionForGamma.cc.

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

Here is the call graph for this function:

void G4AdjointForcedInteractionForGamma::BuildPhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 75 of file G4AdjointForcedInteractionForGamma.cc.

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

Here is the call graph for this function:

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

Implements G4VContinuousDiscreteProcess.

Definition at line 262 of file G4AdjointForcedInteractionForGamma.cc.

266 {return DBL_MAX;
267 }
#define DBL_MAX
Definition: templates.hh:83
G4double G4AdjointForcedInteractionForGamma::GetMeanFreePath ( const G4Track track,
G4double  previousStepSize,
G4ForceCondition condition 
)
protectedvirtual

Implements G4VContinuousDiscreteProcess.

Definition at line 270 of file G4AdjointForcedInteractionForGamma.cc.

273 { return 0.;
274 }
G4VParticleChange * G4AdjointForcedInteractionForGamma::PostStepDoIt ( const G4Track track,
const G4Step  
)
virtual

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 98 of file G4AdjointForcedInteractionForGamma.cc.

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

Here is the call graph for this function:

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

Reimplemented from G4VContinuousDiscreteProcess.

Definition at line 215 of file G4AdjointForcedInteractionForGamma.cc.

219 { G4int step_id = track.GetCurrentStepNumber();
220  *condition = NotForced;
221  copy_gamma_for_forced_interaction = false;
222  G4int track_id=track.GetTrackID();
223  is_free_flight_gamma = (track_id != last_free_flight_trackid+1 || continue_gamma_as_new_free_flight);
224  if (is_free_flight_gamma) {
225  if (step_id == 1 || continue_gamma_as_new_free_flight) {
226  *condition=Forced;
227  //A gamma with same conditions will be generate at next post_step do it for the forced interaction
228  copy_gamma_for_forced_interaction = true;
229  last_free_flight_trackid = track_id;
230  acc_track_length=0.;
231  total_acc_nb_adj_interaction_length=0.;
232  total_acc_nb_fwd_interaction_length=0.;
233  continue_gamma_as_new_free_flight=false;
234  return 1.e-90;
235  }
236  else {
237  //Computation of accumulated length for
238  return DBL_MAX;
239  }
240  }
241  else { //compute the interaction length for forced interaction
242  if (step_id ==1) {
243  G4double min_val= std::exp(-total_acc_nb_adj_interaction_length);
244  theNumberOfInteractionLengthLeft = -std::log( min_val+G4UniformRand()*(1.-min_val));
246  acc_nb_adj_interaction_length=0.;
247  acc_nb_fwd_interaction_length=0.;
248  }
249  G4VPhysicalVolume* thePostPhysVolume = track.GetStep()->GetPreStepPoint()->GetPhysicalVolume();
250  G4double ekin =track.GetKineticEnergy();
251  G4double postCS=0.;
252  if (thePostPhysVolume){
254  ekin,thePostPhysVolume->GetLogicalVolume()->GetMaterialCutsCouple());
255  }
256  if (postCS>0.) return theNumberOfInteractionLengthLeft/postCS;
257  else return DBL_MAX;
258  }
259 }
G4double condition(const G4ErrorSymMatrix &m)
static G4AdjointGamma * AdjointGamma()
G4double GetTotalAdjointCS(G4ParticleDefinition *aPartDef, G4double Ekin, const G4MaterialCutsCouple *aCouple)
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
const G4Step * GetStep() const
int G4int
Definition: G4Types.hh:78
G4StepPoint * GetPreStepPoint() const
G4double GetKineticEnergy() const
#define G4UniformRand()
Definition: Randomize.hh:97
G4int GetCurrentStepNumber() const
G4VPhysicalVolume * GetPhysicalVolume() const
G4int GetTrackID() const
G4LogicalVolume * GetLogicalVolume() const
const G4MaterialCutsCouple * GetMaterialCutsCouple() const
G4double theInitialNumberOfInteractionLength
Definition: G4VProcess.hh:300
double G4double
Definition: G4Types.hh:76
#define DBL_MAX
Definition: templates.hh:83
static G4AdjointCSManager * GetAdjointCSManager()

Here is the call graph for this function:

void G4AdjointForcedInteractionForGamma::PreparePhysicsTable ( const G4ParticleDefinition )
virtual

Reimplemented from G4VProcess.

Definition at line 70 of file G4AdjointForcedInteractionForGamma.cc.

71 {;
72 }
void G4AdjointForcedInteractionForGamma::RegisterAdjointBremModel ( G4VEmAdjointModel aAdjointBremModel)
inline

Definition at line 92 of file G4AdjointForcedInteractionForGamma.hh.

92 {theAdjointBremModel = aAdjointBremModel;}
void G4AdjointForcedInteractionForGamma::RegisterAdjointComptonModel ( G4VEmAdjointModel aAdjointComptonModel)
inline

Definition at line 91 of file G4AdjointForcedInteractionForGamma.hh.

91 {theAdjointComptonModel = aAdjointComptonModel;}

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