MINT2
Public Member Functions | Protected Attributes | List of all members
HistogramBase Class Reference

#include <HistogramBase.h>

Inheritance diagram for HistogramBase:
HyperHistogram

Public Member Functions

 HistogramBase (int nBins)
 
virtual ~HistogramBase ()
 
int checkBinNumber (int bin) const
 
void resetBinContents (int nBins)
 
void clear ()
 
void fillBase (int binNum, double weight)
 
void setBinContent (int bin, double val)
 
void setBinError (int bin, double val)
 
virtual void merge (const HistogramBase &other)
 
double getBinContent (int bin) const
 
double getBinError (int bin) const
 
int getNBins () const
 
void divide (const HistogramBase &other)
 
void multiply (const HistogramBase &other)
 
void add (const HistogramBase &other)
 
void minus (const HistogramBase &other)
 
void pulls (const HistogramBase &other)
 
void pulls (const HistogramBase &other1, const HistogramBase &other2)
 
void asymmetry (const HistogramBase &other)
 
void asymmetry (const HistogramBase &other1, const HistogramBase &other2)
 
void drawPullHistogram (const HistogramBase &other, TString name, int nBins=50, double pmLimits=3.5) const
 
double chi2 (const HistogramBase &other) const
 
double pvalue (const HistogramBase &other, int ndof=-1) const
 
double chi2sig (const HistogramBase &other, int ndof=-1) const
 
double integral () const
 
double integralError () const
 
void randomiseWithinErrors (int seed)
 
double getMin () const
 
double getMax () const
 
void setMin (double min)
 
void setMax (double max)
 
double getMinDensity () const
 
double getMaxDensity () const
 
void setMinDensity (double min)
 
void setMaxDensity (double max)
 
void saveBase ()
 
void saveBase (TString filename)
 
void loadBase (TString filename)
 
void normalise (double area=1.0)
 
virtual double getBinVolume (int bin) const
 
double getFrequencyDensity (int bin) const
 
void reserveCapacity (int nElements)
 
void makeFrequencyDensity ()
 
void print ()
 

Protected Attributes

int _nBins
 
std::vector< double > _binContents
 
std::vector< double > _sumW2
 
double _min
 
double _max
 
double _minDensity
 
double _maxDensity
 

Detailed Description

HyperPlot, Author: Sam Harnew, sam.h.nosp@m.arne.nosp@m.w@gma.nosp@m.il.c.nosp@m.om , Date: Dec 2015

The base class for any histogram object (my version of a TH1)

Definition at line 30 of file HistogramBase.h.

Constructor & Destructor Documentation

◆ HistogramBase()

HistogramBase::HistogramBase ( int  nBins)

Construct a histogram base with a specified number of bins

Definition at line 8 of file HistogramBase.cpp.

