#include <TClonesArray.h>
#include <TH1F.h>
#include <TH2F.h>
+#include <TH3F.h>
#include <TList.h>
#include <TLorentzVector.h>
#include <TRandom3.h>
fHistEmbJetsPtArea[i] = 0;
fHistEmbJetsCorrPtArea[i] = 0;
fHistEmbPartPtvsJetPt[i] = 0;
- fHistLeadPartPtvsArea[i] = 0;
+ fHistEmbPartPtvsJetCorrPt[i] = 0;
fHistDistLeadPart2JetAxis[i] = 0;
fHistEmbBkgArea[i] = 0;
fHistRhoVSEmbBkg[i] = 0;
fHistEmbJetsPtArea[i] = 0;
fHistEmbJetsCorrPtArea[i] = 0;
fHistEmbPartPtvsJetPt[i] = 0;
- fHistLeadPartPtvsArea[i] = 0;
+ fHistEmbPartPtvsJetCorrPt[i] = 0;
fHistDistLeadPart2JetAxis[i] = 0;
fHistEmbBkgArea[i] = 0;
fHistRhoVSEmbBkg[i] = 0;
TString histname;
+ const Int_t nbinsZ = 12;
+ Float_t binsZ[nbinsZ+1] = {0,1,2,3,4,5,6,7,8,9,10,20,1000};
+
+ Float_t *binsPt = GenerateFixedBinArray(fNbins, fMinBinPt, fMaxBinPt);
+ Float_t *binsCorrPt = GenerateFixedBinArray(fNbins*2, -fMaxBinPt, fMaxBinPt);
+ Float_t *binsArea = GenerateFixedBinArray(40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3);
+
for (Int_t i = 0; i < fNcentBins; i++) {
if (!fTracksName.IsNull() || !fCaloName.IsNull()) {
histname = "fHistRCPt_";
if (!fEmbJetsName.IsNull()) {
histname = "fHistEmbJetsPtArea_";
histname += i;
- fHistEmbJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
+ fHistEmbJetsPtArea[i] = new TH3F(histname.Data(), histname.Data(), 40, binsArea, fNbins, binsPt, nbinsZ, binsZ);
fHistEmbJetsPtArea[i]->GetXaxis()->SetTitle("area");
fHistEmbJetsPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,raw} (GeV/#it{c})");
fOutput->Add(fHistEmbJetsPtArea[i]);
histname = "fHistEmbJetsCorrPtArea_";
histname += i;
- fHistEmbJetsCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins * 2, -fMaxBinPt, fMaxBinPt);
+ fHistEmbJetsCorrPtArea[i] = new TH3F(histname.Data(), histname.Data(), 40, binsArea, fNbins * 2, binsCorrPt, nbinsZ, binsZ);
fHistEmbJetsCorrPtArea[i]->GetXaxis()->SetTitle("area");
fHistEmbJetsCorrPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,corr} (GeV/#it{c})");
fOutput->Add(fHistEmbJetsCorrPtArea[i]);
fHistEmbPartPtvsJetPt[i]->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistEmbPartPtvsJetPt[i]);
- histname = "fHistLeadPartPtvsArea_";
+ histname = "fHistEmbPartPtvsJetCorrPt_";
histname += i;
- fHistLeadPartPtvsArea[i] = new TH2F(histname.Data(), histname.Data(), 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
- fHistLeadPartPtvsArea[i]->GetXaxis()->SetTitle("area");
- fHistLeadPartPtvsArea[i]->GetYaxis()->SetTitle("#it{p}_{T,const}^{leading} (GeV/#it{c})");
- fOutput->Add(fHistLeadPartPtvsArea[i]);
+ fHistEmbPartPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
+ fHistEmbPartPtvsJetCorrPt[i]->GetXaxis()->SetTitle("#sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
+ fHistEmbPartPtvsJetCorrPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - A#rho (GeV/#it{c})");
+ fHistEmbPartPtvsJetCorrPt[i]->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistEmbPartPtvsJetCorrPt[i]);
histname = "fHistDistLeadPart2JetAxis_";
histname += i;
if (fJets) {
// Random cones far from leading jet
- Int_t *sortedJets = GetSortedArray(fJets);
+ static Int_t sortedJets[9999] = {-1};
+ GetSortedArray(sortedJets, fJets);
AliEmcalJet* jet = 0;
- if (sortedJets && sortedJets[0] > 0)
+ if (sortedJets[0] >= 0)
jet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[0]));
RCpt = 0;
Int_t countEmbJets = 0;
while (embJet != 0) {
-
+ AliDebug(2,Form("Elaborating embedded jet n. %d", countEmbJets));
countEmbJets++;
Double_t maxClusterPt = 0;
Double_t maxPartEta = 0;
Double_t maxPartPhi = 0;
- AliVCluster *cluster = embJet->GetLeadingCluster(fEmbCaloClusters);
- if (cluster) {
- TLorentzVector nPart;
- cluster->GetMomentum(nPart, fVertex);
-
- maxClusterEta = nPart.Eta();
- maxClusterPhi = nPart.Phi();
- maxClusterPt = nPart.Pt();
- }
-
- AliVParticle *track = embJet->GetLeadingTrack(fEmbTracks);
- if (track) {
- maxTrackEta = track->Eta();
- maxTrackPhi = track->Phi();
- maxTrackPt = track->Pt();
+ if (fLeadingHadronType == 1 || fLeadingHadronType == 2) {
+ AliVCluster *cluster = embJet->GetLeadingCluster(fEmbCaloClusters);
+ if (cluster) {
+ TLorentzVector nPart;
+ cluster->GetMomentum(nPart, fVertex);
+
+ maxClusterEta = nPart.Eta();
+ maxClusterPhi = nPart.Phi();
+ maxClusterPt = nPart.Pt();
+ }
}
- if (!track && !cluster) {
- AliWarning(Form("%s - Embedded jet found but no leading particle was found (?) !", GetName()));
- return kTRUE;
+ if (fLeadingHadronType == 0 || fLeadingHadronType == 2) {
+ AliVParticle *track = embJet->GetLeadingTrack(fEmbTracks);
+ if (track) {
+ maxTrackEta = track->Eta();
+ maxTrackPhi = track->Phi();
+ maxTrackPt = track->Pt();
+ }
}
if (maxTrackPt > maxClusterPt) {
Double_t distLeading2Jet = TMath::Sqrt((embJet->Eta() - maxPartEta) * (embJet->Eta() - maxPartEta) + (embJet->Phi() - maxPartPhi) * (embJet->Phi() - maxPartPhi));
fHistEmbPartPtvsJetPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt());
+ fHistEmbPartPtvsJetCorrPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt() - embJet->Area() * fRhoVal);
fHistLeadPartPhiEta->Fill(maxPartEta, maxPartPhi);
- fHistLeadPartPtvsArea[fCentBin]->Fill(maxPartPt, embJet->Area());
fHistDistLeadPart2JetAxis[fCentBin]->Fill(distLeading2Jet);
- fHistEmbJetsPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt());
- fHistEmbJetsCorrPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - fRhoVal * embJet->Area());
+ fHistEmbJetsPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt(), maxPartPt);
+ fHistEmbJetsCorrPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - fRhoVal * embJet->Area(), maxPartPt);
fHistEmbJetsPhiEta->Fill(embJet->Eta(), embJet->Phi());
fHistEmbBkgArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - embJet->MCPt());
}
if (countEmbJets==0) {
+ AliDebug(1,"No embedded jets found!");
if (fEmbTracks) {
DoEmbTrackLoop();
for (Int_t i = 0; i < fEmbeddedTrackNIds; i++) {
+ AliDebug(2,Form("Embedded track %d found!",i));
AliVParticle *track2 = static_cast<AliVParticle*>(fEmbTracks->At(fEmbeddedTrackIds[i]));
if (!track2) continue;
fHistEmbNotFoundPhiEta[fCentBin]->Fill(track2->Eta(), track2->Phi());
if (fEmbCaloClusters) {
DoEmbClusterLoop();
for (Int_t i = 0; i < fEmbeddedClusterNIds; i++) {
+ AliDebug(2,Form("Embedded cluster %d found!",i));
AliVCluster *cluster2 = static_cast<AliVCluster*>(fEmbCaloClusters->At(fEmbeddedClusterIds[i]));
TLorentzVector nPart;
cluster2->GetMomentum(nPart, fVertex);
AliError(Form("Could not receive jet %d", iJet));
continue;
}
+
+ if (jet->MCPt() < fMCJetPtThreshold)
+ continue;
if (!AcceptJet(jet))
continue;
- if (jet->MCPt() >= fMCJetPtThreshold)
- return jet;
+ return jet;
}
return 0;
class TString;
class TH1;
class TH2;
+class TH3;
#include "AliAnalysisTaskEmcalJet.h"
// Jet embedding
TH2 *fHistEmbNotFoundPhiEta[4]; //!Phi-Eta of "not found" embedded particles
TH1 *fHistEmbNotFoundPt[4]; //!Pt of "not found" embedded particles
- TH2 *fHistEmbJetsPtArea[4]; //!Pt vs. area of embedded jets
- TH2 *fHistEmbJetsCorrPtArea[4]; //!Pt-rho*A vs. area of embedded jets
- TH2 *fHistEmbPartPtvsJetPt[4]; //!MC pt vs total pt of embedded jets
- TH2 *fHistEmbJetsPhiEta; //!Phi-Eta distribution of embedded jets
+ TH3 *fHistEmbJetsPtArea[4]; //!Pt vs. area of embedded jets
+ TH3 *fHistEmbJetsCorrPtArea[4]; //!Pt-rho*A vs. area of embedded jets
+ TH2 *fHistEmbPartPtvsJetPt[4]; //!MC jet pt total jet pt
+ TH2 *fHistEmbPartPtvsJetCorrPt[4];//!MC jet pt total jet pt - rho*A
+ TH2 *fHistEmbJetsPhiEta; //!Phi-Eta distribution of embedded jets<
TH2 *fHistLeadPartPhiEta; //!Phi-Eta distribution of the leading particle of embedded jets
- TH2 *fHistLeadPartPtvsArea[4]; //!Pt of the leading particle vs. area of embedded jets
TH1 *fHistDistLeadPart2JetAxis[4];//!Distance between leading particle and jet axis
TH2 *fHistEmbBkgArea[4]; //!Pt(embjet) - Pt(embtrack) vs. area of embedded jets
TH2 *fHistRhoVSEmbBkg[4]; //!Area(embjet) * rho vs. Pt(embjet) - Pt(embtrack)
fMaxClusterPt(100),
fMaxTrackPt(100),
fLeadingHadronType(0),
+ fNLeadingJets(2),
fJets(0),
fRho(0),
fRhoVal(0)
fMaxClusterPt(100),
fMaxTrackPt(100),
fLeadingHadronType(0),
+ fNLeadingJets(2),
fJets(0),
fRho(0),
fRhoVal(0)
{
// Return true if jet is accepted.
- if (jet->Pt() <= fJetPtCut)
+ if (jet->Pt() <= fJetPtCut)
return kFALSE;
- if (jet->Area() <= fJetAreaCut)
+
+ if (jet->Area() <= fJetAreaCut)
return kFALSE;
+
if (jet->AreaEmc()<fAreaEmcCut)
return kFALSE;
+
if (!AcceptBiasJet(jet))
return kFALSE;
+
if (jet->MaxTrackPt() > fMaxTrackPt || jet->MaxClusterPt() > fMaxClusterPt)
return kFALSE;
}
//________________________________________________________________________
-Int_t* AliAnalysisTaskEmcalJet::GetSortedArray(TClonesArray *array) const
+Bool_t AliAnalysisTaskEmcalJet::GetSortedArray(Int_t indexes[], TClonesArray *array, Double_t rho) const
{
// Get the leading jets.
- static Float_t pt[9999];
- static Int_t indexes[9999];
+ static Float_t pt[9999] = {0};
if (!array)
return 0;
const Int_t n = array->GetEntriesFast();
if (n < 1)
- return 0;
+ return kFALSE;
if (array->GetClass()->GetBaseClass("AliEmcalJet")) {
if (!AcceptJet(jet))
continue;
- pt[i] = jet->Pt() - fRhoVal * jet->Area();
+ pt[i] = jet->Pt() - rho * jet->Area();
}
}
if (pt[indexes[0]] == -FLT_MAX)
return 0;
- return indexes;
+ return kTRUE;
}
//________________________________________________________________________
void SetPtBiasJetClus(Float_t b) { fPtBiasJetClus = b ; }
void SetPtBiasJetTrack(Float_t b) { fPtBiasJetTrack = b ; }
void SetLeadingHadronType(Int_t t) { fLeadingHadronType = t ; }
+ void SetNLeadingJets(Int_t t) { fNLeadingJets = t ; }
protected:
Float_t* GenerateFixedBinArray(Int_t n, Float_t min, Float_t max) const;
Double_t GetLeadingHadronPt(AliEmcalJet* jet) const;
void ExecOnce() ;
AliRhoParameter *GetRhoFromEvent(const char *name) ;
- Int_t *GetSortedArray(TClonesArray *array) const;
+ Bool_t GetSortedArray(Int_t indexes[], TClonesArray *array, Double_t rho=0) const;
Bool_t IsJetTrack(AliEmcalJet* jet, Int_t itrack, Bool_t sorted = kTRUE) const;
Bool_t IsJetCluster(AliEmcalJet* jet, Int_t iclus, Bool_t sorted = kTRUE) const;
Bool_t RetrieveEventObjects() ;
Float_t fMaxClusterPt; // maximum cluster constituent pt to accept the jet
Float_t fMaxTrackPt; // maximum track constituent pt to accept the jet
Int_t fLeadingHadronType; // 0 = charged, 1 = neutral, 2 = both
+ Int_t fNLeadingJets; // how many jets are to be considered the leading jet(s)
TClonesArray *fJets; //!jets
- AliRhoParameter *fRho; //!Event rho
- Double_t fRhoVal; //!Event rho value
+ AliRhoParameter *fRho; //!event rho
+ Double_t fRhoVal; //!event rho value
private:
AliAnalysisTaskEmcalJet(const AliAnalysisTaskEmcalJet&); // not implemented
AliAnalysisTaskEmcalJet &operator=(const AliAnalysisTaskEmcalJet&); // not implemented
- ClassDef(AliAnalysisTaskEmcalJet, 5) // EMCAL Jet base analysis task
+ ClassDef(AliAnalysisTaskEmcalJet, 6) // EMCAL Jet base analysis task
};
#endif
}
if (fJets) {
- Int_t *sortedJets = GetSortedArray(fJets);
+ static Int_t sortedJets[9999] = {-1};
+ GetSortedArray(sortedJets, fJets);
- if (sortedJets) {
+ if (sortedJets[0]>=0) {
AliEmcalJet* leadJet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[0]));
if (leadJet)
fHistLeadingJetPt[fCentBin]->Fill(leadJet->Pt());
// Authors: C.Loizides, S.Aiola
#include <vector>
-#include <algorithm>
#include "AliEmcalJetTask.h"
#include <TChain.h>
#include <TClonesArray.h>
#include <TList.h>
#include <TLorentzVector.h>
+#include <TMath.h>
#include "AliAnalysisManager.h"
#include "AliCentrality.h"
ClassImp(AliEmcalJetTask)
-//________________________________________________________________________
-inline bool AliEmcalJetTask::ComparePseudoJets(fastjet::PseudoJet a, fastjet::PseudoJet b)
-{
- if (a.perp() > b.perp())
- return true;
- else
- return false;
-}
-
//________________________________________________________________________
AliEmcalJetTask::AliEmcalJetTask() :
AliAnalysisTaskSE("AliEmcalJetTask"),
if ((fJetType & kAKT) != 0) {
name = "antikt";
jalgo = fastjet::antikt_algorithm;
+ AliDebug(1,"Using AKT algorithm");
+ }
+ else {
+ AliDebug(1,"Using KT algorithm");
}
if ((fJetType & kR020Jet) != 0)
Double_t vertex[3] = {0, 0, 0};
fEvent->GetPrimaryVertex()->GetXYZ(vertex);
- AliDebug(2,Form("Jet type = %d)", fJetType));
+ AliDebug(2,Form("Jet type = %d", fJetType));
if ((fIsMcPart || ((fJetType & kFullJet) != 0) || ((fJetType & kChargedJet) != 0)) && fTracks) {
const Int_t Ntracks = fTracks->GetEntries();
continue;
}
if (fIsMcPart || t->GetLabel() != 0) {
- if (fMCConstSel != kAllJets && t->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
- AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fMCConstSel, t->TestBits(TObject::kBitMask)));
+ if (fMCConstSel == kNone) {
+ AliDebug(2,Form("Skipping track %d because bit mask is 0.", iTracks));
continue;
}
+ if (fMCConstSel != kAllJets) {
+ if (t->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
+ AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fMCConstSel, t->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Track %d matches the bit mask (%d, %d)", iTracks, fMCConstSel, t->TestBits(TObject::kBitMask)));
+ }
+ }
}
else {
- if (fConstSel != kAllJets && t->TestBits(fConstSel) != (Int_t)fConstSel) {
- AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fConstSel, t->TestBits(TObject::kBitMask)));
+ if (fConstSel == kNone) {
+ AliDebug(2,Form("Skipping track %d because bit mask is 0.", iTracks));
continue;
}
- }
+ if (fConstSel != kAllJets) {
+ if (t->TestBits(fConstSel) != (Int_t)fConstSel) {
+ AliDebug(2,Form("Skipping track %d because it does not match the bit mask (%d, %d)", iTracks, fConstSel, t->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Track %d matches the bit mask (%d, %d)", iTracks, fConstSel, t->TestBits(TObject::kBitMask)));
+ }
+ }
+ }
if (t->Pt() < fMinJetTrackPt)
continue;
Double_t eta = t->Eta();
continue;
// offset of 100 for consistency with cluster ids
+ AliDebug(2,Form("Track %d accepted (label = %d, pt = %f)", iTracks, t->GetLabel(), t->Pt()));
fjw.AddInputVector(t->Px(), t->Py(), t->Pz(), t->P(), iTracks + 100);
}
}
continue;
if (c->GetLabel() > 0) {
- if (fMCConstSel != kAllJets && ep->TestBits(fMCConstSel) != (Int_t)fMCConstSel)
+ if (fMCConstSel == kNone) {
+ AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
continue;
+ }
+ if (fMCConstSel != kAllJets) {
+ if (ep->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
+ AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fMCConstSel, ep->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fMCConstSel, ep->TestBits(TObject::kBitMask)));
+ }
+ }
}
else {
- if (fConstSel != kAllJets && ep->TestBits(fConstSel) != (Int_t)fConstSel)
+ if (fConstSel == kNone) {
+ AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
continue;
+ }
+ if (fConstSel != kAllJets) {
+ if (ep->TestBits(fConstSel) != (Int_t)fConstSel) {
+ AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fConstSel, ep->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fConstSel, ep->TestBits(TObject::kBitMask)));
+ }
+ }
}
cEta = ep->Eta();
continue;
if (c->GetLabel() > 0) {
- if (fMCConstSel != kAllJets && c->TestBits(fMCConstSel) != (Int_t)fMCConstSel)
+ if (fMCConstSel == kNone) {
+ AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
continue;
+ }
+ if (fMCConstSel != kAllJets) {
+ if (c->TestBits(fMCConstSel) != (Int_t)fMCConstSel) {
+ AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fMCConstSel, c->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fMCConstSel, c->TestBits(TObject::kBitMask)));
+ }
+ }
}
else {
- if (fConstSel != kAllJets && c->TestBits(fConstSel) != (Int_t)fConstSel)
+ if (fConstSel == kNone) {
+ AliDebug(2,Form("Skipping cluster %d because bit mask is 0.", iClus));
continue;
+ }
+ if (fConstSel != kAllJets) {
+ if (c->TestBits(fConstSel) != (Int_t)fConstSel) {
+ AliDebug(2,Form("Skipping cluster %d because it does not match the bit mask (%d, %d)", iClus, fConstSel, c->TestBits(TObject::kBitMask)));
+ continue;
+ }
+ else {
+ AliDebug(2,Form("Cluster %d matches the bit mask (%d, %d)", iClus, fConstSel, c->TestBits(TObject::kBitMask)));
+ }
+ }
}
TLorentzVector nP;
(cPhi<fPhiMin) || (cPhi>fPhiMax))
continue;
// offset of 100 to skip ghost particles uid = -1
+ AliDebug(2,Form("Cluster %d accepted (label = %d)", iClus, c->GetLabel()));
fjw.AddInputVector(cPx, cPy, cPz, TMath::Sqrt(cPx*cPx+cPy*cPy+cPz*cPz), -iClus - 100);
}
}
// loop over fastjet jets
std::vector<fastjet::PseudoJet> jets_incl = fjw.GetInclusiveJets();
- if (fMarkConst > 0)
- std::sort(jets_incl.begin(), jets_incl.end(), ComparePseudoJets);
- for (UInt_t ij=0, jetCount=0; ij<jets_incl.size(); ++ij) {
+ // sort jets according to jet pt
+ static Int_t indexes[9999] = {-1};
+ GetSortedArray(indexes, jets_incl);
+
+ AliDebug(1,Form("%d jets found", (Int_t)jets_incl.size()));
+ for (UInt_t ijet=0, jetCount=0; ijet<jets_incl.size(); ++ijet) {
+ Int_t ij = indexes[ijet];
+ AliDebug(3,Form("Jet pt = %f, area = %f", jets_incl[ij].perp(), fjw.GetJetArea(ij)));
+
if (jets_incl[ij].perp()<fMinJetPt)
continue;
if (fjw.GetJetArea(ij)<fMinJetArea)
for(UInt_t ic = 0; ic < constituents.size(); ++ic) {
Int_t uid = constituents[ic].user_index();
-
+ AliDebug(2,Form("Processing constituent %d", uid));
if ((uid == -1) /*&& (constituents[ic].kt2() < 1e-25)*/) { //ghost particle
++gall;
Double_t gphi = constituents[ic].phi();
} else if ((uid > 0) && fTracks) { // track constituent
Int_t tid = uid - 100;
AliVParticle *t = static_cast<AliVParticle*>(fTracks->At(tid));
- if (!t)
+ if (!t) {
+ AliError(Form("Could not find track %d",tid));
continue;
- if (jetCount < fMarkConst)
+ }
+ if (jetCount < fMarkConst) {
+ AliDebug(2,Form("Marking track %d with bit map %d", tid, fJetType));
t->SetBit(fJetType);
+ }
Double_t cEta = t->Eta();
Double_t cPhi = t->Phi();
Double_t cPt = t->Pt();
(jet->Eta() > geom->GetArm1EtaMin()) &&
(jet->Eta() < geom->GetArm1EtaMax()))
jet->SetAxisInEmcal(kTRUE);
+
+ AliDebug(2,Form("Added jet n. %d, pt = %f, area = %f, constituents = %d", jetCount, jet->Pt(), jet->Area(), (Int_t)constituents.size()));
jetCount++;
}
//fJets->Sort();
}
+//________________________________________________________________________
+Bool_t AliEmcalJetTask::GetSortedArray(Int_t indexes[], std::vector<fastjet::PseudoJet> array) const
+{
+ // Get the leading jets.
+
+ static Float_t pt[9999] = {0};
+
+ const Int_t n = (Int_t)array.size();
+
+ if (n < 1)
+ return kFALSE;
+
+ for (Int_t i = 0; i < n; i++)
+ pt[i] = array[i].perp();
+
+ TMath::Sort(n, pt, indexes);
+
+ return kTRUE;
+}
+
//________________________________________________________________________
Bool_t AliEmcalJetTask::DoInit()
{
UInt_t GetJetType() { return fJetType; }
- static bool ComparePseudoJets(fastjet::PseudoJet a, fastjet::PseudoJet b);
-
protected:
void FindJets();
Bool_t DoInit();
+ Bool_t GetSortedArray(Int_t indexes[], std::vector<fastjet::PseudoJet> array) const;
TString fTracksName; // name of track collection
TString fCaloName; // name of calo cluster collection
if (!part->IsPhysicalPrimary())
continue;
+ AliDebug(3, Form("Embedding MC particle with pT = %f, eta = %f, phi = %f", part->Pt(), part->Eta(), part->Phi()));
AddMCParticle(part, i);
}
}
label = 99999;
}
- AliDebug(3, Form("Embedding track with pT = %f, eta = %f, phi = %f", track->Pt(), track->Eta(), track->Phi()));
+ AliDebug(3, Form("Embedding track with pT = %f, eta = %f, phi = %f, label = %d", track->Pt(), track->Eta(), track->Phi(), label));
AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc, label);
}
}
#include <TMath.h>
#include <TString.h>
#include <TRandom.h>
+#include <TParameter.h>
+#include <TSystem.h>
+#include <TH1F.h>
#include "AliVEvent.h"
#include "AliLog.h"
//________________________________________________________________________
AliJetEmbeddingFromPYTHIATask::AliJetEmbeddingFromPYTHIATask() :
- AliJetEmbeddingFromAODTask("AliJetEmbeddingFromPYTHIATask", kFALSE),
+ AliJetEmbeddingFromAODTask("AliJetEmbeddingFromPYTHIATask"),
fPYTHIAPath(),
fPtHardBinScaling(),
fLHC11hAnchorRun(kTRUE),
fAnchorRun(-1),
- fCurrentPtHardBin(-1)
+ fCurrentPtHardBin(-1),
+ fPtHardBinParam(0),
+ fHistPtHardBins(0)
{
// Default constructor.
SetSuffix("PYTHIAEmbedding");
}
//________________________________________________________________________
-AliJetEmbeddingFromPYTHIATask::AliJetEmbeddingFromPYTHIATask(const char *name) :
- AliJetEmbeddingFromAODTask(name, kFALSE),
+AliJetEmbeddingFromPYTHIATask::AliJetEmbeddingFromPYTHIATask(const char *name, Bool_t drawqa) :
+ AliJetEmbeddingFromAODTask(name, drawqa),
fPYTHIAPath("/alice/sim/2012/LHC12a15e"),
fPtHardBinScaling(),
fLHC11hAnchorRun(kTRUE),
fAnchorRun(-1),
- fCurrentPtHardBin(-1)
+ fCurrentPtHardBin(-1),
+ fPtHardBinParam(0),
+ fHistPtHardBins(0)
{
// Standard constructor.
SetSuffix("PYTHIAEmbedding");
// Destructor
}
+//________________________________________________________________________
+void AliJetEmbeddingFromPYTHIATask::UserCreateOutputObjects()
+{
+ if (!fQAhistos)
+ return;
+
+ AliJetModelBaseTask::UserCreateOutputObjects();
+
+ fHistPtHardBins = new TH1F("fHistPtHardBins", "fHistPtHardBins", 11, 0, 11);
+ fHistPtHardBins->GetXaxis()->SetTitle("p_{T} hard bin");
+ fHistPtHardBins->GetYaxis()->SetTitle("total events");
+ fOutput->Add(fHistPtHardBins);
+
+ const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
+ const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
+
+ for (Int_t i = 1; i < 12; i++)
+ fHistPtHardBins->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
+
+ PostData(1, fOutput);
+}
+
//________________________________________________________________________
Bool_t AliJetEmbeddingFromPYTHIATask::ExecOnce()
{
fPtHardBinScaling[i] /= sum;
}
+ fPtHardBinParam = static_cast<TParameter<int>*>(InputEvent()->FindListObject("PYTHIAPtHardBin"));
+ if (!fPtHardBinParam) {
+ fPtHardBinParam = new TParameter<int>("PYTHIAPtHardBin", 0);
+ AliDebug(2,"Adding pt hard bin param object to the event list...");
+ InputEvent()->AddObject(fPtHardBinParam);
+ }
+
return AliJetEmbeddingFromAODTask::ExecOnce();
}
{
Int_t newPtHard = GetRandomPtHardBin();
+ fPtHardBinParam->SetVal(newPtHard);
+
+ if (fHistPtHardBins)
+ fHistPtHardBins->Fill(newPtHard+1);
+
if (newPtHard != fCurrentPtHardBin) {
fCurrentPtHardBin = newPtHard;
if (!OpenNextFile()) return kFALSE;
else if (fCurrentAODFileID < 100)
fileName += "0";
fileName += fCurrentAODFileID;
- fileName += "/AliAOD.root";
+
+ if (!gSystem->AccessPathName(Form("%s/AliAOD.root")))
+ fileName += "/AliAOD.root";
+ else if (!gSystem->AccessPathName(Form("%s/aod_archive.zip")))
+ fileName += "/aod_archive.zip#AliAOD.root";
+ else
+ fileName += "/root_archive.zip#AliAOD.root";
return fileName;
}
#include "AliJetEmbeddingFromAODTask.h"
#include <TArrayD.h>
+template<class T>
+class TParameter;
+
class TString;
class AliJetEmbeddingFromPYTHIATask : public AliJetEmbeddingFromAODTask {
public:
AliJetEmbeddingFromPYTHIATask();
- AliJetEmbeddingFromPYTHIATask(const char *name);
+ AliJetEmbeddingFromPYTHIATask(const char *name, Bool_t drawqa=kFALSE);
virtual ~AliJetEmbeddingFromPYTHIATask();
Bool_t UserNotify();
+ void UserCreateOutputObjects();
void SetPYTHIAPath(const char* p) { fPYTHIAPath = p ; }
void SetPtHardBinScaling(Int_t n, Double_t *scaling) { new (&fPtHardBinScaling) TArrayD(n, scaling) ; }
Bool_t fLHC11hAnchorRun ;// LHC11h anchor runs
Int_t fAnchorRun ;// Anchor run
Int_t fCurrentPtHardBin ;//!Pt hard bin of the current open file
+ TParameter<int> *fPtHardBinParam ;//!Pt hard bin param
+ TH1 *fHistPtHardBins ;//!Embeded pt hard bin distribution
private:
AliJetEmbeddingFromPYTHIATask(const AliJetEmbeddingFromPYTHIATask&); // not implemented
#include <TH2F.h>
#include <TH3F.h>
#include <TLorentzVector.h>
+#include <TSystem.h>
+#include <TFile.h>
+#include <TKey.h>
+#include <TProfile.h>
+#include "AliAnalysisManager.h"
#include "AliVCluster.h"
#include "AliVTrack.h"
#include "AliEmcalJet.h"
fRho2(0),
fRho2Val(0),
fTracks2Map(0),
- fHistNTrials(0),
+ fHistTrialsAfterSel(0),
+ fHistEventsAfterSel(0),
+ fHistTrials(0),
+ fHistXsection(0),
fHistEvents(0),
fHistJets1PhiEta(0),
fHistJets1PtArea(0),
fHistJets1CorrPtArea(0),
+ fHistLeadingJets1PtArea(0),
+ fHistLeadingJets1CorrPtArea(0),
fHistJets2PhiEta(0),
fHistJets2PtArea(0),
fHistJets2CorrPtArea(0),
+ fHistLeadingJets2PtArea(0),
+ fHistLeadingJets2CorrPtArea(0),
fHistJets2PhiEtaAcceptance(0),
fHistJets2PtAreaAcceptance(0),
fHistJets2CorrPtAreaAcceptance(0),
- fHistMatchingLevelvsJet2Pt(0),
+ fHistLeadingJets2PtAreaAcceptance(0),
+ fHistLeadingJets2CorrPtAreaAcceptance(0),
+ fHistCommonEnergy1vsJet1Pt(0),
+ fHistCommonEnergy2vsJet2Pt(0),
+ fHistDistancevsJet1Pt(0),
+ fHistDistancevsJet2Pt(0),
fHistDistancevsCommonEnergy1(0),
fHistDistancevsCommonEnergy2(0),
fHistJet2PtOverJet1PtvsJet2Pt(0),
fHistJet1PtOverJet2PtvsJet1Pt(0),
- fHistDeltaEtaPhivsJet2Pt(0),
+ fHistDeltaEtaPhi(0),
fHistDeltaPtvsJet1Pt(0),
fHistDeltaPtvsJet2Pt(0),
fHistDeltaPtvsMatchingLevel(0),
fRho2(0),
fRho2Val(0),
fTracks2Map(0),
- fHistNTrials(0),
+ fHistTrialsAfterSel(0),
+ fHistEventsAfterSel(0),
+ fHistTrials(0),
+ fHistXsection(0),
fHistEvents(0),
fHistJets1PhiEta(0),
fHistJets1PtArea(0),
fHistJets1CorrPtArea(0),
+ fHistLeadingJets1PtArea(0),
+ fHistLeadingJets1CorrPtArea(0),
fHistJets2PhiEta(0),
fHistJets2PtArea(0),
fHistJets2CorrPtArea(0),
+ fHistLeadingJets2PtArea(0),
+ fHistLeadingJets2CorrPtArea(0),
fHistJets2PhiEtaAcceptance(0),
fHistJets2PtAreaAcceptance(0),
fHistJets2CorrPtAreaAcceptance(0),
- fHistMatchingLevelvsJet2Pt(0),
+ fHistLeadingJets2PtAreaAcceptance(0),
+ fHistLeadingJets2CorrPtAreaAcceptance(0),
+ fHistCommonEnergy1vsJet1Pt(0),
+ fHistCommonEnergy2vsJet2Pt(0),
+ fHistDistancevsJet1Pt(0),
+ fHistDistancevsJet2Pt(0),
fHistDistancevsCommonEnergy1(0),
fHistDistancevsCommonEnergy2(0),
fHistJet2PtOverJet1PtvsJet2Pt(0),
fHistJet1PtOverJet2PtvsJet1Pt(0),
- fHistDeltaEtaPhivsJet2Pt(0),
+ fHistDeltaEtaPhi(0),
fHistDeltaPtvsJet1Pt(0),
fHistDeltaPtvsJet2Pt(0),
fHistDeltaPtvsMatchingLevel(0),
// Destructor
}
+//________________________________________________________________________
+Bool_t AliJetResponseMaker::PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
+{
+ //
+ // Get the cross section and the trails either from pyxsec.root or from pysec_hists.root
+ // Get the pt hard bin from the file path
+ // This is to called in Notify and should provide the path to the AOD/ESD file
+ // (Partially copied from AliAnalysisHelperJetTasks)
+
+ TString file(currFile);
+ fXsec = 0;
+ fTrials = 1;
+
+ if(file.Contains("root_archive.zip#")){
+ Ssiz_t pos1 = file.Index("root_archive",12,0,TString::kExact);
+ Ssiz_t pos = file.Index("#",1,pos1,TString::kExact);
+ Ssiz_t pos2 = file.Index(".root",5,TString::kExact);
+ file.Replace(pos+1,pos2-pos1,"");
+ }
+ else {
+ // not an archive take the basename....
+ file.ReplaceAll(gSystem->BaseName(file.Data()),"");
+ }
+ Printf("%s",file.Data());
+
+ // Get the pt hard bin
+ TString strPthard(file);
+ strPthard.Remove(strPthard.Last('/'));
+ strPthard.Remove(strPthard.Last('/'));
+ strPthard.Remove(0,strPthard.Last('/')+1);
+ if (strPthard.IsDec())
+ pthard = strPthard.Atoi();
+ else
+ AliWarning(Form("Could not extract file number from path %s", strPthard.Data()));
+
+ TFile *fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec.root")); // problem that we cannot really test the existance of a file in a archive so we have to lvie with open error message from root
+ if(!fxsec){
+ // next trial fetch the histgram file
+ fxsec = TFile::Open(Form("%s%s",file.Data(),"pyxsec_hists.root"));
+ if(!fxsec){
+ // not a severe condition but inciate that we have no information
+ return kFALSE;
+ }
+ else{
+ // find the tlist we want to be independtent of the name so use the Tkey
+ TKey* key = (TKey*)fxsec->GetListOfKeys()->At(0);
+ if(!key){
+ fxsec->Close();
+ return kFALSE;
+ }
+ TList *list = dynamic_cast<TList*>(key->ReadObj());
+ if(!list){
+ fxsec->Close();
+ return kFALSE;
+ }
+ fXsec = ((TProfile*)list->FindObject("h1Xsec"))->GetBinContent(1);
+ fTrials = ((TH1F*)list->FindObject("h1Trials"))->GetBinContent(1);
+ fxsec->Close();
+ }
+ } // no tree pyxsec.root
+ else {
+ TTree *xtree = (TTree*)fxsec->Get("Xsection");
+ if(!xtree){
+ fxsec->Close();
+ return kFALSE;
+ }
+ UInt_t ntrials = 0;
+ Double_t xsection = 0;
+ xtree->SetBranchAddress("xsection",&xsection);
+ xtree->SetBranchAddress("ntrials",&ntrials);
+ xtree->GetEntry(0);
+ fTrials = ntrials;
+ fXsec = xsection;
+ fxsec->Close();
+ }
+ return kTRUE;
+}
+
+//________________________________________________________________________
+Bool_t AliJetResponseMaker::UserNotify()
+{
+ TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
+ if (!tree) {
+ AliError(Form("%s - UserNotify: No current tree!",GetName()));
+ return kFALSE;
+ }
+
+ Float_t xsection = 0;
+ Float_t trials = 0;
+ Int_t pthard = 0;
+
+ TFile *curfile = tree->GetCurrentFile();
+ if (!curfile) {
+ AliError(Form("%s - UserNotify: No current file!",GetName()));
+ return kFALSE;
+ }
+
+ PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthard);
+
+ fHistTrials->SetBinContent(pthard + 1, fHistTrials->GetBinContent(pthard + 1) + trials);
+ fHistXsection->SetBinContent(pthard + 1, fHistXsection->GetBinContent(pthard + 1) + xsection);
+ fHistEvents->SetBinContent(pthard + 1, fHistEvents->GetBinContent(pthard + 1) + 1);
+
+ return kTRUE;
+}
+
//________________________________________________________________________
void AliJetResponseMaker::UserCreateOutputObjects()
{
AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
- const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
- const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
+ fHistTrialsAfterSel = new TH1F("fHistTrialsAfterSel", "fHistTrialsAfterSel", 11, 0, 11);
+ fHistTrialsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
+ fHistTrialsAfterSel->GetYaxis()->SetTitle("trials");
+ fOutput->Add(fHistTrialsAfterSel);
+
+ fHistEventsAfterSel = new TH1F("fHistEventsAfterSel", "fHistEventsAfterSel", 11, 0, 11);
+ fHistEventsAfterSel->GetXaxis()->SetTitle("p_{T} hard bin");
+ fHistEventsAfterSel->GetYaxis()->SetTitle("total events");
+ fOutput->Add(fHistEventsAfterSel);
- fHistNTrials = new TH1F("fHistNTrials", "fHistNTrials", 11, 0, 11);
- fHistNTrials->GetXaxis()->SetTitle("p_{T} hard bin");
- fHistNTrials->GetYaxis()->SetTitle("trials");
- fOutput->Add(fHistNTrials);
+ fHistTrials = new TH1F("fHistTrials", "fHistTrials", 11, 0, 11);
+ fHistTrials->GetXaxis()->SetTitle("p_{T} hard bin");
+ fHistTrials->GetYaxis()->SetTitle("trials");
+ fOutput->Add(fHistTrials);
+
+ fHistXsection = new TH1F("fHistXsection", "fHistXsection", 11, 0, 11);
+ fHistXsection->GetXaxis()->SetTitle("p_{T} hard bin");
+ fHistXsection->GetYaxis()->SetTitle("xsection");
+ fOutput->Add(fHistXsection);
fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
fHistEvents->GetYaxis()->SetTitle("total events");
fOutput->Add(fHistEvents);
+ const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
+ const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
+
for (Int_t i = 1; i < 12; i++) {
- fHistNTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
+ fHistTrialsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
+ fHistEventsAfterSel->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
+
+ fHistTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
+ fHistXsection->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
}
fHistJets1PtArea->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets1PtArea);
+ fHistLeadingJets1PtArea = new TH2F("fHistLeadingJets1PtArea", "fHistLeadingJets1PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
+ fHistLeadingJets1PtArea->GetXaxis()->SetTitle("area");
+ fHistLeadingJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
+ fHistLeadingJets1PtArea->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets1PtArea);
+
if (!fRhoName.IsNull()) {
fHistJets1CorrPtArea = new TH2F("fHistJets1CorrPtArea", "fHistJets1CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
fHistJets1CorrPtArea->GetXaxis()->SetTitle("area");
fHistJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
fHistJets1CorrPtArea->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets1CorrPtArea);
+
+ fHistLeadingJets1CorrPtArea = new TH2F("fHistLeadingJets1CorrPtArea", "fHistLeadingJets1CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
+ fHistLeadingJets1CorrPtArea->GetXaxis()->SetTitle("area");
+ fHistLeadingJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
+ fHistLeadingJets1CorrPtArea->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets1CorrPtArea);
}
fHistJets2PhiEta = new TH2F("fHistJets2PhiEta", "fHistJets2PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
fHistJets2PtArea->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets2PtArea);
+ fHistLeadingJets2PtArea = new TH2F("fHistLeadingJets2PtArea", "fHistLeadingJets2PtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
+ fHistLeadingJets2PtArea->GetXaxis()->SetTitle("area");
+ fHistLeadingJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
+ fHistLeadingJets2PtArea->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets2PtArea);
+
fHistJets2PtAreaAcceptance = new TH2F("fHistJets2PtAreaAcceptance", "fHistJets2PtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
fHistJets2PtAreaAcceptance->GetXaxis()->SetTitle("area");
fHistJets2PtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
fHistJets2PtAreaAcceptance->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets2PtAreaAcceptance);
+ fHistLeadingJets2PtAreaAcceptance = new TH2F("fHistLeadingJets2PtAreaAcceptance", "fHistLeadingJets2PtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, fNbins, fMinBinPt, fMaxBinPt);
+ fHistLeadingJets2PtAreaAcceptance->GetXaxis()->SetTitle("area");
+ fHistLeadingJets2PtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
+ fHistLeadingJets2PtAreaAcceptance->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets2PtAreaAcceptance);
+
if (!fRho2Name.IsNull()) {
fHistJets2CorrPtArea = new TH2F("fHistJets2CorrPtArea", "fHistJets2CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
fHistJets2CorrPtArea->GetXaxis()->SetTitle("area");
fHistJets2CorrPtArea->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets2CorrPtArea);
+ fHistLeadingJets2CorrPtArea = new TH2F("fHistLeadingJets2CorrPtArea", "fHistLeadingJets2CorrPtArea", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
+ fHistLeadingJets2CorrPtArea->GetXaxis()->SetTitle("area");
+ fHistLeadingJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
+ fHistLeadingJets2CorrPtArea->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets2CorrPtArea);
+
fHistJets2CorrPtAreaAcceptance = new TH2F("fHistJets2CorrPtAreaAcceptance", "fHistJets2CorrPtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
fHistJets2CorrPtAreaAcceptance->GetXaxis()->SetTitle("area");
fHistJets2CorrPtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
fHistJets2CorrPtAreaAcceptance->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJets2CorrPtAreaAcceptance);
+
+ fHistLeadingJets2CorrPtAreaAcceptance = new TH2F("fHistLeadingJets2CorrPtAreaAcceptance", "fHistLeadingJets2CorrPtAreaAcceptance", 40, 0, fJetRadius * fJetRadius * TMath::Pi() * 3, 2*fNbins, -fMaxBinPt, fMaxBinPt);
+ fHistLeadingJets2CorrPtAreaAcceptance->GetXaxis()->SetTitle("area");
+ fHistLeadingJets2CorrPtAreaAcceptance->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
+ fHistLeadingJets2CorrPtAreaAcceptance->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistLeadingJets2CorrPtAreaAcceptance);
}
- fHistMatchingLevelvsJet2Pt = new TH2F("fHistMatchingLevelvsJet2Pt", "fHistMatchingLevelvsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
- fHistMatchingLevelvsJet2Pt->GetXaxis()->SetTitle("Matching level");
- fHistMatchingLevelvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
- fHistMatchingLevelvsJet2Pt->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistMatchingLevelvsJet2Pt);
+ fHistCommonEnergy1vsJet1Pt = new TH2F("fHistCommonEnergy1vsJet1Pt", "fHistCommonEnergy1vsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
+ fHistCommonEnergy1vsJet1Pt->GetXaxis()->SetTitle("Common energy 1 (%)");
+ fHistCommonEnergy1vsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
+ fHistCommonEnergy1vsJet1Pt->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistCommonEnergy1vsJet1Pt);
+
+ fHistCommonEnergy2vsJet2Pt = new TH2F("fHistCommonEnergy2vsJet2Pt", "fHistCommonEnergy2vsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
+ fHistCommonEnergy2vsJet2Pt->GetXaxis()->SetTitle("Common energy 2 (%)");
+ fHistCommonEnergy2vsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
+ fHistCommonEnergy2vsJet2Pt->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistCommonEnergy2vsJet2Pt);
+
+ fHistDistancevsJet1Pt = new TH2F("fHistDistancevsJet1Pt", "fHistDistancevsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
+ fHistDistancevsJet1Pt->GetXaxis()->SetTitle("Distance");
+ fHistDistancevsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
+ fHistDistancevsJet1Pt->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistDistancevsJet1Pt);
+
+ fHistDistancevsJet2Pt = new TH2F("fHistDistancevsJet2Pt", "fHistDistancevsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
+ fHistDistancevsJet2Pt->GetXaxis()->SetTitle("Distance");
+ fHistDistancevsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
+ fHistDistancevsJet2Pt->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistDistancevsJet2Pt);
fHistDistancevsCommonEnergy1 = new TH2F("fHistDistancevsCommonEnergy1", "fHistDistancevsCommonEnergy1", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
fHistDistancevsCommonEnergy1->GetXaxis()->SetTitle("Distance");
fHistJet1PtOverJet2PtvsJet1Pt->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJet1PtOverJet2PtvsJet1Pt);
- fHistDeltaEtaPhivsJet2Pt = new TH3F("fHistDeltaEtaPhivsJet2Pt", "fHistDeltaEtaPhivsJet2Pt", 40, -1, 1, 128, -1.6, 4.8, fNbins/2, fMinBinPt, fMaxBinPt);
- fHistDeltaEtaPhivsJet2Pt->GetXaxis()->SetTitle("#Delta#eta");
- fHistDeltaEtaPhivsJet2Pt->GetYaxis()->SetTitle("#Delta#phi");
- fHistDeltaEtaPhivsJet2Pt->GetZaxis()->SetTitle("p_{T,2}");
- fOutput->Add(fHistDeltaEtaPhivsJet2Pt);
+ fHistDeltaEtaPhi = new TH2F("fHistDeltaEtaPhi", "fHistDeltaEtaPhi", 200, -1, 1, 250, -1.6, 4.8);
+ fHistDeltaEtaPhi->GetXaxis()->SetTitle("#Delta#eta");
+ fHistDeltaEtaPhi->GetYaxis()->SetTitle("#Delta#phi");
+ fHistDeltaEtaPhi->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistDeltaEtaPhi);
fHistDeltaPtvsJet1Pt = new TH2F("fHistDeltaPtvsJet1Pt", "fHistDeltaPtvsJet1Pt", fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
fHistDeltaPtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
{
// Fill histograms.
- fHistEvents->SetBinContent(fPtHardBin + 1, fHistEvents->GetBinContent(fPtHardBin + 1) + 1);
- fHistNTrials->SetBinContent(fPtHardBin + 1, fHistNTrials->GetBinContent(fPtHardBin + 1) + fNTrials);
+ static Int_t indexes[9999] = {-1};
+
+ fHistEventsAfterSel->SetBinContent(fPtHardBin + 1, fHistEventsAfterSel->GetBinContent(fPtHardBin + 1) + 1);
+ fHistTrialsAfterSel->SetBinContent(fPtHardBin + 1, fHistTrialsAfterSel->GetBinContent(fPtHardBin + 1) + fNTrials);
+
+ GetSortedArray(indexes, fJets2, fRho2Val);
const Int_t nJets2 = fJets2->GetEntriesFast();
+ Int_t naccJets2 = 0;
+ Int_t naccJets2Acceptance = 0;
+
for (Int_t i = 0; i < nJets2; i++) {
- AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(fJets2->At(i));
+ AliEmcalJet* jet2 = static_cast<AliEmcalJet*>(fJets2->At(indexes[i]));
if (!jet2) {
AliError(Form("Could not receive jet2 %d", i));
fHistJets2PtAreaAcceptance->Fill(jet2->Area(), jet2->Pt());
fHistJets2PhiEtaAcceptance->Fill(jet2->Eta(), jet2->Phi());
- if (!fRho2Name.IsNull())
+ if (naccJets2Acceptance < fNLeadingJets)
+ fHistLeadingJets2PtAreaAcceptance->Fill(jet2->Area(), jet2->Pt());
+
+ if (!fRho2Name.IsNull()) {
fHistJets2CorrPtAreaAcceptance->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
+ if (naccJets2Acceptance < fNLeadingJets)
+ fHistLeadingJets2CorrPtAreaAcceptance->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
+ }
+
+ naccJets2Acceptance++;
}
if (!AcceptBiasJet2(jet2))
if (jet2->Eta() < fJet2MinEta || jet2->Eta() > fJet2MaxEta || jet2->Phi() < fJet2MinPhi || jet2->Phi() > fJet2MaxPhi)
continue;
-
+
fHistJets2PtArea->Fill(jet2->Area(), jet2->Pt());
fHistJets2PhiEta->Fill(jet2->Eta(), jet2->Phi());
+
+ if (naccJets2 < fNLeadingJets)
+ fHistLeadingJets2PtArea->Fill(jet2->Area(), jet2->Pt());
- if (!fRho2Name.IsNull())
+ if (!fRho2Name.IsNull()) {
fHistJets2CorrPtArea->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
+ if (naccJets2 < fNLeadingJets)
+ fHistLeadingJets2CorrPtArea->Fill(jet2->Area(), jet2->Pt() - fRho2Val * jet2->Area());
+ }
+
+ naccJets2++;
if (jet2->MatchedJet()) {
if (!AcceptBiasJet(jet2->MatchedJet()) ||
- jet2->MatchedJet()->MaxTrackPt() > fMaxTrackPt || jet2->MatchedJet()->MaxClusterPt() > fMaxClusterPt ||
- jet2->MatchedJet()->Pt() > fMaxBinPt) {
+ jet2->MatchedJet()->MaxTrackPt() > fMaxTrackPt || jet2->MatchedJet()->MaxClusterPt() > fMaxClusterPt) {
fHistMissedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
}
else {
+ if (jet2->MatchedJet()->Pt() > fMaxBinPt)
+ fHistMissedJets2PtArea->Fill(jet2->Area(), jet2->Pt());
+
Double_t d1=-1, d2=-1;
if (jet2->GetMatchingType() == kGeometrical) {
fHistDistancevsCommonEnergy1->Fill(jet2->ClosestJetDistance(), d1);
fHistDistancevsCommonEnergy2->Fill(jet2->ClosestJetDistance(), d2);
+
+ fHistDistancevsJet1Pt->Fill(jet2->ClosestJetDistance(), jet2->MatchedJet()->Pt());
+ fHistDistancevsJet2Pt->Fill(jet2->ClosestJetDistance(), jet2->Pt());
+
+ fHistCommonEnergy1vsJet1Pt->Fill(d1, jet2->MatchedJet()->Pt());
+ fHistCommonEnergy2vsJet2Pt->Fill(d2, jet2->Pt());
}
else if (jet2->GetMatchingType() == kMCLabel || jet2->GetMatchingType() == kSameCollections) {
GetGeometricalMatchingLevel(jet2->MatchedJet(), jet2, d1);
+
fHistDistancevsCommonEnergy1->Fill(d1, jet2->MatchedJet()->ClosestJetDistance());
fHistDistancevsCommonEnergy2->Fill(d1, jet2->ClosestJetDistance());
+
+ fHistDistancevsJet1Pt->Fill(d1, jet2->MatchedJet()->Pt());
+ fHistDistancevsJet2Pt->Fill(d1, jet2->Pt());
+
+ fHistCommonEnergy1vsJet1Pt->Fill(jet2->MatchedJet()->ClosestJetDistance(), jet2->MatchedJet()->Pt());
+ fHistCommonEnergy2vsJet2Pt->Fill(jet2->ClosestJetDistance(), jet2->Pt());
}
-
- fHistMatchingLevelvsJet2Pt->Fill(jet2->ClosestJetDistance(), jet2->Pt());
Double_t deta = jet2->MatchedJet()->Eta() - jet2->Eta();
Double_t dphi = jet2->MatchedJet()->Phi() - jet2->Phi();
- fHistDeltaEtaPhivsJet2Pt->Fill(deta, dphi, jet2->Pt());
+ fHistDeltaEtaPhi->Fill(deta, dphi, jet2->Pt());
Double_t dpt = jet2->MatchedJet()->Pt() - jet2->Pt();
fHistDeltaPtvsJet1Pt->Fill(jet2->MatchedJet()->Pt(), dpt);
}
}
+ GetSortedArray(indexes, fJets, fRhoVal);
+
const Int_t nJets1 = fJets->GetEntriesFast();
+ Int_t naccJets1 = 0;
for (Int_t i = 0; i < nJets1; i++) {
- AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(fJets->At(i));
+ AliDebug(2,Form("Processing jet %d", i));
+ AliEmcalJet* jet1 = static_cast<AliEmcalJet*>(fJets->At(indexes[i]));
if (!jet1) {
AliError(Form("Could not receive jet %d", i));
fHistJets1PtArea->Fill(jet1->Area(), jet1->Pt());
fHistJets1PhiEta->Fill(jet1->Eta(), jet1->Phi());
- if (!fRhoName.IsNull())
+ if (naccJets1 < fNLeadingJets)
+ fHistLeadingJets1PtArea->Fill(jet1->Area(), jet1->Pt());
+
+ if (!fRhoName.IsNull()) {
fHistJets1CorrPtArea->Fill(jet1->Area(), jet1->Pt() - fRhoVal * jet1->Area());
+
+ if (naccJets1 < fNLeadingJets)
+ fHistLeadingJets1CorrPtArea->Fill(jet1->Area(), jet1->Pt() - fRhoVal * jet1->Area());
+ }
+
+ naccJets1++;
}
return kTRUE;
class TClonesArray;
class TH1;
class TH2;
-class TH3;
#include "AliEmcalJet.h"
#include "AliAnalysisTaskEmcalJet.h"
};
void UserCreateOutputObjects();
+ Bool_t UserNotify();
void SetJets2Name(const char *n) { fJets2Name = n ; }
void SetTracks2Name(const char *n) { fTracks2Name = n ; }
void SetAreMCCollections(Bool_t f1, Bool_t f2) { fAreCollections1MC = f1; fAreCollections2MC = f2; }
protected:
+ Bool_t PythiaInfoFromFile(const char* currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard);
Bool_t IsEventSelected();
Bool_t AcceptJet(AliEmcalJet* jet) const;
Bool_t AcceptBiasJet2(AliEmcalJet *jet) const;
Bool_t fAreCollections1MC; // collections 1 MC
Bool_t fAreCollections2MC; // collections 1 MC
MatchingType fMatching; // matching type
- Double_t fMatchingPar1; // matching parameter for generated-reconstructed matching
- Double_t fMatchingPar2; // matching parameter for reconstructed-generated matching
+ Double_t fMatchingPar1; // matching parameter for jet1-jet2 matching
+ Double_t fMatchingPar2; // matching parameter for jet2-jet1 matching
Float_t fJet2MinEta; // minimum eta jet 2 acceptance
Float_t fJet2MaxEta; // maximum eta jet 2 acceptance
Float_t fJet2MinPhi; // minimum phi jet 2 acceptance
Double_t fRho2Val; //!Event rho 2 value
TH1 *fTracks2Map; //!MC particle map
// General histograms
- TH1 *fHistNTrials; //!total number of trials per pt hard bin
+ TH1 *fHistTrialsAfterSel; //!total number of trials per pt hard bin after selection
+ TH1 *fHistEventsAfterSel; //!total number of events per pt hard bin after selection
+ TH1 *fHistTrials; //!trials from pyxsec.root
+ TH1 *fHistXsection; //!x section from pyxsec.root
TH1 *fHistEvents; //!total number of events per pt hard bin
// Jets 1
TH2 *fHistJets1PhiEta; //!phi-eta distribution of jets 1
TH2 *fHistJets1PtArea; //!inclusive jet pt vs area histogram 1
TH2 *fHistJets1CorrPtArea; //!inclusive jet pt vs. area histogram 1
+ TH2 *fHistLeadingJets1PtArea; //!leading jet pt vs area histogram 1
+ TH2 *fHistLeadingJets1CorrPtArea; //!leading jet pt vs. area histogram 1
// Jets 2
TH2 *fHistJets2PhiEta; //!phi-eta distribution of jets 2
TH2 *fHistJets2PtArea; //!inclusive jet pt vs. area histogram 2
TH2 *fHistJets2CorrPtArea; //!inclusive jet pt vs. area histogram 2
+ TH2 *fHistLeadingJets2PtArea; //!leading jet pt vs. area histogram 2
+ TH2 *fHistLeadingJets2CorrPtArea; //!leading jet pt vs. area histogram 2
TH2 *fHistJets2PhiEtaAcceptance; //!phi-eta distribution of jets 2 using jet 1 cuts (acceptance, leading hadron bias, ...)
TH2 *fHistJets2PtAreaAcceptance; //!inclusive jet pt vs. area histogram 2 using jet 1 cuts (acceptance, leading hadron bias, ...)
TH2 *fHistJets2CorrPtAreaAcceptance; //!inclusive jet pt vs. area histogram 2 using jet 1 cuts (acceptance, leading hadron bias, ...)
+ TH2 *fHistLeadingJets2PtAreaAcceptance; //!leading jet pt vs. area histogram 2 using jet 1 cuts (acceptance, leading hadron bias, ...)
+ TH2 *fHistLeadingJets2CorrPtAreaAcceptance; //!leading jet pt vs. area histogram 2 using jet 1 cuts (acceptance, leading hadron bias, ...)
// Jet1-Jet2 matching
- TH2 *fHistMatchingLevelvsJet2Pt; //!matching level vs jet 2 pt
+ TH2 *fHistCommonEnergy1vsJet1Pt; //!common energy 1 (%) vs jet 1 pt
+ TH2 *fHistCommonEnergy2vsJet2Pt; //!common energy 2 (%) vs jet 2 pt
+ TH2 *fHistDistancevsJet1Pt; //!distance vs jet 1 pt
+ TH2 *fHistDistancevsJet2Pt; //!distance vs jet 2 pt
TH2 *fHistDistancevsCommonEnergy1; //!distance vs common energy 1 (%)
TH2 *fHistDistancevsCommonEnergy2; //!distance vs common energy 2 (%)
TH2 *fHistJet2PtOverJet1PtvsJet2Pt; //!jet 2 pt over jet 1 pt vs jet 2 pt
TH2 *fHistJet1PtOverJet2PtvsJet1Pt; //!jet 1 pt over jet 2 pt vs jet 1 pt
- TH3 *fHistDeltaEtaPhivsJet2Pt; //!delta eta-phi between matched jets vs jet 2 pt
+ TH2 *fHistDeltaEtaPhi; //!delta eta-phi between matched jets
TH2 *fHistDeltaPtvsJet1Pt; //!delta pt between matched jets vs jet 1 pt
TH2 *fHistDeltaPtvsJet2Pt; //!delta pt between matched jets vs jet 2 pt
TH2 *fHistDeltaPtvsMatchingLevel; //!delta pt between matched jets vs matching level
AliJetResponseMaker(const AliJetResponseMaker&); // not implemented
AliJetResponseMaker &operator=(const AliJetResponseMaker&); // not implemented
- ClassDef(AliJetResponseMaker, 11) // Jet response matrix producing task
+ ClassDef(AliJetResponseMaker, 12) // Jet response matrix producing task
};
#endif
for (Int_t i = 0; i < 4; i++) {
fHistEvents[i] = 0;
fHistLeadingJetPt[i] = 0;
- fHist2LeadingJetPt[i] = 0;
fHistLeadingJetCorrPt[i] = 0;
fHistRhoVSleadJetPt[i] = 0;
fHistJetPhiEta[i] = 0;
for (Int_t i = 0; i < 4; i++) {
fHistEvents[i] = 0;
fHistLeadingJetPt[i] = 0;
- fHist2LeadingJetPt[i] = 0;
fHistLeadingJetCorrPt[i] = 0;
fHistRhoVSleadJetPt[i] = 0;
fHistJetPhiEta[i] = 0;
fHistLeadingJetPt[i]->GetYaxis()->SetTitle("counts");
fOutput->Add(fHistLeadingJetPt[i]);
- histname = "fHist2LeadingJetPt_";
- histname += i;
- fHist2LeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
- fHist2LeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
- fHist2LeadingJetPt[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHist2LeadingJetPt[i]);
-
histname = "fHistLeadingJetCorrPt_";
histname += i;
fHistLeadingJetCorrPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
return kTRUE;
}
- Int_t *sortedJets = GetSortedArray(fJets);
+ static Int_t sortedJets[9999] = {-1};
+ GetSortedArray(sortedJets, fJets, fRhoVal);
- if (!sortedJets || sortedJets[0] < 0) { // no accepted jets, skipping
+ if (sortedJets[0] < 0) { // no accepted jets, skipping
fHistEvents[fCentBin]->Fill("No jets", 1);
return kTRUE;
}
- AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[0]));
- if (!jet) { // error, I cannot get the leading jet from collection (should never happen), skipping
- fHistEvents[fCentBin]->Fill("Max jet not found", 1);
- return kTRUE;
- }
-
// OK, event accepted
- Float_t maxJetCorrPt = jet->Pt() - fRhoVal * jet->Area();
-
if (fRhoVal == 0)
fHistEvents[fCentBin]->Fill("Rho == 0", 1);
- else if (maxJetCorrPt <= 0)
- fHistEvents[fCentBin]->Fill("Max jet <= 0", 1);
-
else
fHistEvents[fCentBin]->Fill("OK", 1);
- if (jet) {
- fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
- fHistRhoVSleadJetPt[fCentBin]->Fill(fRhoVal, jet->Pt());
- fHistLeadingJetCorrPt[fCentBin]->Fill(maxJetCorrPt);
- }
+ for (Int_t i = 0; i < fNLeadingJets && i < fJets->GetEntriesFast(); i++) {
+ AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[i]));
+
+ if (!jet) {
+ AliError(Form("Could not receive jet %d", sortedJets[i]));
+ continue;
+ }
+
+ if (!AcceptJet(jet))
+ continue;
- AliEmcalJet* jet2 = 0;
- if (sortedJets[1] >= 0)
- jet2 = static_cast<AliEmcalJet*>(fJets->At(sortedJets[1]));
+ Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
- if (jet2)
- fHist2LeadingJetPt[fCentBin]->Fill(jet2->Pt());
+ fHistLeadingJetCorrPt[fCentBin]->Fill(corrPt);
+ fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
+
+ if (i==0)
+ fHistRhoVSleadJetPt[fCentBin]->Fill(fRhoVal, jet->Pt());
+ }
Int_t njets = DoJetLoop();
// General histograms
TH1F *fHistEvents[4]; //!Events accepted/rejected
TH1F *fHistLeadingJetPt[4]; //!Leading jet pt spectrum
- TH1F *fHist2LeadingJetPt[4]; //!Second leading jet pt spectrum
TH1F *fHistLeadingJetCorrPt[4]; //!Corrected leading jet pt spectrum
TH2F *fHistRhoVSleadJetPt[4]; //!Area(leadjet) * rho vs. leading jet pt
TH2F *fNjetsVsCent; //!No. of jets vs. centrality
//-------------------------------------------------------
// Init the task and do settings
//-------------------------------------------------------
- TString name(Form("%s_%s_%s_%s_R0%d_",taskname,ntracks,nclusters,nrho,(Int_t)floor(jetradius*100+0.5)));
+ TString name;
+ if (strcmp(ntracks, "") == 0 && strcmp(nclusters, "") == 0)
+ name = Form("%s_%s_R0%d",taskname,nrho,(Int_t)floor(jetradius*100+0.5));
+ else if (strcmp(ntracks, "") == 0)
+ name = Form("%s_%s_%s_R0%d",taskname,nclusters,nrho,(Int_t)floor(jetradius*100+0.5));
+ else if (strcmp(nclusters, "") == 0)
+ name = Form("%s_%s_%s_R0%d",taskname,ntracks,nrho,(Int_t)floor(jetradius*100+0.5));
+ else
+ name = Form("%s_%s_%s_%s_R0%d",taskname,ntracks,nclusters,nrho,(Int_t)floor(jetradius*100+0.5));
+
if (type == AliAnalysisTaskEmcal::kTPC)
- name += "TPC";
+ name += "_TPC";
else if (type == AliAnalysisTaskEmcal::kEMCAL)
- name += "EMCAL";
+ name += "_EMCAL";
else if (type == AliAnalysisTaskEmcal::kUser)
- name += "USER";
+ name += "_USER";
AliAnalysisTaskDeltaPt* jetTask = new AliAnalysisTaskDeltaPt(name);
jetTask->SetAnaType(type);
AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
TList::Class(),AliAnalysisManager::kOutputContainer,
Form("%s", AliAnalysisManager::GetCommonFileName()));
- mgr->ConnectInput (jetTask, 0, cinput1 );
- mgr->ConnectOutput (jetTask, 1, coutput1 );
+ mgr->ConnectInput(jetTask, 0, cinput1);
+ mgr->ConnectOutput(jetTask, 1, coutput1);
return jetTask;
}
const Int_t nCells = 1234567890,
const Bool_t copyArray = kTRUE,
const Int_t nFiles = 1234567890,
+ const Bool_t makeQA = kFALSE,
const Double_t minPt = 0,
const Double_t maxPt = 1000,
const Double_t minEta = -0.9,
// Init the task and do settings
//-------------------------------------------------------
- AliJetEmbeddingFromPYTHIATask *jetEmb = new AliJetEmbeddingFromPYTHIATask(taskName);
+ AliJetEmbeddingFromPYTHIATask *jetEmb = new AliJetEmbeddingFromPYTHIATask(taskName,makeQA);
jetEmb->SetTracksName(tracksName);
jetEmb->SetClusName(clusName);
jetEmb->SetCellsName(cellsName);
// Create containers for input/output
mgr->ConnectInput(jetEmb, 0, mgr->GetCommonInputContainer());
+ if (makeQA) {
+ TString contName = taskName;
+ contName += "_histos";
+ AliAnalysisDataContainer *outc = mgr->CreateContainer(contName,
+ TList::Class(),
+ AliAnalysisManager::kOutputContainer,
+ "AnalysisResults.root");
+ mgr->ConnectOutput(jetEmb, 1, outc);
+ }
+
return jetEmb;
}
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff