ROOT logo
/*****************************************************************************
 * Project: RooFit                                                           *
 * Package: RooFitCore                                                       *
 * @(#)root/roofitcore:$Id: RooHist.cxx 29139 2009-06-22 14:31:47Z brun $
 * Authors:                                                                  *
 *   WV, Wouter Verkerke, UC Santa Barbara, verkerke@slac.stanford.edu       *
 *   DK, David Kirkby,    UC Irvine,         dkirkby@uci.edu                 *
 *                                                                           *
 * Copyright (c) 2000-2005, Regents of the University of California          *
 *                          and Stanford University. All rights reserved.    *
 *                                                                           *
 * Redistribution and use in source and binary forms,                        *
 * with or without modification, are permitted according to the terms        *
 * listed in LICENSE (http://roofit.sourceforge.net/license.txt)             *
 *****************************************************************************/

//////////////////////////////////////////////////////////////////////////////
//
// BEGIN_HTML
// A RooHist is a graphical representation of binned data based on the
// TGraphAsymmErrors class. Error bars are calculated using either Poisson
// or Binomial statistics. A RooHist is used to represent histograms in
// a RooPlot.
// END_HTML
//

#include "RooFit.h"

#include "RooHist.h"
#include "RooHist.h"
#include "RooHistError.h"
#include "RooCurve.h"
#include "RooMsgService.h"

#include "TH1.h"
#include "TClass.h"
#include "Riostream.h"
#include <iomanip>
#include <math.h>

ClassImp(RooHist)
  ;


//_____________________________________________________________________________
RooHist::RooHist() :
  _nominalBinWidth(1),
  _nSigma(1),
  _entries(0),
  _rawEntries(0)
{
  // Default constructor
}



//_____________________________________________________________________________
  RooHist::RooHist(Double_t nominalBinWidth, Double_t nSigma, Double_t /*xErrorFrac*/, Double_t /*scaleFactor*/) :
    TGraphAsymmErrors(), _nominalBinWidth(nominalBinWidth), _nSigma(nSigma), _rawEntries(-1)
{
  // Create an empty histogram that can be filled with the addBin()
  // and addAsymmetryBin() methods. Use the optional parameter to
  // specify the confidence level in units of sigma to use for
  // calculating error bars. The nominal bin width specifies the
  // default used by addBin(), and is used to set the relative
  // normalization of bins with different widths.

  initialize();
}


//_____________________________________________________________________________
RooHist::RooHist(const TH1 &data, Double_t nominalBinWidth, Double_t nSigma, RooAbsData::ErrorType etype, Double_t xErrorFrac, 
		 Bool_t correctForBinWidth, Double_t scaleFactor) :
  TGraphAsymmErrors(), _nominalBinWidth(nominalBinWidth), _nSigma(nSigma), _rawEntries(-1)
{
  // Create a histogram from the contents of the specified TH1 object
  // which may have fixed or variable bin widths. Error bars are
  // calculated using Poisson statistics. Prints a warning and rounds
  // any bins with non-integer contents. Use the optional parameter to
  // specify the confidence level in units of sigma to use for
  // calculating error bars. The nominal bin width specifies the
  // default used by addBin(), and is used to set the relative
  // normalization of bins with different widths. If not set, the
  // nominal bin width is calculated as range/nbins.

  initialize();
  // copy the input histogram's name and title
  SetName(data.GetName());
  SetTitle(data.GetTitle());
  // calculate our nominal bin width if necessary
  if(_nominalBinWidth == 0) {
    const TAxis *axis= ((TH1&)data).GetXaxis();
    if(axis->GetNbins() > 0) _nominalBinWidth= (axis->GetXmax() - axis->GetXmin())/axis->GetNbins();
  }
  // TH1::GetYaxis() is not const (why!?)
  setYAxisLabel(const_cast<TH1&>(data).GetYaxis()->GetTitle());
  
  // initialize our contents from the input histogram's contents
  Int_t nbin= data.GetNbinsX();
  for(Int_t bin= 1; bin <= nbin; bin++) {
    Axis_t x= data.GetBinCenter(bin);
    Stat_t y= data.GetBinContent(bin);
    Stat_t dy = data.GetBinError(bin) ;
    if (etype==RooAbsData::Poisson) {
      addBin(x,roundBin(y),data.GetBinWidth(bin),xErrorFrac,scaleFactor);
    } else if (etype==RooAbsData::SumW2) {
      addBinWithError(x,y,dy,dy,data.GetBinWidth(bin),xErrorFrac,correctForBinWidth,scaleFactor);
    } else {
      addBinWithError(x,y,0,0,data.GetBinWidth(bin),xErrorFrac,correctForBinWidth,scaleFactor);
    }
  }
  // add over/underflow bins to our event count
  _entries+= data.GetBinContent(0) + data.GetBinContent(nbin+1);
}