8  :
9  _nBins(nBins),
10  _binContents(nBins+1,0.0),
11  _sumW2 (nBins+1,0.0),
12  _min(-999.9),
13  _max(-999.9),
14  _minDensity(-999.9),
15  _maxDensity(-999.9)
16 {
17  WELCOME_LOG << "Good day from the HistogramBase() Constructor";
18 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
double _minDensity
Definition: HistogramBase.h:40
#define WELCOME_LOG
double _maxDensity
Definition: HistogramBase.h:41

◆ ~HistogramBase()

HistogramBase::~HistogramBase ( )
virtual

Destructor

Definition at line 148 of file HistogramBase.cpp.

148  {
149  GOODBYE_LOG << "Goodbye from the HistogramBase() Constructor";
150 }
#define GOODBYE_LOG

Member Function Documentation

◆ add()

void HistogramBase::add ( const HistogramBase other)

Add this hitogram to another. Histograms must have the same number of bins

Definition at line 286 of file HistogramBase.cpp.

286  {
287 
288  if (other._binContents.size() != _binContents.size()){
289  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
290  return;
291  }
292 
293  for(unsigned int i = 0; i < _binContents.size(); i++){
294 
295  double binCont1 = _binContents.at(i);
296  double binCont2 = other._binContents.at(i);
297  double var1 = _sumW2.at(i);
298  double var2 = other._sumW2.at(i);
299 
300  double content = binCont1 + binCont2;
301  double var = var1 + var2;
302 
303  _binContents.at(i) = content;
304  _sumW2.at(i) = var;
305 
306  }
307 
308  //reset min and max to default values
309  _min = -999.9;
310  _max = -999.9;
311 
312 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ asymmetry() [1/2]

void HistogramBase::asymmetry ( const HistogramBase other)

Definition at line 404 of file HistogramBase.cpp.

404  {
405 
406  if (other._binContents.size() != _binContents.size()){
407  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
408  return;
409  }
410 
411  for( int i = 0; i < getNBins(); i++){
412 
413  double binCont1 = _binContents.at(i);
414  double binCont2 = other._binContents.at(i);
415  double var1 = _sumW2.at(i);
416  double var2 = other._sumW2.at(i);
417 
418  double content = (binCont1 - binCont2)/(binCont1 + binCont2);
419 
420  double dzdx = (2.0*binCont2)/((binCont1 + binCont2)*(binCont1 + binCont2));
421  double dzdy = (2.0*binCont1)/((binCont1 + binCont2)*(binCont1 + binCont2));
422  double varsq = dzdx*dzdx*var1*var1 + dzdy*dzdy*var2*var2;
423 
424  _binContents.at(i) = content;
425  _sumW2.at(i) = varsq;
426 
427  }
428 
429  //reset min and max to default values
430  _min = -999.9;
431  _max = -999.9;
432 
433 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ asymmetry() [2/2]

void HistogramBase::asymmetry ( const HistogramBase other1,
const HistogramBase other2 
)

Definition at line 435 of file HistogramBase.cpp.

435  {
436 
437  if (other1._binContents.size() != _binContents.size() || other2._binContents.size() != _binContents.size()){
438  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
439  return;
440  }
441 
442  for( int i = 0; i < getNBins(); i++){
443 
444  double binCont1 = other1._binContents.at(i);
445  double binCont2 = other2._binContents.at(i);
446  double var1 = other1._sumW2.at(i);
447  double var2 = other2._sumW2.at(i);
448 
449  double content = (binCont1 - binCont2)/(binCont1 + binCont2);
450 
451  double dzdx = (2.0*binCont2)/((binCont1 + binCont2)*(binCont1 + binCont2));
452  double dzdy = (2.0*binCont1)/((binCont1 + binCont2)*(binCont1 + binCont2));
453  double varsq = dzdx*dzdx*var1*var1 + dzdy*dzdy*var2*var2;
454 
455  _binContents.at(i) = content;
456  _sumW2.at(i) = varsq;
457 
458  }
459 
460  //reset min and max to default values
461  _min = -999.9;
462  _max = -999.9;
463 
464 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ checkBinNumber()

int HistogramBase::checkBinNumber ( int  bin) const

Check if it's a valid bin number, if not return the overflow/underflow bin

Definition at line 163 of file HistogramBase.cpp.

163  {
164  if(bin == -1){
165  //This is also a valid way to add to the underflow/overflow bin with no error message
166  return _nBins;
167  }
168  if(bin > _nBins || bin < 0) {
169  ERROR_LOG << "Bin " << bin << " does not exist! Adding to underflow/overflow bin " << _nBins;
170  return _nBins;
171  }
172  return bin;
173 }
#define ERROR_LOG

◆ chi2()

double HistogramBase::chi2 ( const HistogramBase other) const

Calculate the chi2 between this histogram and another. Histograms must have the same number of bins

Definition at line 515 of file HistogramBase.cpp.

515  {
516 
517  double chi2 = 0.0;
518 
519  for( int i = 0; i < getNBins(); i++){
520 
521  double binCont1 = _binContents.at(i);
522  double binCont2 = other._binContents.at(i);
523  double var1 = _sumW2.at(i);
524  double var2 = other._sumW2.at(i);
525 
526  double content = binCont1 - binCont2;
527  double var = var1 + var2;
528 
529  double thisTerm = (content * content) / var;
530 
531  chi2 += thisTerm;
532 
533  }
534 
535  return chi2;
536 
537 }
double chi2(const HistogramBase &other) const
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ chi2sig()

double HistogramBase::chi2sig ( const HistogramBase other,
int  ndof = -1 
) const

Calculate the significance of the chi2 value obtained (in #sigma) given that it follows a chi2 distribution for ndof degrees of freedom. If ndof is not given it is assumed that ndof == nbins

Definition at line 504 of file HistogramBase.cpp.

504  {
505 
506  double pval = pvalue(other, ndof);
507  double sigma = TMath::NormQuantile( pval/2.0 );
508  return fabs(sigma);
509 
510 }
double pvalue(const HistogramBase &other, int ndof=-1) const

◆ clear()

void HistogramBase::clear ( )

Definition at line 28 of file HistogramBase.cpp.

28  {
29  for (unsigned i = 0; i < _binContents.size(); i++){
30  _binContents.at(i) = 0.0;
31  _sumW2 .at(i) = 0.0;
32  }
33 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ divide()

void HistogramBase::divide ( const HistogramBase other)

Divide this hitogram by another. Histograms must have the same number of bins

Definition at line 221 of file HistogramBase.cpp.

221  {
222 
223  if (other._binContents.size() != _binContents.size()){
224  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
225  return;
226  }
227 
228  for(unsigned int i = 0; i < _binContents.size(); i++){
229 
230  double binCont1 = _binContents.at(i);
231  double binCont2 = other._binContents.at(i);
232  double var1 = _sumW2.at(i);
233  double var2 = other._sumW2.at(i);
234  double frac1sq = var1/(binCont1*binCont1);
235  double frac2sq = var2/(binCont2*binCont2);
236 
237  double content = binCont1 / binCont2;
238  double var = content*content*(frac1sq + frac2sq);
239 
240  if (binCont2 == 0.0) content = 0.0;
241 
242  _binContents.at(i) = content;
243  _sumW2.at(i) = var;
244 
245  }
246 
247  //reset min and max to default values
248  _min = -999.9;
249  _max = -999.9;
250 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ drawPullHistogram()

void HistogramBase::drawPullHistogram ( const HistogramBase other,
TString  name,
int  nBins = 50,
double  pmLimits = 3.5 
) const

Draw the distribution of pulls between two histograms. They are only drawn between +- pmLimits sigma, using nBins. The plot is saved with the name 'name'.

Definition at line 542 of file HistogramBase.cpp.

542  {
543 
544  TH1D pulls ("pulls", "pulls", nBins, -pmLimits, pmLimits);
545  pulls.GetXaxis()->SetTitle("pull");
546  pulls.GetYaxis()->SetTitle("frequency");
547  for( int i = 0; i < getNBins(); i++){
548 
549  double binCont1 = _binContents.at(i);
550  double binCont2 = other._binContents.at(i);
551  double var1 = _sumW2.at(i);
552  double var2 = other._sumW2.at(i);
553 
554  double content = binCont1 - binCont2;
555  double var = var1 + var2;
556 
557  double pull = (content) / sqrt(var);
558  pulls.Fill(pull);
559 
560  }
561 
562  RootPlotter1D plotter(&pulls);
563  plotter.plot(name);
564 
565 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
void pulls(const HistogramBase &other)

◆ fillBase()

void HistogramBase::fillBase ( int  binNum,
double  weight 
)

Fill the specifed bin with a specifed weight. Note that the bin numbers run from 0 to nBins-1 inclusive so that nBins is overflow/underflow

Definition at line 155 of file HistogramBase.cpp.

155  {
156  binNum = checkBinNumber(binNum);
157  _binContents[binNum] += weight;
158  _sumW2[binNum] += weight*weight;
159 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
int checkBinNumber(int bin) const

◆ getBinContent()

double HistogramBase::getBinContent ( int  bin) const

Get the content of a bin

Definition at line 617 of file HistogramBase.cpp.

617  {
618  bin = checkBinNumber(bin);
619  return _binContents[bin];
620 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
int checkBinNumber(int bin) const

◆ getBinError()

double HistogramBase::getBinError ( int  bin) const

Get the error of a bin

Definition at line 624 of file HistogramBase.cpp.

624  {
625  bin = checkBinNumber(bin);
626  return sqrt(_sumW2[bin]);
627 }
std::vector< double > _sumW2
Definition: HistogramBase.h:36
int checkBinNumber(int bin) const

◆ getBinVolume()

double HistogramBase::getBinVolume ( int  bin) const
virtual

Get the bin volume. This is a virual function, as the bin volume depends on the binning used (which isn't defined here).

Reimplemented in HyperHistogram.

Definition at line 632 of file HistogramBase.cpp.

632  {
633  bin = checkBinNumber(bin);
634  return 1.0;
635 }
int checkBinNumber(int bin) const

◆ getFrequencyDensity()

double HistogramBase::getFrequencyDensity ( int  bin) const

Get the frequency density in a bin - this is just the bin content divided by the bin volume

Definition at line 639 of file HistogramBase.cpp.

639  {
640 
641  return getBinContent(bin)/getBinVolume(bin);
642 
643 }
virtual double getBinVolume(int bin) const
double getBinContent(int bin) const

◆ getMax()

double HistogramBase::getMax ( ) const

If the user hasn't specified a minimum using setMax, then loop over the bins and find the maximum

Definition at line 188 of file HistogramBase.cpp.

188  {
189  if (_max != -999.9) return _max;
190  MinMaxFinder stats;
191  for (int i = 0; i < _nBins; i++){
192  stats.add( _binContents[i] );
193  }
194  return stats.getMax();
195 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
void add(const double &x, const double &weight=1.0)
const double & getMax() const

◆ getMaxDensity()

double HistogramBase::getMaxDensity ( ) const

If the user hasn't specified a maximum using setMaxDensity, then loop over the bins and find the maximum density

Definition at line 210 of file HistogramBase.cpp.

210  {
211  if (_maxDensity != -999.9) return _maxDensity;
212  MinMaxFinder stats;
213  for (int i = 0; i < _nBins; i++){
214  stats.add( getFrequencyDensity(i) );
215  }
216  return stats.getMax();
217 }
void add(const double &x, const double &weight=1.0)
double getFrequencyDensity(int bin) const
const double & getMax() const
double _maxDensity
Definition: HistogramBase.h:41

◆ getMin()

double HistogramBase::getMin ( ) const

If the user hasn't specified a minimum using setMin, then loop over the bins and find the minimum

Definition at line 177 of file HistogramBase.cpp.

177  {
178  if (_min != -999.9) return _min;
179  MinMaxFinder stats;
180  for (int i = 0; i < _nBins; i++){
181  stats.add( _binContents[i] );
182  }
183  return stats.getMin();
184 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
void add(const double &x, const double &weight=1.0)
const double & getMin() const

◆ getMinDensity()

double HistogramBase::getMinDensity ( ) const

If the user hasn't specified a minimum using setMinDensity, then loop over the bins and find the minimum density

Definition at line 199 of file HistogramBase.cpp.

199  {
200  if (_minDensity != -999.9) return _minDensity;
201  MinMaxFinder stats;
202  for (int i = 0; i < _nBins; i++){
203  stats.add( getFrequencyDensity(i) );
204  }
205  return stats.getMin();
206 }
void add(const double &x, const double &weight=1.0)
const double & getMin() const
double getFrequencyDensity(int bin) const
double _minDensity
Definition: HistogramBase.h:40

◆ getNBins()

int HistogramBase::getNBins ( ) const
inline

Get the number of bins in the histogram

Definition at line 63 of file HistogramBase.h.

◆ integral()

double HistogramBase::integral ( ) const

Calcualte the integral (sum of bin contents, excluding the undeflow/overflow)

Definition at line 569 of file HistogramBase.cpp.

569  {
570  double integral = 0.0;
571  for(int i = 0; i < getNBins() ; i++){
572  integral += _binContents.at(i);
573  }
574  return integral;
575 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
double integral() const

◆ integralError()

double HistogramBase::integralError ( ) const

Calcualte the error on the integral (excluding the undeflow/overflow)

Definition at line 579 of file HistogramBase.cpp.

579  {
580  double var = 0.0;
581  for(int i = 0; i < getNBins() ; i++){
582  var += _sumW2.at(i);
583  }
584  return sqrt(var);
585 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ loadBase()

void HistogramBase::loadBase ( TString  filename)

Load the contents, sumw2, and bin numbers from a TTree in the ROOT file specified (opened using READ)

Definition at line 115 of file HistogramBase.cpp.

115  {
116 
117  TFile file(filename, "READ");
118  TTree* tree = (TTree*)file.Get("HistogramBase");
119 
120  int binNumber = -1;
121  double binContent = 0.0;
122  double sumW2 = 0.0;
123 
124  int nEntries = tree->GetEntries();
125 
126  tree->SetBranchAddress("binNumber" , &binNumber );
127  tree->SetBranchAddress("binContent", &binContent);
128  tree->SetBranchAddress("sumW2" , &sumW2 );
129 
130  this->resetBinContents(nEntries - 1);
131 
132 
133  for(int ent = 0; ent < nEntries; ent++){
134  tree->GetEntry(ent);
135 
136  VERBOSE_LOG << "Bin = " << binNumber << " Content = " << binContent << " SumW2 = " << sumW2;
137 
138  _binContents.at((int)binNumber) = binContent;
139  _sumW2 .at((int)binNumber) = sumW2;
140  }
141 
142  file.Close();
143 
144 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define VERBOSE_LOG
void resetBinContents(int nBins)

◆ makeFrequencyDensity()

void HistogramBase::makeFrequencyDensity ( )

Replace all the bin contents with frequency density.

Definition at line 647 of file HistogramBase.cpp.

647  {
648 
649  for(int i = 0; i < _nBins; i++){
650  double binVolume = this->getBinVolume(i);
651  _binContents[i] = _binContents[i] / binVolume;
652  _sumW2[i] = _sumW2[i] / (binVolume*binVolume);
653  }
654 
655 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
virtual double getBinVolume(int bin) const

◆ merge()

void HistogramBase::merge ( const HistogramBase other)
virtual

Merge one HistogramBase with another

Reimplemented in HyperHistogram.

Definition at line 46 of file HistogramBase.cpp.

46  {
47 
48  double overflowCont = _binContents.at(_nBins);
49  overflowCont += other._binContents.at(other._nBins);
50 
51  double overflowSumW2 = _sumW2.at(_nBins);
52  overflowSumW2 += other._sumW2.at(other._nBins);
53 
54  int capacityNeeded = _nBins + other._nBins + 1;
55 
56  reserveCapacity(capacityNeeded);
57 
58  _binContents.at(_nBins) = other._binContents.at(0);
59  _sumW2 .at(_nBins) = other._sumW2 .at(0);
60 
61  for (int i = 1; i < other._nBins; i++){
62  _binContents.push_back( other._binContents.at(i) );
63  _sumW2 .push_back( other._sumW2 .at(i) );
64  }
65 
66  _binContents.push_back( overflowCont );
67  _sumW2 .push_back( overflowSumW2 );
68 
69  _nBins += other._nBins;
70 
71 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
void reserveCapacity(int nElements)

◆ minus()

void HistogramBase::minus ( const HistogramBase other)

Subtract other histogram from this one. Histograms must have the same number of bins

Definition at line 316 of file HistogramBase.cpp.

316  {
317 
318  if (other._binContents.size() != _binContents.size()){
319  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
320  return;
321  }
322 
323  for(unsigned int i = 0; i < _binContents.size(); i++){
324 
325  double binCont1 = _binContents.at(i);
326  double binCont2 = other._binContents.at(i);
327  double var1 = _sumW2.at(i);
328  double var2 = other._sumW2.at(i);
329 
330  double content = binCont1 - binCont2;
331  double var = var1 + var2;
332 
333  _binContents.at(i) = content;
334  _sumW2.at(i) = var;
335 
336  }
337 
338  //reset min and max to default values
339  _min = -999.9;
340  _max = -999.9;
341 
342 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ multiply()

void HistogramBase::multiply ( const HistogramBase other)

Multiply this hitogram by another. Histograms must have the same number of bins

Definition at line 254 of file HistogramBase.cpp.

254  {
255 
256  if (other._binContents.size() != _binContents.size()){
257  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
258  return;
259  }
260 
261  for(unsigned int i = 0; i < _binContents.size(); i++){
262 
263  double binCont1 = _binContents.at(i);
264  double binCont2 = other._binContents.at(i);
265  double var1 = _sumW2.at(i);
266  double var2 = other._sumW2.at(i);
267  double frac1sq = var1/(binCont1*binCont1);
268  double frac2sq = var2/(binCont2*binCont2);
269 
270  double content = binCont1 * binCont2;
271  double var = content*content*(frac1sq + frac2sq);
272 
273  _binContents.at(i) = content;
274  _sumW2.at(i) = var;
275 
276  }
277 
278  //reset min and max to default values
279  _min = -999.9;
280  _max = -999.9;
281 
282 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ normalise()

void HistogramBase::normalise ( double  area = 1.0)

Normalise sum of bin contents to specified area. Errors are also scaled.

Definition at line 589 of file HistogramBase.cpp.

589  {
590 
591  double divisor = integral()/area;
592  double divisor2 = divisor*divisor;
593 
594  for( int i = 0; i < getNBins(); i++){
595  _binContents.at(i) = _binContents.at(i)/divisor;
596  _sumW2 .at(i) = _sumW2 .at(i)/divisor2;
597  }
598 
599 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
double integral() const

◆ print()

void HistogramBase::print ( )

Print the contents and errors of all the bins to the screen.

Definition at line 659 of file HistogramBase.cpp.

659  {
660 
661  for(int i = 0; i < _nBins; i++){
662  INFO_LOG << "Bin Content " << i << ": " << _binContents[i] << " SumW2: " << _sumW2[i];
663  }
664  INFO_LOG << "Overflow: " << _binContents[_nBins];
665 
666 }
#define INFO_LOG
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ pulls() [1/2]

void HistogramBase::pulls ( const HistogramBase other)

Find pulls between this histogram and another. Histograms must have the same number of bins

Definition at line 346 of file HistogramBase.cpp.

346  {
347 
348  if (other._binContents.size() != _binContents.size()){
349  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
350  return;
351  }
352 
353  for( int i = 0; i < getNBins(); i++){
354 
355  double binCont1 = _binContents.at(i);
356  double binCont2 = other._binContents.at(i);
357  double var1 = _sumW2.at(i);
358  double var2 = other._sumW2.at(i);
359 
360  double content = binCont1 - binCont2;
361  double var = var1 + var2;
362 
363  _binContents.at(i) = content/sqrt(var);
364  _sumW2.at(i) = 0.0;
365 
366  }
367 
368  //reset min and max to default values
369  _min = -999.9;
370  _max = -999.9;
371 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ pulls() [2/2]

void HistogramBase::pulls ( const HistogramBase other1,
const HistogramBase other2 
)

Find pulls between two histograms. Histograms must have the same number of bins

Definition at line 376 of file HistogramBase.cpp.

376  {
377 
378  if (other1._binContents.size() != _binContents.size() || other2._binContents.size() != _binContents.size()){
379  ERROR_LOG << "Trying to divide histograms with different numbers of bins. Doing nothing.";
380  return;
381  }
382 
383  for( int i = 0; i < getNBins(); i++){
384 
385  double binCont1 = other1._binContents.at(i);
386  double binCont2 = other2._binContents.at(i);
387  double var1 = other1._sumW2.at(i);
388  double var2 = other2._sumW2.at(i);
389 
390  double content = binCont1 - binCont2;
391  double var = var1 + var2;
392 
393  _binContents.at(i) = content/sqrt(var);
394  _sumW2.at(i) = 0.0;
395 
396  }
397 
398  //reset min and max to default values
399  _min = -999.9;
400  _max = -999.9;
401 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36
#define ERROR_LOG

◆ pvalue()

double HistogramBase::pvalue ( const HistogramBase other,
int  ndof = -1 
) const

Calculate the chi2 between this histogram and another, and use this to find a p-value using the specified degrees of freedom. If ndof isn't specified, nBins is used. Histograms must have the same number of bins

Definition at line 492 of file HistogramBase.cpp.

492  {
493 
494  if (ndof < 0.0) ndof = getNBins();
495  double chi2 = this->chi2(other);
496 
497  return TMath::Prob(chi2, ndof);
498 
499 }
double chi2(const HistogramBase &other) const
int getNBins() const
Definition: HistogramBase.h:63

◆ randomiseWithinErrors()

void HistogramBase::randomiseWithinErrors ( int  seed)

Randomise the histogram within it's Gaussian errors using the given seed. If value is < 0.0 then content is set to 0.0. errors remain the same before and after randomisation.

Definition at line 472 of file HistogramBase.cpp.

472  {
473 
474  TRandom3 random(seed);
475 
476  for( int i = 0; i < getNBins(); i++){
477 
478  double mean = _binContents.at(i);
479  double sigma = sqrt(_sumW2.at(i));
480 
481  double newContent = random.Gaus(mean, sigma);
482  if (newContent < 0.0) newContent = 0.0;
483  _binContents.at(i) = newContent;
484  }
485 
486 }
int getNBins() const
Definition: HistogramBase.h:63
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ reserveCapacity()

void HistogramBase::reserveCapacity ( int  nElements)

Definition at line 38 of file HistogramBase.cpp.

38  {
39  _binContents.reserve(nElements+1);
40  _sumW2 .reserve(nElements+1);
41 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ resetBinContents()

void HistogramBase::resetBinContents ( int  nBins)

Update the number of bins and set all contents to zero

Definition at line 22 of file HistogramBase.cpp.

22  {
23  _nBins = nBins;
24  _binContents = std::vector<double>(nBins+1,0.0);
25  _sumW2 = std::vector<double>(nBins+1,0.0);
26 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ saveBase() [1/2]

void HistogramBase::saveBase ( )

Save the contents, sumw2, and bin numbers in a TTree to an open TFile

Definition at line 77 of file HistogramBase.cpp.

77  {
78 
79  TTree tree("HistogramBase", "HistogramBase");
80 
81  int binNumber = -1;
82  double binContent = 0.0;
83  double sumW2 = 0.0;
84 
85  tree.Branch("binNumber" , &binNumber );
86  tree.Branch("binContent", &binContent);
87  tree.Branch("sumW2" , &sumW2 );
88 
89  for(unsigned int bin = 0; bin < _binContents.size(); bin++ ){
90  binNumber = bin;
91  binContent = _binContents.at(bin);
92  sumW2 = _sumW2 .at(bin);
93  tree.Fill();
94  }
95 
96  tree.Write();
97 
98 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
std::vector< double > _sumW2
Definition: HistogramBase.h:36

◆ saveBase() [2/2]

void HistogramBase::saveBase ( TString  filename)

Save the contents, sumw2, and bin numbers in a TTree into the ROOT file specified (opened using RECREATE)

Definition at line 102 of file HistogramBase.cpp.

102  {
103 
104  TFile file(filename, "RECREATE");
105 
106  saveBase();
107 
108  file.Write();
109  file.Close();
110 
111 }

◆ setBinContent()

void HistogramBase::setBinContent ( int  bin,
double  val 
)

Set the content of a bin (leaves error unchanged)

Definition at line 603 of file HistogramBase.cpp.

603  {
604  bin = checkBinNumber(bin);
605  _binContents[bin] = val;
606 }
std::vector< double > _binContents
Definition: HistogramBase.h:35
int checkBinNumber(int bin) const

◆ setBinError()

void HistogramBase::setBinError ( int  bin,
double  val 
)

Set the error of a bin

Definition at line 610 of file HistogramBase.cpp.

610  {
611  bin = checkBinNumber(bin);
612  _sumW2[bin] = val*val;
613 }
std::vector< double > _sumW2
Definition: HistogramBase.h:36
int checkBinNumber(int bin) const

◆ setMax()

void HistogramBase::setMax ( double  max)
inline

Set the maximum bin content for plotting

Definition at line 91 of file HistogramBase.h.

◆ setMaxDensity()

void HistogramBase::setMaxDensity ( double  max)
inline

Set the maximum frequency density for plotting

Definition at line 96 of file HistogramBase.h.

◆ setMin()

void HistogramBase::setMin ( double  min)
inline

Set the minimum bin content for plotting

Definition at line 90 of file HistogramBase.h.

◆ setMinDensity()

void HistogramBase::setMinDensity ( double  min)
inline

Set the minimum frequency density for plotting

Definition at line 95 of file HistogramBase.h.

Member Data Documentation

◆ _binContents

std::vector<double> HistogramBase::_binContents
protected

Bin contents (note that bin nBins is underflow/overflow)

Definition at line 35 of file HistogramBase.h.

◆ _max

double HistogramBase::_max
protected

Maximum bin content (for plotting)

Definition at line 39 of file HistogramBase.h.

◆ _maxDensity

double HistogramBase::_maxDensity
protected

Maximum bin density (bin content / bin volume) (for plotting)

Definition at line 41 of file HistogramBase.h.

◆ _min

double HistogramBase::_min
protected

Minimum bin content (for plotting)

Definition at line 38 of file HistogramBase.h.

◆ _minDensity

double HistogramBase::_minDensity
protected

Minimum bin density (bin content / bin volume) (for plotting)

Definition at line 40 of file HistogramBase.h.

◆ _nBins

int HistogramBase::_nBins
protected

Number of bins in the histogram

Definition at line 34 of file HistogramBase.h.

◆ _sumW2

std::vector<double> HistogramBase::_sumW2
protected

Sum of weights^2 for each bin (note that bin nBins is underflow/overflow)

Definition at line 36 of file HistogramBase.h.


The documentation for this class was generated from the following files: