Geant4_10
exrdmHisto.cc
Go to the documentation of this file.
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: exrdmHisto.cc 68007 2013-03-13 11:28:03Z gcosmo $
27 //
30 //
31 #ifdef G4ANALYSIS_USE
32 #include <AIDA/AIDA.h>
33 #endif
34 //
35 #ifdef G4ANALYSIS_USE_ROOT
36 #include "TROOT.h"
37 #include "TApplication.h"
38 #include "TGClient.h"
39 #include "TCanvas.h"
40 #include "TSystem.h"
41 #include "TTree.h"
42 #include "TBranch.h"
43 #include "TFile.h"
44 #include "TH1D.h"
45 #include "TNtuple.h"
46 #endif
47 
48 #include "exrdmHisto.hh"
49 #include "exrdmHistoMessenger.hh"
50 #include "G4ParticleTable.hh"
51 
52 #include "G4Tokenizer.hh"
53 
54 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
55 exrdmHisto::exrdmHisto()
56 :fHistName("exrdm"), fHistType("root"),
57  fNHisto(0), fNTuple(0), fVerbose(0),
58  fDefaultAct(1)
59 {
60 #ifdef G4ANALYSIS_USE
61  fAida = 0;
62  fTree = 0;
63 #endif
64 
65 #ifdef G4ANALYSIS_USE_ROOT
66  fROOThisto.clear();
67  fROOTntup.clear();
68  fRarray.clear();
69  fRcol.clear();
70 #endif
71 
72  fActive.clear();
73  fBins.clear();
74  fXmin.clear();
75  fXmax.clear();
76  fUnit.clear();
77  fIds.clear();
78  fTitles.clear();
79  fTupleName.clear();
80  fTupleId.clear();
81  fTupleList.clear();
82  fTupleListROOT.clear();
83 
84  fMessenger = new exrdmHistoMessenger(this);
85 }
86 
87 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
88 
89 exrdmHisto::~exrdmHisto()
90 {
91 #ifdef G4ANALYSIS_USE
92  fHisto.clear();
93  fNtup.clear();
94 #endif
95 #ifdef G4ANALYSIS_USE_ROOT
96  //FIXME : G.Barrand : the below is crashy.
97  // In principle the TH are deleted
98  // when doing the TFile::Close !
99  // In fact the fHfileROOT should
100  // be deleted in Save(). And I am pretty
101  // sure that the TApplication is not needed.
102  //
103  // removed by F.Lei
104  // for(G4int i=0; i<fNHisto; i++) {
105  // if(fROOThisto[i]) delete fROOThisto[i];
106  // }
107 #endif
108  delete fMessenger;
109 }
110 
111 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
112 
113 void exrdmHisto::Book()
114 {
115 #ifdef G4ANALYSIS_USE
116  G4cout << "### exrdmHisto books " << fNHisto << " histograms " << G4endl;
117  // Creating the analysis factory
118  fAida = AIDA_createAnalysisFactory();
119  if(!fAida) {
120  G4cout << "ERROR: can't get AIDA." << G4endl;
121  return;
122  }
123  // Creating the fTree factory
124  {AIDA::ITreeFactory* tf = fAida->createTreeFactory();
125  // Creating a fTree mapped to a new fAida file.
126  G4String fileName = fHistName + "." + fHistType;
127  if (fHistType == "root") fileName = fHistName + "_aida." + fHistType;
128  fTree = tf->create(fileName,fHistType,false,true,"");
129  delete tf;
130  if(!fTree) {
131  G4cout << "ERROR: Tree store " << fHistName << " is not created!" << G4endl;
132  return;
133  }
134  G4cout << "Tree store : " << fTree->storeName() << G4endl;}
135  // Creating a histogram factory, whose histograms will be handled by the fTree
136  {AIDA::IHistogramFactory* hf = fAida->createHistogramFactory(*fTree);
137  // Creating an 1-dimensional histograms in the root directory of the fTree
138  for(G4int i=0; i<fNHisto; i++) {
139  if(fActive[i]) {
140  if(fVerbose>1)
141  G4cout<<"Book: histogram "<< i << " id= " << fIds[i] <<G4endl;
142  G4String tit = fIds[i];
143  if(fHistType == "root") tit = "h" + fIds[i];
144  fHisto[i] = hf->createHistogram1D(tit, fTitles[i], fBins[i], fXmin[i],
145  fXmax[i]);
146  }
147  }
148  delete hf;
149  G4cout << "AIDA histograms are booked" << G4endl;}
150 
151  // Creating a tuple factory, whose tuples will be handled by the fTree
152  {AIDA::ITupleFactory* tpf = fAida->createTupleFactory( *fTree );
153  G4cout << "AIDA will Book " << fNTuple << " ntuples" << G4endl;
154  for(G4int i=0; i<fNTuple; i++) {
155  if(fTupleList[i] != "") {
156  G4cout << "Creating Ntuple: " << fTupleName[i] <<":" <<fTupleList[i]
157  << G4endl;
158  fNtup[i] = tpf->create(fTupleId[i], fTupleName[i], fTupleList[i],"");
159  }
160  }
161  delete tpf;
162  G4cout << "AIDA ntuples are booked" << G4endl;}
163 #endif
164 
165 #ifdef G4ANALYSIS_USE_ROOT
166 // new TApplication("App", ((int *)0), ((char **)0));
167  G4String fileNameROOT = fHistName + G4String(".root");
168  fHfileROOT = new TFile(fileNameROOT.c_str() ,"RECREATE","ROOT file for exRDM");
169  G4cout << "Root file: " << fileNameROOT << G4endl;
170  // Creating an 1-dimensional histograms in the root directory of the fTree
171  for(G4int i=0; i<fNHisto; i++) {
172  if(fActive[i]) {
173  G4String id = G4String("h")+fIds[i];
174  fROOThisto[i] = new TH1D(id, fTitles[i], fBins[i], fXmin[i], fXmax[i]);
175  G4cout << "ROOT Histo " << fIds[i] << " " << fTitles[i] << " booked "
176  << G4endl;
177  }
178  }
179  // Now the ntuples
180  for(G4int i=0; i<fNTuple; i++) {
181  if(fTupleListROOT[i] != "") {
182  G4String id = G4String("t")+fTupleId[i];
183  G4cout << "Creating Ntuple "<<fTupleId[i] << " in ROOT file: "
184  << fTupleName[i] << G4endl;
185  fROOTntup[i] = new TNtuple(id, fTupleName[i], fTupleListROOT[i]);
186  G4cout << "ROOT Ntuple " << id << " " << fTupleName[i] <<" "
187  << fTupleListROOT[i]<< " booked " << G4endl;
188  }
189  }
190 #endif
191 
192 }
193 
194 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
195 
196 void exrdmHisto::Save()
197 {
198 #ifdef G4ANALYSIS_USE
199  // Write histogram file
200  fTree->commit();
201  G4cout << "Closing the AIDA fTree..." << G4endl;
202  fTree->close();
203  G4cout << "Histograms and Ntuples are saved" << G4endl;
204  delete fTree;
205  fTree = 0;
206  delete fAida;
207  fAida = 0;
208  {for(G4int i=0; i<fNHisto; i++) fHisto[i] = 0;}
209  {for(G4int i=0; i<fNTuple; i++) fNtup[i] = 0;}
210 #endif
211 #ifdef G4ANALYSIS_USE_ROOT
212  G4cout << "ROOT: files writing..." << G4endl;
213  fHfileROOT->Write();
214  G4cout << "ROOT: files closing..." << G4endl;
215  fHfileROOT->Close();
216  //
217  // F.Lei added following Guy's suggestion!
218  delete fHfileROOT;
219 
220 #endif
221 }
222 
223 
224 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
225 
226 void exrdmHisto::Add1D(const G4String& id, const G4String& name, G4int nb,
227  G4double x1, G4double x2, G4double u)
228 {
229  if(fVerbose > 0) {
230  G4cout << "New histogram will be booked: #" << id << " <" << name
231  << " " << nb << " " << x1 << " " << x2 << " " << u
232  << G4endl;
233  }
234  fNHisto++;
235  x1 /= u;
236  x2 /= u;
237  fActive.push_back(fDefaultAct);
238  fBins.push_back(nb);
239  fXmin.push_back(x1);
240  fXmax.push_back(x2);
241  fUnit.push_back(u);
242  fIds.push_back(id);
243  fTitles.push_back(name);
244 #ifdef G4ANALYSIS_USE
245  fHisto.push_back(0);
246 #endif
247 #ifdef G4ANALYSIS_USE_ROOT
248  fROOThisto.push_back(0);
249 #endif
250 }
251 
252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
253 
254 void exrdmHisto::SetHisto1D(G4int i, G4int nb, G4double x1, G4double x2, G4double u)
255 {
256  if(i>=0 && i<fNHisto) {
257  if(fVerbose > 0) {
258  G4cout << "Update histogram: #" << i
259  << " " << nb << " " << x1 << " " << x2 << " " << u
260  << G4endl;
261  }
262  fBins[i] = nb;
263  fXmin[i] = x1;
264  fXmax[i] = x2;
265  fUnit[i] = u;
266  } else {
267  G4cout << "exrdmHisto::setexrdmHisto1D: WARNING! wrong histogram index "
268  << i << G4endl;
269  }
270 }
271 
272 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
273 
274 void exrdmHisto::FillHisto(G4int i, G4double x, G4double w)
275 {
276  if(fVerbose > 1) {
277  G4cout << "fill histogram: #" << i << " at x= " << x
278  << " weight= " << w
279  << G4endl;
280  }
281 #ifdef G4ANALYSIS_USE
282  if(i>=0 && i<fNHisto) {
283  fHisto[i]->fill(x/fUnit[i], w);
284  } else {
285  G4cout << "exrdmHisto::fill: WARNING! wrong AIDA histogram index "
286  << i << G4endl;
287  }
288 #endif
289 #ifdef G4ANALYSIS_USE_ROOT
290  if(i>=0 && i<fNHisto) {
291  fROOThisto[i]->Fill(x/fUnit[i],w);
292  } else {
293  G4cout << "exrdmHisto::fill: WARNING! wrong ROOT histogram index "
294  << i << G4endl;
295  }
296 #endif
297 }
298 
299 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
300 
301 void exrdmHisto::ScaleHisto(G4int i, G4double x)
302 {
303  if(fVerbose > 0) {
304  G4cout << "Scale histogram: #" << i << " by factor " << x << G4endl;
305  }
306 #ifdef G4ANALYSIS_USE
307  if(i>=0 && i<fNHisto) {
308  fHisto[i]->scale(x);
309  G4cout << "exrdmHisto::scale: WARNING! wrong AIDA histogram index "
310  << i << G4endl;
311  }
312 #endif
313 #ifdef G4ANALYSIS_USE_ROOT
314  if(i>=0 && i<fNHisto) {
315  fROOThisto[i]->Scale(x);
316  } else {
317  G4cout << "exrdmHisto::scale: WARNING! wrong ROOT histogram index "
318  << i << G4endl;
319  }
320 #endif
321 }
322 
323 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
324 
325 #ifdef G4ANALYSIS_USE
326 void exrdmHisto::AddTuple(const G4String& w1, const G4String& w2,
327  const G4String& w3 )
328 #else
329 #ifdef G4ANALYSIS_USE_ROOT
330 void exrdmHisto::AddTuple(const G4String& w1, const G4String& w2,
331  const G4String& w3 )
332 #else
333 void exrdmHisto::AddTuple(const G4String& w1, const G4String& w2,
334  const G4String& )
335 #endif
336 #endif
337 
338 {
339  //G4cout << w1 << " " << w2 << " " << w3 << G4endl;
340  fNTuple++;
341  fTupleId.push_back(w1);
342  fTupleName.push_back(w2) ;
343 #ifdef G4ANALYSIS_USE
344  fTupleList.push_back(w3);
345  fNtup.push_back(0);
346 #endif
347 
348 #ifdef G4ANALYSIS_USE_ROOT
349  std::vector<float> ar;
350  ar.clear();
351  for (size_t i = 0; i < 20; i++) ar.push_back(0.);
352  fRarray.push_back(ar);
353  // convert AIDA header to ROOT header for ntuple
354  G4Tokenizer next(w3);
355  G4String token = next();
356  G4String ROOTList1 = "" ;
357  G4int col = 0;
358  while ( token != "") {
359  token = next();
360  if (token == ",") token = next();
361  if (token.contains(",")) token.remove(token.size()-1);
362  ROOTList1 = ROOTList1 + token + G4String(":");
363  col++;
364  }
365  G4String ROOTList = ROOTList1.substr(0,ROOTList1.length()-2);
366 // G4cout << ROOTList << G4endl;
367  fTupleListROOT.push_back(ROOTList);
368  fROOTntup.push_back(0);
369  fRcol.push_back(col-1);
370 #endif
371 }
372 
373 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
374 
375 void exrdmHisto::FillTuple(G4int i, const G4String& parname, G4double x)
376 {
377  if(fVerbose > 1)
378  G4cout << "fill tuple # " << i
379  <<" with parameter <" << parname << "> = " << x << G4endl;
380 #ifdef G4ANALYSIS_USE
381  if(fNtup[i]) fNtup[i]->fill(fNtup[i]->findColumn(parname), x);
382 #endif
383 }
384 
385 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
386 
387 void exrdmHisto::FillTuple(G4int i, G4int col, G4double x)
388 {
389  if(fVerbose > 1) {
390  G4cout << "fill tuple # " << i
391  <<" in column < " << col << "> = " << x << G4endl;
392  }
393 #ifdef G4ANALYSIS_USE
394  if(fNtup[i]) fNtup[i]->fill(col,double(x));
395 #endif
396 
397 #ifdef G4ANALYSIS_USE_ROOT
398  if(fROOTntup[i]) (fRarray[i])[col] = float(x);
399 #endif
400 
401 }
402 
403 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
404 
405 void exrdmHisto::FillTuple(G4int i, const G4String& parname, G4String& x)
406 {
407  if(fVerbose > 1) {
408  G4cout << "fill tuple # " << i
409  <<" with parameter <" << parname << "> = " << x << G4endl;
410  }
411 #ifdef G4ANALYSIS_USE
412  if(fNtup[i]) fNtup[i]->fill(fNtup[i]->findColumn(parname), x);
413 #endif
414 
415 }
416 
417 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
418 
419 void exrdmHisto::AddRow(G4int i)
420 {
421  if(fVerbose > 1) G4cout << "Added a raw #" << i << " to tuple" << G4endl;
422 #ifdef G4ANALYSIS_USE
423  if(fNtup[i]) fNtup[i]->addRow();
424 #endif
425 
426 #ifdef G4ANALYSIS_USE_ROOT
427  float ar[4];
428  for (G4int j=0; j < fRcol[i]; j++) {
429 // G4cout << i << " " << fRarray[i][j] << G4endl;
430  ar[j] = fRarray[i][j];
431  }
432  if(fROOTntup[i]) fROOTntup[i]->Fill(ar);
433 #endif
434 
435 }
436 
437 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
438 
439 void exrdmHisto::SetFileName(const G4String& nam)
440 {
441  fHistName = nam;
442 }
443 
444 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
445 
446 const G4String& exrdmHisto::GetFileName() const
447 {
448  return fHistName;
449 }
450 
451 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
452 
453 void exrdmHisto::SetFileType(const G4String& nam)
454 {
455  fHistType = nam;
456 }
457 
458 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
459 
460 const G4String& exrdmHisto::FileType() const
461 {
462  return fHistType;
463 }
464 
465 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
466 
exrdmHisto::GetFileName
const G4String & GetFileName() const
Definition: exrdmHisto.cc:446
x2
Double_t x2[nxs]
Definition: Style.C:19
exrdmHisto::FillTuple
void FillTuple(G4int, const G4String &, G4double)
Definition: exrdmHisto.cc:375
G4String::remove
G4String & remove(str_size)
exrdmHisto.hh
Definition of the exrdmHisto class.
exrdmHisto::FileType
const G4String & FileType() const
Definition: exrdmHisto.cc:460
exrdmHisto::SetFileType
void SetFileType(const G4String &)
Definition: exrdmHisto.cc:453
exrdmHisto::FillHisto
void FillHisto(G4int, G4double, G4double)
Definition: exrdmHisto.cc:274
name
const XML_Char * name
Definition: expat.h:151
exrdmHistoMessenger.hh
Definition of the exrdmHistoMessenger class.
exrdmHisto::AddTuple
void AddTuple(const G4String &, const G4String &, const G4String &)
Definition: exrdmHisto.cc:333
test.x
tuple x
Definition: test.py:50
exrdmHisto::SetHisto1D
void SetHisto1D(G4int, G4int, G4double, G4double, G4double)
Definition: exrdmHisto.cc:254
G4int
int G4int
Definition: G4Types.hh:78
exrdmHisto::SetFileName
void SetFileName(const G4String &)
Definition: exrdmHisto.cc:439
G4Tokenizer.hh
exrdmHisto::AddRow
void AddRow(G4int)
Definition: exrdmHisto.cc:419
G4cout
G4GLOB_DLL std::ostream G4cout
G4ParticleTable.hh
col
Int_t col[ntarg]
Definition: Style.C:29
x1
Double_t x1[nxs]
Definition: Style.C:18
G4String::contains
G4bool contains(const std::string &) const
exrdmHisto::Book
void Book()
Definition: exrdmHisto.cc:113
G4endl
#define G4endl
Definition: G4ios.hh:61
G4double
double G4double
Definition: G4Types.hh:76
exrdmHisto::ScaleHisto
void ScaleHisto(G4int, G4double)
Definition: exrdmHisto.cc:301
exrdmHisto::Add1D
void Add1D(const G4String &, const G4String &, G4int nb=100, G4double x1=0., G4double x2=1., G4double u=1.)
Definition: exrdmHisto.cc:226
exrdmHisto::Save
void Save()
Definition: exrdmHisto.cc:196
gammaraytel.hf
tuple hf
Definition: gammaraytel.py:8
Generated on Sat Dec 14 2013 14:33:56 for Geant4_10 by   doxygen 1.8.5