Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4eSingleCoulombScatteringModel.cc
이 파일의 문서화 페이지로 가기
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // G4eSingleCoulombScatteringModel.cc
27 // -------------------------------------------------------------------
28 //
29 // GEANT4 Class header file
30 //
31 // File name: G4eSingleCoulombScatteringModel
32 //
33 // Author: Cristina Consolandi
34 //
35 // Creation date: 20.10.2012
36 //
37 // Class Description:
38 // Single Scattering model for electron-nuclei interaction.
39 // Suitable for high energy electrons and low scattering angles.
40 //
41 //
42 // Reference:
43 // M.J. Boschini et al. "Non Ionizing Energy Loss induced by Electrons
44 // in the Space Environment" Proc. of the 13th International Conference
45 // on Particle Physics and Advanced Technology
46 //
47 // (13th ICPPAT, Como 3-7/10/2011), World Scientific (Singapore).
48 // Available at: http://arxiv.org/abs/1111.4042v4
49 //
50 //
51 // -------------------------------------------------------------------
52 //
53 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
54 
55 
57 #include "G4PhysicalConstants.hh"
58 #include "G4SystemOfUnits.hh"
59 #include "Randomize.hh"
61 #include "G4Proton.hh"
62 #include "G4ProductionCutsTable.hh"
63 #include "G4NucleiProperties.hh"
64 #include "G4NistManager.hh"
65 #include "G4ParticleTable.hh"
66 #include "G4IonTable.hh"
67 
68 #include "G4UnitsTable.hh"
69 #include "G4EmParameters.hh"
70 
71 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
72 
73 using namespace std;
74 
76  : G4VEmModel(nam),
77  cosThetaMin(1.0)
78 {
81  fParticleChange = nullptr;
82 
83  pCuts=nullptr;
84  currentMaterial = nullptr;
85  currentElement = nullptr;
86  currentCouple = nullptr;
87 
88  lowEnergyLimit = 0*keV;
89  recoilThreshold = 0.*eV;
90  XSectionModel = 1;
91  FormFactor = 0;
92  particle = nullptr;
93  mass=0.0;
95 
97 }
98 
99 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
100 
102 {
103  delete Mottcross;
104 }
105 
106 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
107 
109  const G4DataVector& cuts)
110 {
112 
113  SetupParticle(p);
114  currentCouple = nullptr;
116  //cosThetaMin = cos(PolarAngleLimit());
118 
119  pCuts = &cuts;
120  //G4ProductionCutsTable::GetProductionCutsTable()->GetEnergyCutsVector(3);
121 
122  /*
123  G4cout << "!!! G4eSingleCoulombScatteringModel::Initialise for "
124  << part->GetParticleName() << " cos(TetMin)= " << cosThetaMin
125  << " cos(TetMax)= " << cosThetaMax <<G4endl;
126  G4cout << "cut= " << (*pCuts)[0] << " cut1= " << (*pCuts)[1] << G4endl;
127  */
128 
129  if(!fParticleChange) {
131  }
132 
133  if(IsMaster()) {
135  }
136 
138 
139  //G4cout<<"NUCLEAR FORM FACTOR: "<<FormFactor<<G4endl;
140 }
141 
143  G4VEmModel* masterModel)
144 {
146 }
147 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
148 
150  const G4ParticleDefinition* p,
151  G4double kinEnergy,
152  G4double Z,
153  G4double ,
154  G4double,
155  G4double )
156 {
157  SetupParticle(p);
158 
159  G4double cross =0.0;
160  if(kinEnergy < lowEnergyLimit) return cross;
161 
163 
164  //Total Cross section
165  Mottcross->SetupKinematic(kinEnergy, Z);
167 
168  //cout<< "Compute Cross Section....cross "<<G4BestUnit(cross,"Surface") << " cm2 "<< cross/cm2 <<" Z: "<<Z<<" kinEnergy: "<<kinEnergy<<endl;
169 
170  //G4cout<<"Energy: "<<kinEnergy/MeV<<" Total Cross: "<<cross<<G4endl;
171  return cross;
172 }
173 
174 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
175 
177  std::vector<G4DynamicParticle*>* fvect,
178  const G4MaterialCutsCouple* couple,
179  const G4DynamicParticle* dp,
180  G4double cutEnergy,
181  G4double)
182 {
183  G4double kinEnergy = dp->GetKineticEnergy();
184  //cout<<"--- kinEnergy "<<kinEnergy<<endl;
185 
186  if(kinEnergy < lowEnergyLimit) return;
187 
188  DefineMaterial(couple);
190 
191  // Choose nucleus
192  //last two :cutEnergy= min e kinEnergy=max
194  kinEnergy,cutEnergy,kinEnergy);
195 
196  G4double Z = currentElement->GetZ();
197  G4int iz = G4int(Z);
198  G4int ia = SelectIsotopeNumber(currentElement);
200 
201  //G4cout<<"..Z: "<<Z<<" ..iz: "<<iz<<" ..ia: "<<ia<<" ..mass2: "<<mass2<<G4endl;
202 
203  Mottcross->SetupKinematic(kinEnergy, Z);
205  if(cross == 0.0) { return; }
206  //cout<< "Energy: "<<kinEnergy/MeV<<" Z: "<<Z<<"....cross "<<G4BestUnit(cross,"Surface") << " cm2 "<< cross/cm2 <<endl;
207 
209  G4double sint = sin(z1);
210  G4double cost = sqrt(1.0 - sint*sint);
211  G4double phi = twopi* G4UniformRand();
212 
213  // kinematics in the Lab system
214  G4double ptot = dp->GetTotalMomentum();
215  G4double e1 = dp->GetTotalEnergy();
216  // Lab. system kinematics along projectile direction
217  G4LorentzVector v0 = G4LorentzVector(0, 0, ptot, e1);
218  G4double bet = ptot/(v0.e() + mass2);
219  G4double gam = 1.0/sqrt((1.0 - bet)*(1.0 + bet));
220 
221  //CM Projectile
222  G4double momCM = gam*(ptot - bet*e1);
223  G4double eCM = gam*(e1 - bet*ptot);
224  //energy & momentum after scattering of incident particle
225  G4double pxCM = momCM*sint*cos(phi);
226  G4double pyCM = momCM*sint*sin(phi);
227  G4double pzCM = momCM*cost;
228 
229  //CM--->Lab
230  G4LorentzVector v1(pxCM , pyCM, gam*(pzCM + bet*eCM), gam*(eCM + bet*pzCM));
231 
232  // Rotate to global system
234  G4ThreeVector newDirection = v1.vect().unit();
235  newDirection.rotateUz(dir);
236 
238 
239  // recoil
240  v0 -= v1;
241  G4double trec = v0.e();
242  G4double edep = 0.0;
243 
244  G4double tcut = recoilThreshold;
245 
246  //G4cout<<" Energy Transfered: "<<trec/eV<<G4endl;
247 
248  if(pCuts) {
249  tcut= std::max(tcut,(*pCuts)[currentMaterialIndex]);
250  //G4cout<<"Cuts: "<<(*pCuts)[currentMaterialIndex]/eV<<" eV"<<G4endl;
251  //G4cout<<"Threshold: "<<tcut/eV<<" eV"<<G4endl;
252  }
253 
254  if(trec > tcut) {
256  newDirection = v0.vect().unit();
257  newDirection.rotateUz(dir);
258  G4DynamicParticle* newdp = new G4DynamicParticle(ion, newDirection, trec);
259  fvect->push_back(newdp);
260  } else if(trec > 0.0) {
261  edep = trec;
263  }
264 
265  // finelize primary energy and energy balance
266  G4double finalT = v1.e() - mass;
267  //G4cout<<"Final Energy: "<<finalT/eV<<G4endl;
268  if(finalT <= lowEnergyLimit) {
269  edep += finalT;
270  finalT = 0.0;
271  }
272  edep = std::max(edep, 0.0);
275 
276 }
277 
278 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
void SetElementSelectors(std::vector< G4EmElementSelector * > *)
Definition: G4VEmModel.hh:792
T max(const T t1, const T t2)
brief Return the largest of the two arguments
virtual void SampleSecondaries(std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) final
static G4ParticleTable * GetParticleTable()
static constexpr double keV
Definition: G4SIunits.hh:216
G4int SelectIsotopeNumber(const G4Element *)
Definition: G4VEmModel.cc:335
const G4ThreeVector & GetMomentumDirection() const
const char * p
Definition: xmltok.h:285
G4IonTable * GetIonTable() const
void SetProposedKineticEnergy(G4double proposedKinEnergy)
void SetupParticle(const G4ParticleDefinition *)
Hep3Vector & rotateUz(const Hep3Vector &)
Definition: ThreeVector.cc:38
const G4Element * SelectRandomAtom(const G4MaterialCutsCouple *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy=0.0, G4double maxEnergy=DBL_MAX)
Definition: G4VEmModel.hh:552
std::vector< G4EmElementSelector * > * GetElementSelectors()
Definition: G4VEmModel.hh:784
virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A, G4double cut, G4double emax) final
void ProposeNonIonizingEnergyDeposit(G4double anEnergyPart)
Float_t Z
double G4double
Definition: G4Types.hh:76
G4ParticleDefinition * GetDefinition() const
void Initialise(const G4ParticleDefinition *, G4double cosThetaLim)
void InitialiseElementSelectors(const G4ParticleDefinition *, const G4DataVector &)
Definition: G4VEmModel.cc:146
void DefineMaterial(const G4MaterialCutsCouple *)
void ProposeMomentumDirection(G4double Px, G4double Py, G4double Pz)
#define G4UniformRand()
Definition: Randomize.hh:53
static constexpr double twopi
Definition: G4SIunits.hh:76
static G4double GetNuclearMass(const G4double A, const G4double Z)
static constexpr double eV
Definition: G4SIunits.hh:215
Double_t edep
G4eSingleCoulombScatteringModel(const G4String &nam="eSingleCoulombScat")
const G4MaterialCutsCouple * CurrentCouple() const
Definition: G4VEmModel.hh:462
Hep3Vector unit() const
int G4int
Definition: G4Types.hh:78
G4NuclearFormfactorType NuclearFormfactorType() const
TDirectory * dir
G4double GetKineticEnergy() const
G4bool IsMaster() const
Definition: G4VEmModel.hh:700
Hep3Vector vect() const
G4double GetTotalMomentum() const
CLHEP::HepLorentzVector G4LorentzVector
G4double GetTotalEnergy() const
G4ParticleDefinition * GetIon(G4int Z, G4int A, G4int lvl=0)
Definition: G4IonTable.cc:504
G4ParticleChangeForGamma * GetParticleChangeForGamma()
Definition: G4VEmModel.cc:131
void SetupKinematic(G4double kinEnergy, G4double Z)
void ProposeLocalEnergyDeposit(G4double anEnergyPart)
static G4EmParameters * Instance()
virtual void Initialise(const G4ParticleDefinition *, const G4DataVector &) final
virtual void InitialiseLocal(const G4ParticleDefinition *, G4VEmModel *masterModel) final
static G4NistManager * Instance()