Merged tasks DStar and DStarSpectra (Alessandro, Yifei)

[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / AliAnalysisTaskSECompareHF.cxx

diff --git a/PWG3/vertexingHF/AliAnalysisTaskSECompareHF.cxx b/PWG3/vertexingHF/AliAnalysisTaskSECompareHF.cxx
index 2601aca43cc59cc025e810b21c3c912101a2b951..cbfdba42500f758e73bc8f01ca824fc221a41d64 100644 (file)
--- a/PWG3/vertexingHF/AliAnalysisTaskSECompareHF.cxx
+++ b/PWG3/vertexingHF/AliAnalysisTaskSECompareHF.cxx
@@ -26,6 +26,8 @@
 #include <TList.h>
 #include <TH1F.h>
 
+#include "AliAnalysisManager.h"
+#include "AliAODHandler.h"
 #include "AliAODEvent.h"
 #include "AliAODVertex.h"
 #include "AliESDtrack.h"
@@ -33,6 +35,8 @@
 #include "AliAODMCHeader.h"
 #include "AliAODMCParticle.h"
 #include "AliAODRecoDecayHF2Prong.h"
+#include "AliAODRecoDecayHF3Prong.h"
+#include "AliAODRecoDecayHF4Prong.h"
 #include "AliAODRecoCascadeHF.h"
 #include "AliAnalysisVertexingHF.h"
 #include "AliAnalysisTaskSE.h"
@@ -122,6 +126,7 @@ void AliAnalysisTaskSECompareHF::UserCreateOutputObjects()
   fHistNEvents->SetMinimum(0);
   fOutput->Add(fHistNEvents);
 
+  OpenFile(2); // 2 is the slot number of the ntuple
   fNtupleCmp = new TNtuple("fNtupleCmp","Charm comparison","pdg:nprongs:VxRec:VxTrue:ErrVx:VyRec:VyTrue:ErrVy:VzRec:VzTrue:ErrVz:Chi2toNDF:PtRec:Mrec:CPta:Prodd0");
 
   return;
@@ -135,37 +140,57 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
 
   
   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
-  
-  fHistNEvents->Fill(0); // count event
-  // Post the data already here
-  PostData(1,fOutput);
 
-  // load HF vertices                     
-  TClonesArray *inputArrayVertices =
-    (TClonesArray*)aod->GetList()->FindObject("VerticesHF");
+  TClonesArray *inputArrayVertices = 0;
+  TClonesArray *inputArrayD0toKpi = 0;
+  TClonesArray *inputArrayDstar = 0;
+
+  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.    
+    aod = dynamic_cast<AliAODEvent*> (AODEvent());
+    // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)
+    // have to taken from the AOD event hold by the AliAODExtension
+    AliAODHandler* aodHandler = (AliAODHandler*) 
+      ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
+    if(aodHandler->GetExtensions()) {
+      AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");
+      AliAODEvent *aodFromExt = ext->GetAOD();
+      // load HF vertices                
+      inputArrayVertices = (TClonesArray*)aodFromExt->GetList()->FindObject("VerticesHF");
+      // load D0->Kpi candidates
+      inputArrayD0toKpi = (TClonesArray*)aodFromExt->GetList()->FindObject("D0toKpi");
+      // load D*+ candidates                                                   
+      inputArrayDstar = (TClonesArray*)aodFromExt->GetList()->FindObject("Dstar");
+    }
+  } else {
+    // load HF vertices                
+    inputArrayVertices = (TClonesArray*)aod->GetList()->FindObject("VerticesHF");
+    // load D0->Kpi candidates                                                 
+    inputArrayD0toKpi = (TClonesArray*)aod->GetList()->FindObject("D0toKpi");
+    // load D*+ candidates                                                   
+    inputArrayDstar = (TClonesArray*)aod->GetList()->FindObject("Dstar");
+  }
+
+
   if(!inputArrayVertices) {
     printf("AliAnalysisTaskSECompareHF::UserExec: Vertices branch not found!\n");
     return;
   }
