Hi all,
I have a histogram and a fitting function that I use on a sub-range of the histogram (let’s call it sidebands; in the example that follows it is the 6-10 region). I then use the fit results to evaluate the integral of the same function in the rest of the histogram range (in this example, it is the 10-14 region) and its uncertainty.
The result that I get is that the integral in the 10-14 region has a much larger uncertainty than the integral in the 6-10 region. Is this because the fitting function is only defined between 6 and 10?
Many thanks,
–Christos
PS.
$ which root
/afs/cern.ch/sw/lcg/app/releases/ROOT/5.27.06/slc4_ia32_gcc34/root//bin/root
#include <TROOT.h>
#include <TH1F.h>
#include <TF1.h>
#include <TMath.h>
#include <TStyle.h>
#include <TRandom.h>
#include <TVirtualFitter.h>
#include <iostream>
#include <fstream>
#include <vector>
#include <set>
#include <map>
float BIN_SIZE = 0;
inline Double_t my_gauss(Double_t * x, Double_t * par)
{
// par[0]: # of events, par[1]: mean, par[2]: sigma
Double_t arg = 0;
if (par[2]<0) par[2]=-par[2]; // make sure sigma > 0
if (par[2] != 0) arg = (x[0] - par[1])/par[2];
return par[0]*BIN_SIZE*TMath::Exp(-0.5*arg*arg)/
(TMath::Sqrt(2*TMath::Pi())*par[2]);
}
using namespace std;
int testIntegralError(void)
{
gStyle->Reset();
// mean of gaussian
const float mean = 10;
// upper, lower limits of distribution
const float LOW_LIMIT = mean -4;
const float HI_LIMIT = mean + 4;
const float fit_limit = mean;
// this histogram is the "data"
TH1F * my_hist = new TH1F("my_hist", "Gauss distribution", 200,
LOW_LIMIT, HI_LIMIT);
BIN_SIZE = my_hist->GetBinWidth(0);
// This function will be used for "sampling" values into a histogram
TF1 * fgen = new TF1("gauss1", my_gauss, LOW_LIMIT, HI_LIMIT, 3);
fgen->SetParameters(1, mean, 1); // # of events, mean, sigma
const int NENTRIES = 15000;
// this function will be used to fit the histogram
TF1 * fhist = new TF1("gauss2", my_gauss, LOW_LIMIT, fit_limit, 3);
fhist->SetParNames("Evt count", "Mean", "Sigma");
// initial "guess" for the parameters
fhist->SetParameters(NENTRIES, mean, 0.8); // # of events, mean, sigma
fhist->SetLineColor(kRed);
gStyle->SetOptFit(1111);
gStyle->SetStatH(0.0);
gStyle->SetStatW(0.0);
#if 0
fhist->SetParLimits(1, mean, mean);
#endif
for(int i = 0; i != NENTRIES; ++i)
my_hist->Fill(fgen->GetRandom());
my_hist->Fit(fhist, "ILR");
TVirtualFitter * fitter = TVirtualFitter::Fitter(0,3);
Double_t par[3] = {fitter->GetParameter(0), fitter->GetParameter(1),
fitter->GetParameter(2)};
cout << " Histogram integral: " << my_hist->Integral() << endl;
cout << " Fit-function integral between " << LOW_LIMIT << " and "
<< fit_limit << ": " << fhist->Integral(LOW_LIMIT, fit_limit)
/BIN_SIZE
<< " +- " << fhist->IntegralError(LOW_LIMIT, fit_limit)/BIN_SIZE
<< endl;
cout << " Fit-function integral between " << fit_limit << " and "
<< HI_LIMIT << ": "
<< fhist->Integral(fit_limit, HI_LIMIT)/BIN_SIZE
<< " +- "
<< fhist->IntegralError(fit_limit, HI_LIMIT, par,
fitter->GetCovarianceMatrix())/BIN_SIZE
<< endl;
return 0;
}