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

#include <G4AdjointIonIonisationModel.hh>

Inheritance diagram for G4AdjointIonIonisationModel:
Collaboration diagram for G4AdjointIonIonisationModel:

Public Member Functions

 G4AdjointIonIonisationModel ()
 
virtual ~G4AdjointIonIonisationModel ()
 
virtual void SampleSecondaries (const G4Track &aTrack, G4bool IsScatProjToProjCase, G4ParticleChange *fParticleChange)
 
virtual G4double DiffCrossSectionPerAtomPrimToSecond (G4double kinEnergyProj, G4double kinEnergyProd, G4double Z, G4double A=0.)
 
virtual void CorrectPostStepWeight (G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double GetSecondAdjEnergyMaxForScatProjToProjCase (G4double PrimAdjEnergy)
 
virtual G4double GetSecondAdjEnergyMinForScatProjToProjCase (G4double PrimAdjEnergy, G4double Tcut=0)
 
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase (G4double PrimAdjEnergy)
 
virtual G4double GetSecondAdjEnergyMinForProdToProjCase (G4double PrimAdjEnergy)
 
void SetUseOnlyBragg (G4bool aBool)
 
void SetIon (G4ParticleDefinition *adj_ion, G4ParticleDefinition *fwd_ion)
 
- Public Member Functions inherited from G4VEmAdjointModel
 G4VEmAdjointModel (const G4String &nam)
 
virtual ~G4VEmAdjointModel ()
 
virtual G4double AdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double GetAdjointCrossSection (const G4MaterialCutsCouple *aCouple, G4double primEnergy, G4bool IsScatProjToProjCase)
 
virtual G4double DiffCrossSectionPerAtomPrimToScatPrim (G4double kinEnergyProj, G4double kinEnergyScatProj, G4double Z, G4double A=0.)
 
virtual G4double DiffCrossSectionPerVolumePrimToSecond (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyProd)
 
virtual G4double DiffCrossSectionPerVolumePrimToScatPrim (const G4Material *aMaterial, G4double kinEnergyProj, G4double kinEnergyScatProj)
 
void DefineCurrentMaterial (const G4MaterialCutsCouple *couple)
 
std::vector< std::vector
< double > * > 
ComputeAdjointCrossSectionVectorPerAtomForSecond (G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
 
std::vector< std::vector
< double > * > 
ComputeAdjointCrossSectionVectorPerAtomForScatProj (G4double kinEnergyProd, G4double Z, G4double A=0., G4int nbin_pro_decade=10)
 
std::vector< std::vector
< double > * > 
ComputeAdjointCrossSectionVectorPerVolumeForSecond (G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
 
std::vector< std::vector
< double > * > 
ComputeAdjointCrossSectionVectorPerVolumeForScatProj (G4Material *aMaterial, G4double kinEnergyProd, G4int nbin_pro_decade=10)
 
void SetCSMatrices (std::vector< G4AdjointCSMatrix * > *Vec1CSMatrix, std::vector< G4AdjointCSMatrix * > *Vec2CSMatrix)
 
G4ParticleDefinitionGetAdjointEquivalentOfDirectPrimaryParticleDefinition ()
 
G4ParticleDefinitionGetAdjointEquivalentOfDirectSecondaryParticleDefinition ()
 
G4double GetHighEnergyLimit ()
 
G4double GetLowEnergyLimit ()
 
void SetHighEnergyLimit (G4double aVal)
 
void SetLowEnergyLimit (G4double aVal)
 
void DefineDirectEMModel (G4VEmModel *aModel)
 
void SetAdjointEquivalentOfDirectPrimaryParticleDefinition (G4ParticleDefinition *aPart)
 
void SetAdjointEquivalentOfDirectSecondaryParticleDefinition (G4ParticleDefinition *aPart)
 
void SetSecondPartOfSameType (G4bool aBool)
 
G4bool GetSecondPartOfSameType ()
 
void SetUseMatrix (G4bool aBool)
 
void SetUseMatrixPerElement (G4bool aBool)
 
void SetUseOnlyOneMatrixForAllElements (G4bool aBool)
 
void SetApplyCutInRange (G4bool aBool)
 
G4bool GetUseMatrix ()
 
G4bool GetUseMatrixPerElement ()
 
G4bool GetUseOnlyOneMatrixForAllElements ()
 
G4bool GetApplyCutInRange ()
 
G4String GetName ()
 
virtual void SetCSBiasingFactor (G4double aVal)
 
void SetCorrectWeightForPostStepInModel (G4bool aBool)
 
void SetAdditionalWeightCorrectionFactorForPostStepOutsideModel (G4double factor)
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmAdjointModel
G4double DiffCrossSectionFunction1 (G4double kinEnergyProj)
 
G4double DiffCrossSectionFunction2 (G4double kinEnergyProj)
 
G4double DiffCrossSectionPerVolumeFunctionForIntegrationOverEkinProj (G4double EkinProd)
 
G4double SampleAdjSecEnergyFromCSMatrix (size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
 
G4double SampleAdjSecEnergyFromCSMatrix (G4double prim_energy, G4bool IsScatProjToProjCase)
 
void SelectCSMatrix (G4bool IsScatProjToProjCase)
 
virtual G4double SampleAdjSecEnergyFromDiffCrossSectionPerAtom (G4double prim_energy, G4bool IsScatProjToProjCase)
 
- Protected Attributes inherited from G4VEmAdjointModel
G4VEmModeltheDirectEMModel
 
G4VParticleChangepParticleChange
 
const G4String name
 
G4int ASelectedNucleus
 
G4int ZSelectedNucleus
 
G4MaterialSelectedMaterial
 
G4double kinEnergyProdForIntegration
 
G4double kinEnergyScatProjForIntegration
 
G4double kinEnergyProjForIntegration
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForProdToProjBackwardScattering
 
std::vector< G4AdjointCSMatrix * > * pOnCSMatrixForScatProjToProjBackwardScattering
 
std::vector< G4doubleCS_Vs_ElementForScatProjToProjCase
 
std::vector< G4doubleCS_Vs_ElementForProdToProjCase
 
G4double lastCS
 
G4double lastAdjointCSForScatProjToProjCase
 
G4double lastAdjointCSForProdToProjCase
 
G4ParticleDefinitiontheAdjEquivOfDirectPrimPartDef
 
G4ParticleDefinitiontheAdjEquivOfDirectSecondPartDef
 
G4ParticleDefinitiontheDirectPrimaryPartDef
 
G4bool second_part_of_same_type
 
G4double preStepEnergy
 
G4MaterialcurrentMaterial
 
G4MaterialCutsCouplecurrentCouple
 
size_t currentMaterialIndex
 
size_t currentCoupleIndex
 
G4double currentTcutForDirectPrim
 
G4double currentTcutForDirectSecond
 
G4bool ApplyCutInRange
 
G4double mass_ratio_product
 
G4double mass_ratio_projectile
 
G4double HighEnergyLimit
 
G4double LowEnergyLimit
 
G4double CS_biasing_factor
 
G4bool UseMatrix
 
G4bool UseMatrixPerElement
 
G4bool UseOnlyOneMatrixForAllElements
 
size_t indexOfUsedCrossSectionMatrix
 
size_t model_index
 
G4bool correct_weight_for_post_step_in_model
 
G4double additional_weight_correction_factor_for_post_step_outside_model
 

Detailed Description

Definition at line 71 of file G4AdjointIonIonisationModel.hh.

Constructor & Destructor Documentation

G4AdjointIonIonisationModel::G4AdjointIonIonisationModel ( )

Definition at line 47 of file G4AdjointIonIonisationModel.cc.

47  :
48  G4VEmAdjointModel("Adjoint_IonIonisation")
49 {
50 
51 
52  UseMatrix =true;
53  UseMatrixPerElement = true;
54  ApplyCutInRange = true;
58  use_only_bragg = false; // for the Ion ionisation using the parametrised table model the cross sections and the sample of secondaries is done
59  // as in the BraggIonModel, Therefore the use of this flag;
60 
61  //The direct EM Model is taken has BetheBloch it is only used for the computation
62  // of the differential cross section.
63  //The Bragg model could be used as an alternative as it offers the same differential cross section
64 
65  theBetheBlochDirectEMModel = new G4BetheBlochModel(G4GenericIon::GenericIon());
66  theBraggIonDirectEMModel = new G4BraggIonModel(G4GenericIon::GenericIon());
70  /* theDirectPrimaryPartDef =fwd_ion;
71  theAdjEquivOfDirectPrimPartDef =adj_ion;
72 
73  DefineProjectileProperty();*/
74 
75 
76 
77 
78 }
G4ParticleDefinition * theDirectPrimaryPartDef
static G4AdjointElectron * AdjointElectron()
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
G4VEmAdjointModel(const G4String &nam)
G4bool UseOnlyOneMatrixForAllElements
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef

Here is the call graph for this function:

G4AdjointIonIonisationModel::~G4AdjointIonIonisationModel ( )
virtual

Definition at line 81 of file G4AdjointIonIonisationModel.cc.

82 {;}

Member Function Documentation

void G4AdjointIonIonisationModel::CorrectPostStepWeight ( G4ParticleChange fParticleChange,
G4double  old_weight,
G4double  adjointPrimKinEnergy,
G4double  projectileKinEnergy,
G4bool  IsScatProjToProjCase 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 266 of file G4AdjointIonIonisationModel.cc.

268 {
269  //It is needed because the direct cross section used to compute the differential cross section is not the one used in
270  // the direct model where the GenericIon stuff is considered with correction of effective charge. In the direct model the samnepl of secondaries does
271  // not reflect the integral cross section. The integral fwd cross section that we used to compute the differential CS
272  // match the sample of secondaries in the forward case despite the fact that its is not the same total CS than in the FWD case. For this reasion an extra
273  // weight correction is needed at the end.
274 
275 
276  G4double new_weight=old_weight;
277 
278  //the correction of CS due to the problem explained above
279  G4double kinEnergyProjScaled = massRatio*projectileKinEnergy;
280  theDirectEMModel =theBraggIonDirectEMModel;
281  if (kinEnergyProjScaled >2.*MeV && !use_only_bragg) theDirectEMModel = theBetheBlochDirectEMModel; //Bethe Bloch Model
283  G4double chargeSqRatio =1.;
284  if (chargeSquare>1.) chargeSqRatio = theDirectEMModel->GetChargeSquareRatio(theDirectPrimaryPartDef,currentMaterial,projectileKinEnergy);
285  G4double CorrectFwdCS = chargeSqRatio*theDirectEMModel->ComputeCrossSectionPerAtom(G4GenericIon::GenericIon(),kinEnergyProjScaled,1,1 ,currentTcutForDirectSecond,1.e20);
286  if (UsedFwdCS >0) new_weight*= CorrectFwdCS/UsedFwdCS;//May be some check is needed if UsedFwdCS ==0 probably that then we should avoid a secondary to be produced,
287 
288 
289  //additional CS crorrection needed for cross section biasing in general.
290  //May be wrong for ions!!! Most of the time not used!
291  G4double w_corr =1./CS_biasing_factor;
293  new_weight*=w_corr;
294 
295  new_weight*=projectileKinEnergy/adjointPrimKinEnergy;
296 
297  fParticleChange->SetParentWeightByProcess(false);
298  fParticleChange->SetSecondaryWeightByProcess(false);
299  fParticleChange->ProposeParentWeight(new_weight);
300 }
G4double GetPostStepWeightCorrection()
G4ParticleDefinition * theDirectPrimaryPartDef
void ProposeParentWeight(G4double finalWeight)
G4Material * currentMaterial
void SetSecondaryWeightByProcess(G4bool)
void SetParentWeightByProcess(G4bool)
G4VEmModel * theDirectEMModel
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:321
virtual G4double GetChargeSquareRatio(const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
Definition: G4VEmModel.cc:353
static G4GenericIon * GenericIon()
Definition: G4GenericIon.cc:93
static constexpr double MeV
Definition: G4SIunits.hh:214
G4double currentTcutForDirectSecond
double G4double
Definition: G4Types.hh:76
static G4AdjointCSManager * GetAdjointCSManager()

Here is the call graph for this function:

Here is the caller graph for this function:

G4double G4AdjointIonIonisationModel::DiffCrossSectionPerAtomPrimToSecond ( G4double  kinEnergyProj,
G4double  kinEnergyProd,
G4double  Z,
G4double  A = 0. 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 159 of file G4AdjointIonIonisationModel.cc.

164 {//Probably that here the Bragg Model should be also used for kinEnergyProj/nuc<2MeV
165 
166 
167 
168  G4double dSigmadEprod=0;
169  G4double Emax_proj = GetSecondAdjEnergyMaxForProdToProjCase(kinEnergyProd);
170  G4double Emin_proj = GetSecondAdjEnergyMinForProdToProjCase(kinEnergyProd);
171 
172  G4double kinEnergyProjScaled = massRatio*kinEnergyProj;
173 
174 
175  if (kinEnergyProj>Emin_proj && kinEnergyProj<=Emax_proj){ //the produced particle should have a kinetic energy smaller than the projectile
176  G4double Tmax=kinEnergyProj;
177 
178  G4double E1=kinEnergyProd;
179  G4double E2=kinEnergyProd*1.000001;
180  G4double dE=(E2-E1);
181  G4double sigma1,sigma2;
182  theDirectEMModel =theBraggIonDirectEMModel;
183  if (kinEnergyProjScaled >2.*MeV && !use_only_bragg) theDirectEMModel = theBetheBlochDirectEMModel; //Bethe Bloch Model
184  sigma1=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,Z,A ,E1,1.e20);
185  sigma2=theDirectEMModel->ComputeCrossSectionPerAtom(theDirectPrimaryPartDef,kinEnergyProj,Z,A ,E2,1.e20);
186 
187  dSigmadEprod=(sigma1-sigma2)/dE;
188 
189  //G4double chargeSqRatio = currentModel->GetChargeSquareRatio(theDirectPrimaryPartDef,currentMaterial,E);
190 
191 
192 
193  if (dSigmadEprod>1.) {
194  G4cout<<"sigma1 "<<kinEnergyProj/MeV<<'\t'<<kinEnergyProd/MeV<<'\t'<<sigma1<<G4endl;
195  G4cout<<"sigma2 "<<kinEnergyProj/MeV<<'\t'<<kinEnergyProd/MeV<<'\t'<<sigma2<<G4endl;
196  G4cout<<"dsigma "<<kinEnergyProj/MeV<<'\t'<<kinEnergyProd/MeV<<'\t'<<dSigmadEprod<<G4endl;
197 
198  }
199 
200 
201 
202 
203 
204 
205  if (theDirectEMModel == theBetheBlochDirectEMModel ){
206  //correction of differential cross section at high energy to correct for the suppression of particle at secondary at high
207  //energy used in the Bethe Bloch Model. This correction consist to multiply by g the probability function used
208  //to test the rejection of a secondary
209  //-------------------------
210 
211  //Source code taken from G4BetheBlochModel::SampleSecondaries
212 
213  G4double deltaKinEnergy = kinEnergyProd;
214 
215  //Part of the taken code
216  //----------------------
217 
218 
219 
220  // projectile formfactor - suppresion of high energy
221  // delta-electron production at high energy
222 
223 
224  G4double x = formfact*deltaKinEnergy;
225  if(x > 1.e-6) {
226  G4double totEnergy = kinEnergyProj + mass;
227  G4double etot2 = totEnergy*totEnergy;
228  G4double beta2 = kinEnergyProj*(kinEnergyProj + 2.0*mass)/etot2;
229  G4double f;
230  G4double f1 = 0.0;
231  f = 1.0 - beta2*deltaKinEnergy/Tmax;
232  if( 0.5 == spin ) {
233  f1 = 0.5*deltaKinEnergy*deltaKinEnergy/etot2;
234  f += f1;
235  }
236  G4double x1 = 1.0 + x;
237  G4double gg = 1.0/(x1*x1);
238  if( 0.5 == spin ) {
239  G4double x2 = 0.5*electron_mass_c2*deltaKinEnergy/(mass*mass);
240  gg *= (1.0 + magMoment2*(x2 - f1/f)/(1.0 + x2));
241  }
242  if(gg > 1.0) {
243  G4cout << "### G4BetheBlochModel in Adjoint Sim WARNING: gg= " << gg
244  << G4endl;
245  gg=1.;
246  }
247  //G4cout<<"gg"<<gg<<G4endl;
248  dSigmadEprod*=gg;
249  }
250  }
251 
252  }
253 
254  return dSigmadEprod;
255 }
G4ParticleDefinition * theDirectPrimaryPartDef
virtual G4double GetSecondAdjEnergyMinForProdToProjCase(G4double PrimAdjEnergy)
static constexpr double electron_mass_c2
static const G4double dE
G4VEmModel * theDirectEMModel
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.cc:321
virtual G4double GetSecondAdjEnergyMaxForProdToProjCase(G4double PrimAdjEnergy)
#define G4endl
Definition: G4ios.hh:61
static constexpr double MeV
Definition: G4SIunits.hh:214
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

G4double G4AdjointIonIonisationModel::GetSecondAdjEnergyMaxForProdToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 356 of file G4AdjointIonIonisationModel.cc.

357 { return HighEnergyLimit;
358 }

Here is the caller graph for this function:

G4double G4AdjointIonIonisationModel::GetSecondAdjEnergyMaxForScatProjToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 344 of file G4AdjointIonIonisationModel.cc.

345 {
346  G4double Tmax=PrimAdjEnergy*one_plus_ratio_2/(one_minus_ratio_2-2.*ratio*PrimAdjEnergy/mass);
347  return Tmax;
348 }
double G4double
Definition: G4Types.hh:76
G4double G4AdjointIonIonisationModel::GetSecondAdjEnergyMinForProdToProjCase ( G4double  PrimAdjEnergy)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 361 of file G4AdjointIonIonisationModel.cc.

362 { G4double Tmin= (2*PrimAdjEnergy-4*mass + std::sqrt(4.*PrimAdjEnergy*PrimAdjEnergy +16.*mass*mass + 8.*PrimAdjEnergy*mass*(1/ratio +ratio)))/4.;
363  return Tmin;
364 }
double G4double
Definition: G4Types.hh:76

Here is the caller graph for this function:

G4double G4AdjointIonIonisationModel::GetSecondAdjEnergyMinForScatProjToProjCase ( G4double  PrimAdjEnergy,
G4double  Tcut = 0 
)
virtual

Reimplemented from G4VEmAdjointModel.

Definition at line 351 of file G4AdjointIonIonisationModel.cc.

352 { return PrimAdjEnergy+Tcut;
353 }
void G4AdjointIonIonisationModel::SampleSecondaries ( const G4Track aTrack,
G4bool  IsScatProjToProjCase,
G4ParticleChange fParticleChange 
)
virtual

Implements G4VEmAdjointModel.

Definition at line 85 of file G4AdjointIonIonisationModel.cc.

88 {
89  const G4DynamicParticle* theAdjointPrimary =aTrack.GetDynamicParticle();
90 
91  //Elastic inverse scattering
92  //---------------------------------------------------------
93  G4double adjointPrimKinEnergy = theAdjointPrimary->GetKineticEnergy();
94  G4double adjointPrimP =theAdjointPrimary->GetTotalMomentum();
95 
96  if (adjointPrimKinEnergy>HighEnergyLimit*0.999){
97  return;
98  }
99 
100  //Sample secondary energy
101  //-----------------------
102  G4double projectileKinEnergy = SampleAdjSecEnergyFromCSMatrix(adjointPrimKinEnergy, IsScatProjToProjCase);
103  CorrectPostStepWeight(fParticleChange,
104  aTrack.GetWeight(),
105  adjointPrimKinEnergy,
106  projectileKinEnergy,
107  IsScatProjToProjCase); //Caution !!!this weight correction should be always applied
108 
109 
110  //Kinematic:
111  //we consider a two body elastic scattering for the forward processes where the projectile knock on an e- at rest and gives
112  // him part of its energy
113  //----------------------------------------------------------------------------------------
114 
116  G4double projectileTotalEnergy = projectileM0+projectileKinEnergy;
117  G4double projectileP2 = projectileTotalEnergy*projectileTotalEnergy - projectileM0*projectileM0;
118 
119 
120 
121  //Companion
122  //-----------
124  if (IsScatProjToProjCase) {
126  }
127  G4double companionTotalEnergy =companionM0+ projectileKinEnergy-adjointPrimKinEnergy;
128  G4double companionP2 = companionTotalEnergy*companionTotalEnergy - companionM0*companionM0;
129 
130 
131  //Projectile momentum
132  //--------------------
133  G4double P_parallel = (adjointPrimP*adjointPrimP + projectileP2 - companionP2)/(2.*adjointPrimP);
134  G4double P_perp = std::sqrt( projectileP2 - P_parallel*P_parallel);
135  G4ThreeVector dir_parallel=theAdjointPrimary->GetMomentumDirection();
136  G4double phi =G4UniformRand()*2.*3.1415926;
137  G4ThreeVector projectileMomentum = G4ThreeVector(P_perp*std::cos(phi),P_perp*std::sin(phi),P_parallel);
138  projectileMomentum.rotateUz(dir_parallel);
139 
140 
141 
142  if (!IsScatProjToProjCase ){ //kill the primary and add a secondary
143  fParticleChange->ProposeTrackStatus(fStopAndKill);
144  fParticleChange->AddSecondary(new G4DynamicParticle(theAdjEquivOfDirectPrimPartDef,projectileMomentum));
145  //G4cout<<"projectileMomentum "<<projectileMomentum<<G4endl;
146  }
147  else {
148  fParticleChange->ProposeEnergy(projectileKinEnergy);
149  fParticleChange->ProposeMomentumDirection(projectileMomentum.unit());
150  }
151 
152 
153 
154 
155 }
G4double GetKineticEnergy() const
CLHEP::Hep3Vector G4ThreeVector
const G4DynamicParticle * GetDynamicParticle() const
G4double SampleAdjSecEnergyFromCSMatrix(size_t MatrixIndex, G4double prim_energy, G4bool IsScatProjToProjCase)
virtual void CorrectPostStepWeight(G4ParticleChange *fParticleChange, G4double old_weight, G4double adjointPrimKinEnergy, G4double projectileKinEnergy, G4bool IsScatProjToProjCase)
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
G4double GetTotalMomentum() const
#define G4UniformRand()
Definition: Randomize.hh:97
const G4ThreeVector & GetMomentumDirection() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:38
G4double GetPDGMass() const
Hep3Vector unit() const
void ProposeEnergy(G4double finalEnergy)
void AddSecondary(G4Track *aSecondary)
G4double GetWeight() const
G4ParticleDefinition * theAdjEquivOfDirectSecondPartDef
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
double G4double
Definition: G4Types.hh:76
void ProposeTrackStatus(G4TrackStatus status)

Here is the call graph for this function:

void G4AdjointIonIonisationModel::SetIon ( G4ParticleDefinition adj_ion,
G4ParticleDefinition fwd_ion 
)

Definition at line 258 of file G4AdjointIonIonisationModel.cc.

259 { theDirectPrimaryPartDef =fwd_ion;
261 
262  DefineProjectileProperty();
263 }
G4ParticleDefinition * theDirectPrimaryPartDef
G4ParticleDefinition * theAdjEquivOfDirectPrimPartDef
void G4AdjointIonIonisationModel::SetUseOnlyBragg ( G4bool  aBool)
inline

Definition at line 106 of file G4AdjointIonIonisationModel.hh.

106 {use_only_bragg =aBool;}

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