cl->fIdMax = id;
cl->fEmax = emax;
cl->fTmax = cells->GetCellTime(id);
+ cl->fE2max = GetSecondMaxCell(clus);
if (clus->GetDistanceToBadChannel()<10000)
cl->fDbc = clus->GetDistanceToBadChannel();
if (!TMath::IsNaN(clus->GetDispersion()))
GetSigma(clus,maxAxis,minAxis);
clus->SetTOF(maxAxis); // store sigma in TOF
cl->fSig = maxAxis;
+ Double_t sEtaEta = 0;
+ Double_t sPhiPhi = 0;
+ GetSigmaEtaEta(clus, sEtaEta, sPhiPhi);
+ cl->fSigEtaEta = sEtaEta;
+ cl->fSigPhiPhi = sPhiPhi;
Double_t clusterEcc = 0;
if (maxAxis > 0)
clusterEcc = TMath::Sqrt(1.0 - minAxis*minAxis/(maxAxis*maxAxis));
cl->fTrIso = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist);
cl->fTrIso1 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist, 1);
cl->fTrIso2 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist, 2);
+ cl->fTrIsoD1 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist-0.1);
+ cl->fTrIso1D1 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist-0.1, 1);
+ cl->fTrIso2D1 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist-0.1, 2);
+ cl->fTrIsoD3 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist+0.1);
+ cl->fTrIso1D3 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist+0.1, 1);
+ cl->fTrIso2D3 = GetTrackIsolation(clusterVec.Eta(),clusterVec.Phi(),fIsoDist+0.1, 2);
+ cl->fTrIsoStrip = GetTrackIsoStrip(clusterVec.Eta(), clusterVec.Phi());
cl->fCeCore = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),0.05);
cl->fCeIso = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),fIsoDist);
+ cl->fCeIso1 = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),0.10);
+ cl->fCeIso3 = GetCellIsolation(clsVec.Eta(),clsVec.Phi(),0.30);
+ cl->fCeIso4x4 = GetCellIsoNxM(clsVec.Eta(),clsVec.Phi(), 4, 4);
+ cl->fCeIso5x5 = GetCellIsoNxM(clsVec.Eta(),clsVec.Phi(), 5, 5);
+ cl->fCeIso3x22 = GetCellIsoNxM(clsVec.Eta(),clsVec.Phi(), 3, 22);
if (fAmpInTrigger) { // fill trigger info if present
Double_t trigpen = 0;
}
Float_t tmpR=-1, tmpZ=-1;
+ Double_t dedx = 0;
if (!fDoTrMatGeom) {
AliExternalTrackParam *tParam = 0;
if (fEsdEv) {
AliESDtrack *esdTrack = static_cast<AliESDtrack*>(track);
tParam = new AliExternalTrackParam(*esdTrack->GetTPCInnerParam());
+ dedx = esdTrack->GetTPCsignal();
} else
tParam = new AliExternalTrackParam(track);
cl->fTrDz = tmpZ;
cl->fTrDr = TMath::Sqrt(tmpR*tmpR-tmpZ*tmpZ);
cl->fTrEp = clus->E()/track->P();
+ cl->fTrDedx = dedx;
cl->fIsTrackM = 1;
}
}
}
return cellIsolation;
}
+//________________________________________________________________________
+Double_t AliAnalysisTaskEMCALPi0PbPb::GetCellIsoNxM(Double_t cEta, Double_t cPhi, Int_t N, Int_t M) const
+{
+ // Compute isolation based on cell content, in a NxM rectangle.
+
+ AliVCaloCells *cells = fEsdCells;
+ if (!cells)
+ cells = fAodCells;
+ if (!cells)
+ return 0;
+
+ Double_t cellIsolation = 0;
+ Int_t ncells = cells->GetNumberOfCells();
+ for (Int_t i = 0; i<ncells; ++i) {
+ Int_t absID = TMath::Abs(cells->GetCellNumber(i));
+ Float_t eta=-1, phi=-1;
+ fGeom->EtaPhiFromIndex(absID,eta,phi);
+ Double_t phidiff = TVector2::Phi_mpi_pi(phi-cPhi);
+ Double_t etadiff = (eta-cEta)*(eta-cEta)+phidiff*phidiff;
+ if(TMath::Abs(etadiff)/0.014>N)
+ continue;
+ if(TMath::Abs(phidiff)/0.014>M)
+ continue;
+ Double_t cellE = cells->GetAmplitude(i);
+ cellIsolation += cellE;
+ }
+ return cellIsolation;
+}
+
//________________________________________________________________________
Double_t AliAnalysisTaskEMCALPi0PbPb::GetCellEnergy(const AliVCluster *cluster) const
sigmaMin = TMath::Abs(Sxx + Syy - TMath::Sqrt(TMath::Abs((Sxx-Syy)*(Sxx-Syy)+4.0*Sxy*Sxy)))/2.0;
sigmaMin = TMath::Sqrt(TMath::Abs(sigmaMin));
}
+//________________________________________________________________________
+void AliAnalysisTaskEMCALPi0PbPb::GetSigmaEtaEta(const AliVCluster *c, Double_t& sEtaEta, Double_t &sPhiPhi) const
+{
+ // Calculate the (E) weighted variance along the pseudorapidity.
+
+ sEtaEta = 0;
+ sPhiPhi = 0;
+
+ Double_t Ec = c->E(); // cluster energy
+ if(Ec<=0)
+ return;
+
+ const Int_t ncells = c->GetNCells();
+
+ Double_t EtaC = 0; // cluster first moment along eta
+ Double_t PhiC = 0; // cluster first moment along phi
+ Double_t Setaeta = 0; // cluster second central moment along eta
+ Double_t Sphiphi = 0; // cluster second central moment along phi
+ Double_t w[ncells]; // weight max(0,4.5*log(E_i/Ec))
+ Double_t sumw = 0;
+ Int_t id[ncells];
+
+ AliVCaloCells *cells = fEsdCells;
+ if (!cells)
+ cells = fAodCells;
+
+ if (!cells)
+ return;
+
+ if (ncells==1)
+ return;
+
+ for(Int_t j=0; j<ncells; ++j) {
+ id[j] = TMath::Abs(c->GetCellAbsId(j));
+ Double_t cellen = cells->GetCellAmplitude(id[j]);
+ w[j] = TMath::Max(0., 4.5+TMath::Log(cellen/Ec));
+ TVector3 pos;
+ fGeom->GetGlobal(id[j],pos);
+ EtaC += w[j]*pos.Eta();
+ PhiC += w[j]*pos.Phi();
+ sumw += w[j];
+ }
+ EtaC /= sumw;
+ PhiC /= sumw;
+
+ for(Int_t j=0; j<ncells; ++j) {
+ TVector3 pos;
+ fGeom->GetGlobal(id[j],pos);
+ Setaeta = w[j]*(pos.Eta() - EtaC)*(pos.Eta() - EtaC);
+ Sphiphi = w[j]*(pos.Phi() - PhiC)*(pos.Phi() - PhiC);
+ }
+ Setaeta /= sumw;
+ sEtaEta = TMath::Sqrt(Setaeta);
+ Sphiphi /= sumw;
+ sPhiPhi = TMath::Sqrt(Sphiphi);
+}
+
//________________________________________________________________________
Int_t AliAnalysisTaskEMCALPi0PbPb::GetNCells(const AliVCluster *c, Double_t emin) const
}
return trkIsolation;
}
+//________________________________________________________________________
+Double_t AliAnalysisTaskEMCALPi0PbPb::GetTrackIsoStrip(Double_t cEta, Double_t cPhi, Double_t dEta, Double_t dPhi, Double_t pt) const
+{
+ // Compute isolation based on tracks.
+
+ Double_t trkIsolation = 0;
+ Int_t ntrks = fSelPrimTracks->GetEntries();
+ for(Int_t j = 0; j<ntrks; ++j) {
+ AliVTrack *track = static_cast<AliVTrack*>(fSelTracks->At(j));
+ if (!track)
+ continue;
+ if (track->Pt()<pt)
+ continue;
+ Float_t eta = track->Eta();
+ Float_t phi = track->Phi();
+ Double_t phidiff = TVector2::Phi_mpi_pi(phi-cPhi);
+ Double_t etadiff = (eta-cEta);
+ if(TMath::Abs(etadiff)>dEta || TMath::Abs(phidiff)>dPhi)
+ continue;
+ trkIsolation += track->Pt();
+ }
+ return trkIsolation;
+}
+
//________________________________________________________________________
Bool_t AliAnalysisTaskEMCALPi0PbPb::IsShared(const AliVCluster *c) const
ProcessDaughters(dmc,i,arr);
}
}
+
+//________________________________________________________________________
+Double_t AliAnalysisTaskEMCALPi0PbPb::GetSecondMaxCell(AliVCluster *clus)
+{
+ // Get second maximum cell.
+
+ AliVCaloCells *cells = fEsdCells;
+ if (!cells)
+ cells = fAodCells;
+ if (!cells)
+ return -1;
+
+ Double_t secondEmax=0, firstEmax=0;
+ Double_t cellen;
+ for(Int_t iCell=0;iCell<clus->GetNCells();iCell++){
+ Int_t absId = clus->GetCellAbsId(iCell);
+ cellen = cells->GetCellAmplitude(absId);
+ if(cellen > firstEmax)
+ firstEmax = cellen;
+ }
+ for(Int_t iCell=0;iCell<clus->GetNCells();iCell++){
+ Int_t absId = clus->GetCellAbsId(iCell);
+ cellen = cells->GetCellAmplitude(absId);
+ if(cellen < firstEmax && cellen > secondEmax)
+ secondEmax = cellen;
+ }
+ return secondEmax;
+}
void FillVertex(AliStaVertex *v, const AliESDVertex *esdv);
void FillVertex(AliStaVertex *v, const AliAODVertex *aodv);
Double_t GetCellIsolation(Double_t cEta, Double_t cPhi, Double_t radius=0.2) const;
+ Double_t GetCellIsoNxM(Double_t cEta, Double_t cPhi, Int_t N, Int_t M) const;
Double_t GetCellEnergy(const AliVCluster *c) const;
Double_t GetMaxCellEnergy(const AliVCluster *c) const { Short_t id=-1; return GetMaxCellEnergy(c,id); }
Double_t GetMaxCellEnergy(const AliVCluster *c, Short_t &id) const;
Int_t GetNCells(const AliVCluster *c, Double_t emin=0.) const;
void GetSigma(const AliVCluster *c, Double_t &sigmaMax, Double_t &sigmaMin) const;
+ void GetSigmaEtaEta(const AliVCluster *c, Double_t &sigmaEtaEta, Double_t &sigmaPhiPhi) const;
Double_t GetTrackIsolation(Double_t cEta, Double_t cPhi, Double_t radius=0.2, Double_t pt=0.) const;
+ Double_t GetTrackIsoStrip(Double_t cEta, Double_t cPhi, Double_t dEta=0.015, Double_t dPhi=0.3, Double_t pt=0.) const;
Double_t GetTrigEnergy(const AliVCluster *c) const;
Bool_t IsShared(const AliVCluster *c) const;
void PrintDaughters(const AliVParticle *p, const TObjArray *arr, Int_t level=0) const;
void PrintTrackRefs(AliMCParticle *p) const;
void ProcessDaughters(AliVParticle *p, Int_t index, const TObjArray *arr);
void ProcessDaughters(AliMCParticle *p, Int_t index, const AliMCEvent *arr);
+ Double_t GetSecondMaxCell(AliVCluster *clus);
// input members
TString fCentVar; // variable for centrality determination
{
public:
AliStaCluster() : TObject(), fE(0), fR(0), fEta(0), fPhi(0), fN(0), fN1(0), fN3(0), fIdMax(0), fEmax(0), fTmax(0),
- fDbc(-1), fDisp(-1), fM20(0), fM02(0), fEcc(0), fSig(0), fIsTrackM(0), fTrDz(0), fTrDr(-1),
- fTrEp(0), fTrIso(0), fTrIso1(0), fTrIso2(0), fCeIso(0), fCeCore(0), fIsTrigM(0), fTrigE(-1),
- fTrigMaskE(-1), fIsShared(0), fMcLabel(-1), fEmbE(0) {;}
+ fE2max(0),fDbc(-1), fDisp(-1), fM20(0), fM02(0), fEcc(0), fSig(0), fSigEtaEta(0), fSigPhiPhi(0),
+ fIsTrackM(0), fTrDz(0), fTrDr(-1), fTrEp(0), fTrDedx(0), fTrIso(0), fTrIso1(0), fTrIso2(0),
+ fTrIsoD1(0), fTrIso1D1(0), fTrIso2D1(0), fTrIsoD3(0), fTrIso1D3(0), fTrIso2D3(0),fTrIsoStrip(0),
+ fCeIso(0), fCeIso1(0), fCeIso3(0), fCeIso4x4(0), fCeIso5x5(0), fCeCore(0), fCeIso3x22(0),
+ fIsTrigM(0), fTrigE(-1), fTrigMaskE(-1), fIsShared(0), fMcLabel(-1), fEmbE(0) {;}
public:
Double32_t fE; //[0,0,16] energy
UShort_t fIdMax; // id maximum cell
Double32_t fEmax; //[0,0,16] energy of maximum cell
Double32_t fTmax; //[0,0,16] time of maximum cell
+ Double32_t fE2max; //[0,0,16] energy of second maximum cell
Double32_t fDbc; //[0,0,16] distance to nearest bad channel
Double32_t fDisp; //[0,0,16] cluster dispersion, for shape analysis
Double32_t fM20; //[0,0,16] 2-nd moment along the main eigen axis
Double32_t fM02; //[0,0,16] 2-nd moment along the second eigen axis
Double32_t fEcc; //[0,0,16] eccentricity
Double32_t fSig; //[0,0,16] sigma
+ Double32_t fSigEtaEta; //[0,0,16] sigma eta-eta
+ Double32_t fSigPhiPhi; //[0,0,16] sigma phi-phi
Bool_t fIsTrackM; // if true then track values are set
Double32_t fTrDz; //[0,0,16] dZ to nearest track
Double32_t fTrDr; //[0,0,16] dR to nearest track (in x,y)
Double32_t fTrEp; //[0,0,16] E/P to nearest track
+ Double32_t fTrDedx; //[0,0,16] dE/dx (TPC signal) to nearest track
Double32_t fTrIso; //[0,0,16] track isolation
Double32_t fTrIso1; //[0,0,16] track isolation (pt>1GeV/c)
Double32_t fTrIso2; //[0,0,16] track isolation (pt>2GeV/c)
- Double32_t fCeIso; //[0,0,16] cell isolation
- Double32_t fCeCore; //[0,0,16] cell content in R=0.025
+ Double32_t fTrIsoD1; //[0,0,16] track isolation, iso dist 0.25
+ Double32_t fTrIso1D1; //[0,0,16] track isolation (pt>1GeV/c), iso dist 0.1
+ Double32_t fTrIso2D1; //[0,0,16] track isolation (pt>2GeV/c), iso dist 0.1
+ Double32_t fTrIsoD3; //[0,0,16] track isolation, iso dist 0.3
+ Double32_t fTrIso1D3; //[0,0,16] track isolation (pt>1GeV/c), iso dist 0.3
+ Double32_t fTrIso2D3; //[0,0,16] track isolation (pt>2GeV/c), iso dist 0.3
+ Double32_t fTrIsoStrip; //[0,0,16] track isolation strip, dEtaXdPhi=0.015x0.3
+ Double32_t fCeIso; //[0,0,16] cell isolation in R=0.20
+ Double32_t fCeIso1; //[0,0,16] cell isolation in R=0.10
+ Double32_t fCeIso3; //[0,0,16] cell isolation in R=0.30
+ Double32_t fCeIso4x4; //[0,0,16] cell isolation in 4x4 cells
+ Double32_t fCeIso5x5; //[0,0,16] cell isolation in 5x5 cells
+ Double32_t fCeCore; //[0,0,16] cell content in R=0.05
+ Double32_t fCeIso3x22; //[0,0,16] cell isolation in rectangular strip of dEtaXdPhi=0.042x0.308
Bool_t fIsTrigM; // if true then trigger values are set
Double32_t fTrigE; //[0,0,16] trigger tower energy
Double32_t fTrigMaskE; //[0,0,16] masked trigger tower energy
{
public:
AliStaPart() : TObject(), fPt(0), fEta(0), fPhi(0), fVR(0), fVEta(0), fVPhi(0), fPid(0), fMo(-1), fDet(-2),
- fLab(-1), fNs(0) { memset(fDs,-1,sizeof(Short_t)*9); }
+ fLab(-1), fNs(0) { memset(fDs,-1,sizeof(Short_t)*99); }
Int_t OnEmcal() const { return (fDet==8); }
Int_t IsSim() const { return (fDet!=-2); }
// the following must be filled before first usage
Short_t fLab; //! label (index in array)
Short_t fNs; //! number of daughters
- Short_t fDs[9]; //! daughters
+ Short_t fDs[99]; //! daughters
ClassDef(AliStaPart,1) // Particle class
};