PWGHF/vertexingHF/charmFlow/AliAnalysisTaskFlowD2H.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 // FlowD2H main task:
  18 // >> Make flowEvent with RPcuts and POIcuts given in constructor and passes it
  19 //    to the daughter tasks.
  20 // >> The POIcuts are polymorphic based on the AliRDHFCuts class allowing the
  21 //    use of all charmed candidates reconstructed in the central barrel.
  22 // Author: Carlos Perez (cperez@cern.ch)
  23 //==============================================================================
  24
  25 /* $Id$ */
  26
  27 #include "TChain.h"
  28 #include "TList.h"
  29 #include "TH2D.h"
  30 #include "TProfile.h"
  31
  32 #include "AliAnalysisTaskSE.h"
  33 #include "AliAnalysisManager.h"
  34
  35 #include "AliCentrality.h"
  36
  37 #include "AliAODEvent.h"
  38 #include "AliAODInputHandler.h"
  39 #include "AliAODTrack.h"
  40
  41 #include "TMath.h"
  42 #include "TObjArray.h"
  43
  44 #include "AliFlowCandidateTrack.h"
  45 #include "AliFlowTrackCuts.h"
  46 #include "AliFlowEvent.h"
  47 #include "AliFlowCommonConstants.h"
  48
  49 #include "AliRDHFCuts.h"
  50 #include "AliRDHFCutsD0toKpi.h"
  51 #include "AliRDHFCutsDStartoKpipi.h"
  52 #include "AliRDHFCutsDplustoKpipi.h"
  53
  54 #include "AliAODRecoDecayHF2Prong.h"
  55 #include "AliAODRecoDecayHF3Prong.h"
  56 #include "AliAODRecoCascadeHF.h"
  57
  58 #include "AliAnalysisTaskFlowD2H.h"
  59
  60 ClassImp(AliAnalysisTaskFlowD2H)
  61
  62 //_____________________________________________________________________________
  63 AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H() :
  64   AliAnalysisTaskSE(), fCutsRP(NULL), fCutsPOI(NULL), fSource(0),
  65   fDebugV2(kFALSE), fHList(NULL), fAnaCuts(NULL)
  66 {
  67 // Default constructor
  68   for(int i=0; i!=2; ++i)
  69     fEvent[i] = NULL;
  70   for(int i=0; i!=2; ++i)
  71     fMass[i] = NULL;
  72   for(int i=0; i!=2; ++i)
  73     fPOIEta[i] = 0;
  74   for(int i=0; i!=4; ++i)
  75     fFlowEta[i] = 0;
  76   for(int x=0; x!=2; ++x)
  77     fFlowPts[x] = 0;
  78   for(int x=0; x!=2; ++x)
  79     for(int m=0; m!=5; ++m)
  80       fFlowBands[x][m] = 0;
  81 }
  82
  83 //_____________________________________________________________________________
  84 AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H(const char *name,
  85     AliFlowTrackCuts *cutsRPs, AliRDHFCuts *cutsPOIs, Int_t specie) :
  86   AliAnalysisTaskSE(name), fCutsRP(cutsRPs),
  87   fCutsPOI(cutsPOIs), fSource(specie),
  88   fDebugV2(kFALSE), fHList(NULL), fAnaCuts(NULL)
  89 {
  90 // Standard constructor
  91   for(int i=0; i!=2; ++i)
  92     fEvent[i] = NULL;
  93   for(int i=0; i!=2; ++i)
  94     fMass[i] = NULL;
  95   for(int i=0; i!=2; ++i)
  96     fPOIEta[i] = 0;
  97   for(int i=0; i!=4; ++i)
  98     fFlowEta[i] = 0;
  99   for(int x=0; x!=2; ++x)
 100     fFlowPts[x] = 0;
 101   for(int x=0; x!=2; ++x)
 102     for(int m=0; m!=5; ++m)
 103       fFlowBands[x][m] = 0;
 104
 105   DefineInput( 0,TChain::Class());
 106   DefineOutput(1,TList::Class());
 107   DefineOutput(2,AliFlowEventSimple::Class()); // first band
 108   DefineOutput(3,AliFlowEventSimple::Class()); // second band
 109   DefineOutput(4,AliFlowEventSimple::Class()); // third band
 110   DefineOutput(5,AliFlowEventSimple::Class()); // fourth band
 111   DefineOutput(6,AliFlowEventSimple::Class()); // fifth band
 112 }
 113
 114 //_____________________________________________________________________________
 115 AliAnalysisTaskFlowD2H::~AliAnalysisTaskFlowD2H(){
 116   // delete objects
 117   if(fHList)
 118     delete fHList;
 119 }
 120
 121 //_____________________________________________________________________________
 122 void AliAnalysisTaskFlowD2H::UserCreateOutputObjects(){
 123   // Define output objects + parameters
 124   if(fDebugV2)
 125     printf("DEBUGGER: Creating output\n");
 126   fHList = new TList();
 127   fHList->SetOwner();
 128   AddHistograms();
 129   PostData(1,fHList);
 130
 131   AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();
 132   cc->SetNbinsMult(10000);
 133   cc->SetMultMin(0);
 134   cc->SetMultMax(10000);
 135
 136   cc->SetNbinsPt(fFlowPts[1]-fFlowPts[0]);
 137   cc->SetPtMin(fFlowPts[0]);
 138   cc->SetPtMax(fFlowPts[1]);
 139
 140   cc->SetNbinsPhi(180);
 141   cc->SetPhiMin(0.0);
 142   cc->SetPhiMax(TMath::TwoPi());
 143
 144   cc->SetNbinsEta(200);
 145   cc->SetEtaMin(-5.0);
 146   cc->SetEtaMax(+5.0);
 147
 148   cc->SetNbinsQ(500);
 149   cc->SetQMin(0.0);
 150   cc->SetQMax(3.0);
 151 }
 152 //_____________________________________________________________________________
 153 void AliAnalysisTaskFlowD2H::AddHistograms()
 154 {
 155   // Create histograms
 156   TList *tEvents = new TList();
 157     tEvents->SetName("Events");
 158     tEvents->SetOwner();
 159     for(int i=0; i!=2; ++i) {
 160       fEvent[i] = new TH2D(Form("Event%d",i),"Events;V0M;Arb",10,0,100,10,0,12000);
 161       tEvents->Add(fEvent[i]);
 162     }
 163     fAnaCuts=new TProfile("Cuts","Analysis Cuts", 10,0,10);
 164       fAnaCuts->Fill(0.5,fPOIEta[0],1);fAnaCuts->GetXaxis()->SetBinLabel(1,"ETAm");
 165       fAnaCuts->Fill(1.5,fPOIEta[1],1);fAnaCuts->GetXaxis()->SetBinLabel(2,"ETAM");
 166     tEvents->Add(fAnaCuts);
 167     fHList->Add(tEvents);
 168   TList *tCandidates;
 169     tCandidates = new TList();
 170     tCandidates->SetOwner();
 171     tCandidates->SetName(Form("Candidates%d",fSource));
 172     Double_t dBinMin=0, dBinMax=3;
 173     Int_t nBins=3;
 174     switch (fSource) {
 175       case (AliRDHFCuts::kD0toKpiCuts): // 360/72 (bw=5MeV)
 176         dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
 177       case (AliRDHFCuts::kDstarCuts): // 36/72 (bw=0.5MeV)
 178         dBinMin = 0.137; dBinMax = 0.173; nBins=72; break;
 179       case (AliRDHFCuts::kDplusCuts): // 360/72 (bw=5MeV)
 180         dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
 181     }
 182     for(int i=0; i!=2; ++i) {
 183       fMass[i] = new TH2D( Form("Mass%d",i),
 184                            Form("Mass%d;Mass [GeV];Pt [GeV]",i),
 185                            nBins, dBinMin, dBinMax,
 186                            fFlowPts[1]-fFlowPts[0], fFlowPts[0], fFlowPts[1] );
 187       tCandidates->Add(fMass[i]);
 188     }
 189     fHList->Add(tCandidates);
 190   if (fDebugV2) printf("DEBUGGER: Created histos for DMeson %d\n", fSource );
 191 }
 192 //_____________________________________________________________________________
 193 void AliAnalysisTaskFlowD2H::NotifyRun()
 194 {
 195   // Do nothing
 196 }
 197 //_____________________________________________________________________________
 198 void AliAnalysisTaskFlowD2H::UserExec(Option_t *)
 199 {
 200   // Do analysis + fIll histograms
 201   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 202
 203   if(!fAOD)
 204     return;
 205   Int_t iMulti = fAOD->GetNTracks();  /// TO DO
 206   fEvent[0]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );
 207   if(!fCutsPOI->IsEventSelected(fAOD)) return;
 208   if(fCutsPOI->IsEventSelectedInCentrality(fAOD)>0) return;
 209   fEvent[1]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );
 210
 211   if (fDebugV2) printf("Event selected\n");
 212   fCutsRP->SetEvent(fAOD,MCEvent());
 213   AliFlowTrackCuts* dummy = new AliFlowTrackCuts("null_cuts");
 214   dummy->SetParamType(AliFlowTrackCuts::kGlobal);
 215   dummy->SetPtRange(+1,-1); // select nothing QUICK
 216   dummy->SetEtaRange(+1,-1); // select nothing VZERO
 217   dummy->SetEvent(fAOD,MCEvent());
 218   AliFlowEvent *flowEvent[5];
 219   for(int r=0; r!=5; ++r) {
 220     flowEvent[r] = new AliFlowEvent(fCutsRP,dummy);
 221     flowEvent[r]->SetReferenceMultiplicity( iMulti );
 222     flowEvent[r]->DefineDeadZone(0,0,0,0);
 223     flowEvent[r]->TagSubeventsInEta( fFlowEta[0], fFlowEta[1],
 224                                      fFlowEta[2], fFlowEta[3] );
 225   }
 226   delete dummy;
 227   if (fDebugV2) printf(" ᶫFlowEvent has %d RPs\n", flowEvent[0]->NumberOfTracks() );
 228
 229   switch (fSource) {
 230     case (AliRDHFCuts::kD0toKpiCuts):
 231       FillD0toKpi(fAOD,flowEvent); break;
 232     case (AliRDHFCuts::kDstarCuts):
 233       FillDStartoKpipi(fAOD,flowEvent); break;
 234     case (AliRDHFCuts::kDplusCuts):
 235       FillDplustoKpipi(fAOD,flowEvent); break;
 236   }
 237
 238   for(int m=0; m!=5; ++m)
 239     PostData(2+m,flowEvent[m]);
 240   PostData(1,fHList);
 241
 242 }
 243 //______________________________________________________________________________
 244 void AliAnalysisTaskFlowD2H::FillD0toKpi(AliAODEvent *theAOD,
 245                                          AliFlowEvent *theMB[5] )
 246 {
 247   // Fill D0->Kpi histos
 248   TList *listHF = (TList*) theAOD->GetList();
 249   if(!listHF) return;
 250   TClonesArray *listDzero = (TClonesArray*) listHF->FindObject("D0toKpi");
 251   if(!listDzero) return;
 252   int nEntries = listDzero->GetEntriesFast();
 253   if( !nEntries ) return;
 254   AliRDHFCutsD0toKpi *fCutsD0toKpi = (AliRDHFCutsD0toKpi*) fCutsPOI;
 255   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 256   const Int_t ndght=2;
 257   Int_t nIDs[ndght];
 258   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 259     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*) listDzero->UncheckedAt( iEntry );
 260     if( !d0cand ) continue;
 261     if( !d0cand->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts) ) continue;
 262     if( !fCutsD0toKpi->IsInFiducialAcceptance(d0cand->Pt(),d0cand->Y(421)) )continue;
 263     int topCut  = fCutsD0toKpi->IsSelected( d0cand, AliRDHFCuts::kAll, theAOD );
 264     int nLevel=topCut>0?1:0;
 265     if( (d0cand->Eta()<fPOIEta[0])||(d0cand->Eta()>fPOIEta[1]) )
 266       nLevel=0;
 267     Double_t massD0=topCut>1?d0cand->InvMassD0bar():d0cand->InvMassD0();
 268     for(int h=0; h!=nLevel+1; ++h)
 269       fMass[h]->Fill(massD0,d0cand->Pt());
 270     if( (d0cand->Pt()<fFlowPts[0]) || (d0cand->Pt()>fFlowPts[1]) ) continue;
 271     AliAODTrack* iT;
 272     for(Int_t i=0; i!=ndght; ++i) {
 273       iT = (AliAODTrack*)d0cand->GetDaughter(i);
 274       nIDs[i] = iT->GetID();
 275     }
 276     // Candidates Insertion (done in filling method: faster)
 277     if(nLevel)
 278       for(Int_t r=0; r!=5; ++r)
 279         if( (massD0>=fFlowBands[0][r]) && (massD0<fFlowBands[1][r]) ) {
 280           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 281                                  MakeTrack(massD0, d0cand->Pt(), d0cand->Phi(),
 282                                            d0cand->Eta(), ndght, nIDs);
 283           if(fDebugV2) printf("   ᶫInjecting D0 candidate on band %d \n", r);
 284           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 285             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 286               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 287               if(!iRP->InRPSelection()) continue;
 288               if( TMath::Abs(sTrack->GetIDDaughter(iDau)) == TMath::Abs(iRP->GetID()) ) {
 289                 sTrack->SetDaughter(iDau,iRP);
 290                 iRP->SetForRPSelection(kFALSE);
 291                 if(fDebugV2) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 292               }
 293             }
 294           theMB[r]->AddTrack(sTrack);
 295         }
 296     // END of Candidates Insertion
 297   }
 298 }
 299 //______________________________________________________________________________
 300 void AliAnalysisTaskFlowD2H::FillDStartoKpipi(const AliAODEvent *theAOD,
 301                                               AliFlowEvent *theMB[5] )
 302 {
 303   // Fills D* histos
 304   TList *listHF = (TList*) theAOD->GetList();
 305   if(!listHF) return;
 306   TClonesArray *listDstar = (TClonesArray*) listHF->FindObject("Dstar");
 307   if(!listDstar) return;
 308   int nEntries = listDstar->GetEntriesFast();
 309   if( !nEntries ) return;
 310   AliRDHFCutsDStartoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDStartoKpipi*) fCutsPOI;
 311   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 312   const Int_t ndght=3;
 313   Int_t nIDs[ndght];
 314   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 315     AliAODRecoCascadeHF *dst = (AliAODRecoCascadeHF*) listDstar->UncheckedAt( iEntry );
 316     if( !dst ) continue;
 317     if( !dst->HasSelectionBit(AliRDHFCuts::kDstarCuts) ) continue;
 318     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dst->Pt(),dst->YDstar()) )continue;
 319     int topCut = fCutsDStartoKpipi->IsSelected( dst, AliRDHFCuts::kAll );
 320     int nLevel=topCut>0?1:0;
 321     if( (dst->Eta()<fPOIEta[0])||(dst->Eta()>fPOIEta[1]) )
 322       nLevel=0;
 323     Double_t massDS=dst->DeltaInvMass();
 324     for(int h=0; h!=nLevel+1; ++h)
 325       fMass[h]->Fill(massDS,dst->Pt());
 326     if( (dst->Pt()<fFlowPts[0]) || (dst->Pt()>fFlowPts[1]) ) continue;
 327     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*)dst->Get2Prong();
 328     nIDs[0] = ((AliAODTrack*)d0cand->GetDaughter(0))->GetID();
 329     nIDs[1] = ((AliAODTrack*)d0cand->GetDaughter(1))->GetID();
 330     nIDs[2] = ((AliAODTrack*)dst->GetBachelor() )->GetID();
 331     // Candidates Insertion (done in filling method: faster)
 332     if(nLevel)
 333       for(Int_t r=0; r!=5; ++r)
 334         if( (massDS>=fFlowBands[0][r]) && (massDS<fFlowBands[1][r]) ) {
 335           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 336                                  MakeTrack(massDS, dst->Pt(), dst->Phi(),
 337                                            dst->Eta(), ndght, nIDs);
 338           if(fDebugV2) printf("   ᶫInjecting DStar candidate on band %d \n", r);
 339           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 340             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 341               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 342               if(!iRP->InRPSelection()) continue;
 343               if( TMath::Abs(sTrack->GetIDDaughter(iDau)) == TMath::Abs(iRP->GetID()) ) {
 344                 sTrack->SetDaughter(iDau,iRP);
 345                 iRP->SetForRPSelection(kFALSE);
 346                 if(fDebugV2) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 347               }
 348             }
 349           theMB[r]->AddTrack(sTrack);
 350         }
 351     // END of Candidates Insertion
 352   }
 353 }
 354 //______________________________________________________________________________
 355 void AliAnalysisTaskFlowD2H::FillDplustoKpipi(const AliAODEvent *theAOD,
 356                                               AliFlowEvent *theMB[5] )
 357 {
 358   // Fill D+ histos
 359   TList *listHF = (TList*) theAOD->GetList();
 360   if(!listHF) return;
 361   TClonesArray *listDplus = (TClonesArray*) listHF->FindObject("Charm3Prong");
 362   if(!listDplus) return;
 363   int nEntries = listDplus->GetEntriesFast();
 364   if( !nEntries ) return;
 365   AliRDHFCutsDplustoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDplustoKpipi*) fCutsPOI;
 366   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 367   const Int_t ndght=3;
 368   Int_t nIDs[ndght];
 369   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 370     AliAODRecoDecayHF3Prong *dplu = (AliAODRecoDecayHF3Prong*) listDplus->UncheckedAt( iEntry );
 371     if( !dplu ) continue;
 372     if( !dplu->HasSelectionBit(AliRDHFCuts::kDplusCuts) ) continue;
 373     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dplu->Pt(),dplu->YDplus()) )continue;
 374     int topCut = fCutsDStartoKpipi->IsSelected( dplu, AliRDHFCuts::kAll );
 375     int nLevel=topCut>0?1:0;
 376     if( (dplu->Eta()<fPOIEta[0])||(dplu->Eta()>fPOIEta[1]) )
 377       nLevel=0;
 378     Double_t massDp=dplu->InvMassDplus();
 379     for(int h=0; h!=nLevel+1; ++h)
 380       fMass[h]->Fill(massDp,dplu->Pt());
 381     if( (dplu->Pt()<fFlowPts[0]) || (dplu->Pt()>fFlowPts[1]) ) continue;
 382     nIDs[0] = ((AliAODTrack*)dplu->GetDaughter(0))->GetID();
 383     nIDs[1] = ((AliAODTrack*)dplu->GetDaughter(1))->GetID();
 384     nIDs[2] = ((AliAODTrack*)dplu->GetDaughter(2))->GetID();
 385     // Candidates Insertion (done in filling method: faster)
 386     if(nLevel)
 387       for(Int_t r=0; r!=5; ++r)
 388         if( (massDp>=fFlowBands[0][r]) && (massDp<fFlowBands[1][r]) ) {
 389           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 390                                  MakeTrack(massDp, dplu->Pt(), dplu->Phi(),
 391                                            dplu->Eta(), ndght, nIDs);
 392           if(fDebugV2) printf("   ᶫInjecting DStar candidate on band %d \n", r);
 393           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 394             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 395               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 396               if(!iRP->InRPSelection()) continue;
 397               if( TMath::Abs(sTrack->GetIDDaughter(iDau)) == TMath::Abs(iRP->GetID()) ) {
 398                 sTrack->SetDaughter(iDau,iRP);
 399                 iRP->SetForRPSelection(kFALSE);
 400                 if(fDebugV2) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 401               }
 402             }
 403           theMB[r]->AddTrack(sTrack);
 404         }
 405     // END of Candidates Insertion
 406   }
 407 }
 408 //______________________________________________________________________________
 409 AliFlowCandidateTrack* AliAnalysisTaskFlowD2H::MakeTrack( Double_t mass,
 410                           Double_t pt, Double_t phi, Double_t eta,
 411                           Int_t nDau, const Int_t *iID ) {
 412   // create track for flow tasks
 413   AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
 414   sTrack->SetMass(mass);
 415   sTrack->SetPt(pt);
 416   sTrack->SetPhi(phi);
 417   sTrack->SetEta(eta);
 418   for(Int_t iDau=0; iDau!=nDau; ++iDau)
 419     sTrack->AddDaughter(iID[iDau]);
 420   sTrack->SetForPOISelection(kTRUE);
 421   sTrack->SetForRPSelection(kFALSE);
 422   return sTrack;
 423 }
 424 //_____________________________________________________________________________
 425 void AliAnalysisTaskFlowD2H::Terminate(Option_t *)
 426 {
 427   // Close the analysis, empty for now
 428 }
 429