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