AliPhysics  8bb951a (8bb951a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliEmcalJet.h
Go to the documentation of this file.
1 #ifndef AliEmcalJet_H
2 #define AliEmcalJet_H
3 
4 #include <Riosfwd.h>
5 #include <vector>
6 #include <algorithm>
7 #include <utility>
8 
9 #include <TArrayS.h>
10 #include <TMath.h>
11 #include <TClonesArray.h>
12 #include <TVector2.h>
13 #include <TLorentzVector.h>
14 #include <TString.h>
15 
16 #include "AliVParticle.h"
17 #include "AliVCluster.h"
18 #include "AliVEvent.h"
19 
20 class AliEmcalJet : public AliVParticle
21 {
22  public:
24  kDStar = 1<<0,
25  kD0 = 1<<1,
26  kSig1 = 1<<2,
27  kSig2 = 1<<3,
28  kBckgrd1 = 1<<4,
29  kBckgrd2 = 1<<5,
30  kBckgrd3 = 1<<6
31  //.....
32  };
33 
34  AliEmcalJet();
35  AliEmcalJet(Double_t px, Double_t py, Double_t pz);
36  AliEmcalJet(Double_t pt, Double_t eta, Double_t phi, Double_t m);
37  AliEmcalJet(const AliEmcalJet &jet);
38  AliEmcalJet& operator=(const AliEmcalJet &jet);
39  friend std::ostream &operator<<(std::ostream &in, const AliEmcalJet &jet);
40 
41  std::ostream &Print(std::ostream &in) const;
42  TString toString() const;
43 
44  Double_t Px() const { return fPt*TMath::Cos(fPhi); }
45  Double_t Py() const { return fPt*TMath::Sin(fPhi); }
46  Double_t Pz() const { return fPt*TMath::SinH(fEta); }
47  Double_t Pt() const { return fPt; }
48  Double_t P() const { return fPt*TMath::CosH(fEta); }
49  Bool_t PxPyPz(Double_t p[3]) const { p[0]=Px();p[1]=Py();p[2]=Pz(); return 1; }
50  Double_t Xv() const { return 0.; }
51  Double_t Yv() const { return 0.; }
52  Double_t Zv() const { return 0.; }
53  Bool_t XvYvZv(Double_t x[3]) const { x[0]=0;x[1]=0;x[2]=0; return 1; }
54  Double_t OneOverPt() const { return 1./fPt; }
55  Double_t Phi() const { return fPhi; }
56  Double_t Phi_0_2pi() const { return TVector2::Phi_0_2pi(fPhi); }
57  Double_t Theta() const { return 2*TMath::ATan(TMath::Exp(-fEta)); }
58  Double_t E() const { Double_t p=P(); return TMath::Sqrt(fM*fM+p*p); }
59  Double_t M() const { return fM; }
60  Double_t Eta() const { return fEta; }
61  Double_t Y() const { Double_t e = E(); Double_t pz = Pz(); return 0.5*TMath::Log((e+pz)/(e-pz)); }
62  Short_t Charge() const { return 0; }
63  Int_t GetLabel() const { return fLabel; }
64  Int_t PdgCode() const { return 0; }
65  const Double_t *PID() const { return 0; }
66  void GetMom(TLorentzVector &vec) const {GetMomentum(vec);}
67  void GetMomentum(TLorentzVector &vec) const;
68 
69  Double_t Area() const { return fArea; }
70  Double_t AreaPt() const { return fArea; }
71  Double_t AreaEta() const { return fAreaEta; }
72  Double_t AreaPhi() const { return fAreaPhi; }
73  Double_t AreaE() const { return fAreaE; }
74  Double_t AreaEmc() const { return fAreaEmc; }
75  Bool_t AxisInEmcal() const { return fAxisInEmcal; }
76  Int_t Compare(const TObject* obj) const;
77  Short_t ClusterAt(Int_t idx) const { return fClusterIDs.At(idx); }
78  AliVCluster *ClusterAt(Int_t idx, TClonesArray *ca) const { if (!ca) return 0; return dynamic_cast<AliVCluster*>(ca->At(ClusterAt(idx))); }
79  Int_t ContainsCluster(AliVCluster* cluster, TClonesArray* clusters) const { return clusters==NULL||cluster==NULL ? 0 : ContainsCluster(clusters->IndexOf(cluster)); }
80  Int_t ContainsCluster(Int_t ic) const;
81  AliVCluster *GetLeadingCluster(TClonesArray *clusters) const;
82  UShort_t GetNumberOfClusters() const { return fClusterIDs.GetSize(); }
83  UShort_t GetNumberOfTracks() const { return fTrackIDs.GetSize(); }
85  Double_t FracEmcalArea() const { return fAreaEmc/fArea; }
86  Bool_t IsInsideEmcal() const { return (fAreaEmc/fArea>0.999); }
87  Bool_t IsInEmcal() const { return (fAreaEmc>0); }
88  Bool_t IsMC() const { return (Bool_t)(MCPt() > 0); }
89  Bool_t IsSortable() const { return kTRUE; }
90  Double_t MaxNeutralPt() const { return fMaxNPt; }
91  Double_t MaxChargedPt() const { return fMaxCPt; }
92  Double_t NEF() const { return fNEF; }
93  UShort_t Nn() const { return fNn; }
94  UShort_t Nch() const { return fNch; }
95  UShort_t N() const { return Nch()+Nn(); }
96  Int_t NEmc() const { return fNEmc; }
97  Double_t MCPt() const { return fMCPt; }
98  Double_t MaxClusterPt() const { return MaxNeutralPt(); }
99  Double_t MaxTrackPt() const { return MaxChargedPt(); }
100  Double_t MaxPartPt() const { return fMaxCPt < fMaxNPt ? fMaxNPt : fMaxCPt; }
101  Double_t PtEmc() const { return fPtEmc; }
102  Double_t PtSub() const { return fPtSub; }
103  Double_t PtSubVect() const { return fPtSubVect; }
104  Double_t PtSub(Double_t rho, Bool_t save = kFALSE);
105  Double_t PtSubVect(Double_t rho, Bool_t save = kFALSE);
106  TLorentzVector SubtractRhoVect(Double_t rho, Bool_t save = kFALSE);
107  Short_t TrackAt(Int_t idx) const { return fTrackIDs.At(idx); }
108  AliVParticle *TrackAt(Int_t idx, TClonesArray *ta) const { if (!ta) return 0; return dynamic_cast<AliVParticle*>(ta->At(TrackAt(idx))); }
109  Int_t ContainsTrack(AliVParticle* track, TClonesArray* tracks) const { return tracks==NULL||track==NULL ? 0 : ContainsTrack(tracks->IndexOf(track)); }
110  Int_t ContainsTrack(Int_t it) const;
111  AliVParticle *GetLeadingTrack(TClonesArray *tracks) const;
112 
113 
114  void AddClusterAt(Int_t clus, Int_t idx){ fClusterIDs.AddAt(clus, idx); }
115  void AddTrackAt(Int_t track, Int_t idx) { fTrackIDs.AddAt(track, idx); }
116  void Clear(Option_t */*option*/="") { fClusterIDs.Set(0); fTrackIDs.Set(0); fClosestJets[0] = 0; fClosestJets[1] = 0;
117  fClosestJetsDist[0] = 0; fClosestJetsDist[1] = 0; fMatched = 0; fPtSub = 0;
118  fGhosts.clear(); fHasGhost = kFALSE; }
119  Double_t DeltaR(const AliVParticle* part) const;
120  Double_t GetZ ( const Double_t trkPx, const Double_t trkPy, const Double_t trkPz ) const; // Get Z of constituent trk
121  Double_t GetZ ( const AliVParticle* trk ) const; // Get Z of constituent trk
122  Double_t GetXi ( const AliVParticle* trk ) const { return TMath::Log ( 1/GetZ (trk) ); } // Get Xi of constituent trk
123  Double_t GetXi ( const Double_t trkPx, const Double_t trkPy, const Double_t trkPz ) const { return TMath::Log ( 1/GetZ (trkPx, trkPy, trkPz ) ); } // Get Xi of constituent trk
124 
125  void SetLabel(Int_t l) { fLabel = l; }
126  void SetArea(Double_t a) { fArea = a; }
127  void SetAreaEta(Double_t a) { fAreaEta = a; }
128  void SetAreaPhi(Double_t a) { fAreaPhi = TVector2::Phi_0_2pi(a); }
129  void SetAreaE(Double_t a) { fAreaE = a; }
130  void SetAreaEmc(Double_t a) { fAreaEmc = a; }
131  void SetAxisInEmcal(Bool_t b) { fAxisInEmcal = b; }
132  void SetMaxNeutralPt(Double32_t t) { fMaxNPt = t; }
133  void SetMaxChargedPt(Double32_t t) { fMaxCPt = t; }
134  void SetNEF(Double_t nef) { fNEF = nef; }
135  void SetNumberOfClusters(Int_t n) { fClusterIDs.Set(n); }
136  void SetNumberOfTracks(Int_t n) { fTrackIDs.Set(n); }
137  void SetNumberOfCharged(Int_t n) { fNch = n; }
138  void SetNumberOfNeutrals(Int_t n) { fNn = n; }
139  void SetMCPt(Double_t p) { fMCPt = p; }
140  void SortConstituents();
141  std::vector<int> SortConstituentsPt(TClonesArray *tracks) const;
142  void SetNEmc(Int_t n) { fNEmc = n; }
143  void SetPtEmc(Double_t pt) { fPtEmc = pt; }
144  void SetPtSub(Double_t ps) { fPtSub = ps; }
145  void SetPtSubVect(Double_t ps) { fPtSubVect = ps; }
146 
147  // Trigger
148  Bool_t IsTriggerJet(UInt_t trigger=AliVEvent::kEMCEJE) const { return (Bool_t)((fTriggers & trigger) != 0); }
149  void SetTrigger(UInt_t trigger) { fTriggers = trigger; }
150  void AddTrigger(UInt_t trigger) { fTriggers |= trigger; }
151 
152  // Matching
153  void SetClosestJet(AliEmcalJet *j, Double_t d) { fClosestJets[0] = j; fClosestJetsDist[0] = d ; }
154  void SetSecondClosestJet(AliEmcalJet *j, Double_t d) { fClosestJets[1] = j; fClosestJetsDist[1] = d ; }
155  void SetMatchedToClosest(UShort_t m) { fMatched = 0; fMatchingType = m ; }
156  void SetMatchedToSecondClosest(UShort_t m) { fMatched = 1; fMatchingType = m ; }
157  void ResetMatching();
158  AliEmcalJet* ClosestJet() const { return fClosestJets[0] ; }
159  Double_t ClosestJetDistance() const { return fClosestJetsDist[0] ; }
160  AliEmcalJet* SecondClosestJet() const { return fClosestJets[1] ; }
161  Double_t SecondClosestJetDistance() const { return fClosestJetsDist[1] ; }
162  AliEmcalJet* MatchedJet() const { return fMatched < 2 ? fClosestJets[fMatched] : 0; }
163  UShort_t GetMatchingType() const { return fMatchingType ; }
164 
166  void SetTagStatus(Int_t i) { fTagStatus = i ; }
167  AliEmcalJet* GetTaggedJet() const { return fTaggedJet ; }
168  Int_t GetTagStatus() const { return fTagStatus ; }
169 
170  //jet shape derivatives
171  //jet mass
172  void SetFirstDerivative(Double_t d) { fJetShapeMassFirstDer = d ; }
173  void SetSecondDerivative(Double_t d) { fJetShapeMassSecondDer = d ; }
176  Double_t GetFirstDerivative() const { return fJetShapeMassFirstDer ; }
177  Double_t GetSecondDerivative() const { return fJetShapeMassSecondDer ; }
178  Double_t GetFirstOrderSubtracted() const { return fJetShapeMassFirstSub ; }
179  Double_t GetSecondOrderSubtracted() const { return fJetShapeMassSecondSub ; }
180 
181  //jet structure function
182  TArrayF GetGRNumerator() const { return fGRNumerator ; }
183  TArrayF GetGRDenominator() const { return fGRDenominator ; }
184  TArrayF GetGRNumeratorSub() const { return fGRNumeratorSub ; }
185  TArrayF GetGRDenominatorSub() const { return fGRDenominatorSub ; }
186  void AddGRNumAt(Float_t num, Int_t idx) { fGRNumerator.AddAt(num, idx) ; }
187  void AddGRDenAt(Float_t den, Int_t idx) { fGRDenominator.AddAt(den, idx) ; }
188  void SetGRNumSize(UInt_t s) { fGRNumerator.Set(s) ; }
189  void SetGRDenSize(UInt_t s) { fGRDenominator.Set(s) ; }
190 
191  void AddGRNumSubAt(Float_t num, Int_t idx) { fGRNumeratorSub.AddAt(num, idx) ; }
192  void AddGRDenSubAt(Float_t den, Int_t idx) { fGRDenominatorSub.AddAt(den, idx) ; }
193  void SetGRNumSubSize(UInt_t s) { fGRNumeratorSub.Set(s) ; }
194  void SetGRDenSubSize(UInt_t s) { fGRDenominatorSub.Set(s) ; }
195  void PrintGR();
196 
197  //Angularity
206 
207  //pTD
208  void SetFirstDerivativepTD(Double_t d) { fJetShapepTDFirstDer = d ; }
212  Double_t GetFirstDerivativepTD() const { return fJetShapepTDFirstDer ; }
213  Double_t GetSecondDerivativepTD() const { return fJetShapepTDSecondDer ; }
214  Double_t GetFirstOrderSubtractedpTD() const { return fJetShapepTDFirstSub ; }
216 
217  //Circularity
226 
227  //Sigma2
232  Double_t GetFirstDerivativeSigma2() const { return fJetShapeSigma2FirstDer ; }
236 
237 
238  //number of contituents
247 
248  //leading minus subleading constituent
253  Double_t GetFirstDerivativeLeSub() const { return fJetShapeLeSubFirstDer ; }
254  Double_t GetSecondDerivativeLeSub() const { return fJetShapeLeSubSecondDer ; }
257 
258  //heavy-flavor jets
259  Int_t GetFlavour() const { return fFlavourTagging; }
260  void AddFlavourTag(Int_t tag) { fFlavourTagging |= tag; }
261  void SetFlavour(Int_t flavour) { fFlavourTagging = flavour; }
262  Bool_t TestFlavourTag(Int_t tag) const { return (Bool_t)((tag & fFlavourTagging) !=0); }
263  void AddFlavourTrack(AliVParticle* hftrack){ if (!fFlavourTracks) fFlavourTracks = new TObjArray(); fFlavourTracks->Add(hftrack); }
264  AliVParticle *GetFlavourTrack(Int_t i=0) const { return fFlavourTracks && i >= 0 && i < fFlavourTracks->GetEntriesFast() ? static_cast<AliVParticle*>(fFlavourTracks->At(i)) : 0; }
265  Double_t GetFlavourTrackZ(Int_t i=0) const;
266  AliVParticle *RemoveFlavourTrack(Int_t i=0) { return fFlavourTracks && i >= 0 && i < fFlavourTracks->GetEntriesFast() ? static_cast<AliVParticle*>(fFlavourTracks->RemoveAt(i)) : 0; }
268 
269  void AddGhost(const Double_t dPx, const Double_t dPy, const Double_t dPz, const Double_t dE) {
270  TLorentzVector ghost(dPx, dPy, dPz, dE);
271  fGhosts.push_back(ghost);
272  if (!fHasGhost) fHasGhost = kTRUE;
273  return;
274  }
275 
276  Bool_t HasGhost() const { return fHasGhost; }
277  const std::vector<TLorentzVector> GetGhosts() const { return fGhosts; }
278 
279  void Print(Option_t* /*opt*/ = "") const;
280  void PrintConstituents(TClonesArray* tracks, TClonesArray* clusters) const;
281 
282  protected:
283  Double32_t fPt; //[0,0,12] pt
284  Double32_t fEta; //[-1,1,12] eta
285  Double32_t fPhi; //[0,6.3,12] phi
286  Double32_t fM; //[0,0,8] mass
287  Double32_t fNEF; //[0,1,8] neutral energy fraction
288  Double32_t fArea; //[0,0,12] area (transverse)
289  Double32_t fAreaEta; //[0,0,12] area eta
290  Double32_t fAreaPhi; //[0,0,12] area phi
291  Double32_t fAreaE; //[0,0,12] temporal area component
292  Double32_t fAreaEmc; //[0,0,12] area on EMCAL surface (determined from ghosts)
293  Bool_t fAxisInEmcal; // =true if jet axis inside EMCAL acceptance
294  Int_t fFlavourTagging; // tag jet with a flavour, bit 0 = no tag; bit 1= Dstar; bit 2 = D0
295  TObjArray *fFlavourTracks;
296  Double32_t fMaxCPt; //[0,0,12] pt of maximum charged constituent
297  Double32_t fMaxNPt; //[0,0,12] pt of maximum neutral constituent
298  Double32_t fMCPt; // pt from MC particles contributing to the jet
299  Int_t fNn; // number of neutral constituents
300  Int_t fNch; // number of charged constituents
301  Double32_t fPtEmc; //[0,0,12] pt in EMCAL acceptance
302  Int_t fNEmc; // number of constituents in EMCAL acceptance
303  TArrayI fClusterIDs; // array containing ids of cluster constituents
304  TArrayI fTrackIDs; // array containing ids of track constituents
306  Double32_t fClosestJetsDist[2];
307  UShort_t fMatched;
308  UShort_t fMatchingType;
310  Int_t fTagStatus;
311  Double_t fPtSub;
312  Double_t fPtSubVect;
313  UInt_t fTriggers;
314 
319  Int_t fLabel; // label to inclusive jet for constituent subtracted jet
320 
321  TArrayF fGRNumerator;
322  TArrayF fGRDenominator;
323  TArrayF fGRNumeratorSub;
325 
330 
335 
340 
345 
350 
355 
356  Bool_t fHasGhost;
357  std::vector<TLorentzVector> fGhosts;
358 
359  private:
361  { // sort in decreasing order
362  // first value of the pair is Pt and the second is entry index
363  bool operator () (const std::pair<Double_t, Int_t>& p1, const std::pair<Double_t, Int_t>& p2) { return p1.first > p2.first ; }
364  };
365 
366  ClassDef(AliEmcalJet,17) // Emcal jet class in cylindrical coordinates
367 };
368 
369 std::ostream &operator<<(std::ostream &in, const AliEmcalJet &jet);
370 #endif
TArrayI fTrackIDs
Definition: AliEmcalJet.h:304
Double_t AreaEmc() const
Definition: AliEmcalJet.h:74
void SetMaxNeutralPt(Double32_t t)
Definition: AliEmcalJet.h:132
void SetSecondDerivativeConstituent(Double_t d)
Definition: AliEmcalJet.h:240
Int_t PdgCode() const
Definition: AliEmcalJet.h:64
void SetSecondClosestJet(AliEmcalJet *j, Double_t d)
Definition: AliEmcalJet.h:154
TArrayF fGRDenominatorSub
array with angular structure function numerator
Definition: AliEmcalJet.h:324
void SetSecondDerivativepTD(Double_t d)
Definition: AliEmcalJet.h:209
void SetTagStatus(Int_t i)
Definition: AliEmcalJet.h:166
Double_t fPtSubVect
background subtracted pt (not stored set from outside)
Definition: AliEmcalJet.h:312
Int_t NEmc() const
Definition: AliEmcalJet.h:96
Short_t Charge() const
Definition: AliEmcalJet.h:62
void AddFlavourTag(Int_t tag)
Definition: AliEmcalJet.h:260
Double_t Area() const
Definition: AliEmcalJet.h:69
void AddTrackAt(Int_t track, Int_t idx)
Definition: AliEmcalJet.h:115
Bool_t XvYvZv(Double_t x[3]) const
Definition: AliEmcalJet.h:53
Double_t GetSecondDerivativepTD() const
Definition: AliEmcalJet.h:213
AliEmcalJet * GetTaggedJet() const
Definition: AliEmcalJet.h:167
Double_t Xv() const
Definition: AliEmcalJet.h:50
Int_t fTagStatus
jet tagged to this jet
Definition: AliEmcalJet.h:310
TArrayF GetGRDenominatorSub() const
Definition: AliEmcalJet.h:185
Double_t fJetShapepTDFirstDer
result from shape derivatives for jet Angularity: 2nd order subtracted
Definition: AliEmcalJet.h:331
Double_t fJetShapeConstituentSecondDer
result from shape derivatives for jet const: 1st derivative
Definition: AliEmcalJet.h:347
Bool_t IsMC() const
Definition: AliEmcalJet.h:88
void SetGRDenSubSize(UInt_t s)
Definition: AliEmcalJet.h:194
Double_t fJetShapepTDSecondDer
result from shape derivatives for jet pTD: 1st derivative
Definition: AliEmcalJet.h:332
void SetTaggedJet(AliEmcalJet *j)
Definition: AliEmcalJet.h:165
void SetSecondDerivativeCircularity(Double_t d)
Definition: AliEmcalJet.h:219
Double_t GetXi(const AliVParticle *trk) const
Definition: AliEmcalJet.h:122
Double_t MCPt() const
Definition: AliEmcalJet.h:97
Double32_t fAreaE
Definition: AliEmcalJet.h:291
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:158
TArrayF GetGRDenominator() const
Definition: AliEmcalJet.h:183
Int_t GetTagStatus() const
Definition: AliEmcalJet.h:168
Double_t ClosestJetDistance() const
Definition: AliEmcalJet.h:159
TArrayI fClusterIDs
Definition: AliEmcalJet.h:303
Bool_t AxisInEmcal() const
Definition: AliEmcalJet.h:75
AliVParticle * GetFlavourTrack(Int_t i=0) const
Definition: AliEmcalJet.h:264
void ClearFlavourTracks()
Definition: AliEmcalJet.h:267
Double_t Eta() const
Definition: AliEmcalJet.h:60
Int_t fLabel
result from shape derivatives for jet mass: 2nd order subtracted
Definition: AliEmcalJet.h:319
Double_t fPtSub
status of tagging -1: NA 0: not tagged 1: tagged
Definition: AliEmcalJet.h:311
Int_t GetFlavour() const
Definition: AliEmcalJet.h:259
Double_t Py() const
Definition: AliEmcalJet.h:45
void Clear(Option_t *="")
Definition: AliEmcalJet.h:116
Double_t GetSecondDerivativeCircularity() const
Definition: AliEmcalJet.h:223
Double_t GetFirstDerivativeConstituent() const
Definition: AliEmcalJet.h:243
void SetSecondOrderSubtractedSigma2(Double_t d)
Definition: AliEmcalJet.h:231
Double_t Phi() const
Definition: AliEmcalJet.h:55
Double_t GetSecondOrderSubtractedLeSub() const
Definition: AliEmcalJet.h:256
Bool_t fHasGhost
result from shape derivatives for jet LeSub: 2nd order subtracted
Definition: AliEmcalJet.h:356
void AddClusterAt(Int_t clus, Int_t idx)
Definition: AliEmcalJet.h:114
Double_t GetFirstOrderSubtractedConstituent() const
Definition: AliEmcalJet.h:245
Double_t fJetShapeLeSubFirstSub
result from shape derivatives for jet LeSub: 2nd derivative
Definition: AliEmcalJet.h:353
void SetFirstOrderSubtractedAngularity(Double_t d)
Definition: AliEmcalJet.h:200
AliEmcalJet * MatchedJet() const
Definition: AliEmcalJet.h:162
void SetFirstOrderSubtractedSigma2(Double_t d)
Definition: AliEmcalJet.h:230
Double_t GetSecondDerivativeConstituent() const
Definition: AliEmcalJet.h:244
void SetFirstOrderSubtractedLeSub(Double_t d)
Definition: AliEmcalJet.h:251
Double_t AreaE() const
Definition: AliEmcalJet.h:73
Double_t fJetShapeMassSecondDer
result from shape derivatives for jet mass: 1st derivative
Definition: AliEmcalJet.h:316
Bool_t IsSortable() const
Definition: AliEmcalJet.h:89
Int_t GetLabel() const
Definition: AliEmcalJet.h:63
void SetGRDenSize(UInt_t s)
Definition: AliEmcalJet.h:189
void SetArea(Double_t a)
Definition: AliEmcalJet.h:126
Double_t FracEmcalArea() const
Definition: AliEmcalJet.h:85
Double_t MaxChargedPt() const
Definition: AliEmcalJet.h:91
AliEmcalJet & operator=(const AliEmcalJet &jet)
Double_t GetSecondDerivativeSigma2() const
Definition: AliEmcalJet.h:233
void GetMom(TLorentzVector &vec) const
Definition: AliEmcalJet.h:66
Double_t GetFlavourTrackZ(Int_t i=0) const
void SetAreaE(Double_t a)
Definition: AliEmcalJet.h:129
Double_t fJetShapeAngularitySecondDer
result from shape derivatives for jet Angularity: 1st derivative
Definition: AliEmcalJet.h:327
friend std::ostream & operator<<(std::ostream &in, const AliEmcalJet &jet)
void SetMCPt(Double_t p)
Definition: AliEmcalJet.h:139
Double_t E() const
Definition: AliEmcalJet.h:58
Double_t fJetShapeLeSubSecondSub
result from shape derivatives for jet LeSub: 1st order subtracted
Definition: AliEmcalJet.h:354
void SetPtSubVect(Double_t ps)
Definition: AliEmcalJet.h:145
Double_t GetFirstOrderSubtracted() const
Definition: AliEmcalJet.h:178
UShort_t Nch() const
Definition: AliEmcalJet.h:94
void SetLabel(Int_t l)
Definition: AliEmcalJet.h:125
void SetFirstDerivativeAngularity(Double_t d)
Definition: AliEmcalJet.h:198
Double_t AreaPt() const
Definition: AliEmcalJet.h:70
UShort_t GetNumberOfConstituents() const
Definition: AliEmcalJet.h:84
Int_t ContainsCluster(AliVCluster *cluster, TClonesArray *clusters) const
Definition: AliEmcalJet.h:79
void SetMatchedToSecondClosest(UShort_t m)
Definition: AliEmcalJet.h:156
Double_t GetFirstOrderSubtractedLeSub() const
Definition: AliEmcalJet.h:255
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:83
void SetSecondDerivativeLeSub(Double_t d)
Definition: AliEmcalJet.h:250
void SetGRNumSize(UInt_t s)
Definition: AliEmcalJet.h:188
Double_t fJetShapeCircularityFirstDer
result from shape derivatives for jet pTD: 2nd order subtracted
Definition: AliEmcalJet.h:336
Double_t Px() const
Definition: AliEmcalJet.h:44
std::ostream & Print(std::ostream &in) const
void SetSecondOrderSubtractedLeSub(Double_t d)
Definition: AliEmcalJet.h:252
void ResetMatching()
void GetMomentum(TLorentzVector &vec) const
Double32_t fMCPt
Definition: AliEmcalJet.h:298
Double_t GetFirstDerivativeAngularity() const
Definition: AliEmcalJet.h:202
Double32_t fMaxNPt
Definition: AliEmcalJet.h:297
Double32_t fM
Definition: AliEmcalJet.h:286
Double_t fJetShapepTDFirstSub
result from shape derivatives for jet pTD: 2nd derivative
Definition: AliEmcalJet.h:333
Double_t AreaEta() const
Definition: AliEmcalJet.h:71
void SetFirstOrderSubtractedConstituent(Double_t d)
Definition: AliEmcalJet.h:241
Short_t ClusterAt(Int_t idx) const
Definition: AliEmcalJet.h:77
void SetFirstOrderSubtracted(Double_t d)
Definition: AliEmcalJet.h:174
Bool_t IsInEmcal() const
Definition: AliEmcalJet.h:87
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:82
AliVParticle * GetLeadingTrack(TClonesArray *tracks) const
Double_t GetSecondDerivativeLeSub() const
Definition: AliEmcalJet.h:254
Double_t GetXi(const Double_t trkPx, const Double_t trkPy, const Double_t trkPz) const
Definition: AliEmcalJet.h:123
TArrayF fGRNumeratorSub
array with angular structure function denominator
Definition: AliEmcalJet.h:323
UShort_t GetMatchingType() const
Definition: AliEmcalJet.h:163
Double32_t fAreaPhi
Definition: AliEmcalJet.h:290
Double_t fJetShapeConstituentSecondSub
result from shape derivatives for jet const: 1st order subtracted
Definition: AliEmcalJet.h:349
Double_t Yv() const
Definition: AliEmcalJet.h:51
void AddFlavourTrack(AliVParticle *hftrack)
Definition: AliEmcalJet.h:263
Double_t GetSecondOrderSubtractedSigma2() const
Definition: AliEmcalJet.h:235
void AddGRDenAt(Float_t den, Int_t idx)
Definition: AliEmcalJet.h:187
UShort_t fMatched
distance to closest jets (see above)
Definition: AliEmcalJet.h:307
std::vector< int > SortConstituentsPt(TClonesArray *tracks) const
std::vector< TLorentzVector > fGhosts
Definition: AliEmcalJet.h:357
Double_t OneOverPt() const
Definition: AliEmcalJet.h:54
TString toString() const
void SetFirstDerivativeSigma2(Double_t d)
Definition: AliEmcalJet.h:228
void AddGRNumAt(Float_t num, Int_t idx)
Definition: AliEmcalJet.h:186
Double_t PtSubVect() const
Definition: AliEmcalJet.h:103
Double_t GetSecondDerivativeAngularity() const
Definition: AliEmcalJet.h:203
void SetSecondOrderSubtractedAngularity(Double_t d)
Definition: AliEmcalJet.h:201
Double32_t fArea
Definition: AliEmcalJet.h:288
void SetSecondOrderSubtracted(Double_t d)
Definition: AliEmcalJet.h:175
Double_t GetFirstDerivativepTD() const
Definition: AliEmcalJet.h:212
Double_t Theta() const
Definition: AliEmcalJet.h:57
Double_t MaxTrackPt() const
Definition: AliEmcalJet.h:99
AliVCluster * GetLeadingCluster(TClonesArray *clusters) const
AliEmcalJet * fTaggedJet
matching type
Definition: AliEmcalJet.h:309
Double_t GetSecondDerivative() const
Definition: AliEmcalJet.h:177
Int_t ContainsTrack(AliVParticle *track, TClonesArray *tracks) const
Definition: AliEmcalJet.h:109
Double_t Phi_0_2pi() const
Definition: AliEmcalJet.h:56
void SetFirstDerivativepTD(Double_t d)
Definition: AliEmcalJet.h:208
TLorentzVector SubtractRhoVect(Double_t rho, Bool_t save=kFALSE)
Double_t PtSub() const
Definition: AliEmcalJet.h:102
Double_t GetFirstDerivativeCircularity() const
Definition: AliEmcalJet.h:222
void SetFirstDerivativeCircularity(Double_t d)
Definition: AliEmcalJet.h:218
Double_t Y() const
Definition: AliEmcalJet.h:61
Double_t fJetShapeCircularityFirstSub
result from shape derivatives for jet circularity: 2nd derivative
Definition: AliEmcalJet.h:338
Double32_t fAreaEta
Definition: AliEmcalJet.h:289
void SetSecondOrderSubtractedConstituent(Double_t d)
Definition: AliEmcalJet.h:242
void SetPtSub(Double_t ps)
Definition: AliEmcalJet.h:144
void SetMaxChargedPt(Double32_t t)
Definition: AliEmcalJet.h:133
void SetNEF(Double_t nef)
Definition: AliEmcalJet.h:134
const Double_t * PID() const
Definition: AliEmcalJet.h:65
Double_t MaxNeutralPt() const
Definition: AliEmcalJet.h:90
void AddGhost(const Double_t dPx, const Double_t dPy, const Double_t dPz, const Double_t dE)
Definition: AliEmcalJet.h:269
Double_t GetFirstDerivativeSigma2() const
Definition: AliEmcalJet.h:232
Double_t GetZ(const Double_t trkPx, const Double_t trkPy, const Double_t trkPz) const
void SetNEmc(Int_t n)
Definition: AliEmcalJet.h:142
Double_t fJetShapeMassFirstDer
triggers that the jet might have fired (AliVEvent::EOfflineTriggerTypes)
Definition: AliEmcalJet.h:315
Double_t PtEmc() const
Definition: AliEmcalJet.h:101
Int_t fFlavourTagging
Definition: AliEmcalJet.h:294
AliEmcalJet * SecondClosestJet() const
Definition: AliEmcalJet.h:160
Double_t DeltaR(const AliVParticle *part) const
Double_t fJetShapeSigma2SecondSub
result from shape derivatives for jet sigma2: 1st order subtracted
Definition: AliEmcalJet.h:344
UShort_t N() const
Definition: AliEmcalJet.h:95
Double_t fJetShapeLeSubFirstDer
result from shape derivatives for jet const: 2nd order subtracted
Definition: AliEmcalJet.h:351
void SetFirstDerivativeLeSub(Double_t d)
Definition: AliEmcalJet.h:249
bool operator()(const std::pair< Double_t, Int_t > &p1, const std::pair< Double_t, Int_t > &p2)
Definition: AliEmcalJet.h:363
Double_t AreaPhi() const
Definition: AliEmcalJet.h:72
Double_t GetSecondOrderSubtracted() const
Definition: AliEmcalJet.h:179
Double_t Pt() const
Definition: AliEmcalJet.h:47
Double_t GetSecondOrderSubtractedpTD() const
Definition: AliEmcalJet.h:215
Double32_t fPtEmc
Definition: AliEmcalJet.h:301
void SetSecondDerivativeAngularity(Double_t d)
Definition: AliEmcalJet.h:199
Bool_t HasGhost() const
Definition: AliEmcalJet.h:276
Double_t GetFirstDerivativeLeSub() const
Definition: AliEmcalJet.h:253
void SetNumberOfCharged(Int_t n)
Definition: AliEmcalJet.h:137
Double_t fJetShapeCircularitySecondDer
result from shape derivatives for jet circularity: 1st derivative
Definition: AliEmcalJet.h:337
void SetGRNumSubSize(UInt_t s)
Definition: AliEmcalJet.h:193
void SetAreaEta(Double_t a)
Definition: AliEmcalJet.h:127
Double_t fJetShapeAngularityFirstSub
result from shape derivatives for jet Angularity: 2nd derivative
Definition: AliEmcalJet.h:328
Double_t GetFirstOrderSubtractedpTD() const
Definition: AliEmcalJet.h:214
Double32_t fPt
Definition: AliEmcalJet.h:283
void AddGRNumSubAt(Float_t num, Int_t idx)
Definition: AliEmcalJet.h:191
void SetSecondOrderSubtractedpTD(Double_t d)
Definition: AliEmcalJet.h:211
void AddGRDenSubAt(Float_t den, Int_t idx)
Definition: AliEmcalJet.h:192
Double_t P() const
Definition: AliEmcalJet.h:48
std::ostream & operator<<(std::ostream &in, const AliEmcalJet &jet)
Bool_t IsInsideEmcal() const
Definition: AliEmcalJet.h:86
void SetClosestJet(AliEmcalJet *j, Double_t d)
Definition: AliEmcalJet.h:153
Double32_t fEta
Definition: AliEmcalJet.h:284
void SetFirstOrderSubtractedpTD(Double_t d)
Definition: AliEmcalJet.h:210
Double32_t fMaxCPt
heavy flavour candidate tracks matched to the jet
Definition: AliEmcalJet.h:296
AliVParticle * RemoveFlavourTrack(Int_t i=0)
Definition: AliEmcalJet.h:266
Double_t GetFirstOrderSubtractedAngularity() const
Definition: AliEmcalJet.h:204
Short_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:107
void SetFirstDerivative(Double_t d)
Definition: AliEmcalJet.h:172
Double_t SecondClosestJetDistance() const
Definition: AliEmcalJet.h:161
Double_t GetFirstDerivative() const
Definition: AliEmcalJet.h:176
Bool_t TestFlavourTag(Int_t tag) const
Definition: AliEmcalJet.h:262
TArrayF fGRDenominator
array with angular structure function numerator
Definition: AliEmcalJet.h:322
Double_t fJetShapeConstituentFirstSub
result from shape derivatives for jet const: 2nd derivative
Definition: AliEmcalJet.h:348
Double32_t fAreaEmc
Definition: AliEmcalJet.h:292
Double_t fJetShapeSigma2FirstDer
result from shape derivatives for jetcircularity: 2nd order subtracted
Definition: AliEmcalJet.h:341
void SetSecondOrderSubtractedCircularity(Double_t d)
Definition: AliEmcalJet.h:221
UShort_t Nn() const
Definition: AliEmcalJet.h:93
Double_t fJetShapeSigma2SecondDer
result from shape derivatives for jet sigma2: 1st derivative
Definition: AliEmcalJet.h:342
Double_t fJetShapeConstituentFirstDer
result from shape derivatives for jetsigma2: 2nd order subtracted
Definition: AliEmcalJet.h:346
Double_t fJetShapeSigma2FirstSub
result from shape derivatives for jet sigma2: 2nd derivative
Definition: AliEmcalJet.h:343
void SetSecondDerivativeSigma2(Double_t d)
Definition: AliEmcalJet.h:229
Double_t GetFirstOrderSubtractedCircularity() const
Definition: AliEmcalJet.h:224
Double_t fJetShapeCircularitySecondSub
result from shape derivatives for jet circularity: 1st order subtracted
Definition: AliEmcalJet.h:339
Double_t MaxClusterPt() const
Definition: AliEmcalJet.h:98
Double_t fJetShapeMassSecondSub
result from shape derivatives for jet mass: 1st order subtracted
Definition: AliEmcalJet.h:318
Double_t fJetShapeLeSubSecondDer
result from shape derivatives for jet LeSub: 1st derivative
Definition: AliEmcalJet.h:352
void SetNumberOfTracks(Int_t n)
Definition: AliEmcalJet.h:136
Bool_t PxPyPz(Double_t p[3]) const
Definition: AliEmcalJet.h:49
Double_t Pz() const
Definition: AliEmcalJet.h:46
void SetFirstOrderSubtractedCircularity(Double_t d)
Definition: AliEmcalJet.h:220
const std::vector< TLorentzVector > GetGhosts() const
Definition: AliEmcalJet.h:277
void SetTrigger(UInt_t trigger)
Definition: AliEmcalJet.h:149
Double_t fJetShapeAngularitySecondSub
result from shape derivatives for jet Angularity: 1st order subtracted
Definition: AliEmcalJet.h:329
void SetFlavour(Int_t flavour)
Definition: AliEmcalJet.h:261
UInt_t fTriggers
background vector subtracted pt (not stored set from outside)
Definition: AliEmcalJet.h:313
Double_t GetSecondOrderSubtractedAngularity() const
Definition: AliEmcalJet.h:205
Bool_t fAxisInEmcal
Definition: AliEmcalJet.h:293
void PrintGR()
void SortConstituents()
TArrayF fGRNumerator
Definition: AliEmcalJet.h:321
Double32_t fPhi
Definition: AliEmcalJet.h:285
Double_t GetSecondOrderSubtractedCircularity() const
Definition: AliEmcalJet.h:225
AliEmcalJet * fClosestJets[2]
Definition: AliEmcalJet.h:305
TArrayF GetGRNumeratorSub() const
Definition: AliEmcalJet.h:184
void AddTrigger(UInt_t trigger)
Definition: AliEmcalJet.h:150
Double_t NEF() const
Definition: AliEmcalJet.h:92
void SetAreaPhi(Double_t a)
Definition: AliEmcalJet.h:128
void SetPtEmc(Double_t pt)
Definition: AliEmcalJet.h:143
Double_t fJetShapeMassFirstSub
result from shape derivatives for jet mass: 2nd derivative
Definition: AliEmcalJet.h:317
Bool_t IsTriggerJet(UInt_t trigger=AliVEvent::kEMCEJE) const
Definition: AliEmcalJet.h:148
void SetSecondDerivative(Double_t d)
Definition: AliEmcalJet.h:173
AliVParticle * TrackAt(Int_t idx, TClonesArray *ta) const
Definition: AliEmcalJet.h:108
void SetFirstDerivativeConstituent(Double_t d)
Definition: AliEmcalJet.h:239
Double32_t fNEF
Definition: AliEmcalJet.h:287
Double32_t fClosestJetsDist[2]
if this is MC it contains the two closest detector level jets in order of distance and viceversa ...
Definition: AliEmcalJet.h:306
UShort_t fMatchingType
0,1 if it is matched with one of the closest jets; 2 if it is not matched
Definition: AliEmcalJet.h:308
Double_t M() const
Definition: AliEmcalJet.h:59
Int_t Compare(const TObject *obj) const
AliVCluster * ClusterAt(Int_t idx, TClonesArray *ca) const
Definition: AliEmcalJet.h:78
Double_t MaxPartPt() const
Definition: AliEmcalJet.h:100
TArrayF GetGRNumerator() const
Definition: AliEmcalJet.h:182
Double_t GetSecondOrderSubtractedConstituent() const
Definition: AliEmcalJet.h:246
void SetMatchedToClosest(UShort_t m)
Definition: AliEmcalJet.h:155
Double_t fJetShapepTDSecondSub
result from shape derivatives for jet pTD: 1st order subtracted
Definition: AliEmcalJet.h:334
Double_t fJetShapeAngularityFirstDer
array with angular structure function denominator
Definition: AliEmcalJet.h:326
void SetNumberOfClusters(Int_t n)
Definition: AliEmcalJet.h:135
Double_t GetFirstOrderSubtractedSigma2() const
Definition: AliEmcalJet.h:234
Double_t Zv() const
Definition: AliEmcalJet.h:52
void SetAxisInEmcal(Bool_t b)
Definition: AliEmcalJet.h:131
TObjArray * fFlavourTracks
Definition: AliEmcalJet.h:295
void SetAreaEmc(Double_t a)
Definition: AliEmcalJet.h:130
void PrintConstituents(TClonesArray *tracks, TClonesArray *clusters) const
void SetNumberOfNeutrals(Int_t n)
Definition: AliEmcalJet.h:138