#include <TH2C.h>
#include <TList.h>
#include <TStreamerInfo.h>
+#include <TRandom.h>
#include "AliVEvent.h"
#include "AliAODTrack.h"
#include "AliVTrack.h"
#include "AliVCluster.h"
#include "AliVCaloCells.h"
+#include "AliAODMCParticle.h"
#include "AliCentrality.h"
#include "AliVHeader.h"
#include "AliVVertex.h"
AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask() :
AliJetModelBaseTask("AliJetEmbeddingFromAODTask"),
fFileList(0),
+ fRandomAccess(kFALSE),
fAODTreeName(),
fAODHeaderName(),
fAODVertexName(),
fAODTrackName(),
fAODClusName(),
fAODCellsName(),
+ fAODMCParticlesName(),
fMinCentrality(0),
fMaxCentrality(10),
fTriggerMask(AliVEvent::kAny),
fZVertexCut(10),
fIncludeNoITS(kTRUE),
- fTotalFiles(2200),
+ fUseNegativeLabels(kTRUE),
+ fTrackEfficiency(1),
+ fIsMC(kFALSE),
+ fTotalFiles(2050),
+ fAttempts(5),
fEsdTreeMode(kFALSE),
fCurrentFileID(0),
- fCurrentAODFileID(-1),
+ fCurrentAODFileID(0),
fCurrentAODFile(0),
fPicoTrackVersion(0),
+ fCurrentAODTree(0),
fAODHeader(0),
fAODVertex(0),
fAODTracks(0),
fAODClusters(0),
fAODCaloCells(0),
- fHistFileIDs(0)
+ fAODMCParticles(0),
+ fCurrentAODEntry(-1),
+ fHistFileMatching(0),
+ fHistAODFileError(0),
+ fHistNotEmbedded(0)
{
// Default constructor.
SetSuffix("AODEmbedding");
SetMarkMC(0);
fAODfilterBits[0] = -1;
fAODfilterBits[1] = -1;
+ fEtaMin = -1;
+ fEtaMax = 1;
+ fPhiMin = -10;
+ fPhiMax = 10;
+ fPtMin = 0;
+ fPtMax = 1000;
}
//________________________________________________________________________
AliJetEmbeddingFromAODTask::AliJetEmbeddingFromAODTask(const char *name, Bool_t drawqa) :
AliJetModelBaseTask(name, drawqa),
fFileList(0),
+ fRandomAccess(kFALSE),
fAODTreeName("aodTree"),
fAODHeaderName("header"),
fAODVertexName("vertices"),
fAODTrackName("tracks"),
- fAODClusName("caloClusters"),
+ fAODClusName(""),
fAODCellsName("emcalCells"),
+ fAODMCParticlesName(AliAODMCParticle::StdBranchName()),
fMinCentrality(0),
fMaxCentrality(10),
fTriggerMask(AliVEvent::kAny),
fZVertexCut(10),
fIncludeNoITS(kTRUE),
- fTotalFiles(2200),
+ fUseNegativeLabels(kTRUE),
+ fTrackEfficiency(1),
+ fIsMC(kFALSE),
+ fTotalFiles(2050),
+ fAttempts(5),
fEsdTreeMode(kFALSE),
fCurrentFileID(0),
- fCurrentAODFileID(-1),
+ fCurrentAODFileID(0),
fCurrentAODFile(0),
fPicoTrackVersion(0),
+ fCurrentAODTree(0),
fAODHeader(0),
fAODVertex(0),
fAODTracks(0),
fAODClusters(0),
- fAODCaloCells(0),
- fHistFileIDs(0)
+ fAODCaloCells(0),
+ fAODMCParticles(0),
+ fCurrentAODEntry(-1),
+ fHistFileMatching(0),
+ fHistAODFileError(0),
+ fHistNotEmbedded(0)
{
// Standard constructor.
SetSuffix("AODEmbedding");
SetMarkMC(0);
fAODfilterBits[0] = -1;
fAODfilterBits[1] = -1;
+ fEtaMin = -1;
+ fEtaMax = 1;
+ fPhiMin = -10;
+ fPhiMax = 10;
+ fPtMin = 0;
+ fPtMax = 1000;
}
//________________________________________________________________________
AliJetModelBaseTask::UserCreateOutputObjects();
- fHistFileIDs = new TH2C("fHistFileIDs", "fHistFileIDs", fTotalFiles, -0.5, fTotalFiles - 0.5, fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
- fHistFileIDs->GetXaxis()->SetTitle("File no. (PYTHIA)");
- fHistFileIDs->GetYaxis()->SetTitle("File no. (Embedded AOD)");
- fOutput->Add(fHistFileIDs);
+ fHistFileMatching = new TH2C("fHistFileMatching", "fHistFileMatching", fTotalFiles, -0.5, fTotalFiles - 0.5, fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
+ fHistFileMatching->GetXaxis()->SetTitle("File no. (PYTHIA)");
+ fHistFileMatching->GetYaxis()->SetTitle("File no. (Embedded AOD)");
+ fHistFileMatching->GetZaxis()->SetTitle("counts");
+ fOutput->Add(fHistFileMatching);
+
+ fHistAODFileError = new TH1C("fHistAODFileError", "fHistAODFileError", fFileList->GetEntriesFast(), -0.5, fFileList->GetEntriesFast() -0.5);
+ fHistAODFileError->GetXaxis()->SetTitle("File no. (Embedded AOD)");
+ fHistAODFileError->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistAODFileError);
+
+ fHistNotEmbedded = new TH1C("fHistNotEmbedded", "fHistNotEmbedded", fTotalFiles, -0.5, fTotalFiles - 0.5);
+ fHistNotEmbedded->GetXaxis()->SetTitle("File no. (PYTHIA)");
+ fHistNotEmbedded->GetYaxis()->SetTitle("counts");
+ fOutput->Add(fHistNotEmbedded);
PostData(1, fOutput);
}
return kTRUE;
}
fCurrentFileID = path.Atoi();
- fCurrentAODFileID = fFileList->GetEntriesFast() * fCurrentFileID / fTotalFiles;
- AliInfo(Form("Start embedding from file ID %d", fCurrentAODFileID));
+ if (!fRandomAccess) {
+ fCurrentAODFileID = fFileList->GetEntriesFast() * fCurrentFileID / fTotalFiles-1;
+ AliInfo(Form("Start embedding from file ID %d", fCurrentAODFileID));
+ }
return kTRUE;
}
fCurrentAODFile = 0;
}
- if (fCurrentAODFileID >= fFileList->GetEntriesFast())
- return kFALSE;
+ Int_t i = 0;
+ TString fileName;
- TObjString *objFileName = static_cast<TObjString*>(fFileList->At(fCurrentAODFileID));
- TString fileName = objFileName->GetString();
-
- if (fileName.BeginsWith("alien://") && !gGrid) {
+ do {
+ if (i>0) {
+ AliDebug(3,Form("Failed to open file %s...", fileName.Data()));
+ if (fHistAODFileError)
+ fHistAODFileError->Fill(fCurrentAODFileID);
+ }
+
+ fileName = GetNextFileName();
+
+ if (fileName.IsNull())
+ break;
+
+ if (fileName.BeginsWith("alien://") && !gGrid) {
AliInfo("Trying to connect to AliEn ...");
TGrid::Connect("alien://");
- }
+ }
- fCurrentAODFile = TFile::Open(fileName);
+ AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
+ fCurrentAODFile = TFile::Open(fileName);
- if (!fCurrentAODFile || fCurrentAODFile->IsZombie()) {
+ i++;
+
+ } while ((!fCurrentAODFile || fCurrentAODFile->IsZombie()) && i < fAttempts);
+
+ if (!fCurrentAODFile || fCurrentAODFile->IsZombie())
return kFALSE;
- }
const TList *clist = fCurrentAODFile->GetStreamerInfoCache();
if(clist) {
fPicoTrackVersion = 0;
}
- if (fQAhistos)
- fHistFileIDs->Fill(fCurrentFileID, fCurrentAODFileID);
+ fCurrentAODTree = static_cast<TTree*>(fCurrentAODFile->Get(fAODTreeName));
+ if (!fCurrentAODTree)
+ return kFALSE;
+
+ if (!fAODHeaderName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODHeaderName, &fAODHeader);
+
+ if (!fAODVertexName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODVertexName, &fAODVertex);
+
+ if (!fAODTrackName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODTrackName, &fAODTracks);
+
+ if (!fAODClusName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODClusName, &fAODClusters);
+
+ if (!fAODCellsName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODCellsName, &fAODCaloCells);
+
+ if (!fAODMCParticlesName.IsNull())
+ fCurrentAODTree->SetBranchAddress(fAODMCParticlesName, &fAODMCParticles);
+
+ if (fRandomAccess)
+ fCurrentAODEntry = TMath::Nint(gRandom->Rndm()*fCurrentAODTree->GetEntries());
+ else
+ fCurrentAODEntry = 0;
+
+ AliDebug(2,Form("Will start embedding from entry %d", fCurrentAODEntry));
+
+ if (fHistFileMatching)
+ fHistFileMatching->Fill(fCurrentFileID, fCurrentAODFileID-1);
- fCurrentAODFileID++;
return kTRUE;
}
//________________________________________________________________________
-Bool_t AliJetEmbeddingFromAODTask::GetNextEntry()
+TString AliJetEmbeddingFromAODTask::GetNextFileName()
{
- static TTree *tree = 0;
- static Int_t entry = 0;
+ if (fRandomAccess)
+ fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*fFileList->GetEntriesFast());
+ else
+ fCurrentAODFileID++;
+ if (fCurrentAODFileID >= fFileList->GetEntriesFast())
+ return "";
+
+ TObjString *objFileName = static_cast<TObjString*>(fFileList->At(fCurrentAODFileID));
+ return objFileName->GetString();
+}
+
+//________________________________________________________________________
+Bool_t AliJetEmbeddingFromAODTask::GetNextEntry()
+{
Int_t attempts = 0;
do {
-
- if (!fCurrentAODFile || !tree || entry >= tree->GetEntries()) {
+ if (!fCurrentAODFile || !fCurrentAODTree || fCurrentAODEntry >= fCurrentAODTree->GetEntries()) {
if (!OpenNextFile())
return kFALSE;
-
- tree = static_cast<TTree*>(fCurrentAODFile->Get(fAODTreeName));
- if (!tree)
- return kFALSE;
-
- if (!fAODHeaderName.IsNull())
- tree->SetBranchAddress(fAODHeaderName, &fAODHeader);
-
- if (!fAODVertexName.IsNull())
- tree->SetBranchAddress(fAODVertexName, &fAODVertex);
-
- if (!fAODTrackName.IsNull())
- tree->SetBranchAddress(fAODTrackName, &fAODTracks);
-
- if (!fAODClusName.IsNull())
- tree->SetBranchAddress(fAODClusName, &fAODClusters);
-
- if (!fAODCellsName.IsNull())
- tree->SetBranchAddress(fAODCellsName, &fAODCaloCells);
-
- entry = 0;
}
- tree->GetEntry(entry);
- entry++;
+ fCurrentAODTree->GetEntry(fCurrentAODEntry);
+ fCurrentAODEntry++;
attempts++;
if (attempts == 1000)
AliWarning("After 1000 attempts no event has been accepted by the event selection (trigger, centrality...)!");
- } while (!IsAODEventSelected() && tree);
+ } while (!IsAODEventSelected() && fCurrentAODTree);
- return (tree!=0);
+ return (fCurrentAODTree!=0);
}
//________________________________________________________________________
if (!fEsdTreeMode && fAODHeader) {
AliAODHeader *aodHeader = static_cast<AliAODHeader*>(fAODHeader);
- 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));
- return kFALSE;
+ 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));
+ return kFALSE;
+ }
}
- 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));
- return kFALSE;
+ 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));
+ return kFALSE;
+ }
}
}
void AliJetEmbeddingFromAODTask::Run()
{
if (!GetNextEntry()) {
- AliError("Unable to get AOD event to embed. Nothing will be embedded.");
+ if (fHistNotEmbedded)
+ fHistNotEmbedded->Fill(fCurrentFileID);
+ AliError("Unable to get the AOD event to embed. Nothing will be embedded.");
return;
}
- if (fTracks && fOutTracks) {
+ if (fOutMCParticles) {
+
+ if (fCopyArray && fMCParticles)
+ CopyMCParticles();
+
+ if (fAODMCParticles) {
+ AliDebug(2, 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) {
+ AliError(Form("Could not find MC particle %d in branch %s of tree %s!", i, fAODMCParticlesName.Data(), fAODTreeName.Data()));
+ continue;
+ }
+
+ AliDebug(3, Form("Processing MC particle %d with pT = %f, eta = %f, phi = %f", i, part->Pt(), part->Eta(), part->Phi()));
+
+ if (!part->IsPhysicalPrimary())
+ continue;
+
+ if (part->Pt() < fPtMin || part->Pt() > fPtMax ||
+ part->Eta() < fEtaMin || part->Eta() > fEtaMax ||
+ part->Phi() < fPhiMin || part->Phi() > fPhiMax)
+ continue;
+
+ AddMCParticle(part, i);
+ AliDebug(3, "Embedded!");
+ }
+ }
+ }
+
+ if (fOutTracks) {
- if (fCopyArray)
+ if (fCopyArray && fTracks)
CopyTracks();
+ AliDebug(2, Form("Start embedding with %d tracks.", fOutTracks->GetEntriesFast()));
+
if (fAODTracks) {
AliDebug(2, Form("%d tracks will be processed for embedding.", fAODTracks->GetEntriesFast()));
for (Int_t i = 0; i < fAODTracks->GetEntriesFast(); i++) {
continue;
}
+ AliDebug(3, Form("Processing track %d with pT = %f, eta = %f, phi = %f, label = %d", i, track->Pt(), track->Eta(), track->Phi(), track->GetLabel()));
+
Int_t type = 0;
Bool_t isEmc = kFALSE;
if (!fEsdTreeMode) {
if ((aodtrack->GetStatus()&AliESDtrack::kITSrefit)==0) {
if (fIncludeNoITS)
type = 2;
- else
+ else {
+ AliDebug(3, "Track not embedded because ITS refit failed.");
continue;
}
+ }
else {
type = 1;
}
}
else { /*not a good track*/
+ AliDebug(3, "Track not embedded because not an hybrid track.");
continue;
}
}
}
- AliDebug(3, Form("Embedding track with pT = %f, eta = %f, phi = %f", track->Pt(), track->Eta(), track->Phi()));
- AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc);
+ if (track->Pt() < fPtMin || track->Pt() > fPtMax ||
+ track->Eta() < fEtaMin || track->Eta() > fEtaMax ||
+ track->Phi() < fPhiMin || track->Phi() > fPhiMax) {
+ AliDebug(3, "Track not embedded because out of limits.");
+ continue;
+ }
+
+ if (fTrackEfficiency < 1) {
+ Double_t r = gRandom->Rndm();
+ if (fTrackEfficiency < r) {
+ AliDebug(3, "Track not embedded because of artificial inefiiciency.");
+ continue;
+ }
+ }
+
+ Int_t label = 0;
+ if (fIsMC) {
+ 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;
+ }
+ }
+
+ AddTrack(track->Pt(), track->Eta(), track->Phi(), type, track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), isEmc, label);
+ AliDebug(3, "Track embedded!");
}
}
}
- if (fClusters && fOutClusters) {
+ if (fOutClusters) {
- if (fCopyArray)
+ if (fCopyArray && fClusters)
CopyClusters();
if (fAODClusters) {
TLorentzVector vect;
Double_t vert[3] = {0,0,0};
clus->GetMomentum(vect,vert);
- AddCluster(clus->E(), vect.Eta(), vect.Phi());
+
+ if (vect.Pt() < fPtMin || vect.Pt() > fPtMax ||
+ vect.Eta() < fEtaMin || vect.Eta() > fEtaMax ||
+ vect.Phi() < fPhiMin || vect.Phi() > fPhiMax)
+ continue;
+
+ AddCluster(clus->E(), vect.Eta(), vect.Phi(), clus->GetLabel());
}
}
}
- if (fCaloCells && fOutCaloCells) {
+ if (fOutCaloCells) {
+
+ Int_t totalCells = 0;
- Int_t totalCells = fCaloCells->GetNumberOfCells();
+ if (fCaloCells)
+ totalCells += fCaloCells->GetNumberOfCells();
if (fAODCaloCells)
totalCells += fAODCaloCells->GetNumberOfCells();
SetNumberOfOutCells(totalCells);
- CopyCells();
+ if (fCaloCells)
+ CopyCells();
if (fAODCaloCells) {
for (Short_t i = 0; i < fAODCaloCells->GetNumberOfCells(); i++) {
- Short_t mclabel = 0;
+ Int_t mclabel = 0;
Double_t efrac = 0.;
Double_t time = -1;
Short_t cellNum = -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));
- AddCell(amp, cellNum, time);
+ AddCell(amp, cellNum, time, mclabel);
}
}