-// $Id$
-
-/**************************************************************************
- * 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. *
- **************************************************************************/
-
-/////////////////////////////////////////////////////////////
-//
-// AliAnalysisTaskSE for Dmesons - jet correlations analysis
-//
-// Author: Xiaoming Zhang, xmzhang@lbl.gov
-///////////////////////////////////////////////////////////////
-
-#include <iostream>
-
-#include <TMath.h>
-#include <TVector3.h>
-#include <TLorentzVector.h>
-#include <TH1D.h>
-#include <TH2D.h>
-#include <TList.h>
-#include <TClonesArray.h>
-#include <TDatabasePDG.h>
-
-#include "AliVTrack.h"
-#include "AliVCluster.h"
-#include "AliEmcalJet.h"
-#include "AliAnalysisTaskEmcal.h"
-#include "AliAnalysisTaskEmcalJet.h"
-#include "AliAODRecoDecayHF2Prong.h"
-#include "AliAODRecoCascadeHF.h"
-#include "AliAnalysisTaskFlavourJetCorrelations.h"
-
-ClassImp(AliAnalysisTaskFlavourJetCorrelations)
-
-//_____________________________________________________________________________
-AliAnalysisTaskFlavourJetCorrelations::AliAnalysisTaskFlavourJetCorrelations() :
-AliAnalysisTaskEmcalJet(),
-fUsedDzeros(0),
-fUsedD0bars(0),
-fUsedDstars(0),
-fListControlHistos(0),
-fListAnDzeroHistos(0),
-fListAnDstarHistos(0)
-{
-//
-// Default constructor
-//
-}
-
-//_____________________________________________________________________________
-AliAnalysisTaskFlavourJetCorrelations::AliAnalysisTaskFlavourJetCorrelations(const char *name, Bool_t bIsHisto) :
-AliAnalysisTaskEmcalJet(name, bIsHisto),
-fUsedDzeros(0),
-fUsedD0bars(0),
-fUsedDstars(0),
-fListControlHistos(0),
-fListAnDzeroHistos(0),
-fListAnDstarHistos(0)
-{
-//
-// Constructor
-//
-
- SetMakeGeneralHistograms(bIsHisto);
- DefineOutput(2, TList::Class());
- DefineOutput(3, TList::Class());
- DefineOutput(4, TList::Class());
-}
-
-//_____________________________________________________________________________
-AliAnalysisTaskFlavourJetCorrelations::~AliAnalysisTaskFlavourJetCorrelations()
-{
-//
-// Default destructor
-//
-
- if (fUsedDzeros) { delete fUsedDzeros; fUsedDzeros=0; }
- if (fUsedD0bars) { delete fUsedD0bars; fUsedD0bars=0; }
- if (fUsedDstars) { delete fUsedDstars; fUsedDstars=0; }
-
- if (fListControlHistos) { delete fListControlHistos; fListControlHistos=0; }
- if (fListAnDzeroHistos) { delete fListAnDzeroHistos; fListAnDzeroHistos=0; }
- if (fListAnDstarHistos) { delete fListAnDstarHistos; fListAnDstarHistos=0; }
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::UserCreateOutputObjects()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::UserCreateOutputObjects
-//
-
-
- AliAnalysisTaskEmcal::UserCreateOutputObjects();
-
- if (!fListControlHistos) fListControlHistos = new TList(); fListControlHistos->SetOwner();
- if (!fListAnDzeroHistos) fListAnDzeroHistos = new TList(); fListAnDzeroHistos->SetOwner();
- if (!fListAnDstarHistos) fListAnDstarHistos = new TList(); fListAnDstarHistos->SetOwner();
-
- MakeControlHistograms();
- CreateDzeroHistograms();
- CreateDstarHistograms();
-
- PostData(2, fListControlHistos);
- PostData(3, fListAnDzeroHistos);
- PostData(4, fListAnDstarHistos);
- return;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::Run()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::Run
-//
-
- const Int_t nJetPtBin = 7;
- const Double_t dJetPtBin[] = { 0., 5., 10., 20., 30., 50., 80., 120. };
-//=============================================================================
-
- AliEmcalJet *pJet = 0;
- for (Int_t iJet=0; iJet<fJets->GetEntriesFast(); iJet++) {
- pJet = static_cast<AliEmcalJet*>(fJets->At(iJet)); if (!pJet) continue;
- if (!AcceptJet(pJet)) { pJet=0; continue; }
-
- Int_t iJetPtBin = TMath::BinarySearch(nJetPtBin,dJetPtBin,pJet->Pt());
- if (iJetPtBin<0 || iJetPtBin>=nJetPtBin) continue;
-
- if (fUsedDzeros) RunDzeroJet(pJet, iJetPtBin, kTRUE);
- if (fUsedD0bars) RunDzeroJet(pJet, iJetPtBin, kFALSE);
- if (fUsedDstars) RunDstarJet(pJet, iJetPtBin);
- pJet = 0;
- }
-
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::RunDzeroJet(AliEmcalJet const *pJet, const Int_t iJetPtBin, const Bool_t bIsD0)
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::RunDzeroJet
-//
-
- Int_t nCands = 0;
- if (bIsD0) nCands = fUsedDzeros->GetEntriesFast();
- else nCands = fUsedD0bars->GetEntriesFast(); if (nCands==0) return;
-//=============================================================================
-
- const Int_t nms = kDzeroMatchType;
- const TString sMatched[nms] = { "MatchConeCandi", "MatchAreaCandi", "MatchConeProng", "MatchAreaProng" };
-//=============================================================================
-
- Double_t dJetPt = pJet->Pt();
- TVector3 vJet, vRec, vec0, vec1;
- vJet.SetPtEtaPhi(dJetPt, pJet->Eta(), pJet->Phi());
-
- Double_t dMag1 = vJet.Mag();
- Double_t dMag2 = vJet.Mag2();
- Double_t dJetArea = pJet->Area();
-//=============================================================================
-
- Bool_t bIsMatch[nms], bIsUsed[nms];
- for (Int_t i=0; i<nms; i++) { bIsMatch[i]=kFALSE; bIsUsed[i]=kFALSE; }
-
- AliAODRecoDecayHF2Prong *candi = 0;
- for (Int_t iCand=0; iCand<nCands; iCand++) {
- if (bIsD0)
- candi = dynamic_cast<AliAODRecoDecayHF2Prong*>(fUsedDzeros->At(iCand));
- else
- candi = dynamic_cast<AliAODRecoDecayHF2Prong*>(fUsedD0bars->At(iCand));
- if (!candi) continue;
- if (candi->Pt()>dJetPt) continue;
-//=============================================================================
-
- Double_t dMassDzero = candi->InvMassD0();
- Double_t dMassD0bar = candi->InvMassD0bar();
- vRec.SetXYZ(candi->Px(), candi->Pt(), candi->Pz());
- vec0.SetXYZ(candi->PxProng(0), candi->PyProng(0), candi->PzProng(0));
- vec1.SetXYZ(candi->PxProng(1), candi->PyProng(1), candi->PzProng(1));
-
- Double_t delRecoR = vRec.DeltaR(vJet);
- Double_t delTrk0R = vec0.DeltaR(vJet);
- Double_t delTrk1R = vec1.DeltaR(vJet);
-
- bIsMatch[kMatchConeCandi] = (delRecoR<fJetRadius);
- bIsMatch[kMatchAreaCandi] = (delRecoR<dJetArea);
- bIsMatch[kMatchConeProng] = ((delTrk0R<fJetRadius) && (delTrk1R<fJetRadius));
- bIsMatch[kMatchAreaProng] = ((delTrk0R<dJetArea) && (delTrk1R<dJetArea));
-
- Bool_t bNoMatch = kTRUE;
- for (Int_t i=0; i<nms; i++) if (bIsMatch[i]) { bIsUsed[i]=kTRUE; bNoMatch=kFALSE; } if (bNoMatch) continue;
-//=============================================================================
-
- Double_t dRecoFFs = vRec.Dot(vJet) / dMag2;
- Double_t dRecoRel = vRec.Cross(vJet).Mag() / dMag1;
- Double_t dRecoRho = 0.5 / TMath::Pi() / ((delRecoR<1e-12) ? 1e-12 : delRecoR);
-
- Double_t dTrk0FFs = vec0.Dot(vJet) / dMag2;
- Double_t dTrk0Rel = vec0.Cross(vJet).Mag() / dMag1;
- Double_t dTrk0Rho = 0.5 / TMath::Pi() / ((delTrk0R<1e-12) ? 1e-12 : delTrk0R);
-
- Double_t dTrk1FFs = vec1.Dot(vJet) / dMag2;
- Double_t dTrk1Rel = vec1.Cross(vJet).Mag() / dMag1;
- Double_t dTrk1Rho = 0.5 / TMath::Pi() / ((delTrk1R<1e-12) ? 1e-12 : delTrk1R);
-
- for (Int_t i=0; i<nms; i++) if (bIsMatch[i]) {
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_InvM_JetPt_%s",sMatched[i].Data())))->Fill(dMassDzero,dJetPt);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Reco_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoFFs,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Reco_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoRel,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Reco_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delRecoR,dMassDzero,dRecoRho);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0FFs,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0Rel,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk0R,dMassDzero,dTrk0Rho);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1FFs,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1Rel,dMassDzero);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk1R,dMassDzero,dTrk1Rho);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hDzero_2ProngsRho_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk0R,delTrk1R);
-//=============================================================================
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_InvM_JetPt_%s",sMatched[i].Data())))->Fill(dMassD0bar,dJetPt);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Reco_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoFFs,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Reco_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoRel,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Reco_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delRecoR,dMassD0bar,dRecoRho);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0FFs,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0Rel,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk0R,dMassD0bar,dTrk0Rho);
-
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1FFs,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1Rel,dMassD0bar);
- ((TH2D*)fListAnDzeroHistos->FindObject(Form("hD0bar_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk1R,dMassD0bar,dTrk1Rho);
- }
-
- candi = 0;
- }
-
- for (Int_t i=0; i<nms; i++) if (bIsUsed[i]) ((TH1D*)fListAnDzeroHistos->FindObject(Form("hDzero_UsedJetPt_%s",sMatched[i].Data())))->Fill(dJetPt);
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::RunDstarJet(AliEmcalJet const *pJet, const Int_t iJetPtBin)
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::RunDstarJet
-//
-
- Int_t nCands = fUsedDstars->GetEntriesFast(); if (nCands==0) return;
-//=============================================================================
-
- const Int_t nms = kDzeroMatchType;
- const TString sMatched[nms] = { "MatchConeCandi", "MatchAreaCandi", "MatchConeProng", "MatchAreaProng" };
-//=============================================================================
-
- Double_t dJetPt = pJet->Pt();
- TVector3 vJet, vRec, vec0, vec1;
- vJet.SetPtEtaPhi(dJetPt, pJet->Eta(), pJet->Phi());
-
- Double_t dMag1 = vJet.Mag();
- Double_t dMag2 = vJet.Mag2();
- Double_t dJetArea = pJet->Area();
-//=============================================================================
-
- Bool_t bIsMatch[nms], bIsUsed[nms];
- for (Int_t i=0; i<nms; i++) { bIsMatch[i]=kFALSE; bIsUsed[i]=kFALSE; }
-
-//AliAODTrack *bache = 0
- AliAODRecoCascadeHF *candi = 0;
- AliAODRecoDecayHF2Prong *dZero = 0;
- for (Int_t iCand=0; iCand<nCands; iCand++) {
- candi = dynamic_cast<AliAODRecoCascadeHF*>(fUsedDstars->At(iCand)); if (!candi) continue;
- dZero = candi->Get2Prong(); if (!dZero) continue;
-// bache = candi->GetBachelor(); if (!bache) continue;
- if (candi->Pt()>dJetPt) continue;
-//=============================================================================
-
- Double_t dInvM = candi->DeltaInvMass();
- vRec.SetXYZ(candi->Px(), candi->Pt(), candi->Pz());
- vec0.SetXYZ(dZero->PxProng(0), dZero->PyProng(0), dZero->PzProng(0));
- vec1.SetXYZ(dZero->PxProng(1), dZero->PyProng(1), dZero->PzProng(1));
-
- Double_t delRecoR = vRec.DeltaR(vJet);
- Double_t delTrk0R = vec0.DeltaR(vJet);
- Double_t delTrk1R = vec1.DeltaR(vJet);
-
- bIsMatch[kMatchConeCandi] = (delRecoR<fJetRadius);
- bIsMatch[kMatchAreaCandi] = (delRecoR<dJetArea);
- bIsMatch[kMatchConeProng] = ((delTrk0R<fJetRadius) && (delTrk1R<fJetRadius));
- bIsMatch[kMatchAreaProng] = ((delTrk0R<dJetArea) && (delTrk1R<dJetArea));
-
- Bool_t bNoMatch = kTRUE;
- for (Int_t i=0; i<nms; i++) if (bIsMatch[i]) { bIsUsed[i]=kTRUE; bNoMatch=kFALSE; } if (bNoMatch) continue;
-//=============================================================================
-
- Double_t dRecoFFs = vRec.Dot(vJet) / dMag2;
- Double_t dRecoRel = vRec.Cross(vJet).Mag() / dMag1;
- Double_t dRecoRho = 0.5 / TMath::Pi() / ((delRecoR<1e-12) ? 1e-12 : delRecoR);
-
- Double_t dTrk0FFs = vec0.Dot(vJet) / dMag2;
- Double_t dTrk0Rel = vec0.Cross(vJet).Mag() / dMag1;
- Double_t dTrk0Rho = 0.5 / TMath::Pi() / ((delTrk0R<1e-12) ? 1e-12 : delTrk0R);
-
- Double_t dTrk1FFs = vec1.Dot(vJet) / dMag2;
- Double_t dTrk1Rel = vec1.Cross(vJet).Mag() / dMag1;
- Double_t dTrk1Rho = 0.5 / TMath::Pi() / ((delTrk1R<1e-12) ? 1e-12 : delTrk1R);
-
- for (Int_t i=0; i<nms; i++) if (bIsMatch[i]) {
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_InvM_JetPt_%s",sMatched[i].Data())))->Fill(dInvM,dJetPt);
-
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Reco_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoFFs,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Reco_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dRecoRel,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Reco_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delRecoR,dInvM,dRecoRho);
-
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0FFs,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk0Rel,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk0R,dInvM,dTrk0Rho);
-
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_FFs_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1FFs,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_Rel_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(dTrk1Rel,dInvM);
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_Trks_Rho_InvM_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk1R,dInvM,dTrk1Rho);
-
- ((TH2D*)fListAnDstarHistos->FindObject(Form("hDstar_2ProngsRho_%s_JetPtBin%d",sMatched[i].Data(),iJetPtBin)))->Fill(delTrk0R,delTrk1R);
- }
-
- dZero = 0;
- candi = 0;
- }
-
- for (Int_t i=0; i<nms; i++) if (bIsUsed[i]) ((TH1D*)fListAnDstarHistos->FindObject(Form("hDstar_UsedJetPt_%s",sMatched[i].Data())))->Fill(dJetPt);
-
- return;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::RetrieveEventObjects()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::RetrieveEventObjects
-//
-
- if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects()) return kFALSE;
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::IsEventSelected()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::IsEventSelected
-//
-
- if (!AliAnalysisTaskEmcal::IsEventSelected()) return kFALSE;
-
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::FillHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::FillHistograms
-//
-
- if (!AliAnalysisTaskEmcal::FillHistograms()) return kFALSE;
-
- if (!FillControlHistograms()) return kFALSE;
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::FillGeneralHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::FillGeneralHistograms
-//
-
- if (!AliAnalysisTaskEmcal::FillGeneralHistograms()) return kFALSE;
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::ExecOnce()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::ExecOnce
-//
-
- AliAnalysisTaskEmcalJet::ExecOnce();
-
- if (!fUsedDzeros) fUsedDzeros = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AnaUsedDzero"));
- if (!fUsedD0bars) fUsedD0bars = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AnaUsedD0bar"));
- if (!fUsedDstars) fUsedDstars = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("AnaUsedDstar"));
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::CreateDzeroHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::CreateDzeroHistograms
-//
-
- if (!fListAnDzeroHistos) return;
- Bool_t bStatusTmpH = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
- const Int_t nms = kDzeroMatchType;
- const TString sMatched[nms] = { "MatchConeCandi", "MatchAreaCandi", "MatchConeProng", "MatchAreaProng" };
- const Double_t dMass = TDatabasePDG::Instance()->GetParticle(421)->Mass();
-
- const Int_t nMassBin = 200;
- const Double_t dMassMin = dMass - 0.15;
- const Double_t dMassMax = dMass + 0.15;
-
- const Int_t nJetPtBin = 7;
- const Double_t aJetPtBin[] = { 0., 5., 10., 20., 30., 50., 80., 120. };
-
- const Int_t nFFsBin = 21;
- const Double_t dFFsBin[] = { 0.0030, 0.0040, 0.0052, 0.0069, 0.0092, 0.0121, 0.0160,
- 0.0212, 0.0280, 0.0371, 0.0490, 0.0648, 0.0856, 0.1132,
- 0.1497, 0.1980, 0.2618, 0.3461, 0.4576, 0.6050, 0.8000, 1. };
-
- const Int_t nRelBin = 12;
- const Double_t dRelBin[] = { 0., 0.06, 0.12, 0.2084, 0.2802, 0.3769, 0.5069, 0.6817, 0.9169, 1.2331, 1.6585, 2.2306, 3. };
-
- const Int_t nRhoBin = 20;
- const Double_t dRhoMin = 0.;
- const Double_t dRhoMax = 1.;
-//=============================================================================
-
- TH1D *h1 = 0;
- TH2D *h2 = 0;
- for (Int_t i=0; i<nms; i++) {
- h1 = new TH1D(Form("hDzero_UsedJetPt_%s",sMatched[i].Data()), "", nJetPtBin, aJetPtBin);
- h1->Sumw2(); fListAnDzeroHistos->Add(h1); h1=0;
-
- h2 = new TH2D(Form("hDzero_InvM_JetPt_%s",sMatched[i].Data()), "", nMassBin, dMassMin, dMassMax, nJetPtBin, aJetPtBin);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- h2 = new TH2D(Form("hD0bar_InvM_JetPt_%s",sMatched[i].Data()), "", nMassBin, dMassMin, dMassMax, nJetPtBin, aJetPtBin);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- TString sName;
- for (Int_t j=0; j<nJetPtBin; j++) {
- sName = "hDzero_Reco_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_Reco_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_Reco_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_Trks_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_Trks_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_Trks_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hDzero_2ProngsRho";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nRhoBin, dRhoMin, dRhoMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-//=============================================================================
-
- sName = "hD0bar_Reco_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hD0bar_Reco_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hD0bar_Reco_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hD0bar_Trks_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hD0bar_Trks_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
-
- sName = "hD0bar_Trks_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDzeroHistos->Add(h2); h2=0;
- }
- }
-
- TH1::AddDirectory(bStatusTmpH);
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::CreateDstarHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::CreateDstarHistograms
-//
-
- if (!fListAnDstarHistos) return;
- Bool_t bStatusTmpH = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
- const Int_t nms = kDzeroMatchType;
- const TString sMatched[nms] = { "MatchConeCandi", "MatchAreaCandi", "MatchConeProng", "MatchAreaProng" };
- Double_t dMassDzero = TDatabasePDG::Instance()->GetParticle(421)->Mass();
- Double_t dMassDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
- Double_t dMassDelta = dMassDstar - dMassDzero;
-
- const Int_t nMassBin = 200;
- const Double_t dMassMin = dMassDelta - 0.015;
- const Double_t dMassMax = dMassDelta + 0.015;
-
- const Int_t nJetPtBin = 7;
- const Double_t aJetPtBin[] = { 0., 5., 10., 20., 30., 50., 80., 120. };
-
- const Int_t nFFsBin = 21;
- const Double_t dFFsBin[] = { 0.0030, 0.0040, 0.0052, 0.0069, 0.0092, 0.0121, 0.0160,
- 0.0212, 0.0280, 0.0371, 0.0490, 0.0648, 0.0856, 0.1132,
- 0.1497, 0.1980, 0.2618, 0.3461, 0.4576, 0.6050, 0.8000, 1. };
-
- const Int_t nRelBin = 12;
- const Double_t dRelBin[] = { 0., 0.06, 0.12, 0.2084, 0.2802, 0.3769, 0.5069, 0.6817, 0.9169, 1.2331, 1.6585, 2.2306, 3. };
-
- const Int_t nRhoBin = 20;
- const Double_t dRhoMin = 0.;
- const Double_t dRhoMax = 1.;
-//=============================================================================
-
- TH1D *h1 = 0;
- TH2D *h2 = 0;
- for (Int_t i=0; i<nms; i++) {
- h1 = new TH1D(Form("hDstar_UsedJetPt_%s",sMatched[i].Data()), "", nJetPtBin, aJetPtBin);
- h1->Sumw2(); fListAnDstarHistos->Add(h1); h1=0;
-
- h2 = new TH2D(Form("hDstar_InvM_JetPt_%s",sMatched[i].Data()), "", nMassBin, dMassMin, dMassMax, nJetPtBin, aJetPtBin);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- TString sName;
- for (Int_t j=0; j<nJetPtBin; j++) {
- sName = "hDstar_Reco_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_Reco_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_Reco_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_Trks_FFs_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nFFsBin, dFFsBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_Trks_Rel_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRelBin, dRelBin, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_Trks_Rho_InvM";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nMassBin, dMassMin, dMassMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
-
- sName = "hDstar_2ProngsRho";
- h2 = new TH2D(Form("%s_%s_JetPtBin%d",sName.Data(),sMatched[i].Data(),j), "", nRhoBin, dRhoMin, dRhoMax, nRhoBin, dRhoMin, dRhoMax);
- h2->Sumw2(); fListAnDstarHistos->Add(h2); h2=0;
- }
- }
-
- TH1::AddDirectory(bStatusTmpH);
- return;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskFlavourJetCorrelations::FillControlHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::FillControlHistograms
-//
-
- if (!fListControlHistos) return kFALSE;
-
- if (fTracks) {
- AliVTrack *pTrk = 0;
- for (Int_t itrk=0; itrk<fTracks->GetEntriesFast(); itrk++) {
- pTrk = static_cast<AliVTrack*>(fTracks->At(itrk)); if (!pTrk) continue;
- if (!AcceptTrack(pTrk)) { pTrk = 0; continue; }
-
- ((TH1D*)fListControlHistos->FindObject("hTrksPt"))->Fill(pTrk->Pt()); pTrk=0;
- }
- }
-
- if (fCaloClusters) {
- AliVCluster *pCls = 0;
- TLorentzVector vLorentz;
- for (Int_t iClus=0; iClus<fCaloClusters->GetEntriesFast(); iClus++) {
- pCls = static_cast<AliVCluster*>(fCaloClusters->At(iClus)); if (!pCls) continue;
-
- pCls->GetMomentum(vLorentz, fVertex);
- ((TH1D*)fListControlHistos->FindObject("hClusPt"))->Fill(vLorentz.Pt()); pCls=0;
- }
- }
-
- if (fJets) {
- AliEmcalJet* pJet = 0;
- static Int_t sortedJets[9999] = { -1 };
- if (GetSortedArray(sortedJets,fJets)) if (sortedJets[0]>=0) {
- pJet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[0]));
- if (pJet) ((TH1D*)fListControlHistos->FindObject("hLeadingJets_Pt"))->Fill(pJet->Pt());
- } pJet = 0;
-
- for (Int_t iJet=0; iJet<fJets->GetEntriesFast(); iJet++) {
- pJet = static_cast<AliEmcalJet*>(fJets->At(iJet)); if (!pJet) continue;
- if (!AcceptJet(pJet)) continue;
-
- Double_t dPt = pJet->Pt();
- Double_t dArea = pJet->Area();
- ((TH2D*)fListControlHistos->FindObject("hJets_Pt_Area"))->Fill(dPt,dArea);
- ((TH2D*)fListControlHistos->FindObject("hJets_Eta_Phi"))->Fill(pJet->Eta(),pJet->Phi());
- ((TH2D*)fListControlHistos->FindObject("hJets_LeadingPt"))->Fill(dPt,GetLeadingHadronPt(pJet));
-
- if (fRho) {
- ((TH2D*)fListControlHistos->FindObject("hJets_CorrPt_Area"))->Fill(dPt-fRhoVal*dArea,dArea);
- } pJet = 0;
- }
- }
-
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskFlavourJetCorrelations::MakeControlHistograms()
-{
-//
-// AliAnalysisTaskFlavourJetCorrelations::MakeControlHistograms
-//
-
- if (!fListControlHistos) return;
- Bool_t bStatusTmpH = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
- TH1D *h1 = 0;
- TH2D *h2 = 0;
-
- if (!fTracksName.IsNull()) {
- h1 = new TH1D("hTrksPt", "", fNbins/2, fMinBinPt, fMaxBinPt/2.);
- h1->Sumw2(); fListControlHistos->Add(h1); h1=0;
- }
-
- if (!fCaloName.IsNull()) {
- h1 = new TH1D("hClusPt", "", fNbins/2, fMinBinPt, fMaxBinPt/2.);
- h1->Sumw2(); fListControlHistos->Add(h1); h1=0;
- }
-
- if (!fJetsName.IsNull()) {
- h1 = new TH1D("hLeadingJets_Pt", "", fNbins, fMinBinPt, fMaxBinPt);
- h1->Sumw2(); fListControlHistos->Add(h1); h1=0;
-
- h2 = new TH2D("hJets_Pt_Area", "", fNbins, fMinBinPt, fMaxBinPt, 30, 0., 3.);
- h2->Sumw2(); fListControlHistos->Add(h2); h2=0;
-
- h2 = new TH2D("hJets_Eta_Phi", "", 50, -1., 1., 101, 0., 2.*TMath::Pi() + TMath::Pi()/200.);
- h2->Sumw2(); fListControlHistos->Add(h2); h2=0;
-
- h2 = new TH2D("hJets_LeadingPt", "", fNbins, fMinBinPt, fMaxBinPt, fNbins/2, fMinBinPt, fMaxBinPt/2.);
- h2->Sumw2(); fListControlHistos->Add(h2); h2=0;
-
- if (!fRhoName.IsNull()) {
- h2 = new TH2D("hJets_CorrPt_Area", "", 2*fNbins, -1.*fMaxBinPt, fMaxBinPt, 30, 0., 3.);
- }
- }
-
- TH1::AddDirectory(bStatusTmpH);
- return;
-}
+// $Id$\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: The ALICE Off-line Project. *\r
+ * Contributors are mentioned in the code where appropriate. *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+//\r
+//\r
+// Analysis Taks for reconstructed particle correlation \r
+// (first implementation done for D mesons) with jets \r
+// (use the so called Emcal framework)\r
+//\r
+//-----------------------------------------------------------------------\r
+// Authors:\r
+// C. Bianchin (Utrecht University) chiara.bianchin@cern.ch\r
+// A. Grelli (Utrecht University) a.grelli@uu.nl\r
+// X. Zhang (LBNL) XMZhang@lbl.gov\r
+//-----------------------------------------------------------------------\r
+\r
+#include <TDatabasePDG.h>\r
+#include <TParticle.h>\r
+#include <TVector3.h>\r
+#include "TROOT.h"\r
+#include <TH3F.h>\r
+\r
+#include "AliAnalysisTaskFlavourJetCorrelations.h"\r
+#include "AliAODMCHeader.h"\r
+#include "AliAODHandler.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliLog.h"\r
+#include "AliEmcalJet.h"\r
+#include "AliAODRecoDecay.h"\r
+#include "AliAODRecoCascadeHF.h"\r
+#include "AliAODRecoDecayHF2Prong.h"\r
+#include "AliESDtrack.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliPicoTrack.h"\r
+#include "AliRDHFCutsD0toKpi.h"\r
+#include "AliRDHFCutsDStartoKpipi.h"\r
+\r
+ClassImp(AliAnalysisTaskFlavourJetCorrelations)\r
+\r
+//__________________________________________________________________________\r
+AliAnalysisTaskFlavourJetCorrelations::AliAnalysisTaskFlavourJetCorrelations() :\r
+ AliAnalysisTaskEmcalJet(),\r
+ fUseMCInfo(kTRUE), \r
+ fCandidateType(),\r
+ fPDGmother(),\r
+ fNProngs(),\r
+ fPDGdaughters(),\r
+ fBranchName(),\r
+ fOutput(0),\r
+ fCuts(0),\r
+ fMinMass(),\r
+ fMaxMass(), \r
+ fJetArrName(0),\r
+ fCandArrName(0),\r
+ fLeadingJetOnly(kFALSE)\r
+{\r
+ //\r
+ // Default ctor\r
+ //\r
+}\r
+//___________________________________________________________________________\r
+AliAnalysisTaskFlavourJetCorrelations::AliAnalysisTaskFlavourJetCorrelations(const Char_t* name, AliRDHFCuts* cuts,ECandidateType candtype, TString jetArrName) :\r
+ AliAnalysisTaskEmcalJet(name),\r
+ fUseMCInfo(kTRUE),\r
+ fCandidateType(),\r
+ fPDGmother(),\r
+ fNProngs(),\r
+ fPDGdaughters(),\r
+ fBranchName(),\r
+ fOutput(0),\r
+ fCuts(0),\r
+ fMinMass(),\r
+ fMaxMass(), \r
+ fJetArrName(0),\r
+ fCandArrName(0),\r
+ fLeadingJetOnly(kFALSE)\r
+{\r
+ //\r
+ // Constructor. Initialization of Inputs and Outputs\r
+ //\r
+ \r
+ Info("AliAnalysisTaskFlavourJetCorrelations","Calling Constructor");\r
+ fCuts=cuts;\r
+ fCandidateType=candtype;\r
+ const Int_t nptbins=fCuts->GetNPtBins();\r
+ Float_t defaultSigmaD013[13]={0.012, 0.012, 0.012, 0.015, 0.015,0.018,0.018,0.020,0.020,0.030,0.030,0.037,0.040};\r
+\r
+ switch(fCandidateType){\r
+ case 0:\r
+ fPDGmother=421;\r
+ fNProngs=2;\r
+ fPDGdaughters[0]=211;//pi \r
+ fPDGdaughters[1]=321;//K\r
+ fPDGdaughters[2]=0; //empty\r
+ fPDGdaughters[3]=0; //empty\r
+ fBranchName="D0toKpi";\r
+ fCandArrName="D0";\r
+ break;\r
+ case 1: \r
+ fPDGmother=413;\r
+ fNProngs=3;\r
+ fPDGdaughters[1]=211;//pi soft\r
+ fPDGdaughters[0]=421;//D0\r
+ fPDGdaughters[2]=211;//pi fromD0\r
+ fPDGdaughters[3]=321; // K from D0\r
+ fBranchName="Dstar";\r
+ fCandArrName="Dstar";\r
+ //fSigmaD0=new Float_t[nptbins];\r
+ if(nptbins<=13){\r
+ for(Int_t ipt=0;ipt<nptbins;ipt++) fSigmaD0[ipt]= defaultSigmaD013[ipt];\r
+ } else {\r
+ AliFatal(Form("Default sigma D0 not enough for %d pt bins, use SetSigmaD0ForDStar to set them",nptbins));\r
+ }\r
+ break;\r
+ default:\r
+ printf("%d not accepted!!\n",fCandidateType);\r
+ break;\r
+ }\r
+\r
+ if(fCandidateType==kD0toKpi)SetMassLimits(0.15,fPDGmother);\r
+ if(fCandidateType==kDstartoKpipi) SetMassLimits(0.015, fPDGmother);\r
+ fJetArrName=jetArrName;\r
+ Printf("Jet read is %s",fJetArrName.Data());\r
+\r
+\r
+\r
+ DefineOutput(1,TList::Class()); // histos\r
+ DefineOutput(2,AliRDHFCuts::Class()); // my cuts\r
+}\r
+//___________________________________________________________________________\r
+AliAnalysisTaskFlavourJetCorrelations::~AliAnalysisTaskFlavourJetCorrelations() {\r
+ //\r
+ // destructor\r
+ //\r
+\r
+ Info("~AliAnalysisTaskFlavourJetCorrelations","Calling Destructor"); \r
+ \r
+ delete fOutput;\r
+ delete fCuts;\r
+ \r
+}\r
+\r
+//___________________________________________________________\r
+void AliAnalysisTaskFlavourJetCorrelations::Init(){\r
+ //\r
+ // Initialization\r
+ //\r
+ if(fDebug > 1) printf("AnalysisTaskRecoJetCorrelations::Init() \n");\r
+ switch(fCandidateType){\r
+ case 0:\r
+ {\r
+ AliRDHFCutsD0toKpi* copyfCuts=new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));\r
+ copyfCuts->SetName("AnalysisCutsDzero");\r
+ // Post the data\r
+ PostData(2,copyfCuts);\r
+ }\r
+ break;\r
+ case 1:\r
+ {\r
+ AliRDHFCutsDStartoKpipi* copyfCuts=new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));\r
+ copyfCuts->SetName("AnalysisCutsDStar");\r
+ // Post the cuts\r
+ PostData(2,copyfCuts);\r
+ }\r
+ break;\r
+ default:\r
+ return;\r
+ }\r
+\r
+ return;\r
+}\r
+\r
+//_________________________________________________\r
+void AliAnalysisTaskFlavourJetCorrelations::UserExec(Option_t *)\r
+{\r
+ // user exec\r
+\r
+ // Load the event\r
+ AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);\r
+ \r
+ TClonesArray *arrayDStartoD0pi=0;\r
+ TClonesArray *trackArr = 0;\r
+ TClonesArray *clusArr = 0;\r
+ TClonesArray *jetArr = 0;\r
+ TClonesArray *candidatesArr = 0;\r
+//TClonesArray *isselArr = 0;\r
+\r
+ if (!aodEvent && AODEvent() && IsStandardAOD()) {\r
+\r
+ // In case there is an AOD handler writing a standard AOD, use the AOD \r
+ // event in memory rather than the input (ESD) event. \r
+ aodEvent = dynamic_cast<AliAODEvent*> (AODEvent());\r
+\r
+ // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)\r
+ // have to taken from the AOD event hold by the AliAODExtension\r
+ AliAODHandler* aodHandler = (AliAODHandler*) \r
+ ((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());\r
+ if(aodHandler->GetExtensions()) {\r
+ AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");\r
+ AliAODEvent *aodFromExt = ext->GetAOD();\r
+ arrayDStartoD0pi=(TClonesArray*)aodFromExt->GetList()->FindObject(fBranchName.Data());\r
+ }\r
+ } else {\r
+ arrayDStartoD0pi=(TClonesArray*)aodEvent->GetList()->FindObject(fBranchName.Data());\r
+ }\r
+\r
+ if (!arrayDStartoD0pi) {\r
+ AliInfo(Form("Could not find array %s, skipping the event",fBranchName.Data()));\r
+// return;\r
+ } else AliDebug(2, Form("Found %d vertices",arrayDStartoD0pi->GetEntriesFast())); \r
+\r
+ TClonesArray* mcArray = 0x0;\r
+ if (fUseMCInfo) {\r
+ mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));\r
+ if (!mcArray) {\r
+ printf("AliAnalysisTaskSEDStarSpectra::UserExec: MC particles not found!\n");\r
+ return;\r
+ }\r
+ }\r
+\r
+ //retrieve jets\r
+ trackArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("PicoTracks"));\r
+ clusArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject("CaloClustersCorr"));\r
+ jetArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fJetArrName));\r
+\r
+ if(!trackArr){\r
+ AliInfo(Form("Could not find the track array\n"));\r
+ return;\r
+ }\r
+\r
+ if(!jetArr){\r
+ Printf("JET array not found");\r
+ return;\r
+ }\r
+\r
+ //retrieve reconstructed particles selected\r
+ /*\r
+ TString listname="AliAnalysisTaskSEDmesonsForJetCorrelations";\r
+ TList *l=dynamic_cast<TList*>(InputEvent()->FindListObject(listname));\r
+ TClonesArray *cla=dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(listname));\r
+ // if(l){\r
+ // l->ls();\r
+ // } else{\r
+ // Printf("list not found!!!!!!!!!!!");\r
+ // return;\r
+ // } \r
+ if(!cla){\r
+ Printf("cla not found!!!!!!!!!!!");\r
+ return;\r
+ } else {\r
+ cla->ls();\r
+ }\r
+ */\r
+\r
+ TString arrname="fCandidateArray";\r
+ candidatesArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s%s",arrname.Data(),fCandArrName.Data())));\r
+ if (!candidatesArr) {\r
+ Printf("%s%s array not found",arrname.Data(),fCandArrName.Data());\r
+ InputEvent()->GetList()->ls();\r
+ return;\r
+ }\r
+\r
+ /*\r
+ arrname="fIsSelectedArray";\r
+ isselArr = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(Form("%s%s",arrname.Data(),fCandArrName.Data())));\r
+ if(!isselArr){\r
+ Printf("%s%s array not found",arrname.Data(),fCandArrName.Data());\r
+ InputEvent()->ls();\r
+ return;\r
+ }\r
+ */\r
+\r
+ //Histograms\r
+ TH1I* hstat=(TH1I*)fOutput->FindObject("hstat");\r
+ TH1F* hPtJetTrks=(TH1F*)fOutput->FindObject("hPtJetTrks");\r
+ TH1F* hPhiJetTrks=(TH1F*)fOutput->FindObject("hPhiJetTrks");\r
+ TH1F* hEtaJetTrks=(TH1F*)fOutput->FindObject("hEtaJetTrks");\r
+ TH1F* hEjetTrks=(TH1F*)fOutput->FindObject("hEjetTrks");\r
+ TH1F* hPtJet=(TH1F*)fOutput->FindObject("hPtJet");\r
+ TH1F* hPhiJet=(TH1F*)fOutput->FindObject("hPhiJet");\r
+ TH1F* hEtaJet=(TH1F*)fOutput->FindObject("hEtaJet");\r
+ TH1F* hEjet=(TH1F*)fOutput->FindObject("hEjet");\r
+ TH1F* hNtrArr=(TH1F*)fOutput->FindObject("hNtrArr");\r
+ TH1F* hNJetPerEv=(TH1F*)fOutput->FindObject("hNJetPerEv");\r
+ TH1F* hdeltaRJetTracks=((TH1F*)fOutput->FindObject("hdeltaRJetTracks"));\r
+\r
+ hstat->Fill(0);\r
+\r
+ // fix for temporary bug in ESDfilter \r
+ // the AODs with null vertex pointer didn't pass the PhysSel\r
+ if(!aodEvent->GetPrimaryVertex() || TMath::Abs(aodEvent->GetMagneticField())<0.001) return;\r
+\r
+ //Event selection\r
+ Bool_t iseventselected=fCuts->IsEventSelected(aodEvent);\r
+ TString firedTriggerClasses=((AliAODEvent*)aodEvent)->GetFiredTriggerClasses();\r
+ if(!iseventselected) return;\r
+ \r
+ hstat->Fill(1);\r
+ \r
+ //trigger on jets\r
+ Int_t njets=jetArr->GetEntriesFast();\r
+ if(njets == 0) return;\r
+\r
+ const Int_t nD=arrayDStartoD0pi->GetEntriesFast();\r
+ hstat->Fill(2,nD);\r
+ \r
+ // counters for efficiencies\r
+ Int_t icountReco = 0;\r
+ \r
+ //D* and D0 prongs needed to MatchToMC method\r
+ // Int_t pdgDgDStartoD0pi[2]={421,211};\r
+ // Int_t pdgDgD0toKpi[2]={321,211};\r
+ //D0 from D0 bar\r
+ \r
+ // we start with jets\r
+ Double_t ejet = 0;\r
+ Double_t phiJet = 0;\r
+ Double_t etaJet = 0;\r
+ Double_t ptjet = 0;\r
+ Double_t leadingJet =0;\r
+\r
+ Int_t ntrarr=trackArr->GetEntriesFast();\r
+ hNtrArr->Fill(ntrarr);\r
+ \r
+ for(Int_t i=0;i<ntrarr;i++){\r
+ AliVTrack *jtrack=static_cast<AliVTrack*>(trackArr->At(i));\r
+ //reusing histograms\r
+ hPtJetTrks->Fill(jtrack->Pt());\r
+ hPhiJetTrks->Fill(jtrack->Phi());\r
+ hEtaJetTrks->Fill(jtrack->Eta());\r
+ hEjetTrks->Fill(jtrack->E());\r
+ }\r
+\r
+ \r
+ //option to use only the leading jet\r
+ if(fLeadingJetOnly){\r
+ for (Int_t iJetsL = 0; iJetsL<njets; iJetsL++) { \r
+ AliEmcalJet* jetL = (AliEmcalJet*)jetArr->At(iJetsL);\r
+ ptjet = jetL->Pt();\r
+ if(ptjet>leadingJet ) leadingJet = ptjet;\r
+\r
+ }\r
+ }\r
+\r
+ Int_t cntjet=0;\r
+ //loop over jets and consider only the leading jet in the event\r
+ for (Int_t iJets = 0; iJets<njets; iJets++) { \r
+ AliEmcalJet* jet = (AliEmcalJet*)jetArr->At(iJets);\r
+ if(!AcceptJet(jet)) continue;\r
+\r
+ vector<int> DmesonInJetLabels(10);\r
+ //Int_t iD=0;\r
+ //jets variables\r
+ ejet = jet->E();\r
+ phiJet = jet->Phi();\r
+ etaJet = jet->Eta();\r
+ ptjet = jet->Pt();\r
+ \r
+ // choose the leading jet\r
+ if(fLeadingJetOnly && (ejet<leadingJet)) continue;\r
+ //used for normalization\r
+ hstat->Fill(3);\r
+ cntjet++;\r
+ // fill energy, eta and phi of the jet\r
+ hEjet ->Fill(ejet);\r
+ hPhiJet ->Fill(phiJet);\r
+ hEtaJet ->Fill(etaJet);\r
+ hPtJet ->Fill(ptjet);\r
+ \r
+ //loop on jet particles\r
+ Int_t ntrjet= jet->GetNumberOfTracks(); \r
+ for(Int_t itrk=0;itrk<ntrjet;itrk++){\r
+ \r
+ AliPicoTrack* jetTrk=(AliPicoTrack*)jet->TrackAt(itrk,trackArr); \r
+ hdeltaRJetTracks->Fill(DeltaR(jet,jetTrk));\r
+ \r
+ // check MC in for the traks within the jet, look at D mesons\r
+ // within the jet area\r
+ //if(fUseMCInfo) FillMCDJetInfo(jetTrk,jet,mcArray,ptjet);\r
+ \r
+ }//end loop on jet tracks\r
+ \r
+ //retrieve charm candidates selected\r
+ Int_t candidates = candidatesArr->GetEntriesFast();\r
+ for(Int_t ic = 0; ic < candidates; ic++) {\r
+ // D* candidates\r
+ AliAODRecoDecayHF* charm = 0x0;\r
+ charm=(AliAODRecoDecayHF*)candidatesArr->At(ic);\r
+ \r
+ FillHistogramsRecoJetCorr(charm, jet); \r
+\r
+ }\r
+ } // end of jet loop \r
+\r
+ hNJetPerEv->Fill(cntjet);\r
+\r
+ AliDebug(2, Form("Found %i Reco particles that are D*!!",icountReco));\r
+ \r
+ PostData(1,fOutput);\r
+\r
+}\r
+\r
+//________________________________________ terminate ___________________________\r
+void AliAnalysisTaskFlavourJetCorrelations::Terminate(Option_t*)\r
+{ \r
+ // The Terminate() function is the last function to be called during\r
+ // a query. It always runs on the client, it can be used to present\r
+ // the results graphically or save the results to file.\r
+ \r
+ Info("Terminate"," terminate");\r
+ AliAnalysisTaskSE::Terminate();\r
+\r
+ fOutput = dynamic_cast<TList*> (GetOutputData(1));\r
+ if (!fOutput) { \r
+ printf("ERROR: fOutput not available\n");\r
+ return;\r
+ }\r
+}\r
+//___________________________________________________________________________\r
+\r
+void AliAnalysisTaskFlavourJetCorrelations::UserCreateOutputObjects() { \r
+ // output \r
+ Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());\r
+ \r
+ //slot #1 \r
+ OpenFile(1);\r
+ fOutput = new TList();\r
+ fOutput->SetOwner();\r
+ // define histograms\r
+ DefineHistoForAnalysis();\r
+\r
+ PostData(1,fOutput);\r
+ return;\r
+}\r
+//_________________________________________________________________\r
+void AliAnalysisTaskFlavourJetCorrelations::SetMassLimits(Double_t range, Int_t pdg){\r
+ Float_t mass=0;\r
+ Int_t abspdg=TMath::Abs(pdg);\r
+ \r
+ mass=TDatabasePDG::Instance()->GetParticle(abspdg)->Mass();\r
+ // compute the Delta mass for the D*\r
+ if(fCandidateType==kDstartoKpipi){\r
+ Float_t mass1=0;\r
+ mass1=TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
+ mass = mass-mass1;\r
+ }\r
+ cout<<"mass ---------------"<<endl;\r
+ fMinMass = mass-range;\r
+ fMaxMass = mass+range;\r
+ \r
+ AliInfo(Form("Setting mass limits to %f, %f",fMinMass,fMaxMass));\r
+ if (fMinMass<0 || fMaxMass<=0 || fMaxMass<fMinMass) AliFatal("Wrong mass limits!\n");\r
+}\r
+//_________________________________________________________________\r
+void AliAnalysisTaskFlavourJetCorrelations::SetMassLimits(Double_t lowlimit, Double_t uplimit){\r
+ if(uplimit>lowlimit)\r
+ {\r
+ fMinMass = lowlimit;\r
+ fMaxMass = uplimit;\r
+ }\r
+ else{\r
+ printf("Error! Lower limit larger than upper limit!\n");\r
+ fMinMass = uplimit - uplimit*0.2;\r
+ fMaxMass = uplimit;\r
+ }\r
+}\r
+\r
+//__________________________________________________________________\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::SetD0WidthForDStar(Int_t nptbins,Float_t *width){\r
+ if(nptbins>30) {\r
+ AliInfo("Maximum number of bins allowed is 30!");\r
+ return kFALSE;\r
+ }\r
+ if(!width) return kFALSE;\r
+ for(Int_t ipt=0;ipt<nptbins;ipt++) fSigmaD0[ipt]=width[ipt];\r
+ return kTRUE;\r
+}\r
+\r
+//__________________________________________________________________\r
+Double_t AliAnalysisTaskFlavourJetCorrelations::Z(AliVParticle* part,AliEmcalJet* jet) const{\r
+ if(!part || !jet){\r
+ printf("Particle or jet do not exist!\n");\r
+ return -99;\r
+ }\r
+ Double_t p[3],pj[3];\r
+ Bool_t okpp=part->PxPyPz(p);\r
+ Bool_t okpj=jet->PxPyPz(pj);\r
+ if(!okpp || !okpj){\r
+ printf("Problems getting momenta\n");\r
+ return -999;\r
+ }\r
+ Double_t pjet=jet->P();\r
+ Double_t z=(p[0]*pj[0]+p[1]*pj[1]+p[2]*pj[2])/(pjet*pjet);\r
+ return z;\r
+}\r
+//___________________________________ histograms _______________________________________\r
+\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::DefineHistoForAnalysis(){\r
+ \r
+ // Statistics \r
+ TH1I* hstat=new TH1I("hstat","Statistics",5,-0.5,4.5);\r
+ hstat->GetXaxis()->SetBinLabel(1,"N ev anal");\r
+ hstat->GetXaxis()->SetBinLabel(2,"N ev sel");\r
+ hstat->GetXaxis()->SetBinLabel(3,"N cand sel cuts");\r
+ hstat->GetXaxis()->SetBinLabel(4,"N jets");\r
+ hstat->GetXaxis()->SetBinLabel(5,"N cand in jet");\r
+ // if(fUseMCInfo) {\r
+ // hstat->GetXaxis()->SetBinLabel(7,"N D");\r
+ // hstat->GetXaxis()->SetBinLabel(8,"N D in jet");\r
+\r
+ // }\r
+\r
+ hstat->SetNdivisions(1);\r
+ fOutput->Add(hstat);\r
+\r
+ const Int_t nbinsmass=200;\r
+\r
+ \r
+ if(fCandidateType==kDstartoKpipi) \r
+ {\r
+ // TH2F *hDiff = new TH2F("hDiff","M(kpipi)-M(kpi)",500,fMinMass,fMaxMass,100, 0.,30.);\r
+ // hDiff->SetStats(kTRUE);\r
+ // hDiff->GetXaxis()->SetTitle("M(kpipi)-M(kpi) GeV/c^{2}");\r
+ // hDiff->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ // fOutput->Add(hDiff);\r
+ \r
+ TH2F* hDiffSideBand = new TH2F("hDiffSideBand","M(kpipi)-M(kpi) Side Band Background",nbinsmass,fMinMass,fMaxMass,100, 0.,30.);\r
+ hDiffSideBand->SetStats(kTRUE);\r
+ hDiffSideBand->GetXaxis()->SetTitle("M(kpipi)-M(Kpi) GeV/c^{2}");\r
+ hDiffSideBand->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ fOutput->Add(hDiffSideBand); \r
+ \r
+ //correlation histograms\r
+ // fPhi = new TH1F("phi","Delta phi between Jet axis and DStar ",25,-1.57,4.72);\r
+ // fPhi->SetStats(kTRUE);\r
+ // fPhi->GetXaxis()->SetTitle("#Delta #phi (rad)");\r
+ // fPhi->GetYaxis()->SetTitle("Entries");\r
+ // fOutput->Add(fPhi);\r
+\r
+ // fDphiD0Dstar = new TH1F("phiD0Dstar","Delta phi between D0 and DStar ",1000,-6.5,6.5);\r
+ // fOutput->Add(fDphiD0Dstar);\r
+\r
+ // fPhiBkg = new TH1F("phiBkg","Delta phi between Jet axis and DStar background ",25,-1.57,4.72);\r
+ // fPhiBkg->SetStats(kTRUE);\r
+ // fPhiBkg->GetXaxis()->SetTitle("#Delta #phi (rad)");\r
+ // fPhiBkg->GetYaxis()->SetTitle("Entries");\r
+ // fOutput->Add(fPhiBkg);\r
+\r
+ // TH1F* hRECOPtDStar = new TH1F("hRECODStar","RECO DStar pt distribution",30,0,30);\r
+ // hRECOPtDStar->SetStats(kTRUE);\r
+ // hRECOPtDStar->SetLineColor(2);\r
+ // hRECOPtDStar->GetXaxis()->SetTitle("GeV/c");\r
+ // hRECOPtDStar->GetYaxis()->SetTitle("Entries");\r
+ // fOutput->Add(hRECOPtDStar);\r
+ \r
+ // TH1F* hRECOPtBkg = new TH1F("hRECOptBkg","RECO pt distribution side bands",30,0,30);\r
+ // hRECOPtBkg->SetStats(kTRUE);\r
+ // hRECOPtBkg->SetLineColor(2);\r
+ // hRECOPtBkg->GetXaxis()->SetTitle("GeV/c");\r
+ // hRECOPtBkg->GetYaxis()->SetTitle("Entries");\r
+ // fOutput->Add(hRECOPtBkg);\r
+\r
+ TH1F* hPtPion = new TH1F("hPtPion","Primary pions candidates pt ",500,0,10);\r
+ hPtPion->SetStats(kTRUE);\r
+ hPtPion->GetXaxis()->SetTitle("GeV/c");\r
+ hPtPion->GetYaxis()->SetTitle("Entries");\r
+ fOutput->Add(hPtPion);\r
+\r
+ }\r
+\r
+ // jet related fistograms\r
+ \r
+ TH1F* hEjetTrks = new TH1F("hEjetTrks", "Jet tracks energy distribution;Energy (GeV)",500,0,200);\r
+ TH1F* hPhiJetTrks = new TH1F("hPhiJetTrks","Jet tracks #phi distribution; #phi (rad)", 200,0,6.30);\r
+ TH1F* hEtaJetTrks = new TH1F("hEtaJetTrks","Jet tracks #eta distribution; #eta", 100,-1.5,1.5);\r
+ TH1F* hPtJetTrks = new TH1F("hPtJetTrks", "Jet tracks Pt distribution; p_{T} (GeV/c^{2})",500,0,200);\r
+ \r
+ TH1F* hEjet = new TH1F("hEjet", "Jet energy distribution;Energy (GeV)",500,0,200);\r
+ TH1F* hPhiJet = new TH1F("hPhiJet","Jet #phi distribution; #phi (rad)", 200,0,6.30);\r
+ TH1F* hEtaJet = new TH1F("hEtaJet","Jet #eta distribution; #eta", 100,-1.5,1.5);\r
+ TH1F* hPtJet = new TH1F("hPtJet", "Jet Pt distribution; p_{T} (GeV/c^{2})",500,0,200);\r
+\r
+ TH1F* hPtJetWithD=new TH1F("hPtJetWithD","D-Jet Pt distribution; p_{T} (GeV/c^{2})",500,0,200);\r
+ TH1F* hdeltaRJetTracks=new TH1F("hdeltaRJetTracks","Delta R of tracks in the jets",200, 0.,10.);\r
+\r
+ TH1F* hNtrArr= new TH1F("hNtrArr", "Number of tracks in the array of jets; number of tracks",500,0,1000);\r
+ TH1F *hNJetPerEv=new TH1F("hNJetPerEv","Number of jets used per event; number of jets/ev",10,-0.5,9.5);\r
+\r
+ fOutput->Add(hEjetTrks);\r
+ fOutput->Add(hPhiJetTrks);\r
+ fOutput->Add(hEtaJetTrks);\r
+ fOutput->Add(hPtJetTrks);\r
+ fOutput->Add(hEjet);\r
+ fOutput->Add(hPhiJet);\r
+ fOutput->Add(hEtaJet);\r
+ fOutput->Add(hPtJet);\r
+ fOutput->Add(hPtJetWithD);\r
+ fOutput->Add(hdeltaRJetTracks);\r
+ fOutput->Add(hNtrArr);\r
+ fOutput->Add(hNJetPerEv);\r
+\r
+ /*\r
+ if(fUseMCInfo){\r
+ fhMomjetpartPdg=new TH1F("fhMomjetpartPdg","Jet particles' mothers distribution;PDG code;Counts;",1100,-550,550);\r
+ fOutput->Add(fhMomjetpartPdg);\r
+ fptDinjetallvsptjMC=new TH2F("fptDinjetallvsptjMC","p_{t} distribution of D within a jet (all DeltaR) vs p_{t}^{jet};p_{t}^{D};p_{t}^{jet}",100, 0.,30.,500,0.,200.);\r
+ fOutput->Add(fptDinjetallvsptjMC);\r
+ fptJwithDMC=new TH1F("fptJwithDMC","p_{t}^{jet} with a D meson (all DeltaR);p_{t}^{jet};Counts",500,0.,200.);\r
+ fOutput->Add(fptJwithDMC);\r
+\r
+ fptDinjet04vsptjMC=new TH2F("fptDinjet04vsptjMC","p_{t} distribution of D within a jet (DeltaR 0.4) vs p_{t}^{jet};p_{t}^{D};p_{t}^{jet}",100, 0.,30.,500,0.,200.);\r
+ fOutput->Add(fptDinjet04vsptjMC);\r
+ TH1F* hdeltaRDMC=new TH1F("hdeltaRDMC","#Delta R for MC tagged D mesons; #Delta R_{D}^{MC}",200, 0.,10.);\r
+ fOutput->Add(hdeltaRDMC);\r
+ }\r
+ */\r
+ TH1F* hDeltaRD=new TH1F("hDeltaRD","#Delta R distribution of D candidates selected;#Delta R",200, 0.,10.);\r
+ fOutput->Add(hDeltaRD);\r
+\r
+ TH3F* hdeltaPhiDja=new TH3F("hdeltaPhiDja", "#Delta#phi D-jet (jet p_{T} > threshold)",nbinsmass,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ hdeltaPhiDja->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ hdeltaPhiDja->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c)");\r
+ hdeltaPhiDja->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ // TH3F* hdeltaPhiDjl=new TH3F("hdeltaPhiDjl", Form("#Delta#phi D-jet (jet p_{T} < %.0f (GeV/c^{2}))",fJetPtThr),500,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ // hdeltaPhiDjl->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ // hdeltaPhiDjl->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ // hdeltaPhiDjl->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ // TH3F* hdeltaPhiDjh=new TH3F("hdeltaPhiDjh", Form("#Delta#phi D-jet (jet p_{T} > %.0f (GeV/c^{2}))",fJetPtThr),500,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ // hdeltaPhiDjh->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ // hdeltaPhiDjh->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ // hdeltaPhiDjh->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ fOutput->Add(hdeltaPhiDja);\r
+ // fOutput->Add(hdeltaPhiDjl);\r
+ // fOutput->Add(hdeltaPhiDjh);\r
+\r
+ //background (side bands for the Dstar and like sign for D0)\r
+\r
+ TH3F* hdeltaPhiDjaB=new TH3F("hdeltaPhiDjaB", "#Delta#phi D-jet (all jet p_{T})",nbinsmass,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ hdeltaPhiDjaB->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ hdeltaPhiDjaB->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c)");\r
+ hdeltaPhiDjaB->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ // TH3F* hdeltaPhiDjlB=new TH3F("hdeltaPhiDjlB", Form("#Delta#phi D-jet (jet p_{T} < %.0f (GeV/c^{2}))",fJetPtThr),1500,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ // hdeltaPhiDjlB->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ // hdeltaPhiDjlB->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ // hdeltaPhiDjlB->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ // TH3F* hdeltaPhiDjhB=new TH3F("hdeltaPhiDjhB", Form("#Delta#phi D-jet (jet p_{T} > %.0f (GeV/c^{2}))",fJetPtThr),1500,fMinMass,fMaxMass,100, 0.,30.,50,(-1)*TMath::Pi()/2.,3./2.*TMath::Pi());\r
+ // hdeltaPhiDjhB->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ // hdeltaPhiDjhB->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c^{2})");\r
+ // hdeltaPhiDjhB->GetYaxis()->SetTitle("#Delta#phi (rad)");\r
+ fOutput->Add(hdeltaPhiDjaB);\r
+ // fOutput->Add(hdeltaPhiDjlB);\r
+ // fOutput->Add(hdeltaPhiDjhB);\r
+\r
+ TH2F* hInvMassptD = new TH2F("hInvMassptD","D (Delta R < 0.4) invariant mass distribution p_{T}^{j} > threshold",nbinsmass,fMinMass,fMaxMass,100,0.,50.);\r
+ hInvMassptD->SetStats(kTRUE);\r
+ hInvMassptD->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ hInvMassptD->GetYaxis()->SetTitle("p_{t}^{D} (GeV/c)");\r
+\r
+ fOutput->Add(hInvMassptD);\r
+ // fMasspjDeltaR=new TH3F("fMasspjDeltaR","Mass vs p^{jet} vs #Delta R;Mass (Gev/c);p^{jet}(Gev/c^{2});#Delta R",1500,fMinMass,fMaxMass,100,0.,50.,100,0.,1.);\r
+ // fMasspjDeltaR->SetStats(kFALSE);\r
+ // fOutput->Add(fMasspjDeltaR);\r
+\r
+ TH3F* hzptD=new TH3F("hzptD","Fragmentation function (DeltaR < 0.4)",100,0.,1.2,nbinsmass,fMinMass,fMaxMass,100,0,50);\r
+ hzptD->SetStats(kTRUE);\r
+ hzptD->GetXaxis()->SetTitle("z=p_{D} #cdot p_{j}/|p_{j}|^{2}");\r
+ hzptD->GetYaxis()->SetTitle("mass (GeV/c)");\r
+ hzptD->GetZaxis()->SetTitle("p_{t}^{D} (GeV/c)");\r
+ fOutput->Add(hzptD);\r
+\r
+ TH3F* hzptDB=new TH3F("hzptDB","Fragmentation function (DeltaR < 0.4) - Side Bands",100,0.,1.2,nbinsmass,fMinMass,fMaxMass,100,0.,50.);\r
+ hzptDB->SetStats(kTRUE);\r
+ hzptDB->GetXaxis()->SetTitle("z=p_{D} #cdot p_{j}/|p_{j}|^{2}");\r
+ hzptDB->GetYaxis()->SetTitle("mass (GeV/c)");\r
+ hzptDB->GetZaxis()->SetTitle("p_{t}^{D} (GeV/c)");\r
+ fOutput->Add(hzptDB);\r
+ //TH1F* hResZ = new TH1F("hResZ","Fragmentation function ",50,0,1);\r
+ // TH1F* hResZBkg = new TH1F("hResZBkg","Fragmentation function background",50,0,1); \r
+\r
+ //fOutput->Add(hResZ);\r
+ //fOutput->Add(hResZBkg);\r
+\r
+\r
+ return kTRUE; \r
+}\r
+\r
+void AliAnalysisTaskFlavourJetCorrelations::FillHistogramsRecoJetCorr(AliAODRecoDecayHF* candidate, AliEmcalJet *jet){\r
+\r
+ Double_t ptD=candidate->Pt();\r
+ Double_t ptjet=jet->Pt();\r
+ Double_t deltaR=DeltaR(candidate,jet);\r
+ Double_t deltaphi = jet->Phi()-candidate->Phi();\r
+ if(deltaphi<=-(TMath::Pi())/2) deltaphi = deltaphi+2*(TMath::Pi());\r
+ if(deltaphi>(3*(TMath::Pi()))/2) deltaphi = deltaphi-2*(TMath::Pi());\r
+ Double_t z=Z(candidate,jet);\r
+\r
+ TH1F* hDeltaRD=(TH1F*)fOutput->FindObject("hDeltaRD");\r
+ hDeltaRD->Fill(deltaR); \r
+ TH1I* hstat=(TH1I*)fOutput->FindObject("hstat");\r
+ hstat->Fill(4);\r
+ TH1F* hPtJetWithD=(TH1F*)fOutput->FindObject("hPtJetWithD");\r
+ hPtJetWithD->Fill(ptjet);\r
+\r
+ if(fCandidateType==kD0toKpi) {\r
+\r
+ FillHistogramsD0JetCorr(candidate,deltaphi,z,ptD,deltaR, AODEvent());\r
+\r
+ }\r
+\r
+ if(fCandidateType==kDstartoKpipi) {\r
+ AliAODRecoCascadeHF* dstar = (AliAODRecoCascadeHF*)candidate;\r
+ FillHistogramsDstarJetCorr(dstar,deltaphi,z,ptD,deltaR);\r
+\r
+ }\r
+\r
+}\r
+\r
+void AliAnalysisTaskFlavourJetCorrelations::FillHistogramsD0JetCorr(AliAODRecoDecayHF* candidate, Double_t dPhi, Double_t z, Double_t ptD, Double_t deltaR, AliAODEvent* aodEvent){\r
+ Double_t masses[2]={0.,0.};\r
+ Int_t pdgdaughtersD0[2]={211,321};//pi,K \r
+ Int_t pdgdaughtersD0bar[2]={321,211};//K,pi \r
+\r
+ masses[0]=candidate->InvMass(fNProngs,(UInt_t*)pdgdaughtersD0); //D0\r
+ masses[1]=candidate->InvMass(fNProngs,(UInt_t*)pdgdaughtersD0bar); //D0bar\r
+ \r
+ \r
+ Int_t isselected=fCuts->IsSelected(candidate,AliRDHFCuts::kAll,aodEvent);\r
+ if(isselected==1 || isselected==3) {\r
+\r
+ FillHistogramsD(masses[0],dPhi,z, ptD, deltaR);\r
+ }\r
+ if(isselected>=2) {\r
+\r
+ FillHistogramsD(masses[1],dPhi,z, ptD, deltaR);\r
+ }\r
+\r
+}\r
+\r
+void AliAnalysisTaskFlavourJetCorrelations::FillHistogramsDstarJetCorr(AliAODRecoCascadeHF* dstar, Double_t dPhi, Double_t z, Double_t ptD, Double_t deltaR){\r
+ AliAODTrack *softpi = (AliAODTrack*)dstar->GetBachelor();\r
+ Double_t deltamass= dstar->DeltaInvMass(); \r
+ Double_t massD0= dstar->InvMassD0();\r
+\r
+ TH1F* hPtPion=(TH1F*)fOutput->FindObject("hPtPion");\r
+ hPtPion->Fill(softpi->Pt());\r
+ FillHistogramsD(deltamass,dPhi,z, ptD, deltaR);\r
+ // evaluate side band background\r
+ TH2F* hDiffSideBand=(TH2F*)fOutput->FindObject("hDiffSideBand");\r
+\r
+ TH3F* hdeltaPhiDjaB=(TH3F*)fOutput->FindObject("hdeltaPhiDjaB");\r
+\r
+ TH3F* hzptDB=(TH3F*)fOutput->FindObject("hzptDB");\r
+\r
+ Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
+\r
+ Int_t bin = fCuts->PtBin(ptD);\r
+ Float_t fourSigmal= mPDGD0-3.5*fSigmaD0[bin] , sixSigmal= mPDGD0-5.*fSigmaD0[bin];\r
+ Float_t fourSigmar= mPDGD0+3.5*fSigmaD0[bin] , sixSigmar= mPDGD0+5.*fSigmaD0[bin];\r
+\r
+ if((massD0>sixSigmal && massD0<fourSigmal) || (massD0>fourSigmar && massD0<=sixSigmar)){ \r
+ hDiffSideBand->Fill(deltamass,ptD); // M(Kpipi)-M(Kpi) side band background \r
+ hdeltaPhiDjaB->Fill(deltamass,ptD,dPhi);\r
+\r
+ if(deltaR<0.4){ // evaluate in the near side \r
+ hzptDB->Fill(z,deltamass,ptD); \r
+ }\r
+\r
+ } // SideBandBackground(dstar, dPhi, z, ptD, deltaR);\r
+ \r
+}\r
+\r
+void AliAnalysisTaskFlavourJetCorrelations::FillHistogramsD(Double_t mass,Double_t dphi, Double_t z,Double_t ptD, Double_t deltaR){\r
+ TH3F* hdeltaPhiDja=((TH3F*)fOutput->FindObject("hdeltaPhiDja"));\r
+ hdeltaPhiDja->Fill(mass,ptD,dphi);\r
+\r
+ if(deltaR<0.4) {\r
+ TH3F* hzptD=(TH3F*)fOutput->FindObject("hzptD");\r
+ hzptD->Fill(z,mass,ptD);\r
+\r
+ TH2F* hInvMassptD=(TH2F*)fOutput->FindObject("hInvMassptD");\r
+ hInvMassptD->Fill(mass,ptD);\r
+ }\r
+}\r
+//______________________________ side band background for D*___________________________________\r
+\r
+// void AliAnalysisTaskFlavourJetCorrelations::SideBandBackground(AliAODRecoCascadeHF *candDstar, Double_t dPhi, Double_t z, Double_t ptD, Double_t deltaR){\r
+\r
+// // D* side band background method. Two side bands, in M(Kpi) are taken at ~6 sigmas \r
+// // (expected detector resolution) on the left and right frm the D0 mass. Each band\r
+// // has a width of ~5 sigmas. Two band needed for opening angle considerations \r
+// TH2F* hDiffSideBand=(TH2F*)fOutput->FindObject("hDiffSideBand");\r
+\r
+// TH3F* hdeltaPhiDjaB=(TH3F*)fOutput->FindObject("hdeltaPhiDjaB");\r
+\r
+// TH3F* hzptDB=(TH3F*)fOutput->FindObject("hzptDB");\r
+\r
+// Double_t mPDGD0=TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
+\r
+// Int_t bin = fCuts->PtBin(candDstar->Pt());\r
+// Float_t fourSigmal= mPDGD0-3.5*fSigmaD0[bin] , sixSigmal= mPDGD0-5.*fSigmaD0[bin];\r
+// Float_t fourSigmar= mPDGD0+3.5*fSigmaD0[bin] , sixSigmar= mPDGD0+5.*fSigmaD0[bin];\r
+\r
+// Double_t invM=candDstar->InvMassD0(), deltaM=candDstar->DeltaInvMass(); //invMDstar=candDstar->InvMassDstarKpipi(),\r
+// Printf("Inv mass = %f between %f and %f or %f and %f?",invM, sixSigmal,fourSigmal,fourSigmar,sixSigmar);\r
+// Double_t ptD=candDstar->Pt();\r
+// //Double_t ptjet=jet->Pt();\r
+// Double_t dPhi=jet->Phi()-candDstar->Phi();\r
+// Double_t deltaR=DeltaR(candDstar,jet);\r
+// if(dPhi>TMath::Pi())dPhi = dPhi - 2.*TMath::Pi();\r
+// if(dPhi<(-1.*TMath::Pi()))dPhi = dPhi + 2.*TMath::Pi();\r
+\r
+// if((invM>sixSigmal && invM<fourSigmal) || (invM>fourSigmar && invM<=sixSigmar)){ \r
+// hDiffSideBand->Fill(deltaM,ptD); // M(Kpipi)-M(Kpi) side band background \r
+// hdeltaPhiDjaB->Fill(deltaM,ptD,dPhi);\r
+\r
+// if(deltaR<0.4){ // evaluate in the near side \r
+// hzptDB->Fill(Z(candDstar,jet),deltaM,ptD); \r
+// }\r
+\r
+// }\r
+// }\r
+\r
+Float_t AliAnalysisTaskFlavourJetCorrelations::DeltaR(AliVParticle *p1, AliVParticle *p2) const {\r
+ //Calculate DeltaR between p1 and p2: DeltaR=sqrt(Delataphi^2+DeltaEta^2)\r
+\r
+ if(!p1 || !p2) return -1;\r
+ Double_t phi1=p1->Phi(),eta1=p1->Eta();\r
+ Double_t phi2 = p2->Phi(),eta2 = p2->Eta() ;\r
+\r
+ Double_t dPhi=phi1-phi2;\r
+ if(dPhi<=-(TMath::Pi())/2) dPhi = dPhi+2*(TMath::Pi());\r
+ if(dPhi>(3*(TMath::Pi()))/2) dPhi = dPhi-2*(TMath::Pi());\r
+\r
+ Double_t dEta=eta1-eta2;\r
+ Double_t deltaR=TMath::Sqrt(dEta*dEta + dPhi*dPhi );\r
+ return deltaR;\r
+ \r
+}\r
+\r
+/*\r
+//_____________________________________________________________________\r
+\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::IsD(Int_t pdg) const{\r
+ Int_t abspdg=TMath::Abs(pdg);\r
+ if(abspdg>400 && abspdg<500) return kTRUE;\r
+ else return kFALSE;\r
+}\r
+\r
+//_____________________________________________________________________\r
+\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::IsD(Int_t pdg,Int_t abspdgD) const{\r
+ Int_t abspdg=TMath::Abs(pdg);\r
+ if(abspdg==abspdgD) return kTRUE;\r
+ else return kFALSE;\r
+}\r
+\r
+//_____________________________________________________________________\r
+\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::PartFromC(AliMCParticle* mother) const{\r
+ Int_t pdgmoth=mother->PdgCode();\r
+ if(TMath::Abs(pdgmoth)==4) return kTRUE;\r
+ else return kFALSE;\r
+}\r
+\r
+Int_t AliAnalysisTaskFlavourJetCorrelations::GetFirstMother(Int_t labpart, TClonesArray *mcArray) const{\r
+ AliAODMCParticle* partMC=(AliAODMCParticle*)mcArray->UncheckedAt(labpart);\r
+ if (!partMC) return -2;\r
+ Int_t labmom=labpart;\r
+ Printf("Starting from %d",labmom);\r
+ while(labmom>-1){\r
+ labmom=labpart;\r
+ partMC=(AliAODMCParticle*)mcArray->UncheckedAt(labmom);\r
+ if (!partMC) return -2;\r
+ labmom= partMC->GetMother();\r
+ Printf("Lab mom %d",labmom);\r
+ }\r
+ Printf("Return labpart %d", labpart);\r
+ return labpart;\r
+}\r
+\r
+\r
+// -------------------------------------- check the PDG -----------------------------------------\r
+\r
+Int_t AliAnalysisTaskFlavourJetCorrelations::FindPDGInFamily(Int_t labpart,Int_t pdgcode, TClonesArray *mcArray) const{\r
+\r
+ //return the label of the particle which is a "pdgcode" in the family\r
+ Printf("FindPDGInFamily label %d, pdg %d, mcarray %p",labpart,pdgcode,mcArray);\r
+ AliAODMCParticle* partMC=(AliAODMCParticle*)mcArray->UncheckedAt(labpart);\r
+ if (!partMC) return -2;\r
+ Int_t labmom=labpart;\r
+ Printf("Starting from %d",labmom);\r
+ while(labmom>-1){\r
+\r
+ partMC=(AliAODMCParticle*)mcArray->UncheckedAt(labmom);\r
+ if (!partMC) return -2;\r
+ Int_t pdgmom=partMC->GetPdgCode();\r
+ if(pdgmom==pdgcode) return labmom;\r
+ labmom= partMC->GetMother();\r
+ Printf("Lab mom %d",labmom);\r
+ }\r
+\r
+ return -1;\r
+\r
+}\r
+\r
+// ------------------------- check on MC the distance between D meson and jet ----------------------\r
+\r
+Bool_t AliAnalysisTaskFlavourJetCorrelations::FillMCDJetInfo(AliPicoTrack *jetTrk, AliEmcalJet* jet, TClonesArray *mcArray, Double_t ptjet){\r
+ \r
+ Bool_t foundD = kFALSE;\r
+ vector<int> DmesonInJetLabels(10);\r
+ \r
+ Int_t jtlabel=jetTrk->GetLabel();\r
+ if(jtlabel<=0) return foundD;\r
+ AliAODMCParticle* jetMCpart=(AliAODMCParticle*)mcArray->UncheckedAt(jtlabel);\r
+ if(!jetMCpart) return foundD;\r
+ printf("AliMCParticle %d, %p\n",1,jetMCpart);\r
+ \r
+ Int_t labDmeson=FindPDGInFamily(jtlabel,fPDGmother, mcArray);\r
+ if(labDmeson>0){\r
+ AliAODMCParticle *partDmeson=(AliAODMCParticle*)mcArray->UncheckedAt(labDmeson);\r
+ fhMomjetpartPdg->Fill(partDmeson->GetPdgCode());\r
+ \r
+ //tmp\r
+ Int_t momjetpartlabel=labDmeson;\r
+ \r
+ Int_t iD=5;\r
+ Bool_t exists=kFALSE;\r
+ for(Int_t k=0;k<iD;k++){\r
+ if(momjetpartlabel==DmesonInJetLabels[k]) {//mother already found\r
+ exists=kTRUE;\r
+ break;\r
+ }\r
+ }\r
+ if(!exists) {\r
+ DmesonInJetLabels[iD]=momjetpartlabel;\r
+ AliDebug(2,Form("D meson number %d found: label %d\n",iD,DmesonInJetLabels[iD]));\r
+ hstat->Fill(6);\r
+ iD++;\r
+ foundD=kTRUE;\r
+ \r
+ }\r
+ }\r
+ \r
+ if(fUseMCInfo && foundD) {\r
+ fptJwithDMC->Fill(ptjet); //filled only once per jet, not per each D meson\r
+ Int_t iD=5;\r
+ // loop over the D within the jet \r
+ for(Int_t kD=0;kD<iD;kD++){\r
+ \r
+ AliAODMCParticle* momjetpart=(AliAODMCParticle*)mcArray->At(DmesonInJetLabels[kD]);\r
+ Double_t ptD=momjetpart->Pt(),etaD=momjetpart->Eta(), phiD=momjetpart->Phi();\r
+ fptDinjetallvsptjMC->Fill(ptD,ptjet);\r
+ \r
+ Double_t deltaRD=DeltaR(jet,momjetpart);\r
+ \r
+ ((TH1F*)fOutput->FindObject("hdeltaRDMC"))->Fill(deltaRD); //Delta R of D mesons (MC)\r
+ \r
+ Double_t z=Z(momjetpart,jet);\r
+ \r
+ if(deltaRD<0.4) {\r
+ hstat->Fill(7);\r
+ //comment if you prefer to ask for DeltaR of the daughters < 0.4 and uncomment below\r
+ fptDinjet04vsptjMC->Fill(ptD,ptjet);\r
+ fzptDptj->Fill(z,ptD,ptjet);\r
+ }\r
+ \r
+ \r
+ }//end loop on MC D\r
+ \r
+ return foundD;\r
+ \r
+ }\r
+}\r
+*/ \r
-#ifndef ALIANALYSISTASKSERECOJETCORRELATIONS_H
-#define ALIANALYSISTASKSERECOJetCORRELATIONS_H
-
-// $Id$
-
-/* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-
-//*************************************************************************
-// Class AliAnalysisTaskFlavourJetCorrelations
-// AliAnalysisTaskSE for Dmesons - jet correlations analysis
-// Author: Xiaoming Zhang, xmzhang@lbl.gov
-//*************************************************************************
-
-#include "AliAnalysisTaskEmcalJet.h"
-
-class TList;
-class TClonesArray;
-
-class AliAnalysisTaskFlavourJetCorrelations : public AliAnalysisTaskEmcalJet {
-
- public :
-
- enum {
- kMatchConeCandi,
- kMatchAreaCandi,
- kMatchConeProng,
- kMatchAreaProng,
- kDzeroMatchType
- };
-
- AliAnalysisTaskFlavourJetCorrelations();
- AliAnalysisTaskFlavourJetCorrelations(const char *name, Bool_t bIsHisto=kTRUE);
- virtual ~AliAnalysisTaskFlavourJetCorrelations();
-
- virtual void UserCreateOutputObjects();
-
- private :
-
- AliAnalysisTaskFlavourJetCorrelations(const AliAnalysisTaskFlavourJetCorrelations &);
- AliAnalysisTaskFlavourJetCorrelations& operator=(const AliAnalysisTaskFlavourJetCorrelations &);
-
- virtual void ExecOnce();
- virtual Bool_t FillGeneralHistograms();
- virtual Bool_t FillHistograms();
- virtual Bool_t IsEventSelected();
- virtual Bool_t RetrieveEventObjects();
- virtual Bool_t Run();
-
- void MakeControlHistograms();
- Bool_t FillControlHistograms();
-
- void CreateDzeroHistograms();
- void CreateDstarHistograms();
-
- void RunDzeroJet(AliEmcalJet const *pJet, const Int_t iJetPtBin, const Bool_t bIsD0);
- void RunDstarJet(AliEmcalJet const *pJet, const Int_t iJetPtBin);
-
- TClonesArray *fUsedDzeros; //! input Dzero candidates array
- TClonesArray *fUsedD0bars; //! input D0bar candidates array
- TClonesArray *fUsedDstars; //! input Dstar candidates array
-
- TList *fListControlHistos; //! list of output contral histograms
- TList *fListAnDzeroHistos; //! list of output Dzero - jet correlation histograms
- TList *fListAnDstarHistos; //! list of output Dstar - jet correlation histograms
-
- ClassDef(AliAnalysisTaskFlavourJetCorrelations, 1);
-};
-
-#endif
+#ifndef ALIANALYSISTASKFLAVOURJETCORRELATIONS_H\r
+#define ALIANALYSISTASKFLAVOURJETCORRELATIONS_H\r
+// $Id$\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: The ALICE Off-line Project. *\r
+ * Contributors are mentioned in the code where appropriate. *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+//-----------------------------------------------------------------------\r
+// Author : A. Grelli, Utrecht University\r
+// C. Bianchin, Utrecht University\r
+// X. Zhang, LBNL\r
+//-----------------------------------------------------------------------\r
+\r
+\r
+#include <TH2F.h>\r
+#include "AliAODEvent.h"\r
+#include "AliPicoTrack.h"\r
+#include "AliAnalysisTaskEmcalJet.h"\r
+\r
+class TH3F;\r
+class TParticle ;\r
+class TClonesArray ;\r
+class AliMCParticle;\r
+class AliAODMCParticle;\r
+class AliRDHFCuts;\r
+class AliEmcalJet;\r
+class AliAODRecoDecayHF;\r
+class AliAODRecoCascadeHF;\r
+class AliAODEvent;\r
+\r
+class AliAnalysisTaskFlavourJetCorrelations : public AliAnalysisTaskEmcalJet \r
+{\r
+\r
+ public:\r
+\r
+ enum ECandidateType{ kD0toKpi, kDstartoKpipi };\r
+\r
+ AliAnalysisTaskFlavourJetCorrelations();\r
+ AliAnalysisTaskFlavourJetCorrelations(const Char_t* name,AliRDHFCuts* cuts, ECandidateType candtype, TString jetArrName);\r
+ virtual ~AliAnalysisTaskFlavourJetCorrelations();\r
+\r
+ virtual void UserCreateOutputObjects();\r
+ virtual void UserExec(Option_t *option);\r
+ virtual void Terminate(Option_t *);\r
+ virtual void Init();\r
+ virtual void LocalInit() {Init();}\r
+\r
+ // inizializations\r
+ Bool_t DefineHistoForAnalysis(); \r
+\r
+ // set MC usage\r
+ void SetMC(Bool_t theMCon) {fUseMCInfo = theMCon;}\r
+ Bool_t GetMC() const {return fUseMCInfo;}\r
+ \r
+ void SetMassLimits(Double_t range, Int_t pdg);\r
+ void SetMassLimits(Double_t lowlimit, Double_t uplimit);\r
+\r
+ //jet reconstruction algorithm\r
+ void SetJetArrayName(TString jetArrName) {fJetArrName=jetArrName;};\r
+ TString GetJetArrayName() const {return fJetArrName;};\r
+\r
+ // trigger on jet events\r
+ void SetTriggerOnLeadingJet(Bool_t triggerOnLeadingJet) {fLeadingJetOnly=triggerOnLeadingJet;};\r
+ Bool_t GetTriggerOnLeadingJet() const {return fLeadingJetOnly;}\r
+\r
+\r
+ // Array of D0 width for the Dstar\r
+ Bool_t SetD0WidthForDStar(Int_t nptbins,Float_t* width);\r
+\r
+ //Bool_t FillMCDJetInfo(AliPicoTrack *jetTrk,AliEmcalJet* jet, TClonesArray *mcArray,Double_t ptjet);\r
+ void FillHistogramsRecoJetCorr(AliAODRecoDecayHF* candidate, AliEmcalJet *jet);\r
+ void FillHistogramsD0JetCorr(AliAODRecoDecayHF* candidate, Double_t dPhi, Double_t z, Double_t ptD, Double_t deltaR, AliAODEvent* aodEvent);\r
+\r
+ void FillHistogramsDstarJetCorr(AliAODRecoCascadeHF* dstar, Double_t dPhi, Double_t z, Double_t ptD, Double_t deltaR);\r
+\r
+ void FillHistogramsD(Double_t mass,Double_t dphi, Double_t z,Double_t ptD, Double_t deltaR);\r
+\r
+ private:\r
+ \r
+ AliAnalysisTaskFlavourJetCorrelations(const AliAnalysisTaskFlavourJetCorrelations &source);\r
+ AliAnalysisTaskFlavourJetCorrelations& operator=(const AliAnalysisTaskFlavourJetCorrelations& source); \r
+\r
+ Double_t Z(AliVParticle* part,AliEmcalJet* jet) const;\r
+ Float_t DeltaR(AliVParticle *p1, AliVParticle *p2) const;\r
+\r
+ /*\r
+ Bool_t IsD(Int_t pdg) const;\r
+ Bool_t IsD(Int_t pdg,Int_t abspdgD) const;\r
+ Bool_t PartFromC(AliMCParticle* mother) const;\r
+ Int_t GetFirstMother(Int_t lab,TClonesArray* mcarr) const; \r
+ Int_t FindPDGInFamily(Int_t labpart,Int_t pdgcode, TClonesArray *mcArray) const;\r
+ */\r
+\r
+ Bool_t fUseMCInfo; // Use MC info\r
+ Int_t fCandidateType; // Dstar or D0\r
+ Int_t fPDGmother; // PDG code of D meson\r
+ Int_t fNProngs; // number of prong of the decay channel \r
+ Int_t fPDGdaughters[4]; // PDG codes of daughters\r
+ Float_t fSigmaD0[30]; //\r
+ TString fBranchName; // AOD branch name\r
+ TList *fOutput; //! user output\r
+ AliRDHFCuts *fCuts; // Cuts \r
+\r
+ Double_t fMinMass; // mass lower limit histogram\r
+ Double_t fMaxMass; // mass upper limit histogram\r
+\r
+ TString fJetArrName; // name of the jet array, taken from the task running the jet finder\r
+ TString fCandArrName; // string which correspond to the candidate type\r
+ Bool_t fLeadingJetOnly; // use only the leading jet in the event to make the correlations\r
+\r
+ ClassDef(AliAnalysisTaskFlavourJetCorrelations,1); // class for charm-jet correlations\r
+};\r
+\r
+#endif\r
-// $Id$
-
-/**************************************************************************
- * 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. *
- **************************************************************************/
-
-/////////////////////////////////////////////////////////////
-//
-// AliAnalysisTaskSE for Dmesons - jet correlations analysis
-//
-// Author: Xiaoming Zhang, xmzhang@lbl.gov
-///////////////////////////////////////////////////////////////
-
-#include <TMath.h>
-#include <TH1D.h>
-#include <TH2D.h>
-#include <TList.h>
-#include <TClonesArray.h>
-#include <TDatabasePDG.h>
-
-#include "AliAnalysisManager.h"
-#include "AliAODTrack.h"
-#include "AliAODEvent.h"
-#include "AliAODHandler.h"
-#include "AliAODExtension.h"
-#include "AliAODRecoDecayHF.h"
-#include "AliAODRecoDecayHF2Prong.h"
-#include "AliAODRecoCascadeHF.h"
-#include "AliRDHFCuts.h"
-#include "AliRDHFCutsD0toKpi.h"
-#include "AliRDHFCutsDStartoKpipi.h"
-#include "AliAnalysisTaskSEDmesonsFilterCJ.h"
-
-ClassImp(AliAnalysisTaskSEDmesonsFilterCJ)
-
-//_____________________________________________________________________________
-AliAnalysisTaskSEDmesonsFilterCJ::AliAnalysisTaskSEDmesonsFilterCJ() :
-AliAnalysisTaskSE(),
-fEventAOD(0),
-fCutDzero(0),
-fCutDstar(0),
-fDzeroClArr(0),
-fDstarClArr(0),
-fUsedDzeros(0),
-fUsedD0bars(0),
-fUsedDstars(0),
-fIsNotExecOnce(kTRUE),
-fListCutsDmesons(0),
-fListDzeroHistos(0),
-fListDstarHistos(0)
-{
-//
-// Default constructor
-//
-}
-
-//_____________________________________________________________________________
-AliAnalysisTaskSEDmesonsFilterCJ::AliAnalysisTaskSEDmesonsFilterCJ(const char *name) :
-AliAnalysisTaskSE(name),
-fEventAOD(0),
-fCutDzero(0),
-fCutDstar(0),
-fDzeroClArr(0),
-fDstarClArr(0),
-fUsedDzeros(0),
-fUsedD0bars(0),
-fUsedDstars(0),
-fIsNotExecOnce(kTRUE),
-fListCutsDmesons(0),
-fListDzeroHistos(0),
-fListDstarHistos(0)
-{
-//
-// Constructor
-//
-
- DefineOutput(1, TList::Class());
- DefineOutput(2, TList::Class());
- DefineOutput(3, TList::Class());
-}
-
-//_____________________________________________________________________________
-AliAnalysisTaskSEDmesonsFilterCJ::~AliAnalysisTaskSEDmesonsFilterCJ()
-{
-//
-// Default destructor
-//
-
- if (fEventAOD) { delete fEventAOD; fEventAOD =0; }
-
- if (fCutDzero) { delete fCutDzero; fCutDzero =0; }
- if (fCutDstar) { delete fCutDstar; fCutDstar =0; }
-
- if (fDzeroClArr) { delete fDzeroClArr; fDzeroClArr=0; }
- if (fDstarClArr) { delete fDstarClArr; fDstarClArr=0; }
-
- if (fUsedDzeros) { delete fUsedDzeros; fUsedDzeros=0; }
- if (fUsedD0bars) { delete fUsedD0bars; fUsedD0bars=0; }
- if (fUsedDstars) { delete fUsedDstars; fUsedDstars=0; }
-
- if (fListCutsDmesons) { delete fListCutsDmesons; fListCutsDmesons=0; }
- if (fListDzeroHistos) { delete fListDzeroHistos; fListDzeroHistos=0; }
- if (fListDstarHistos) { delete fListDstarHistos; fListDstarHistos=0; }
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::Init()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::Init()
-//
-
- if (fListCutsDmesons) return;
- fListCutsDmesons = new TList(); fListCutsDmesons->SetOwner();
-
- if (fCutDzero) {
- AliRDHFCutsD0toKpi *cutDzero = new AliRDHFCutsD0toKpi(*fCutDzero);
- cutDzero->SetName("AnalysisCutsDzero");
- fListCutsDmesons->Add(cutDzero);
- }
-
- if (fCutDstar) {
- AliRDHFCutsDStartoKpipi *cutDstar = new AliRDHFCutsDStartoKpipi(*fCutDstar);
- cutDstar->SetName("AnalysisCutsDstar");
- fListCutsDmesons->Add(cutDstar);
- }
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::UserCreateOutputObjects()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::UserCreateOutputObjects
-//
-
- if (fCutDzero) { fUsedDzeros = new TClonesArray("AliAODRecoDecayHF2Prong",0); fUsedDzeros->SetName("AnaUsedDzero");}
- if (fCutDzero) { fUsedD0bars = new TClonesArray("AliAODRecoDecayHF2Prong",0); fUsedD0bars->SetName("AnaUsedD0bar");}
- if (fCutDstar) { fUsedDstars = new TClonesArray("AliAODRecoCascadeHF", 0); fUsedDstars->SetName("AnaUsedDstar"); }
-
- if (!fListDzeroHistos) fListDzeroHistos = new TList(); fListDzeroHistos->SetOwner();
- if (!fListDstarHistos) fListDstarHistos = new TList(); fListDstarHistos->SetOwner();
-
- if (fCutDzero) CreateDzeroHistograms();
- if (fCutDstar) CreateDstarHistograms();
-
- PostData(1, fListCutsDmesons);
- PostData(2, fListDzeroHistos);
- PostData(3, fListDstarHistos);
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::UserExec(Option_t *)
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::UserExec
-//
-
- if (fIsNotExecOnce) {
- fIsNotExecOnce=ExecOnce();
- if (fIsNotExecOnce) return;
- }
-
- if (ExecEach()) ExecAnas();
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::Terminate(Option_t *)
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::Terminate
-//
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::NotifyRun()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::NotifyRun
-//
-
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::ExecAnas()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::ExecAnas
-//
-
- if (fUsedDzeros) fUsedDzeros->Delete();
- if (fUsedD0bars) fUsedD0bars->Delete();
- if (fUsedDstars) fUsedDstars->Delete();
-
- if (fDzeroClArr) ExecAnasDzero();
- if (fDstarClArr) ExecAnasDstar();
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::ExecAnasDzero()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::ExecAnasDzero
-//
-
- if (!fCutDzero->IsEventSelected(fEventAOD)) return;
- const Int_t nCands = fDzeroClArr->GetEntriesFast(); if (nCands==0) return;
-
- ((TH1D*)fListDzeroHistos->FindObject("hDzero_candi"))->Fill(nCands);
-
- Int_t countN = 0;
- AliAODRecoDecayHF2Prong *candi = 0;
- for (Int_t iCand=0; iCand<nCands; iCand++) {
- candi = dynamic_cast<AliAODRecoDecayHF2Prong*>(fDzeroClArr->At(iCand)); if (!candi) continue;
-
- Double_t dPt = candi->Pt();
- if (!fCutDzero->IsInFiducialAcceptance(dPt,candi->YD0())) continue;
- Int_t mask = fCutDzero->IsSelected(candi,AliRDHFCuts::kAll,fEventAOD); if (mask==0) continue;
-
- Double_t dMass = 0.;
- if ((mask==1) || (mask==3)) dMass = candi->InvMassD0();
- if (mask==2) dMass = candi->InvMassD0bar();
- ((TH2D*)fListDzeroHistos->FindObject("hDzero_Mass_Pt"))->Fill(dMass,dPt);
-
- if ((mask==1) || (mask==3)) new((*fUsedDzeros)[countN++]) AliAODRecoDecayHF2Prong(*candi);
- if (mask==2) new((*fUsedD0bars)[countN++]) AliAODRecoDecayHF2Prong(*candi);
- candi = 0;
- }
-
- ((TH1D*)fListDzeroHistos->FindObject("hDzero_seled"))->Fill(countN);
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::ExecAnasDstar()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::ExecAnasDstar
-//
-
- if (!fCutDstar->IsEventSelected(fEventAOD)) return;
- const Int_t nCands = fDstarClArr->GetEntriesFast(); if (nCands==0) return;
-
- Double_t dMassDzero = TDatabasePDG::Instance()->GetParticle(421)->Mass();
- Double_t dMassDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
- Double_t dMassDelta = dMassDstar - dMassDzero;
- Double_t dSigmaD0Pt[13]={ 0.012, 0.012, 0.012, 0.015, 0.015, 0.018, 0.018, 0.020, 0.020, 0.030, 0.030, 0.037, 0.040 };
-
- ((TH1D*)fListDstarHistos->FindObject("hDstar_candi"))->Fill(nCands);
-
- Int_t countN = 0;
- AliAODRecoCascadeHF *candi = 0;
- for (Int_t iCand=0; iCand<nCands; iCand++) {
- candi = dynamic_cast<AliAODRecoCascadeHF*>(fDstarClArr->At(iCand)); if (!candi) continue;
-
- Double_t dPt = candi->Pt();
- if (!fCutDstar->IsInFiducialAcceptance(dPt,candi->YDstar())) continue;
- if (!fCutDstar->IsSelected(candi,AliRDHFCuts::kAll)) continue;
-
- Int_t iBin = fCutDzero->PtBin(dPt);
- if ((iBin<0) || (iBin>=fCutDzero->GetNPtBins())) { AliWarning(Form("pT(D*)=%f out of bounds!!!",dPt)); continue; }
-
- Double_t dMassResD0 = TMath::Abs(candi->InvMassD0() - dMassDzero); if (dMassResD0>3.*dSigmaD0Pt[iBin]) continue;
-
- Double_t deltaM = candi->DeltaInvMass();
- ((TH2D*)fListDstarHistos->FindObject("hDstar_deltaMass_Pt"))->Fill(deltaM,dPt);
-
- if (TMath::Abs(deltaM-dMassDelta)) ((TH1D*)fListDstarHistos->FindObject("hDstar_SoftPiPt"))->Fill(((AliAODTrack*)candi->GetBachelor())->Pt());
-
- new((*fUsedDstars)[countN++]) AliAODRecoCascadeHF(*candi);
- candi = 0;
- }
-
- ((TH1D*)fListDstarHistos->FindObject("hDstar_seled"))->Fill(countN);
- return;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskSEDmesonsFilterCJ::ExecOnce()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::ExecOnce
-//
-
- if (!InputEvent()) return kTRUE;
- fEventAOD = dynamic_cast<AliAODEvent*>(InputEvent());
-
- if (fEventAOD) {
- if (fCutDzero) fDzeroClArr = dynamic_cast<TClonesArray*>(fEventAOD->FindListObject("D0toKpi"));
- if (fCutDstar) fDstarClArr = dynamic_cast<TClonesArray*>(fEventAOD->FindListObject("Dstar"));
- } else if (AODEvent() && IsStandardAOD()) {
- fEventAOD = dynamic_cast<AliAODEvent*>(AODEvent());
- AliAODHandler* aodH = dynamic_cast<AliAODHandler*>((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());
-
- if(aodH->GetExtensions()) {
- AliAODExtension *aodExt = dynamic_cast<AliAODExtension*>(aodH->GetExtensions()->FindObject("AliAOD.VertexingHF.root"));
- if (fCutDzero) fDzeroClArr = dynamic_cast<TClonesArray*>(aodExt->GetAOD()->FindListObject("D0toKpi"));
- if (fCutDstar) fDstarClArr = dynamic_cast<TClonesArray*>(aodExt->GetAOD()->FindListObject("Dstar"));
- }
- }
-
- if (!fEventAOD) { AliError("No AOD Event!!!"); return kTRUE; }
- if (!fDzeroClArr && fCutDzero) { AliError("No Dzero Branch!!!"); return kTRUE; }
- if (!fDstarClArr && fCutDstar) { AliError("No Dstar Branch!!!"); return kTRUE; }
-
- if (fDzeroClArr && fUsedDzeros) {
- fUsedDzeros->Delete();
- if (!(InputEvent()->FindListObject("AnaUsedDzero"))) InputEvent()->AddObject(fUsedDzeros);
- }
-
- if (fDzeroClArr && fUsedD0bars) {
- fUsedD0bars->Delete();
- if (!(InputEvent()->FindListObject("AnaUsedD0bar"))) InputEvent()->AddObject(fUsedD0bars);
- }
-
- if (fDstarClArr && fUsedDstars) {
- fUsedDstars->Delete();
- if (!(InputEvent()->FindListObject("AnaUsedDstar"))) InputEvent()->AddObject(fUsedDstars);
- }
-
- return kFALSE;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskSEDmesonsFilterCJ::ExecEach()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::ExecEach
-//
-
- if (!fEventAOD) {
- AliWarning("No Input Event, Skip Event!!!");
- return kFALSE;
- }
-
- if ((!fDzeroClArr) && (!fDstarClArr)) {
- AliWarning("No D meson Branch, Skip Event!!!");
- return kFALSE;
- }
-
- if ((!fEventAOD->GetPrimaryVertex()) || (TMath::Abs(fEventAOD->GetMagneticField())<1e-3)) return kFALSE;
-
- Bool_t bIsDmeson = kTRUE;
- if (fDzeroClArr) bIsDmeson = ((fDzeroClArr->GetEntriesFast()==0) && bIsDmeson);
- if (fDstarClArr) bIsDmeson = ((fDstarClArr->GetEntriesFast()==0) && bIsDmeson);
- if (bIsDmeson) return kFALSE;
-
- return kTRUE;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::CreateDzeroHistograms()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::CreateDzeroHistograms
-//
-
- if (!fListDzeroHistos) return;
- Bool_t bStatusTmpH = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
- TH1D *h1D = 0;
- h1D = new TH1D("hDzero_candi", "", 100, -0.5, 99.5); h1D->Sumw2(); fListDzeroHistos->Add(h1D); h1D=0;
- h1D = new TH1D("hDzero_seled", "", 100, -0.5, 99.5); h1D->Sumw2(); fListDzeroHistos->Add(h1D); h1D=0;
-
- TH2D *h2D = 0;
- Double_t dMass = TDatabasePDG::Instance()->GetParticle(421)->Mass();
- h2D = new TH2D("hDzero_Mass_Pt", "", 200, dMass-0.15, dMass+0.15, 100, 0., 50.);
- h2D->Sumw2(); fListDzeroHistos->Add(h2D); h2D=0;
-
- TH1::AddDirectory(bStatusTmpH);
- return;
-}
-
-//_____________________________________________________________________________
-void AliAnalysisTaskSEDmesonsFilterCJ::CreateDstarHistograms()
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::CreateDstarHistograms
-//
-
- if (!fListDstarHistos) return;
- Bool_t bStatusTmpH = TH1::AddDirectoryStatus();
- TH1::AddDirectory(kFALSE);
-
- TH1D *h1D = 0;
- h1D = new TH1D("hDstar_candi", "", 100, -0.5, 99.5); h1D->Sumw2(); fListDstarHistos->Add(h1D); h1D=0;
- h1D = new TH1D("hDstar_seled", "", 100, -0.5, 99.5); h1D->Sumw2(); fListDstarHistos->Add(h1D); h1D=0;
- h1D = new TH1D("hDstar_SoftPiPt", "", 500, 0., 10.0); h1D->Sumw2(); fListDstarHistos->Add(h1D); h1D=0;
-
- TH2D *h2D = 0;
- Double_t dMassDzero = TDatabasePDG::Instance()->GetParticle(421)->Mass();
- Double_t dMassDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();
- Double_t dMassDelta = dMassDstar - dMassDzero;
- h2D = new TH2D("hDstar_deltaMass_Pt", "", 200, dMassDelta-0.015, dMassDelta+0.015, 100, 0., 50.);
- h2D->Sumw2(); fListDstarHistos->Add(h2D); h2D=0;
-
- TH1::AddDirectory(bStatusTmpH);
- return;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskSEDmesonsFilterCJ::SetCutDzero(AliRDHFCutsD0toKpi *cut)
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::SetCutDzero
-//
-
- if (!cut) return kTRUE;
-
- if (fCutDzero) { delete fCutDzero; fCutDzero=0; }
- fCutDzero = new AliRDHFCutsD0toKpi(*cut);
- if (!fCutDzero) { AliError("No Dzero Cuts!!!"); return kTRUE; }
-
- return kFALSE;
-}
-
-//_____________________________________________________________________________
-Bool_t AliAnalysisTaskSEDmesonsFilterCJ::SetCutDstar(AliRDHFCutsDStartoKpipi *cut)
-{
-//
-// AliAnalysisTaskSEDmesonsFilterCJ::SetCutDstar
-//
-
- if (!cut) return kTRUE;
-
- if (fCutDstar) { delete fCutDstar; fCutDstar=0; }
- fCutDstar = new AliRDHFCutsDStartoKpipi(*cut);
- if (!fCutDstar) { AliError("No Dstar Cuts!!!"); return kTRUE; }
-
- return kFALSE;
-}
+// $Id$\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: The ALICE Off-line Project. *\r
+ * Contributors are mentioned in the code where appropriate. *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+//\r
+//\r
+// Task for selecting D mesons to be used as an input for D-Jet correlations\r
+//\r
+//-----------------------------------------------------------------------\r
+// Authors:\r
+// C. Bianchin (Utrecht University) chiara.bianchin@cern.ch\r
+// A.Grelli (Utrecht University) a.grelli@uu.nl\r
+// Xiaoming Zhang (LBNL) XMZhang@lbl.gov\r
+//-----------------------------------------------------------------------\r
+\r
+#include <TDatabasePDG.h>\r
+#include <TParticle.h>\r
+#include <TVector3.h>\r
+#include "TROOT.h"\r
+#include <TH3F.h>\r
+\r
+#include "AliRDHFCutsDStartoKpipi.h"\r
+#include "AliRDHFCutsD0toKpi.h"\r
+#include "AliAODMCHeader.h"\r
+#include "AliAODHandler.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliLog.h"\r
+#include "AliAODVertex.h"\r
+#include "AliAODRecoDecay.h"\r
+#include "AliAODRecoCascadeHF.h"\r
+#include "AliAODRecoDecayHF2Prong.h"\r
+#include "AliESDtrack.h"\r
+#include "AliAODMCParticle.h"\r
+#include "AliAnalysisTaskSEDmesonsFilterCJ.h"\r
+\r
+ClassImp(AliAnalysisTaskSEDmesonsFilterCJ)\r
+\r
+//_____________________________________________________________________________\r
+AliAnalysisTaskSEDmesonsFilterCJ::AliAnalysisTaskSEDmesonsFilterCJ() :\r
+ AliAnalysisTaskSE(),\r
+ fUseMCInfo(kFALSE),\r
+ fCandidateType(0),\r
+ fCandidateName(""),\r
+ fPDGmother(0),\r
+ fNProngs(0),\r
+ fBranchName(""),\r
+ fOutput(0),\r
+//fOutputCandidates(0),\r
+ fCuts(0),\r
+ fMinMass(0.),\r
+ fMaxMass(0.),\r
+ fCandidateArray(0)\r
+//fIsSelectedArray(0)\r
+{\r
+//\r
+// Default constructor\r
+//\r
+\r
+ for (Int_t i=4; i--;) fPDGdaughters[i] = 0;\r
+ for (Int_t i=30; i--;) fSigmaD0[i] = 0.;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliAnalysisTaskSEDmesonsFilterCJ::AliAnalysisTaskSEDmesonsFilterCJ(const char *name, AliRDHFCuts *cuts, ECandidateType candtype) :\r
+ AliAnalysisTaskSE(name),\r
+ fUseMCInfo(kFALSE),\r
+ fCandidateType(candtype),\r
+ fCandidateName(""),\r
+ fPDGmother(0),\r
+ fNProngs(0),\r
+//fPDGdaughters(),\r
+ fBranchName(""),\r
+ fOutput(0),\r
+//fOutputCandidates(0),\r
+ fCuts(cuts),\r
+ fMinMass(0.),\r
+ fMaxMass(0.),\r
+ fCandidateArray(0)\r
+//fIsSelectedArray(0)\r
+{\r
+//\r
+// Constructor. Initialization of Inputs and Outputs\r
+//\r
+ \r
+ Info("AliAnalysisTaskSEDmesonsFilterCJ","Calling Constructor");\r
+\r
+ for (Int_t i=4; i--;) fPDGdaughters[i] = 0;\r
+ for (Int_t i=30; i--;) fSigmaD0[i] = 0.;\r
+\r
+ const Int_t nptbins = fCuts->GetNPtBins();\r
+ Float_t defaultSigmaD013[13] = { 0.012, 0.012, 0.012, 0.015, 0.015, 0.018, 0.018, 0.020, 0.020, 0.030, 0.030, 0.037, 0.040 };\r
+\r
+ switch (fCandidateType) {\r
+ case 0 :\r
+ fCandidateName = "D0";\r
+ fPDGmother = 421;\r
+ fNProngs = 2;\r
+ fPDGdaughters[0] = 211; // pi \r
+ fPDGdaughters[1] = 321; // K\r
+ fPDGdaughters[2] = 0; // empty\r
+ fPDGdaughters[3] = 0; // empty\r
+ fBranchName = "D0toKpi";\r
+ break;\r
+ case 1 :\r
+ fCandidateName = "Dstar";\r
+ fPDGmother = 413;\r
+ fNProngs = 3;\r
+ fPDGdaughters[1] = 211; // pi soft\r
+ fPDGdaughters[0] = 421; // D0\r
+ fPDGdaughters[2] = 211; // pi fromD0\r
+ fPDGdaughters[3] = 321; // K from D0\r
+ fBranchName = "Dstar";\r
+\r
+ if (nptbins<=13) {\r
+ for (Int_t ipt=0; ipt<nptbins;ipt++) fSigmaD0[ipt] = defaultSigmaD013[ipt];\r
+ } else {\r
+ AliFatal(Form("Default sigma D0 not enough for %d pt bins, use SetSigmaD0ForDStar to set them",nptbins));\r
+ }\r
+ break;\r
+ default :\r
+ printf("%d not accepted!!\n",fCandidateType);\r
+ break;\r
+ }\r
+\r
+ if (fCandidateType==kD0toKpi) SetMassLimits(0.15, fPDGmother);\r
+ if (fCandidateType==kDstartoKpipi) SetMassLimits(0.015, fPDGmother);\r
+\r
+ DefineOutput(1, TList::Class()); // histos\r
+ DefineOutput(2, AliRDHFCuts::Class()); // my cuts\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AliAnalysisTaskSEDmesonsFilterCJ::~AliAnalysisTaskSEDmesonsFilterCJ()\r
+{\r
+//\r
+// destructor\r
+//\r
+\r
+ Info("~AliAnalysisTaskSEDmesonsFilterCJ","Calling Destructor"); \r
+ \r
+ if (fOutput) { delete fOutput; fOutput = 0; }\r
+ if (fCuts) { delete fCuts; fCuts = 0; }\r
+ if (fCandidateArray) { delete fCandidateArray; fCandidateArray = 0; }\r
+//if (fIsSelectedArray) { delete fIsSelectedArray; fIsSelectedArray = 0; }\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::Init()\r
+{\r
+//\r
+// Initialization\r
+//\r
+\r
+ if(fDebug>1) printf("AnalysisTaskSEDmesonsForJetCorrelations::Init() \n");\r
+\r
+ switch (fCandidateType) {\r
+ case 0: {\r
+ AliRDHFCutsD0toKpi* copyfCutsDzero = new AliRDHFCutsD0toKpi(*(static_cast<AliRDHFCutsD0toKpi*>(fCuts)));\r
+ copyfCutsDzero->SetName("AnalysisCutsDzero");\r
+ PostData(2, copyfCutsDzero); // Post the data\r
+ } break;\r
+ case 1: {\r
+ AliRDHFCutsDStartoKpipi* copyfCutsDstar = new AliRDHFCutsDStartoKpipi(*(static_cast<AliRDHFCutsDStartoKpipi*>(fCuts)));\r
+ copyfCutsDstar->SetName("AnalysisCutsDStar");\r
+ PostData(2, copyfCutsDstar); // Post the cuts\r
+ } break;\r
+ default: return;\r
+ }\r
+\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::UserCreateOutputObjects()\r
+{ \r
+//\r
+// output \r
+//\r
+\r
+ Info("UserCreateOutputObjects","CreateOutputObjects of task %s\n", GetName());\r
+ \r
+//OpenFile(1);\r
+ fOutput = new TList(); fOutput->SetOwner();\r
+ DefineHistoForAnalysis(); // define histograms\r
+ \r
+//fOutputCandidates= new TList();\r
+//fOutputCandidates->SetOwner();\r
+\r
+ fCandidateArray = new TClonesArray("AliAODRecoDecayHF",0);\r
+ fCandidateArray->SetName(Form("fCandidateArray%s",fCandidateName.Data()));\r
+\r
+//fOutputCandidates->Add(fCandidateArray);\r
+//fIsSelectedArray = new TClonesArray("Int_t",0);\r
+//fIsSelectedArray->SetName(Form("fIsSelectedArray%s",suffix.Data()));\r
+//fOutputCandidates->Add(fIsSelectedArray);\r
+\r
+ PostData(1, fOutput);\r
+//PostData(3, fOutputCandidates);\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::UserExec(Option_t *)\r
+{\r
+//\r
+// user exec\r
+//\r
+\r
+ // add cadidate branch\r
+ fCandidateArray->Delete();\r
+ if (!(InputEvent()->FindListObject(Form("fCandidateArray%s",fCandidateName.Data())))) InputEvent()->AddObject(fCandidateArray);\r
+\r
+ // Load the event\r
+ AliAODEvent *aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);\r
+ \r
+ TClonesArray *arrayDStartoD0pi = 0;\r
+ if (!aodEvent && AODEvent() && IsStandardAOD()) {\r
+\r
+ // In case there is an AOD handler writing a standard AOD, use the AOD \r
+ // event in memory rather than the input (ESD) event. \r
+ aodEvent = dynamic_cast<AliAODEvent*>(AODEvent());\r
+\r
+ // in this case the braches in the deltaAOD (AliAOD.VertexingHF.root)\r
+ // have to taken from the AOD event hold by the AliAODExtension\r
+ AliAODHandler *aodHandler = (AliAODHandler*)((AliAnalysisManager::GetAnalysisManager())->GetOutputEventHandler());\r
+ if(aodHandler->GetExtensions()) {\r
+ AliAODExtension *ext = (AliAODExtension*)aodHandler->GetExtensions()->FindObject("AliAOD.VertexingHF.root");\r
+ AliAODEvent *aodFromExt = ext->GetAOD();\r
+ arrayDStartoD0pi = (TClonesArray*)aodFromExt->GetList()->FindObject(fBranchName.Data());\r
+ }\r
+ } else {\r
+ arrayDStartoD0pi = (TClonesArray*)aodEvent->GetList()->FindObject(fBranchName.Data());\r
+ }\r
+\r
+ if (!arrayDStartoD0pi) {\r
+ AliInfo(Form("Could not find array %s, skipping the event",fBranchName.Data()));\r
+ return;\r
+ } else {\r
+ AliDebug(2, Form("Found %d vertices",arrayDStartoD0pi->GetEntriesFast())); \r
+ }\r
+ \r
+ TClonesArray* mcArray = 0x0;\r
+ if (fUseMCInfo) {\r
+ mcArray = dynamic_cast<TClonesArray*>(aodEvent->FindListObject(AliAODMCParticle::StdBranchName()));\r
+ if (!mcArray) {\r
+ printf("AliAnalysisTaskSEDStarSpectra::UserExec: MC particles not found!\n");\r
+ return;\r
+ }\r
+ }\r
+\r
+ //Histograms\r
+ TH1I* hstat = (TH1I*)fOutput->FindObject("hstat");\r
+ TH2F* hInvMassptD = (TH2F*)fOutput->FindObject("hInvMassptD");\r
+\r
+ TH1F* hPtPion=0x0;\r
+ if (fCandidateType==kDstartoKpipi) hPtPion = (TH1F*)fOutput->FindObject("hPtPion");\r
+ hstat->Fill(0);\r
+\r
+ // fix for temporary bug in ESDfilter \r
+ // the AODs with null vertex pointer didn't pass the PhysSel\r
+ if(!aodEvent->GetPrimaryVertex() || TMath::Abs(aodEvent->GetMagneticField())<0.001) return;\r
+\r
+ //Event selection\r
+ Bool_t iseventselected=fCuts->IsEventSelected(aodEvent);\r
+//TString firedTriggerClasses=((AliAODEvent*)aodEvent)->GetFiredTriggerClasses();\r
+ if (!iseventselected) return;\r
+ hstat->Fill(1);\r
+\r
+ const Int_t nD = arrayDStartoD0pi->GetEntriesFast();\r
+ hstat->Fill(2, nD);\r
+\r
+//Int_t icountReco = 0; // counters for efficiencies\r
+\r
+ //D* and D0 prongs needed to MatchToMC method\r
+ Int_t pdgDgDStartoD0pi[2] = { 421, 211 }; // D0,pi\r
+ Int_t pdgDgD0toKpi[2] = { 321, 211 }; // K, pi\r
+\r
+ //D0 from D0 bar\r
+ Int_t pdgdaughtersD0[2] = { 211, 321 }; // pi,K \r
+ Int_t pdgdaughtersD0bar[2] = { 321, 211 }; // K,pi \r
+\r
+ Int_t iCand =0;\r
+ Int_t isSelected = 0;\r
+ AliAODRecoDecayHF *charmCand = 0;\r
+\r
+ Int_t mcLabel = -9999;\r
+ Int_t pdgMeson = 413;\r
+ if (fCandidateType==kD0toKpi) pdgMeson = 421;\r
+\r
+ for (Int_t icharm=0; icharm<nD; icharm++) { //loop over D candidates\r
+ charmCand = (AliAODRecoDecayHF*)arrayDStartoD0pi->At(icharm); // D candidates\r
+ if (!charmCand) return;\r
+\r
+\r
+ if (fUseMCInfo) { // Look in MC, try to simulate the z\r
+ if (fCandidateType==kDstartoKpipi) {\r
+ AliAODRecoCascadeHF *temp = (AliAODRecoCascadeHF*)charmCand;\r
+ mcLabel = temp->MatchToMC(413,421,pdgDgDStartoD0pi,pdgDgD0toKpi,mcArray);\r
+ }\r
+\r
+ if (fCandidateType==kD0toKpi) mcLabel = charmCand->MatchToMC(421,mcArray,fNProngs,fPDGdaughters);\r
+\r
+ if (mcLabel<=0) continue; \r
+ }\r
+\r
+ // region of interest + cuts\r
+ if (!fCuts->IsInFiducialAcceptance(charmCand->Pt(),charmCand->Y(pdgMeson))) continue; \r
+\r
+ isSelected = fCuts->IsSelected(charmCand,AliRDHFCuts::kAll,aodEvent); //selected\r
+ if (!isSelected) continue;\r
+\r
+ new ((*fCandidateArray)[iCand]) AliAODRecoDecayHF(*charmCand);\r
+// new ((*fIsSelectedArray)[iCand]) Int_t(isSelected);\r
+ iCand++;\r
+\r
+ Double_t masses[2];\r
+ Double_t ptD = charmCand->Pt();\r
+ if (fCandidateType==kDstartoKpipi) { //D*->D0pi->Kpipi\r
+\r
+ //softpion from D* decay\r
+ AliAODRecoCascadeHF *temp = (AliAODRecoCascadeHF*)charmCand;\r
+ AliAODTrack *track2 = (AliAODTrack*)temp->GetBachelor(); \r
+ hstat->Fill(3);\r
+\r
+ // select D* in the D0 window.\r
+ // In the cut object window is loose to allow side bands\r
+ Double_t mPDGD0 = TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
+\r
+ // retrieve the corresponding pt bin for the candidate\r
+ // and set the expected D0 width (x3)\r
+// static const Int_t n = fCuts->GetNPtBins();\r
+ Int_t bin = fCuts->PtBin(ptD);\r
+ if(bin<0 || bin>=fCuts->GetNPtBins()) {\r
+ AliError(Form("Pt %.3f out of bounds",ptD));\r
+ continue;\r
+ }\r
+\r
+ AliInfo(Form("Pt bin %d and sigma D0 %.4f",bin,fSigmaD0[bin]));\r
+ if ((temp->InvMassD0()>=(mPDGD0-3.*fSigmaD0[bin])) && (temp->InvMassD0()<=(mPDGD0+3.*fSigmaD0[bin]))) { \r
+ masses[0] = temp->DeltaInvMass(); //D*\r
+ masses[1] = 0.; //dummy for D*\r
+\r
+ //D* delta mass\r
+ hInvMassptD->Fill(masses[0], ptD); // 2 D slice for pt bins\r
+\r
+ // D* pt and soft pion pt for good candidates \r
+ Double_t mPDGDstar = TDatabasePDG::Instance()->GetParticle(413)->Mass();\r
+ Double_t invmassDelta = temp->DeltaInvMass();\r
+ if (TMath::Abs(invmassDelta-(mPDGDstar-mPDGD0))<0.0021) hPtPion->Fill(track2->Pt());\r
+ }\r
+ } //Dstar specific\r
+\r
+ if (fCandidateType==kD0toKpi) { //D0->Kpi\r
+\r
+ //needed quantities\r
+ masses[0] = charmCand->InvMass(fNProngs, (UInt_t*)pdgdaughtersD0); //D0\r
+ masses[1] = charmCand->InvMass(fNProngs, (UInt_t*)pdgdaughtersD0bar); //D0bar\r
+ hstat->Fill(3);\r
+\r
+ // mass vs pt\r
+ if (isSelected==1 || isSelected==3) hInvMassptD->Fill(masses[0],ptD);\r
+ if (isSelected>=2) hInvMassptD->Fill(masses[1],ptD);\r
+ } //D0 specific\r
+\r
+ charmCand = 0; \r
+ } // end of D cand loop\r
+\r
+//AliDebug(2, Form("Found %i Reco particles that are D*!!",icountReco));\r
+\r
+//PostData(1,fOutput);\r
+//PostData(3,fOutputCandidates);\r
+\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::Terminate(Option_t*)\r
+{\r
+// The Terminate() function is the last function to be called during\r
+// a query. It always runs on the client, it can be used to present\r
+// the results graphically or save the results to file.\r
+ \r
+ Info("Terminate"," terminate");\r
+ AliAnalysisTaskSE::Terminate();\r
+\r
+ fOutput = dynamic_cast<TList*>(GetOutputData(1));\r
+ if (!fOutput) {\r
+ printf("ERROR: fOutput not available\n");\r
+ return;\r
+ }\r
+\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::SetMassLimits(Double_t range, Int_t pdg)\r
+{\r
+//\r
+// AliAnalysisTaskSEDmesonsFilterCJ::SetMassLimits\r
+//\r
+\r
+ Float_t mass = TDatabasePDG::Instance()->GetParticle(TMath::Abs(pdg))->Mass();\r
+\r
+ // compute the Delta mass for the D*\r
+ if (fCandidateType==kDstartoKpipi) mass -= TDatabasePDG::Instance()->GetParticle(421)->Mass();\r
+//cout << "mass ---------------" << endl;\r
+\r
+ fMinMass = mass - range;\r
+ fMaxMass = mass + range;\r
+\r
+ AliInfo(Form("Setting mass limits to %f, %f",fMinMass,fMaxMass));\r
+ if ((fMinMass<0.) || (fMaxMass<=0.) || (fMaxMass<fMinMass)) AliFatal("Wrong mass limits!\n");\r
+\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AliAnalysisTaskSEDmesonsFilterCJ::SetMassLimits(Double_t lowlimit, Double_t uplimit)\r
+{\r
+//\r
+// AliAnalysisTaskSEDmesonsFilterCJ::SetMassLimits\r
+//\r
+\r
+ if (uplimit>lowlimit) {\r
+ fMinMass = lowlimit;\r
+ fMaxMass = uplimit;\r
+ } else {\r
+ printf("Error! Lower limit larger than upper limit!\n");\r
+ fMinMass = uplimit - uplimit*0.2;\r
+ fMaxMass = uplimit;\r
+ }\r
+\r
+ return;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Bool_t AliAnalysisTaskSEDmesonsFilterCJ::SetD0WidthForDStar(Int_t nptbins, Float_t *width)\r
+{\r
+//\r
+// AliAnalysisTaskSEDmesonsFilterCJ::SetD0WidthForDStar\r
+//\r
+\r
+ if (nptbins>30) {\r
+ AliInfo("Maximum number of bins allowed is 30!");\r
+ return kFALSE;\r
+ }\r
+\r
+ if (!width) return kFALSE;\r
+ for (Int_t ipt=0; ipt<nptbins; ipt++) fSigmaD0[ipt]=width[ipt];\r
+\r
+ return kTRUE;\r
+}\r
+\r
+\r
+//_____________________________________________________________________________\r
+Bool_t AliAnalysisTaskSEDmesonsFilterCJ::DefineHistoForAnalysis()\r
+{\r
+//\r
+// AliAnalysisTaskSEDmesonsFilterCJ::DefineHistoForAnalysis\r
+//\r
+\r
+ // Statistics \r
+ TH1I* hstat = new TH1I("hstat","Statistics",6,-0.5,5.5);\r
+ hstat->GetXaxis()->SetBinLabel(1, "N ev anal");\r
+ hstat->GetXaxis()->SetBinLabel(2, "N ev sel");\r
+ hstat->GetXaxis()->SetBinLabel(3, "N cand");\r
+ hstat->GetXaxis()->SetBinLabel(4, "N cand sel cuts");\r
+ hstat->GetXaxis()->SetBinLabel(5, "N jets");\r
+ hstat->GetXaxis()->SetBinLabel(6, "N cand in jet");\r
+/*if(fUseMCInfo) {\r
+ hstat->GetXaxis()->SetBinLabel(7,"N D");\r
+ hstat->GetXaxis()->SetBinLabel(8,"N D in jet");\r
+ }*/\r
+ hstat->SetNdivisions(1);\r
+ fOutput->Add(hstat);\r
+\r
+ // Invariant mass related histograms\r
+ const Int_t nbinsmass = 200;\r
+ TH2F *hInvMass = new TH2F("hInvMassptD", "D invariant mass distribution", nbinsmass, fMinMass, fMaxMass, 100, 0., 50.);\r
+ hInvMass->SetStats(kTRUE);\r
+ hInvMass->GetXaxis()->SetTitle("mass (GeV/c)");\r
+ hInvMass->GetYaxis()->SetTitle("p_{T} (GeV/c)");\r
+ fOutput->Add(hInvMass);\r
+ \r
+ if (fCandidateType==kDstartoKpipi) {\r
+ TH1F* hPtPion = new TH1F("hPtPion", "Primary pions candidates pt", 500, 0., 10.);\r
+ hPtPion->SetStats(kTRUE);\r
+ hPtPion->GetXaxis()->SetTitle("p_{T} (GeV/c)");\r
+ hPtPion->GetYaxis()->SetTitle("entries");\r
+ fOutput->Add(hPtPion);\r
+ }\r
+\r
+ return kTRUE; \r
+}\r
-#ifndef ALIANALYSISTASKSEDMESONSFILTERCJ_H
-#define ALIANALYSISTASKSEDMESONSFILTERCJ_H
-
-// $Id$
-
-/* Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice */
-
-//*************************************************************************
-// Class AliAnalysisTaskSEDmesonsCJ
-// AliAnalysisTaskSE for Dmesons - jet correlations analysis
-// Author: Xiaoming Zhang, xmzhang@lbl.gov
-//*************************************************************************
-
-#include "AliAnalysisTaskSE.h"
-
-class TList;
-class AliRDHFCutsD0toKpi;
-class AliRDHFCutsDStartoKpipi;
-
-class AliAnalysisTaskSEDmesonsFilterCJ : public AliAnalysisTaskSE {
-
- public :
-
- AliAnalysisTaskSEDmesonsFilterCJ();
- AliAnalysisTaskSEDmesonsFilterCJ(const char *name);
- virtual ~AliAnalysisTaskSEDmesonsFilterCJ();
-
- virtual void Init();
- virtual void LocalInit() { Init(); }
- virtual void UserCreateOutputObjects();
- virtual void UserExec(Option_t *opt);
- virtual void Terminate(Option_t *opt);
- virtual void NotifyRun();
-
- Bool_t SetCutDzero(AliRDHFCutsD0toKpi *cut);
- Bool_t SetCutDstar(AliRDHFCutsDStartoKpipi *cut);
-
- private :
-
- AliAnalysisTaskSEDmesonsFilterCJ(const AliAnalysisTaskSEDmesonsFilterCJ &);
- AliAnalysisTaskSEDmesonsFilterCJ& operator=(const AliAnalysisTaskSEDmesonsFilterCJ &);
-
- Bool_t ExecOnce();
- Bool_t ExecEach();
-
- void ExecAnas();
- void ExecAnasDzero();
- void ExecAnasDstar();
-
- void CreateDzeroHistograms();
- void CreateDstarHistograms();
-//=============================================================================
-
- AliAODEvent *fEventAOD; //! in put AOD event
-
- AliRDHFCutsD0toKpi *fCutDzero; //! Dzero selection cut
- AliRDHFCutsDStartoKpipi *fCutDstar; //! Dstar selection cut
-
- TClonesArray *fDzeroClArr; //! intput Dzero candidates array
- TClonesArray *fDstarClArr; //! intput Dstar candidates array
-
- TClonesArray *fUsedDzeros; //! selected Dzero candidates array
- TClonesArray *fUsedD0bars; //! selected D0bar candidates array
- TClonesArray *fUsedDstars; //! selected Dstar candidates array
-
- Bool_t fIsNotExecOnce; // flag for ExecOnce
-
- TList *fListCutsDmesons; //! list of Dmeson selection cuts
- TList *fListDzeroHistos; //! list of Dzero output control histrograms
- TList *fListDstarHistos; //! list of Dstar output control histrograms
-
- ClassDef(AliAnalysisTaskSEDmesonsFilterCJ, 1);
-};
-
-#endif
+#ifndef ALIANALYSISTASKSEDMESONSFILTERCJ_H\r
+#define ALIANALYSISTASKSEDMESONSFILTERCJ_H\r
+// $Id$\r
+/**************************************************************************\r
+ * Copyright(c) 1998-2009, ALICE Experiment at CERN, All rights reserved. *\r
+ * *\r
+ * Author: The ALICE Off-line Project. *\r
+ * Contributors are mentioned in the code where appropriate. *\r
+ * *\r
+ * Permission to use, copy, modify and distribute this software and its *\r
+ * documentation strictly for non-commercial purposes is hereby granted *\r
+ * without fee, provided that the above copyright notice appears in all *\r
+ * copies and that both the copyright notice and this permission notice *\r
+ * appear in the supporting documentation. The authors make no claims *\r
+ * about the suitability of this software for any purpose. It is *\r
+ * provided "as is" without express or implied warranty. *\r
+ **************************************************************************/\r
+\r
+//-----------------------------------------------------------------------\r
+// Author : A. Grelli, Utrecht University\r
+// C. Bianchin, Utrecht University\r
+// X. Zhang, LBNL\r
+//-----------------------------------------------------------------------\r
+\r
+\r
+#include <TH2F.h>\r
+#include "AliAODEvent.h"\r
+#include "AliPicoTrack.h"\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+class TH3F;\r
+class TString;\r
+class TParticle ;\r
+class TClonesArray ;\r
+class AliMCParticle;\r
+class AliAODMCParticle;\r
+class AliRDHFCuts;\r
+class AliAODRecoCascadeHF;\r
+\r
+class AliAnalysisTaskSEDmesonsFilterCJ : public AliAnalysisTaskSE \r
+{\r
+\r
+ public :\r
+\r
+ enum ECandidateType{ kD0toKpi, kDstartoKpipi };\r
+ \r
+ AliAnalysisTaskSEDmesonsFilterCJ();\r
+ AliAnalysisTaskSEDmesonsFilterCJ(const Char_t* name,AliRDHFCuts* cuts,ECandidateType candtype);\r
+ virtual ~AliAnalysisTaskSEDmesonsFilterCJ();\r
+\r
+ virtual void UserCreateOutputObjects();\r
+ virtual void UserExec(Option_t *option);\r
+ virtual void Terminate(Option_t *);\r
+ virtual void Init();\r
+ virtual void LocalInit() { Init(); }\r
+\r
+ // inizializations\r
+ Bool_t DefineHistoForAnalysis();\r
+\r
+ // set MC usage\r
+ void SetMC(Bool_t theMCon) { fUseMCInfo = theMCon; }\r
+ Bool_t GetMC() const { return fUseMCInfo; }\r
+ \r
+ void SetMassLimits(Double_t range, Int_t pdg);\r
+ void SetMassLimits(Double_t lowlimit, Double_t uplimit);\r
+\r
+ // Array of D0 width for the Dstar\r
+ Bool_t SetD0WidthForDStar(Int_t nptbins, Float_t *width);\r
+\r
+ private :\r
+ \r
+ AliAnalysisTaskSEDmesonsFilterCJ(const AliAnalysisTaskSEDmesonsFilterCJ &source);\r
+ AliAnalysisTaskSEDmesonsFilterCJ& operator=(const AliAnalysisTaskSEDmesonsFilterCJ& source); \r
+\r
+ Bool_t fUseMCInfo; // Use MC info\r
+\r
+ UInt_t fCandidateType; // Dstar or D0\r
+ TString fCandidateName; // Dstar or D0\r
+\r
+ Int_t fPDGmother; // PDG code of D meson\r
+ Int_t fNProngs; // number of prong of the decay channel \r
+ Int_t fPDGdaughters[4]; // PDG codes of daughters\r
+ Float_t fSigmaD0[30]; // D0 sigma for Dstar\r
+\r
+ TString fBranchName; // AOD branch name\r
+ TList *fOutput; //! user output\r
+//TList *fOutputCandidates; //! output of array of candidates (kExchange)\r
+\r
+ AliRDHFCuts *fCuts; //! Cuts \r
+ Double_t fMinMass; // mass lower limit histogram\r
+ Double_t fMaxMass; // mass upper limit histogram\r
+\r
+ TClonesArray *fCandidateArray; //! contains candidates selected by AliRDHFCuts\r
+//TClonesArray *fIsSelectedArray; //! contains result of IsSelected for candidates which pass the cuts (needed for D0)\r
+\r
+ ClassDef(AliAnalysisTaskSEDmesonsFilterCJ,1); // class for charm-jet correlations\r
+};\r
+\r
+#endif\r
-// $Id$
-
-AliAnalysisTaskFlavourJetCorrelations *AddTaskFlavourJetCorrelations(
- TString sUsedTrks = "",
- TString sUsedClus = "",
- TString sUsedJets = "",
- TString sUsedRho = "",
- TString sTaskName = "AliAnalysisTaskFlavourJetCorrelations",
- Bool_t bIsMakeHisto = kTRUE)
-{
- AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
-
- if (!mgr) {
- ::Error("AddTaskFlavourJetCorrelations.C::AddTaskFlavourJetCorrelations", "No analysis manager to connect to.");
- return NULL;
- }
-
- if (!mgr->GetInputEventHandler()) {
- ::Error("AddTaskFlavourJetCorrelations.C::AddTaskFlavourJetCorrelations", "This task requires an input event handler");
- return NULL;
- }
-
- TString type = mgr->GetInputEventHandler()->GetDataType();
- if (!type.Contains("ESD") && !type.Contains("AOD")) {
- ::Error("AddTaskFlavourJetCorrelations.C::AddTaskFlavourJetCorrelations", "Task manager to have an ESD or AOD input handler.");
- return NULL;
- }
-//=============================================================================
-
- AliAnalysisTaskFlavourJetCorrelations *taskFlavourCJ = new AliAnalysisTaskFlavourJetCorrelations(sTaskName.Data(),bIsMakeHisto);
- taskFlavourCJ->SetClusName(sUsedClus.Data());
- taskFlavourCJ->SetTracksName(sUsedTrks.Data());
- taskFlavourCJ->SetJetsName(sUsedJets.Data());
- taskFlavourCJ->SetRhoName(sUsedRho.Data());
-
- mgr->AddTask(taskFlavourCJ);
- mgr->ConnectInput( taskFlavourCJ, 0, mgr->GetCommonInputContainer());
-//mgr->ConnectOutput(taskFlavourCJ, 0, mgr->GetCommonOutputContainer());
-
- TString sOutput0 = sTaskName + "_" + sUsedJets;
- TString sOutput1 = sOutput0 + "_GeneralHistograms";
- TString sOutput2 = sOutput0 + "_ControlHistograms";
- TString sOutput3 = sOutput0 + "_AnDzeroHistograms";
- TString sOutput4 = sOutput0 + "_AnDstarHistograms";
- TString sCommon = AliAnalysisManager::GetCommonFileName();
- mgr->ConnectOutput(taskFlavourCJ, 1, mgr->CreateContainer(sOutput1.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,sCommon.Data()));
- mgr->ConnectOutput(taskFlavourCJ, 2, mgr->CreateContainer(sOutput2.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,sCommon.Data()));
- mgr->ConnectOutput(taskFlavourCJ, 3, mgr->CreateContainer(sOutput3.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,sCommon.Data()));
- mgr->ConnectOutput(taskFlavourCJ, 4, mgr->CreateContainer(sOutput4.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,sCommon.Data()));
-
- return taskFlavourCJ;
-}
+// $Id$\r
+\r
+AliAnalysisTaskFlavourJetCorrelations *AddTaskFlavourJetCorrelations(\r
+ AliAnalysisTaskRecoJetCorrelations::ECandidateType cand = AliAnalysisTaskRecoJetCorrelations::kDstartoKpipi,\r
+ TString filename = "DStartoKpipiCuts.root",\r
+ Bool_t theMCon = kFALSE,\r
+ TString jetArrname = "",\r
+ TString suffix = "",\r
+ Bool_t triggerOnLeadingJet = kFALSE,\r
+ Float_t R = 0.4,\r
+ Float_t jptcut = 10.,\r
+ Int_t acceptance = 1 /*1= 0-2pi, 2=emcal cut*/,\r
+ Double_t areaCut = 0.)\r
+{\r
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
+ if (!mgr) {\r
+ ::Error("AddTaskFlavourJetCorrelations::AddTaskFlavourJetCorrelations", "No analysis manager to connect to.");\r
+ return NULL;\r
+ } \r
+\r
+ Bool_t useStdC = kFALSE;\r
+ TFile* filecuts=TFile::Open(filename);\r
+ if (!filecuts || (filecuts && !filecuts->IsOpen())) {\r
+ cout<<"Input file not found: use std cuts"<<endl;\r
+ useStdC = kTRUE;\r
+ }\r
+\r
+ AliRDHFCuts *analysiscuts = 0x0;\r
+ switch (cand) {\r
+ case 0 :\r
+ if (useStdC) {\r
+ analysiscuts = new AliRDHFCutsD0toKpi();\r
+ analysiscuts->SetStandardCutsPP2010();\r
+ } else\r
+ analysiscuts = (AliRDHFCutsD0toKpi*)filecuts->Get("D0toKpiCuts");\r
+ break;\r
+ case 1 :\r
+ if(useStdC) {\r
+ analysiscuts = new AliRDHFCutsDStartoKpipi();\r
+ analysiscuts->SetStandardCutsPP2010();\r
+ } else\r
+ analysiscuts = (AliRDHFCutsDStartoKpipi*)filecuts->Get("DStartoKpipiCuts");\r
+ analysiscuts->SetName("DStartoKpipiCuts");\r
+ break;\r
+ }\r
+ \r
+ if (!analysiscuts) { // mm let's see if everything is ok\r
+ AliFatal("Specific AliRDHFCuts not found");\r
+ return;\r
+ }\r
+\r
+ printf("CREATE TASK\n"); //CREATE THE TASK\r
+\r
+ // create the task\r
+ AliAnalysisTaskFlavourJetCorrelations *task = new AliAnalysisTaskFlavourJetCorrelations("AnaTaskFlavourJetCorrelations", analysiscuts, cand, jetArrname);\r
+\r
+ //D meson settings\r
+ task->SetMC(theMCon);\r
+ task->SetTriggerOnLeadingJet(triggerOnLeadingJet);\r
+\r
+ //jet settings\r
+ task->SetJetRadius(R);\r
+ task->SetJetPtCut(jptcut);\r
+\r
+ Float_t etaCov=0.9;\r
+ if (acceptance==2) etaCov=0.7; //EMCal\r
+\r
+ Float_t minEta = -etaCov+R;\r
+ Float_t maxEta = etaCov-R;\r
+ if (acceptance) task->SetJetEtaLimits(minEta, maxEta);\r
+\r
+ Float_t minPhi = 0.;\r
+ Float_t maxPhi = 2.*TMath::Pi();\r
+//if (acceptance==2) { /*80-180 degree*/ }\r
+ if (acceptance) task->SetJetPhiLimits(minPhi, maxPhi);\r
+\r
+ //Float_t area=0.6*TMath::Pi()*R*R;\r
+ task->SetJetAreaCut(areaCut);\r
+ mgr->AddTask(task);\r
+\r
+ // Create and connect containers for input/output\r
+ TString outputfile = AliAnalysisManager::GetCommonFileName();\r
+ outputfile += ":PWG3_D2H_DEmcalJet";\r
+ outputfile += suffix;\r
+\r
+ TString nameContainer0="histogramsCorrelations";\r
+ TString nameContainer1="cuts";\r
+\r
+ nameContainer0 += suffix;\r
+ nameContainer1 += suffix;\r
+\r
+ // ------ input data ------\r
+ AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer();\r
+\r
+ // ----- output data -----\r
+\r
+ // output TH1I for event counting\r
+ AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(nameContainer0, TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());\r
+ AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(nameContainer1, AliRDHFCuts::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data());\r
+\r
+ mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());\r
+ mgr->ConnectOutput(task,1,coutput1);\r
+ mgr->ConnectOutput(task,2,coutput2);\r
+\r
+ Printf("Input and Output connected to the manager");\r
+ return task ;\r
+}\r
-// $Id$
-
-AliAnalysisTaskSEDmesonsFilterCJ *AddTaskSEDmesonsFilterCJ(TString sCutDzero = "", TString sCutDstar = "")
-{
- AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
-
- if (!mgr) {
- ::Error("AddTaskSEDmesonsFilterCJ.C::AddTaskSEDmesonsFilterCJ", "No analysis manager to connect to.");
- return NULL;
- }
-
- TString type = mgr->GetInputEventHandler()->GetDataType();
- if (!type.Contains("ESD") && !type.Contains("AOD")) {
- ::Error("AddTaskSEDmesonsFilterCJ.C::AddTaskSEDmesonsFilterCJ", "Task manager to have an ESD or AOD input handler.");
- return NULL;
- }
-
- if (!mgr->GetInputEventHandler()) {
- ::Error("AddTaskSEEAddTaskSEDmesonsFilterCJ.C::AddTaskSEEAddTaskSEDmesonsFilterCJ", "This task requires an input event handler");
- return NULL;
- }
-//=============================================================================
-
- TFile *file = 0;
- AliRDHFCutsD0toKpi *cutDzero = 0;
- AliRDHFCutsDStartoKpipi *cutDstar = 0;
-
- if (!sCutDzero.IsNull()) {
- file = TFile::Open(sCutDzero.Data(), "READ");
- cutDzero = dynamic_cast<AliRDHFCutsD0toKpi*>(file->Get("D0toKpiCuts"));
- file->Close();
- }
-
- if (!sCutDstar.IsNull()) {
- file = TFile::Open(sCutDstar.Data(), "READ");
- cutDstar = dynamic_cast<AliRDHFCutsDStartoKpipi*>(file->Get("DStartoKpipiCuts"));
- file->Close();
- }
-
- AliAnalysisTaskSEDmesonsFilterCJ *taskDemsonFilterCJ = new AliAnalysisTaskSEDmesonsFilterCJ("AliAnalysisTaskSEDmesonsFilterCJ");
- taskDemsonFilterCJ->SetCutDzero(cutDzero);
- taskDemsonFilterCJ->SetCutDstar(cutDstar);
- mgr->AddTask(taskDemsonFilterCJ);
-
- mgr->ConnectInput( taskDemsonFilterCJ, 0, mgr->GetCommonInputContainer());
-//mgr->ConnectOutput(taskDemsonFilterCJ, 0, mgr->GetCommonOutputContainer());
-
- TString sCommon = AliAnalysisManager::GetCommonFileName();
- TString sCutName = "AnalysisCutsDmeson.root";
- TString sOutput1 = "AliAnalysisTaskSEDmesonsFilterCJ_ListDmesonCuts";
- TString sOutput2 = "AliAnalysisTaskSEDmesonsFilterCJ_Dzero_ControlHistograms";
- TString sOutput3 = "AliAnalysisTaskSEDmesonsFilterCJ_Dstar_ControlHistograms";
- mgr->ConnectOutput(taskDemsonFilterCJ, 1, mgr->CreateContainer(sOutput1.Data(),TList::Class(),AliAnalysisManager::kOutputContainer,sCutName.Data()));
- mgr->ConnectOutput(taskDemsonFilterCJ, 2, mgr->CreateContainer(sOutput2.Data(),TList::Class(),AliAnalysisManager::kOutputContainer, sCommon.Data()));
- mgr->ConnectOutput(taskDemsonFilterCJ, 3, mgr->CreateContainer(sOutput3.Data(),TList::Class(),AliAnalysisManager::kOutputContainer, sCommon.Data()));
- return taskDemsonFilterCJ;
-}
+// $Id$\r
+\r
+AliAnalysisTaskSEDmesonsFilterCJ *AddTaskSEDmesonsFilterCJ(\r
+ AliAnalysisTaskSEDmesonsForJetCorrelations::ECandidateType cand = AliAnalysisTaskSEDmesonsForJetCorrelations::kDstartoKpipi,\r
+ TString filename = "DStartoKpipiCuts.root",\r
+ Bool_t theMCon = kFALSE,\r
+ TString suffix = "")\r
+{\r
+ AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();\r
+ if (!mgr) {\r
+ ::Error("AddTaskSEDmesonsFilterCJ", "No analysis manager to connect to.");\r
+ return NULL;\r
+ } \r
+\r
+ Bool_t useStdC = kFALSE;\r
+ TFile* filecuts=TFile::Open(filename);\r
+ if(!filecuts || (filecuts && !filecuts->IsOpen())) {\r
+ cout<<"Input file not found: use std cuts"<<endl;\r
+ useStdC = kTRUE;\r
+ }\r
+\r
+ AliRDHFCuts *analysiscuts=0x0;\r
+ switch (cand) {\r
+ case 0 :\r
+ if(useStdC) {\r
+ analysiscuts = new AliRDHFCutsD0toKpi();\r
+ analysiscuts->SetStandardCutsPP2010();\r
+ } else\r
+ analysiscuts = (AliRDHFCutsD0toKpi*)filecuts->Get("D0toKpiCuts");\r
+ break;\r
+ case 1 :\r
+ if(useStdC) {\r
+ analysiscuts = new AliRDHFCutsDStartoKpipi();\r
+ analysiscuts->SetStandardCutsPP2010();\r
+ } else\r
+ analysiscuts = (AliRDHFCutsDStartoKpipi*)filecuts->Get("DStartoKpipiCuts");\r
+ analysiscuts->SetName("DStartoKpipiCuts");\r
+ break;\r
+ }\r
+ \r
+ if (!analysiscuts) { // mm let's see if everything is ok\r
+ AliFatal("Specific AliRDHFCuts not found");\r
+ return;\r
+ } \r
+\r
+ printf("CREATE TASK\n"); //CREATE THE TASK\r
+\r
+ // create the task\r
+ AliAnalysisTaskSEDmesonsFilterCJ *task = new AliAnalysisTaskSEDmesonsFilterCJ("AnaTaskSEDmesonsFilterCJ",analysiscuts,cand);\r
+ task->SetMC(theMCon); //D meson settings\r
+ mgr->AddTask(task);\r
+\r
+ // Create and connect containers for input/output\r
+ TString outputfile = AliAnalysisManager::GetCommonFileName();\r
+ outputfile += ":PWG3_D2H_DmesonsForJetCorrelations";\r
+ outputfile += suffix;\r
+\r
+ TString nameContainer0="histograms";\r
+ TString nameContainer1="cuts";\r
+//TString nameContainer2="DcandidatesSel";\r
+\r
+ nameContainer0 += suffix;\r
+ nameContainer1 += suffix;\r
+//nameContainer2 += suffix;\r
+\r
+ // ------ input data ------\r
+ AliAnalysisDataContainer *cinput0 = mgr->GetCommonInputContainer();\r
+ \r
+ // ----- output data -----\r
+ \r
+ // output TH1I for event counting\r
+ AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(nameContainer0, TList::Class(),AliAnalysisManager::kOutputContainer,outputfile.Data());\r
+ AliAnalysisDataContainer *coutput2 = mgr->CreateContainer(nameContainer1, AliRDHFCuts::Class(),AliAnalysisManager::kOutputContainer, outputfile.Data());\r
+//AliAnalysisDataContainer *coutput3 = mgr->CreateContainer(nameContainer2, TList::Class(),AliAnalysisManager::kExchangeContainer, outputfile.Data());\r
+ \r
+ mgr->ConnectInput(task,0,mgr->GetCommonInputContainer());\r
+ mgr->ConnectOutput(task,1,coutput1);\r
+ mgr->ConnectOutput(task,2,coutput2);\r
+//mgr->ConnectOutput(task,3,coutput3);\r
+\r
+ Printf("Input and Output connected to the manager");\r
+ return task ;\r
+}\r
+\r
// $Id$
-void AddTasksFlavourJet()
+void AddTasksFlavourJet(const Int_t iCandType = 1 /*0 = D0, 1=Dstar...*/,
+ const TString sCutFile = "cutsHF/D0toKpiCutsppRecVtxNoPileupRejNoEMCAL.root",
+ const Double_t dJetPtCut = 1.,
+ const Double_t dJetAreaCut = 0.,
+ const Int_t iAccCut = 1,
+ const TString sRunPeriod = "LHC10b",
+ const Int_t uBeamType = 0,
+ const UInt_t uTriggerMask = AliVEvent::kMB, /*for jets; the D mesons trigger is defined in the cut object*/
+ const Bool_t bIsMC = kFALSE,
+ TString sText=""/*completes the name of the candidate task lists*/
+)
{
- const UInt_t uTriggerMask = AliVEvent::kMB;
-
- const TString sCutDzero = "cutsHF/D0toKpiCutsppRecVtxNoPileupRejNoEMCAL.root";
- const TString sCutDstar = "cutsHF/DStartoKpipiCuts.root";
-
const TString sInputTrk = "tracks";
const TString sUsedTrks = "PicoTracks";
const TString sUsedClus = "";
- const TString sRunPeriod = "-1:-1";
const Int_t iJetAlgo = 1;
const Int_t iJetType = 1;
- const Int_t uBeamType = 0;
- const Double_t dJetPtCut = 0.;
- const Double_t dJetAreaCut = 0.;
-
const Int_t nRadius = 3;
const Double_t aRadius[] = { 0.2, 0.4, 0.6 };
const TString sRadius[] = { "R02", "R04", "R06" };
//=============================================================================
+
AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
if (!mgr) {
::Error("AddTasksFlavourJet.C::AddTasksFlavourJet", "No analysis manager to connect to.");
- return NULL;
+ return;
}
TString type = mgr->GetInputEventHandler()->GetDataType();
if (!type.Contains("ESD") && !type.Contains("AOD")) {
::Error("AddTasksFlavourJet.C::AddTasksFlavourJet", "Task manager to have an ESD or AOD input handler.");
- return NULL;
+ return;
}
if (!mgr->GetInputEventHandler()) {
::Error("AddTasksFlavourJet.C::AddTasksFlavourJet", "This task requires an input event handler");
- return NULL;
+ return;
}
//=============================================================================
UInt_t uAnaType = (((iJetType==0) || (iJetType==2)) ? 1 : 0);
Int_t iLeading = ((iJetType==0) ? 3 : ((iJetType==1) ? 0 : 1));
+ //D mesons -- PID
gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
AliAnalysisTaskSE *taskRespPID = AddTaskPIDResponse(bIsMC);
- taskRespPID->SelectCollisionCandidates(uTriggerMask);
- gROOT->LoadMacro("AddTaskSEDmesonsFilterCJ.C");
- AliAnalysisTaskSEDmesonsFilterCJ *taskDmesonsFilter = AddTaskSEDmesonsFilterCJ(sCutDzero,sCutDstar);
- taskDmesonsFilter->SelectCollisionCandidates(uTriggerMask);
+ // -- D meson selection
+ gROOT->LoadMacro("AliAnalysisTaskSEDmesonsForJetCorrelations.cxx++g");
+ gROOT->LoadMacro("AddTaskDmesonsForJetCorrelations.C");
+ AliAnalysisTaskSEDmesonsForJetCorrelations *taskDmesonsFilter = AddTaskDmesonsForJetCorrelations(iCandType,sCutFile,bIsMC,sText);
+ // EMCal framework
+ // -- Physics selection task
+ gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalPhysicsSelection.C");
+ AliPhysicsSelectionTask *physSelTask = AddTaskEmcalPhysicsSelection(kTRUE, kTRUE, uTriggerMask, 5, 5, 10, kTRUE, -1, -1, -1, -1);
+
+ if (!physSelTask) {
+ cout << "no physSelTask";
+ return;
+ }
+
+ // --
gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalSetup.C");
AliEmcalSetupTask *taskSetupEMCal = AddTaskEmcalSetup();
taskSetupEMCal->SetGeoPath("$ALICE_ROOT/OADB/EMCAL");
taskSetupEMCal->SelectCollisionCandidates(uTriggerMask);
- gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalPicoTrackMaker.C");
+ // Jet preparation
+ gROOT->LoadMacro("/data/Work/jets/testEMCalJetFramework/code/v4/AddTaskJetPreparation.C");
+ AddTaskJetPreparation(type,bIsMC,sRunPeriod);
+
+/*gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalPicoTrackMaker.C");
AliEmcalPicoTrackMaker *taskPicoTrack = AddTaskEmcalPicoTrackMaker(sUsedTrks.Data(),sInputTrk.Data(),sRunPeriod.Data());
- taskPicoTrack->SelectCollisionCandidates(uTriggerMask);
+ taskPicoTrack->SelectCollisionCandidates(uTriggerMask);*/
gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskEmcalJet.C");
- gROOT->LoadMacro("AddTaskFlavourJetCorrelations.C");
+ gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskEmcalJetSample.C");
+ gROOT->LoadMacro("AliAnalysisTaskRecoJetCorrelations.cxx++g");
+ gROOT->LoadMacro("AddTaskRecoJetCorrelations.C");
for (Int_t i=0; i<nRadius; i++) {
AliEmcalJetTask *taskFJ = AddTaskEmcalJet(sUsedTrks.Data(),sUsedClus.Data(),iJetAlgo,aRadius[i],iJetType);
taskFJ->SelectCollisionCandidates(uTriggerMask);
- AliAnalysisTaskFlavourJetCorrelations *taskFlavourCJ = AddTaskFlavourJetCorrelations(sUsedTrks,sUsedClus,taskFJ->GetName());
- taskFlavourCJ->SetName(Form("AliAnalysisTaskFlavourJetCorrelations_%s",sRadius[i].Data()));
- taskFlavourCJ->SetForceBeamType(uBeamType);
- taskFlavourCJ->SetAnaType(uAnaType);
- taskFlavourCJ->SetLeadingHadronType(iLeading);
- taskFlavourCJ->SetJetRadius(aRadius[i]);
- taskFlavourCJ->SetJetPtCut(dJetPtCut);
- taskFlavourCJ->SetPercAreaCut(dJetAreaCut);
- taskFlavourCJ->SelectCollisionCandidates(uTriggerMask);
+ AliAnalysisTaskFlavourJetCorrelations *taskDmesonCJ = AddTaskFlavourJetCorrelations(iCandType,
+ sCutFile,
+ bIsMC,
+ taskFJ->GetName(),
+ Form("JetR%s",sRadius[i].Data()),
+ iLeading,
+ aRadius[i],
+ dJetPtCut,
+ iAccCut,
+ dJetAreaCut);
+
+ taskDmesonCJ->SetName(Form("AliAnalysisTaskSEEmcalJetDmesonsCJ_%s",sRadius[i].Data()));
+ taskDmesonCJ->SetForceBeamType(uBeamType);
+ taskDmesonCJ->SetAnaType(uAnaType);
+ taskDmesonCJ->SetLeadingHadronType(iLeading);
+// taskDmesonCJ->SelectCollisionCandidates(uTriggerMask);
}
return;
// $Id$
-const Int_t mode = 0; // 0-->Local, 1-->Grid interactive
-const Bool_t bIsAOD = kTRUE;
-const Bool_t bIsMC = kFALSE;
-
-const Bool_t bIsPhysSel = kFALSE;
-const Bool_t bIsCentSel = kFALSE;
-const Bool_t bIsEvnPSel = kFALSE;
-const Bool_t bIsRespPID = kFALSE;
-
-const TString setData = "dataset.txt";
-//=============================================================================
-
-const Bool_t bIsEMCalAna = kFALSE;
-
-const Bool_t bFastOnly = kTRUE;
-const Bool_t bHistoW = kTRUE;
-const Double_t dMinEclus = 5.;
-const Double_t dMinPtClus = 5.;
-const Double_t dVz = 10.;
-const Bool_t bVzDiff = kTRUE;
-const Double_t dCentMin = -1.;
-const Double_t dCentMax = -1.;
-const Double_t dMinScaleCT = -1.;
-const Double_t dMaxScaleCT = -1.;
-//=============================================================================
-
-void AnalysisTrainCorrJetsLocal()
+// runEMCalJetAnalysis.C
+// =====================
+// This macro can be used to run a jet analysis within the EMCal Jet Framework.
+//
+// Examples:
+// -> Analyze ESDs from the pA pilot run on the AliEn grid with your task in AnaClass.cxx/.h
+// dataType = "ESD", useGrid = kTRUE, pattern = "*ESDs/pass2/*ESDs.root", addCXXs = "AnaClass.cxx",
+// addHs = "AnaClass.h", gridDir = "/alice/data/2012/LHC12g", gridMode = "full", runNumbers = "188359 188362"
+//
+// -> Analyze AODs (up to 96 files) locally given in files_aod.txt
+// dataType = "AOD", useGrid = kFALSE, numLocalFiles = 96
+//
+// MERGING ON ALIEN
+// ++++++++++++++++
+// If you run on the grid, you can monitor the jobs with alimonitor.cern.ch. When enough of them are in DONE state,
+// you have to merge the output. This can be done automatically, if you just change the gridMode to "terminate" and
+// give the EXACT name of the task whose output should be merged in uniqueName.
+//
+//
+// Authors: R. Haake, S. Aiola
+
+#include <ctime>
+#include "TGrid.h"
+
+AliAnalysisGrid* CreateAlienHandler(const char* uniqueName, const char* gridDir, const char* gridMode, const char* runNumbers, const Int_t nrunspermaster,
+ const char* pattern, TString additionalCode, TString additionalHeaders, Int_t maxFilesPerWorker,
+ Int_t workerTTL, Bool_t isMC);
+
+//______________________________________________________________________________
+void AnalysisTrainCorrJetsLocal (
+ const char* dataType = "AOD", // set the analysis type, AOD, ESD or sESD
+ Bool_t useGrid = kTRUE, // local or grid
+ TString localfilename = "file_aodlhc10d.txt",
+ const char* gridMode = "test", // set the grid run mode (can be "full", "test", "offline", "submit" or "terminate")
+ const char* pattern = "*ESDs/pass2/AOD099/*AliAOD.root", //"*/*/AliAOD.root" "*ESDs/pass1/AOD106/*AliAOD.root", // file pattern (here one can specify subdirs like passX etc.) (used on grid)
+ const char* gridDir = "/alice/data/2010/LHC10d", // /alice/data/2011/LHC11d /alice/sim/2012/LHC12f2b dir on alien, where the files live (used on grid)
+ const char* runNumbers = /*126437*/" 126432 126425 126424 126422 126409 126408 126407 126406 126405 126404 126403 126359 126352 126351 126350 126285 126284 126283 126168 126167 126160 126158 126097 126090 126088 126082 126081 126078 126073 126008 126007 126004 125855 125851 125850 125849 125848 125847 125844 125843 125842 125633 125632 125630 125296 125134 125101 125100 125097 125085 125023 124751 122375 122374", // considered run numbers (used on grid) /*LHC12g 188359 188362, LHC11a 146860 146859*/ /*LHC12f2b 158285 159582 */
+ const Int_t nrunspermaster= 100,
+ UInt_t numLocalFiles = 3, // number of files analyzed locally
+ const char* runPeriod = "LHC10d", // set the run period (used on grid)
+ const char* uniqueName = "DJetNewCode", // sets base string for the name of the task on the grid
+ UInt_t pSel = AliVEvent::kAny, // used event selection for every task except for the analysis tasks
+ Bool_t useTender = kFALSE, // trigger, if tender task should be used
+ Bool_t isMC = kFALSE, // trigger, if MC handler should be used
+
+ // Here you have to specify additional code files you want to use but that are not in aliroot
+ const char* addCXXs = "AliAnalysisTaskRecoJetCorrelations.cxx AliAnalysisTaskSEDmesonsForJetCorrelations.cxx",
+ const char* addHs = "AliAnalysisTaskRecoJetCorrelations.h AliAnalysisTaskSEDmesonsForJetCorrelations.h",
+
+ // These two settings depend on the dataset and your quotas on the AliEN services
+ Int_t maxFilesPerWorker = 20,
+ Int_t workerTTL = 86000,
+ Int_t nfiletestmode = 1
+ )
{
- if (mode==1 && !TGrid::Connect("alien://")) {
- ::Error("AnalysisTrainCorrJetsLocal.C::AnalysisTrainCorrJetsLocal", "Can not connect to the Grid!");
- return;
- }
- if (LoadLibraries()) {
- ::Error("AnalysisTrainCorrJetsLocal.C::AnalysisTrainCorrJetsLocal", "Load libraries failed!");
- return;
- }
-//=============================================================================
- TChain *chain = CreateAODFriendChain(setData);
- if (!chain) {
- ::Error("AnalysisTrainCorrJetsLocal.C::AnalysisTrainCorrJetsLocal", "Creating input chain failed!");
- return;
+ // Some pre-settings and constants
+
+ enum AlgoType {kKT, kANTIKT};
+ enum JetType {kFULLJETS, kCHARGEDJETS, kNEUTRALJETS};
+ gSystem->SetFPEMask();
+ gSystem->Setenv("ETRAIN_ROOT", ".");
+ gSystem->Setenv("ETRAIN_PERIOD", runPeriod);
+ // change this objects to strings
+ TString usedData(dataType);
+ TString additionalCXXs(addCXXs);
+ TString additionalHs(addHs);
+ cout << dataType << " analysis chosen" << endl;
+ if (useGrid)
+ {
+ cout << "-- using AliEn grid.\n";
+ if (usedData == "sESD")
+ {
+ cout << "Skimmed ESD analysis not available on the grid!" << endl;
+ return;
+ }
}
-//=============================================================================
+ else
+ cout << "-- using local analysis.\n";
+
+
+ // Load necessary libraries
+ LoadLibs();
- const UInt_t triggerMask = AliVEvent::kMB;
- AliAnalysisManager *mgr = new AliAnalysisManager("AnalysisTrainCorrJetsLocal", "Analysis Train Jet Correlation Local");
+ // Create analysis manager
+ AliAnalysisManager* mgr = new AliAnalysisManager(uniqueName);
- if (bIsAOD) {
+ // Check type of input and create handler for it
+ TString localFiles("-1");
+ if(usedData == "AOD")
+ {
+ localFiles = localfilename;
gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/train/AddAODHandler.C");
- AliAODInputHandler *aodIH = AddAODHandler();
- } else {
+ AliAODInputHandler* aodH = AddAODHandler();
+ }
+ else if((usedData == "ESD") || (usedData == "sESD"))
+ {
+ if (usedData == "ESD")
+ localFiles = "files_esd.txt";
+ else
+ localFiles = "files_sesd.txt";
+
gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/train/AddESDHandler.C");
- AliESDInputHandler *esdIH = AddESDHandler();
-// esdIH->SetReadFriends(kFALSE);
+ AliESDInputHandler* esdH = AddESDHandler();
}
- if (bIsMC && !bIsAOD) {
- AliMCEventHandler *mctEH = new AliMCEventHandler();
- mcH->SetPreReadMode(AliMCEventHandler::kLmPreRead);
- mcH->SetReadTR(kTRUE);
- mgr->SetMCtruthEventHandler(mcH);
+ else
+ {
+ cout << "Data type not recognized! You have to specify ESD, AOD, or sESD!\n";
}
-//=============================================================================
-
- if (bIsPhysSel && !bIsAOD && !bIsEMCalAna) {
- gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
- AliPhysicsSelectionTask *taskPhysSel = AddTaskPhysicsSelection(bIsMC);
+
+ if(!useGrid)
+ cout << "Using " << localFiles.Data() << " as input file list.\n";
+
+ gROOT->LoadMacro("AddTasksCorrJets.C");
+ //List of arguments:
+
+// const Int_t iCandtype=1 /*0 = D0, 1=Dstar...*/,
+// const TString sCutFile = "cutsHF/D0toKpiCutsppRecVtxNoPileupRejNoEMCAL.root",
+// const Double_t dJetPtCut = 1.,
+// const Double_t dJetAreaCut = 0.,
+// const Int_t iAccCut = 1,
+// const TString sRunPeriod = "LHC10b",
+// const Int_t uBeamType = 0,
+// const UInt_t uTriggerMask = AliVEvent::kMB, /*for jets and phys sel; the D mesons trigger is defined in the cut object*/
+// const Bool_t bIsMC=kFALSE,
+// TString sText=""/*completes the name of the candidate task lists*/
+
+ AddTasksCorrJets(1,"/data/Work/jets/testEMCalJetFramework/CutFilesMB/DStartoKpipiCuts_new.root",10.,0.,1,runPeriod,0,pSel,isMC,"");
+
+ // Set the physics selection for all given tasks
+ TObjArray *toptasks = mgr->GetTasks();
+ for (Int_t i=0; i<toptasks->GetEntries(); ++i)
+ {
+ AliAnalysisTaskSE *task = dynamic_cast<AliAnalysisTaskSE*>(toptasks->At(i));
+ if (!task)
+ continue;
+ if (task->InheritsFrom("AliPhysicsSelectionTask"))
+ continue;
+ ::Info("setPSel", "Set physics selection for %s (%s)", task->GetName(), task->ClassName());
+ task->SelectCollisionCandidates(pSel);
}
- if (bIsCentSel && !bIsAOD) {
- gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskCentrality.C");
- AliCentralitySelectionTask *taskCentrality = AddTaskCentrality();
- if (bIsMC) taskCentrality->SetMCInput();
- }
+ if(gridMode=="full") mgr->SetUseProgressBar(1, 25);
+
+ if (!mgr->InitAnalysis())
+ return;
+ mgr->PrintStatus();
- if (bIsEvnPSel) {
- gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskEventplane.C");
- AliEPSelectionTask *taskEventPlane = AddTaskEventplane();
- if (bIsMC) taskEventPlane->SetUseMCRP();
- }
+ if (useGrid)
+ { // GRID CALCULATION
- if (bIsRespPID) {
- gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPIDResponse.C");
- AliAnalysisTaskSE *taskRespPID = AddTaskPIDResponse(bIsMC);
- }
-//=============================================================================
-
- if (bIsEMCalAna) {
- if (bIsEvnPSel) {
- gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalPhysicsSelection.C");
- AliPhysicsSelectionTask *taskPhysSel = AddTaskEmcalPhysicsSelection(bFastOnly, bHistoW, triggerMask,
- dMinEclus, dMinPtClus,
- dVz, bVzDiff,
- bIsCentSel ? dCentMin : -1.,
- bIsCentSel ? dCentMax : -1.,
- dMinScaleCT, dMaxScaleCT);
- }
+ AliAnalysisGrid *plugin = CreateAlienHandler(uniqueName, gridDir, gridMode, runNumbers, nrunspermaster, pattern, additionalCXXs, additionalHs, maxFilesPerWorker, workerTTL, isMC);
+ plugin->SetNtestFiles(nfiletestmode);
- gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/AddTaskEmcalSetup.C");
- AliEmcalSetupTask *taskSetupEMCal = AddTaskEmcalSetup();
- taskSetupEMCal->SetGeoPath("$ALICE_ROOT/OADB/EMCAL");
+ mgr->SetGridHandler(plugin);
- gROOT->LoadMacro("$ALICE_ROOT/PWGJE/EMCALJetTasks/macros/AddTaskJetPreparation.C");
- AddTaskJetPreparation(bIsAOD ? "AOD" : "ESD"); //TODO
+ // start analysis
+ cout << "Starting GRID Analysis...";
+ if(gridMode=="test") mgr->SetDebugLevel(10);
+ else mgr->SetDebugLevel(0);
+ mgr->StartAnalysis("grid");
+ }
+ else
+ { // LOCAL CALCULATION
+
+ TChain* chain = 0;
+ if (usedData == "AOD")
+ {
+ Printf("Run Create AOD Chain");
+ gROOT->LoadMacro("/data/Work/jets/testEMCalJetFramework/AODchainWithFriend/CreateAODChain.C");
+ chain = CreateAODChain(localFiles.Data(), numLocalFiles,0,kTRUE,kTRUE);
+ //Printf("Chain Friend has %d entries", ((TTree*)chain->GetFriend())->GetEntriesFast());
+ }
+ else
+ { // ESD or skimmed ESD
+ gROOT->LoadMacro("$ALICE_ROOT/PWG/EMCAL/macros/CreateESDChain.C");
+ chain = CreateESDChain(localFiles.Data(), numLocalFiles);
+ }
+
+ // start analysis
+ cout << "Starting LOCAL Analysis...";
+ mgr->SetDebugLevel(10);
+ mgr->StartAnalysis("local", chain);
}
-//=============================================================================
-
- gROOT->LoadMacro("AddTasksFlavorJet.C"); AddTasksFlavorJet();
- if (mgr->InitAnalysis()) { mgr->PrintStatus(); mgr->StartAnalysis("local",chain); }
- return;
}
-//=============================================================================
-TChain *CreateAODFriendChain(TString setData)
-{
- if (!setData.EndsWith(".txt")) return 0;
-
- TChain *chain = new TChain("aodTree");
- TChain *cFrid = new TChain("aodTree");
-
- TString dataFile;
- ifstream dataList(setData.Data(), ios::in);
- while (!dataList.eof()) {
- dataFile.ReadLine(dataList,kFALSE);
- if (!dataFile.EndsWith("AliAOD.root")) continue;
- if (!gSystem->AccessPathName(dataFile.Data())) chain->Add(dataFile.Data());
-
- dataFile.ReplaceAll("AliAOD.root","AliAOD.VertexingHF.root");
- if (!gSystem->AccessPathName(dataFile.Data())) cFrid->Add(dataFile.Data());
- } dataList.close();
-
- chain->AddFriend(cFrid);
- return chain;
+//______________________________________________________________________________
+void LoadLibs()
+{
+ // Load common libraries (better too many than too few)
+ gSystem->Load("libTree");
+ gSystem->Load("libVMC");
+ gSystem->Load("libGeom");
+ gSystem->Load("libGui");
+ gSystem->Load("libXMLParser");
+ gSystem->Load("libMinuit");
+ gSystem->Load("libMinuit2");
+ gSystem->Load("libProof");
+ gSystem->Load("libPhysics");
+ gSystem->Load("libSTEERBase");
+ gSystem->Load("libESD");
+ gSystem->Load("libAOD");
+ gSystem->Load("libOADB");
+ gSystem->Load("libANALYSIS");
+ gSystem->Load("libCDB");
+ gSystem->Load("libRAWDatabase");
+ gSystem->Load("libSTEER");
+ gSystem->Load("libEVGEN");
+ gSystem->Load("libANALYSISalice");
+ gSystem->Load("libCORRFW");
+ gSystem->Load("libTOFbase");
+ //gSystem->Load("libTOFrec");
+ gSystem->Load("libRAWDatabase.so");
+ gSystem->Load("libRAWDatarec.so");
+ gSystem->Load("libTPCbase.so");
+ gSystem->Load("libTPCrec.so");
+ gSystem->Load("libITSbase.so");
+ gSystem->Load("libITSrec.so");
+ gSystem->Load("libTRDbase.so");
+ gSystem->Load("libTENDER.so");
+ gSystem->Load("libSTAT.so");
+ gSystem->Load("libTRDrec.so");
+ gSystem->Load("libHMPIDbase.so");
+ gSystem->Load("libPWGPP.so");
+ gSystem->Load("libPWGHFbase");
+ gSystem->Load("libPWGDQdielectron");
+ gSystem->Load("libPWGHFhfe");
+ gSystem->Load("libPWGflowBase.so");
+ gSystem->Load("libPWGflowTasks.so");
+ gSystem->Load("libPWGHFvertexingHF");
+ gSystem->Load("libEMCALUtils");
+ gSystem->Load("libPHOSUtils");
+ gSystem->Load("libPWGCaloTrackCorrBase");
+ gSystem->Load("libEMCALraw");
+ gSystem->Load("libEMCALbase");
+ gSystem->Load("libEMCALrec");
+ gSystem->Load("libTRDbase");
+ gSystem->Load("libVZERObase");
+ gSystem->Load("libVZEROrec");
+ gSystem->Load("libTENDER");
+ gSystem->Load("libTENDERSupplies");
+ gSystem->Load("libPWGEMCAL");
+ gSystem->Load("libPWGGAEMCALTasks");
+ gSystem->Load("libPWGTools");
+ gSystem->Load("libPWGCFCorrelationsBase");
+ gSystem->Load("libPWGCFCorrelationsDPhi");
+
+ //load CGAL, Fastjet and SISCone
+ //gSystem->Load("libCGAL");
+ gSystem->Load("libfastjet");
+ gSystem->Load("libSISConePlugin");
+
+ gSystem->Load("libJETAN");
+ gSystem->Load("libFASTJETAN");
+ gSystem->Load("libPWGJEEMCALJetTasks");
+
+
+ // include paths
+ gSystem->AddIncludePath("-Wno-deprecated");
+ gSystem->AddIncludePath("-I$ALICE_ROOT -I$ALICE_ROOT/include -I$ALICE_ROOT/EMCAL");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/PWGDQ/dielectron -I$ALICE_ROOT/PWGHF/hfe");
+ gSystem->AddIncludePath("-I$ALICE_ROOT/JETAN -I$ALICE_ROOT/JETAN/fastjet");
}
-//=============================================================================
-Bool_t LoadLibraries()
+AliAnalysisGrid* CreateAlienHandler(const char* uniqueName, const char* gridDir, const char* gridMode, const char* runNumbers, const Int_t nrunspermaster,
+ const char* pattern, TString additionalCode, TString additionalHeaders, Int_t maxFilesPerWorker,
+ Int_t workerTTL, Bool_t isMC)
{
- if (gSystem->Load("libTree") <0) return kTRUE;
- if (gSystem->Load("libGeom") <0) return kTRUE;
- if (gSystem->Load("libPhysics") <0) return kTRUE;
- if (gSystem->Load("libVMC") <0) return kTRUE;
- if (gSystem->Load("libMinuit") <0) return kTRUE;
- if (gSystem->Load("libMinuit2") <0) return kTRUE;
-
- if (gSystem->Load("libCore") <0) return kTRUE;
- if (gSystem->Load("libXMLIO") <0) return kTRUE;
- if (gSystem->Load("libXMLParser") <0) return kTRUE;
- if (gSystem->Load("libProof") <0) return kTRUE;
- if (gSystem->Load("libProofPlayer")<0) return kTRUE;
- if (gSystem->Load("libGui") <0) return kTRUE;
-//=============================================================================
-
- gSystem->AddIncludePath("-Wno-deprecated");
- gSystem->AddIncludePath("-I$ALICE_ROOT/include");
- gSystem->AddIncludePath("-I$ALICE_ROOT/PWGHF/vertexingHF");
- gSystem->AddIncludePath("-I$ALICE_ROOT/EMCAL");
- gSystem->AddIncludePath("-I$ALICE_ROOT/JETAN");
- gSystem->AddIncludePath("-I$ALICE_ROOT/JETAN/fastjet");
-
- if (gSystem->Load("libSTEERBase") <0) return kTRUE;
- if (gSystem->Load("libESD") <0) return kTRUE;
- if (gSystem->Load("libAOD") <0) return kTRUE;
- if (gSystem->Load("libANALYSIS") <0) return kTRUE;
- if (gSystem->Load("libOADB") <0) return kTRUE;
- if (gSystem->Load("libANALYSISalice") <0) return kTRUE;
- if (gSystem->Load("libCORRFW") <0) return kTRUE;
-
- if (gSystem->Load("libPWGTools") <0) return kTRUE;
- if (gSystem->Load("libPWGflowBase") <0) return kTRUE;
- if (gSystem->Load("libPWGflowTasks") <0) return kTRUE;
- if (gSystem->Load("libPWGHFbase") <0) return kTRUE;
- if (gSystem->Load("libPWGHFvertexingHF") <0) return kTRUE;
-
- if (gSystem->Load("libSTAT") <0) return kTRUE;
- if (gSystem->Load("libEMCALUtils") <0) return kTRUE;
-//if (gSystem->Load("libPHOSUtils") <0) return kTRUE;
-
- if (gSystem->Load("libCDB") <0) return kTRUE;
- if (gSystem->Load("libRAWDatabase") <0) return kTRUE;
- if (gSystem->Load("libRAWDatarec") <0) return kTRUE;
- if (gSystem->Load("libSTEER") <0) return kTRUE;
- if (gSystem->Load("libITSbase") <0) return kTRUE;
- if (gSystem->Load("libITSrec") <0) return kTRUE;
- if (gSystem->Load("libTPCbase") <0) return kTRUE;
- if (gSystem->Load("libTPCrec") <0) return kTRUE;
- if (gSystem->Load("libTRDbase") <0) return kTRUE;
- if (gSystem->Load("libTRDrec") <0) return kTRUE;
- if (gSystem->Load("libTOFbase") <0) return kTRUE;
-//if (gSystem->Load("libTOFrec") <0) return kTRUE;
- if (gSystem->Load("libHMPIDbase") <0) return kTRUE;
- if (gSystem->Load("libEMCALraw") <0) return kTRUE;
- if (gSystem->Load("libEMCALbase") <0) return kTRUE;
- if (gSystem->Load("libEMCALrec") <0) return kTRUE;
- if (gSystem->Load("libVZERObase") <0) return kTRUE;
- if (gSystem->Load("libVZEROrec") <0) return kTRUE;
- if (gSystem->Load("libTENDER") <0) return kTRUE;
- if (gSystem->Load("libTENDERSupplies") <0) return kTRUE;
-
- if (gSystem->Load("libCGAL") <0) return kTRUE;
- if (gSystem->Load("libfastjet") <0) return kTRUE;
- if (gSystem->Load("libsiscone") <0) return kTRUE;
- if (gSystem->Load("libSISConePlugin") <0) return kTRUE;
-
- if (gSystem->Load("libJETAN") <0) return kTRUE;
- if (gSystem->Load("libFASTJETAN") <0) return kTRUE;
- if (gSystem->Load("libPWGEMCAL") <0) return kTRUE;
- if (gSystem->Load("libPWGGAEMCALTasks") <0) return kTRUE;
- if (gSystem->Load("libPWGJEEMCALJetTasks")<0) return kTRUE;
-
- if (gROOT->LoadMacro("AliAnalysisTaskSEDmesonsFilterCJ.cxx+") <0) return kTRUE;
- if (gROOT->LoadMacro("AliAnalysisTaskFlavorJetCorrelations.cxx+")<0) return kTRUE;
- return kFALSE;
+ TDatime currentTime;
+ TString tmpName(uniqueName);
+ /*
+ // Only add current date and time when not in terminate mode! In this case the exact name has to be supplied by the user
+ if(strcmp(gridMode, "terminate"))
+ {
+ tmpName += "_";
+ tmpName += currentTime.GetDate();
+ tmpName += "_";
+ tmpName += currentTime.GetTime();
+ }else tmpName +="_20130412_122423";
+ */
+ TString tmpAdditionalLibs("");
+ tmpAdditionalLibs = Form("libTree.so libVMC.so libGeom.so libGui.so libXMLParser.so libMinuit.so libMinuit2.so libProof.so libPhysics.so libSTEERBase.so libESD.so libAOD.so libOADB.so libANALYSIS.so libCDB.so libRAWDatabase.so libSTEER.so libANALYSISalice.so libCORRFW.so libTOFbase.so libRAWDatabase.so libRAWDatarec.so libTPCbase.so libTPCrec.so libITSbase.so libITSrec.so libTRDbase.so libTENDER.so libSTAT.so libTRDrec.so libHMPIDbase.so libPWGPP.so libPWGHFbase.so libPWGDQdielectron.so libPWGHFhfe.so libPWGflowBase.so libPWGflowTasks.so libPWGHFvertexingHF.so libEMCALUtils.so libPHOSUtils.so libPWGCaloTrackCorrBase.so libEMCALraw.so libEMCALbase.so libEMCALrec.so libTRDbase.so libVZERObase.so libVZEROrec.so libTENDER.so libTENDERSupplies.so libPWGEMCAL.so libPWGGAEMCALTasks.so libPWGTools.so libPWGCFCorrelationsBase.so libPWGCFCorrelationsDPhi.so libCGAL.so libJETAN.so libfastjet.so libSISConePlugin.so libFASTJETAN.so libPWGJE.so libPWGmuon.so libPWGJEEMCALJetTasks.so %s %s",additionalCode.Data(),additionalHeaders.Data());
+
+
+ TString macroName("");
+ TString execName("");
+ TString jdlName("");
+ macroName = Form("%s.C", tmpName.Data());
+ execName = Form("%s.sh", tmpName.Data());
+ jdlName = Form("%s.jdl", tmpName.Data());
+
+ AliAnalysisAlien *plugin = new AliAnalysisAlien();
+ plugin->SetOverwriteMode();
+ plugin->SetRunMode(gridMode);
+
+ // Here you can set the (Ali)ROOT version you want to use
+ plugin->SetAPIVersion("V1.1x");
+ plugin->SetROOTVersion("v5-34-05");
+ plugin->SetAliROOTVersion("v5-04-38-AN");
+
+ plugin->SetGridDataDir(gridDir); // e.g. "/alice/sim/LHC10a6"
+ plugin->SetDataPattern(pattern); //dir structure in run directory
+ plugin->SetFriendChainName("./AliAOD.VertexingHF.root");
+ if (!isMC)
+ plugin->SetRunPrefix("000");
+
+ plugin->AddRunList(runNumbers);
+ plugin->SetNrunsPerMaster(nrunspermaster);
+
+ plugin->SetGridWorkingDir(Form("%s",tmpName.Data()));
+ plugin->SetGridOutputDir("output");
+
+ plugin->SetAnalysisSource(additionalCode.Data());
+ plugin->SetAdditionalLibs(tmpAdditionalLibs.Data());
+ plugin->AddIncludePath("-I. -I$ROOTSYS/include -I$ALICE_ROOT -I$ALICE_ROOT/include -I$ALICE_ROOT/STEER/STEER -I$ALICE_ROOT/STEER/STEERBase -I$ALICE_ROOT/STEER/ESD -I$ALICE_ROOT/STEER/AOD -I$ALICE_ROOT/ANALYSIS -I$ALICE_ROOT/OADB -I$ALICE_ROOT/PWGHF -I$ALICE_ROOT/PWGHF/base -I$ALICE_ROOT/PWGHF/vertexingHF -I$ALICE_ROOT/PWG/FLOW/Base -I$ALICE_ROOT/PWG/FLOW/Tasks -I$ALICE_ROOT/PWGJE -I$ALICE_ROOT/JETAN -I$ALICE_ROOT/PWGJE/EMCALJetTasks -g");
+ plugin->AddExternalPackage("boost::v1_43_0");
+ plugin->AddExternalPackage("cgal::v3.6");
+ plugin->AddExternalPackage("fastjet::v2.4.2");
+
+ plugin->SetDefaultOutputs(kTRUE);
+ // merging via jdl
+ plugin->SetMergeViaJDL(kTRUE);
+ plugin->SetOneStageMerging(kFALSE);
+ plugin->SetMaxMergeStages(2);
+
+ //plugin->SetMergeExcludes("");
+ plugin->SetAnalysisMacro(macroName.Data());
+ plugin->SetSplitMaxInputFileNumber(maxFilesPerWorker);
+ plugin->SetExecutable(execName.Data());
+ plugin->SetTTL(workerTTL);
+ plugin->SetInputFormat("xml-single");
+ plugin->SetJDLName(jdlName.Data());
+ plugin->SetPrice(1);
+ plugin->SetSplitMode("se");
+
+ return plugin;
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff