ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
TSelHist.cxx
Go to the documentation of this file.
1 // @(#)root/proof:$Id$
2 // Author: Sangsu Ryu 22/06/2010
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2005, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 
13 /** \class TSelHist
14 \ingroup proofbench
15 
16 PROOF selector for CPU-intensive benchmark test.
17 Events are generated and 1-D, 2-D, and/or 3-D histograms are filled.
18 
19 */
20 
21 #define TSelHist_cxx
22 
23 #include "TSelHist.h"
24 #include "TProofBenchTypes.h"
25 #include <TCanvas.h>
26 #include <TPaveText.h>
27 #include <TFormula.h>
28 #include <TF1.h>
29 #include <TH1F.h>
30 #include <TH2F.h>
31 #include <TH3F.h>
32 #include <TMath.h>
33 #include <TRandom3.h>
34 #include <TString.h>
35 #include <TStyle.h>
36 #include <TSystem.h>
37 #include <TParameter.h>
38 #include <TROOT.h>
39 
41 
42 ////////////////////////////////////////////////////////////////////////////////
43 ///Constructor
44 
46  : fHistType(0), fNHists(16), fDraw(0), fHist1D(0), fHist2D(0), fHist3D(0),
47  fRandom(0), fCHist1D(0), fCHist2D(0), fCHist3D(0)
48 {
49 }
50 
51 ////////////////////////////////////////////////////////////////////////////////
52 /// Destructor
53 
55 {
56  //if (fRandom) delete fRandom;
58 
59  // Info("TSelHist","destroying ...");
60 
61  if (!fDraw){
62  for (Int_t i=0; i < fNHists; i++) {
63  if (fHist1D && fHist1D[i] && !fOutput->FindObject(fHist1D[i])) {
64  SafeDelete(fHist1D[i]);
65  }
66  if (fHist2D && fHist2D[i] && !fOutput->FindObject(fHist2D[i])) {
67  SafeDelete(fHist2D[i]);
68  }
69  if (fHist3D && fHist3D[i] && !fOutput->FindObject(fHist3D[i])) {
70  SafeDelete(fHist3D[i]);
71  }
72  }
73  }
77 }
78 
79 ////////////////////////////////////////////////////////////////////////////////
80 /// The Begin() function is called at the start of the query.
81 /// When running with PROOF Begin() is only called on the client.
82 /// The tree argument is deprecated (on PROOF 0 is passed).
83 
84 void TSelHist::Begin(TTree * /*tree*/)
85 {
86  TString option = GetOption();
87 
88  Bool_t found_histtype=kFALSE;
89  Bool_t found_nhists=kFALSE;
90  Bool_t found_draw=kFALSE;
91 
92  TIter nxt(fInput);
93  TString sinput;
94  TObject *obj;
95 
96 
97  while ((obj = nxt())){
98  sinput=obj->GetName();
99  //Info("Begin", "object name=%s", sinput.Data());
100  if (sinput.Contains("PROOF_Benchmark_HistType")){
101  if ((fHistType = dynamic_cast<TPBHistType *>(obj))) found_histtype = kTRUE;
102  continue;
103  }
104  if (sinput.Contains("PROOF_BenchmarkNHists")){
105  TParameter<Int_t>* a=dynamic_cast<TParameter<Int_t>*>(obj);
106  if (a){
107  fNHists= a->GetVal();
108  found_nhists=kTRUE;
109  //Info("Begin", "PROOF_BenchmarkNHists=%d", fNHists);
110  }
111  else{
112  Error("Begin", "PROOF_BenchmarkNHists not type TParameter<Int_t>*");
113  }
114  continue;
115  }
116  if (sinput.Contains("PROOF_BenchmarkDraw")){
117  TParameter<Int_t>* a=dynamic_cast<TParameter<Int_t>*>(obj);
118  if (a){
119  fDraw= a->GetVal();
120  found_draw=kTRUE;
121  //Info("Begin", "PROOF_BenchmarkDraw=%d", fDraw);
122  }
123  else{
124  Error("Begin", "PROOF_BenchmarkDraw not type TParameter<Int_t>*");
125  }
126  continue;
127  }
128  }
129 
130  if (!found_histtype){
132  Warning("Begin", "PROOF_Benchmark_HistType not found; using default: %d",
133  (Int_t) fHistType->GetType());
134  }
135  if (!found_nhists){
136  Warning("Begin", "PROOF_BenchmarkNHists not found; using default: %d",
137  fNHists);
138  }
139  if (!found_draw){
140  Warning("Begin", "PROOF_BenchmarkDraw not found; using default: %d",
141  fDraw);
142  }
143 
144  if (fDraw) {
148  }
149 }
150 
151 ////////////////////////////////////////////////////////////////////////////////
152 /// The SlaveBegin() function is called after the Begin() function.
153 /// When running with PROOF SlaveBegin() is called on each slave server.
154 /// The tree argument is deprecated (on PROOF 0 is passed).
155 
156 void TSelHist::SlaveBegin(TTree * /*tree*/)
157 {
158  TString option = GetOption();
159 
160  Bool_t found_histtype=kFALSE;
161  Bool_t found_nhists=kFALSE;
162  Bool_t found_draw=kFALSE;
163 
164  TIter nxt(fInput);
165  TString sinput;
166  TObject *obj;
167 
168  while ((obj = nxt())){
169  sinput=obj->GetName();
170  //Info("SlaveBegin", "object name=%s", sinput.Data());
171  if (sinput.Contains("PROOF_Benchmark_HistType")){
172  if ((fHistType = dynamic_cast<TPBHistType *>(obj))) found_histtype = kTRUE;
173  continue;
174  }
175  if (sinput.Contains("PROOF_BenchmarkNHists")){
176  TParameter<Int_t>* a=dynamic_cast<TParameter<Int_t>*>(obj);
177  if (a){
178  fNHists=a->GetVal();
179  found_nhists=kTRUE;
180  //Info("SlaveBegin", "PROOF_BenchmarkNHists=%d", fNHists);
181  }
182  else{
183  Error("SlaveBegin", "PROOF_BenchmarkNHists not type TParameter"
184  "<Int_t>*");
185  }
186  continue;
187  }
188  if (sinput.Contains("PROOF_BenchmarkDraw")){
189  TParameter<Int_t>* a=dynamic_cast<TParameter<Int_t>*>(obj);
190  if (a){
191  fDraw=a->GetVal();
192  found_draw=kTRUE;
193  //Info("SlaveBegin", "PROOF_BenchmarkDraw=%d", fDraw);
194  }
195  else{
196  Error("SlaveBegin", "PROOF_BenchmarkDraw not type TParameter"
197  "<Int_t>*");
198  }
199  continue;
200  }
201  }
202 
203  if (!found_histtype){
205  Warning("SlaveBegin", "PROOF_Benchmark_HistType not found; using default: %d",
206  fHistType->GetType());
207  }
208  if (!found_nhists){
209  Warning("SlaveBegin", "PROOF_BenchmarkNHists not found; using default: %d",
210  fNHists);
211  }
212  if (!found_draw){
213  Warning("SlaveBegin", "PROOF_BenchmarkDraw not found; using default: %d",
214  fDraw);
215  }
216 
217  // Create the histogram
219  fHist1D = new TH1F*[fNHists];
220  for (Int_t i=0; i < fNHists; i++) {
221  fHist1D[i] = new TH1F(Form("h1d_%d",i), Form("h1d_%d",i), 100, -3., 3.);
223  if (fDraw) fOutput->Add(fHist1D[i]);
224  }
225  }
227  fHist2D = new TH2F*[fNHists];
228  for (Int_t i=0; i < fNHists; i++) {
229  fHist2D[i] = new TH2F(Form("h2d_%d",i), Form("h2d_%d",i), 100, -3., 3.,
230  100, -3., 3.);
232  if (fDraw) fOutput->Add(fHist2D[i]);
233  }
234  }
236  fHist3D = new TH3F*[fNHists];
237  for (Int_t i=0; i < fNHists; i++) {
238  fHist3D[i] = new TH3F(Form("h3d_%d",i), Form("h3d_%d",i), 100, -3., 3.,
239  100, -3., 3., 100, -3., 3.);
241  if (fDraw) fOutput->Add(fHist3D[i]);
242  }
243  }
244  // Set random seed
245  fRandom = new TRandom3(0);
246 }
247 
248 ////////////////////////////////////////////////////////////////////////////////
249 /// The Process() function is called for each entry in the tree (or possibly
250 /// keyed object in the case of PROOF) to be processed. The entry argument
251 /// specifies which entry in the currently loaded tree is to be processed.
252 /// It can be passed to either TSelHist::GetEntry() or TBranch::GetEntry()
253 /// to read either all or the required parts of the data. When processing
254 /// keyed objects with PROOF, the object is already loaded and is available
255 /// via the fObject pointer.
256 ///
257 /// This function should contain the "body" of the analysis. It can contain
258 /// simple or elaborate selection criteria, run algorithms on the data
259 /// of the event and typically fill histograms.
260 ///
261 /// The processing can be stopped by calling Abort().
262 ///
263 /// Use fStatus to set the return value of TTree::Process().
264 ///
265 /// The return value is currently not used.
266 
268 {
269  Double_t x, y, z;
271  for (Int_t i=0; i < fNHists; i++) {
272  if (fRandom && fHist1D[i]) {
273  x = fRandom->Gaus(0.,1.);
274  fHist1D[i]->Fill(x);
275  }
276  }
277  }
279  for (Int_t i=0; i < fNHists; i++) {
280  if (fRandom && fHist2D[i]) {
281  x = fRandom->Gaus(0.,1.);
282  y = fRandom->Gaus(0.,1.);
283  fHist2D[i]->Fill(x, y);
284  }
285  }
286  }
288  for (Int_t i=0; i < fNHists; i++) {
289  if (fRandom && fHist3D[i]) {
290  x = fRandom->Gaus(0.,1.);
291  y = fRandom->Gaus(0.,1.);
292  z = fRandom->Gaus(0.,1.);
293  fHist3D[i]->Fill(x, y, z);
294  }
295  }
296  }
297 
298  return kTRUE;
299 }
300 
301 ////////////////////////////////////////////////////////////////////////////////
302 /// The SlaveTerminate() function is called after all entries or objects
303 /// have been processed. When running with PROOF SlaveTerminate() is called
304 /// on each slave server.
305 
307 {
308 }
309 
310 ////////////////////////////////////////////////////////////////////////////////
311 /// The Terminate() function is the last function to be called during
312 /// a query. It always runs on the client, it can be used to present
313 /// the results graphically or save the results to file.
314 
316 {
317  //
318  // Create a canvas, with 100 pads
319  //
320 
321  if (!fDraw || gROOT->IsBatch()){
322  return;
323  }
324 
326  fCHist1D=dynamic_cast<TCanvas*>(gROOT->FindObject("CHist1D"));
327  if (!fCHist1D){
328  fCHist1D = new TCanvas("CHist1D","Proof TSelHist Canvas (1D)", 200, 10,
329  700,700);
331  nside = (nside*nside < fNHists) ? nside+1 : nside;
332  fCHist1D->Divide(nside,nside,0,0);
333  }
334 
335  for (Int_t i=0; i < fNHists; i++) {
336  fHist1D[i] = dynamic_cast<TH1F *>
337  (fOutput->FindObject(Form("h1d_%d",i)));
338  fCHist1D->cd(i+1);
339  if (fHist1D[i]) fHist1D[i]->Draw();
340  }
341  // Final update
342  fCHist1D->cd();
343  fCHist1D->Update();
344  }
346  fCHist2D=dynamic_cast<TCanvas*>(gROOT->FindObject("CHist2D"));
347  if (!fCHist2D){
348  fCHist2D = new TCanvas("CHist2D","Proof TSelHist Canvas (2D)", 200, 10,
349  700,700);
351  nside = (nside*nside < fNHists) ? nside+1 : nside;
352  fCHist2D->Divide(nside,nside,0,0);
353  }
354  for (Int_t i=0; i < fNHists; i++) {
355  fHist2D[i] = dynamic_cast<TH2F *>
356  (fOutput->FindObject(Form("h2d_%d",i)));
357  fCHist2D->cd(i+1);
358  if (fHist2D[i]) fHist2D[i]->Draw("SURF");
359  }
360  // Final update
361  fCHist2D->cd();
362  fCHist2D->Update();
363  }
364 
366  fCHist3D=dynamic_cast<TCanvas*>(gROOT->FindObject("CHist3D"));
367  if (!fCHist3D){
368  fCHist3D = new TCanvas("CHist3D","Proof TSelHist Canvas (3D)", 200, 10,
369  700,700);
371  nside = (nside*nside < fNHists) ? nside+1 : nside;
372  fCHist3D->Divide(nside,nside,0,0);
373  }
374 
375  fOutput->Print("a");
376  for (Int_t i=0; i < fNHists; i++) {
377  fHist3D[i] = dynamic_cast<TH3F *>
378  (fOutput->FindObject(Form("h3d_%d",i)));
379  fCHist3D->cd(i+1);
380  if (fHist3D[i]) printf("fHist3D[%d] found\n", i);
381  if (fHist3D[i]) fHist3D[i]->Draw();
382  }
383  // Final update
384  fCHist3D->cd();
385  fCHist3D->Update();
386  }
387 
388 }
virtual const char * GetOption() const
Definition: TSelector.h:65
virtual void Begin(TTree *tree)
The Begin() function is called at the start of the query.
Definition: TSelHist.cxx:84
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3159
TCanvas * fCHist3D
Definition: TSelHist.h:47
Random number generator class based on M.
Definition: TRandom3.h:29
TSelectorList * fOutput
Definition: TSelector.h:50
long long Long64_t
Definition: RtypesCore.h:69
float Float_t
Definition: RtypesCore.h:53
3-D histogram with a float per channel (see TH1 documentation)}
Definition: TH3.h:272
virtual Double_t Gaus(Double_t mean=0, Double_t sigma=1)
Samples a random number from the standard Normal (Gaussian) Distribution with the given mean and sigm...
Definition: TRandom.cxx:235
Definition: Rtypes.h:61
TObject * FindObject(const char *name) const
Find object using its name.
Definition: THashList.cxx:213
TVirtualPad * cd(Int_t subpadnumber=0)
Set current canvas & pad.
Definition: TCanvas.cxx:659
#define gROOT
Definition: TROOT.h:344
TCanvas * fCHist1D
Definition: TSelHist.h:45
Basic string class.
Definition: TString.h:137
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:570
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
virtual void SlaveTerminate()
The SlaveTerminate() function is called after all entries or objects have been processed.
Definition: TSelHist.cxx:306
#define SafeDelete(p)
Definition: RConfig.h:436
Double_t x[n]
Definition: legend1.C:17
Bool_t fDraw
Definition: TSelHist.h:40
virtual void SlaveBegin(TTree *tree)
The SlaveBegin() function is called after the Begin() function.
Definition: TSelHist.cxx:156
TCanvas * fCHist2D
Definition: TSelHist.h:46
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
Float_t z[5]
Definition: Ifit.C:16
Named parameter, streamable and storable.
Definition: TParameter.h:49
virtual void Terminate()
The Terminate() function is the last function to be called during a query.
Definition: TSelHist.cxx:315
EHistType GetType() const
TH3F ** fHist3D
Definition: TSelHist.h:43
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2878
virtual void SetFillColor(Color_t fcolor)
Definition: TAttFill.h:50
2-D histogram with a float per channel (see TH1 documentation)}
Definition: TH2.h:256
ClassImp(TSelHist) TSelHist
Constructor.
Definition: TSelHist.cxx:40
TH2F ** fHist2D
Definition: TSelHist.h:42
char * Form(const char *fmt,...)
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH3.cxx:275
virtual void Print(Option_t *option="") const
Default print for collections, calls Print(option, 1).
virtual Bool_t Process(Long64_t entry)
The Process() function is called for each entry in the tree (or possibly keyed object in the case of ...
Definition: TSelHist.cxx:267
The Canvas class.
Definition: TCanvas.h:48
virtual const char * GetName() const
Returns name of object.
Definition: TObject.cxx:415
double Double_t
Definition: RtypesCore.h:55
ClassImp(TMCParticle) void TMCParticle printf(": p=(%7.3f,%7.3f,%9.3f) ;", fPx, fPy, fPz)
Double_t y[n]
Definition: legend1.C:17
virtual ~TSelHist()
Destructor.
Definition: TSelHist.cxx:54
Mother of all ROOT objects.
Definition: TObject.h:58
Int_t fNHists
Definition: TSelHist.h:39
TList * fInput
Current object if processing object (vs. TTree)
Definition: TSelector.h:49
TPBHistType * fHistType
Definition: TSelHist.h:38
virtual void Divide(Int_t nx=1, Int_t ny=1, Float_t xmargin=0.01, Float_t ymargin=0.01, Int_t color=0)
Automatic pad generation by division.
Definition: TPad.cxx:1073
virtual void Add(TObject *obj)
Definition: TList.h:81
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:567
TH1F ** fHist1D
Definition: TSelHist.h:41
A TTree object has a header with a name and a title.
Definition: TTree.h:98
PROOF selector for CPU-intensive benchmark test.
Definition: TSelHist.h:34
const AParamType & GetVal() const
Definition: TParameter.h:77
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
virtual void Update()
Update canvas pad buffers.
Definition: TCanvas.cxx:2179
const Bool_t kTRUE
Definition: Rtypes.h:91
Int_t Fill(Double_t)
Invalid Fill method.
Definition: TH2.cxx:287
TObject * obj
TRandom3 * fRandom
Definition: TSelHist.h:44
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:904