Updates (Chiara)

diff --git a/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx
index 26abd0c69ca78c58455681bd61ebbae6e5b38664..9e681e96ef10ee7c9c4d560daf18ac4b59e647ee 100644 (file)
--- a/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx
+++ b/PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx
@@ -80,6 +80,7 @@ fVHFmycuts(0),
 fArray(0),
 fReadMC(0),
 fLsNormalization(1.)
+
 {
   // Default constructor
   for(Int_t i=0;i<5;i++) {fTotPosPairs[i]=0; fTotNegPairs[i]=0;}
@@ -311,8 +312,9 @@ void AliAnalysisTaskSED0Mass::UserCreateOutputObjects()
     tmpS27l->Sumw2();
  
     //distribution w/o cuts
-    TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,0,3.73);
-    TH1F *tmpMB=(TH1F*)tmpMt->Clone();
+    //   TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,0.7,3.);
+    TH1F* tmpMS = new TH1F(nameMassNocutsS.Data(),"D^{0} invariant mass; M [GeV]; Entries",300,1.56484,2.16484); //range (MD0-300MeV, mD0 + 300MeV)
+    TH1F *tmpMB=(TH1F*)tmpMS->Clone();
     tmpMB->SetName(nameMassNocutsB.Data());
     tmpMS->Sumw2();
     tmpMB->Sumw2();
@@ -704,7 +706,7 @@ void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
 
          }
          
-      }
+       }
 
        if(pt>0. && pt<=1.){
          ((TH1F*)fDistr->FindObject("hd0d0S_1"))->Fill(d->Prodd0d0());
@@ -864,10 +866,12 @@ void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
 //     printf("    cosThetaPoint    > %f\n",fD0toKpiCuts[8]);
   
     Int_t ptbin=0;
+    Bool_t isInRange=kFALSE;
     
     //cout<<"P_t = "<<pt<<endl;
     if (pt>0. && pt<=1.) {
       ptbin=1;
+      isInRange=kTRUE;
       /*
       //test d0 cut
       fVHFPPR->SetD0toKpiCuts(0.7,0.04,0.8,0.5,0.5,0.05,0.1,-0.0002,0.5);
@@ -886,6 +890,7 @@ void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
     if(pt>1. && pt<=3.) {
       if(pt>1. && pt<=2.) ptbin=2;  
       if(pt>2. && pt<=3.) ptbin=3;  
+      isInRange=kTRUE;
       /*
       //test d0 cut
       fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.0002,0.6);
@@ -899,20 +904,22 @@ void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
     }
  
     if(pt>3. && pt<=5.){
-       ptbin=4;  
-       /*
-       //test d0 cut
-       fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.0001,0.8);
-       fVHFmycuts->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.00015,0.8);
-       */
-       //original
-       fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.05,-0.0001,0.8);
-       fVHFmycuts->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.05,-0.00015,0.8);
-       
-       //printf("I'm in the bin %d\n",ptbin);
+      ptbin=4;  
+      isInRange=kTRUE;
+      /*
+      //test d0 cut
+      fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.0001,0.8);
+      fVHFmycuts->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.00015,0.8);
+      */
+      //original
+      fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.05,-0.0001,0.8);
+      fVHFmycuts->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.05,-0.00015,0.8);
+      
+      //printf("I'm in the bin %d\n",ptbin);
     }
