Combined (simultaneous) fit of two histogram with separate functions and some common parameters
See http://root.cern/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
Author: Lorenzo Moneta
This notebook tutorial was automatically generated with ROOTBOOK-izer from the macro found in the ROOT repository on Tuesday, March 19, 2024 at 07:07 PM.
%%cpp -d
#include <Fit/Fitter.h>
#include <Fit/BinData.h>
#include <Fit/Chi2FCN.h>
#include <TH1.h>
#include <Math/WrappedMultiTF1.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
};
%%cpp -d
input_line_44:6:1: error: expected expression %%cpp -d ^ input_line_44:6:2: error: expected expression %%cpp -d ^ input_line_44:6:3: error: use of undeclared identifier 'cpp' %%cpp -d ^ input_line_44:6:8: error: use of undeclared identifier 'd' %%cpp -d ^
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
};
%%cpp -d
input_line_53:9:1: error: expected expression %%cpp -d ^ input_line_53:9:2: error: expected expression %%cpp -d ^ input_line_53:9:3: error: use of undeclared identifier 'cpp' %%cpp -d ^ input_line_53:9:8: error: use of undeclared identifier 'd' %%cpp -d ^
Create the GlobalCHi2 structure
struct GlobalChi2 {
GlobalChi2(ROOT::Math::IMultiGenFunction &f1, ROOT::Math::IMultiGenFunction &f2) : 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);
}
const ROOT::Math::IMultiGenFunction *fChi2_1;
const ROOT::Math::IMultiGenFunction *fChi2_2;
};
input_line_54:10:22: error: use of undeclared identifier 'iparB' p1[i] = par[iparB[i]]; ^ input_line_54:14:22: error: use of undeclared identifier 'iparSB' p2[i] = par[iparSB[i]]; ^
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);
ROOT::Math::WrappedMultiTF1 wfB(*fB, 1);
ROOT::Math::WrappedMultiTF1 wfSB(*fSB, 1);
ROOT::Fit::DataOptions opt;
ROOT::Fit::DataRange rangeB;
set the data range
rangeB.SetRange(10, 90);
ROOT::Fit::BinData dataB(opt, rangeB);
ROOT::Fit::FillData(dataB, hB);
ROOT::Fit::DataRange rangeSB;
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);
ROOT::Fit::Fitter fitter;
const int Npar = 6;
double par0[Npar] = {5, 5, -0.1, 100, 30, 10};
input_line_66:14:1: error: unknown type name 'GlobalChi2' GlobalChi2 globalChi2(chi2_B, chi2_SB); ^
create before the parameter settings in order to fix or set range on them
fitter.Config().SetParamsSettings(6, par0);
input_line_68:2:3: error: use of undeclared identifier 'fitter' (fitter.Config().SetParamsSettings(6, par0)) ^ input_line_68:2:40: error: use of undeclared identifier 'par0' (fitter.Config().SetParamsSettings(6, par0)) ^ Error in <HandleInterpreterException>: Error evaluating expression (fitter.Config().SetParamsSettings(6, par0)) Execution of your code was aborted.
fix 5-th parameter
fitter.Config().ParSettings(4).Fix();
input_line_70:2:3: error: use of undeclared identifier 'fitter' (fitter.Config().ParSettings(4).Fix()) ^ Error in <HandleInterpreterException>: Error evaluating expression (fitter.Config().ParSettings(4).Fix()) Execution of your code was aborted.
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().MinimizerOptions().SetPrintLevel(0);
fitter.Config().SetMinimizer("Minuit2", "Migrad");
input_line_72:2:3: error: use of undeclared identifier 'fitter' (fitter.Config().ParSettings(2).SetLimits(-10, -1.0E-4)) ^ Error in <HandleInterpreterException>: Error evaluating expression (fitter.Config().ParSettings(2).SetLimits(-10, -1.0E-4)) Execution of your code was aborted.
fit FCN function directly (specify optionally data size and flag to indicate that is a chi2 fit)
fitter.FitFCN(6, globalChi2, nullptr, dataB.Size() + dataSB.Size(), true);
ROOT::Fit::FitResult result = fitter.Result();
result.Print(std::cout);
std::cout << "Combined fit Chi2 = " << result.Chi2() << std::endl;
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);
fB->SetLineColor(kBlue);
hB->GetListOfFunctions()->Add(fB);
hB->Draw();
c1->cd(2);
fSB->SetFitResult(result, iparSB);
fSB->SetRange(rangeSB().first, rangeSB().second);
fSB->SetLineColor(kRed);
hSB->GetListOfFunctions()->Add(fSB);
hSB->Draw();
[runStaticInitializersOnce]: Failed to materialize symbols: { (main, { _ZNK4ROOT3Fit9FitResult4Chi2Ev, _Z31__fd_init_order__cling_Un1Qu312v, __vd_init_order__cling_Un1Qu313, __cxx_global_var_initcling_module_346_.6, _GLOBAL__sub_I_cling_module_346, _ZN12__cling_N52916__cling_Un1Qu330EPv, __cxx_global_var_initcling_module_346_.7, __vd_init_order__cling_Un1Qu311, _Z31__fd_init_order__cling_Un1Qu310v, _ZN12__cling_N5296resultE, __orc_init_func.cling-module-346, __cxx_global_var_initcling_module_346_.4, __cxx_global_var_initcling_module_346_, _ZN12__cling_N5292c1E, __cxx_global_var_initcling_module_346_.1, _ZNK3TH118GetListOfFunctionsEv, $.cling-module-346.__inits.0, _ZN12__cling_N52924__dynamic__cling_Un1Qu30E, _ZNK5cling7runtime8internal15LifetimeHandler9getMemoryEv }) } IncrementalExecutor::executeFunction: symbol '_ZN5cling7runtime8internal15LifetimeHandlerC1EPNS1_15DynamicExprInfoEPN5clang11DeclContextEPKcPNS_11InterpreterE' unresolved while linking [cling interface function]! You are probably missing the definition of cling::runtime::internal::LifetimeHandler::LifetimeHandler(cling::runtime::internal::DynamicExprInfo*, clang::DeclContext*, char const*, cling::Interpreter*) Maybe you need to load the corresponding shared library? IncrementalExecutor::executeFunction: symbol '_ZN5cling7runtime8internal15LifetimeHandlerD1Ev' unresolved while linking [cling interface function]! You are probably missing the definition of cling::runtime::internal::LifetimeHandler::~LifetimeHandler() Maybe you need to load the corresponding shared library?
Draw all canvases
gROOT->GetListOfCanvases()->Draw()
[runStaticInitializersOnce]: Failed to materialize symbols: { (main, { __orc_init_func.cling-module-346 }) }