//_____________________________________________________________________________
RooHist::RooHist(const TH1 &data1, const TH1 &data2, Double_t nominalBinWidth, Double_t nSigma, 
		 Double_t xErrorFrac, Bool_t efficiency, Double_t scaleFactor) :
  TGraphAsymmErrors(), _nominalBinWidth(nominalBinWidth), _nSigma(nSigma), _rawEntries(-1)
{
  // Create a histogram from the asymmetry between the specified TH1 objects
  // which may have fixed or variable bin widths, but which must both have
  // the same binning. The asymmetry is calculated as (1-2)/(1+2). Error bars are
  // calculated using Binomial statistics. Prints a warning and rounds
  // any bins with non-integer contents. Use the optional parameter to
  // specify the confidence level in units of sigma to use for
  // calculating error bars. The nominal bin width specifies the
  // default used by addAsymmetryBin(), and is used to set the relative
  // normalization of bins with different widths. If not set, the
  // nominal bin width is calculated as range/nbins.

  initialize();
  // copy the first input histogram's name and title
  SetName(data1.GetName());
  SetTitle(data1.GetTitle());
  // calculate our nominal bin width if necessary
  if(_nominalBinWidth == 0) {
    const TAxis *axis= ((TH1&)data1).GetXaxis();
    if(axis->GetNbins() > 0) _nominalBinWidth= (axis->GetXmax() - axis->GetXmin())/axis->GetNbins();
  }

  if (!efficiency) {
    setYAxisLabel(Form("Asymmetry (%s - %s)/(%s + %s)",
		     data1.GetName(),data2.GetName(),data1.GetName(),data2.GetName()));
  } else {
    setYAxisLabel(Form("Efficiency (%s)/(%s + %s)",
		     data1.GetName(),data1.GetName(),data2.GetName()));
  }
  // initialize our contents from the input histogram contents
  Int_t nbin= data1.GetNbinsX();
  if(data2.GetNbinsX() != nbin) {
    coutE(InputArguments) << "RooHist::RooHist: histograms have different number of bins" << endl;
    return;
  }
  for(Int_t bin= 1; bin <= nbin; bin++) {
    Axis_t x= data1.GetBinCenter(bin);
    if(fabs(data2.GetBinCenter(bin)-x)>1e-10) {
      coutW(InputArguments) << "RooHist::RooHist: histograms have different centers for bin " << bin << endl;
    }
    Stat_t y1= data1.GetBinContent(bin);
    Stat_t y2= data2.GetBinContent(bin);
    if (!efficiency) {
      addAsymmetryBin(x,roundBin(y1),roundBin(y2),data1.GetBinWidth(bin),xErrorFrac,scaleFactor);
    } else {
      addEfficiencyBin(x,roundBin(y1),roundBin(y2),data1.GetBinWidth(bin),xErrorFrac,scaleFactor);
    }
  }
  // we do not have a meaningful number of entries
  _entries= -1;
}



//_____________________________________________________________________________
RooHist::RooHist(const RooHist& hist1, const RooHist& hist2, Double_t wgt1, Double_t wgt2, 
		 RooAbsData::ErrorType etype, Double_t xErrorFrac) : _rawEntries(-1)
{
  // Create histogram as sum of two existing histograms. If Poisson errors are selected the histograms are
  // added and Poisson confidence intervals are calculated for the summed content. If wgt1 and wgt2 are not
  // 1 in this mode, a warning message is printed. If SumW2 errors are selectd the histograms are added
  // and the histograms errors are added in quadrature, taking the weights into account.

  // Initialize the histogram
  initialize() ;
     
  // Copy all non-content properties from hist1
  SetName(hist1.GetName()) ;
  SetTitle(hist1.GetTitle()) ;  
  _nominalBinWidth=hist1._nominalBinWidth ;
  _nSigma=hist1._nSigma ;
  setYAxisLabel(hist1.getYAxisLabel()) ;

  if (!hist1.hasIdenticalBinning(hist2)) {
    coutE(InputArguments) << "RooHist::RooHist input histograms have incompatible binning, combined histogram will remain empty" << endl ;
    return ;
  }

  if (etype==RooAbsData::Poisson) {
    // Add histograms with Poisson errors

    // Issue warning if weights are not 1
    if (wgt1!=1.0 || wgt2 != 1.0) {
      coutW(InputArguments) << "RooHist::RooHist: WARNING: Poisson errors of weighted sum of two histograms is not well defined! " << endl
			    << "                  Summed histogram bins will rounded to nearest integer for Poisson confidence interval calculation" << endl ;
    }

    // Add histograms, calculate Poisson confidence interval on sum value
    Int_t i,n=hist1.GetN() ;
    for(i=0 ; i<n ; i++) {
      Double_t x1,y1,x2,y2,dx1 ;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
      hist1.GetPoint(i,x1,y1) ;
#else
      const_cast<RooHist&>(hist1).GetPoint(i,x1,y1) ;
#endif
      dx1 = hist1.GetErrorX(i) ;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
      hist2.GetPoint(i,x2,y2) ;
#else
      const_cast<RooHist&>(hist2).GetPoint(i,x2,y2) ;
#endif
      addBin(x1,roundBin(wgt1*y1+wgt2*y2),2*dx1/xErrorFrac,xErrorFrac) ;
    }    

  } else {
    // Add histograms with SumW2 errors

    // Add histograms, calculate combined sum-of-weights error
    Int_t i,n=hist1.GetN() ;
    for(i=0 ; i<n ; i++) {
      Double_t x1,y1,x2,y2,dx1,dy1,dy2 ;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
      hist1.GetPoint(i,x1,y1) ;
#else
      const_cast<RooHist&>(hist1).GetPoint(i,x1,y1) ;
#endif
      dx1 = hist1.GetErrorX(i) ;
      dy1 = hist1.GetErrorY(i) ;
      dy2 = hist2.GetErrorY(i) ;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
      hist2.GetPoint(i,x2,y2) ;
#else
      const_cast<RooHist&>(hist2).GetPoint(i,x2,y2) ;
#endif
      Double_t dy = sqrt(wgt1*wgt1*dy1*dy1+wgt2*wgt2*dy2*dy2) ;
      addBinWithError(x1,wgt1*y1+wgt2*y2,dy,dy,2*dx1/xErrorFrac,xErrorFrac) ;
    }       
  }

}