-    if(pt>5.){
+    if(pt>5.&& pt<=10.){ //additional upper limit to compare with Correction Framework
       ptbin=5;
+      isInRange=kTRUE;
       /*
       //test d0 cut
       fVHFPPR->SetD0toKpiCuts(0.7,0.02,0.8,0.7,0.7,0.05,0.1,-0.00005,0.8);
@@ -926,7 +933,7 @@ void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
     //printf("I'm in the bin %d\n",ptbin);
     //old
     //fVHF->SetD0toKpiCuts(0.7,0.03,0.8,0.06,0.06,0.05,0.05,-0.0002,0.6); //2.p-p vertex reconstructed    
-    if(prongsinacc){
+    if(prongsinacc && isInRange){
       FillHists(ptbin,d,mcArray,fVHFPPR,fOutputPPR);
       FillHists(ptbin,d,mcArray,fVHFmycuts,fOutputmycuts);
     }
@@ -972,7 +979,7 @@ void AliAnalysisTaskSED0Mass::FillHists(Int_t ptbin, AliAODRecoDecayHF2Prong *pa
     //count candidates selected by cuts
     fNentries->Fill(2);
     //count true D0 selected by cuts
-    if (labD0>=0) fNentries->Fill(3);
+    if (fReadMC && labD0>=0) fNentries->Fill(3);
     PostData(3,fNentries);
 
     if (okD0==1) {
@@ -1056,6 +1063,7 @@ void AliAnalysisTaskSED0Mass::Terminate(Option_t */*option*/)
   //
   if(fDebug > 1) printf("AnalysisTaskSED0Mass: Terminate() \n");
 
+
   fOutputPPR = dynamic_cast<TList*> (GetOutputData(1));
   if (!fOutputPPR) {     
     printf("ERROR: fOutputthight not available\n");

diff --git a/PWG3/vertexingHF/AliHFMassFitter.cxx b/PWG3/vertexingHF/AliHFMassFitter.cxx
index 5b97517c4b2f3753e9cee6829d0002f38c432944..3c0618c7027be9f95f94b79548895e4fb17e8710 100644 (file)
--- a/PWG3/vertexingHF/AliHFMassFitter.cxx
+++ b/PWG3/vertexingHF/AliHFMassFitter.cxx
@@ -86,6 +86,7 @@ AliHFMassFitter::AliHFMassFitter (TH1F *histoToFit, Double_t minvalue, Double_t
   // standard constructor
 
   fhistoInvMass= (TH1F*)histoToFit->Clone("fhistoInvMass");
+  fhistoInvMass->SetDirectory(0);
   fminMass=minvalue; 
   fmaxMass=maxvalue;
   if(rebin!=1) RebinMass(rebin); 
@@ -234,6 +235,7 @@ void AliHFMassFitter::ComputeParSize() {
 void AliHFMassFitter::SetHisto(TH1F *histoToFit){
   //fhistoInvMass = (TH1F*)histoToFit->Clone();
   fhistoInvMass = new TH1F(*histoToFit);
+  fhistoInvMass->SetDirectory(0);
   cout<<"SetHisto pointer "<<fhistoInvMass<<endl;
 }
 
@@ -256,7 +258,7 @@ void AliHFMassFitter::Reset() {
   cout<<"Reset "<<fhistoInvMass<<endl;
   if(fhistoInvMass) {
     //cout<<"esiste"<<endl;
-    //delete fhistoInvMass;
+    delete fhistoInvMass;
     fhistoInvMass=NULL;
     cout<<fhistoInvMass<<endl;
   }
@@ -747,7 +749,7 @@ 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<<"Primo fit: \nintbkg1 = "<<intbkg1<<"\t(Compare with par0 = "<<funcbkg->GetParameter(0)<<")\nslope1= "<<slope1<<"\nconc1 = "<<conc1<<endl;
+    cout<<"First fit: \nintbkg1 = "<<intbkg1<<"\t(Compare with par0 = "<<funcbkg->GetParameter(0)<<")\nslope1= "<<slope1<<"\nconc1 = "<<conc1<<endl;
   } 
   else cout<<"\t\t//"<<endl;
   
@@ -768,8 +770,9 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
       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;
+
+      //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
        {
@@ -784,8 +787,12 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
          funcbkg1->SetParameters(0.5*(totInt-intbkg1),fMass,ffactor*fSigmaSgn,intbkg1,slope1);
        }
     }
-    fhistoInvMass->Fit(bkg1name.Data(),"R,L,E,+,0");
-  
+    Int_t status=fhistoInvMass->Fit(bkg1name.Data(),"R,L,E,+,0");
+    if (status != 0){
+      cout<<"Minuit returned "<<status<<endl;
+      return kFALSE;
+    }
+
     for(Int_t i=0;i<bkgPar;i++){
       fFitPars[bkgPar-3+i]=funcbkg1->GetParameter(i);
       //cout<<bkgPar-3+i<<"\t"<<funcbkg1->GetParameter(i);
@@ -829,7 +836,10 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
   Double_t sgnInt;
   sgnInt = totInt-bkgInt;
   cout<<"------TotInt = "<<totInt<<"\tsgnInt = "<<sgnInt<<endl;
-
+  if (sgnInt <= 0){
+    cout<<"Setting sgnInt = - sgnInt"<<endl;
+    sgnInt=- sgnInt;
+  }
   /*Fit All Mass distribution with exponential + gaussian (+gaussiam braodened) */
   TF1 *funcmass = new TF1(massname.Data(),this,&AliHFMassFitter::FitFunction4MassDistr,fminMass,fmaxMass,nFitPars,"AliHFMassFitter","FitFunction4MassDistr");
   cout<<"Function name = "<<funcmass->GetName()<<endl<<endl;
@@ -842,14 +852,14 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     funcmass->SetParNames("TotInt","Slope","SgnInt","Mean","Sigma");
     funcmass->SetParameters(totInt,slope1,sgnInt,fMass,fSigmaSgn);
     //cout<<"Parameters set to: "<<totInt<<"\t"<<slope1<<"\t"<<sgnInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
-//     cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<nFitPars<<"\tgsidebands = "<<fSideBands<<endl;
+    //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<nFitPars<<"\tgsidebands = "<<fSideBands<<endl;
     funcmass->FixParameter(0,totInt);
   }
   if (nFitPars==6){
     funcmass->SetParNames("TotInt","Coef1","Coef2","SgnInt","Mean","Sigma");
     funcmass->SetParameters(totInt,slope1,conc1,sgnInt,fMass,fSigmaSgn);
-//     cout<<"Parameters set to: "<<totInt<<"\t"<<slope1<<"\t"<<conc1<<"\t"<<sgnInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
-//     cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<nFitPars<<"\tgsidebands = "<<fSideBands<<endl;
+    //cout<<"Parameters set to: "<<totInt<<"\t"<<slope1<<"\t"<<conc1<<"\t"<<sgnInt<<"\t"<<fMass<<"\t"<<fSigmaSgn<<"\t"<<endl;
+    //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<nFitPars<<"\tgsidebands = "<<fSideBands<<endl;
     funcmass->FixParameter(0,totInt);
   }
   if(nFitPars==4){
@@ -861,8 +871,15 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
     //cout<<"Limits: ("<<fminMass<<","<<fmaxMass<<")\tnPar = "<<nFitPars<<"\tgsidebands = "<<fSideBands<<endl;
 
   }
-  
-  fhistoInvMass->Fit(massname.Data(),"R,L,E,+,0");
+
+  Int_t status;
+
+  status = fhistoInvMass->Fit(massname.Data(),"R,L,E,+,0");
+  if (status != 0){
+    cout<<"Minuit returned "<<status<<endl;
+    return kFALSE;
+  }
+
   cout<<"fit done"<<endl;
   //reset value of fMass and fSigmaSgn to those found from fit
   fMass=funcmass->GetParameter(nFitPars-2);
@@ -913,6 +930,7 @@ Bool_t AliHFMassFitter::MassFitter(Bool_t draw){
        for (Int_t i=0; i<2; i++) {
          gMinuit->SetErrorDef(errDef[i]);
          cont[i] = (TGraph*)gMinuit->Contour(80,kpar,jpar);
+         cout<<"Minuit Status = "<<gMinuit->GetStatus()<<endl;
        }
        
        if(!cont[0] || !cont[1]){
@@ -994,6 +1012,15 @@ void  AliHFMassFitter::GetFitPars(Float_t *vector) const {
 }
 
 
+//_________________________________________________________________________
+void AliHFMassFitter::IntS(Float_t *valuewitherror){
+  Int_t index=fParsSize/2 - 3;
+  valuewitherror[0]=fFitPars[index];
+  index=fParsSize - 3;
+  valuewitherror[1]=fFitPars[index];
+  }
+
+
 //_________________________________________________________________________
 void AliHFMassFitter::AddFunctionsToHisto(){
 

diff --git a/PWG3/vertexingHF/AliHFMassFitter.h b/PWG3/vertexingHF/AliHFMassFitter.h
index 6b3c5f32751f937da79cda01ac30027e46dd827d..56918215260a54eb3345f4c6fb35dfa443a67e47 100644 (file)
--- a/PWG3/vertexingHF/AliHFMassFitter.h
+++ b/PWG3/vertexingHF/AliHFMassFitter.h
@@ -40,12 +40,12 @@ class AliHFMassFitter : public TNamed {
   void     SetBinN(Int_t newbinN){fNbin=newbinN;}
   void     SetType(Int_t fittypeb, Int_t fittypes);
   void     SetReflectionSigmaFactor(Int_t constant) {ffactor=constant;}
-  void     SetInitialGaussianMean(Double_t mean) {fMass=mean;}
-  void     SetInitialGaussianSigma(Double_t sigma) {fSigmaSgn=sigma;}
-  void     SetSideBands(Bool_t onlysidebands=kTRUE) {fSideBands=onlysidebands;}
+  void     SetInitialGaussianMean(Double_t mean) {fMass=mean;} // change the default value of the mean
+  void     SetInitialGaussianSigma(Double_t sigma) {fSigmaSgn=sigma;} // change the default value of the sigma
+  void     SetSideBands(Bool_t onlysidebands=kTRUE) {fSideBands=onlysidebands;} // consider only side bands
 
   //getters
-  TH1F*    GetHistoClone() const;
+  TH1F*    GetHistoClone() const; //return the histogram
   void     GetRangeFit(Double_t &minvalue, Double_t &maxvalue) const {minvalue=fminMass; maxvalue=fmaxMass;}
   Double_t GetMinRangeFit()const {return fminMass;}
   Double_t GetMaxRangeFit()const {return fmaxMass;}
@@ -61,21 +61,23 @@ class AliHFMassFitter : public TNamed {
 
   void     PrintParTitles() const;
 
-  void     InitNtuParam(char *ntuname="ntupar");
-  void     FillNtuParam();
-  TNtuple* GetNtuParam() const {return fntuParam;}
-  TNtuple* NtuParamOneShot(char *ntuname="ntupar");
-  void     WriteHisto(TString path="./");
-  void     WriteNtuple(TString path="./") const;
+  void     InitNtuParam(char *ntuname="ntupar"); // initialize TNtuple to store the parameters
+  void     FillNtuParam(); //Fill the TNtuple with the current parameters
+  TNtuple* GetNtuParam() const {return fntuParam;} // return the TNtuple
+  TNtuple* NtuParamOneShot(char *ntuname="ntupar"); // the three functions above all together
+  void     WriteHisto(TString path="./"); // write the histogram
+  void     WriteNtuple(TString path="./") const; // write the TNtuple
   void     DrawFit() const;
   void     Reset();
 
-  void     Signal(Double_t nOfSigma,Double_t &signal,Double_t &errsignal) const; 
-  void     Signal(Double_t min,Double_t max,Double_t &signal,Double_t &errsignal) const; 
-  void     Background(Double_t nOfSigma,Double_t &background,Double_t &errbackground) const; 
-  void     Background(Double_t min,Double_t max,Double_t &background,Double_t &errbackground) const; 
-  void     Significance(Double_t nOfSigma,Double_t &significance,Double_t &errsignificance) const;
-  void     Significance(Double_t min,Double_t max,Double_t &significance,Double_t &errsignificance) const;
+  void     IntS(Float_t *valuewitherror);    // integral of signal given my the fit with error
+  Double_t IntTot(){return fhistoInvMass->Integral("width");}  // return total integral of the histogram
+  void     Signal(Double_t nOfSigma,Double_t &signal,Double_t &errsignal) const; // signal in nsigma with error 
+  void     Signal(Double_t min,Double_t max,Double_t &signal,Double_t &errsignal) const; // signal in (min, max) with error 
+  void     Background(Double_t nOfSigma,Double_t &background,Double_t &errbackground) const; // backgournd in nsigma with error 
+  void     Background(Double_t min,Double_t max,Double_t &background,Double_t &errbackground) const; // backgournd in (min, max) with error 
+  void     Significance(Double_t nOfSigma,Double_t &significance,Double_t &errsignificance) const; // significance in nsigma with error 
+  void     Significance(Double_t min,Double_t max,Double_t &significance,Double_t &errsignificance) const; // significance in (min, max) with error 
 
   Double_t FitFunction4MassDistr (Double_t*, Double_t*);
   Double_t FitFunction4Sgn (Double_t*, Double_t*);