--- /dev/null
+// Jet v2 task using QA method, based on jet sample task (S.Aiola).
+//
+// Authors: Jason Mueller (CERN summer student 2014) & Alice Ohlson
+
+
+#include <TClonesArray.h>
+#include <TH1F.h>
+#include <TH2F.h>
+#include <TH3F.h>
+#include <TList.h>
+#include <TLorentzVector.h>
+#include <TF1.h>
+
+#include "AliVCluster.h"
+#include "AliAODCaloCluster.h"
+#include "AliESDCaloCluster.h"
+#include "AliVTrack.h"
+#include "AliEmcalJet.h"
+#include "AliRhoParameter.h"
+#include "AliLog.h"
+#include "AliJetContainer.h"
+#include "AliParticleContainer.h"
+#include "AliClusterContainer.h"
+#include "AliPicoTrack.h"
+
+#include "AliAnalysisTaskEmcalJetv2QA.h"
+
+ClassImp(AliAnalysisTaskEmcalJetv2QA)
+
+//________________________________________________________________________
+AliAnalysisTaskEmcalJetv2QA::AliAnalysisTaskEmcalJetv2QA() :
+AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetv2QA", kTRUE),
+ nCentBins(0),
+ nJetPtBins(0),
+ centBins(0x0),
+ jetPtBins(0x0),
+ fJetv2(0x0),
+ doPtWeight(kFALSE),
+ fHistTracksPt(0),
+ fHistClustersPt(0),
+ fHistLeadingJetPt(0),
+ fHistLeadingJetPtCorr(0),
+ fHistJetsPhiEta(0),
+ fHistJetsPtArea(0),
+ fHistJetsPtLeadHad(0),
+ fHistJetsCorrPtArea(0),
+ fHistPtDEtaDPhiTrackClus(0),
+ fHistPtDEtaDPhiClusTrack(0),
+ fDPhiJet(0),
+ fDPhiJetPythia(0),
+ fDPhiEP(0),
+ hGx(0),
+ hGy2(0),
+ hGxGy2(0),
+ hGy4(0),
+ hGx2(0),
+ hGx2Gy2(0),
+ hGxGy4(0),
+ hGy6(0),
+ hGx2Gy4(0),
+ hGxGy6(0),
+ hGy8(0),
+ hGy(0),
+ hN(0),
+ htv2std(0),
+ htjv2std(0),
+ htj2v2std(0),
+ hV0jv2std(0),
+ htdPsi(0),
+ htjdPsi(0),
+ htj2dPsi(0),
+ hV0jdPsi(0),
+ hAx(0),
+ hAxDijet(0),
+ hQx(0),
+ hQy(0),
+ hEventData(0),
+ hNTracks(0),
+ hNTracksCent(0),
+ hGxTracks(0),
+ hGyTracks(0),
+ hGy2Tracks(0),
+ hGxGy2Tracks(0),
+ hGy4Tracks(0),
+ htEPRes(0),
+ htjEPRes(0),
+ htj2EPRes(0),
+ fJetsCont(0),
+ fTracksCont(0),
+ fCaloClustersCont(0)
+{
+ // Default constructor.
+
+ cout << endl;
+ cout << "*******************************************************" << endl;
+ cout << "*** AliAnalysisTaskEmcalJetv2QA constructor! ***" << endl;
+ cout << "*******************************************************" << endl;
+ cout << endl;
+
+
+
+ SetMakeGeneralHistograms(kTRUE);
+}
+
+//________________________________________________________________________
+AliAnalysisTaskEmcalJetv2QA::AliAnalysisTaskEmcalJetv2QA(const char *name) :
+ AliAnalysisTaskEmcalJet(name, kTRUE),
+ nCentBins(0),
+ nJetPtBins(0),
+ centBins(0x0),
+ jetPtBins(0x0),
+ fJetv2(0x0),
+ doPtWeight(kFALSE),
+ fHistTracksPt(0),
+ fHistClustersPt(0),
+ fHistLeadingJetPt(0),
+ fHistLeadingJetPtCorr(0),
+ fHistJetsPhiEta(0),
+ fHistJetsPtArea(0),
+ fHistJetsPtLeadHad(0),
+ fHistJetsCorrPtArea(0),
+ fHistPtDEtaDPhiTrackClus(0),
+ fHistPtDEtaDPhiClusTrack(0),
+ fDPhiJet(0),
+ fDPhiJetPythia(0),
+ fDPhiEP(0),
+ hGx(0),
+ hGy2(0),
+ hGxGy2(0),
+ hGy4(0),
+ hGx2(0),
+ hGx2Gy2(0),
+ hGxGy4(0),
+ hGy6(0),
+ hGx2Gy4(0),
+ hGxGy6(0),
+ hGy8(0),
+ hGy(0),
+ hN(0),
+ htv2std(0),
+ htjv2std(0),
+ htj2v2std(0),
+ hV0jv2std(0),
+ htdPsi(0),
+ htjdPsi(0),
+ htj2dPsi(0),
+ hV0jdPsi(0),
+ hAx(0),
+ hAxDijet(0),
+ hQx(0),
+ hQy(0),
+ hEventData(0),
+ hNTracks(0),
+ hNTracksCent(0),
+ hGxTracks(0),
+ hGyTracks(0),
+ hGy2Tracks(0),
+ hGxGy2Tracks(0),
+ hGy4Tracks(0),
+ htEPRes(0),
+ htjEPRes(0),
+ htj2EPRes(0),
+ fJetsCont(0),
+ fTracksCont(0),
+ fCaloClustersCont(0)
+{
+ // Standard constructor.
+
+
+ SetMakeGeneralHistograms(kTRUE);
+}
+
+//________________________________________________________________________
+AliAnalysisTaskEmcalJetv2QA::~AliAnalysisTaskEmcalJetv2QA()
+{
+ // Destructor.
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetv2QA::UserCreateOutputObjects()
+{
+ // Create user output.
+
+ AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
+
+ fJetsCont = GetJetContainer(0);
+ if(fJetsCont) { //get particles and clusters connected to jets
+ fTracksCont = fJetsCont->GetParticleContainer();
+ fCaloClustersCont = fJetsCont->GetClusterContainer();
+ } else { //no jets, just analysis tracks and clusters
+ fTracksCont = GetParticleContainer(0);
+ fCaloClustersCont = GetClusterContainer(0);
+ }
+
+ nCentBins = 6; // set number of centrality bins
+ nJetPtBins = 4; // set number of jetPtBins
+ centBins = new Double_t[nCentBins+1];
+ jetPtBins = new Double_t[nJetPtBins+1];
+ centBins[0] = 0.; centBins[1] = 5.; centBins[2] = 10.; centBins[3] = 20.; centBins[4] = 30.; centBins[5] = 40.; centBins[6] = 50.; // set edges of bins
+ jetPtBins[0] = 50.; jetPtBins[1] = 70.; jetPtBins[2] = 90.; jetPtBins[3] = 120.; jetPtBins[4] = 200.; // set edges of bins
+
+ fTracksCont->SetClassName("AliVTrack");
+ fCaloClustersCont->SetClassName("AliVCluster");
+
+ fHistTracksPt = new TH1F("fHistTracksPt","fHistTracksPt", fNbins / 2, fMinBinPt, fMaxBinPt / 2);
+ fHistTracksPt->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
+ fHistTracksPt->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTracksPt);
+
+ fHistClustersPt = new TH1F("fHistClustersPt", "fHistClustersPt", fNbins / 2, fMinBinPt, fMaxBinPt / 2);
+ fHistClustersPt->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
+ fHistClustersPt->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistClustersPt);
+
+ fHistLeadingJetPt = new TH1F("fHistLeadingJetPt", "fHistLeadingJetPt", fNbins, fMinBinPt, fMaxBinPt);
+ fHistLeadingJetPt->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
+ fHistLeadingJetPt->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJetPt);
+
+ fHistLeadingJetPtCorr = new TH1F("fHistLeadingJetPtCorr", "fHistLeadingJetPtCorr", fNbins, fMinBinPt, fMaxBinPt);
+ fHistLeadingJetPtCorr->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
+ fHistLeadingJetPtCorr->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJetPtCorr);
+
+ fHistJetsPhiEta = new TH2F("fHistJetsPhiEta", "fHistJetsPhiEta", 50, -1, 1, 101, 0, TMath::Pi()*2 + TMath::Pi()/200);
+ fHistJetsPhiEta->GetXaxis()->SetTitle("#eta");
+ fHistJetsPhiEta->GetYaxis()->SetTitle("#phi");
+ fOutput->Add(fHistJetsPhiEta);
+
+ fHistJetsPtArea = new TH2F("fHistJetsPtArea", "fHistJetsPtArea", fNbins, fMinBinPt, fMaxBinPt, 30, 0, 3);
+ fHistJetsPtArea->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
+ fHistJetsPtArea->GetYaxis()->SetTitle("area");
+ fOutput->Add(fHistJetsPtArea);
+
+ fHistJetsPtLeadHad = new TH2F("fHistJetsPtLeadHad", "fHistJetsPtLeadHad", fNbins, fMinBinPt, fMaxBinPt, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
+ fHistJetsPtLeadHad->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
+ fHistJetsPtLeadHad->GetYaxis()->SetTitle("p_{T,lead} (GeV/c)");
+ fHistJetsPtLeadHad->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistJetsPtLeadHad);
+
+ fHistJetsCorrPtArea = new TH2F("fHistJetsCorrPtArea", "fHistJetsCorrPtArea", fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3);
+ fHistJetsCorrPtArea->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]");
+ fHistJetsCorrPtArea->GetYaxis()->SetTitle("area");
+ fOutput->Add(fHistJetsCorrPtArea);
+
+ fHistPtDEtaDPhiTrackClus = new TH3F("fHistPtDEtaDPhiTrackClus","fHistPtDEtaDPhiTrackClus;#it{p}_{T}^{track};#Delta#eta;#Delta#varphi",100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+ fOutput->Add(fHistPtDEtaDPhiTrackClus);
+
+ fHistPtDEtaDPhiClusTrack = new TH3F("fHistPtDEtaDPhiClusTrack","fHistPtDEtaDPhiClusTrack;#it{p}_{T}^{clus};#Delta#eta;#Delta#varphi",100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+ fOutput->Add(fHistPtDEtaDPhiClusTrack);
+
+ fDPhiJet = new TH1F("fDPhiJet","fDPhiJet",90, -TMath::Pi()/3, 5*TMath::Pi()/3);
+ fDPhiJet->GetXaxis()->SetTitle("#Delta#varphi");
+ fDPhiJet->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fDPhiJet);
+
+ fDPhiJetPythia = new TH1F("fDPhiJetPythia","fDPhiJetPythia",90, -TMath::Pi()/3, 5*TMath::Pi()/3);
+ fDPhiJetPythia->GetXaxis()->SetTitle("#Delta#varphi");
+ fDPhiJetPythia->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fDPhiJetPythia);
+
+ fDPhiEP = new TH1F("fDPhiEP","fDPhiEP",90, 0, 2*TMath::Pi());
+ fDPhiEP->GetXaxis()->SetTitle("#Delta#varphi");
+ fDPhiEP->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fDPhiEP);
+
+ hGx = new TH2D("hGx", "Gx v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGx->GetXaxis()->SetTitle("Centrality (%)");
+ hGx->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGx);
+
+ hGy2 = new TH2D("hGy2", "Gy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGy2->GetXaxis()->SetTitle("Centrality (%)");
+ hGy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGy2);
+
+ hGxGy2 = new TH2D("hGxGy2", "GxGy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGxGy2->GetXaxis()->SetTitle("Centrality (%)");
+ hGxGy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGxGy2);
+
+ hGy4 = new TH2D("hGy4", "Gy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGy4->GetXaxis()->SetTitle("Centrality (%)");
+ hGy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGy4);
+
+ hGx2 = new TH2D("hGx2", "Gx2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGx2->GetXaxis()->SetTitle("Centrality (%)");
+ hGx2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGx2);
+
+ hGx2Gy2 = new TH2D("hGx2Gy2", "Gx2Gy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGx2Gy2->GetXaxis()->SetTitle("Centrality (%)");
+ hGx2Gy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGx2Gy2);
+
+ hGxGy4 = new TH2D("hGxGy4", "GxGy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGxGy4->GetXaxis()->SetTitle("Centrality (%)");
+ hGxGy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGxGy4);
+
+ hGy6 = new TH2D("hGy6", "Gy6 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGy6->GetXaxis()->SetTitle("Centrality (%)");
+ hGy6->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGy6);
+
+ hGx2Gy4 = new TH2D("hGx2Gy4", "Gx2Gy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGx2Gy4->GetXaxis()->SetTitle("Centrality (%)");
+ hGx2Gy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGx2Gy4);
+
+ hGxGy6 = new TH2D("hGxGy6", "GxGy6 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGxGy6->GetXaxis()->SetTitle("Centrality (%)");
+ hGxGy6->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGxGy6);
+
+ hGy8 = new TH2D("hGy8", "Gy8 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGy8->GetXaxis()->SetTitle("Centrality (%)");
+ hGy8->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGy8);
+
+ hGy = new TH2D("hGy", "Gy v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hGy->GetXaxis()->SetTitle("Centrality (%)");
+ hGy->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hGy);
+
+ hN = new TH2D("hN", "N v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hN->GetXaxis()->SetTitle("Centrality (%)");
+ hN->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hN);
+
+ htv2std = new TH2D("htv2std", "v2std v Centrality v JetPt w/o jets", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htv2std->GetXaxis()->SetTitle("Centrality (%)");
+ htv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htv2std);
+
+ htjv2std = new TH2D("htjv2std", "v2std v Centrality v JetPt w/ jets", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htjv2std->GetXaxis()->SetTitle("Centrality (%)");
+ htjv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htjv2std);
+
+ htj2v2std = new TH2D("htj2v2std", "v2std v Centrality v JetPt w/ trackPt < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htj2v2std->GetXaxis()->SetTitle("Centrality (%)");
+ htj2v2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htj2v2std);
+
+ hV0jv2std = new TH2D("hV0jv2std", "v2std v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
+ hV0jv2std->GetXaxis()->SetTitle("Centrality (%)");
+ hV0jv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hV0jv2std);
+
+ Int_t ndpsibins = 100;
+ Double_t dpsibins[101];
+ for(Int_t t = 0; t < 101; t++) dpsibins[t] = TMath::Pi()*t/50.;
+
+ htdPsi = new TH3D("htdPsi", "JetAxis - EventPlane w/o jets", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
+ htdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
+ fOutput->Add(htdPsi);
+
+ htjdPsi = new TH3D("htjdPsi", "JetAxis - EventPlane w/ jets", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
+ htjdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
+ fOutput->Add(htjdPsi);
+
+ htj2dPsi = new TH3D("htj2dPsi", "JetAxis - EventPlane w/ trackPt < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
+ htj2dPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
+ fOutput->Add(htj2dPsi);
+
+ hV0jdPsi = new TH3D("hV0jdPsi", "JetAxis - EventPlane", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
+ hV0jdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
+ fOutput->Add(hV0jdPsi);
+
+ hQx = new TH2D("hQx", "Qx Distribution in EP frame", 100, -0.3, 0.3, nCentBins, centBins);
+ hQx->GetXaxis()->SetTitle("Qx");
+ hQx->GetYaxis()->SetTitle("Centrality (%)");
+ fOutput->Add(hQx);
+
+ hQy = new TH2D("hQy", "Qy Distribution in EP frame", 100, -0.3, 0.3, nCentBins, centBins);
+ hQy->GetXaxis()->SetTitle("Qy");
+ hQy->GetYaxis()->SetTitle("Centrality (%)");
+ fOutput->Add(hQy);
+
+ hAx = new TH2D("hAx", "Ax Distribution in EP frame w/o Dijets", 100, -35, 70, nJetPtBins, jetPtBins);
+ hAx->GetXaxis()->SetTitle("Ax");
+ hAx->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hAx);
+
+ hAxDijet = new TH2D("hAxDijet", "Ax Distribution in EP frame w/ Dijets", 100, -35, 70, nJetPtBins, jetPtBins);
+ hAxDijet->GetXaxis()->SetTitle("Ax");
+ hAxDijet->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(hAxDijet);
+
+ hEventData = new TH1F("hEventData","Events Kept and Discarded", 9, 0, 9);
+ hEventData->GetYaxis()->SetTitle("counts");
+ fOutput->Add(hEventData);
+
+ hNTracks = new TH2F("hNTracks","Number of Tracks Per Event", 100, 0, 3000, 3, 0, 3);
+ hNTracks->GetXaxis()->SetTitle("# tracks");
+ fOutput->Add(hNTracks);
+
+ hNTracksCent = new TH2D("hNTracksCent", "NTracks by centrality", 100, 0, 3000, (Int_t)centBins[nCentBins], centBins[0], centBins[nCentBins]);
+ hNTracksCent->GetXaxis()->SetTitle("# tracks");
+ hNTracksCent->GetYaxis()->SetTitle("Centrality (%)");
+ fOutput->Add(hNTracksCent);
+
+ hGxTracks = new TH2D("hGxTracks", "Gx by NTracks", 200, -200, 200, 100, 0, 3000);
+ hGxTracks->GetXaxis()->SetTitle("Gx");
+ hGxTracks->GetYaxis()->SetTitle("# tracks");
+ fOutput->Add(hGxTracks);
+
+ hGyTracks = new TH2D("hGyTracks", "Gy by NTracks", 200, -200, 200, 100, 0, 3000);
+ hGyTracks->GetXaxis()->SetTitle("Gy");
+ hGyTracks->GetYaxis()->SetTitle("# tracks");
+ fOutput->Add(hGyTracks);
+
+ hGy2Tracks = new TH2D("hGy2Tracks", "Gy2 by NTracks", 100, 0, 20000, 100, 0, 3000);
+ hGy2Tracks->GetXaxis()->SetTitle("Gy2");
+ hGy2Tracks->GetYaxis()->SetTitle("# tracks");
+ fOutput->Add(hGy2Tracks);
+
+ hGxGy2Tracks = new TH2D("hGxGy2Tracks", "GxGy2 by NTracks", 100, -100000, 100000, 100, 0, 3000);
+ hGxGy2Tracks->GetXaxis()->SetTitle("GxGy2");
+ hGxGy2Tracks->GetYaxis()->SetTitle("# tracks");
+ fOutput->Add(hGxGy2Tracks);
+
+ hGy4Tracks = new TH2D("hGy4Tracks", "Gy4 by NTracks", 100, 0, 100000000, 100, 0, 3000);
+ hGy4Tracks->GetXaxis()->SetTitle("Gy4");
+ hGy4Tracks->GetYaxis()->SetTitle("# tracks");
+ fOutput->Add(hGy4Tracks);
+
+ htEPRes = new TH2D("htEPRes", "EP Resolution w/o Jets", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htEPRes->GetXaxis()->SetTitle("Centrality (%)");
+ htEPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htEPRes);
+
+ htjEPRes = new TH2D("htjEPRes", "EP Resolution w/ Jets", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htjEPRes->GetXaxis()->SetTitle("Centrality (%)");
+ htjEPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htjEPRes);
+
+ htj2EPRes = new TH2D("htj2EPRes", "EP Resolution w/ trackPT < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins);
+ htj2EPRes->GetXaxis()->SetTitle("Centrality (%)");
+ htj2EPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
+ fOutput->Add(htj2EPRes);
+
+ PostData(1, fOutput); // Post data for ALL output slots > 0 here.
+}
+
+//________________________________________________________________________
+Bool_t AliAnalysisTaskEmcalJetv2QA::FillHistograms()
+{
+ // Fill histograms.
+
+ return kTRUE;
+}
+
+
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetv2QA::ExecOnce() {
+
+ AliAnalysisTaskEmcalJet::ExecOnce();
+
+ if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
+ if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
+ if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
+
+}
+
+//________________________________________________________________________
+Bool_t AliAnalysisTaskEmcalJetv2QA::Run() // this part loops over each event
+{
+ // Run analysis code here, if needed. It will be executed before FillHistograms().
+
+ Double_t jetPhi = -999; // angle of leading jet
+ Double_t jetPt = -999; // pt of leading jet
+ Double_t jetArea = -999;
+ Double_t trackPt = -999;
+ Double_t phi = -999; // track phi
+ Double_t dPhi = -999; // track phi - jet axis
+ Double_t dPhiQA = -999; // track phi - EP
+ Double_t tSin = 0; // used for std EP calc
+ Double_t tCos = 0; // used for std EP calc
+ Double_t jSin = 0; // used for std EP calc
+ Double_t jCos = 0; // used for std EP calc
+ Double_t tSin2 = 0; // used for std EP calc with trackPt < 2 GeV
+ Double_t tCos2 = 0; // used for std EP calc with trackPt < 2 GeV
+ Double_t qx = 0; // used for Qx distribution
+ Double_t qy = 0; // used for Qy distribution
+ Double_t ax = 0; // used for Ax distribution
+ Double_t tEP = -999; // EP w/o jets
+ Double_t tjEP = -999; // EP w/ jets
+ Double_t tjEP2 = -999; // EP w/ jets
+ Double_t dPsi = -999; // jet axis - EP
+ Double_t gx = 0; // used for G moment calc
+ Double_t gy = 0; // used for G moment calc
+ Int_t isDijet = 0; // if 0, no dijet. if 1, dijet.
+ Int_t nTracksBkgnd = 0; // used to keep track of number of tracks in background
+ Int_t nTracksJet =0; // used to keep track of number of tracks in jets
+
+ if(fJetsCont)
+ {
+
+ AliEmcalJet *jettest = fJetsCont->GetNextAcceptJet(0);
+ while(jettest) {
+
+ fHistJetsPtArea->Fill(jettest->Pt(), jettest->Area());
+ fHistJetsPhiEta->Fill(jettest->Eta(), jettest->Phi());
+
+ Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jettest);
+ fHistJetsPtLeadHad->Fill(jettest->Pt(), ptLeading);
+
+ if (fHistJetsCorrPtArea) {
+ Float_t corrPt = jettest->Pt() - fJetsCont->GetRhoVal() * jettest->Area();
+ fHistJetsCorrPtArea->Fill(corrPt, jettest->Area());
+ }
+ jettest = fJetsCont->GetNextAcceptJet();
+ }
+
+
+ AliEmcalJet *jet = fJetsCont->GetLeadingJet();
+ if(jet)
+ {
+ jetPhi = jet->Phi(); // get leading jet phi (jet axis)
+ jetPt = jet->Pt(); // get leading jet pT to filter out low pT jet events
+ jetArea = jet->Area();
+ }
+ }
+
+ // event cuts
+ if(!fTracksCont)
+ {
+ hEventData->Fill("!fTracksCont",1);
+ return kFALSE;
+ }
+ if(!fJetsCont)
+ {
+ hEventData->Fill("!fJetsCont",1);
+ return kFALSE;
+ }
+ if(jetPt == -999)
+ {
+ hEventData->Fill("jetPt=-999",1);
+ return kFALSE;
+ }
+ if(jetPt < jetPtBins[0])
+ {
+ hEventData->Fill("leadingJetPt<jetPtMin",1);
+ return kFALSE;
+ }
+ if(jetPt > jetPtBins[nJetPtBins])
+ {
+ hEventData->Fill("leadingJetPt>jetPtMax",1);
+ return kFALSE;
+ }
+ if(fCent < centBins[0])
+ {
+ hEventData->Fill("cent<centMin",1);
+ return kFALSE;
+ }
+ if(fCent > centBins[nCentBins])
+ {
+ hEventData->Fill("cent>centMax",1);
+ return kFALSE;
+ }
+ hEventData->Fill("good event",1);
+
+ AliEmcalJet *dijet = fJetsCont->GetNextAcceptJet(0); // check for dijet events
+ while(dijet)
+ {
+ if(dijet->Pt() > 50 && fabs(jetPhi-dijet->Phi()-TMath::Pi()) < 0.4) // loop over jets with pT>50 and exclude leading jet and check that angular separation is < 0.4
+ isDijet = 1;
+ dijet = fJetsCont->GetNextAcceptJet();
+ }
+
+ if (fCaloClustersCont)
+ {
+ AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
+ while(cluster) {
+ TLorentzVector nPart;
+ cluster->GetMomentum(nPart, fVertex);
+ fHistClustersPt->Fill(nPart.Pt());
+ cluster = fCaloClustersCont->GetNextAcceptCluster();
+ }
+ }
+
+ fHistLeadingJetPt->Fill(jetPt);
+ fHistLeadingJetPtCorr->Fill(jetPt-fJetsCont->GetRhoVal()*jetArea);
+
+ AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0));
+ while(track)
+ { // loop over all particles (including jet tracks)
+ trackPt = track->Pt();
+ fHistTracksPt->Fill(trackPt);
+ phi = track->Phi(); // get track phi
+
+ dPhi = phi-jetPhi; // get track phi - jet axis
+ if(dPhi < 0) dPhi += TMath::TwoPi();
+ if(dPhi > TMath::TwoPi()) dPhi -= TMath::TwoPi();
+
+ dPhiQA = phi-fEPV0; // get track phi - EP
+ if(dPhiQA < 0) dPhiQA += TMath::TwoPi();
+ if(dPhiQA > TMath::TwoPi()) dPhiQA -= TMath::TwoPi();
+ fDPhiEP->Fill(dPhiQA);
+
+ // fill jet-hadron correlation just to check if track labels make sense...
+ if(track->Pt()>1.)
+ {
+ Double_t dphiJet = dPhi;
+ if(dphiJet > 5*TMath::Pi()/3) dphiJet -= 2*TMath::Pi();
+ if(dphiJet < -TMath::Pi()/3) dphiJet += 2*TMath::Pi();
+ if(track->GetLabel() == 0)
+ fDPhiJet->Fill(dphiJet);
+ else
+ fDPhiJetPythia->Fill(dphiJet);
+ }
+
+ Double_t weight = 1.;
+ if(doPtWeight) weight = trackPt;
+
+ gx += weight*cos(2*dPhi);
+ gy += weight*sin(2*dPhi);
+
+ if(track->GetLabel() == 0)
+ { // sum for std EP method
+ tSin += weight*sin(2*phi); // bkgnd has label = 0
+ tCos += weight*cos(2*phi);
+ qx += weight*cos(2*dPhiQA);
+ qy += weight*sin(2*dPhiQA);
+ nTracksBkgnd += 1;
+ }
+ else
+ {
+ jSin += weight*sin(2*phi); // jets have label =/= 0
+ jCos += weight*cos(2*phi);
+ ax += weight*cos(2*dPhi);
+ nTracksJet += 1;
+ }
+
+ if(track->Pt() < 2.)
+ { // sum for std EP method w/ trackPt < 2 GeV
+ tSin2 += weight*sin(2*phi);
+ tCos2 += weight*cos(2*phi);
+ }
+
+
+ track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()); // increment to next track
+ } // close loop over particles
+
+ hNTracks->Fill(nTracksBkgnd,"Bkgnd Tracks",1);
+ hNTracks->Fill(nTracksJet,"Jet Tracks",1);
+ hNTracks->Fill(nTracksBkgnd+nTracksJet,"Total Tracks",1);
+ hNTracksCent->Fill(nTracksBkgnd+nTracksJet,fCent,1);
+
+ if(nTracksBkgnd == 0)
+ {
+ hEventData->Fill("no tracks",1);
+ hEventData->Fill("good event",-1);
+ return kFALSE;
+ }
+
+ Double_t v2weight = 1+2*fJetv2*cos(2*(jetPhi-fEPV0)); // set v2 weight for event
+
+ hGx->Fill(fCent, jetPt, v2weight*gx); // fill histograms for QA method
+ hGy2->Fill(fCent, jetPt, v2weight*gy*gy);
+ hGxGy2->Fill(fCent, jetPt, v2weight*gx*gy*gy);
+ hGy4->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy);
+ hGx2->Fill(fCent, jetPt, v2weight*gx*gx);
+ hGx2Gy2->Fill(fCent, jetPt, v2weight*gx*gx*gy*gy);
+ hGxGy4->Fill(fCent, jetPt, v2weight*gx*gy*gy*gy*gy);
+ hGy6->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy*gy*gy);
+ hGx2Gy4->Fill(fCent, jetPt, v2weight*gx*gx*gy*gy*gy*gy);
+ hGxGy6->Fill(fCent, jetPt, v2weight*gx*gy*gy*gy*gy*gy*gy);
+ hGy8->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy*gy*gy*gy*gy);
+ hGy->Fill(fCent, jetPt, v2weight*gy);
+ hN->Fill(fCent, jetPt, v2weight);
+
+ hGxTracks->Fill(gx,nTracksBkgnd+nTracksJet);
+ hGyTracks->Fill(gy,nTracksBkgnd+nTracksJet);
+ hGy2Tracks->Fill(gy*gy,nTracksBkgnd+nTracksJet);
+ hGxGy2Tracks->Fill(gx*gy*gy,nTracksBkgnd+nTracksJet);
+ hGy4Tracks->Fill(gy*gy*gy*gy,nTracksBkgnd+nTracksJet);
+ hQx->Fill(qx/nTracksBkgnd,fCent);
+ hQy->Fill(qy/nTracksBkgnd,fCent);
+
+ if(isDijet == 0)
+ hAx->Fill(ax,jetPt);
+ if(isDijet == 1)
+ hAxDijet->Fill(ax,jetPt);
+
+ tEP = 0.5*atan2(tSin,tCos); // calculate EP w/o jets
+ tjEP = 0.5*atan2((tSin+jSin),(tCos+jCos)); // calculate EP w/ jets
+ tjEP2 = 0.5*atan2(tSin2,tCos2); // calculate EP w/ trackPt < 2 GeV
+
+ htEPRes->Fill(fCent, jetPt, v2weight*cos(2*(tEP-fEPV0)));
+ htjEPRes->Fill(fCent, jetPt, v2weight*cos(2*(tjEP-fEPV0)));
+ htj2EPRes->Fill(fCent, jetPt, v2weight*cos(2*(tjEP2-fEPV0)));
+
+
+ dPsi = jetPhi-tEP;
+ if(dPsi < 0) dPsi += TMath::TwoPi();
+ if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
+
+ htv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/o jets
+ htdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/o jets
+
+ dPsi = jetPhi-tjEP;
+ if(dPsi < 0) dPsi += TMath::TwoPi();
+ if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
+
+ htjv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/ jets
+ htjdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/ jets
+
+ dPsi = jetPhi-tjEP2;
+ if(dPsi < 0) dPsi += TMath::TwoPi();
+ if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
+
+ htj2v2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/ trackPt < 2 GeV
+ htj2dPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/ trackPt < 2 GeV
+
+ dPsi = jetPhi-fEPV0;
+ if(dPsi < 0) dPsi += TMath::TwoPi();
+ if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
+
+ hV0jv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data
+ hV0jdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EPV0
+
+
+ return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetv2QA::Terminate(Option_t *)
+{
+ // Called once at the end of the analysis.
+ if(centBins) delete [] centBins;
+ if(jetPtBins) delete [] jetPtBins;
+}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff