AliPhysics  vAN-20150427 (e6e7aad)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliAnalysisTaskMuonPerformance.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2007, 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 //-----------------------------------------------------------------------------
22 //-----------------------------------------------------------------------------
23 
24 // ROOT includes
25 #include "TROOT.h"
26 #include "TH1.h"
27 #include "TH2.h"
28 #include "TH3.h"
29 #include "TF1.h"
30 #include "TLegend.h"
31 #include "TGraphErrors.h"
32 #include "TGraphAsymmErrors.h"
33 #include "TAxis.h"
34 #include "TObjArray.h"
35 #include "TCanvas.h"
36 #include "TMath.h"
37 #include "TMCProcess.h"
38 #include "TGeoGlobalMagField.h"
39 #include "TGeoManager.h"
40 
41 // STEER includes
42 #include "AliMCEvent.h"
43 #include "AliMCParticle.h"
44 #include "AliMCEventHandler.h"
45 #include "AliESDEvent.h"
46 #include "AliESDMuonTrack.h"
47 #include "AliCDBManager.h"
48 #include "AliGeomManager.h"
49 
50 // ANALYSIS includes
51 #include "AliAnalysisManager.h"
52 #include "AliAnalysisDataSlot.h"
53 #include "AliAnalysisDataContainer.h"
54 #include "AliCentrality.h"
55 
56 // CORRFW includes
57 #include "AliCFContainer.h"
58 #include "AliCFGridSparse.h"
59 #include "AliCFEffGrid.h"
60 
61 // MUON includes
62 #include "AliMUONConstants.h"
63 #include "AliMUONTrack.h"
64 #include "AliMUONTriggerTrack.h"
65 #include "AliMUONLocalTrigger.h"
66 #include "AliMUONVTrackStore.h"
67 #include "AliMUONVTriggerTrackStore.h"
68 #include "AliMUONESDInterface.h"
69 #include "AliMUONRecoParam.h"
70 #include "AliMUONCDB.h"
71 #include "AliMUONTrackExtrap.h"
72 #include "AliMUONTrackParam.h"
73 #include "AliMUONRecoCheck.h"
74 #include "AliMUONVCluster.h"
75 #include "AliMUONVTrackReconstructor.h"
76 
77 // MUON mapping includes
78 #include "AliMpSegmentation.h"
79 #include "AliMpDEIterator.h"
80 
82 
83 
84 using std::cout;
85 using std::endl;
86 using std::flush;
87 
89 ClassImp(AliAnalysisTaskMuonPerformance) // Class implementation in ROOT context
91 
92 
93 //________________________________________________________________________
95 AliAnalysisTaskSE(),
96 fDefaultStorage(""),
97 fAlignOCDBpath(""),
98 fRecoParamOCDBpath(""),
99 fNPBins(30),
100 fCorrectForSystematics(kFALSE),
101 fFitResiduals(kFALSE),
102 fEnforceTrkCriteria(kFALSE),
103 fUseMCKinematics(kFALSE),
104 fMCTrigLevelFromMatchTrk(kFALSE),
105 fRequestedStationMask(0),
106 fRequest2ChInSameSt45(0),
107 fSigmaCutTrig(-1.),
108 fNDE(0),
109 fCFContainer(0x0),
110 fEfficiencyList(0x0),
111 fTriggerList(0x0),
112 fTrackerList(0x0),
113 fPAtVtxList(0x0),
114 fSlopeAtVtxList(0x0),
115 fEtaAtVtxList(0x0),
116 fPhiAtVtxList(0x0),
117 fPAt1stClList(0x0),
118 fSlopeAt1stClList(0x0),
119 fDCAList(0x0),
120 fClusterList(0x0)
121 {
122  //
124  //
125  fPRange[0] = 0.;
126  fPRange[1] = 0.;
127  fClusterMaxRes[0] = 0.;
128  fClusterMaxRes[1] = 0.;
129  for (Int_t i = 0; i < 1100; i++) fDEIndices[i] = 0;
130  for (Int_t i = 0; i < 200; i++) fDEIds[i] = 0;
131 }
132 
133 
134 //________________________________________________________________________
136 AliAnalysisTaskSE(name),
137 fDefaultStorage("raw://"),
138 fAlignOCDBpath(""),
139 fRecoParamOCDBpath(""),
140 fNPBins(30),
141 fCorrectForSystematics(kFALSE),
142 fFitResiduals(kFALSE),
143 fEnforceTrkCriteria(kFALSE),
144 fUseMCKinematics(kFALSE),
145 fMCTrigLevelFromMatchTrk(kFALSE),
146 fRequestedStationMask(0),
147 fRequest2ChInSameSt45(0),
148 fSigmaCutTrig(-1.),
149 fNDE(0),
150 fCFContainer(0x0),
151 fEfficiencyList(0x0),
152 fTriggerList(0x0),
153 fTrackerList(0x0),
154 fPAtVtxList(0x0),
155 fSlopeAtVtxList(0x0),
156 fEtaAtVtxList(0x0),
157 fPhiAtVtxList(0x0),
158 fPAt1stClList(0x0),
159 fSlopeAt1stClList(0x0),
160 fDCAList(0x0),
161 fClusterList(0x0)
162 {
163  //
165  //
166  fPRange[0] = 0.;
167  fPRange[1] = 300.;
168  fClusterMaxRes[0] = 0.;
169  fClusterMaxRes[1] = 0.;
170  for (Int_t i = 0; i < 1100; i++) fDEIndices[i] = 0;
171  for (Int_t i = 0; i < 200; i++) fDEIds[i] = 0;
172 
173  DefineOutput(1, AliCFContainer::Class());
174  DefineOutput(2, TObjArray::Class());
175  DefineOutput(3, TObjArray::Class());
176  DefineOutput(4, TObjArray::Class());
177  DefineOutput(5, TObjArray::Class());
178  DefineOutput(6, TObjArray::Class());
179  DefineOutput(7, TObjArray::Class());
180  DefineOutput(8, TObjArray::Class());
181  DefineOutput(9, TObjArray::Class());
182  DefineOutput(10, TObjArray::Class());
183  DefineOutput(11, TObjArray::Class());
184  DefineOutput(12, TObjArray::Class());
185 }
186 
187 
188 //________________________________________________________________________
190 {
191  //
193  //
194  if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
195  delete fCFContainer;
196  delete fEfficiencyList;
197  delete fTriggerList;
198  delete fTrackerList;
199  delete fPAtVtxList;
200  delete fSlopeAtVtxList;
201  delete fEtaAtVtxList;
202  delete fPhiAtVtxList;
203  delete fPAt1stClList;
204  delete fSlopeAt1stClList;
205  delete fDCAList;
206  delete fClusterList;
207  }
208 }
209 
210 
211 //___________________________________________________________________________
213 {
214  //
217  //
218 
219  // load OCDB objects only once
220  if (fSigmaCutTrig > 0) return;
221 
222  // set OCDB location
223  AliCDBManager* cdbm = AliCDBManager::Instance();
224  if (cdbm->IsDefaultStorageSet()) printf("PerformanceTask: CDB default storage already set!\n");
225  else {
226  cdbm->SetDefaultStorage(fDefaultStorage.Data());
227  if (!fAlignOCDBpath.IsNull()) cdbm->SetSpecificStorage("MUON/Align/Data",fAlignOCDBpath.Data());
228  if (!fRecoParamOCDBpath.IsNull()) cdbm->SetSpecificStorage("MUON/Calib/RecoParam",fRecoParamOCDBpath.Data());
229  }
230  if (cdbm->GetRun() > -1) printf("PerformanceTask: run number already set!\n");
231  else cdbm->SetRun(fCurrentRunNumber);
232 
233  // load magnetic field for track extrapolation
234  if (!TGeoGlobalMagField::Instance()->GetField()) {
235  if (!AliMUONCDB::LoadField()) return;
236  }
237 
238  // load mapping
239  if (!AliMpSegmentation::Instance(kFALSE)) {
240  if (!AliMUONCDB::LoadMapping(kTRUE)) return;
241  }
242 
243  // load geometry for track extrapolation to vertex and for checking hits are under pads in reconstructible tracks
244  if (!AliGeomManager::GetGeometry()) {
245  AliGeomManager::LoadGeometry();
246  if (!AliGeomManager::GetGeometry()) return;
247  if (!AliGeomManager::ApplyAlignObjsFromCDB("MUON")) return;
248  }
249 
250  // load recoParam
251  const AliMUONRecoParam* recoParam = (AliMUONESDInterface::GetTracker())
252  ? AliMUONESDInterface::GetTracker()->GetRecoParam()
253  : AliMUONCDB::LoadRecoParam();
254 
255  if (!recoParam) {
257  fRequest2ChInSameSt45 = kFALSE;
258  fSigmaCutTrig = -1.;
259  AliError("--> skip this run");
260  return;
261  }
262 
263  // compute the mask of requested stations from recoParam
265  for (Int_t i = 0; i < 5; i++) if (recoParam->RequestStation(i)) fRequestedStationMask |= ( 1 << i );
266 
267  // get from recoParam whether a track need 2 chambers hit in the same station (4 or 5) or not to be reconstructible
268  fRequest2ChInSameSt45 = !recoParam->MakeMoreTrackCandidates();
269 
270  // get sigma cut from recoParam to associate trigger track to triggerable track
271  fSigmaCutTrig = recoParam->GetSigmaCutForTrigger();
272 
273  if (!AliMUONESDInterface::GetTracker()) AliMUONESDInterface::ResetTracker(recoParam);
274 
275  for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
276 
277  // find the highest chamber resolution and set histogram bins
278  if (recoParam->GetDefaultNonBendingReso(i) > fClusterMaxRes[0]) fClusterMaxRes[0] = recoParam->GetDefaultNonBendingReso(i);
279  if (recoParam->GetDefaultBendingReso(i) > fClusterMaxRes[1]) fClusterMaxRes[1] = recoParam->GetDefaultBendingReso(i);
280 
281  // fill correspondence between DEId and indices for histo (starting from 1)
282  AliMpDEIterator it;
283  it.First(i);
284  while (!it.IsDone()) {
285  fNDE++;
286  fDEIndices[it.CurrentDEId()] = fNDE;
287  fDEIds[fNDE] = it.CurrentDEId();
288  it.Next();
289  }
290 
291  }
292 
294 }
295 
296 
297 //___________________________________________________________________________
299 {
300  //
302  //
303 
304  // do it once the OCDB has been loaded (i.e. from NotifyRun())
305  if (fSigmaCutTrig < 0) return;
306 
307  // ------ CFContainer ------
308 
309  // define axes of particle container
310  Int_t nPtBins = 30;
311  Double_t ptMin = 0., ptMax = 15.;
312  TString ptTitle("p_{t}"), ptUnits("GeV/c");
313 
314  Int_t nEtaBins = 15;
315  Double_t etaMin = -4., etaMax = -2.5;
316  TString etaTitle("#eta"), etaUnits("a.u.");
317 
318  Int_t nPhiBins = 15;
319  Double_t phiMin = 0.; Double_t phiMax = 2.*TMath::Pi();
320  TString phiTitle("#phi"), phiUnits("rad");
321 
322  Int_t nThetaAbsEndBins = 4;
323  Double_t thetaAbsEndMin = -0.5, thetaAbsEndMax = 3.5;
324  TString thetaAbsEndTitle("#theta_{abs}"), thetaAbsEndUnits("a.u.");
325 
326  Int_t nChargeBins = 2;
327  Double_t chargeMin = -2, chargeMax = 2.;
328  TString chargeTitle("charge"), chargeUnits("e");
329 
330  Int_t nHasTrackerBins = 2;
331  Double_t hasTrackerMin = -0.5, hasTrackerMax = (Double_t)nHasTrackerBins - 0.5;
332  TString hasTrackerTitle("Has tracker"), hasTrackerUnits("");
333 
334  Int_t nTriggerBins = kNtrigCuts;
335  Double_t triggerMin = -0.5, triggerMax = (Double_t)nTriggerBins - 0.5;
336  TString triggerTitle("Trigger match"), triggerUnits("");
337 
338  Int_t nMotherTypeBins = kNtrackSources;
339  Double_t motherTypeMin = -0.5, motherTypeMax = (Double_t)kNtrackSources - 0.5;
340  TString motherTypeTitle("motherType"), motherTypeUnits("");
341 
342  Int_t nMatchMCBins = kNMatchMC;
343  Double_t matchMCMin = -0.5, matchMCMax = (Double_t)kNMatchMC - 0.5;
344  TString matchMCTitle("MatchMC"), matchMCUnits("");
345 
346  Int_t nMCTriggerBins = kNtrigCuts;
347  Double_t mcTriggerMin = -0.5, mcTriggerMax = (Double_t)nMCTriggerBins - 0.5;
348  TString mcTriggerTitle("MC Trigger match"), mcTriggerUnits("");
349 
350  Int_t nCentBins = 22;
351  Double_t centMin = -5., centMax = 105.;
352  TString centTitle("centrality"), centUnits("%");
353 
354  Int_t nDupliTrgBins = 2;
355  Double_t dupliTrgMin = -0.5, dupliTrgMax = 1.5;
356  TString dupliTrgTitle("duplicate trigger"), dupliTrgUnits("");
357 
358  Int_t nbins[kNvars] = {nPtBins, nEtaBins, nPhiBins, nThetaAbsEndBins, nChargeBins, nHasTrackerBins, nTriggerBins, nMotherTypeBins, nMatchMCBins, nMCTriggerBins, nCentBins, nDupliTrgBins};
359  Double_t xmin[kNvars] = {ptMin, etaMin, phiMin, thetaAbsEndMin, chargeMin, hasTrackerMin, triggerMin, motherTypeMin, matchMCMin, mcTriggerMin, centMin, dupliTrgMin};
360  Double_t xmax[kNvars] = {ptMax, etaMax, phiMax, thetaAbsEndMax, chargeMax, hasTrackerMax, triggerMax, motherTypeMax, matchMCMax, mcTriggerMax, centMax, dupliTrgMax};
361  TString axisTitle[kNvars] = {ptTitle, etaTitle, phiTitle, thetaAbsEndTitle, chargeTitle, hasTrackerTitle, triggerTitle, motherTypeTitle, matchMCTitle, mcTriggerTitle, centTitle, dupliTrgTitle};
362  TString axisUnits[kNvars] = {ptUnits, etaUnits, phiUnits, thetaAbsEndUnits, chargeUnits, hasTrackerUnits, triggerUnits, motherTypeUnits, matchMCUnits, mcTriggerUnits, centUnits, dupliTrgUnits};
363 
364  // create particle container
365  fCFContainer = new AliCFContainer(GetOutputSlot(1)->GetContainer()->GetName(),"container for tracks",kNsteps,kNvars,nbins);
366 
367  // set axes title and limits
368  for ( Int_t idim = 0; idim<kNvars; idim++) {
369  TString histoTitle = Form("%s (%s)", axisTitle[idim].Data(), axisUnits[idim].Data());
370  histoTitle.ReplaceAll("()","");
371  fCFContainer->SetVarTitle(idim, histoTitle.Data());
372  fCFContainer->SetBinLimits(idim, xmin[idim], xmax[idim]);
373  }
374 
375  // define histo names
376  TString stepTitle[kNsteps] = {"reconstructed", "generated"};
377 
378  // define axes labels if any
379  TString trigName[kNtrigCuts];
380  trigName[kNoMatchTrig] = "NoMatch";
381  trigName[kOtherTrig] = "Other";
382  trigName[kAllPtTrig] = "AllPt";
383  trigName[kLowPtTrig] = "LowPt";
384  trigName[kHighPtTrig] = "HighPt";
385 
386  TString srcName[kNtrackSources];
387  srcName[kCharmMu] = "Charm";
388  srcName[kBeautyMu] = "Beauty";
389  srcName[kPrimaryMu] = "Decay";
390  srcName[kSecondaryMu] = "Secondary";
391  srcName[kRecoHadron] = "Hadrons";
392  srcName[kUnknownPart] = "Fakes";
393 
394  TString mMCName[kNMatchMC];
395  mMCName[kNoMatch] = "NoMatch";
396  mMCName[kTrackerOnly] = "TrackerOnly";
397  mMCName[kMatchedSame] = "MatchedSame";
398  mMCName[kMatchedDiff] = "MatchedDiff";
399  mMCName[kTriggerOnly] = "TriggerOnly";
400 
401  // set histo name and axis labels if any
402  for (Int_t istep=0; istep<kNsteps; istep++) {
403 
404  // histo name
405  fCFContainer->SetStepTitle(istep, stepTitle[istep].Data());
406  AliCFGridSparse* gridSparse = fCFContainer->GetGrid(istep);
407 
408  // trigger info
409  TAxis* triggerAxis = gridSparse->GetAxis(kVarTrigger);
410  for ( Int_t ibin=0; ibin<kNtrigCuts; ibin++ ) {
411  triggerAxis->SetBinLabel(ibin+1,trigName[ibin]);
412  }
413 
414  // mother type
415  TAxis* motherTypeAxis = gridSparse->GetAxis(kVarMotherType);
416  for ( Int_t ibin=0; ibin<kNtrackSources; ibin++ ) {
417  motherTypeAxis->SetBinLabel(ibin+1,srcName[ibin]);
418  }
419 
420  // MC matching flag
421  TAxis* matchMCAxis = gridSparse->GetAxis(kVarMatchMC);
422  for ( Int_t ibin=0; ibin<kNMatchMC; ibin++ ) {
423  matchMCAxis->SetBinLabel(ibin+1,mMCName[ibin]);
424  }
425 
426  // MC trigger info
427  TAxis* mcTriggerAxis = gridSparse->GetAxis(kVarMCTrigger);
428  for ( Int_t ibin=0; ibin<kNtrigCuts; ibin++ ) {
429  mcTriggerAxis->SetBinLabel(ibin+1,trigName[ibin]);
430  }
431 
432  }
433 
434  // ------ trigger histograms ------
435 
436  fTriggerList = new TObjArray(100);
437  fTriggerList->SetOwner();
438 
439  TH1F* h1 = new TH1F("hResTrigX11", "Residual X11;X11_{reco} - X11_{MC} (cm)", 100, -10., 10.);
440  fTriggerList->AddAt(h1, kResTrigX11);
441  h1 = new TH1F("hResTrigY11", "Residual Y11;Y11_{reco} - Y11_{MC} (cm)", 100, -10., 10.);
442  fTriggerList->AddAt(h1, kResTrigY11);
443  h1 = new TH1F("hResTrigSlopeY", "Residual slope y;ySlope_{reco} - ySlope_{MC} (rad)", 100, -0.1, 0.1);
444  fTriggerList->AddAt(h1, kResTrigSlopeY);
445 
446  // ------ tracker histograms ------
447 
448  fTrackerList = new TObjArray(100);
449  fTrackerList->SetOwner();
450 
451  // momentum resolution at vertex
452  const Int_t deltaPAtVtxNBins = 250;
453  Double_t deltaPAtVtxEdges[2];
454  deltaPAtVtxEdges[0] = -20. - 0.05 * fPRange[1];
455  deltaPAtVtxEdges[1] = 5. + 0.05 * fPRange[1];
456 
457  h1 = new TH1F("hResPAtVtx"," delta P at vertex;#Delta_{p} (GeV/c)",deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
458  fTrackerList->AddAt(h1, kResPAtVtx);
459  TH2F *h2 = new TH2F("hResPAtVtxVsP","#Delta_{p} at vertex versus p;p (GeV/c);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
460  fTrackerList->AddAt(h2, kResPAtVtxVsP);
461  h2 = new TH2F("hResPAtVtxVsPIn23deg","#Delta_{p} at vertex versus p for tracks between 2 and 3 degrees at absorber end;p (GeV/c);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
463  h2 = new TH2F("hResPAtVtxVsPIn310deg","#Delta_{p} at vertex versus p for tracks between 3 and 10 degrees at absorber end;p (GeV/c);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
465  h2 = new TH2F("hResPAtVtxVsPIn02degMC","#Delta_{p} at vertex versus p for tracks with MC angle below 2 degrees;p (GeV/c);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],deltaPAtVtxNBins/10,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
467  h2 = new TH2F("hResPAtVtxVsPosAbsEndIn02degMC","#Delta_{p} at vertex versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);#Delta_{p} (GeV/c)",1000,0.,100.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
469  h2 = new TH2F("hResPAtVtxVsPosAbsEndIn23degMC","#Delta_{p} at vertex versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);#Delta_{p} (GeV/c)",1000,0.,100.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
471  h2 = new TH2F("hResPAtVtxVsPosAbsEndIn310degMC","#Delta_{p} at vertex versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);#Delta_{p} (GeV/c)",1000,0.,100.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
473  h2 = new TH2F("hResPAtVtxVsAngleAtAbsEnd","#Delta_{p} at vertex versus track position at absorber end converted to degrees;angle (Deg);#Delta_{p} (GeV/c)",10,0.,10.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
475  h2 = new TH2F("hResPAtVtxVsMCAngle","#Delta_{p} at vertex versus MC angle;MC angle (Deg);#Delta_{p} (GeV/c)",10,0.,10.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
476  fTrackerList->AddAt(h2, kResPAtVtxVsMCAngle);
477  TH3F *h3 = new TH3F("hResPAtVtxVsAngleAtAbsEndVsP","#Delta_{p} at vertex versus track position at absorber end converted to degrees versus momentum;p (GeV/c);angle (Deg);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],100,0.,10.,deltaPAtVtxNBins,deltaPAtVtxEdges[0],deltaPAtVtxEdges[1]);
479 
480  // transverse momentum resolution at vertex
481  h2 = new TH2F("hResPtAtVtxVsPt","#Delta_{p_{t}} at vertex versus p_{t};p_{t} (GeV/c);#Delta_{p_{t}} (GeV/c)",2*fNPBins,fPRange[0]/10.,fPRange[1]/10.,deltaPAtVtxNBins,deltaPAtVtxEdges[0]/10.,deltaPAtVtxEdges[1]/10.);
482  fTrackerList->AddAt(h2, kResPtAtVtxVsPt);
483 
484  // momentum resolution at first cluster
485  const Int_t deltaPAtFirstClNBins = 500;
486  Double_t deltaPAtFirstClEdges[2];
487  deltaPAtFirstClEdges[0] = -5. - 0.05 * fPRange[1];
488  deltaPAtFirstClEdges[1] = 5. + 0.05 * fPRange[1];
489 
490  h1 = new TH1F("hResPAt1stCl"," delta P at first cluster;#Delta_{p} (GeV/c)",deltaPAtFirstClNBins,deltaPAtFirstClEdges[0],deltaPAtFirstClEdges[1]);
491  fTrackerList->AddAt(h1, kResPAt1stCl);
492  h2 = new TH2F("hResPAt1stClVsP","#Delta_{p} at first cluster versus p;p (GeV/c);#Delta_{p} (GeV/c)",2*fNPBins,fPRange[0],fPRange[1],deltaPAtFirstClNBins,deltaPAtFirstClEdges[0],deltaPAtFirstClEdges[1]);
493  fTrackerList->AddAt(h2, kResPAt1stClVsP);
494 
495  // transverse momentum resolution at first cluster
496  h2 = new TH2F("hResPtAt1stClVsPt","#Delta_{p_{t}} at first cluster versus p_{t};p_{t} (GeV/c);#Delta_{p_{t}} (GeV/c)",2*fNPBins,fPRange[0]/10.,fPRange[1]/10.,deltaPAtFirstClNBins,deltaPAtFirstClEdges[0]/10.,deltaPAtFirstClEdges[1]/10.);
497  fTrackerList->AddAt(h2, kResPtAt1stClVsPt);
498 
499  // angular resolution at vertex
500  const Int_t deltaSlopeAtVtxNBins = 500;
501  const Double_t deltaSlopeAtVtxEdges[2] = {-0.05, 0.05};
502 
503  h1 = new TH1F("hResSlopeXAtVtx","#Delta_{slope_{X}} at vertex;#Delta_{slope_{X}}", deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
504  fTrackerList->AddAt(h1, kResSlopeXAtVtx);
505  h1 = new TH1F("hResSlopeYAtVtx","#Delta_{slope_{Y}} at vertex;#Delta_{slope_{Y}}", deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
506  fTrackerList->AddAt(h1, kResSlopeYAtVtx);
507  h2 = new TH2F("hResSlopeXAtVtxVsP","#Delta_{slope_{X}} at vertex versus p;p (GeV/c);#Delta_{slope_{X}}",2*fNPBins,fPRange[0],fPRange[1], deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
508  fTrackerList->AddAt(h2, kResSlopeXAtVtxVsP);
509  h2 = new TH2F("hResSlopeYAtVtxVsP","#Delta_{slope_{Y}} at vertex versus p;p (GeV/c);#Delta_{slope_{Y}}",2*fNPBins,fPRange[0],fPRange[1], deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
510  fTrackerList->AddAt(h2, kResSlopeYAtVtxVsP);
511  h2 = new TH2F("hResSlopeXAtVtxVsPosAbsEndIn02degMC","#Delta_{slope_{X}} at vertex versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);#Delta_{slope_{X}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
513  h2 = new TH2F("hResSlopeYAtVtxVsPosAbsEndIn02degMC","#Delta_{slope_{Y}} at vertex versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);#Delta_{slope_{Y}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
515  h2 = new TH2F("hResSlopeXAtVtxVsPosAbsEndIn23degMC","#Delta_{slope_{X}} at vertex versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);#Delta_{slope_{X}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
517  h2 = new TH2F("hResSlopeYAtVtxVsPosAbsEndIn23degMC","#Delta_{slope_{Y}} at vertex versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);#Delta_{slope_{Y}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
519  h2 = new TH2F("hResSlopeXAtVtxVsPosAbsEndIn310degMC","#Delta_{slope_{X}} at vertex versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);#Delta_{slope_{X}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
521  h2 = new TH2F("hResSlopeYAtVtxVsPosAbsEndIn310degMC","#Delta_{slope_{Y}} at vertex versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);#Delta_{slope_{Y}}",1000,0.,100.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
523  h2 = new TH2F("hResSlopeXAtVtxVsAngleAtAbsEnd","#Delta_{slope_{X}} at vertex versus track position at absorber end converted to degrees;angle (Deg);#Delta_{slope_{X}}",10,0.,10.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
525  h2 = new TH2F("hResSlopeYAtVtxVsAngleAtAbsEnd","#Delta_{slope_{Y}} at vertex versus track position at absorber end converted to degrees;angle (Deg);#Delta_{slope_{Y}}",10,0.,10.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
527  h2 = new TH2F("hResSlopeXAtVtxVsMCAngle","#Delta_{slope_{X}} at vertex versus MC angle;MC angle (Deg);#Delta_{slope_{X}}",10,0.,10.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
529  h2 = new TH2F("hResSlopeYAtVtxVsMCAngle","#Delta_{slope_{Y}} at vertex versus MC angle;MC angle (Deg);#Delta_{slope_{Y}}",10,0.,10.,deltaSlopeAtVtxNBins, deltaSlopeAtVtxEdges[0], deltaSlopeAtVtxEdges[1]);
531 
532  // angular resolution at first cluster
533  const Int_t deltaSlopeAtFirstClNBins = 500;
534  const Double_t deltaSlopeAtFirstClEdges[2] = {-0.01, 0.01};
535 
536  h1 = new TH1F("hResSlopeXAt1stCl","#Delta_{slope_{X}} at first cluster;#Delta_{slope_{X}}", deltaSlopeAtFirstClNBins, deltaSlopeAtFirstClEdges[0], deltaSlopeAtFirstClEdges[1]);
537  fTrackerList->AddAt(h1, kResSlopeXAt1stCl);
538  h1 = new TH1F("hResSlopeYAt1stCl","#Delta_{slope_{Y}} at first cluster;#Delta_{slope_{Y}}", deltaSlopeAtFirstClNBins, deltaSlopeAtFirstClEdges[0], deltaSlopeAtFirstClEdges[1]);
539  fTrackerList->AddAt(h1, kResSlopeYAt1stCl);
540  h2 = new TH2F("hResSlopeXAt1stClVsP","#Delta_{slope_{X}} at first cluster versus p;p (GeV/c);#Delta_{slope_{X}}",2*fNPBins,fPRange[0],fPRange[1], deltaSlopeAtFirstClNBins, deltaSlopeAtFirstClEdges[0], deltaSlopeAtFirstClEdges[1]);
542  h2 = new TH2F("hResSlopeYAt1stClVsP","#Delta_{slope_{Y}} at first cluster versus p;p (GeV/c);#Delta_{slope_{Y}}",2*fNPBins,fPRange[0],fPRange[1], deltaSlopeAtFirstClNBins, deltaSlopeAtFirstClEdges[0], deltaSlopeAtFirstClEdges[1]);
544 
545  // eta resolution at vertex
546  const Int_t deltaEtaAtVtxNBins = 500;
547  const Double_t deltaEtaAtVtxEdges[2] = {-0.5, 0.5};
548 
549  h1 = new TH1F("hResEtaAtVtx","#Delta_{eta} at vertex;#Delta_{eta}", deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
550  fTrackerList->AddAt(h1, kResEtaAtVtx);
551  h2 = new TH2F("hResEtaAtVtxVsP","#Delta_{eta} at vertex versus p;p (GeV/c);#Delta_{eta}",2*fNPBins,fPRange[0],fPRange[1], deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
552  fTrackerList->AddAt(h2, kResEtaAtVtxVsP);
553  h2 = new TH2F("hResEtaAtVtxVsPosAbsEndIn02degMC","#Delta_{eta} at vertex versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);#Delta_{eta}",1000,0.,100.,deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
555  h2 = new TH2F("hResEtaAtVtxVsPosAbsEndIn23degMC","#Delta_{eta} at vertex versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);#Delta_{eta}",1000,0.,100.,deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
557  h2 = new TH2F("hResEtaAtVtxVsPosAbsEndIn310degMC","#Delta_{eta} at vertex versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);#Delta_{eta}",1000,0.,100.,deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
559  h2 = new TH2F("hResEtaAtVtxVsAngleAtAbsEnd","#Delta_{eta} at vertex versus track position at absorber end converted to degrees;angle (Deg);#Delta_{eta}",10,0.,10.,deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
561  h2 = new TH2F("hResEtaAtVtxVsMCAngle","#Delta_{eta} at vertex versus MC angle;MC angle (Deg);#Delta_{eta}",10,0.,10.,deltaEtaAtVtxNBins, deltaEtaAtVtxEdges[0], deltaEtaAtVtxEdges[1]);
563 
564  // phi resolution at vertex
565  const Int_t deltaPhiAtVtxNBins = 500;
566  const Double_t deltaPhiAtVtxEdges[2] = {-0.5, 0.5};
567 
568  h1 = new TH1F("hResPhiAtVtx","#Delta_{phi} at vertex;#Delta_{phi}", deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
569  fTrackerList->AddAt(h1, kResPhiAtVtx);
570  h2 = new TH2F("hResPhiAtVtxVsP","#Delta_{phi} at vertex versus p;p (GeV/c);#Delta_{phi}",2*fNPBins,fPRange[0],fPRange[1], deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
571  fTrackerList->AddAt(h2, kResPhiAtVtxVsP);
572  h2 = new TH2F("hResPhiAtVtxVsPosAbsEndIn02degMC","#Delta_{phi} at vertex versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);#Delta_{phi}",1000,0.,100.,deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
574  h2 = new TH2F("hResPhiAtVtxVsPosAbsEndIn23degMC","#Delta_{phi} at vertex versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);#Delta_{phi}",1000,0.,100.,deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
576  h2 = new TH2F("hResPhiAtVtxVsPosAbsEndIn310degMC","#Delta_{phi} at vertex versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);#Delta_{phi}",1000,0.,100.,deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
578  h2 = new TH2F("hResPhiAtVtxVsAngleAtAbsEnd","#Delta_{phi} at vertex versus track position at absorber end converted to degrees;angle (Deg);#Delta_{phi}",10,0.,10.,deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
580  h2 = new TH2F("hResPhiAtVtxVsMCAngle","#Delta_{phi} at vertex versus MC angle;MC angle (Deg);#Delta_{phi}",10,0.,10.,deltaPhiAtVtxNBins, deltaPhiAtVtxEdges[0], deltaPhiAtVtxEdges[1]);
582 
583  // DCA resolution
584  const Int_t deltaPDCANBins = 500;
585  const Double_t deltaPDCAEdges[2] = {0., 1000.};
586  const Double_t deltaPMCSAngEdges[2] = {-1.5, 1.5};
587 
588  h1 = new TH1F("hPDCA","p #times DCA at vertex;p #times DCA (GeV #times cm)", deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
589  fTrackerList->AddAt(h1, kPDCA);
590  h2 = new TH2F("hPDCAVsPIn23deg","p #times DCA versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);p #times DCA (GeV #times cm)",2*fNPBins,fPRange[0],fPRange[1], deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
591  fTrackerList->AddAt(h2, kPDCAVsPIn23deg);
592  h2 = new TH2F("hPDCAVsPIn310deg","p #times DCA versus p for tracks within [3,10[ degrees at absorber end;p (GeV/c);p #times DCA (GeV #times cm)",2*fNPBins,fPRange[0],fPRange[1], deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
593  fTrackerList->AddAt(h2, kPDCAVsPIn310deg);
594  h2 = new TH2F("hPDCAVsPosAbsEndIn02degMC","p #times DCA versus track position at absorber end for tracks with MC angle < 2 degrees;position (cm);p #times DCA (GeV #times cm)",1000,0.,100.,deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
596  h2 = new TH2F("hPDCAVsPosAbsEndIn23degMC","p #times DCA}versus track position at absorber end for tracks with MC angle in [2,3[ degrees;position (cm);p #times DCA (GeV #times cm)",1000,0.,100.,deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
598  h2 = new TH2F("hPDCAVsPosAbsEndIn310degMC","p #times DCA versus track position at absorber end for tracks with MC angle in [3,10[ degrees;position (cm);p #times DCA (GeV #times cm)",1000,0.,100.,deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
600  h2 = new TH2F("hPDCAVsAngleAtAbsEnd","p #times DCA versus track position at absorber end converted to degrees;angle (Deg);p #times DCA (GeV #times cm)",10,0.,10.,deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
602  h2 = new TH2F("hPDCAVsMCAngle","p #times DCA versus MC angle;MC angle (Deg);p #times DCA (GeV #times cm)",10,0.,10.,deltaPDCANBins, deltaPDCAEdges[0], deltaPDCAEdges[1]);
603  fTrackerList->AddAt(h2, kPDCAVsMCAngle);
604 
605  // MCS angular dispersion
606  h2 = new TH2F("hPMCSAngVsPIn23deg","p #times #Delta#theta_{MCS} versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);p #times #Delta#theta_{MCS} (GeV)",2*fNPBins,fPRange[0],fPRange[1], deltaPDCANBins, deltaPMCSAngEdges[0], deltaPMCSAngEdges[1]);
607  fTrackerList->AddAt(h2, kPMCSAngVsPIn23deg);
608  h2 = new TH2F("hPMCSAngVsPIn310deg","p #times #Delta#theta_{MCS} versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);p #times #Delta#theta_{MCS} (GeV)",2*fNPBins,fPRange[0],fPRange[1], deltaPDCANBins, deltaPMCSAngEdges[0], deltaPMCSAngEdges[1]);
609  fTrackerList->AddAt(h2, kPMCSAngVsPIn310deg);
610 
611  // cluster resolution
612  Int_t nCh = AliMUONConstants::NTrackingCh();
613  const Int_t clusterResNBins = 5000;
614  Double_t clusterResMaxX = 10.*fClusterMaxRes[0];
615  Double_t clusterResMaxY = 10.*fClusterMaxRes[1];
616 
617  h2 = new TH2F("hResClXVsCh", "cluster-track residual-X distribution per chamber;chamber ID;#Delta_{X} (cm)", nCh, 0.5, nCh+0.5, clusterResNBins, -clusterResMaxX, clusterResMaxX);
618  fTrackerList->AddAt(h2, kResClXVsCh);
619  h2 = new TH2F("hResClYVsCh", "cluster-track residual-Y distribution per chamber;chamber ID;#Delta_{Y} (cm)", nCh, 0.5, nCh+0.5, clusterResNBins, -clusterResMaxY, clusterResMaxY);
620  fTrackerList->AddAt(h2, kResClYVsCh);
621  h2 = new TH2F("hResClXVsDE", "cluster-track residual-X distribution per DE;DE ID;#Delta_{X} (cm)", fNDE, 0.5, fNDE+0.5, clusterResNBins, -clusterResMaxX, clusterResMaxX);
622  for (Int_t i = 1; i <= fNDE; i++) h2->GetXaxis()->SetBinLabel(i, Form("%d",fDEIds[i]));
623  fTrackerList->AddAt(h2, kResClXVsDE);
624  h2 = new TH2F("hResClYVsDE", "cluster-track residual-Y distribution per DE;DE ID;#Delta_{Y} (cm)", fNDE, 0.5, fNDE+0.5, clusterResNBins, -clusterResMaxY, clusterResMaxY);
625  for (Int_t i = 1; i <= fNDE; i++) h2->GetXaxis()->SetBinLabel(i, Form("%d",fDEIds[i]));
626  fTrackerList->AddAt(h2, kResClYVsDE);
627 
628  // Disable printout of AliMCEvent
629  AliLog::SetClassDebugLevel("AliMCEvent",-1);
630 
631  PostData(1, fCFContainer);
632  PostData(2, fTriggerList);
633  PostData(3, fTrackerList);
634 }
635 
636 //________________________________________________________________________
637 void AliAnalysisTaskMuonPerformance::UserExec(Option_t * /*option*/)
638 {
639  //
642  //
643 
644  // check that OCDB objects have been properly set
645  if (fSigmaCutTrig < 0) return;
646 
647  // Load ESD event
648  AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
649  if ( ! esd ) {
650  AliError ("ESD event not found. Nothing done!");
651  return;
652  }
653 
654  // Load MC event
655  AliMCEventHandler* mcH = static_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
656  if ( ! mcH ) {
657  AliError ("MCH event handler not found. Nothing done!");
658  return;
659  }
660 
661  // get centrality
662  Double_t centrality = esd->GetCentrality()->GetCentralityPercentileUnchecked("V0M");
663 
664  // Get reference tracks
665  AliMUONRecoCheck rc(esd,mcH);
666  AliMUONVTriggerTrackStore* triggerTrackRefStore = rc.TriggerableTracks(-1);
667  AliMUONVTrackStore* trackRefStore = rc.TrackRefs(-1);
668  AliMUONVTrackStore* reconstructibleStore = rc.ReconstructibleTracks(-1, fRequestedStationMask, fRequest2ChInSameSt45);
669 
670  Double_t containerInput[kNvars];
671  containerInput[kVarCent] = centrality;
672  AliMUONTrackParam *trackParam;
673  Double_t x1,y1,z1,slopex1,slopey1,pX1,pY1,pZ1,p1,pT1,eta1,phi1;
674  Double_t x2,y2,z2,slopex2,slopey2,pX2,pY2,pZ2,p2,pT2,eta2,phi2;
675  Double_t dPhi;
676  Double_t xAbs,yAbs,dAbs,aAbs,aMCS,aMC;
677  Double_t xDCA,yDCA,dca,pU;
678 
679  // ------ Loop over reconstructed tracks ------
680  AliESDMuonTrack *esdTrack = 0x0;
681  Int_t nMuTracks = esd->GetNumberOfMuonTracks();
682  Int_t *loCircuit = new Int_t[nMuTracks];
683  Int_t nTrgTracks = 0;
684  for (Int_t iMuTrack = 0; iMuTrack < nMuTracks; ++iMuTrack) {
685 
686  esdTrack = esd->GetMuonTrack(iMuTrack);
687 
688  // Tracker
689  AliMUONTrack* matchedTrackRef = 0x0;
690  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
691  Bool_t isValid = kFALSE;
692  if (esdTrack->ContainTrackerData()) {
693 
694  // convert ESD track to MUON track (without recomputing track parameters at each clusters)
695  AliMUONTrack muonTrack;
696  AliMUONESDInterface::ESDToMUON(*esdTrack, muonTrack, kFALSE);
697 
698  // decide whether the track is valid for efficiency calculations
699  isValid = (!fEnforceTrkCriteria || muonTrack.IsValid(fRequestedStationMask, fRequest2ChInSameSt45));
700 
701  // get the associated simulated track (discard decays)
702  Int_t mcLabel = esdTrack->GetLabel();
703  if (mcLabel >= 0 && !esdTrack->TestBit(BIT(22)))
704  matchedTrackRef = static_cast<AliMUONTrack*>(trackRefStore->FindObject(mcLabel));
705  if (matchedTrackRef && isValid) containerInput[kVarMatchMC] = static_cast<Double_t>(kTrackerOnly);
706 
707  // compute track resolution (discard not-reconstructible trackRef)
708  if (matchedTrackRef && !esdTrack->TestBit(BIT(23))) {
709 
710  // simulated track parameters at vertex
711  trackParam = matchedTrackRef->GetTrackParamAtVertex();
712  x1 = trackParam->GetNonBendingCoor();
713  y1 = trackParam->GetBendingCoor();
714  z1 = trackParam->GetZ();
715  slopex1 = trackParam->GetNonBendingSlope();
716  slopey1 = trackParam->GetBendingSlope();
717  pX1 = trackParam->Px();
718  pY1 = trackParam->Py();
719  pZ1 = trackParam->Pz();
720  p1 = trackParam->P();
721  pT1 = TMath::Sqrt(pX1*pX1 + pY1*pY1);
722  aMC = TMath::ATan(-pT1/pZ1) * TMath::RadToDeg();
723  eta1 = TMath::Log(TMath::Tan(0.5*TMath::ATan(-pT1/pZ1)));
724  phi1 = TMath::Pi()+TMath::ATan2(-pY1, -pX1);
725 
726  // reconstructed track parameters at the end of the absorber
727  AliMUONTrackParam trackParamAtAbsEnd(*((AliMUONTrackParam*)muonTrack.GetTrackParamAtCluster()->First()));
728  AliMUONTrackExtrap::ExtrapToZ(&trackParamAtAbsEnd, AliMUONConstants::AbsZEnd());
729  xAbs = trackParamAtAbsEnd.GetNonBendingCoor();
730  yAbs = trackParamAtAbsEnd.GetBendingCoor();
731  dAbs = TMath::Sqrt(xAbs*xAbs + yAbs*yAbs);
732  aAbs = TMath::ATan(-dAbs/AliMUONConstants::AbsZEnd()) * TMath::RadToDeg();
733  pX2 = trackParamAtAbsEnd.Px();
734  pY2 = trackParamAtAbsEnd.Py();
735  pZ2 = trackParamAtAbsEnd.Pz();
736  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
737  aMCS = TMath::ATan(-pT2/pZ2) * TMath::RadToDeg();
738 
739  // reconstructed track parameters at vertex
740  trackParam = muonTrack.GetTrackParamAtVertex();
741  x2 = trackParam->GetNonBendingCoor();
742  y2 = trackParam->GetBendingCoor();
743  z2 = trackParam->GetZ();
744  slopex2 = trackParam->GetNonBendingSlope();
745  slopey2 = trackParam->GetBendingSlope();
746  pX2 = trackParam->Px();
747  pY2 = trackParam->Py();
748  pZ2 = trackParam->Pz();
749  p2 = trackParam->P();
750  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
751  eta2 = TMath::Log(TMath::Tan(0.5*TMath::ATan(-pT2/pZ2)));
752  phi2 = TMath::Pi()+TMath::ATan2(-pY2, -pX2);
753 
754  // reconstructed track parameters at DCA
755  AliMUONTrackParam trackParamAtDCA(*((AliMUONTrackParam*) muonTrack.GetTrackParamAtCluster()->First()));
756  pU = trackParamAtDCA.P();
757  AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&trackParamAtDCA, z2);
758  xDCA = trackParamAtDCA.GetNonBendingCoor();
759  yDCA = trackParamAtDCA.GetBendingCoor();
760  dca = TMath::Sqrt(xDCA*xDCA + yDCA*yDCA);
761 
762  dPhi = phi2-phi1;
763  if (dPhi < -TMath::Pi()) dPhi += 2.*TMath::Pi();
764  else if (dPhi > TMath::Pi()) dPhi -= 2.*TMath::Pi();
765 
766  // fill histograms
767  static_cast<TH1*>(fTrackerList->At(kResPAtVtx))->Fill(p2-p1);
768  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsP))->Fill(p1,p2-p1);
769  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsAngleAtAbsEnd))->Fill(aAbs,p2-p1);
770  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsMCAngle))->Fill(aMC,p2-p1);
771  static_cast<TH3*>(fTrackerList->At(kResPAtVtxVsAngleAtAbsEndVsP))->Fill(p1,aAbs,p2-p1);
772  static_cast<TH2*>(fTrackerList->At(kResPtAtVtxVsPt))->Fill(pT1,pT2-pT1);
773 
774  static_cast<TH1*>(fTrackerList->At(kResSlopeXAtVtx))->Fill(slopex2-slopex1);
775  static_cast<TH1*>(fTrackerList->At(kResSlopeYAtVtx))->Fill(slopey2-slopey1);
776  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsP))->Fill(p1,slopex2-slopex1);
777  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsP))->Fill(p1,slopey2-slopey1);
778  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsAngleAtAbsEnd))->Fill(aAbs,slopex2-slopex1);
779  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsAngleAtAbsEnd))->Fill(aAbs,slopey2-slopey1);
780  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsMCAngle))->Fill(aMC,slopex2-slopex1);
781  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsMCAngle))->Fill(aMC,slopey2-slopey1);
782 
783  static_cast<TH1*>(fTrackerList->At(kResEtaAtVtx))->Fill(eta2-eta1);
784  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsP))->Fill(p1,eta2-eta1);
785  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsAngleAtAbsEnd))->Fill(aAbs,eta2-eta1);
786  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsMCAngle))->Fill(aMC,eta2-eta1);
787 
788  static_cast<TH1*>(fTrackerList->At(kResPhiAtVtx))->Fill(dPhi);
789  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsP))->Fill(p1,dPhi);
790  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsAngleAtAbsEnd))->Fill(aAbs,dPhi);
791  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsMCAngle))->Fill(aMC,dPhi);
792 
793  static_cast<TH1*>(fTrackerList->At(kPDCA))->Fill(0.5*(p2+pU)*dca);
794  static_cast<TH2*>(fTrackerList->At(kPDCAVsAngleAtAbsEnd))->Fill(aAbs,0.5*(p2+pU)*dca);
795  static_cast<TH2*>(fTrackerList->At(kPDCAVsMCAngle))->Fill(aMC,0.5*(p2+pU)*dca);
796 
797  if (aAbs > 2. && aAbs < 3.) {
798 
799  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn23deg))->Fill(p1,p2-p1);
800  static_cast<TH2*>(fTrackerList->At(kPDCAVsPIn23deg))->Fill(p1,0.5*(p2+pU)*dca);
801  static_cast<TH2*>(fTrackerList->At(kPMCSAngVsPIn23deg))->Fill(p1,0.5*(p2+pU)*(aMCS-aMC)*TMath::DegToRad());
802 
803  } else if (aAbs >= 3. && aAbs < 10.) {
804 
805  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn310deg))->Fill(p1,p2-p1);
806  static_cast<TH2*>(fTrackerList->At(kPDCAVsPIn310deg))->Fill(p1,0.5*(p2+pU)*dca);
807  static_cast<TH2*>(fTrackerList->At(kPMCSAngVsPIn310deg))->Fill(p1,0.5*(p2+pU)*(aMCS-aMC)*TMath::DegToRad());
808 
809  }
810 
811  if (aMC < 2.) {
812 
813  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn02degMC))->Fill(p1,p2-p1);
814  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,p2-p1);
815  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,slopex2-slopex1);
816  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,slopey2-slopey1);
817  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,eta2-eta1);
818  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,dPhi);
819  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn02degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
820 
821  } else if (aMC >= 2. && aMC < 3) {
822 
823  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,p2-p1);
824  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,slopex2-slopex1);
825  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,slopey2-slopey1);
826  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,eta2-eta1);
827  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,dPhi);
828  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn23degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
829 
830  } else if (aMC >= 3. && aMC < 10.) {
831 
832  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,p2-p1);
833  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,slopex2-slopex1);
834  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,slopey2-slopey1);
835  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,eta2-eta1);
836  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,dPhi);
837  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn310degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
838 
839  }
840 
841  // simulated track parameters at first cluster
842  trackParam = (AliMUONTrackParam*) matchedTrackRef->GetTrackParamAtCluster()->First();
843  x1 = trackParam->GetNonBendingCoor();
844  y1 = trackParam->GetBendingCoor();
845  z1 = trackParam->GetZ();
846  slopex1 = trackParam->GetNonBendingSlope();
847  slopey1 = trackParam->GetBendingSlope();
848  pX1 = trackParam->Px();
849  pY1 = trackParam->Py();
850  pZ1 = trackParam->Pz();
851  p1 = trackParam->P();
852  pT1 = TMath::Sqrt(pX1*pX1 + pY1*pY1);
853 
854  // reconstructed track parameters at first cluster
855  trackParam = (AliMUONTrackParam*) muonTrack.GetTrackParamAtCluster()->First();
856  x2 = trackParam->GetNonBendingCoor();
857  y2 = trackParam->GetBendingCoor();
858  z2 = trackParam->GetZ();
859  slopex2 = trackParam->GetNonBendingSlope();
860  slopey2 = trackParam->GetBendingSlope();
861  pX2 = trackParam->Px();
862  pY2 = trackParam->Py();
863  pZ2 = trackParam->Pz();
864  p2 = trackParam->P();
865  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
866 
867  // fill histograms
868  static_cast<TH1*>(fTrackerList->At(kResPAt1stCl))->Fill(p2-p1);
869  static_cast<TH2*>(fTrackerList->At(kResPAt1stClVsP))->Fill(p1,p2-p1);
870  static_cast<TH2*>(fTrackerList->At(kResPtAt1stClVsPt))->Fill(pT1,pT2-pT1);
871  static_cast<TH1*>(fTrackerList->At(kResSlopeXAt1stCl))->Fill(slopex2-slopex1);
872  static_cast<TH1*>(fTrackerList->At(kResSlopeYAt1stCl))->Fill(slopey2-slopey1);
873  static_cast<TH2*>(fTrackerList->At(kResSlopeXAt1stClVsP))->Fill(p1,slopex2-slopex1);
874  static_cast<TH2*>(fTrackerList->At(kResSlopeYAt1stClVsP))->Fill(p1,slopey2-slopey1);
875 
876  // clusters residuals
877  for (Int_t iCl1 = 0; iCl1 < muonTrack.GetNClusters(); iCl1++) {
878 
879  AliMUONVCluster* cluster1 = static_cast<AliMUONTrackParam*>(muonTrack.GetTrackParamAtCluster()->UncheckedAt(iCl1))->GetClusterPtr();
880  Int_t chId = cluster1->GetChamberId();
881  Int_t deId = cluster1->GetDetElemId();
882 
883  for (Int_t iCl2 = 0; iCl2 < matchedTrackRef->GetNClusters(); iCl2++) {
884 
885  AliMUONVCluster* cluster2 = static_cast<AliMUONTrackParam*>(matchedTrackRef->GetTrackParamAtCluster()->UncheckedAt(iCl2))->GetClusterPtr();
886 
887  if (cluster2->GetDetElemId() == deId) {
888 
889  // fill histograms
890  static_cast<TH2*>(fTrackerList->At(kResClXVsCh))->Fill(chId+1, cluster1->GetX()-cluster2->GetX());
891  static_cast<TH2*>(fTrackerList->At(kResClYVsCh))->Fill(chId+1, cluster1->GetY()-cluster2->GetY());
892  static_cast<TH2*>(fTrackerList->At(kResClXVsDE))->Fill(fDEIndices[deId], cluster1->GetX()-cluster2->GetX());
893  static_cast<TH2*>(fTrackerList->At(kResClYVsDE))->Fill(fDEIndices[deId], cluster1->GetY()-cluster2->GetY());
894 
895  break;
896 
897  }
898 
899  }
900 
901  }
902 
903  } // end if (matchedTrackRef)
904 
905  } // end if (esdTrack->ContainTrackerData())
906 
907  // look for MC trigger associated to MC track
908  AliMUONTriggerTrack* triggerTrackRef = (isValid && matchedTrackRef)
909  ? static_cast<AliMUONTriggerTrack*>(triggerTrackRefStore->FindObject(matchedTrackRef->GetUniqueID()))
910  : 0x0;
911 
912  // Trigger
913  AliMUONTriggerTrack* matchedTrigTrackRef = 0x0;
914  containerInput[kVarDupliTrg] = 0.;
915  if (esdTrack->ContainTriggerData()) {
916 
917  // check if this track is not already accounted for and record it if not
918  Bool_t trackExist = kFALSE;
919  for (Int_t i=0; i<nTrgTracks; i++) if (esdTrack->LoCircuit() == loCircuit[i]) trackExist = kTRUE;
920  if (trackExist) containerInput[kVarDupliTrg] = 1.;
921  else loCircuit[nTrgTracks++] = esdTrack->LoCircuit();
922 
923  // Convert ESD track to trigger track
924  AliMUONLocalTrigger locTrg;
925  AliMUONESDInterface::ESDToMUON(*esdTrack, locTrg);
926  AliMUONTriggerTrack trigTrack;
927  rc.TriggerToTrack(locTrg, trigTrack);
928 
929  // try to match the trigger track with the same MC track if any
930  if (triggerTrackRef && trigTrack.Match(*triggerTrackRef, fSigmaCutTrig)) matchedTrigTrackRef = triggerTrackRef;
931  if (matchedTrigTrackRef) containerInput[kVarMatchMC] = static_cast<Double_t>(kMatchedSame);
932  else { // or try to match with any triggerable track
933  matchedTrigTrackRef = rc.FindCompatibleTrack(trigTrack, *triggerTrackRefStore, fSigmaCutTrig);
934  if (matchedTrigTrackRef) {
935  if (isValid && matchedTrackRef) {
936  containerInput[kVarMatchMC] = static_cast<Double_t>(kMatchedDiff);
937  if (!fMCTrigLevelFromMatchTrk) triggerTrackRef = matchedTrigTrackRef;
938  } else {
939  containerInput[kVarMatchMC] = static_cast<Double_t>(kTriggerOnly);
940  triggerTrackRef = matchedTrigTrackRef;
941  }
942  }
943  }
944 
945  // fill histograms
946  if (matchedTrigTrackRef) {
947  static_cast<TH1*>(fTriggerList->At(kResTrigX11))->Fill(trigTrack.GetX11() - matchedTrigTrackRef->GetX11());
948  static_cast<TH1*>(fTriggerList->At(kResTrigY11))->Fill(trigTrack.GetY11() - matchedTrigTrackRef->GetY11());
949  static_cast<TH1*>(fTriggerList->At(kResTrigSlopeY))->Fill(trigTrack.GetSlopeY() - matchedTrigTrackRef->GetSlopeY());
950  }
951 
952  }
953 
954  // fill container
955  if (triggerTrackRef) {
956  if (triggerTrackRef->GetPtCutLevel() == 0) containerInput[kVarMCTrigger] = static_cast<Double_t>(kOtherTrig);
957  else if (triggerTrackRef->GetPtCutLevel() == 1) containerInput[kVarMCTrigger] = static_cast<Double_t>(kAllPtTrig);
958  else if (triggerTrackRef->GetPtCutLevel() == 2) containerInput[kVarMCTrigger] = static_cast<Double_t>(kLowPtTrig);
959  else if (triggerTrackRef->GetPtCutLevel() == 3) containerInput[kVarMCTrigger] = static_cast<Double_t>(kHighPtTrig);
960  } else containerInput[kVarMCTrigger] = kNoMatchTrig;
961 
962  // get MC particle ID
963  Int_t mcID = -1;
964  if (isValid && matchedTrackRef) mcID = static_cast<Int_t>(matchedTrackRef->GetUniqueID());
965  else if (matchedTrigTrackRef) mcID = static_cast<Int_t>(matchedTrigTrackRef->GetUniqueID());
966 
967  // fill particle container
968  FillContainerInfoReco(containerInput, esdTrack, isValid, mcID);
969  fCFContainer->Fill(containerInput, kStepReconstructed);
970 
971  }
972 
973  // clean memory
974  delete[] loCircuit;
975  containerInput[kVarDupliTrg] = 0.;
976 
977  // ------ Loop over reconstructible tracks ------
978  AliMUONTrack* trackRef = 0x0;
979  TIter next(reconstructibleStore->CreateIterator());
980  while ((trackRef = static_cast<AliMUONTrack*>(next()))) {
981 
982  // find the corresponding triggerable track if any
983  UInt_t mcID = trackRef->GetUniqueID();
984  AliMUONTriggerTrack* trigTrackRef = static_cast<AliMUONTriggerTrack*>(triggerTrackRefStore->FindObject(mcID));
985 
986  // fill particle container
987  FillContainerInfoMC(containerInput, static_cast<AliMCParticle*>(fMCEvent->GetTrack(static_cast<Int_t>(mcID))));
988  containerInput[kVarHasTracker] = 1.;
989  if (trigTrackRef) {
990  if (trigTrackRef->GetPtCutLevel() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kOtherTrig);
991  else if (trigTrackRef->GetPtCutLevel() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
992  else if (trigTrackRef->GetPtCutLevel() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
993  else if (trigTrackRef->GetPtCutLevel() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
994  } else containerInput[kVarTrigger] = static_cast<Double_t>(kNoMatchTrig);
995  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
996  containerInput[kVarMCTrigger] = static_cast<Double_t>(kNoMatchTrig);
997  fCFContainer->Fill(containerInput, kStepGeneratedMC);
998 
999  }
1000 
1001  // ------ Loop over triggerable ghosts ------
1002  AliMUONTriggerTrack* trigTrackRef = 0x0;
1003  TIter nextTrig(triggerTrackRefStore->CreateIterator());
1004  while ((trigTrackRef = static_cast<AliMUONTriggerTrack*>(nextTrig()))) {
1005 
1006  // discard tracks also reconstructible
1007  UInt_t mcID = trigTrackRef->GetUniqueID();
1008  if (reconstructibleStore->FindObject(mcID)) continue;
1009 
1010  // fill particle container
1011  FillContainerInfoMC(containerInput, static_cast<AliMCParticle*>(fMCEvent->GetTrack(static_cast<Int_t>(mcID))));
1012  containerInput[kVarHasTracker] = 0.;
1013  if (trigTrackRef->GetPtCutLevel() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kOtherTrig);
1014  else if (trigTrackRef->GetPtCutLevel() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
1015  else if (trigTrackRef->GetPtCutLevel() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
1016  else if (trigTrackRef->GetPtCutLevel() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
1017  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
1018  containerInput[kVarMCTrigger] = static_cast<Double_t>(kNoMatchTrig);
1019  fCFContainer->Fill(containerInput, kStepGeneratedMC);
1020 
1021  }
1022 
1023  // Post final data
1024  PostData(1, fCFContainer);
1025  PostData(2, fTriggerList);
1026  PostData(3, fTrackerList);
1027 }
1028 
1029 //________________________________________________________________________
1031 {
1032  //
1034  //
1035 
1036  // do it only locally
1037  //if (gROOT->IsBatch()) return;
1038 
1039  // get output containers
1040  fCFContainer = dynamic_cast<AliCFContainer*>(GetOutputData(1));
1041  fTriggerList = dynamic_cast<TObjArray*>(GetOutputData(2));
1042  fTrackerList = dynamic_cast<TObjArray*>(GetOutputData(3));
1043  if (!fCFContainer || !fTriggerList || !fTrackerList) {
1044  AliWarning("Output containers not found: summary histograms are not created");
1045  return;
1046  }
1047 
1048  // --- compute efficiencies ---
1049  fEfficiencyList = new TObjArray(100);
1050  fEfficiencyList->SetOwner();
1051 
1052  TObjArray* effAnyPt = new TObjArray(100);
1053  effAnyPt->SetName("effAnyPt");
1054  effAnyPt->SetOwner();
1055  fEfficiencyList->AddLast(effAnyPt);
1056 
1057  TObjArray* effAllPt = new TObjArray(100);
1058  effAllPt->SetName("effAllPt");
1059  effAllPt->SetOwner();
1060  fEfficiencyList->AddLast(effAllPt);
1061 
1062  TObjArray* effLowPt = new TObjArray(100);
1063  effLowPt->SetName("effLowPt");
1064  effLowPt->SetOwner();
1065  fEfficiencyList->AddLast(effLowPt);
1066 
1067  TObjArray* effHighPt = new TObjArray(100);
1068  effHighPt->SetName("effHighPt");
1069  effHighPt->SetOwner();
1070  fEfficiencyList->AddLast(effHighPt);
1071 
1072  TObjArray* notTrgable = new TObjArray(100);
1073  notTrgable->SetName("notTrgable");
1074  notTrgable->SetOwner();
1075  fEfficiencyList->AddLast(notTrgable);
1076 
1077  TObjArray* trgableNoPtOnly = new TObjArray(100);
1078  trgableNoPtOnly->SetName("trgableNoPtOnly");
1079  trgableNoPtOnly->SetOwner();
1080  fEfficiencyList->AddLast(trgableNoPtOnly);
1081 
1082  TObjArray* trgableAPtOnly = new TObjArray(100);
1083  trgableAPtOnly->SetName("trgableAPtOnly");
1084  trgableAPtOnly->SetOwner();
1085  fEfficiencyList->AddLast(trgableAPtOnly);
1086 
1087  TObjArray* trgableLPtOnly = new TObjArray(100);
1088  trgableLPtOnly->SetName("trgableLPtOnly");
1089  trgableLPtOnly->SetOwner();
1090  fEfficiencyList->AddLast(trgableLPtOnly);
1091 
1092  TObjArray* trgableHPtOnly = new TObjArray(100);
1093  trgableHPtOnly->SetName("trgableHPtOnly");
1094  trgableHPtOnly->SetOwner();
1095  fEfficiencyList->AddLast(trgableHPtOnly);
1096 
1097  AliCFEffGrid* efficiency = new AliCFEffGrid("eff","",*fCFContainer);
1098  efficiency->CalculateEfficiency(kStepReconstructed, kStepGeneratedMC);
1099  Double_t totalEff = 0., totalEffErr = 0.;
1100  Int_t sumEffBin = 0;
1101 
1102  // add histogram summarizing global efficiencies
1103  TH1D* effSummary = new TH1D("effSummary", "Efficiency summary", 1, 0., 0.);
1104  effSummary->GetYaxis()->SetTitle("Efficiency");
1105  fEfficiencyList->AddLast(effSummary);
1106 
1107  // ------ Tracker only ------
1108 
1109  // Tracker efficiency using all reconstructed tracks
1110  efficiency->GetNum()->SetRangeUser(kVarHasTracker, 1., 1.);
1111  efficiency->GetDen()->SetRangeUser(kVarHasTracker, 1., 1.);
1112  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1113  efficiency->GetDen()->GetAxis(kVarTrigger)->SetRange();
1114  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1115  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1116  FillEffHistos(efficiency, "trackerTracks", fEfficiencyList);
1117  GetEfficiency(efficiency, totalEff, totalEffErr);
1118  effSummary->Fill("Tracker_all", totalEff);
1119  effSummary->SetBinError(++sumEffBin, totalEffErr);
1120  printf("Tracker efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1121 
1122  // Tracker efficiency using tracks matched with reconstructible ones
1123  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1124  FillEffHistos(efficiency, "trackerTracksMatchMC", fEfficiencyList);
1125  GetEfficiency(efficiency, totalEff, totalEffErr);
1126  effSummary->Fill("Tracker_MCId", totalEff);
1127  effSummary->SetBinError(++sumEffBin, totalEffErr);
1128  printf("Tracker efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1129 
1130  // ------ Tracker matched with any trigger ------
1131 
1132  // Matched efficiency using all reconstructed tracks
1133  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1134  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kHighPtTrig);
1135  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1136  FillEffHistos(efficiency, "matchedTracks", effAnyPt);
1137  GetEfficiency(efficiency, totalEff, totalEffErr);
1138  effSummary->Fill("Matched_all", totalEff);
1139  effSummary->SetBinError(++sumEffBin, totalEffErr);
1140  printf("Matched efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1141 
1142  // Matched efficiency using tracks matched with reconstructible ones, triggerable or not
1143  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1144  FillEffHistos(efficiency, "matchedTracksMatchMC", effAnyPt);
1145  GetEfficiency(efficiency, totalEff, totalEffErr);
1146  effSummary->Fill("Matched_MCId", totalEff);
1147  effSummary->SetBinError(++sumEffBin, totalEffErr);
1148  printf("Matched efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1149 
1150  // Matched efficiency using tracks matched with reconstructible ones triggerable
1151  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kHighPtTrig);
1152  FillEffHistos(efficiency, "matchedTracksMatchMCAnypt", effAnyPt);
1153  GetEfficiency(efficiency, totalEff, totalEffErr);
1154  effSummary->Fill("Matched_MCIdAnypt", totalEff);
1155  effSummary->SetBinError(++sumEffBin, totalEffErr);
1156  printf("Matched efficiency using reconstructed tracks matching MC-anyPt = %f +- %f\n", totalEff, totalEffErr);
1157 
1158  // Matched efficiency using tracks matched with reconstructible ones not triggerable
1159  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1160  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effAnyPt);
1161  GetEfficiency(efficiency, totalEff, totalEffErr);
1162  effSummary->Fill("Matched_MCIdNoTrig", totalEff);
1163  effSummary->SetBinError(++sumEffBin, totalEffErr);
1164  printf("Matched efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1165 
1166  // Matched efficiency using tracks matched with same reconstructible & triggerable ones
1167  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1168  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1169  FillEffHistos(efficiency, "matchedTracksMatchSameMC", effAnyPt);
1170  GetEfficiency(efficiency, totalEff, totalEffErr);
1171  effSummary->Fill("Matched_SameMCId", totalEff);
1172  effSummary->SetBinError(++sumEffBin, totalEffErr);
1173  printf("Matched efficiency using reconstructed tracks matching same MC = %f +- %f\n", totalEff, totalEffErr);
1174 
1175  // Matched efficiency using tracks matched with different reconstructible & triggerable ones
1176  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1177  FillEffHistos(efficiency, "matchedTracksMatchDiffMC", effAnyPt);
1178  GetEfficiency(efficiency, totalEff, totalEffErr);
1179  effSummary->Fill("Matched_DiffMCId", totalEff);
1180  effSummary->SetBinError(++sumEffBin, totalEffErr);
1181  printf("Matched efficiency using reconstructed tracks matching different MC = %f +- %f\n", totalEff, totalEffErr);
1182 
1183  // Matched efficiency using tracks matched with reconstructible ones only
1184  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1185  FillEffHistos(efficiency, "matchedTracksMatchTrkMC", effAnyPt);
1186  GetEfficiency(efficiency, totalEff, totalEffErr);
1187  effSummary->Fill("Matched_TrkMCId", totalEff);
1188  effSummary->SetBinError(++sumEffBin, totalEffErr);
1189  printf("Matched efficiency using reconstructed tracks matching tracker MC = %f +- %f\n", totalEff, totalEffErr);
1190 
1191  // ------ Tracker matched with all pt trigger ------
1192 
1193  // Matched all pt efficiency using all reconstructed tracks
1194  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1195  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1196  FillEffHistos(efficiency, "matchedTracks", effAllPt);
1197  GetEfficiency(efficiency, totalEff, totalEffErr);
1198  printf("Matched Apt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1199 
1200  // Matched all pt efficiency using tracks matched with reconstructible ones, triggerable or not
1201  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1202  FillEffHistos(efficiency, "matchedTracksMatchMC", effAllPt);
1203  GetEfficiency(efficiency, totalEff, totalEffErr);
1204  printf("Matched Apt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1205 
1206  // Matched all pt efficiency using tracks matched with reconstructible ones triggerable Apt
1207  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1208  FillEffHistos(efficiency, "matchedTracksMatchMCApt", effAllPt);
1209  GetEfficiency(efficiency, totalEff, totalEffErr);
1210  printf("Matched Apt efficiency using reconstructed tracks matching MC-Apt = %f +- %f\n", totalEff, totalEffErr);
1211 
1212  // Matched all pt efficiency using tracks matched with reconstructible ones triggerable other pt
1213  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1214  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effAllPt);
1215  GetEfficiency(efficiency, totalEff, totalEffErr);
1216  printf("Matched Apt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1217 
1218  // Matched all pt efficiency using tracks matched with reconstructible ones not triggerable
1219  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1220  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effAllPt);
1221  GetEfficiency(efficiency, totalEff, totalEffErr);
1222  printf("Matched Apt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1223 
1224  // Matched all pt efficiency using tracks matched with same reconstructible & triggerable ones (all pt MC trig)
1225  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1226  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1227  FillEffHistos(efficiency, "matchedTracksMatchSameMCApt", effAllPt);
1228 
1229  // Matched all pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1230  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1231  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effAllPt);
1232 
1233  // Matched all pt efficiency using tracks matched with different reconstructible & triggerable ones (all pt MC trig)
1234  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1235  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1236  FillEffHistos(efficiency, "matchedTracksMatchDiffMCApt", effAllPt);
1237 
1238  // Matched all pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1239  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1240  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effAllPt);
1241 
1242  // Matched efficiency using tracks matched with reconstructible ones only (all pt MC trig)
1243  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1244  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1245  FillEffHistos(efficiency, "matchedTracksMatchTrkMCApt", effAllPt);
1246 
1247  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1248  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kOtherTrig);
1249  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effAllPt);
1250 
1251  // ------ Tracker matched with low pt trigger ------
1252 
1253  // Matched low pt efficiency using all reconstructed tracks
1254  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1255  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1256  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1257  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1258  FillEffHistos(efficiency, "matchedTracks", effLowPt);
1259  GetEfficiency(efficiency, totalEff, totalEffErr);
1260  printf("Matched Lpt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1261 
1262  // Matched low pt efficiency using tracks matched with reconstructible ones, triggerable or not
1263  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1264  FillEffHistos(efficiency, "matchedTracksMatchMC", effLowPt);
1265  GetEfficiency(efficiency, totalEff, totalEffErr);
1266  printf("Matched Lpt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1267 
1268  // Matched low pt efficiency using tracks matched with reconstructible ones triggerable Lpt
1269  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1270  FillEffHistos(efficiency, "matchedTracksMatchMCLpt", effLowPt);
1271  GetEfficiency(efficiency, totalEff, totalEffErr);
1272  printf("Matched Lpt efficiency using reconstructed tracks matching MC-Lpt = %f +- %f\n", totalEff, totalEffErr);
1273 
1274  // Matched low pt efficiency using tracks matched with reconstructible ones triggerable other pt
1275  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1276  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effLowPt);
1277  GetEfficiency(efficiency, totalEff, totalEffErr);
1278  printf("Matched Lpt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1279 
1280  // Matched low pt efficiency using tracks matched with reconstructible ones not triggerable
1281  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1282  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effLowPt);
1283  GetEfficiency(efficiency, totalEff, totalEffErr);
1284  printf("Matched Lpt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1285 
1286  // Matched low pt efficiency using tracks matched with same reconstructible & triggerable ones (low pt MC trig)
1287  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1288  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1289  FillEffHistos(efficiency, "matchedTracksMatchSameMCLpt", effLowPt);
1290 
1291  // Matched low pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1292  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1293  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effLowPt);
1294 
1295  // Matched low pt efficiency using tracks matched with different reconstructible & triggerable ones (low pt MC trig)
1296  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1297  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1298  FillEffHistos(efficiency, "matchedTracksMatchDiffMCLpt", effLowPt);
1299 
1300  // Matched low pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1301  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1302  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effLowPt);
1303 
1304  // Matched efficiency using tracks matched with reconstructible ones only (low pt MC trig)
1305  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1306  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1307  FillEffHistos(efficiency, "matchedTracksMatchTrkMCLpt", effLowPt);
1308 
1309  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1310  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kAllPtTrig);
1311  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effLowPt);
1312 
1313  // ------ Tracker matched with high pt trigger ------
1314 
1315  // Matched high pt efficiency using all reconstructed tracks
1316  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1317  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1318  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1319  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1320  FillEffHistos(efficiency, "matchedTracks", effHighPt);
1321  GetEfficiency(efficiency, totalEff, totalEffErr);
1322  printf("Matched Hpt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1323 
1324  // Matched high pt efficiency using tracks matched with reconstructible ones, triggerable or not
1325  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1326  FillEffHistos(efficiency, "matchedTracksMatchMC", effHighPt);
1327  GetEfficiency(efficiency, totalEff, totalEffErr);
1328  printf("Matched Hpt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1329 
1330  // Matched high pt efficiency using tracks matched with reconstructible ones triggerable Hpt
1331  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1332  FillEffHistos(efficiency, "matchedTracksMatchMCHpt", effHighPt);
1333  GetEfficiency(efficiency, totalEff, totalEffErr);
1334  printf("Matched Hpt efficiency using reconstructed tracks matching MC-Hpt = %f +- %f\n", totalEff, totalEffErr);
1335 
1336  // Matched high pt efficiency using tracks matched with reconstructible ones triggerable other pt
1337  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1338  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effHighPt);
1339  GetEfficiency(efficiency, totalEff, totalEffErr);
1340  printf("Matched Hpt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1341 
1342  // Matched high pt efficiency using tracks matched with reconstructible ones not triggerable
1343  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1344  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effHighPt);
1345  GetEfficiency(efficiency, totalEff, totalEffErr);
1346  printf("Matched Hpt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1347 
1348  // Matched high pt efficiency using tracks matched with same reconstructible & triggerable ones (high pt MC trig)
1349  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1350  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1351  FillEffHistos(efficiency, "matchedTracksMatchSameMCHpt", effHighPt);
1352 
1353  // Matched high pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1354  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1355  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effHighPt);
1356 
1357  // Matched high pt efficiency using tracks matched with different reconstructible & triggerable ones (high pt MC trig)
1358  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1359  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1360  FillEffHistos(efficiency, "matchedTracksMatchDiffMCHpt", effHighPt);
1361 
1362  // Matched high pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1363  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1364  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effHighPt);
1365 
1366  // Matched efficiency using tracks matched with reconstructible ones only (high pt MC trig)
1367  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1368  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1369  FillEffHistos(efficiency, "matchedTracksMatchTrkMCHpt", effHighPt);
1370 
1371  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1372  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kLowPtTrig);
1373  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effHighPt);
1374 
1375  // ------ Trigger only ------
1376 
1377  // Trigger efficiency using all reconstructed tracks
1378  efficiency->GetNum()->GetAxis(kVarHasTracker)->SetRange();
1379  efficiency->GetDen()->GetAxis(kVarHasTracker)->SetRange();
1380  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1381  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kHighPtTrig);
1382  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1383  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1384  efficiency->GetNum()->SetRangeUser(kVarDupliTrg, 0., 0.);
1385  FillEffHistos(efficiency, "triggerTracks", effAnyPt);
1386  GetEfficiency(efficiency, totalEff, totalEffErr);
1387  effSummary->Fill("Trigger_all", totalEff);
1388  effSummary->SetBinError(++sumEffBin, totalEffErr);
1389  printf("Trigger efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1390 
1391  // Trigger efficiency using tracks matched with triggerable ones
1392  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1393  FillEffHistos(efficiency, "triggerTracksMatchMC", effAnyPt);
1394  GetEfficiency(efficiency, totalEff, totalEffErr);
1395  effSummary->Fill("Trigger_MCId", totalEff);
1396  effSummary->SetBinError(++sumEffBin, totalEffErr);
1397  printf("Trigger efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1398 
1399  // ------ All pt trigger only ------
1400 
1401  // All pt trigger efficiency using all reconstructed tracks
1402  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1403  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1404  FillEffHistos(efficiency, "triggerTracks", effAllPt);
1405 
1406  // All pt trigger efficiency using tracks matched with all pt triggerable ones
1407  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1408  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1409  FillEffHistos(efficiency, "triggerTracksMatchMCApt", effAllPt);
1410 
1411  // All pt trigger efficiency using tracks matched with other triggerable ones
1412  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1413  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effAllPt);
1414 
1415  // ------ Low pt trigger only ------
1416 
1417  // Low pt trigger efficiency using all reconstructed tracks
1418  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1419  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1420  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1421  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1422  FillEffHistos(efficiency, "triggerTracks", effLowPt);
1423 
1424  // Low pt trigger efficiency using tracks matched with Low pt triggerable ones
1425  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1426  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1427  FillEffHistos(efficiency, "triggerTracksMatchMCLpt", effLowPt);
1428 
1429  // Low pt trigger efficiency using tracks matched with other triggerable ones
1430  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1431  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effLowPt);
1432 
1433  // ------ High pt trigger only ------
1434 
1435  // High pt trigger efficiency using all reconstructed tracks
1436  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1437  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1438  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1439  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1440  FillEffHistos(efficiency, "triggerTracks", effHighPt);
1441 
1442  // High pt trigger efficiency using tracks matched with High pt triggerable ones
1443  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1444  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1445  FillEffHistos(efficiency, "triggerTracksMatchMCHpt", effHighPt);
1446 
1447  // All pt trigger efficiency using tracks matched with other triggerable ones
1448  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1449  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effHighPt);
1450 
1451  // ------ Tracker reconstructible not triggerable ------
1452 
1453  // all tracker tracks
1454  efficiency->GetNum()->SetRangeUser(kVarHasTracker, 1., 1.);
1455  efficiency->GetDen()->SetRangeUser(kVarHasTracker, 1., 1.);
1456  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1457  efficiency->GetDen()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1458  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1459  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1460  efficiency->GetNum()->GetAxis(kVarDupliTrg)->SetRange();
1461  FillEffHistos(efficiency, "allTracks", notTrgable);
1462 
1463  // tracker not matched
1464  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1465  FillEffHistos(efficiency, "notMatched", notTrgable);
1466 
1467  // tracker matched with all pt trigger
1468  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1469  FillEffHistos(efficiency, "MatchedApt", notTrgable);
1470 
1471  // tracker matched with low pt trigger
1472  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1473  FillEffHistos(efficiency, "MatchedLpt", notTrgable);
1474 
1475  // tracker matched with high pt trigger
1476  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1477  FillEffHistos(efficiency, "MatchedHpt", notTrgable);
1478 
1479  // ------ Tracker reconstructible triggerable no pt ------
1480 
1481  // all tracker tracks
1482  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1483  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kOtherTrig);
1484  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1485  FillEffHistos(efficiency, "allTracks", trgableNoPtOnly);
1486 
1487  // tracker not matched
1488  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1489  FillEffHistos(efficiency, "notMatched", trgableNoPtOnly);
1490 
1491  // tracker matched with all pt trigger
1492  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1493  FillEffHistos(efficiency, "MatchedApt", trgableNoPtOnly);
1494 
1495  // tracker matched with low pt trigger
1496  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1497  FillEffHistos(efficiency, "MatchedLpt", trgableNoPtOnly);
1498 
1499  // tracker matched with high pt trigger
1500  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1501  FillEffHistos(efficiency, "MatchedHpt", trgableNoPtOnly);
1502 
1503  // ------ Tracker reconstructible triggerable Apt ------
1504 
1505  // all tracker tracks
1506  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1507  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1508  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kAllPtTrig);
1509  FillEffHistos(efficiency, "allTracks", trgableAPtOnly);
1510 
1511  // tracker not matched
1512  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1513  FillEffHistos(efficiency, "notMatched", trgableAPtOnly);
1514 
1515  // tracker matched with all pt trigger
1516  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1517  FillEffHistos(efficiency, "MatchedApt", trgableAPtOnly);
1518 
1519  // tracker matched with low pt trigger
1520  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1521  FillEffHistos(efficiency, "MatchedLpt", trgableAPtOnly);
1522 
1523  // tracker matched with high pt trigger
1524  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1525  FillEffHistos(efficiency, "MatchedHpt", trgableAPtOnly);
1526 
1527  // ------ Tracker reconstructible triggerable Lpt ------
1528 
1529  // all tracker tracks
1530  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1531  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1532  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kLowPtTrig);
1533  FillEffHistos(efficiency, "allTracks", trgableLPtOnly);
1534 
1535  // tracker not matched
1536  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1537  FillEffHistos(efficiency, "notMatched", trgableLPtOnly);
1538 
1539  // tracker matched with all pt trigger
1540  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1541  FillEffHistos(efficiency, "MatchedApt", trgableLPtOnly);
1542 
1543  // tracker matched with low pt trigger
1544  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1545  FillEffHistos(efficiency, "MatchedLpt", trgableLPtOnly);
1546 
1547  // tracker matched with high pt trigger
1548  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1549  FillEffHistos(efficiency, "MatchedHpt", trgableLPtOnly);
1550 
1551  // ------ Tracker reconstructible triggerable Hpt ------
1552 
1553  // all tracker tracks
1554  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1555  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1556  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1557  FillEffHistos(efficiency, "allTracks", trgableHPtOnly);
1558 
1559  // tracker not matched
1560  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1561  FillEffHistos(efficiency, "notMatched", trgableHPtOnly);
1562 
1563  // tracker matched with all pt trigger
1564  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1565  FillEffHistos(efficiency, "MatchedApt", trgableHPtOnly);
1566 
1567  // tracker matched with low pt trigger
1568  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1569  FillEffHistos(efficiency, "MatchedLpt", trgableHPtOnly);
1570 
1571  // tracker matched with high pt trigger
1572  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1573  FillEffHistos(efficiency, "MatchedHpt", trgableHPtOnly);
1574 
1575  // ------ reset ranges before saving CF containers ------
1576  efficiency->GetNum()->GetAxis(kVarHasTracker)->SetRange();
1577  efficiency->GetDen()->GetAxis(kVarHasTracker)->SetRange();
1578  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1579  efficiency->GetDen()->GetAxis(kVarTrigger)->SetRange();
1580  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1581  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1582  efficiency->GetNum()->GetAxis(kVarDupliTrg)->SetRange();
1583 
1584  // ------ Plot summary ------
1585 
1586  // plot histogram summarizing global efficiencies
1587  TCanvas* cEffSummary = new TCanvas("cEffSummary","Efficiency summary",20,20,310,310);
1588  cEffSummary->SetFillColor(10); cEffSummary->SetHighLightColor(10);
1589  cEffSummary->SetLeftMargin(0.15); cEffSummary->SetBottomMargin(0.15);
1590  effSummary->DrawCopy("etext");
1591 
1592  // --- plot trigger resolution ---
1593  TCanvas* cTriggerResolution = new TCanvas("cTriggerResolution","Trigger resolution",10,10,310,310);
1594  cTriggerResolution->SetFillColor(10); cTriggerResolution->SetHighLightColor(10);
1595  cTriggerResolution->SetLeftMargin(0.15); cTriggerResolution->SetBottomMargin(0.15);
1596  cTriggerResolution->Divide(2,2);
1597  cTriggerResolution->cd(1);
1598  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigX11))->DrawCopy();
1599  cTriggerResolution->cd(2);
1600  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigY11))->DrawCopy();
1601  cTriggerResolution->cd(3);
1602  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigSlopeY))->DrawCopy();
1603 
1604  // --- compute momentum resolution at vertex ---
1605  fPAtVtxList = new TObjArray(100);
1606  fPAtVtxList->SetOwner();
1607 
1608  // define graphs
1609  TGraphAsymmErrors* gMeanResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1610  gMeanResPAtVtxVsP->SetName("gMeanResPAtVtxVsP");
1611  gMeanResPAtVtxVsP->SetTitle("<#Delta_{p}> at vertex versus p;p (GeV/c);<#Delta_{p}> (GeV/c)");
1612  fPAtVtxList->AddAt(gMeanResPAtVtxVsP, kMeanResPAtVtxVsP);
1613  TGraphAsymmErrors* gMostProbResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1614  gMostProbResPAtVtxVsP->SetName("gMostProbResPAtVtxVsP");
1615  gMostProbResPAtVtxVsP->SetTitle("Most probable #Delta_{p} at vertex versus p;p (GeV/c);Most prob. #Delta_{p} (GeV/c)");
1616  fPAtVtxList->AddAt(gMostProbResPAtVtxVsP, kMostProbResPAtVtxVsP);
1617  TGraphAsymmErrors* gSigmaResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1618  gSigmaResPAtVtxVsP->SetName("gSigmaResPAtVtxVsP");
1619  gSigmaResPAtVtxVsP->SetTitle("#sigma_{p}/p at vertex versus p;p (GeV/c);#sigma_{p}/p (%)");
1620  fPAtVtxList->AddAt(gSigmaResPAtVtxVsP, kSigmaResPAtVtxVsP);
1621 
1622  // fit histo and fill graphs
1623  TH2* h = static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsP));
1624  FitLandauGausResVsP(h, "momentum residuals at vertex", gMeanResPAtVtxVsP, gMostProbResPAtVtxVsP, gSigmaResPAtVtxVsP);
1625 
1626  // convert resolution into relative resolution
1627  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
1628  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
1629  Double_t x,y;
1630  gSigmaResPAtVtxVsP->GetPoint(i/rebinFactorX-1, x, y);
1631  gSigmaResPAtVtxVsP->SetPoint(i/rebinFactorX-1, x, 100.*y/x);
1632  gSigmaResPAtVtxVsP->SetPointEYlow(i/rebinFactorX-1, 100.*gSigmaResPAtVtxVsP->GetErrorYlow(i/rebinFactorX-1)/x);
1633  gSigmaResPAtVtxVsP->SetPointEYhigh(i/rebinFactorX-1, 100.*gSigmaResPAtVtxVsP->GetErrorYhigh(i/rebinFactorX-1)/x);
1634  }
1635 
1636  // --- compute momentum resolution at first cluster ---
1637  fPAt1stClList = new TObjArray(100);
1638  fPAt1stClList->SetOwner();
1639 
1640  // define graphs
1641  TGraphAsymmErrors* gMeanResPAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1642  gMeanResPAt1stClVsP->SetName("gMeanResPAt1stClVsP");
1643  gMeanResPAt1stClVsP->SetTitle("<#Delta_{p}> at first cluster versus p;p (GeV/c);<#Delta_{p}> (GeV/c)");
1644  fPAt1stClList->AddAt(gMeanResPAt1stClVsP, kMeanResPAt1stClVsP);
1645  TGraphAsymmErrors* gSigmaResPAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1646  gSigmaResPAt1stClVsP->SetName("gSigmaResPAt1stClVsP");
1647  gSigmaResPAt1stClVsP->SetTitle("#sigma_{p}/p at first cluster versus p;p (GeV/c);#sigma_{p}/p (%)");
1648  fPAt1stClList->AddAt(gSigmaResPAt1stClVsP, kSigmaResPAt1stClVsP);
1649 
1650  // fit histo and fill graphs
1651  h = static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAt1stClVsP));
1652  FitGausResVsMom(h, 0., 1., "momentum residuals at first cluster", gMeanResPAt1stClVsP, gSigmaResPAt1stClVsP);
1653 
1654  // convert resolution into relative resolution
1655  rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
1656  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
1657  Double_t x,y;
1658  gSigmaResPAt1stClVsP->GetPoint(i/rebinFactorX-1, x, y);
1659  gSigmaResPAt1stClVsP->SetPoint(i/rebinFactorX-1, x, 100.*y/x);
1660  gSigmaResPAt1stClVsP->SetPointEYlow(i/rebinFactorX-1, 100.*gSigmaResPAt1stClVsP->GetErrorYlow(i/rebinFactorX-1)/x);
1661  gSigmaResPAt1stClVsP->SetPointEYhigh(i/rebinFactorX-1, 100.*gSigmaResPAt1stClVsP->GetErrorYhigh(i/rebinFactorX-1)/x);
1662  }
1663 
1664  // --- compute slope resolution at vertex ---
1665  fSlopeAtVtxList = new TObjArray(100);
1666  fSlopeAtVtxList->SetOwner();
1667 
1668  // define graphs
1669  TGraphAsymmErrors* gMeanResSlopeXAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1670  gMeanResSlopeXAtVtxVsP->SetName("gMeanResSlopeXAtVtxVsP");
1671  gMeanResSlopeXAtVtxVsP->SetTitle("<#Delta_{slope_{X}}> at vertex versus p;p (GeV/c);<#Delta_{slope_{X}}>");
1672  fSlopeAtVtxList->AddAt(gMeanResSlopeXAtVtxVsP, kMeanResSlopeXAtVtxVsP);
1673  TGraphAsymmErrors* gMeanResSlopeYAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1674  gMeanResSlopeYAtVtxVsP->SetName("gMeanResSlopeYAtVtxVsP");
1675  gMeanResSlopeYAtVtxVsP->SetTitle("<#Delta_{slope_{Y}}> at vertex versus p;p (GeV/c);<#Delta_{slope_{Y}}>");
1676  fSlopeAtVtxList->AddAt(gMeanResSlopeYAtVtxVsP, kMeanResSlopeYAtVtxVsP);
1677  TGraphAsymmErrors* gSigmaResSlopeXAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1678  gSigmaResSlopeXAtVtxVsP->SetName("gSigmaResSlopeXAtVtxVsP");
1679  gSigmaResSlopeXAtVtxVsP->SetTitle("#sigma_{slope_{X}} at vertex versus p;p (GeV/c);#sigma_{slope_{X}}");
1680  fSlopeAtVtxList->AddAt(gSigmaResSlopeXAtVtxVsP, kSigmaResSlopeXAtVtxVsP);
1681  TGraphAsymmErrors* gSigmaResSlopeYAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1682  gSigmaResSlopeYAtVtxVsP->SetName("gSigmaResSlopeYAtVtxVsP");
1683  gSigmaResSlopeYAtVtxVsP->SetTitle("#sigma_{slope_{Y}} at vertex versus p;p (GeV/c);#sigma_{slope_{Y}}");
1684  fSlopeAtVtxList->AddAt(gSigmaResSlopeYAtVtxVsP, kSigmaResSlopeYAtVtxVsP);
1685 
1686  // fit histo and fill graphs
1687  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsP)), 0., 2.e-3,
1688  "slopeX residuals at vertex", gMeanResSlopeXAtVtxVsP, gSigmaResSlopeXAtVtxVsP);
1689  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsP)), 0., 2.e-3,
1690  "slopeY residuals at vertex", gMeanResSlopeYAtVtxVsP, gSigmaResSlopeYAtVtxVsP);
1691 
1692  // --- compute slope resolution at first cluster ---
1693  fSlopeAt1stClList = new TObjArray(100);
1694  fSlopeAt1stClList->SetOwner();
1695 
1696  // define graphs
1697  TGraphAsymmErrors* gMeanResSlopeXAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1698  gMeanResSlopeXAt1stClVsP->SetName("gMeanResSlopeXAt1stClVsP");
1699  gMeanResSlopeXAt1stClVsP->SetTitle("<#Delta_{slope_{X}}> at first cluster versus p;p (GeV/c);<#Delta_{slope_{X}}>");
1700  fSlopeAt1stClList->AddAt(gMeanResSlopeXAt1stClVsP, kMeanResSlopeXAt1stClVsP);
1701  TGraphAsymmErrors* gMeanResSlopeYAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1702  gMeanResSlopeYAt1stClVsP->SetName("gMeanResSlopeYAt1stClVsP");
1703  gMeanResSlopeYAt1stClVsP->SetTitle("<#Delta_{slope_{Y}}> at first cluster versus p;p (GeV/c);<#Delta_{slope_{Y}}>");
1704  fSlopeAt1stClList->AddAt(gMeanResSlopeYAt1stClVsP, kMeanResSlopeYAt1stClVsP);
1705  TGraphAsymmErrors* gSigmaResSlopeXAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1706  gSigmaResSlopeXAt1stClVsP->SetName("gSigmaResSlopeXAt1stClVsP");
1707  gSigmaResSlopeXAt1stClVsP->SetTitle("#sigma_{slope_{X}} at first cluster versus p;p (GeV/c);#sigma_{slope_{X}}");
1708  fSlopeAt1stClList->AddAt(gSigmaResSlopeXAt1stClVsP, kSigmaResSlopeXAt1stClVsP);
1709  TGraphAsymmErrors* gSigmaResSlopeYAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1710  gSigmaResSlopeYAt1stClVsP->SetName("gSigmaResSlopeYAt1stClVsP");
1711  gSigmaResSlopeYAt1stClVsP->SetTitle("#sigma_{slope_{Y}} at first cluster versus p;p (GeV/c);#sigma_{slope_{Y}}");
1712  fSlopeAt1stClList->AddAt(gSigmaResSlopeYAt1stClVsP, kSigmaResSlopeYAt1stClVsP);
1713 
1714  // fit histo and fill graphs
1715  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAt1stClVsP)), 0., 3.e-4,
1716  "slopeX residuals at first cluster", gMeanResSlopeXAt1stClVsP, gSigmaResSlopeXAt1stClVsP);
1717  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAt1stClVsP)), 0., 3.e-4,
1718  "slopeY residuals at first cluster", gMeanResSlopeYAt1stClVsP, gSigmaResSlopeYAt1stClVsP);
1719 
1720  // --- compute eta resolution at vertex ---
1721  fEtaAtVtxList = new TObjArray(100);
1722  fEtaAtVtxList->SetOwner();
1723 
1724  // define graphs
1725  TGraphAsymmErrors* gMeanResEtaAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1726  gMeanResEtaAtVtxVsP->SetName("gMeanResEtaAtVtxVsP");
1727  gMeanResEtaAtVtxVsP->SetTitle("<#Delta_{eta}> at vertex versus p;p (GeV/c);<#Delta_{eta}>");
1728  fEtaAtVtxList->AddAt(gMeanResEtaAtVtxVsP, kMeanResEtaAtVtxVsP);
1729  TGraphAsymmErrors* gSigmaResEtaAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1730  gSigmaResEtaAtVtxVsP->SetName("gSigmaResEtaAtVtxVsP");
1731  gSigmaResEtaAtVtxVsP->SetTitle("#sigma_{eta} at vertex versus p;p (GeV/c);#sigma_{eta}");
1732  fEtaAtVtxList->AddAt(gSigmaResEtaAtVtxVsP, kSigmaResEtaAtVtxVsP);
1733 
1734  // fit histo and fill graphs
1735  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsP)), 0., 0.1,
1736  "eta residuals at vertex", gMeanResEtaAtVtxVsP, gSigmaResEtaAtVtxVsP);
1737 
1738  // --- compute phi resolution at vertex ---
1739  fPhiAtVtxList = new TObjArray(100);
1740  fPhiAtVtxList->SetOwner();
1741 
1742  // define graphs
1743  TGraphAsymmErrors* gMeanResPhiAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1744  gMeanResPhiAtVtxVsP->SetName("gMeanResPhiAtVtxVsP");
1745  gMeanResPhiAtVtxVsP->SetTitle("<#Delta_{phi}> at vertex versus p;p (GeV/c);<#Delta_{phi}>");
1746  fPhiAtVtxList->AddAt(gMeanResPhiAtVtxVsP, kMeanResPhiAtVtxVsP);
1747  TGraphAsymmErrors* gSigmaResPhiAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1748  gSigmaResPhiAtVtxVsP->SetName("gSigmaResPhiAtVtxVsP");
1749  gSigmaResPhiAtVtxVsP->SetTitle("#sigma_{phi} at vertex versus p;p (GeV/c);#sigma_{phi}");
1750  fPhiAtVtxList->AddAt(gSigmaResPhiAtVtxVsP, kSigmaResPhiAtVtxVsP);
1751 
1752  // fit histo and fill graphs
1753  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsP)), 0., 0.01,
1754  "phi residuals at vertex", gMeanResPhiAtVtxVsP, gSigmaResPhiAtVtxVsP);
1755 
1756  // --- compute DCA resolution and MCS dispersion ---
1757  fDCAList = new TObjArray(100);
1758  fDCAList->SetOwner();
1759 
1760  // define graphs
1761  TGraphAsymmErrors* gMeanPDCAVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1762  gMeanPDCAVsPIn23deg->SetName("gMeanPDCAVsPIn23deg");
1763  gMeanPDCAVsPIn23deg->SetTitle("<p #times DCA> versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);<p #times DCA> (GeV #times cm)");
1764  fDCAList->AddAt(gMeanPDCAVsPIn23deg, kMeanPDCAVsPIn23deg);
1765  TGraphAsymmErrors* gSigmaPDCAVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1766  gSigmaPDCAVsPIn23deg->SetName("gSigmaPDCAVsPIn23deg");
1767  gSigmaPDCAVsPIn23deg->SetTitle("#sigma_{p #times DCA} versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);#sigma_{p #times DCA} (GeV #times cm)");
1768  fDCAList->AddAt(gSigmaPDCAVsPIn23deg, kSigmaPDCAVsPIn23deg);
1769  TGraphAsymmErrors* gMeanPDCAVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1770  gMeanPDCAVsPIn310deg->SetName("gMeanPDCAVsPIn310deg");
1771  gMeanPDCAVsPIn310deg->SetTitle("<p #times DCA> versus p for tracks within [3,10[ degrees at absorber end;p (GeV/c);<p #times DCA> (GeV #times cm)");
1772  fDCAList->AddAt(gMeanPDCAVsPIn310deg, kMeanPDCAVsPIn310deg);
1773  TGraphAsymmErrors* gSigmaPDCAVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1774  gSigmaPDCAVsPIn310deg->SetName("gSigmaPDCAVsPIn310deg");
1775  gSigmaPDCAVsPIn310deg->SetTitle("#sigma_{p #times DCA} versus p for tracks within [3,10[ degrees at absorber end;p (GeV/c);#sigma_{p #times DCA} (GeV #times cm)");
1776  fDCAList->AddAt(gSigmaPDCAVsPIn310deg, kSigmaPDCAVsPIn310deg);
1777 
1778  TGraphAsymmErrors* gMeanPMCSAngVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1779  gMeanPMCSAngVsPIn23deg->SetName("gMeanPMCSAngVsPIn23deg");
1780  gMeanPMCSAngVsPIn23deg->SetTitle("<p #times #Delta#theta_{MCS}> versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);<p #times #Delta#theta_{MCS}> (GeV)");
1781  fDCAList->AddAt(gMeanPMCSAngVsPIn23deg, kMeanPMCSAngVsPIn23deg);
1782  TGraphAsymmErrors* gSigmaPMCSAngVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1783  gSigmaPMCSAngVsPIn23deg->SetName("gSigmaPMCSAngVsPIn23deg");
1784  gSigmaPMCSAngVsPIn23deg->SetTitle("#sigma_{p #times #Delta#theta_{MCS}} versus p for tracks within [2,3[ degrees at absorber end;p (GeV/c);#sigma_{p #times #Delta#theta_{MCS}} (GeV)");
1785  fDCAList->AddAt(gSigmaPMCSAngVsPIn23deg, kSigmaPMCSAngVsPIn23deg);
1786  TGraphAsymmErrors* gMeanPMCSAngVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1787  gMeanPMCSAngVsPIn310deg->SetName("gMeanPMCSAngVsPIn310deg");
1788  gMeanPMCSAngVsPIn310deg->SetTitle("<p #times #Delta#theta_{MCS}> versus p for tracks within [3,10[ degrees at absorber end;p (GeV/c);<p #times #Delta#theta_{MCS}> (GeV)");
1789  fDCAList->AddAt(gMeanPMCSAngVsPIn310deg, kMeanPMCSAngVsPIn310deg);
1790  TGraphAsymmErrors* gSigmaPMCSAngVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1791  gSigmaPMCSAngVsPIn310deg->SetName("gSigmaPMCSAngVsPIn310deg");
1792  gSigmaPMCSAngVsPIn310deg->SetTitle("#sigma_{p #times #Delta#theta_{MCS}} versus p for tracks within [3,10[ degrees at absorber end;p (GeV/c);#sigma_{p #times #Delta#theta_{MCS}} (GeV)");
1793  fDCAList->AddAt(gSigmaPMCSAngVsPIn310deg, kSigmaPMCSAngVsPIn310deg);
1794 
1795  // fit histo and fill graphs
1796  FitPDCAVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPIn23deg)),
1797  "p*DCA (tracks in [2,3] deg.)", gMeanPDCAVsPIn23deg, gSigmaPDCAVsPIn23deg);
1798  FitPDCAVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPIn310deg)),
1799  "p*DCA (tracks in [3,10] deg.)", gMeanPDCAVsPIn310deg, gSigmaPDCAVsPIn310deg);
1800  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPMCSAngVsPIn23deg)), 0., 2.e-3,
1801  "p*MCSAngle (tracks in [2,3] deg.)", gMeanPMCSAngVsPIn23deg, gSigmaPMCSAngVsPIn23deg);
1802  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPMCSAngVsPIn310deg)), 0., 2.e-3,
1803  "p*MCSAngle (tracks in [3,10] deg.)", gMeanPMCSAngVsPIn310deg, gSigmaPMCSAngVsPIn310deg);
1804 
1805  // --- compute cluster resolution ---
1806  fClusterList = new TObjArray(100);
1807  fClusterList->SetOwner();
1808 
1809  // define graphs per chamber
1810  TGraphErrors* gMeanResClXVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1811  gMeanResClXVsCh->SetName("gMeanResClXVsCh");
1812  gMeanResClXVsCh->SetTitle("cluster-trackRef residual-X per Ch: mean;chamber ID;<#Delta_{X}> (cm)");
1813  gMeanResClXVsCh->SetMarkerStyle(kFullDotLarge);
1814  fClusterList->AddAt(gMeanResClXVsCh, kMeanResClXVsCh);
1815  TGraphErrors* gMeanResClYVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1816  gMeanResClYVsCh->SetName("gMeanResClYVsCh");
1817  gMeanResClYVsCh->SetTitle("cluster-trackRef residual-Y per Ch: mean;chamber ID;<#Delta_{Y}> (cm)");
1818  gMeanResClYVsCh->SetMarkerStyle(kFullDotLarge);
1819  fClusterList->AddAt(gMeanResClYVsCh, kMeanResClYVsCh);
1820  TGraphErrors* gSigmaResClXVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1821  gSigmaResClXVsCh->SetName("gSigmaResClXVsCh");
1822  gSigmaResClXVsCh->SetTitle("cluster-trackRef residual-X per Ch: sigma;chamber ID;#sigma_{X} (cm)");
1823  gSigmaResClXVsCh->SetMarkerStyle(kFullDotLarge);
1824  fClusterList->AddAt(gSigmaResClXVsCh, kSigmaResClXVsCh);
1825  TGraphErrors* gSigmaResClYVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1826  gSigmaResClYVsCh->SetName("gSigmaResClYVsCh");
1827  gSigmaResClYVsCh->SetTitle("cluster-trackRef residual-Y per Ch: sigma;chamber ID;#sigma_{Y} (cm)");
1828  gSigmaResClYVsCh->SetMarkerStyle(kFullDotLarge);
1829  fClusterList->AddAt(gSigmaResClYVsCh, kSigmaResClYVsCh);
1830 
1831  // define graphs per DE
1832  TGraphErrors* gMeanResClXVsDE = new TGraphErrors(fNDE);
1833  gMeanResClXVsDE->SetName("gMeanResClXVsDE");
1834  gMeanResClXVsDE->SetTitle("cluster-trackRef residual-X per DE: mean;DE ID;<#Delta_{X}> (cm)");
1835  gMeanResClXVsDE->SetMarkerStyle(kFullDotLarge);
1836  fClusterList->AddAt(gMeanResClXVsDE, kMeanResClXVsDE);
1837  TGraphErrors* gMeanResClYVsDE = new TGraphErrors(fNDE);
1838  gMeanResClYVsDE->SetName("gMeanResClYVsDE");
1839  gMeanResClYVsDE->SetTitle("cluster-trackRef residual-Y per dE: mean;DE ID;<#Delta_{Y}> (cm)");
1840  gMeanResClYVsDE->SetMarkerStyle(kFullDotLarge);
1841  fClusterList->AddAt(gMeanResClYVsDE, kMeanResClYVsDE);
1842  TGraphErrors* gSigmaResClXVsDE = new TGraphErrors(fNDE);
1843  gSigmaResClXVsDE->SetName("gSigmaResClXVsDE");
1844  gSigmaResClXVsDE->SetTitle("cluster-trackRef residual-X per DE: sigma;DE ID;#sigma_{X} (cm)");
1845  gSigmaResClXVsDE->SetMarkerStyle(kFullDotLarge);
1846  fClusterList->AddAt(gSigmaResClXVsDE, kSigmaResClXVsDE);
1847  TGraphErrors* gSigmaResClYVsDE = new TGraphErrors(fNDE);
1848  gSigmaResClYVsDE->SetName("gSigmaResClYVsDE");
1849  gSigmaResClYVsDE->SetTitle("cluster-trackRef residual-Y per DE: sigma;DE ID;#sigma_{Y} (cm)");
1850  gSigmaResClYVsDE->SetMarkerStyle(kFullDotLarge);
1851  fClusterList->AddAt(gSigmaResClYVsDE, kSigmaResClYVsDE);
1852 
1853  // fit histo and fill graphs per chamber
1854  Double_t clusterResPerCh[10][2];
1855  for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
1856  TH1D *tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsCh))->ProjectionY("tmp",i+1,i+1,"e");
1857  FitClusterResidual(tmp, i, clusterResPerCh[i][0], gMeanResClXVsCh, gSigmaResClXVsCh);
1858  delete tmp;
1859  tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClYVsCh))->ProjectionY("tmp",i+1,i+1,"e");
1860  FitClusterResidual(tmp, i, clusterResPerCh[i][1], gMeanResClYVsCh, gSigmaResClYVsCh);
1861  delete tmp;
1862  }
1863 
1864  // fit histo and fill graphs per DE
1865  Double_t clusterResPerDE[200][2];
1866  for (Int_t i = 0; i < fNDE; i++) {
1867  TH1D *tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsDE))->ProjectionY("tmp",i+1,i+1,"e");
1868  FitClusterResidual(tmp, i, clusterResPerDE[i][0], gMeanResClXVsDE, gSigmaResClXVsDE);
1869  delete tmp;
1870  tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClYVsDE))->ProjectionY("tmp",i+1,i+1,"e");
1871  FitClusterResidual(tmp, i, clusterResPerDE[i][1], gMeanResClYVsDE, gSigmaResClYVsDE);
1872  delete tmp;
1873  }
1874 
1875  // set DE graph labels
1876  TAxis* xAxis = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsDE))->GetXaxis();
1877  gMeanResClXVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1878  gMeanResClYVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1879  gSigmaResClXVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1880  gSigmaResClYVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1881  for (Int_t i = 1; i <= fNDE; i++) {
1882  const char* label = xAxis->GetBinLabel(i);
1883  gMeanResClXVsDE->GetXaxis()->SetBinLabel(i, label);
1884  gMeanResClYVsDE->GetXaxis()->SetBinLabel(i, label);
1885  gSigmaResClXVsDE->GetXaxis()->SetBinLabel(i, label);
1886  gSigmaResClYVsDE->GetXaxis()->SetBinLabel(i, label);
1887  }
1888 
1889  // --- diplay momentum residuals at vertex ---
1890  TCanvas* cResPAtVtx = DrawVsAng("cResPAtVtx", "momentum residual at vertex in 3 angular regions",
1891  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPAtVtx)),
1892  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsAngleAtAbsEnd)));
1893  fPAtVtxList->AddAt(cResPAtVtx, kcResPAtVtx);
1894  TCanvas* cResPAtVtxMC = DrawVsAng("cResPAtVtxMC", "momentum residual at vertex in 3 MC angular regions",
1895  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPAtVtx)),
1896  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsMCAngle)));
1897  fPAtVtxList->AddAt(cResPAtVtxMC, kcResPAtVtxMC);
1898  TCanvas* cResPAtVtxVsPosAbsEndMC = DrawVsPos("cResPAtVtxVsPosAbsEndMC", "momentum residual at vertex versus position at absorber end in 3 MC angular regions",
1899  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn02degMC)),
1900  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn23degMC)),
1901  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn310degMC)));
1902  fPAtVtxList->AddAt(cResPAtVtxVsPosAbsEndMC, kcResPAtVtxVsPosAbsEndMC);
1903  TCanvas* cResPAtVtxVsPIn23deg = DrawFitLandauGausResPVsP("cResPAtVtxVsPIn23deg", "momentum residual for tracks between 2 and 3 degrees",
1904  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn23deg)),
1905  10, "momentum residuals at vertex (tracks in [2,3] deg.)");
1906  fPAtVtxList->AddAt(cResPAtVtxVsPIn23deg, kcResPAtVtxVsPIn23deg);
1907  TCanvas* cResPAtVtxVsPIn310deg = DrawFitLandauGausResPVsP("cResPAtVtxVsPIn310deg", "momentum residual for tracks between 3 and 10 degrees",
1908  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn310deg)),
1909  10, "momentum residuals at vertex (tracks in [3,10] deg.)");
1910  fPAtVtxList->AddAt(cResPAtVtxVsPIn310deg, kcResPAtVtxVsPIn310deg);
1911  TCanvas* cResPAtVtxVsPIn02degMC = DrawResPVsP("cResPAtVtxVsPIn02degMC", "momentum residuals for tracks with MC angle < 2 degrees",
1912  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn02degMC)), 5);
1913  fPAtVtxList->AddAt(cResPAtVtxVsPIn02degMC, kcResPAtVtxVsPIn02degMC);
1914 
1915  // --- diplay slopeX residuals at vertex ---
1916  TCanvas* cResSlopeXAtVtx = DrawVsAng("cResSlopeXAtVtx", "slope_{X} residual at vertex in 3 angular regions",
1917  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeXAtVtx)),
1918  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsAngleAtAbsEnd)));
1919  fSlopeAtVtxList->AddAt(cResSlopeXAtVtx, kcResSlopeXAtVtx);
1920  TCanvas* cResSlopeYAtVtx = DrawVsAng("cResSlopeYAtVtx", "slope_{Y} residual at vertex in 3 angular regions",
1921  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeYAtVtx)),
1922  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsAngleAtAbsEnd)));
1923  fSlopeAtVtxList->AddAt(cResSlopeYAtVtx, kcResSlopeYAtVtx);
1924  TCanvas* cResSlopeXAtVtxMC = DrawVsAng("cResSlopeXAtVtxMC", "slope_{X} residual at vertex in 3 MC angular regions",
1925  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeXAtVtx)),
1926  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsMCAngle)));
1927  fSlopeAtVtxList->AddAt(cResSlopeXAtVtxMC, kcResSlopeXAtVtxMC);
1928  TCanvas* cResSlopeYAtVtxMC = DrawVsAng("cResSlopeYAtVtxMC", "slope_{Y} residual at vertex in 3 MC angular regions",
1929  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeYAtVtx)),
1930  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsMCAngle)));
1931  fSlopeAtVtxList->AddAt(cResSlopeYAtVtxMC, kcResSlopeYAtVtxMC);
1932  TCanvas* cResSlopeXAtVtxVsPosAbsEndMC = DrawVsPos("cResSlopeXAtVtxVsPosAbsEndMC", "slope_{X} residual at vertex versus position at absorber end in 3 MC angular regions",
1933  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn02degMC)),
1934  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn23degMC)),
1935  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn310degMC)));
1936  fSlopeAtVtxList->AddAt(cResSlopeXAtVtxVsPosAbsEndMC, kcResSlopeXAtVtxVsPosAbsEndMC);
1937  TCanvas* cResSlopeYAtVtxVsPosAbsEndMC = DrawVsPos("cResSlopeYAtVtxVsPosAbsEndMC", "slope_{Y} residual at vertex versus position at absorber end in 3 MC angular regions",
1938  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn02degMC)),
1939  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn23degMC)),
1940  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn310degMC)));
1941  fSlopeAtVtxList->AddAt(cResSlopeYAtVtxVsPosAbsEndMC, kcResSlopeYAtVtxVsPosAbsEndMC);
1942 
1943  // --- diplay eta residuals at vertex ---
1944  TCanvas* cResEtaAtVtx = DrawVsAng("cResEtaAtVtx", "eta residual at vertex in 3 angular regions",
1945  static_cast<TH1*>(fTrackerList->UncheckedAt(kResEtaAtVtx)),
1946  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsAngleAtAbsEnd)));
1947  fEtaAtVtxList->AddAt(cResEtaAtVtx, kcResEtaAtVtx);
1948  TCanvas* cResEtaAtVtxMC = DrawVsAng("cResEtaAtVtxMC", "eta residual at vertex in 3 MC angular regions",
1949  static_cast<TH1*>(fTrackerList->UncheckedAt(kResEtaAtVtx)),
1950  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsMCAngle)));
1951  fEtaAtVtxList->AddAt(cResEtaAtVtxMC, kcResEtaAtVtxMC);
1952  TCanvas* cResEtaAtVtxVsPosAbsEndMC = DrawVsPos("cResEtaAtVtxVsPosAbsEndMC", "eta residual at vertex versus position at absorber end in 3 MC angular regions",
1953  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn02degMC)),
1954  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn23degMC)),
1955  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn310degMC)));
1956  fEtaAtVtxList->AddAt(cResEtaAtVtxVsPosAbsEndMC, kcResEtaAtVtxVsPosAbsEndMC);
1957 
1958  // --- diplay phi residuals at vertex ---
1959  TCanvas* cResPhiAtVtx = DrawVsAng("cResPhiAtVtx", "phi residual at vertex in 3 angular regions",
1960  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPhiAtVtx)),
1961  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsAngleAtAbsEnd)));
1962  fPhiAtVtxList->AddAt(cResPhiAtVtx, kcResPhiAtVtx);
1963  TCanvas* cResPhiAtVtxMC = DrawVsAng("cResPhiAtVtxMC", "phi residual at vertex in 3 MC angular regions",
1964  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPhiAtVtx)),
1965  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsMCAngle)));
1966  fPhiAtVtxList->AddAt(cResPhiAtVtxMC, kcResPhiAtVtxMC);
1967  TCanvas* cResPhiAtVtxVsPosAbsEndMC = DrawVsPos("cResPhiAtVtxVsPosAbsEndMC", "phi residual at vertex versus position at absorber end in 3 MC angular regions",
1968  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn02degMC)),
1969  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn23degMC)),
1970  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn310degMC)));
1971  fPhiAtVtxList->AddAt(cResPhiAtVtxVsPosAbsEndMC, kcResPhiAtVtxVsPosAbsEndMC);
1972 
1973  // --- diplay P*DCA ---
1974  TCanvas* cPDCA = DrawVsAng("cPDCA", "p #times DCA in 3 angular regions",
1975  static_cast<TH1*>(fTrackerList->UncheckedAt(kPDCA)),
1976  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsAngleAtAbsEnd)));
1977  fDCAList->AddAt(cPDCA, kcPDCA);
1978  TCanvas* cPDCAMC = DrawVsAng("cPDCAMC", "p #times DCA in 3 MC angular regions",
1979  static_cast<TH1*>(fTrackerList->UncheckedAt(kPDCA)),
1980  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsMCAngle)));
1981  fDCAList->AddAt(cPDCAMC, kcPDCAMC);
1982  TCanvas* cPDCAVsPosAbsEndMC = DrawVsPos("cPDCAVsPosAbsEndMC", "p #times DCA versus position at absorber end in 3 MC angular regions",
1983  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn02degMC)),
1984  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn23degMC)),
1985  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn310degMC)));
1986  fDCAList->AddAt(cPDCAVsPosAbsEndMC, kcPDCAVsPosAbsEndMC);
1987 
1988  // post param containers
1989  PostData(4, fEfficiencyList);
1990  PostData(5, fPAtVtxList);
1991  PostData(6, fSlopeAtVtxList);
1992  PostData(7, fEtaAtVtxList);
1993  PostData(8, fPhiAtVtxList);
1994  PostData(9, fPAt1stClList);
1995  PostData(10, fSlopeAt1stClList);
1996  PostData(11, fDCAList);
1997  PostData(12, fClusterList);
1998 }
1999 
2000 
2001 //________________________________________________________________________
2002 Bool_t AliAnalysisTaskMuonPerformance::GetEfficiency(AliCFEffGrid* efficiency, Double_t& calcEff, Double_t& calcEffErr)
2003 {
2004  //
2006  //
2007 
2008  Bool_t isGood = kTRUE;
2009 
2010  TH1* histo = 0x0;
2011  Double_t sum[2] = {0., 0.};
2012  for ( Int_t ihisto=0; ihisto<2; ihisto++ ) {
2013  histo = ( ihisto==0 ) ? efficiency->GetNum()->Project(kVarCharge) : efficiency->GetDen()->Project(kVarCharge);
2014  sum[ihisto] = histo->Integral();
2015  delete histo;
2016  }
2017 
2018  if ( sum[1] == 0. ) isGood = kFALSE;
2019 
2020  calcEff = ( isGood ) ? sum[0]/sum[1] : 0.;
2021  if ( calcEff > 1. ) isGood = kFALSE;
2022 
2023  calcEffErr = ( isGood ) ? TMath::Sqrt(calcEff*(1-calcEff)/sum[1]) : 0.;
2024 
2025  return isGood;
2026 }
2027 
2028 //________________________________________________________________________
2029 Int_t AliAnalysisTaskMuonPerformance::RecoTrackMother(AliMCParticle* mcParticle)
2030 {
2031  //
2033  //
2034 
2035  if ( ! mcParticle )
2036  return kUnknownPart;
2037 
2038  Int_t recoPdg = mcParticle->PdgCode();
2039 
2040  // Track is not a muon
2041  if ( TMath::Abs(recoPdg) != 13 ) return kRecoHadron;
2042 
2043  Int_t imother = mcParticle->GetMother();
2044 
2045  Bool_t isFirstMotherHF = kFALSE;
2046  Int_t step = 0;
2047 
2048  //printf("\n"); // REMEMBER TO CUT
2049  while ( imother >= 0 ) {
2050  TParticle* part = static_cast<AliMCParticle*>(fMCEvent->GetTrack(imother))->Particle();
2051  //printf("Particle %s -> %s\n", part->GetName(), TMCProcessName[part->GetUniqueID()]); // REMEMBER TO CUT
2052 
2053  if ( part->GetUniqueID() == kPHadronic )
2054  return kSecondaryMu;
2055 
2056  Int_t absPdg = TMath::Abs(part->GetPdgCode());
2057 
2058  step++;
2059  if ( step == 1 ) // only for first mother
2060  isFirstMotherHF = ( ( absPdg >= 400 && absPdg < 600 ) ||
2061  ( absPdg >= 4000 && absPdg < 6000 ) );
2062 
2063  if ( absPdg < 4 )
2064  return kPrimaryMu;
2065 
2066  if ( isFirstMotherHF) {
2067  if ( absPdg == 4 )
2068  return kCharmMu;
2069  else if ( absPdg == 5 )
2070  return kBeautyMu;
2071  }
2072 
2073  imother = part->GetFirstMother();
2074  }
2075 
2076  return kPrimaryMu;
2077 
2078 }
2079 
2080 //________________________________________________________________________
2081 Float_t AliAnalysisTaskMuonPerformance::GetBinThetaAbsEnd(Float_t RAtAbsEnd, Bool_t isTheta)
2082 {
2083  //
2085  //
2086  Float_t thetaDeg = ( isTheta ) ? RAtAbsEnd : TMath::ATan( RAtAbsEnd / 505. );
2087  thetaDeg *= TMath::RadToDeg();
2088  if ( thetaDeg < 2. )
2089  return 0.;
2090  else if ( thetaDeg < 3. )
2091  return 1.;
2092  else if ( thetaDeg < 10. )
2093  return 2.;
2094 
2095  return 3.;
2096 }
2097 
2098 //________________________________________________________________________
2099 void AliAnalysisTaskMuonPerformance::FillContainerInfoReco(Double_t* containerInput, AliESDMuonTrack* esdTrack,
2100  Bool_t isValid, Int_t mcID)
2101 {
2102  //
2104  //
2105 
2106  AliMCParticle* mcPart = (mcID >= 0) ? static_cast<AliMCParticle*>(fMCEvent->GetTrack(mcID)) : 0x0;
2107 
2108  if (fUseMCKinematics && mcPart) {
2109  containerInput[kVarPt] = mcPart->Pt();
2110  containerInput[kVarEta] = mcPart->Eta();
2111  containerInput[kVarPhi] = mcPart->Phi();
2112  } else {
2113  containerInput[kVarPt] = esdTrack->Pt();
2114  containerInput[kVarEta] = esdTrack->Eta();
2115  containerInput[kVarPhi] = esdTrack->Phi();
2116  }
2117  containerInput[kVarThetaZones] = GetBinThetaAbsEnd(esdTrack->GetRAtAbsorberEnd());
2118  containerInput[kVarCharge] = static_cast<Double_t>(esdTrack->Charge());
2119  containerInput[kVarHasTracker] = static_cast<Double_t>(esdTrack->ContainTrackerData() && isValid);
2120  if (esdTrack->GetMatchTrigger() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kNoMatchTrig);
2121  else if (esdTrack->GetMatchTrigger() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
2122  else if (esdTrack->GetMatchTrigger() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
2123  else if (esdTrack->GetMatchTrigger() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
2124  containerInput[kVarMotherType] = static_cast<Double_t>(RecoTrackMother(mcPart));
2125 
2126 }
2127 
2128 //________________________________________________________________________
2129 void AliAnalysisTaskMuonPerformance::FillContainerInfoMC(Double_t* containerInput, AliMCParticle* mcPart)
2130 {
2131  //
2133  //
2134 
2135  containerInput[kVarPt] = mcPart->Pt();
2136  containerInput[kVarEta] = mcPart->Eta();
2137  containerInput[kVarPhi] = mcPart->Phi();
2138  containerInput[kVarThetaZones] = GetBinThetaAbsEnd(TMath::Pi()-mcPart->Theta(),kTRUE);
2139  containerInput[kVarCharge] = static_cast<Double_t>(mcPart->Charge())/3.;
2140  containerInput[kVarMotherType] = static_cast<Double_t>(RecoTrackMother(mcPart));
2141 
2142 }
2143 
2144 //________________________________________________________________________
2145 Double_t langaufun(Double_t *x, Double_t *par)
2146 {
2147  //Fit parameters:
2148  //par[0]=Width (scale) parameter of Landau density
2149  //par[1]=Most Probable (MP, location) parameter of Landau density
2150  //par[2]=Total area (integral -inf to inf, normalization constant)
2151  //par[3]=Width (sigma) of convoluted Gaussian function
2152  //
2153  //In the Landau distribution (represented by the CERNLIB approximation),
2154  //the maximum is located at x=-0.22278298 with the location parameter=0.
2155  //This shift is corrected within this function, so that the actual
2156  //maximum is identical to the MP parameter.
2157 
2158  // Numeric constants
2159  Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
2160  Double_t mpshift = -0.22278298; // Landau maximum location
2161 
2162  // MP shift correction
2163  Double_t mpc = par[1] - mpshift * par[0];
2164 
2165  // Range of convolution integral
2166  Double_t xlow = x[0] - 5. * par[3];
2167  Double_t xupp = x[0] + 5. * par[3];
2168  Double_t step = TMath::Min(0.2*par[0],0.1*par[3]);
2169 
2170  // Convolution integral of Landau and Gaussian by sum
2171  Double_t xx = xlow + 0.5 * step;
2172  Double_t sum = 0.0;
2173  while (xx < xupp) {
2174  sum += TMath::Landau(-xx,mpc,par[0]) * TMath::Gaus(x[0],xx,par[3]);
2175  xx += step;
2176  }
2177 
2178  return (par[2] * step * sum * invsq2pi / par[3] / par[0]);
2179 
2180 }
2181 
2182 //________________________________________________________________________
2183 void AliAnalysisTaskMuonPerformance::FitLandauGausResVsP(TH2* h, const char* fitting, TGraphAsymmErrors* gMean,
2184  TGraphAsymmErrors* gMostProb, TGraphAsymmErrors* gSigma)
2185 {
2187 
2188  static TF1 *fGaus2 = 0x0;
2189  if (!fGaus2) fGaus2 = new TF1("fGaus2","gaus");
2190  static TF1* fLandauGaus = 0x0;
2191  if (!fLandauGaus) {
2192  fLandauGaus = new TF1("fLandauGaus",langaufun,h->GetYaxis()->GetBinLowEdge(1),h->GetYaxis()->GetBinLowEdge(h->GetNbinsY()+1),4);
2193  fLandauGaus->SetNpx(10000);
2194  }
2195 
2196  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2197  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2198 
2199  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2200 
2201  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2202 
2203  // first fit
2204  fGaus2->SetParameters((Double_t)tmp->GetEntries(), 0., 1.);
2205  tmp->Fit("fGaus2", "WWNQ");
2206 
2207  // rebin histo
2208  Double_t sigma = fGaus2->GetParameter(2);
2209  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*tmp->GetNbinsX(),TMath::Max(0.5*sigma/tmp->GetBinWidth(1),1.)));
2210  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2211  tmp->Rebin(rebin);
2212 
2213  // second fit
2214  Double_t mean = fGaus2->GetParameter(1);
2215  fLandauGaus->SetParameters(0.25*sigma*TMath::Sqrt(8.*log(2.)), mean, tmp->GetEntries()*tmp->GetXaxis()->GetBinWidth(1), 0.5*sigma);
2216  fLandauGaus->SetParLimits(0, 0.0025*sigma*TMath::Sqrt(8.*log(2.)), 1000.);
2217  fLandauGaus->SetParLimits(3, 0., 2.*sigma);
2218  Double_t xMin = TMath::Max(mean-50.*sigma, tmp->GetXaxis()->GetXmin());
2219  Double_t xMax = TMath::Min(mean+10.*sigma, tmp->GetXaxis()->GetXmax());
2220  if (xMin < tmp->GetXaxis()->GetXmax() && xMax > tmp->GetXaxis()->GetXmin()) fLandauGaus->SetRange(xMin, xMax);
2221  tmp->Fit("fLandauGaus","RNQ");
2222 
2223  // get fit results and fill graph
2224  Double_t fwhm = 4.*fLandauGaus->GetParameter(0);
2225  sigma = fLandauGaus->GetParameter(3);
2226  Double_t sigmaP = TMath::Sqrt(sigma*sigma + fwhm*fwhm/(8.*log(2.)));
2227  Double_t fwhmErr = fLandauGaus->GetParError(0);
2228  Double_t sigmaErr = fLandauGaus->GetParError(3);
2229  Double_t sigmaPErr = TMath::Sqrt(sigma*sigma*sigmaErr*sigmaErr + fwhm*fwhm*fwhmErr*fwhmErr/(64.*log(2.)*log(2.))) / sigmaP;
2230  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2231  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2232  h->GetXaxis()->SetRange();
2233  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2234  gMean->SetPoint(i/rebinFactorX-1, p, tmp->GetMean());
2235  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], tmp->GetMeanError(), tmp->GetMeanError());
2236  gMostProb->SetPoint(i/rebinFactorX-1, p, -fLandauGaus->GetParameter(1));
2237  gMostProb->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fLandauGaus->GetParError(1), fLandauGaus->GetParError(1));
2238  gSigma->SetPoint(i/rebinFactorX-1, p, sigmaP);
2239  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], sigmaPErr, sigmaPErr);
2240 
2241  // clean memory
2242  delete tmp;
2243  }
2244 
2245  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2246 
2247 }
2248 
2249 //________________________________________________________________________
2250 void AliAnalysisTaskMuonPerformance::FitGausResVsMom(TH2* h, const Double_t mean0,
2251  const Double_t sigma0, const char* fitting,
2252  TGraphAsymmErrors* gMean, TGraphAsymmErrors* gSigma)
2253 {
2255 
2256  static TF1* fGaus = 0x0;
2257  if (!fGaus) fGaus = new TF1("fGaus","gaus");
2258 
2259  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2260  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2261 
2262  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2263 
2264  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2265 
2266  // first fit
2267  fGaus->SetParameters(tmp->GetEntries(), mean0, sigma0);
2268  tmp->Fit("fGaus","WWNQ");
2269 
2270  // rebin histo
2271  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*tmp->GetNbinsX(),TMath::Max(0.5*fGaus->GetParameter(2)/tmp->GetBinWidth(1),1.)));
2272  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2273  tmp->Rebin(rebin);
2274 
2275  // second fit
2276  tmp->Fit("fGaus","NQ");
2277 
2278  // get fit results and fill histograms
2279  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2280  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2281  h->GetXaxis()->SetRange();
2282  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2283  gMean->SetPoint(i/rebinFactorX-1, p, fGaus->GetParameter(1));
2284  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fGaus->GetParError(1), fGaus->GetParError(1));
2285  gSigma->SetPoint(i/rebinFactorX-1, p, fGaus->GetParameter(2));
2286  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fGaus->GetParError(2), fGaus->GetParError(2));
2287 
2288  // clean memory
2289  delete tmp;
2290  }
2291 
2292  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2293 
2294 }
2295 
2296 //________________________________________________________________________
2297 void AliAnalysisTaskMuonPerformance::FitPDCAVsMom(TH2* h, const char* fitting,
2298  TGraphAsymmErrors* gMean, TGraphAsymmErrors* gSigma)
2299 {
2301 
2302  static TF1* fPGaus = 0x0;
2303  if (!fPGaus) fPGaus = new TF1("fPGaus","x*gaus");
2304 
2305  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2306  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2307 
2308  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2309 
2310  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2311 
2312  // rebin histo
2313  Int_t rebin = static_cast<Int_t>(25*(tmp->GetNbinsX()/(tmp->GetBinLowEdge(tmp->GetNbinsX()+1)-tmp->GetBinLowEdge(1))));
2314  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2315  tmp->Rebin(rebin);
2316 
2317  // fit
2318  fPGaus->SetParameters(1.,0.,80.);
2319  tmp->Fit("fPGaus","NQ");
2320 
2321  // get fit results and fill histograms
2322  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2323  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2324  h->GetXaxis()->SetRange();
2325  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2326  gMean->SetPoint(i/rebinFactorX-1, p, fPGaus->GetParameter(1));
2327  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fPGaus->GetParError(1), fPGaus->GetParError(1));
2328  gSigma->SetPoint(i/rebinFactorX-1, p, fPGaus->GetParameter(2));
2329  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fPGaus->GetParError(2), fPGaus->GetParError(2));
2330 
2331  // clean memory
2332  delete tmp;
2333  }
2334 
2335  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2336 
2337 }
2338 
2339 //________________________________________________________________________
2340 void AliAnalysisTaskMuonPerformance::FitClusterResidual(TH1* h, Int_t i, Double_t& sigma,
2341  TGraphErrors* gMean, TGraphErrors* gSigma)
2342 {
2344 
2345  static TF1* fRGaus = 0x0;
2346  Double_t mean, meanErr, sigmaErr;
2347 
2348  if (fFitResiduals) {
2349 
2350  if (!fRGaus) fRGaus = new TF1("fRGaus","gaus");
2351 
2352  // first fit
2353  Double_t xMin = h->GetXaxis()->GetXmin();
2354  Double_t xMax = h->GetXaxis()->GetXmax();
2355  fRGaus->SetRange(xMin, xMax);
2356  fRGaus->SetParameters(h->GetEntries(), 0., 0.1);
2357  fRGaus->SetParLimits(1, xMin, xMax);
2358  h->Fit("fRGaus", "WWNQ");
2359 
2360  // rebin histo
2361  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*h->GetNbinsX(),TMath::Max(0.3*fRGaus->GetParameter(2)/h->GetBinWidth(1),1.)));
2362  while (h->GetNbinsX()%rebin!=0) rebin--;
2363  h->Rebin(rebin);
2364 
2365  // second fit
2366  xMin = TMath::Max(fRGaus->GetParameter(1)-10.*fRGaus->GetParameter(2), h->GetXaxis()->GetXmin());
2367  xMax = TMath::Min(fRGaus->GetParameter(1)+10.*fRGaus->GetParameter(2), h->GetXaxis()->GetXmax());
2368  fRGaus->SetRange(xMin, xMax);
2369  fRGaus->SetParLimits(1, xMin, xMax);
2370  h->Fit("fRGaus","NQR");
2371 
2372  mean = fRGaus->GetParameter(1);
2373  meanErr = fRGaus->GetParError(1);
2374  sigma = fRGaus->GetParameter(2);
2375  sigmaErr = fRGaus->GetParError(2);
2376 
2377  } else {
2378 
2379  Zoom(h);
2380  mean = h->GetMean();
2381  meanErr = h->GetMeanError();
2382  sigma = h->GetRMS();
2383  sigmaErr = h->GetRMSError();
2384  h->GetXaxis()->SetRange(0,0);
2385 
2386  }
2387 
2388  gMean->SetPoint(i, i+1, mean);
2389  gMean->SetPointError(i, 0., meanErr);
2390 
2391  if (fCorrectForSystematics) {
2392  Double_t s = TMath::Sqrt(mean*mean + sigma*sigma);
2393  sigmaErr = (s>0.) ? TMath::Sqrt(sigma*sigma*sigmaErr*sigmaErr + mean*mean*meanErr*meanErr) / s : 0.;
2394  sigma = s;
2395  }
2396 
2397  gSigma->SetPoint(i, i+1, sigma);
2398  gSigma->SetPointError(i, 0., sigmaErr);
2399 
2400 }
2401 
2402 //________________________________________________________________________
2403 TCanvas* AliAnalysisTaskMuonPerformance::DrawVsAng(const char* name, const char* title, TH1* h1, TH2* h2)
2404 {
2406  TCanvas* c = new TCanvas(name, title);
2407  c->cd();
2408  h1->DrawCopy();
2409  TH1D *proj1 = h2->ProjectionY(Form("%s_proj_0_2",h2->GetName()),1,2);
2410  proj1->SetLineColor(2);
2411  proj1->Draw("sames");
2412  TH1D *proj2 = h2->ProjectionY(Form("%s_proj_2_3",h2->GetName()),3,3);
2413  proj2->SetLineColor(4);
2414  proj2->Draw("sames");
2415  TH1D *proj3 = h2->ProjectionY(Form("%s__proj_3_10",h2->GetName()),4,10);
2416  proj3->SetLineColor(3);
2417  proj3->Draw("sames");
2418  return c;
2419 }
2420 
2421 //________________________________________________________________________
2422 TCanvas* AliAnalysisTaskMuonPerformance::DrawVsPos(const char* name, const char* title, TH2* h1, TH2* h2, TH2* h3)
2423 {
2425  TCanvas* c = new TCanvas(name, title);
2426  c->cd();
2427  h1->SetMarkerColor(2);
2428  h1->DrawCopy();
2429  h2->SetMarkerColor(4);
2430  h2->DrawCopy("sames");
2431  h3->SetMarkerColor(3);
2432  h3->DrawCopy("sames");
2433  return c;
2434 }
2435 
2436 //________________________________________________________________________
2437 TCanvas* AliAnalysisTaskMuonPerformance::DrawFitLandauGausResPVsP(const char* name, const char* title,
2438  TH2* h, const Int_t nBins, const char* fitting)
2439 {
2441 
2442  static TF1 *fGaus3 = 0x0;
2443  if (!fGaus3) fGaus3 = new TF1("fGaus3","gaus");
2444  static TF1* fLandauGaus2 = 0x0;
2445  if (!fLandauGaus2) {
2446  fLandauGaus2 = new TF1("fLandauGaus2",langaufun,h->GetYaxis()->GetBinLowEdge(1),h->GetYaxis()->GetBinLowEdge(h->GetNbinsY()+1),4);
2447  fLandauGaus2->SetNpx(10000);
2448  }
2449 
2450  TLegend* l = new TLegend(0.15,0.25,0.3,0.85);
2451  TCanvas* c = new TCanvas(name, title);
2452  c->cd();
2453 
2454  h->Sumw2();
2455 
2456  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/nBins, 1);
2457  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2458 
2459  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,nBins)<<flush;
2460 
2461  // draw projection
2462  TH1D* proj = h->ProjectionY(Form("%s_%d",h->GetName(),i/rebinFactorX),i-rebinFactorX+1,i);
2463  if (proj->GetEntries() > 0) proj->Scale(1./proj->GetEntries());
2464  proj->Draw((i==rebinFactorX)?"hist":"histsames");
2465  proj->SetLineColor(i/rebinFactorX);
2466 
2467  if (proj->GetEntries() > 10) {
2468 
2469  // first fit
2470  fGaus3->SetParameters(1., 0., 1.);
2471  proj->Fit("fGaus3", "WWNQ");
2472 
2473  // rebin histo
2474  Double_t sigma = fGaus3->GetParameter(2);
2475  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*proj->GetNbinsX(),TMath::Max(0.5*sigma/proj->GetBinWidth(1),1.)));
2476  while (proj->GetNbinsX()%rebin!=0) rebin--;
2477  proj->Rebin(rebin);
2478  proj->Scale(1./rebin);
2479 
2480  // second fit
2481  Double_t mean = fGaus3->GetParameter(1);
2482  fLandauGaus2->SetParameters(0.25*sigma*TMath::Sqrt(8.*log(2.)), mean, proj->GetXaxis()->GetBinWidth(1), 0.5*sigma);
2483  fLandauGaus2->SetParLimits(0, 0.0025*sigma*TMath::Sqrt(8.*log(2.)), 1000.);
2484  fLandauGaus2->SetParLimits(3, 0., 2.*sigma);
2485  Double_t xMin = TMath::Max(mean-50.*sigma, proj->GetXaxis()->GetXmin());
2486  Double_t xMax = TMath::Min(mean+10.*sigma, proj->GetXaxis()->GetXmax());
2487  if (xMin < proj->GetXaxis()->GetXmax() && xMax > proj->GetXaxis()->GetXmin()) fLandauGaus2->SetRange(xMin, xMax);
2488  fLandauGaus2->SetLineColor(i/rebinFactorX);
2489  proj->Fit("fLandauGaus2","RQ","sames");
2490 
2491  }
2492 
2493  // set label
2494  Double_t p = 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2495  l->AddEntry(proj,Form("%5.1f GeV",p));
2496 
2497  }
2498 
2499  cout<<Form("\rFitting %s... %d/%d",fitting,nBins,nBins)<<endl;
2500 
2501  l->Draw("same");
2502 
2503  return c;
2504 
2505 }
2506 
2507 //________________________________________________________________________
2508 TCanvas* AliAnalysisTaskMuonPerformance::DrawResPVsP(const char* name, const char* title, TH2* h, const Int_t nBins)
2509 {
2511 
2512  TLegend* l = new TLegend(0.15,0.25,0.3,0.85);
2513  TCanvas* c = new TCanvas(name, title);
2514  c->cd();
2515 
2516  h->Sumw2();
2517 
2518  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/nBins, 1);
2519  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2520 
2521  // draw projection
2522  TH1D* proj = h->ProjectionY(Form("%s_%d",h->GetName(),i/rebinFactorX),i-rebinFactorX+1,i);
2523  if (proj->GetEntries() > 0) proj->Scale(1./proj->GetEntries());
2524  proj->Draw((i==rebinFactorX)?"hist":"histsames");
2525  proj->SetLineColor(i/rebinFactorX);
2526  proj->SetLineWidth(2);
2527 
2528  // set label
2529  Double_t p = 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2530  l->AddEntry(proj,Form("%5.1f GeV",p));
2531 
2532  }
2533 
2534  l->Draw("same");
2535 
2536  return c;
2537 }
2538 
2539 //________________________________________________________________________
2540 void AliAnalysisTaskMuonPerformance::Zoom(TH1* h, Double_t fractionCut)
2541 {
2543 
2544  Double_t maxEventsCut = fractionCut * h->GetEntries();
2545  Int_t nBins = h->GetNbinsX();
2546 
2547  // set low edge
2548  Int_t minBin;
2549  Double_t eventsCut = 0.;
2550  for (minBin = 1; minBin <= nBins; minBin++) {
2551  eventsCut += h->GetBinContent(minBin);
2552  if (eventsCut > maxEventsCut) break;
2553  }
2554 
2555  // set high edge
2556  Int_t maxBin;
2557  eventsCut = 0.;
2558  for (maxBin = nBins; maxBin >= 1; maxBin--) {
2559  eventsCut += h->GetBinContent(maxBin);
2560  if (eventsCut > maxEventsCut) break;
2561  }
2562 
2563  // set new axis range
2564  h->GetXaxis()->SetRange(--minBin, ++maxBin);
2565 }
2566 
2567 //________________________________________________________________________
2568 void AliAnalysisTaskMuonPerformance::FillEffHistos(AliCFEffGrid* efficiency, const char* suffix, TObjArray* list)
2569 {
2571 
2572  TH1* auxHisto = efficiency->Project(kVarPt);
2573  auxHisto->SetName(Form("effVsPt_%s",suffix));
2574  list->AddLast(auxHisto);
2575 
2576  auxHisto = efficiency->Project(kVarEta);
2577  auxHisto->SetName(Form("effVsEta_%s",suffix));
2578  list->AddLast(auxHisto);
2579 
2580  auxHisto = efficiency->Project(kVarPhi);
2581  auxHisto->SetName(Form("effVsPhi_%s",suffix));
2582  list->AddLast(auxHisto);
2583 
2584  auxHisto = efficiency->Project(kVarCent);
2585  auxHisto->SetName(Form("effVsCent_%s",suffix));
2586  list->AddLast(auxHisto);
2587 
2588 }
2589 
momentum residual at first cluster versus P
eta residual at vertex in 3 MC angular regions
TCanvas * DrawVsAng(const char *name, const char *title, TH1 *h1, TH2 *h2)
mean P * MCS deviation angle versus P for tracks in ]3,10[ degrees at absorber end ...
momentum residual at vertex versus P for tracks in ]3,10[ degrees at absorber end ...
ClassImp(AliAnalysisTaskTriggerRates) AliAnalysisTaskTriggerRates
momentum residual at vertex versus position at absorber end in 3 MC angular regions ...
phi residual at vertex in 3 MC angular regions
slope-Y residual at vertex in 3 angular regions
eta residual at vertex versus position at absorber end in 3 MC angular regions
P * DCA resolution versus P for tracks in ]2,3] degrees at absorber end.
void Zoom(TH1 *h, Double_t fractionCut=0.01)
Int_t RecoTrackMother(AliMCParticle *mcParticle)
slope-Y residual at vertex versus position at absorber end in 3 MC angular regions ...
phi residual at vertex versus angle at absorber end
TObjArray * fSlopeAtVtxList
List of graph and canvas about slope resolution at vertex.
Bool_t fFitResiduals
fit or not the cluster residuals to extract means and sigmas
P * DCA resolution versus P for tracks in ]3,10[ degrees at absorber end.
slope-X residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
slope-X residual at vertex in 3 MC angular regions
Int_t fDEIndices[1100]
index of DE in histograms refered by ID
Bool_t fUseMCKinematics
use kinematics from MC track when available
eta residual at vertex versus position at absorber end for tracks with MC angle in ]3...
transverse momentum residual at vertex versus pT
P * MCS deviation angle distribution versus P for tracks in ]2,3] degrees at absorber end...
momentum residual at vertex in 3 MC angular regions
TObjArray * fSlopeAt1stClList
List of graph and canvas about slope resolution at first cluster.
P * DCA distribution versus position at absorber end for tracks with MC angle in ]3,10[ degrees.
TObjArray * fClusterList
List of graph and canvas about cluster resolution.
momentum residual at vertex versus angle at absorber end versus momentum
void FillEffHistos(AliCFEffGrid *efficiency, const char *suffix, TObjArray *list)
slope-X residual at vertex in 3 angular regions
TList * list
void FitClusterResidual(TH1 *h, Int_t i, Double_t &sigma, TGraphErrors *gMean, TGraphErrors *gSigma)
eta residual at vertex versus position at absorber end for tracks with MC angle in ]2...
matched with reconstructible track and triggerable track of same ID
P * DCA versus position at absorber end in 3 MC angular regions.
TObjArray * fTriggerList
List of histograms for trigger resolution.
Double_t ptMin
matched with reconstructible track and triggerable track of different ID
slope-Y residual at vertex versus position at absorber end for tracks with MC angle in ]3...
eta residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
void FitPDCAVsMom(TH2 *h, const char *fitting, TGraphAsymmErrors *gMean, TGraphAsymmErrors *gSigma)
phi residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
P * DCA distribution versus position at absorber end for tracks with MC angle in ]2,3] degrees.
mean P * DCA versus P for tracks in ]2,3] degrees at absorber end
UInt_t fRequestedStationMask
mask of requested stations
mean P * DCA versus P for tracks in ]3,10[ degrees at absorber end
momentum residual at vertex versus position at absorber end for tracks with MC angle in ]3...
TCanvas * DrawResPVsP(const char *name, const char *title, TH2 *h, const Int_t nBins)
relative momentum resolution at first cluster versus P
void FillContainerInfoMC(Double_t *containerInput, AliMCParticle *mcPart)
P * DCA distribution versus angle at absorber end.
mean P * MCS deviation angle versus P for tracks in ]2,3] degrees at absorber end ...
phi residual at vertex versus position at absorber end in 3 MC angular regions
TObjArray * fPhiAtVtxList
List of graph and canvas about phi resolution at vertex.
phi residual at vertex versus position at absorber end for tracks with MC angle in ]3...
momentum residual at vertex in 3 angular regions
Double_t fSigmaCutTrig
sigma cut to associate trigger track to triggerable track
P * DCA distribution versus position at absorber end for tracks with MC angle <= 2 degrees...
void FitLandauGausResVsP(TH2 *h, const char *fitting, TGraphAsymmErrors *gMean, TGraphAsymmErrors *gMostProb, TGraphAsymmErrors *gSigma)
TCanvas * DrawFitLandauGausResPVsP(const char *name, const char *title, TH2 *h, const Int_t nBins, const char *fitting)
TObjArray * fEtaAtVtxList
List of graph and canvas about eta resolution at vertex.
TCanvas * DrawVsPos(const char *name, const char *title, TH2 *h1, TH2 *h2, TH2 *h3)
P * DCA distribution versus P for tracks in ]3,10[ degrees at absorber end.
TObjArray * fEfficiencyList
List of histograms for tracker/trigger efficiencies.
phi residual at vertex in 3 angular regions
Bool_t fMCTrigLevelFromMatchTrk
set the triggerable level from the MC matching the trk part
momentum residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
transverse momentum residual at first cluster versus pT
Bool_t fRequest2ChInSameSt45
2 fired chambers requested in the same station (4 or 5) or not
slope-X residual at vertex versus position at absorber end for tracks with MC angle in ]3...
TObjArray * fPAtVtxList
List of graph and canvas about momentum resolution at vertex.
momentum residual at vertex versus P for tracks in ]2,3] degrees at absorber end
TObjArray * fDCAList
List of graph and canvas about DCA.
mean momentum residual at first cluster versus P
Residual of x position in first trigger chamber.
P * DCA distribution versus P for tracks in ]2,3] degrees at absorber end.
TObjArray * fTrackerList
List of histograms for tracker resolution.
TString fDefaultStorage
location of the default OCDB storage
TObjArray * fPAt1stClList
List of graph and canvas about momentum resolution at first cluster.
eta residual at vertex in 3 angular regions
Double_t centMax
slope-X residual at vertex versus position at absorber end in 3 MC angular regions ...
Int_t fDEIds[200]
ID of DE refered by index in histograms.
slope-Y residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
slope-X residual at vertex versus angle at absorber end
momentum residual at vertex versus position at absorber end for tracks with MC angle in ]2...
slope-Y residual at vertex versus angle at absorber end
AliCFContainer * fCFContainer
Pointer to the CF container.
momentum residual for tracks between 3 and 10 degrees
not matched with either reconstructible track or triggerable track
void FitGausResVsMom(TH2 *h, const Double_t mean0, const Double_t sigma0, const char *fitting, TGraphAsymmErrors *gMean, TGraphAsymmErrors *gSigma)
Bool_t fEnforceTrkCriteria
select reconstructed tracks still reconstructible
momentum residual at vertex versus P for tracks with MC angle < 2 degrees
momentum residual at vertex versus angle at absorber end
phi residual at vertex versus position at absorber end for tracks with MC angle in ]2...
TString fRecoParamOCDBpath
OCDB path to the recoParam file.
Float_t GetBinThetaAbsEnd(Float_t RAtAbsEnd, Bool_t isTheta=kFALSE)
Residual of y position in first trigger chamber.
Double_t ptMax
P * MCS deviation angle dispersion versus P for tracks in ]2,3] degrees at absorber end...
Bool_t fCorrectForSystematics
add or not the systematic shifts of the residuals to the resolution
slope-X residual at vertex versus position at absorber end for tracks with MC angle in ]2...
momentum residual for tracks between 2 and 3 degrees
Double_t langaufun(Double_t *x, Double_t *par)
slope-Y residual at vertex in 3 MC angular regions
P * MCS deviation angle distribution versus P for tracks in ]3,10[ degrees at absorber end...
void FillContainerInfoReco(Double_t *containerInput, AliESDMuonTrack *esdTrack, Bool_t isValid, Int_t mcID)
eta residual at vertex versus angle at absorber end
momentum residuals for tracks with MC angle < 2 degrees
relative momentum resolution at vertex versus P
TString fAlignOCDBpath
OCDB path to the alignment file.
Match trigger not passing any Pt threshold (MC)
slope-Y residual at vertex versus position at absorber end for tracks with MC angle in ]2...
most probable momentum residual at vertex versus P
Double_t centMin
P * MCS deviation angle dispersion versus P for tracks in ]3,10[ degrees at absorber end...
Double_t fClusterMaxRes[2]
highest chamber resolution in both directions
Bool_t GetEfficiency(AliCFEffGrid *efficiency, Double_t &calcEff, Double_t &calcEffErr)