Update
authordainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 7 Feb 2011 10:49:32 +0000 (10:49 +0000)
committerdainese <dainese@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 7 Feb 2011 10:49:32 +0000 (10:49 +0000)
PWG3/vertexingHF/AliAnalysisTaskSEHFQA.cxx
PWG3/vertexingHF/AliAnalysisTaskSEHFQA.h
PWG3/vertexingHF/macros/AddTaskHFQA.C

index 60ad19386dadd2a2cb1ddab5546f8ecb541b4df2..bba57b05809864668face742f28470a60a1993c0 100644 (file)
 #include "AliAODVertex.h"
 #include "AliAODTrack.h"
 #include "AliAODMCParticle.h"
+#include "AliAODMCHeader.h"
 #include "AliAODRecoDecayHF2Prong.h"
 #include "AliAODRecoCascadeHF.h"
 #include "AliAnalysisVertexingHF.h"
 #include "AliAnalysisTaskSE.h"
+#include "AliCounterCollection.h"
 #include "AliRDHFCuts.h"
 #include "AliRDHFCutsDplustoKpipi.h"
 #include "AliRDHFCutsD0toKpipipi.h"
@@ -63,20 +65,28 @@ AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA():AliAnalysisTaskSE(),
   fNEntries(0x0),
   fOutputPID(0x0),
   fOutputTrack(0x0),
+  fOutputCounters(0x0),
+  fOutputCheckCentrality(0x0),
   fDecayChannel(AliAnalysisTaskSEHFQA::kD0toKpi),
-  fCuts(0x0)
+  fCuts(0x0),
+  fEstimator(AliRDHFCuts::kCentTRK),
+  fReadMC(kFALSE)
 {
   //default constructor
 }
 
 //____________________________________________________________________________
 AliAnalysisTaskSEHFQA::AliAnalysisTaskSEHFQA(const char *name, AliAnalysisTaskSEHFQA::DecChannel ch,AliRDHFCuts* cuts):
-AliAnalysisTaskSE(name),
-fNEntries(0x0),
-fOutputPID(0x0),
-fOutputTrack(0x0),
-fDecayChannel(ch),
-fCuts(0x0)
+  AliAnalysisTaskSE(name),
+  fNEntries(0x0),
+  fOutputPID(0x0),
+  fOutputTrack(0x0),
+  fOutputCounters(0x0),
+  fOutputCheckCentrality(0x0),
+  fDecayChannel(ch),
+  fCuts(0x0),
+  fEstimator(AliRDHFCuts::kCentTRK),
+  fReadMC(kFALSE)
 {
   //constructor
 
@@ -106,28 +116,30 @@ fCuts(0x0)
   case 4:
     DefineOutput(4,AliRDHFCutsD0toKpipipi::Class());  //My private output
     break;
-   case 5:
+  case 5:
     DefineOutput(4,AliRDHFCutsLctopKpi::Class());  //My private output
     break;
- }
+  }
+  // Output slot #5 writes into a TList container (AliCounterCollection)
+  DefineOutput(5,TList::Class());  //My private output
+  // Output slot #6 writes into a TList container (TH1F)
+  DefineOutput(6,TList::Class());  //My private output
 }
 
 //___________________________________________________________________________
 AliAnalysisTaskSEHFQA::~AliAnalysisTaskSEHFQA()
 {
   //destructor
-  if(fNEntries){
-    delete fNEntries;
-    fNEntries=0;
-  }
-  if(fOutputPID){
-    delete fOutputPID;
-    fOutputPID=0;
-  }
-  if(fOutputTrack){
-    delete fOutputTrack;
-    fOutputTrack=0;
-  }
+  delete fNEntries;
+
+  delete fOutputPID;
+
+  delete fOutputTrack;
+
+  delete fOutputCounters;
+
+  delete fOutputCheckCentrality;
+
 }
 
 //___________________________________________________________________________
@@ -197,13 +209,21 @@ void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
 
   //count events
 
-  fNEntries=new TH1F(GetOutputSlot(1)->GetContainer()->GetName(), "Counts the number of events", 6,-0.5,5.5);
+  fNEntries=new TH1F(GetOutputSlot(1)->GetContainer()->GetName(), "Counts the number of events", 9,-0.5,8.5);
   fNEntries->GetXaxis()->SetBinLabel(1,"nEventsAnal");
   fNEntries->GetXaxis()->SetBinLabel(2,"Pile-up Rej");
   fNEntries->GetXaxis()->SetBinLabel(3,"No VertexingHF");
   fNEntries->GetXaxis()->SetBinLabel(4,"nCandidates(AnCuts)");
   fNEntries->GetXaxis()->SetBinLabel(5,"EventsWithGoodVtx");
   fNEntries->GetXaxis()->SetBinLabel(6,"N. of 0SMH");
