Hadr04/src/Run.cc

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 /// \file Run.cc
0027 /// \brief Implementation of the Run class
0028
0029 #include "Run.hh"
0030
0031 #include "DetectorConstruction.hh"
0032 #include "HistoManager.hh"
0033 #include "PrimaryGeneratorAction.hh"
0034
0035 #include "G4SystemOfUnits.hh"
0036 #include "G4UnitsTable.hh"
0037
0038 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0039
0040 Run::Run(DetectorConstruction* det) : fDetector(det) {}
0041
0042 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0043
0044 void Run::SetPrimary(G4ParticleDefinition* particle, G4double energy)
0045 {
0046   fParticle = particle;
0047   fEkin = energy;
0048 }
0049
0050 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0051
0052 void Run::CountProcesses(const G4VProcess* process)
0053 {
0054   if (process == nullptr) return;
0055   G4String procName = process->GetProcessName();
0056   std::map<G4String, G4int>::iterator it = fProcCounter.find(procName);
0057   if (it == fProcCounter.end()) {
0058     fProcCounter[procName] = 1;
0059   }
0060   else {
0061     fProcCounter[procName]++;
0062   }
0063 }
0064
0065 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0066
0067 void Run::ParticleCount(G4String name, G4double Ekin)
0068 {
0069   std::map<G4String, ParticleData>::iterator it = fParticleDataMap.find(name);
0070   if (it == fParticleDataMap.end()) {
0071     fParticleDataMap[name] = ParticleData(1, Ekin, Ekin, Ekin);
0072   }
0073   else {
0074     ParticleData& data = it->second;
0075     data.fCount++;
0076     data.fEmean += Ekin;
0077     // update min max
0078     G4double emin = data.fEmin;
0079     if (Ekin < emin) data.fEmin = Ekin;
0080     G4double emax = data.fEmax;
0081     if (Ekin > emax) data.fEmax = Ekin;
0082   }
0083 }
0084
0085 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0086
0087 void Run::SumTrackLength(G4int nstep1, G4int nstep2, G4double trackl1, G4double trackl2,
0088                          G4double time1, G4double time2)
0089 {
0090   fNbStep1 += nstep1;
0091   fNbStep2 += nstep2;
0092   fTrackLen1 += trackl1;
0093   fTrackLen2 += trackl2;
0094   fTime1 += time1;
0095   fTime2 += time2;
0096 }
0097
0098 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0099
0100 void Run::Merge(const G4Run* run)
0101 {
0102   const Run* localRun = static_cast<const Run*>(run);
0103
0104   // primary particle info
0105   //
0106   fParticle = localRun->fParticle;
0107   fEkin = localRun->fEkin;
0108
0109   // accumulate sums
0110   //
0111   fNbStep1 += localRun->fNbStep1;
0112   fNbStep2 += localRun->fNbStep2;
0113   fTrackLen1 += localRun->fTrackLen1;
0114   fTrackLen2 += localRun->fTrackLen2;
0115   fTime1 += localRun->fTime1;
0116   fTime2 += localRun->fTime2;
0117
0118   // map: processes count
0119   std::map<G4String, G4int>::const_iterator itp;
0120   for (itp = localRun->fProcCounter.begin(); itp != localRun->fProcCounter.end(); ++itp) {
0121     G4String procName = itp->first;
0122     G4int localCount = itp->second;
0123     if (fProcCounter.find(procName) == fProcCounter.end()) {
0124       fProcCounter[procName] = localCount;
0125     }
0126     else {
0127       fProcCounter[procName] += localCount;
0128     }
0129   }
0130
0131   // map: created particles count
0132   std::map<G4String, ParticleData>::const_iterator itn;
0133   for (itn = localRun->fParticleDataMap.begin(); itn != localRun->fParticleDataMap.end(); ++itn) {
0134     G4String name = itn->first;
0135     const ParticleData& localData = itn->second;
0136     if (fParticleDataMap.find(name) == fParticleDataMap.end()) {
0137       fParticleDataMap[name] =
0138         ParticleData(localData.fCount, localData.fEmean, localData.fEmin, localData.fEmax);
0139     }
0140     else {
0141       ParticleData& data = fParticleDataMap[name];
0142       data.fCount += localData.fCount;
0143       data.fEmean += localData.fEmean;
0144       G4double emin = localData.fEmin;
0145       if (emin < data.fEmin) data.fEmin = emin;
0146       G4double emax = localData.fEmax;
0147       if (emax > data.fEmax) data.fEmax = emax;
0148     }
0149   }
0150
0151   G4Run::Merge(run);
0152 }
0153
0154 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0155
0156 void Run::EndOfRun()
0157 {
0158   G4int prec = 5, wid = prec + 2;
0159   G4int dfprec = G4cout.precision(prec);
0160
0161   // run condition
0162   //
0163   G4Material* material = fDetector->GetMaterial();
0164   G4double density = material->GetDensity();
0165
0166   G4String Particle = fParticle->GetParticleName();
0167   G4cout << "\n The run is " << numberOfEvent << " " << Particle << " of "
0168          << G4BestUnit(fEkin, "Energy") << " through "
0169          << G4BestUnit(0.5 * (fDetector->GetSize()), "Length") << " of " << material->GetName()
0170          << " (density: " << G4BestUnit(density, "Volumic Mass") << ")" << G4endl;
0171
0172   if (numberOfEvent == 0) {
0173     G4cout.precision(dfprec);
0174     return;
0175   }
0176
0177   // frequency of processes
0178   //
0179   G4cout << "\n Process calls frequency :" << G4endl;
0180   G4int survive = 0;
0181   std::map<G4String, G4int>::iterator it;
0182   for (it = fProcCounter.begin(); it != fProcCounter.end(); it++) {
0183     G4String procName = it->first;
0184     G4int count = it->second;
0185     G4cout << "\t" << procName << "= " << count;
0186     if (procName == "Transportation") survive = count;
0187   }
0188   G4cout << G4endl;
0189
0190   if (survive > 0) {
0191     G4cout << "\n Nb of incident particles surviving after "
0192            << G4BestUnit(0.5 * (fDetector->GetSize()), "Length") << " of "
0193            << fDetector->GetMaterial()->GetName() << " : " << survive << G4endl;
0194   }
0195
0196   // total track length of incident neutron
0197   //
0198   G4cout << "\n Parcours of incident neutron:";
0199
0200   G4double meanCollision1 = (G4double)fNbStep1 / numberOfEvent;
0201   G4double meanCollision2 = (G4double)fNbStep2 / numberOfEvent;
0202   G4double meanCollisTota = meanCollision1 + meanCollision2;
0203
0204   G4cout << "\n   nb of collisions    E>1*eV= " << meanCollision1
0205          << "      E<1*eV= " << meanCollision2 << "       total= " << meanCollisTota;
0206
0207   G4double meanTrackLen1 = fTrackLen1 / numberOfEvent;
0208   G4double meanTrackLen2 = fTrackLen2 / numberOfEvent;
0209   G4double meanTrackLtot = meanTrackLen1 + meanTrackLen2;
0210
0211   G4cout << "\n   track length        E>1*eV= " << G4BestUnit(meanTrackLen1, "Length")
0212          << "  E<1*eV= " << G4BestUnit(meanTrackLen2, "Length")
0213          << "   total= " << G4BestUnit(meanTrackLtot, "Length");
0214
0215   G4double meanTime1 = fTime1 / numberOfEvent;
0216   G4double meanTime2 = fTime2 / numberOfEvent;
0217   G4double meanTimeTo = meanTime1 + meanTime2;
0218
0219   G4cout << "\n   time of flight      E>1*eV= " << G4BestUnit(meanTime1, "Time")
0220          << "  E<1*eV= " << G4BestUnit(meanTime2, "Time")
0221          << "   total= " << G4BestUnit(meanTimeTo, "Time") << G4endl;
0222
0223   // particles count
0224   //
0225   G4cout << "\n List of generated particles:" << G4endl;
0226
0227   std::map<G4String, ParticleData>::iterator itn;
0228   for (itn = fParticleDataMap.begin(); itn != fParticleDataMap.end(); itn++) {
0229     G4String name = itn->first;
0230     ParticleData data = itn->second;
0231     G4int count = data.fCount;
0232     G4double eMean = data.fEmean / count;
0233     G4double eMin = data.fEmin;
0234     G4double eMax = data.fEmax;
0235
0236     G4cout << "  " << std::setw(13) << name << ": " << std::setw(7) << count
0237            << "  Emean = " << std::setw(wid) << G4BestUnit(eMean, "Energy") << "\t( "
0238            << G4BestUnit(eMin, "Energy") << " --> " << G4BestUnit(eMax, "Energy") << ")" << G4endl;
0239   }
0240
0241   // normalize histograms
0242   ////G4AnalysisManager* analysisManager = G4AnalysisManager::Instance();
0243   ////G4double factor = 1./numberOfEvent;
0244   ////analysisManager->ScaleH1(3,factor);
0245
0246   // remove all contents in fProcCounter, fCount
0247   fProcCounter.clear();
0248   fParticleDataMap.clear();
0249
0250   // restore default format
0251   G4cout.precision(dfprec);
0252 }
0253
0254 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......