-// $Id: AliJetEmbeddingFromAODTask.cxx $
+// $Id$
//
// Jet embedding from AOD task.
//
#include "AliJetEmbeddingFromAODTask.h"
+// C++ standard library
+#include <vector>
+
+// ROOT
#include <TFile.h>
#include <TTree.h>
#include <TClonesArray.h>
#include <TList.h>
#include <TStreamerInfo.h>
#include <TRandom.h>
+#include <TSystem.h>
+#include <TLorentzVector.h>
+// AliRoot
#include "AliVEvent.h"
#include "AliAODTrack.h"
#include "AliESDtrack.h"
#include "AliVHeader.h"
#include "AliVVertex.h"
#include "AliAODHeader.h"
+#include "AliFJWrapper.h"
#include "AliLog.h"
#include "AliInputEventHandler.h"
+#include "AliNamedString.h"
ClassImp(AliJetEmbeddingFromAODTask)
fAODClusName(),
fAODCellsName(),
fAODMCParticlesName(),
- fMinCentrality(0),
- fMaxCentrality(10),
+ fMinCentrality(-1),
+ fMaxCentrality(-1),
fTriggerMask(AliVEvent::kAny),
fZVertexCut(10),
- fIncludeNoITS(kTRUE),
+ fMaxVertexDist(999),
+ fJetMinPt(0),
+ fJetMinEta(-0.5),
+ fJetMaxEta(0.5),
+ fJetMinPhi(-999),
+ fJetMaxPhi(999),
+ fJetConstituentMinPt(0),
+ fJetRadius(0.4),
+ fJetType(0),
+ fJetAlgo(1),
+ fJetParticleLevel(kTRUE),
+ fParticleMinPt(0),
+ fParticleMaxPt(0),
+ fParticleSelection(0),
+ fIncludeNoITS(kFALSE),
+ fCutMaxFractionSharedTPCClusters(0.4),
fUseNegativeLabels(kTRUE),
fTrackEfficiency(1),
- fIsMC(kFALSE),
+ fIsAODMC(kFALSE),
fTotalFiles(2050),
fAttempts(5),
fEsdTreeMode(kFALSE),
fAODCaloCells(0),
fAODMCParticles(0),
fCurrentAODEntry(-1),
+ fAODFilePath(0),
fHistFileMatching(0),
fHistAODFileError(0),
- fHistNotEmbedded(0)
+ fHistNotEmbedded(0),
+ fHistEmbeddingQA(0),
+ fHistRejectedEvents(0),
+ fEmbeddingCount(0)
{
// Default constructor.
SetSuffix("AODEmbedding");
- SetMarkMC(0);
fAODfilterBits[0] = -1;
fAODfilterBits[1] = -1;
fEtaMin = -1;
fAODClusName(""),
fAODCellsName("emcalCells"),
fAODMCParticlesName(AliAODMCParticle::StdBranchName()),
- fMinCentrality(0),
- fMaxCentrality(10),
+ fMinCentrality(-1),
+ fMaxCentrality(-1),
fTriggerMask(AliVEvent::kAny),
fZVertexCut(10),
- fIncludeNoITS(kTRUE),
+ fMaxVertexDist(999),
+ fJetMinPt(0),
+ fJetMinEta(-0.5),
+ fJetMaxEta(0.5),
+ fJetMinPhi(-999),
+ fJetMaxPhi(999),
+ fJetConstituentMinPt(0),
+ fJetRadius(0.4),
+ fJetType(0),
+ fJetAlgo(1),
+ fJetParticleLevel(kTRUE),
+ fParticleMinPt(0),
+ fParticleMaxPt(0),
+ fParticleSelection(0),
+ fIncludeNoITS(kFALSE),
+ fCutMaxFractionSharedTPCClusters(0.4),
fUseNegativeLabels(kTRUE),
fTrackEfficiency(1),
- fIsMC(kFALSE),
+ fIsAODMC(kFALSE),
fTotalFiles(2050),
fAttempts(5),
fEsdTreeMode(kFALSE),
fAODCaloCells(0),
fAODMCParticles(0),
fCurrentAODEntry(-1),
+ fAODFilePath(0),
fHistFileMatching(0),
fHistAODFileError(0),
- fHistNotEmbedded(0)
+ fHistNotEmbedded(0),
+ fHistEmbeddingQA(0),
+ fHistRejectedEvents(0),
+ fEmbeddingCount(0)
{
// Standard constructor.
SetSuffix("AODEmbedding");
- SetMarkMC(0);
fAODfilterBits[0] = -1;
fAODfilterBits[1] = -1;
fEtaMin = -1;
fHistNotEmbedded->GetYaxis()->SetTitle("counts");
fOutput->Add(fHistNotEmbedded);
+ fHistEmbeddingQA = new TH1F("fHistEmbeddingQA", "fHistEmbeddingQA", 2, 0, 2);
+ fHistEmbeddingQA->GetXaxis()->SetTitle("Event state");
+ fHistEmbeddingQA->GetYaxis()->SetTitle("counts");
+ fHistEmbeddingQA->GetXaxis()->SetBinLabel(1, "OK");
+ fHistEmbeddingQA->GetXaxis()->SetBinLabel(2, "Not embedded");
+ fOutput->Add(fHistEmbeddingQA);
+
+ fHistRejectedEvents = new TH1F("fHistRejectedEvents", "fHistRejectedEvents", 500, -0.5, 499.5);
+ fHistRejectedEvents->GetXaxis()->SetTitle("# of rejected events");
+ fHistRejectedEvents->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistRejectedEvents);
+
PostData(1, fOutput);
}
else
fEsdTreeMode = kTRUE;
+ fAODFilePath = static_cast<AliNamedString*>(InputEvent()->FindListObject("AODEmbeddingFile"));
+ if (!fAODFilePath) {
+ fAODFilePath = new AliNamedString("AODEmbeddingFile", "");
+ AliDebug(3,"Adding AOD embedding file path object to the event list...");
+ InputEvent()->AddObject(fAODFilePath);
+ }
+
return AliJetModelBaseTask::ExecOnce();
}
}
Int_t i = 0;
- TString fileName;
- do {
- if (i>0) {
- AliDebug(3,Form("Failed to open file %s...", fileName.Data()));
- if (fHistAODFileError)
+ while ((!fCurrentAODFile || fCurrentAODFile->IsZombie()) && i < fAttempts) {
+ if (i > 0 && fHistAODFileError) {
fHistAODFileError->Fill(fCurrentAODFileID);
}
- fileName = GetNextFileName();
-
- if (fileName.IsNull())
- break;
-
- if (fileName.BeginsWith("alien://") && !gGrid) {
- AliInfo("Trying to connect to AliEn ...");
- TGrid::Connect("alien://");
- }
-
- AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
- fCurrentAODFile = TFile::Open(fileName);
-
+ fCurrentAODFile = GetNextFile();
i++;
+ }
- } while ((!fCurrentAODFile || fCurrentAODFile->IsZombie()) && i < fAttempts);
-
- if (!fCurrentAODFile || fCurrentAODFile->IsZombie())
+ if (!fCurrentAODFile || fCurrentAODFile->IsZombie()) {
+ AliError("Could not open AOD file to embed!");
return kFALSE;
+ }
const TList *clist = fCurrentAODFile->GetStreamerInfoCache();
if(clist) {
}
fCurrentAODTree = static_cast<TTree*>(fCurrentAODFile->Get(fAODTreeName));
- if (!fCurrentAODTree)
+ if (!fCurrentAODTree) {
+ AliError(Form("Could not get tree %s from file %s", fAODTreeName.Data(), fCurrentAODFile->GetName()));
return kFALSE;
+ }
if (!fAODHeaderName.IsNull())
fCurrentAODTree->SetBranchAddress(fAODHeaderName, &fAODHeader);
if (fRandomAccess)
fCurrentAODEntry = TMath::Nint(gRandom->Rndm()*fCurrentAODTree->GetEntries());
else
- fCurrentAODEntry = 0;
+ fCurrentAODEntry = -1;
- AliDebug(2,Form("Will start embedding from entry %d", fCurrentAODEntry));
+ AliDebug(3,Form("Will start embedding from entry %d", fCurrentAODEntry+1));
if (fHistFileMatching)
fHistFileMatching->Fill(fCurrentFileID, fCurrentAODFileID-1);
+
+ fEmbeddingCount = 0;
return kTRUE;
}
//________________________________________________________________________
-TString AliJetEmbeddingFromAODTask::GetNextFileName()
+TFile* AliJetEmbeddingFromAODTask::GetNextFile()
{
- if (fRandomAccess)
- fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*fFileList->GetEntriesFast());
- else
- fCurrentAODFileID++;
+ if (fRandomAccess)
+ fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*fFileList->GetEntriesFast());
+ else
+ fCurrentAODFileID++;
+
+ if (fCurrentAODFileID >= fFileList->GetEntriesFast()) {
+ AliError("No more file in the list!");
+ return 0;
+ }
+
+ TObjString *objFileName = static_cast<TObjString*>(fFileList->At(fCurrentAODFileID));
+ TString fileName(objFileName->GetString());
- if (fCurrentAODFileID >= fFileList->GetEntriesFast())
- return "";
-
- TObjString *objFileName = static_cast<TObjString*>(fFileList->At(fCurrentAODFileID));
- return objFileName->GetString();
+ if (fileName.BeginsWith("alien://") && !gGrid) {
+ AliInfo("Trying to connect to AliEn ...");
+ TGrid::Connect("alien://");
+ }
+
+ TString baseFileName(fileName);
+ if (baseFileName.Contains(".zip#")) {
+ Ssiz_t pos = baseFileName.Last('#');
+ baseFileName.Remove(pos);
+ }
+
+ if (gSystem->AccessPathName(baseFileName)) {
+ AliError(Form("File %s does not exist!", baseFileName.Data()));
+ return 0;
+ }
+
+ AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
+ TFile *file = TFile::Open(fileName);
+
+ if (!file || file->IsZombie()) {
+ AliError(Form("Unable to open file: %s!", fileName.Data()));
+ return 0;
+ }
+
+ return file;
}
//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::GetNextEntry()
{
- Int_t attempts = 0;
+ Int_t attempts = -1;
do {
- if (!fCurrentAODFile || !fCurrentAODTree || fCurrentAODEntry >= fCurrentAODTree->GetEntries()) {
- if (!OpenNextFile())
+ if (!fCurrentAODFile || !fCurrentAODTree || fCurrentAODEntry+1 >= fCurrentAODTree->GetEntries()) {
+ if (!OpenNextFile()) {
+ AliError("Could not open the next file!");
return kFALSE;
+ }
+ }
+
+ if (!fCurrentAODTree) {
+ AliError("Could not get the tree!");
+ return kFALSE;
}
- fCurrentAODTree->GetEntry(fCurrentAODEntry);
fCurrentAODEntry++;
- attempts++;
+ fCurrentAODTree->GetEntry(fCurrentAODEntry);
+ attempts++;
if (attempts == 1000)
AliWarning("After 1000 attempts no event has been accepted by the event selection (trigger, centrality...)!");
- } while (!IsAODEventSelected() && fCurrentAODTree);
+ } while (!IsAODEventSelected());
+
+ if (fHistRejectedEvents)
+ fHistRejectedEvents->Fill(attempts);
- return (fCurrentAODTree!=0);
+ if (!fCurrentAODTree)
+ return kFALSE;
+
+ fEmbeddingCount++;
+
+ return kTRUE;
}
//________________________________________________________________________
Bool_t AliJetEmbeddingFromAODTask::IsAODEventSelected()
{
+ // AOD event selection.
+
if (!fEsdTreeMode && fAODHeader) {
AliAODHeader *aodHeader = static_cast<AliAODHeader*>(fAODHeader);
+ // Trigger selection
if (fTriggerMask != AliVEvent::kAny) {
UInt_t offlineTrigger = aodHeader->GetOfflineTrigger();
if ((offlineTrigger & fTriggerMask) == 0) {
- AliDebug(2,Form("Event rejected due to physics selection. Event trigger mask: %d, trigger mask selection: %d.", offlineTrigger, fTriggerMask));
+ AliDebug(2,Form("Event rejected due to physics selection. Event trigger mask: %d, trigger mask selection: %d.",
+ offlineTrigger, fTriggerMask));
return kFALSE;
}
}
+ // Centrality selection
if (fMinCentrality >= 0) {
AliCentrality *cent = aodHeader->GetCentralityP();
Float_t centVal = cent->GetCentralityPercentile("V0M");
if (centVal < fMinCentrality || centVal >= fMaxCentrality) {
- AliDebug(2,Form("Event rejected due to centrality selection. Event centrality: %f, centrality range selection: %f to %f", centVal, fMinCentrality, fMaxCentrality));
+ AliDebug(2,Form("Event rejected due to centrality selection. Event centrality: %f, centrality range selection: %f to %f",
+ centVal, fMinCentrality, fMaxCentrality));
return kFALSE;
}
}
}
+ // Vertex selection
if (fAODVertex) {
Double_t vert[3]={0};
((AliVVertex*)fAODVertex->At(0))->GetXYZ(vert);
if (TMath::Abs(vert[2]) > fZVertexCut) {
- AliDebug(2,Form("Event rejected due to Z vertex selection. Event Z vertex: %f, Z vertex cut: %f", vert[2], fZVertexCut));
+ AliDebug(2,Form("Event rejected due to Z vertex selection. Event Z vertex: %f, Z vertex cut: %f",
+ vert[2], fZVertexCut));
+ return kFALSE;
+ }
+ Double_t dist = TMath::Sqrt((vert[0]-fVertex[0])*(vert[0]-fVertex[0])+(vert[1]-fVertex[1])*(vert[1]-fVertex[1])+(vert[2]-fVertex[2])*(vert[2]-fVertex[2]));
+ if (dist > fMaxVertexDist) {
+ AliDebug(2,Form("Event rejected because the distance between the current and embedded verteces is > %f. "
+ "Current event vertex (%f, %f, %f), embedded event vertex (%f, %f, %f). Distance = %f",
+ fMaxVertexDist, fVertex[0], fVertex[1], fVertex[2], vert[0], vert[1], vert[2], dist));
return kFALSE;
}
+
+ }
+
+ // Particle selection
+ if ((fParticleSelection == 1 && FindParticleInRange(fAODTracks)==kFALSE) ||
+ (fParticleSelection == 2 && FindParticleInRange(fAODClusters)==kFALSE) ||
+ (fParticleSelection == 3 && FindParticleInRange(fAODMCParticles)==kFALSE))
+ return kFALSE;
+
+ // Jet selection
+ if (fJetMinPt > 0) {
+ TLorentzVector jet;
+
+ if (fJetParticleLevel) {
+ if (fAODMCParticles)
+ jet = GetLeadingJet(fAODMCParticles);
+ else {
+ AliWarning("Particle level jets selected, but not MC particles found. The jet event selection will be skipped.");
+ return kTRUE;
+ }
+ }
+ else {
+ if (fAODTracks || fAODClusters)
+ jet = GetLeadingJet(fAODTracks, fAODClusters);
+ else {
+ AliWarning("Detector level jets selected, but not tracks or clusters found. The jet event selection will be skipped.");
+ return kTRUE;
+ }
+ }
+
+ AliDebug(1, Form("Leading jet pt = %f", jet.Pt()));
+ if (jet.Pt() < fJetMinPt)
+ return kFALSE;
}
return kTRUE;
}
+//________________________________________________________________________
+Bool_t AliJetEmbeddingFromAODTask::FindParticleInRange(TClonesArray *array)
+{
+ if (!array) return kFALSE;
+
+ if (array->GetClass()->InheritsFrom("AliVParticle")) {
+ const Int_t nentries = array->GetEntriesFast();
+ for (Int_t i = 0; i < nentries; i++) {
+ AliVParticle *part = static_cast<AliVParticle*>(array->At(i));
+ if (!part) continue;
+ if (part->Pt() > fParticleMinPt && part->Pt() < fParticleMaxPt) return kTRUE;
+ }
+ }
+ else if (array->GetClass()->InheritsFrom("AliVCluster")) {
+ const Int_t nentries = array->GetEntriesFast();
+ for (Int_t i = 0; i < nentries; i++) {
+ AliVCluster *clus = static_cast<AliVCluster*>(array->At(i));
+ if (!clus) continue;
+
+ TLorentzVector vect;
+ Double_t vert[3] = {0,0,0};
+ clus->GetMomentum(vect,vert);
+
+ if (vect.Pt() > fParticleMinPt && vect.Pt() < fParticleMaxPt) return kTRUE;
+ }
+ }
+ else {
+ AliWarning(Form("Unable to do event selection based on particle pT: %s class type not recognized.",array->GetClass()->GetName()));
+ }
+
+ return kFALSE;
+}
+
//________________________________________________________________________
void AliJetEmbeddingFromAODTask::Run()
{
if (!GetNextEntry()) {
if (fHistNotEmbedded)
fHistNotEmbedded->Fill(fCurrentFileID);
+ if (fHistEmbeddingQA)
+ fHistEmbeddingQA->Fill("Not embedded", 1);
AliError("Unable to get the AOD event to embed. Nothing will be embedded.");
return;
}
+ if (fHistEmbeddingQA)
+ fHistEmbeddingQA->Fill("OK", 1);
+
+ fAODFilePath->SetString(fCurrentAODFile->GetName());
+
if (fOutMCParticles) {
if (fCopyArray && fMCParticles)
CopyMCParticles();
if (fAODMCParticles) {
- AliDebug(2, Form("%d MC particles will be processed for embedding.", fAODMCParticles->GetEntriesFast()));
+ AliDebug(3, Form("%d MC particles will be processed for embedding.", fAODMCParticles->GetEntriesFast()));
for (Int_t i = 0; i < fAODMCParticles->GetEntriesFast(); i++) {
AliAODMCParticle *part = static_cast<AliAODMCParticle*>(fAODMCParticles->At(i));
if (!part) {
if (fCopyArray && fTracks)
CopyTracks();
- AliDebug(2, Form("Start embedding with %d tracks.", fOutTracks->GetEntriesFast()));
+ AliDebug(3, Form("Start embedding with %d tracks.", fOutTracks->GetEntriesFast()));
if (fAODTracks) {
- AliDebug(2, Form("%d tracks will be processed for embedding.", fAODTracks->GetEntriesFast()));
+ AliDebug(3, Form("%d tracks will be processed for embedding.", fAODTracks->GetEntriesFast()));
for (Int_t i = 0; i < fAODTracks->GetEntriesFast(); i++) {
AliVTrack *track = static_cast<AliVTrack*>(fAODTracks->At(i));
if (!track) {
else {
AliDebug(3, "Track not embedded because ITS refit failed.");
continue;
- }
+ }
}
else {
type = 1;
AliDebug(3, "Track not embedded because not an hybrid track.");
continue;
}
-
+ if (fCutMaxFractionSharedTPCClusters > 0) {
+ Double_t frac = Double_t(aodtrack->GetTPCnclsS()) / Double_t(aodtrack->GetTPCncls());
+ if (frac > fCutMaxFractionSharedTPCClusters)
+ continue;
+ }
if (TMath::Abs(aodtrack->GetTrackEtaOnEMCal()) < 0.75 &&
aodtrack->GetTrackPhiOnEMCal() > 70 * TMath::DegToRad() &&
aodtrack->GetTrackPhiOnEMCal() < 190 * TMath::DegToRad())
else
type = ptrack->GetLabel();
isEmc = ptrack->IsEMCAL();
+
+ if (!fIncludeNoITS && type==2) {
+ AliDebug(3, "Track not embedded because ITS refit failed.");
+ continue;
+ }
}
}
if (fTrackEfficiency < r) {
AliDebug(3, "Track not embedded because of artificial inefiiciency.");
continue;
- }
+ }
}
Int_t label = 0;
- if (fIsMC) {
+ if (fIsAODMC) {
if (fUseNegativeLabels)
label = track->GetLabel();
else
label = TMath::Abs(track->GetLabel());
- if (label == 0) {
- AliDebug(2,Form("%s: Track %d with label==0", GetName(), i));
- label = 99999;
- }
+ if (label == 0)
+ AliDebug(1,Form("%s: Track %d with label==0", GetName(), i));
}
- AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc, label);
+ AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), track->GetTrackPtOnEMCal(), isEmc, label, track->Charge());
AliDebug(3, "Track embedded!");
}
}
AliError(Form("Could not find cluster %d in branch %s of tree %s!", i, fAODClusName.Data(), fAODTreeName.Data()));
continue;
}
+
+ if (!clus->IsEMCAL())
+ continue;
+
TLorentzVector vect;
Double_t vert[3] = {0,0,0};
clus->GetMomentum(vect,vert);
vect.Eta() < fEtaMin || vect.Eta() > fEtaMax ||
vect.Phi() < fPhiMin || vect.Phi() > fPhiMax)
continue;
-
- AddCluster(clus->E(), vect.Eta(), vect.Phi(), clus->GetLabel());
+
+ Int_t label = 0;
+ if (fIsAODMC) {
+ label = clus->GetLabel();
+ if (label <= 0)
+ AliDebug(1,Form("%s: Clus %d with label<=0", GetName(), i));
+ }
+ AddCluster(clus);
}
}
}
if (fOutCaloCells) {
+ Double_t totalEnergy = 0;
Int_t totalCells = 0;
if (fCaloCells)
Double_t amp = -1;
fAODCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
- AliDebug(3,Form("Adding cell with amplitude %f, absolute ID %d, time %f", amp, cellNum, time));
+
+ if (fIsAODMC) {
+ if (mclabel <= 0)
+ AliDebug(1,Form("%s: Cell %d with label<=0", GetName(), i));
+ }
+ else {
+ mclabel = 0;
+ }
+
+ AliDebug(2,Form("Adding cell with amplitude %f, absolute ID %d, time %f, mc label %d", amp, cellNum, time, mclabel));
AddCell(amp, cellNum, time, mclabel);
+ totalEnergy += amp;
}
}
+ AliDebug(2,Form("Added cells = %d (energy = %f), total cells = %d", fAddedCells, totalEnergy, totalCells));
+ }
+}
+
+//________________________________________________________________________
+TLorentzVector AliJetEmbeddingFromAODTask::GetLeadingJet(TClonesArray *tracks, TClonesArray *clusters)
+{
+ TString name("kt");
+ fastjet::JetAlgorithm jalgo(fastjet::kt_algorithm);
+ if (fJetAlgo == 1) {
+ name = "antikt";
+ jalgo = fastjet::antikt_algorithm;
+ }
+
+ // setup fj wrapper
+ AliFJWrapper fjw(name, name);
+ fjw.SetAreaType(fastjet::active_area_explicit_ghosts);
+ fjw.SetGhostArea(1); // set a very large ghost area to speed up jet finding
+ fjw.SetR(fJetRadius);
+ fjw.SetAlgorithm(jalgo);
+ fjw.SetMaxRap(1);
+ fjw.Clear();
+
+ if (tracks) {
+ const Int_t Ntracks = tracks->GetEntries();
+ for (Int_t iTracks = 0; iTracks < Ntracks; ++iTracks) {
+ AliVParticle *t = static_cast<AliVParticle*>(tracks->At(iTracks));
+ if (!t)
+ continue;
+
+ AliAODMCParticle *aodmcpart = dynamic_cast<AliAODMCParticle*>(t);
+ if (aodmcpart && !aodmcpart->IsPhysicalPrimary())
+ continue;
+
+ if ((fJetType == 1 && t->Charge() == 0) ||
+ (fJetType == 2 && t->Charge() != 0))
+ continue;
+
+ if (t->Pt() < fJetConstituentMinPt)
+ continue;
+
+ fjw.AddInputVector(t->Px(), t->Py(), t->Pz(), t->P(), iTracks + 100);
+ }
+ }
+
+ if (clusters && fJetType != 1) {
+ Double_t vert[3]={0};
+ if (fAODVertex)
+ ((AliVVertex*)fAODVertex->At(0))->GetXYZ(vert);
+
+ const Int_t Nclus = clusters->GetEntries();
+ for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
+ AliVCluster *c = static_cast<AliVCluster*>(clusters->At(iClus));
+ if (!c)
+ continue;
+
+ if (!c->IsEMCAL())
+ continue;
+
+ TLorentzVector nP;
+ c->GetMomentum(nP, vert);
+
+ if (nP.Pt() < fJetConstituentMinPt)
+ continue;
+
+ fjw.AddInputVector(nP.Px(), nP.Py(), nP.Pz(), nP.P(), -iClus - 100);
+ }
+ }
+
+ // run jet finder
+ fjw.Run();
+
+ std::vector<fastjet::PseudoJet> jets_incl = fjw.GetInclusiveJets();
+ AliDebug(1,Form("%d jets found", (Int_t)jets_incl.size()));
+
+ TLorentzVector jet;
- AliDebug(2,Form("Added cells = %d, total cells = %d", fAddedCells, totalCells));
+ Int_t njets = jets_incl.size();
+
+ if (njets > 0) {
+ //std::vector<fastjet::PseudoJet> jets_incl_sorted = fastjet::sorted_by_pt(jets_incl);
+ for (Int_t i = 0; i < njets; i++) {
+ if (jet.Pt() >= jets_incl[i].perp())
+ continue;
+ if (jets_incl[i].eta() < fJetMinEta || jets_incl[i].eta() > fJetMaxEta || jets_incl[i].phi() < fJetMinPhi || jets_incl[i].phi() > fJetMaxPhi)
+ continue;
+ jet.SetPxPyPzE(jets_incl[i].px(), jets_incl[i].py(), jets_incl[i].pz(), jets_incl[i].E());
+ }
}
+
+ return jet;
}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff