source: trunk/source/processes/electromagnetic/lowenergy/test/hTest/src/hTestHisto.cc @ 1350

Last change on this file since 1350 was 1350, checked in by garnier, 13 years ago

update to last version 4.9.4
File size: 8.0 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//---------------------------------------------------------------------------
27//
28// ClassName:   hTestHisto
29// 
30//
31// Author:      V.Ivanchenko 30/01/01
32//
33//----------------------------------------------------------------------------
34//
35
36//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
37//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
38
39#include "hTestHisto.hh"
40#include <iomanip>
41
42#include <memory> // for the auto_ptr(T>
43
44#include "AIDA/AIDA.h"
45
46//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47
48hTestHisto* hTestHisto::fManager = 0;
49
50//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
51
52hTestHisto* hTestHisto::GetPointer()
53{
54  if(!fManager) {
55    fManager = new hTestHisto();
56  }
57  return fManager;
58}
59
60//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
61
62hTestHisto::hTestHisto() 
63{
64  verbose = 1;
65  histName = G4String("histo.hbook");
66  ntup = 0;
67  nHisto = 1;
68  nAbsSaved = 0;
69  maxEnergy = 0.0;
70  nTuple = false;
71}
72
73//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
74
75hTestHisto::~hTestHisto() 
76{
77  histo.clear(); 
78  G4cout << "hTestHisto: Histograms are deleted for " << theName << G4endl;
79}
80
81//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
82
83void hTestHisto::BeginOfHisto(G4int num)
84{ 
85  if(0 < verbose) G4cout << "hTestHisto # " << num << " started " << G4endl;
86  zend     = 0.0;
87  zend2    = 0.0;
88  zEvt     = 0.0;
89 
90  if(0 < nHisto) bookHisto();
91
92  if(verbose > 0) {
93    G4cout << "hTestHisto: Histograms are booked and run has been started" 
94           << G4endl;
95  }
96}
97
98//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
99
100void hTestHisto::EndOfHisto()
101{
102
103  G4cout << "hTestHisto: End of run actions are started" << G4endl;
104
105  // Zend average
106  if(zEvt > 0.0) {
107    zend  /= zEvt;
108    zend2 /= zEvt;
109    zend2 -= zend*zend;
110    G4double sig = 0.0;
111    if(zend2 > 0.) sig = std::sqrt(zend2);
112    zend2 = sig / std::sqrt(zEvt);
113    G4cout<<"========================================================"<<G4endl;
114    G4cout << setprecision(4) << "Range(mm)= " << zend/mm
115           << "; Stragling(mm)= " << sig/mm
116           << setprecision(2) << " +- " << zend2/mm
117           << "    " << zEvt << " events for range" << G4endl;
118    G4cout<<"========================================================"<<G4endl;
119  } 
120
121   // Write histogram file
122  if(0 < nHisto || ntup) {
123    tree->commit();
124    std::cout << "Closing the tree..." << std::endl;
125    tree->close();
126    G4cout << "Histograms and Ntuples are saved" << G4endl;
127  }
128}
129
130//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
131
132void hTestHisto::SaveEvent()
133{
134  if(ntup) {
135    ntup->addRow();
136  }                       
137}
138
139//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
140
141void hTestHisto::SaveToTuple(const G4String& parname, G4double val)
142{
143  if(ntup) ntup->fill( ntup->findColumn(parname), (float)val);
144}
145
146//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
147
148void hTestHisto::SaveToTuple(const G4String& parname,G4double val,G4double defval)
149{
150  if(ntup) ntup->fill( ntup->findColumn(parname), (float)val);
151}
152
153//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
154
155void hTestHisto::bookHisto()
156{
157  G4double zmax = (AbsorberThickness + gap) * NumberOfAbsorbers / mm;
158  G4cout << "hTestHisto: Histograms will be saved to the file <" 
159         << histName << ">"
160         << " AbsThick(mm)= " << AbsorberThickness/mm
161         << " Nabs= " << NumberOfAbsorbers
162         << " zmax= " << zmax
163         << " emax= " << maxEnergy
164         << " nHisto= " << nHisto
165         << G4endl;
166
167  // Creating the analysis factory
168  std::auto_ptr< AIDA::IAnalysisFactory > af( AIDA_createAnalysisFactory() );
169
170  // Creating the tree factory
171  std::auto_ptr< AIDA::ITreeFactory > tf( af->createTreeFactory() );
172
173  // Creating a tree mapped to a new hbook file.
174  tree = tf->create(histName,false,false,"hbook");
175  G4cout << "Tree store : " << tree->storeName() << G4endl;
176 
177  histo.resize(nHisto);
178
179  // Creating a histogram factory, whose histograms will be handled by the tree
180  std::auto_ptr< AIDA::IHistogramFactory > hf(af->createHistogramFactory( *tree ));
181
182  // Creating an 1-dimensional histograms in the root directory of the tree
183
184  if(0 < nHisto) histo[0] = hf->createHistogram1D("10",
185    "Energy deposit (MeV) in absorber (mm)",NumberOfAbsorbers,0.0,zmax);
186
187  if(1 < nHisto) histo[1] = hf->createHistogram1D("11",
188    "Energy (MeV) of delta-electrons",50,0.0,maxEnergy/MeV);
189
190  if(2 < nHisto) histo[2] = hf->createHistogram1D("12",
191    "Theta (degrees) of delta-electrons",36,0.0,180.);
192
193  if(3 < nHisto) histo[3] = hf->createHistogram1D("13",
194    "Energy (MeV) of secondary gamma",50,0.0,maxEnergy/MeV);
195
196  if(4 < nHisto) histo[4] = hf->createHistogram1D("14",
197    "Theta (degrees) of secondary gamma",36,0.0,180.);
198
199  // Creating a tuple factory, whose tuples will be handled by the tree
200  std::auto_ptr< AIDA::ITupleFactory > tpf( af->createTupleFactory( *tree ) );
201
202  // If using Anaphe HBOOK implementation, there is a limitation on the
203  // length of the variable names in a ntuple
204  if(nTuple) ntup = tpf->createITuple( "100", "Range/Energy", 
205  "float xend, yend, zend, ltpk, tend, teta, loss, dedx, back, leak, edep" );
206
207}
208
209//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
210
211void hTestHisto::AddEnergy(G4double edep, G4double z)
212{
213  if(0 < nHisto) histo[0]->fill((float)z/mm, (float)edep/MeV);
214}
215
216//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
217
218void hTestHisto::AddEndPoint(G4double z)
219{
220  zend  += z;
221  zend2 += z*z;
222  zEvt  += 1.0;
223}
224
225//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
226
227void hTestHisto::AddDeltaElectron(const G4DynamicParticle* elec)
228{
229  if(1 < nHisto) histo[1]->fill(elec->GetKineticEnergy()/MeV,1.0);
230  if(2 < nHisto)
231     histo[2]->fill((elec->GetMomentumDirection()).theta()/deg,1.0);
232}
233
234//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
235
236void hTestHisto::AddPhoton(const G4DynamicParticle* ph)
237{
238  if(3 < nHisto) histo[3]->fill(ph->GetKineticEnergy()/MeV,1.0);
239  if(4 < nHisto)
240     histo[4]->fill((ph->GetMomentumDirection()).theta()/deg,1.0);
241}
242
243//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
244
245
Note: See TracBrowser for help on using the repository browser.