source: trunk/examples/extended/electromagnetic/TestEm7/src/RunAction.cc@ 1353

Last change on this file since 1353 was 1342, checked in by garnier, 15 years ago

update ti head
File size: 6.5 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// $Id: RunAction.cc,v 1.25 2010/09/17 18:45:43 maire Exp $
27// GEANT4 tag $Name: examples-V09-03-09 $
28//
29//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
30//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
31
32#include "RunAction.hh"
33#include "DetectorConstruction.hh"
34#include "PhysicsList.hh"
35#include "StepMax.hh"
36#include "HistoManager.hh"
37#include "PrimaryGeneratorAction.hh"
38
39#include "G4Run.hh"
40#include "G4RunManager.hh"
41#include "G4UnitsTable.hh"
42#include "G4ios.hh"
43
44#include "Randomize.hh"
45
46//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
47
48RunAction::RunAction(DetectorConstruction* det, PhysicsList* phys,
49 HistoManager* histo, PrimaryGeneratorAction* kin)
50:detector(det), physics(phys),histoManager(histo), kinematic(kin)
51{
52 tallyEdep = new G4double[MaxTally];
53}
54
55//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
56
57RunAction::~RunAction()
58{
59 delete [] tallyEdep;
60}
61
62//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
63
64void RunAction::BeginOfRunAction(const G4Run* aRun)
65{
66 G4cout << "### Run " << aRun->GetRunID() << " start." << G4endl;
67
68 // save Rndm status
69 G4RunManager::GetRunManager()->SetRandomNumberStore(true);
70 CLHEP::HepRandom::showEngineStatus();
71
72 //initialize projected range, tallies, Ebeam, and book histograms
73 //
74 nPrimarySteps = 0;
75 nRange = 0;
76 projRange = projRange2 = 0.;
77 edeptot = eniel = 0.;
78 for (G4int j=0; j<MaxTally; j++) tallyEdep[j] = 0.;
79 kinematic->ResetEbeamCumul();
80
81 // define "1" histogram binning
82 // histogram "1" is defined by the length of the target
83 // zoomed histograms are defined by UI command
84 G4double length = detector->GetAbsorSizeX();
85 G4double stepMax = physics->GetStepMaxProcess()->GetMaxStep();
86 G4int nbmin = 100;
87 G4int nbBins = (G4int)(0.5 + length/stepMax);
88 if (nbBins < nbmin) nbBins = nbmin;
89 histoManager->SetHisto(1, nbBins, 0., length, "mm");
90
91 histoManager->book();
92}
93
94//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
95
96void RunAction::EndOfRunAction(const G4Run* aRun)
97{
98 G4int NbofEvents = aRun->GetNumberOfEvent();
99 if (NbofEvents == 0) return;
100
101 //run conditions
102 //
103 G4Material* material = detector->GetAbsorMaterial();
104 G4double density = material->GetDensity();
105
106 G4String particle = kinematic->GetParticleGun()->GetParticleDefinition()
107 ->GetParticleName();
108 G4double energy = kinematic->GetParticleGun()->GetParticleEnergy();
109 G4cout << "\n The run consists of " << NbofEvents << " "<< particle << " of "
110 << G4BestUnit(energy,"Energy") << " through "
111 << G4BestUnit(detector->GetAbsorSizeX(),"Length") << " of "
112 << material->GetName() << " (density: "
113 << G4BestUnit(density,"Volumic Mass") << ")" << G4endl;
114
115 //compute projected range and straggling
116 //
117 if(nRange > 0) {
118 projRange /= nRange;
119 projRange2 /= nRange;
120 }
121 G4double rms = projRange2 - projRange*projRange;
122 if (rms>0.) rms = std::sqrt(rms); else rms = 0.;
123
124 G4double nstep = G4double(nPrimarySteps)/G4double(NbofEvents);
125
126 G4cout.precision(6);
127 G4cout << "\n Projected Range= "<< G4BestUnit(projRange,"Length")
128 << " rms= " << G4BestUnit( rms,"Length")
129 << G4endl;
130 G4cout << " Mean number of primary steps = "<< nstep << G4endl;
131
132 //compute energy deposition and niel
133 //
134 edeptot /= NbofEvents;
135 G4cout << " Total energy deposit= "<< G4BestUnit(edeptot,"Energy")
136 << G4endl;
137 eniel /= NbofEvents;
138 G4cout << " niel energy deposit = "<< G4BestUnit(eniel,"Energy")
139 << G4endl;
140
141 //print dose in tallies
142 //
143 G4int tallyNumber = detector->GetTallyNumber();
144 if (tallyNumber > 0) {
145 G4double tallyMass = detector->GetTallyMass();
146 G4double Ebeam = kinematic->GetEbeamCumul();
147 G4cout << "\n---------------------------------------------------------\n";
148 G4cout << " Cumulated Doses : \tEdep \tEdep/Ebeam \tDose" << G4endl;
149 for (G4int j=0; j<tallyNumber; j++) {
150 G4double Edep = tallyEdep[j], ratio = 100*Edep/Ebeam;
151 G4double Dose = Edep/tallyMass;
152 G4cout << " tally " << j << ": \t \t"
153 << G4BestUnit(Edep,"Energy") << "\t"
154 << ratio << " % \t"
155 << G4BestUnit(Dose,"Dose") << G4endl;
156 }
157 G4cout << "\n---------------------------------------------------------\n";
158 G4cout << G4endl;
159 }
160
161 // normalize histograms
162 for (G4int j=1; j<3; j++) {
163 G4double binWidth = histoManager->GetBinWidth(j);
164 G4double fac = (mm/MeV)/(NbofEvents * binWidth);
165 histoManager->Scale(j, fac);
166 }
167 histoManager->Scale(3, 1./NbofEvents);
168
169 // save and clean histo
170 histoManager->save();
171
172 // show Rndm status
173 CLHEP::HepRandom::showEngineStatus();
174}
175
176//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
Note: See TracBrowser for help on using the repository browser.