Re: Storing multiple TH1::Fit functions with different fit ranges in same histo

From: Rene Brun <Rene.Brun_at_cern.ch>
Date: Fri, 2 Dec 2005 09:09:00 +0100 (MET)


Hi Remi,

See the simple solution in the code below.

Rene Brun

{
  TH1F *trigPt = new TH1F("trigPt", "Pt efficiency", 100, 0, 50);   trigPt->FillRandom("landau");

  const Float_t ptbins[5] = {1.5, 3, 5, 10, 50};   Float_t eff[4];
  Float_t effSigma[4];
  for (Int_t i = 0; i < 4; ++i) {

     trigPt->Fit("pol0","+","",ptbins[i],ptbins[i+1]);

    TF1 *f1 = (TF1*)trigPt->GetListOfFunctions()->At(i);     eff[i] = f1->GetParameter(0);
    effSigma[i] = f1->GetParError(0);
    std::cout << "Ptbin " << i << ": " << eff[i] << " +/- " << effSigma[i] <<std::endl;
  }
}

On Thu, 1 Dec 2005, Remi Mommsen
wrote:

> Dear ROOTers,
>
> I came across a problem trying to fit the same function multiple times to the
> same hist but using different fit ranges. I like to access the fitted
> function after each fit for a given range and store their parameters. This
> works well if I only store the last fitted function, i.e. I do not use the
> option "+" for TH1::Fit. If I use it, I always get the first stored function
> back. Is there a solution w/o using a different function in each iteration
> with an unique name?
>
>
> The following macro illustrates the problem:
>
> TH1F *trigPt = new TH1F("trigPt", "Pt efficiency", 100, 0, 50);
> trigPt->FillRandom("landau");
>
> const Float_t ptbins[5] = {1.5, 3, 5, 10, 50};
> Float_t eff[4];
> Float_t effSigma[4];
> for (Int_t i = 0; i < 4; ++i) {
> // case 1: this gets the correct function
> trigPt->Fit("pol0","","",ptbins[i],ptbins[i+1]);
>
> // case 2: this does get always the first pol0 function
> // trigPt->Fit("pol0","+","",ptbins[i],ptbins[i+1]);
>
> eff[i] = trigPt->GetFunction("pol0")->GetParameter(0);
> effSigma[i] = trigPt->GetFunction("pol0")->GetParError(0);
> std::cout << "Ptbin " << i << ": " << eff[i] << " +/- " << effSigma[i] <<
> std::endl;
> }
>
> Output in case 1:
> root [] .x rootbug.C
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 342.285725 10.681537
> Ptbin 0: 342.286 +/- 10.6815
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 172.255836 6.562314
> Ptbin 1: 172.256 +/- 6.56231
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 69.657010 2.639261
> Ptbin 2: 69.657 +/- 2.63926
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 3.242070 0.205195
> Ptbin 3: 3.24207 +/- 0.205195
>
>
> Output in case 2:
> note that the cout line always prints the values for function 1
> root [] .x rootbug.C
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 355.334138 10.883231
> Ptbin 0: 355.334 +/- 10.8832
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 193.946272 6.963230
> Ptbin 1: 355.334 +/- 10.8832
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 67.710358 2.602121
> Ptbin 2: 355.334 +/- 10.8832
> Fitting results:
> Parameters:
> NO. VALUE ERROR
> 0 3.305274 0.208544
> Ptbin 3: 355.334 +/- 10.8832
>
>
> If it matters: I'm using root 5.06/00 on Mac OS X 10.4.3
>
>
> Cheers,
> Remi
>
>
>
> --
> Truth decays into beauty, while beauty soon becomes merely charm.
> Charm ends up as strangeness, and even that doesn't last, but up and
> down are forever. (Anonymous)
>
> *********************************************************************
> Remigius K. Mommsen e-mail: mommsen_at_fnal.gov
> University of Manchester URL: http://cern.ch/mommsen
> Fermilab, MS 357 voice: ++1 (630) 840-8321
> P.O. Box 500 fax: ++1 (630) 840-2649
> Batavia, Il 60510, US home: ++1 (630) 236-0932
> *********************************************************************
>
>
Received on Fri Dec 02 2005 - 09:09:06 MET

This archive was generated by hypermail 2.2.0 : Tue Jan 02 2007 - 14:45:14 MET