| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 // 25 //
26 // 26 //
27 // ------------------------------------------- 27 // -------------------------------------------------------------------
28 // 28 //
29 // GEANT4 Class file 29 // GEANT4 Class file
30 // 30 //
31 // 31 //
32 // File name: G4eeToHadronsMultiModel 32 // File name: G4eeToHadronsMultiModel
33 // 33 //
34 // Author: Vladimir Ivanchenko 34 // Author: Vladimir Ivanchenko
35 // 35 //
36 // Creation date: 02.08.2004 36 // Creation date: 02.08.2004
37 // 37 //
38 // Modifications: 38 // Modifications:
39 // 08-11-04 Migration to new interface of Stor 39 // 08-11-04 Migration to new interface of Store/Retrieve tables (V.Ivantchenko)
40 // 08-04-05 Major optimisation of internal int 40 // 08-04-05 Major optimisation of internal interfaces (V.Ivantchenko)
41 // 41 //
42 42
43 // 43 //
44 // ------------------------------------------- 44 // -------------------------------------------------------------------
45 // 45 //
46 //....oooOO0OOooo........oooOO0OOooo........oo 46 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47 //....oooOO0OOooo........oooOO0OOooo........oo 47 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
48 48
49 #include "G4eeToHadronsMultiModel.hh" 49 #include "G4eeToHadronsMultiModel.hh"
50 #include "G4PhysicalConstants.hh" 50 #include "G4PhysicalConstants.hh"
51 #include "G4SystemOfUnits.hh" 51 #include "G4SystemOfUnits.hh"
52 #include "G4eeToTwoPiModel.hh" 52 #include "G4eeToTwoPiModel.hh"
53 #include "G4eeTo3PiModel.hh" 53 #include "G4eeTo3PiModel.hh"
54 #include "G4eeToPGammaModel.hh" 54 #include "G4eeToPGammaModel.hh"
55 #include "G4ee2KNeutralModel.hh" 55 #include "G4ee2KNeutralModel.hh"
56 #include "G4ee2KChargedModel.hh" 56 #include "G4ee2KChargedModel.hh"
57 #include "G4eeCrossSections.hh" 57 #include "G4eeCrossSections.hh"
58 #include "G4Vee2hadrons.hh" 58 #include "G4Vee2hadrons.hh"
59 #include "G4Positron.hh" 59 #include "G4Positron.hh"
60 60
61 //....oooOO0OOooo........oooOO0OOooo........oo 61 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
62 62
63 using namespace std; 63 using namespace std;
64 64
65 G4eeToHadronsMultiModel::G4eeToHadronsMultiMod 65 G4eeToHadronsMultiModel::G4eeToHadronsMultiModel(G4int ver,
66 const G4String& mname) : G4VEmModel(mname), 66 const G4String& mname) : G4VEmModel(mname), verbose(ver)
67 { 67 {
68 maxKineticEnergy = 4.521*CLHEP::GeV; //crre 68 maxKineticEnergy = 4.521*CLHEP::GeV; //crresponding to 10TeV in lab
69 delta = 1.0*CLHEP::MeV; //for bin width 69 delta = 1.0*CLHEP::MeV; //for bin width
70 } 70 }
71 71
72 //....oooOO0OOooo........oooOO0OOooo........oo 72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
73 73
74 G4eeToHadronsMultiModel::~G4eeToHadronsMultiMo 74 G4eeToHadronsMultiModel::~G4eeToHadronsMultiModel()
75 { 75 {
76 delete cross; 76 delete cross;
77 } 77 }
78 78
79 //....oooOO0OOooo........oooOO0OOooo........oo 79 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
80 80
81 void G4eeToHadronsMultiModel::Initialise(const 81 void G4eeToHadronsMultiModel::Initialise(const G4ParticleDefinition*,
82 const G4DataVector& cuts) 82 const G4DataVector& cuts)
83 { 83 {
84 if(!isInitialised) { 84 if(!isInitialised) {
85 isInitialised = true; 85 isInitialised = true;
86 86
87 //G4cout<<"###Initialise in HadronMultiMod 87 //G4cout<<"###Initialise in HadronMultiModel###"<<G4endl;
88 88
89 cross = new G4eeCrossSections(); 89 cross = new G4eeCrossSections();
90 90
91 G4eeToTwoPiModel* m2pi = 91 G4eeToTwoPiModel* m2pi =
92 new G4eeToTwoPiModel(cross,maxKineticEne 92 new G4eeToTwoPiModel(cross,maxKineticEnergy,delta);
93 AddEEModel(m2pi,cuts); 93 AddEEModel(m2pi,cuts);
94 94
95 G4eeTo3PiModel* m3pi = 95 G4eeTo3PiModel* m3pi =
96 new G4eeTo3PiModel(cross,maxKineticEnerg 96 new G4eeTo3PiModel(cross,maxKineticEnergy,delta);
97 AddEEModel(m3pi,cuts); 97 AddEEModel(m3pi,cuts);
98 98
99 G4ee2KChargedModel* m2kc = 99 G4ee2KChargedModel* m2kc =
100 new G4ee2KChargedModel(cross,maxKineticE 100 new G4ee2KChargedModel(cross,maxKineticEnergy,delta);
101 AddEEModel(m2kc,cuts); 101 AddEEModel(m2kc,cuts);
102 102
103 G4ee2KNeutralModel* m2kn = 103 G4ee2KNeutralModel* m2kn =
104 new G4ee2KNeutralModel(cross,maxKineticE 104 new G4ee2KNeutralModel(cross,maxKineticEnergy,delta);
105 AddEEModel(m2kn,cuts); 105 AddEEModel(m2kn,cuts);
106 106
107 G4eeToPGammaModel* mpg1 = 107 G4eeToPGammaModel* mpg1 =
108 new G4eeToPGammaModel(cross,"pi0",maxKin 108 new G4eeToPGammaModel(cross,"pi0",maxKineticEnergy,delta);
109 AddEEModel(mpg1,cuts); 109 AddEEModel(mpg1,cuts);
110 110
111 G4eeToPGammaModel* mpg2 = 111 G4eeToPGammaModel* mpg2 =
112 new G4eeToPGammaModel(cross,"eta",maxKin 112 new G4eeToPGammaModel(cross,"eta",maxKineticEnergy,delta);
113 AddEEModel(mpg2,cuts); 113 AddEEModel(mpg2,cuts);
114 114
115 nModels = (G4int)models.size(); 115 nModels = (G4int)models.size();
116 116
117 fParticleChange = GetParticleChangeForGamm 117 fParticleChange = GetParticleChangeForGamma();
118 } 118 }
119 } 119 }
120 120
121 //....oooOO0OOooo........oooOO0OOooo........oo 121 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
122 122
123 void G4eeToHadronsMultiModel::AddEEModel(G4Vee 123 void G4eeToHadronsMultiModel::AddEEModel(G4Vee2hadrons* mod,
124 const G4DataVector& cuts) 124 const G4DataVector& cuts)
125 { 125 {
126 G4eeToHadronsModel* model = new G4eeToHadron 126 G4eeToHadronsModel* model = new G4eeToHadronsModel(mod, verbose);
127 models.push_back(model); 127 models.push_back(model);
128 G4double elow = mod->LowEnergy(); 128 G4double elow = mod->LowEnergy();
129 ekinMin.push_back(elow); 129 ekinMin.push_back(elow);
130 if(thKineticEnergy > elow) { thKineticEnergy 130 if(thKineticEnergy > elow) { thKineticEnergy = elow; }
131 ekinMax.push_back(mod->HighEnergy()); 131 ekinMax.push_back(mod->HighEnergy());
132 ekinPeak.push_back(mod->PeakEnergy()); 132 ekinPeak.push_back(mod->PeakEnergy());
133 cumSum.push_back(0.0); 133 cumSum.push_back(0.0);
134 134
135 const G4ParticleDefinition* positron = G4Pos 135 const G4ParticleDefinition* positron = G4Positron::Positron();
136 model->Initialise(positron,cuts); 136 model->Initialise(positron,cuts);
137 } 137 }
138 138
139 //....oooOO0OOooo........oooOO0OOooo........oo 139 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
140 140
141 G4double G4eeToHadronsMultiModel::CrossSection 141 G4double G4eeToHadronsMultiModel::CrossSectionPerVolume(
142 const G4Material* mat, 142 const G4Material* mat,
143 const G4ParticleDefinition* p, 143 const G4ParticleDefinition* p,
144 G4double kineticEnergy, 144 G4double kineticEnergy,
145 G4double, G4double) 145 G4double, G4double)
146 { 146 {
147 return mat->GetElectronDensity()* 147 return mat->GetElectronDensity()*
148 ComputeCrossSectionPerElectron(p, kineticE 148 ComputeCrossSectionPerElectron(p, kineticEnergy);
149 } 149 }
150 150
151 //....oooOO0OOooo........oooOO0OOooo........oo 151 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
152 152
153 G4double G4eeToHadronsMultiModel::ComputeCross 153 G4double G4eeToHadronsMultiModel::ComputeCrossSectionPerAtom(
154 const G4 154 const G4ParticleDefinition* p,
155 G4double kineticEnergy, 155 G4double kineticEnergy,
156 G4double Z, G4double, 156 G4double Z, G4double,
157 G4double, G4double) 157 G4double, G4double)
158 { 158 {
159 return Z*ComputeCrossSectionPerElectron(p, k 159 return Z*ComputeCrossSectionPerElectron(p, kineticEnergy);
160 } 160 }
161 161
162 //....oooOO0OOooo........oooOO0OOooo........oo 162 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
163 163
164 G4double G4eeToHadronsMultiModel::ComputeCross 164 G4double G4eeToHadronsMultiModel::ComputeCrossSectionPerElectron(const G4ParticleDefinition*,
165 G4double kineticEnergy, 165 G4double kineticEnergy,
166 G4double, G4double) 166 G4double, G4double)
167 { 167 {
168 G4double res = 0.0; 168 G4double res = 0.0;
169 169
170 G4double energy = LabToCM(kineticEnergy); 170 G4double energy = LabToCM(kineticEnergy);
171 171
172 if (energy > thKineticEnergy) { 172 if (energy > thKineticEnergy) {
173 for(G4int i=0; i<nModels; i++) { 173 for(G4int i=0; i<nModels; i++) {
174 if(energy >= ekinMin[i] && energy <= eki 174 if(energy >= ekinMin[i] && energy <= ekinMax[i]){
175 res += (models[i])->ComputeCrossSectio 175 res += (models[i])->ComputeCrossSectionPerElectron(0,energy);
176 } 176 }
177 cumSum[i] = res; 177 cumSum[i] = res;
178 } 178 }
179 } 179 }
180 return res*csFactor; 180 return res*csFactor;
181 } 181 }
182 182
183 //....oooOO0OOooo........oooOO0OOooo........oo 183 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
184 184
185 void G4eeToHadronsMultiModel::SampleSecondarie 185 void G4eeToHadronsMultiModel::SampleSecondaries(
186 std::vector<G4DynamicParticle*>* ne 186 std::vector<G4DynamicParticle*>* newp,
187 const G4MaterialCutsCouple* couple, 187 const G4MaterialCutsCouple* couple,
188 const G4DynamicParticle* dp, 188 const G4DynamicParticle* dp,
189 G4double, G4double) 189 G4double, G4double)
190 { 190 {
191 G4double kinEnergy = dp->GetKineticEnergy(); 191 G4double kinEnergy = dp->GetKineticEnergy();
192 G4double energy = LabToCM(kinEnergy); 192 G4double energy = LabToCM(kinEnergy);
193 if (energy > thKineticEnergy) { 193 if (energy > thKineticEnergy) {
194 G4double q = cumSum[nModels-1]*G4UniformRa 194 G4double q = cumSum[nModels-1]*G4UniformRand();
195 for(G4int i=0; i<nModels; ++i) { 195 for(G4int i=0; i<nModels; ++i) {
196 if(q <= cumSum[i]) { 196 if(q <= cumSum[i]) {
197 (models[i])->SampleSecondaries(newp, c 197 (models[i])->SampleSecondaries(newp, couple,dp);
198 if(newp->size() > 0) { 198 if(newp->size() > 0) {
199 fParticleChange->ProposeTrackStatus(fStopA 199 fParticleChange->ProposeTrackStatus(fStopAndKill);
200 } 200 }
201 break; 201 break;
202 } 202 }
203 } 203 }
204 } 204 }
205 } 205 }
206 206
207 //....oooOO0OOooo........oooOO0OOooo........oo 207 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
208 208
209 void G4eeToHadronsMultiModel::ModelDescription 209 void G4eeToHadronsMultiModel::ModelDescription(std::ostream& outFile) const
210 { 210 {
211 if(verbose > 0) { 211 if(verbose > 0) {
212 G4double e1 = 0.5*thKineticEnergy*thKineti 212 G4double e1 = 0.5*thKineticEnergy*thKineticEnergy/electron_mass_c2
213 - 2.0*electron_mass_c2; 213 - 2.0*electron_mass_c2;
214 G4double e2 = 0.5*maxKineticEnergy*maxKine 214 G4double e2 = 0.5*maxKineticEnergy*maxKineticEnergy/electron_mass_c2
215 - 2.0*electron_mass_c2; 215 - 2.0*electron_mass_c2;
216 outFile << " e+ annihilation into had 216 outFile << " e+ annihilation into hadrons active from "
217 << e1/GeV << " GeV to " << e2/GeV << " G 217 << e1/GeV << " GeV to " << e2/GeV << " GeV" << G4endl;
218 } 218 }
219 } 219 }
220 220
221 //....oooOO0OOooo........oooOO0OOooo........oo 221 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
222 222
223 void G4eeToHadronsMultiModel::SetCrossSecFacto 223 void G4eeToHadronsMultiModel::SetCrossSecFactor(G4double fac)
224 { 224 {
225 if(fac > 1.0) { 225 if(fac > 1.0) {
226 csFactor = fac; 226 csFactor = fac;
227 if(verbose > 0) { 227 if(verbose > 0) {
228 G4cout << "### G4eeToHadronsMultiModel: 228 G4cout << "### G4eeToHadronsMultiModel: The cross section for "
229 << "G4eeToHadronsMultiModel is increase 229 << "G4eeToHadronsMultiModel is increased by "
230 << csFactor << " times" << G4endl; 230 << csFactor << " times" << G4endl;
231 } 231 }
232 } 232 }
233 } 233 }
234 234
235 //....oooOO0OOooo........oooOO0OOooo........oo 235 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
236 236