PWG3/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   fDebug(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   fDebug(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 {
 117   if(fHList)
 118     delete fHList;
 119 }
 120
 121 //_____________________________________________________________________________
 122 void AliAnalysisTaskFlowD2H::UserCreateOutputObjects()
 123 {
 124   if(fDebug)
 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   TList *tEvents = new TList();
 156     tEvents->SetName("Events");
 157     tEvents->SetOwner();
 158     for(int i=0; i!=2; ++i) {
 159       fEvent[i] = new TH2D(Form("Event%d",i),"Events;V0M;Arb",10,0,100,10,0,12000);
 160       tEvents->Add(fEvent[i]);
 161     }
 162     fAnaCuts=new TProfile("Cuts","Analysis Cuts", 10,0,10);
 163       fAnaCuts->Fill(0.5,fPOIEta[0],1);fAnaCuts->GetXaxis()->SetBinLabel(1,"ETAm");
 164       fAnaCuts->Fill(1.5,fPOIEta[1],1);fAnaCuts->GetXaxis()->SetBinLabel(2,"ETAM");
 165     tEvents->Add(fAnaCuts);
 166     fHList->Add(tEvents);
 167   TList *tCandidates;
 168     tCandidates = new TList();
 169     tCandidates->SetOwner();
 170     tCandidates->SetName(Form("Candidates%d",fSource));
 171     Double_t dBinMin=0, dBinMax=3;
 172     Int_t nBins=3;
 173     switch (fSource) {
 174       case (AliRDHFCuts::kD0toKpiCuts): // 360/72 (bw=5MeV)
 175         dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
 176       case (AliRDHFCuts::kDstarCuts): // 36/72 (bw=0.5MeV)
 177         dBinMin = 0.137; dBinMax = 0.173; nBins=72; break;
 178       case (AliRDHFCuts::kDplusCuts): // 360/72 (bw=5MeV)
 179         dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
 180     }
 181     for(int i=0; i!=2; ++i) {
 182       fMass[i] = new TH2D( Form("Mass%d",i),
 183                            Form("Mass%d;Mass [GeV];Pt [GeV]",i),
 184                            nBins, dBinMin, dBinMax,
 185                            fFlowPts[1]-fFlowPts[0], fFlowPts[0], fFlowPts[1] );
 186       tCandidates->Add(fMass[i]);
 187     }
 188     fHList->Add(tCandidates);
 189   if (fDebug) printf("DEBUGGER: Created histos for DMeson %d\n", fSource );
 190 }
 191 //_____________________________________________________________________________
 192 void AliAnalysisTaskFlowD2H::NotifyRun()
 193 {
 194 }
 195 //_____________________________________________________________________________
 196 void AliAnalysisTaskFlowD2H::UserExec(Option_t *)
 197 {
 198   AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
 199
 200   if(!fAOD)
 201     return;
 202   Int_t iMulti = fAOD->GetNTracks();  /// TO DO
 203   fEvent[0]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );
 204   if(!fCutsPOI->IsEventSelected(fAOD)) return;
 205   if(fCutsPOI->IsEventSelectedInCentrality(fAOD)>0) return;
 206   fEvent[1]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );
 207
 208   if (fDebug) printf("Event selected\n");
 209   fCutsRP->SetEvent(fAOD,MCEvent());
 210   AliFlowTrackCuts* dummy = new AliFlowTrackCuts("null_cuts");
 211   dummy->SetParamType(AliFlowTrackCuts::kGlobal);
 212   dummy->SetPtRange(+1,-1); // select nothing QUICK
 213   dummy->SetEtaRange(+1,-1); // select nothing VZERO
 214   dummy->SetEvent(fAOD,MCEvent());
 215   AliFlowEvent *flowEvent[5];
 216   for(int r=0; r!=5; ++r) {
 217     flowEvent[r] = new AliFlowEvent(fCutsRP,dummy);
 218     flowEvent[r]->SetReferenceMultiplicity( iMulti );
 219     flowEvent[r]->DefineDeadZone(0,0,0,0);
 220     flowEvent[r]->TagSubeventsInEta( fFlowEta[0], fFlowEta[1],
 221                                      fFlowEta[2], fFlowEta[3] );
 222   }
 223   delete dummy;
 224   if (fDebug) printf(" ᶫFlowEvent has %d RPs\n", flowEvent[0]->NumberOfTracks() );
 225
 226   switch (fSource) {
 227     case (AliRDHFCuts::kD0toKpiCuts):
 228       FillD0toKpi(fAOD,flowEvent); break;
 229     case (AliRDHFCuts::kDstarCuts):
 230       FillDStartoKpipi(fAOD,flowEvent); break;
 231     case (AliRDHFCuts::kDplusCuts):
 232       FillDplustoKpipi(fAOD,flowEvent); break;
 233   }
 234
 235   for(int m=0; m!=5; ++m)
 236     PostData(2+m,flowEvent[m]);
 237   PostData(1,fHList);
 238
 239 }
 240 //______________________________________________________________________________
 241 void AliAnalysisTaskFlowD2H::FillD0toKpi(AliAODEvent *theAOD,
 242                                          AliFlowEvent *theMB[5] )
 243 {
 244   TList *listHF = (TList*) theAOD->GetList();
 245   if(!listHF) return;
 246   TClonesArray *listDzero = (TClonesArray*) listHF->FindObject("D0toKpi");
 247   if(!listDzero) return;
 248   int nEntries = listDzero->GetEntriesFast();
 249   if( !nEntries ) return;
 250   AliRDHFCutsD0toKpi *fCutsD0toKpi = (AliRDHFCutsD0toKpi*) fCutsPOI;
 251   if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 252   const Int_t ndght=2;
 253   Int_t nIDs[ndght];
 254   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 255     AliAODRecoDecayHF2Prong *D0 = (AliAODRecoDecayHF2Prong*) listDzero->UncheckedAt( iEntry );
 256     if( !D0 ) continue;
 257     if( !D0->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts) ) continue;
 258     if( !fCutsD0toKpi->IsInFiducialAcceptance(D0->Pt(),D0->Y(421)) )continue;
 259     int topCut  = fCutsD0toKpi->IsSelected( D0, AliRDHFCuts::kAll, theAOD );
 260     int nLevel=topCut>0?1:0;
 261     if( (D0->Eta()<fPOIEta[0])||(D0->Eta()>fPOIEta[1]) )
 262       nLevel=0;
 263     Double_t MassD0=topCut>1?D0->InvMassD0bar():D0->InvMassD0();
 264     for(int h=0; h!=nLevel+1; ++h)
 265       fMass[h]->Fill(MassD0,D0->Pt());
 266     if( (D0->Pt()<fFlowPts[0]) || (D0->Pt()>fFlowPts[1]) ) continue;
 267     AliAODTrack* iT;
 268     for(Int_t i=0; i!=ndght; ++i) {
 269       iT = (AliAODTrack*)D0->GetDaughter(i);
 270       nIDs[i] = iT->GetID();
 271     }
 272     // Candidates Insertion (done in filling method: faster)
 273     if(nLevel)
 274       for(Int_t r=0; r!=5; ++r)
 275         if( (MassD0>=fFlowBands[0][r]) && (MassD0<fFlowBands[1][r]) ) {
 276           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 277                                  MakeTrack(MassD0, D0->Pt(), D0->Phi(),
 278                                            D0->Eta(), ndght, nIDs);
 279           if(fDebug) printf("   ᶫInjecting D0 candidate on band %d \n", r);
 280           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 281             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 282               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 283               if(!iRP->InRPSelection()) continue;
 284               if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
 285                 sTrack->SetDaughter(iDau,iRP);
 286                 iRP->SetForRPSelection(kFALSE);
 287                 if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 288               }
 289             }
 290           theMB[r]->AddTrack(sTrack);
 291         }
 292     // END of Candidates Insertion
 293   }
 294 }
 295 //______________________________________________________________________________
 296 void AliAnalysisTaskFlowD2H::FillDStartoKpipi(AliAODEvent *theAOD,
 297                                               AliFlowEvent *theMB[5] )
 298 {
 299   TList *listHF = (TList*) theAOD->GetList();
 300   if(!listHF) return;
 301   TClonesArray *listDstar = (TClonesArray*) listHF->FindObject("Dstar");
 302   if(!listDstar) return;
 303   int nEntries = listDstar->GetEntriesFast();
 304   if( !nEntries ) return;
 305   AliRDHFCutsDStartoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDStartoKpipi*) fCutsPOI;
 306   if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 307   const Int_t ndght=3;
 308   Int_t nIDs[ndght];
 309   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 310     AliAODRecoCascadeHF *DS = (AliAODRecoCascadeHF*) listDstar->UncheckedAt( iEntry );
 311     if( !DS ) continue;
 312     if( !DS->HasSelectionBit(AliRDHFCuts::kDstarCuts) ) continue;
 313     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(DS->Pt(),DS->YDstar()) )continue;
 314     int topCut = fCutsDStartoKpipi->IsSelected( DS, AliRDHFCuts::kAll );
 315     int nLevel=topCut>0?1:0;
 316     if( (DS->Eta()<fPOIEta[0])||(DS->Eta()>fPOIEta[1]) )
 317       nLevel=0;
 318     Double_t MassDS=DS->DeltaInvMass();
 319     for(int h=0; h!=nLevel+1; ++h)
 320       fMass[h]->Fill(MassDS,DS->Pt());
 321     if( (DS->Pt()<fFlowPts[0]) || (DS->Pt()>fFlowPts[1]) ) continue;
 322     AliAODRecoDecayHF2Prong *D0 = (AliAODRecoDecayHF2Prong*)DS->Get2Prong();
 323     nIDs[0] = ((AliAODTrack*)D0->GetDaughter(0))->GetID();
 324     nIDs[1] = ((AliAODTrack*)D0->GetDaughter(1))->GetID();
 325     nIDs[2] = ((AliAODTrack*)DS->GetBachelor() )->GetID();
 326     // Candidates Insertion (done in filling method: faster)
 327     if(nLevel)
 328       for(Int_t r=0; r!=5; ++r)
 329         if( (MassDS>=fFlowBands[0][r]) && (MassDS<fFlowBands[1][r]) ) {
 330           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 331                                  MakeTrack(MassDS, DS->Pt(), DS->Phi(),
 332                                            DS->Eta(), ndght, nIDs);
 333           if(fDebug) printf("   ᶫInjecting DStar candidate on band %d \n", r);
 334           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 335             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 336               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 337               if(!iRP->InRPSelection()) continue;
 338               if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
 339                 sTrack->SetDaughter(iDau,iRP);
 340                 iRP->SetForRPSelection(kFALSE);
 341                 if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 342               }
 343             }
 344           theMB[r]->AddTrack(sTrack);
 345         }
 346     // END of Candidates Insertion
 347   }
 348 }
 349 //______________________________________________________________________________
 350 void AliAnalysisTaskFlowD2H::FillDplustoKpipi(AliAODEvent *theAOD,
 351                                               AliFlowEvent *theMB[5] )
 352 {
 353   TList *listHF = (TList*) theAOD->GetList();
 354   if(!listHF) return;
 355   TClonesArray *listDplus = (TClonesArray*) listHF->FindObject("Charm3Prong");
 356   if(!listDplus) return;
 357   int nEntries = listDplus->GetEntriesFast();
 358   if( !nEntries ) return;
 359   AliRDHFCutsDplustoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDplustoKpipi*) fCutsPOI;
 360   if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
 361   const Int_t ndght=3;
 362   Int_t nIDs[ndght];
 363   for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
 364     AliAODRecoDecayHF3Prong *Dp = (AliAODRecoDecayHF3Prong*) listDplus->UncheckedAt( iEntry );
 365     if( !Dp ) continue;
 366     if( !Dp->HasSelectionBit(AliRDHFCuts::kDplusCuts) ) continue;
 367     if( !fCutsDStartoKpipi->IsInFiducialAcceptance(Dp->Pt(),Dp->YDplus()) )continue;
 368     int topCut = fCutsDStartoKpipi->IsSelected( Dp, AliRDHFCuts::kAll );
 369     int nLevel=topCut>0?1:0;
 370     if( (Dp->Eta()<fPOIEta[0])||(Dp->Eta()>fPOIEta[1]) )
 371       nLevel=0;
 372     Double_t MassDp=Dp->InvMassDplus();
 373     for(int h=0; h!=nLevel+1; ++h)
 374       fMass[h]->Fill(MassDp,Dp->Pt());
 375     if( (Dp->Pt()<fFlowPts[0]) || (Dp->Pt()>fFlowPts[1]) ) continue;
 376     nIDs[0] = ((AliAODTrack*)Dp->GetDaughter(0))->GetID();
 377     nIDs[1] = ((AliAODTrack*)Dp->GetDaughter(1))->GetID();
 378     nIDs[2] = ((AliAODTrack*)Dp->GetDaughter(2))->GetID();
 379     // Candidates Insertion (done in filling method: faster)
 380     if(nLevel)
 381       for(Int_t r=0; r!=5; ++r)
 382         if( (MassDp>=fFlowBands[0][r]) && (MassDp<fFlowBands[1][r]) ) {
 383           AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
 384                                  MakeTrack(MassDp, Dp->Pt(), Dp->Phi(),
 385                                            Dp->Eta(), ndght, nIDs);
 386           if(fDebug) printf("   ᶫInjecting DStar candidate on band %d \n", r);
 387           for(Int_t iDau=0; iDau!=ndght; ++iDau)
 388             for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
 389               AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
 390               if(!iRP->InRPSelection()) continue;
 391               if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
 392                 sTrack->SetDaughter(iDau,iRP);
 393                 iRP->SetForRPSelection(kFALSE);
 394                 if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
 395               }
 396             }
 397           theMB[r]->AddTrack(sTrack);
 398         }
 399     // END of Candidates Insertion
 400   }
 401 }
 402 //______________________________________________________________________________
 403 AliFlowCandidateTrack* AliAnalysisTaskFlowD2H::MakeTrack( Double_t mass,
 404                           Double_t pt, Double_t phi, Double_t eta,
 405                           Int_t nDau, Int_t *iID ) {
 406   AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
 407   sTrack->SetMass(mass);
 408   sTrack->SetPt(pt);
 409   sTrack->SetPhi(phi);
 410   sTrack->SetEta(eta);
 411   for(Int_t iDau=0; iDau!=nDau; ++iDau)
 412     sTrack->AddDaughter(iID[iDau]);
 413   sTrack->SetForPOISelection(kTRUE);
 414   sTrack->SetForRPSelection(kFALSE);
 415   return sTrack;
 416 }
 417 //_____________________________________________________________________________
 418 void AliAnalysisTaskFlowD2H::Terminate(Option_t *)
 419 {
 420
 421 }
 422