[u/mrichter/AliRoot.git] / PWG3 / vertexingHF / charmFlow / AliAnalysisTaskFlowD2H.cxx

/**************************************************************************
* Copyright(c) 1998-2008,ALICE Experiment at CERN,All rights reserved.    *
*                                                                         *
* Author: The ALICE Off-line Project.                                     *
* Contributors are mentioned in the code where appropriate.               *
*                                                                         *
* Permission to use,copy,modify and distribute this software and its      *
* documentation strictly for non-commercial purposes is hereby granted    *
* without fee,provided that the above copyright notice appears in all     *
* copies and that both the copyright notice and this permission notice    *
* appear in the supporting documentation. The authors make no claims      *
* about the suitability of this software for any purpose. It is           *
* provided "as is" without express or implied warranty.                   *
**************************************************************************/

//==============================================================================
// FlowD2H main task:
// >> Make flowEvent with RPcuts and POIcuts given in constructor and passes it
//    to the daughter tasks.
// >> The POIcuts are polymorphic based on the AliRDHFCuts class allowing the 
//    use of all charmed candidates reconstructed in the central barrel.
// Author: Carlos Perez (cperez@cern.ch)
//==============================================================================

#include "TChain.h"
#include "TList.h"
#include "TH2D.h"
#include "TProfile.h"

#include "AliAnalysisTaskSE.h"
#include "AliAnalysisManager.h"

#include "AliCentrality.h"

#include "AliAODEvent.h"
#include "AliAODInputHandler.h"
#include "AliAODTrack.h"

#include "TMath.h"
#include "TObjArray.h"

#include "AliFlowCandidateTrack.h"
#include "AliFlowTrackCuts.h"
#include "AliFlowEvent.h"
#include "AliFlowCommonConstants.h"

#include "AliRDHFCuts.h"
#include "AliRDHFCutsD0toKpi.h"
#include "AliRDHFCutsDStartoKpipi.h"
#include "AliRDHFCutsDplustoKpipi.h"

#include "AliAODRecoDecayHF2Prong.h"
#include "AliAODRecoDecayHF3Prong.h"
#include "AliAODRecoCascadeHF.h"

#include "AliAnalysisTaskFlowD2H.h"

ClassImp(AliAnalysisTaskFlowD2H)

//_____________________________________________________________________________
AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H() :
  AliAnalysisTaskSE(), fCutsRP(NULL), fCutsPOI(NULL), fSource(0),
  fDebug(kFALSE), fHList(NULL), fAnaCuts(NULL)
{
// Default constructor
  for(int i=0; i!=2; ++i)
    fEvent[i] = NULL;
  for(int i=0; i!=2; ++i)
    fMass[i] = NULL;
  for(int i=0; i!=2; ++i)
    fPOIEta[i] = 0;
  for(int i=0; i!=4; ++i)
    fFlowEta[i] = 0;
  for(int x=0; x!=2; ++x)
    fFlowPts[x] = 0;
  for(int x=0; x!=2; ++x)
    for(int m=0; m!=5; ++m)
      fFlowBands[x][m] = 0;
}

//_____________________________________________________________________________
AliAnalysisTaskFlowD2H::AliAnalysisTaskFlowD2H(const char *name,
    AliFlowTrackCuts *cutsRPs, AliRDHFCuts *cutsPOIs, Int_t specie) :
  AliAnalysisTaskSE(name), fCutsRP(cutsRPs), 
  fCutsPOI(cutsPOIs), fSource(specie),
  fDebug(kFALSE), fHList(NULL), fAnaCuts(NULL)
{
// Standard constructor
  for(int i=0; i!=2; ++i)
    fEvent[i] = NULL;
  for(int i=0; i!=2; ++i)
    fMass[i] = NULL;
  for(int i=0; i!=2; ++i)
    fPOIEta[i] = 0;
  for(int i=0; i!=4; ++i)
    fFlowEta[i] = 0;
  for(int x=0; x!=2; ++x)
    fFlowPts[x] = 0;
  for(int x=0; x!=2; ++x)
    for(int m=0; m!=5; ++m)
      fFlowBands[x][m] = 0;

  DefineInput( 0,TChain::Class());
  DefineOutput(1,TList::Class());
  DefineOutput(2,AliFlowEventSimple::Class()); // first band
  DefineOutput(3,AliFlowEventSimple::Class()); // second band
  DefineOutput(4,AliFlowEventSimple::Class()); // third band
  DefineOutput(5,AliFlowEventSimple::Class()); // fourth band
  DefineOutput(6,AliFlowEventSimple::Class()); // fifth band
}