+  if(fReadMC){
+    fNEntries->GetXaxis()->SetBinLabel(7,"MC Cand from c");
+    fNEntries->GetXaxis()->SetBinLabel(8,"MC Cand from b");
+    fNEntries->GetXaxis()->SetBinLabel(9,"N fakes Trks");
+  } else{
+    fNEntries->GetXaxis()->SetBinLabel(7,"N candidates");
+  }
+
   fNEntries->GetXaxis()->SetNdivisions(1,kFALSE);
 
   //PID
@@ -234,14 +254,24 @@ void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
   hname="hTOFsigmaProton160";
   TH2F* hTOFsigmaProton160=new TH2F(hname.Data(),"(TOFsignal-timep)/160 ps;p[GeV/c];(TOFsignal-time p)/160 ps",200,0.,4.,100,-5,5);
 
-  hname="hTOFsigmaK120";
-  TH2F* hTOFsigmaK120=new TH2F(hname.Data(),"(TOFsignal-timeK)/120 ps;p[GeV/c];(TOFsignal-timeK)/120 ps",200,0.,4.,100,-5,5);
+  hname="hTOFsigmaKSigPid";
+  TH2F* hTOFsigmaKSigPid=new TH2F(hname.Data(),"(TOFsignal-timeK)/tofSigPid;p[GeV/c];(TOFsignal-timeK)/tofSigPid",200,0.,4.,100,-5,5);
+
+  hname="hTOFsigmaPionSigPid";
+  TH2F* hTOFsigmaPionSigPid=new TH2F(hname.Data(),"(TOFsignal-time#pi)/tofSigPid;p[GeV/c];(TOFsignal-time#pi)/tofSigPid",200,0.,4.,100,-5,5);
+
+  hname="hTOFsigmaProtonSigPid";
+  TH2F* hTOFsigmaProtonSigPid=new TH2F(hname.Data(),"(TOFsignal-timep)/tofSigPid;p[GeV/c];(TOFsignal-time p)/tofSigPid",200,0.,4.,100,-5,5);
+
+  hname="hTOFsigPid3sigPion";
+  TH1F* hTOFsigPid3sigPion=new TH1F(hname.Data(),"TOF PID resolution (#pi) [ps]",500,0.,1000.);
 
-  hname="hTOFsigmaPion120";
-  TH2F* hTOFsigmaPion120=new TH2F(hname.Data(),"(TOFsignal-time#pi)/120 ps;p[GeV/c];(TOFsignal-time#pi)/120 ps",200,0.,4.,100,-5,5);
+  hname="hTOFsigPid3sigKaon";
+  TH1F* hTOFsigPid3sigKaon=new TH1F(hname.Data(),"TOF PID resolution (K) [ps]",500,0.,1000.);
+
+  hname="hTOFsigPid3sigProton";
+  TH1F* hTOFsigPid3sigProton=new TH1F(hname.Data(),"TOF PID resolution (p) [ps]",500,0.,1000.);
 
-  hname="hTOFsigmaProton120";
-  TH2F* hTOFsigmaProton120=new TH2F(hname.Data(),"(TOFsignal-timep)/120 ps;p[GeV/c];(TOFsignal-time p)/120 ps",200,0.,4.,100,-5,5);
 
   //TPC pid
   hname="hTPCsig";
@@ -270,10 +300,13 @@ void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
   fOutputPID->Add(hTOFsigmaK160);
   fOutputPID->Add(hTOFsigmaPion160);
   fOutputPID->Add(hTOFsigmaProton160);
-  fOutputPID->Add(hTOFsigmaK120);
-  fOutputPID->Add(hTOFsigmaPion120);
-  fOutputPID->Add(hTOFsigmaProton120);
-  fOutputPID->Add(hTPCsigvsp);
+  fOutputPID->Add(hTOFsigmaKSigPid);
+  fOutputPID->Add(hTOFsigmaPionSigPid);
+  fOutputPID->Add(hTOFsigmaProtonSigPid);
+  fOutputPID->Add(hTOFsigPid3sigPion);
+  fOutputPID->Add(hTOFsigPid3sigKaon);
+  fOutputPID->Add(hTOFsigPid3sigProton);
+ fOutputPID->Add(hTPCsigvsp);
   fOutputPID->Add(hTPCsigvspAC);
   fOutputPID->Add(hTPCsigmaK);
   fOutputPID->Add(hTPCsigmaPion);