//_____________________________________________________________________________
void RooHist::initialize() 
{
  // Perform common initialization for all constructors.

  SetMarkerStyle(8);
  _entries= 0;
}


//_____________________________________________________________________________
Double_t RooHist::getFitRangeNEvt() const 
{
  // Return the number of events of the dataset associated with this RooHist.
  // This is the number of events in the RooHist itself, unless a different
  // value was specified through setRawEntries()

  return (_rawEntries==-1 ? _entries : _rawEntries) ;
}


//_____________________________________________________________________________
Double_t RooHist::getFitRangeNEvt(Double_t xlo, Double_t xhi) const 
{
  // Calculate integral of histogram in given range 

  Double_t sum(0) ;
  for (int i=0 ; i<GetN() ; i++) {
    Double_t x,y ;

#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
    GetPoint(i,x,y) ;
#else
    const_cast<RooHist*>(this)->GetPoint(i,x,y) ;
#endif

    if (x>=xlo && x<=xhi) {
      sum += y ;
    }
  }
  
  if (_rawEntries!=-1) {
    coutW(Plotting) << "RooHist::getFitRangeNEvt() WARNING: Number of normalization events associated to histogram is not equal to number of events in histogram" << endl
		    << "                           due cut made in RooAbsData::plotOn() call. Automatic normalization over sub-range of plot variable assumes"    << endl
		    << "                           that the effect of that cut is uniform across the plot, which may be an incorrect assumption. To be sure of"   << endl 
		    << "                           correct normalization explicit pass normalization information to RooAbsPdf::plotOn() call using Normalization()" << endl ;
    sum *= _rawEntries / _entries ;
  }

  return sum ;
}



//_____________________________________________________________________________
Double_t RooHist::getFitRangeBinW() const 
{
  // Return (average) bin width of this RooHist
  return _nominalBinWidth ;
}



//_____________________________________________________________________________
Int_t RooHist::roundBin(Double_t y) 
{
  // Return the nearest positive integer to the input value
  // and print a warning if an adjustment is required.

  if(y < 0) {
    coutW(Plotting) << fName << "::roundBin: rounding negative bin contents to zero: " << y << endl;
    return 0;
  }
  Int_t n= (Int_t)(y+0.5);
  if(fabs(y-n)>1e-6) {
    coutW(Plotting) << fName << "::roundBin: rounding non-integer bin contents: " << y << endl;
  }
  return n;
}


//_____________________________________________________________________________
void RooHist::addBin(Axis_t binCenter, Int_t n, Double_t binWidth, Double_t xErrorFrac, Double_t scaleFactor) 
{
  // Add a bin to this histogram with the specified integer bin contents
  // and using an error bar calculated with Poisson statistics. The bin width
  // is used to set the relative scale of bins with different widths.

  Double_t scale= 1;
  if(binWidth > 0) {
    scale= _nominalBinWidth/binWidth;
  }  
  _entries+= n;
  Int_t index= GetN();

  // calculate Poisson errors for this bin
  Double_t ym,yp,dx(0.5*binWidth);
  if(!RooHistError::instance().getPoissonInterval(n,ym,yp,_nSigma)) {
    coutE(Plotting) << "RooHist::addBin: unable to add bin with " << n << " events" << endl;
    return;
  }

  SetPoint(index,binCenter,n*scale*scaleFactor);
  SetPointError(index,dx*xErrorFrac,dx*xErrorFrac,scale*(n-ym)*scaleFactor,scale*(yp-n)*scaleFactor);
  updateYAxisLimits(scale*yp);
  updateYAxisLimits(scale*ym);
}



//_____________________________________________________________________________
void RooHist::addBinWithError(Axis_t binCenter, Double_t n, Double_t elow, Double_t ehigh, Double_t binWidth, 
			      Double_t xErrorFrac, Bool_t correctForBinWidth, Double_t scaleFactor) 
{
  // Add a bin to this histogram with the specified bin contents
  // and error. The bin width is used to set the relative scale of 
  // bins with different widths.

  Double_t scale= 1;
  if(binWidth > 0 && correctForBinWidth) {
    scale= _nominalBinWidth/binWidth;
  }  
  _entries+= n;
  Int_t index= GetN();

  Double_t dx(0.5*binWidth) ;
  SetPoint(index,binCenter,n*scale*scaleFactor);
  SetPointError(index,dx*xErrorFrac,dx*xErrorFrac,elow*scale*scaleFactor,ehigh*scale*scaleFactor);
  updateYAxisLimits(scale*(n-elow));
  updateYAxisLimits(scale*(n+ehigh));
}




//_____________________________________________________________________________
void RooHist::addBinWithXYError(Axis_t binCenter, Double_t n, Double_t exlow, Double_t exhigh, Double_t eylow, Double_t eyhigh, 
				Double_t scaleFactor)
{
  // Add a bin to this histogram with the specified bin contents
  // and error. The bin width is used to set the relative scale of 
  // bins with different widths.

  _entries+= n;
  Int_t index= GetN();

  SetPoint(index,binCenter,n*scaleFactor);
  SetPointError(index,exlow,exhigh,eylow*scaleFactor,eyhigh*scaleFactor);
  updateYAxisLimits(scaleFactor*(n-eylow));
  updateYAxisLimits(scaleFactor*(n+eyhigh));
}





