]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWG3/vertexingHF/AliHFMassFitter.cxx
Added LHC10h run list for flow analysis (Giacomo)
[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliHFMassFitter.cxx
index 9913fa63742ce3fd676c88886084683cbadb46ff..f10fe80f06fb523566bfb7d88eeb880354ebf8c1 100644 (file)
@@ -13,6 +13,8 @@
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
+/* $Id$ */
+
 /////////////////////////////////////////////////////////////
 //
 // AliHFMassFitter for the fit of invariant mass distribution
@@ -35,6 +37,7 @@
 #include <TMinuit.h>
 #include <TStyle.h>
 #include <TPaveText.h>
+#include <TDatabasePDG.h>
 
 #include "AliHFMassFitter.h"
 
@@ -113,7 +116,7 @@ AliHFMassFitter::AliHFMassFitter (const TH1F *histoToFit, Double_t minvalue, Dou
   CheckRangeFit();
   ftypeOfFit4Bkg=fittypeb;
   ftypeOfFit4Sgn=fittypes;
-  if(ftypeOfFit4Bkg!=0 && ftypeOfFit4Bkg!=1 && ftypeOfFit4Bkg!=2) fWithBkg=kFALSE;
+  if(ftypeOfFit4Bkg!=0 && ftypeOfFit4Bkg!=1 && ftypeOfFit4Bkg!=2 && ftypeOfFit4Bkg!=4 && ftypeOfFit4Bkg!=5) fWithBkg=kFALSE;
   else fWithBkg=kTRUE;
   if (!fWithBkg) cout<<"Fit Histogram of Signal only"<<endl;
   else  cout<<"Type of fit For Background = "<<ftypeOfFit4Bkg<<endl;
@@ -172,30 +175,15 @@ AliHFMassFitter::~AliHFMassFitter() {
   //destructor
 
   cout<<"AliHFMassFitter destructor called"<<endl;
-  if(fhistoInvMass) {
-    cout<<"deleting histogram..."<<endl;
-    delete fhistoInvMass;
-    fhistoInvMass=NULL;
-  }
-  if(fntuParam){
-    cout<<"deleting ntuple..."<<endl;
-    delete fntuParam;
 
-    fntuParam=NULL;
-  }
+  delete fhistoInvMass;
 
-  if(fFitPars) {
-    delete[] fFitPars;
-    cout<<"deleting parameter array..."<<endl;
-    fFitPars=NULL;
-  }
-
-  if(fFixPar) {
-    delete[] fFixPar;
-    cout<<"deleting bool array..."<<endl;
-    fFixPar=NULL;
-  }
+  delete fntuParam;
 
+  delete[] fFitPars;
+  
+  delete[] fFixPar;
+  
   fcounter = 0;
 }
 
@@ -225,17 +213,13 @@ AliHFMassFitter& AliHFMassFitter::operator=(const AliHFMassFitter &mfit){
   fContourGraph= mfit.fContourGraph;
 
   if(mfit.fParsSize > 0){
-    if(fFitPars) {
-      delete[] fFitPars;
-      fFitPars=NULL;
-    }
+    delete[] fFitPars;
+    
     fFitPars=new Float_t[fParsSize];
     memcpy(fFitPars,mfit.fFitPars,mfit.fParsSize*sizeof(Float_t));
 
-    if(fFixPar) {
-      delete[] fFixPar;
-      fFixPar=NULL;
-    }
+    delete[] fFixPar;
+    
     fFixPar=new Bool_t[fNFinalPars];
     memcpy(fFixPar,mfit.fFixPar,mfit.fNFinalPars*sizeof(Float_t));
   }
@@ -264,6 +248,12 @@ void AliHFMassFitter::ComputeNFinalPars() {
   case 3:
     fNFinalPars=1;
     break;
+  case 4:
+    fNFinalPars=2;     
+    break;
+  case 5:
+    fNFinalPars=3;     
+    break;
   default:
     cout<<"Error in computing fNFinalPars: check ftypeOfFit4Bkg"<<endl;
     break;
@@ -291,6 +281,12 @@ void AliHFMassFitter::ComputeParSize() {
   case 3:
     fParsSize = 1*3;
     break;
+  case 4:
+    fParsSize = 2*3; 
+    break;
+  case 5:
+    fParsSize = 3*3; 
+    break;
   default:
     cout<<"Error in computing fParsSize: check ftypeOfFit4Bkg"<<endl;
     break;
@@ -362,7 +358,7 @@ void AliHFMassFitter::SetHisto(const TH1F *histoToFit){
 
   fhistoInvMass = new TH1F(*histoToFit);
   fhistoInvMass->SetDirectory(0);
-  cout<<"SetHisto pointer "<<fhistoInvMass<<endl;
+  //cout<<"SetHisto pointer "<<fhistoInvMass<<endl;
 }
 
 //___________________________________________________________________________
@@ -390,14 +386,7 @@ void AliHFMassFitter::Reset() {
   fMass=1.85;
   fSigmaSgn=0.012;
   cout<<"Reset "<<fhistoInvMass<<endl;
-  if(fhistoInvMass) {
-    delete fhistoInvMass;
-    fhistoInvMass=NULL;
-    cout<<fhistoInvMass<<endl;
-  }
-  else cout<<"histogram doesn't exist, do not delete"<<endl;
-  
-
+  delete fhistoInvMass;
 }
 
 //_________________________________________________________________________
@@ -632,6 +621,52 @@ Double_t AliHFMassFitter::FitFunction4MassDistr (Double_t *x, Double_t *par){
     }
   }
 
+  //Power fit
+
+    // par[0] = tot integral
+    // par[1] = coef1
+    // par[2] = gaussian integral
+    // par[3] = gaussian mean
+    // par[4] = gaussian sigma
+
+  if (ftypeOfFit4Bkg==4) {
+    
+    if(ftypeOfFit4Sgn == 0) {
+      
+      Double_t parbkg[2] = {par[0]-par[2], par[1]}; 
+      bkg = FitFunction4Bkg(x,parbkg);
+    }
+    if(ftypeOfFit4Sgn == 1) {
+      
+      Double_t parbkg[5] = {par[2],par[3],ffactor*par[4],par[0]-par[2], par[1]}; 
+      bkg = FitFunction4Bkg(x,parbkg);
+    }
+    sgn = FitFunction4Sgn(x,&par[2]);
+  }
+
+
+  //Power and exponential fit
+
+    // par[0] = tot integral
+    // par[1] = coef1
+    // par[2] = coef2
+    // par[3] = gaussian integral
+    // par[4] = gaussian mean
+    // par[5] = gaussian sigma
+
+  if (ftypeOfFit4Bkg==5) {      
+   
+    if(ftypeOfFit4Sgn == 0) {
+       Double_t parbkg[3] = {par[0]-par[3],par[1],par[2]}; 
+       bkg = FitFunction4Bkg(x,parbkg); 
+    }
+    if(ftypeOfFit4Sgn == 1) {
+     Double_t parbkg[6] = {par[3],par[4],ffactor*par[5],par[0]-par[3], par[1], par[2]}; 
+     bkg = FitFunction4Bkg(x,parbkg);
+    }
+    sgn = FitFunction4Sgn(x,&par[3]);
+  }                            
+
   total = bkg + sgn;
   
   return  total;
@@ -707,12 +742,37 @@ Double_t AliHFMassFitter::FitFunction4Bkg (Double_t *x, Double_t *par){
   case 3:
     total=par[0+firstPar];
     break;
+  case 4:  
+    //power function 
+    //y=a(x-m_pi)^b -> integral = a/(b+1)*((max-m_pi)^(b+1)-(min-m_pi)^(b+1))
+    //
+    //a = integral*(b+1)/((max-m_pi)^(b+1)-(min-m_pi)^(b+1))
+    // * [0] = integralBkg;
+    // * [1] = b;
+    // a(power function) = [0]*([1]+1)/((max-m_pi)^([1]+1)-(min-m_pi)^([1]+1))*(x-m_pi)^[1]
+    {
+    Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
+
+    total = par[0+firstPar]*(par[1+firstPar]+1.)/(TMath::Power(fmaxMass-mpi,par[1+firstPar]+1.)-TMath::Power(fminMass-mpi,par[1+firstPar]+1.))*TMath::Power(x[0]-mpi,par[1+firstPar]);
+    }
+    break;
+  case 5:
+   //power function wit exponential
+    //y=a*Sqrt(x-m_pi)*exp(-b*(x-m_pi))  
+    { 
+    Double_t mpi = TDatabasePDG::Instance()->GetParticle(211)->Mass();
+
+    total = par[1+firstPar]*TMath::Sqrt(x[0] - mpi)*TMath::Exp(-1.*par[2+firstPar]*(x[0]-mpi));
+    } 
+    break;
 //   default:
 //     Types of Fit Functions for Background:
 //     * 0 = exponential;
 //     * 1 = linear;
 //     * 2 = polynomial 2nd order
 //     * 3 = no background"<<endl;
+//     * 4 = Power function 
+//     * 5 = Power function with exponential 
 
   }
   return total+gaus2;
@@ -913,10 +973,10 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
 
   //Total integral
   Double_t totInt = fhistoInvMass->Integral(fminBinMass,fmaxBinMass, "width");
-  cout<<"Here tot integral is = "<<totInt<<"; integral in whole range is "<<fhistoInvMass->Integral("width")<<endl;
+  //cout<<"Here tot integral is = "<<totInt<<"; integral in whole range is "<<fhistoInvMass->Integral("width")<<endl;
   fSideBands = kTRUE;
   Double_t width=fhistoInvMass->GetBinWidth(8);
-  cout<<"fNbin"<<fNbin<<endl;
+  //cout<<"fNbin = "<<fNbin<<endl;
   if (fNbin==0) fNbin=fhistoInvMass->GetNbinsX();
 
   Bool_t ok=SideBandsBounds();
@@ -961,6 +1021,14 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
       funcbkg->FixParameter(0,0.);
     }
     break;
+  case 4:
+    funcbkg->SetParNames("BkgInt","Coef2");  
+    funcbkg->SetParameters(sideBandsInt,0.5);
+    break;
+ case 5:
+    funcbkg->SetParNames("BkgInt","Coef1","Coef2");
+    funcbkg->SetParameters(sideBandsInt, -10., 5.);
+    break;
   default:
     cout<<"Wrong choise of ftypeOfFit4Bkg ("<<ftypeOfFit4Bkg<<")"<<endl;
     return kFALSE;
@@ -987,7 +1055,10 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     intbkg1 = funcbkg->Integral(fminMass,fmaxMass);
     if(ftypeOfFit4Bkg!=3) slope1 = funcbkg->GetParameter(1);
     if(ftypeOfFit4Bkg==2) conc1 = funcbkg->GetParameter(2);
-    cout<<"First fit: \nintbkg1 = "<<intbkg1<<"\t(Compare with par0 = "<<funcbkg->GetParameter(0)<<")\nslope1= "<<slope1<<"\nconc1 = "<<conc1<<endl;
+    if(ftypeOfFit4Bkg==5) conc1 = funcbkg->GetParameter(2); 
+
+    //cout<<"First fit: \nintbkg1 = "<<intbkg1<<"\t(Compare with par0 = "<<funcbkg->GetParameter(0)<<")\nslope1= "<<slope1<<"\nconc1 = "<<conc1<<endl;
   } 
   else cout<<"\t\t//"<<endl;
   
@@ -997,34 +1068,62 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     cout<<"\nBACKGROUND FIT WITH REFLECTION"<<endl;
     bkgPar+=3;
     
-    cout<<"fNFinalPars = "<<fNFinalPars<<"\tbkgPar = "<<bkgPar<<endl;
+    //cout<<"fNFinalPars = "<<fNFinalPars<<"\tbkgPar = "<<bkgPar<<endl;
 
     funcbkg1 = new TF1(bkg1name.Data(),this,&AliHFMassFitter::FitFunction4Bkg,fminMass,fmaxMass,bkgPar,"AliHFMassFitter","FitFunction4Bkg");
     cout<<"Function name = "<<funcbkg1->GetName()<<endl;
 
     funcbkg1->SetLineColor(2); //red
 
-    if(ftypeOfFit4Bkg==2){
-      cout<<"*** Polynomial Fit ***"<<endl;
-      funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Coef1","Coef2");
-      funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1,conc1);
-
-      //cout<<"Parameters set to: "<<0.5*(totInt-intbkg1)<<"\t"<<fMass<<"\t"<<ffactor*fSigmaSgn<<"\t"<<intbkg1<<"\t"<<slope1<<"\t"<<conc1<<"\t"<<endl;
-      //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<bkgPar<<"\tgsidebands = "<<fSideBands<<endl;
-    } else{
-      if(ftypeOfFit4Bkg==3) //no background: gaus sign+ gaus broadened
+ switch (ftypeOfFit4Bkg) {
+    case 0:
+       {
+        cout<<"*** Exponential Fit ***"<<endl;
+        funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Slope");
+       funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1);
+       }
+        break;
+     case 1: 
+       {
+       cout<<"*** Linear Fit ***"<<endl;
+        funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Slope");
+       funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1);
+       }
+        break;    
+     case 2:
+        {
+        cout<<"*** Polynomial Fit ***"<<endl;
+        funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Coef1","Coef2");
+        funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1,conc1);
+        }
+        break;
+    case 3:
+       //no background: gaus sign+ gaus broadened
        {
-         cout<<"*** No background Fit ***"<<endl;
-         funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","Const");
-         funcbkg1->SetParameters(0.5*totInt,fMass,ffactor*fSigmaSgn,0.); 
-         funcbkg1->FixParameter(3,0.);
-       } else{ //expo or linear
-         if(ftypeOfFit4Bkg==0) cout<<"*** Exponential Fit ***"<<endl;
-         if(ftypeOfFit4Bkg==1) cout<<"*** Linear Fit ***"<<endl;
-         funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Slope");
-         funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1);
+       cout<<"*** No background Fit ***"<<endl;
+       funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","Const");
+       funcbkg1->SetParameters(0.5*totInt,fMass,ffactor*fSigmaSgn,0.); 
+       funcbkg1->FixParameter(3,0.);
+       } 
+        break;
+     case 4:
+       {       
+       cout<<"*** Power function Fit ***"<<endl;
+       funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Coef2");
+        funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1);
+               }
+        break;
+      case 5:
+       { 
+       cout<<"*** Power function conv. with exponential Fit ***"<<endl;
+        funcbkg1->SetParNames("IntGB","MeanGB","SigmaGB","BkgInt","Coef1","Coef2");
+        funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1,conc1);
        }
