#include "AliJetModelBaseTask.h"
#include <TClonesArray.h>
-#include <TH1.h>
+#include <TH1I.h>
#include <TLorentzVector.h>
#include <TRandom3.h>
#include <TList.h>
+#include <TF2.h>
#include "AliVEvent.h"
#include "AliAODCaloCluster.h"
#include "AliEMCALRecPoint.h"
#include "AliESDCaloCells.h"
#include "AliAODCaloCells.h"
+#include "AliAODMCParticle.h"
#include "AliVCaloCells.h"
#include "AliPicoTrack.h"
#include "AliEMCALGeometry.h"
#include "AliLog.h"
+#include "AliNamedArrayI.h"
ClassImp(AliJetModelBaseTask)
fOutCaloName(),
fCellsName(),
fOutCellsName(),
+ fMCParticlesName(),
+ fOutMCParticlesName(),
+ fPythiaInfoName(""),
+ fIsMC(kFALSE),
fSuffix(),
fEtaMin(-1),
fEtaMax(1),
fNClusters(0),
fNCells(0),
fNTracks(0),
- fMarkMC(kTRUE),
+ fMarkMC(99999),
fPtSpectrum(0),
+ fPtPhiEvPlDistribution(0),
+ fDensitySpectrum(0),
+ fDifferentialV2(0),
+ fAddV2(kFALSE),
+ fFlowFluctuations(kFALSE),
fQAhistos(kFALSE),
+ fPsi(0),
fIsInit(0),
fGeom(0),
fClusters(0),
fCaloCells(0),
fOutCaloCells(0),
fAddedCells(0),
+ fMCParticles(0),
+ fMCParticlesMap(0),
+ fOutMCParticles(0),
+ fOutMCParticlesMap(0),
+ fMCLabelShift(0),
fEsdMode(kFALSE),
- fOutput(0)
+ fOutput(0),
+ fPythiaInfo(0x0)
{
// Default constructor.
+
+ fVertex[0] = 0;
+ fVertex[1] = 0;
+ fVertex[2] = 0;
}
//________________________________________________________________________
fOutCaloName("CaloClustersCorrEmbedded"),
fCellsName(""),
fOutCellsName(""),
+ fMCParticlesName(""),
+ fOutMCParticlesName(""),
+ fPythiaInfoName(""),
+ fIsMC(kFALSE),
fSuffix("Processed"),
fEtaMin(-1),
fEtaMax(1),
fNClusters(0),
fNCells(0),
fNTracks(1),
- fMarkMC(kTRUE),
+ fMarkMC(99999),
fPtSpectrum(0),
+ fPtPhiEvPlDistribution(0),
+ fDensitySpectrum(0),
+ fDifferentialV2(0),
+ fAddV2(kFALSE),
+ fFlowFluctuations(kFALSE),
fQAhistos(drawqa),
+ fPsi(0),
fIsInit(0),
fGeom(0),
fClusters(0),
fCaloCells(0),
fOutCaloCells(0),
fAddedCells(0),
+ fMCParticles(0),
+ fMCParticlesMap(0),
+ fOutMCParticles(0),
+ fOutMCParticlesMap(0),
+ fMCLabelShift(0),
fEsdMode(kFALSE),
- fOutput(0)
+ fOutput(0),
+ fPythiaInfo(0x0)
{
// Standard constructor.
if (fQAhistos) {
DefineOutput(1, TList::Class());
}
+
+ fVertex[0] = 0;
+ fVertex[1] = 0;
+ fVertex[2] = 0;
}
//________________________________________________________________________
if (!fIsInit)
return;
+ fVertex[0] = 0;
+ fVertex[1] = 0;
+ fVertex[2] = 0;
+
+ const AliVVertex *vert = InputEvent()->GetPrimaryVertex();
+ if (vert)
+ vert->GetXYZ(fVertex);
+
if (fCopyArray) {
if (fOutTracks)
fOutTracks->Delete();
if (fOutClusters)
fOutClusters->Delete();
+ if (fOutMCParticles)
+ fOutMCParticles->Delete();
}
+ if (fDensitySpectrum) {
+ fNTracks = TMath::Nint(fDensitySpectrum->GetRandom());
+ fNCells = TMath::Nint(fDensitySpectrum->GetRandom());
+ fNClusters = TMath::Nint(fDensitySpectrum->GetRandom());
+ }
+
+ // Clear map
+ if (fOutMCParticlesMap)
+ fOutMCParticlesMap->Clear();
+
AliVCaloCells *tempCaloCells = 0;
if (fCaloCells) {
}
}
+ if (fPtPhiEvPlDistribution || fAddV2)
+ fPsi = gRandom->Rndm() * TMath::Pi();
+
Run();
-
- if (fCaloCells) {
+
+ if (fCaloCells && !fCopyArray) {
delete fCaloCells;
fCaloCells = tempCaloCells;
}
fPhiMax = fPhiMin;
}
- if (!fCaloCells && !fCellsName.IsNull()) {
+ if (!fCellsName.IsNull()) {
fCaloCells = dynamic_cast<AliVCaloCells*>(InputEvent()->FindListObject(fCellsName));
if (!fCaloCells) {
- AliError(Form("%s: Couldn't retrieve calo cells with name %s!", GetName(), fCellsName.Data()));
- return kFALSE;
+ AliWarning(Form("%s: Couldn't retrieve calo cells with name %s!", GetName(), fCellsName.Data()));
}
-
- if (!fCaloCells->InheritsFrom("AliVCaloCells")) {
+ else if (!fCaloCells->InheritsFrom("AliVCaloCells")) {
AliError(Form("%s: Collection %s does not contain a AliVCaloCells object!", GetName(), fCellsName.Data()));
fCaloCells = 0;
return kFALSE;
if (!fOutCaloCells) {
fOutCellsName = fCellsName;
- if (fCopyArray) {
+ if (fCopyArray)
fOutCellsName += fSuffix;
-
+ if (fCopyArray || !fCaloCells) {
if (fEsdMode)
fOutCaloCells = new AliESDCaloCells(fOutCellsName,fOutCellsName);
else
}
}
- if (fNTracks > 0 && !fTracks && !fTracksName.IsNull()) {
+ if (!fTracksName.IsNull()) {
fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
if (!fTracks) {
- AliError(Form("%s: Couldn't retrieve tracks with name %s!", GetName(), fTracksName.Data()));
- return kFALSE;
+ AliWarning(Form("%s: Couldn't retrieve tracks with name %s!", GetName(), fTracksName.Data()));
}
-
- if (!fTracks->GetClass()->GetBaseClass("AliPicoTrack")) {
+ else if (!fTracks->GetClass()->GetBaseClass("AliPicoTrack")) {
AliError(Form("%s: Collection %s does not contain AliPicoTrack objects!", GetName(), fTracksName.Data()));
fTracks = 0;
return kFALSE;
if (!fOutTracks) {
fOutTracksName = fTracksName;
- if (fCopyArray) {
+ if (fCopyArray)
fOutTracksName += fSuffix;
- fOutTracks = new TClonesArray("AliPicoTrack", fTracks->GetSize());
+ if (fCopyArray || !fTracks) {
+ fOutTracks = new TClonesArray("AliPicoTrack");
fOutTracks->SetName(fOutTracksName);
if (InputEvent()->FindListObject(fOutTracksName)) {
AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutTracksName.Data()));
}
}
- if (fNClusters > 0 && !fClusters && !fCaloName.IsNull()) {
+ if(fAddV2 && (!fDifferentialV2)) {
+ AliWarning(Form("%s: Cannot add v2 without diffential v2!", GetName()));
+ }
+
+ if (!fCaloName.IsNull()) {
fClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
if (!fClusters) {
- AliError(Form("%s: Couldn't retrieve clusters with name %s!", GetName(), fCaloName.Data()));
- return kFALSE;
+ AliWarning(Form("%s: Couldn't retrieve clusters with name %s!", GetName(), fCaloName.Data()));
}
-
- if (!fClusters->GetClass()->GetBaseClass("AliVCluster")) {
+ else if (!fClusters->GetClass()->GetBaseClass("AliVCluster")) {
AliError(Form("%s: Collection %s does not contain AliVCluster objects!", GetName(), fCaloName.Data()));
fClusters = 0;
return kFALSE;
if (!fOutClusters) {
fOutCaloName = fCaloName;
- if (fCopyArray) {
+ if (fCopyArray)
fOutCaloName += fSuffix;
- fOutClusters = new TClonesArray(fClusters->GetClass()->GetName(), fClusters->GetSize());
+ TString className;
+ if (fClusters)
+ className = fClusters->GetClass()->GetName();
+ else if (fEsdMode)
+ className = "AliESDCaloCluster";
+ else
+ className = "AliAODCaloCluster";
+
+ if (fCopyArray || !fClusters) {
+ fOutClusters = new TClonesArray(className.Data());
fOutClusters->SetName(fOutCaloName);
if (InputEvent()->FindListObject(fOutCaloName)) {
AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutCaloName.Data()));
}
}
- if (!fCaloName.IsNull() || !fCellsName.IsNull()) {
- if (!fGeom) {
- if (fGeomName.Length() > 0) {
- fGeom = AliEMCALGeometry::GetInstance(fGeomName);
- if (!fGeom)
- AliError(Form("Could not get geometry with name %s!", fGeomName.Data()));
- } else {
- fGeom = AliEMCALGeometry::GetInstance();
- if (!fGeom)
- AliError("Could not get default geometry!");
+ if (!fMCParticlesName.IsNull()) {
+ fMCParticles = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fMCParticlesName));
+ if (!fMCParticles) {
+ AliWarning(Form("%s: Couldn't retrieve MC particles with name %s!", GetName(), fMCParticlesName.Data()));
+ }
+ else {
+ if (!fMCParticles->GetClass()->GetBaseClass("AliAODMCParticle")) {
+ AliError(Form("%s: Collection %s does not contain AliAODMCParticle objects!", GetName(), fMCParticlesName.Data()));
+ fMCParticles = 0;
+ return kFALSE;
+ }
+
+ fMCParticlesMap = dynamic_cast<AliNamedArrayI*>(InputEvent()->FindListObject(fMCParticlesName + "_Map"));
+
+ if (!fMCParticlesMap) {
+ AliWarning(Form("%s: Could not retrieve map for MC particles %s! Will assume MC labels consistent with indexes...", GetName(), fMCParticlesName.Data()));
+ fMCParticlesMap = new AliNamedArrayI(fMCParticlesName + "_Map", 99999);
+ for (Int_t i = 0; i < 99999; i++) {
+ fMCParticlesMap->AddAt(i,i);
+ }
+ }
+ }
+
+ if (!fOutMCParticles) {
+ fOutMCParticlesName = fMCParticlesName;
+ if (fCopyArray)
+ fOutMCParticlesName += fSuffix;
+ if (fCopyArray || !fMCParticles) {
+
+ fOutMCParticles = new TClonesArray("AliAODMCParticle");
+ fOutMCParticles->SetName(fOutMCParticlesName);
+ if (InputEvent()->FindListObject(fOutMCParticlesName)) {
+ AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fOutMCParticlesName.Data()));
+ return kFALSE;
+ }
+ else {
+ InputEvent()->AddObject(fOutMCParticles);
+ }
+
+ fOutMCParticlesMap = new AliNamedArrayI(fOutMCParticlesName + "_Map",99999);
+ if (InputEvent()->FindListObject(fOutMCParticlesName + "_Map")) {
+ AliFatal(Form("%s: Map %s_Map is already present in the event!", GetName(), fOutMCParticlesName.Data()));
+ return kFALSE;
+ }
+ else {
+ InputEvent()->AddObject(fOutMCParticlesMap);
+ }
+ }
+ else {
+ fOutMCParticles = fMCParticles;
+ fOutMCParticlesMap = fMCParticlesMap;
+ }
+ }
+ }
+
+ if (!fGeom) {
+ if (fGeomName.Length() > 0) {
+ fGeom = AliEMCALGeometry::GetInstance(fGeomName);
+ if (!fGeom) {
+ AliFatal(Form("Could not get geometry with name %s!", fGeomName.Data()));
+ return kFALSE;
+ }
+ } else {
+ fGeom = AliEMCALGeometry::GetInstance();
+ if (!fGeom) {
+ AliFatal("Could not get default geometry!");
+ return kFALSE;
}
}
-
- const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
- const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
- const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
- const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
-
- if (fEtaMax > EmcalMaxEta) fEtaMax = EmcalMaxEta;
- if (fEtaMax < EmcalMinEta) fEtaMax = EmcalMinEta;
- if (fEtaMin > EmcalMaxEta) fEtaMin = EmcalMaxEta;
- if (fEtaMin < EmcalMinEta) fEtaMin = EmcalMinEta;
-
- if (fPhiMax > EmcalMaxPhi) fPhiMax = EmcalMaxPhi;
- if (fPhiMax < EmcalMinPhi) fPhiMax = EmcalMinPhi;
- if (fPhiMin > EmcalMaxPhi) fPhiMin = EmcalMaxPhi;
- if (fPhiMin < EmcalMinPhi) fPhiMin = EmcalMinPhi;
}
return kTRUE;
return n;
}
-//________________________________________________________________________
-void AliJetModelBaseTask::CopyCells()
-{
- for (Short_t i = 0; i < fCaloCells->GetNumberOfCells(); i++) {
- Short_t mclabel = -1;
- Double_t efrac = 0.;
- Double_t time = -1;
- Short_t cellNum = -1;
- Double_t amp = -1;
-
- fCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
- fOutCaloCells->SetCell(i, cellNum, amp, time, mclabel, efrac);
- }
-
- fAddedCells = fCaloCells->GetNumberOfCells();
-
- AliDebug(2, Form("%d cells from the PYTHIA event", fAddedCells));
-}
-
//________________________________________________________________________
Int_t AliJetModelBaseTask::AddCell(Double_t e, Double_t eta, Double_t phi)
{
}
//________________________________________________________________________
-Int_t AliJetModelBaseTask::AddCell(Double_t e, Int_t absId, Double_t time)
+Int_t AliJetModelBaseTask::AddCell(Double_t e, Int_t absId, Double_t time, Int_t label)
{
// Add a cell to the event.
- Int_t label = fMarkMC ? 100 : -1;
+ if (label < 0)
+ label = 0;
+
+ label += fMarkMC + fMCLabelShift;
+
+ Short_t pos = -1;
+ if (fCaloCells)
+ pos = fCaloCells->GetCellPosition(absId);
+
+ Double_t efrac = 1;
+ Bool_t increaseOnSuccess = kFALSE;
+
+ if (pos < 0) {
+ increaseOnSuccess = kTRUE;
+ pos = fAddedCells;
+ }
+ else {
+ Short_t cellNumber = -1;
+ Double_t old_e = 0;
+ Double_t old_time = 0;
+ Int_t old_label = 0;
+ Double_t old_efrac = 0;
+ fOutCaloCells->GetCell(pos, cellNumber, old_e, old_time, old_label, old_efrac);
+
+ efrac = e / (old_e + e);
+
+ if (old_label > 0 && e < old_efrac * old_e) {
+ label = old_label;
+ efrac = old_efrac;
+ time = old_time;
+ }
+
+ e += old_e;
+ }
- Bool_t r = fOutCaloCells->SetCell(fAddedCells, absId, e, time, label, 0);
+ Bool_t r = fOutCaloCells->SetCell(pos, absId, e, time, label, efrac);
if (r) {
- fAddedCells++;
+ if (increaseOnSuccess)
+ fAddedCells++;
return fAddedCells;
}
- else {
+ else
return 0;
- }
}
+//________________________________________________________________________
+AliVCluster* AliJetModelBaseTask::AddCluster(AliVCluster *oc)
+{
+ // Add a cluster to the event.
+
+ const Int_t nClusters = fOutClusters->GetEntriesFast();
+ AliVCluster *dc = static_cast<AliVCluster*>(fOutClusters->New(nClusters));
+ dc->SetType(AliVCluster::kEMCALClusterv1);
+ dc->SetE(oc->E());
+ Float_t pos[3] = {0};
+ oc->GetPosition(pos);
+ dc->SetPosition(pos);
+ dc->SetNCells(oc->GetNCells());
+ dc->SetCellsAbsId(oc->GetCellsAbsId());
+ dc->SetCellsAmplitudeFraction(oc->GetCellsAmplitudeFraction());
+ dc->SetID(oc->GetID());
+ dc->SetDispersion(oc->GetDispersion());
+ dc->SetEmcCpvDistance(-1);
+ dc->SetChi2(-1);
+ dc->SetTOF(oc->GetTOF()); //time-of-flight
+ dc->SetNExMax(oc->GetNExMax()); //number of local maxima
+ dc->SetM02(oc->GetM02());
+ dc->SetM20(oc->GetM20());
+ dc->SetDistanceToBadChannel(oc->GetDistanceToBadChannel());
+
+ //MC
+ UInt_t nlabels = oc->GetNLabels();
+ Int_t *labels = oc->GetLabels();
+ TArrayI parents;
+ if (nlabels != 0 && labels)
+ parents.Set(nlabels, labels);
+ else {
+ nlabels = 1;
+ parents.Set(1);
+ parents[0] = 0;
+ }
+
+ if (fMarkMC+fMCLabelShift != 0) {
+ for (UInt_t i = 0; i < nlabels; i++) {
+ parents[i] += fMarkMC+fMCLabelShift;
+ }
+ }
+
+ AliESDCaloCluster *esdClus = dynamic_cast<AliESDCaloCluster*>(dc);
+ if (esdClus) {
+ esdClus->AddLabels(parents);
+ }
+ else {
+ AliAODCaloCluster *aodClus = dynamic_cast<AliAODCaloCluster*>(dc);
+ if (aodClus)
+ aodClus->SetLabel(parents.GetArray(), nlabels);
+ }
+
+ return dc;
+}
//________________________________________________________________________
-AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Double_t eta, Double_t phi)
+AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Double_t eta, Double_t phi, Int_t label)
{
// Add a cluster to the event.
e = nPart.E();
}
- return AddCluster(e, absId);
+ return AddCluster(e, absId, label);
}
//________________________________________________________________________
-AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Int_t absId)
+AliVCluster* AliJetModelBaseTask::AddCluster(Double_t e, Int_t absId, Int_t label)
{
// Add a cluster to the event.
cluster->SetCellsAmplitudeFraction(&fract);
cluster->SetID(nClusters);
cluster->SetEmcCpvDistance(-1);
- if (fMarkMC)
- cluster->SetChi2(100); // MC flag!
+ //MC label
+ if (label < 0)
+ label = 0;
+
+ label += fMarkMC+fMCLabelShift;
+
+ if (fEsdMode) {
+ AliESDCaloCluster *esdClus = static_cast<AliESDCaloCluster*>(cluster);
+ TArrayI parents(1, &label);
+ esdClus->AddLabels(parents);
+ }
+ else {
+ AliAODCaloCluster *aodClus = static_cast<AliAODCaloCluster*>(cluster);
+ aodClus->SetLabel(&label, 1);
+ }
+
return cluster;
}
//________________________________________________________________________
-AliPicoTrack* AliJetModelBaseTask::AddTrack(Double_t pt, Double_t eta, Double_t phi, Byte_t type, Double_t etaemc, Double_t phiemc, Bool_t ise)
+AliPicoTrack* AliJetModelBaseTask::AddTrack(Double_t pt, Double_t eta, Double_t phi, Byte_t type, Double_t etaemc, Double_t phiemc, Double_t ptemc, Bool_t ise, Int_t label, Short_t charge, Double_t mass)
{
// Add a track to the event.
+
+ if (pt < 0 && eta < -100 && phi < 0) {
+ GetRandomParticle(pt,eta,phi);
+ }
+ else {
+ if (pt < 0)
+ pt = GetRandomPt();
+ if (eta < -100)
+ eta = GetRandomEta();
+ if (phi < 0)
+ phi = GetRandomPhi(pt);
+ }
- const Int_t nTracks = fOutTracks->GetEntriesFast();
+ if (label >= 0)
+ label += fMarkMC+fMCLabelShift;
+ else if (label < 0)
+ label -= fMarkMC+fMCLabelShift;
- if (pt < 0)
- pt = GetRandomPt();
- if (eta < -100)
- eta = GetRandomEta();
- if (phi < 0)
- phi = GetRandomPhi();
+ if(fAddV2) AddV2(phi, pt);
- Int_t label = fMarkMC ? 100 : -1;
+ const Int_t nTracks = fOutTracks->GetEntriesFast();
AliPicoTrack *track = new ((*fOutTracks)[nTracks]) AliPicoTrack(pt,
eta,
phi,
- 1,
- label, // MC flag!
+ charge,
+ label,
type,
etaemc,
- phiemc,
- ise);
+ phiemc,
+ ptemc,
+ ise,
+ mass);
return track;
}
+//________________________________________________________________________
+AliAODMCParticle* AliJetModelBaseTask::AddMCParticle(AliAODMCParticle *part, Int_t origIndex)
+{
+ const Int_t nPart = fOutMCParticles->GetEntriesFast();
+
+ AliAODMCParticle *aodpart = new ((*fOutMCParticles)[nPart]) AliAODMCParticle(*part);
+
+ if (origIndex + fMCLabelShift >= fOutMCParticlesMap->GetSize())
+ fOutMCParticlesMap->Set((origIndex + fMCLabelShift)*2);
+
+ fOutMCParticlesMap->AddAt(nPart, origIndex + fMCLabelShift);
+ AliDebug(2, Form("Setting bin %d to %d (fMCLabelShift=%d, origIndex=%d)",
+ origIndex + fMCLabelShift, fOutMCParticlesMap->At(origIndex + fMCLabelShift), fMCLabelShift, origIndex));
+
+ return aodpart;
+}
+
+//_____________________________________________________________________________
+void AliJetModelBaseTask::AddV2(Double_t &phi, Double_t &pt) const
+{
+ // similar to AliFlowTrackSimple::AddV2, except for the flow fluctuations
+ Double_t phi0(phi), v2(0.), f(0.), fp(0.), phiprev(0.);
+ if(fDifferentialV2) v2 = fDifferentialV2->Eval(pt);
+ if(TMath::AreEqualAbs(v2, 0, 1e-5)) return;
+ // introduce flow fluctuations (gaussian)
+ if(fFlowFluctuations) v2 += TMath::Sqrt(2*(v2*.25)*(v2*.25))*TMath::ErfInverse(2*(gRandom->Uniform(0, fFlowFluctuations))-1);
+ for (Int_t i(0); i < 100; i++) {
+ phiprev=phi; //store last value for comparison
+ f = phi-phi0+v2*TMath::Sin(2.*(phi-fPsi));
+ fp = 1.0+2.0*v2*TMath::Cos(2.*(phi-fPsi)); //first derivative
+ phi -= f/fp;
+ if (TMath::AreEqualAbs(phiprev, phi, 1e-10)) break;
+ }
+}
+
+//_____________________________________________________________________________
+void AliJetModelBaseTask::CopyCells()
+{
+ if (!fCaloCells)
+ return;
+
+ fAddedCells = 0;
+ fCaloCells->Sort();
+ for (Short_t i = 0; i < fCaloCells->GetNumberOfCells(); i++) {
+ Int_t mclabel = 0;
+ Double_t efrac = 0.;
+ Double_t time = -1;
+ Short_t cellNum = -1;
+ Double_t amp = -1;
+
+ fCaloCells->GetCell(i, cellNum, amp, time, mclabel, efrac);
+
+ if (!fIsMC)
+ mclabel = 0;
+
+ fOutCaloCells->SetCell(i, cellNum, amp, time, mclabel, efrac);
+ fAddedCells++;
+ }
+
+ AliDebug(2, Form("%d cells from the current event", fAddedCells));
+}
+
//________________________________________________________________________
void AliJetModelBaseTask::CopyClusters()
{
// Copy all the clusters in the new collection
+ if (!fClusters)
+ return;
const Int_t nClusters = fClusters->GetEntriesFast();
Int_t nCopiedClusters = 0;
AliESDCaloCluster *esdcluster = static_cast<AliESDCaloCluster*>(fClusters->At(i));
if (!esdcluster || !esdcluster->IsEMCAL())
continue;
- new ((*fOutClusters)[nCopiedClusters]) AliESDCaloCluster(*esdcluster);
+ AliESDCaloCluster *clus = new ((*fOutClusters)[nCopiedClusters]) AliESDCaloCluster(*esdcluster);
+ if (!fIsMC) {
+ TArrayI *labels = clus->GetLabelsArray();
+ if (labels)
+ labels->Reset();
+ }
nCopiedClusters++;
}
}
AliAODCaloCluster *aodcluster = static_cast<AliAODCaloCluster*>(fClusters->At(i));
if (!aodcluster || !aodcluster->IsEMCAL())
continue;
- new ((*fOutClusters)[nCopiedClusters]) AliAODCaloCluster(*aodcluster);
+ AliAODCaloCluster *clus = new ((*fOutClusters)[nCopiedClusters]) AliAODCaloCluster(*aodcluster);
+ if (!fIsMC)
+ clus->SetLabel(0,0);
nCopiedClusters++;
}
}
//________________________________________________________________________
void AliJetModelBaseTask::CopyTracks()
{
+ // Copy all the tracks in the new collection
+
+ if (!fTracks)
+ return;
+
const Int_t nTracks = fTracks->GetEntriesFast();
Int_t nCopiedTracks = 0;
for (Int_t i = 0; i < nTracks; ++i) {
- AliPicoTrack *track = static_cast<AliPicoTrack*>(fTracks->At(i));
- if (!track)
+ AliPicoTrack *
picotrack = static_cast<AliPicoTrack*>(fTracks->At(i));
+ if (!picotrack)
continue;
- new ((*fOutTracks)[nCopiedTracks]) AliPicoTrack(*track);
+ AliPicoTrack *track = new ((*fOutTracks)[nCopiedTracks]) AliPicoTrack(*picotrack);
+ if (!fIsMC && track->GetLabel() != 0)
+ track->SetLabel(0);
nCopiedTracks++;
}
}
+//________________________________________________________________________
+void AliJetModelBaseTask::CopyMCParticles()
+{
+ // Copy all the MC particles in the new collection
+
+ if (!fMCParticles)
+ return;
+
+ const Int_t nPart = fMCParticles->GetEntriesFast();
+ Int_t nCopiedPart = 0;
+ for (Int_t i = 0; i < nPart; ++i) {
+ AliAODMCParticle *part = static_cast<AliAODMCParticle*>(fMCParticles->At(i));
+ if (!part)
+ continue;
+ new ((*fOutMCParticles)[nCopiedPart]) AliAODMCParticle(*part);
+
+ nCopiedPart++;
+ }
+
+ if (!fMCParticlesMap || !fOutMCParticlesMap)
+ return;
+
+ if (fOutMCParticlesMap->GetSize() < fMCParticlesMap->GetSize())
+ fOutMCParticlesMap->Set(fMCParticlesMap->GetSize() * 2);
+
+ for (Int_t i = 0; i < fMCParticlesMap->GetSize(); i++) {
+ fOutMCParticlesMap->AddAt(fMCParticlesMap->At(i), i);
+ if (fMCParticlesMap->At(i) >= 0)
+ fMCLabelShift = i;
+ }
+
+ AliDebug(2,Form("MC particles copied. fMCLabelShift=%d",fMCLabelShift));
+}
+
//________________________________________________________________________
void AliJetModelBaseTask::GetRandomCell(Double_t &eta, Double_t &phi, Int_t &absId)
{
Double_t rndPhi = phi;
do {
if (eta < -100)
- rndEta = GetRandomEta();
+ rndEta = GetRandomEta(kTRUE);
if (phi < 0)
- rndPhi = GetRandomPhi();
+ rndPhi = GetRandomPhi(kTRUE);
fGeom->GetAbsCellIdFromEtaPhi(rndEta, rndPhi, absId);
repeats++;
} while (absId == -1 && repeats < 100);
}
//________________________________________________________________________
-Double_t AliJetModelBaseTask::GetRandomEta()
+Double_t AliJetModelBaseTask::GetRandomEta(Bool_t emcal)
{
// Get random eta.
- return gRandom->Rndm() * (fEtaMax - fEtaMin) + fEtaMin;
+ Double_t etamax = fEtaMax;
+ Double_t etamin = fEtaMin;
+
+ if (emcal) {
+ const Double_t EmcalMinEta = fGeom->GetArm1EtaMin();
+ const Double_t EmcalMaxEta = fGeom->GetArm1EtaMax();
+
+ if (etamax > EmcalMaxEta) etamax = EmcalMaxEta;
+ if (etamax < EmcalMinEta) etamax = EmcalMinEta;
+ if (etamin > EmcalMaxEta) etamin = EmcalMaxEta;
+ if (etamin < EmcalMinEta) etamin = EmcalMinEta;
+ }
+
+ return gRandom->Rndm() * (etamax - etamin) + etamin;
}
//________________________________________________________________________
-Double_t AliJetModelBaseTask::GetRandomPhi()
+Double_t AliJetModelBaseTask::GetRandomPhi(Bool_t emcal)
{
// Get random phi.
+
+ Double_t phimax = fPhiMax;
+ Double_t phimin = fPhiMin;
+
+ if (emcal) {
+ const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
+ const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
+
+ if (phimax > EmcalMaxPhi) phimax = EmcalMaxPhi;
+ if (phimax < EmcalMinPhi) phimax = EmcalMinPhi;
+ if (phimin > EmcalMaxPhi) phimin = EmcalMaxPhi;
+ if (phimin < EmcalMinPhi) phimin = EmcalMinPhi;
+ }
+
+ Double_t result = gRandom->Rndm() * (phimax - phimin) + phimin;
- return gRandom->Rndm() * (fPhiMax - fPhiMin) + fPhiMin;
+ return result;
}
//________________________________________________________________________
return gRandom->Rndm() * (fPtMax - fPtMin) + fPtMin;
}
+//________________________________________________________________________
+void AliJetModelBaseTask::GetRandomParticle(Double_t &pt, Double_t &eta, Double_t &phi, Bool_t emcal)
+{
+ // Get a random particle.
+
+ eta = GetRandomEta(emcal);
+
+ if (fPtPhiEvPlDistribution) {
+ Double_t phimax = fPhiMax;
+ Double_t phimin = fPhiMin;
+
+ if (emcal) {
+ const Double_t EmcalMinPhi = fGeom->GetArm1PhiMin() * TMath::DegToRad();
+ const Double_t EmcalMaxPhi = fGeom->GetArm1PhiMax() * TMath::DegToRad();
+
+ if (phimax > EmcalMaxPhi) phimax = EmcalMaxPhi;
+ if (phimax < EmcalMinPhi) phimax = EmcalMinPhi;
+ if (phimin > EmcalMaxPhi) phimin = EmcalMaxPhi;
+ if (phimin < EmcalMinPhi) phimin = EmcalMinPhi;
+ }
+
+ if (fPtPhiEvPlDistribution->GetXmin() > phimax || fPtPhiEvPlDistribution->GetXmax() < phimin) {
+ AliWarning(Form("The hisogram %s does not overlap with the EMCal acceptance limits. It will be ignored.",fPtPhiEvPlDistribution->GetName()));
+ pt = GetRandomPt();
+ phi = GetRandomPhi(emcal);
+ }
+ else {
+ do {
+ fPtPhiEvPlDistribution->GetRandom2(pt,phi);
+ phi += fPsi;
+ if (phi > TMath::Pi() * 2) phi -= TMath::Pi() * 2;
+ } while (phi > phimax || phi < phimin);
+ }
+ }
+ else {
+ pt = GetRandomPt();
+ phi = GetRandomPhi(emcal);
+ }
+}
+
//________________________________________________________________________
void AliJetModelBaseTask::Run()
{
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff