[ROOT] Problem retrieving fitparameters from self-defined functions...

From: Marc-Andre Pleier (pleier@fnal.gov)
Date: Sat Nov 22 2003 - 20:35:45 MET 

Dear Rooters,

I encounter the following problem (see attached sourcecode example
for reproducing this...): 

Using a fitfunction I created myself which is used in a subroutine
for a TF1 instantiation a la TF1 * myfit = new TF1("myfit",function,...);
I am not able to extract the fit parameters that are printed on the
screen by Minuit, if this subroutine is called more than once.

For example fitting the same histogram more than once, I get the same
Minuit output on the screen, but when I try to get the fit parameters
from the second call on, the values returned are the parameter initial 
values instead of the Minuit values.
This problem disappears when I call the TF1 destructor at the end
of the subroutine.

If I use the same subroutine with a predefined function like gaus,
this problem doesn't show up.
But if I don't put the TF1 destructor call in that routine as well,
any call to my self-defined function after will again only
return the parameter initial values instead of the Minuit values.

Could somebody explain me why ROOT is behaving that way?

To reproduce this effect, you can simply compile the attached code
and switch the destructor calls on and off.

I am using ROOT Version 3.05/00  23 February 2003.

Cheers,
Marc-Andre.

-----------------------------------------------------------------------------
  Marc-Andre Pleier                            Phone : +1 630 840 8575
  Fermilab D0                                   Fax  : +1 630 840 8886
  P.O. Box 500, MS 352                          Email: pleier@fnal.gov
  Batavia, IL 60510-0500, USA
-----------------------------------------------------------------------------


Here is my code:

#include <TF1.h>
#include <TH1.h>
#include <TFile.h>
#include <iostream>
using namespace std;

double mygauss(Double_t *x, Double_t *par)
{
  Double_t fitval = par[0]*exp(-.5*((x[0]-par[1])/par[2])*((x[0]-par[1])/par[2]));
  return fitval;
}

void fitmegood(TH1D * loose)
{
  TF1 * myfit = new TF1("myfit", "gaus",50.,130.);
  myfit->SetParameters(7.,90.,4.);
  loose->Fit("myfit");
  cout << "fitpars : " << myfit->GetParameter(1) << " | " << myfit->GetParameter(2) << endl;
  // if we want fitmebad to work for all calls in the second block...
  //myfit->~TF1();
  return;
}

void fitmebad(TH1D * loose)
{
  TF1 * myfit = new TF1("myfit", mygauss,50., 130., 3);
  myfit->SetParameters(7.,90.,4.);
  loose->Fit("myfit");
  cout << "fitpars : " << myfit->GetParameter(1) << " | " << myfit->GetParameter(2) << endl;
  // if we want fitmebad to work for second++ calls in the first block...
  //myfit->~TF1();
  return;
}

int main(void)
{
  // Create a plot to be fitted...
  TF1 * willy = new TF1("willy", "gaus", 40., 140.);
  willy->SetParameters(7.42, 92., 3.5);
  TH1D* loose = new TH1D("willyplot", "test the fitting",100, 40., 140.);
  loose->FillRandom("willy", 150.);

  cout << "Fit using self-defined gaussian, block 1:" << endl;
  fitmebad(loose);
  fitmebad(loose);
  fitmebad(loose);

  cout << "Fit using predefined gaussian:" << endl;
  fitmegood(loose);
  fitmegood(loose);
  fitmegood(loose);

  cout << "Fit using self-defined gaussian, block 2:" << endl;
  fitmebad(loose);
  fitmebad(loose);
  fitmebad(loose);

  return 0;
}

This archive was generated by hypermail 2b29 : Thu Jan 01 2004 - 17:50:17 MET