//_____________________________________________________________________________
void RooHist::addAsymmetryBin(Axis_t binCenter, Int_t n1, Int_t n2, Double_t binWidth, Double_t xErrorFrac, Double_t scaleFactor) 
{
  // Add a bin to this histogram with the value (n1-n2)/(n1+n2)
  // using an error bar calculated with Binomial statistics.
  
  Double_t scale= 1;
  if(binWidth > 0) scale= _nominalBinWidth/binWidth;
  Int_t index= GetN();

  // calculate Binomial errors for this bin
  Double_t ym,yp,dx(0.5*binWidth);
  if(!RooHistError::instance().getBinomialIntervalAsym(n1,n2,ym,yp,_nSigma)) {
    coutE(Plotting) << "RooHist::addAsymmetryBin: unable to calculate binomial error for bin with " << n1 << "," << n2 << " events" << endl;
    return;
  }

  Double_t a= (Double_t)(n1-n2)/(n1+n2);
  SetPoint(index,binCenter,a*scaleFactor);
  SetPointError(index,dx*xErrorFrac,dx*xErrorFrac,(a-ym)*scaleFactor,(yp-a)*scaleFactor);
  updateYAxisLimits(scale*yp);
  updateYAxisLimits(scale*ym);
}



//_____________________________________________________________________________
void RooHist::addEfficiencyBin(Axis_t binCenter, Int_t n1, Int_t n2, Double_t binWidth, Double_t xErrorFrac, Double_t scaleFactor) 
{
  // Add a bin to this histogram with the value n1/(n1+n2)
  // using an error bar calculated with Binomial statistics.
  
  Double_t scale= 1;
  if(binWidth > 0) scale= _nominalBinWidth/binWidth;
  Int_t index= GetN();

  // calculate Binomial errors for this bin
  Double_t ym,yp,dx(0.5*binWidth);
  if(!RooHistError::instance().getBinomialIntervalEff(n1,n2,ym,yp,_nSigma)) {
    coutE(Plotting) << "RooHist::addEfficiencyBin: unable to calculate binomial error for bin with " << n1 << "," << n2 << " events" << endl;
    return;
  }

  Double_t a= (Double_t)(n1)/(n1+n2);
  SetPoint(index,binCenter,a*scaleFactor);
  SetPointError(index,dx*xErrorFrac,dx*xErrorFrac,(a-ym)*scaleFactor,(yp-a)*scaleFactor);
  updateYAxisLimits(scale*yp);
  updateYAxisLimits(scale*ym);
}



//_____________________________________________________________________________
RooHist::~RooHist() 
{ 
  // Destructor
}



//_____________________________________________________________________________
Bool_t RooHist::hasIdenticalBinning(const RooHist& other) const 
{
  // Return kTRUE if binning of this RooHist is identical to that of 'other'

  // First check if number of bins is the same
  if (GetN() != other.GetN()) {
    return kFALSE ;
  }

  // Next require that all bin centers are the same
  Int_t i ;
  for (i=0 ; i<GetN() ; i++) {
    Double_t x1,x2,y1,y2 ;
    
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
    GetPoint(i,x1,y1) ;
    other.GetPoint(i,x2,y2) ;
#else
    const_cast<RooHist&>(*this).GetPoint(i,x1,y1) ;
    const_cast<RooHist&>(other).GetPoint(i,x2,y2) ;
#endif

    if (fabs(x1-x2)>1e-10) {
      return kFALSE ;
    }

  }

  return kTRUE ;
}



//_____________________________________________________________________________
Bool_t RooHist::isIdentical(const RooHist& other, Double_t tol) const 
{
  // Return kTRUE if contents of this RooHIst is identical within given
  // relative tolerance to that of 'other'

  // Make temporary TH1s output of RooHists to perform kolmogorov test
  TH1::AddDirectory(kFALSE) ;
  TH1F h_self("h_self","h_self",GetN(),0,1) ;
  TH1F h_other("h_other","h_other",GetN(),0,1) ;
  TH1::AddDirectory(kTRUE) ;

  for (Int_t i=0 ; i<GetN() ; i++) {
    h_self.SetBinContent(i+1,GetY()[i]) ;
    h_other.SetBinContent(i+1,other.GetY()[i]) ;
  }  

  Double_t M = h_self.KolmogorovTest(&h_other,"M") ;
  if (M>tol) {
    Double_t kprob = h_self.KolmogorovTest(&h_other) ;
    cout << "RooHist::isIdentical() tolerance exceeded M=" << M << " (tol=" << tol << "), corresponding prob = " << kprob << endl ;
    return kFALSE ;
  }

  return kTRUE ;
}



//_____________________________________________________________________________
void RooHist::printMultiline(ostream& os, Int_t contents, Bool_t verbose, TString indent) const 
{
  // Print info about this histogram to the specified output stream.
  //
  //   Standard: number of entries
  //      Shape: error CL and maximum value
  //    Verbose: print our bin contents and errors
  
  RooPlotable::printMultiline(os,contents,verbose,indent);
  os << indent << "--- RooHist ---" << endl;
  Int_t n= GetN();
  os << indent << "  Contains " << n << " bins" << endl;
  if(verbose) {
    os << indent << "  Errors calculated at" << _nSigma << "-sigma CL" << endl;
    os << indent << "  Bin Contents:" << endl;
    for(Int_t i= 0; i < n; i++) {
      os << indent << setw(3) << i << ") x= " <<  fX[i];
      if(fEXhigh[i] > 0 || fEXlow[i] > 0) {
	os << " +" << fEXhigh[i] << " -" << fEXlow[i];
      }
      os << " , y = " << fY[i] << " +" << fEYhigh[i] << " -" << fEYlow[i] << endl;
    }
  }
}



