]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
-rm obsolete fit functionality
authorjbook <jbook@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 1 Nov 2012 13:28:49 +0000 (13:28 +0000)
committerjbook <jbook@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 1 Nov 2012 13:28:49 +0000 (13:28 +0000)
-add easy way to combine a peak + any bgrd function

PWGDQ/dielectron/AliDielectronSignalFunc.cxx
PWGDQ/dielectron/AliDielectronSignalFunc.h

index cfa9cbc8a1e62f7927f2ea6765f4c46baad4368b..0ab9e0957bcc59b0b0ffe73e0982f74c6791cc05 100644 (file)
@@ -45,17 +45,21 @@ can be used.
 ClassImp(AliDielectronSignalFunc)
 
 TH1F* AliDielectronSignalFunc::fgHistSimPM=0x0;
 ClassImp(AliDielectronSignalFunc)
 
 TH1F* AliDielectronSignalFunc::fgHistSimPM=0x0;
+TF1*  AliDielectronSignalFunc::fFuncSignal=0x0;
+TF1*  AliDielectronSignalFunc::fFuncBackground=0x0;
+Int_t AliDielectronSignalFunc::fNparPeak=0;
+Int_t AliDielectronSignalFunc::fNparBgnd=0;
+
 
 AliDielectronSignalFunc::AliDielectronSignalFunc() :
 AliDielectronSignalBase(),
 
 AliDielectronSignalFunc::AliDielectronSignalFunc() :
 AliDielectronSignalBase(),
-fFuncSignal(0x0),
-fFuncBackground(0x0),
 fFuncSigBack(0x0),
 fParMass(1),
 fParMassWidth(2),
 fFitOpt("SMNQE"),
 fUseIntegral(kFALSE),
 fPolDeg(0),
 fFuncSigBack(0x0),
 fParMass(1),
 fParMassWidth(2),
 fFitOpt("SMNQE"),
 fUseIntegral(kFALSE),
 fPolDeg(0),
+fDof(0),
 fChi2Dof(0.0)
 {
   //
 fChi2Dof(0.0)
 {
   //
@@ -67,14 +71,13 @@ fChi2Dof(0.0)
 //______________________________________________
 AliDielectronSignalFunc::AliDielectronSignalFunc(const char* name, const char* title) :
 AliDielectronSignalBase(name, title),
 //______________________________________________
 AliDielectronSignalFunc::AliDielectronSignalFunc(const char* name, const char* title) :
 AliDielectronSignalBase(name, title),
-fFuncSignal(0x0),
-fFuncBackground(0x0),
 fFuncSigBack(0x0),
 fParMass(1),
 fParMassWidth(2),
 fFitOpt("SMNQE"),
 fUseIntegral(kFALSE),
 fPolDeg(0),
 fFuncSigBack(0x0),
 fParMass(1),
 fParMassWidth(2),
 fFitOpt("SMNQE"),
 fUseIntegral(kFALSE),
 fPolDeg(0),
+fDof(0),
 fChi2Dof(0.0)
 {
   //
 fChi2Dof(0.0)
 {
   //
@@ -88,8 +91,6 @@ AliDielectronSignalFunc::~AliDielectronSignalFunc()
   //
   // Default Destructor
   //
   //
   // Default Destructor
   //
-  if(fFuncSignal) delete fFuncSignal;
-  if(fFuncBackground) delete fFuncBackground;
   if(fFuncSigBack) delete fFuncSigBack;
 }
 
   if(fFuncSigBack) delete fFuncSigBack;
 }
 
@@ -101,10 +102,6 @@ void AliDielectronSignalFunc::Process(TObjArray * const arrhist)
   // Fit the invariant mass histograms and retrieve the signal and background
   //
   switch(fMethod) {
   // Fit the invariant mass histograms and retrieve the signal and background
   //
   switch(fMethod) {
-  case kFittedMC :
-       ProcessFitIKF(arrhist);
-       break;
-                 
   case kFitted :
     ProcessFit(arrhist);
     break;
   case kFitted :
     ProcessFit(arrhist);
     break;
@@ -121,237 +118,8 @@ void AliDielectronSignalFunc::Process(TObjArray * const arrhist)
     AliError("Background substraction method not known!");
   }
 }
     AliError("Background substraction method not known!");
   }
 }