-
-  // load D0->Kpi candidates                                                   
-  TClonesArray *inputArrayD0toKpi =
-    (TClonesArray*)aod->GetList()->FindObject("D0toKpi");
   if(!inputArrayD0toKpi) {
     printf("AliAnalysisTaskSECompareHF::UserExec: D0toKpi branch not found!\n");
     return;
   }
-
- 
-  // load D*+ candidates                                                   
-  TClonesArray *inputArrayDstar =
-    (TClonesArray*)aod->GetList()->FindObject("Dstar");
   if(!inputArrayDstar) {
     printf("AliAnalysisTaskSECompareHF::UserExec: Dstar branch not found!\n");
     return;
   }
   
 
+  fHistNEvents->Fill(0); // count event
+  // Post the data already here
+  PostData(1,fOutput);
+
   // AOD primary vertex
   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
   //vtx1->Print();
@@ -193,9 +218,13 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
   AliAODRecoDecayHF3Prong *d3=0;
   AliAODRecoDecayHF4Prong *d4=0;
 
+  Int_t pdgDgD0toKpi[2]={321,211};
+  Int_t pdgDgDplustoKpipi[3]={321,211,211};
+  Int_t pdgDgD0toKpipipi[4]={321,211,211,211};
+
   // loop over vertices
   Int_t nVertices = inputArrayVertices->GetEntriesFast();
-  printf("Number of vertices: %d\n",nVertices);
+  if(fDebug>1) printf("Number of vertices: %d\n",nVertices);
 
   for (Int_t iVtx = 0; iVtx < nVertices; iVtx++) {
     AliAODVertex *vtx = (AliAODVertex*)inputArrayVertices->UncheckedAt(iVtx);
@@ -225,7 +254,7 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
       d2 = (AliAODRecoDecayHF2Prong*)vtx->GetParent();
       if(d2->IsLikeSign()) continue;
       if(d2->Charge() != 0) continue; // these are D* 
-      lab = d2->MatchToMC(421,mcArray);
+      lab = d2->MatchToMC(421,mcArray,2,pdgDgD0toKpi);
       if(lab>=0) {
        unsetvtx=kFALSE;
        if(!d2->GetOwnPrimaryVtx()) {
@@ -237,12 +266,12 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
          AliAODMCParticle *dMC = (AliAODMCParticle*)mcArray->At(lab);
          pdg = dMC->GetPdgCode();
          invmass = (pdg==421 ? d2->InvMassD0() : d2->InvMassD0bar());
-         fHistMass->Fill(invmass);
-         // Post the data already here
-         PostData(1,fOutput);
          // get a daughter for true pos of decay vertex
          AliAODMCParticle *dg0MC = (AliAODMCParticle*)mcArray->At(dMC->GetDaughter(0));
          dg0MC->XvYvZv(posTrue);
+         fHistMass->Fill(invmass);
+         // Post the data already here
+         PostData(1,fOutput);
          Float_t tmp[16]={(Float_t)pdg,(Float_t)nprongs,
                           (Float_t)posRec[0],(Float_t)posTrue[0],(Float_t)errx,
                           (Float_t)posRec[1],(Float_t)posTrue[1],(Float_t)erry,
@@ -258,7 +287,7 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
     case 3: // look for D+
       d3 = (AliAODRecoDecayHF3Prong*)vtx->GetParent();
       if(d3->IsLikeSign()) continue;
-      lab = d3->MatchToMC(411,mcArray);
+      lab = d3->MatchToMC(411,mcArray,3,pdgDgDplustoKpipi);
       if(lab>=0) {
        unsetvtx=kFALSE;
        if(!d3->GetOwnPrimaryVtx()) {
@@ -287,7 +316,7 @@ void AliAnalysisTaskSECompareHF::UserExec(Option_t */*option*/)
     case 4: // look for D0->Kpipipi
       d4 = (AliAODRecoDecayHF4Prong*)vtx->GetParent();
       if(d4->IsLikeSign()) continue;
-      lab = d4->MatchToMC(421,mcArray);
+      lab = d4->MatchToMC(421,mcArray,4,pdgDgD0toKpipipi);
       if(lab>=0) {
        unsetvtx=kFALSE;
        if(!d4->GetOwnPrimaryVtx()) {