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

Detailed Description

View in nbviewer Open in SWAN Combined (simultaneous) fit of two histogram with separate functions and some common parameters

See http://root.cern.ch/phpBB3//viewtopic.php?f=3&t=11740#p50908 for a modified version working with Fumili or GSLMultiFit

N.B. this macro must be compiled with ACliC

****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 131.104
NDf = 115
Edm = 2.11602e-08
NCalls = 225
Par_0 = 5.5396 +/- 0.0354094
Par_1 = 4.66089 +/- 0.050106
Par_2 = -0.0514037 +/- 0.00108539 (limited)
Par_3 = 77.2733 +/- 3.93105 (limited)
Par_4 = 30 (fixed)
Par_5 = 4.864 +/- 0.243005
#include "Fit/Fitter.h"
#include "Fit/BinData.h"
#include "Fit/Chi2FCN.h"
#include "TH1.h"
#include "TList.h"
#include "HFitInterface.h"
#include "TCanvas.h"
#include "TStyle.h"
// definition of shared parameter
// background function
int iparB[2] = { 0, // exp amplitude in B histo
2 // exp common parameter
};
// signal + background function
int iparSB[5] = { 1, // exp amplitude in S+B histo
2, // exp common parameter
3, // gaussian amplitude
4, // gaussian mean
5 // gaussian sigma
};
// Create the GlobalCHi2 structure
struct GlobalChi2 {
fChi2_1(&f1), fChi2_2(&f2) {}
// parameter vector is first background (in common 1 and 2)
// and then is signal (only in 2)
double operator() (const double *par) const {
double p1[2];
for (int i = 0; i < 2; ++i) p1[i] = par[iparB[i] ];
double p2[5];
for (int i = 0; i < 5; ++i) p2[i] = par[iparSB[i] ];
return (*fChi2_1)(p1) + (*fChi2_2)(p2);
}
};
void combinedFit() {
TH1D * hB = new TH1D("hB","histo B",100,0,100);
TH1D * hSB = new TH1D("hSB","histo S+B",100, 0,100);
TF1 * fB = new TF1("fB","expo",0,100);
fB->SetParameters(1,-0.05);
hB->FillRandom("fB");
TF1 * fS = new TF1("fS","gaus",0,100);
fS->SetParameters(1,30,5);
hSB->FillRandom("fB",2000);
hSB->FillRandom("fS",1000);
// perform now global fit
TF1 * fSB = new TF1("fSB","expo + gaus(2)",0,100);
// set the data range
rangeB.SetRange(10,90);
ROOT::Fit::BinData dataB(opt,rangeB);
ROOT::Fit::FillData(dataB, hB);
rangeSB.SetRange(10,50);
ROOT::Fit::BinData dataSB(opt,rangeSB);
ROOT::Fit::FillData(dataSB, hSB);
ROOT::Fit::Chi2Function chi2_B(dataB, wfB);
ROOT::Fit::Chi2Function chi2_SB(dataSB, wfSB);
GlobalChi2 globalChi2(chi2_B, chi2_SB);
const int Npar = 6;
double par0[Npar] = { 5,5,-0.1,100, 30,10};
// create before the parameter settings in order to fix or set range on them
fitter.Config().SetParamsSettings(6,par0);
// fix 5-th parameter
fitter.Config().ParSettings(4).Fix();
// set limits on the third and 4-th parameter
fitter.Config().ParSettings(2).SetLimits(-10,-1.E-4);
fitter.Config().ParSettings(3).SetLimits(0,10000);
fitter.Config().ParSettings(3).SetStepSize(5);
fitter.Config().SetMinimizer("Minuit2","Migrad");
// fit FCN function directly
// (specify optionally data size and flag to indicate that is a chi2 fit)
fitter.FitFCN(6,globalChi2,0,dataB.Size()+dataSB.Size(),true);
ROOT::Fit::FitResult result = fitter.Result();
result.Print(std::cout);
TCanvas * c1 = new TCanvas("Simfit","Simultaneous fit of two histograms",
10,10,700,700);
c1->Divide(1,2);
c1->cd(1);
gStyle->SetOptFit(1111);
fB->SetFitResult( result, iparB);
fB->SetRange(rangeB().first, rangeB().second);
hB->GetListOfFunctions()->Add(fB);
hB->Draw();
c1->cd(2);
fSB->SetFitResult( result, iparSB);
fSB->SetRange(rangeSB().first, rangeSB().second);
hSB->GetListOfFunctions()->Add(fSB);
hSB->Draw();
}
@ kRed
Definition Rtypes.h:66
@ kBlue
Definition Rtypes.h:66
TRObject operator()(const T1 &t1) const
R__EXTERN TStyle * gStyle
Definition TStyle.h:413
Class describing the binned data sets : vectors of x coordinates, y values and optionally error on y ...
Definition BinData.h:52
Chi2FCN class for binnned fits using the least square methods.
Definition Chi2FCN.h:46
class describing the range in the coordinates it supports multiple range in a coordinate.
Definition DataRange.h:35
void SetRange(unsigned int icoord, double xmin, double xmax)
set a range [xmin,xmax] for the new coordinate icoord If more range exists for other coordinates,...
void SetParamsSettings(unsigned int npar, const double *params, const double *vstep=0)
set the parameter settings from number of parameters and a vector of values and optionally step value...
void SetMinimizer(const char *type, const char *algo=0)
set minimizer type
Definition FitConfig.h:181
const ParameterSettings & ParSettings(unsigned int i) const
get the parameter settings for the i-th parameter (const method)
Definition FitConfig.h:76
ROOT::Math::MinimizerOptions & MinimizerOptions()
access to the minimizer control parameter (non const method)
Definition FitConfig.h:167
class containg the result of the fit and all the related information (fitted parameter values,...
Definition FitResult.h:47
void Print(std::ostream &os, bool covmat=false) const
print the result and optionaly covariance matrix and correlations
Fitter class, entry point for performing all type of fits.
Definition Fitter.h:77
bool FitFCN(unsigned int npar, Function &fcn, const double *params=0, unsigned int dataSize=0, bool chi2fit=false)
Fit using the a generic FCN function as a C++ callable object implementing double () (const double *)...
Definition Fitter.h:610
const FitResult & Result() const
get fit result
Definition Fitter.h:384
const FitConfig & Config() const
access to the fit configuration (const method)
Definition Fitter.h:412
void SetStepSize(double err)
set the step size
void SetLimits(double low, double up)
set a double side limit, if low == up the parameter is fixed if low > up the limits are removed The c...
void Fix()
fix the parameter
Documentation for the abstract class IBaseFunctionMultiDim.
Definition IFunction.h:62
void SetPrintLevel(int level)
set print level
Class to Wrap a ROOT Function class (like TF1) in a IParamMultiFunction interface of multi-dimensions...
virtual void SetLineColor(Color_t lcolor)
Set the line color.
Definition TAttLine.h:40
The Canvas class.
Definition TCanvas.h:23
1-Dim function class
Definition TF1.h:213
virtual void SetFitResult(const ROOT::Fit::FitResult &result, const Int_t *indpar=0)
Set the result from the fit parameter values, errors, chi2, etc... Optionally a pointer to a vector (...
Definition TF1.cxx:3375
virtual void SetRange(Double_t xmin, Double_t xmax)
Initialize the upper and lower bounds to draw the function.
Definition TF1.cxx:3539
virtual void SetParameters(const Double_t *params)
Definition TF1.h:644
1-D histogram with a double per channel (see TH1 documentation)}
Definition TH1.h:618
virtual void FillRandom(const char *fname, Int_t ntimes=5000, TRandom *rng=nullptr)
Fill histogram following distribution in function fname.
Definition TH1.cxx:3526
TList * GetListOfFunctions() const
Definition TH1.h:243
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition TH1.cxx:3074
virtual void Add(TObject *obj)
Definition TList.h:81
void SetOptFit(Int_t fit=1)
The type of information about fit parameters printed in the histogram statistics box can be selected ...
Definition TStyle.cxx:1541
return c1
Definition legend1.C:41
TF1 * f1
Definition legend1.C:11
void FillData(BinData &dv, const TH1 *hist, TF1 *func=0)
fill the data vector from a TH1.
Definition first.py:1
DataOptions : simple structure holding the options on how the data are filled.
Definition DataOptions.h:28
Author
Lorenzo Moneta

Definition in file combinedFit.C.