-
-//______________________________________________
-void AliDielectronSignalFunc::ProcessFitIKF(TObjArray * const arrhist) {
-  //
-  // Fit the +- invariant mass distribution only
-  // 
-  //
-  
-  const Double_t bigNumber = 100000.;
-  Double_t chi2ndfm1 = bigNumber;
-  Double_t ratiom1   = bigNumber;
-  Double_t chi2ndf   = bigNumber;
-  Int_t nDP =0;
-  
-  Int_t maxPolDeg = 8;
-      
-  fHistDataPM = (TH1F*)(arrhist->At(1))->Clone("histPM");  // +-    SE
-    if(fRebin>1) fHistDataPM->Rebin(fRebin);
-  
-  fgHistSimPM = (TH1F*)(arrhist->At(3))->Clone("histPMsim");  // +- mc shape
-  if (!fgHistSimPM) {
-       AliFatal("No mc peak shape found at idx 3.");
-       return;
-  }
-  if(fRebin>1) fgHistSimPM->Rebin(fRebin);
-  
-  // try out the polynomial degrees
-  for (Int_t iPD=0; iPD<=maxPolDeg; iPD++) {
-    TH1F *hf1 = (TH1F *) fHistDataPM->Clone(Form("hf1_PD%d",iPD));
-       
-       FitOneMinv(hf1, fgHistSimPM, iPD);
-       if (fChi2Dof > 0) chi2ndf = fChi2Dof;
-    AliInfo(Form("nDP: %d, iPD: %d, chi2ndf: %f", nDP, iPD, chi2ndf));
-       
-    ratiom1 = TMath::Abs(fChi2Dof - 1);
-    if (chi2ndfm1 > ratiom1) { // search for the closest to 1.
-      chi2ndfm1 = ratiom1;
-      nDP       = iPD;
-    }
-  }
-    
-  
-  // fit again with the best polynomial degree
-  TH1F *h2 = (TH1F *) fHistDataPM->Clone(Form("h2_PD%d",nDP));
-  
-  FitOneMinv(h2, fgHistSimPM, nDP);
-  AliInfo(Form("Best Fit: PD %d, chi^2/ndf %.3f, S/B %.3f",nDP,fChi2Dof,fValues(3)));
-    
-}
-//______________________________________________
-void AliDielectronSignalFunc::FitOneMinv(TH1F *hMinv, TH1F *hSim, Int_t pod) {
-  //
-  // main function to fit an inv mass spectrum
-  //
-  
-  TObjArray *arrResults = new TObjArray;
-  arrResults->SetOwner();
-  arrResults->AddAt(hMinv,0);
-  
-  // Degree of polynomial
-  fPolDeg = pod;
-    
-  // inclusion and exclusion areas (values)
-  const Double_t kJPMass    = 3.096916;
-  // inclusion and exclusion areas (bin numbers)
-  const Int_t binIntLo = hMinv->FindBin(fIntMin);
-  const Int_t binIntHi = hMinv->FindBin(fIntMax);
-  // for error calculation
-  Double_t intAreaEdgeLo = hMinv->GetBinLowEdge(binIntLo);
-  Double_t intAreaEdgeHi = hMinv->GetBinLowEdge(binIntHi)+hMinv->GetBinWidth(binIntHi);
-  Double_t norm = (binIntHi-binIntLo)/(fIntMax - fIntMin); 
-  
-  TH1F  *hBfit = (TH1F *) hMinv->Clone(); // for bg only fit (excluding peak region)
-  TH1F  *hSigF = (TH1F *) hMinv->Clone(); // signal with subtracted bg
-  TH1F  *hBgrd = (TH1F *) hMinv->Clone(); // bg histogram
-  
-  hBfit->Reset();
-  hSigF->Reset();
-  hBgrd->Reset();
-  
-  // extract start parameter for the MC signal fit
-  Double_t bgvalAv = (hMinv->Integral(1,hMinv->GetNbinsX()+1) - hMinv->Integral(binIntLo,binIntHi)) / (hMinv->GetNbinsX()+1 - (binIntHi-binIntLo));
-  Double_t pkval     = hMinv->GetBinContent(hMinv->FindBin(kJPMass)) - bgvalAv;
-  Double_t heightMC  = hSim->GetBinContent(hSim->FindBin(kJPMass));
-  Double_t peakScale = (heightMC > 0. ? pkval/heightMC : 0.0);
-  
-  Int_t nBgnd = 2 + fPolDeg; // degree + c1st oefficient + higher coefficients
-  Int_t nMinv = nBgnd + 1;  // bgrd + peakscale
-  
-  // Create the spectra without peak region
-  for (Int_t iBin = 0; iBin <= hMinv->GetNbinsX(); iBin++) {
-    if ((iBin < binIntLo) || (iBin > binIntHi)) {
-      hBfit->SetBinContent(iBin,hMinv->GetBinContent(iBin));
-      hBfit->SetBinError(iBin,hMinv->GetBinError(iBin));
-       }
-  }
-    
-  
-  // =======
-  //   1.
-  // =======
-  // Do the fit to the background spectrum
-  TF1 *fBo = new TF1("bgrd_fit",BgndFun,fFitMin,fFitMax,nBgnd);
-  for (Int_t iPar=0; iPar<nBgnd; iPar++) {
-    if (iPar == 0) fBo->FixParameter(0, fPolDeg);
-       if (iPar == 1) fBo->SetParameter(iPar, bgvalAv);
-    if (iPar >= 2) fBo->SetParameter(iPar, 0.);
-  }
-  hBfit->Fit(fBo,"0qR");
-  //hBfit->SetNameTitle("bgrd_fit");
-  arrResults->AddAt(fBo,1);
-  
-  
-  // =======
-  //   2.
-  // =======
-  // Fit the whole spectrum with peak and background
-  TF1 *fSB = new TF1("bgrd_peak_fit",MinvFun,fFitMin,fFitMax,nMinv);
-  fSB->FixParameter(0, fPolDeg);
-  fSB->SetParameter(1, peakScale);
-  // copy the polynomial parameters
-  for (Int_t iPar=0; iPar<=fPolDeg; iPar++) 
-       fSB->SetParameter(2+iPar, fBo->GetParameter(iPar+1));
-  hMinv->Fit(fSB,"0qR");
-  arrResults->AddAt(fSB,2);
-  
-  
-  // =======
-  //   3.
-  // =======
-  // Create the background function
-  TF1 *fB = new TF1("bgrdOnly_fkt",BgndFun,fFitMin,fFitMax,nBgnd);
-  fB->FixParameter(0,fPolDeg);
-  for (Int_t iDeg=0; iDeg<=fPolDeg; iDeg++) {
-       fB->SetParameter(1+iDeg,fSB->GetParameter(2+iDeg));
-       fB->SetParError(1+iDeg,fSB->GetParError(2+iDeg));
-  }
-  // create background histogram from background function
-  hBgrd->Eval(fB);
-  hBgrd->Fit(fB,"0qR");
-  // calculate the integral and integral error from fit function
-  Double_t intc = fB->Integral(intAreaEdgeLo, intAreaEdgeHi) * norm;
-  Double_t inte = fB->IntegralError(intAreaEdgeLo, intAreaEdgeHi) * norm;
-  arrResults->AddAt(fB,3);
-
-  // Fill the background spectrum erros. Use the error from the fit function for the background fB
-  for (Int_t iBin = 0; iBin <= hBgrd->GetNbinsX(); iBin++) {
-    Double_t x = hBgrd->GetBinCenter(iBin);
-    if ((x >= fFitMin) && (x <= fFitMax)) {
-         Double_t binte = inte / TMath::Sqrt((binIntHi-binIntLo)+1);
-         hBgrd->SetBinError(iBin,binte); 
-    }
-  }
-  arrResults->AddAt(hBgrd,4);
-  
-  // =======
-  //   4.
-  // =======
-  // Subtract the background  
-  hSigF->Add(hMinv,hBgrd,1.0,-1.0);
-  for (Int_t iBin = 0; iBin <= hSigF->GetNbinsX(); iBin++) {
-    Double_t x = hSigF->GetBinCenter(iBin);
-    if ((x < fFitMin) || (x > fFitMax)) {
-      hSigF->SetBinContent(iBin,0.0);
-      hSigF->SetBinError(iBin,0.0);
-    }
-  }
-  hSigF->SetNameTitle("peak_only","");
-  arrResults->AddAt(hSigF,5);
-  
-  // =======
-  //   5.
-  // =======
-  // Fit the background-subtracted spectrum
-  TF1 *fS  = new TF1("peakOnly_fit",PeakFunCB,fFitMin,fFitMax,5);
-  fS->SetParameters(-.05,1,kJPMass,.003,700);
-  fS->SetParNames("alpha","n","meanx","sigma","N");
-  hSigF->Fit(fS,"0qR");
-  arrResults->AddAt(fS,6);
-  
-  
-  // connect data members 
-  fFuncSignal     = (TF1*) arrResults->At(6)->Clone();
-  fFuncBackground = (TF1*) arrResults->At(3)->Clone();
-  fFuncSigBack    = (TF1*) arrResults->At(2)->Clone();
-  fHistSignal     = (TH1F*)arrResults->At(5)->Clone();
-  fHistBackground = (TH1F*)arrResults->At(4)->Clone();
-  
-  // fit results
-  Double_t chi2 = fSB->GetChisquare();
-  Int_t    ndf  = fSB->GetNDF();
-  fChi2Dof = chi2/ndf;
-  
-  // signal + signal error
-  fValues(0) = hSigF->IntegralAndError(binIntLo, binIntHi, fErrors(0));
-  fValues(1) = intc; // background
-  fErrors(1) = inte; // background error
-  // S/B (2) and significance (3)
-  SetSignificanceAndSOB();
-  fValues(4) = fS->GetParameter(2); // mass
-  fErrors(4) = fS->GetParError(2);  // mass error
-  fValues(5) = fS->GetParameter(3); // mass wdth
-  fErrors(5) = fS->GetParError(3);  // mass wdth error
-  
-    
-  delete arrResults;
-  
-}
-//______________________________________________________________________________
-Double_t AliDielectronSignalFunc::BgndFun(const Double_t *x, const Double_t *par) {
-  // parameters
-  // [0]:   degree of polynomial
-  // [1]:   constant polynomial coefficient
-  // [2]..: higher polynomial coefficients
-  
-  Int_t    deg   = ((Int_t) par[0]);
-  
-  Double_t f     = 0.0;
-  Double_t yy    = 1.0;
-  Double_t xx    = x[0];
-
-  for (Int_t i = 0; i <= deg; i++) {
-    f  += par[i+1] * yy;
-    yy *= xx;
-  } 
-  
-    
-  return f;
-}
 //______________________________________________________________________________
 //______________________________________________________________________________
