| 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 // $Id: Run.cc 71376 2013-06-14 07:44:50Z maire $
>> 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 "HistoManager.hh" <<
36 #include "PrimaryGeneratorAction.hh" 35 #include "PrimaryGeneratorAction.hh"
>> 36 #include "HistoManager.hh"
37 37
38 #include "G4PhysicalConstants.hh" <<
39 #include "G4SystemOfUnits.hh" 38 #include "G4SystemOfUnits.hh"
40 #include "G4UnitsTable.hh" 39 #include "G4UnitsTable.hh"
>> 40 #include "G4PhysicalConstants.hh"
41 41
42 //....oooOO0OOooo........oooOO0OOooo........oo 42 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
43 43
44 void Run::SetPrimary(G4ParticleDefinition* par << 44 Run::Run()
>> 45 : G4Run(),
>> 46 fParticle(0), fEkin(0.),
>> 47 fDecayCount(0), fTimeCount(0), fPrimaryTime(0.),
>> 48 fTimeWindow1(0.), fTimeWindow2(0.)
45 { 49 {
>> 50 fEkinTot[0] = fPbalance[0] = fEventTime[0] = fEvisEvent[0] = 0. ;
>> 51 fEkinTot[1] = fPbalance[1] = fEventTime[1] = fEvisEvent[1] = DBL_MAX;
>> 52 fEkinTot[2] = fPbalance[2] = fEventTime[2] = fEvisEvent[2] = 0. ;
>> 53 }
>> 54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 55
>> 56 Run::~Run()
>> 57 { }
>> 58
>> 59 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 60
>> 61 void Run::SetPrimary(G4ParticleDefinition* particle, G4double energy)
>> 62 {
46 fParticle = particle; 63 fParticle = particle;
47 fEkin = energy; 64 fEkin = energy;
48 } << 65 }
49 66
50 //....oooOO0OOooo........oooOO0OOooo........oo 67 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
51 68
52 void Run::ParticleCount(G4String name, G4doubl 69 void Run::ParticleCount(G4String name, G4double Ekin, G4double meanLife)
53 { 70 {
54 std::map<G4String, ParticleData>::iterator i 71 std::map<G4String, ParticleData>::iterator it = fParticleDataMap.find(name);
55 if (it == fParticleDataMap.end()) { << 72 if ( it == fParticleDataMap.end()) {
56 fParticleDataMap[name] = ParticleData(1, E 73 fParticleDataMap[name] = ParticleData(1, Ekin, Ekin, Ekin, meanLife);
57 } 74 }
58 else { 75 else {
59 ParticleData& data = it->second; 76 ParticleData& data = it->second;
60 data.fCount++; 77 data.fCount++;
61 data.fEmean += Ekin; 78 data.fEmean += Ekin;
62 // update min max << 79 //update min max
63 G4double emin = data.fEmin; 80 G4double emin = data.fEmin;
64 if (Ekin < emin) data.fEmin = Ekin; 81 if (Ekin < emin) data.fEmin = Ekin;
65 G4double emax = data.fEmax; 82 G4double emax = data.fEmax;
66 if (Ekin > emax) data.fEmax = Ekin; 83 if (Ekin > emax) data.fEmax = Ekin;
67 data.fTmean = meanLife; 84 data.fTmean = meanLife;
68 } << 85 }
69 } 86 }
70 87
71 //....oooOO0OOooo........oooOO0OOooo........oo 88 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
72 89
73 void Run::SetTimeWindow(G4double t1, G4double 90 void Run::SetTimeWindow(G4double t1, G4double t2)
74 { 91 {
75 fTimeWindow1 = t1; 92 fTimeWindow1 = t1;
76 fTimeWindow2 = t2; 93 fTimeWindow2 = t2;
77 } 94 }
78 95
79 //....oooOO0OOooo........oooOO0OOooo........oo 96 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
80 97
81 void Run::CountInTimeWindow(G4String name, G4b << 98 void Run::CountInTimeWindow(G4String name, G4bool life1,
>> 99 G4bool life2, G4bool decay)
82 { 100 {
83 std::map<G4String, ActivityData>::iterator i 101 std::map<G4String, ActivityData>::iterator it = fActivityMap.find(name);
84 if (it == fActivityMap.end()) { << 102 if ( it == fActivityMap.end()) {
85 G4int n1(0), n2(0), nd(0); 103 G4int n1(0), n2(0), nd(0);
86 if (life1) n1 = 1; << 104 if(life1) n1 = 1;
87 if (life2) n2 = 1; << 105 if(life2) n2 = 1;
88 if (decay) nd = 1; << 106 if(decay) nd = 1;
89 fActivityMap[name] = ActivityData(n1, n2, 107 fActivityMap[name] = ActivityData(n1, n2, nd);
90 } 108 }
91 else { 109 else {
92 ActivityData& data = it->second; 110 ActivityData& data = it->second;
93 if (life1) data.fNlife_t1++; << 111 if(life1) data.fNlife_t1++;
94 if (life2) data.fNlife_t2++; << 112 if(life2) data.fNlife_t2++;
95 if (decay) data.fNdecay_t1t2++; << 113 if(decay) data.fNdecay_t1t2++;
96 } 114 }
97 } 115 }
98 116
99 //....oooOO0OOooo........oooOO0OOooo........oo 117 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
100 118
101 void Run::Balance(G4double Ekin, G4double Pbal 119 void Run::Balance(G4double Ekin, G4double Pbal)
102 { 120 {
103 fDecayCount++; 121 fDecayCount++;
104 fEkinTot[0] += Ekin; 122 fEkinTot[0] += Ekin;
105 // update min max << 123 //update min max
106 if (fDecayCount == 1) fEkinTot[1] = fEkinTot 124 if (fDecayCount == 1) fEkinTot[1] = fEkinTot[2] = Ekin;
107 if (Ekin < fEkinTot[1]) fEkinTot[1] = Ekin; 125 if (Ekin < fEkinTot[1]) fEkinTot[1] = Ekin;
108 if (Ekin > fEkinTot[2]) fEkinTot[2] = Ekin; 126 if (Ekin > fEkinTot[2]) fEkinTot[2] = Ekin;
109 << 127
110 fPbalance[0] += Pbal; 128 fPbalance[0] += Pbal;
111 // update min max << 129 //update min max
112 if (fDecayCount == 1) fPbalance[1] = fPbalan << 130 if (fDecayCount == 1) fPbalance[1] = fPbalance[2] = Pbal;
113 if (Pbal < fPbalance[1]) fPbalance[1] = Pbal 131 if (Pbal < fPbalance[1]) fPbalance[1] = Pbal;
114 if (Pbal > fPbalance[2]) fPbalance[2] = Pbal << 132 if (Pbal > fPbalance[2]) fPbalance[2] = Pbal;
115 } 133 }
116 134
117 //....oooOO0OOooo........oooOO0OOooo........oo 135 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
118 136
119 void Run::EventTiming(G4double time) 137 void Run::EventTiming(G4double time)
120 { 138 {
121 fTimeCount++; << 139 fTimeCount++;
122 fEventTime[0] += time; 140 fEventTime[0] += time;
123 if (fTimeCount == 1) fEventTime[1] = fEventT << 141 if (fTimeCount == 1) fEventTime[1] = fEventTime[2] = time;
124 if (time < fEventTime[1]) fEventTime[1] = ti 142 if (time < fEventTime[1]) fEventTime[1] = time;
125 if (time > fEventTime[2]) fEventTime[2] = ti << 143 if (time > fEventTime[2]) fEventTime[2] = time;
126 } 144 }
127 145
128 //....oooOO0OOooo........oooOO0OOooo........oo 146 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
129 147
130 void Run::PrimaryTiming(G4double ptime) 148 void Run::PrimaryTiming(G4double ptime)
131 { 149 {
132 fPrimaryTime += ptime; 150 fPrimaryTime += ptime;
133 } 151 }
134 152
135 //....oooOO0OOooo........oooOO0OOooo........oo 153 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
136 154
137 void Run::EvisEvent(G4double Evis) 155 void Run::EvisEvent(G4double Evis)
138 { 156 {
139 fEvisEvent[0] += Evis; 157 fEvisEvent[0] += Evis;
140 if (fTimeCount == 1) fEvisEvent[1] = fEvisEv << 158 if (fTimeCount == 1) fEvisEvent[1] = fEvisEvent[2] = Evis;
141 if (Evis < fEvisEvent[1]) fEvisEvent[1] = Ev 159 if (Evis < fEvisEvent[1]) fEvisEvent[1] = Evis;
142 if (Evis > fEvisEvent[2]) fEvisEvent[2] = Ev << 160 if (Evis > fEvisEvent[2]) fEvisEvent[2] = Evis;
143 } 161 }
144 162
145 //....oooOO0OOooo........oooOO0OOooo........oo 163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
146 164
147 void Run::Merge(const G4Run* run) 165 void Run::Merge(const G4Run* run)
148 { 166 {
149 const Run* localRun = static_cast<const Run* 167 const Run* localRun = static_cast<const Run*>(run);
150 168
151 // primary particle info << 169 //primary particle info
152 // 170 //
153 fParticle = localRun->fParticle; 171 fParticle = localRun->fParticle;
154 fEkin = localRun->fEkin; << 172 fEkin = localRun->fEkin;
155 << 173
156 // accumulate sums 174 // accumulate sums
157 // 175 //
158 fDecayCount += localRun->fDecayCount; << 176 fDecayCount += localRun->fDecayCount;
159 fTimeCount += localRun->fTimeCount; << 177 fTimeCount += localRun->fTimeCount;
160 fPrimaryTime += localRun->fPrimaryTime; 178 fPrimaryTime += localRun->fPrimaryTime;
161 179
162 fEkinTot[0] += localRun->fEkinTot[0]; << 180 fEkinTot[0] += localRun->fEkinTot[0];
163 fPbalance[0] += localRun->fPbalance[0]; << 181 fPbalance[0] += localRun->fPbalance[0];
164 fEventTime[0] += localRun->fEventTime[0]; 182 fEventTime[0] += localRun->fEventTime[0];
165 fEvisEvent[0] += localRun->fEvisEvent[0]; << 183 fEvisEvent[0] += localRun->fEvisEvent[0];
166 << 184
167 G4double min, max; << 185 G4double min,max;
168 min = localRun->fEkinTot[1]; << 186 min = localRun->fEkinTot[1]; max = localRun->fEkinTot[2];
169 max = localRun->fEkinTot[2]; <<
170 if (fEkinTot[1] > min) fEkinTot[1] = min; 187 if (fEkinTot[1] > min) fEkinTot[1] = min;
171 if (fEkinTot[2] < max) fEkinTot[2] = max; 188 if (fEkinTot[2] < max) fEkinTot[2] = max;
172 // 189 //
173 min = localRun->fPbalance[1]; << 190 min = localRun->fPbalance[1]; max = localRun->fPbalance[2];
174 max = localRun->fPbalance[2]; <<
175 if (fPbalance[1] > min) fPbalance[1] = min; 191 if (fPbalance[1] > min) fPbalance[1] = min;
176 if (fPbalance[2] < max) fPbalance[2] = max; 192 if (fPbalance[2] < max) fPbalance[2] = max;
177 // 193 //
178 min = localRun->fEventTime[1]; << 194 min = localRun->fEventTime[1]; max = localRun->fEventTime[2];
179 max = localRun->fEventTime[2]; <<
180 if (fEventTime[1] > min) fEventTime[1] = min 195 if (fEventTime[1] > min) fEventTime[1] = min;
181 if (fEventTime[2] < max) fEventTime[2] = max 196 if (fEventTime[2] < max) fEventTime[2] = max;
182 // 197 //
183 min = localRun->fEvisEvent[1]; << 198 min = localRun->fEvisEvent[1]; max = localRun->fEvisEvent[2];
184 max = localRun->fEvisEvent[2]; <<
185 if (fEvisEvent[1] > min) fEvisEvent[1] = min 199 if (fEvisEvent[1] > min) fEvisEvent[1] = min;
186 if (fEvisEvent[2] < max) fEvisEvent[2] = max 200 if (fEvisEvent[2] < max) fEvisEvent[2] = max;
187 << 201
188 // maps << 202 //maps
189 std::map<G4String, ParticleData>::const_iter << 203 std::map<G4String,ParticleData>::const_iterator itn;
190 for (itn = localRun->fParticleDataMap.begin( << 204 for (itn = localRun->fParticleDataMap.begin();
>> 205 itn != localRun->fParticleDataMap.end(); ++itn) {
>> 206
191 G4String name = itn->first; 207 G4String name = itn->first;
192 const ParticleData& localData = itn->secon << 208 const ParticleData& localData = itn->second;
193 if (fParticleDataMap.find(name) == fPartic << 209 if ( fParticleDataMap.find(name) == fParticleDataMap.end()) {
194 fParticleDataMap[name] = ParticleData(lo << 210 fParticleDataMap[name]
195 lo << 211 = ParticleData(localData.fCount,
>> 212 localData.fEmean,
>> 213 localData.fEmin,
>> 214 localData.fEmax,
>> 215 localData.fTmean);
196 } 216 }
197 else { 217 else {
198 ParticleData& data = fParticleDataMap[na << 218 ParticleData& data = fParticleDataMap[name];
199 data.fCount += localData.fCount; 219 data.fCount += localData.fCount;
200 data.fEmean += localData.fEmean; 220 data.fEmean += localData.fEmean;
201 G4double emin = localData.fEmin; 221 G4double emin = localData.fEmin;
202 if (emin < data.fEmin) data.fEmin = emin 222 if (emin < data.fEmin) data.fEmin = emin;
203 G4double emax = localData.fEmax; 223 G4double emax = localData.fEmax;
204 if (emax > data.fEmax) data.fEmax = emax 224 if (emax > data.fEmax) data.fEmax = emax;
205 data.fTmean = localData.fTmean; 225 data.fTmean = localData.fTmean;
206 } << 226 }
207 } 227 }
208 << 228
209 // activity << 229 //activity
210 fTimeWindow1 = localRun->fTimeWindow1; 230 fTimeWindow1 = localRun->fTimeWindow1;
211 fTimeWindow2 = localRun->fTimeWindow2; 231 fTimeWindow2 = localRun->fTimeWindow2;
212 << 232
213 std::map<G4String, ActivityData>::const_iter << 233 std::map<G4String,ActivityData>::const_iterator ita;
214 for (ita = localRun->fActivityMap.begin(); i << 234 for (ita = localRun->fActivityMap.begin();
>> 235 ita != localRun->fActivityMap.end(); ++ita) {
>> 236
215 G4String name = ita->first; 237 G4String name = ita->first;
216 const ActivityData& localData = ita->secon << 238 const ActivityData& localData = ita->second;
217 if (fActivityMap.find(name) == fActivityMa << 239 if ( fActivityMap.find(name) == fActivityMap.end()) {
218 fActivityMap[name] = << 240 fActivityMap[name]
219 ActivityData(localData.fNlife_t1, loca << 241 = ActivityData(localData.fNlife_t1,
220 } << 242 localData.fNlife_t2,
221 else { << 243 localData.fNdecay_t1t2);
222 ActivityData& data = fActivityMap[name]; << 244 } else {
>> 245 ActivityData& data = fActivityMap[name];
223 data.fNlife_t1 += localData.fNlife_t1; 246 data.fNlife_t1 += localData.fNlife_t1;
224 data.fNlife_t2 += localData.fNlife_t2; 247 data.fNlife_t2 += localData.fNlife_t2;
225 data.fNdecay_t1t2 += localData.fNdecay_t 248 data.fNdecay_t1t2 += localData.fNdecay_t1t2;
226 } 249 }
227 } 250 }
228 << 251
229 G4Run::Merge(run); << 252 G4Run::Merge(run);
230 } << 253 }
231 << 254
232 //....oooOO0OOooo........oooOO0OOooo........oo 255 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
233 256
234 void Run::EndOfRun() << 257 void Run::EndOfRun()
235 { 258 {
236 G4int nbEvents = numberOfEvent; << 259 G4int nbEvents = numberOfEvent;
237 G4String partName = fParticle->GetParticleNa << 260 G4String partName = fParticle->GetParticleName();
238 << 261
239 G4cout << "\n ======================== run s << 262 G4cout << "\n ======================== run summary ======================";
240 G4cout << "\n The run was " << nbEvents << " << 263 G4cout << "\n The run was " << nbEvents << " " << partName << " of "
241 << G4BestUnit(fEkin, "Energy"); << 264 << G4BestUnit(fEkin,"Energy");
242 G4cout << "\n ============================== << 265 G4cout << "\n ===========================================================\n";
243 G4cout << G4endl; << 266 G4cout << G4endl;
244 if (nbEvents == 0) { << 267 if (nbEvents == 0) { return; }
245 return; << 268
246 } << 269 G4int prec = 4, wid = prec + 2;
247 << 270 G4int dfprec = G4cout.precision(prec);
248 G4int prec = 4, wid = prec + 2; << 271
249 G4int dfprec = G4cout.precision(prec); << 272 //particle count
250 << 273 //
251 // particle count << 274 G4cout << " Nb of generated particles: \n" << G4endl;
252 // << 275
253 G4cout << " Nb of generated particles: \n" < << 276 std::map<G4String,ParticleData>::iterator it;
254 << 277 for (it = fParticleDataMap.begin(); it != fParticleDataMap.end(); it++) {
255 std::map<G4String, ParticleData>::iterator i << 278 G4String name = it->first;
256 for (it = fParticleDataMap.begin(); it != fP <<
257 G4String name = it->first; <<
258 ParticleData data = it->second; 279 ParticleData data = it->second;
259 G4int count = data.fCount; << 280 G4int count = data.fCount;
260 G4double eMean = data.fEmean / count; << 281 G4double eMean = data.fEmean/count;
261 G4double eMin = data.fEmin; << 282 G4double eMin = data.fEmin;
262 G4double eMax = data.fEmax; << 283 G4double eMax = data.fEmax;
263 G4double meanLife = data.fTmean; 284 G4double meanLife = data.fTmean;
264 << 285
265 G4cout << " " << std::setw(15) << name << 286 G4cout << " " << std::setw(15) << name << ": " << std::setw(7) << count
266 << " Emean = " << std::setw(wid) < << 287 << " Emean = " << std::setw(wid) << G4BestUnit(eMean, "Energy")
267 << G4BestUnit(eMin, "Energy") << " << 288 << "\t( " << G4BestUnit(eMin, "Energy")
>> 289 << " --> " << G4BestUnit(eMax, "Energy") << ")";
268 if (meanLife > 0.) 290 if (meanLife > 0.)
269 G4cout << "\tmean life = " << G4BestUnit << 291 G4cout << "\tmean life = " << G4BestUnit(meanLife, "Time") << G4endl;
270 else if (meanLife < 0.) << 292 else if (meanLife < 0.) G4cout << "\tstable" << G4endl;
271 G4cout << "\tstable" << G4endl; << 293 else G4cout << G4endl;
272 else << 294 }
273 G4cout << G4endl; << 295
274 } << 296 //energy momentum balance
275 << 297 //
276 // energy momentum balance << 298
277 // << 299 if (fDecayCount > 0) {
278 << 300 G4double Ebmean = fEkinTot[0]/fDecayCount;
279 if (fDecayCount > 0) { << 301 G4double Pbmean = fPbalance[0]/fDecayCount;
280 G4double Ebmean = fEkinTot[0] / fDecayCoun << 302
281 G4double Pbmean = fPbalance[0] / fDecayCou << 303 G4cout << "\n Ekin Total (Q single decay): mean = "
282 << 304 << std::setw(wid) << G4BestUnit(Ebmean, "Energy")
283 G4cout << "\n Ekin Total (Q single decay << 305 << "\t( " << G4BestUnit(fEkinTot[1], "Energy")
284 << G4BestUnit(Ebmean, "Energy") << << 306 << " --> " << G4BestUnit(fEkinTot[2], "Energy")
285 << G4BestUnit(fEkinTot[2], "Energy" << 307 << ")" << G4endl;
286 << 308
287 G4cout << "\n Momentum balance (excludin << 309 G4cout << "\n Momentum balance (excluding gamma desexcitation): mean = "
288 << G4BestUnit(Pbmean, "Energy") << << 310 << std::setw(wid) << G4BestUnit(Pbmean, "Energy")
289 << " --> " << G4BestUnit(fPbalance[ << 311 << "\t( " << G4BestUnit(fPbalance[1], "Energy")
290 } << 312 << " --> " << G4BestUnit(fPbalance[2], "Energy")
291 << 313 << ")" << G4endl;
292 // total time of life << 314 }
293 // << 315
294 if (fTimeCount > 0) { << 316 //total time of life
295 G4double Tmean = fEventTime[0] / fTimeCoun << 317 //
296 G4double halfLife = Tmean * std::log(2.); << 318 if (fTimeCount > 0) {
297 << 319 G4double Tmean = fEventTime[0]/fTimeCount;
298 G4cout << "\n Total time of life (full c << 320 G4double halfLife = Tmean*std::log(2.);
299 << G4BestUnit(Tmean, "Time") << " << 321
300 << G4BestUnit(halfLife, "Time") << << 322 G4cout << "\n Total time of life (full chain): mean = "
301 << " --> " << G4BestUnit(fEventTime << 323 << std::setw(wid) << G4BestUnit(Tmean, "Time")
302 } << 324 << " half-life = "
303 << 325 << std::setw(wid) << G4BestUnit(halfLife, "Time")
304 // total visible Energy << 326 << " ( " << G4BestUnit(fEventTime[1], "Time")
305 // << 327 << " --> " << G4BestUnit(fEventTime[2], "Time")
306 if (fTimeCount > 0) { << 328 << ")" << G4endl;
307 G4double Evmean = fEvisEvent[0] / fTimeCou << 329 }
308 << 330
309 G4cout << "\n Total visible energy (full << 331 //total visible Energy
310 << G4BestUnit(Evmean, "Energy") << << 332 //
311 << " --> " << G4BestUnit(fEvisEvent << 333 if (fTimeCount > 0) {
312 } << 334 G4double Evmean = fEvisEvent[0]/fTimeCount;
313 << 335
314 // activity of primary ion << 336 G4cout << "\n Total visible energy (full chain) : mean = "
315 // << 337 << std::setw(wid) << G4BestUnit(Evmean, "Energy")
316 G4double pTimeMean = fPrimaryTime / nbEvents << 338 << " ( " << G4BestUnit(fEvisEvent[1], "Energy")
317 G4double molMass = fParticle->GetAtomicMass( << 339 << " --> " << G4BestUnit(fEvisEvent[2], "Energy")
318 G4double nAtoms_perUnitOfMass = Avogadro / m << 340 << ")" << G4endl;
319 G4double Activity_perUnitOfMass = 0.0; << 341 }
320 if (pTimeMean > 0.0) { << 342
321 Activity_perUnitOfMass = nAtoms_perUnitOfM << 343 //activity of primary ion
322 } << 344 //
323 << 345 G4double pTimeMean = fPrimaryTime/nbEvents;
324 G4cout << "\n Activity of " << partName << << 346 G4double molMass = fParticle->GetAtomicMass()*g/mole;
325 << Activity_perUnitOfMass * g / becqu << 347 G4double nAtoms_perUnitOfMass = Avogadro/molMass;
326 << Activity_perUnitOfMass * g / curie << 348 G4double Activity_perUnitOfMass = 0.0;
327 << G4endl; << 349 if (pTimeMean > 0.0)
328 << 350 { Activity_perUnitOfMass = nAtoms_perUnitOfMass/pTimeMean; }
329 // activities in time window << 351
330 // << 352 G4cout << "\n Activity of " << partName << " = "
331 if (fTimeWindow2 > 0.) { << 353 << std::setw(wid) << Activity_perUnitOfMass*g/becquerel
332 G4double dt = fTimeWindow2 - fTimeWindow1; << 354 << " Bq/g (" << Activity_perUnitOfMass*g/curie
333 G4cout << " Activities in time window [t << 355 << " Ci/g) \n"
334 << ", " << G4BestUnit(fTimeWindow2, << 356 << G4endl;
335 << "] (delta time = " << G4BestUni << 357
336 << G4endl; << 358
337 << 359 //activities in time window
338 std::map<G4String, ActivityData>::iterator << 360 //
339 for (ita = fActivityMap.begin(); ita != fA << 361 if (fTimeWindow2 > 0.) {
340 G4String name = ita->first; << 362 G4double dt = fTimeWindow2 - fTimeWindow1;
>> 363 G4cout << " Activities in time window [t1, t2] = ["
>> 364 << G4BestUnit(fTimeWindow1, "Time") << ", "
>> 365 << G4BestUnit(fTimeWindow2, "Time") << "] (delta time = "
>> 366 << G4BestUnit(dt, "Time") << ") : \n" << G4endl;
>> 367
>> 368 std::map<G4String,ActivityData>::iterator ita;
>> 369 for (ita = fActivityMap.begin(); ita != fActivityMap.end(); ita++) {
>> 370 G4String name = ita->first;
341 ActivityData data = ita->second; 371 ActivityData data = ita->second;
342 G4int n1 = data.fNlife_t1; << 372 G4int n1 = data.fNlife_t1;
343 G4int n2 = data.fNlife_t2; << 373 G4int n2 = data.fNlife_t2;
344 G4int ndecay = data.fNdecay_t1t2; 374 G4int ndecay = data.fNdecay_t1t2;
345 G4double actv = ndecay / dt; << 375 G4double actv = ndecay/dt;
346 376
347 G4cout << " " << std::setw(15) << name 377 G4cout << " " << std::setw(15) << name << ": "
348 << " n(t1) = " << std::setw(7) < << 378 << " n(t1) = " << std::setw(7) << n1
349 << "\t decays = " << std::setw( << 379 << "\tn(t2) = " << std::setw(7) << n2
350 << " ---> <actv> = " << G4BestU << 380 << "\t decays = " << std::setw(7) << ndecay
351 } << 381 << " ---> <actv> = " << G4BestUnit(actv, "Activity") << "\n";
352 } << 382 }
353 G4cout << G4endl; << 383 }
354 << 384 G4cout << G4endl;
355 // normalize histograms << 385
356 // << 386 //normalize histograms
357 G4AnalysisManager* analysisManager = G4Analy << 387 //
358 G4double factor = 100. / nbEvents; << 388 G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
359 analysisManager->ScaleH1(1, factor); << 389 G4double factor = 100./nbEvents;
360 analysisManager->ScaleH1(2, factor); << 390 analysisManager->ScaleH1(1,factor);
361 analysisManager->ScaleH1(3, factor); << 391 analysisManager->ScaleH1(2,factor);
362 analysisManager->ScaleH1(4, factor); << 392 analysisManager->ScaleH1(3,factor);
363 analysisManager->ScaleH1(5, factor); << 393 analysisManager->ScaleH1(4,factor);
364 << 394 analysisManager->ScaleH1(5,factor);
365 // remove all contents in fParticleDataMap << 395
366 // << 396 // remove all contents in fParticleDataMap
367 fParticleDataMap.clear(); << 397 //
368 fActivityMap.clear(); << 398 fParticleDataMap.clear();
369 << 399 fActivityMap.clear();
370 // restore default precision << 400
371 // << 401 // restore default precision
372 G4cout.precision(dfprec); << 402 //
>> 403 G4cout.precision(dfprec);
373 } 404 }
374 405
375 //....oooOO0OOooo........oooOO0OOooo........oo 406 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
376 407