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