Logo ROOT   6.08/07
Reference Guide
HybridPlot.h
Go to the documentation of this file.
1 // @(#)root/roostats:$Id$
2 
3 /*************************************************************************
4  * Project: RooStats *
5  * Package: RooFit/RooStats *
6  * Authors: *
7  * Kyle Cranmer, Lorenzo Moneta, Gregory Schott, Wouter Verkerke *
8  *************************************************************************
9  * Copyright (C) 1995-2008, Rene Brun and Fons Rademakers. *
10  * All rights reserved. *
11  * *
12  * For the licensing terms see $ROOTSYS/LICENSE. *
13  * For the list of contributors see $ROOTSYS/README/CREDITS. *
14  *************************************************************************/
15 
16 #ifndef ROOSTATS_HybridPlot
17 #define ROOSTATS_HybridPlot
18 
19 #include <vector>
20 #include <iostream>
21 
22 #ifndef ROOT_TNamed
23 #include "TNamed.h"
24 #endif
25 
26 // these should be maybe forward decleared
27 // by moving implementations in source file
28 #include "TH1.h"
29 
30 
31 class TLine;
32 class TLegend;
33 class TH1;
34 class TVirtualPad;
35 
36 
37 namespace RooStats {
38 
39 
40  /**
41  This class provides the plots for the result of a study performed with the
42  HybridCalculatorOriginal class.
43 
44 
45  Authors: D. Piparo, G. Schott - Universitaet Karlsruhe
46 
47 
48  An example plot is available here:
49  http://www-ekp.physik.uni-karlsruhe.de/~schott/roostats/hybridplot_example.png
50 */
51 
52 
53  class HybridPlot : public TNamed {
54 
55  public:
56 
57  /// Constructor
58  HybridPlot(const char* name,
59  const char* title,
60  const std::vector<double> & sb_values,
61  const std::vector<double> & b_values,
62  double testStat_data,
63  int n_bins,
64  bool verbosity=true);
65 
66  /// Destructor
67  ~HybridPlot();
68 
69  /// Draw on current pad
70  void Draw (const char* options="");
71 
72  /// All the objects are written to rootfile
73  void DumpToFile (const char* RootFileName, const char* options);
74 
75  /// Get B histo mean
76  double GetBmean(){return fB_histo->GetMean();};
77 
78  /// Get B histo RMS
79  double GetBrms(){return fB_histo->GetRMS();};
80 
81  /// Get B histo
82  TH1F * GetBhisto(){return fB_histo;}
83 
84  /// Get B histo center
85  double GetBCenter(double n_sigmas=1, bool display=false)
86  {return GetHistoCenter(fB_histo,n_sigmas,display);};
87 
88  /// Get B histo integration extremes to obtain the requested area fraction
89  double* GetBIntExtremes(double frac)
90  {return GetHistoPvals(fB_histo,frac);};
91 
92  /// Get SB histo mean
93  double GetSBmean(){return fSb_histo->GetMean();};
94 
95  /// Get SB histo center
96  double GetSBCenter(double n_sigmas=1, bool display=false)
97  {return GetHistoCenter(fSb_histo,n_sigmas,display);};
98 
99  /// Get SB histo RMS
100  double GetSBrms(){return fSb_histo->GetRMS();};
101 
102  /// Get SB histo integration extremes to obtain the requested area fraction
103  double* GetSBIntExtremes(double frac)
104  {return GetHistoPvals(fSb_histo,frac);};
105 
106  /// Get B histo
108 
109  /// Get the pad (or canvas) where it has been drawn
110  TVirtualPad * GetCanvas() { return fPad; }
111 
112  /// Write an image on disk
113  void DumpToImage (const char* filename);
114 
115 
116  /// Get the center of the histo
117  double GetHistoCenter(TH1* histo, double n_rms=1,bool display_result=false);
118 
119  /// Get the "effective sigmas" of the histo
120  double* GetHistoPvals (TH1* histo, double percentage);
121 
122  /// Get the median of an histogram
123  double GetMedian(TH1* histo);
124 
125  private:
126 
127  TH1F* fSb_histo; // The sb Histo
128  TH1F* fSb_histo_shaded; // The sb Histo shaded
129  TH1F* fB_histo; // The b Histo
130  TH1F* fB_histo_shaded; // The b Histo shaded
131  TLine* fData_testStat_line; // The line for the data value of the test statistic
132  TLegend* fLegend; // The legend of the plot
133  TVirtualPad * fPad; // The pad where it has been drawn
134  bool fVerbose; // verbosity flag
135 
136  ClassDef(HybridPlot,1) // Provides the plots for an HybridResult
137  };
138 }
139 
140 #endif
double * GetBIntExtremes(double frac)
Get B histo integration extremes to obtain the requested area fraction.
Definition: HybridPlot.h:89
This class displays a legend box (TPaveText) containing several legend entries.
Definition: TLegend.h:27
double * GetHistoPvals(TH1 *histo, double percentage)
Get the "effective sigmas" of the histo.
Definition: HybridPlot.cxx:302
~HybridPlot()
Destructor.
Definition: HybridPlot.cxx:117
double GetBCenter(double n_sigmas=1, bool display=false)
Get B histo center.
Definition: HybridPlot.h:85
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
Definition: TH1.cxx:6760
tomato 1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:575
void Draw(const char *options="")
Draw on current pad.
Definition: HybridPlot.cxx:131
Double_t GetRMS(Int_t axis=1) const
Definition: TH1.h:320
void DumpToImage(const char *filename)
Write an image on disk.
Definition: HybridPlot.cxx:220
#define ClassDef(name, id)
Definition: Rtypes.h:254
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
This class provides the plots for the result of a study performed with the HybridCalculatorOriginal c...
Definition: HybridPlot.h:53
double GetSBCenter(double n_sigmas=1, bool display=false)
Get SB histo center.
Definition: HybridPlot.h:96
TVirtualPad is an abstract base class for the Pad and Canvas classes.
Definition: TVirtualPad.h:59
TVirtualPad * fPad
Definition: HybridPlot.h:133
double GetSBrms()
Get SB histo RMS.
Definition: HybridPlot.h:100
HybridPlot(const char *name, const char *title, const std::vector< double > &sb_values, const std::vector< double > &b_values, double testStat_data, int n_bins, bool verbosity=true)
Constructor.
Definition: HybridPlot.cxx:33
void DumpToFile(const char *RootFileName, const char *options)
All the objects are written to rootfile.
Definition: HybridPlot.cxx:193
A simple line.
Definition: TLine.h:33
TLine * fData_testStat_line
Definition: HybridPlot.h:131
double GetSBmean()
Get SB histo mean.
Definition: HybridPlot.h:93
Namespace for the RooStats classes.
Definition: Asimov.h:20
double * GetSBIntExtremes(double frac)
Get SB histo integration extremes to obtain the requested area fraction.
Definition: HybridPlot.h:103
double GetMedian(TH1 *histo)
Get the median of an histogram.
Definition: HybridPlot.cxx:353
The TH1 histogram class.
Definition: TH1.h:80
double GetHistoCenter(TH1 *histo, double n_rms=1, bool display_result=false)
Get the center of the histo.
Definition: HybridPlot.cxx:238
double GetBrms()
Get B histo RMS.
Definition: HybridPlot.h:79
double GetBmean()
Get B histo mean.
Definition: HybridPlot.h:76
TH1F * GetSBhisto()
Get B histo.
Definition: HybridPlot.h:107
TVirtualPad * GetCanvas()
Get the pad (or canvas) where it has been drawn.
Definition: HybridPlot.h:110
TH1F * GetBhisto()
Get B histo.
Definition: HybridPlot.h:82
char name[80]
Definition: TGX11.cxx:109