@@ -302,12 +335,6 @@ void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
   TH1F* hdistrGoodTr=new TH1F(hname.Data(),"Distribution of number of 'good' tracks per event;no.good-tracks/ev;Entries",4000,-0.5,3999.5);
   hdistrGoodTr->SetTitleOffset(1.3,"Y");
 
-  hname="hNtracklets";
-  TH1F* hNtracklets=new TH1F(hname.Data(),"Number of tracklets;ntracklets;Entries",5000,-0.5,4999.5);
-
-  hname="hMult";
-  TH1F* hMult=new TH1F(hname.Data(),"Multiplicity;multiplicity;Entries",10000,-0.5,9999.5);
-
   hname="hd0";
   TH1F* hd0=new TH1F(hname.Data(),"Impact parameter distribution of 'good' tracks;d_{0}[cm];Entries/10^{3} cm",200,-0.1,0.1);
 
@@ -316,15 +343,65 @@ void AliAnalysisTaskSEHFQA::UserCreateOutputObjects()
   fOutputTrack->Add(hnClsSPD);
   fOutputTrack->Add(hptGoodTr);
   fOutputTrack->Add(hdistrGoodTr);
-  fOutputTrack->Add(hNtracklets);
-  fOutputTrack->Add(hMult);
   fOutputTrack->Add(hd0);
   
+  if(fCuts->GetUseCentrality()){
+
+    //Centrality (Counters)
+    fOutputCounters=new TList();
+    fOutputCounters->SetOwner();
+    fOutputCounters->SetName(GetOutputSlot(5)->GetContainer()->GetName());
+
+    AliCounterCollection *stdEstimator=new AliCounterCollection("stdEstimator");
+    stdEstimator->AddRubric("run",500000);
+    stdEstimator->AddRubric("centralityclass","0_10/10_20/20_30/30_40/40_50/50_60/60_70/70_80/80_90/90_100");
+    stdEstimator->Init();
+    AliCounterCollection *secondEstimator=new AliCounterCollection("secondEstimator");
+    secondEstimator->AddRubric("run",500000);
+    secondEstimator->AddRubric("centralityclass","0_10/10_20/20_30/30_40/40_50/50_60/60_70/70_80/80_90/90_100/-990_-980");
+    secondEstimator->Init();
+
+    fOutputCounters->Add(stdEstimator);
+    fOutputCounters->Add(secondEstimator);
+
+    //Centrality (Checks)
+    fOutputCheckCentrality=new TList();
+    fOutputCheckCentrality->SetOwner();
+    fOutputCheckCentrality->SetName(GetOutputSlot(6)->GetContainer()->GetName());
+
+    hname="hNtrackletsIn";
+    TH1F* hNtrackletsIn=new TH1F(hname.Data(),"Number of tracklets in Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
+
+    hname="hMultIn";
+    TH1F* hMultIn=new TH1F(hname.Data(),"Multiplicity;multiplicity in Centrality range;Entries",10000,-0.5,9999.5);
+
+    hname="hNtrackletsOut";
+    TH1F* hNtrackletsOut=new TH1F(hname.Data(),"Number of tracklets out of Centrality range;ntracklets;Entries",5000,-0.5,4999.5);
+
+    hname="hMultOut";
+    TH1F* hMultOut=new TH1F(hname.Data(),"Multiplicity out of Centrality range;multiplicity;Entries",10000,-0.5,9999.5);
+
+    fOutputCheckCentrality->Add(hNtrackletsIn);
+    fOutputCheckCentrality->Add(hNtrackletsOut);
+    fOutputCheckCentrality->Add(hMultIn);
+    fOutputCheckCentrality->Add(hMultOut);
+  } else{
+    hname="hNtracklets";
+    TH1F* hNtracklets=new TH1F(hname.Data(),"Number of tracklets;ntracklets;Entries",5000,-0.5,4999.5);
+
+    hname="hMult";
+    TH1F* hMult=new TH1F(hname.Data(),"Multiplicity;multiplicity;Entries",10000,-0.5,9999.5);
+    fOutputTrack->Add(hNtracklets);
+    fOutputTrack->Add(hMult);
+  }
+
   // Post the data
   PostData(1,fNEntries);
   PostData(2,fOutputPID);
   PostData(3,fOutputTrack);
   PostData(4,fCuts);
+  PostData(5,fOutputCounters);
+  PostData(6,fOutputCheckCentrality);
 }
 
 //___________________________________________________________________________
