Geant4  10.02.p03
G4eSingleCoulombScatteringModel Class Reference

#include <G4eSingleCoulombScatteringModel.hh>

Inheritance diagram for G4eSingleCoulombScatteringModel:
Collaboration diagram for G4eSingleCoulombScatteringModel:

Public Member Functions

 G4eSingleCoulombScatteringModel (const G4String &nam="eSingleCoulombScat")
 
virtual ~G4eSingleCoulombScatteringModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A, G4double cut, G4double emax)
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle *> *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy)
 
void SetRecoilThreshold (G4double eth)
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void InitialiseForMaterial (const G4ParticleDefinition *, const G4Material *)
 
virtual void InitialiseForElement (const G4ParticleDefinition *, G4int Z)
 
virtual G4double ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
virtual G4double CrossSectionPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double GetPartialCrossSection (const G4Material *, G4int, const G4ParticleDefinition *, G4double)
 
virtual G4double ComputeCrossSectionPerShell (const G4ParticleDefinition *, G4int Z, G4int shellIdx, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
virtual G4double ChargeSquareRatio (const G4Track &)
 
virtual G4double GetChargeSquareRatio (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual G4double GetParticleCharge (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void StartTracking (G4Track *)
 
virtual void CorrectionsAlongStep (const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double &eloss, G4double &niel, G4double length)
 
virtual G4double Value (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double cut=0.0)
 
virtual G4double MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *)
 
virtual void SetupForMaterial (const G4ParticleDefinition *, const G4Material *, G4double kineticEnergy)
 
virtual void DefineForRegion (const G4Region *)
 
virtual void ModelDescription (std::ostream &outFile) const
 
void InitialiseElementSelectors (const G4ParticleDefinition *, const G4DataVector &)
 
std::vector< G4EmElementSelector * > * GetElementSelectors ()
 
void SetElementSelectors (std::vector< G4EmElementSelector *> *)
 
G4double ComputeDEDX (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=DBL_MAX)
 
G4double CrossSection (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeMeanFreePath (const G4ParticleDefinition *, G4double kineticEnergy, const G4Material *, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, const G4Element *, G4double kinEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectIsotopeNumber (const G4Element *)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectRandomAtom (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
G4int SelectRandomAtomNumber (const G4Material *)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=0)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
void SetAngularDistribution (G4VEmAngularDistribution *)
 
G4double HighEnergyLimit () const
 
G4double LowEnergyLimit () const
 
G4double HighEnergyActivationLimit () const
 
G4double LowEnergyActivationLimit () const
 
G4double PolarAngleLimit () const
 
G4double SecondaryThreshold () const
 
G4bool LPMFlag () const
 
G4bool DeexcitationFlag () const
 
G4bool ForceBuildTableFlag () const
 
G4bool UseAngularGeneratorFlag () const
 
void SetAngularGeneratorFlag (G4bool)
 
void SetHighEnergyLimit (G4double)
 
void SetLowEnergyLimit (G4double)
 
void SetActivationHighEnergyLimit (G4double)
 
void SetActivationLowEnergyLimit (G4double)
 
G4bool IsActive (G4double kinEnergy)
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
G4double MaxSecondaryKinEnergy (const G4DynamicParticle *dynParticle)
 
const G4StringGetName () const
 
void SetCurrentCouple (const G4MaterialCutsCouple *)
 
const G4ElementGetCurrentElement () const
 
const G4IsotopeGetCurrentIsotope () const
 
G4bool IsLocked () const
 
void SetLocked (G4bool)
 

Private Member Functions

void DefineMaterial (const G4MaterialCutsCouple *)
 
void SetupParticle (const G4ParticleDefinition *)
 
G4eSingleCoulombScatteringModeloperator= (const G4eSingleCoulombScatteringModel &right)
 
 G4eSingleCoulombScatteringModel (const G4eSingleCoulombScatteringModel &)
 

Private Attributes

G4IonTabletheIonTable
 
G4ParticleChangeForGamma * fParticleChange
 
G4NistManagerfNistManager
 
G4ScreeningMottCrossSectionMottcross
 
const std::vector< G4double > * pCuts
 
const G4MaterialCutsCouplecurrentCouple
 
const G4MaterialcurrentMaterial
 
const G4ElementcurrentElement
 
G4int currentMaterialIndex
 
G4double cosThetaMin
 
G4double recoilThreshold
 
const G4ParticleDefinitionparticle
 
G4double mass
 
G4double lowEnergyLimit
 

Additional Inherited Members

- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLoss * GetParticleChangeForLoss ()
 
G4ParticleChangeForGamma * GetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 
- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChange * pParticleChange
 
G4PhysicsTablexSectionTable
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 
- Static Protected Attributes inherited from G4VEmModel
static const G4double inveplus = 1.0/CLHEP::eplus
 

Detailed Description

Definition at line 66 of file G4eSingleCoulombScatteringModel.hh.

Constructor & Destructor Documentation

◆ G4eSingleCoulombScatteringModel() [1/2]

G4eSingleCoulombScatteringModel::G4eSingleCoulombScatteringModel ( const G4String nam = "eSingleCoulombScat")

Definition at line 75 of file G4eSingleCoulombScatteringModel.cc.

76  : G4VEmModel(nam),
77  cosThetaMin(1.0)
78 {
81  fParticleChange = 0;
82 
83  pCuts=0;
84  currentMaterial = 0;
85  currentElement = 0;
86  currentCouple = 0;
87 
88  lowEnergyLimit = 0*keV;
89  recoilThreshold = 0.*eV;
90  particle = 0;
91  mass=0;
93 
95 }
static G4NistManager * Instance()
G4VEmModel(const G4String &nam)
Definition: G4VEmModel.cc:69
G4IonTable * GetIonTable() const
static const double eV
Definition: G4SIunits.hh:212
static G4ParticleTable * GetParticleTable()
static const double keV
Definition: G4SIunits.hh:213
Here is the call graph for this function:

◆ ~G4eSingleCoulombScatteringModel()

G4eSingleCoulombScatteringModel::~G4eSingleCoulombScatteringModel ( )
virtual

Definition at line 99 of file G4eSingleCoulombScatteringModel.cc.

100 {
101  delete Mottcross;
102 }

◆ G4eSingleCoulombScatteringModel() [2/2]

G4eSingleCoulombScatteringModel::G4eSingleCoulombScatteringModel ( const G4eSingleCoulombScatteringModel )
private

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4eSingleCoulombScatteringModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition p,
G4double  kinEnergy,
G4double  Z,
G4double  A,
G4double  cut,
G4double  emax 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 141 of file G4eSingleCoulombScatteringModel.cc.

148 {
149  SetupParticle(p);
150 
151  G4double cross =0.0;
152  if(kinEnergy < lowEnergyLimit) return cross;
153 
155 
156  //Total Cross section
157  Mottcross->SetupKinematic(kinEnergy, Z);
158  cross = Mottcross->NuclearCrossSection();
159 
160  //cout<< "Compute Cross Section....cross "<<G4BestUnit(cross,"Surface") << " cm2 "<< cross/cm2 <<" Z: "<<Z<<" kinEnergy: "<<kinEnergy<<endl;
161  return cross;
162 }
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:452
Float_t Z
void DefineMaterial(const G4MaterialCutsCouple *)
void SetupParticle(const G4ParticleDefinition *)
double G4double
Definition: G4Types.hh:76
void SetupKinematic(G4double kinEnergy, G4double Z)
Here is the call graph for this function:

◆ DefineMaterial()

void G4eSingleCoulombScatteringModel::DefineMaterial ( const G4MaterialCutsCouple cup)
inlineprivate

Definition at line 129 of file G4eSingleCoulombScatteringModel.hh.

130 {
131  if(cup != currentCouple) {
132  currentCouple = cup;
133  currentMaterial = cup->GetMaterial();
135  }
136 }
const G4Material * GetMaterial() const
Here is the call graph for this function:
Here is the caller graph for this function:

◆ Initialise()

void G4eSingleCoulombScatteringModel::Initialise ( const G4ParticleDefinition p,
const G4DataVector cuts 
)
virtual

Implements G4VEmModel.

Definition at line 106 of file G4eSingleCoulombScatteringModel.cc.

108 {
109  SetupParticle(p);
110  currentCouple = 0;
112  //cosThetaMin = cos(PolarAngleLimit());
114 
115  pCuts = &cuts;
116  //G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3);
117 
118  /*
119  G4cout << "!!! G4eSingleCoulombScatteringModel::Initialise for "
120  << part->GetParticleName() << " cos(TetMin)= " << cosThetaMin
121  << " cos(TetMax)= " << cosThetaMax <<G4endl;
122  G4cout << "cut= " << (*pCuts)[0] << " cut1= " << (*pCuts)[1] << G4endl;
123  */
124 
125  if(!fParticleChange) {
127  }
128 
129  if(IsMaster()) {
131  }
132 }
void Initialise(const G4ParticleDefinition *, G4double cosThetaLim)
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4VEmModel.cc:148
G4bool IsMaster() const
Definition: G4VEmModel.hh:718
void SetupParticle(const G4ParticleDefinition *)
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:134
Here is the call graph for this function:

◆ InitialiseLocal()

void G4eSingleCoulombScatteringModel::InitialiseLocal ( const G4ParticleDefinition ,
G4VEmModel masterModel 
)
virtual

Reimplemented from G4VEmModel.

Definition at line 134 of file G4eSingleCoulombScatteringModel.cc.

136 {
138 }
void SetElementSelectors(std::vector< G4EmElementSelector *> *)
Definition: G4VEmModel.hh:810
std::vector< G4EmElementSelector * > * GetElementSelectors()
Definition: G4VEmModel.hh:802
Here is the call graph for this function:

◆ operator=()

G4eSingleCoulombScatteringModel& G4eSingleCoulombScatteringModel::operator= ( const G4eSingleCoulombScatteringModel right)
private

◆ SampleSecondaries()

void G4eSingleCoulombScatteringModel::SampleSecondaries ( std::vector< G4DynamicParticle *> *  fvect,
const G4MaterialCutsCouple couple,
const G4DynamicParticle dp,
G4double  tmin,
G4double  maxEnergy 
)
virtual

Implements G4VEmModel.

Definition at line 166 of file G4eSingleCoulombScatteringModel.cc.

172 {
173  G4double kinEnergy = dp->GetKineticEnergy();
174  //cout<<"--- kinEnergy "<<kinEnergy<<endl;
175 
176  if(kinEnergy < lowEnergyLimit) return;
177 
178  DefineMaterial(couple);
180 
181  // Choose nucleus
182  //last two :cutEnergy= min e kinEnergy=max
184  kinEnergy,cutEnergy,kinEnergy);
185 
186  G4double Z = currentElement->GetZ();
187  G4int iz = G4int(Z);
188  G4int ia = SelectIsotopeNumber(currentElement);
190 
191  //G4cout<<"..Z: "<<Z<<" ..iz: "<<iz<<" ..ia: "<<ia<<" ..mass2: "<<mass2<<G4endl;
192 
193  Mottcross->SetupKinematic(kinEnergy, Z);
194  G4double cross= Mottcross->NuclearCrossSection(); //MODIFY TO LOAD TABLE
195  if(cross == 0.0) { return; }
196  //cout<< "Energy: "<<kinEnergy/MeV<<" Z: "<<Z<<"....cross "<<G4BestUnit(cross,"Surface") << " cm2 "<< cross/cm2 <<endl;
197 
199  G4double sint = sin(z1);
200  G4double cost = sqrt(1.0 - sint*sint);
201  G4double phi = twopi* G4UniformRand();
202 
203  // kinematics in the Lab system
204  G4double ptot = dp->GetTotalMomentum();
205  G4double e1 = dp->GetTotalEnergy();
206  // Lab. system kinematics along projectile direction
207  G4LorentzVector v0 = G4LorentzVector(0, 0, ptot, e1);
208  G4double bet = ptot/(v0.e() + mass2);
209  G4double gam = 1.0/sqrt((1.0 - bet)*(1.0 + bet));
210 
211  //CM Projectile
212  G4double momCM = gam*(ptot - bet*e1);
213  G4double eCM = gam*(e1 - bet*ptot);
214  //energy & momentum after scattering of incident particle
215  G4double pxCM = momCM*sint*cos(phi);
216  G4double pyCM = momCM*sint*sin(phi);
217  G4double pzCM = momCM*cost;
218 
219  //CM--->Lab
220  G4LorentzVector v1(pxCM , pyCM, gam*(pzCM + bet*eCM), gam*(eCM + bet*pzCM));
221 
222  // Rotate to global system
224  G4ThreeVector newDirection = v1.vect().unit();
225  newDirection.rotateUz(dir);
226 
227  fParticleChange->ProposeMomentumDirection(newDirection);
228 
229  // recoil
230  v0 -= v1;
231  G4double trec = v0.e();
232  G4double edep = 0.0;
233 
234  G4double tcut = recoilThreshold;
235 
236  //G4cout<<" Energy Transfered: "<<trec/eV<<G4endl;
237 
238  if(pCuts) {
239  tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]);
240  }
241 
242  if(trec > tcut) {
244  newDirection = v0.vect().unit();
245  newDirection.rotateUz(dir);
246  G4DynamicParticle* newdp = new G4DynamicParticle(ion, newDirection, trec);
247  fvect->push_back(newdp);
248  } else if(trec > 0.0) {
249  edep = trec;
250  fParticleChange->ProposeNonIonizingEnergyDeposit(edep);
251  }
252 
253  // finelize primary energy and energy balance
254  G4double finalT = v1.e() - mass;
255  //G4cout<<"Final Energy: "<<finalT/eV<<G4endl;
256  if(finalT <= lowEnergyLimit) {
257  edep += finalT;
258  finalT = 0.0;
259  }
260  edep = std::max(edep, 0.0);
261  fParticleChange->SetProposedKineticEnergy(finalT);
262  fParticleChange->ProposeLocalEnergyDeposit(edep);
263 
264 }
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:491
G4double GetTotalMomentum() const
TDirectory * dir
int G4int
Definition: G4Types.hh:78
Double_t edep
G4double GetTotalEnergy() const
Hep3Vector vect() const
G4double GetKineticEnergy() const
#define G4UniformRand()
Definition: Randomize.hh:97
Float_t Z
G4double iz
Definition: TRTMaterials.hh:39
Hep3Vector unit() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:38
static const double twopi
Definition: G4SIunits.hh:75
G4int SelectIsotopeNumber(const G4Element *)
Definition: G4VEmModel.hh:585
static const G4double e1
void DefineMaterial(const G4MaterialCutsCouple *)
void SetupParticle(const G4ParticleDefinition *)
const G4ThreeVector & GetMomentumDirection() const
G4ParticleDefinition * GetDefinition() const
double G4double
Definition: G4Types.hh:76
void SetupKinematic(G4double kinEnergy, G4double Z)
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:544
CLHEP::HepLorentzVector G4LorentzVector
Here is the call graph for this function:

◆ SetRecoilThreshold()

void G4eSingleCoulombScatteringModel::SetRecoilThreshold ( G4double  eth)
inline

Definition at line 151 of file G4eSingleCoulombScatteringModel.hh.

◆ SetupParticle()

void G4eSingleCoulombScatteringModel::SetupParticle ( const G4ParticleDefinition p)
inlineprivate

Definition at line 141 of file G4eSingleCoulombScatteringModel.hh.

142 {
143  if(p != particle) {
144  particle = p;
145  mass = particle->GetPDGMass();
147  }
148 }
void SetupParticle(const G4ParticleDefinition *)
Here is the call graph for this function:
Here is the caller graph for this function:

Member Data Documentation

◆ cosThetaMin

G4double G4eSingleCoulombScatteringModel::cosThetaMin
private

Definition at line 116 of file G4eSingleCoulombScatteringModel.hh.

◆ currentCouple

const G4MaterialCutsCouple* G4eSingleCoulombScatteringModel::currentCouple
private

Definition at line 111 of file G4eSingleCoulombScatteringModel.hh.

◆ currentElement

const G4Element* G4eSingleCoulombScatteringModel::currentElement
private

Definition at line 113 of file G4eSingleCoulombScatteringModel.hh.

◆ currentMaterial

const G4Material* G4eSingleCoulombScatteringModel::currentMaterial
private

Definition at line 112 of file G4eSingleCoulombScatteringModel.hh.

◆ currentMaterialIndex

G4int G4eSingleCoulombScatteringModel::currentMaterialIndex
private

Definition at line 114 of file G4eSingleCoulombScatteringModel.hh.

◆ fNistManager

G4NistManager* G4eSingleCoulombScatteringModel::fNistManager
private

Definition at line 107 of file G4eSingleCoulombScatteringModel.hh.

◆ fParticleChange

G4ParticleChangeForGamma* G4eSingleCoulombScatteringModel::fParticleChange
private

Definition at line 106 of file G4eSingleCoulombScatteringModel.hh.

◆ lowEnergyLimit

G4double G4eSingleCoulombScatteringModel::lowEnergyLimit
private

Definition at line 122 of file G4eSingleCoulombScatteringModel.hh.

◆ mass

G4double G4eSingleCoulombScatteringModel::mass
private

Definition at line 121 of file G4eSingleCoulombScatteringModel.hh.

◆ Mottcross

G4ScreeningMottCrossSection* G4eSingleCoulombScatteringModel::Mottcross
private

Definition at line 108 of file G4eSingleCoulombScatteringModel.hh.

◆ particle

const G4ParticleDefinition* G4eSingleCoulombScatteringModel::particle
private

Definition at line 120 of file G4eSingleCoulombScatteringModel.hh.

◆ pCuts

const std::vector<G4double>* G4eSingleCoulombScatteringModel::pCuts
private

Definition at line 110 of file G4eSingleCoulombScatteringModel.hh.

◆ recoilThreshold

G4double G4eSingleCoulombScatteringModel::recoilThreshold
private

Definition at line 117 of file G4eSingleCoulombScatteringModel.hh.

◆ theIonTable

G4IonTable* G4eSingleCoulombScatteringModel::theIonTable
private

Definition at line 105 of file G4eSingleCoulombScatteringModel.hh.


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