-Double_t AliDielectronSignalFunc::PeakFun(const Double_t *x, const Double_t *par) {
+Double_t AliDielectronSignalFunc::PeakFunMC(const Double_t *x, const Double_t *par) {
   // Fit MC signal shape
   // parameters
   // [0]:   scale for simpeak
   // Fit MC signal shape
   // parameters
   // [0]:   scale for simpeak
@@ -370,42 +138,19 @@ Double_t AliDielectronSignalFunc::PeakFun(const Double_t *x, const Double_t *par
     return 0.0;
   }
   
     return 0.0;
   }
   
-  return (par[0] * hPeak->GetBinContent(idx));
-  
-}
-
-//______________________________________________________________________________
-Double_t AliDielectronSignalFunc::MinvFun(const Double_t *x, const Double_t *par) {
-  // parameters
-  // [0]:   degree of polynomial             -> [0]   for BgndFun
-  // [1]:   scale for simpeak                -> [0]   for PeakFun
-  // [2]:   constant polynomial coefficient  -> [1]   for BgndFun
-  // [3]..: higher polynomial coefficients   -> [2].. for BgndFun
-  
-  Int_t    deg = ((Int_t) par[0]);
-  Double_t parPK[25], parBG[25];
-  
-  parBG[0] = par[0]; // degree of polynomial
-  
-  parPK[0] = par[1]; // MC minv scale
-  for (Int_t i = 0; i <= deg; i++) parBG[i+1] = par[i+2]; // polynomial coefficients
-  
-  Double_t peak = PeakFun(x,parPK);
-  Double_t bgnd = BgndFun(x,parBG);
-  Double_t f    = peak + bgnd;
+  return (par[0+fNparBgnd] * hPeak->GetBinContent(idx));
   
   
-  return f;
 }
 
 //______________________________________________________________________________
 Double_t AliDielectronSignalFunc::PeakFunCB(const Double_t *x, const Double_t *par) {
 }
 
 //______________________________________________________________________________
 Double_t AliDielectronSignalFunc::PeakFunCB(const Double_t *x, const Double_t *par) {
-  // Crystal Ball function fit
+  // Crystal Ball fit function
   
   
-  Double_t alpha = par[0];
-  Double_t     n = par[1];
-  Double_t meanx = par[2];
-  Double_t sigma = par[3];
-  Double_t    nn = par[4];
+  Double_t alpha = par[0+fNparBgnd];
+  Double_t     n = par[1+fNparBgnd];
+  Double_t meanx = par[2+fNparBgnd];
+  Double_t sigma = par[3+fNparBgnd];
+  Double_t    nn = par[4+fNparBgnd];
    
   Double_t a = TMath::Power((n/TMath::Abs(alpha)), n) * TMath::Exp(-.5*alpha*alpha);
   Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
    
   Double_t a = TMath::Power((n/TMath::Abs(alpha)), n) * TMath::Exp(-.5*alpha*alpha);
   Double_t b = n/TMath::Abs(alpha) - TMath::Abs(alpha);
@@ -422,6 +167,17 @@ Double_t AliDielectronSignalFunc::PeakFunCB(const Double_t *x, const Double_t *p
   return fitval;
 }
 
   return fitval;
 }
 