+        break;
     }
+      //cout<<"Parameters set to: "<<0.5*(totInt-intbkg1)<<"\t"<<fMass<<"\t"<<ffactor*fSigmaSgn<<"\t"<<intbkg1<<"\t"<<slope1<<"\t"<<conc1<<"\t"<<endl;
+      //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<bkgPar<<"\tgsidebands = "<<fSideBands<<endl;
+
     Int_t status=fhistoInvMass->Fit(bkg1name.Data(),"R,L,E,+,0");
     if (status != 0){
       cout<<"Minuit returned "<<status<<endl;
@@ -1041,6 +1140,8 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     intbkg1=funcbkg1->GetParameter(3);
     if(ftypeOfFit4Bkg!=3) slope1 = funcbkg1->GetParameter(4);
     if(ftypeOfFit4Bkg==2) conc1 = funcbkg1->GetParameter(5);
+    if(ftypeOfFit4Bkg==5) conc1 = funcbkg1->GetParameter(5); 
+
 
   } else {
     bkgPar+=3;
@@ -1065,18 +1166,18 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
   //sidebands integral - second approx (from fit)
   fSideBands = kFALSE;
   Double_t bkgInt;
-  cout<<"Compare intbkg1 = "<<intbkg1<<" and integral = ";
+  //cout<<"Compare intbkg1 = "<<intbkg1<<" and integral = ";
   if(ftypeOfFit4Sgn == 1) bkgInt=funcbkg1->Integral(fminMass,fmaxMass);
   else bkgInt=funcbkg->Integral(fminMass,fmaxMass);
-  cout<</*"------BkgInt(Fit) = "<<*/bkgInt<<endl;
+  //cout<</*"------BkgInt(Fit) = "<<*/bkgInt<<endl;
 
   //Signal integral - first approx
   Double_t sgnInt;
   sgnInt = totInt-bkgInt;
-  cout<<"------TotInt = "<<totInt<<"\tsgnInt = "<<sgnInt<<endl;
+  //cout<<"------TotInt = "<<totInt<<"\tsgnInt = "<<sgnInt<<endl;
   if (sgnInt <= 0){
     cout<<"Setting sgnInt = - sgnInt"<<endl;
-    sgnInt=sgnInt;
+    sgnInt=(-1)*sgnInt;
   }
   /*Fit All Mass distribution with exponential + gaussian (+gaussian braodened) */
   TF1 *funcmass = new TF1(massname.Data(),this,&AliHFMassFitter::FitFunction4MassDistr,fminMass,fmaxMass,fNFinalPars,"AliHFMassFitter","FitFunction4MassDistr");
@@ -1093,23 +1194,18 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     //cout<<"Parameters set to: "<<totInt<<"\t"<<slope1<<"\t"<<sgnInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
     //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<fNFinalPars<<"\tgsidebands = "<<fSideBands<<endl;
     if(fFixPar[0]){
-      cout<<"fix1"<<endl;
       funcmass->FixParameter(0,totInt);
     }
     if(fFixPar[1]){
-      cout<<"fix2"<<endl;
       funcmass->FixParameter(1,slope1);
     }
     if(fFixPar[2]){
-      cout<<"fix3"<<endl;
       funcmass->FixParameter(2,sgnInt);
     }
     if(fFixPar[3]){
-      cout<<"fix4"<<endl;
       funcmass->FixParameter(3,fMass);
     }
     if(fFixPar[4]){
-      cout<<"fix5"<<endl;
       funcmass->FixParameter(4,fSigmaSgn);
     }
   }