//_____________________________________________________________________________
AliAnalysisTaskFlowD2H::~AliAnalysisTaskFlowD2H()
{
  if(fHList)
    delete fHList;
}

//_____________________________________________________________________________
void AliAnalysisTaskFlowD2H::UserCreateOutputObjects()
{
  if(fDebug)
    printf("DEBUGGER: Creating output\n");
  fHList = new TList();
  fHList->SetOwner();
  AddHistograms();
  PostData(1,fHList);

  AliFlowCommonConstants* cc = AliFlowCommonConstants::GetMaster();
  cc->SetNbinsMult(10000);
  cc->SetMultMin(0);
  cc->SetMultMax(10000);

  cc->SetNbinsPt(fFlowPts[1]-fFlowPts[0]);
  cc->SetPtMin(fFlowPts[0]);
  cc->SetPtMax(fFlowPts[1]);

  cc->SetNbinsPhi(180);
  cc->SetPhiMin(0.0);
  cc->SetPhiMax(TMath::TwoPi());

  cc->SetNbinsEta(200);
  cc->SetEtaMin(-5.0);
  cc->SetEtaMax(+5.0);

  cc->SetNbinsQ(500);
  cc->SetQMin(0.0);
  cc->SetQMax(3.0);
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowD2H::AddHistograms()
{
  TList *tEvents = new TList();
    tEvents->SetName("Events");
    tEvents->SetOwner();
    for(int i=0; i!=2; ++i) {
      fEvent[i] = new TH2D(Form("Event%d",i),"Events;V0M;Arb",10,0,100,10,0,12000);
      tEvents->Add(fEvent[i]);
    }
    fAnaCuts=new TProfile("Cuts","Analysis Cuts", 10,0,10);
      fAnaCuts->Fill(0.5,fPOIEta[0],1);fAnaCuts->GetXaxis()->SetBinLabel(1,"ETAm");
      fAnaCuts->Fill(1.5,fPOIEta[1],1);fAnaCuts->GetXaxis()->SetBinLabel(2,"ETAM");
    tEvents->Add(fAnaCuts);
    fHList->Add(tEvents);
  TList *tCandidates;
    tCandidates = new TList();
    tCandidates->SetOwner();
    tCandidates->SetName(Form("Candidates%d",fSource));
    Double_t dBinMin=0, dBinMax=3;
    Int_t nBins=3;
    switch (fSource) {
      case (AliRDHFCuts::kD0toKpiCuts): // 360/72 (bw=5MeV)
        dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
      case (AliRDHFCuts::kDstarCuts): // 36/72 (bw=0.5MeV)
        dBinMin = 0.137; dBinMax = 0.173; nBins=72; break;
      case (AliRDHFCuts::kDplusCuts): // 360/72 (bw=5MeV)
        dBinMin = 1.695; dBinMax = 2.055; nBins=72; break;
    }
    for(int i=0; i!=2; ++i) {
      fMass[i] = new TH2D( Form("Mass%d",i),
                           Form("Mass%d;Mass [GeV];Pt [GeV]",i),
                           nBins, dBinMin, dBinMax,
                           fFlowPts[1]-fFlowPts[0], fFlowPts[0], fFlowPts[1] );
      tCandidates->Add(fMass[i]);
    }
    fHList->Add(tCandidates);
  if (fDebug) printf("DEBUGGER: Created histos for DMeson %d\n", fSource );
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowD2H::NotifyRun()
{
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowD2H::UserExec(Option_t *)
{
  AliAODEvent *fAOD = dynamic_cast<AliAODEvent*>(InputEvent());

  if(!fAOD)
    return;
  Int_t iMulti = fAOD->GetNTracks();  /// TO DO
  fEvent[0]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );
  if(!fCutsPOI->IsEventSelected(fAOD)) return;
  if(fCutsPOI->IsEventSelectedInCentrality(fAOD)>0) return;
  fEvent[1]->Fill( fCutsPOI->GetCentrality(fAOD), iMulti );

  if (fDebug) printf("Event selected\n");
  fCutsRP->SetEvent(fAOD,MCEvent());
  AliFlowTrackCuts* dummy = new AliFlowTrackCuts("null_cuts");
  dummy->SetParamType(AliFlowTrackCuts::kGlobal);
  dummy->SetPtRange(+1,-1); // select nothing QUICK
  dummy->SetEtaRange(+1,-1); // select nothing VZERO
  dummy->SetEvent(fAOD,MCEvent());
  AliFlowEvent *flowEvent[5];
  for(int r=0; r!=5; ++r) {
    flowEvent[r] = new AliFlowEvent(fCutsRP,dummy);
    flowEvent[r]->SetReferenceMultiplicity( iMulti );
    flowEvent[r]->DefineDeadZone(0,0,0,0);
    flowEvent[r]->TagSubeventsInEta( fFlowEta[0], fFlowEta[1],
                                     fFlowEta[2], fFlowEta[3] );
  }
  delete dummy;
  if (fDebug) printf(" ᶫFlowEvent has %d RPs\n", flowEvent[0]->NumberOfTracks() );

  switch (fSource) {
    case (AliRDHFCuts::kD0toKpiCuts):
      FillD0toKpi(fAOD,flowEvent); break;
    case (AliRDHFCuts::kDstarCuts):
      FillDStartoKpipi(fAOD,flowEvent); break;
    case (AliRDHFCuts::kDplusCuts):
      FillDplustoKpipi(fAOD,flowEvent); break;
  }

  for(int m=0; m!=5; ++m)
    PostData(2+m,flowEvent[m]);
  PostData(1,fHList);

}
//______________________________________________________________________________
void AliAnalysisTaskFlowD2H::FillD0toKpi(AliAODEvent *theAOD, 
                                         AliFlowEvent *theMB[5] )
{
  TList *listHF = (TList*) theAOD->GetList();
  if(!listHF) return;
  TClonesArray *listDzero = (TClonesArray*) listHF->FindObject("D0toKpi");
  if(!listDzero) return;
  int nEntries = listDzero->GetEntriesFast();
  if( !nEntries ) return;
  AliRDHFCutsD0toKpi *fCutsD0toKpi = (AliRDHFCutsD0toKpi*) fCutsPOI;
  if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
  const Int_t ndght=2;
  Int_t nIDs[ndght];
  for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
    AliAODRecoDecayHF2Prong *D0 = (AliAODRecoDecayHF2Prong*) listDzero->UncheckedAt( iEntry );
    if( !D0 ) continue;
    if( !D0->HasSelectionBit(AliRDHFCuts::kD0toKpiCuts) ) continue;
    if( !fCutsD0toKpi->IsInFiducialAcceptance(D0->Pt(),D0->Y(421)) )continue;
    int topCut  = fCutsD0toKpi->IsSelected( D0, AliRDHFCuts::kAll, theAOD );
    int nLevel=topCut>0?1:0;
    if( (D0->Eta()<fPOIEta[0])||(D0->Eta()>fPOIEta[1]) )
      nLevel=0;
    Double_t MassD0=topCut>1?D0->InvMassD0bar():D0->InvMassD0();
    for(int h=0; h!=nLevel+1; ++h)
      fMass[h]->Fill(MassD0,D0->Pt());
    if( (D0->Pt()<fFlowPts[0]) || (D0->Pt()>fFlowPts[1]) ) continue;
    AliAODTrack* iT;
    for(Int_t i=0; i!=ndght; ++i) {
      iT = (AliAODTrack*)D0->GetDaughter(i);
      nIDs[i] = iT->GetID();
    }
    // Candidates Insertion (done in filling method: faster)
    if(nLevel)
      for(Int_t r=0; r!=5; ++r)
        if( (MassD0>=fFlowBands[0][r]) && (MassD0<fFlowBands[1][r]) ) {
          AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
                                 MakeTrack(MassD0, D0->Pt(), D0->Phi(),
                                           D0->Eta(), ndght, nIDs);
          if(fDebug) printf("   ᶫInjecting D0 candidate on band %d \n", r);
          for(Int_t iDau=0; iDau!=ndght; ++iDau)
            for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
              AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
              if(!iRP->InRPSelection()) continue;
              if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
                sTrack->SetDaughter(iDau,iRP);
                iRP->SetForRPSelection(kFALSE);
                if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
              }
            }
          theMB[r]->AddTrack(sTrack);
        }
    // END of Candidates Insertion
  }
}
//______________________________________________________________________________
void AliAnalysisTaskFlowD2H::FillDStartoKpipi(AliAODEvent *theAOD, 
                                              AliFlowEvent *theMB[5] )
{
  TList *listHF = (TList*) theAOD->GetList();
  if(!listHF) return;
  TClonesArray *listDstar = (TClonesArray*) listHF->FindObject("Dstar");
  if(!listDstar) return;
  int nEntries = listDstar->GetEntriesFast();
  if( !nEntries ) return;
  AliRDHFCutsDStartoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDStartoKpipi*) fCutsPOI;
  if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
  const Int_t ndght=3;
  Int_t nIDs[ndght];
  for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
    AliAODRecoCascadeHF *DS = (AliAODRecoCascadeHF*) listDstar->UncheckedAt( iEntry );
    if( !DS ) continue;
    if( !DS->HasSelectionBit(AliRDHFCuts::kDstarCuts) ) continue;
    if( !fCutsDStartoKpipi->IsInFiducialAcceptance(DS->Pt(),DS->YDstar()) )continue;
    int topCut = fCutsDStartoKpipi->IsSelected( DS, AliRDHFCuts::kAll );
    int nLevel=topCut>0?1:0;
    if( (DS->Eta()<fPOIEta[0])||(DS->Eta()>fPOIEta[1]) )
      nLevel=0;
    Double_t MassDS=DS->DeltaInvMass();
    for(int h=0; h!=nLevel+1; ++h)
      fMass[h]->Fill(MassDS,DS->Pt());
    if( (DS->Pt()<fFlowPts[0]) || (DS->Pt()>fFlowPts[1]) ) continue;
    AliAODRecoDecayHF2Prong *D0 = (AliAODRecoDecayHF2Prong*)DS->Get2Prong();
    nIDs[0] = ((AliAODTrack*)D0->GetDaughter(0))->GetID();
    nIDs[1] = ((AliAODTrack*)D0->GetDaughter(1))->GetID();
    nIDs[2] = ((AliAODTrack*)DS->GetBachelor() )->GetID();
    // Candidates Insertion (done in filling method: faster)
    if(nLevel)
      for(Int_t r=0; r!=5; ++r)
        if( (MassDS>=fFlowBands[0][r]) && (MassDS<fFlowBands[1][r]) ) {
          AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
                                 MakeTrack(MassDS, DS->Pt(), DS->Phi(),
                                           DS->Eta(), ndght, nIDs);
          if(fDebug) printf("   ᶫInjecting DStar candidate on band %d \n", r);
          for(Int_t iDau=0; iDau!=ndght; ++iDau)
            for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
              AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
              if(!iRP->InRPSelection()) continue;
              if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
                sTrack->SetDaughter(iDau,iRP);
                iRP->SetForRPSelection(kFALSE);
                if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
              }
            }
          theMB[r]->AddTrack(sTrack);
        }
    // END of Candidates Insertion
  }
}
//______________________________________________________________________________
void AliAnalysisTaskFlowD2H::FillDplustoKpipi(AliAODEvent *theAOD, 
                                              AliFlowEvent *theMB[5] )
{
  TList *listHF = (TList*) theAOD->GetList();
  if(!listHF) return;
  TClonesArray *listDplus = (TClonesArray*) listHF->FindObject("Charm3Prong");
  if(!listDplus) return;
  int nEntries = listDplus->GetEntriesFast();
  if( !nEntries ) return;
  AliRDHFCutsDplustoKpipi *fCutsDStartoKpipi = (AliRDHFCutsDplustoKpipi*) fCutsPOI;
  if (fDebug) printf("  ᶫ%d candidates found. Looping...\n", nEntries);
  const Int_t ndght=3;
  Int_t nIDs[ndght];
  for( int iEntry=0; iEntry!=nEntries; ++iEntry ) {
    AliAODRecoDecayHF3Prong *Dp = (AliAODRecoDecayHF3Prong*) listDplus->UncheckedAt( iEntry );
    if( !Dp ) continue;
    if( !Dp->HasSelectionBit(AliRDHFCuts::kDplusCuts) ) continue;
    if( !fCutsDStartoKpipi->IsInFiducialAcceptance(Dp->Pt(),Dp->YDplus()) )continue;
    int topCut = fCutsDStartoKpipi->IsSelected( Dp, AliRDHFCuts::kAll );
    int nLevel=topCut>0?1:0;
    if( (Dp->Eta()<fPOIEta[0])||(Dp->Eta()>fPOIEta[1]) )
      nLevel=0;
    Double_t MassDp=Dp->InvMassDplus();
    for(int h=0; h!=nLevel+1; ++h)
      fMass[h]->Fill(MassDp,Dp->Pt());
    if( (Dp->Pt()<fFlowPts[0]) || (Dp->Pt()>fFlowPts[1]) ) continue;
    nIDs[0] = ((AliAODTrack*)Dp->GetDaughter(0))->GetID();
    nIDs[1] = ((AliAODTrack*)Dp->GetDaughter(1))->GetID();
    nIDs[2] = ((AliAODTrack*)Dp->GetDaughter(2))->GetID();
    // Candidates Insertion (done in filling method: faster)
    if(nLevel)
      for(Int_t r=0; r!=5; ++r)
        if( (MassDp>=fFlowBands[0][r]) && (MassDp<fFlowBands[1][r]) ) {
          AliFlowCandidateTrack *sTrack = (AliFlowCandidateTrack*)
                                 MakeTrack(MassDp, Dp->Pt(), Dp->Phi(),
                                           Dp->Eta(), ndght, nIDs);
          if(fDebug) printf("   ᶫInjecting DStar candidate on band %d \n", r);
          for(Int_t iDau=0; iDau!=ndght; ++iDau)
            for(Int_t iRPs=0; iRPs!=theMB[r]->NumberOfTracks(); ++iRPs) {
              AliFlowTrack *iRP = (AliFlowTrack*) (theMB[r]->GetTrack(iRPs));
              if(!iRP->InRPSelection()) continue;
              if( fabs(sTrack->GetIDDaughter(iDau)) == fabs(iRP->GetID()) ) {
                sTrack->SetDaughter(iDau,iRP);
                iRP->SetForRPSelection(kFALSE);
                if(fDebug) printf("    ᶫdaughter%d with fID %d was removed from this RP set\n", iDau, sTrack->GetIDDaughter(iDau));
              }
            }
          theMB[r]->AddTrack(sTrack);
        }
    // END of Candidates Insertion
  }
}
//______________________________________________________________________________
AliFlowCandidateTrack* AliAnalysisTaskFlowD2H::MakeTrack( Double_t mass, 
                          Double_t pt, Double_t phi, Double_t eta, 
                          Int_t nDau, Int_t *iID ) {
  AliFlowCandidateTrack *sTrack = new AliFlowCandidateTrack();
  sTrack->SetMass(mass);
  sTrack->SetPt(pt);
  sTrack->SetPhi(phi);
  sTrack->SetEta(eta);
  for(Int_t iDau=0; iDau!=nDau; ++iDau)
    sTrack->AddDaughter(iID[iDau]);
  sTrack->SetForPOISelection(kTRUE);
  sTrack->SetForRPSelection(kFALSE);
  return sTrack;
}
//_____________________________________________________________________________
void AliAnalysisTaskFlowD2H::Terminate(Option_t *)
{

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