| 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 /// \file medical/fanoCavity/src/Run.cc 26 /// \file medical/fanoCavity/src/Run.cc
27 /// \brief Implementation of the Run class 27 /// \brief Implementation of the Run class
28 // 28 //
29 // << 29 // $Id: Run.cc 71035 2013-06-10 09:17:35Z gcosmo $
>> 30 //
30 //....oooOO0OOooo........oooOO0OOooo........oo 31 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
31 //....oooOO0OOooo........oooOO0OOooo........oo 32 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
32 33
33 #include "Run.hh" 34 #include "Run.hh"
34 <<
35 #include "DetectorConstruction.hh" 35 #include "DetectorConstruction.hh"
36 #include "HistoManager.hh" <<
37 #include "PrimaryGeneratorAction.hh" 36 #include "PrimaryGeneratorAction.hh"
>> 37 #include "HistoManager.hh"
38 38
39 #include "G4Electron.hh" <<
40 #include "G4EmCalculator.hh" <<
41 #include "G4Run.hh" 39 #include "G4Run.hh"
42 #include "G4RunManager.hh" 40 #include "G4RunManager.hh"
43 #include "G4SystemOfUnits.hh" <<
44 #include "G4UnitsTable.hh" 41 #include "G4UnitsTable.hh"
45 #include "Randomize.hh" << 42 #include "G4EmCalculator.hh"
>> 43 #include "G4Electron.hh"
46 44
>> 45 #include "G4SystemOfUnits.hh"
>> 46 #include "Randomize.hh"
47 #include <iomanip> 47 #include <iomanip>
48 48
49 //....oooOO0OOooo........oooOO0OOooo........oo 49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50 50
51 Run::Run(DetectorConstruction* det, PrimaryGen << 51 Run::Run(DetectorConstruction* det,PrimaryGeneratorAction *kin )
52 : fDetector(det), fKinematic(kin), fProcCoun << 52 :fDetector(det),fKinematic(kin), fProcCounter(0), fMateWall(0),fMateCavity(0)
53 53
54 { 54 {
55 // geometry <<
56 // <<
57 fWallThickness = fDetector->GetWallThickness <<
58 fWallRadius = fDetector->GetWallRadius(); <<
59 fMateWall = fDetector->GetWallMaterial(); <<
60 fDensityWall = fMateWall->GetDensity(); <<
61 <<
62 fCavityThickness = fDetector->GetCavityThick <<
63 fCavityRadius = fDetector->GetCavityRadius() <<
64 fSurfaceCavity = CLHEP::pi * fCavityRadius * <<
65 fVolumeCavity = fSurfaceCavity * fCavityThic <<
66 fMateCavity = fDetector->GetCavityMaterial() <<
67 fDensityCavity = fMateCavity->GetDensity(); <<
68 fMassCavity = fVolumeCavity * fDensityCavity <<
69 <<
70 // process counter <<
71 // <<
72 fProcCounter = new ProcessesCount; <<
73 <<
74 // kinetic energy of charged secondary a cre <<
75 // <<
76 fEsecondary = fEsecondary2 = 0.; <<
77 fNbSec = 0; <<
78 <<
79 // charged particles and energy flow in cavi <<
80 // <<
81 fPartFlowCavity[0] = fPartFlowCavity[1] = 0; <<
82 fEnerFlowCavity[0] = fEnerFlowCavity[1] = 0. <<
83 <<
84 // total energy deposit and charged track se <<
85 // <<
86 fEdepCavity = fEdepCavity2 = fTrkSegmCavity <<
87 fNbEventCavity = 0; <<
88 55
89 // survey convergence << 56 //geometry
90 // << 57 //
91 fOldEmean = fOldDose = 0.; << 58 fWallThickness = fDetector->GetWallThickness();
92 << 59 fWallRadius = fDetector->GetWallRadius();
93 // stepLenth of charged particles << 60 fMateWall = fDetector->GetWallMaterial();
94 // << 61 fDensityWall = fMateWall->GetDensity();
95 fStepWall = fStepWall2 = fStepCavity = fStep << 62
96 fNbStepWall = fNbStepCavity = 0; << 63 fCavityThickness = fDetector->GetCavityThickness();
97 << 64 fCavityRadius = fDetector->GetCavityRadius();
98 // histograms << 65 fSurfaceCavity = CLHEP::pi*fCavityRadius*fCavityRadius;
99 // << 66 fVolumeCavity = fSurfaceCavity*fCavityThickness;
100 G4AnalysisManager* analysisManager = G4Analy << 67 fMateCavity = fDetector->GetCavityMaterial();
101 if (analysisManager->IsActive()) { << 68 fDensityCavity = fMateCavity->GetDensity();
102 analysisManager->OpenFile(); << 69 fMassCavity = fVolumeCavity*fDensityCavity;
103 } << 70
>> 71 //process counter
>> 72 //
>> 73 fProcCounter = new ProcessesCount;
>> 74
>> 75 //kinetic energy of charged secondary a creation
>> 76 //
>> 77 fEsecondary = fEsecondary2 = 0.;
>> 78 fNbSec = 0;
>> 79
>> 80 //charged particles and energy flow in cavity
>> 81 //
>> 82 fPartFlowCavity[0] = fPartFlowCavity[1] = 0;
>> 83 fEnerFlowCavity[0] = fEnerFlowCavity[1] = 0.;
>> 84
>> 85 //total energy deposit and charged track segment in cavity
>> 86 //
>> 87 fEdepCavity = fEdepCavity2 = fTrkSegmCavity = 0.;
>> 88 fNbEventCavity = 0;
>> 89
>> 90 //survey convergence
>> 91 //
>> 92 fOldEmean = fOldDose = 0.;
>> 93
>> 94 //stepLenth of charged particles
>> 95 //
>> 96 fStepWall = fStepWall2 = fStepCavity = fStepCavity2 =0.;
>> 97 fNbStepWall = fNbStepCavity = 0;
>> 98
>> 99 //histograms
>> 100 //
>> 101 G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
>> 102 if ( analysisManager->IsActive() ) {
>> 103 analysisManager->OpenFile();
>> 104 }
104 } 105 }
105 106
106 //....oooOO0OOooo........oooOO0OOooo........oo 107 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
107 108
108 Run::~Run() {} << 109 Run::~Run()
>> 110 {
>> 111 }
109 112
110 //....oooOO0OOooo........oooOO0OOooo........oo 113 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
111 114
112 void Run::EndOfRun() 115 void Run::EndOfRun()
113 { // Only call by Master thread << 116 { // Only call by Master thread
114 std::ios::fmtflags mode = G4cout.flags(); << 117 std::ios::fmtflags mode = G4cout.flags();
115 G4cout.setf(std::ios::fixed, std::ios::float << 118 G4cout.setf(std::ios::fixed,std::ios::floatfield);
116 << 119
117 if (numberOfEvent == 0) return; << 120 if (numberOfEvent == 0) return;
118 << 121
119 // run conditions << 122 //run conditions
120 // << 123 //
121 G4ParticleDefinition* particle = fKinematic- << 124 G4ParticleDefinition* particle = fKinematic->GetParticleGun()
122 G4String partName = particle->GetParticleNam << 125 ->GetParticleDefinition();
123 G4double energy = fKinematic->GetParticleGun << 126 G4String partName = particle->GetParticleName();
124 << 127 G4double energy = fKinematic->GetParticleGun()->GetParticleEnergy();
125 G4cout << "\n ======================== run s << 128
126 << 129 G4cout <<"\n ======================== run summary ======================\n";
127 G4int prec = G4cout.precision(3); << 130
128 << 131 G4int prec = G4cout.precision(3);
129 G4cout << "\n The run consists of " << numbe << 132
130 << G4BestUnit(energy, "Energy") << " << 133 G4cout <<"\n The run consists of "<<numberOfEvent<<" "<< partName << " of "
131 << " of " << fMateWall->GetName() << 134 << G4BestUnit(energy,"Energy") << " through 2*"
132 << " (density: " << G4BestUnit(fDensi << 135 << G4BestUnit(fWallThickness,"Length") << " of "
133 << 136 << fMateWall->GetName() << " (density: "
134 G4cout << "\n the cavity is " << G4BestUnit( << 137 << G4BestUnit(fDensityWall,"Volumic Mass") << ")" << G4endl;
135 << fMateCavity->GetName() << " (densi << 138
136 << "); Mass = " << G4BestUnit(fMassCa << 139 G4cout << "\n the cavity is "
137 << 140 << G4BestUnit(fCavityThickness,"Length") << " of "
138 G4cout << "\n ============================== << 141 << fMateCavity->GetName() << " (density: "
139 << 142 << G4BestUnit(fDensityCavity,"Volumic Mass") << "); Mass = "
140 // frequency of processes << 143 << G4BestUnit(fMassCavity,"Mass") << G4endl;
141 // << 144
142 G4cout << "\n Process calls frequency --->"; << 145 G4cout<<"\n ============================================================\n";
143 for (size_t i = 0; i < fProcCounter->size(); << 146
144 G4String procName = (*fProcCounter)[i]->Ge << 147 //frequency of processes
145 G4int count = (*fProcCounter)[i]->GetCount << 148 //
146 G4cout << " " << procName << "= " << coun << 149 G4cout << "\n Process calls frequency --->";
147 } << 150 for (size_t i=0; i< fProcCounter->size();i++) {
148 G4cout << G4endl; << 151 G4String procName = (*fProcCounter)[i]->GetName();
149 << 152 G4int count = (*fProcCounter)[i]->GetCounter();
150 // extract cross sections with G4EmCalculato << 153 G4cout << " " << procName << "= " << count;
151 // <<
152 G4EmCalculator emCalculator; <<
153 G4cout << "\n Gamma crossSections in wall ma <<
154 G4double sumc = 0.0; <<
155 for (size_t i = 0; i < fProcCounter->size(); <<
156 G4String procName = (*fProcCounter)[i]->Ge <<
157 G4double massSigma = <<
158 emCalculator.ComputeCrossSectionPerVolum <<
159 / fDensityWall; <<
160 if (massSigma > 0.) { <<
161 sumc += massSigma; <<
162 G4cout << " " << procName << "= " << G4 <<
163 } 154 }
164 } << 155 G4cout << G4endl;
165 G4cout << " --> total= " << G4BestUnit(sum <<
166 <<
167 // mean kinetic energy of secondary electron <<
168 // <<
169 if (fNbSec == 0) return; <<
170 G4double meanEsecond = fEsecondary / fNbSec, <<
171 G4double varianceEsec = meanEsecond2 - meanE <<
172 G4double dToverT = 0.; <<
173 if (varianceEsec > 0.) dToverT = std::sqrt(v <<
174 G4double csdaRange = emCalculator.GetCSDARan <<
175 <<
176 G4cout.precision(4); <<
177 G4cout << "\n Mean energy of secondary e- = <<
178 << 100 * dToverT << " %" <<
179 << " (--> range in wall material = " <<
180 <<
181 // compute mass energy transfer coefficient <<
182 // <<
183 G4double massTransfCoef = sumc * meanEsecond <<
184 <<
185 G4cout << " Mass_energy_transfer coef: " << <<
186 156
187 // stopping power from EmCalculator << 157 //extract cross sections with G4EmCalculator
188 // << 158 //
189 G4double dedxWall = emCalculator.GetDEDX(mea << 159 G4EmCalculator emCalculator;
190 dedxWall /= fDensityWall; << 160 G4cout << "\n Gamma crossSections in wall material :";
191 G4double dedxCavity = emCalculator.GetDEDX(m << 161 G4double sumc = 0.0;
192 dedxCavity /= fDensityCavity; << 162 for (size_t i=0; i< fProcCounter->size();i++) {
193 << 163 G4String procName = (*fProcCounter)[i]->GetName();
194 G4cout << "\n StoppingPower in wall = " << << 164 G4double massSigma =
195 << "\n in cavity = " << << 165 emCalculator.ComputeCrossSectionPerVolume(energy,particle,
196 << G4endl; << 166 procName,fMateWall)/fDensityWall;
197 << 167 if (massSigma > 0.) {
198 // charged particle flow in cavity << 168 sumc += massSigma;
199 // << 169 G4cout << " " << procName << "= "
200 G4cout << "\n Charged particle flow in cavit << 170 << G4BestUnit(massSigma, "Surface/Mass");
201 << "\n Enter --> nbParticles = " << 171 }
202 << "\t Energy = " << G4BestUnit(fEner << 172 }
203 << "\n Exit --> nbParticles = " << 173 G4cout << " --> total= "
204 << "\t Energy = " << G4BestUnit(fEner << 174 << G4BestUnit(sumc, "Surface/Mass") << G4endl;
205 <<
206 if (fPartFlowCavity[0] == 0) return; <<
207 <<
208 // beam energy fluence <<
209 // <<
210 G4double rBeam = fWallRadius * (fKinematic-> <<
211 G4double surfaceBeam = CLHEP::pi * rBeam * r <<
212 <<
213 // error on Edep in cavity <<
214 // <<
215 if (fNbEventCavity == 0) return; <<
216 G4double meanEdep = fEdepCavity / fNbEventCa <<
217 G4double meanEdep2 = fEdepCavity2 / fNbEvent <<
218 G4double varianceEdep = meanEdep2 - meanEdep <<
219 G4double dEoverE = 0.; <<
220 if (varianceEdep > 0.) dEoverE = std::sqrt(v <<
221 <<
222 // total dose in cavity <<
223 // <<
224 G4double doseCavity = fEdepCavity / fMassCav <<
225 G4double doseOverBeam = doseCavity * surface <<
226 <<
227 // track length in cavity <<
228 G4double meantrack = fTrkSegmCavity / fPartF <<
229 <<
230 G4cout.precision(4); <<
231 G4cout << "\n Total edep in cavity = " << G4 <<
232 << 100 * dEoverE << " %" <<
233 << "\t Total charged trackLength = " <<
234 << " (mean value = " << G4BestUnit( <<
235 << "\n Total dose in cavity = " << do <<
236 << "\n Dose/EnergyFluence = " << G4 <<
237 <<
238 // ratio simulation/theory <<
239 // <<
240 G4double ratio = doseOverBeam / massTransfCo <<
241 G4double error = ratio * std::sqrt(dEoverE * <<
242 <<
243 G4cout.precision(5); <<
244 G4cout << "\n (Dose/EnergyFluence)/Mass_ener <<
245 <<
246 // compute mean step size of charged particl <<
247 // <<
248 fStepWall /= fNbStepWall; <<
249 fStepWall2 /= fNbStepWall; <<
250 G4double rms = fStepWall2 - fStepWall * fSte <<
251 if (rms > 0.) <<
252 rms = std::sqrt(rms); <<
253 else <<
254 rms = 0.; <<
255 <<
256 G4cout.precision(4); <<
257 G4cout << "\n StepSize of ch. tracks in wall <<
258 << G4BestUnit(rms, "Length") << "\t ( <<
259 << ")"; <<
260 <<
261 fStepCavity /= fNbStepCavity; <<
262 fStepCavity2 /= fNbStepCavity; <<
263 rms = fStepCavity2 - fStepCavity * fStepCavi <<
264 if (rms > 0.) <<
265 rms = std::sqrt(rms); <<
266 else <<
267 rms = 0.; <<
268 <<
269 G4cout << "\n StepSize of ch. tracks in cavi <<
270 << G4BestUnit(rms, "Length") <<
271 << "\t (nbSteps/track = " << double(f <<
272 <<
273 G4cout << G4endl; <<
274 175
275 // reset default formats << 176 //mean kinetic energy of secondary electrons
276 G4cout.setf(mode, std::ios::floatfield); << 177 //
277 G4cout.precision(prec); << 178 if (fNbSec == 0) return;
>> 179 G4double meanEsecond = fEsecondary/fNbSec,meanEsecond2= fEsecondary2/fNbSec;
>> 180 G4double varianceEsec = meanEsecond2 - meanEsecond*meanEsecond;
>> 181 G4double dToverT = 0.;
>> 182 if (varianceEsec>0.) dToverT = std::sqrt(varianceEsec/fNbSec)/meanEsecond;
>> 183 G4double csdaRange =
>> 184 emCalculator.GetCSDARange(meanEsecond,G4Electron::Electron(),fMateWall);
>> 185
>> 186 G4cout.precision(4);
>> 187 G4cout
>> 188 << "\n Mean energy of secondary e- = " << G4BestUnit(meanEsecond,"Energy")
>> 189 << " +- " << 100*dToverT << " %"
>> 190 << " (--> range in wall material = " << G4BestUnit(csdaRange,"Length")
>> 191 << ")"
>> 192 << G4endl;
>> 193
>> 194 //compute mass energy transfer coefficient
>> 195 //
>> 196 G4double massTransfCoef = sumc*meanEsecond/energy;
>> 197
>> 198 G4cout << " Mass_energy_transfer coef: "
>> 199 << G4BestUnit(massTransfCoef, "Surface/Mass")
>> 200 << G4endl;
>> 201
>> 202 //stopping power from EmCalculator
>> 203 //
>> 204 G4double dedxWall =
>> 205 emCalculator.GetDEDX(meanEsecond,G4Electron::Electron(),fMateWall);
>> 206 dedxWall /= fDensityWall;
>> 207 G4double dedxCavity =
>> 208 emCalculator.GetDEDX(meanEsecond,G4Electron::Electron(),fMateCavity);
>> 209 dedxCavity /= fDensityCavity;
>> 210
>> 211 G4cout
>> 212 << "\n StoppingPower in wall = "
>> 213 << G4BestUnit(dedxWall,"Energy*Surface/Mass")
>> 214 << "\n in cavity = "
>> 215 << G4BestUnit(dedxCavity,"Energy*Surface/Mass")
>> 216 << G4endl;
>> 217
>> 218 //charged particle flow in cavity
>> 219 //
>> 220 G4cout
>> 221 << "\n Charged particle flow in cavity :"
>> 222 << "\n Enter --> nbParticles = " << fPartFlowCavity[0]
>> 223 << "\t Energy = " << G4BestUnit (fEnerFlowCavity[0], "Energy")
>> 224 << "\n Exit --> nbParticles = " << fPartFlowCavity[1]
>> 225 << "\t Energy = " << G4BestUnit (fEnerFlowCavity[1], "Energy")
>> 226 << G4endl;
>> 227
>> 228 if (fPartFlowCavity[0] == 0) return;
>> 229
>> 230 //beam energy fluence
>> 231 //
>> 232 G4double rBeam = fWallRadius*(fKinematic->GetBeamRadius());
>> 233 G4double surfaceBeam = CLHEP::pi*rBeam*rBeam;
>> 234
>> 235 //error on Edep in cavity
>> 236 //
>> 237 if (fNbEventCavity == 0) return;
>> 238 G4double meanEdep = fEdepCavity/fNbEventCavity;
>> 239 G4double meanEdep2 = fEdepCavity2/fNbEventCavity;
>> 240 G4double varianceEdep = meanEdep2 - meanEdep*meanEdep;
>> 241 G4double dEoverE = 0.;
>> 242 if(varianceEdep>0.) dEoverE=std::sqrt(varianceEdep/fNbEventCavity)/meanEdep;
>> 243
>> 244 //total dose in cavity
>> 245 //
>> 246 G4double doseCavity = fEdepCavity/fMassCavity;
>> 247 G4double doseOverBeam = doseCavity*surfaceBeam/(numberOfEvent*energy);
>> 248
>> 249 //track length in cavity
>> 250 G4double meantrack = fTrkSegmCavity/fPartFlowCavity[0];
>> 251
>> 252 G4cout.precision(4);
>> 253 G4cout
>> 254 << "\n Total edep in cavity = " << G4BestUnit(fEdepCavity,"Energy")
>> 255 << " +- " << 100*dEoverE << " %"
>> 256 << "\t Total charged trackLength = " <<G4BestUnit(fTrkSegmCavity,"Length")
>> 257 << " (mean value = " << G4BestUnit(meantrack,"Length") << ")"
>> 258 << "\n Total dose in cavity = " << doseCavity/(MeV/mg) << " MeV/mg"
>> 259 << "\n Dose/EnergyFluence = " << G4BestUnit(doseOverBeam,"Surface/Mass")
>> 260 << G4endl;
>> 261
>> 262 //ratio simulation/theory
>> 263 //
>> 264 G4double ratio = doseOverBeam/massTransfCoef;
>> 265 G4double error = ratio*std::sqrt(dEoverE*dEoverE + dToverT*dToverT);
>> 266
>> 267 G4cout.precision(5);
>> 268 G4cout
>> 269 << "\n (Dose/EnergyFluence)/Mass_energy_transfer = " << ratio
>> 270 << " +- " << error << G4endl;
>> 271
>> 272 //compute mean step size of charged particles
>> 273 //
>> 274 fStepWall /= fNbStepWall; fStepWall2 /= fNbStepWall;
>> 275 G4double rms = fStepWall2 - fStepWall*fStepWall;
>> 276 if (rms>0.) rms = std::sqrt(rms); else rms = 0.;
>> 277
>> 278 G4cout.precision(4);
>> 279 G4cout
>> 280 << "\n StepSize of ch. tracks in wall = "
>> 281 << G4BestUnit(fStepWall,"Length") << " +- " << G4BestUnit( rms,"Length")
>> 282 << "\t (nbSteps/track = " << double(fNbStepWall)/fNbSec << ")";
>> 283
>> 284 fStepCavity /= fNbStepCavity; fStepCavity2 /= fNbStepCavity;
>> 285 rms = fStepCavity2 - fStepCavity*fStepCavity;
>> 286 if (rms>0.) rms = std::sqrt(rms); else rms = 0.;
>> 287
>> 288 G4cout
>> 289 << "\n StepSize of ch. tracks in cavity = "
>> 290 << G4BestUnit(fStepCavity,"Length") << " +- " << G4BestUnit( rms,"Length")
>> 291 << "\t (nbSteps/track = "<<double(fNbStepCavity)/fPartFlowCavity[0] << ")";
>> 292
>> 293 G4cout << G4endl;
>> 294
>> 295 // reset default formats
>> 296 G4cout.setf(mode,std::ios::floatfield);
>> 297 G4cout.precision(prec);
>> 298
>> 299 // delete and remove all contents in fProcCounter
>> 300 while (fProcCounter->size()>0){
>> 301 OneProcessCount* aProcCount=fProcCounter->back();
>> 302 fProcCounter->pop_back();
>> 303 delete aProcCount;
>> 304 }
>> 305 delete fProcCounter;
278 306
279 // delete and remove all contents in fProcCo <<
280 while (fProcCounter->size() > 0) { <<
281 OneProcessCount* aProcCount = fProcCounter <<
282 fProcCounter->pop_back(); <<
283 delete aProcCount; <<
284 } <<
285 delete fProcCounter; <<
286 } 307 }
287 308
288 //....oooOO0OOooo........oooOO0OOooo........oo 309 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
289 <<
290 void Run::SurveyConvergence(G4int NbofEvents) 310 void Run::SurveyConvergence(G4int NbofEvents)
291 { 311 {
292 if (NbofEvents == 0) return; 312 if (NbofEvents == 0) return;
293 313
294 // mean kinetic energy of secondary electron << 314 //mean kinetic energy of secondary electrons
295 // 315 //
296 G4double meanEsecond = 0.; 316 G4double meanEsecond = 0.;
297 if (fNbSec > 0) meanEsecond = fEsecondary / << 317 if (fNbSec > 0) meanEsecond = fEsecondary/fNbSec;
298 G4double rateEmean = 0.; 318 G4double rateEmean = 0.;
299 // compute variation rate (%), iteration to 319 // compute variation rate (%), iteration to iteration
300 if (fOldEmean > 0.) rateEmean = 100 * (meanE << 320 if (fOldEmean > 0.) rateEmean = 100*(meanEsecond/fOldEmean - 1.);
301 fOldEmean = meanEsecond; 321 fOldEmean = meanEsecond;
302 322
303 // beam energy fluence << 323 //beam energy fluence
304 // 324 //
305 G4double rBeam = fWallRadius * (fKinematic-> << 325 G4double rBeam = fWallRadius*(fKinematic->GetBeamRadius());
306 G4double surfaceBeam = CLHEP::pi * rBeam * r << 326 G4double surfaceBeam = CLHEP::pi*rBeam*rBeam;
307 G4double beamEnergy = fKinematic->GetParticl 327 G4double beamEnergy = fKinematic->GetParticleGun()->GetParticleEnergy();
308 328
309 // total dose in cavity << 329 //total dose in cavity
310 // 330 //
311 G4double doseCavity = fEdepCavity / fMassCav << 331 G4double doseCavity = fEdepCavity/fMassCavity;
312 G4double doseOverBeam = doseCavity * surface << 332 G4double doseOverBeam = doseCavity*surfaceBeam/(NbofEvents*beamEnergy);
313 G4double rateDose = 0.; 333 G4double rateDose = 0.;
314 // compute variation rate (%), iteration to 334 // compute variation rate (%), iteration to iteration
315 if (fOldDose > 0.) rateDose = 100 * (doseOve << 335 if (fOldDose > 0.) rateDose = 100*(doseOverBeam/fOldDose - 1.);
316 fOldDose = doseOverBeam; 336 fOldDose = doseOverBeam;
317 337
318 std::ios::fmtflags mode = G4cout.flags(); 338 std::ios::fmtflags mode = G4cout.flags();
319 G4cout.setf(std::ios::fixed, std::ios::float << 339 G4cout.setf(std::ios::fixed,std::ios::floatfield);
320 G4int prec = G4cout.precision(3); 340 G4int prec = G4cout.precision(3);
321 341
322 G4cout << " ---> NbofEvents= " << NbofEvents << 342 G4cout << "\n ---> NbofEvents= " << NbofEvents
323 << " Tkin= " << G4BestUnit(meanEsec << 343 << " NbOfelectr= " << fNbSec
>> 344 << " Tkin= " << G4BestUnit(meanEsecond,"Energy")
>> 345 << " (" << rateEmean << " %)"
324 << " NbOfelec in cav= " << fPartFlo 346 << " NbOfelec in cav= " << fPartFlowCavity[0]
325 << " Dose/EnFluence= " << G4BestUni << 347 << " Dose/EnFluence= " << G4BestUnit(doseOverBeam,"Surface/Mass")
326 << " %) \n" << 348 << " (" << rateDose << " %)"
327 << G4endl; 349 << G4endl;
328 350
329 // reset default formats 351 // reset default formats
330 G4cout.setf(mode, std::ios::floatfield); << 352 G4cout.setf(mode,std::ios::floatfield);
331 G4cout.precision(prec); 353 G4cout.precision(prec);
332 } 354 }
333 355
334 //....oooOO0OOooo........oooOO0OOooo........oo 356 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
335 357
336 void Run::Merge(const G4Run* run) << 358 void Run::Merge(const G4Run* run) {
337 { <<
338 const Run* localRun = static_cast<const Run* <<
339 359
340 fPartFlowCavity[0] += localRun->fPartFlowCav << 360 const Run* localRun = static_cast<const Run*>(run);
341 fPartFlowCavity[1] += localRun->fPartFlowCav <<
342 fEnerFlowCavity[0] += localRun->fEnerFlowCav <<
343 fEnerFlowCavity[1] += localRun->fEnerFlowCav <<
344 fEdepCavity += localRun->fEdepCavity; <<
345 fEdepCavity2 += localRun->fEdepCavity2; <<
346 fTrkSegmCavity += localRun->fTrkSegmCavity; <<
347 fNbEventCavity += localRun->fNbEventCavity; <<
348 <<
349 fStepWall += localRun->fStepWall; <<
350 fStepWall2 += localRun->fStepWall2; <<
351 fStepCavity += localRun->fStepCavity; <<
352 fStepCavity2 += localRun->fStepCavity2; <<
353 fNbStepWall += localRun->fNbStepWall; <<
354 fNbStepCavity += localRun->fNbStepCavity; <<
355 361
356 fEsecondary += localRun->fEsecondary; << 362 fPartFlowCavity[0]+= localRun->fPartFlowCavity[0];
357 fEsecondary2 += localRun->fEsecondary2; << 363 fPartFlowCavity[1]+= localRun->fPartFlowCavity[1];
>> 364 fEnerFlowCavity[0]+= localRun->fEnerFlowCavity[0];
>> 365 fEnerFlowCavity[1]+= localRun->fEnerFlowCavity[1];
>> 366 fEdepCavity += localRun->fEdepCavity;
>> 367 fEdepCavity2 += localRun->fEdepCavity2;
>> 368 fTrkSegmCavity += localRun->fTrkSegmCavity;
>> 369 fNbEventCavity += localRun->fNbEventCavity;
>> 370
>> 371 fStepWall += localRun->fStepWall;
>> 372 fStepWall2 += localRun->fStepWall2;
>> 373 fStepCavity += localRun->fStepCavity;
>> 374 fStepCavity2 += localRun->fStepCavity2;
>> 375 fNbStepWall += localRun->fNbStepWall;
>> 376 fNbStepCavity += localRun->fNbStepCavity;
358 377
359 fNbSec += localRun->fNbSec; << 378 fEsecondary += localRun->fEsecondary;
>> 379 fEsecondary2 += localRun->fEsecondary2;
360 380
361 // ??? G4double fOldEmean << 381 fNbSec += localRun->fNbSec;
>> 382
>> 383 // ??? G4double fOldEmean
362 // ??? G4Double fOldDose; 384 // ??? G4Double fOldDose;
363 385
364 // Merge ProcessCount varaibles 386 // Merge ProcessCount varaibles
365 std::vector<OneProcessCount*>::iterator it; 387 std::vector<OneProcessCount*>::iterator it;
366 for (it = localRun->fProcCounter->begin(); i << 388 for ( it = localRun->fProcCounter->begin();it !=localRun->fProcCounter->end();
>> 389 it++ )
>> 390 {
367 OneProcessCount* process = *it; 391 OneProcessCount* process = *it;
368 for (G4int i = 0; i < process->GetCounter( << 392 for ( G4int i = 0; i < process->GetCounter() ; i++)
369 this->CountProcesses(process->GetName()) 393 this->CountProcesses(process->GetName());
370 } 394 }
371 395
372 G4Run::Merge(run); 396 G4Run::Merge(run);
>> 397
373 } 398 }
374 399
375 //....oooOO0OOooo........oooOO0OOooo........oo 400 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
376 401
377 void Run::CountProcesses(G4String procName) 402 void Run::CountProcesses(G4String procName)
378 { 403 {
379 // does the process already encounted ? << 404 //does the process already encounted ?
380 size_t nbProc = fProcCounter->size(); << 405 size_t nbProc = fProcCounter->size();
381 size_t i = 0; << 406 size_t i = 0;
382 while ((i < nbProc) && ((*fProcCounter)[i]-> << 407 while ((i<nbProc)&&((*fProcCounter)[i]->GetName()!=procName)) i++;
383 i++; << 408 if (i == nbProc) fProcCounter->push_back( new OneProcessCount(procName));
384 if (i == nbProc) fProcCounter->push_back(new <<
385 409
386 (*fProcCounter)[i]->Count(); << 410 (*fProcCounter)[i]->Count();
387 } 411 }
388 412
389 //....oooOO0OOooo........oooOO0OOooo........oo 413 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
390 414