SensorThickness was defined twice. Set inner chip thickness to 250mum to bypass bug...
[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
0267cfa6 62using std::cerr;
63using std::endl;
fe17d4cb 64ClassImp(AliEMCALTracker)
8ba062b1 65
fe17d4cb 66//
67//------------------------------------------------------------------------------
68//
7ead3142 69AliEMCALTracker::AliEMCALTracker() :
70 AliTracker(),
8fc351e3 71 fCutPt(0),
72 fCutNITS(0),
73 fCutNTPC(50),
da34fafe 74 fStep(20),
8fc351e3 75 fTrackCorrMode(kTrackCorrMMB),
a29b2a8a 76 fEMCalSurfaceDistance(440),
8fc351e3 77 fClusterWindow(50),
78 fCutEta(0.025),
79 fCutPhi(0.05),
42ceff04 80 fITSTrackSA(kFALSE),
8fc351e3 81 fTracks(0),
82 fClusters(0),
83 fGeom(0)
fe17d4cb 84{
5970dfe2 85 //
86 // Default constructor.
87 // Initializes all simple data members to default values,
88 // and all collections to NULL.
89 // Output file name is set to a default value.
90 //
0832a2bf 91 InitParameters();
fe17d4cb 92}
93//
94//------------------------------------------------------------------------------
95//
7ead3142 96AliEMCALTracker::AliEMCALTracker(const AliEMCALTracker& copy) :
97 AliTracker(),
98 fCutPt(copy.fCutPt),
99 fCutNITS(copy.fCutNITS),
100 fCutNTPC(copy.fCutNTPC),
101 fStep(copy.fStep),
102 fTrackCorrMode(copy.fTrackCorrMode),
103 fEMCalSurfaceDistance(copy.fEMCalSurfaceDistance),
104 fClusterWindow(copy.fClusterWindow),
105 fCutEta(copy.fCutEta),
106 fCutPhi(copy.fCutPhi),
107 fITSTrackSA(copy.fITSTrackSA),
108 fTracks((TObjArray*)copy.fTracks->Clone()),
109 fClusters((TObjArray*)copy.fClusters->Clone()),
110 fGeom(copy.fGeom)
fe17d4cb 111{
5970dfe2 112 //
113 // Copy constructor
114 // Besides copying all parameters, duplicates all collections.
115 //
fe17d4cb 116}
117//
118//------------------------------------------------------------------------------
119//
f1d9131f 120AliEMCALTracker& AliEMCALTracker::operator=(const AliEMCALTracker& source)
121{ // assignment operator; use copy ctor
122 if (&source == this) return *this;
5970dfe2 123
f1d9131f 124 new (this) AliEMCALTracker(source);
125 return *this;
fe17d4cb 126}
127//
128//------------------------------------------------------------------------------
129//
8ba062b1 130void AliEMCALTracker::InitParameters()
131{
5970dfe2 132 //
133 // Retrieve initialization parameters
134 //
8ba062b1 135
136 // Check if the instance of AliEMCALRecParam exists,
3e3faf55 137 const AliEMCALRecParam* recParam = AliEMCALReconstructor::GetRecParam();
3a2a23e1 138
a4974ace 139 if (!recParam) {
ba6de5ea 140 AliFatal("Reconstruction parameters for EMCAL not set!");
7ead3142 141 } else {
5970dfe2 142 fCutEta = recParam->GetMthCutEta();
143 fCutPhi = recParam->GetMthCutPhi();
8fc351e3 144 fStep = recParam->GetExtrapolateStep();
5970dfe2 145 fCutPt = recParam->GetTrkCutPt();
8fc351e3 146 fCutNITS = recParam->GetTrkCutNITS();
147 fCutNTPC = recParam->GetTrkCutNTPC();
41f05b8c 148 }
8ba062b1 149}
c61f0e70 150//
151//------------------------------------------------------------------------------
152//
fe17d4cb 153void AliEMCALTracker::Clear(Option_t* option)
154{
7ead3142 155 //
156 // Clearing method
157 // Deletes all objects in arrays and the arrays themselves
158 //
fe17d4cb 159
7ead3142 160 TString opt(option);
161 Bool_t clearTracks = opt.Contains("TRACKS");
162 Bool_t clearClusters = opt.Contains("CLUSTERS");
163 if (opt.Contains("ALL")) {
164 clearTracks = kTRUE;
165 clearClusters = kTRUE;
166 }
fe17d4cb 167
7ead3142 168 //fTracks is a collection of esdTrack
169 //When clearing this array, the linked objects should not be deleted
170 if (fTracks != 0x0 && clearTracks) {
171 fTracks->Clear();
172 delete fTracks;
173 fTracks = 0;
174 }
175 if (fClusters != 0x0 && clearClusters) {
176 fClusters->Delete();
177 delete fClusters;
178 fClusters = 0;
179 }
fe17d4cb 180}
181//
182//------------------------------------------------------------------------------
183//
184Int_t AliEMCALTracker::LoadClusters(TTree *cTree)
185{
7ead3142 186 //
187 // Load EMCAL clusters in the form of AliEMCALRecPoint,
188 // from simulation temporary files.
189 // (When included in reconstruction chain, this method is used automatically)
190 //
fe17d4cb 191
7ead3142 192 Clear("CLUSTERS");
193
194 cTree->SetBranchStatus("*",0); //disable all branches
195 cTree->SetBranchStatus("EMCALECARP",1); //Enable only the branch we need
bce21ea7 196
7ead3142 197 TBranch *branch = cTree->GetBranch("EMCALECARP");
198 if (!branch) {
199 AliError("Can't get the branch with the EMCAL clusters");
200 return 1;
201 }
fe17d4cb 202
7ead3142 203 TClonesArray *clusters = new TClonesArray("AliEMCALRecPoint", 1000);
204 branch->SetAddress(&clusters);
fe17d4cb 205
7ead3142 206 //cTree->GetEvent(0);
207 branch->GetEntry(0);
208 Int_t nClusters = (Int_t)clusters->GetEntries();
a4974ace 209 if (fClusters) fClusters->Delete();
7ead3142 210 else fClusters = new TObjArray(0);
211 for (Int_t i = 0; i < nClusters; i++) {
212 AliEMCALRecPoint *cluster = (AliEMCALRecPoint*)clusters->At(i);
213 if (!cluster) continue;
214 AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
215 fClusters->AddLast(matchCluster);
216 }
2ad4424e 217
7ead3142 218 branch->SetAddress(0);
219 clusters->Delete();
220 delete clusters;
fe17d4cb 221
7ead3142 222 AliInfo(Form("Collected %d RecPoints from Tree", fClusters->GetEntries()));
223
224 return 0;
fe17d4cb 225}
226//
227//------------------------------------------------------------------------------
228//
af885e0f 229Int_t AliEMCALTracker::LoadClusters(AliESDEvent *esd)
fe17d4cb 230{
5970dfe2 231 //
232 // Load EMCAL clusters in the form of AliESDCaloClusters,
233 // from an AliESD object.
234 //
235
236 // make sure that tracks/clusters collections are empty
237 Clear("CLUSTERS");
238 fClusters = new TObjArray(0);
239
240 Int_t nClusters = esd->GetNumberOfCaloClusters();
7ead3142 241 for (Int_t i=0; i<nClusters; i++) {
242 AliESDCaloCluster *cluster = esd->GetCaloCluster(i);
243 if (!cluster || !cluster->IsEMCAL()) continue ;
244 AliEMCALMatchCluster *matchCluster = new AliEMCALMatchCluster(i, cluster);
245 fClusters->AddLast(matchCluster);
246 }
5970dfe2 247
248 AliInfo(Form("Collected %d clusters from ESD", fClusters->GetEntries()));
249 return 0;
fe17d4cb 250}
251//
252//------------------------------------------------------------------------------
253//
af885e0f 254Int_t AliEMCALTracker::LoadTracks(AliESDEvent *esd)
fe17d4cb 255{
5970dfe2 256 //
257 // Load ESD tracks.
258 //
42ceff04 259
26050574 260 UInt_t mask1 = esd->GetESDRun()->GetDetectorsInDAQ();
261 UInt_t mask2 = esd->GetESDRun()->GetDetectorsInReco();
262 Bool_t desc1 = (mask1 >> 3) & 0x1;
263 Bool_t desc2 = (mask2 >> 3) & 0x1;
264 if (desc1==0 || desc2==0) {
30e29b2a 265// AliError(Form("TPC not in DAQ/RECO: %u (%u)/%u (%u)",
266// mask1, esd->GetESDRun()->GetDetectorsInReco(),
267// mask2, esd->GetESDRun()->GetDetectorsInDAQ()));
26050574 268 fITSTrackSA = kTRUE;
42ceff04 269 }
270
5970dfe2 271 Clear("TRACKS");
272 fTracks = new TObjArray(0);
fe17d4cb 273
5970dfe2 274 Int_t nTracks = esd->GetNumberOfTracks();
8fc351e3 275 //Bool_t isKink=kFALSE;
7ead3142 276 for (Int_t i = 0; i < nTracks; i++) {
277 AliESDtrack *esdTrack = esd->GetTrack(i);
278 // set by default the value corresponding to "no match"
279 esdTrack->SetEMCALcluster(kUnmatched);
280 esdTrack->ResetStatus(AliESDtrack::kEMCALmatch);
281
282 //Select good quaulity tracks
a4974ace 283 if (esdTrack->Pt()<fCutPt) continue;
284 if (!fITSTrackSA)
285 if (esdTrack->GetNcls(1)<fCutNTPC) continue;
7ead3142 286
287 //Loose geometric cut
288 Double_t phi = esdTrack->Phi()*TMath::RadToDeg();
a4974ace 289 if (TMath::Abs(esdTrack->Eta())>0.9 || phi <= 10 || phi >= 250) continue;
7ead3142 290 fTracks->AddLast(esdTrack);
291 }
292 AliInfo(Form("Collected %d tracks", fTracks->GetEntries()));
293 return 0;
5970dfe2 294}
295//
296//------------------------------------------------------------------------------
297//
298void AliEMCALTracker::SetTrackCorrectionMode(Option_t *option)
299{
300 //
301 // Set track correction mode
302 // gest the choice in string format and converts into
303 // internal enum
304 //
305
306 TString opt(option);
307 opt.ToUpper();
308
7ead3142 309 if (!opt.CompareTo("NONE")) {
310 fTrackCorrMode = kTrackCorrNone;
311 } else if (!opt.CompareTo("MMB")) {
312 fTrackCorrMode = kTrackCorrMMB;
313 } else {
314 cerr << "E-AliEMCALTracker::SetTrackCorrectionMode '" << option << "': Unrecognized option" << endl;
315 }
fe17d4cb 316}
317//
318//------------------------------------------------------------------------------
319//
af885e0f 320Int_t AliEMCALTracker::PropagateBack(AliESDEvent* esd)
fe17d4cb 321{
7ead3142 322 //
323 // Main operation method.
324 // Gets external AliESD containing tracks to be matched.
325 // After executing match finding, stores in the same ESD object all infos
326 // and releases the object for further reconstruction steps.
327 //
328 //
329 // Note: should always return 0=OK, because otherwise all tracking
330 // is aborted for this event
5970dfe2 331
852a34a0 332 if (!esd)
333 {
7ead3142 334 AliError("NULL ESD passed");
335 return 1;
336 }
fe17d4cb 337
7ead3142 338 // step 1: collect clusters
339 Int_t okLoadClusters, nClusters;
852a34a0 340
341 if (!fClusters || (fClusters && fClusters->IsEmpty()))
7ead3142 342 okLoadClusters = LoadClusters(esd);
852a34a0 343
7ead3142 344 nClusters = fClusters->GetEntries();
5970dfe2 345
7ead3142 346 // step 2: collect ESD tracks
347 Int_t nTracks, okLoadTracks;
348 okLoadTracks = LoadTracks(esd);
349 nTracks = fTracks->GetEntries();
350
852a34a0 351 AliDebug(5,Form("Propagate back %d tracks ok %d, for %d clusters ok %d",
352 nTracks,okLoadTracks,nClusters,okLoadClusters));
353
7ead3142 354 // step 3: for each track, find the closest cluster as matched within residual cuts
355 Int_t index=-1;
852a34a0 356 for (Int_t it = 0; it < nTracks; it++)
357 {
7ead3142 358 AliESDtrack *track = (AliESDtrack*)fTracks->At(it);
359 index = FindMatchedCluster(track);
852a34a0 360 if (index>-1)
361 {
7ead3142 362 AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(index);
363 track->SetEMCALcluster(cluster->Index());
364 track->SetStatus(AliESDtrack::kEMCALmatch);
365 }
366 }
fe17d4cb 367
7ead3142 368 return 0;
fe17d4cb 369}
fe17d4cb 370
371//
372//------------------------------------------------------------------------------
373//
5970dfe2 374Int_t AliEMCALTracker::FindMatchedCluster(AliESDtrack *track)
375{
376 //
377 // For each track, extrapolate it to all the clusters
378 // Find the closest one as matched if the residuals (dEta, dPhi) satisfy the cuts
379 //
17773e2e 380
ee602376 381 Float_t maxEta=fCutEta;
382 Float_t maxPhi=fCutPhi;
5970dfe2 383 Int_t index = -1;
384
385 // If the esdFriend is available, use the TPCOuter point as the starting point of extrapolation
386 // Otherwise use the TPCInner point
8108ab9e 387 AliExternalTrackParam *trkParam = 0;
42ceff04 388
852a34a0 389 if (!fITSTrackSA)
390 {
42ceff04 391 const AliESDfriendTrack* friendTrack = track->GetFriendTrack();
852a34a0 392
a4974ace 393 if (friendTrack && friendTrack->GetTPCOut())
42ceff04 394 trkParam = const_cast<AliExternalTrackParam*>(friendTrack->GetTPCOut());
a4974ace 395 else if (track->GetInnerParam())
42ceff04 396 trkParam = const_cast<AliExternalTrackParam*>(track->GetInnerParam());
397 }
5970dfe2 398 else
42ceff04 399 trkParam = new AliExternalTrackParam(*track);
400
a4974ace 401 if (!trkParam) return index;
42ceff04 402
8fc351e3 403 AliExternalTrackParam trkParamTmp(*trkParam);
a29b2a8a 404 Float_t eta, phi, pt;
852a34a0 405 if (!AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&trkParamTmp, fEMCalSurfaceDistance, track->GetMass(kTRUE), fStep, eta, phi, pt))
406 {
a4974ace 407 if (fITSTrackSA) delete trkParam;
408 return index;
30e29b2a 409 }
852a34a0 410
a29b2a8a 411 track->SetTrackPhiEtaPtOnEMCal(phi,eta,pt);
852a34a0 412
413 if (TMath::Abs(eta)>0.75 || (phi) < 70*TMath::DegToRad() || (phi) > 190*TMath::DegToRad())
414 {
a4974ace 415 if (fITSTrackSA) delete trkParam;
416 return index;
30e29b2a 417 }
8fc351e3 418
5970dfe2 419 //Perform extrapolation
420 Double_t trkPos[3];
8fc351e3 421 trkParamTmp.GetXYZ(trkPos);
5970dfe2 422 Int_t nclusters = fClusters->GetEntries();
852a34a0 423 for (Int_t ic=0; ic<nclusters; ic++)
424 {
a4974ace 425 AliEMCALMatchCluster *cluster = (AliEMCALMatchCluster*)fClusters->At(ic);
852a34a0 426
427 Float_t clsPos[3] = {static_cast<Float_t>(cluster->X()),
428 static_cast<Float_t>(cluster->Y()),
429 static_cast<Float_t>(cluster->Z())};
430
a4974ace 431 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));
432 //printf("\n dR=%f,wind=%f\n",dR,fClusterWindow); //MARCEL
852a34a0 433
a4974ace 434 if (dR > fClusterWindow) continue;
8fc351e3 435
a4974ace 436 AliExternalTrackParam trkParTmp(trkParamTmp);
8fc351e3 437
a4974ace 438 Float_t tmpEta, tmpPhi;
439 if (!AliEMCALRecoUtils::ExtrapolateTrackToPosition(&trkParTmp, clsPos,track->GetMass(kTRUE), 5, tmpEta, tmpPhi)) continue;
852a34a0 440
441 if (TMath::Abs(tmpPhi)<TMath::Abs(maxPhi) && TMath::Abs(tmpEta)<TMath::Abs(maxEta))
442 {
a4974ace 443 maxPhi=tmpPhi;
444 maxEta=tmpEta;
445 index=ic;
446 }
447 }
30e29b2a 448
a4974ace 449 if (fITSTrackSA) delete trkParam;
852a34a0 450
5970dfe2 451 return index;
fe17d4cb 452}
5970dfe2 453
fe17d4cb 454//
455//------------------------------------------------------------------------------
456//
457void AliEMCALTracker::UnloadClusters()
458{
7ead3142 459 //
460 // Free memory from all arrays
461 // This method is called after the local tracking step
462 // so we can safely delete everything
463 //
fe17d4cb 464
7ead3142 465 Clear();
fe17d4cb 466}
04475328 467
468//
469//------------------------------------------------------------------------------
470//
7ead3142 471AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliEMCALRecPoint *recPoint) :
472 fIndex(index),
473 fX(0.),
474 fY(0.),
475 fZ(0.)
fe17d4cb 476{
7ead3142 477 //
478 // Translates an AliEMCALRecPoint object into the internal format.
479 // Index of passed cluster in its native array must be specified.
480 //
481 TVector3 clpos;
482 recPoint->GetGlobalPosition(clpos);
483
484 fX = clpos.X();
485 fY = clpos.Y();
486 fZ = clpos.Z();
fe17d4cb 487}
488//
489//------------------------------------------------------------------------------
490//
7ead3142 491AliEMCALTracker::AliEMCALMatchCluster::AliEMCALMatchCluster(Int_t index, AliESDCaloCluster *caloCluster) :
492 fIndex(index),
493 fX(0.),
494 fY(0.),
495 fZ(0.)
fe17d4cb 496{
7ead3142 497 //
498 // Translates an AliESDCaloCluster object into the internal format.
499 // Index of passed cluster in its native array must be specified.
500 //
501 Float_t clpos[3]= {0., 0., 0.};
502 caloCluster->GetPosition(clpos);
fe17d4cb 503
7ead3142 504 fX = (Double_t)clpos[0];
505 fY = (Double_t)clpos[1];
506 fZ = (Double_t)clpos[2];
fe17d4cb 507}