@@ -339,9 +416,13 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
   PostData(2,fOutputPID);
   PostData(3,fOutputTrack);
   PostData(4,fCuts);
+  PostData(5,fOutputCounters);
+  PostData(6,fOutputCheckCentrality);
 
   TClonesArray *arrayProng =0;
   Int_t pdg=0;
+  Int_t *pdgdaughters=0x0;
+
   if(!aod && AODEvent() && IsStandardAOD()) { 
     // In case there is an AOD handler writing a standard AOD, use the AOD 
     // event in memory rather than the input (ESD) event.    
@@ -361,26 +442,61 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
       case 0:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
        pdg=411;
+       if(fReadMC){
+         pdgdaughters =new Int_t[3];
+         pdgdaughters[0]=211;//pi
+         pdgdaughters[1]=321;//K
+         pdgdaughters[2]=211;//pi
+       }
        break; 
       case 1:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
        pdg=421;
+       if(fReadMC){
+         pdgdaughters =new Int_t[2];
+         pdgdaughters[0]=211;//pi 
+         pdgdaughters[1]=321;//K
+       }
        break; 
       case 2:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
        pdg=413;
+       if(fReadMC){
+         pdgdaughters[1]=211;//pi
+         pdgdaughters[0]=321;//K
+         pdgdaughters[2]=211;//pi (soft?)
+       }
        break; 
       case 3:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
        pdg=431;
+       if(fReadMC){
+         pdgdaughters =new Int_t[3];
+         pdgdaughters[0]=321;//K
+         pdgdaughters[1]=321;//K
+         pdgdaughters[2]=211;//pi
+       }
        break; 
       case 4:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm4Prong");
        pdg=421;
+       if(fReadMC){
+         pdgdaughters =new Int_t[4];
+         pdgdaughters[0]=321;
+         pdgdaughters[1]=211;
+         pdgdaughters[2]=211;
+         pdgdaughters[3]=211;
+       }
        break; 
       case 5:
        arrayProng=(TClonesArray*)aodFromExt->GetList()->FindObject("Charm3Prong");
        pdg=4122;
+       if(fReadMC){
+         pdgdaughters =new Int_t[3];
+         pdgdaughters[0]=2212;//p
+         pdgdaughters[1]=321;//K
+         pdgdaughters[2]=211;//pi
+       }
        break; 
       }
     }
@@ -389,26 +505,61 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
     case 0:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
       pdg=411;
+      if(fReadMC){
+       pdgdaughters =new Int_t[3];
+       pdgdaughters[0]=211;//pi
+       pdgdaughters[1]=321;//K
+       pdgdaughters[2]=211;//pi
+      }
       break; 
     case 1:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("D0toKpi");
       pdg=421;
+      if(fReadMC){
+       pdgdaughters =new Int_t[2];
+       pdgdaughters[0]=211;//pi 
+       pdgdaughters[1]=321;//K
+      }
       break; 
     case 2:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Dstar");
       pdg=413;
+      if(fReadMC){
+       pdgdaughters[1]=211;//pi
+       pdgdaughters[0]=321;//K
+       pdgdaughters[2]=211;//pi (soft?)
+      }
       break; 
     case 3:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
       pdg=431;
+      if(fReadMC){
+       pdgdaughters =new Int_t[3];
+       pdgdaughters[0]=321;//K
+       pdgdaughters[1]=321;//K
+       pdgdaughters[2]=211;//pi
+      }
       break; 
     case 4:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm4Prong");
       pdg=421;
+      if(fReadMC){
+       pdgdaughters =new Int_t[4];
+       pdgdaughters[0]=321;
+       pdgdaughters[1]=211;
+       pdgdaughters[2]=211;
+       pdgdaughters[3]=211;
+      }
       break; 
     case 5:
       arrayProng=(TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
       pdg=4122;
+      if(fReadMC){
+       pdgdaughters =new Int_t[3];
+       pdgdaughters[0]=2212;//p
+       pdgdaughters[1]=321;//K
+       pdgdaughters[2]=211;//pi
+      }
       break; 
     }
   }
@@ -419,6 +570,25 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
     fNEntries->Fill(2);
   }
   
