coverity
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALTracker.cxx
CommitLineData
fe17d4cb 1//========================================================================
2// Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved.
3//
4// Author: The ALICE Off-line Project.
5// Contributors are mentioned in the code where appropriate.
6//
7// Permission to use, copy, modify and distribute this software and its
8// documentation strictly for non-commercial purposes is hereby granted
9// without fee, provided that the above copyright notice appears in all
10// copies and that both the copyright notice and this permission notice
11// appear in the supporting documentation. The authors make no claims
12// about the suitability of this software for any purpose. It is
13// provided "as is" without express or implied warranty.
14//========================================================================
15//
16// Class AliEMCALTracker
17// -----------------------
18// Implementation of the track matching method between barrel tracks and
19// EMCAL clusters.
20// Besides algorithm implementation, some cuts are required to be set
21// in order to define, for each track, an acceptance window where clusters
22// are searched to find best match (if any).
23// The class accepts as input an ESD container, and works directly on it,
24// simply setting, for each of its tracks, the fEMCALindex flag, for each
25// track which is matched to a cluster.
26// In order to use method, one must launch PropagateBack().
27//
28// ------------------------------------------------------------------------
29// author: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
5970dfe2 30// Revised by Rongrong 2010-05-31 (rongrong.ma@cern.ch)
fe17d4cb 31//=========================================================================
32
33#include <Riostream.h>
34#include <iomanip>
35
36#include <TFile.h>
37#include <TTree.h>
fe17d4cb 38#include <TList.h>
39#include <TString.h>
40#include <TVector3.h>
41#include <TClonesArray.h>
c61f0e70 42#include <TGeoMatrix.h>
fe17d4cb 43
44#include "AliLog.h"
af885e0f 45#include "AliESDEvent.h"
fe17d4cb 46#include "AliESDtrack.h"
89ffc0b0 47#include "AliESDCaloCluster.h"
fe17d4cb 48#include "AliEMCALRecPoint.h"
49#include "AliRunLoader.h"
50#include "AliEMCALTrack.h"
51#include "AliEMCALLoader.h"
c61f0e70 52#include "AliEMCALGeometry.h"
8ba062b1 53#include "AliEMCALReconstructor.h"
54#include "AliEMCALRecParam.h"
55#include "AliCDBEntry.h"
56#include "AliCDBManager.h"
3e3faf55 57#include "AliEMCALReconstructor.h"
ee602376 58#include "AliEMCALRecoUtils.h"
c61f0e70 59
fe17d4cb 60#include "AliEMCALTracker.h"
61
62ClassImp(AliEMCALTracker)
8ba062b1 63
fe17d4cb 64//
65//------------------------------------------------------------------------------
66//
67AliEMCALTracker::AliEMCALTracker()
8fc351e3 68: AliTracker(),
69 fCutPt(0),
70 fCutNITS(0),
71 fCutNTPC(50),
da34fafe 72 fStep(20),
8fc351e3 73 fTrackCorrMode(kTrackCorrMMB),
74 fClusterWindow(50),
75 fCutEta(0.025),
76 fCutPhi(0.05),
77 fTracks(0),
78 fClusters(0),
79 fGeom(0)
fe17d4cb 80{
5970dfe2 81 //
82 // Default constructor.
83 // Initializes all simple data members to default values,
84 // and all collections to NULL.
85 // Output file name is set to a default value.
86 //
0832a2bf 87 InitParameters();
fe17d4cb 88}
89//
90//------------------------------------------------------------------------------
91//
92AliEMCALTracker::AliEMCALTracker(const AliEMCALTracker& copy)
93 : AliTracker(),
5970dfe2 94 fCutPt(copy.fCutPt),
dcd86c5d 95 fCutNITS(copy.fCutNITS),
96 fCutNTPC(copy.fCutNTPC),
5970dfe2 97 fStep(copy.fStep),
98 fTrackCorrMode(copy.fTrackCorrMode),
8fc351e3 99 fClusterWindow(copy.fClusterWindow),
5970dfe2 100 fCutEta(copy.fCutEta),
101 fCutPhi(copy.fCutPhi),
fe17d4cb 102 fTracks((TObjArray*)copy.fTracks->Clone()),
103 fClusters((TObjArray*)copy.fClusters->Clone()),
c61f0e70 104 fGeom(copy.fGeom)
fe17d4cb 105{
5970dfe2 106 //
107 // Copy constructor
108 // Besides copying all parameters, duplicates all collections.
109 //
fe17d4cb 110}
111//
112//------------------------------------------------------------------------------
113//
f1d9131f 114AliEMCALTracker& AliEMCALTracker::operator=(const AliEMCALTracker& source)
115{ // assignment operator; use copy ctor
116 if (&source == this) return *this;
5970dfe2 117
f1d9131f 118 new (this) AliEMCALTracker(source);
119 return *this;
fe17d4cb 120}
121//
122//------------------------------------------------------------------------------
123//
8ba062b1 124void AliEMCALTracker::InitParameters()
125{
5970dfe2 126 //
127 // Retrieve initialization parameters
128 //
8ba062b1 129
130 // Check if the instance of AliEMCALRecParam exists,
3e3faf55 131 const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
3a2a23e1 132
ba6de5ea 133 if(!recParam){
134 AliFatal("Reconstruction parameters for EMCAL not set!");
8ba062b1 135 }
41f05b8c 136 else{
5970dfe2 137
138 fCutEta = recParam->GetMthCutEta();
139 fCutPhi = recParam->GetMthCutPhi();
8fc351e3 140 fStep = recParam->GetExtrapolateStep();
5970dfe2 141 fCutPt = recParam->GetTrkCutPt();
8fc351e3 142 fCutNITS = recParam->GetTrkCutNITS();
143 fCutNTPC = recParam->GetTrkCutNTPC();
41f05b8c 144 }
3e3faf55 145
8ba062b1 146}
5970dfe2 147
c61f0e70 148//
149//------------------------------------------------------------------------------
150//
fe17d4cb 151void AliEMCALTracker::Clear(Option_t* option)
152{
153 //
154 // Clearing method
044225d9 155 // Deletes all objects in arrays and the arrays themselves
fe17d4cb 156 //
157
158 TString opt(option);
159 Bool_t clearTracks = opt.Contains("TRACKS");
160 Bool_t clearClusters = opt.Contains("CLUSTERS");
fe17d4cb 161 if (opt.Contains("ALL")) {
162 clearTracks = kTRUE;
163 clearClusters = kTRUE;
fe17d4cb 164 }
165
5970dfe2 166 //fTracks is a collection of esdTrack
167 //When clearing this array, the linked objects should not be deleted
fe17d4cb 168 if (fTracks != 0x0 && clearTracks) {
5970dfe2 169 fTracks->Clear();
044225d9 170 delete fTracks;
171 fTracks = 0;
fe17d4cb 172 }
173 if (fClusters != 0x0 && clearClusters) {
044225d9 174 fClusters->Delete();
175 delete fClusters;
176 fClusters = 0;
fe17d4cb 177 }
fe17d4cb 178}
179//
180//------------------------------------------------------------------------------
181//
182Int_t AliEMCALTracker::LoadClusters(TTree *cTree)
183{
184 //
185 // Load EMCAL clusters in the form of AliEMCALRecPoint,
186 // from simulation temporary files.
187 // (When included in reconstruction chain, this method is used automatically)
188 //
189
190 Clear("CLUSTERS");
191
bce21ea7 192 cTree->SetBranchStatus("*",0); //disable all branches
193 cTree->SetBranchStatus("EMCALECARP",1); //Enable only the branch we need
194
fe17d4cb 195 TBranch *branch = cTree->GetBranch("EMCALECARP");
196 if (!branch) {
c61f0e70 197 AliError("Can't get the branch with the EMCAL clusters");
fe17d4cb 198 return 1;
199 }
200
9596d957 201 TClonesArray *clusters = new TClonesArray("AliEMCALRecPoint", 1000);
fe17d4cb 202 branch->SetAddress(&clusters);
fe17d4cb 203
bce21ea7 204 //cTree->GetEvent(0);
205 branch->GetEntry(0);
9596d957 206 Int_t nClusters = (Int_t)clusters->GetEntries();
bce21ea7 207 if(fClusters) fClusters->Delete();
208 else fClusters = new TObjArray(0);
fe17d4cb 209 for (Int_t i = 0; i < nClusters; i++) {
210 AliEMCALRecPoint *cluster = (AliEMCALRecPoint*)clusters->At(i);
211 if (!cluster) continue;
fe17d4cb 212 AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
213 fClusters->AddLast(matchCluster);
214 }
3d9e8b15 215
2ad4424e 216 branch->SetAddress(0);
217 clusters->Delete();
218 delete clusters;
219
5970dfe2 220 AliInfo(Form("Collected %d RecPoints from Tree", fClusters->GetEntries()));
fe17d4cb 221
222 return 0;
223}
224//
225//------------------------------------------------------------------------------
226//
af885e0f 227Int_t AliEMCALTracker::LoadClusters(AliESDEvent *esd)
fe17d4cb 228{
5970dfe2 229 //
230 // Load EMCAL clusters in the form of AliESDCaloClusters,
231 // from an AliESD object.
232 //
233
234 // make sure that tracks/clusters collections are empty
235 Clear("CLUSTERS");
236 fClusters = new TObjArray(0);
237
238 Int_t nClusters = esd->GetNumberOfCaloClusters();
239 for (Int_t i=0; i<nClusters; i++)
240 {
241 AliESDCaloCluster *cluster = esd->GetCaloCluster(i);
242 if (!cluster || !cluster->IsEMCAL()) continue ;
243 AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
244 fClusters->AddLast(matchCluster);
245 }
246
247 AliInfo(Form("Collected %d clusters from ESD", fClusters->GetEntries()));
248 return 0;
fe17d4cb 249}
250//
251//------------------------------------------------------------------------------
252//
af885e0f 253Int_t AliEMCALTracker::LoadTracks(AliESDEvent *esd)
fe17d4cb 254{
5970dfe2 255 //
256 // Load ESD tracks.
257 //
fe17d4cb 258
5970dfe2 259 Clear("TRACKS");
260 fTracks = new TObjArray(0);
fe17d4cb 261
5970dfe2 262 Int_t nTracks = esd->GetNumberOfTracks();
8fc351e3 263 //Bool_t isKink=kFALSE;
5970dfe2 264 for (Int_t i = 0; i < nTracks; i++)
265 {
266 AliESDtrack *esdTrack = esd->GetTrack(i);
267 // set by default the value corresponding to "no match"
268 esdTrack->SetEMCALcluster(kUnmatched);
8fc351e3 269 esdTrack->ResetStatus(AliESDtrack::kEMCALmatch);
5970dfe2 270
271 //Select good quaulity tracks
272 if(esdTrack->Pt()<fCutPt) continue;
273 if(esdTrack->GetNcls(1)<fCutNTPC)continue;
274
5970dfe2 275 //Loose geometric cut
276 Double_t phi = esdTrack->Phi()*TMath::RadToDeg();
d874f457 277 if(TMath::Abs(esdTrack->Eta())>0.8 || phi <= 20 || phi >= 240 ) continue;
5970dfe2 278
279 fTracks->AddLast(esdTrack);
280 }
281
282 AliInfo(Form("Collected %d tracks", fTracks->GetEntries()));
283 return 0;
284}
285//
286//------------------------------------------------------------------------------
287//
288void AliEMCALTracker::SetTrackCorrectionMode(Option_t *option)
289{
290 //
291 // Set track correction mode
292 // gest the choice in string format and converts into
293 // internal enum
294 //
295
296 TString opt(option);
297 opt.ToUpper();
298
299 if (!opt.CompareTo("NONE"))
300 {
301 fTrackCorrMode = kTrackCorrNone;
302 }
303 else if (!opt.CompareTo("MMB"))
304 {
305 fTrackCorrMode = kTrackCorrMMB;
306 }
307 else
308 {
309 cerr << "E-AliEMCALTracker::SetTrackCorrectionMode '" << option << "': Unrecognized option" << endl;
310 }
fe17d4cb 311}
312//
313//------------------------------------------------------------------------------
314//
af885e0f 315Int_t AliEMCALTracker::PropagateBack(AliESDEvent* esd)
fe17d4cb 316{
317 //
318 // Main operation method.
319 // Gets external AliESD containing tracks to be matched.
320 // After executing match finding, stores in the same ESD object all infos
321 // and releases the object for further reconstruction steps.
322 //
3d9e8b15 323 //
324 // Note: should always return 0=OK, because otherwise all tracking
325 // is aborted for this event
5970dfe2 326
fe17d4cb 327 if (!esd) {
328 AliError("NULL ESD passed");
329 return 1;
330 }
331
5970dfe2 332 // step 1: collect clusters
c61f0e70 333 Int_t okLoadClusters, nClusters;
fe17d4cb 334 if (!fClusters || (fClusters && fClusters->IsEmpty())) {
fe17d4cb 335 okLoadClusters = LoadClusters(esd);
fe17d4cb 336 }
c61f0e70 337 nClusters = fClusters->GetEntries();
5970dfe2 338
339 // step 2: collect ESD tracks
340 Int_t nTracks, okLoadTracks;
341 okLoadTracks = LoadTracks(esd);
c61f0e70 342 nTracks = fTracks->GetEntries();
fe17d4cb 343
5970dfe2 344 // step 3: for each track, find the closest cluster as matched within residual cuts
345 Int_t index=-1;
346 for (Int_t it = 0; it < nTracks; it++)
347 {
348 AliESDtrack *track = (AliESDtrack*)fTracks->At(it);
349 index = FindMatchedCluster(track);
350 if (index>-1)
351 {
352 AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(index);
353 track->SetEMCALcluster(cluster->Index());
952d023a 354 track->SetStatus(AliESDtrack::kEMCALmatch);
5970dfe2 355 }
356 }
fe17d4cb 357
358 return 0;
359}
fe17d4cb 360
361//
362//------------------------------------------------------------------------------
363//
5970dfe2 364Int_t AliEMCALTracker::FindMatchedCluster(AliESDtrack *track)
365{
366 //
367 // For each track, extrapolate it to all the clusters
368 // Find the closest one as matched if the residuals (dEta, dPhi) satisfy the cuts
369 //
17773e2e 370
ee602376 371 Float_t maxEta=fCutEta;
372 Float_t maxPhi=fCutPhi;
5970dfe2 373 Int_t index = -1;
374
375 // If the esdFriend is available, use the TPCOuter point as the starting point of extrapolation
376 // Otherwise use the TPCInner point
8108ab9e 377 AliExternalTrackParam *trkParam = 0;
5970dfe2 378 const AliESDfriendTrack* friendTrack = track->GetFriendTrack();
379 if(friendTrack && friendTrack->GetTPCOut())
380 trkParam = const_cast<AliExternalTrackParam*>(friendTrack->GetTPCOut());
381 else
382 trkParam = const_cast<AliExternalTrackParam*>(track->GetInnerParam());
383 if(!trkParam) return index;
384
8fc351e3 385
386 AliExternalTrackParam trkParamTmp(*trkParam);
ee602376 387 Float_t eta, phi;
388 if(!AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&trkParamTmp, 430., track->GetMass(), fStep, eta, phi)) return index;
8fc351e3 389 if(TMath::Abs(eta)>0.75 || (phi) < 70*TMath::DegToRad() || (phi) > 190*TMath::DegToRad()) return index;
390
5970dfe2 391 //Perform extrapolation
392 Double_t trkPos[3];
8fc351e3 393 trkParamTmp.GetXYZ(trkPos);
5970dfe2 394 Int_t nclusters = fClusters->GetEntries();
395 for(Int_t ic=0; ic<nclusters; ic++)
396 {
5970dfe2 397 AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(ic);
8fc351e3 398 Float_t clsPos[3] = {cluster->X(),cluster->Y(),cluster->Z()};
399 Double_t dR = TMath::Sqrt(TMath::Power(trkPos[0]-clsPos[0],2)+TMath::Power(trkPos[1]-clsPos[1],2)+TMath::Power(trkPos[2]-clsPos[2],2));
400 if(dR > fClusterWindow) continue;
401
402 AliExternalTrackParam trkParTmp(trkParamTmp);
403
ee602376 404 Float_t tmpEta, tmpPhi;
da34fafe 405 if(!AliEMCALRecoUtils::ExtrapolateTrackToPosition(&trkParTmp, clsPos,track->GetMass(), 5, tmpEta, tmpPhi)) continue;
5970dfe2 406 if(TMath::Abs(tmpPhi)<TMath::Abs(maxPhi) && TMath::Abs(tmpEta)<TMath::Abs(maxEta))
407 {
408 maxPhi=tmpPhi;
409 maxEta=tmpEta;
410 index=ic;
411 }
5970dfe2 412 }
413 return index;
fe17d4cb 414}
5970dfe2 415
fe17d4cb 416//
417//------------------------------------------------------------------------------
418//
419void AliEMCALTracker::UnloadClusters()
420{
421 //
044225d9 422 // Free memory from all arrays
423 // This method is called after the local tracking step
424 // so we can safely delete everything
fe17d4cb 425 //
426
044225d9 427 Clear();
fe17d4cb 428}
04475328 429
430//
431//------------------------------------------------------------------------------
432//
fe17d4cb 433AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliEMCALRecPoint *recPoint)
434 : fIndex(index),
fe17d4cb 435 fX(0.),
436 fY(0.),
437 fZ(0.)
438{
439 //
440 // Translates an AliEMCALRecPoint object into the internal format.
441 // Index of passed cluster in its native array must be specified.
442 //
443 TVector3 clpos;
444 recPoint->GetGlobalPosition(clpos);
445
446 fX = clpos.X();
447 fY = clpos.Y();
448 fZ = clpos.Z();
449}
450//
451//------------------------------------------------------------------------------
452//
453AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliESDCaloCluster *caloCluster)
454 : fIndex(index),
fe17d4cb 455 fX(0.),
456 fY(0.),
457 fZ(0.)
458{
459 //
460 // Translates an AliESDCaloCluster object into the internal format.
461 // Index of passed cluster in its native array must be specified.
462 //
53e430a3 463 Float_t clpos[3]= {0., 0., 0.};
7592dfc4 464 caloCluster->GetPosition(clpos);
fe17d4cb 465
466 fX = (Double_t)clpos[0];
467 fY = (Double_t)clpos[1];
468 fZ = (Double_t)clpos[2];
469}