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

#include <G4hCoulombScatteringModel.hh>

+ Inheritance diagram for G4hCoulombScatteringModel:

Public Member Functions

 G4hCoulombScatteringModel (G4bool combined=true)
 
virtual ~G4hCoulombScatteringModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &) override
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel) override
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A, G4double cut, G4double emax) override
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) override
 
virtual G4double MinPrimaryEnergy (const G4Material *, const G4ParticleDefinition *, G4double) final
 
void SetLowEnergyThreshold (G4double val)
 
void SetRecoilThreshold (G4double eth)
 
void SetFixedCut (G4double)
 
G4double GetFixedCut () const
 
- Public Member Functions inherited from G4VEmModel
 G4VEmModel (const G4String &nam)
 
virtual ~G4VEmModel ()
 
virtual void Initialise (const G4ParticleDefinition *, const G4DataVector &)=0
 
virtual void SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin=0.0, G4double tmax=DBL_MAX)=0
 
virtual void InitialiseLocal (const G4ParticleDefinition *, G4VEmModel *masterModel)
 
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 level, const G4ParticleDefinition *, G4double kineticEnergy)
 
virtual G4double ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0., G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
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 * > *)
 
virtual 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)
 
const G4ElementSelectRandomAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
 
const G4ElementSelectTargetAtom (const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double logKineticEnergy, 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 *)
 
G4int SelectIsotopeNumber (const G4Element *)
 
void SetParticleChange (G4VParticleChange *, G4VEmFluctuationModel *f=nullptr)
 
void SetCrossSectionTable (G4PhysicsTable *, G4bool isLocal)
 
G4ElementDataGetElementData ()
 
G4PhysicsTableGetCrossSectionTable ()
 
G4VEmFluctuationModelGetModelOfFluctuations ()
 
G4VEmAngularDistributionGetAngularDistribution ()
 
G4VEmModelGetTripletModel ()
 
void SetTripletModel (G4VEmModel *)
 
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) const
 
void SetPolarAngleLimit (G4double)
 
void SetSecondaryThreshold (G4double)
 
void SetLPMFlag (G4bool val)
 
void SetDeexcitationFlag (G4bool val)
 
void SetForceBuildTable (G4bool val)
 
void SetFluctuationFlag (G4bool val)
 
void SetMasterThread (G4bool val)
 
G4bool IsMaster () const
 
void SetUseBaseMaterials (G4bool val)
 
G4bool UseBaseMaterials () 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)
 
G4VEmModeloperator= (const G4VEmModel &right)=delete
 
 G4VEmModel (const G4VEmModel &)=delete
 

Protected Member Functions

void DefineMaterial (const G4MaterialCutsCouple *)
 
void SetupParticle (const G4ParticleDefinition *)
 
- Protected Member Functions inherited from G4VEmModel
G4ParticleChangeForLossGetParticleChangeForLoss ()
 
G4ParticleChangeForGammaGetParticleChangeForGamma ()
 
virtual G4double MaxSecondaryEnergy (const G4ParticleDefinition *, G4double kineticEnergy)
 
const G4MaterialCutsCoupleCurrentCouple () const
 
void SetCurrentElement (const G4Element *)
 

Additional Inherited Members

- Protected Attributes inherited from G4VEmModel
G4ElementDatafElementData
 
G4VParticleChangepParticleChange
 
G4PhysicsTablexSectionTable
 
const G4MaterialpBaseMaterial
 
const std::vector< G4double > * theDensityFactor
 
const std::vector< G4int > * theDensityIdx
 
size_t idxTable
 
G4bool lossFlucFlag
 
G4double inveplus
 
G4double pFactor
 

Detailed Description

Definition at line 66 of file G4hCoulombScatteringModel.hh.

Constructor & Destructor Documentation

◆ G4hCoulombScatteringModel()

G4hCoulombScatteringModel::G4hCoulombScatteringModel ( G4bool  combined = true)
explicit

Definition at line 66 of file G4hCoulombScatteringModel.cc.

67 : G4VEmModel("hCoulombScattering"),
68 cosThetaMin(1.0),
69 cosThetaMax(-1.0),
70 isCombined(combined)
71{
72 fParticleChange = nullptr;
73 fNistManager = G4NistManager::Instance();
75 theProton = G4Proton::Proton();
76 currentMaterial = nullptr;
77 fixedCut = -1.0;
78
79 pCuts = nullptr;
80
81 recoilThreshold = 0.0; // by default does not work
82
83 particle = nullptr;
84 currentCouple = nullptr;
85 wokvi = new G4WentzelVIRelXSection();
86
87 currentMaterialIndex = 0;
88 mass = CLHEP::proton_mass_c2;
89 elecRatio = 0.0;
90}
static G4NistManager * Instance()
G4IonTable * GetIonTable() const
static G4ParticleTable * GetParticleTable()
static G4Proton * Proton()
Definition: G4Proton.cc:92