+  TClonesArray *mcArray = 0;
+  AliAODMCHeader *mcHeader = 0;
+
+  //check if MC
+  if(fReadMC) {
+    // load MC particles
+    mcArray = (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
+    if(!mcArray) {
+      printf("AliAnalysisTaskSEHFQA::UserExec: MC particles branch not found!\n");
+      return;
+    }
+    
+    // load MC header
+    mcHeader = (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
+    if(!mcHeader) {
+      printf("AliAnalysisTaskSEHFQA::UserExec: MC header branch not found!\n");
+      return;
+    }
+  }
   if(!aod) return;
   // fix for temporary bug in ESDfilter 
   // the AODs with null vertex pointer didn't pass the PhysSel
@@ -445,16 +615,32 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
     return;
   }
 
+  if(fCuts->GetUseCentrality()){
+    Int_t runNumber = aod->GetRunNumber();
+    Int_t stdCent = (Int_t)(fCuts->GetCentrality(aod)+0.5);
+    Int_t secondCent = (Int_t)(fCuts->GetCentrality(aod,fEstimator)+0.5);
+    Int_t mincent=stdCent-stdCent%10;
+    ((AliCounterCollection*)fOutputCounters->FindObject("stdEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
+    mincent=secondCent-secondCent%10;
+    ((AliCounterCollection*)fOutputCounters->FindObject("secondEstimator"))->Count(Form("centralityclass:%d_%d/Run:%d",mincent,mincent+10,runNumber));
+
+    if(stdCent<fCuts->GetMinCentrality() || stdCent>fCuts->GetMaxCentrality()){
+      ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsOut"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
+      ((TH1F*)fOutputCheckCentrality->FindObject("hMultOut"))->Fill(aod->GetHeader()->GetRefMultiplicity());
+    }else{
+      ((TH1F*)fOutputCheckCentrality->FindObject("hNtrackletsIn"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
+      ((TH1F*)fOutputCheckCentrality->FindObject("hMultIn"))->Fill(aod->GetHeader()->GetRefMultiplicity());
+    }
+  } else{
+    ((TH1F*)fOutputTrack->FindObject("hNtracklets"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
+    ((TH1F*)fOutputTrack->FindObject("hMult"))->Fill(aod->GetHeader()->GetRefMultiplicity());
+  }
+
   Int_t ntracks=0;
   Int_t isGoodTrack=0;
 
   if(aod) ntracks=aod->GetNTracks();
 
-
-  ((TH1F*)fOutputTrack->FindObject("hNtracklets"))->Fill(aod->GetTracklets()->GetNumberOfTracklets());
-  ((TH1F*)fOutputTrack->FindObject("hMult"))->Fill(aod->GetHeader()->GetRefMultiplicity());
-
-
   //loop on tracks in the event
   for (Int_t k=0;k<ntracks;k++){
     AliAODTrack* track=aod->GetTrack(k);
@@ -465,6 +651,9 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
     Double_t times[AliPID::kSPECIES];
     pid->GetIntegratedTimes(times);
     
+    Double_t tofRes[AliPID::kSPECIES];
+    pid->GetTOFpidResolution(tofRes);
+
     //check TOF
     if(pidHF && pidHF->CheckStatus(track,"TOF")){
       ((TH1F*)fOutputPID->FindObject("hTOFtime"))->Fill(times[AliPID::kProton]);
@@ -472,13 +661,34 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
       ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(pid->GetTOFsignal());
       if (pid->GetTOFsignal()< 0) ((TH1F*)fOutputPID->FindObject("hTOFsig"))->Fill(-1);
 
+      // test a "simple" 160 ps TOF sigma PID
       ((TH2F*)fOutputPID->FindObject("hTOFsigmaK160"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kKaon])/160);
       ((TH2F*)fOutputPID->FindObject("hTOFsigmaPion160"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kPion])/160);
       ((TH2F*)fOutputPID->FindObject("hTOFsigmaProton160"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kProton])/160);
-      ((TH2F*)fOutputPID->FindObject("hTOFsigmaK120"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kKaon])/120);
-      ((TH2F*)fOutputPID->FindObject("hTOFsigmaPion120"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kPion])/120);
-      ((TH2F*)fOutputPID->FindObject("hTOFsigmaProton120"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kProton])/120);
-
+      
+      // test TOF sigma PID
+      if (tofRes[2] != 0.) {   // protection against 'old' AODs...
+       ((TH2F*)fOutputPID->FindObject("hTOFsigmaKSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kKaon])/tofRes[3]);
+       ((TH2F*)fOutputPID->FindObject("hTOFsigmaPionSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kPion])/tofRes[2]);
+       ((TH2F*)fOutputPID->FindObject("hTOFsigmaProtonSigPid"))->Fill(track->P(),(pid->GetTOFsignal()-times[AliPID::kProton])/tofRes[4]);
+       for (Int_t iS=2; iS<5; iS++){ //we plot TOF Pid resolution for 3-sigma identified particles
+         if ( (TMath::Abs(times[iS]-pid->GetTOFsignal())/tofRes[iS])<3.){
+           switch (iS) {
+           case AliPID::kPion:
+             ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigPion"))->Fill(tofRes[iS]);
+             break;
+           case AliPID::kKaon:
+             ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigKaon"))->Fill(tofRes[iS]);
+             break;
+           case AliPID::kProton:
+             ((TH1F*)fOutputPID->FindObject("hTOFsigPid3sigProton"))->Fill(tofRes[iS]);
+             break;
+           default:
+             break;
+           }
+         }
+       }
+      }
     }//if TOF status
 
     if(pidHF && pidHF->CheckStatus(track,"TPC")){ 
@@ -503,8 +713,11 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
       Double_t TPCsignal=pid->GetTPCsignal();
       ((TH1F*)fOutputPID->FindObject("hTPCsig"))->Fill(TPCsignal);
       ((TH1F*)fOutputPID->FindObject("hTPCsigvsp"))->Fill(TPCp,TPCsignal);
+      //if (pidHF->IsKaonRaw(track, "TOF"))
       ((TH2F*)fOutputPID->FindObject("hTPCsigmaK"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kKaon));
-      ((TH2F*)fOutputPID->FindObject("hTPCsigmaPion"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kPion));
+      //if (pidHF->IsPionRaw(track, "TOF"))
+         ((TH2F*)fOutputPID->FindObject("hTPCsigmaPion"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kPion));
+      //if (pidHF->IsProtonRaw(track,"TOF"))
       ((TH2F*)fOutputPID->FindObject("hTPCsigmaProton"))->Fill(TPCp,tpcres->GetNumberOfSigmas(TPCp,TPCsignal,track->GetTPCNcls(),AliPID::kProton));
       delete tpcres;
     }//if TPC status
@@ -543,49 +756,77 @@ void AliAnalysisTaskSEHFQA::UserExec(Option_t */*option*/)
   } //end loop on tracks
 
   if(!isSimpleMode){
-  // loop over candidates
-  Int_t nCand = arrayProng->GetEntriesFast();
-  Int_t ndaugh=3;
-  if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi) ndaugh=2;
-  if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpipipi) ndaugh=4;
-
-  for (Int_t iCand = 0; iCand < nCand; iCand++) {
-    AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayProng->UncheckedAt(iCand);
+    // loop over candidates
+    Int_t nCand = arrayProng->GetEntriesFast();
+    Int_t ndaugh=3;
+    if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpi) ndaugh=2;
+    if(fDecayChannel==AliAnalysisTaskSEHFQA::kD0toKpipipi) ndaugh=4;
+
+    for (Int_t iCand = 0; iCand < nCand; iCand++) {
+      AliAODRecoDecayHF *d = (AliAODRecoDecayHF*)arrayProng->UncheckedAt(iCand);
+
+      if(fReadMC){ 
+       Int_t labD = d->MatchToMC(pdg,mcArray,ndaugh,pdgdaughters);
+       if(labD>=0){
+         AliAODMCParticle *partD = (AliAODMCParticle*)mcArray->At(labD);
+         Int_t label=partD->GetMother();
+         AliAODMCParticle *mot = (AliAODMCParticle*)mcArray->At(label);
+         while(label>=0){//get first mother
+           mot = (AliAODMCParticle*)mcArray->At(label);
+           label=mot->GetMother();
+         }
+         Int_t pdgMotCode = mot->GetPdgCode();
+       
+         if(TMath::Abs(pdgMotCode)==4) fNEntries->Fill(6); //from primary charm
+         if(TMath::Abs(pdgMotCode)==5) fNEntries->Fill(7); //from beauty
 
-    for(Int_t id=0;id<ndaugh;id++){
+       }
+      }//end MC
+      else fNEntries->Fill(6); //count the candidates (data)
 
-      //other histograms to be filled when the cut object is given
-      AliAODTrack* track=(AliAODTrack*)d->GetDaughter(id);
+      for(Int_t id=0;id<ndaugh;id++){
 
-      //track quality
+       //other histograms to be filled when the cut object is given
+       AliAODTrack* track=(AliAODTrack*)d->GetDaughter(id);
+       if(fReadMC){
+         Int_t label=track->GetLabel();
+         if (label<0)fNEntries->Fill(8);
+       }
+       //track quality
 
-      if (fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(pdg)) && fCuts->IsSelected(d,AliRDHFCuts::kTracks,aod)) {
+       if (fCuts->IsInFiducialAcceptance(d->Pt(),d->Y(pdg)) && fCuts->IsSelected(d,AliRDHFCuts::kTracks,aod)) {
        
-       ((TH1F*)fOutputTrack->FindObject("hptGoodTr"))->Fill(track->Pt());
-       isGoodTrack++;
+         ((TH1F*)fOutputTrack->FindObject("hptGoodTr"))->Fill(track->Pt());
+         isGoodTrack++;
       
-       ((TH1F*)fOutputTrack->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
+         ((TH1F*)fOutputTrack->FindObject("hdistrGoodTr"))->Fill(isGoodTrack);
       
-       ((TH1F*)fOutputTrack->FindObject("hd0"))->Fill(d->Getd0Prong(id));
+         ((TH1F*)fOutputTrack->FindObject("hd0"))->Fill(d->Getd0Prong(id));
   
-       if (fCuts->IsSelected(d,AliRDHFCuts::kAll,aod)){
+         if (fCuts->IsSelected(d,AliRDHFCuts::kAll,aod)){
          
-         AliAODPid *pid = track->GetDetPid();
-         AliAODPidHF* pidHF=fCuts->GetPidHF();
-         Double_t times[5];
-         pid->GetIntegratedTimes(times);
-         if(pidHF && pidHF->CheckStatus(track,"TOF")) ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptimeAC"))->Fill(track->P(),pid->GetTOFsignal()-times[AliPID::kKaon]);
-         if(pidHF && pidHF->CheckStatus(track,"TPC")) ((TH2F*)fOutputPID->FindObject("hTPCsigvspAC"))->Fill(pid->GetTPCmomentum(),pid->GetTPCsignal());
-
-         fNEntries->Fill(3);
-       } //end analysis cuts
-      } //end acceptance and track cuts
-    } //end loop on tracks in the candidate
-  } //end loop on candidates  
+           AliAODPid *pid = track->GetDetPid();
+           AliAODPidHF* pidHF=fCuts->GetPidHF();
+           Double_t times[5];
+           pid->GetIntegratedTimes(times);
+           if(pidHF && pidHF->CheckStatus(track,"TOF")) ((TH2F*)fOutputPID->FindObject("hTOFtimeKaonHyptimeAC"))->Fill(track->P(),pid->GetTOFsignal()-times[AliPID::kKaon]);
+           if(pidHF && pidHF->CheckStatus(track,"TPC")) ((TH2F*)fOutputPID->FindObject("hTPCsigvspAC"))->Fill(pid->GetTPCmomentum(),pid->GetTPCsignal());
+
+           fNEntries->Fill(3);
+         } //end analysis cuts
+       } //end acceptance and track cuts
+      } //end loop on tracks in the candidate
+    } //end loop on candidates  
   }
+  PostData(1,fNEntries);
+  PostData(2,fOutputPID);
+  PostData(3,fOutputTrack);
+  PostData(4,fCuts);
+  PostData(5,fOutputCounters);
+  PostData(6,fOutputCheckCentrality);
 }
 
-  //____________________________________________________________________________
+//____________________________________________________________________________
 void AliAnalysisTaskSEHFQA::Terminate(Option_t */*option*/){
   //terminate analysis
 
index 699fad763f0704c0ba4d9382f6305e3a537888aa..938edae6955e5ea591858694fdc557340a5649ae 100644 (file)
@@ -37,7 +37,7 @@ class AliAnalysisTaskSEHFQA : public AliAnalysisTaskSE
   virtual void Terminate(Option_t *option);
 
   //setters
-
+  void SetReadMC(Bool_t mcflag){fReadMC=mcflag;}
 
   //getters
   AliRDHFCuts* GetCutObject() const {return fCuts;}
@@ -50,10 +50,14 @@ class AliAnalysisTaskSEHFQA : public AliAnalysisTaskSE
  TH1F*  fNEntries;         //! histogram with number of events on output slot 1
  TList* fOutputPID;        //! list sent on output slot 2
  TList* fOutputTrack;      //! list sent on output slot 3
+ TList* fOutputCounters;   //! list sent on output slot 5
+ TList* fOutputCheckCentrality;   //! list sent on output slot 6
  DecChannel fDecayChannel; //identify the decay channel
  AliRDHFCuts* fCuts;       // object containing cuts
+ AliRDHFCuts::ECentrality fEstimator; //2nd estimator for centrality
+ Bool_t fReadMC;           // flag to read MC
 
- ClassDef(AliAnalysisTaskSEHFQA,1); //AnalysisTaskSE for the quality assurance of HF in hadrons
+ ClassDef(AliAnalysisTaskSEHFQA,3); //AnalysisTaskSE for the quality assurance of HF in hadrons
 
 };
 
index e2057cdd557fa0ff0372befa8f13f51f96dfad0b..a9116cf5c147beb076885463740c51a5f9ed2bcf 100644 (file)
@@ -1,4 +1,4 @@
-AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString filecutsname="D0toKpiCuts.root"){
+AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString filecutsname="D0toKpiCuts.root",Bool_t readMC=kFALSE){
   //
   // Test macro for the AliAnalysisTaskSE for HF mesons quality assurance
   //Author: C.Bianchin chiara.bianchin@pd.infn.it
@@ -18,7 +18,7 @@ AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString
 
   AliRDHFCuts *analysiscuts=0x0;
 
-  TString filename="",out1name="nEntriesQA",out2name="outputPid",out3name="outputTrack",out4name="cuts",inname="input",suffix="",cutsobjname="";
+  TString filename="",out1name="nEntriesQA",out2name="outputPid",out3name="outputTrack",out4name="cuts",out5name="countersCentrality",out6name="outputCentrCheck",inname="input",suffix="",cutsobjname="",centr="";
   filename = AliAnalysisManager::GetCommonFileName();
   filename += ":PWG3_D2H_QA";
 
@@ -27,7 +27,7 @@ AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString
     cutsobjname="AnalysisCuts";
     if(stdcuts) {
       analysiscuts = new AliRDHFCutsDplustoKpipi();
-      analysiscuts->SetStandardCutsPP2010();
+      analysiscuts->SetStandardCutsPbPb2010();
     }
     else analysiscuts = (AliRDHFCutsDplustoKpipi*)filecuts->Get(cutsobjname);
     suffix="Dplus";
@@ -36,7 +36,7 @@ AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString
     cutsobjname="D0toKpiCuts";
     if(stdcuts) {
       analysiscuts = new AliRDHFCutsD0toKpi();
-      analysiscuts->SetStandardCutsPP2010();
+      analysiscuts->SetStandardCutsPbPb2010();
     }
     else analysiscuts = (AliRDHFCutsD0toKpi*)filecuts->Get(cutsobjname);
     suffix="D0";
@@ -84,14 +84,27 @@ AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString
   out2name+=suffix;
   out3name+=suffix;
   out4name=cutsobjname;
+  out5name+=suffix;
+  out6name+=suffix;
 
   if(!analysiscuts && filecutsname!="none"){
     cout<<"Specific AliRDHFCuts not found"<<endl;
     return;
   }
+
+  centr=Form("%.0f%.0f",analysiscuts->GetMinCentrality(),analysiscuts->GetMaxCentrality());
+  inname+=centr;
+  out1name+=centr;
+  out2name+=centr;
+  out3name+=centr;
+  out4name+=centr;
+  out5name+=centr;
+  out6name+=centr;
+
  
   AliAnalysisTaskSEHFQA* taskQA=new AliAnalysisTaskSEHFQA(Form("QA%s",suffix.Data()),ch,analysiscuts);
 
+  taskQA->SetReadMC(readMC);
   mgr->AddTask(taskQA);
 
   //
@@ -111,5 +124,11 @@ AliAnalysisTaskSEHFQA* AddTaskHFQA(AliAnalysisTaskSEHFQA::DecChannel ch,TString
   AliAnalysisDataContainer *coutput4 = mgr->CreateContainer(out4name,AliRDHFCuts::Class(),AliAnalysisManager::kOutputContainer, filename.Data()); //cuts
   mgr->ConnectOutput(taskQA,4,coutput4);
 
+  AliAnalysisDataContainer *coutput5 = mgr->CreateContainer(out5name,TList::Class(),AliAnalysisManager::kOutputContainer, filename.Data()); //quality of tracks
+  mgr->ConnectOutput(taskQA,5,coutput5);
+
+  AliAnalysisDataContainer *coutput6 = mgr->CreateContainer(out6name,TList::Class(),AliAnalysisManager::kOutputContainer, filename.Data()); //quality of tracks
+  mgr->ConnectOutput(taskQA,6,coutput6);
+
  return taskQA;
 }