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), fSwap(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), fSwap(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(500);
 129   cc->SetNbinsPt(24/fPtBinWidth);
 130   cc->SetNbinsPhi(360);
 131   cc->SetNbinsEta(90);
 132   cc->SetNbinsQ(100);
 133   cc->SetNbinsMass( fMassBins );
 134   cc->SetMultMin(0);
 135   cc->SetMultMax(500);
 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     fTPCEvent->SetNumberOfPOIs(fTPCEvent->GetNumberOfPOIs()+1);
 227     fVZEEvent->SetNumberOfPOIs(fVZEEvent->GetNumberOfPOIs()+1);
 228   }
 229   if(fDebugV2) printf("TPCevent %d | VZEevent %d\n", fTPCEvent->NumberOfTracks(), fVZEEvent->NumberOfTracks() );
 230
 231   PostData(1,fHList);
 232   PostData(2,fTPCEvent);
 233   PostData(3,fVZEEvent);
 234 }
 235 //______________________________________________________________________________
 236 void AliAnalysisTaskFlowD2H::FillD0toKpi(const AliAODEvent *theAOD)
 237 {
 238   // Fill D0->Kpi histos
 239   TList *listHF = (TList*) theAOD->GetList();
 240   if(!listHF) return;
 241   TClonesArray *listDzero = (TClonesArray*) listHF->FindObject("D0toKpi");
 242   if(!listDzero) return;
 243   int nEntries = listDzero->GetEntriesFast();
 244   if( !nEntries ) return;
 245   fEvent->Fill( 2 ); // EVERYTHING OKAY
 246
 247   AliRDHFCutsD0toKpi *fCutsD0toKpi = (AliRDHFCutsD0toKpi*) fCutsPOI;
 248   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 249   Int_t nIDs[2];
 250   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 251     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*) listDzero->UncheckedAt( iEntry );
 252     if( !d0cand ) continue;
 253     // APPLYING CUTS
 254     if( !d0cand->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts) ) continue;
 255     if( !fCutsD0toKpi->IsInFiducialAcceptance(d0cand->Pt(),d0cand->Y(421)) )continue;
 256     Int_t topCut = fCutsD0toKpi->IsSelected( d0cand, AliRDHFCuts::kAll, NULL );
 257     if(!topCut) continue;
 258     Double_t massD0=topCut>1?d0cand->InvMassD0bar():d0cand->InvMassD0();
 259     if((fSwap)&&(topCut==3)) {
 260       massD0=d0cand->InvMassD0();
 261     }
 262
 263     // TO HANDLE AUTOCORRELATIONS
 264     nIDs[0] = ( (AliAODTrack*)d0cand->GetDaughter(0) )->GetID();
 265     nIDs[1] = ( (AliAODTrack*)d0cand->GetDaughter(1) )->GetID();
 266     // ADDING TRACK
 267     MakeTrack(massD0, d0cand->Pt(), d0cand->Phi(), d0cand->Eta(), 2, nIDs);
 268     if(fDebugV2) printf("   ᶫInjecting D0 candidate\n");
 269   }
 270 }
 271 //______________________________________________________________________________
 272 void AliAnalysisTaskFlowD2H::FillDStartoKpipi(const AliAODEvent *theAOD )
 273 {
 274   // Fills D* to Kpipi
 275   TList *listHF = (TList*) theAOD->GetList();
 276   if(!listHF) return;
 277   TClonesArray *listDstar = (TClonesArray*) listHF->FindObject("Dstar");
 278   if(!listDstar) return;
 279   int nEntries = listDstar->GetEntriesFast();
 280   if( !nEntries ) return;
 281   fEvent->Fill( 2 ); // EVERYTHING OKAY
 282
 283   AliRDHFCutsDStartoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDStartoKpipi*) fCutsPOI;
 284   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 285   Int_t nIDs[3];
 286   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 287     AliAODRecoCascadeHF *dst = (AliAODRecoCascadeHF*) listDstar->UncheckedAt( iEntry );
 288     if( !dst ) continue;
 289     AliAODRecoDecayHF2Prong *d0cand = (AliAODRecoDecayHF2Prong*)dst->Get2Prong();
 290     if(!d0cand) continue;
 291     if( !d0cand->GetDaughter(0) ) continue;
 292     if( !d0cand->GetDaughter(1) ) continue;
 293     if( !dst->GetBachelor() ) continue;
 294     // APPLYING CUTS
 295     if( !dst->HasSelectionBit(AliRDHFCuts::kDstarCuts) ) continue;
 296     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dst->Pt(),dst->YDstar()) )continue;
 297     Int_t topCut=0;
 298     topCut = fCutsDStartoKpipi->IsSelected( dst, AliRDHFCuts::kAll );
 299     if(!topCut) continue;
 300     Double_t massDS=dst->DeltaInvMass();
 301     // TO HANDLE AUTOCORRELATIONS
 302     nIDs[0] = ((AliAODTrack*)d0cand->GetDaughter(0))->GetID();
 303     nIDs[1] = ((AliAODTrack*)d0cand->GetDaughter(1))->GetID();
 304     nIDs[2] = ((AliAODTrack*)dst->GetBachelor() )->GetID();
 305     // ADDING TRACK
 306     MakeTrack(massDS, dst->Pt(), dst->Phi(), dst->Eta(), 3, nIDs);
 307     if(fDebugV2) printf("   ᶫInjecting DStar candidate %d\n",iEntry);
 308   }
 309 }
 310 //______________________________________________________________________________
 311 void AliAnalysisTaskFlowD2H::FillDplustoKpipi(const AliAODEvent *theAOD )
 312 {
 313   // Fill D+ to Kpipi
 314   TList *listHF = (TList*) theAOD->GetList();
 315   if(!listHF) return;
 316   TClonesArray *listDplus = (TClonesArray*) listHF->FindObject("Charm3Prong");
 317   if(!listDplus) return;
 318   int nEntries = listDplus->GetEntriesFast();
 319   if( !nEntries ) return;
 320   fEvent->Fill( 2 ); // EVERYTHING OKAY
 321
 322   AliRDHFCutsDplustoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDplustoKpipi*) fCutsPOI;
 323   if (fDebugV2) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 324   Int_t nIDs[3];
 325   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 326     AliAODRecoDecayHF3Prong *dplu = (AliAODRecoDecayHF3Prong*) listDplus->UncheckedAt( iEntry );
 327     if( !dplu ) continue;
 328     // APPLYING CUTS
 329     if( !dplu->HasSelectionBit(AliRDHFCuts::kDplusCuts) ) continue;
 330     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(dplu->Pt(),dplu->YDplus()) )continue;
 331     Int_t topCut = fCutsDStartoKpipi->IsSelected( dplu, AliRDHFCuts::kAll );
 332     if(!topCut) continue;
 333     Double_t massDp=dplu->InvMassDplus();
 334     // TO HANDLE AUTOCORRELATIONS
 335     nIDs[0] = ((AliAODTrack*)dplu->GetDaughter(0))->GetID();
 336     nIDs[1] = ((AliAODTrack*)dplu->GetDaughter(1))->GetID();
 337     nIDs[2] = ((AliAODTrack*)dplu->GetDaughter(2))->GetID();
 338     // ADDING TRACK
 339     MakeTrack(massDp, dplu->Pt(), dplu->Phi(), dplu->Eta(), 3, nIDs);
 340     if(fDebugV2) printf("   ᶫInjecting Dplus candidate %d\n",iEntry);
 341   }
 342 }
 343 //_______________________________________________________________________________
 344 void AliAnalysisTaskFlowD2H::MakeTrack( Double_t mass, Double_t pt, Double_t phi,
 345                                         Double_t eta, Int_t nDau, const Int_t *iID ) {
 346   // create track for flow tasks
 347   Bool_t overwrite = kTRUE;
 348   AliFlowCandidateTrack *oTrack = (static_cast<AliFlowCandidateTrack*> (fCandidates->At( fCandidates->GetLast()+1 )));
 349   if( !oTrack ) { // creates new
 350     oTrack = new AliFlowCandidateTrack();
 351     overwrite = kFALSE;
 352   } else { // overwrites
 353     oTrack->ClearMe();
 354   }
 355   oTrack->SetMass(mass);
 356   oTrack->SetPt(pt);
 357   oTrack->SetPhi(phi);
 358   oTrack->SetEta(eta);
 359   for(Int_t iDau=0; iDau!=nDau; ++iDau)
 360     oTrack->AddDaughter(iID[iDau]);
 361   oTrack->SetForPOISelection(kTRUE);
 362   oTrack->SetForRPSelection(kFALSE);
 363   if(overwrite) {
 364     fCandidates->SetLast( fCandidates->GetLast()+1 );
 365   } else {
 366     fCandidates->AddLast(oTrack);
 367   }
 368   return;
 369 }
 370
 371 void AliAnalysisTaskFlowD2H::SetCommonConstants(Int_t massBins, Double_t minMass, Double_t maxMass, Int_t ptWidth) {
 372   fMassBins = massBins;
 373   fMinMass = minMass;
 374   fMaxMass = maxMass;
 375   fPtBinWidth = ptWidth;
 376 }
 377
 378 //=======================================================================
 379 void AliAnalysisTaskFlowD2H::Terminate(Option_t *)
 380 {
 381   //terminate
 382 }