+//______________________________________________________________________________
+Double_t AliDielectronSignalFunc::PeakFunGaus(const Double_t *x, const Double_t *par) {
+  // Gaussian fit function
+
+  Double_t     n = par[0+fNparBgnd];
+  Double_t  mean = par[1+fNparBgnd];
+  Double_t sigma = par[2+fNparBgnd];
+  Double_t    xx = x[0];
+
+  return ( n*exp(-0.5*TMath::Power((xx-mean)/sigma,2)) );
+}
 
 //______________________________________________
 void AliDielectronSignalFunc::ProcessFit(TObjArray * const arrhist) {
 
 //______________________________________________
 void AliDielectronSignalFunc::ProcessFit(TObjArray * const arrhist) {
@@ -448,8 +204,9 @@ void AliDielectronSignalFunc::ProcessFit(TObjArray * const arrhist) {
   fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
   TFitResultPtr pmFitPtr = fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
   //TFitResult *pmFitResult = pmFitPtr.Get(); // used only with TF1Helper
   fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
   TFitResultPtr pmFitPtr = fHistDataPM->Fit(fFuncSigBack, fFitOpt.Data(), "", fFitMin, fFitMax);
   //TFitResult *pmFitResult = pmFitPtr.Get(); // used only with TF1Helper
-  fFuncSignal->SetParameters(fFuncSigBack->GetParameters());
-  fFuncBackground->SetParameters(fFuncSigBack->GetParameters()+fFuncSignal->GetNpar());
+  fFuncBackground->SetParameters(fFuncSigBack->GetParameters());
+  fFuncSignal->SetParameters(fFuncSigBack->GetParameters()+fFuncBackground->GetNpar());
+  //  fFuncBackground->SetParameters(fFuncSigBack->GetParameters()+fFuncSignal->GetNpar());
   
   for(Int_t iBin=1; iBin<=fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
     Double_t m = fHistDataPM->GetBinCenter(iBin);
   
   for(Int_t iBin=1; iBin<=fHistDataPM->GetXaxis()->GetNbins(); iBin++) {
     Double_t m = fHistDataPM->GetBinCenter(iBin);
@@ -529,6 +286,9 @@ void AliDielectronSignalFunc::ProcessFit(TObjArray * const arrhist) {
   fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
   
   fProcessed = kTRUE;
   fErrors(5) = fFuncSigBack->GetParError(fParMassWidth);
   
   fProcessed = kTRUE;
+  
+  fHistBackground->GetListOfFunctions()->Add(fFuncBackground);
+
 }
 
 //______________________________________________
 }
 
 //______________________________________________
@@ -664,6 +424,7 @@ void AliDielectronSignalFunc::SetFunctions(TF1 * const combined, TF1 * const sig
   fFuncSigBack=combined;
   fParMass=parM;
   fParMassWidth=parMres;
   fFuncSigBack=combined;
   fParMass=parM;
   fParMassWidth=parMres;
+
 }
 
 //______________________________________________
 }
 
 //______________________________________________
@@ -777,3 +538,33 @@ void AliDielectronSignalFunc::Draw(const Option_t* option)
   if (optStat) DrawStats();
   
 }
   if (optStat) DrawStats();
   
 }
+
+
+//______________________________________________________________________________
+void AliDielectronSignalFunc::CombineFunc(TF1 * const peak, TF1 * const bgnd) {
+  //
+  // combine the bgnd and the peak function
+  //
+
+  if (!peak||!bgnd) {
+    AliError("Both, signal and background function need to be set!");
+    return;
+  }
+  fFuncSignal=peak;
+  fFuncBackground=bgnd;
+
+  fNparPeak     = fFuncSignal->GetNpar();
+  fNparBgnd     = fFuncBackground->GetNpar();
+
+  fFuncSigBack = new TF1("BgndPeak",AliDielectronSignalFunc::PeakBgndFun, 0.,5.,fNparPeak+fNparBgnd);
+  return;
+}
+
+//______________________________________________________________________________
+Double_t AliDielectronSignalFunc::PeakBgndFun(const Double_t *x, const Double_t *par) {
+  //
+  // merge peak and bgnd functions
+  //
+  return (fFuncBackground->EvalPar(x,par) + fFuncSignal->EvalPar(x,par));
+}
+
index c8efe5a4bc0d4dbc9a2379011cffc7a576da979c..2a63a9e2954cdc52d0789f9c93c2e0e8c2a59c90 100644 (file)
@@ -82,15 +82,13 @@ public:
   virtual ~AliDielectronSignalFunc();
 
   virtual void Process(TObjArray * const arrhist);
   virtual ~AliDielectronSignalFunc();
 
   virtual void Process(TObjArray * const arrhist);
-  void ProcessFitIKF(TObjArray * const arrhist);      // fit the SE +- distribution with mc shape and best chi2/dof
-  void FitOneMinv(TH1F *hMinv, TH1F *hSim, Int_t pod); // fit procedure for one single minv spectrum
-    
   void ProcessFit(TObjArray * const arrhist);      // fit the SE +- distribution
   void ProcessLS(TObjArray * const arrhist);       // substract the fitted SE like-sign background
   void ProcessEM(TObjArray * const arrhist);       // substract the fitted SE+ME like-sign background
   void ProcessFit(TObjArray * const arrhist);      // fit the SE +- distribution
   void ProcessLS(TObjArray * const arrhist);       // substract the fitted SE like-sign background
   void ProcessEM(TObjArray * const arrhist);       // substract the fitted SE+ME like-sign background
-  
+
   void SetUseIntegral(Bool_t flag=kTRUE) {fUseIntegral = flag;};
   void SetFunctions(TF1 * const combined, TF1 * const sig=0, TF1 * const back=0, Int_t parM=1, Int_t parMres=2);
   void SetUseIntegral(Bool_t flag=kTRUE) {fUseIntegral = flag;};
   void SetFunctions(TF1 * const combined, TF1 * const sig=0, TF1 * const back=0, Int_t parM=1, Int_t parMres=2);
+  void CombineFunc(TF1 * const peak=0, TF1 * const bgnd=0);
   void SetFitOption(const char* opt) {
     fFitOpt=opt; 
     fFitOpt.ToLower(); 
   void SetFitOption(const char* opt) {
     fFitOpt=opt; 
     fFitOpt.ToLower(); 
@@ -103,20 +101,22 @@ public:
   TF1*  GetCombinedFunction()   const { return fFuncSigBack;       }
   
   Int_t GetPolDeg()             const { return fPolDeg;            }
   TF1*  GetCombinedFunction()   const { return fFuncSigBack;       }
   
   Int_t GetPolDeg()             const { return fPolDeg;            }
+  Int_t GetDof()                const { return fDof;            }
   Double_t GetChi2Dof()         const { return fChi2Dof;           }
   
   virtual void Draw(const Option_t* option = "");
   
 private:
   Double_t GetChi2Dof()         const { return fChi2Dof;           }
   
   virtual void Draw(const Option_t* option = "");
   
 private:
-  
-  static Double_t BgndFun(const Double_t *x, const Double_t *par); // polynomial fit function
-  static Double_t PeakFun(const Double_t *x, const Double_t *par); // create peak function from mc shape
-  static Double_t MinvFun(const Double_t *x, const Double_t *par); // combined function of BgndFun + PeakFun
-  static Double_t PeakFunCB(const Double_t *x, const Double_t *par); // crystal ball function for background-subtracted peak fit
-  
 
 
-  TF1 *fFuncSignal;                // Function for the signal description
-  TF1 *fFuncBackground;            // Function for the background description
+  // peak functions
+  static Double_t PeakFunMC(const Double_t *x, const Double_t *par); // peak function from a mc histo
+  static Double_t PeakFunCB(const Double_t *x, const Double_t *par); // crystal ball function
+  static Double_t PeakFunGaus(const Double_t *x, const Double_t *par); // gaussian
+
+  static Double_t PeakBgndFun(const Double_t *x, const Double_t *par); // combine any bgrd and any peak function
+
+  static TF1 *fFuncSignal;                // Function for the signal description
+  static TF1 *fFuncBackground;            // Function for the background description
   TF1 *fFuncSigBack;               // Combined function signal plus background
   Int_t fParMass;                  // the index of the parameter corresponding to the resonance mass
   Int_t fParMassWidth;             // the index of the parameter corresponding to the resonance mass width
   TF1 *fFuncSigBack;               // Combined function signal plus background
   Int_t fParMass;                  // the index of the parameter corresponding to the resonance mass
   Int_t fParMassWidth;             // the index of the parameter corresponding to the resonance mass width
@@ -125,7 +125,12 @@ private:
   Bool_t fUseIntegral;         // use the integral of the fitted functions to extract signal and background
 
   Int_t    fPolDeg;                // polynomial degree of the background function
   Bool_t fUseIntegral;         // use the integral of the fitted functions to extract signal and background
 
   Int_t    fPolDeg;                // polynomial degree of the background function
+  Int_t    fDof;                   // degrees of freedom
   Double_t fChi2Dof;               // chi2/dof of the fitted inv mass spectra
   Double_t fChi2Dof;               // chi2/dof of the fitted inv mass spectra
+  
+  static Int_t    fNparPeak;              // number of parameters for peak function
+  static Int_t    fNparBgnd;              // number of parameters for background function
+
   static TH1F* fgHistSimPM;         // simulated peak shape 
   
   ClassDef(AliDielectronSignalFunc,2)         // Dielectron SignalFunc
   static TH1F* fgHistSimPM;         // simulated peak shape 
   
   ClassDef(AliDielectronSignalFunc,2)         // Dielectron SignalFunc