@@ -1124,7 +1220,7 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     if(fFixPar[2])funcmass->FixParameter(2,conc1);
     if(fFixPar[3])funcmass->FixParameter(3,sgnInt);
     if(fFixPar[4])funcmass->FixParameter(4,fMass);
-    if(fFixPar[5])funcmass->FixParameter(5,fSigmaSgn);
+    if(fFixPar[5])funcmass->FixParameter(5,fSigmaSgn); 
     //
     //funcmass->FixParameter(2,sgnInt);
   }
@@ -1133,7 +1229,10 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     if(ftypeOfFit4Sgn == 1) funcmass->SetParameters(0.,0.5*totInt,fMass,fSigmaSgn);
     else funcmass->SetParameters(0.,totInt,fMass,fSigmaSgn);
     if(fFixPar[0]) funcmass->FixParameter(0,0.);
-    //cout<<"Parameters set to: "<<0.5*totInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
+    if(fFixPar[1])funcmass->FixParameter(1,sgnInt);
+    if(fFixPar[2])funcmass->FixParameter(2,fMass);
+    if(fFixPar[3])funcmass->FixParameter(3,fSigmaSgn);
+   //cout<<"Parameters set to: "<<0.5*totInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
     //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<fNFinalPars<<"\tgsidebands = "<<fSideBands<<endl;
 
   }
