// Jet v2 task using QA method, based on jet sample task (S.Aiola). // // Authors: Jason Mueller (CERN summer student 2014) & Alice Ohlson #include #include #include #include #include #include #include #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 jetPtBins[nJetPtBins]) { hEventData->Fill("leadingJetPt>jetPtMax",1); return kFALSE; } if(fCent < centBins[0]) { hEventData->Fill("cent 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(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(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; }