ROOT logo

From $ROOTSYS/tutorials/fit/fitExclude.C

// Illustrate how to fit excluding points in a given range
// Author: Rene Brun
#include "TH1.h"
#include "TF1.h"
#include "TList.h"

Bool_t reject;
Double_t fline(Double_t *x, Double_t *par)
{
    if (reject && x[0] > 2.5 && x[0] < 3.5) {
      TF1::RejectPoint();
      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 = kTRUE;
   h->Fit(fl,"0");
   reject = kFALSE;
   //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();
}
   
 fitExclude.C:1
 fitExclude.C:2
 fitExclude.C:3
 fitExclude.C:4
 fitExclude.C:5
 fitExclude.C:6
 fitExclude.C:7
 fitExclude.C:8
 fitExclude.C:9
 fitExclude.C:10
 fitExclude.C:11
 fitExclude.C:12
 fitExclude.C:13
 fitExclude.C:14
 fitExclude.C:15
 fitExclude.C:16
 fitExclude.C:17
 fitExclude.C:18
 fitExclude.C:19
 fitExclude.C:20
 fitExclude.C:21
 fitExclude.C:22
 fitExclude.C:23
 fitExclude.C:24
 fitExclude.C:25
 fitExclude.C:26
 fitExclude.C:27
 fitExclude.C:28
 fitExclude.C:29
 fitExclude.C:30
 fitExclude.C:31
 fitExclude.C:32
 fitExclude.C:33
 fitExclude.C:34
 fitExclude.C:35
 fitExclude.C:36
 fitExclude.C:37
 fitExclude.C:38
 fitExclude.C:39
 fitExclude.C:40
 fitExclude.C:41
 fitExclude.C:42