PWG3/vertexingHF/AliAnalysisTaskSEDplus.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /////////////////////////////////////////////////////////////
  17 //
  18 // AliAnalysisTaskSE for the extraction of signal(e.g D+) of heavy flavor
  19 // decay candidates with the MC truth.
  20
  21 /////////////////////////////////////////////////////////////
  22
  23 #include <TClonesArray.h>
  24 #include <TNtuple.h>
  25 #include <TList.h>
  26 #include <TH1F.h>
  27 #include "AliAODEvent.h"
  28 #include "AliAODVertex.h"
  29 #include "AliAODTrack.h"
  30 #include "AliAODMCHeader.h"
  31 #include "AliAODMCParticle.h"
  32 #include "AliAODRecoDecayHF3Prong.h"
  33 #include "AliAnalysisVertexingHF.h"
  34 #include "AliAnalysisTaskSE.h"
  35 #include "AliAnalysisTaskSEDplus.h"
  36
  37 ClassImp(AliAnalysisTaskSEDplus)
  38
  39
  40 //________________________________________________________________________
  41 AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus():
  42 AliAnalysisTaskSE(),
  43 fOutput(0),
  44 fNtupleDplus(0),
  45 fNtupleDplusbackg(0),
  46 fHistMass(0),
  47 fHistSignal(0),
  48 fHistBackground(0),
  49 fVHF(0)
  50 {
  51   // Default constructor
  52 }
  53
  54 //________________________________________________________________________
  55 AliAnalysisTaskSEDplus::AliAnalysisTaskSEDplus(const char *name):
  56 AliAnalysisTaskSE(name),
  57 fOutput(0),
  58 fNtupleDplus(0),
  59 fNtupleDplusbackg(0),
  60 fHistMass(0),
  61 fHistSignal(0),
  62 fHistBackground(0),
  63 fVHF(0)
  64 {
  65   // Default constructor
  66
  67   // Output slot #1 writes into a TList container
  68   DefineOutput(1,TList::Class());  //My private output
  69 }
  70
  71 //________________________________________________________________________
  72 AliAnalysisTaskSEDplus::~AliAnalysisTaskSEDplus()
  73 {
  74   // Destructor
  75   if (fOutput) {
  76     delete fOutput;
  77     fOutput = 0;
  78   }
  79   if (fVHF) {
  80     delete fVHF;
  81     fVHF = 0;
  82   }
  83 }
  84
  85 //________________________________________________________________________
  86 void AliAnalysisTaskSEDplus::Init()
  87 {
  88   // Initialization
  89
  90   if(fDebug > 1) printf("AnalysisTaskSEDplus::Init() \n");
  91
  92   gROOT->LoadMacro("ConfigVertexingHF.C");
  93
  94   fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()");
  95   fVHF->PrintStatus();
  96
  97   return;
  98 }
  99
 100 //________________________________________________________________________
 101 void AliAnalysisTaskSEDplus::UserCreateOutputObjects()
 102 {
 103   // Create the output container
 104   //
 105   if(fDebug > 1) printf("AnalysisTaskSEDplus::UserCreateOutputObjects() \n");
 106
 107   // Several histograms are more conveniently managed in a TList
 108   fOutput = new TList();
 109   fOutput->SetOwner();
 110
 111   fHistMass = new TH1F("fHistMass", "D^{+} invariant mass; M [GeV]; Entries",100,1.765,1.965);
 112   fHistSignal = new TH1F("fHistSignal", "D^{+} invariant mass - MC; M [GeV]; Entries",100,1.765,1.965);
 113   fHistBackground =new TH1F("fHistBackground", "Background invariant mass - MC; M [GeV]; Entries",100,1.765,1.965);
 114
 115
 116
 117
 118   fHistMass->Sumw2();
 119   fHistSignal->Sumw2();
 120   fHistBackground->Sumw2();
 121
 122
 123   //fHistMass->SetMinimum(0);
 124   //fHistSignal->SetMinimum(0);
 125   //fHistBackground->SetMinimum(0);
 126   fOutput->Add(fHistMass);
 127   fOutput->Add(fHistSignal);
 128   fOutput->Add(fHistBackground);
 129
 130
 131   fNtupleDplus = new TNtuple("fNtupleDplus","D +","pdg:Px:Py:Pz:Ptpi:Ptpi2:PtK:PtRec:PtTrue:PointingAngle:DecLeng:VxTrue:VxRec:InvMass:sigvert");
 132   fNtupleDplusbackg = new TNtuple("fNtupleDplusbackg","D + backg","Ptpibkg:Ptpi2bkg:PtKbkg:PtRecbkg:PointingAnglebkg:DLbkg:VxRecbkg:InvMassbkg:sigvertbkg");
 133
 134   fOutput->Add(fNtupleDplus);
 135   fOutput->Add(fNtupleDplusbackg);
 136
 137   return;
 138 }
 139
 140 //________________________________________________________________________
 141 void AliAnalysisTaskSEDplus::UserExec(Option_t */*option*/)
 142 {
 143   // Execute analysis for current event:
 144   // heavy flavor candidates association to MC truth
 145
 146
 147
 148   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
 149
 150   // load Dplus->Kpipi candidates
 151   TClonesArray *array3Prong =
 152     (TClonesArray*)aod->GetList()->FindObject("Charm3Prong");
 153   if(!array3Prong) {
 154     printf("AliAnalysisTaskSEDplus::UserExec: Charm3Prong branch not found!\n");
 155     return;
 156   }
 157
 158   // AOD primary vertex
 159     AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
 160     //    vtx1->Print();
 161
 162   // load MC particles
 163   TClonesArray *arrayMC =
 164     (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 165   if(!arrayMC) {
 166     printf("AliAnalysisTaskSEDplus::UserExec: MC particles branch not found!\n");
 167     return;
 168   }
 169
 170   // load MC header
 171   AliAODMCHeader *mcHeader =
 172     (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 173   if(!mcHeader) {
 174     printf("AliAnalysisTaskSEDplus::UserExec: MC header branch not found!\n");
 175     return;
 176   }
 177     Int_t n3Prong = array3Prong->GetEntriesFast();
 178     printf("Number of D+->Kpipi: %d\n",n3Prong);
 179
 180
 181     for (Int_t i3Prong = 0; i3Prong < n3Prong; i3Prong++) {
 182       AliAODRecoDecayHF3Prong *d = (AliAODRecoDecayHF3Prong*)array3Prong->UncheckedAt(i3Prong);
 183
 184
 185       Bool_t unsetvtx=kFALSE;
 186       if(!d->GetOwnPrimaryVtx()){
 187         d->SetOwnPrimaryVtx(vtx1);
 188         unsetvtx=kTRUE;
 189       }
 190       if(d->SelectDplus(fVHF->GetDplusCuts())) {
 191
 192
 193         Int_t labDp = d->MatchToMC(411,arrayMC);
 194
 195         if(labDp>=0) {
 196           AliAODMCParticle *partDp = (AliAODMCParticle*)arrayMC->At(labDp);
 197           Int_t pdgDp = TMath::Abs(partDp->GetPdgCode());
 198           if(pdgDp==411){
 199
 200             fHistSignal->Fill(d->InvMassDplus());
 201             fHistMass->Fill(d->InvMassDplus());
 202
 203             // Post the data already here
 204             PostData(1,fOutput);
 205
 206             fNtupleDplus->Fill(pdgDp,partDp->Px()-d->Px(),partDp->Py()-d->Py(),partDp->Pz()-d->Pz(),d->PtProng(0),d->PtProng(2),d->PtProng(1),d->Pt(),partDp->Pt(),d->CosPointingAngle(),d->DecayLength(),partDp->Xv(),d->Xv(),d->InvMassDplus(),d->GetSigmaVert());
 207           }
 208         }
 209         else {
 210           fHistBackground->Fill(d->InvMassDplus());
 211           fNtupleDplusbackg->Fill(d->PtProng(0),d->PtProng(2),d->PtProng(1),d->Pt(),d->CosPointingAngle(),d->DecayLength(),d->Xv(),d->InvMassDplus(),d->GetSigmaVert());
 212
 213           fHistMass->Fill(d->InvMassDplus());
 214         }
 215
 216
 217       }
 218
 219
 220
 221       if(unsetvtx) d->UnsetOwnPrimaryVtx();
 222
 223     }
 224      // end loop on D+->Kpipi
 225
 226
 227   return;
 228 }
 229
 230 //________________________________________________________________________
 231 void AliAnalysisTaskSEDplus::Terminate(Option_t */*option*/)
 232 {
 233   // Terminate analysis
 234   //
 235   if(fDebug > 1) printf("AnalysisTaskSEDplus: Terminate() \n");
 236
 237   fOutput = dynamic_cast<TList*> (GetOutputData(1));
 238   if (!fOutput) {
 239     printf("ERROR: fOutput not available\n");
 240     return;
 241   }
 242
 243   fNtupleDplus = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleDplus"));
 244   fHistMass = dynamic_cast<TH1F*>(fOutput->FindObject("fHistMass"));
 245   fHistSignal = dynamic_cast<TH1F*>(fOutput->FindObject("fHistSignal"));
 246   fHistBackground = dynamic_cast<TH1F*>(fOutput->FindObject("fHistBackground"));
 247   fNtupleDplusbackg = dynamic_cast<TNtuple*>(fOutput->FindObject("fNtupleDplusbackg"));
 248     return;
 249 }
 250