PWG3/vertexingHF/AliAnalysisTaskSED0Mass.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-2009, 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 D0 candidates invariant mass histogram
  19 // and comparison with the MC truth.
  20 //
  21 // Authors: A.Dainese, andrea.dainese@lnl.infn.it
  22 // and Chiara Bianchin, chiara.bianchin@pd.infn.it
  23 /////////////////////////////////////////////////////////////
  24
  25 #include <Riostream.h>
  26 #include <TClonesArray.h>
  27 #include <TNtuple.h>
  28 #include <TList.h>
  29 #include <TH1F.h>
  30
  31
  32 #include "AliAODEvent.h"
  33 #include "AliAODVertex.h"
  34 #include "AliAODTrack.h"
  35 #include "AliAODMCHeader.h"
  36 #include "AliAODMCParticle.h"
  37 #include "AliAODRecoDecayHF2Prong.h"
  38 #include "AliAODRecoCascadeHF.h"
  39 #include "AliAnalysisVertexingHF.h"
  40 #include "AliAnalysisTaskSE.h"
  41 #include "AliAnalysisTaskSED0Mass.h"
  42
  43 ClassImp(AliAnalysisTaskSED0Mass)
  44
  45
  46 //________________________________________________________________________
  47 AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass():
  48 AliAnalysisTaskSE(),
  49 fOutput(0),
  50 //fNtupleD0Cmp(0),
  51 fhistMass(0),
  52 fhistSgn(0),
  53 fhistBkg(0),
  54 fVHF(0)
  55 {
  56   // Default constructor
  57 }
  58
  59 //________________________________________________________________________
  60 AliAnalysisTaskSED0Mass::AliAnalysisTaskSED0Mass(const char *name):
  61 AliAnalysisTaskSE(name),
  62 fOutput(0),
  63 //fNtupleD0Cmp(0),
  64 fhistMass(0),
  65 fhistSgn(0),
  66 fhistBkg(0),
  67 fVHF(0)
  68 {
  69   // Default constructor
  70
  71   // Output slot #1 writes into a TList container
  72   DefineOutput(1,TList::Class());  //My private output
  73 }
  74
  75 //________________________________________________________________________
  76 AliAnalysisTaskSED0Mass::~AliAnalysisTaskSED0Mass()
  77 {
  78   // Destructor
  79
  80   if (fhistMass){
  81     delete fhistMass;
  82     fhistMass=0;
  83   }
  84
  85   if (fhistSgn){
  86     delete fhistSgn;
  87     fhistSgn=0;
  88   }
  89
  90   if (fhistBkg){
  91     delete fhistBkg;
  92     fhistBkg=0;
  93   }
  94
  95   if (fOutput) {
  96     delete fOutput;
  97     fOutput = 0;
  98   }
  99   if (fVHF) {
 100     delete fVHF;
 101     fVHF = 0;
 102   }
 103
 104 }
 105
 106 //________________________________________________________________________
 107 void AliAnalysisTaskSED0Mass::Init()
 108 {
 109   // Initialization
 110
 111   if(fDebug > 1) printf("AnalysisTaskSED0Mass::Init() \n");
 112
 113   gROOT->LoadMacro("$ALICE_ROOT/PWG3/vertexingHF/ConfigVertexingHF.C");
 114
 115   fVHF = (AliAnalysisVertexingHF*)gROOT->ProcessLine("ConfigVertexingHF()");
 116   // set dedidcated cuts
 117   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
 118   fVHF->PrintStatus();
 119
 120   return;
 121 }
 122
 123 //________________________________________________________________________
 124 void AliAnalysisTaskSED0Mass::UserCreateOutputObjects()
 125 {
 126
 127   // Create the output container
 128   //
 129   if(fDebug > 1) printf("AnalysisTaskSED0Mass::UserCreateOutputObjects() \n");
 130
 131   // Several histograms are more conveniently managed in a TList
 132   fOutput = new TList();
 133   fOutput->SetOwner();
 134
 135   const Int_t nhist=3;
 136   fhistMass=new TClonesArray("TH1F",nhist);
 137   fhistMass->SetName("fhistMass");
 138   fhistSgn =new TClonesArray("TH1F",nhist);
 139   fhistSgn->SetName("fhistSgn");
 140   fhistBkg =new TClonesArray("TH1F",nhist);
 141   fhistBkg->SetName("fhistBkg");
 142
 143   TString nameMass, nameSgn, nameBkg;
 144
 145   for(Int_t i=0;i<nhist;i++){
 146     nameMass="fhistMass_";
 147     nameMass+=i+1;
 148     cout<<"Creating TH1F with name "<<nameMass<<endl;
 149     new ((*fhistMass)[i]) TH1F(nameMass.Data(),"D^{0} invariant mass; M [GeV]; Entries",200,1.765,1.965);
 150     nameSgn="fhistSgn_";
 151     nameSgn+=i+1;
 152     new ((*fhistSgn)[i]) TH1F(nameSgn.Data(), "D^{0} invariant mass - MC; M [GeV]; Entries",200,1.765,1.965);
 153     nameBkg="fhistBkg_";
 154     nameBkg+=i+1;
 155     new ((*fhistBkg)[i]) TH1F(nameBkg.Data(), "Background invariant mass - MC; M [GeV]; Entries",200,1.765,1.965);
 156     printf("Created histograms %s\t%s\t%s\n",fhistMass->UncheckedAt(i)->GetName(),fhistSgn->UncheckedAt(i)->GetName(),fhistBkg->UncheckedAt(i)->GetName());
 157   }
 158
 159   fOutput->Add(fhistMass);
 160   fOutput->Add(fhistSgn);
 161   fOutput->Add(fhistBkg);
 162
 163   return;
 164 }
 165
 166 //________________________________________________________________________
 167 void AliAnalysisTaskSED0Mass::UserExec(Option_t */*option*/)
 168 {
 169   // Execute analysis for current event:
 170   // heavy flavor candidates association to MC truth
 171
 172   AliAODEvent *aod = dynamic_cast<AliAODEvent*> (InputEvent());
 173
 174   // load D0->Kpi candidates
 175   TClonesArray *inputArrayD0toKpi =
 176     (TClonesArray*)aod->GetList()->FindObject("D0toKpi");
 177   if(!inputArrayD0toKpi) {
 178     printf("AliAnalysisTaskSECompareHFpt::UserExec: D0toKpi branch not found!\n");
 179     return;
 180   }
 181
 182   // AOD primary vertex
 183   AliAODVertex *vtx1 = (AliAODVertex*)aod->GetPrimaryVertex();
 184   //vtx1->Print();
 185
 186   // load MC particles
 187   TClonesArray *mcArray =
 188     (TClonesArray*)aod->GetList()->FindObject(AliAODMCParticle::StdBranchName());
 189   if(!mcArray) {
 190     printf("AliAnalysisTaskSED0Mass::UserExec: MC particles branch not found!\n");
 191     return;
 192   }
 193
 194   // load MC header
 195   AliAODMCHeader *mcHeader =
 196     (AliAODMCHeader*)aod->GetList()->FindObject(AliAODMCHeader::StdBranchName());
 197   if(!mcHeader) {
 198     printf("AliAnalysisTaskSED0Mass::UserExec: MC header branch not found!\n");
 199     return;
 200   }
 201
 202
 203   // loop over D0->Kpi candidates
 204   Int_t nInD0toKpi = inputArrayD0toKpi->GetEntriesFast();
 205   printf("Number of D0->Kpi: %d\n",nInD0toKpi);
 206
 207   for (Int_t iD0toKpi = 0; iD0toKpi < nInD0toKpi; iD0toKpi++) {
 208
 209     AliAODRecoDecayHF2Prong *d = (AliAODRecoDecayHF2Prong*)inputArrayD0toKpi->UncheckedAt(iD0toKpi);
 210     Bool_t unsetvtx=kFALSE;
 211     if(!d->GetOwnPrimaryVtx()) {
 212       d->SetOwnPrimaryVtx(vtx1); // needed to compute all variables
 213       unsetvtx=kTRUE;
 214     }
 215
 216     //  cout<<"Cuts applied: ";
 217     //     Double_t *cutsapplied=(Double_t*)fVHF->GetD0toKpiCuts();
 218     //     for (Int_t i=0;i<9;i++){
 219     //       printf(cutsapplied[i]\t);
 220     //     }
 221     //    cout<<endl;
 222
 223
 224     Double_t pt = d->Pt();
 225     Int_t ptbin=0;
 226
 227     //cout<<"P_t = "<<pt<<endl;
 228     if (pt>0. && pt<=1.) {
 229       ptbin=1;
 230       //printf("I'm in the bin %d\n",ptbin);
 231     }
 232     else {
 233       if(pt>1. && pt<=3.) {
 234         ptbin=2;
 235         //printf("I'm in the bin %d\n",ptbin);
 236       }
 237       else {
 238         ptbin=3;
 239         //printf("I'm in the bin %d\n",ptbin);
 240       }//if(pt>3)
 241     }
 242
 243     FillHists(ptbin,d,mcArray);
 244
 245     if(unsetvtx) d->UnsetOwnPrimaryVtx();
 246
 247   }
 248
 249
 250
 251   // Post the data
 252   PostData(1,fOutput);
 253
 254
 255   return;
 256 }
 257 //____________________________________________________________________________*
 258 void AliAnalysisTaskSED0Mass::FillHists(Int_t ptbin, AliAODRecoDecayHF2Prong *part, TClonesArray *arrMC){
 259   //
 260   // function used in UserExec:
 261   //
 262   Int_t okD0=0,okD0bar=0;
 263
 264   if(part->SelectD0(fVHF->GetD0toKpiCuts(),okD0,okD0bar)) {//selected
 265     Double_t invmassD0 = part->InvMassD0(), invmassD0bar = part->InvMassD0bar();
 266     //printf("SELECTED\n");
 267     Int_t labD0 = part->MatchToMC(421,arrMC); //return MC particle label if the array corresponds to a D0, -1 if not (cf. AliAODRecoDecay.cxx)
 268     //printf("labD0 %d",labD0);
 269
 270     if(labD0>=0) {
 271
 272       AliAODMCParticle *partD0 = (AliAODMCParticle*)arrMC->At(labD0);
 273
 274       Int_t pdgD0 = partD0->GetPdgCode();
 275
 276       //printf(" pdgD0 %d\n",pdgD0);
 277
 278       if (pdgD0==421){ //D0
 279         //cout<<"Fill S with D0"<<endl;
 280         ((TH1F*)(fhistSgn->UncheckedAt(ptbin-1)))->Fill(invmassD0);
 281       }
 282       else {//D0bar
 283         //printf("Fill S with D0bar");
 284         ((TH1F*)(fhistSgn->UncheckedAt(ptbin-1)))->Fill(invmassD0bar);
 285       }
 286
 287     }
 288     else {//background
 289           //printf("Fill background");
 290       ((TH1F*)(fhistBkg->UncheckedAt(ptbin-1)))->Fill(invmassD0);
 291       ((TH1F*)(fhistBkg->UncheckedAt(ptbin-1)))->Fill(invmassD0bar);
 292     }
 293
 294     //no MC info, just cut selection
 295     if (okD0==1) {
 296       //printf("Fill mass with D0");
 297       ((TH1F*)(fhistMass->UncheckedAt(ptbin-1)))->Fill(invmassD0);
 298
 299     }
 300     if (okD0bar==1) {
 301       ((TH1F*)fhistMass->UncheckedAt(ptbin-1))->Fill(invmassD0bar);
 302       //printf("Fill mass with D0bar");
 303     }
 304
 305   } //else cout<<"NOT SELECTED"<<endl;
 306
 307 }
 308 //________________________________________________________________________
 309 void AliAnalysisTaskSED0Mass::Terminate(Option_t */*option*/)
 310 {
 311   // Terminate analysis
 312   //
 313   if(fDebug > 1) printf("AnalysisTaskSED0Mass: Terminate() \n");
 314
 315   fOutput = dynamic_cast<TList*> (GetOutputData(1));
 316   if (!fOutput) {
 317     printf("ERROR: fOutput not available\n");
 318     return;
 319   }
 320
 321
 322   fhistMass = dynamic_cast<TClonesArray*>(fOutput->FindObject("fhistMass"));
 323   fhistSgn = dynamic_cast<TClonesArray*>(fOutput->FindObject("fhistSgn"));
 324   fhistBkg = dynamic_cast<TClonesArray*>(fOutput->FindObject("fhistBkg"));
 325
 326   return;
 327 }
 328