| 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 Run.cc 26 /// \file Run.cc
27 /// \brief Implementation of the Run class 27 /// \brief Implementation of the Run class
28 // 28 //
29 // << 29 //
30 30
31 #include "Run.hh" 31 #include "Run.hh"
32 <<
33 #include "DetectorConstruction.hh" 32 #include "DetectorConstruction.hh"
34 #include "HistoManager.hh" <<
35 #include "PrimaryGeneratorAction.hh" 33 #include "PrimaryGeneratorAction.hh"
>> 34 #include "HistoManager.hh"
36 35
37 #include "G4ProcessTable.hh" 36 #include "G4ProcessTable.hh"
38 #include "G4Radioactivation.hh" 37 #include "G4Radioactivation.hh"
39 #include "G4SystemOfUnits.hh" <<
40 #include "G4TwoVector.hh" 38 #include "G4TwoVector.hh"
41 #include "G4UnitsTable.hh" 39 #include "G4UnitsTable.hh"
>> 40 #include "G4SystemOfUnits.hh"
42 41
43 #include <fstream> 42 #include <fstream>
44 43
45 //....oooOO0OOooo........oooOO0OOooo........oo 44 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
46 45
47 Run::Run(DetectorConstruction* det) : fDetecto << 46 Run::Run(DetectorConstruction* det)
>> 47 : G4Run(),
>> 48 fDetector(det), fParticle(0), fEkin(0.)
>> 49 {
>> 50 fEdepTarget = fEdepTarget2 = 0.;
>> 51 fEdepDetect = fEdepDetect2 = 0.;
>> 52 }
>> 53
>> 54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 55
>> 56 Run::~Run()
>> 57 { }
48 58
49 //....oooOO0OOooo........oooOO0OOooo........oo 59 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
50 60
51 void Run::SetPrimary(G4ParticleDefinition* par 61 void Run::SetPrimary(G4ParticleDefinition* particle, G4double energy)
52 { << 62 {
53 fParticle = particle; 63 fParticle = particle;
54 fEkin = energy; 64 fEkin = energy;
55 } 65 }
56 << 66
57 //....oooOO0OOooo........oooOO0OOooo........oo 67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
58 68
59 void Run::CountProcesses(const G4VProcess* pro << 69 void Run::CountProcesses(const G4VProcess* process, G4int iVol)
60 { 70 {
61 if (process == nullptr) return; << 71 if (process == nullptr) return;
62 G4String procName = process->GetProcessName( 72 G4String procName = process->GetProcessName();
63 if (iVol == 1) { 73 if (iVol == 1) {
64 std::map<G4String, G4int>::iterator it1 = << 74 std::map<G4String,G4int>::iterator it1 = fProcCounter1.find(procName);
65 if (it1 == fProcCounter1.end()) { << 75 if ( it1 == fProcCounter1.end()) {
66 fProcCounter1[procName] = 1; 76 fProcCounter1[procName] = 1;
67 } 77 }
68 else { 78 else {
69 fProcCounter1[procName]++; << 79 fProcCounter1[procName]++;
70 } 80 }
71 } 81 }
72 << 82
73 if (iVol == 2) { 83 if (iVol == 2) {
74 std::map<G4String, G4int>::iterator it2 = << 84 std::map<G4String,G4int>::iterator it2 = fProcCounter2.find(procName);
75 if (it2 == fProcCounter2.end()) { << 85 if ( it2 == fProcCounter2.end()) {
76 fProcCounter2[procName] = 1; 86 fProcCounter2[procName] = 1;
77 } 87 }
78 else { 88 else {
79 fProcCounter2[procName]++; << 89 fProcCounter2[procName]++;
80 } << 90 }
81 } 91 }
82 } 92 }
83 << 93
84 //....oooOO0OOooo........oooOO0OOooo........oo 94 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
85 95
86 void Run::ParticleCount(G4String name, G4doubl 96 void Run::ParticleCount(G4String name, G4double Ekin, G4int iVol)
87 { 97 {
88 if (iVol == 1) { 98 if (iVol == 1) {
89 std::map<G4String, ParticleData>::iterator << 99 std::map<G4String, ParticleData>::iterator it = fParticleDataMap1.find(name);
90 if (it == fParticleDataMap1.end()) { << 100 if ( it == fParticleDataMap1.end()) {
91 fParticleDataMap1[name] = ParticleData(1 << 101 fParticleDataMap1[name] = ParticleData(1, Ekin, Ekin, Ekin);
92 } << 102 }
93 else { << 103 else {
94 ParticleData& data = it->second; << 104 ParticleData& data = it->second;
95 data.fCount++; << 105 data.fCount++;
96 data.fEmean += Ekin; << 106 data.fEmean += Ekin;
97 // update min max << 107 //update min max
98 G4double emin = data.fEmin; << 108 G4double emin = data.fEmin;
99 if (Ekin < emin) data.fEmin = Ekin; << 109 if (Ekin < emin) data.fEmin = Ekin;
100 G4double emax = data.fEmax; << 110 G4double emax = data.fEmax;
101 if (Ekin > emax) data.fEmax = Ekin; << 111 if (Ekin > emax) data.fEmax = Ekin;
102 } << 112 }
103 } 113 }
104 << 114
105 if (iVol == 2) { 115 if (iVol == 2) {
106 std::map<G4String, ParticleData>::iterator << 116 std::map<G4String, ParticleData>::iterator it = fParticleDataMap2.find(name);
107 if (it == fParticleDataMap2.end()) { << 117 if ( it == fParticleDataMap2.end()) {
108 fParticleDataMap2[name] = ParticleData(1 << 118 fParticleDataMap2[name] = ParticleData(1, Ekin, Ekin, Ekin);
109 } << 119 }
110 else { << 120 else {
111 ParticleData& data = it->second; << 121 ParticleData& data = it->second;
112 data.fCount++; << 122 data.fCount++;
113 data.fEmean += Ekin; << 123 data.fEmean += Ekin;
114 // update min max << 124 //update min max
115 G4double emin = data.fEmin; << 125 G4double emin = data.fEmin;
116 if (Ekin < emin) data.fEmin = Ekin; << 126 if (Ekin < emin) data.fEmin = Ekin;
117 G4double emax = data.fEmax; << 127 G4double emax = data.fEmax;
118 if (Ekin > emax) data.fEmax = Ekin; << 128 if (Ekin > emax) data.fEmax = Ekin;
119 } << 129 }
120 } << 130 }
121 } 131 }
122 << 132
123 //....oooOO0OOooo........oooOO0OOooo........oo 133 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
124 134
125 void Run::AddEdep(G4double edep1, G4double ede 135 void Run::AddEdep(G4double edep1, G4double edep2)
126 { << 136 {
127 fEdepTarget += edep1; << 137 fEdepTarget += edep1;
128 fEdepTarget2 += edep1 * edep1; << 138 fEdepTarget2 += edep1*edep1;
129 fEdepDetect += edep2; << 139 fEdepDetect += edep2;
130 fEdepDetect2 += edep2 * edep2; << 140 fEdepDetect2 += edep2*edep2;
131 } 141 }
132 << 142
133 //....oooOO0OOooo........oooOO0OOooo........oo 143 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
134 144
135 void Run::Merge(const G4Run* run) 145 void Run::Merge(const G4Run* run)
136 { 146 {
137 const Run* localRun = static_cast<const Run* 147 const Run* localRun = static_cast<const Run*>(run);
138 << 148
139 // primary particle info << 149 //primary particle info
140 // 150 //
141 fParticle = localRun->fParticle; 151 fParticle = localRun->fParticle;
142 fEkin = localRun->fEkin; << 152 fEkin = localRun->fEkin;
143 << 153
144 // accumulate sums 154 // accumulate sums
145 // 155 //
146 fEdepTarget += localRun->fEdepTarget; << 156 fEdepTarget += localRun->fEdepTarget;
147 fEdepTarget2 += localRun->fEdepTarget2; << 157 fEdepTarget2 += localRun->fEdepTarget2;
148 fEdepDetect += localRun->fEdepDetect; << 158 fEdepDetect += localRun->fEdepDetect;
149 fEdepDetect2 += localRun->fEdepDetect2; << 159 fEdepDetect2 += localRun->fEdepDetect2;
>> 160
>> 161 //map: processes count in target
>> 162
>> 163 std::map<G4String,G4int>::const_iterator itp1;
>> 164 for ( itp1 = localRun->fProcCounter1.begin();
>> 165 itp1 != localRun->fProcCounter1.end(); ++itp1 ) {
150 166
151 // map: processes count in target <<
152 <<
153 std::map<G4String, G4int>::const_iterator it <<
154 for (itp1 = localRun->fProcCounter1.begin(); <<
155 G4String procName = itp1->first; 167 G4String procName = itp1->first;
156 G4int localCount = itp1->second; 168 G4int localCount = itp1->second;
157 if (fProcCounter1.find(procName) == fProcC << 169 if ( fProcCounter1.find(procName) == fProcCounter1.end()) {
158 fProcCounter1[procName] = localCount; 170 fProcCounter1[procName] = localCount;
159 } 171 }
160 else { 172 else {
161 fProcCounter1[procName] += localCount; 173 fProcCounter1[procName] += localCount;
162 } << 174 }
163 } 175 }
>> 176
>> 177 //map: processes count in detector
>> 178
>> 179 std::map<G4String,G4int>::const_iterator itp2;
>> 180 for ( itp2 = localRun->fProcCounter2.begin();
>> 181 itp2 != localRun->fProcCounter2.end(); ++itp2 ) {
164 182
165 // map: processes count in detector <<
166 <<
167 std::map<G4String, G4int>::const_iterator it <<
168 for (itp2 = localRun->fProcCounter2.begin(); <<
169 G4String procName = itp2->first; 183 G4String procName = itp2->first;
170 G4int localCount = itp2->second; 184 G4int localCount = itp2->second;
171 if (fProcCounter2.find(procName) == fProcC << 185 if ( fProcCounter2.find(procName) == fProcCounter2.end()) {
172 fProcCounter2[procName] = localCount; 186 fProcCounter2[procName] = localCount;
173 } 187 }
174 else { 188 else {
175 fProcCounter2[procName] += localCount; 189 fProcCounter2[procName] += localCount;
176 } << 190 }
177 } 191 }
178 << 192
179 // map: created particles in target << 193 //map: created particles in target
180 std::map<G4String, ParticleData>::const_iter << 194 std::map<G4String,ParticleData>::const_iterator itc;
181 for (itc = localRun->fParticleDataMap1.begin << 195 for (itc = localRun->fParticleDataMap1.begin();
>> 196 itc != localRun->fParticleDataMap1.end(); ++itc) {
>> 197
182 G4String name = itc->first; 198 G4String name = itc->first;
183 const ParticleData& localData = itc->secon << 199 const ParticleData& localData = itc->second;
184 if (fParticleDataMap1.find(name) == fParti << 200 if ( fParticleDataMap1.find(name) == fParticleDataMap1.end()) {
185 fParticleDataMap1[name] = << 201 fParticleDataMap1[name]
186 ParticleData(localData.fCount, localDa << 202 = ParticleData(localData.fCount,
>> 203 localData.fEmean,
>> 204 localData.fEmin,
>> 205 localData.fEmax);
187 } 206 }
188 else { 207 else {
189 ParticleData& data = fParticleDataMap1[n << 208 ParticleData& data = fParticleDataMap1[name];
190 data.fCount += localData.fCount; 209 data.fCount += localData.fCount;
191 data.fEmean += localData.fEmean; 210 data.fEmean += localData.fEmean;
192 G4double emin = localData.fEmin; 211 G4double emin = localData.fEmin;
193 if (emin < data.fEmin) data.fEmin = emin 212 if (emin < data.fEmin) data.fEmin = emin;
194 G4double emax = localData.fEmax; 213 G4double emax = localData.fEmax;
195 if (emax > data.fEmax) data.fEmax = emax << 214 if (emax > data.fEmax) data.fEmax = emax;
196 } << 215 }
197 } 216 }
198 << 217
199 // map: created particle in detector << 218 //map: created particle in detector
200 std::map<G4String, ParticleData>::const_iter << 219 std::map<G4String,ParticleData>::const_iterator itn;
201 for (itn = localRun->fParticleDataMap2.begin << 220 for (itn = localRun->fParticleDataMap2.begin();
>> 221 itn != localRun->fParticleDataMap2.end(); ++itn) {
>> 222
202 G4String name = itn->first; 223 G4String name = itn->first;
203 const ParticleData& localData = itn->secon << 224 const ParticleData& localData = itn->second;
204 if (fParticleDataMap2.find(name) == fParti << 225 if ( fParticleDataMap2.find(name) == fParticleDataMap2.end()) {
205 fParticleDataMap2[name] = << 226 fParticleDataMap2[name]
206 ParticleData(localData.fCount, localDa << 227 = ParticleData(localData.fCount,
>> 228 localData.fEmean,
>> 229 localData.fEmin,
>> 230 localData.fEmax);
207 } 231 }
208 else { 232 else {
209 ParticleData& data = fParticleDataMap2[n << 233 ParticleData& data = fParticleDataMap2[name];
210 data.fCount += localData.fCount; 234 data.fCount += localData.fCount;
211 data.fEmean += localData.fEmean; 235 data.fEmean += localData.fEmean;
212 G4double emin = localData.fEmin; 236 G4double emin = localData.fEmin;
213 if (emin < data.fEmin) data.fEmin = emin 237 if (emin < data.fEmin) data.fEmin = emin;
214 G4double emax = localData.fEmax; 238 G4double emax = localData.fEmax;
215 if (emax > data.fEmax) data.fEmax = emax << 239 if (emax > data.fEmax) data.fEmax = emax;
216 } << 240 }
217 } 241 }
218 242
219 G4Run::Merge(run); << 243 G4Run::Merge(run);
220 } << 244 }
221 245
222 //....oooOO0OOooo........oooOO0OOooo........oo 246 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
223 247
224 void Run::EndOfRun() << 248 void Run::EndOfRun()
225 { 249 {
226 G4int prec = 5, wid = prec + 2; << 250 G4int prec = 5, wid = prec + 2;
227 G4int dfprec = G4cout.precision(prec); 251 G4int dfprec = G4cout.precision(prec);
228 << 252
229 // run condition 253 // run condition
230 // << 254 //
231 G4String Particle = fParticle->GetParticleNa << 255 G4String Particle = fParticle->GetParticleName();
232 G4cout << "\n The run is " << numberOfEvent << 256 G4cout << "\n The run is " << numberOfEvent << " "<< Particle << " of "
233 << G4BestUnit(fEkin, "Energy") << " t << 257 << G4BestUnit(fEkin,"Energy") << " through : ";
234 << 258
235 G4cout << "\n Target : Length = " << G4Bes << 259 G4cout << "\n Target : Length = "
236 << " Radius = " << G4BestUnit(fDet << 260 << G4BestUnit(fDetector->GetTargetLength(),"Length")
237 << " Material = " << fDetector->GetTa << 261 << " Radius = "
238 G4cout << "\n Detector : Length = " << G4Bes << 262 << G4BestUnit(fDetector->GetTargetRadius(),"Length")
239 << " Thickness = " << G4BestUnit(fDet << 263 << " Material = "
240 << " Material = " << fDetector->GetDe << 264 << fDetector->GetTargetMaterial()->GetName();
241 << 265 G4cout << "\n Detector : Length = "
242 if (numberOfEvent == 0) { << 266 << G4BestUnit(fDetector->GetDetectorLength(),"Length")
243 G4cout.precision(dfprec); << 267 << " Thickness = "
244 return; << 268 << G4BestUnit(fDetector->GetDetectorThickness(),"Length")
245 } << 269 << " Material = "
>> 270 << fDetector->GetDetectorMaterial()->GetName() << G4endl;
246 271
>> 272 if (numberOfEvent == 0) { G4cout.precision(dfprec); return;}
>> 273
247 // compute mean Energy deposited and rms in 274 // compute mean Energy deposited and rms in target
248 // 275 //
249 G4int TotNbofEvents = numberOfEvent; 276 G4int TotNbofEvents = numberOfEvent;
250 fEdepTarget /= TotNbofEvents; << 277 fEdepTarget /= TotNbofEvents; fEdepTarget2 /= TotNbofEvents;
251 fEdepTarget2 /= TotNbofEvents; << 278 G4double rmsEdep = fEdepTarget2 - fEdepTarget*fEdepTarget;
252 G4double rmsEdep = fEdepTarget2 - fEdepTarge << 279 if (rmsEdep>0.) rmsEdep = std::sqrt(rmsEdep);
253 if (rmsEdep > 0.) << 280 else rmsEdep = 0.;
254 rmsEdep = std::sqrt(rmsEdep); << 281
255 else <<
256 rmsEdep = 0.; <<
257 <<
258 G4cout << "\n Mean energy deposit in target, 282 G4cout << "\n Mean energy deposit in target, in time window = "
259 << G4BestUnit(fEdepTarget, "Energy") << 283 << G4BestUnit(fEdepTarget,"Energy") << "; rms = "
>> 284 << G4BestUnit(rmsEdep, "Energy")
260 << G4endl; 285 << G4endl;
261 286
262 // compute mean Energy deposited and rms in 287 // compute mean Energy deposited and rms in detector
263 // 288 //
264 fEdepDetect /= TotNbofEvents; << 289 fEdepDetect /= TotNbofEvents; fEdepDetect2 /= TotNbofEvents;
265 fEdepDetect2 /= TotNbofEvents; << 290 rmsEdep = fEdepDetect2 - fEdepDetect*fEdepDetect;
266 rmsEdep = fEdepDetect2 - fEdepDetect * fEdep << 291 if (rmsEdep>0.) rmsEdep = std::sqrt(rmsEdep);
267 if (rmsEdep > 0.) << 292 else rmsEdep = 0.;
268 rmsEdep = std::sqrt(rmsEdep); << 293
269 else <<
270 rmsEdep = 0.; <<
271 <<
272 G4cout << " Mean energy deposit in detector, 294 G4cout << " Mean energy deposit in detector, in time window = "
273 << G4BestUnit(fEdepDetect, "Energy") << 295 << G4BestUnit(fEdepDetect,"Energy") << "; rms = "
>> 296 << G4BestUnit(rmsEdep, "Energy")
274 << G4endl; 297 << G4endl;
275 298
276 // frequency of processes in target 299 // frequency of processes in target
277 // 300 //
278 G4cout << "\n Process calls frequency in tar 301 G4cout << "\n Process calls frequency in target :" << G4endl;
279 G4int index = 0; 302 G4int index = 0;
280 std::map<G4String, G4int>::iterator it1; << 303 std::map<G4String,G4int>::iterator it1;
281 for (it1 = fProcCounter1.begin(); it1 != fPr 304 for (it1 = fProcCounter1.begin(); it1 != fProcCounter1.end(); it1++) {
282 G4String procName = it1->first; << 305 G4String procName = it1->first;
283 G4int count = it1->second; << 306 G4int count = it1->second;
284 G4String space = " "; << 307 G4String space = " "; if (++index%3 == 0) space = "\n";
285 if (++index % 3 == 0) space = "\n"; << 308 G4cout << " " << std::setw(20) << procName << "="<< std::setw(7) << count
286 G4cout << " " << std::setw(20) << procName << 309 << space;
287 } 310 }
288 G4cout << G4endl; 311 G4cout << G4endl;
289 << 312
290 // frequency of processes in detector 313 // frequency of processes in detector
291 // 314 //
292 G4cout << "\n Process calls frequency in det 315 G4cout << "\n Process calls frequency in detector:" << G4endl;
293 index = 0; 316 index = 0;
294 std::map<G4String, G4int>::iterator it2; << 317 std::map<G4String,G4int>::iterator it2;
295 for (it2 = fProcCounter2.begin(); it2 != fPr 318 for (it2 = fProcCounter2.begin(); it2 != fProcCounter2.end(); it2++) {
296 G4String procName = it2->first; << 319 G4String procName = it2->first;
297 G4int count = it2->second; << 320 G4int count = it2->second;
298 G4String space = " "; << 321 G4String space = " "; if (++index%3 == 0) space = "\n";
299 if (++index % 3 == 0) space = "\n"; << 322 G4cout << " " << std::setw(20) << procName << "="<< std::setw(7) << count
300 G4cout << " " << std::setw(20) << procName << 323 << space;
301 } 324 }
302 G4cout << G4endl; 325 G4cout << G4endl;
303 << 326
304 // particles count in target 327 // particles count in target
305 // 328 //
306 G4cout << "\n List of generated particles in 329 G4cout << "\n List of generated particles in target:" << G4endl;
307 << 330
308 std::map<G4String, ParticleData>::iterator i << 331 std::map<G4String,ParticleData>::iterator itc;
309 for (itc = fParticleDataMap1.begin(); itc != 332 for (itc = fParticleDataMap1.begin(); itc != fParticleDataMap1.end(); itc++) {
310 G4String name = itc->first; << 333 G4String name = itc->first;
311 ParticleData data = itc->second; << 334 ParticleData data = itc->second;
312 G4int count = data.fCount; << 335 G4int count = data.fCount;
313 G4double eMean = data.fEmean / count; << 336 G4double eMean = data.fEmean/count;
314 G4double eMin = data.fEmin; << 337 G4double eMin = data.fEmin;
315 G4double eMax = data.fEmax; << 338 G4double eMax = data.fEmax;
316 << 339
317 G4cout << " " << std::setw(13) << name << 340 G4cout << " " << std::setw(13) << name << ": " << std::setw(7) << count
318 << " Emean = " << std::setw(wid) < << 341 << " Emean = " << std::setw(wid) << G4BestUnit(eMean, "Energy")
319 << G4BestUnit(eMin, "Energy") << " << 342 << "\t( " << G4BestUnit(eMin, "Energy")
320 } << 343 << " --> " << G4BestUnit(eMax, "Energy")
321 << 344 << ")" << G4endl;
322 // particles count in detector << 345 }
323 // << 346
324 G4cout << "\n List of generated particles in << 347 // particles count in detector
325 << 348 //
326 std::map<G4String, ParticleData>::iterator i << 349 G4cout << "\n List of generated particles in detector:" << G4endl;
327 for (itn = fParticleDataMap2.begin(); itn != << 350
>> 351 std::map<G4String,ParticleData>::iterator itn;
>> 352 for (itn = fParticleDataMap2.begin(); itn != fParticleDataMap2.end(); itn++) {
328 G4String name = itn->first; 353 G4String name = itn->first;
329 ParticleData data = itn->second; 354 ParticleData data = itn->second;
330 G4int count = data.fCount; 355 G4int count = data.fCount;
331 G4double eMean = data.fEmean / count; << 356 G4double eMean = data.fEmean/count;
332 G4double eMin = data.fEmin; 357 G4double eMin = data.fEmin;
333 G4double eMax = data.fEmax; 358 G4double eMax = data.fEmax;
334 << 359
335 G4cout << " " << std::setw(13) << name << 360 G4cout << " " << std::setw(13) << name << ": " << std::setw(7) << count
336 << " Emean = " << std::setw(wid) < << 361 << " Emean = " << std::setw(wid) << G4BestUnit(eMean, "Energy")
337 << G4BestUnit(eMin, "Energy") << " << 362 << "\t( " << G4BestUnit(eMin, "Energy")
>> 363 << " --> " << G4BestUnit(eMax, "Energy") << ")" << G4endl;
338 } 364 }
339 G4cout << G4endl; 365 G4cout << G4endl;
340 << 366
341 // activities in VR mode 367 // activities in VR mode
342 // 368 //
343 WriteActivity(numberOfEvent); 369 WriteActivity(numberOfEvent);
344 << 370
345 // remove all contents in fProcCounter, fCou << 371 //remove all contents in fProcCounter, fCount
346 fProcCounter1.clear(); 372 fProcCounter1.clear();
347 fProcCounter2.clear(); 373 fProcCounter2.clear();
348 fParticleDataMap1.clear(); << 374 fParticleDataMap1.clear();
349 fParticleDataMap2.clear(); 375 fParticleDataMap2.clear();
350 << 376
351 // restore default format << 377 //restore default format
352 G4cout.precision(dfprec); << 378 G4cout.precision(dfprec);
353 } 379 }
354 380
355 //....oooOO0OOooo........oooOO0OOooo........oo 381 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
356 382
357 void Run::WriteActivity(G4int nevent) 383 void Run::WriteActivity(G4int nevent)
358 { 384 {
359 G4ProcessTable* pTable = G4ProcessTable::Get << 385 G4ProcessTable *pTable = G4ProcessTable::GetProcessTable();
360 G4Radioactivation* rDecay = << 386 G4Radioactivation* rDecay = (G4Radioactivation *)
361 (G4Radioactivation*)pTable->FindProcess("R << 387 pTable->FindProcess("Radioactivation", "GenericIon");
362 << 388
363 // output the induced radioactivities (in VR << 389 // output the induced radioactivities (in VR mode only)
364 // << 390 //
365 if ((rDecay == 0) || (rDecay->IsAnalogueMont << 391 if ((rDecay == 0) || (rDecay->IsAnalogueMonteCarlo())) return;
>> 392
>> 393 G4String fileName = G4AnalysisManager::Instance()->GetFileName() + ".activity";
>> 394 std::ofstream outfile (fileName, std::ios::out );
>> 395
>> 396 std::vector<G4RadioactivityTable*> theTables =
>> 397 rDecay->GetTheRadioactivityTables();
366 398
367 G4String fileName = G4AnalysisManager::Insta << 399 for (size_t i = 0 ; i < theTables.size(); i++) {
368 std::ofstream outfile(fileName, std::ios::ou <<
369 <<
370 std::vector<G4RadioactivityTable*> theTables <<
371 <<
372 for (size_t i = 0; i < theTables.size(); i++ <<
373 G4double rate, error; 400 G4double rate, error;
374 outfile << "Radioactivities in decay windo 401 outfile << "Radioactivities in decay window no. " << i << G4endl;
375 outfile << "Z \tA \tE \tActivity (decays/w << 402 outfile << "Z \tA \tE \tActivity (decays/window) \tError (decays/window) "
>> 403 << G4endl;
376 404
377 map<G4ThreeVector, G4TwoVector>* aMap = th << 405 map<G4ThreeVector,G4TwoVector> *aMap = theTables[i]->GetTheMap();
378 map<G4ThreeVector, G4TwoVector>::iterator << 406 map<G4ThreeVector,G4TwoVector>::iterator iter;
379 for (iter = aMap->begin(); iter != aMap->e << 407 for (iter=aMap->begin(); iter != aMap->end(); iter++) {
380 rate = iter->second.x() / nevent; << 408 rate = iter->second.x()/nevent;
381 error = std::sqrt(iter->second.y()) / ne << 409 error = std::sqrt(iter->second.y())/nevent;
382 if (rate < 0.) rate = 0.; // statically << 410 if (rate < 0.) rate = 0.; // statically it can be < 0.
383 outfile << iter->first.x() << "\t" << it << 411 outfile << iter->first.x() <<"\t"<< iter->first.y() <<"\t"
384 << rate << "\t" << error << G4en << 412 << iter->first.z() << "\t" << rate <<"\t" << error << G4endl;
385 } 413 }
386 outfile << G4endl; 414 outfile << G4endl;
387 } << 415 }
388 outfile.close(); << 416 outfile.close();
389 } 417 }
390 418
391 //....oooOO0OOooo........oooOO0OOooo........oo 419 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
392 420