Hi
could someone please tell me if there is a way to simultaneously fit multiple functions to multiple graphs?
To be more explicit I have several functions f_i which have some common parameters shared amongst them and I would like to fit these to different graphs, is there a way to fit these? I am not an expert at C so if you could try to be explicit in how I would implement this that would be much appreciated.
thanks
Mark
So I have an updated status to this post. I have managed to perform a simultaneous fit for multiple graphs,
the issue I am now facing is that the fitted parameters have very large errors even though the fit seems to be very good, and I would like to know how to resolve this issue. The way I have performed this fit is quite simple it is a modified version of the tutorial combinedFit.C it is as follows:
//+ 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
//
//Author: L. Moneta - Dec 2010
#include "Fit/Fitter.h"
#include "Fit/BinData.h"
#include "Fit/Chi2FCN.h"
#include "TH1.h"
#include "TF1.h"
#include "TList.h"
#include "Math/WrappedMultiTF1.h"
#include "HFitInterface.h"
#include "TCanvas.h"
#include "TStyle.h"
#include "TMath.h"
#include "TGraph.h"
#include "TRandom.h"
// definition of shared parameter
// background function
int iparB1[2] = { 0, // exp amplitude in B histo
2 // exp common parameter
};
// signal + background function
int iparB2[2] = { 1, // exp amplitude in S+B histo
2, // exp common parameter
};
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[iparB1[i] ];
double p2[2];
for (int i = 0; i < 2; ++i) p2[i] = par[iparB2[i] ];
return (*fChi2_1)(p1) + (*fChi2_2)(p2);
}
const ROOT::Math::IMultiGenFunction * fChi2_1;
const ROOT::Math::IMultiGenFunction * fChi2_2;
};
double myfunc1(double *x, double *par){
float xx =x[0];
double f = TMath::Exp(par[0]+par[1]*xx);
return f;
}
double myfunc2(double *x, double *par){
float xx =x[0];
double f = TMath::Exp(par[0]+par[1]*xx);
return f;
}
double myfunction1(double *x, double *par){
double constant = 2.0;
float xx =x[0];
double f = TMath::Exp(constant*par[0]+par[1]*xx);
return f;
}
double myfunction2(double *x, double *par){
float constant = 2.0;
float xx =x[0];
double f = TMath::Exp(constant*par[0]+par[1]*xx);
return f;
}
void combinedFit(){
int nd = 200;
double rnd1,x, rnd2,xp[nd],y1[nd],y2[nd];
double e = 0.1;
TRandom r1,r2;
x = (double) 100.0/nd;
//TH1D * hB1 = new TH1D("hB1","histo B1",100,0,100);
//TH1D * hB2 = new TH1D("hB2","histo B2",100, 0,100);
TF1 * fB1 = new TF1("fB1",myfunc1,0,100,2);
TF1 * fB2 = new TF1("fB2",myfunc2,0,100,2);
fB1->SetParameters(1,-0.05);
fB2->SetParameters(1,-0.05);
//hB1->FillRandom("fB1",2000);
//hB2->FillRandom("fB1",2000);
for (Int_t i=0; i<nd; i++) {
rnd1 = r1.Uniform(-e,e); // Generate a random number in [-e,e]
rnd2 = r2.Uniform(-e,e); // Generate a random number in [-e,e]
y1[i] = fB1->Eval(x*i)*(1+rnd1);
y2[i] = fB2->Eval(x*i)*(1+rnd2);
xp[i] = x*i;
//cout<< rnd2<<" "<<xp[i] <<" " << y2[i]<<endl;
}
TGraph *gb1 = new TGraph(nd,xp,y1);
TGraph *gb2 = new TGraph(nd,xp,y2);
TF1 * ffitB1 = new TF1("ffitB1",myfunction1,0,100,2);
TF1 * ffitB2 = new TF1("ffitB2",myfunction2,0,100,2);
ROOT::Math::WrappedMultiTF1 wfB1(*ffitB1,1);
ROOT::Math::WrappedMultiTF1 wfB2(*ffitB2,1);
ROOT::Fit::DataOptions opt;
ROOT::Fit::DataRange rangeB1;
ROOT::Fit::DataRange rangeB2;
// set the data range
rangeB1.SetRange(0,100);
rangeB2.SetRange(0,100);
ROOT::Fit::BinData dataB1(opt,rangeB1);
ROOT::Fit::BinData dataB2(opt,rangeB2);
ROOT::Fit::FillData(dataB1, gb1);
ROOT::Fit::FillData(dataB2, gb2);
ROOT::Fit::Chi2Function chi2_B1(dataB1, wfB1);
ROOT::Fit::Chi2Function chi2_B2(dataB2, wfB2);
GlobalChi2 globalChi2(chi2_B1, chi2_B2);
ROOT::Fit::Fitter fitter;
const int Npar = 3;
double par0[Npar] = {50,50,-0.1};
// create before the parameter settings in order to fix or set range on them
fitter.Config().SetParamsSettings(3,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(2).SetLimits(0,10000);
//fitter.Config().ParSettings(3).SetStepSize(5);
fitter.Config().MinimizerOptions().SetPrintLevel(0);
fitter.Config().SetMinimizer("Minuit","Minimize");
// fit FCN function directly
// (specify optionally data size and flag to indicate that is a chi2 fit)
fitter.FitFCN(3,globalChi2,0,dataB1.Size()+dataB2.Size(),true);
ROOT::Fit::FitResult result = fitter.Result();
result.Print(std::cout);
TCanvas * c1 = new TCanvas("Simfit","Simultaneous fit of two graphs",
10,10,700,700);
c1->Divide(1,2);
c1->cd(1);
gStyle->SetOptFit(1111);
ffitB1->SetFitResult(result, iparB1);
ffitB1->SetRange(rangeB1().first, rangeB1().second);
ffitB1->SetLineColor(kBlue);
gb1->GetListOfFunctions()->Add(ffitB1);
gb1->Draw("alp");
//
c1->cd(2);
ffitB2->SetFitResult( result, iparB2);
ffitB2->SetRange(rangeB2().first, rangeB2().second);
ffitB2->SetLineColor(kRed);
gb2->GetListOfFunctions()->Add(ffitB2);
gb2->Draw("alp");
}This is essentially how the fit is performed in my actual code though the functions I am fitting are slightly more complicated unfortunately this piece of code does not have the problem of large errors I can post my actual code but it is more complicated and not very well written…
Is there a reason why the errors minuit give are large and is there any way of minimizing them?