AliPhysics  8d00e07 (8d00e07)
AliEmcalClusTrackMatcherTask.cxx
Go to the documentation of this file.
1 //
2 // Track/cluster matcher
3 //
4 // Author: C.Loizides, S.Aiola
5 
7 
8 #include <TClonesArray.h>
9 #include <TClass.h>
10 
11 #include <AliAODCaloCluster.h>
12 #include <AliESDCaloCluster.h>
13 #include <AliLog.h>
14 #include <AliVCluster.h>
15 #include <AliVTrack.h>
16 #include <AliEMCALRecoUtils.h>
17 
18 #include "AliEmcalParticle.h"
19 #include "AliParticleContainer.h"
20 #include "AliClusterContainer.h"
21 
23 
24 //________________________________________________________________________
27  fPropDist(440),
28  fDoPropagation(kFALSE),
29  fAttemptProp(kTRUE),
30  fAttemptPropMatch(kFALSE),
31  fMaxDistance(0.1),
32  fAttachEmcalParticles(kFALSE),
33  fUpdateTracks(kTRUE),
34  fUpdateClusters(kTRUE),
35  fEmcalTracks(0),
36  fEmcalClusters(0),
37  fNEmcalTracks(0),
38  fNEmcalClusters(0),
39  fHistMatchEtaAll(0),
40  fHistMatchPhiAll(0)
41 {
42  // Constructor.
43 
44  for(Int_t icent=0; icent<8; ++icent) {
45  for(Int_t ipt=0; ipt<9; ++ipt) {
46  for(Int_t ieta=0; ieta<2; ++ieta) {
47  fHistMatchEta[icent][ipt][ieta] = 0;
48  fHistMatchPhi[icent][ipt][ieta] = 0;
49  }
50  }
51  }
52 }
53 
54 //________________________________________________________________________
56  AliAnalysisTaskEmcal(name, histo),
57  fPropDist(440),
58  fDoPropagation(kFALSE),
59  fAttemptProp(kTRUE),
60  fAttemptPropMatch(kFALSE),
61  fMaxDistance(0.1),
62  fAttachEmcalParticles(kFALSE),
63  fUpdateTracks(kTRUE),
64  fUpdateClusters(kTRUE),
65  fEmcalTracks(0),
66  fEmcalClusters(0),
67  fNEmcalTracks(0),
68  fNEmcalClusters(0),
69  fHistMatchEtaAll(0),
70  fHistMatchPhiAll(0)
71 {
72  // Standard constructor.
73 
74  for(Int_t icent=0; icent<8; ++icent) {
75  for(Int_t ipt=0; ipt<9; ++ipt) {
76  for(Int_t ieta=0; ieta<2; ++ieta) {
77  fHistMatchEta[icent][ipt][ieta] = 0;
78  fHistMatchPhi[icent][ipt][ieta] = 0;
79  }
80  }
81  }
82 }
83 
84 //________________________________________________________________________
86 {
87  // Destructor.
88 }
89 
90 //________________________________________________________________________
92 {
93  // Initialize the analysis.
94 
96  if (tracks) {
97  TClass trackClass(tracks->GetClassName());
98  if (!trackClass.InheritsFrom("AliVTrack")) {
99  tracks->SetClassName("AliVTrack"); // enforce only AliVTrack and derived classes
100  }
101  }
102 
104 
106  if (!fLocalInitialized) return;
107 
108  TString emcalTracksName(Form("EmcalTracks_%s", tracks->GetArrayName().Data()));
109  TString emcalClustersName(Form("EmcalClusters_%s", clusters->GetArrayName().Data()));
110 
111  fEmcalTracks = new TClonesArray("AliEmcalParticle");
112  fEmcalTracks->SetName(emcalTracksName);
113  fEmcalClusters = new TClonesArray("AliEmcalParticle");
114  fEmcalClusters->SetName(emcalClustersName);
115 
116  if (fAttachEmcalParticles) {
119  }
120 }
121 
122 //________________________________________________________________________
124 {
125  // Create my user objects.
126 
127  if (!fCreateHisto) return;
128 
130 
131  const Int_t nCentChBins = fNcentBins * 2;
132 
133  fHistMatchEtaAll = new TH1F("fHistMatchEtaAll", "fHistMatchEtaAll", 400, -0.2, 0.2);
134  fHistMatchPhiAll = new TH1F("fHistMatchPhiAll", "fHistMatchPhiAll", 400, -0.2, 0.2);
137 
138  for(Int_t icent=0; icent<nCentChBins; ++icent) {
139  for(Int_t ipt=0; ipt<9; ++ipt) {
140  for(Int_t ieta=0; ieta<2; ++ieta) {
141  TString nameEta(Form("fHistMatchEta_%i_%i_%i",icent,ipt,ieta));
142  fHistMatchEta[icent][ipt][ieta] = new TH1F(nameEta, nameEta, 400, -0.2, 0.2);
143  fHistMatchEta[icent][ipt][ieta]->SetXTitle("#Delta#eta");
144  TString namePhi(Form("fHistMatchPhi_%i_%i_%i",icent,ipt,ieta));
145  fHistMatchPhi[icent][ipt][ieta] = new TH1F(namePhi, namePhi, 400, -0.2, 0.2);
146  fHistMatchPhi[icent][ipt][ieta]->SetXTitle("#Delta#phi");
147  fOutput->Add(fHistMatchEta[icent][ipt][ieta]);
148  fOutput->Add(fHistMatchPhi[icent][ipt][ieta]);
149  }
150  }
151  }
152 
153  PostData(1, fOutput);
154 }
155 
156 //________________________________________________________________________
158 {
159  // Get momenum bin.
160 
161  Int_t pbin=-1;
162  if (p<0.5)
163  pbin=0;
164  else if (p>=0.5 && p<1.0)
165  pbin=1;
166  else if (p>=1.0 && p<1.5)
167  pbin=2;
168  else if (p>=1.5 && p<2.)
169  pbin=3;
170  else if (p>=2. && p<3.)
171  pbin=4;
172  else if (p>=3. && p<4.)
173  pbin=5;
174  else if (p>=4. && p<5.)
175  pbin=6;
176  else if (p>=5. && p<8.)
177  pbin=7;
178  else if (p>=8.)
179  pbin=8;
180 
181  return pbin;
182 }
183 
184 //________________________________________________________________________
186 {
187  // Run the matching.
188 
190  DoMatching();
193 
194  return kTRUE;
195 }
196 
197 //________________________________________________________________________
199 {
200  // Create AliEmcalParticle collections to handle the matching efficiently.
201  // At the same time propagates tracks, if requested.
202 
205 
206  fEmcalTracks->Delete();
207  fEmcalClusters->Delete();
208 
209  fNEmcalTracks = 0;
210  fNEmcalClusters = 0;
211 
212  AliVCluster* cluster = 0;
213  AliVTrack* track = 0;
214 
215  clusters->ResetCurrentID();
216  while ((cluster = static_cast<AliVCluster*>(clusters->GetNextAcceptCluster()))) {
217 
218  // Clears the matching info
219  cluster->SetEmcCpvDistance(-1);
220  cluster->SetTrackDistance(1024, 1024);
221  AliAODCaloCluster *ac = dynamic_cast<AliAODCaloCluster*>(cluster);
222  AliESDCaloCluster *ec = 0;
223  if (ac) {
224  const Int_t N = ac->GetNTracksMatched();
225  for (Int_t i = N - 1; i >= 0; i--) {
226  TObject *ptr = ac->GetTrackMatched(i);
227  ac->RemoveTrackMatched(ptr);
228  }
229  }
230  else {
231  ec = dynamic_cast<AliESDCaloCluster*>(cluster);
232  TArrayI *arr = ec->GetTracksMatched();
233  if (arr) arr->Set(0);
234  }
235 
236  // Create AliEmcalParticle objects to handle the matching
237  AliEmcalParticle* emcalCluster = new ((*fEmcalClusters)[fNEmcalClusters])
238  AliEmcalParticle(cluster, clusters->GetCurrentID(), fVertex[0], fVertex[1], fVertex[2], AliVCluster::kNonLinCorr);
239  emcalCluster->SetMatchedPtr(fEmcalTracks);
240 
241  fNEmcalClusters++;
242  }
243 
244  tracks->ResetCurrentID();
245  while ((track = static_cast<AliVTrack*>(tracks->GetNextAcceptParticle()))) {
246 
247  // Clears the matching info
248  track->ResetStatus(AliVTrack::kEMCALmatch);
249  track->SetEMCALcluster(-1);
250 
251  // Propagate tracks if requested
252  Bool_t propthistrack = kFALSE;
253  if (fDoPropagation) {
254  propthistrack = kTRUE;
255  }
256  else if (!track->IsExtrapolatedToEMCAL()) {
257  if (fAttemptProp) {
258  propthistrack = kTRUE;
259  }
260  else if (fAttemptPropMatch && IsTrackInEmcalAcceptance(track)) {
261  propthistrack = kTRUE;
262  }
263  }
264  if (propthistrack) AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(track, fPropDist);
265 
266  // Create AliEmcalParticle objects to handle the matching
267  AliEmcalParticle* emcalTrack = new ((*fEmcalTracks)[fNEmcalTracks]) AliEmcalParticle(track, tracks->GetCurrentID());
268  emcalTrack->SetMatchedPtr(fEmcalClusters);
269 
270  AliDebug(2, Form("Now adding track (pT = %.3f, eta = %.3f, phi = %.3f)"
271  "Phi, Eta on EMCal = %.3f, %.3f",
272  emcalTrack->Pt(), emcalTrack->Eta(), emcalTrack->Phi(),
273  track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal()));
274 
275  fNEmcalTracks++;
276  }
277 }
278 
279 //________________________________________________________________________
281 {
282  // Set the links between tracks and clusters.
283 
284  const Double_t maxd2 = fMaxDistance*fMaxDistance;
285 
286  for (Int_t itrack = 0; itrack < fNEmcalTracks; itrack++) {
287  AliEmcalParticle* emcalTrack = static_cast<AliEmcalParticle*>(fEmcalTracks->At(itrack));
288  AliVTrack* track = emcalTrack->GetTrack();
289 
290  for (Int_t icluster = 0; icluster < fNEmcalClusters; icluster++) {
291  AliEmcalParticle* emcalCluster = static_cast<AliEmcalParticle*>(fEmcalClusters->At(icluster));
292  AliVCluster* cluster = emcalCluster->GetCluster();
293 
294  Double_t deta = 999;
295  Double_t dphi = 999;
296  GetEtaPhiDiff(track, cluster, dphi, deta);
297  Double_t d2 = deta * deta + dphi * dphi;
298  if (d2 > maxd2) continue;
299 
300  Double_t d = TMath::Sqrt(d2);
301  emcalCluster->AddMatchedObj(itrack, d);
302  emcalTrack->AddMatchedObj(icluster, d);
303  AliDebug(2, Form("Now matching cluster E = %.3f, pT = %.3f, eta = %.3f, phi = %.3f "
304  "with track pT = %.3f, eta = %.3f, phi = %.3f"
305  "Track eta, phi on EMCal = %.3f, %.3f, d = %.3f",
306  cluster->GetNonLinCorrEnergy(), emcalCluster->Pt(), emcalCluster->Eta(), emcalCluster->Phi(),
307  emcalTrack->Pt(), emcalTrack->Eta(), emcalTrack->Phi(),
308  track->GetTrackEtaOnEMCal(), track->GetTrackPhiOnEMCal(), d));
309 
310  if (fCreateHisto) {
311  Int_t mombin = GetMomBin(track->P());
312  Int_t centbinch = fCentBin;
313  if (track->Charge() < 0) centbinch += fNcentBins;
314  Int_t etabin = 0;
315  if(track->Eta() > 0) etabin = 1;
316 
317  fHistMatchEta[centbinch][mombin][etabin]->Fill(deta);
318  fHistMatchPhi[centbinch][mombin][etabin]->Fill(dphi);
319  fHistMatchEtaAll->Fill(deta);
320  fHistMatchPhiAll->Fill(dphi);
321  }
322  }
323  }
324 }
325 
326 //________________________________________________________________________
328 {
329  // Update clusters with matching info.
330 
331  for (Int_t icluster = 0; icluster < fNEmcalClusters; icluster++) {
332  AliEmcalParticle* emcalCluster = static_cast<AliEmcalParticle*>(fEmcalClusters->At(icluster));
333  const Int_t N = emcalCluster->GetNumberOfMatchedObj();
334  AliVCluster* cluster = emcalCluster->GetCluster();
335  AliDebug(3, Form("Cluster E = %.2f, eta = %.2f, phi = %.2f, Nmatch = %d", cluster->GetNonLinCorrEnergy(), emcalCluster->Eta(), emcalCluster->Phi(), N));
336 
337  if (N <= 0) continue;
338 
339  // Set the first match distance
340  const UInt_t firstMatchId = emcalCluster->GetMatchedObjId();
341  AliEmcalParticle* emcalTrackFirstMatch = static_cast<AliEmcalParticle*>(fEmcalTracks->At(firstMatchId));
342  AliVTrack* trackFirstMatch = emcalTrackFirstMatch->GetTrack();
343  Double_t deta = 999;
344  Double_t dphi = 999;
345  GetEtaPhiDiff(trackFirstMatch, cluster, dphi, deta);
346  cluster->SetTrackDistance(dphi, deta);
347 
348  // Cast into ESD/AOD objects
349  AliAODCaloCluster *ac = dynamic_cast<AliAODCaloCluster*>(cluster);
350  AliESDCaloCluster *ec = 0;
351  if (!ac) ec = dynamic_cast<AliESDCaloCluster*>(cluster);
352 
353  // Copy the matched tracks in the cluster. Note: different methods for ESD/AOD
354  if (ac) {
355  for (Int_t i=0; i < N; ++i) {
356  Int_t id = emcalCluster->GetMatchedObjId(i);
357  AliEmcalParticle* emcalTrack = static_cast<AliEmcalParticle*>(fEmcalTracks->At(id));
358 
359  AliDebug(3, Form("Pt = %.2f, eta = %.2f, phi = %.2f", emcalTrack->Pt(), emcalTrack->Eta(), emcalTrack->Phi()));
360 
361  TObject *obj = emcalTrack->GetTrack();
362  ac->AddTrackMatched(obj);
363  }
364  }
365  else {
366  TArrayI arr(N);
367  for (Int_t i = 0; i < N; ++i) {
368  Int_t id = emcalCluster->GetMatchedObjId(i);
369  AliEmcalParticle* emcalTrack = static_cast<AliEmcalParticle*>(fEmcalTracks->At(id));
370  arr.AddAt(emcalTrack->IdInCollection(), i);
371  }
372  ec->AddTracksMatched(arr);
373  }
374  }
375 }
376 
377 //________________________________________________________________________
379 {
380  // Update tracks with matching info.
381 
382  for (Int_t itrack = 0; itrack < fNEmcalTracks; itrack++) {
383  AliEmcalParticle* emcalTrack = static_cast<AliEmcalParticle*>(fEmcalTracks->At(itrack));
384  if (emcalTrack->GetNumberOfMatchedObj() <= 0) continue;
385  AliEmcalParticle* emcalCluster = static_cast<AliEmcalParticle*>(fEmcalClusters->At(emcalTrack->GetMatchedObjId()));
386 
387  AliVTrack* track = emcalTrack->GetTrack();
388  track->SetEMCALcluster(emcalCluster->IdInCollection());
389  track->SetStatus(AliVTrack::kEMCALmatch);
390  }
391 }
virtual AliVParticle * GetNextAcceptParticle()
double Double_t
Definition: External.C:58
Int_t GetMatchedObjId(UShort_t i=0) const
Double_t Phi() const
UShort_t GetNumberOfMatchedObj() const
Base task in the EMCAL framework.
Bool_t fLocalInitialized
whether or not the task has been already initialized
void SetMatchedPtr(TObjArray *arr)
static void GetEtaPhiDiff(const AliVTrack *t, const AliVCluster *v, Double_t &phidiff, Double_t &etadiff)
Calculate and difference between a track (t) and a cluster (c).
TClonesArray * fEmcalClusters
emcal tracks
Int_t fCentBin
!event centrality bin
AliVTrack * GetTrack() const
Container for particles within the EMCAL framework.
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
Bool_t IsTrackInEmcalAcceptance(AliVParticle *part, Double_t edges=0.9) const
Determines if a track is inside the EMCal acceptance.
int Int_t
Definition: External.C:63
Double_t Eta() const
Double_t Pt() const
unsigned int UInt_t
Definition: External.C:33
Int_t fNcentBins
how many centrality bins
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
Int_t IdInCollection() const
void ExecOnce()
Perform steps needed to initialize the analysis.
TH1 * fHistMatchPhi[10][9][2]
deta distribution
Bool_t Run()
Run function. This is the core function of the analysis and contains the user code. Therefore users have to implement this function.
Int_t fNEmcalClusters
number of emcal tracks
TH1 * fHistMatchEta[10][9][2]
dphi distribution
AliEmcalList * fOutput
!output list
AliVCluster * GetCluster() const
Double_t fVertex[3]
!event vertex
Bool_t fCreateHisto
whether or not create histograms
void AddMatchedObj(Int_t id, Double_t d)
virtual void ExecOnce()
Perform steps needed to initialize the analysis.
void AddObjectToEvent(TObject *obj, Bool_t attempt=kFALSE)
Add object to event.
void UserCreateOutputObjects()
Main initialization function on the worker.
bool Bool_t
Definition: External.C:53
Container structure for EMCAL clusters.
AliVCluster * GetNextAcceptCluster()
TH1 * fHistMatchEtaAll
number of emcal clusters