File indexing completed on 2025-04-10 08:06:48
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029
0030
0031
0032
0033
0034
0035
0036
0037
0038
0039 #include "Run.hh"
0040
0041 #include "DetectorConstruction.hh"
0042 #include "PrimaryGeneratorAction.hh"
0043
0044 #include "G4Material.hh"
0045 #include "G4SystemOfUnits.hh"
0046 #include "G4UnitsTable.hh"
0047
0048
0049
0050 Run::Run(const DetectorConstruction* detector)
0051 : G4Run(), fDetector(detector), fParticle(0), fEkin(0.), fSP(0.), fSP2(0.)
0052 {}
0053
0054
0055
0056 Run::~Run() {}
0057
0058
0059
0060 void Run::SetPrimary(G4ParticleDefinition* particle, G4double energy)
0061 {
0062 fParticle = particle;
0063 fEkin = energy;
0064 }
0065
0066
0067
0068 void Run::AddSP(G4double t)
0069 {
0070 fSP += t;
0071 fSP2 += t * t;
0072 }
0073
0074
0075
0076 void Run::Merge(const G4Run* run)
0077 {
0078 const Run* localRun = static_cast<const Run*>(run);
0079
0080
0081 fParticle = localRun->fParticle;
0082 fEkin = localRun->fEkin;
0083
0084
0085 fSP += localRun->fSP;
0086 fSP2 += localRun->fSP2;
0087
0088 G4Run::Merge(run);
0089 }
0090
0091
0092
0093 void Run::EndOfRun()
0094 {
0095 std::ios::fmtflags mode = G4cout.flags();
0096 G4cout.setf(std::ios::fixed, std::ios::floatfield);
0097 G4int prec = G4cout.precision(2);
0098
0099
0100 G4Material* material = fDetector->GetAbsorMaterial();
0101 G4double density = material->GetDensity();
0102 G4String partName = fParticle->GetParticleName();
0103
0104 G4cout << "\n ======================== run summary =====================\n";
0105 G4cout << "\n The run is " << numberOfEvent << " " << partName << " of "
0106 << G4BestUnit(fEkin, "Energy") << " through a sphere of radius "
0107 << G4BestUnit(fDetector->GetAbsorRadius(), "Length") << "of " << material->GetName()
0108 << " (density: " << G4BestUnit(density, "Volumic Mass") << ")" << G4endl;
0109
0110 if (numberOfEvent == 0) {
0111 G4cout.setf(mode, std::ios::floatfield);
0112 G4cout.precision(prec);
0113 return;
0114 }
0115
0116
0117 fSP /= numberOfEvent;
0118 fSP2 /= numberOfEvent;
0119 G4double rms = fSP2 - fSP * fSP;
0120
0121 if (rms > 0.)
0122 rms = std::sqrt(rms);
0123 else
0124 rms = 0.;
0125
0126 G4cout << "\n total Stopping Power (keV/um) = " << fSP / (keV / um) << " +- "
0127 << rms / (keV / um) << G4endl;
0128
0129
0130 FILE* myFile;
0131 myFile = fopen("spower.txt", "a");
0132 fprintf(myFile, "%e %e %e \n", fEkin / eV, fSP / (keV / um), rms / (keV / um));
0133 fclose(myFile);
0134
0135
0136 G4cout.setf(mode, std::ios::floatfield);
0137 G4cout.precision(prec);
0138 }