//_____________________________________________________________________________
void RooHist::printName(ostream& os) const 
{
  // Print name of RooHist

  os << GetName() ;
}



//_____________________________________________________________________________
void RooHist::printTitle(ostream& os) const 
{
  // Print title of RooHist

  os << GetTitle() ;
}



//_____________________________________________________________________________
void RooHist::printClassName(ostream& os) const 
{
  // Print class name of RooHist

  os << IsA()->GetName() ;
}



//_____________________________________________________________________________
RooHist* RooHist::makeResidHist(const RooCurve& curve,bool normalize) const 
{
  // Create and return RooHist containing  residuals w.r.t to given curve.
  // If normalize is true, the residuals are normalized by the histogram
  // errors creating a RooHist with pull values


  // Copy all non-content properties from hist1
  RooHist* hist = new RooHist(_nominalBinWidth) ;
  hist->SetName(Form(normalize?"pull_%s_s":"resid_%s_s",GetName(),curve.GetName())) ;
  hist->SetTitle(Form(normalize?"Pull of %s and %s":"Residual of %s and %s",GetTitle(),curve.GetTitle())) ;  

  // Determine range of curve 
  Double_t xstart,xstop,y ;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
  curve.GetPoint(0,xstart,y) ;
  curve.GetPoint(curve.GetN()-1,xstop,y) ;
#else
  const_cast<RooCurve&>(curve).GetPoint(0,xstart,y) ;
  const_cast<RooCurve&>(curve).GetPoint(curve.GetN()-1,xstop,y) ;
#endif
  
  // Add histograms, calculate Poisson confidence interval on sum value
  for(Int_t i=0 ; i<GetN() ; i++) {    
    Double_t x,point;
#if ROOT_VERSION_CODE >= ROOT_VERSION(4,0,1)
    GetPoint(i,x,point) ;
#else
    const_cast<RooHist&>(*this).GetPoint(i,x,point) ;
#endif

    // Only calculate pull for bins inside curve range
    if (x<xstart || x>xstop) continue ;

    Double_t yy = point - curve.interpolate(x) ;
    Double_t dyl = GetErrorYlow(i) ;
    Double_t dyh = GetErrorYhigh(i) ;
    if (normalize) {
        Double_t norm = (yy>0?dyl:dyh);
	if (norm==0.) {
	  coutW(Plotting) << "RooHist::makeResisHist(" << GetName() << ") WARNING: point " << i << " has zero error, setting residual to zero" << endl ;
	  yy=0 ;
	  dyh=0 ;
	  dyl=0 ;
	} else {
	  yy   /= norm;
	  dyh /= norm;
	  dyl /= norm;
	}
    }
    hist->addBinWithError(x,yy,dyl,dyh);
  }
  return hist ;
}
 RooHist.cxx:1
 RooHist.cxx:2
 RooHist.cxx:3
 RooHist.cxx:4
 RooHist.cxx:5
 RooHist.cxx:6
 RooHist.cxx:7
 RooHist.cxx:8
 RooHist.cxx:9
 RooHist.cxx:10
 RooHist.cxx:11
 RooHist.cxx:12
 RooHist.cxx:13
 RooHist.cxx:14
 RooHist.cxx:15
 RooHist.cxx:16
 RooHist.cxx:17
 RooHist.cxx:18
 RooHist.cxx:19
 RooHist.cxx:20
 RooHist.cxx:21
 RooHist.cxx:22
 RooHist.cxx:23
 RooHist.cxx:24
 RooHist.cxx:25
 RooHist.cxx:26
 RooHist.cxx:27
 RooHist.cxx:28
 RooHist.cxx:29
 RooHist.cxx:30
 RooHist.cxx:31
 RooHist.cxx:32
 RooHist.cxx:33
 RooHist.cxx:34
 RooHist.cxx:35
 RooHist.cxx:36
 RooHist.cxx:37
 RooHist.cxx:38
 RooHist.cxx:39
 RooHist.cxx:40
 RooHist.cxx:41
 RooHist.cxx:42
 RooHist.cxx:43
 RooHist.cxx:44
 RooHist.cxx:45
 RooHist.cxx:46
 RooHist.cxx:47
 RooHist.cxx:48
 RooHist.cxx:49
 RooHist.cxx:50
 RooHist.cxx:51
 RooHist.cxx:52
 RooHist.cxx:53
 RooHist.cxx:54
 RooHist.cxx:55
 RooHist.cxx:56
 RooHist.cxx:57
 RooHist.cxx:58
 RooHist.cxx:59
 RooHist.cxx:60
 RooHist.cxx:61
 RooHist.cxx:62
 RooHist.cxx:63
 RooHist.cxx:64
 RooHist.cxx:65
 RooHist.cxx:66
 RooHist.cxx:67
 RooHist.cxx:68
 RooHist.cxx:69
 RooHist.cxx:70
 RooHist.cxx:71
 RooHist.cxx:72
 RooHist.cxx:73
 RooHist.cxx:74
 RooHist.cxx:75
 RooHist.cxx:76
 RooHist.cxx:77
 RooHist.cxx:78
 RooHist.cxx:79
 RooHist.cxx:80
 RooHist.cxx:81
 RooHist.cxx:82
 RooHist.cxx:83
 RooHist.cxx:84
 RooHist.cxx:85
 RooHist.cxx:86
 RooHist.cxx:87
 RooHist.cxx:88
 RooHist.cxx:89
 RooHist.cxx:90
 RooHist.cxx:91
 RooHist.cxx:92
 RooHist.cxx:93
 RooHist.cxx:94
 RooHist.cxx:95
 RooHist.cxx:96
 RooHist.cxx:97
 RooHist.cxx:98
 RooHist.cxx:99
 RooHist.cxx:100
 RooHist.cxx:101
 RooHist.cxx:102
 RooHist.cxx:103
 RooHist.cxx:104
 RooHist.cxx:105
 RooHist.cxx:106
 RooHist.cxx:107
 RooHist.cxx:108
 RooHist.cxx:109
 RooHist.cxx:110
 RooHist.cxx:111
 RooHist.cxx:112
 RooHist.cxx:113
 RooHist.cxx:114
 RooHist.cxx:115
 RooHist.cxx:116
 RooHist.cxx:117
 RooHist.cxx:118
 RooHist.cxx:119
 RooHist.cxx:120
 RooHist.cxx:121
 RooHist.cxx:122
 RooHist.cxx:123
 RooHist.cxx:124
 RooHist.cxx:125
 RooHist.cxx:126
 RooHist.cxx:127
 RooHist.cxx:128
 RooHist.cxx:129
 RooHist.cxx:130
 RooHist.cxx:131
 RooHist.cxx:132
 RooHist.cxx:133
 RooHist.cxx:134
 RooHist.cxx:135
 RooHist.cxx:136
 RooHist.cxx:137
 RooHist.cxx:138
 RooHist.cxx:139
 RooHist.cxx:140
 RooHist.cxx:141
 RooHist.cxx:142
 RooHist.cxx:143
 RooHist.cxx:144
 RooHist.cxx:145
 RooHist.cxx:146
 RooHist.cxx:147
 RooHist.cxx:148
 RooHist.cxx:149
 RooHist.cxx:150
 RooHist.cxx:151
 RooHist.cxx:152
 RooHist.cxx:153
 RooHist.cxx:154
 RooHist.cxx:155
 RooHist.cxx:156
 RooHist.cxx:157
 RooHist.cxx:158
 RooHist.cxx:159
 RooHist.cxx:160
 RooHist.cxx:161
 RooHist.cxx:162
 RooHist.cxx:163
 RooHist.cxx:164
 RooHist.cxx:165
 RooHist.cxx:166
 RooHist.cxx:167
 RooHist.cxx:168
 RooHist.cxx:169
 RooHist.cxx:170
 RooHist.cxx:171
 RooHist.cxx:172
 RooHist.cxx:173
 RooHist.cxx:174
 RooHist.cxx:175
 RooHist.cxx:176
 RooHist.cxx:177
 RooHist.cxx:178
 RooHist.cxx:179
 RooHist.cxx:180
 RooHist.cxx:181
 RooHist.cxx:182
 RooHist.cxx:183
 RooHist.cxx:184
 RooHist.cxx:185
 RooHist.cxx:186
 RooHist.cxx:187
 RooHist.cxx:188
 RooHist.cxx:189
 RooHist.cxx:190
 RooHist.cxx:191
 RooHist.cxx:192
 RooHist.cxx:193
 RooHist.cxx:194
 RooHist.cxx:195
 RooHist.cxx:196
 RooHist.cxx:197
 RooHist.cxx:198
 RooHist.cxx:199
 RooHist.cxx:200
 RooHist.cxx:201
 RooHist.cxx:202
 RooHist.cxx:203
 RooHist.cxx:204
 RooHist.cxx:205
 RooHist.cxx:206
 RooHist.cxx:207
 RooHist.cxx:208
 RooHist.cxx:209
 RooHist.cxx:210
 RooHist.cxx:211
 RooHist.cxx:212
 RooHist.cxx:213
 RooHist.cxx:214
 RooHist.cxx:215
 RooHist.cxx:216
 RooHist.cxx:217
 RooHist.cxx:218
 RooHist.cxx:219
 RooHist.cxx:220
 RooHist.cxx:221
 RooHist.cxx:222
 RooHist.cxx:223
 RooHist.cxx:224
 RooHist.cxx:225
 RooHist.cxx:226
 RooHist.cxx:227
 RooHist.cxx:228
 RooHist.cxx:229
 RooHist.cxx:230
 RooHist.cxx:231
 RooHist.cxx:232
 RooHist.cxx:233
 RooHist.cxx:234
 RooHist.cxx:235
 RooHist.cxx:236
 RooHist.cxx:237
 RooHist.cxx:238
 RooHist.cxx:239
 RooHist.cxx:240
 RooHist.cxx:241
 RooHist.cxx:242
 RooHist.cxx:243
 RooHist.cxx:244
 RooHist.cxx:245
 RooHist.cxx:246
 RooHist.cxx:247
 RooHist.cxx:248
 RooHist.cxx:249
 RooHist.cxx:250
 RooHist.cxx:251
 RooHist.cxx:252
 RooHist.cxx:253
 RooHist.cxx:254
 RooHist.cxx:255
 RooHist.cxx:256
 RooHist.cxx:257
 RooHist.cxx:258
 RooHist.cxx:259
 RooHist.cxx:260
 RooHist.cxx:261
 RooHist.cxx:262
 RooHist.cxx:263
 RooHist.cxx:264
 RooHist.cxx:265
 RooHist.cxx:266
 RooHist.cxx:267
 RooHist.cxx:268
 RooHist.cxx:269
 RooHist.cxx:270
 RooHist.cxx:271
 RooHist.cxx:272
 RooHist.cxx:273
 RooHist.cxx:274
 RooHist.cxx:275
 RooHist.cxx:276
 RooHist.cxx:277
 RooHist.cxx:278
 RooHist.cxx:279
 RooHist.cxx:280
 RooHist.cxx:281
 RooHist.cxx:282
 RooHist.cxx:283
 RooHist.cxx:284
 RooHist.cxx:285
 RooHist.cxx:286
 RooHist.cxx:287
 RooHist.cxx:288
 RooHist.cxx:289
 RooHist.cxx:290
 RooHist.cxx:291
 RooHist.cxx:292
 RooHist.cxx:293
 RooHist.cxx:294
 RooHist.cxx:295
 RooHist.cxx:296
 RooHist.cxx:297
 RooHist.cxx:298
 RooHist.cxx:299
 RooHist.cxx:300
 RooHist.cxx:301
 RooHist.cxx:302
 RooHist.cxx:303
 RooHist.cxx:304
 RooHist.cxx:305
 RooHist.cxx:306
 RooHist.cxx:307
 RooHist.cxx:308
 RooHist.cxx:309
 RooHist.cxx:310
 RooHist.cxx:311
 RooHist.cxx:312
 RooHist.cxx:313
 RooHist.cxx:314
 RooHist.cxx:315
 RooHist.cxx:316
 RooHist.cxx:317
 RooHist.cxx:318
 RooHist.cxx:319
 RooHist.cxx:320
 RooHist.cxx:321
 RooHist.cxx:322
 RooHist.cxx:323
 RooHist.cxx:324
 RooHist.cxx:325
 RooHist.cxx:326
 RooHist.cxx:327
 RooHist.cxx:328
 RooHist.cxx:329
 RooHist.cxx:330
 RooHist.cxx:331
 RooHist.cxx:332
 RooHist.cxx:333
 RooHist.cxx:334
 RooHist.cxx:335
 RooHist.cxx:336
 RooHist.cxx:337
 RooHist.cxx:338
 RooHist.cxx:339
 RooHist.cxx:340
 RooHist.cxx:341
 RooHist.cxx:342
 RooHist.cxx:343
 RooHist.cxx:344
 RooHist.cxx:345
 RooHist.cxx:346
 RooHist.cxx:347
 RooHist.cxx:348
 RooHist.cxx:349
 RooHist.cxx:350
 RooHist.cxx:351
 RooHist.cxx:352
 RooHist.cxx:353
 RooHist.cxx:354
 RooHist.cxx:355
 RooHist.cxx:356
 RooHist.cxx:357
 RooHist.cxx:358
 RooHist.cxx:359
 RooHist.cxx:360
 RooHist.cxx:361
 RooHist.cxx:362
 RooHist.cxx:363
 RooHist.cxx:364
 RooHist.cxx:365
 RooHist.cxx:366
 RooHist.cxx:367
 RooHist.cxx:368
 RooHist.cxx:369
 RooHist.cxx:370
 RooHist.cxx:371
 RooHist.cxx:372
 RooHist.cxx:373
 RooHist.cxx:374
 RooHist.cxx:375
 RooHist.cxx:376
 RooHist.cxx:377
 RooHist.cxx:378
 RooHist.cxx:379
 RooHist.cxx:380
 RooHist.cxx:381
 RooHist.cxx:382
 RooHist.cxx:383
 RooHist.cxx:384
 RooHist.cxx:385
 RooHist.cxx:386
 RooHist.cxx:387
 RooHist.cxx:388
 RooHist.cxx:389
 RooHist.cxx:390
 RooHist.cxx:391
 RooHist.cxx:392
 RooHist.cxx:393
 RooHist.cxx:394
 RooHist.cxx:395
 RooHist.cxx:396
 RooHist.cxx:397
 RooHist.cxx:398
 RooHist.cxx:399
 RooHist.cxx:400
 RooHist.cxx:401
 RooHist.cxx:402
 RooHist.cxx:403
 RooHist.cxx:404
 RooHist.cxx:405
 RooHist.cxx:406
 RooHist.cxx:407
 RooHist.cxx:408
 RooHist.cxx:409
 RooHist.cxx:410
 RooHist.cxx:411
 RooHist.cxx:412
 RooHist.cxx:413
 RooHist.cxx:414
 RooHist.cxx:415
 RooHist.cxx:416
 RooHist.cxx:417
 RooHist.cxx:418
 RooHist.cxx:419
 RooHist.cxx:420
 RooHist.cxx:421
 RooHist.cxx:422
 RooHist.cxx:423
 RooHist.cxx:424
 RooHist.cxx:425
 RooHist.cxx:426
 RooHist.cxx:427
 RooHist.cxx:428
 RooHist.cxx:429
 RooHist.cxx:430
 RooHist.cxx:431
 RooHist.cxx:432
 RooHist.cxx:433
 RooHist.cxx:434
 RooHist.cxx:435
 RooHist.cxx:436
 RooHist.cxx:437
 RooHist.cxx:438
 RooHist.cxx:439
 RooHist.cxx:440
 RooHist.cxx:441
 RooHist.cxx:442
 RooHist.cxx:443
 RooHist.cxx:444
 RooHist.cxx:445
 RooHist.cxx:446
 RooHist.cxx:447
 RooHist.cxx:448
 RooHist.cxx:449
 RooHist.cxx:450
 RooHist.cxx:451
 RooHist.cxx:452
 RooHist.cxx:453
 RooHist.cxx:454
 RooHist.cxx:455
 RooHist.cxx:456
 RooHist.cxx:457
 RooHist.cxx:458
 RooHist.cxx:459
 RooHist.cxx:460
 RooHist.cxx:461
 RooHist.cxx:462
 RooHist.cxx:463
 RooHist.cxx:464
 RooHist.cxx:465
 RooHist.cxx:466
 RooHist.cxx:467
 RooHist.cxx:468
 RooHist.cxx:469
 RooHist.cxx:470
 RooHist.cxx:471
 RooHist.cxx:472
 RooHist.cxx:473
 RooHist.cxx:474
 RooHist.cxx:475
 RooHist.cxx:476
 RooHist.cxx:477
 RooHist.cxx:478
 RooHist.cxx:479
 RooHist.cxx:480
 RooHist.cxx:481
 RooHist.cxx:482
 RooHist.cxx:483
 RooHist.cxx:484
 RooHist.cxx:485
 RooHist.cxx:486
 RooHist.cxx:487
 RooHist.cxx:488
 RooHist.cxx:489
 RooHist.cxx:490
 RooHist.cxx:491
 RooHist.cxx:492
 RooHist.cxx:493
 RooHist.cxx:494
 RooHist.cxx:495
 RooHist.cxx:496
 RooHist.cxx:497
 RooHist.cxx:498
 RooHist.cxx:499
 RooHist.cxx:500
 RooHist.cxx:501
 RooHist.cxx:502
 RooHist.cxx:503
 RooHist.cxx:504
 RooHist.cxx:505
 RooHist.cxx:506
 RooHist.cxx:507
 RooHist.cxx:508
 RooHist.cxx:509
 RooHist.cxx:510
 RooHist.cxx:511
 RooHist.cxx:512
 RooHist.cxx:513
 RooHist.cxx:514
 RooHist.cxx:515
 RooHist.cxx:516
 RooHist.cxx:517
 RooHist.cxx:518
 RooHist.cxx:519
 RooHist.cxx:520
 RooHist.cxx:521
 RooHist.cxx:522
 RooHist.cxx:523
 RooHist.cxx:524
 RooHist.cxx:525
 RooHist.cxx:526
 RooHist.cxx:527
 RooHist.cxx:528
 RooHist.cxx:529
 RooHist.cxx:530
 RooHist.cxx:531
 RooHist.cxx:532
 RooHist.cxx:533
 RooHist.cxx:534
 RooHist.cxx:535
 RooHist.cxx:536
 RooHist.cxx:537
 RooHist.cxx:538
 RooHist.cxx:539
 RooHist.cxx:540
 RooHist.cxx:541
 RooHist.cxx:542
 RooHist.cxx:543
 RooHist.cxx:544
 RooHist.cxx:545
 RooHist.cxx:546
 RooHist.cxx:547
 RooHist.cxx:548
 RooHist.cxx:549
 RooHist.cxx:550
 RooHist.cxx:551
 RooHist.cxx:552
 RooHist.cxx:553
 RooHist.cxx:554
 RooHist.cxx:555
 RooHist.cxx:556
 RooHist.cxx:557
 RooHist.cxx:558
 RooHist.cxx:559
 RooHist.cxx:560
 RooHist.cxx:561
 RooHist.cxx:562
 RooHist.cxx:563
 RooHist.cxx:564
 RooHist.cxx:565
 RooHist.cxx:566
 RooHist.cxx:567
 RooHist.cxx:568
 RooHist.cxx:569
 RooHist.cxx:570
 RooHist.cxx:571
 RooHist.cxx:572
 RooHist.cxx:573
 RooHist.cxx:574
 RooHist.cxx:575
 RooHist.cxx:576
 RooHist.cxx:577
 RooHist.cxx:578
 RooHist.cxx:579
 RooHist.cxx:580
 RooHist.cxx:581
 RooHist.cxx:582
 RooHist.cxx:583
 RooHist.cxx:584
 RooHist.cxx:585
 RooHist.cxx:586
 RooHist.cxx:587
 RooHist.cxx:588
 RooHist.cxx:589
 RooHist.cxx:590
 RooHist.cxx:591
 RooHist.cxx:592
 RooHist.cxx:593
 RooHist.cxx:594
 RooHist.cxx:595
 RooHist.cxx:596
 RooHist.cxx:597
 RooHist.cxx:598
 RooHist.cxx:599
 RooHist.cxx:600
 RooHist.cxx:601
 RooHist.cxx:602
 RooHist.cxx:603
 RooHist.cxx:604
 RooHist.cxx:605
 RooHist.cxx:606
 RooHist.cxx:607
 RooHist.cxx:608
 RooHist.cxx:609
 RooHist.cxx:610
 RooHist.cxx:611
 RooHist.cxx:612
 RooHist.cxx:613
 RooHist.cxx:614
 RooHist.cxx:615
 RooHist.cxx:616
 RooHist.cxx:617
 RooHist.cxx:618
 RooHist.cxx:619
 RooHist.cxx:620
 RooHist.cxx:621
 RooHist.cxx:622
 RooHist.cxx:623
 RooHist.cxx:624
 RooHist.cxx:625
 RooHist.cxx:626
 RooHist.cxx:627
 RooHist.cxx:628
 RooHist.cxx:629
 RooHist.cxx:630
 RooHist.cxx:631
 RooHist.cxx:632
 RooHist.cxx:633
 RooHist.cxx:634
 RooHist.cxx:635
 RooHist.cxx:636
 RooHist.cxx:637
 RooHist.cxx:638
 RooHist.cxx:639
 RooHist.cxx:640
 RooHist.cxx:641
 RooHist.cxx:642
 RooHist.cxx:643
 RooHist.cxx:644
 RooHist.cxx:645
 RooHist.cxx:646
 RooHist.cxx:647