How to Use Unbin Fit in 1D Histogram(TH1)

[ZuyiChen]

Hi
I am learning the Unbin Fit of root in 1D Histogram(TH1F). When I was following the instruction of unbin fit from the root web , I was confused with the code it provided. See following original code I copied:

double * buffer = histogram->GetBuffer();
// Number of entry is first entry in the buffer.
int n = buffer[0];
// When creating the data object,
// it is important to create it with the size of the data.
ROOT::Fit::UnBinData data(n);
for (int i = 0; i < n; ++i)
   data.add(buffer[2*i+1]);

However, when I run “int n = buffer[0];”, my code broke down. What’s more, I found that no matter what kind of histogram I created(like “gaus”), the command “histogram->GetBufferLength()” or “histogram->GetBufferSize()” always return result “0”. Like following simple code:

TH1F *histo = new TH1F("histo", "istogramma", 100, -3, 3);
histo->FillRandom("gaus",5000);
cout << histo->GetBufferSize() << '\n';
cout << histo->GetBufferLength() << '\n';  

Could any one show me where the problem is?

[bellenot]

Maybe @moneta can help and explain why GetBuffer() returns a nullptr and GetBufferSize() returns 0

[jonas]

Hi, welcome to the ROOT forum!

A TH1 is a histogram, so it is meant to reduce the data into bins and not keep the original unbinned dataset! It’s true that there is this buffer, but the buffer might or might not be cleared already. It is an implementation detail of the histogram.

Showing this example on the website with this hack of using the TH1 buffer in an attempt to get back the unbinned dataset is bad, and I suggest to remove it from the website.

So what you need to do is to fill your unbinned data directly from e.g. random numbers, a TTree, or whatever other collections. But not from a TH1.

For example, where we fill a TH1 and a ROOT::Fit::UnBinData with the same random numbers:

int nEntries = 5000;

TH1F *histo = new TH1F("histo", "istogramma", 100, -3, 3);
ROOT::Fit::UnBinData data(nEntries);

auto gaus = static_cast<TF1*>(gROOT->GetFunction("gaus"));


for(int i = 0; i < nEntries; ++i) {
   double val = gaus->GetRandom(-3, 3);
   histo->Fill(val);
   data.Add(val);
}

I hope this explanation helps!

Cheers,
Jonas

[moneta]

Hello,

It is correct using the buffer of the histogram to retrieve its unbinned data is a low-level (expert) functionality which is not fully documented.
The buffer of the histogram is designed for computing automatically the histogram axis limits in case they are knot known in advance (see ROOT: TH1 Class Reference ).
Your example will work, but you need to call first:

histo->SetBuffer(nmax_events)

to create the histogram buffer. An alternative is to provide in the constructor a xmin axis value larger than xmax:

TH1F *histo = new TH1F("histo", "istogramma", 100, 1, 0);

Note that in this case the default buffer size will be used (1000) events. If you fill the histogram with more events the buffer will be deleted.

Lorenzo

[ZuyiChen]

Thank you so much. I understand now. But here comes another problem(I’m afraid that I have to bother you again).
I used a TTree Class to fill my unbinned data. But I cannot get the correct fitting result even if I set initial parameter values and limits, or let us say, it returned the boundary of the limits I set. Could you help me with this new question? Here is my main code:

TTree *tree = new TTree("TreeName","");// I tried to create a tree with only one branch.
tree->Branch("mq1",&mq1,"mq1/D");
in.open("mvsq.txt");// original data, the forum forbidden new user from uploading files
while(in>>mq1)
{
    if(mq1>xrange1&&mq1<xrange2)
     {
        tree->Fill();
     }
}
in.close();
tree->Print();

//Unbindata defination
int nevt = tree->Draw("mq1","","goff");
std::vector<double> dataX(  tree->GetV1(),  tree->GetV1() + nevt); 
ROOT::Fit::UnBinData data(nevt,dataX.data());

//Function defination and parameter settings

TF1 *func=new TF1("func","[0]*exp(-(x-[1])*(x-[1])/(2*[2]*[2]))+[3]*exp(-(x-[4])*(x-[4])/(2*[5]*[5]))+[6]",xrange1 ,xrange2 );
double par[7]={80,1.750425 ,2e-5 ,150,1.750525 ,2e-5 ,0.00000001};
// func->SetParameters(par);
// func->SetParLimits(0,50,150);
// func->SetParLimits(1,(1.7504 ) ,(1.75045 ) );
// func->SetParLimits(2,0.000001 ,0.00005 );
// func->SetParLimits(3,110,190);
// func->SetParLimits(4,(1.7505 ) ,(1.75055 ) );
// func->SetParLimits(5,0.000001 ,0.00005 );
// func->SetParLimits(6,0,10);
ROOT::Math::WrappedMultiTF1 fitFunction(*func,func->GetNdim());
ROOT::Fit::Fitter fitter;
fitter.SetFunction(fitFunction, false);
fitter.Config().SetParamsSettings(7,par);
fitter.Config().ParSettings(0).SetLimits(0,200);
fitter.Config().ParSettings(1).SetLimits(1.7504,1.75045 );
fitter.Config().ParSettings(2).SetLimits(0.000001 ,0.00005 );
fitter.Config().ParSettings(2).SetStepSize(1e-7 );
fitter.Config().ParSettings(3).SetLimits(0,300);
fitter.Config().ParSettings(4).SetLimits(1.7505 ,1.75055 );
fitter.Config().ParSettings(5).SetLimits(0.000001 ,0.00005 );
fitter.Config().ParSettings(5).SetStepSize(1e-7);
fitter.Config().ParSettings(6).SetLimits(0,10);
fitter.Config().SetMinimizer("Minuit", "Migrad");
fitter.LikelihoodFit(data);
ROOT::Fit::FitResult result=fitter.Result();
result.Print(std::cout);

[moneta]

Hi,
Your fit probably fails because your function is not normalised. When doing an unbinned likelihood fit, you need to provide a PDF as function to fit, and its integral must be equal to 1. So you need to re-define your function accordingly.

You can also perform an unbinned likelihood fit, in this case you need to pass the extended flag to true in
Fitter::LikelihoodFit, see here, changing the line to:

fitter.LikelihoodFit(data, true))

In this case winternallythe fit will try to rnoalize your function using numerical integration. However, this might not work very well, due to some difficulty that can be encountered in computing the integral.
If you still will have issue, please post your data as well so we can try to reproduce your problems

Lorenzo

[ZuyiChen]

Thank you for your explanation . But I am still confused with what you said " PDF function". I used the command “fitter.LikelihoodFit(data, true))” but it still doesn’t work. I think it is more efficient to upload my original data.
Mg21O142.cpp (2.7 KB)
mvsq.txt (207.4 KB)

[ZuyiChen]

Thank you for your explanation . But I am still confused with what you said " PDF function". I used the command “fitter.LikelihoodFit(data, true))” but it still didn’t work. I think it is more efficient to upload my original data.
Mg21O142.cpp (2.7 KB)
mvsq.txt (207.4 KB)