AliPhysics  5364b50 (5364b50)
AliCaloTrackMatcher.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: Daniel Mühlheim *
5  * Version 1.0 *
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 
17 //---------------------------------------------
18 // Basic Track Matching Class
19 //---------------------------------------------
21 
22 
23 #include "AliAnalysisManager.h"
24 #include "AliAODEvent.h"
25 #include "AliCaloTrackMatcher.h"
26 #include "AliEMCALRecoUtils.h"
27 #include "AliESDEvent.h"
28 #include "AliESDtrack.h"
29 #include "AliESDtrackCuts.h"
30 #include "AliInputEventHandler.h"
31 #include "AliTrackerBase.h"
32 #include "AliV0ReaderV1.h"
33 
34 #include "TAxis.h"
35 #include "TChain.h"
36 #include "TH1F.h"
37 #include "TF1.h"
38 
39 #include <vector>
40 #include <map>
41 #include <utility>
42 
43 class iostream;
44 
45 using namespace std;
46 
47 
48 ClassImp(AliCaloTrackMatcher)
49 
50 //________________________________________________________________________
51 AliCaloTrackMatcher::AliCaloTrackMatcher(const char *name, Int_t clusterType) : AliAnalysisTaskSE(name),
52  fClusterType(clusterType),
53  fV0ReaderName(""),
54  fCorrTaskSetting(""),
55  fAnalysisTrainMode("Grid"),
56  fMatchingWindow(200),
57  fMatchingResidual(0.2),
58  fRunNumber(-1),
59  fGeomEMCAL(NULL),
60  fGeomPHOS(NULL),
61  fMapTrackToCluster(),
62  fMapClusterToTrack(),
63  fNEntries(1),
64  fVectorDeltaEtaDeltaPhi(0),
65  fMap_TrID_ClID_ToIndex(),
66  fSecMapTrackToCluster(),
67  fSecMapClusterToTrack(),
68  fSecNEntries(1),
69  fSecVectorDeltaEtaDeltaPhi(0),
70  fSecMap_TrID_ClID_ToIndex(),
71  fSecMap_TrID_ClID_AlreadyTried(),
72  fListHistos(NULL),
73  fHistControlMatches(NULL),
74  fSecHistControlMatches(NULL)
75 {
76  // Default constructor
77  DefineInput(0, TChain::Class());
78 }
79 
80 //________________________________________________________________________
82  // default deconstructor
83  fMapTrackToCluster.clear();
84  fMapClusterToTrack.clear();
85  fVectorDeltaEtaDeltaPhi.clear();
86  fMap_TrID_ClID_ToIndex.clear();
87 
88  fSecMapTrackToCluster.clear();
89  fSecMapClusterToTrack.clear();
90  fSecVectorDeltaEtaDeltaPhi.clear();
91  fSecMap_TrID_ClID_ToIndex.clear();
92  fSecMap_TrID_ClID_AlreadyTried.clear();
93 
94  if(fHistControlMatches) delete fHistControlMatches;
95  if(fSecHistControlMatches) delete fSecHistControlMatches;
96  if(fAnalysisTrainMode.EqualTo("Grid")){
97  if(fListHistos != NULL){
98  delete fListHistos;
99  }
100  }
101 }
102 
103 //________________________________________________________________________
105  fMapTrackToCluster.clear();
106  fMapClusterToTrack.clear();
107  fVectorDeltaEtaDeltaPhi.clear();
108  fMap_TrID_ClID_ToIndex.clear();
109 
110  fSecMapTrackToCluster.clear();
111  fSecMapClusterToTrack.clear();
112  fSecVectorDeltaEtaDeltaPhi.clear();
113  fSecMap_TrID_ClID_ToIndex.clear();
114  fSecMap_TrID_ClID_AlreadyTried.clear();
115 }
116 
117 //________________________________________________________________________
119  if(fListHistos != NULL){
120  delete fListHistos;
121  fListHistos = NULL;
122  }
123  if(fListHistos == NULL){
124  fListHistos = new TList();
125  fListHistos->SetOwner(kTRUE);
126  fListHistos->SetName(Form("CaloTrackMatcher_%i",fClusterType));
127  }
128 
129  // Create User Output Objects
130  fHistControlMatches = new TH2F(Form("ControlMatches_%i",fClusterType),Form("ControlMatches_%i",fClusterType),7,-0.5,6.5,50,0.15,200);
131  SetLogBinningYTH2(fHistControlMatches);
132  fHistControlMatches->GetYaxis()->SetTitle("track pT (GeV/c)");
133  fHistControlMatches->GetXaxis()->SetBinLabel(1,"nTr in");
134  fHistControlMatches->GetXaxis()->SetBinLabel(2,"no inner Tr params || track not in right direction");
135  fHistControlMatches->GetXaxis()->SetBinLabel(3,"failed 1st pro-step");
136  fHistControlMatches->GetXaxis()->SetBinLabel(4,"Tr not in EMCal acc");
137  fHistControlMatches->GetXaxis()->SetBinLabel(5,"failed 2nd pro-step");
138  fHistControlMatches->GetXaxis()->SetBinLabel(6,"w/o match to cluster");
139  fHistControlMatches->GetXaxis()->SetBinLabel(7,"nTr out, w/ match");
140  fListHistos->Add(fHistControlMatches);
141 
142  fSecHistControlMatches = new TH2F(Form("ControlSecMatches_%i",fClusterType),Form("ControlSecMatches_%i",fClusterType),7,-0.5,6.5,50,0.15,200);
143  SetLogBinningYTH2(fSecHistControlMatches);
144  fSecHistControlMatches->GetYaxis()->SetTitle("track pT (GeV/c)");
145  fSecHistControlMatches->GetXaxis()->SetBinLabel(1,"nTr in");
146  fSecHistControlMatches->GetXaxis()->SetBinLabel(2,"no inner Tr params");
147  fSecHistControlMatches->GetXaxis()->SetBinLabel(3,"failed 1st pro-step");
148  fSecHistControlMatches->GetXaxis()->SetBinLabel(4,"Tr not in EMCal acc");
149  fSecHistControlMatches->GetXaxis()->SetBinLabel(5,"failed 2nd pro-step");
150  fSecHistControlMatches->GetXaxis()->SetBinLabel(6,"w/o match to cluster");
151  fSecHistControlMatches->GetXaxis()->SetBinLabel(7,"nTr out, w/ match");
152  fListHistos->Add(fSecHistControlMatches);
153 }
154 
155 //________________________________________________________________________
157  // Initialize function to be called once before analysis
158  fMapTrackToCluster.clear();
159  fMapClusterToTrack.clear();
160  fNEntries = 1;
161  fVectorDeltaEtaDeltaPhi.clear();
162  fMap_TrID_ClID_ToIndex.clear();
163 
164  fSecMapTrackToCluster.clear();
165  fSecMapClusterToTrack.clear();
166  fSecNEntries = 1;
167  fSecVectorDeltaEtaDeltaPhi.clear();
168  fSecMap_TrID_ClID_ToIndex.clear();
169  fSecMap_TrID_ClID_AlreadyTried.clear();
170 
171  if(fRunNumber == -1 || fRunNumber != runNumber){
172  if(fClusterType == 1 || fClusterType == 3){
173  fGeomEMCAL = AliEMCALGeometry::GetInstanceFromRunNumber(runNumber);
174  if(!fGeomEMCAL){ AliFatal("EMCal geometry not initialized!");}
175  }else if(fClusterType == 2){
176  fGeomPHOS = AliPHOSGeometry::GetInstance();
177  if(!fGeomPHOS) AliFatal("PHOS geometry not initialized!");
178  }
179  fRunNumber = runNumber;
180  }
181 }
182 
183 //________________________________________________________________________
185 
186  // main method of AliCaloTrackMatcher, first initialize and then process event
187 
188  //DebugV0Matching();
189 
190  // do processing only for EMCal (1), DCal (3) or PHOS (2) clusters, otherwise do nothing
191  if(fClusterType == 1 || fClusterType == 2 || fClusterType == 3){
192  Initialize(fInputEvent->GetRunNumber());
193  ProcessEvent(fInputEvent);
194  }
195 
196  //DebugMatching();
197 
198  return;
199 }
200 
201 //________________________________________________________________________
202 void AliCaloTrackMatcher::ProcessEvent(AliVEvent *event){
203  Int_t nClus = 0;
204  TClonesArray * arrClusters = NULL;
205  if(!fCorrTaskSetting.CompareTo("")){
206  nClus = event->GetNumberOfCaloClusters();
207  } else {
208  arrClusters = dynamic_cast<TClonesArray*>(event->FindListObject(Form("%sClustersBranch",fCorrTaskSetting.Data())));
209  nClus = arrClusters->GetEntries();
210  }
211  Int_t nModules = 0;
212  if(fClusterType == 1 || fClusterType == 3) nModules = fGeomEMCAL->GetNumberOfSuperModules();
213  else if(fClusterType == 2) nModules = fGeomPHOS->GetNModules();
214 
215  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
216  AliAODEvent *aodev = 0;
217  if (!esdev) {
218  aodev = dynamic_cast<AliAODEvent*>(event);
219  if (!aodev) {
220  AliError("Task needs AOD or ESD event, returning");
221  return;
222  }
223  }
224  static AliESDtrackCuts *EsdTrackCuts = 0x0;
225  static int prevRun = -1;
226  // Using standard function for setting Cuts
227  Int_t runNumber = fInputEvent->GetRunNumber();
228  if (esdev){
229  if (prevRun!=runNumber) {
230  delete EsdTrackCuts;
231  EsdTrackCuts = 0;
232  prevRun = runNumber;
233  }
234  if (!EsdTrackCuts) {
235  // if LHC11a or earlier or if LHC13g or if LHC12a-i -> use 2010 cuts
236  if( (runNumber<=146860) || (runNumber>=197470 && runNumber<=197692) || (runNumber>=172440 && runNumber<=193766) ){
237  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010();
238  // else if run2 data use 2015 PbPb cuts
239  } else if (runNumber>=209122){
240  //EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb();
241  // hard coded track cuts for the moment, because AliESDtrackCuts::GetStandardITSTPCTrackCuts2015PbPb() gives spams warnings
242  EsdTrackCuts = new AliESDtrackCuts();
243  // TPC; clusterCut = 1, cutAcceptanceEdges = kTRUE, removeDistortedRegions = kFALSE
244  EsdTrackCuts->AliESDtrackCuts::SetMinNCrossedRowsTPC(70);
245  EsdTrackCuts->AliESDtrackCuts::SetMinRatioCrossedRowsOverFindableClustersTPC(0.8);
246  EsdTrackCuts->SetCutGeoNcrNcl(2., 130., 1.5, 0.0, 0.0); // only dead zone and not clusters per length
247  //EsdTrackCuts->AliESDtrackCuts::SetCutOutDistortedRegionsTPC(kTRUE);
248  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterTPC(4);
249  EsdTrackCuts->AliESDtrackCuts::SetAcceptKinkDaughters(kFALSE);
250  EsdTrackCuts->AliESDtrackCuts::SetRequireTPCRefit(kTRUE);
251  // ITS; selPrimaries = 1
252  EsdTrackCuts->AliESDtrackCuts::SetRequireITSRefit(kTRUE);
253  EsdTrackCuts->AliESDtrackCuts::SetClusterRequirementITS(AliESDtrackCuts::kSPD,
254  AliESDtrackCuts::kAny);
255  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexXYPtDep("0.0105+0.0350/pt^1.1");
256  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2TPCConstrainedGlobal(36);
257  EsdTrackCuts->AliESDtrackCuts::SetMaxDCAToVertexZ(2);
258  EsdTrackCuts->AliESDtrackCuts::SetDCAToVertex2D(kFALSE);
259  EsdTrackCuts->AliESDtrackCuts::SetRequireSigmaToVertex(kFALSE);
260  EsdTrackCuts->AliESDtrackCuts::SetMaxChi2PerClusterITS(36);
261  // else use 2011 version of track cuts
262  }else{
263  EsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011();
264  }
265  EsdTrackCuts->SetMaxDCAToVertexZ(2);
266  }
267  }
268 
269  for (Int_t itr=0;itr<event->GetNumberOfTracks();itr++){
270  AliExternalTrackParam *trackParam = 0;
271  AliVTrack *inTrack = 0x0;
272  if(esdev){
273  inTrack = esdev->GetTrack(itr);
274  if(!inTrack) continue;
275  if(TMath::Abs(inTrack->Eta())>0.8 && (fClusterType==1 || fClusterType==3 )) continue;
276  if(TMath::Abs(inTrack->Eta())>0.3 && fClusterType==2 ) continue;
277  if(inTrack->Pt()<0.5) continue;
278  AliESDtrack *esdt = dynamic_cast<AliESDtrack*>(inTrack);
279  if(!EsdTrackCuts->AcceptTrack(esdt)) continue;
280  fHistControlMatches->Fill(0.,inTrack->Pt());
281 
282  const AliExternalTrackParam *in = esdt->GetInnerParam();
283  if (!in){AliDebug(2, "Could not get InnerParam of Track, continue"); fHistControlMatches->Fill(1.,inTrack->Pt()); continue;}
284  trackParam = new AliExternalTrackParam(*in);
285  } else if(aodev) {
286  inTrack = dynamic_cast<AliVTrack*>(aodev->GetTrack(itr));
287  if(!inTrack) continue;
288  if(inTrack->GetID()<0) continue; // Avoid double counting of tracks
289 
290  if(TMath::Abs(inTrack->Eta())>0.8 && (fClusterType==1 || fClusterType==3 )) continue;
291  if(TMath::Abs(inTrack->Eta())>0.3 && fClusterType==2 ) continue;
292  if(inTrack->Pt()<0.5) continue;
293 
294 
295  AliAODTrack *aodt = dynamic_cast<AliAODTrack*>(inTrack);
296  if(!aodt->IsHybridGlobalConstrainedGlobal()) continue;
297  Double_t xyz[3] = {0}, pxpypz[3] = {0}, cv[21] = {0};
298  aodt->GetPxPyPz(pxpypz);
299  aodt->GetXYZ(xyz);
300  aodt->GetCovarianceXYZPxPyPz(cv);
301  fHistControlMatches->Fill(0.,inTrack->Pt());
302 
303  // check for EMC tracks already propagated tracks are out of bounds
304  if (fClusterType == 1 || fClusterType == 3){
305  if( TMath::Abs(aodt->GetTrackEtaOnEMCal()) > 0.75 )
306  continue;
307 
308  // conditions for run1
309  if( fClusterType == 1 && nModules < 13 && ( aodt->GetTrackPhiOnEMCal() < 70*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 190*TMath::DegToRad()))
310  continue;
311 
312  // conditions for run2
313  if( nModules > 12 ){
314  if (fClusterType == 3 && ( aodt->GetTrackPhiOnEMCal() < 250*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 340*TMath::DegToRad()))
315  continue;
316  if( fClusterType == 1 && ( aodt->GetTrackPhiOnEMCal() < 70*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 190*TMath::DegToRad()))
317  continue;
318  }
319  }
320  trackParam = new AliExternalTrackParam(xyz,pxpypz,cv,aodt->Charge());
321  }
322 
323  if (!trackParam) {AliError("Could not get TrackParameters, continue"); fHistControlMatches->Fill(1.,inTrack->Pt()); continue;}
324  AliExternalTrackParam emcParam(*trackParam);
325  Float_t eta, phi, pt;
326 
327  //propagate tracks to emc surfaces
328  if(fClusterType == 1 || fClusterType == 3){
329  if (!AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&emcParam, 440., 0.139, 20., eta, phi, pt)) {
330  delete trackParam;
331  fHistControlMatches->Fill(2.,inTrack->Pt());
332  continue;
333  }
334  if( TMath::Abs(eta) > 0.75 ) {
335  delete trackParam;
336  fHistControlMatches->Fill(3.,inTrack->Pt());
337  continue;
338  }
339  // Save some time and memory in case of no DCal present
340  if( fClusterType == 1 && nModules < 13 && ( phi < 70*TMath::DegToRad() || phi > 190*TMath::DegToRad())){
341  delete trackParam;
342  fHistControlMatches->Fill(3.,inTrack->Pt());
343  continue;
344  }
345  // Save some time and memory in case of run2
346  if( nModules > 12 ){
347  if (fClusterType == 3 && ( phi < 250*TMath::DegToRad() || phi > 340*TMath::DegToRad())){
348  delete trackParam;
349  fHistControlMatches->Fill(3.,inTrack->Pt());
350  continue;
351  }
352  if( fClusterType == 1 && ( phi < 70*TMath::DegToRad() || phi > 190*TMath::DegToRad())){
353  delete trackParam;
354  fHistControlMatches->Fill(3.,inTrack->Pt());
355  continue;
356  }
357  }
358  }else if(fClusterType == 2){
359  if( !AliTrackerBase::PropagateTrackToBxByBz(&emcParam, 460., 0.139, 20, kTRUE, 0.8, -1)){
360  delete trackParam;
361  fHistControlMatches->Fill(2.,inTrack->Pt());
362  continue;
363  }
364 
365  if (TMath::Abs(eta) > 0.25 && (fClusterType == 2)){
366  delete trackParam;
367  fHistControlMatches->Fill(3.,inTrack->Pt());
368  continue;
369  }
370 
371 //to do: implement of distance checks
372  }
373 
374  Float_t dEta=-999, dPhi=-999;
375  Float_t clsPos[3] = {0.,0.,0.};
376  Double_t exPos[3] = {0.,0.,0.};
377  if (!emcParam.GetXYZ(exPos)){
378  fHistControlMatches->Fill(2.,inTrack->Pt());
379  delete trackParam;
380  continue;
381  }
382 
383 // cout << inTrack->GetID() << " - " << trackParam << endl;
384 // cout << "eta/phi: " << eta << ", " << phi << endl;
385 // cout << "nClus: " << nClus << endl;
386  Int_t nClusterMatchesToTrack = 0;
387  for(Int_t iclus=0;iclus < nClus;iclus++){
388  AliVCluster* cluster = NULL;
389  if(arrClusters){
390  if(esdev){
391  if(arrClusters)
392  cluster = new AliESDCaloCluster(*(AliESDCaloCluster*)arrClusters->At(iclus));
393  } else if(aodev){
394  if(arrClusters)
395  cluster = new AliAODCaloCluster(*(AliAODCaloCluster*)arrClusters->At(iclus));
396  }
397  }
398  else
399  cluster = event->GetCaloCluster(iclus);
400  if (!cluster){
401  continue;
402  }
403 // cout << "-------------------------LOOPING: " << iclus << ", " << cluster->GetID() << endl;
404  cluster->GetPosition(clsPos);
405  Double_t dR = TMath::Sqrt(TMath::Power(exPos[0]-clsPos[0],2)+TMath::Power(exPos[1]-clsPos[1],2)+TMath::Power(exPos[2]-clsPos[2],2));
406 //cout << "dR: " << dR << endl;
407  if (dR > fMatchingWindow){
408  continue;
409  }
410  Double_t clusterR = TMath::Sqrt( clsPos[0]*clsPos[0] + clsPos[1]*clsPos[1] );
411  AliExternalTrackParam trackParamTmp(emcParam);//Retrieve the starting point every time before the extrapolation
412  if(fClusterType == 1 || fClusterType == 3){
413  if (!cluster->IsEMCAL()){
414  continue;
415  }
416  if(!AliEMCALRecoUtils::ExtrapolateTrackToCluster(&trackParamTmp, cluster, 0.139, 5., dEta, dPhi)){
417  fHistControlMatches->Fill(4.,inTrack->Pt());
418  continue;
419  }
420  }else if(fClusterType == 2){
421  if (!cluster->IsPHOS()){
422  continue;
423  }
424  if(!AliTrackerBase::PropagateTrackToBxByBz(&trackParamTmp, clusterR, 0.139, 5., kTRUE, 0.8, -1)){
425  fHistControlMatches->Fill(4.,inTrack->Pt());
426  continue;
427  }
428  Double_t trkPos[3] = {0,0,0};
429  trackParamTmp.GetXYZ(trkPos);
430  TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
431  TVector3 clsPosVec(clsPos);
432  dPhi = clsPosVec.DeltaPhi(trkPosVec);
433  dEta = clsPosVec.Eta()-trkPosVec.Eta();
434  }
435 
436  Float_t dR2 = dPhi*dPhi + dEta*dEta;
437 
438 //cout << dEta << " - " << dPhi << " - " << dR2 << endl;
439  if(dR2 > fMatchingResidual){
440  continue;
441  }
442  nClusterMatchesToTrack++;
443  if(aodev){
444  fMapTrackToCluster.insert(make_pair(itr,cluster->GetID()));
445  fMapClusterToTrack.insert(make_pair(cluster->GetID(),itr));
446  }else{
447  fMapTrackToCluster.insert(make_pair(inTrack->GetID(),cluster->GetID()));
448  fMapClusterToTrack.insert(make_pair(cluster->GetID(),inTrack->GetID()));
449  }
450  fVectorDeltaEtaDeltaPhi.push_back(make_pair(dEta,dPhi));
451  fMap_TrID_ClID_ToIndex[make_pair(inTrack->GetID(),cluster->GetID())] = fNEntries++;
452  if( (Int_t)fVectorDeltaEtaDeltaPhi.size() != (fNEntries-1)) AliFatal("Fatal error in AliCaloTrackMatcher, vector and map are not in sync!");
453  }
454  if(nClusterMatchesToTrack == 0) fHistControlMatches->Fill(5.,inTrack->Pt());
455  else fHistControlMatches->Fill(6.,inTrack->Pt());
456  delete trackParam;
457  }
458 
459  return;
460 }
461 
462 //________________________________________________________________________
463 Bool_t AliCaloTrackMatcher::PropagateV0TrackToClusterAndGetMatchingResidual(AliVTrack* inSecTrack, AliVCluster* cluster, AliVEvent* event, Float_t &dEta, Float_t &dPhi){
464 
465  //if V0-track to cluster match is already available return stored residuals
466  if(GetSecTrackClusterMatchingResidual(inSecTrack->GetID(),cluster->GetID(), dEta, dPhi)){
467 //cout << "RESIDUAL ALREADY AVAILABLE! - " << dEta << "/" << dPhi << endl;
468  return kTRUE;
469  }
470 
471  if(IsSecTrackClusterAlreadyTried(inSecTrack->GetID(),cluster->GetID())){
472 //cout << "PROPAGATION ALREADY FAILED! - " << inSecTrack->GetID() << "/" << cluster->GetID() << endl;
473  return kFALSE;
474  }
475 
476 //cout << "running matching! - " << inSecTrack->GetID() << "/" << cluster->GetID() << endl;
477  //if match has not yet been computed, go on:
478  Int_t nModules = 0;
479  if(fClusterType == 1 || fClusterType == 3) nModules = fGeomEMCAL->GetNumberOfSuperModules();
480  else if(fClusterType == 2) nModules = fGeomPHOS->GetNModules();
481 
482  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(event);
483  AliAODEvent *aodev = 0;
484  if (!esdev) {
485  aodev = dynamic_cast<AliAODEvent*>(event);
486  if (!aodev) {
487  AliError("Task needs AOD or ESD event, returning");
488  return kFALSE;
489  }
490  }
491 
492  if(!cluster->IsEMCAL() && !cluster->IsPHOS()){AliError("Cluster is neither EMCAL nor PHOS, returning");return kFALSE;}
493 
494  Float_t clusterPosition[3] = {0,0,0};
495  cluster->GetPosition(clusterPosition);
496  Double_t clusterR = TMath::Sqrt( clusterPosition[0]*clusterPosition[0] + clusterPosition[1]*clusterPosition[1] );
497 
498  if(!inSecTrack) return kFALSE;
499  fSecHistControlMatches->Fill(0.,inSecTrack->Pt());
500 
501  AliESDtrack *esdt = dynamic_cast<AliESDtrack*>(inSecTrack);
502  AliAODTrack *aodt = 0;
503  if (!esdt) {
504  aodt = dynamic_cast<AliAODTrack*>(inSecTrack);
505  if (!aodt){
506  AliError("Track is neither ESD nor AOD, continue");
507  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
508  return kFALSE;
509  }
510  }
511 
512  AliExternalTrackParam *trackParam = 0;
513  if (esdt) {
514  const AliExternalTrackParam *in = esdt->GetInnerParam();
515  if (!in){
516  AliDebug(2, "Could not get InnerParam of Track, continue");
517  fSecHistControlMatches->Fill(1.,inSecTrack->Pt());
518  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
519  return kFALSE;
520  }
521  trackParam = new AliExternalTrackParam(*in);
522  } else {
523  // check if tracks should be propagated at all
524  if (fClusterType == 1 || fClusterType == 3){
525  if (TMath::Abs(aodt->GetTrackEtaOnEMCal()) > 0.8)
526  return kFALSE;
527  if( nModules < 13 ){
528  if (( aodt->GetTrackPhiOnEMCal() < 60*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 200*TMath::DegToRad()))
529  return kFALSE;
530  } else if( nModules > 12 ){
531  if (fClusterType == 3 && ( aodt->GetTrackPhiOnEMCal() < 250*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 340*TMath::DegToRad()))
532  return kFALSE;
533  if( fClusterType == 1 && ( aodt->GetTrackPhiOnEMCal() < 60*TMath::DegToRad() || aodt->GetTrackPhiOnEMCal() > 200*TMath::DegToRad()))
534  return kFALSE;
535  }
536  } else {
537  if ( aodt->Phi() < 60*TMath::DegToRad() || aodt->Phi() > 200*TMath::DegToRad())
538  return kFALSE;
539  if (TMath::Abs(aodt->Eta()) > 0.3 )
540  return kFALSE;
541  }
542 
543  Double_t xyz[3] = {0}, pxpypz[3] = {0}, cv[21] = {0};
544  aodt->GetPxPyPz(pxpypz);
545  aodt->GetXYZ(xyz);
546  aodt->GetCovarianceXYZPxPyPz(cv);
547  trackParam = new AliExternalTrackParam(xyz,pxpypz,cv,aodt->Charge());
548  }
549  if(!trackParam){
550  AliError("Could not get TrackParameters, continue");
551  fSecHistControlMatches->Fill(1.,inSecTrack->Pt());
552  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
553  return kFALSE;
554  }
555 
556  Bool_t propagated = kFALSE;
557  AliExternalTrackParam emcParam(*trackParam);
558  Float_t dPhiTemp = 0;
559  Float_t dEtaTemp = 0;
560 
561  if(cluster->IsEMCAL()){
562  if (inSecTrack->Pt() < 0.3 ) {
563  fSecHistControlMatches->Fill(3.,inSecTrack->Pt());
564  }
565  Float_t eta = 0;Float_t phi = 0;Float_t pt = 0;
566  propagated = AliEMCALRecoUtils::ExtrapolateTrackToEMCalSurface(&emcParam, 430, 0.000510999, 20, eta, phi, pt);
567  if(propagated){
568  if( TMath::Abs(eta) > 0.8 ) {
569  delete trackParam;
570  fSecHistControlMatches->Fill(3.,inSecTrack->Pt());
571  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
572  return kFALSE;
573  }
574  // Save some time and memory in case of no DCal present
575  if( nModules < 13 && ( phi < 60*TMath::DegToRad() || phi > 200*TMath::DegToRad())){
576  delete trackParam;
577  fSecHistControlMatches->Fill(3.,inSecTrack->Pt());
578  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
579  return kFALSE;
580  }
581 
582  propagated = AliEMCALRecoUtils::ExtrapolateTrackToCluster(&emcParam, cluster, 0.000510999, 5, dEtaTemp, dPhiTemp);
583  if(!propagated){
584  delete trackParam;
585  fSecHistControlMatches->Fill(4.,inSecTrack->Pt());
586  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
587  return kFALSE;
588  }
589  }else{
590  delete trackParam;
591  fSecHistControlMatches->Fill(2.,inSecTrack->Pt());
592  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
593  return kFALSE;
594  }
595 
596  }else if(cluster->IsPHOS()){
597  propagated = AliTrackerBase::PropagateTrackToBxByBz(&emcParam, clusterR, 0.000510999, 20, kTRUE, 0.8, -1);
598  if (propagated){
599  Double_t trkPos[3] = {0,0,0};
600  emcParam.GetXYZ(trkPos);
601  TVector3 trkPosVec(trkPos[0],trkPos[1],trkPos[2]);
602  TVector3 clsPosVec(clusterPosition);
603  dPhiTemp = clsPosVec.DeltaPhi(trkPosVec);
604  dEtaTemp = clsPosVec.Eta()-trkPosVec.Eta();
605  }else{
606  delete trackParam;
607  fSecHistControlMatches->Fill(2.,inSecTrack->Pt());
608  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
609  return kFALSE;}
610  }
611 
612  if (propagated){
613  Float_t dR2 = dPhiTemp*dPhiTemp + dEtaTemp*dEtaTemp;
614 //cout << dEtaTemp << " - " << dPhiTemp << " - " << dR2 << endl;
615  if(dR2 > fMatchingResidual){
616  fSecHistControlMatches->Fill(5.,inSecTrack->Pt());
617  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
618 //cout << "NO MATCH! - " << inSecTrack->GetID() << "/" << cluster->GetID() << endl;
619  delete trackParam;
620  return kFALSE;
621  }
622 //cout << "MATCHED!!!!!!!" << endl;
623 
624  if(aodev){
625  //need to search for position in case of AOD
626  Int_t TrackPos = -1;
627  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
628  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
629  if(currTrack->GetID() == inSecTrack->GetID()){
630  TrackPos = iTrack;
631  break;
632  }
633  }
634  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: PropagateV0TrackToClusterAndGetMatchingResidual - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",inSecTrack->GetID()));
635  fSecMapTrackToCluster.insert(make_pair(TrackPos,cluster->GetID()));
636  fSecMapClusterToTrack.insert(make_pair(cluster->GetID(),TrackPos));
637  }else{
638  fSecMapTrackToCluster.insert(make_pair(inSecTrack->GetID(),cluster->GetID()));
639  fSecMapClusterToTrack.insert(make_pair(cluster->GetID(),inSecTrack->GetID()));
640  }
641  fSecVectorDeltaEtaDeltaPhi.push_back(make_pair(dEtaTemp,dPhiTemp));
642  fSecMap_TrID_ClID_ToIndex[make_pair(inSecTrack->GetID(),cluster->GetID())] = fSecNEntries++;
643  if( (Int_t)fSecVectorDeltaEtaDeltaPhi.size() != (fSecNEntries-1)) AliFatal("Fatal error in AliCaloTrackMatcher, vector and map are not in sync!");
644 
645  fSecHistControlMatches->Fill(6.,inSecTrack->Pt());
646  dEta = dEtaTemp;
647  dPhi = dPhiTemp;
648  delete trackParam;
649  return kTRUE;
650  }else AliFatal("Fatal error in AliCaloTrackMatcher, track is labeled as sucessfully propagated although this should be impossible!");
651 
652  fSecMap_TrID_ClID_AlreadyTried[make_pair(inSecTrack->GetID(),cluster->GetID())] = 1.;
653  delete trackParam;
654  return kFALSE;
655 }
656 
657 //________________________________________________________________________
658 //________________________________________________________________________
659 //________________________________________________________________________
660 //________________________________________________________________________
661 //________________________________________________________________________
663  Int_t position = fMap_TrID_ClID_ToIndex[make_pair(trackID,clusterID)];
664  if(position == 0) return kFALSE;
665 
666  pairFloat tempEtaPhi = fVectorDeltaEtaDeltaPhi.at(position-1);
667  dEta = tempEtaPhi.first;
668  dPhi = tempEtaPhi.second;
669  return kTRUE;
670 }
671 //________________________________________________________________________
672 Int_t AliCaloTrackMatcher::GetNMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
673  Int_t matched = 0;
674  multimap<Int_t,Int_t>::iterator it;
675  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
676  if(it->first == clusterID){
677  Float_t tempDEta, tempDPhi;
678  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
679  if(!tempTrack) continue;
680  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
681  if(tempTrack->Charge()>0){
682  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) matched++;
683  }else if(tempTrack->Charge()<0){
684  dPhiMin*=-1;
685  dPhiMax*=-1;
686  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) matched++;
687  }
688  }
689  }
690  }
691 
692  return matched;
693 }
694 
695 //________________________________________________________________________
696 Int_t AliCaloTrackMatcher::GetNMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
697  Int_t matched = 0;
698  multimap<Int_t,Int_t>::iterator it;
699  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
700  if(it->first == clusterID){
701  Float_t tempDEta, tempDPhi;
702  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
703  if(!tempTrack) continue;
704  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
705  Bool_t match_dEta = kFALSE;
706  Bool_t match_dPhi = kFALSE;
707  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
708  else match_dEta = kFALSE;
709 
710  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
711  else match_dPhi = kFALSE;
712 
713  if (match_dPhi && match_dEta )matched++;
714  }
715  }
716  }
717  return matched;
718 }
719 
720 //________________________________________________________________________
722  Int_t matched = 0;
723  multimap<Int_t,Int_t>::iterator it;
724  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
725  if(it->first == clusterID){
726  Float_t tempDEta, tempDPhi;
727  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
728  if(!tempTrack) continue;
729  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
730  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) matched++;
731  }
732  }
733  }
734  return matched;
735 }
736 
737 //________________________________________________________________________
738 Int_t AliCaloTrackMatcher::GetNMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
739 
740  Int_t TrackPos = -1;
741  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
742  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
743  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
744  if(currTrack->GetID() == trackID){
745  TrackPos = iTrack;
746  break;
747  }
748  }
749  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
750  }else TrackPos = trackID; // for ESD just take trackID
751 
752  Int_t matched = 0;
753  multimap<Int_t,Int_t>::iterator it;
754  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
755  if(!tempTrack) return matched;
756  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
757  if(it->first == TrackPos){
758  Float_t tempDEta, tempDPhi;
759  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
760  if(tempTrack->Charge()>0){
761  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) matched++;
762  }else if(tempTrack->Charge()<0){
763  dPhiMin*=-1;
764  dPhiMax*=-1;
765  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) matched++;
766  }
767  }
768  }
769  }
770  return matched;
771 }
772 
773 //________________________________________________________________________
774 Int_t AliCaloTrackMatcher::GetNMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
775  Int_t TrackPos = -1;
776  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
777  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
778  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
779  if(currTrack->GetID() == trackID){
780  TrackPos = iTrack;
781  break;
782  }
783  }
784  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
785  }else TrackPos = trackID; // for ESD just take trackID
786 
787  Int_t matched = 0;
788  multimap<Int_t,Int_t>::iterator it;
789  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
790  if(!tempTrack) return matched;
791  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
792  if(it->first == TrackPos){
793  Float_t tempDEta, tempDPhi;
794  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
795  Bool_t match_dEta = kFALSE;
796  Bool_t match_dPhi = kFALSE;
797  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
798  else match_dEta = kFALSE;
799 
800  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
801  else match_dPhi = kFALSE;
802 
803  if (match_dPhi && match_dEta )matched++;
804 
805  }
806  }
807  }
808  return matched;
809 }
810 
811 //________________________________________________________________________
813  Int_t TrackPos = -1;
814  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
815  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
816  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
817  if(currTrack->GetID() == trackID){
818  TrackPos = iTrack;
819  break;
820  }
821  }
822  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
823  }else TrackPos = trackID; // for ESD just take trackID
824 
825  Int_t matched = 0;
826  multimap<Int_t,Int_t>::iterator it;
827  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
828  if(!tempTrack) return matched;
829  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
830  if(it->first == TrackPos){
831  Float_t tempDEta, tempDPhi;
832  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
833  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) matched++;
834  }
835  }
836  }
837  return matched;
838 }
839 
840 //________________________________________________________________________
841 vector<Int_t> AliCaloTrackMatcher::GetMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
842  vector<Int_t> tempMatchedTracks;
843  multimap<Int_t,Int_t>::iterator it;
844  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
845  if(it->first == clusterID){
846  Float_t tempDEta, tempDPhi;
847  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
848  if(!tempTrack) continue;
849  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
850  if(tempTrack->Charge()>0){
851  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) tempMatchedTracks.push_back(it->second);
852  }else if(tempTrack->Charge()<0){
853  dPhiMin*=-1;
854  dPhiMax*=-1;
855  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) tempMatchedTracks.push_back(it->second);
856  }
857  }
858  }
859  }
860  return tempMatchedTracks;
861 }
862 
863 //________________________________________________________________________
864 vector<Int_t> AliCaloTrackMatcher::GetMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
865  vector<Int_t> tempMatchedTracks;
866  multimap<Int_t,Int_t>::iterator it;
867  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
868  if(it->first == clusterID){
869  Float_t tempDEta, tempDPhi;
870  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
871  if(!tempTrack) continue;
872  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
873  Bool_t match_dEta = kFALSE;
874  Bool_t match_dPhi = kFALSE;
875  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
876  else match_dEta = kFALSE;
877 
878  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
879  else match_dPhi = kFALSE;
880 
881  if (match_dPhi && match_dEta )tempMatchedTracks.push_back(it->second);
882 
883  }
884  }
885  }
886  return tempMatchedTracks;
887 }
888 
889 //________________________________________________________________________
890 vector<Int_t> AliCaloTrackMatcher::GetMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dR){
891  vector<Int_t> tempMatchedTracks;
892  multimap<Int_t,Int_t>::iterator it;
893  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it){
894  if(it->first == clusterID){
895  Float_t tempDEta, tempDPhi;
896  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
897  if(!tempTrack) continue;
898  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
899  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) tempMatchedTracks.push_back(it->second);
900  }
901  }
902  }
903  return tempMatchedTracks;
904 }
905 
906 //________________________________________________________________________
907 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
908  Int_t TrackPos = -1;
909  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
910  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
911  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
912  if(currTrack->GetID() == trackID){
913  TrackPos = iTrack;
914  break;
915  }
916  }
917  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
918  }else TrackPos = trackID; // for ESD just take trackID
919 
920  vector<Int_t> tempMatchedClusters;
921  multimap<Int_t,Int_t>::iterator it;
922  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
923  if(!tempTrack) return tempMatchedClusters;
924  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
925  if(it->first == TrackPos){
926  Float_t tempDEta, tempDPhi;
927  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
928  if(tempTrack->Charge()>0){
929  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) tempMatchedClusters.push_back(it->second);
930  }else if(tempTrack->Charge()<0){
931  dPhiMin*=-1;
932  dPhiMax*=-1;
933  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) tempMatchedClusters.push_back(it->second);
934  }
935  }
936  }
937  }
938 
939  return tempMatchedClusters;
940 }
941 
942 //________________________________________________________________________
943 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
944  Int_t TrackPos = -1;
945  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
946  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
947  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
948  if(currTrack->GetID() == trackID){
949  TrackPos = iTrack;
950  break;
951  }
952  }
953  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
954  }else TrackPos = trackID; // for ESD just take trackID
955 
956  vector<Int_t> tempMatchedClusters;
957  multimap<Int_t,Int_t>::iterator it;
958  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
959  if(!tempTrack) return tempMatchedClusters;
960  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
961  if(it->first == TrackPos){
962  Float_t tempDEta, tempDPhi;
963  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
964  Bool_t match_dEta = kFALSE;
965  Bool_t match_dPhi = kFALSE;
966  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
967  else match_dEta = kFALSE;
968 
969  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
970  else match_dPhi = kFALSE;
971 
972  if (match_dPhi && match_dEta )tempMatchedClusters.push_back(it->second);
973  }
974  }
975  }
976  return tempMatchedClusters;
977 }
978 
979 //________________________________________________________________________
980 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, Float_t dR){
981  Int_t TrackPos = -1;
982  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
983  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
984  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
985  if(currTrack->GetID() == trackID){
986  TrackPos = iTrack;
987  break;
988  }
989  }
990  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
991  }else TrackPos = trackID; // for ESD just take trackID
992 
993  vector<Int_t> tempMatchedClusters;
994  multimap<Int_t,Int_t>::iterator it;
995  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
996  if(!tempTrack) return tempMatchedClusters;
997  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it){
998  if(it->first == TrackPos){
999  Float_t tempDEta, tempDPhi;
1000  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1001  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) tempMatchedClusters.push_back(it->second);
1002  }
1003  }
1004  }
1005  return tempMatchedClusters;
1006 }
1007 
1008 //________________________________________________________________________
1009 //________________________________________________________________________
1010 //________________________________________________________________________
1011 //________________________________________________________________________
1012 //________________________________________________________________________
1014  Int_t position = fSecMap_TrID_ClID_ToIndex[make_pair(trackID,clusterID)];
1015  if(position == 0) return kFALSE;
1016 
1017  pairFloat tempEtaPhi = fSecVectorDeltaEtaDeltaPhi.at(position-1);
1018  dEta = tempEtaPhi.first;
1019  dPhi = tempEtaPhi.second;
1020  return kTRUE;
1021 }
1022 //________________________________________________________________________
1024  Int_t position = fSecMap_TrID_ClID_AlreadyTried[make_pair(trackID,clusterID)];
1025  if(position == 0) return kFALSE;
1026  else return kTRUE;
1027 }
1028 //________________________________________________________________________
1029 Int_t AliCaloTrackMatcher::GetNMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
1030  Int_t matched = 0;
1031  multimap<Int_t,Int_t>::iterator it;
1032  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1033  if(it->first == clusterID){
1034  Float_t tempDEta, tempDPhi;
1035  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1036  if(!tempTrack) continue;
1037  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1038  if(tempTrack->Charge()>0){
1039  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) matched++;
1040  }else if(tempTrack->Charge()<0){
1041  dPhiMin*=-1;
1042  dPhiMax*=-1;
1043  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) matched++;
1044  }
1045  }
1046  }
1047  }
1048 
1049  return matched;
1050 }
1051 
1052 //________________________________________________________________________
1053 Int_t AliCaloTrackMatcher::GetNMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
1054  Int_t matched = 0;
1055  multimap<Int_t,Int_t>::iterator it;
1056  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1057  if(it->first == clusterID){
1058  Float_t tempDEta, tempDPhi;
1059  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1060  if(!tempTrack) continue;
1061  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1062  Bool_t match_dEta = kFALSE;
1063  Bool_t match_dPhi = kFALSE;
1064  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
1065  else match_dEta = kFALSE;
1066 
1067  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
1068  else match_dPhi = kFALSE;
1069 
1070  if (match_dPhi && match_dEta )matched++;
1071  }
1072  }
1073  }
1074 
1075  return matched;
1076 }
1077 
1078 //________________________________________________________________________
1080  Int_t matched = 0;
1081  multimap<Int_t,Int_t>::iterator it;
1082  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1083  if(it->first == clusterID){
1084  Float_t tempDEta, tempDPhi;
1085  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1086  if(!tempTrack) continue;
1087  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1088  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) matched++;
1089  }
1090  }
1091  }
1092 
1093  return matched;
1094 }
1095 
1096 //________________________________________________________________________
1097 Int_t AliCaloTrackMatcher::GetNMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
1098  Int_t TrackPos = -1;
1099  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1100  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1101  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1102  if(currTrack->GetID() == trackID){
1103  TrackPos = iTrack;
1104  break;
1105  }
1106  }
1107  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1108  }else TrackPos = trackID; // for ESD just take trackID
1109 
1110  Int_t matched = 0;
1111  multimap<Int_t,Int_t>::iterator it;
1112  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1113  if(!tempTrack) return matched;
1114  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1115  if(it->first == TrackPos){
1116  Float_t tempDEta, tempDPhi;
1117  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1118  if(tempTrack->Charge()>0){
1119  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) matched++;
1120  }else if(tempTrack->Charge()<0){
1121  dPhiMin*=-1;
1122  dPhiMax*=-1;
1123  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) matched++;
1124  }
1125  }
1126  }
1127  }
1128 
1129  return matched;
1130 }
1131 
1132 //________________________________________________________________________
1133 Int_t AliCaloTrackMatcher::GetNMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
1134  Int_t TrackPos = -1;
1135  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1136  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1137  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1138  if(currTrack->GetID() == trackID){
1139  TrackPos = iTrack;
1140  break;
1141  }
1142  }
1143  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1144  }else TrackPos = trackID; // for ESD just take trackID
1145 
1146  Int_t matched = 0;
1147  multimap<Int_t,Int_t>::iterator it;
1148  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1149  if(!tempTrack) return matched;
1150  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1151  if(it->first == TrackPos){
1152  Float_t tempDEta, tempDPhi;
1153  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1154  Bool_t match_dEta = kFALSE;
1155  Bool_t match_dPhi = kFALSE;
1156  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
1157  else match_dEta = kFALSE;
1158 
1159  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
1160  else match_dPhi = kFALSE;
1161 
1162  if (match_dPhi && match_dEta )matched++;
1163 
1164  }
1165  }
1166  }
1167 
1168  return matched;
1169 }
1170 
1171 //________________________________________________________________________
1173  Int_t TrackPos = -1;
1174  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1175  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1176  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1177  if(currTrack->GetID() == trackID){
1178  TrackPos = iTrack;
1179  break;
1180  }
1181  }
1182  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1183  }else TrackPos = trackID; // for ESD just take trackID
1184 
1185  Int_t matched = 0;
1186  multimap<Int_t,Int_t>::iterator it;
1187  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1188  if(!tempTrack) return matched;
1189  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1190  if(it->first == TrackPos){
1191  Float_t tempDEta, tempDPhi;
1192  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1193  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) matched++;
1194  }
1195  }
1196  }
1197 
1198  return matched;
1199 }
1200 
1201 //________________________________________________________________________
1202 vector<Int_t> AliCaloTrackMatcher::GetMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
1203  vector<Int_t> tempMatchedTracks;
1204  multimap<Int_t,Int_t>::iterator it;
1205  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1206  if(it->first == clusterID){
1207  Float_t tempDEta, tempDPhi;
1208  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1209  if(!tempTrack) continue;
1210  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1211  if(tempTrack->Charge()>0){
1212  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) tempMatchedTracks.push_back(it->second);
1213  }else if(tempTrack->Charge()<0){
1214  dPhiMin*=-1;
1215  dPhiMax*=-1;
1216  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) tempMatchedTracks.push_back(it->second);
1217  }
1218  }
1219  }
1220  }
1221 
1222  return tempMatchedTracks;
1223 }
1224 
1225 //________________________________________________________________________
1226 vector<Int_t> AliCaloTrackMatcher::GetMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
1227  vector<Int_t> tempMatchedTracks;
1228  multimap<Int_t,Int_t>::iterator it;
1229  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1230  if(it->first == clusterID){
1231  Float_t tempDEta, tempDPhi;
1232  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1233  if(!tempTrack) continue;
1234  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1235  Bool_t match_dEta = kFALSE;
1236  Bool_t match_dPhi = kFALSE;
1237  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
1238  else match_dEta = kFALSE;
1239 
1240  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
1241  else match_dPhi = kFALSE;
1242 
1243  if (match_dPhi && match_dEta )tempMatchedTracks.push_back(it->second);
1244  }
1245  }
1246  }
1247 
1248  return tempMatchedTracks;
1249 }
1250 
1251 //________________________________________________________________________
1252 vector<Int_t> AliCaloTrackMatcher::GetMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dR){
1253  vector<Int_t> tempMatchedTracks;
1254  multimap<Int_t,Int_t>::iterator it;
1255  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it){
1256  if(it->first == clusterID){
1257  Float_t tempDEta, tempDPhi;
1258  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(it->second));
1259  if(!tempTrack) continue;
1260  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->first,tempDEta,tempDPhi)){
1261  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) tempMatchedTracks.push_back(it->second);
1262  }
1263  }
1264  }
1265 
1266  return tempMatchedTracks;
1267 }
1268 
1269 //________________________________________________________________________
1270 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin){
1271  Int_t TrackPos = -1;
1272  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1273  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1274  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1275  if(currTrack->GetID() == trackID){
1276  TrackPos = iTrack;
1277  break;
1278  }
1279  }
1280  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1281  }else TrackPos = trackID; // for ESD just take trackID
1282 
1283  vector<Int_t> tempMatchedClusters;
1284  multimap<Int_t,Int_t>::iterator it;
1285  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1286  if(!tempTrack) return tempMatchedClusters;
1287  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1288  if(it->first == TrackPos){
1289  Float_t tempDEta, tempDPhi;
1290  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1291  if(tempTrack->Charge()>0){
1292  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin < tempDPhi) && (tempDPhi < dPhiMax) ) tempMatchedClusters.push_back(it->second);
1293  }else if(tempTrack->Charge()<0){
1294  dPhiMin*=-1;
1295  dPhiMax*=-1;
1296  if( (dEtaMin < tempDEta) && (tempDEta < dEtaMax) && (dPhiMin > tempDPhi) && (tempDPhi > dPhiMax) ) tempMatchedClusters.push_back(it->second);
1297  }
1298  }
1299  }
1300  }
1301 
1302  return tempMatchedClusters;
1303 }
1304 
1305 //________________________________________________________________________
1306 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, TF1* fFuncPtDepEta, TF1* fFuncPtDepPhi){
1307  Int_t TrackPos = -1;
1308  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1309  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1310  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1311  if(currTrack->GetID() == trackID){
1312  TrackPos = iTrack;
1313  break;
1314  }
1315  }
1316  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1317  }else TrackPos = trackID; // for ESD just take trackID
1318 
1319  vector<Int_t> tempMatchedClusters;
1320  multimap<Int_t,Int_t>::iterator it;
1321  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1322  if(!tempTrack) return tempMatchedClusters;
1323  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1324  if(it->first == TrackPos){
1325  Float_t tempDEta, tempDPhi;
1326  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1327  Bool_t match_dEta = kFALSE;
1328  Bool_t match_dPhi = kFALSE;
1329  if( TMath::Abs(tempDEta) < fFuncPtDepEta->Eval(tempTrack->Pt())) match_dEta = kTRUE;
1330  else match_dEta = kFALSE;
1331 
1332  if( TMath::Abs(tempDPhi) < fFuncPtDepPhi->Eval(tempTrack->Pt())) match_dPhi = kTRUE;
1333  else match_dPhi = kFALSE;
1334 
1335  if (match_dPhi && match_dEta )tempMatchedClusters.push_back(it->second);
1336  }
1337  }
1338  }
1339 
1340  return tempMatchedClusters;
1341 }
1342 
1343 //________________________________________________________________________
1344 vector<Int_t> AliCaloTrackMatcher::GetMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, Float_t dR){
1345  Int_t TrackPos = -1;
1346  if(event->IsA()==AliAODEvent::Class()){ // for AOD, we have to look for position of track in the event
1347  for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++){
1348  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(iTrack));
1349  if(currTrack->GetID() == trackID){
1350  TrackPos = iTrack;
1351  break;
1352  }
1353  }
1354  if(TrackPos == -1) AliFatal(Form("AliCaloTrackMatcher: GetNMatchedClusterIDsForTrack - track (ID: '%i') cannot be retrieved from event, should be impossible as it has been used in maim task before!",trackID));
1355  }else TrackPos = trackID; // for ESD just take trackID
1356 
1357  vector<Int_t> tempMatchedClusters;
1358  multimap<Int_t,Int_t>::iterator it;
1359  AliVTrack* tempTrack = dynamic_cast<AliVTrack*>(event->GetTrack(TrackPos));
1360  if(!tempTrack) return tempMatchedClusters;
1361  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it){
1362  if(it->first == TrackPos){
1363  Float_t tempDEta, tempDPhi;
1364  if(GetTrackClusterMatchingResidual(tempTrack->GetID(),it->second,tempDEta,tempDPhi)){
1365  if (TMath::Sqrt(tempDEta*tempDEta + tempDPhi*tempDPhi) < dR ) tempMatchedClusters.push_back(it->second);
1366  }
1367  }
1368  }
1369 
1370  return tempMatchedClusters;
1371 }
1372 
1373 //________________________________________________________________________
1375  Float_t sumTrackEt = 0.;
1376  vector<Int_t> labelsMatched = GetMatchedTrackIDsForCluster(event, clusterID, dR);
1377  if((Int_t) labelsMatched.size()<1) return sumTrackEt;
1378 
1379  TLorentzVector vecTrack;
1380  for (UInt_t i = 0; i < labelsMatched.size(); i++){
1381  AliVTrack* currTrack = dynamic_cast<AliVTrack*>(event->GetTrack(labelsMatched.at(i)));
1382  if(!currTrack) continue;
1383  vecTrack.SetPxPyPzE(currTrack->Px(),currTrack->Py(),currTrack->Pz(),currTrack->E());
1384  sumTrackEt += vecTrack.Et();
1385  }
1386 
1387  return sumTrackEt;
1388 }
1389 
1390 //________________________________________________________________________
1392  TAxis *axisafter = histoRebin->GetYaxis();
1393  Int_t bins = axisafter->GetNbins();
1394  Double_t from = axisafter->GetXmin();
1395  Double_t to = axisafter->GetXmax();
1396  Double_t *newbins = new Double_t[bins+1];
1397  newbins[0] = from;
1398  Double_t factor = TMath::Power(to/from, 1./bins);
1399  for(Int_t i=1; i<=bins; ++i) newbins[i] = factor * newbins[i-1];
1400  axisafter->Set(bins, newbins);
1401  delete [] newbins;
1402  return;
1403 }
1404 
1405 //________________________________________________________________________
1407  if(fSecVectorDeltaEtaDeltaPhi.size()>0){
1408  cout << "******************************" << endl;
1409  cout << "******************************" << endl;
1410  cout << "NEW EVENT !" << endl;
1411  cout << "vector etaphi:" << endl;
1412  cout << fSecVectorDeltaEtaDeltaPhi.size() << endl;
1413  cout << "multimap" << endl;
1414  mapT::iterator iter = fSecMap_TrID_ClID_ToIndex.begin();
1415  for (iter = fSecMap_TrID_ClID_ToIndex.begin(); iter != fSecMap_TrID_ClID_ToIndex.end(); ++iter){
1416  Float_t dEta, dPhi = 0;
1417  if(!GetSecTrackClusterMatchingResidual(iter->first.first,iter->first.second,dEta,dPhi)) continue;
1418  cout << " [" << iter->first.first << "/" << iter->first.second << ", " << iter->second << "] - (" << dEta << "/" << dPhi << ")" << endl;
1419  }
1420  cout << "mapTrackToCluster" << endl;
1421  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(fInputEvent);
1422  for (Int_t itr=0;itr<esdev->GetNumberOfTracks();itr++){
1423  AliVTrack *inTrack = esdev->GetTrack(itr);
1424  if(!inTrack) continue;
1425  TString tCharge;
1426  if(inTrack->Charge()>0) tCharge = "+";
1427  else if(inTrack->Charge()<0) tCharge = "-";
1428  cout << itr << " (" << tCharge << ") - " << GetNMatchedClusterIDsForSecTrack(fInputEvent,inTrack->GetID(),5,-5,0.2,-0.4) << "\t\t";
1429  }
1430  cout << endl;
1431  multimap<Int_t,Int_t>::iterator it;
1432  for (it=fSecMapTrackToCluster.begin(); it!=fSecMapTrackToCluster.end(); ++it) cout << it->first << " => " << it->second << '\n';
1433  cout << "mapClusterToTrack" << endl;
1434  Int_t tempClus = it->second;
1435  for (it=fSecMapClusterToTrack.begin(); it!=fSecMapClusterToTrack.end(); ++it) cout << it->first << " => " << it->second << '\n';
1436  vector<Int_t> tempTracks = GetMatchedSecTrackIDsForCluster(fInputEvent,tempClus, 5, -5, 0.2, -0.4);
1437  for(UInt_t iJ=0; iJ<tempTracks.size();iJ++){
1438  cout << tempClus << " - " << tempTracks.at(iJ) << endl;
1439  }
1440  }
1441  return;
1442 }
1443 
1444 //________________________________________________________________________
1446  if(fVectorDeltaEtaDeltaPhi.size()>0){
1447  cout << "******************************" << endl;
1448  cout << "******************************" << endl;
1449  cout << "NEW EVENT !" << endl;
1450  cout << "vector etaphi:" << endl;
1451  cout << fVectorDeltaEtaDeltaPhi.size() << endl;
1452  cout << "multimap" << endl;
1453  mapT::iterator iter = fMap_TrID_ClID_ToIndex.begin();
1454  for (iter = fMap_TrID_ClID_ToIndex.begin(); iter != fMap_TrID_ClID_ToIndex.end(); ++iter){
1455  Float_t dEta, dPhi = 0;
1456  if(!GetTrackClusterMatchingResidual(iter->first.first,iter->first.second,dEta,dPhi)) continue;
1457  cout << " [" << iter->first.first << "/" << iter->first.second << ", " << iter->second << "] - (" << dEta << "/" << dPhi << ")" << endl;
1458  }
1459  cout << "mapTrackToCluster" << endl;
1460  AliESDEvent *esdev = dynamic_cast<AliESDEvent*>(fInputEvent);
1461  for (Int_t itr=0;itr<esdev->GetNumberOfTracks();itr++){
1462  AliVTrack *inTrack = esdev->GetTrack(itr);
1463  if(!inTrack) continue;
1464  TString tCharge;
1465  if(inTrack->Charge()>0) tCharge = "+";
1466  else if(inTrack->Charge()<0) tCharge = "-";
1467  cout << itr << " (" << tCharge << ") - " << GetNMatchedClusterIDsForTrack(fInputEvent,inTrack->GetID(),5,-5,0.2,-0.4) << "\t\t";
1468  }
1469  cout << endl;
1470  multimap<Int_t,Int_t>::iterator it;
1471  for (it=fMapTrackToCluster.begin(); it!=fMapTrackToCluster.end(); ++it) cout << it->first << " => " << it->second << '\n';
1472  cout << "mapClusterToTrack" << endl;
1473  Int_t tempClus = it->second;
1474  for (it=fMapClusterToTrack.begin(); it!=fMapClusterToTrack.end(); ++it) cout << it->first << " => " << it->second << '\n';
1475  vector<Int_t> tempTracks = GetMatchedTrackIDsForCluster(fInputEvent,tempClus, 5, -5, 0.2, -0.4);
1476  for(UInt_t iJ=0; iJ<tempTracks.size();iJ++){
1477  cout << tempClus << " - " << tempTracks.at(iJ) << endl;
1478  }
1479  }
1480  return;
1481 }
double Double_t
Definition: External.C:58
Definition: External.C:236
Int_t GetNMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
void ProcessEvent(AliVEvent *event)
Bool_t GetTrackClusterMatchingResidual(Int_t trackID, Int_t clusterID, Float_t &dEta, Float_t &dPhi)
void Initialize(Int_t runNumber)
vector< Int_t > GetMatchedClusterIDsForSecTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
Int_t GetNMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
float Float_t
Definition: External.C:68
virtual void Terminate(Option_t *)
Bool_t GetSecTrackClusterMatchingResidual(Int_t trackID, Int_t clusterID, Float_t &dEta, Float_t &dPhi)
Float_t SumTrackEtAroundCluster(AliVEvent *event, Int_t clusterID, Float_t dR)
Int_t GetNMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
vector< Int_t > GetMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
virtual void UserExec(Option_t *option)
Definition: External.C:220
Int_t GetNMatchedSecTrackIDsForCluster(AliVEvent *event, Int_t trackID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
Bool_t PropagateV0TrackToClusterAndGetMatchingResidual(AliVTrack *inSecTrack, AliVCluster *cluster, AliVEvent *event, Float_t &dEta, Float_t &dPhi)
pair< Float_t, Float_t > pairFloat
Bool_t IsSecTrackClusterAlreadyTried(Int_t trackID, Int_t clusterID)
const char Option_t
Definition: External.C:48
bool Bool_t
Definition: External.C:53
vector< Int_t > GetMatchedClusterIDsForTrack(AliVEvent *event, Int_t trackID, TF1 *fFuncPtDepEta, TF1 *fFuncPtDepPhi)
vector< Int_t > GetMatchedTrackIDsForCluster(AliVEvent *event, Int_t clusterID, Float_t dEtaMax, Float_t dEtaMin, Float_t dPhiMax, Float_t dPhiMin)
void SetLogBinningYTH2(TH2 *histoRebin)