@@ -1224,19 +1323,12 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     
   }
 
-  if (ftypeOfFit4Sgn == 1 && funcbkg1) {
+  if (ftypeOfFit4Sgn == 1) {
     delete funcbkg1;
-    funcbkg1=NULL;
-  }
-  if (funcbkg) {
-    delete funcbkg;
-    funcbkg=NULL;
   }
-  if (funcmass) {
-    delete funcmass;
-    funcmass=NULL;
-  }
-
+  delete funcbkg;
+  delete funcmass;
+  
   AddFunctionsToHisto();
   if (draw) DrawFit();
  
@@ -1274,13 +1366,21 @@ Bool_t AliHFMassFitter::RefitWithBkgOnly(Bool_t draw){
     funcbkg->SetParameters(integral,-10.,5);
     break;
   case 3:
-    cout<<"Warning! This choice does not have a lot of sense..."<<endl;
+    cout<<"Warning! This choice does not make a lot of sense..."<<endl;
     if(ftypeOfFit4Sgn==0){
       funcbkg->SetParNames("Const");
       funcbkg->SetParameter(0,0.);
       funcbkg->FixParameter(0,0.);
     }
     break;
+  case 4:     
+    funcbkg->SetParNames("BkgInt","Coef1");
+    funcbkg->SetParameters(integral,0.5);
+    break;
+  case 5:    
+    funcbkg->SetParNames("BkgInt","Coef1","Coef2");
+    funcbkg->SetParameters(integral,-10.,5.);
+    break;
   default:
     cout<<"Wrong choise of ftypeOfFit4Bkg ("<<ftypeOfFit4Bkg<<")"<<endl;
     return kFALSE;
@@ -1302,12 +1402,18 @@ Bool_t AliHFMassFitter::RefitWithBkgOnly(Bool_t draw){
 }
 //_________________________________________________________________________
 Double_t AliHFMassFitter::GetChiSquare() const{
+  //Get Chi^2 method
   TF1 *funcmass=(TF1*)fhistoInvMass->GetFunction("funcmass");
+  if(!funcmass) {
+    cout<<"funcmass not found"<<endl;
+    return -1;
+  }
   return funcmass->GetChisquare();
 }
 
 //_________________________________________________________________________
 Double_t AliHFMassFitter::GetReducedChiSquare() const{
+  //Get reduced Chi^2 method
   TF1 *funcmass=(TF1*)fhistoInvMass->GetFunction("funcmass");
   if(!funcmass) {
     cout<<"funcmass not found"<<endl;
@@ -1332,13 +1438,16 @@ void  AliHFMassFitter::GetFitPars(Float_t *vector) const {
 void AliHFMassFitter::IntS(Float_t *valuewitherror) const {
 
   //gives the integral of signal obtained from fit parameters
-  if(!valuewitherror)valuewitherror=new Float_t[2];
+  if(!valuewitherror) {
+    printf("AliHFMassFitter::IntS: got a null pointer\n");
+    return;
+  }
 
   Int_t index=fParsSize/2 - 3;
   valuewitherror[0]=fFitPars[index];
   index=fParsSize - 3;
   valuewitherror[1]=fFitPars[index];
-  }
+}
 
 
 //_________________________________________________________________________
@@ -1346,7 +1455,7 @@ void AliHFMassFitter::AddFunctionsToHisto(){
 
   //Add the background function in the complete range to the list of functions attached to the histogram
 
-  cout<<"AddFunctionsToHisto called"<<endl;
+  //cout<<"AddFunctionsToHisto called"<<endl;
   TString bkgname = "funcbkg";
 
   Bool_t done1=kFALSE,done2=kFALSE;
@@ -1381,10 +1490,7 @@ void AliHFMassFitter::AddFunctionsToHisto(){
     }else{
       bonly->SetLineColor(kBlue+3);
       hlist->Add((TF1*)bonly->Clone());
-      if(bonly) {
-       delete bonly;
-       bonly=NULL;
-      }
+      delete bonly;
     }
 
   }
@@ -1397,10 +1503,9 @@ void AliHFMassFitter::AddFunctionsToHisto(){
     }
 
     bkgname += "FullRange";
-    TF1 *bfullrange=new TF1(bkgname.Data(),this,&AliHFMassFitter::FitFunction4Bkg,fminMass,fmaxMass,fNFinalPars,"AliHFMassFitter","FitFunction4Bkg");
+    TF1 *bfullrange=new TF1(bkgname.Data(),this,&AliHFMassFitter::FitFunction4Bkg,fminMass,fmaxMass,fNFinalPars-3,"AliHFMassFitter","FitFunction4Bkg");
     //cout<<bfullrange->GetName()<<endl;
-    for(Int_t i=0;i<fNFinalPars;i++){
-      //cout<<i<<" di "<<fNFinalPars<<endl;
+    for(Int_t i=0;i<fNFinalPars-3;i++){
       bfullrange->SetParName(i,b->GetParName(i));
       bfullrange->SetParameter(i,b->GetParameter(i));
       bfullrange->SetParError(i,b->GetParError(i));
@@ -1410,7 +1515,7 @@ void AliHFMassFitter::AddFunctionsToHisto(){
 
     bkgnamesave += "Recalc";
 
-    TF1 *blastpar=new TF1(bkgnamesave.Data(),this,&AliHFMassFitter::FitFunction4Bkg,fminMass,fmaxMass,fNFinalPars,"AliHFMassFitter","FitFunction4Bkg");
+    TF1 *blastpar=new TF1(bkgnamesave.Data(),this,&AliHFMassFitter::FitFunction4Bkg,fminMass,fmaxMass,fNFinalPars-3,"AliHFMassFitter","FitFunction4Bkg");
 
     TF1 *mass=fhistoInvMass->GetFunction("funcmass");
 
@@ -1422,11 +1527,11 @@ void AliHFMassFitter::AddFunctionsToHisto(){
     //intBkg=intTot-intS
     blastpar->SetParameter(0,mass->GetParameter(0)-mass->GetParameter(fNFinalPars-3));
     blastpar->SetParError(0,mass->GetParError(fNFinalPars-3));
-    if (fNFinalPars>=2) {
+    if (fNFinalPars>=5) {
       blastpar->SetParameter(1,mass->GetParameter(1));
       blastpar->SetParError(1,mass->GetParError(1));
     }
-    if (fNFinalPars==3) {
+    if (fNFinalPars==6) {
       blastpar->SetParameter(2,mass->GetParameter(2));
       blastpar->SetParError(2,mass->GetParError(2));
     }
@@ -1438,14 +1543,9 @@ void AliHFMassFitter::AddFunctionsToHisto(){
     hlist->Add((TF1*)blastpar->Clone());
     hlist->ls();
   
-    if(bfullrange) {
-      delete bfullrange;
-      bfullrange=NULL;
-    }
-    if(blastpar) {
-      delete blastpar;
-      blastpar=NULL;
-    }
+    delete bfullrange;
+    delete blastpar;
+    
   }
 
 
@@ -1484,11 +1584,8 @@ void AliHFMassFitter::WriteHisto(TString path) const {
 
   cout<<fcounter<<" "<<hget->GetName()<<" written in "<<path<<endl;
 
-  if(output) {
-    delete output;
-    output=NULL;
-  }
-
+  delete output;
+  
 }
 
 //_________________________________________________________________________
@@ -1509,10 +1606,8 @@ void AliHFMassFitter::WriteNtuple(TString path) const{
     nget=NULL;
   }
   */
-  if(output) {
-    delete output;
-    output=NULL;
-  }
+
+  delete output;
 }
 
 //_________________________________________________________________________
@@ -1539,6 +1634,12 @@ void AliHFMassFitter::WriteCanvas(TString userIDstring,TString path,Double_t nsi
   case 3:
     type="noB"; //3+1
     break;
+  case 4:  
+    type="Pow"; //3+3
+    break;
+  case 5:
+    type="PowExp"; //3+3
+    break;
   }
 
   TString filename=Form("%sMassFit.root",type.Data());
@@ -1603,7 +1704,7 @@ void AliHFMassFitter::PlotFit(TVirtualPad* pd,Double_t nsigma,Int_t writeFitInfo
       TF1* f=hdraw->GetFunction("funcbkgonly");
       for (Int_t i=0;i<fNFinalPars-3;i++){
        pinfo->SetTextColor(kBlue+3);
-       TString str=Form("%s = %f #pm %f",f->GetParName(i),f->GetParameter(i),f->GetParError(i));
+       TString str=Form("%s = %.3f #pm %.3f",f->GetParName(i),f->GetParameter(i),f->GetParError(i));
        pinfo->AddText(str);
       }
 
@@ -1622,8 +1723,8 @@ void AliHFMassFitter::PlotFit(TVirtualPad* pd,Double_t nsigma,Int_t writeFitInfo
   pd->cd();
   hdraw->SetMarkerStyle(20);
   hdraw->DrawClone("PE");
-  if(hdraw->GetFunction("funcbkgFullRange")) hdraw->GetFunction("funcbkgFullRange")->DrawClone("same");
-  if(hdraw->GetFunction("funcbkgRecalc")) hdraw->GetFunction("funcbkgRecalc")->DrawClone("same");
+//   if(hdraw->GetFunction("funcbkgFullRange")) hdraw->GetFunction("funcbkgFullRange")->DrawClone("same");
+//   if(hdraw->GetFunction("funcbkgRecalc")) hdraw->GetFunction("funcbkgRecalc")->DrawClone("same");
   if(hdraw->GetFunction("funcmass")) hdraw->GetFunction("funcmass")->DrawClone("same");
 
   if(writeFitInfo > 0){
@@ -1637,7 +1738,7 @@ void AliHFMassFitter::PlotFit(TVirtualPad* pd,Double_t nsigma,Int_t writeFitInfo
 
     for (Int_t i=fNFinalPars-3;i<fNFinalPars;i++){
       pinfom->SetTextColor(kBlue);
-      TString str=Form("%s = %f #pm %f",ff->GetParName(i),ff->GetParameter(i),ff->GetParError(i));
+      TString str=Form("%s = %.3f #pm %.3f",ff->GetParName(i),ff->GetParameter(i),ff->GetParError(i));
       if(!(writeFitInfo==1 && i==fNFinalPars-3)) pinfom->AddText(str);
     }
     pd->cd();
@@ -1906,6 +2007,7 @@ void AliHFMassFitter::Significance(Double_t min, Double_t max, Double_t &signifi
     cout<<"Cannot calculate significance because of div by 0!"<<endl;
     significance=-1;
     errsignificance=0;
+    return;
   }
 
   significance =  signal/TMath::Sqrt(signal+background);