◆ ~G4hCoulombScatteringModel()

G4hCoulombScatteringModel::~G4hCoulombScatteringModel ( )
virtual

Definition at line 94 of file G4hCoulombScatteringModel.cc.

95{
96 delete wokvi;
97}

Member Function Documentation

◆ ComputeCrossSectionPerAtom()

G4double G4hCoulombScatteringModel::ComputeCrossSectionPerAtom ( const G4ParticleDefinition p,
G4double  kinEnergy,
G4double  Z,
G4double  A,
G4double  cut,
G4double  emax 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 182 of file G4hCoulombScatteringModel.cc.

187{
188 //G4cout << "### G4hCoulombScatteringModel::ComputeCrossSectionPerAtom for "
189 //<< p->GetParticleName()<<" Z= "<<Z<<" e(MeV)= "<< kinEnergy/MeV << G4endl;
190 G4double cross = 0.0;
191 elecRatio = 0.0;
192 if(p != particle) { SetupParticle(p); }
193
194 // cross section is set to zero to avoid problems in sample secondary
195 if(kinEnergy <= 0.0) { return cross; }
197
198 G4int iz = G4lrint(Z);
199 G4double tmass = (1 == iz) ? proton_mass_c2 :
200 fNistManager->GetAtomicMassAmu(iz)*amu_c2;
201 wokvi->SetTargetMass(tmass);
202
203 G4double costmin =
204 wokvi->SetupKinematic(kinEnergy, currentMaterial);
205
206 if(cosThetaMax < costmin) {
207 G4double cut = (0.0 < fixedCut) ? fixedCut : cutEnergy;
208 costmin = wokvi->SetupTarget(iz, cut);
209 G4double costmax =
210 (1 == iz && particle == theProton && cosThetaMax < 0.0)
211 ? 0.0 : cosThetaMax;
212 if(costmin > costmax) {
213 cross = wokvi->ComputeNuclearCrossSection(costmin, costmax)
214 + wokvi->ComputeElectronCrossSection(costmin, costmax);
215 }
216 /*
217 if(p->GetParticleName() == "mu+")
218 G4cout << "e(MeV)= " << kinEnergy/MeV << " cross(b)= " << cross/barn
219 << " 1-costmin= " << 1-costmin
220 << " 1-costmax= " << 1-costmax
221 << " 1-cosThetaMax= " << 1-cosThetaMax
222 << G4endl;
223 */
224 }
225 return cross;
226}
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
G4double GetAtomicMassAmu(const G4String &symb) const
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:480
G4double SetupTarget(G4int Z, G4double cut)
G4double ComputeElectronCrossSection(G4double CosThetaMin, G4double CosThetaMax)
G4double ComputeNuclearCrossSection(G4double CosThetaMin, G4double CosThetaMax)
virtual G4double SetupKinematic(G4double kinEnergy, const G4Material *mat)
void SetupParticle(const G4ParticleDefinition *)
void DefineMaterial(const G4MaterialCutsCouple *)
int G4lrint(double ad)
Definition: templates.hh:134

◆ DefineMaterial()

void G4hCoulombScatteringModel::DefineMaterial ( const G4MaterialCutsCouple cup)
inlineprotected

Definition at line 153 of file G4hCoulombScatteringModel.hh.

154{
155 if(cup != currentCouple) {
156 currentCouple = cup;
157 currentMaterial = cup->GetMaterial();
158 currentMaterialIndex = currentCouple->GetIndex();
159 }
160}
const G4Material * GetMaterial() const

Referenced by ComputeCrossSectionPerAtom(), and SampleSecondaries().

◆ GetFixedCut()

G4double G4hCoulombScatteringModel::GetFixedCut ( ) const
inline

Definition at line 191 of file G4hCoulombScatteringModel.hh.

192{
193 return fixedCut;
194}

◆ Initialise()

void G4hCoulombScatteringModel::Initialise ( const G4ParticleDefinition part,
const G4DataVector cuts 
)
overridevirtual

Implements G4VEmModel.

Definition at line 101 of file G4hCoulombScatteringModel.cc.

103{
104 SetupParticle(part);
105 currentCouple = nullptr;
106
107 // defined theta limit between single and multiple scattering
108 isCombined = true;
110
111 if(tet <= 0.0) {
112 cosThetaMin = 1.0;
113 isCombined = false;
114 } else if(tet >= CLHEP::pi) {
115 cosThetaMin = -1.0;
116 } else {
117 cosThetaMin = cos(tet);
118 }
119
120 wokvi->Initialise(part, cosThetaMin);
121 /*
122 G4cout << "G4hCoulombScatteringModel: " << particle->GetParticleName()
123 << " 1-cos(ThetaLimit)= " << 1 - cosThetaMin
124 << " cos(thetaMax)= " << cosThetaMax
125 << G4endl;
126 */
127 pCuts = &cuts;
128 //G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3);
129 /*
130 G4cout << "!!! G4hCoulombScatteringModel::Initialise for "
131 << part->GetParticleName() << " cos(TetMin)= " << cosThetaMin
132 << " cos(TetMax)= " << cosThetaMax <<G4endl;
133 G4cout << "cut= " << (*pCuts)[0] << " cut1= " << (*pCuts)[1] << G4endl;
134 */
135 if(!fParticleChange) {
136 fParticleChange = GetParticleChangeForGamma();
137 }
138 if(IsMaster() && mass < GeV && part->GetParticleName() != "GenericIon") {
139 InitialiseElementSelectors(part, cuts);
140 }
141}
G4double PolarAngleLimit() const
Definition: G4VEmModel.hh:673
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:133
G4bool IsMaster() const
Definition: G4VEmModel.hh:736
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4VEmModel.cc:148
void Initialise(const G4ParticleDefinition *, G4double CosThetaLim)

◆ InitialiseLocal()

void G4hCoulombScatteringModel::InitialiseLocal ( const G4ParticleDefinition ,
G4VEmModel masterModel 
)
overridevirtual

Reimplemented from G4VEmModel.

Definition at line 145 of file G4hCoulombScatteringModel.cc.

147{
149}
void SetElementSelectors(std::vector< G4EmElementSelector * > *)
Definition: G4VEmModel.hh:842
std::vector< G4EmElementSelector * > * GetElementSelectors()
Definition: G4VEmModel.hh:834

◆ MinPrimaryEnergy()

G4double G4hCoulombScatteringModel::MinPrimaryEnergy ( const G4Material material,
const G4ParticleDefinition part,
G4double   
)
finalvirtual

Reimplemented from G4VEmModel.

Definition at line 154 of file G4hCoulombScatteringModel.cc.

157{
158 SetupParticle(part);
159
160 // define cut using cuts for proton
161 G4double cut =
162 std::max(recoilThreshold, (*pCuts)[CurrentCouple()->GetIndex()]);
163
164 // find out lightest element
165 const G4ElementVector* theElementVector = material->GetElementVector();
166 G4int nelm = material->GetNumberOfElements();
167
168 // select lightest element
169 G4int Z = 300;
170 for (G4int j=0; j<nelm; ++j) {
171 Z = std::min(Z,(*theElementVector)[j]->GetZasInt());
172 }
173 G4int A = G4lrint(fNistManager->GetAtomicMassAmu(Z));
175 G4double t = std::max(cut, 0.5*(cut + sqrt(2*cut*targetMass)));
176
177 return t;
178}
double A(double temperature)
std::vector< G4Element * > G4ElementVector
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
static G4double GetNuclearMass(const G4double A, const G4double Z)

◆ SampleSecondaries()

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

Implements G4VEmModel.

Definition at line 230 of file G4hCoulombScatteringModel.cc.

236{
237 G4double kinEnergy = dp->GetKineticEnergy();
239 DefineMaterial(couple);
240
241 // Choose nucleus
242 G4double cut = (0.0 < fixedCut) ? fixedCut : cutEnergy;
243
244 const G4Element* elm = SelectRandomAtom(couple,particle,
245 kinEnergy,cut,kinEnergy);
246
247 G4int iz = elm->GetZasInt();
248 G4int ia = SelectIsotopeNumber(elm);
250
251 wokvi->SetTargetMass(mass2);
252 wokvi->SetupKinematic(kinEnergy, currentMaterial);
253 G4double costmin = wokvi->SetupTarget(iz, cut);
254 G4double costmax = (1 == iz && particle == theProton && cosThetaMax < 0.0)
255 ? 0.0 : cosThetaMax;
256 if(costmin <= costmax) { return; }
257
258 G4double cross = wokvi->ComputeNuclearCrossSection(costmin, costmax);
259 G4double ecross = wokvi->ComputeElectronCrossSection(costmin, costmax);
260 G4double ratio = ecross/(cross + ecross);
261
262 G4ThreeVector newDirection =
263 wokvi->SampleSingleScattering(costmin, costmax, ratio);
264
265 // kinematics in the Lab system
266 G4double ptot = sqrt(kinEnergy*(kinEnergy + 2.0*mass));
267 G4double e1 = mass + kinEnergy;
268
269 // Lab. system kinematics along projectile direction
270 G4LorentzVector v0 = G4LorentzVector(0, 0, ptot, e1+mass2);
271 G4LorentzVector v1 = G4LorentzVector(0, 0, ptot, e1);
272 G4ThreeVector bst = v0.boostVector();
273 v1.boost(-bst);
274 // CM projectile
275 G4double momCM = v1.pz();
276
277 // Momentum after scattering of incident particle
278 v1.setX(momCM*newDirection.x());
279 v1.setY(momCM*newDirection.y());
280 v1.setZ(momCM*newDirection.z());
281
282 // CM--->Lab
283 v1.boost(bst);
284
286 newDirection = v1.vect().unit();
287 newDirection.rotateUz(dir);
288
289 fParticleChange->ProposeMomentumDirection(newDirection);
290
291 // recoil
292 v0 -= v1;
293 G4double trec = std::max(v0.e() - mass2, 0.0);
294 G4double edep = 0.0;
295
296 G4double tcut = recoilThreshold;
297 if(pCuts) { tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]); }
298
299 if(trec > tcut) {
300 G4ParticleDefinition* ion = theIonTable->GetIon(iz, ia, 0);
301 newDirection = v0.vect().unit();
302 newDirection.rotateUz(dir);
303 G4DynamicParticle* newdp = new G4DynamicParticle(ion, newDirection, trec);
304 fvect->push_back(newdp);
305 } else if(trec > 0.0) {
306 edep = trec;
307 fParticleChange->ProposeNonIonizingEnergyDeposit(edep);
308 }
309
310 // finelize primary energy and energy balance
311 G4double finalT = v1.e() - mass;
312 if(finalT < 0.0) {
313 edep += finalT;
314 finalT = 0.0;
315 }
316 edep = std::max(edep, 0.0);
317 fParticleChange->SetProposedKineticEnergy(finalT);
318 fParticleChange->ProposeLocalEnergyDeposit(edep);
319}
CLHEP::HepLorentzVector G4LorentzVector
double z() const
Hep3Vector unit() const
double x() const
double y() const
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:33
Hep3Vector boostVector() const
HepLorentzVector & boost(double, double, double)
Hep3Vector vect() const
const G4ThreeVector & GetMomentumDirection() const
G4ParticleDefinition * GetDefinition() const
G4double GetKineticEnergy() const
G4int GetZasInt() const
Definition: G4Element.hh:131
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:522
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:570
G4int SelectIsotopeNumber(const G4Element *)
Definition: G4VEmModel.cc:337
void ProposeNonIonizingEnergyDeposit(G4double anEnergyPart)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
G4ThreeVector & SampleSingleScattering(G4double CosThetaMin, G4double CosThetaMax, G4double elecRatio)

◆ SetFixedCut()

void G4hCoulombScatteringModel::SetFixedCut ( G4double  val)
inline

Definition at line 184 of file G4hCoulombScatteringModel.hh.

185{
186 fixedCut = val;
187}

◆ SetLowEnergyThreshold()

void G4hCoulombScatteringModel::SetLowEnergyThreshold ( G4double  val)
inline

◆ SetRecoilThreshold()

void G4hCoulombScatteringModel::SetRecoilThreshold ( G4double  eth)
inline

Definition at line 177 of file G4hCoulombScatteringModel.hh.

178{
179 recoilThreshold = eth;
180}

◆ SetupParticle()

void G4hCoulombScatteringModel::SetupParticle ( const G4ParticleDefinition p)
inlineprotected

Definition at line 165 of file G4hCoulombScatteringModel.hh.

166{
167 // Initialise mass and charge
168 if(p != particle) {
169 particle = p;
170 mass = particle->GetPDGMass();
171 wokvi->SetupParticle(p);
172 }
173}
void SetupParticle(const G4ParticleDefinition *)

Referenced by ComputeCrossSectionPerAtom(), Initialise(), MinPrimaryEnergy(), and SampleSecondaries().


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