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 // >> Select candidates and passes flowevents to the daughter tasks.
  19 // >> The POIcuts are polymorphic based on the AliRDHFCuts class allowing the
  20 //    use of all charmed candidates reconstructed in the central barrel.
  21 // Author: Carlos Perez (cperez@cern.ch)
  22 //==============================================================================
  23
  24 /* $Id$ */
  25
  26 #include "TChain.h"
  27
  28 #include "AliAnalysisTaskSE.h"
  29 #include "AliAnalysisManager.h"
  30
  31 #include "AliCentrality.h"
  32
  33 #include "AliAODEvent.h"
  34 #include "AliAODInputHandler.h"
  35 #include "AliAODTrack.h"
  36
  37 #include "AliRDHFCuts.h"
  38
  39 #include "AliRDHFCutsD0toKpi.h"
  40 #include "AliRDHFCutsDStartoKpipi.h"
  41 #include "AliRDHFCutsDplustoKpipi.h"
  42
  43 #include "AliAODRecoDecayHF2Prong.h"
  44 #include "AliAODRecoDecayHF3Prong.h"
  45 #include "AliAODRecoCascadeHF.h"
  46
  47 #include "TMath.h"
  48
  49 #include "AliFlowCommonConstants.h"
  50 #include "AliFlowEvent.h"
  51 #include "AliFlowTrackCuts.h"
  52 #include "AliFlowCandidateTrack.h"
  53
  54 #include "TObjArray.h"
  55 #include "TList.h"
  56 #include "TH1D.h"
  57
  58 #include "AliAnalysisTaskFlowD2H.h"
  59
  60 ClassImp(AliAnalysisTaskFlowD2H)
  61
  62 //_____________________________________________________________________________
  63 AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H() :
  64   AliAnalysisTaskSE(), fTPCEvent(NULL), fVZEEvent(NULL),
  65   fCutsTPC(NULL), fCutsVZE(NULL), fNoPOIs(NULL), fCutsPOI(NULL),
  66   fSource(0), fDebugV2(kFALSE), fMassBins(0), fMinMass(0.),
  67   fMaxMass(0.), fPtBinWidth(0), fHList(NULL), fEvent(NULL),
  68   fCC(NULL), fRFPMTPC(NULL), fRFPPhiTPC(NULL), fCandidates(NULL)
  69 {
  70 // Default constructor
  71 }
  72
  73 //_____________________________________________________________________________
  74 AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H(const char *name,
  75                                                AliFlowTrackCuts *cutsTPC,
  76                                                AliFlowTrackCuts *cutsVZE,
  77                                                AliRDHFCuts *cutsPOIs,
  78                                                Int_t specie) :
  79   AliAnalysisTaskSE(name), fTPCEvent(NULL), fVZEEvent(NULL),
  80   fCutsTPC(cutsTPC), fCutsVZE(cutsVZE), fNoPOIs(NULL), fCutsPOI(cutsPOIs),
  81   fSource(specie), fDebugV2(kFALSE), fMassBins(0), fMinMass(0.),
  82   fMaxMass(0.), fPtBinWidth(0), fHList(NULL), fEvent(NULL),
  83   fCC(NULL), fRFPMTPC(NULL), fRFPPhiTPC(NULL), fCandidates(NULL)
  84 {
  85 // Standard constructor
  86   DefineInput( 0,TChain::Class());
  87   DefineOutput(1,TList::Class());
  88   DefineOutput(2,AliFlowEventSimple::Class());
  89   DefineOutput(3,AliFlowEventSimple::Class());
  90 }
  91
  92 //_____________________________________________________________________________
  93 AliAnalysisTaskFlowD2H::~AliAnalysisTaskFlowD2H(){
  94   // delete objects
  95   if(fTPCEvent) delete fTPCEvent;
  96   if(fVZEEvent) delete fVZEEvent;
  97   if(fCutsTPC) delete fCutsTPC;
  98   if(fCutsVZE) delete fCutsVZE;
  99   if(fNoPOIs) delete fNoPOIs;
 100   if(fCutsPOI) delete fCutsPOI;
 101   if(fHList) delete fHList;
 102   if(fCandidates) delete fCandidates;
 103 }
 104
 105 //_____________________________________________________________________________
 106 void AliAnalysisTaskFlowD2H::UserCreateOutputObjects(){
 107   // Define output objects + parameters
 108   if(fDebugV2)
 109     printf("DEBUGGER: Creating output\n");
 110   fHList = new TList();
 111   fHList->SetOwner();
 112   fEvent = new TH1D("Events","Events",3,0,3);
 113   fEvent->GetXaxis()->SetBinLabel(1,"REACHED");
 114   fEvent->GetXaxis()->SetBinLabel(2,"SELECTED");
 115   fEvent->GetXaxis()->SetBinLabel(3,"DELTA AOD REACHED");
 116   fHList->Add( fEvent );
 117   fCC = new TH1D("CentralityClass","Centrality Class",50,0,100);
 118   fHList->Add( fCC );
 119   fRFPMTPC = new TH1D("RFPMultiplicityTPC","RFP Multiplicity TPC",300,0,3000);
 120   fHList->Add( fRFPMTPC );
 121   fRFPPhiTPC = new TH1D("RFPPhiTPC","RFP Phi TPC",180,0,TMath::TwoPi());
 122   fHList->Add( fRFPPhiTPC );
 123
 124   fCandidates = new TObjArray(100);
 125   fCandidates->SetOwner();
 126
 127   AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();
 128   cc->SetNbinsMult(1);
 129   cc->SetNbinsPt(24/fPtBinWidth);
 130   cc->SetNbinsPhi(1);
 131   cc->SetNbinsEta(15);
 132   cc->SetNbinsQ(1);
 133   cc->SetNbinsMass( fMassBins );
 134   cc->SetMultMin(1);
 135   cc->SetMultMax(2);
 136   cc->SetPtMin(0);
 137   cc->SetPtMax(24);
 138   cc->SetPhiMin(0);
 139   cc->SetPhiMax(TMath::TwoPi());
 140   cc->SetEtaMin(-3.9);
 141   cc->SetEtaMax(+5.1);
 142   cc->SetQMin(0);
 143   cc->SetQMax(1);
 144   cc->SetMassMin( fMinMass );
 145   cc->SetMassMax( fMaxMass );
 146
 147   fTPCEvent = new AliFlowEvent(3000);
 148   fVZEEvent = new AliFlowEvent(170);
 149
 150   fNoPOIs = new AliFlowTrackCuts( "noPOIs" );
 151   fNoPOIs->SetParamType(AliFlowTrackCuts::kGlobal);
 152   fNoPOIs->SetPtRange(+1,-1);
 153
 154   PostData(1,fHList);
 155   PostData(2,fTPCEvent);
 156   PostData(3,fVZEEvent);
 157 }
 158 //_____________________________________________________________________________
 159 void AliAnalysisTaskFlowD2H::UserExec(Option_t *)
 160 {
 161   // Do analysis + fIll histograms
 162   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 163
 164   if(!fAOD) return;
 165   fEvent->Fill( 0 );
 166
 167   // floweventcuts::isselected() and alirdhfcuts::iseventselected() cut on the same
 168   // values in the same way BUT the latter also loads the PIDresponse object from the
 169   // event header!!!
 170   //  if(!fEventCuts->IsSelected(fAOD)) return;
 171   if(!fCutsPOI->IsEventSelected(fAOD)) return;
 172   fEvent->Fill( 1 );
 173
 174   fCC->Fill( fCutsPOI->GetCentrality(fAOD) );
 175
 176   fCutsTPC->SetEvent( fAOD, MCEvent() );
 177   fCutsVZE->SetEvent( fAOD, MCEvent() );
 178   fNoPOIs->SetEvent( fAOD, MCEvent() );
 179   fTPCEvent->Fill( fCutsTPC, fNoPOIs );
 180   fVZEEvent->Fill( fCutsVZE, fNoPOIs );
 181
 182   Int_t rfps = fTPCEvent->GetNumberOfRPs();
 183   fRFPMTPC->Fill( rfps );
 184   for(int iRPs=0; iRPs!=rfps; ++iRPs ) {
 185     AliFlowTrack *iRP = dynamic_cast<AliFlowTrack*>(fTPCEvent->GetTrack( iRPs ));
 186     if (!iRP) continue;
 187     fRFPPhiTPC->Fill( iRP->Phi() );
 188   }
 189
 190   if (fDebugV2) printf("Event selected\n");
 191   fCandidates->SetLast(-1); // resets the array
 192
 193   switch (fSource) {
 194     case (AliRDHFCuts::kD0toKpiCuts):
 195       FillD0toKpi(fAOD); break;
 196     case (AliRDHFCuts::kDstarCuts):
 197       FillDStartoKpipi(fAOD); break;
 198     case (AliRDHFCuts::kDplusCuts):
 199       FillDplustoKpipi(fAOD); break;
 200   }
 201
 202   if(fDebugV2) printf("TPCevent %d | VZEevent %d\n", fTPCEvent->NumberOfTracks(), fVZEEvent->NumberOfTracks() );
 203   //inject candidates
 204   if (fDebugV2)  printf("I received %d candidates\n",fCandidates->GetEntriesFast());
 205   for(int iCand=0; iCand!=fCandidates->GetEntriesFast(); ++iCand ) {
 206     AliFlowCandidateTrack *cand = dynamic_cast<AliFlowCandidateTrack*>(fCandidates->At(iCand));
 207     if (!cand) continue;
 208     if (fDebugV2) printf(" >Checking at candidate %d with %d daughters: mass %f\n",iCand,cand->GetNDaughters(),cand->Mass());
 209     for(int iDau=0; iDau!=cand->GetNDaughters(); ++iDau) {
 210       if(fDebugV2) printf("  >Daughter %d with fID %d", iDau, cand->GetIDDaughter(iDau));
 211       for(int iRPs=0; iRPs!=fTPCEvent->NumberOfTracks(); ++iRPs ) {
 212         AliFlowTrack *iRP = dynamic_cast<AliFlowTrack*>(fTPCEvent->GetTrack( iRPs ));
 213         if (!iRP) continue;
 214         if( !iRP->InRPSelection() ) continue;
 215         if( cand->GetIDDaughter(iDau) == iRP->GetID() ) {
 216           if(fDebugV2) printf(" was in RP set");
 217           iRP->SetForRPSelection(kFALSE);
 218           fTPCEvent->SetNumberOfRPs( fTPCEvent->GetNumberOfRPs() -1 );
 219         }
 220       }
 221       if(fDebugV2) printf("\n");
 222     }
 223     cand->SetForPOISelection(kTRUE);
 224     fTPCEvent->InsertTrack( ((AliFlowTrack*) cand) );
 225     fVZEEvent->InsertTrack( ((AliFlowTrack*) cand) );
 226   }
 227   if(fDebugV2) printf("TPCevent %d | VZEevent %d\n", fTPCEvent->NumberOfTracks(), fVZEEvent->NumberOfTracks() );
 228
 229   PostData(1,fHList);
 230   PostData(2,fTPCEvent);
 231   PostData(3,fVZEEvent);
 232 }
 233 //______________________________________________________________________________
 234 void AliAnalysisTaskFlowD2H::FillD0toKpi(const AliAODEvent *theAOD)
 235 {
 236   // Fill D0->Kpi histos
 237   TList *listHF = (TList*) theAOD->GetList();
 238   if(!listHF) return;
 239   TClonesArray *listDzero = (TClonesArray*) listHF->FindObject("D0toKpi");
 240   if(!listDzero) return;
 241   int nEntries = listDzero->GetEntriesFast();
 242   if( !nEntries ) return;
 243   fEvent->Fill( 2 ); // EVERYTHING OKAY
 244
 245   AliRDHFCutsD0toKpi *fCutsD0toKpi = (AliRDHFCutsD0toKpi*) fCutsPOI;
 246   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 247   Int_t nIDs[2];
 248   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 249     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*) listDzero->UncheckedAt( iEntry );
 250     if( !d0cand ) continue;
 251     // APPLYING CUTS
 252     if( !d0cand->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts) ) continue;
 253     if( !fCutsD0toKpi->IsInFiducialAcceptance(d0cand->Pt(),d0cand->Y(421)) )continue;
 254     Int_t topCut = fCutsD0toKpi->IsSelected( d0cand, AliRDHFCuts::kAll, NULL );
 255     if(!topCut) continue;
 256     Double_t massD0=topCut>1?d0cand->InvMassD0bar():d0cand->InvMassD0();
 257     // TO HANDLE AUTOCORRELATIONS
 258     nIDs[0] = ( (AliAODTrack*)d0cand->GetDaughter(0) )->GetID();
 259     nIDs[1] = ( (AliAODTrack*)d0cand->GetDaughter(1) )->GetID();
 260     // ADDING TRACK
 261     MakeTrack(massD0, d0cand->Pt(), d0cand->Phi(), d0cand->Eta(), 2, nIDs);
 262     if(fDebugV2) printf("   ᶫInjecting D0 candidate\n");
 263   }
 264 }
 265 //______________________________________________________________________________
 266 void AliAnalysisTaskFlowD2H::FillDStartoKpipi(const AliAODEvent *theAOD )
 267 {
 268   // Fills D* to Kpipi
 269   TList *listHF = (TList*) theAOD->GetList();
 270   if(!listHF) return;
 271   TClonesArray *listDstar = (TClonesArray*) listHF->FindObject("Dstar");
 272   if(!listDstar) return;
 273   int nEntries = listDstar->GetEntriesFast();
 274   if( !nEntries ) return;
 275   fEvent->Fill( 2 ); // EVERYTHING OKAY
 276
 277   AliRDHFCutsDStartoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDStartoKpipi*) fCutsPOI;
 278   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 279   Int_t nIDs[3];
 280   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 281     AliAODRecoCascadeHF *dst = (AliAODRecoCascadeHF*) listDstar->UncheckedAt( iEntry );
 282     if( !dst ) continue;
 283     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*)dst->Get2Prong();
 284     if(!d0cand) continue;
 285     if( !d0cand->GetDaughter(0) ) continue;
 286     if( !d0cand->GetDaughter(1) ) continue;
 287     if( !dst->GetBachelor() ) continue;
 288     // APPLYING CUTS
 289     if( !dst->HasSelectionBit(AliRDHFCuts::kDstarCuts) ) continue;
 290     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dst->Pt(),dst->YDstar()) )continue;
 291     Int_t topCut=0;
 292     if(dst->Pt() > 5) {
 293       fCutsDStartoKpipi->SetUsePID(kFALSE);
 294       topCut = fCutsDStartoKpipi->IsSelected( dst, AliRDHFCuts::kAll );
 295       fCutsDStartoKpipi->SetUsePID(kTRUE);
 296     } else topCut = fCutsDStartoKpipi->IsSelected( dst, AliRDHFCuts::kAll );
 297     if(!topCut) continue;
 298     Double_t massDS=dst->DeltaInvMass();
 299     // TO HANDLE AUTOCORRELATIONS
 300     nIDs[0] = ((AliAODTrack*)d0cand->GetDaughter(0))->GetID();
 301     nIDs[1] = ((AliAODTrack*)d0cand->GetDaughter(1))->GetID();
 302     nIDs[2] = ((AliAODTrack*)dst->GetBachelor() )->GetID();
 303     // ADDING TRACK
 304     MakeTrack(massDS, dst->Pt(), dst->Phi(), dst->Eta(), 3, nIDs);
 305     if(fDebugV2) printf("   ᶫInjecting DStar candidate %d\n",iEntry);
 306   }
 307 }
 308 //______________________________________________________________________________
 309 void AliAnalysisTaskFlowD2H::FillDplustoKpipi(const AliAODEvent *theAOD )
 310 {
 311   // Fill D+ to Kpipi
 312   TList *listHF = (TList*) theAOD->GetList();
 313   if(!listHF) return;
 314   TClonesArray *listDplus = (TClonesArray*) listHF->FindObject("Charm3Prong");
 315   if(!listDplus) return;
 316   int nEntries = listDplus->GetEntriesFast();
 317   if( !nEntries ) return;
 318   fEvent->Fill( 2 ); // EVERYTHING OKAY
 319
 320   AliRDHFCutsDplustoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDplustoKpipi*) fCutsPOI;
 321   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 322   Int_t nIDs[3];
 323   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 324     AliAODRecoDecayHF3Prong *dplu = (AliAODRecoDecayHF3Prong*) listDplus->UncheckedAt( iEntry );
 325     if( !dplu ) continue;
 326     // APPLYING CUTS
 327     if( !dplu->HasSelectionBit(AliRDHFCuts::kDplusCuts) ) continue;
 328     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dplu->Pt(),dplu->YDplus()) )continue;
 329     Int_t topCut = fCutsDStartoKpipi->IsSelected( dplu, AliRDHFCuts::kAll );
 330     if(!topCut) continue;
 331     Double_t massDp=dplu->InvMassDplus();
 332     // TO HANDLE AUTOCORRELATIONS
 333     nIDs[0] = ((AliAODTrack*)dplu->GetDaughter(0))->GetID();
 334     nIDs[1] = ((AliAODTrack*)dplu->GetDaughter(1))->GetID();
 335     nIDs[2] = ((AliAODTrack*)dplu->GetDaughter(2))->GetID();
 336     // ADDING TRACK
 337     MakeTrack(massDp, dplu->Pt(), dplu->Phi(), dplu->Eta(), 3, nIDs);
 338     if(fDebugV2) printf("   ᶫInjecting Dplus candidate %d\n",iEntry);
 339   }
 340 }
 341 //_______________________________________________________________________________
 342 void AliAnalysisTaskFlowD2H::MakeTrack( Double_t mass, Double_t pt, Double_t phi,
 343                                         Double_t eta, Int_t nDau, const Int_t *iID ) {
 344   // create track for flow tasks
 345   Bool_t overwrite = kTRUE;
 346   AliFlowCandidateTrack *oTrack = (static_cast<AliFlowCandidateTrack*> (fCandidates->At( fCandidates->GetLast()+1 )));
 347   if( !oTrack ) { // creates new
 348     oTrack = new AliFlowCandidateTrack();
 349     overwrite = kFALSE;
 350   } else { // overwrites
 351     oTrack->ClearMe();
 352   }
 353   oTrack->SetMass(mass);
 354   oTrack->SetPt(pt);
 355   oTrack->SetPhi(phi);
 356   oTrack->SetEta(eta);
 357   for(Int_t iDau=0; iDau!=nDau; ++iDau)
 358     oTrack->AddDaughter(iID[iDau]);
 359   oTrack->SetForPOISelection(kTRUE);
 360   oTrack->SetForRPSelection(kFALSE);
 361   if(overwrite) {
 362     fCandidates->SetLast( fCandidates->GetLast()+1 );
 363   } else {
 364     fCandidates->AddLast(oTrack);
 365   }
 366   return;
 367 }
 368
 369 void AliAnalysisTaskFlowD2H::SetCommonConstants(Int_t massBins, Double_t minMass, Double_t maxMass, Int_t ptWidth) {
 370   fMassBins = massBins;
 371   fMinMass = minMass;
 372   fMaxMass = maxMass;
 373   fPtBinWidth = ptWidth;
 374 }