//
// Author: S.Aiola
-#include <TChain.h>
#include <TClonesArray.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TH3F.h>
+#include <THnSparse.h>
#include <TList.h>
#include <TLorentzVector.h>
#include "AliExternalTrackParam.h"
#include "AliTrackerBase.h"
#include "AliLog.h"
+#include "AliEMCALGeometry.h"
+#include "AliEMCALGeoParams.h"
+#include "AliPicoTrack.h"
+#include "AliVVZERO.h"
+#include "AliESDUtils.h"
#include "AliAnalysisTaskSAQA.h"
AliAnalysisTaskSAQA::AliAnalysisTaskSAQA() :
AliAnalysisTaskEmcalJet("AliAnalysisTaskSAQA", kTRUE),
fCellEnergyCut(0.1),
- fDoTrigger(kFALSE),
- fRepropagateTracks(kFALSE),
- fTrgClusName("ClustersL1GAMMAFEE"),
- fTrgClusters(0),
- fHistCentrality(0),
- fHistZVertex(0),
- fHistTracksCent(0),
- fHistClusCent(0),
- fHistClusTracks(0),
- fHistCellsCent(0),
- fHistCellsTracks(0),
- fHistMaxL1FastORCent(0),
- fHistMaxL1ClusCent(0),
- fHistMaxL1ThrCent(0),
- fHistTracksPt(0),
- fHistTrPhiEta(0),
- fHistTrEmcPhiEta(0),
- fHistTrPhiEtaNonProp(0),
- fHistDeltaEtaPt(0),
- fHistDeltaPhiPt(0),
- fHistDeltaEtaNewProp(0),
- fHistDeltaPhiNewProp(0),
- fHistClusPhiEtaEnergy(0),
+ fParticleLevel(kFALSE),
+ fIsMC(kFALSE),
+ fCentMethod2(""),
+ fCentMethod3(""),
+ fDoV0QA(0),
+ fDoEPQA(0),
+ fDoLeadingObjectPosition(0),
+ fMaxCellsInCluster(30),
+ fCent2(0),
+ fCent3(0),
+ fVZERO(0),
+ fV0ATotMult(0),
+ fV0CTotMult(0),
+ fHistEventQA(0),
+ fHistTrNegativeLabels(0),
+ fHistTrZeroLabels(0),
fHistNCellsEnergy(0),
+ fHistFcrossEnergy(0),
fHistClusTimeEnergy(0),
- fHistCellsEnergy(0),
+ fHistClusEnergyMinusCellEnergy(0),
+ fHistCellsAbsIdEnergy(0),
fHistChVSneCells(0),
fHistChVSneClus(0),
fHistChVSneCorrCells(0)
{
// Default constructor.
- for (Int_t i = 0; i < 5; i++) {
- fHistTrackPhi[i] = 0;
- fHistTrackEta[i] = 0;
- }
-
for (Int_t i = 0; i < 4; i++) {
+ for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
+ fHistTrPhiEtaZeroLab[i] = 0;
+ fHistTrPtZeroLab[i] = 0;
+ fHistTrEmcPhiEta[i] = 0;
+ fHistTrEmcPt[i] = 0;
+ fHistTrPhiEtaNonProp[i] = 0;
+ fHistTrPtNonProp[i] = 0;
+ fHistDeltaEtaPt[i] = 0;
+ fHistDeltaPhiPt[i] = 0;
+ fHistDeltaPtvsPt[i] = 0;
+ fHistClusPhiEtaEnergy[i] = 0;
+ fHistClusDeltaPhiEPEnergy[i] = 0;
+ fHistClusMCEnergyFraction[i] = 0;
fHistJetsPhiEta[i] = 0;
- fHistJetsPtNonBias[i] = 0;
- fHistJetsPtTrack[i] = 0;
- fHistJetsPtClus[i] = 0;
- fHistJetsPt[i] = 0;
fHistJetsPtArea[i] = 0;
- fHistJetsPtAreaNonBias[i] = 0;
}
+
+ SetMakeGeneralHistograms(kTRUE);
}
//________________________________________________________________________
AliAnalysisTaskSAQA::AliAnalysisTaskSAQA(const char *name) :
AliAnalysisTaskEmcalJet(name, kTRUE),
fCellEnergyCut(0.1),
- fDoTrigger(kFALSE),
- fRepropagateTracks(kFALSE),
- fTrgClusName("ClustersL1GAMMAFEE"),
- fTrgClusters(0),
- fHistCentrality(0),
- fHistZVertex(0),
- fHistTracksCent(0),
- fHistClusCent(0),
- fHistClusTracks(0),
- fHistCellsCent(0),
- fHistCellsTracks(0),
- fHistMaxL1FastORCent(0),
- fHistMaxL1ClusCent(0),
- fHistMaxL1ThrCent(0),
- fHistTracksPt(0),
- fHistTrPhiEta(0),
- fHistTrEmcPhiEta(0),
- fHistTrPhiEtaNonProp(0),
- fHistDeltaEtaPt(0),
- fHistDeltaPhiPt(0),
- fHistDeltaEtaNewProp(0),
- fHistDeltaPhiNewProp(0),
- fHistClusPhiEtaEnergy(0),
+ fParticleLevel(kFALSE),
+ fIsMC(kFALSE),
+ fCentMethod2(""),
+ fCentMethod3(""),
+ fDoV0QA(0),
+ fDoEPQA(0),
+ fDoLeadingObjectPosition(0),
+ fMaxCellsInCluster(30),
+ fCent2(0),
+ fCent3(0),
+ fVZERO(0),
+ fV0ATotMult(0),
+ fV0CTotMult(0),
+ fHistEventQA(0),
+ fHistTrNegativeLabels(0),
+ fHistTrZeroLabels(0),
fHistNCellsEnergy(0),
+ fHistFcrossEnergy(0),
fHistClusTimeEnergy(0),
- fHistCellsEnergy(0),
+ fHistClusEnergyMinusCellEnergy(0),
+ fHistCellsAbsIdEnergy(0),
fHistChVSneCells(0),
fHistChVSneClus(0),
fHistChVSneCorrCells(0)
{
// Standard constructor.
- for (Int_t i = 0; i < 5; i++) {
- fHistTrackPhi[i] = 0;
- fHistTrackEta[i] = 0;
- }
-
for (Int_t i = 0; i < 4; i++) {
+ for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
+ fHistTrPhiEtaZeroLab[i] = 0;
+ fHistTrPtZeroLab[i] = 0;
+ fHistTrEmcPhiEta[i] = 0;
+ fHistTrEmcPt[i] = 0;
+ fHistTrPhiEtaNonProp[i] = 0;
+ fHistTrPtNonProp[i] = 0;
+ fHistDeltaEtaPt[i] = 0;
+ fHistDeltaPhiPt[i] = 0;
+ fHistDeltaPtvsPt[i] = 0;
+ fHistClusPhiEtaEnergy[i] = 0;
+ fHistClusDeltaPhiEPEnergy[i] = 0;
+ fHistClusMCEnergyFraction[i] = 0;
fHistJetsPhiEta[i] = 0;
- fHistJetsPtNonBias[i] = 0;
- fHistJetsPtTrack[i] = 0;
- fHistJetsPtClus[i] = 0;
- fHistJetsPt[i] = 0;
fHistJetsPtArea[i] = 0;
- fHistJetsPtAreaNonBias[i] = 0;
}
+
+ SetMakeGeneralHistograms(kTRUE);
}
//________________________________________________________________________
{
// Create histograms
- OpenFile(1);
- fOutput = new TList();
- fOutput->SetOwner();
-
- fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", 100, 0, 100);
- fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
- fHistCentrality->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistCentrality);
-
- fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
- fHistZVertex->GetXaxis()->SetTitle("z");
- fHistZVertex->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistZVertex);
-
- fHistTracksCent = new TH2F("fHistTracksCent","Tracks vs. centrality", 100, 0, 100, fNbins, 0, 4000);
- fHistTracksCent->GetXaxis()->SetTitle("Centrality (%)");
- fHistTracksCent->GetYaxis()->SetTitle("No. of tracks");
- fOutput->Add(fHistTracksCent);
-
- if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
- fHistClusCent = new TH2F("fHistClusCent","Clusters vs. centrality", 100, 0, 100, fNbins, 0, 2000);
- fHistClusCent->GetXaxis()->SetTitle("Centrality (%)");
- fHistClusCent->GetYaxis()->SetTitle("No. of clusters");
- fOutput->Add(fHistClusCent);
-
- fHistClusTracks = new TH2F("fHistClusTracks","Clusters vs. tracks", fNbins, 0, 4000, fNbins, 0, 2000);
- fHistClusTracks->GetXaxis()->SetTitle("No. of tracks");
- fHistClusTracks->GetYaxis()->SetTitle("No. of clusters");
- fOutput->Add(fHistClusTracks);
-
- fHistCellsCent = new TH2F("fHistCellsCent","Cells vs. centrality", 100, 0, 100, fNbins, 0, 6000);
- fHistCellsCent->GetXaxis()->SetTitle("Centrality (%)");
- fHistCellsCent->GetYaxis()->SetTitle("No. of EMCal cells");
- fOutput->Add(fHistCellsCent);
-
- fHistCellsTracks = new TH2F("fHistCellsTracks","Cells vs. tracks", fNbins, 0, 4000, fNbins, 0, 6000);
- fHistCellsTracks->GetXaxis()->SetTitle("No. of tracks");
- fHistCellsTracks->GetYaxis()->SetTitle("No. of EMCal cells");
- fOutput->Add(fHistCellsTracks);
-
- fHistClusTimeEnergy = new TH2F("fHistClusTimeEnergy","Time vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1e-6);
- fHistClusTimeEnergy->GetXaxis()->SetTitle("Energy (GeV)");
- fHistClusTimeEnergy->GetYaxis()->SetTitle("Time");
- fOutput->Add(fHistClusTimeEnergy);
+ AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
+
+ if (fParticleCollArray.GetEntriesFast()>0) {
+ if (!fParticleLevel && fIsMC) {
+ fHistTrNegativeLabels = new TH1F("fHistTrNegativeLabels","fHistTrNegativeLabels", 500, 0, 1);
+ fHistTrNegativeLabels->GetXaxis()->SetTitle("% of negative labels");
+ fHistTrNegativeLabels->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrNegativeLabels);
- if (fDoTrigger) {
- fHistMaxL1FastORCent = new TH2F("fHistMaxL1FastORCent","fHistMaxL1ClusCent", 100, 0, 100, 250, 0, 250);
- fHistMaxL1FastORCent->GetXaxis()->SetTitle("Centrality [%]");
- fHistMaxL1FastORCent->GetYaxis()->SetTitle("Maximum L1 FastOR");
- fOutput->Add(fHistMaxL1FastORCent);
-
- fHistMaxL1ClusCent = new TH2F("fHistMaxL1ClusCent","fHistMaxL1ClusCent", 100, 0, 100, 250, 0, 250);
- fHistMaxL1ClusCent->GetXaxis()->SetTitle("Centrality [%]");
- fHistMaxL1ClusCent->GetYaxis()->SetTitle("Maximum L1 trigger cluster");
- fOutput->Add(fHistMaxL1ClusCent);
-
- fHistMaxL1ThrCent = new TH2F("fHistMaxL1ThrCent","fHistMaxL1ThrCent", 100, 0, 100, 250, 0, 250);
- fHistMaxL1ThrCent->GetXaxis()->SetTitle("Centrality [%]");
- fHistMaxL1ThrCent->GetYaxis()->SetTitle("Maximum L1 threshold");
- fOutput->Add(fHistMaxL1ThrCent);
+ fHistTrZeroLabels = new TH1F("fHistTrZeroLabels","fHistTrZeroLabels", 500, 0, 1);
+ fHistTrZeroLabels->GetXaxis()->SetTitle("% of negative labels");
+ fHistTrZeroLabels->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrZeroLabels);
}
- }
- fHistTracksPt = new TH1F("fHistTracksPt","p_{T} spectrum of reconstructed tracks", fNbins, fMinBinPt, fMaxBinPt);
- fHistTracksPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistTracksPt->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistTracksPt);
-
- fHistTrPhiEta = new TH2F("fHistTrPhiEta","Phi-Eta distribution of tracks", 80, -2, 2, 128, 0, 6.4);
- fHistTrPhiEta->GetXaxis()->SetTitle("#eta");
- fHistTrPhiEta->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistTrPhiEta);
-
- fHistTrEmcPhiEta = new TH2F("fHistTrEmcPhiEta","Phi-Eta emcal distribution of tracks", 80, -2, 2, 128, 0, 6.4);
- fHistTrEmcPhiEta->GetXaxis()->SetTitle("#eta");
- fHistTrEmcPhiEta->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistTrEmcPhiEta);
-
- fHistTrPhiEtaNonProp = new TH2F("fHistTrPhiEtaNonProp","fHistTrPhiEtaNonProp", 80, -2, 2, 128, 0, 6.4);
- fHistTrPhiEtaNonProp->GetXaxis()->SetTitle("#eta");
- fHistTrPhiEtaNonProp->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistTrPhiEtaNonProp);
-
- fHistDeltaEtaPt = new TH2F("fHistDeltaEtaPt","fHistDeltaEtaPt", fNbins, fMinBinPt, fMaxBinPt, 80, -0.5, 0.5);
- fHistDeltaEtaPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistDeltaEtaPt->GetYaxis()->SetTitle("#delta#eta");
- fOutput->Add(fHistDeltaEtaPt);
-
- fHistDeltaPhiPt = new TH2F("fHistDeltaPhiPt","fHistDeltaPhiPt", fNbins, fMinBinPt, fMaxBinPt, 256, -1.6, 4.8);
- fHistDeltaPhiPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistDeltaPhiPt->GetYaxis()->SetTitle("#delta#phi");
- fOutput->Add(fHistDeltaPhiPt);
-
- if (fRepropagateTracks) {
- fHistDeltaEtaNewProp = new TH1F("fHistDeltaEtaNewProp","fHistDeltaEtaNewProp", 800, -0.4, 0.4);
- fHistDeltaEtaNewProp->GetXaxis()->SetTitle("#delta#eta");
- fHistDeltaEtaNewProp->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaEtaNewProp);
-
- fHistDeltaPhiNewProp = new TH1F("fHistDeltaPhiNewProp","fHistDeltaPhiNewProp", 2560, -1.6, 3.2);
- fHistDeltaPhiNewProp->GetXaxis()->SetTitle("#delta#phi");
- fHistDeltaPhiNewProp->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPhiNewProp);
+ TString histname;
+
+ Int_t nlabels = 4;
+ if (fParticleLevel)
+ nlabels = 1;
+
+ for (Int_t i = 0; i < fNcentBins; i++) {
+ for (Int_t j = 0; j < nlabels; j++) {
+ histname = Form("fHistTrPhiEtaPt_%d_%d",i,j);
+ fHistTrPhiEtaPt[i][j] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt);
+ fHistTrPhiEtaPt[i][j]->GetXaxis()->SetTitle("#eta");
+ fHistTrPhiEtaPt[i][j]->GetYaxis()->SetTitle("#phi");
+ fHistTrPhiEtaPt[i][j]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
+ fOutput->Add(fHistTrPhiEtaPt[i][j]);
+ }
+
+ if (!fParticleLevel) {
+ if (fIsMC) {
+ histname = Form("fHistTrPhiEtaZeroLab_%d",i);
+ fHistTrPhiEtaZeroLab[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
+ fHistTrPhiEtaZeroLab[i]->GetXaxis()->SetTitle("#eta");
+ fHistTrPhiEtaZeroLab[i]->GetYaxis()->SetTitle("#phi");
+ fHistTrPhiEtaZeroLab[i]->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrPhiEtaZeroLab[i]);
+
+ histname = Form("fHistTrPtZeroLab_%d",i);
+ fHistTrPtZeroLab[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
+ fHistTrPtZeroLab[i]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistTrPtZeroLab[i]->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrPtZeroLab[i]);
+ }
+
+ histname = Form("fHistTrEmcPhiEta_%d",i);
+ fHistTrEmcPhiEta[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
+ fHistTrEmcPhiEta[i]->GetXaxis()->SetTitle("#eta");
+ fHistTrEmcPhiEta[i]->GetYaxis()->SetTitle("#phi");
+ fOutput->Add(fHistTrEmcPhiEta[i]);
+
+ histname = Form("fHistTrEmcPt_%d",i);
+ fHistTrEmcPt[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
+ fHistTrEmcPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistTrEmcPt[i]->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrEmcPt[i]);
+
+ histname = Form("fHistTrPhiEtaNonProp_%d",i);
+ fHistTrPhiEtaNonProp[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
+ fHistTrPhiEtaNonProp[i]->GetXaxis()->SetTitle("#eta");
+ fHistTrPhiEtaNonProp[i]->GetYaxis()->SetTitle("#phi");
+ fOutput->Add(fHistTrPhiEtaNonProp[i]);
+
+ histname = Form("fHistTrPtNonProp_%d",i);
+ fHistTrPtNonProp[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
+ fHistTrPtNonProp[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistTrPtNonProp[i]->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistTrPtNonProp[i]);
+
+ histname = Form("fHistDeltaEtaPt_%d",i);
+ fHistDeltaEtaPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 50, -0.5, 0.5);
+ fHistDeltaEtaPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistDeltaEtaPt[i]->GetYaxis()->SetTitle("#delta#eta");
+ fOutput->Add(fHistDeltaEtaPt[i]);
+
+ histname = Form("fHistDeltaPhiPt_%d",i);
+ fHistDeltaPhiPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 200, -2, 2);
+ fHistDeltaPhiPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistDeltaPhiPt[i]->GetYaxis()->SetTitle("#delta#phi");
+ fOutput->Add(fHistDeltaPhiPt[i]);
+
+ histname = Form("fHistDeltaPtvsPt_%d",i);
+ fHistDeltaPtvsPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, fNbins, -fMaxBinPt/2, fMaxBinPt/2);
+ fHistDeltaPtvsPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistDeltaPtvsPt[i]->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
+ fHistDeltaPtvsPt[i]->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistDeltaPtvsPt[i]);
+ }
+ }
}
- if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
- fHistClusPhiEtaEnergy = new TH3F("fHistClusPhiEtaEnergy","Phi-Eta-Energy distribution of clusters", fNbins, fMinBinPt, fMaxBinPt, 80, -2, 2, 128, 0, 6.4);
- fHistClusPhiEtaEnergy->GetXaxis()->SetTitle("E [GeV]");
- fHistClusPhiEtaEnergy->GetYaxis()->SetTitle("#eta");
- fHistClusPhiEtaEnergy->GetZaxis()->SetTitle("#phi");
- fOutput->Add(fHistClusPhiEtaEnergy);
+ if (fClusterCollArray.GetEntriesFast()>0) {
+ TString histname;
+
+ for (Int_t i = 0; i < fNcentBins; i++) {
+ histname = "fHistClusPhiEtaEnergy_";
+ histname += i;
+ fHistClusPhiEtaEnergy[i] = new TH3F(histname, histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt);
+ fHistClusPhiEtaEnergy[i]->GetXaxis()->SetTitle("#eta");
+ fHistClusPhiEtaEnergy[i]->GetYaxis()->SetTitle("#phi");
+ fHistClusPhiEtaEnergy[i]->GetZaxis()->SetTitle("E_{cluster} (GeV)");
+ fOutput->Add(fHistClusPhiEtaEnergy[i]);
+
+ histname = "fHistClusDeltaPhiEPEnergy_";
+ histname += i;
+ fHistClusDeltaPhiEPEnergy[i] = new TH2F(histname, histname, fNbins, fMinBinPt, fMaxBinPt, 100, 0, TMath::Pi());
+ fHistClusDeltaPhiEPEnergy[i]->GetXaxis()->SetTitle("E_{cluster} (GeV)");
+ fHistClusDeltaPhiEPEnergy[i]->GetYaxis()->SetTitle("#phi_{cluster} - #psi_{RP}");
+ fOutput->Add(fHistClusDeltaPhiEPEnergy[i]);
+
+ if (fIsEmbedded) {
+ histname = "fHistClusMCEnergyFraction_";
+ histname += i;
+ fHistClusMCEnergyFraction[i] = new TH1F(histname, histname, fNbins, 0, 1.2);
+ fHistClusMCEnergyFraction[i]->GetXaxis()->SetTitle("MC fraction");
+ fHistClusMCEnergyFraction[i]->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistClusMCEnergyFraction[i]);
+ }
+ }
+
+ fHistClusTimeEnergy = new TH2F("fHistClusTimeEnergy","Time vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, fNbins, -1e-6, 1e-6);
+ fHistClusTimeEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
+ fHistClusTimeEnergy->GetYaxis()->SetTitle("Time");
+ fOutput->Add(fHistClusTimeEnergy);
- fHistNCellsEnergy = new TH2F("fHistNCellsEnergy","Number of cells vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, 30, 0, 30);
- fHistNCellsEnergy->GetXaxis()->SetTitle("E [GeV]");
+ Int_t nbins = fMaxCellsInCluster;
+ while (nbins > fNbins) nbins /= 2;
+ fHistNCellsEnergy = new TH2F("fHistNCellsEnergy","Number of cells vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, nbins, -0.5, fMaxCellsInCluster-0.5);
+ fHistNCellsEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
fHistNCellsEnergy->GetYaxis()->SetTitle("N_{cells}");
fOutput->Add(fHistNCellsEnergy);
- }
- if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
-
- fHistCellsEnergy = new TH1F("fHistCellsEnergy","Energy spectrum of cells", fNbins, fMinBinPt, fMaxBinPt);
- fHistCellsEnergy->GetXaxis()->SetTitle("E [GeV]");
- fHistCellsEnergy->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistCellsEnergy);
+ fHistFcrossEnergy = new TH2F("fHistFcrossEnergy","fHistFcrossEnergy", fNbins, fMinBinPt, fMaxBinPt, 200, -3.5, 1.5);
+ fHistFcrossEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
+ fHistFcrossEnergy->GetYaxis()->SetTitle("F_{cross}");
+ fOutput->Add(fHistFcrossEnergy);
+
+ fHistClusEnergyMinusCellEnergy = new TH2F("fHistClusEnergyMinusCellEnergy","fHistClusEnergyMinusCellEnergy",
+ fNbins, fMinBinPt, fMaxBinPt, fNbins, -fMaxBinPt/2, fMaxBinPt/2);
+ fHistClusEnergyMinusCellEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
+ fHistClusEnergyMinusCellEnergy->GetYaxis()->SetTitle("E_{cluster} - #Sigma_{i}E_{cell,i} (GeV)");
+ fOutput->Add(fHistClusEnergyMinusCellEnergy);
+
+ fHistCellsAbsIdEnergy = new TH2F("fHistCellsAbsIdEnergy","fHistCellsAbsIdEnergy", 11600,0,11599,(Int_t)(fNbins / 2), fMinBinPt, fMaxBinPt / 2);
+ fHistCellsAbsIdEnergy->GetXaxis()->SetTitle("cell abs. Id");
+ fHistCellsAbsIdEnergy->GetYaxis()->SetTitle("E_{cluster} (GeV)");
+ fHistCellsAbsIdEnergy->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistCellsAbsIdEnergy);
fHistChVSneCells = new TH2F("fHistChVSneCells","Charged energy vs. neutral (cells) energy",
(Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistChVSneCells->GetXaxis()->SetTitle("E [GeV]");
- fHistChVSneCells->GetYaxis()->SetTitle("P [GeV/c]");
+ fHistChVSneCells->GetXaxis()->SetTitle("Energy (GeV)");
+ fHistChVSneCells->GetYaxis()->SetTitle("Momentum (GeV/c)");
fOutput->Add(fHistChVSneCells);
fHistChVSneClus = new TH2F("fHistChVSneClus","Charged energy vs. neutral (clusters) energy",
(Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistChVSneClus->GetXaxis()->SetTitle("E [GeV]");
- fHistChVSneClus->GetYaxis()->SetTitle("P [GeV/c]");
+ fHistChVSneClus->GetXaxis()->SetTitle("Energy (GeV)");
+ fHistChVSneClus->GetYaxis()->SetTitle("Momentum (GeV/c)");
fOutput->Add(fHistChVSneClus);
fHistChVSneCorrCells = new TH2F("fHistChVSneCorrCells","Charged energy vs. neutral (corrected cells) energy",
(Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt , fMaxBinPt * 2.5);
- fHistChVSneCorrCells->GetXaxis()->SetTitle("E [GeV]");
- fHistChVSneCorrCells->GetYaxis()->SetTitle("P [GeV/c]");
+ fHistChVSneCorrCells->GetXaxis()->SetTitle("Energy (GeV)");
+ fHistChVSneCorrCells->GetYaxis()->SetTitle("Momentum (GeV/c)");
fOutput->Add(fHistChVSneCorrCells);
}
-
- for (Int_t i = 0; i < 5; i++) {
- TString histnamephi("fHistTrackPhi_");
- histnamephi += i;
- fHistTrackPhi[i] = new TH1F(histnamephi.Data(),histnamephi.Data(), 128, 0, 6.4);
- fHistTrackPhi[i]->GetXaxis()->SetTitle("Phi");
- fOutput->Add(fHistTrackPhi[i]);
-
- TString histnameeta("fHistTrackEta_");
- histnameeta += i;
- fHistTrackEta[i] = new TH1F(histnameeta.Data(),histnameeta.Data(), 100, -2, 2);
- fHistTrackEta[i]->GetXaxis()->SetTitle("Eta");
- fOutput->Add(fHistTrackEta[i]);
- }
-
- fHistTrackPhi[0]->SetLineColor(kRed);
- fHistTrackEta[0]->SetLineColor(kRed);
- fHistTrackPhi[1]->SetLineColor(kBlue);
- fHistTrackEta[1]->SetLineColor(kBlue);
- fHistTrackPhi[2]->SetLineColor(kGreen);
- fHistTrackEta[2]->SetLineColor(kGreen);
- fHistTrackPhi[3]->SetLineColor(kOrange);
- fHistTrackEta[3]->SetLineColor(kOrange);
- fHistTrackPhi[4]->SetLineColor(kBlack);
- fHistTrackEta[4]->SetLineColor(kBlack);
-
- if (!fJetsName.IsNull()) {
+
+ if (fJetCollArray.GetEntriesFast()>0) {
TString histname;
- for (Int_t i = 0; i < 4; i++) {
+ for (Int_t i = 0; i < fNcentBins; i++) {
histname = "fHistJetsPhiEta_";
histname += i;
- fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 80, -2, 2, 128, 0, 6.4);
+ fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta");
fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi");
- fHistJetsPhiEta[i]->GetZaxis()->SetTitle("p_{T} [GeV/c]");
+ fHistJetsPhiEta[i]->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJetsPhiEta[i]);
- histname = "fHistJetsPtNonBias_";
+ histname = "fHistJetsPtArea_";
histname += i;
- fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistJetsPtNonBias[i]);
+ fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 50, 0, 1.5);
+ fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
+ fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
+ fOutput->Add(fHistJetsPtArea[i]);
+ }
+ }
+
+ Int_t dim = 0;
+ TString title[20];
+ Int_t nbins[20] = {0};
+ Double_t min[20] = {0};
+ Double_t max[20] = {0};
+
+ if (fForceBeamType != AliAnalysisTaskEmcal::kpp) {
+ title[dim] = "Centrality %";
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = 101;
+ dim++;
+
+ if (!fCentMethod2.IsNull()) {
+ title[dim] = Form("Centrality %s %%", fCentMethod2.Data());
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = 101;
+ dim++;
+ }
- histname = "fHistJetsPtTrack_";
- histname += i;
- fHistJetsPtTrack[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistJetsPtTrack[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPtTrack[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistJetsPtTrack[i]);
+ if (!fCentMethod3.IsNull()) {
+ title[dim] = Form("Centrality %s %%", fCentMethod3.Data());
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = 101;
+ dim++;
+ }
+
+ if (fDoV0QA==1) {
+ title[dim] = "V0A total multiplicity";
+ nbins[dim] = 200;
+ min[dim] = 0;
+ max[dim] = 20000;
+ dim++;
+
+ title[dim] = "V0C total multiplicity";
+ nbins[dim] = 200;
+ min[dim] = 0;
+ max[dim] = 20000;
+ dim++;
+ }
+ else if (fDoV0QA==2) {
+ title[dim] = "V0A+V0C total multiplicity";
+ nbins[dim] = 300;
+ min[dim] = 0;
+ max[dim] = 30000;
+ dim++;
+ }
- if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
- histname = "fHistJetsPtClus_";
- histname += i;
- fHistJetsPtClus[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistJetsPtClus[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPtClus[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistJetsPtClus[i]);
- }
+ if (!fRhoName.IsNull()) {
+ title[dim] = "#rho (GeV/c)";
+ nbins[dim] = fNbins*4;
+ min[dim] = 0;
+ max[dim] = 400;
+ dim++;
+ }
- histname = "fHistJetsPt_";
- histname += i;
- fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
- fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistJetsPt[i]);
+ if (fDoEPQA) {
+ title[dim] = "#psi_{RP}";
+ nbins[dim] = 200;
+ min[dim] = -TMath::Pi();
+ max[dim] = TMath::Pi();
+ dim++;
+ }
+ }
- histname = "fHistJetsPtAreaNonBias_";
- histname += i;
- fHistJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
- fHistJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
- fOutput->Add(fHistJetsPtAreaNonBias[i]);
+ if (fParticleCollArray.GetEntriesFast()>0) {
+ title[dim] = "No. of tracks";
+ nbins[dim] = 3000;
+ min[dim] = -0.5;
+ max[dim] = 6000-0.5;
+ dim++;
+
+ title[dim] = "p_{T,track}^{leading} (GeV/c)";
+ nbins[dim] = fNbins;
+ min[dim] = fMinBinPt;
+ max[dim] = fMaxBinPt;
+ dim++;
+
+ if (fDoLeadingObjectPosition) {
+ title[dim] = "#eta_{track}^{leading}";
+ nbins[dim] = 100;
+ min[dim] = -1;
+ max[dim] = 1;
+ dim++;
+
+ title[dim] = "#phi_{track}^{leading}";
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = TMath::Pi() * 2.02;
+ dim++;
+ }
+ }
- histname = "fHistJetsPtArea_";
- histname += i;
- fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
- fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
- fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
- fOutput->Add(fHistJetsPtArea[i]);
+ if (fClusterCollArray.GetEntriesFast()>0) {
+ title[dim] = "No. of clusters";
+ nbins[dim] = 2000;
+ min[dim] = 0;
+ max[dim] = 4000-0.5;
+ dim++;
+
+ title[dim] = "E_{cluster}^{leading} (GeV)";
+ nbins[dim] = fNbins;
+ min[dim] = fMinBinPt;
+ max[dim] = fMaxBinPt;
+ dim++;
+
+ if (fDoLeadingObjectPosition) {
+ title[dim] = "#eta_{cluster}^{leading}";
+ nbins[dim] = 100;
+ min[dim] = -1;
+ max[dim] = 1;
+ dim++;
+
+ title[dim] = "#phi_{cluster}^{leading}";
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = TMath::Pi() * 2.02;
+ dim++;
}
}
+ if (!fCaloCellsName.IsNull()) {
+ title[dim] = "No. of cells";
+ nbins[dim] = 3000;
+ min[dim] = 0;
+ max[dim] = 6000-0.5;
+ dim++;
+ }
+
+ if (fJetCollArray.GetEntriesFast()>0) {
+ title[dim] = "No. of jets";
+ nbins[dim] = 200;
+ min[dim] = 0;
+ max[dim] = 200-0.5;
+ dim++;
+
+ title[dim] = "p_{T,jet}^{leading} (GeV/c)";
+ nbins[dim] = fNbins;
+ min[dim] = fMinBinPt;
+ max[dim] = fMaxBinPt;
+ dim++;
+
+ if (fDoLeadingObjectPosition) {
+ title[dim] = "#eta_{jet}^{leading}";
+ nbins[dim] = 100;
+ min[dim] = -1;
+ max[dim] = 1;
+ dim++;
+
+ title[dim] = "#phi_{jet}^{leading}";
+ nbins[dim] = 101;
+ min[dim] = 0;
+ max[dim] = TMath::Pi() * 2.02;
+ dim++;
+ }
+ }
+
+ fHistEventQA = new THnSparseF("fHistEventQA","fHistEventQA",dim,nbins,min,max);
+ for (Int_t i = 0; i < dim; i++)
+ fHistEventQA->GetAxis(i)->SetTitle(title[i]);
+ fOutput->Add(fHistEventQA);
+
PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
}
+//________________________________________________________________________
+void AliAnalysisTaskSAQA::ExecOnce()
+{
+ AliAnalysisTaskEmcalJet::ExecOnce();
+
+ if (fDoV0QA) {
+ fVZERO = InputEvent()->GetVZEROData();
+ if (!fVZERO) {
+ AliError("AliVVZERO not available");
+ }
+ }
+}
+
//________________________________________________________________________
Bool_t AliAnalysisTaskSAQA::RetrieveEventObjects()
{
if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
return kFALSE;
- if (!fTrgClusName.IsNull() && fDoTrigger && !fTrgClusters) {
- fTrgClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTrgClusName));
- if (!fTrgClusters) {
- AliError(Form("%s: Could not retrieve trigger clusters %s!", GetName(), fTrgClusName.Data()));
- return kFALSE;
- }
- else {
- TClass *cl = fTrgClusters->GetClass();
- if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
- AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fTrgClusName.Data()));
- fTrgClusters = 0;
- return kFALSE;
+ if (!fCentMethod2.IsNull() || !fCentMethod3.IsNull()) {
+ if (fBeamType == kAA || fBeamType == kpA ) {
+ AliCentrality *aliCent = InputEvent()->GetCentrality();
+ if (aliCent) {
+ if (!fCentMethod2.IsNull())
+ fCent2 = aliCent->GetCentralityPercentile(fCentMethod2);
+ if (!fCentMethod3.IsNull())
+ fCent3 = aliCent->GetCentralityPercentile(fCentMethod3);
}
}
}
+ if (fVZERO) {
+ fV0ATotMult = AliESDUtils::GetCorrV0A(fVZERO->GetMTotV0A(),fVertex[2]);
+ fV0CTotMult = AliESDUtils::GetCorrV0C(fVZERO->GetMTotV0C(),fVertex[2]);
+ }
+
return kTRUE;
}
+
//________________________________________________________________________
Bool_t AliAnalysisTaskSAQA::FillHistograms()
{
// Fill histograms.
- fHistCentrality->Fill(fCent);
- if (fTracks)
- fHistTracksCent->Fill(fCent, fTracks->GetEntriesFast());
- if (fCaloClusters)
- fHistClusCent->Fill(fCent, fCaloClusters->GetEntriesFast());
+ Float_t trackSum = 0;
+ Float_t clusSum = 0;
+ Float_t cellSum = 0;
+ Float_t cellCutSum = 0;
- fHistZVertex->Fill(fVertex[2]);
+ Int_t ntracks = 0;
+ Int_t nclusters = 0;
+ Int_t ncells = 0;
+ Int_t njets = 0;
- Float_t trackSum = DoTrackLoop();
+ Float_t leadingClusE = 0;
+ Float_t leadingClusEta = 0;
+ Float_t leadingClusPhi = 0;
- DoJetLoop();
+ Float_t leadingTrackPt = 0;
+ Float_t leadingTrackEta = 0;
+ Float_t leadingTrackPhi = 0;
- if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
+ Float_t leadingJetPt = 0;
+ Float_t leadingJetEta = 0;
+ Float_t leadingJetPhi = 0;
+
+ if (fTracks) {
+ AliVParticle *leadingTrack = 0;
- if (fTracks && fCaloClusters)
- fHistClusTracks->Fill(fTracks->GetEntriesFast(), fCaloClusters->GetEntriesFast());
+ ntracks = DoTrackLoop(trackSum, leadingTrack);
+ AliDebug(2,Form("%d tracks found in the event", ntracks));
- Float_t clusSum = DoClusterLoop();
+ if (leadingTrack) {
+ leadingTrackPt = leadingTrack->Pt();
+ leadingTrackEta = leadingTrack->Eta();
+ leadingTrackPhi = leadingTrack->Phi();
+ }
+ }
- Float_t cellSum = 0, cellCutSum = 0;
-
- Int_t ncells = DoCellLoop(cellSum, cellCutSum);
+ if (fCaloClusters) {
+ AliVCluster *leadingClus = 0;
+
+ nclusters = DoClusterLoop(clusSum, leadingClus);
+ AliDebug(2,Form("%d clusters found in the event", nclusters));
- if (fTracks)
- fHistCellsTracks->Fill(fTracks->GetEntriesFast(), ncells);
+ if (leadingClus) {
+ TLorentzVector leadingClusVect;
+ leadingClus->GetMomentum(leadingClusVect, fVertex);
+ leadingClusE = leadingClus->E();
+ leadingClusEta = leadingClusVect.Eta();
+ leadingClusPhi = leadingClusVect.Phi();
+ }
- fHistCellsCent->Fill(fCent, ncells);
+ fHistChVSneClus->Fill(clusSum, trackSum);
+ }
+
+ if (fCaloCells) {
+ ncells = DoCellLoop(cellSum, cellCutSum);
+ AliDebug(2,Form("%d cells found in the event", ncells));
fHistChVSneCells->Fill(cellSum, trackSum);
- fHistChVSneClus->Fill(clusSum, trackSum);
fHistChVSneCorrCells->Fill(cellCutSum, trackSum);
+ }
- if (fDoTrigger) {
- Float_t maxTrgClus = DoTriggerClusLoop();
- fHistMaxL1ClusCent->Fill(fCent, maxTrgClus);
-
- Int_t maxL1amp = -1;
- Int_t maxL1thr = -1;
-
- DoTriggerPrimitives(maxL1amp, maxL1thr);
-
- if (maxL1amp > -1)
- fHistMaxL1FastORCent->Fill(fCent, maxL1amp);
-
- if (maxL1thr > -1)
- fHistMaxL1ThrCent->Fill(fCent, maxL1thr);
+ if (fJets) {
+ AliEmcalJet *leadingJet = 0;
+
+ njets = DoJetLoop(leadingJet);
+ AliDebug(2,Form("%d jets found in the event", njets));
+
+ if (leadingJet) {
+ leadingJetPt = leadingJet->Pt();
+ leadingJetEta = leadingJet->Eta();
+ leadingJetPhi = leadingJet->Phi();
}
}
+ FillEventQAHisto(fCent, fCent2, fCent3, fV0ATotMult, fV0CTotMult, fEPV0, fRhoVal,
+ ntracks, nclusters, ncells, njets,
+ leadingTrackPt, leadingTrackEta, leadingTrackPhi,
+ leadingClusE, leadingClusEta, leadingClusPhi,
+ leadingJetPt, leadingJetEta, leadingJetPhi);
+
return kTRUE;
}
+//________________________________________________________________________
+void AliAnalysisTaskSAQA::FillEventQAHisto(Float_t cent, Float_t cent2, Float_t cent3, Float_t v0a, Float_t v0c,
+ Float_t ep, Float_t rho, Int_t ntracks, Int_t nclusters, Int_t ncells, Int_t njets,
+ Float_t maxTrackPt, Float_t maxTrackEta, Float_t maxTrackPhi,
+ Float_t maxClusterE, Float_t maxClusterEta, Float_t maxClusterPhi,
+ Float_t maxJetPt, Float_t maxJetEta, Float_t maxJetPhi)
+{
+ Double_t contents[20]={0};
+
+ for (Int_t i = 0; i < fHistEventQA->GetNdimensions(); i++) {
+ TString title(fHistEventQA->GetAxis(i)->GetTitle());
+ if (title=="Centrality %")
+ contents[i] = cent;
+ else if (title==Form("Centrality %s %%", fCentMethod2.Data()))
+ contents[i] = cent2;
+ else if (title==Form("Centrality %s %%", fCentMethod3.Data()))
+ contents[i] = cent3;
+ else if (title=="V0A total multiplicity")
+ contents[i] = v0a;
+ else if (title=="V0C total multiplicity")
+ contents[i] = v0c;
+ else if (title=="V0A+V0C total multiplicity")
+ contents[i] = v0a+v0c;
+ else if (title=="#psi_{RP}")
+ contents[i] = ep;
+ else if (title=="#rho (GeV/c)")
+ contents[i] = rho;
+ else if (title=="No. of tracks")
+ contents[i] = ntracks;
+ else if (title=="No. of clusters")
+ contents[i] = nclusters;
+ else if (title=="No. of cells")
+ contents[i] = ncells;
+ else if (title=="No. of jets")
+ contents[i] = njets;
+ else if (title=="p_{T,track}^{leading} (GeV/c)")
+ contents[i] = maxTrackPt;
+ else if (title=="#eta_{track}^{leading}")
+ contents[i] = maxTrackEta;
+ else if (title=="#phi_{track}^{leading}")
+ contents[i] = maxTrackPhi;
+ else if (title=="E_{cluster}^{leading} (GeV)")
+ contents[i] = maxClusterE;
+ else if (title=="#eta_{cluster}^{leading}")
+ contents[i] = maxClusterEta;
+ else if (title=="#phi_{cluster}^{leading}")
+ contents[i] = maxClusterPhi;
+ else if (title=="p_{T,jet}^{leading} (GeV/c)")
+ contents[i] = maxJetPt;
+ else if (title=="#eta_{jet}^{leading}")
+ contents[i] = maxJetEta;
+ else if (title=="#phi_{jet}^{leading}")
+ contents[i] = maxJetPhi;
+ else
+ AliWarning(Form("Unable to fill dimension %s!",title.Data()));
+ }
+
+ fHistEventQA->Fill(contents);
+}
+
//________________________________________________________________________
Int_t AliAnalysisTaskSAQA::DoCellLoop(Float_t &sum, Float_t &sum_cut)
{
const Int_t ncells = cells->GetNumberOfCells();
for (Int_t pos = 0; pos < ncells; pos++) {
- Float_t amp = cells->GetAmplitude(pos);
- fHistCellsEnergy->Fill(amp);
+ Float_t amp = cells->GetAmplitude(pos);
+ Int_t absId = cells->GetCellNumber(pos);
+ fHistCellsAbsIdEnergy->Fill(absId,amp);
sum += amp;
if (amp < fCellEnergyCut)
continue;
}
//________________________________________________________________________
-Float_t AliAnalysisTaskSAQA::DoClusterLoop()
+Double_t AliAnalysisTaskSAQA::GetCellEnergySum(AliVCluster *cluster, AliVCaloCells *cells)
+{
+ Double_t sum = 0;
+ for (Int_t i = 0; i < cluster->GetNCells(); i++)
+ sum += cells->GetCellAmplitude(cluster->GetCellAbsId(i));
+ return sum;
+}
+
+//________________________________________________________________________
+Double_t AliAnalysisTaskSAQA::GetFcross(AliVCluster *cluster, AliVCaloCells *cells)
+{
+ Int_t AbsIdseed = -1;
+ Double_t Eseed = 0;
+ for (Int_t i = 0; i < cluster->GetNCells(); i++) {
+ if (cells->GetCellAmplitude(cluster->GetCellAbsId(i)) > AbsIdseed) {
+ Eseed = cells->GetCellAmplitude(cluster->GetCellAbsId(i));
+ AbsIdseed = cluster->GetCellAbsId(i);
+ }
+ }
+
+ if (Eseed < 1e-9)
+ return 100;
+
+ Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
+ fGeom->GetCellIndex(AbsIdseed,imod,iTower,iIphi,iIeta);
+ fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,iphi,ieta);
+
+ //Get close cells index and energy, not in corners
+
+ Int_t absID1 = -1;
+ Int_t absID2 = -1;
+
+ if (iphi < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta);
+ if (iphi > 0) absID2 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta);
+
+ // In case of cell in eta = 0 border, depending on SM shift the cross cell index
+
+ Int_t absID3 = -1;
+ Int_t absID4 = -1;
+
+ if (ieta == AliEMCALGeoParams::fgkEMCALCols-1 && !(imod%2)) {
+ absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod+1, iphi, 0);
+ absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
+ }
+ else if (ieta == 0 && imod%2) {
+ absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
+ absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod-1, iphi, AliEMCALGeoParams::fgkEMCALCols-1);
+ }
+ else {
+ if (ieta < AliEMCALGeoParams::fgkEMCALCols-1)
+ absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
+ if (ieta > 0)
+ absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
+ }
+
+ Double_t ecell1 = cells->GetCellAmplitude(absID1);
+ Double_t ecell2 = cells->GetCellAmplitude(absID2);
+ Double_t ecell3 = cells->GetCellAmplitude(absID3);
+ Double_t ecell4 = cells->GetCellAmplitude(absID4);
+
+ Double_t Ecross = ecell1 + ecell2 + ecell3 + ecell4;
+
+ Double_t Fcross = 1 - Ecross/Eseed;
+
+ return Fcross;
+}
+
+//________________________________________________________________________
+Int_t AliAnalysisTaskSAQA::DoClusterLoop(Float_t &sum, AliVCluster* &leading)
{
// Do cluster loop.
if (!fCaloClusters)
return 0;
- Float_t sum = 0;
+ Int_t nAccClusters = 0;
+
+ AliVCaloCells *cells = InputEvent()->GetEMCALCells();
+
+ sum = 0;
+ leading = 0;
// Cluster loop
- Int_t nclusters = fCaloClusters->GetEntriesFast();
+ const Int_t nclusters = fCaloClusters->GetEntriesFast();
for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
continue;
}
- if (!AcceptCluster(cluster, kTRUE))
+ if (!AcceptCluster(cluster))
continue;
sum += cluster->E();
+ if (!leading || leading->E() < cluster->E()) leading = cluster;
+
TLorentzVector nPart;
cluster->GetMomentum(nPart, fVertex);
- fHistClusPhiEtaEnergy->Fill(cluster->E(), nPart.Eta(), nPart.Phi());
+ fHistClusPhiEtaEnergy[fCentBin]->Fill(nPart.Eta(), nPart.Phi(), cluster->E());
+
+ Double_t ep = nPart.Phi() - fEPV0;
+ while (ep < 0) ep += TMath::Pi();
+ while (ep >= TMath::Pi()) ep -= TMath::Pi();
+ fHistClusDeltaPhiEPEnergy[fCentBin]->Fill(cluster->E(), ep);
+
fHistNCellsEnergy->Fill(cluster->E(), cluster->GetNCells());
fHistClusTimeEnergy->Fill(cluster->E(), cluster->GetTOF());
+
+ if (cells) fHistFcrossEnergy->Fill(cluster->E(), GetFcross(cluster, cells));
+
+ if (cells) fHistClusEnergyMinusCellEnergy->Fill(cluster->E(), cluster->E() - GetCellEnergySum(cluster,cells));
+
+ if (fHistClusMCEnergyFraction[fCentBin])
+ fHistClusMCEnergyFraction[fCentBin]->Fill(cluster->GetMCEnergyFraction());
+
+ nAccClusters++;
}
- return sum;
+ return nAccClusters;
}
//________________________________________________________________________
-Float_t AliAnalysisTaskSAQA::DoTrackLoop()
+Int_t AliAnalysisTaskSAQA::DoTrackLoop(Float_t &sum, AliVParticle* &leading)
{
// Do track loop.
if (!fTracks)
return 0;
- Float_t sum = 0;
+ Int_t nAccTracks = 0;
- Int_t ntracks = fTracks->GetEntriesFast();
- Int_t nclusters = 0;
- if (fCaloClusters)
- nclusters = fCaloClusters->GetEntriesFast();
+ sum = 0;
+ leading = 0;
+
+ const Int_t ntracks = fTracks->GetEntriesFast();
+ Int_t neg = 0;
+ Int_t zero = 0;
for (Int_t i = 0; i < ntracks; i++) {
continue;
}
- AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track);
-
- if (vtrack && !AcceptTrack(vtrack, kTRUE))
+ if (!AcceptTrack(track))
continue;
-
- fHistTracksPt->Fill(track->Pt());
+
+ nAccTracks++;
sum += track->P();
- Int_t label = track->GetLabel();
-
- fHistTrPhiEta->Fill(track->Eta(), track->Phi());
-
- fHistTrackEta[4]->Fill(track->Eta());
- fHistTrackPhi[4]->Fill(track->Phi());
+ if (!leading || leading->Pt() < track->Pt()) leading = track;
- if (label >= 0 && label < 4) {
- fHistTrackEta[label]->Fill(track->Eta());
- fHistTrackPhi[label]->Fill(track->Phi());
+ if (fParticleLevel) {
+ fHistTrPhiEtaPt[fCentBin][0]->Fill(track->Eta(), track->Phi(), track->Pt());
}
+ else {
+ fHistTrPhiEtaPt[fCentBin][3]->Fill(track->Eta(), track->Phi(), track->Pt());
+ if (track->GetLabel() == 0) {
+ zero++;
+ if (fHistTrPhiEtaZeroLab[fCentBin]) {
+ fHistTrPhiEtaZeroLab[fCentBin]->Fill(track->Eta(), track->Phi());
+ fHistTrPtZeroLab[fCentBin]->Fill(track->Pt());
+ }
+ }
- if (!vtrack)
- continue;
+ if (track->GetLabel() < 0)
+ neg++;
- if (vtrack->GetTrackEtaOnEMCal() == -999 || vtrack->GetTrackPhiOnEMCal() == -999)
- fHistTrPhiEtaNonProp->Fill(vtrack->Eta(), vtrack->Phi());
+ Int_t type = 0;
- fHistTrEmcPhiEta->Fill(vtrack->GetTrackEtaOnEMCal(), vtrack->GetTrackPhiOnEMCal());
- fHistDeltaEtaPt->Fill(vtrack->Pt(), vtrack->Eta() - vtrack->GetTrackEtaOnEMCal());
- fHistDeltaPhiPt->Fill(vtrack->Pt(), vtrack->Phi() - vtrack->GetTrackPhiOnEMCal());
+ AliPicoTrack* ptrack = dynamic_cast<AliPicoTrack*>(track);
+ if (ptrack)
+ type = ptrack->GetTrackType();
- if (fRepropagateTracks && vtrack->GetTrackEtaOnEMCal() > -2) {
- Float_t propeta = -999, propphi = -999;
- PropagateTrack(vtrack, propeta, propphi);
- fHistDeltaEtaNewProp->Fill(propeta - vtrack->GetTrackEtaOnEMCal());
- fHistDeltaPhiNewProp->Fill(propphi - vtrack->GetTrackPhiOnEMCal());
+ if (type >= 0 && type < 3)
+ fHistTrPhiEtaPt[fCentBin][type]->Fill(track->Eta(), track->Phi(), track->Pt());
+ else
+ AliDebug(2,Form("%s: track type %d not recognized!", GetName(), type));
}
- }
-
- return sum;
-}
-//____________________________________________________________________________
-void AliAnalysisTaskSAQA::PropagateTrack(AliVTrack *track, Float_t &eta, Float_t &phi)
-{
- eta = -999;
- phi = -999;
-
- if (!track)
- return;
-
- if (track->Pt() == 0)
- return;
-
- // init the magnetic field if not already on
- if(!TGeoGlobalMagField::Instance()->GetField()) {
- AliInfo("Init the magnetic field\n");
- AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(InputEvent());
- if (aodevent) {
- Double_t curSol = 30000*aodevent->GetMagneticField()/5.00668;
- Double_t curDip = 6000 *aodevent->GetMuonMagFieldScale();
- AliMagF *field = AliMagF::CreateFieldMap(curSol,curDip);
- TGeoGlobalMagField::Instance()->SetField(field);
+ AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track);
+
+ if (!vtrack)
+ continue;
+
+ if ((vtrack->GetTrackEtaOnEMCal() == -999 || vtrack->GetTrackPhiOnEMCal() == -999) && fHistTrPhiEtaNonProp[fCentBin]) {
+ fHistTrPhiEtaNonProp[fCentBin]->Fill(vtrack->Eta(), vtrack->Phi());
+ fHistTrPtNonProp[fCentBin]->Fill(vtrack->Pt());
}
+
+ if (fHistTrEmcPhiEta[fCentBin])
+ fHistTrEmcPhiEta[fCentBin]->Fill(vtrack->GetTrackEtaOnEMCal(), vtrack->GetTrackPhiOnEMCal());
+ if (fHistTrEmcPt[fCentBin])
+ fHistTrEmcPt[fCentBin]->Fill(vtrack->GetTrackPtOnEMCal());
+ if (fHistDeltaEtaPt[fCentBin])
+ fHistDeltaEtaPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Eta() - vtrack->GetTrackEtaOnEMCal());
+ if (fHistDeltaPhiPt[fCentBin])
+ fHistDeltaPhiPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Phi() - vtrack->GetTrackPhiOnEMCal());
+ if (fHistDeltaPtvsPt[fCentBin])
+ fHistDeltaPtvsPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Pt() - vtrack->GetTrackPtOnEMCal());
}
-
- Double_t cv[21];
- for (Int_t i = 0; i < 21; i++) cv[i] = 0;
-
- Double_t pos[3], mom[3];
- track->GetXYZ(pos);
- track->GetPxPyPz(mom);
- AliExternalTrackParam *trackParam = new AliExternalTrackParam(pos, mom, cv, track->Charge());
-
- if(!AliTrackerBase::PropagateTrackToBxByBz(trackParam, 430., 0.139, 20, kTRUE, 0.8, -1)) return;
- Double_t trkPos[3] = {0., 0., 0.};
- if(!trackParam->GetXYZ(trkPos)) return;
- TVector3 trkPosVec(trkPos[0], trkPos[1], trkPos[2]);
- eta = trkPosVec.Eta();
- phi = trkPosVec.Phi();
- if(phi < 0)
- phi += 2 * TMath::Pi();
+
+ if (fHistTrNegativeLabels)
+ fHistTrNegativeLabels->Fill(1. * neg / ntracks);
+
+ if (fHistTrZeroLabels)
+ fHistTrZeroLabels->Fill(1. * zero / ntracks);
+
+ return nAccTracks;
}
//________________________________________________________________________
-void AliAnalysisTaskSAQA::DoJetLoop()
+Int_t AliAnalysisTaskSAQA::DoJetLoop(AliEmcalJet* &leading)
{
// Do jet loop.
if (!fJets)
- return;
+ return 0;
+
+ Int_t nAccJets = 0;
+
+ leading = 0;
Int_t njets = fJets->GetEntriesFast();
continue;
}
- if (!AcceptJet(jet, kFALSE))
+ if (!AcceptJet(jet))
continue;
- fHistJetsPtNonBias[fCentBin]->Fill(jet->Pt());
- fHistJetsPtAreaNonBias[fCentBin]->Fill(jet->Pt(), jet->Area());
-
- if (jet->MaxTrackPt() > fPtBiasJetTrack)
- fHistJetsPtTrack[fCentBin]->Fill(jet->Pt());
-
- if (fAnaType == kEMCAL && jet->MaxClusterPt() > fPtBiasJetClus)
- fHistJetsPtClus[fCentBin]->Fill(jet->Pt());
-
- if (!AcceptBiasJet(jet))
- continue;
+ if (!leading || leading->Pt() < jet->Pt()) leading = jet;
- fHistJetsPt[fCentBin]->Fill(jet->Pt());
- fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
+ nAccJets++;
fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
+ fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
}
-}
-
-//________________________________________________________________________
-Float_t AliAnalysisTaskSAQA::DoTriggerClusLoop()
-{
- // Do trigger cluster loop.
-
- if (!fTrgClusters)
- return 0;
-
- Int_t ntrgclusters = fTrgClusters->GetEntriesFast();
- Float_t maxTrgClus = 0;
-
- for (Int_t iClusters = 0; iClusters < ntrgclusters; iClusters++) {
- AliVCluster* cluster = static_cast<AliVCluster*>(fTrgClusters->At(iClusters));
- if (!cluster) {
- AliError(Form("Could not receive cluster %d", iClusters));
- continue;
- }
-
- if (!(cluster->IsEMCAL())) continue;
-
- if (cluster->E() > maxTrgClus)
- maxTrgClus = cluster->E();
-
- }
- return maxTrgClus;
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskSAQA::DoTriggerPrimitives(Int_t &maxL1amp, Int_t &maxL1thr)
-{
- // Do trigger primitives loop.
-
- AliVCaloTrigger *triggers = InputEvent()->GetCaloTrigger("EMCAL");
-
- if (!triggers || triggers->GetEntries() == 0)
- return;
-
- triggers->Reset();
- Int_t L1amp = 0;
- Int_t L1thr = 0;
- maxL1amp = -1;
- maxL1thr = -1;
-
- while (triggers->Next()) {
- triggers->GetL1TimeSum(L1amp);
- if (maxL1amp < L1amp)
- maxL1amp = L1amp;
-
- triggers->GetL1Threshold(L1thr);
- if (maxL1thr < L1thr)
- maxL1thr = L1thr;
- }
-}
-//________________________________________________________________________
-void AliAnalysisTaskSAQA::Terminate(Option_t *)
-{
- // Called once at the end of the analysis.
+ return nAccJets;
}