Logo ROOT  
Reference Guide
 
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
Loading...
Searching...
No Matches
fitExclude.C File Reference

Detailed Description

View in nbviewer Open in SWAN
Illustrates how to fit excluding points in a given range.

****************************************
Minimizer is Minuit2 / Migrad
Chi2 = 84.1529
NDf = 78
Edm = 8.41446e-22
NCalls = 32
p0 = 28.4731 +/- 0.946564
p1 = -4.81581 +/- 0.26533
#include <TH1.h>
#include <TF1.h>
#include <TROOT.h>
bool reject;
double fline(double *x, double *par)
{
if (reject && x[0] > 2.5 && x[0] < 3.5) {
return 0;
}
return par[0] + par[1]*x[0];
}
void fitExclude() {
//Create a source function
TF1 *f1 = new TF1("f1","[0] +[1]*x +gaus(2)",0,5);
f1->SetParameters(6,-1,5,3,0.2);
// create and fill histogram according to the source function
TH1F *h = new TH1F("h","background + signal",100,0,5);
h->FillRandom("f1",2000);
TF1 *fl = new TF1("fl",fline,0,5,2);
fl->SetParameters(2,-1);
//fit only the linear background excluding the signal area
reject = true;
h->Fit(fl,"0");
reject = false;
//store 2 separate functions for visualization
TF1 *fleft = new TF1("fleft",fline,0,2.5,2);
fleft->SetParameters(fl->GetParameters());
h->GetListOfFunctions()->Add(fleft);
gROOT->GetListOfFunctions()->Remove(fleft);
TF1 *fright = new TF1("fright",fline,3.5,5,2);
fright->SetParameters(fl->GetParameters());
h->GetListOfFunctions()->Add(fright);
gROOT->GetListOfFunctions()->Remove(fright);
h->Draw();
}
#define h(i)
Definition RSha256.hxx:106
ROOT::Detail::TRangeCast< T, true > TRangeDynCast
TRangeDynCast is an adapter class that allows the typed iteration through a TCollection.
#define gROOT
Definition TROOT.h:406
1-Dim function class
Definition TF1.h:233
static void RejectPoint(Bool_t reject=kTRUE)
Static function to set the global flag to reject points the fgRejectPoint global flag is tested by al...
Definition TF1.cxx:3684
virtual Double_t * GetParameters() const
Definition TF1.h:548
virtual void SetParameters(const Double_t *params)
Definition TF1.h:677
1-D histogram with a float per channel (see TH1 documentation)
Definition TH1.h:645
Double_t x[n]
Definition legend1.C:17
TF1 * f1
Definition legend1.C:11
Author
Rene Brun

Definition in file fitExclude.C.