source: trunk/examples/extended/hadronic/Hadr01/src/Histo.cc @ 1337

Last change on this file since 1337 was 1230, checked in by garnier, 14 years ago

update to geant4.9.3
File size: 9.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:   Histo - Generic histogram/ntuple manager class
29//
30//
31// Author:      V.Ivanchenko 30.10.03
32//
33//----------------------------------------------------------------------------
34//
35
36//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
37//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
38
39#include "Histo.hh"
40
41#ifdef G4ANALYSIS_USE
42#include <AIDA/AIDA.h>
43#include "HistoMessenger.hh"
44#endif
45
46//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
47//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
48
49Histo::Histo(G4int ver)
50{
51  verbose    = ver;
52  histName   = "histo";
53  histType   = "root";
54  option     = "--noErrors uncompress";
55  nHisto     = 0;
56  defaultAct = 1;
57  tupleName  = "tree.root";
58  tupleId    = "100";
59  tupleList  = "";
60  ntup       = 0;
61  messenger  = 0;
62
63#ifdef G4ANALYSIS_USE
64  messenger = new HistoMessenger(this);
65  tree = 0;
66  af   = 0;
67#endif
68}
69
70//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
71
72Histo::~Histo()
73{
74#ifdef G4ANALYSIS_USE
75  delete messenger;
76  delete af;
77#endif
78}
79
80//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
81
82void Histo::book()
83{
84#ifdef G4ANALYSIS_USE
85  G4cout << "### Histo books " << nHisto << " histograms " << G4endl;
86  // Creating the analysis factory
87  if(!af) af = AIDA_createAnalysisFactory();
88  if(verbose>0)
89    G4cout<<"HIsto books analysis factory ......... "<<G4endl;
90
91  // Creating the tree factory
92  AIDA::ITreeFactory* tf = af->createTreeFactory();
93  if(verbose>0)
94    G4cout<<"Histo books tree factory ......... "<<G4endl;
95
96  G4String histExt = "";
97  char* path = getenv("PHYSLIST");
98  if (path) histExt = "_" + G4String(path);
99
100  G4String histDir = "";
101  path = getenv("HISTODIR");
102  if (path) histDir = G4String(path) + "/";
103
104  // Creating a tree mapped to a new hbook file.
105  G4String name = histDir + histName + histExt + "." + histType;
106  tree = tf->create(name, histType, false, true, option);
107  G4cout << "Histo: tree store : " << tree->storeName() << G4endl;
108  delete tf;
109
110  // Creating a histogram factory, whose histograms will be handled by the tree
111  AIDA::IHistogramFactory* hf = af->createHistogramFactory( *tree );
112
113  // Creating an 1-dimensional histograms in the root directory of the tree
114  for(G4int i=0; i<nHisto; i++) {
115    if(active[i]) {
116      if(verbose>0)
117        G4cout<<" I am in book: histogram "<< i << " id= " << ids[i] <<G4endl;
118
119      G4String idd;
120      if(histType == "root") idd = "h" +  ids[i];
121      else                   idd = ids[i];
122      histo[i] = hf->createHistogram1D(idd, tittles[i], bins[i], xmin[i], xmax[i]);
123    } else {
124      histo[i] = 0;
125    }
126  }
127  delete hf;
128  // Creating a tuple factory, whose tuples will be handled by the tree
129  if(tupleList != "") {
130    if(verbose>0)
131      G4cout<<"Histo books tuple factory for "<<tupleName <<G4endl;
132
133    AIDA::ITupleFactory* tpf = af->createTupleFactory( *tree );
134    ntup = tpf->create(tupleId, tupleName, tupleList);
135    delete tpf;
136  }
137#endif
138} 
139
140//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
141
142void Histo::save()
143{
144#ifdef G4ANALYSIS_USE
145  // Write histogram file
146  tree->commit();
147  G4cout << "Closing the tree..." << G4endl;
148  tree->close();
149  G4cout << "Histograms and Ntuples are saved" << G4endl;
150  delete tree;
151  tree = 0;
152#endif
153}
154
155//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
156
157void Histo::reset()
158{
159#ifdef G4ANALYSIS_USE
160  delete tree;
161  tree = 0;
162#endif
163}
164
165//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
166
167void Histo::setFileType(const G4String& nam) 
168{
169  if(nam == "hbook" || nam == "root" || nam == "aida") histType = nam;
170  else if(nam == "XML" || nam == "xml") histType = "aida";
171}
172
173//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
174
175void Histo::add1D(const G4String& id, const G4String& name, G4int nb,
176                  G4double x1, G4double x2, G4double u)
177{
178  if(nHisto > 0) {
179    for(G4int i=0; i<nHisto; i++) {
180      if(ids[i] == id) return;
181    }
182  }
183
184  if(verbose > 0) 
185    G4cout << "New histogram will be booked: #" << id << "  <" << name
186           << "  " << nb << "  " << x1 << "  " << x2 << "  " << u
187           << G4endl;
188
189  nHisto++;
190  x1 /= u;
191  x2 /= u;
192  active.push_back(defaultAct);
193  bins.push_back(nb);
194  xmin.push_back(x1);
195  xmax.push_back(x2);
196  unit.push_back(u);
197  ids.push_back(id);
198  tittles.push_back(name);
199  histo.push_back(0);
200}
201
202//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
203
204void Histo::setHisto1D(G4int i, G4int nb, G4double x1, G4double x2, G4double u)
205{
206  if(i>=0 && i<nHisto) {
207    if(verbose > 0) 
208    {
209      G4cout << "Update histogram: #" << i 
210             << "  " << nb << "  " << x1 << "  " << x2 << "  " << u
211             << G4endl;
212    }
213    bins[i] = nb;
214    xmin[i] = x1;
215    xmax[i] = x2;
216    unit[i] = u;
217  } else {
218    G4cout << "Histo::setHisto1D: WARNING! wrong histogram index " << i << G4endl;
219    G4cout << "nHisto = " << nHisto << G4endl;
220  }
221}
222
223//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
224
225void Histo::fill(G4int i, G4double x, G4double w)
226{
227  if(verbose > 1) 
228    G4cout << "fill histogram: #" << i << " at x= " << x
229           << "  weight= " << w << " unit= " << unit[i]
230           << G4endl;   
231#ifdef G4ANALYSIS_USE 
232  if(i>=0 && i<nHisto) {
233    histo[i]->fill(x/unit[i], w);
234  } else {
235    G4cout << "Histo::fill: WARNING! wrong histogram index " << i << G4endl;
236  }
237#endif
238}
239
240//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
241
242void Histo::scale(G4int i, G4double x)
243{
244  if(verbose > 0) 
245    G4cout << "Scale histogram: #" << i << " by factor " << x << G4endl;   
246
247#ifdef G4ANALYSIS_USE 
248  if(i>=0 && i<nHisto) {
249    histo[i]->scale(x);
250  } else {
251    G4cout << "Histo::scale: WARNING! wrong histogram index " << i << G4endl;
252  }
253#endif
254}
255
256//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
257
258void Histo::addTuple(const G4String& w1, const G4String& w2, const G4String& w3)
259{
260  tupleId = w1;
261  tupleName = w2;
262  tupleList = w3;
263
264  G4cout<<" addTuple: Id "<<w1<<" Name "<<w2<<" List "<<w3<<G4endl;
265}
266
267//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
268
269void Histo::fillTuple(const G4String& parname, G4double x)
270{
271  if(verbose > 1) 
272    G4cout << "fill tuple by parameter <" << parname << "> = " << x << G4endl; 
273 
274#ifdef G4ANALYSIS_USE 
275  if(ntup) ntup->fill(ntup->findColumn(parname), (float)x);
276#endif
277}
278
279//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
280
281void Histo::addRow()
282{
283#ifdef G4ANALYSIS_USE
284  if(ntup) ntup->addRow();
285#endif
286} 
287
288//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
289void Histo::print(G4int i)
290{
291  G4cout<<"### Histogram  "<<i<<"  ###"<<G4endl;
292#ifdef G4ANALYSIS_USE
293  if(i>=0 && i<nHisto) {
294    G4double step = (xmax[i] - xmin[i])/G4double( bins[i]);
295    G4double x    =  xmin[i] - step*0.5;
296    G4double y, maxX=0, maxY=0;
297    G4int    maxJ=0;
298
299    for(G4int j=0; j<bins[i]; j++) {
300      x += step;
301      y  = histo[i]->binHeight(j);
302      if(maxY < y) {maxY = y; maxX = x; maxJ = j;}
303
304      G4cout<<x<<"  "<<y<<G4endl;
305    }
306    G4cout<<" maxJ  "<<maxJ<<"  maxX  "<<maxX<<"  maxY  "<<maxY<<G4endl;
307  }
308#endif
309}
310//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
311
Note: See TracBrowser for help on using the repository browser.