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