Logo ROOT  
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 #include "TNamed.h"
23 
24 // these should be maybe forward declared
25 // by moving implementations in source file
26 
27 #include "TH1.h"
28 
29 class TLine;
30 class TLegend;
31 class TH1;
32 class TVirtualPad;
33 
34 namespace RooStats {
35 
36  class HybridPlot : public TNamed {
37 
38  public:
39 
40  /// Constructor
41  HybridPlot(const char* name,
42  const char* title,
43  const std::vector<double> & sb_values,
44  const std::vector<double> & b_values,
45  double testStat_data,
46  int n_bins,
47  bool verbosity=true);
48 
49  /// Destructor
50  ~HybridPlot();
51 
52  /// Draw on current pad
53  void Draw (const char* options="");
54 
55  /// All the objects are written to rootfile
56  void DumpToFile (const char* RootFileName, const char* options);
57 
58  /// Get B histo mean
59  double GetBmean(){return fB_histo->GetMean();};
60 
61  /// Get B histo RMS
62  double GetBrms(){return fB_histo->GetRMS();};
63 
64  /// Get B histo
65  TH1F * GetBhisto(){return fB_histo;}
66 
67  /// Get B histo center
68  double GetBCenter(double n_sigmas=1, bool display=false)
69  {return GetHistoCenter(fB_histo,n_sigmas,display);};
70 
71  /// Get B histo integration extremes to obtain the requested area fraction
72  /// call delete [] res to release memory
73  double* GetBIntExtremes(double frac)
74  {return GetHistoPvals(fB_histo,frac);};
75 
76  /// Get SB histo mean
77  double GetSBmean(){return fSb_histo->GetMean();};
78 
79  /// Get SB histo center
80  double GetSBCenter(double n_sigmas=1, bool display=false)
81  {return GetHistoCenter(fSb_histo,n_sigmas,display);};
82 
83  /// Get SB histo RMS
84  double GetSBrms(){return fSb_histo->GetRMS();};
85 
86  /// Get SB histo integration extremes to obtain the requested area fraction
87  /// call delete [] res to release memory
88  double* GetSBIntExtremes(double frac)
89  {return GetHistoPvals(fSb_histo,frac);};
90 
91  /// Get B histo
93 
94  /// Get the pad (or canvas) where it has been drawn
95  TVirtualPad * GetCanvas() { return fPad; }
96 
97  /// Write an image on disk
98  void DumpToImage (const char* filename);
99 
100 
101  /// Get the center of the histo
102  double GetHistoCenter(TH1* histo, double n_rms=1,bool display_result=false);
103 
104  /// Get the "effective sigmas" of the histo, call delete [] res to release memory
105  double* GetHistoPvals (TH1* histo, double percentage);
106 
107  /// Get the median of an histogram
108  double GetMedian(TH1* histo);
109 
110  private:
111 
112  TH1F* fSb_histo; // The sb Histo
113  TH1F* fSb_histo_shaded; // The sb Histo shaded
114  TH1F* fB_histo; // The b Histo
115  TH1F* fB_histo_shaded; // The b Histo shaded
116  TLine* fData_testStat_line; // The line for the data value of the test statistic
117  TLegend* fLegend; // The legend of the plot
118  TVirtualPad * fPad; // The pad where it has been drawn
119  bool fVerbose; // verbosity flag
120 
121  ClassDef(HybridPlot,1) // Provides the plots for an HybridResult
122  };
123 }
124 
125 #endif
RooStats::HybridPlot::fB_histo_shaded
TH1F * fB_histo_shaded
Definition: HybridPlot.h:115
RooStats::HybridPlot::fB_histo
TH1F * fB_histo
Definition: HybridPlot.h:114
RooStats::HybridPlot::GetSBCenter
double GetSBCenter(double n_sigmas=1, bool display=false)
Get SB histo center.
Definition: HybridPlot.h:80
TLine
Definition: TLine.h:22
RooStats::HybridPlot::GetMedian
double GetMedian(TH1 *histo)
Get the median of an histogram.
Definition: HybridPlot.cxx:352
RooStats::HybridPlot::GetSBrms
double GetSBrms()
Get SB histo RMS.
Definition: HybridPlot.h:84
RooStats::HybridPlot::fVerbose
bool fVerbose
Definition: HybridPlot.h:119
RooStats::HybridPlot::GetBrms
double GetBrms()
Get B histo RMS.
Definition: HybridPlot.h:62
TNamed.h
RooStats::HybridPlot::GetBhisto
TH1F * GetBhisto()
Get B histo.
Definition: HybridPlot.h:65
RooStats::HybridPlot::GetSBhisto
TH1F * GetSBhisto()
Get B histo.
Definition: HybridPlot.h:92
RooStats::HybridPlot::HybridPlot
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:42
RooStats::HybridPlot::DumpToFile
void DumpToFile(const char *RootFileName, const char *options)
All the objects are written to rootfile.
Definition: HybridPlot.cxx:198
RooStats::HybridPlot::GetSBIntExtremes
double * GetSBIntExtremes(double frac)
Get SB histo integration extremes to obtain the requested area fraction call delete [] res to release...
Definition: HybridPlot.h:88
RooStats::HybridPlot
Definition: HybridPlot.h:36
RooStats::HybridPlot::GetHistoPvals
double * GetHistoPvals(TH1 *histo, double percentage)
Get the "effective sigmas" of the histo, call delete [] res to release memory.
Definition: HybridPlot.cxx:302
RooStats::HybridPlot::fSb_histo
TH1F * fSb_histo
Definition: HybridPlot.h:112
RooStats::HybridPlot::GetHistoCenter
double GetHistoCenter(TH1 *histo, double n_rms=1, bool display_result=false)
Get the center of the histo.
Definition: HybridPlot.cxx:239
TH1::GetMean
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:7085
TNamed
Definition: TNamed.h:29
RooStats::HybridPlot::GetCanvas
TVirtualPad * GetCanvas()
Get the pad (or canvas) where it has been drawn.
Definition: HybridPlot.h:95
RooStats::HybridPlot::GetBmean
double GetBmean()
Get B histo mean.
Definition: HybridPlot.h:59
RooStats::HybridPlot::GetBIntExtremes
double * GetBIntExtremes(double frac)
Get B histo integration extremes to obtain the requested area fraction call delete [] res to release ...
Definition: HybridPlot.h:73
TH1::GetRMS
Double_t GetRMS(Int_t axis=1) const
Definition: TH1.h:312
TVirtualPad
Definition: TVirtualPad.h:50
RooStats
Definition: Asimov.h:19
RooStats::HybridPlot::fPad
TVirtualPad * fPad
Definition: HybridPlot.h:118
TH1F
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:572
RooStats::HybridPlot::~HybridPlot
~HybridPlot()
Destructor.
Definition: HybridPlot.cxx:124
ClassDef
#define ClassDef(name, id)
Definition: Rtypes.h:325
TH1
Definition: TH1.h:57
name
char name[80]
Definition: TGX11.cxx:110
RooStats::HybridPlot::GetBCenter
double GetBCenter(double n_sigmas=1, bool display=false)
Get B histo center.
Definition: HybridPlot.h:68
RooStats::HybridPlot::Draw
void Draw(const char *options="")
Draw on current pad.
Definition: HybridPlot.cxx:138
RooStats::HybridPlot::fLegend
TLegend * fLegend
Definition: HybridPlot.h:117
RooStats::HybridPlot::DumpToImage
void DumpToImage(const char *filename)
Write an image on disk.
Definition: HybridPlot.cxx:226
TLegend
Definition: TLegend.h:23
TH1.h
RooStats::HybridPlot::GetSBmean
double GetSBmean()
Get SB histo mean.
Definition: HybridPlot.h:77
RooStats::HybridPlot::fSb_histo_shaded
TH1F * fSb_histo_shaded
Definition: HybridPlot.h:113
RooStats::HybridPlot::fData_testStat_line
TLine * fData_testStat_line
Definition: HybridPlot.h:116