Logo ROOT  
Reference Guide
 
Loading...
Searching...
No Matches
fitFraction.C File Reference

Detailed Description

View in nbviewer Open in SWAN
FractionFitter example À la HMCMLL, see R.

Barlow and C. Beeston, Comp. Phys. Comm. 77 (1993) 219-228,

See also
https://indico.in2p3.fr/event/2635/contributions/25070/
Note
An alternative interface is described in https://root.cern.ch/doc/master/rf709__BarlowBeeston_8C_source.html
Minuit2Minimizer: Minimize with max-calls 1345 convergence for edm < 0.01 strategy 1
Minuit2Minimizer : Valid minimum - status = 0
FVAL = -10009.679561972478
Edm = 3.20762774239878909e-07
Nfcn = 153
frac0 = 3.70672e-05 +/- 0.613661 (limited)
frac1 = 0.283223 +/- 0.0555675 (limited)
frac2 = 0.716752 +/- 0.0534245 (limited)
Status: 0
Parameter 0: true 0.010, estim. 0.000 +/- 0.614
Parameter 1: true 0.300, estim. 0.283 +/- 0.056
Parameter 2: true 0.690, estim. 0.717 +/- 0.053
#include "TH1F.h"
#include "TF1.h"
#include "TRandom.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TLatex.h"
#include <iostream>
{
// parameters and functions to generate the data
Int_t Ndata = 1000;
Int_t N0 = 1000;
Int_t N1 = 1000;
Int_t N2 = 1000;
Int_t nBins = 40;
// contribution 0
TF1 f0("f0", "[0]*(1-cos(x))/TMath::Pi()", 0., TMath::Pi());
f0.SetParameter(0, 1.);
f0.SetLineColor(2);
Double_t int0 = f0.Integral(0., TMath::Pi());
// contribution 1
TF1 f1("f1", "[0]*(1-cos(x)*cos(x))*2./TMath::Pi()", 0., TMath::Pi());
f1.SetParameter(0, 1.);
// contribution 2
TF1 f2("f2", "[0]*(1+cos(x))/TMath::Pi()", 0., TMath::Pi());
f2.SetParameter(0, 1.);
f2.SetLineColor(4);
Double_t int2 = f2.Integral(0., TMath::Pi());
// generate data
TH1F data("data", "Data angle distribution", nBins, 0, TMath::Pi()); // data histogram
data.SetXTitle("x");
data.SetMarkerStyle(20);
data.SetMarkerSize(.7);
data.SetMinimum(0);
htruemc0.SetLineColor(2);
htruemc1.SetLineColor(3);
htruemc2.SetLineColor(4);
for (Int_t i = 0; i < Ndata; i++) {
if (p < trueP0) {
x = f0.GetRandom();
htruemc0.Fill(x);
} else if (p < trueP0 + trueP1) {
x = f1.GetRandom();
htruemc1.Fill(x);
} else {
x = f2.GetRandom();
htruemc2.Fill(x);
}
data.Fill(x);
}
// generate MC samples
TH1F mc0("mc0", "MC sample 0 angle distribution", nBins, 0, TMath::Pi()); // first MC histogram
mc0.SetXTitle("x");
mc0.SetLineColor(2);
mc0.SetMarkerColor(2);
mc0.SetMarkerStyle(24);
mc0.SetMarkerSize(.7);
for (Int_t i = 0; i < N0; i++) {
mc0.Fill(f0.GetRandom());
}
TH1F mc1("mc1", "MC sample 1 angle distribution", nBins, 0, TMath::Pi()); // second MC histogram
mc1.SetXTitle("x");
mc1.SetLineColor(3);
mc1.SetMarkerColor(3);
mc1.SetMarkerStyle(24);
mc1.SetMarkerSize(.7);
for (Int_t i = 0; i < N1; i++) {
mc1.Fill(f1.GetRandom());
}
TH1F mc2("mc2", "MC sample 2 angle distribution", nBins, 0, TMath::Pi()); // third MC histogram
mc2.SetXTitle("x");
mc2.SetLineColor(4);
mc2.SetMarkerColor(4);
mc2.SetMarkerStyle(24);
mc2.SetMarkerSize(.7);
for (Int_t i = 0; i < N2; i++) {
mc2.Fill(f2.GetRandom());
}
// FractionFitter
TObjArray mc(3); // MC histograms are put in this array
mc.Add(&mc0);
mc.Add(&mc1);
mc.Add(&mc2);
TFractionFitter fit(&data, &mc); // initialise
fit.Constrain(0, 0.0, 1.0); // constrain fraction 1 to be between 0 and 1
fit.Constrain(1, 0.0, 1.0); // constrain fraction 1 to be between 0 and 1
fit.Constrain(2, 0.0, 1.0); // constrain fraction 1 to be between 0 and 1
// fit.SetRangeX(1,15); // use only the first 15 bins in the fit
Int_t status = fit.Fit(); // perform the fit
std::cout << "Status: " << status << std::endl;
// Display
TCanvas c("c", "FractionFitter example", 700, 700);
c.Divide(2, 2);
c.cd(1);
f0.DrawClone();
f0.GetHistogram()->SetTitle("Original MC distributions");
f1.DrawClone("same");
f2.DrawClone("same");
c.cd(2);
data.SetTitle("Data distribution with true contributions");
data.DrawClone("EP");
htruemc0.Draw("same");
htruemc1.Draw("same");
htruemc2.Draw("same");
c.cd(3);
mc0.SetTitle("MC generated samples with fit predictions");
mc0.Draw("PE");
mc1.Draw("PEsame");
mc2.Draw("PEsame");
if (status == 0) { // check on fit status
auto mcp0 = (TH1F *)fit.GetMCPrediction(0);
mcp0->SetLineColor(2);
mcp0->Draw("same");
auto mcp1 = (TH1F *)fit.GetMCPrediction(1);
mcp1->SetLineColor(3);
mcp1->Draw("same");
auto mcp2 = (TH1F *)fit.GetMCPrediction(2);
mcp2->SetLineColor(4);
mcp2->Draw("same");
}
c.cd(4);
l.SetTextSize(.035);
Char_t text[200];
if (status == 0) { // check on fit status
auto result = (TH1F *)fit.GetPlot();
fit.GetResult(0, p0, errP0);
printf(" Parameter %d: true %.3f, estim. %.3f +/- %.3f\n", 0, trueP0, p0, errP0);
fit.GetResult(1, p1, errP1);
printf(" Parameter %d: true %.3f, estim. %.3f +/- %.3f\n", 1, trueP1, p1, errP1);
fit.GetResult(2, p2, errP2);
printf(" Parameter %d: true %.3f, estim. %.3f +/- %.3f\n", 2, trueP2, p2, errP2);
data.SetTitle("Data distribution with fitted contributions");
data.DrawClone("Ep");
result->Draw("same");
f0.SetParameter(0, Ndata * p0 / int0 * data.GetBinWidth(1));
f0.SetLineStyle(2);
f0.DrawClone("same");
f1.SetParameter(0, Ndata * p1 / int1 * data.GetBinWidth(1));
f1.DrawClone("same");
f2.SetParameter(0, Ndata * p2 / int2 * data.GetBinWidth(1));
f2.SetLineStyle(2);
f2.DrawClone("same");
sprintf(text, "%d: true %.2f, estimated %.2f +/- %.2f\n", 0, trueP0, p0, errP0);
l.DrawTextNDC(.45, .30, text);
sprintf(text, "%d: true %.2f, estimated %.2f +/- %.2f\n", 1, trueP1, p1, errP1);
l.DrawTextNDC(.45, .25, text);
sprintf(text, "%d: true %.2f, estimated %.2f +/- %.2f\n", 2, trueP2, p2, errP2);
l.DrawTextNDC(.45, .20, text);
}
auto cnew = c.DrawClone();
// Cleanup
// delete cnew;
}
#define c(i)
Definition RSha256.hxx:101
int Int_t
Signed integer 4 bytes (int)
Definition RtypesCore.h:59
char Char_t
Character 1 byte (char)
Definition RtypesCore.h:51
double Double_t
Double 8 bytes.
Definition RtypesCore.h:73
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
winID h TVirtualViewer3D TVirtualGLPainter p
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void data
Option_t Option_t TPoint TPoint const char GetTextMagnitude GetFillStyle GetLineColor GetLineWidth GetMarkerStyle GetTextAlign GetTextColor GetTextSize void char Point_t Rectangle_t WindowAttributes_t Float_t Float_t Float_t Int_t Int_t UInt_t UInt_t Rectangle_t result
Option_t Option_t TPoint TPoint const char text
R__EXTERN TRandom * gRandom
Definition TRandom.h:62
R__EXTERN TStyle * gStyle
Definition TStyle.h:442
virtual void SetLineStyle(Style_t lstyle)
Set the line style.
Definition TAttLine.h:44
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:42
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:182
virtual Double_t GetRandom(TRandom *rng=nullptr, Option_t *opt=nullptr)
Return a random number following this function shape.
Definition TF1.cxx:2218
virtual Double_t Integral(Double_t a, Double_t b, Double_t epsrel=1.e-12)
IntegralOneDim or analytical integral.
Definition TF1.cxx:2557
virtual void SetParameter(Int_t param, Double_t value)
Definition TF1.h:623
Fits MC fractions to data histogram.
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:879
To draw Mathematical Formula.
Definition TLatex.h:18
An array of TObjects.
Definition TObjArray.h:31
virtual TObject * DrawClone(Option_t *option="") const
Draw a clone of this object in the current selected pad with: gROOT->SetSelectedPad(c1).
Definition TObject.cxx:319
virtual Double_t Uniform(Double_t x1=1)
Returns a uniform deviate on the interval (0, x1).
Definition TRandom.cxx:682
void SetOptStat(Int_t stat=1)
The type of information printed in the histogram statistics box can be selected via the parameter mod...
Definition TStyle.cxx:1642
Double_t x[n]
Definition legend1.C:17
TF1 * f1
Definition legend1.C:11
constexpr Double_t Pi()
Definition TMath.h:38
TLine l
Definition textangle.C:4
Author

Definition in file fitFraction.C.