AliPhysics  41af4b0 (41af4b0)
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 "AliMultSelection.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 //________________________________________________________________________
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 //________________________________________________________________________
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  AliMultSelection *multSelection = static_cast<AliMultSelection*>(esd->FindListObject("MultSelection"));
663  Double_t centrality = multSelection ? multSelection->GetMultiplicityPercentile("V0M") : -1.;
664 
665  // Get reference tracks
666  AliMUONRecoCheck rc(esd,mcH);
667  AliMUONVTriggerTrackStore* triggerTrackRefStore = rc.TriggerableTracks(-1);
668  AliMUONVTrackStore* trackRefStore = rc.TrackRefs(-1);
669  AliMUONVTrackStore* reconstructibleStore = rc.ReconstructibleTracks(-1, fRequestedStationMask, fRequest2ChInSameSt45);
670 
671  Double_t containerInput[kNvars];
672  containerInput[kVarCent] = centrality;
673  AliMUONTrackParam *trackParam;
674  Double_t x1,y1,z1,slopex1,slopey1,pX1,pY1,pZ1,p1,pT1,eta1,phi1;
675  Double_t x2,y2,z2,slopex2,slopey2,pX2,pY2,pZ2,p2,pT2,eta2,phi2;
676  Double_t dPhi;
677  Double_t xAbs,yAbs,dAbs,aAbs,aMCS,aMC;
678  Double_t xDCA,yDCA,dca,pU;
679 
680  // ------ Loop over reconstructed tracks ------
681  AliESDMuonTrack *esdTrack = 0x0;
682  Int_t nMuTracks = esd->GetNumberOfMuonTracks();
683  Int_t *loCircuit = new Int_t[nMuTracks];
684  Int_t nTrgTracks = 0;
685  for (Int_t iMuTrack = 0; iMuTrack < nMuTracks; ++iMuTrack) {
686 
687  esdTrack = esd->GetMuonTrack(iMuTrack);
688 
689  // Tracker
690  AliMUONTrack* matchedTrackRef = 0x0;
691  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
692  Bool_t isValid = kFALSE;
693  if (esdTrack->ContainTrackerData()) {
694 
695  // convert ESD track to MUON track (without recomputing track parameters at each clusters)
696  AliMUONTrack muonTrack;
697  AliMUONESDInterface::ESDToMUON(*esdTrack, muonTrack, kFALSE);
698 
699  // decide whether the track is valid for efficiency calculations
700  isValid = (!fEnforceTrkCriteria || muonTrack.IsValid(fRequestedStationMask, fRequest2ChInSameSt45));
701 
702  // get the associated simulated track (discard decays)
703  Int_t mcLabel = esdTrack->GetLabel();
704  if (mcLabel >= 0 && !esdTrack->TestBit(BIT(22)))
705  matchedTrackRef = static_cast<AliMUONTrack*>(trackRefStore->FindObject(mcLabel));
706  if (matchedTrackRef && isValid) containerInput[kVarMatchMC] = static_cast<Double_t>(kTrackerOnly);
707 
708  // compute track resolution (discard not-reconstructible trackRef)
709  if (matchedTrackRef && !esdTrack->TestBit(BIT(23))) {
710 
711  // simulated track parameters at vertex
712  trackParam = matchedTrackRef->GetTrackParamAtVertex();
713  x1 = trackParam->GetNonBendingCoor();
714  y1 = trackParam->GetBendingCoor();
715  z1 = trackParam->GetZ();
716  slopex1 = trackParam->GetNonBendingSlope();
717  slopey1 = trackParam->GetBendingSlope();
718  pX1 = trackParam->Px();
719  pY1 = trackParam->Py();
720  pZ1 = trackParam->Pz();
721  p1 = trackParam->P();
722  pT1 = TMath::Sqrt(pX1*pX1 + pY1*pY1);
723  aMC = TMath::ATan(-pT1/pZ1) * TMath::RadToDeg();
724  eta1 = TMath::Log(TMath::Tan(0.5*TMath::ATan(-pT1/pZ1)));
725  phi1 = TMath::Pi()+TMath::ATan2(-pY1, -pX1);
726 
727  // reconstructed track parameters at the end of the absorber
728  AliMUONTrackParam trackParamAtAbsEnd(*((AliMUONTrackParam*)muonTrack.GetTrackParamAtCluster()->First()));
729  AliMUONTrackExtrap::ExtrapToZ(&trackParamAtAbsEnd, AliMUONConstants::AbsZEnd());
730  xAbs = trackParamAtAbsEnd.GetNonBendingCoor();
731  yAbs = trackParamAtAbsEnd.GetBendingCoor();
732  dAbs = TMath::Sqrt(xAbs*xAbs + yAbs*yAbs);
733  aAbs = TMath::ATan(-dAbs/AliMUONConstants::AbsZEnd()) * TMath::RadToDeg();
734  pX2 = trackParamAtAbsEnd.Px();
735  pY2 = trackParamAtAbsEnd.Py();
736  pZ2 = trackParamAtAbsEnd.Pz();
737  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
738  aMCS = TMath::ATan(-pT2/pZ2) * TMath::RadToDeg();
739 
740  // reconstructed track parameters at vertex
741  trackParam = muonTrack.GetTrackParamAtVertex();
742  x2 = trackParam->GetNonBendingCoor();
743  y2 = trackParam->GetBendingCoor();
744  z2 = trackParam->GetZ();
745  slopex2 = trackParam->GetNonBendingSlope();
746  slopey2 = trackParam->GetBendingSlope();
747  pX2 = trackParam->Px();
748  pY2 = trackParam->Py();
749  pZ2 = trackParam->Pz();
750  p2 = trackParam->P();
751  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
752  eta2 = TMath::Log(TMath::Tan(0.5*TMath::ATan(-pT2/pZ2)));
753  phi2 = TMath::Pi()+TMath::ATan2(-pY2, -pX2);
754 
755  // reconstructed track parameters at DCA
756  AliMUONTrackParam trackParamAtDCA(*((AliMUONTrackParam*) muonTrack.GetTrackParamAtCluster()->First()));
757  pU = trackParamAtDCA.P();
758  AliMUONTrackExtrap::ExtrapToVertexWithoutBranson(&trackParamAtDCA, z2);
759  xDCA = trackParamAtDCA.GetNonBendingCoor();
760  yDCA = trackParamAtDCA.GetBendingCoor();
761  dca = TMath::Sqrt(xDCA*xDCA + yDCA*yDCA);
762 
763  dPhi = phi2-phi1;
764  if (dPhi < -TMath::Pi()) dPhi += 2.*TMath::Pi();
765  else if (dPhi > TMath::Pi()) dPhi -= 2.*TMath::Pi();
766 
767  // fill histograms
768  static_cast<TH1*>(fTrackerList->At(kResPAtVtx))->Fill(p2-p1);
769  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsP))->Fill(p1,p2-p1);
770  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsAngleAtAbsEnd))->Fill(aAbs,p2-p1);
771  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsMCAngle))->Fill(aMC,p2-p1);
772  static_cast<TH3*>(fTrackerList->At(kResPAtVtxVsAngleAtAbsEndVsP))->Fill(p1,aAbs,p2-p1);
773  static_cast<TH2*>(fTrackerList->At(kResPtAtVtxVsPt))->Fill(pT1,pT2-pT1);
774 
775  static_cast<TH1*>(fTrackerList->At(kResSlopeXAtVtx))->Fill(slopex2-slopex1);
776  static_cast<TH1*>(fTrackerList->At(kResSlopeYAtVtx))->Fill(slopey2-slopey1);
777  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsP))->Fill(p1,slopex2-slopex1);
778  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsP))->Fill(p1,slopey2-slopey1);
779  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsAngleAtAbsEnd))->Fill(aAbs,slopex2-slopex1);
780  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsAngleAtAbsEnd))->Fill(aAbs,slopey2-slopey1);
781  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsMCAngle))->Fill(aMC,slopex2-slopex1);
782  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsMCAngle))->Fill(aMC,slopey2-slopey1);
783 
784  static_cast<TH1*>(fTrackerList->At(kResEtaAtVtx))->Fill(eta2-eta1);
785  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsP))->Fill(p1,eta2-eta1);
786  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsAngleAtAbsEnd))->Fill(aAbs,eta2-eta1);
787  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsMCAngle))->Fill(aMC,eta2-eta1);
788 
789  static_cast<TH1*>(fTrackerList->At(kResPhiAtVtx))->Fill(dPhi);
790  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsP))->Fill(p1,dPhi);
791  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsAngleAtAbsEnd))->Fill(aAbs,dPhi);
792  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsMCAngle))->Fill(aMC,dPhi);
793 
794  static_cast<TH1*>(fTrackerList->At(kPDCA))->Fill(0.5*(p2+pU)*dca);
795  static_cast<TH2*>(fTrackerList->At(kPDCAVsAngleAtAbsEnd))->Fill(aAbs,0.5*(p2+pU)*dca);
796  static_cast<TH2*>(fTrackerList->At(kPDCAVsMCAngle))->Fill(aMC,0.5*(p2+pU)*dca);
797 
798  if (aAbs > 2. && aAbs < 3.) {
799 
800  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn23deg))->Fill(p1,p2-p1);
801  static_cast<TH2*>(fTrackerList->At(kPDCAVsPIn23deg))->Fill(p1,0.5*(p2+pU)*dca);
802  static_cast<TH2*>(fTrackerList->At(kPMCSAngVsPIn23deg))->Fill(p1,0.5*(p2+pU)*(aMCS-aMC)*TMath::DegToRad());
803 
804  } else if (aAbs >= 3. && aAbs < 10.) {
805 
806  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn310deg))->Fill(p1,p2-p1);
807  static_cast<TH2*>(fTrackerList->At(kPDCAVsPIn310deg))->Fill(p1,0.5*(p2+pU)*dca);
808  static_cast<TH2*>(fTrackerList->At(kPMCSAngVsPIn310deg))->Fill(p1,0.5*(p2+pU)*(aMCS-aMC)*TMath::DegToRad());
809 
810  }
811 
812  if (aMC < 2.) {
813 
814  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPIn02degMC))->Fill(p1,p2-p1);
815  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,p2-p1);
816  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,slopex2-slopex1);
817  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,slopey2-slopey1);
818  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,eta2-eta1);
819  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn02degMC))->Fill(dAbs,dPhi);
820  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn02degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
821 
822  } else if (aMC >= 2. && aMC < 3) {
823 
824  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,p2-p1);
825  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,slopex2-slopex1);
826  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,slopey2-slopey1);
827  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,eta2-eta1);
828  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn23degMC))->Fill(dAbs,dPhi);
829  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn23degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
830 
831  } else if (aMC >= 3. && aMC < 10.) {
832 
833  static_cast<TH2*>(fTrackerList->At(kResPAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,p2-p1);
834  static_cast<TH2*>(fTrackerList->At(kResSlopeXAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,slopex2-slopex1);
835  static_cast<TH2*>(fTrackerList->At(kResSlopeYAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,slopey2-slopey1);
836  static_cast<TH2*>(fTrackerList->At(kResEtaAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,eta2-eta1);
837  static_cast<TH2*>(fTrackerList->At(kResPhiAtVtxVsPosAbsEndIn310degMC))->Fill(dAbs,dPhi);
838  static_cast<TH2*>(fTrackerList->At(kPDCAVsPosAbsEndIn310degMC))->Fill(dAbs,0.5*(p2+pU)*dca);
839 
840  }
841 
842  // simulated track parameters at first cluster
843  trackParam = (AliMUONTrackParam*) matchedTrackRef->GetTrackParamAtCluster()->First();
844  x1 = trackParam->GetNonBendingCoor();
845  y1 = trackParam->GetBendingCoor();
846  z1 = trackParam->GetZ();
847  slopex1 = trackParam->GetNonBendingSlope();
848  slopey1 = trackParam->GetBendingSlope();
849  pX1 = trackParam->Px();
850  pY1 = trackParam->Py();
851  pZ1 = trackParam->Pz();
852  p1 = trackParam->P();
853  pT1 = TMath::Sqrt(pX1*pX1 + pY1*pY1);
854 
855  // reconstructed track parameters at first cluster
856  trackParam = (AliMUONTrackParam*) muonTrack.GetTrackParamAtCluster()->First();
857  x2 = trackParam->GetNonBendingCoor();
858  y2 = trackParam->GetBendingCoor();
859  z2 = trackParam->GetZ();
860  slopex2 = trackParam->GetNonBendingSlope();
861  slopey2 = trackParam->GetBendingSlope();
862  pX2 = trackParam->Px();
863  pY2 = trackParam->Py();
864  pZ2 = trackParam->Pz();
865  p2 = trackParam->P();
866  pT2 = TMath::Sqrt(pX2*pX2 + pY2*pY2);
867 
868  // fill histograms
869  static_cast<TH1*>(fTrackerList->At(kResPAt1stCl))->Fill(p2-p1);
870  static_cast<TH2*>(fTrackerList->At(kResPAt1stClVsP))->Fill(p1,p2-p1);
871  static_cast<TH2*>(fTrackerList->At(kResPtAt1stClVsPt))->Fill(pT1,pT2-pT1);
872  static_cast<TH1*>(fTrackerList->At(kResSlopeXAt1stCl))->Fill(slopex2-slopex1);
873  static_cast<TH1*>(fTrackerList->At(kResSlopeYAt1stCl))->Fill(slopey2-slopey1);
874  static_cast<TH2*>(fTrackerList->At(kResSlopeXAt1stClVsP))->Fill(p1,slopex2-slopex1);
875  static_cast<TH2*>(fTrackerList->At(kResSlopeYAt1stClVsP))->Fill(p1,slopey2-slopey1);
876 
877  // clusters residuals
878  for (Int_t iCl1 = 0; iCl1 < muonTrack.GetNClusters(); iCl1++) {
879 
880  AliMUONVCluster* cluster1 = static_cast<AliMUONTrackParam*>(muonTrack.GetTrackParamAtCluster()->UncheckedAt(iCl1))->GetClusterPtr();
881  Int_t chId = cluster1->GetChamberId();
882  Int_t deId = cluster1->GetDetElemId();
883 
884  for (Int_t iCl2 = 0; iCl2 < matchedTrackRef->GetNClusters(); iCl2++) {
885 
886  AliMUONVCluster* cluster2 = static_cast<AliMUONTrackParam*>(matchedTrackRef->GetTrackParamAtCluster()->UncheckedAt(iCl2))->GetClusterPtr();
887 
888  if (cluster2->GetDetElemId() == deId) {
889 
890  // fill histograms
891  static_cast<TH2*>(fTrackerList->At(kResClXVsCh))->Fill(chId+1, cluster1->GetX()-cluster2->GetX());
892  static_cast<TH2*>(fTrackerList->At(kResClYVsCh))->Fill(chId+1, cluster1->GetY()-cluster2->GetY());
893  static_cast<TH2*>(fTrackerList->At(kResClXVsDE))->Fill(fDEIndices[deId], cluster1->GetX()-cluster2->GetX());
894  static_cast<TH2*>(fTrackerList->At(kResClYVsDE))->Fill(fDEIndices[deId], cluster1->GetY()-cluster2->GetY());
895 
896  break;
897 
898  }
899 
900  }
901 
902  }
903 
904  } // end if (matchedTrackRef)
905 
906  } // end if (esdTrack->ContainTrackerData())
907 
908  // look for MC trigger associated to MC track
909  AliMUONTriggerTrack* triggerTrackRef = (isValid && matchedTrackRef)
910  ? static_cast<AliMUONTriggerTrack*>(triggerTrackRefStore->FindObject(matchedTrackRef->GetUniqueID()))
911  : 0x0;
912 
913  // Trigger
914  AliMUONTriggerTrack* matchedTrigTrackRef = 0x0;
915  containerInput[kVarDupliTrg] = 0.;
916  if (esdTrack->ContainTriggerData()) {
917 
918  // check if this track is not already accounted for and record it if not
919  Bool_t trackExist = kFALSE;
920  for (Int_t i=0; i<nTrgTracks; i++) if (esdTrack->LoCircuit() == loCircuit[i]) trackExist = kTRUE;
921  if (trackExist) containerInput[kVarDupliTrg] = 1.;
922  else loCircuit[nTrgTracks++] = esdTrack->LoCircuit();
923 
924  // Convert ESD track to trigger track
925  AliMUONLocalTrigger locTrg;
926  AliMUONESDInterface::ESDToMUON(*esdTrack, locTrg);
927  AliMUONTriggerTrack trigTrack;
928  rc.TriggerToTrack(locTrg, trigTrack);
929 
930  // try to match the trigger track with the same MC track if any
931  if (triggerTrackRef && trigTrack.Match(*triggerTrackRef, fSigmaCutTrig)) matchedTrigTrackRef = triggerTrackRef;
932  if (matchedTrigTrackRef) containerInput[kVarMatchMC] = static_cast<Double_t>(kMatchedSame);
933  else { // or try to match with any triggerable track
934  matchedTrigTrackRef = rc.FindCompatibleTrack(trigTrack, *triggerTrackRefStore, fSigmaCutTrig);
935  if (matchedTrigTrackRef) {
936  if (isValid && matchedTrackRef) {
937  containerInput[kVarMatchMC] = static_cast<Double_t>(kMatchedDiff);
938  if (!fMCTrigLevelFromMatchTrk) triggerTrackRef = matchedTrigTrackRef;
939  } else {
940  containerInput[kVarMatchMC] = static_cast<Double_t>(kTriggerOnly);
941  triggerTrackRef = matchedTrigTrackRef;
942  }
943  }
944  }
945 
946  // fill histograms
947  if (matchedTrigTrackRef) {
948  static_cast<TH1*>(fTriggerList->At(kResTrigX11))->Fill(trigTrack.GetX11() - matchedTrigTrackRef->GetX11());
949  static_cast<TH1*>(fTriggerList->At(kResTrigY11))->Fill(trigTrack.GetY11() - matchedTrigTrackRef->GetY11());
950  static_cast<TH1*>(fTriggerList->At(kResTrigSlopeY))->Fill(trigTrack.GetSlopeY() - matchedTrigTrackRef->GetSlopeY());
951  }
952 
953  }
954 
955  // fill container
956  if (triggerTrackRef) {
957  if (triggerTrackRef->GetPtCutLevel() == 0) containerInput[kVarMCTrigger] = static_cast<Double_t>(kOtherTrig);
958  else if (triggerTrackRef->GetPtCutLevel() == 1) containerInput[kVarMCTrigger] = static_cast<Double_t>(kAllPtTrig);
959  else if (triggerTrackRef->GetPtCutLevel() == 2) containerInput[kVarMCTrigger] = static_cast<Double_t>(kLowPtTrig);
960  else if (triggerTrackRef->GetPtCutLevel() == 3) containerInput[kVarMCTrigger] = static_cast<Double_t>(kHighPtTrig);
961  } else containerInput[kVarMCTrigger] = kNoMatchTrig;
962 
963  // get MC particle ID
964  Int_t mcID = -1;
965  if (isValid && matchedTrackRef) mcID = static_cast<Int_t>(matchedTrackRef->GetUniqueID());
966  else if (matchedTrigTrackRef) mcID = static_cast<Int_t>(matchedTrigTrackRef->GetUniqueID());
967 
968  // fill particle container
969  FillContainerInfoReco(containerInput, esdTrack, isValid, mcID);
970  fCFContainer->Fill(containerInput, kStepReconstructed);
971 
972  }
973 
974  // clean memory
975  delete[] loCircuit;
976  containerInput[kVarDupliTrg] = 0.;
977 
978  // ------ Loop over reconstructible tracks ------
979  AliMUONTrack* trackRef = 0x0;
980  TIter next(reconstructibleStore->CreateIterator());
981  while ((trackRef = static_cast<AliMUONTrack*>(next()))) {
982 
983  // find the corresponding triggerable track if any
984  UInt_t mcID = trackRef->GetUniqueID();
985  AliMUONTriggerTrack* trigTrackRef = static_cast<AliMUONTriggerTrack*>(triggerTrackRefStore->FindObject(mcID));
986 
987  // fill particle container
988  FillContainerInfoMC(containerInput, static_cast<AliMCParticle*>(fMCEvent->GetTrack(static_cast<Int_t>(mcID))));
989  containerInput[kVarHasTracker] = 1.;
990  if (trigTrackRef) {
991  if (trigTrackRef->GetPtCutLevel() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kOtherTrig);
992  else if (trigTrackRef->GetPtCutLevel() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
993  else if (trigTrackRef->GetPtCutLevel() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
994  else if (trigTrackRef->GetPtCutLevel() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
995  } else containerInput[kVarTrigger] = static_cast<Double_t>(kNoMatchTrig);
996  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
997  containerInput[kVarMCTrigger] = static_cast<Double_t>(kNoMatchTrig);
998  fCFContainer->Fill(containerInput, kStepGeneratedMC);
999 
1000  }
1001 
1002  // ------ Loop over triggerable ghosts ------
1003  AliMUONTriggerTrack* trigTrackRef = 0x0;
1004  TIter nextTrig(triggerTrackRefStore->CreateIterator());
1005  while ((trigTrackRef = static_cast<AliMUONTriggerTrack*>(nextTrig()))) {
1006 
1007  // discard tracks also reconstructible
1008  UInt_t mcID = trigTrackRef->GetUniqueID();
1009  if (reconstructibleStore->FindObject(mcID)) continue;
1010 
1011  // fill particle container
1012  FillContainerInfoMC(containerInput, static_cast<AliMCParticle*>(fMCEvent->GetTrack(static_cast<Int_t>(mcID))));
1013  containerInput[kVarHasTracker] = 0.;
1014  if (trigTrackRef->GetPtCutLevel() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kOtherTrig);
1015  else if (trigTrackRef->GetPtCutLevel() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
1016  else if (trigTrackRef->GetPtCutLevel() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
1017  else if (trigTrackRef->GetPtCutLevel() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
1018  containerInput[kVarMatchMC] = static_cast<Double_t>(kNoMatch);
1019  containerInput[kVarMCTrigger] = static_cast<Double_t>(kNoMatchTrig);
1020  fCFContainer->Fill(containerInput, kStepGeneratedMC);
1021 
1022  }
1023 
1024  // Post final data
1025  PostData(1, fCFContainer);
1026  PostData(2, fTriggerList);
1027  PostData(3, fTrackerList);
1028 }
1029 
1030 //________________________________________________________________________
1032 {
1033  //
1035  //
1036 
1037  // do it only locally
1038  //if (gROOT->IsBatch()) return;
1039 
1040  // get output containers
1041  fCFContainer = dynamic_cast<AliCFContainer*>(GetOutputData(1));
1042  fTriggerList = dynamic_cast<TObjArray*>(GetOutputData(2));
1043  fTrackerList = dynamic_cast<TObjArray*>(GetOutputData(3));
1044  if (!fCFContainer || !fTriggerList || !fTrackerList) {
1045  AliWarning("Output containers not found: summary histograms are not created");
1046  return;
1047  }
1048 
1049  // --- compute efficiencies ---
1050  fEfficiencyList = new TObjArray(100);
1051  fEfficiencyList->SetOwner();
1052 
1053  TObjArray* effAnyPt = new TObjArray(100);
1054  effAnyPt->SetName("effAnyPt");
1055  effAnyPt->SetOwner();
1056  fEfficiencyList->AddLast(effAnyPt);
1057 
1058  TObjArray* effAllPt = new TObjArray(100);
1059  effAllPt->SetName("effAllPt");
1060  effAllPt->SetOwner();
1061  fEfficiencyList->AddLast(effAllPt);
1062 
1063  TObjArray* effLowPt = new TObjArray(100);
1064  effLowPt->SetName("effLowPt");
1065  effLowPt->SetOwner();
1066  fEfficiencyList->AddLast(effLowPt);
1067 
1068  TObjArray* effHighPt = new TObjArray(100);
1069  effHighPt->SetName("effHighPt");
1070  effHighPt->SetOwner();
1071  fEfficiencyList->AddLast(effHighPt);
1072 
1073  TObjArray* notTrgable = new TObjArray(100);
1074  notTrgable->SetName("notTrgable");
1075  notTrgable->SetOwner();
1076  fEfficiencyList->AddLast(notTrgable);
1077 
1078  TObjArray* trgableNoPtOnly = new TObjArray(100);
1079  trgableNoPtOnly->SetName("trgableNoPtOnly");
1080  trgableNoPtOnly->SetOwner();
1081  fEfficiencyList->AddLast(trgableNoPtOnly);
1082 
1083  TObjArray* trgableAPtOnly = new TObjArray(100);
1084  trgableAPtOnly->SetName("trgableAPtOnly");
1085  trgableAPtOnly->SetOwner();
1086  fEfficiencyList->AddLast(trgableAPtOnly);
1087 
1088  TObjArray* trgableLPtOnly = new TObjArray(100);
1089  trgableLPtOnly->SetName("trgableLPtOnly");
1090  trgableLPtOnly->SetOwner();
1091  fEfficiencyList->AddLast(trgableLPtOnly);
1092 
1093  TObjArray* trgableHPtOnly = new TObjArray(100);
1094  trgableHPtOnly->SetName("trgableHPtOnly");
1095  trgableHPtOnly->SetOwner();
1096  fEfficiencyList->AddLast(trgableHPtOnly);
1097 
1098  AliCFEffGrid* efficiency = new AliCFEffGrid("eff","",*fCFContainer);
1099  efficiency->CalculateEfficiency(kStepReconstructed, kStepGeneratedMC);
1100  Double_t totalEff = 0., totalEffErr = 0.;
1101  Int_t sumEffBin = 0;
1102 
1103  // add histogram summarizing global efficiencies
1104  TH1D* effSummary = new TH1D("effSummary", "Efficiency summary", 1, 0., 0.);
1105  effSummary->GetYaxis()->SetTitle("Efficiency");
1106  fEfficiencyList->AddLast(effSummary);
1107 
1108  // ------ Tracker only ------
1109 
1110  // Tracker efficiency using all reconstructed tracks
1111  efficiency->GetNum()->SetRangeUser(kVarHasTracker, 1., 1.);
1112  efficiency->GetDen()->SetRangeUser(kVarHasTracker, 1., 1.);
1113  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1114  efficiency->GetDen()->GetAxis(kVarTrigger)->SetRange();
1115  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1116  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1117  FillEffHistos(efficiency, "trackerTracks", fEfficiencyList);
1118  GetEfficiency(efficiency, totalEff, totalEffErr);
1119  effSummary->Fill("Tracker_all", totalEff);
1120  effSummary->SetBinError(++sumEffBin, totalEffErr);
1121  printf("Tracker efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1122 
1123  // Tracker efficiency using tracks matched with reconstructible ones
1124  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1125  FillEffHistos(efficiency, "trackerTracksMatchMC", fEfficiencyList);
1126  GetEfficiency(efficiency, totalEff, totalEffErr);
1127  effSummary->Fill("Tracker_MCId", totalEff);
1128  effSummary->SetBinError(++sumEffBin, totalEffErr);
1129  printf("Tracker efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1130 
1131  // ------ Tracker matched with any trigger ------
1132 
1133  // Matched efficiency using all reconstructed tracks
1134  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1135  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kHighPtTrig);
1136  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1137  FillEffHistos(efficiency, "matchedTracks", effAnyPt);
1138  GetEfficiency(efficiency, totalEff, totalEffErr);
1139  effSummary->Fill("Matched_all", totalEff);
1140  effSummary->SetBinError(++sumEffBin, totalEffErr);
1141  printf("Matched efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1142 
1143  // Matched efficiency using tracks matched with reconstructible ones, triggerable or not
1144  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1145  FillEffHistos(efficiency, "matchedTracksMatchMC", effAnyPt);
1146  GetEfficiency(efficiency, totalEff, totalEffErr);
1147  effSummary->Fill("Matched_MCId", totalEff);
1148  effSummary->SetBinError(++sumEffBin, totalEffErr);
1149  printf("Matched efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1150 
1151  // Matched efficiency using tracks matched with reconstructible ones triggerable
1152  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kHighPtTrig);
1153  FillEffHistos(efficiency, "matchedTracksMatchMCAnypt", effAnyPt);
1154  GetEfficiency(efficiency, totalEff, totalEffErr);
1155  effSummary->Fill("Matched_MCIdAnypt", totalEff);
1156  effSummary->SetBinError(++sumEffBin, totalEffErr);
1157  printf("Matched efficiency using reconstructed tracks matching MC-anyPt = %f +- %f\n", totalEff, totalEffErr);
1158 
1159  // Matched efficiency using tracks matched with reconstructible ones not triggerable
1160  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1161  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effAnyPt);
1162  GetEfficiency(efficiency, totalEff, totalEffErr);
1163  effSummary->Fill("Matched_MCIdNoTrig", totalEff);
1164  effSummary->SetBinError(++sumEffBin, totalEffErr);
1165  printf("Matched efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1166 
1167  // Matched efficiency using tracks matched with same reconstructible & triggerable ones
1168  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1169  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1170  FillEffHistos(efficiency, "matchedTracksMatchSameMC", effAnyPt);
1171  GetEfficiency(efficiency, totalEff, totalEffErr);
1172  effSummary->Fill("Matched_SameMCId", totalEff);
1173  effSummary->SetBinError(++sumEffBin, totalEffErr);
1174  printf("Matched efficiency using reconstructed tracks matching same MC = %f +- %f\n", totalEff, totalEffErr);
1175 
1176  // Matched efficiency using tracks matched with different reconstructible & triggerable ones
1177  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1178  FillEffHistos(efficiency, "matchedTracksMatchDiffMC", effAnyPt);
1179  GetEfficiency(efficiency, totalEff, totalEffErr);
1180  effSummary->Fill("Matched_DiffMCId", totalEff);
1181  effSummary->SetBinError(++sumEffBin, totalEffErr);
1182  printf("Matched efficiency using reconstructed tracks matching different MC = %f +- %f\n", totalEff, totalEffErr);
1183 
1184  // Matched efficiency using tracks matched with reconstructible ones only
1185  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1186  FillEffHistos(efficiency, "matchedTracksMatchTrkMC", effAnyPt);
1187  GetEfficiency(efficiency, totalEff, totalEffErr);
1188  effSummary->Fill("Matched_TrkMCId", totalEff);
1189  effSummary->SetBinError(++sumEffBin, totalEffErr);
1190  printf("Matched efficiency using reconstructed tracks matching tracker MC = %f +- %f\n", totalEff, totalEffErr);
1191 
1192  // ------ Tracker matched with all pt trigger ------
1193 
1194  // Matched all pt efficiency using all reconstructed tracks
1195  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1196  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1197  FillEffHistos(efficiency, "matchedTracks", effAllPt);
1198  GetEfficiency(efficiency, totalEff, totalEffErr);
1199  printf("Matched Apt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1200 
1201  // Matched all pt efficiency using tracks matched with reconstructible ones, triggerable or not
1202  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1203  FillEffHistos(efficiency, "matchedTracksMatchMC", effAllPt);
1204  GetEfficiency(efficiency, totalEff, totalEffErr);
1205  printf("Matched Apt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1206 
1207  // Matched all pt efficiency using tracks matched with reconstructible ones triggerable Apt
1208  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1209  FillEffHistos(efficiency, "matchedTracksMatchMCApt", effAllPt);
1210  GetEfficiency(efficiency, totalEff, totalEffErr);
1211  printf("Matched Apt efficiency using reconstructed tracks matching MC-Apt = %f +- %f\n", totalEff, totalEffErr);
1212 
1213  // Matched all pt efficiency using tracks matched with reconstructible ones triggerable other pt
1214  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1215  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effAllPt);
1216  GetEfficiency(efficiency, totalEff, totalEffErr);
1217  printf("Matched Apt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1218 
1219  // Matched all pt efficiency using tracks matched with reconstructible ones not triggerable
1220  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1221  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effAllPt);
1222  GetEfficiency(efficiency, totalEff, totalEffErr);
1223  printf("Matched Apt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1224 
1225  // Matched all pt efficiency using tracks matched with same reconstructible & triggerable ones (all pt MC trig)
1226  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1227  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1228  FillEffHistos(efficiency, "matchedTracksMatchSameMCApt", effAllPt);
1229 
1230  // Matched all pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1231  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1232  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effAllPt);
1233 
1234  // Matched all pt efficiency using tracks matched with different reconstructible & triggerable ones (all pt MC trig)
1235  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1236  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1237  FillEffHistos(efficiency, "matchedTracksMatchDiffMCApt", effAllPt);
1238 
1239  // Matched all pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1240  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1241  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effAllPt);
1242 
1243  // Matched efficiency using tracks matched with reconstructible ones only (all pt MC trig)
1244  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1245  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1246  FillEffHistos(efficiency, "matchedTracksMatchTrkMCApt", effAllPt);
1247 
1248  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1249  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kOtherTrig);
1250  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effAllPt);
1251 
1252  // ------ Tracker matched with low pt trigger ------
1253 
1254  // Matched low pt efficiency using all reconstructed tracks
1255  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1256  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1257  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1258  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1259  FillEffHistos(efficiency, "matchedTracks", effLowPt);
1260  GetEfficiency(efficiency, totalEff, totalEffErr);
1261  printf("Matched Lpt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1262 
1263  // Matched low pt efficiency using tracks matched with reconstructible ones, triggerable or not
1264  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1265  FillEffHistos(efficiency, "matchedTracksMatchMC", effLowPt);
1266  GetEfficiency(efficiency, totalEff, totalEffErr);
1267  printf("Matched Lpt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1268 
1269  // Matched low pt efficiency using tracks matched with reconstructible ones triggerable Lpt
1270  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1271  FillEffHistos(efficiency, "matchedTracksMatchMCLpt", effLowPt);
1272  GetEfficiency(efficiency, totalEff, totalEffErr);
1273  printf("Matched Lpt efficiency using reconstructed tracks matching MC-Lpt = %f +- %f\n", totalEff, totalEffErr);
1274 
1275  // Matched low pt efficiency using tracks matched with reconstructible ones triggerable other pt
1276  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1277  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effLowPt);
1278  GetEfficiency(efficiency, totalEff, totalEffErr);
1279  printf("Matched Lpt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1280 
1281  // Matched low pt efficiency using tracks matched with reconstructible ones not triggerable
1282  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1283  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effLowPt);
1284  GetEfficiency(efficiency, totalEff, totalEffErr);
1285  printf("Matched Lpt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1286 
1287  // Matched low pt efficiency using tracks matched with same reconstructible & triggerable ones (low pt MC trig)
1288  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1289  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1290  FillEffHistos(efficiency, "matchedTracksMatchSameMCLpt", effLowPt);
1291 
1292  // Matched low pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1293  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1294  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effLowPt);
1295 
1296  // Matched low pt efficiency using tracks matched with different reconstructible & triggerable ones (low pt MC trig)
1297  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1298  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1299  FillEffHistos(efficiency, "matchedTracksMatchDiffMCLpt", effLowPt);
1300 
1301  // Matched low pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1302  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1303  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effLowPt);
1304 
1305  // Matched efficiency using tracks matched with reconstructible ones only (low pt MC trig)
1306  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1307  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1308  FillEffHistos(efficiency, "matchedTracksMatchTrkMCLpt", effLowPt);
1309 
1310  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1311  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kAllPtTrig);
1312  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effLowPt);
1313 
1314  // ------ Tracker matched with high pt trigger ------
1315 
1316  // Matched high pt efficiency using all reconstructed tracks
1317  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1318  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1319  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1320  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1321  FillEffHistos(efficiency, "matchedTracks", effHighPt);
1322  GetEfficiency(efficiency, totalEff, totalEffErr);
1323  printf("Matched Hpt efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1324 
1325  // Matched high pt efficiency using tracks matched with reconstructible ones, triggerable or not
1326  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1327  FillEffHistos(efficiency, "matchedTracksMatchMC", effHighPt);
1328  GetEfficiency(efficiency, totalEff, totalEffErr);
1329  printf("Matched Hpt efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1330 
1331  // Matched high pt efficiency using tracks matched with reconstructible ones triggerable Hpt
1332  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1333  FillEffHistos(efficiency, "matchedTracksMatchMCHpt", effHighPt);
1334  GetEfficiency(efficiency, totalEff, totalEffErr);
1335  printf("Matched Hpt efficiency using reconstructed tracks matching MC-Hpt = %f +- %f\n", totalEff, totalEffErr);
1336 
1337  // Matched high pt efficiency using tracks matched with reconstructible ones triggerable other pt
1338  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1339  FillEffHistos(efficiency, "matchedTracksMatchMCOther", effHighPt);
1340  GetEfficiency(efficiency, totalEff, totalEffErr);
1341  printf("Matched Hpt efficiency using reconstructed tracks matching MC-other = %f +- %f\n", totalEff, totalEffErr);
1342 
1343  // Matched high pt efficiency using tracks matched with reconstructible ones not triggerable
1344  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1345  FillEffHistos(efficiency, "matchedTracksMatchMCNoTrig", effHighPt);
1346  GetEfficiency(efficiency, totalEff, totalEffErr);
1347  printf("Matched Hpt efficiency using reconstructed tracks matching MC-noTrig = %f +- %f\n", totalEff, totalEffErr);
1348 
1349  // Matched high pt efficiency using tracks matched with same reconstructible & triggerable ones (high pt MC trig)
1350  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kMatchedSame);
1351  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1352  FillEffHistos(efficiency, "matchedTracksMatchSameMCHpt", effHighPt);
1353 
1354  // Matched high pt efficiency using tracks matched with same reconstructible & triggerable ones (other MC trig)
1355  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1356  FillEffHistos(efficiency, "matchedTracksMatchSameMCOther", effHighPt);
1357 
1358  // Matched high pt efficiency using tracks matched with different reconstructible & triggerable ones (high pt MC trig)
1359  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedDiff, kMatchedDiff);
1360  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1361  FillEffHistos(efficiency, "matchedTracksMatchDiffMCHpt", effHighPt);
1362 
1363  // Matched high pt efficiency using tracks matched with different reconstructible & triggerable ones (other MC trig)
1364  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1365  FillEffHistos(efficiency, "matchedTracksMatchDiffMCOther", effHighPt);
1366 
1367  // Matched efficiency using tracks matched with reconstructible ones only (high pt MC trig)
1368  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kTrackerOnly);
1369  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1370  FillEffHistos(efficiency, "matchedTracksMatchTrkMCHpt", effHighPt);
1371 
1372  // Matched efficiency using tracks matched with reconstructible ones only (other MC trig)
1373  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kLowPtTrig);
1374  FillEffHistos(efficiency, "matchedTracksMatchTrkMCOther", effHighPt);
1375 
1376  // ------ Trigger only ------
1377 
1378  // Trigger efficiency using all reconstructed tracks
1379  efficiency->GetNum()->GetAxis(kVarHasTracker)->SetRange();
1380  efficiency->GetDen()->GetAxis(kVarHasTracker)->SetRange();
1381  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1382  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kHighPtTrig);
1383  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1384  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1385  efficiency->GetNum()->SetRangeUser(kVarDupliTrg, 0., 0.);
1386  FillEffHistos(efficiency, "triggerTracks", effAnyPt);
1387  GetEfficiency(efficiency, totalEff, totalEffErr);
1388  effSummary->Fill("Trigger_all", totalEff);
1389  effSummary->SetBinError(++sumEffBin, totalEffErr);
1390  printf("Trigger efficiency using all reconstructed tracks = %f +- %f\n", totalEff, totalEffErr);
1391 
1392  // Trigger efficiency using tracks matched with triggerable ones
1393  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1394  FillEffHistos(efficiency, "triggerTracksMatchMC", effAnyPt);
1395  GetEfficiency(efficiency, totalEff, totalEffErr);
1396  effSummary->Fill("Trigger_MCId", totalEff);
1397  effSummary->SetBinError(++sumEffBin, totalEffErr);
1398  printf("Trigger efficiency using reconstructed tracks matching MC = %f +- %f\n", totalEff, totalEffErr);
1399 
1400  // ------ All pt trigger only ------
1401 
1402  // All pt trigger efficiency using all reconstructed tracks
1403  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kHighPtTrig);
1404  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1405  FillEffHistos(efficiency, "triggerTracks", effAllPt);
1406 
1407  // All pt trigger efficiency using tracks matched with all pt triggerable ones
1408  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1409  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kHighPtTrig);
1410  FillEffHistos(efficiency, "triggerTracksMatchMCApt", effAllPt);
1411 
1412  // All pt trigger efficiency using tracks matched with other triggerable ones
1413  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1414  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effAllPt);
1415 
1416  // ------ Low pt trigger only ------
1417 
1418  // Low pt trigger efficiency using all reconstructed tracks
1419  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1420  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kHighPtTrig);
1421  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1422  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1423  FillEffHistos(efficiency, "triggerTracks", effLowPt);
1424 
1425  // Low pt trigger efficiency using tracks matched with Low pt triggerable ones
1426  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1427  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kHighPtTrig);
1428  FillEffHistos(efficiency, "triggerTracksMatchMCLpt", effLowPt);
1429 
1430  // Low pt trigger efficiency using tracks matched with other triggerable ones
1431  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kAllPtTrig);
1432  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effLowPt);
1433 
1434  // ------ High pt trigger only ------
1435 
1436  // High pt trigger efficiency using all reconstructed tracks
1437  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1438  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1439  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1440  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1441  FillEffHistos(efficiency, "triggerTracks", effHighPt);
1442 
1443  // High pt trigger efficiency using tracks matched with High pt triggerable ones
1444  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kMatchedSame, kTriggerOnly);
1445  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1446  FillEffHistos(efficiency, "triggerTracksMatchMCHpt", effHighPt);
1447 
1448  // All pt trigger efficiency using tracks matched with other triggerable ones
1449  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kLowPtTrig);
1450  FillEffHistos(efficiency, "triggerTracksMatchMCOther", effHighPt);
1451 
1452  // ------ Tracker reconstructible not triggerable ------
1453 
1454  // all tracker tracks
1455  efficiency->GetNum()->SetRangeUser(kVarHasTracker, 1., 1.);
1456  efficiency->GetDen()->SetRangeUser(kVarHasTracker, 1., 1.);
1457  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1458  efficiency->GetDen()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1459  efficiency->GetNum()->SetRangeUser(kVarMatchMC, kTrackerOnly, kMatchedDiff);
1460  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kNoMatchTrig, kNoMatchTrig);
1461  efficiency->GetNum()->GetAxis(kVarDupliTrg)->SetRange();
1462  FillEffHistos(efficiency, "allTracks", notTrgable);
1463 
1464  // tracker not matched
1465  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1466  FillEffHistos(efficiency, "notMatched", notTrgable);
1467 
1468  // tracker matched with all pt trigger
1469  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1470  FillEffHistos(efficiency, "MatchedApt", notTrgable);
1471 
1472  // tracker matched with low pt trigger
1473  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1474  FillEffHistos(efficiency, "MatchedLpt", notTrgable);
1475 
1476  // tracker matched with high pt trigger
1477  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1478  FillEffHistos(efficiency, "MatchedHpt", notTrgable);
1479 
1480  // ------ Tracker reconstructible triggerable no pt ------
1481 
1482  // all tracker tracks
1483  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1484  efficiency->GetDen()->SetRangeUser(kVarTrigger, kOtherTrig, kOtherTrig);
1485  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kOtherTrig, kOtherTrig);
1486  FillEffHistos(efficiency, "allTracks", trgableNoPtOnly);
1487 
1488  // tracker not matched
1489  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1490  FillEffHistos(efficiency, "notMatched", trgableNoPtOnly);
1491 
1492  // tracker matched with all pt trigger
1493  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1494  FillEffHistos(efficiency, "MatchedApt", trgableNoPtOnly);
1495 
1496  // tracker matched with low pt trigger
1497  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1498  FillEffHistos(efficiency, "MatchedLpt", trgableNoPtOnly);
1499 
1500  // tracker matched with high pt trigger
1501  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1502  FillEffHistos(efficiency, "MatchedHpt", trgableNoPtOnly);
1503 
1504  // ------ Tracker reconstructible triggerable Apt ------
1505 
1506  // all tracker tracks
1507  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1508  efficiency->GetDen()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1509  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kAllPtTrig, kAllPtTrig);
1510  FillEffHistos(efficiency, "allTracks", trgableAPtOnly);
1511 
1512  // tracker not matched
1513  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1514  FillEffHistos(efficiency, "notMatched", trgableAPtOnly);
1515 
1516  // tracker matched with all pt trigger
1517  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1518  FillEffHistos(efficiency, "MatchedApt", trgableAPtOnly);
1519 
1520  // tracker matched with low pt trigger
1521  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1522  FillEffHistos(efficiency, "MatchedLpt", trgableAPtOnly);
1523 
1524  // tracker matched with high pt trigger
1525  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1526  FillEffHistos(efficiency, "MatchedHpt", trgableAPtOnly);
1527 
1528  // ------ Tracker reconstructible triggerable Lpt ------
1529 
1530  // all tracker tracks
1531  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1532  efficiency->GetDen()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1533  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kLowPtTrig, kLowPtTrig);
1534  FillEffHistos(efficiency, "allTracks", trgableLPtOnly);
1535 
1536  // tracker not matched
1537  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1538  FillEffHistos(efficiency, "notMatched", trgableLPtOnly);
1539 
1540  // tracker matched with all pt trigger
1541  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1542  FillEffHistos(efficiency, "MatchedApt", trgableLPtOnly);
1543 
1544  // tracker matched with low pt trigger
1545  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1546  FillEffHistos(efficiency, "MatchedLpt", trgableLPtOnly);
1547 
1548  // tracker matched with high pt trigger
1549  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1550  FillEffHistos(efficiency, "MatchedHpt", trgableLPtOnly);
1551 
1552  // ------ Tracker reconstructible triggerable Hpt ------
1553 
1554  // all tracker tracks
1555  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1556  efficiency->GetDen()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1557  efficiency->GetNum()->SetRangeUser(kVarMCTrigger, kHighPtTrig, kHighPtTrig);
1558  FillEffHistos(efficiency, "allTracks", trgableHPtOnly);
1559 
1560  // tracker not matched
1561  efficiency->GetNum()->SetRangeUser(kVarTrigger, kNoMatchTrig, kNoMatchTrig);
1562  FillEffHistos(efficiency, "notMatched", trgableHPtOnly);
1563 
1564  // tracker matched with all pt trigger
1565  efficiency->GetNum()->SetRangeUser(kVarTrigger, kAllPtTrig, kAllPtTrig);
1566  FillEffHistos(efficiency, "MatchedApt", trgableHPtOnly);
1567 
1568  // tracker matched with low pt trigger
1569  efficiency->GetNum()->SetRangeUser(kVarTrigger, kLowPtTrig, kLowPtTrig);
1570  FillEffHistos(efficiency, "MatchedLpt", trgableHPtOnly);
1571 
1572  // tracker matched with high pt trigger
1573  efficiency->GetNum()->SetRangeUser(kVarTrigger, kHighPtTrig, kHighPtTrig);
1574  FillEffHistos(efficiency, "MatchedHpt", trgableHPtOnly);
1575 
1576  // ------ reset ranges before saving CF containers ------
1577  efficiency->GetNum()->GetAxis(kVarHasTracker)->SetRange();
1578  efficiency->GetDen()->GetAxis(kVarHasTracker)->SetRange();
1579  efficiency->GetNum()->GetAxis(kVarTrigger)->SetRange();
1580  efficiency->GetDen()->GetAxis(kVarTrigger)->SetRange();
1581  efficiency->GetNum()->GetAxis(kVarMatchMC)->SetRange();
1582  efficiency->GetNum()->GetAxis(kVarMCTrigger)->SetRange();
1583  efficiency->GetNum()->GetAxis(kVarDupliTrg)->SetRange();
1584 
1585  // ------ Plot summary ------
1586 
1587  // plot histogram summarizing global efficiencies
1588  TCanvas* cEffSummary = new TCanvas("cEffSummary","Efficiency summary",20,20,310,310);
1589  cEffSummary->SetFillColor(10); cEffSummary->SetHighLightColor(10);
1590  cEffSummary->SetLeftMargin(0.15); cEffSummary->SetBottomMargin(0.15);
1591  effSummary->DrawCopy("etext");
1592 
1593  // --- plot trigger resolution ---
1594  TCanvas* cTriggerResolution = new TCanvas("cTriggerResolution","Trigger resolution",10,10,310,310);
1595  cTriggerResolution->SetFillColor(10); cTriggerResolution->SetHighLightColor(10);
1596  cTriggerResolution->SetLeftMargin(0.15); cTriggerResolution->SetBottomMargin(0.15);
1597  cTriggerResolution->Divide(2,2);
1598  cTriggerResolution->cd(1);
1599  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigX11))->DrawCopy();
1600  cTriggerResolution->cd(2);
1601  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigY11))->DrawCopy();
1602  cTriggerResolution->cd(3);
1603  static_cast<TH1*>(fTriggerList->UncheckedAt(kResTrigSlopeY))->DrawCopy();
1604 
1605  // --- compute momentum resolution at vertex ---
1606  fPAtVtxList = new TObjArray(100);
1607  fPAtVtxList->SetOwner();
1608 
1609  // define graphs
1610  TGraphAsymmErrors* gMeanResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1611  gMeanResPAtVtxVsP->SetName("gMeanResPAtVtxVsP");
1612  gMeanResPAtVtxVsP->SetTitle("<#Delta_{p}> at vertex versus p;p (GeV/c);<#Delta_{p}> (GeV/c)");
1613  fPAtVtxList->AddAt(gMeanResPAtVtxVsP, kMeanResPAtVtxVsP);
1614  TGraphAsymmErrors* gMostProbResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1615  gMostProbResPAtVtxVsP->SetName("gMostProbResPAtVtxVsP");
1616  gMostProbResPAtVtxVsP->SetTitle("Most probable #Delta_{p} at vertex versus p;p (GeV/c);Most prob. #Delta_{p} (GeV/c)");
1617  fPAtVtxList->AddAt(gMostProbResPAtVtxVsP, kMostProbResPAtVtxVsP);
1618  TGraphAsymmErrors* gSigmaResPAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1619  gSigmaResPAtVtxVsP->SetName("gSigmaResPAtVtxVsP");
1620  gSigmaResPAtVtxVsP->SetTitle("#sigma_{p}/p at vertex versus p;p (GeV/c);#sigma_{p}/p (%)");
1621  fPAtVtxList->AddAt(gSigmaResPAtVtxVsP, kSigmaResPAtVtxVsP);
1622 
1623  // fit histo and fill graphs
1624  TH2* h = static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsP));
1625  FitLandauGausResVsP(h, "momentum residuals at vertex", gMeanResPAtVtxVsP, gMostProbResPAtVtxVsP, gSigmaResPAtVtxVsP);
1626 
1627  // convert resolution into relative resolution
1628  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
1629  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
1630  Double_t x,y;
1631  gSigmaResPAtVtxVsP->GetPoint(i/rebinFactorX-1, x, y);
1632  gSigmaResPAtVtxVsP->SetPoint(i/rebinFactorX-1, x, 100.*y/x);
1633  gSigmaResPAtVtxVsP->SetPointEYlow(i/rebinFactorX-1, 100.*gSigmaResPAtVtxVsP->GetErrorYlow(i/rebinFactorX-1)/x);
1634  gSigmaResPAtVtxVsP->SetPointEYhigh(i/rebinFactorX-1, 100.*gSigmaResPAtVtxVsP->GetErrorYhigh(i/rebinFactorX-1)/x);
1635  }
1636 
1637  // --- compute momentum resolution at first cluster ---
1638  fPAt1stClList = new TObjArray(100);
1639  fPAt1stClList->SetOwner();
1640 
1641  // define graphs
1642  TGraphAsymmErrors* gMeanResPAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1643  gMeanResPAt1stClVsP->SetName("gMeanResPAt1stClVsP");
1644  gMeanResPAt1stClVsP->SetTitle("<#Delta_{p}> at first cluster versus p;p (GeV/c);<#Delta_{p}> (GeV/c)");
1645  fPAt1stClList->AddAt(gMeanResPAt1stClVsP, kMeanResPAt1stClVsP);
1646  TGraphAsymmErrors* gSigmaResPAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1647  gSigmaResPAt1stClVsP->SetName("gSigmaResPAt1stClVsP");
1648  gSigmaResPAt1stClVsP->SetTitle("#sigma_{p}/p at first cluster versus p;p (GeV/c);#sigma_{p}/p (%)");
1649  fPAt1stClList->AddAt(gSigmaResPAt1stClVsP, kSigmaResPAt1stClVsP);
1650 
1651  // fit histo and fill graphs
1652  h = static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAt1stClVsP));
1653  FitGausResVsMom(h, 0., 1., "momentum residuals at first cluster", gMeanResPAt1stClVsP, gSigmaResPAt1stClVsP);
1654 
1655  // convert resolution into relative resolution
1656  rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
1657  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
1658  Double_t x,y;
1659  gSigmaResPAt1stClVsP->GetPoint(i/rebinFactorX-1, x, y);
1660  gSigmaResPAt1stClVsP->SetPoint(i/rebinFactorX-1, x, 100.*y/x);
1661  gSigmaResPAt1stClVsP->SetPointEYlow(i/rebinFactorX-1, 100.*gSigmaResPAt1stClVsP->GetErrorYlow(i/rebinFactorX-1)/x);
1662  gSigmaResPAt1stClVsP->SetPointEYhigh(i/rebinFactorX-1, 100.*gSigmaResPAt1stClVsP->GetErrorYhigh(i/rebinFactorX-1)/x);
1663  }
1664 
1665  // --- compute slope resolution at vertex ---
1666  fSlopeAtVtxList = new TObjArray(100);
1667  fSlopeAtVtxList->SetOwner();
1668 
1669  // define graphs
1670  TGraphAsymmErrors* gMeanResSlopeXAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1671  gMeanResSlopeXAtVtxVsP->SetName("gMeanResSlopeXAtVtxVsP");
1672  gMeanResSlopeXAtVtxVsP->SetTitle("<#Delta_{slope_{X}}> at vertex versus p;p (GeV/c);<#Delta_{slope_{X}}>");
1673  fSlopeAtVtxList->AddAt(gMeanResSlopeXAtVtxVsP, kMeanResSlopeXAtVtxVsP);
1674  TGraphAsymmErrors* gMeanResSlopeYAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1675  gMeanResSlopeYAtVtxVsP->SetName("gMeanResSlopeYAtVtxVsP");
1676  gMeanResSlopeYAtVtxVsP->SetTitle("<#Delta_{slope_{Y}}> at vertex versus p;p (GeV/c);<#Delta_{slope_{Y}}>");
1677  fSlopeAtVtxList->AddAt(gMeanResSlopeYAtVtxVsP, kMeanResSlopeYAtVtxVsP);
1678  TGraphAsymmErrors* gSigmaResSlopeXAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1679  gSigmaResSlopeXAtVtxVsP->SetName("gSigmaResSlopeXAtVtxVsP");
1680  gSigmaResSlopeXAtVtxVsP->SetTitle("#sigma_{slope_{X}} at vertex versus p;p (GeV/c);#sigma_{slope_{X}}");
1681  fSlopeAtVtxList->AddAt(gSigmaResSlopeXAtVtxVsP, kSigmaResSlopeXAtVtxVsP);
1682  TGraphAsymmErrors* gSigmaResSlopeYAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1683  gSigmaResSlopeYAtVtxVsP->SetName("gSigmaResSlopeYAtVtxVsP");
1684  gSigmaResSlopeYAtVtxVsP->SetTitle("#sigma_{slope_{Y}} at vertex versus p;p (GeV/c);#sigma_{slope_{Y}}");
1685  fSlopeAtVtxList->AddAt(gSigmaResSlopeYAtVtxVsP, kSigmaResSlopeYAtVtxVsP);
1686 
1687  // fit histo and fill graphs
1688  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsP)), 0., 2.e-3,
1689  "slopeX residuals at vertex", gMeanResSlopeXAtVtxVsP, gSigmaResSlopeXAtVtxVsP);
1690  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsP)), 0., 2.e-3,
1691  "slopeY residuals at vertex", gMeanResSlopeYAtVtxVsP, gSigmaResSlopeYAtVtxVsP);
1692 
1693  // --- compute slope resolution at first cluster ---
1694  fSlopeAt1stClList = new TObjArray(100);
1695  fSlopeAt1stClList->SetOwner();
1696 
1697  // define graphs
1698  TGraphAsymmErrors* gMeanResSlopeXAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1699  gMeanResSlopeXAt1stClVsP->SetName("gMeanResSlopeXAt1stClVsP");
1700  gMeanResSlopeXAt1stClVsP->SetTitle("<#Delta_{slope_{X}}> at first cluster versus p;p (GeV/c);<#Delta_{slope_{X}}>");
1701  fSlopeAt1stClList->AddAt(gMeanResSlopeXAt1stClVsP, kMeanResSlopeXAt1stClVsP);
1702  TGraphAsymmErrors* gMeanResSlopeYAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1703  gMeanResSlopeYAt1stClVsP->SetName("gMeanResSlopeYAt1stClVsP");
1704  gMeanResSlopeYAt1stClVsP->SetTitle("<#Delta_{slope_{Y}}> at first cluster versus p;p (GeV/c);<#Delta_{slope_{Y}}>");
1705  fSlopeAt1stClList->AddAt(gMeanResSlopeYAt1stClVsP, kMeanResSlopeYAt1stClVsP);
1706  TGraphAsymmErrors* gSigmaResSlopeXAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1707  gSigmaResSlopeXAt1stClVsP->SetName("gSigmaResSlopeXAt1stClVsP");
1708  gSigmaResSlopeXAt1stClVsP->SetTitle("#sigma_{slope_{X}} at first cluster versus p;p (GeV/c);#sigma_{slope_{X}}");
1709  fSlopeAt1stClList->AddAt(gSigmaResSlopeXAt1stClVsP, kSigmaResSlopeXAt1stClVsP);
1710  TGraphAsymmErrors* gSigmaResSlopeYAt1stClVsP = new TGraphAsymmErrors(fNPBins);
1711  gSigmaResSlopeYAt1stClVsP->SetName("gSigmaResSlopeYAt1stClVsP");
1712  gSigmaResSlopeYAt1stClVsP->SetTitle("#sigma_{slope_{Y}} at first cluster versus p;p (GeV/c);#sigma_{slope_{Y}}");
1713  fSlopeAt1stClList->AddAt(gSigmaResSlopeYAt1stClVsP, kSigmaResSlopeYAt1stClVsP);
1714 
1715  // fit histo and fill graphs
1716  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAt1stClVsP)), 0., 3.e-4,
1717  "slopeX residuals at first cluster", gMeanResSlopeXAt1stClVsP, gSigmaResSlopeXAt1stClVsP);
1718  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAt1stClVsP)), 0., 3.e-4,
1719  "slopeY residuals at first cluster", gMeanResSlopeYAt1stClVsP, gSigmaResSlopeYAt1stClVsP);
1720 
1721  // --- compute eta resolution at vertex ---
1722  fEtaAtVtxList = new TObjArray(100);
1723  fEtaAtVtxList->SetOwner();
1724 
1725  // define graphs
1726  TGraphAsymmErrors* gMeanResEtaAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1727  gMeanResEtaAtVtxVsP->SetName("gMeanResEtaAtVtxVsP");
1728  gMeanResEtaAtVtxVsP->SetTitle("<#Delta_{eta}> at vertex versus p;p (GeV/c);<#Delta_{eta}>");
1729  fEtaAtVtxList->AddAt(gMeanResEtaAtVtxVsP, kMeanResEtaAtVtxVsP);
1730  TGraphAsymmErrors* gSigmaResEtaAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1731  gSigmaResEtaAtVtxVsP->SetName("gSigmaResEtaAtVtxVsP");
1732  gSigmaResEtaAtVtxVsP->SetTitle("#sigma_{eta} at vertex versus p;p (GeV/c);#sigma_{eta}");
1733  fEtaAtVtxList->AddAt(gSigmaResEtaAtVtxVsP, kSigmaResEtaAtVtxVsP);
1734 
1735  // fit histo and fill graphs
1736  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsP)), 0., 0.1,
1737  "eta residuals at vertex", gMeanResEtaAtVtxVsP, gSigmaResEtaAtVtxVsP);
1738 
1739  // --- compute phi resolution at vertex ---
1740  fPhiAtVtxList = new TObjArray(100);
1741  fPhiAtVtxList->SetOwner();
1742 
1743  // define graphs
1744  TGraphAsymmErrors* gMeanResPhiAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1745  gMeanResPhiAtVtxVsP->SetName("gMeanResPhiAtVtxVsP");
1746  gMeanResPhiAtVtxVsP->SetTitle("<#Delta_{phi}> at vertex versus p;p (GeV/c);<#Delta_{phi}>");
1747  fPhiAtVtxList->AddAt(gMeanResPhiAtVtxVsP, kMeanResPhiAtVtxVsP);
1748  TGraphAsymmErrors* gSigmaResPhiAtVtxVsP = new TGraphAsymmErrors(fNPBins);
1749  gSigmaResPhiAtVtxVsP->SetName("gSigmaResPhiAtVtxVsP");
1750  gSigmaResPhiAtVtxVsP->SetTitle("#sigma_{phi} at vertex versus p;p (GeV/c);#sigma_{phi}");
1751  fPhiAtVtxList->AddAt(gSigmaResPhiAtVtxVsP, kSigmaResPhiAtVtxVsP);
1752 
1753  // fit histo and fill graphs
1754  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsP)), 0., 0.01,
1755  "phi residuals at vertex", gMeanResPhiAtVtxVsP, gSigmaResPhiAtVtxVsP);
1756 
1757  // --- compute DCA resolution and MCS dispersion ---
1758  fDCAList = new TObjArray(100);
1759  fDCAList->SetOwner();
1760 
1761  // define graphs
1762  TGraphAsymmErrors* gMeanPDCAVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1763  gMeanPDCAVsPIn23deg->SetName("gMeanPDCAVsPIn23deg");
1764  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)");
1765  fDCAList->AddAt(gMeanPDCAVsPIn23deg, kMeanPDCAVsPIn23deg);
1766  TGraphAsymmErrors* gSigmaPDCAVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1767  gSigmaPDCAVsPIn23deg->SetName("gSigmaPDCAVsPIn23deg");
1768  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)");
1769  fDCAList->AddAt(gSigmaPDCAVsPIn23deg, kSigmaPDCAVsPIn23deg);
1770  TGraphAsymmErrors* gMeanPDCAVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1771  gMeanPDCAVsPIn310deg->SetName("gMeanPDCAVsPIn310deg");
1772  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)");
1773  fDCAList->AddAt(gMeanPDCAVsPIn310deg, kMeanPDCAVsPIn310deg);
1774  TGraphAsymmErrors* gSigmaPDCAVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1775  gSigmaPDCAVsPIn310deg->SetName("gSigmaPDCAVsPIn310deg");
1776  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)");
1777  fDCAList->AddAt(gSigmaPDCAVsPIn310deg, kSigmaPDCAVsPIn310deg);
1778 
1779  TGraphAsymmErrors* gMeanPMCSAngVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1780  gMeanPMCSAngVsPIn23deg->SetName("gMeanPMCSAngVsPIn23deg");
1781  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)");
1782  fDCAList->AddAt(gMeanPMCSAngVsPIn23deg, kMeanPMCSAngVsPIn23deg);
1783  TGraphAsymmErrors* gSigmaPMCSAngVsPIn23deg = new TGraphAsymmErrors(fNPBins);
1784  gSigmaPMCSAngVsPIn23deg->SetName("gSigmaPMCSAngVsPIn23deg");
1785  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)");
1786  fDCAList->AddAt(gSigmaPMCSAngVsPIn23deg, kSigmaPMCSAngVsPIn23deg);
1787  TGraphAsymmErrors* gMeanPMCSAngVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1788  gMeanPMCSAngVsPIn310deg->SetName("gMeanPMCSAngVsPIn310deg");
1789  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)");
1790  fDCAList->AddAt(gMeanPMCSAngVsPIn310deg, kMeanPMCSAngVsPIn310deg);
1791  TGraphAsymmErrors* gSigmaPMCSAngVsPIn310deg = new TGraphAsymmErrors(fNPBins);
1792  gSigmaPMCSAngVsPIn310deg->SetName("gSigmaPMCSAngVsPIn310deg");
1793  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)");
1794  fDCAList->AddAt(gSigmaPMCSAngVsPIn310deg, kSigmaPMCSAngVsPIn310deg);
1795 
1796  // fit histo and fill graphs
1797  FitPDCAVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPIn23deg)),
1798  "p*DCA (tracks in [2,3] deg.)", gMeanPDCAVsPIn23deg, gSigmaPDCAVsPIn23deg);
1799  FitPDCAVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPIn310deg)),
1800  "p*DCA (tracks in [3,10] deg.)", gMeanPDCAVsPIn310deg, gSigmaPDCAVsPIn310deg);
1801  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPMCSAngVsPIn23deg)), 0., 2.e-3,
1802  "p*MCSAngle (tracks in [2,3] deg.)", gMeanPMCSAngVsPIn23deg, gSigmaPMCSAngVsPIn23deg);
1803  FitGausResVsMom(static_cast<TH2*>(fTrackerList->UncheckedAt(kPMCSAngVsPIn310deg)), 0., 2.e-3,
1804  "p*MCSAngle (tracks in [3,10] deg.)", gMeanPMCSAngVsPIn310deg, gSigmaPMCSAngVsPIn310deg);
1805 
1806  // --- compute cluster resolution ---
1807  fClusterList = new TObjArray(100);
1808  fClusterList->SetOwner();
1809 
1810  // define graphs per chamber
1811  TGraphErrors* gMeanResClXVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1812  gMeanResClXVsCh->SetName("gMeanResClXVsCh");
1813  gMeanResClXVsCh->SetTitle("cluster-trackRef residual-X per Ch: mean;chamber ID;<#Delta_{X}> (cm)");
1814  gMeanResClXVsCh->SetMarkerStyle(kFullDotLarge);
1815  fClusterList->AddAt(gMeanResClXVsCh, kMeanResClXVsCh);
1816  TGraphErrors* gMeanResClYVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1817  gMeanResClYVsCh->SetName("gMeanResClYVsCh");
1818  gMeanResClYVsCh->SetTitle("cluster-trackRef residual-Y per Ch: mean;chamber ID;<#Delta_{Y}> (cm)");
1819  gMeanResClYVsCh->SetMarkerStyle(kFullDotLarge);
1820  fClusterList->AddAt(gMeanResClYVsCh, kMeanResClYVsCh);
1821  TGraphErrors* gSigmaResClXVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1822  gSigmaResClXVsCh->SetName("gSigmaResClXVsCh");
1823  gSigmaResClXVsCh->SetTitle("cluster-trackRef residual-X per Ch: sigma;chamber ID;#sigma_{X} (cm)");
1824  gSigmaResClXVsCh->SetMarkerStyle(kFullDotLarge);
1825  fClusterList->AddAt(gSigmaResClXVsCh, kSigmaResClXVsCh);
1826  TGraphErrors* gSigmaResClYVsCh = new TGraphErrors(AliMUONConstants::NTrackingCh());
1827  gSigmaResClYVsCh->SetName("gSigmaResClYVsCh");
1828  gSigmaResClYVsCh->SetTitle("cluster-trackRef residual-Y per Ch: sigma;chamber ID;#sigma_{Y} (cm)");
1829  gSigmaResClYVsCh->SetMarkerStyle(kFullDotLarge);
1830  fClusterList->AddAt(gSigmaResClYVsCh, kSigmaResClYVsCh);
1831 
1832  // define graphs per DE
1833  TGraphErrors* gMeanResClXVsDE = new TGraphErrors(fNDE);
1834  gMeanResClXVsDE->SetName("gMeanResClXVsDE");
1835  gMeanResClXVsDE->SetTitle("cluster-trackRef residual-X per DE: mean;DE ID;<#Delta_{X}> (cm)");
1836  gMeanResClXVsDE->SetMarkerStyle(kFullDotLarge);
1837  fClusterList->AddAt(gMeanResClXVsDE, kMeanResClXVsDE);
1838  TGraphErrors* gMeanResClYVsDE = new TGraphErrors(fNDE);
1839  gMeanResClYVsDE->SetName("gMeanResClYVsDE");
1840  gMeanResClYVsDE->SetTitle("cluster-trackRef residual-Y per dE: mean;DE ID;<#Delta_{Y}> (cm)");
1841  gMeanResClYVsDE->SetMarkerStyle(kFullDotLarge);
1842  fClusterList->AddAt(gMeanResClYVsDE, kMeanResClYVsDE);
1843  TGraphErrors* gSigmaResClXVsDE = new TGraphErrors(fNDE);
1844  gSigmaResClXVsDE->SetName("gSigmaResClXVsDE");
1845  gSigmaResClXVsDE->SetTitle("cluster-trackRef residual-X per DE: sigma;DE ID;#sigma_{X} (cm)");
1846  gSigmaResClXVsDE->SetMarkerStyle(kFullDotLarge);
1847  fClusterList->AddAt(gSigmaResClXVsDE, kSigmaResClXVsDE);
1848  TGraphErrors* gSigmaResClYVsDE = new TGraphErrors(fNDE);
1849  gSigmaResClYVsDE->SetName("gSigmaResClYVsDE");
1850  gSigmaResClYVsDE->SetTitle("cluster-trackRef residual-Y per DE: sigma;DE ID;#sigma_{Y} (cm)");
1851  gSigmaResClYVsDE->SetMarkerStyle(kFullDotLarge);
1852  fClusterList->AddAt(gSigmaResClYVsDE, kSigmaResClYVsDE);
1853 
1854  // fit histo and fill graphs per chamber
1855  Double_t clusterResPerCh[10][2];
1856  for (Int_t i = 0; i < AliMUONConstants::NTrackingCh(); i++) {
1857  TH1D *tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsCh))->ProjectionY("tmp",i+1,i+1,"e");
1858  FitClusterResidual(tmp, i, clusterResPerCh[i][0], gMeanResClXVsCh, gSigmaResClXVsCh);
1859  delete tmp;
1860  tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClYVsCh))->ProjectionY("tmp",i+1,i+1,"e");
1861  FitClusterResidual(tmp, i, clusterResPerCh[i][1], gMeanResClYVsCh, gSigmaResClYVsCh);
1862  delete tmp;
1863  }
1864 
1865  // fit histo and fill graphs per DE
1866  Double_t clusterResPerDE[200][2];
1867  for (Int_t i = 0; i < fNDE; i++) {
1868  TH1D *tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsDE))->ProjectionY("tmp",i+1,i+1,"e");
1869  FitClusterResidual(tmp, i, clusterResPerDE[i][0], gMeanResClXVsDE, gSigmaResClXVsDE);
1870  delete tmp;
1871  tmp = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClYVsDE))->ProjectionY("tmp",i+1,i+1,"e");
1872  FitClusterResidual(tmp, i, clusterResPerDE[i][1], gMeanResClYVsDE, gSigmaResClYVsDE);
1873  delete tmp;
1874  }
1875 
1876  // set DE graph labels
1877  TAxis* xAxis = static_cast<TH2*>(fTrackerList->UncheckedAt(kResClXVsDE))->GetXaxis();
1878  gMeanResClXVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1879  gMeanResClYVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1880  gSigmaResClXVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1881  gSigmaResClYVsDE->GetXaxis()->Set(fNDE, 0.5, fNDE+0.5);
1882  for (Int_t i = 1; i <= fNDE; i++) {
1883  const char* label = xAxis->GetBinLabel(i);
1884  gMeanResClXVsDE->GetXaxis()->SetBinLabel(i, label);
1885  gMeanResClYVsDE->GetXaxis()->SetBinLabel(i, label);
1886  gSigmaResClXVsDE->GetXaxis()->SetBinLabel(i, label);
1887  gSigmaResClYVsDE->GetXaxis()->SetBinLabel(i, label);
1888  }
1889 
1890  // --- diplay momentum residuals at vertex ---
1891  TCanvas* cResPAtVtx = DrawVsAng("cResPAtVtx", "momentum residual at vertex in 3 angular regions",
1892  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPAtVtx)),
1893  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsAngleAtAbsEnd)));
1894  fPAtVtxList->AddAt(cResPAtVtx, kcResPAtVtx);
1895  TCanvas* cResPAtVtxMC = DrawVsAng("cResPAtVtxMC", "momentum residual at vertex in 3 MC angular regions",
1896  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPAtVtx)),
1897  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsMCAngle)));
1898  fPAtVtxList->AddAt(cResPAtVtxMC, kcResPAtVtxMC);
1899  TCanvas* cResPAtVtxVsPosAbsEndMC = DrawVsPos("cResPAtVtxVsPosAbsEndMC", "momentum residual at vertex versus position at absorber end in 3 MC angular regions",
1900  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn02degMC)),
1901  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn23degMC)),
1902  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPosAbsEndIn310degMC)));
1903  fPAtVtxList->AddAt(cResPAtVtxVsPosAbsEndMC, kcResPAtVtxVsPosAbsEndMC);
1904  TCanvas* cResPAtVtxVsPIn23deg = DrawFitLandauGausResPVsP("cResPAtVtxVsPIn23deg", "momentum residual for tracks between 2 and 3 degrees",
1905  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn23deg)),
1906  10, "momentum residuals at vertex (tracks in [2,3] deg.)");
1907  fPAtVtxList->AddAt(cResPAtVtxVsPIn23deg, kcResPAtVtxVsPIn23deg);
1908  TCanvas* cResPAtVtxVsPIn310deg = DrawFitLandauGausResPVsP("cResPAtVtxVsPIn310deg", "momentum residual for tracks between 3 and 10 degrees",
1909  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn310deg)),
1910  10, "momentum residuals at vertex (tracks in [3,10] deg.)");
1911  fPAtVtxList->AddAt(cResPAtVtxVsPIn310deg, kcResPAtVtxVsPIn310deg);
1912  TCanvas* cResPAtVtxVsPIn02degMC = DrawResPVsP("cResPAtVtxVsPIn02degMC", "momentum residuals for tracks with MC angle < 2 degrees",
1913  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPAtVtxVsPIn02degMC)), 5);
1914  fPAtVtxList->AddAt(cResPAtVtxVsPIn02degMC, kcResPAtVtxVsPIn02degMC);
1915 
1916  // --- diplay slopeX residuals at vertex ---
1917  TCanvas* cResSlopeXAtVtx = DrawVsAng("cResSlopeXAtVtx", "slope_{X} residual at vertex in 3 angular regions",
1918  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeXAtVtx)),
1919  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsAngleAtAbsEnd)));
1920  fSlopeAtVtxList->AddAt(cResSlopeXAtVtx, kcResSlopeXAtVtx);
1921  TCanvas* cResSlopeYAtVtx = DrawVsAng("cResSlopeYAtVtx", "slope_{Y} residual at vertex in 3 angular regions",
1922  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeYAtVtx)),
1923  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsAngleAtAbsEnd)));
1924  fSlopeAtVtxList->AddAt(cResSlopeYAtVtx, kcResSlopeYAtVtx);
1925  TCanvas* cResSlopeXAtVtxMC = DrawVsAng("cResSlopeXAtVtxMC", "slope_{X} residual at vertex in 3 MC angular regions",
1926  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeXAtVtx)),
1927  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsMCAngle)));
1928  fSlopeAtVtxList->AddAt(cResSlopeXAtVtxMC, kcResSlopeXAtVtxMC);
1929  TCanvas* cResSlopeYAtVtxMC = DrawVsAng("cResSlopeYAtVtxMC", "slope_{Y} residual at vertex in 3 MC angular regions",
1930  static_cast<TH1*>(fTrackerList->UncheckedAt(kResSlopeYAtVtx)),
1931  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsMCAngle)));
1932  fSlopeAtVtxList->AddAt(cResSlopeYAtVtxMC, kcResSlopeYAtVtxMC);
1933  TCanvas* cResSlopeXAtVtxVsPosAbsEndMC = DrawVsPos("cResSlopeXAtVtxVsPosAbsEndMC", "slope_{X} residual at vertex versus position at absorber end in 3 MC angular regions",
1934  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn02degMC)),
1935  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn23degMC)),
1936  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeXAtVtxVsPosAbsEndIn310degMC)));
1937  fSlopeAtVtxList->AddAt(cResSlopeXAtVtxVsPosAbsEndMC, kcResSlopeXAtVtxVsPosAbsEndMC);
1938  TCanvas* cResSlopeYAtVtxVsPosAbsEndMC = DrawVsPos("cResSlopeYAtVtxVsPosAbsEndMC", "slope_{Y} residual at vertex versus position at absorber end in 3 MC angular regions",
1939  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn02degMC)),
1940  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn23degMC)),
1941  static_cast<TH2*>(fTrackerList->UncheckedAt(kResSlopeYAtVtxVsPosAbsEndIn310degMC)));
1942  fSlopeAtVtxList->AddAt(cResSlopeYAtVtxVsPosAbsEndMC, kcResSlopeYAtVtxVsPosAbsEndMC);
1943 
1944  // --- diplay eta residuals at vertex ---
1945  TCanvas* cResEtaAtVtx = DrawVsAng("cResEtaAtVtx", "eta residual at vertex in 3 angular regions",
1946  static_cast<TH1*>(fTrackerList->UncheckedAt(kResEtaAtVtx)),
1947  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsAngleAtAbsEnd)));
1948  fEtaAtVtxList->AddAt(cResEtaAtVtx, kcResEtaAtVtx);
1949  TCanvas* cResEtaAtVtxMC = DrawVsAng("cResEtaAtVtxMC", "eta residual at vertex in 3 MC angular regions",
1950  static_cast<TH1*>(fTrackerList->UncheckedAt(kResEtaAtVtx)),
1951  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsMCAngle)));
1952  fEtaAtVtxList->AddAt(cResEtaAtVtxMC, kcResEtaAtVtxMC);
1953  TCanvas* cResEtaAtVtxVsPosAbsEndMC = DrawVsPos("cResEtaAtVtxVsPosAbsEndMC", "eta residual at vertex versus position at absorber end in 3 MC angular regions",
1954  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn02degMC)),
1955  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn23degMC)),
1956  static_cast<TH2*>(fTrackerList->UncheckedAt(kResEtaAtVtxVsPosAbsEndIn310degMC)));
1957  fEtaAtVtxList->AddAt(cResEtaAtVtxVsPosAbsEndMC, kcResEtaAtVtxVsPosAbsEndMC);
1958 
1959  // --- diplay phi residuals at vertex ---
1960  TCanvas* cResPhiAtVtx = DrawVsAng("cResPhiAtVtx", "phi residual at vertex in 3 angular regions",
1961  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPhiAtVtx)),
1962  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsAngleAtAbsEnd)));
1963  fPhiAtVtxList->AddAt(cResPhiAtVtx, kcResPhiAtVtx);
1964  TCanvas* cResPhiAtVtxMC = DrawVsAng("cResPhiAtVtxMC", "phi residual at vertex in 3 MC angular regions",
1965  static_cast<TH1*>(fTrackerList->UncheckedAt(kResPhiAtVtx)),
1966  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsMCAngle)));
1967  fPhiAtVtxList->AddAt(cResPhiAtVtxMC, kcResPhiAtVtxMC);
1968  TCanvas* cResPhiAtVtxVsPosAbsEndMC = DrawVsPos("cResPhiAtVtxVsPosAbsEndMC", "phi residual at vertex versus position at absorber end in 3 MC angular regions",
1969  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn02degMC)),
1970  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn23degMC)),
1971  static_cast<TH2*>(fTrackerList->UncheckedAt(kResPhiAtVtxVsPosAbsEndIn310degMC)));
1972  fPhiAtVtxList->AddAt(cResPhiAtVtxVsPosAbsEndMC, kcResPhiAtVtxVsPosAbsEndMC);
1973 
1974  // --- diplay P*DCA ---
1975  TCanvas* cPDCA = DrawVsAng("cPDCA", "p #times DCA in 3 angular regions",
1976  static_cast<TH1*>(fTrackerList->UncheckedAt(kPDCA)),
1977  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsAngleAtAbsEnd)));
1978  fDCAList->AddAt(cPDCA, kcPDCA);
1979  TCanvas* cPDCAMC = DrawVsAng("cPDCAMC", "p #times DCA in 3 MC angular regions",
1980  static_cast<TH1*>(fTrackerList->UncheckedAt(kPDCA)),
1981  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsMCAngle)));
1982  fDCAList->AddAt(cPDCAMC, kcPDCAMC);
1983  TCanvas* cPDCAVsPosAbsEndMC = DrawVsPos("cPDCAVsPosAbsEndMC", "p #times DCA versus position at absorber end in 3 MC angular regions",
1984  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn02degMC)),
1985  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn23degMC)),
1986  static_cast<TH2*>(fTrackerList->UncheckedAt(kPDCAVsPosAbsEndIn310degMC)));
1987  fDCAList->AddAt(cPDCAVsPosAbsEndMC, kcPDCAVsPosAbsEndMC);
1988 
1989  // post param containers
1990  PostData(4, fEfficiencyList);
1991  PostData(5, fPAtVtxList);
1992  PostData(6, fSlopeAtVtxList);
1993  PostData(7, fEtaAtVtxList);
1994  PostData(8, fPhiAtVtxList);
1995  PostData(9, fPAt1stClList);
1996  PostData(10, fSlopeAt1stClList);
1997  PostData(11, fDCAList);
1998  PostData(12, fClusterList);
1999 }
2000 
2001 
2002 //________________________________________________________________________
2003 Bool_t AliAnalysisTaskMuonPerformance::GetEfficiency(AliCFEffGrid* efficiency, Double_t& calcEff, Double_t& calcEffErr)
2004 {
2005  //
2007  //
2008 
2009  Bool_t isGood = kTRUE;
2010 
2011  TH1* histo = 0x0;
2012  Double_t sum[2] = {0., 0.};
2013  for ( Int_t ihisto=0; ihisto<2; ihisto++ ) {
2014  histo = ( ihisto==0 ) ? efficiency->GetNum()->Project(kVarCharge) : efficiency->GetDen()->Project(kVarCharge);
2015  sum[ihisto] = histo->Integral();
2016  delete histo;
2017  }
2018 
2019  if ( sum[1] == 0. ) isGood = kFALSE;
2020 
2021  calcEff = ( isGood ) ? sum[0]/sum[1] : 0.;
2022  if ( calcEff > 1. ) isGood = kFALSE;
2023 
2024  calcEffErr = ( isGood ) ? TMath::Sqrt(calcEff*(1-calcEff)/sum[1]) : 0.;
2025 
2026  return isGood;
2027 }
2028 
2029 //________________________________________________________________________
2031 {
2032  //
2034  //
2035 
2036  if ( ! mcParticle )
2037  return kUnknownPart;
2038 
2039  Int_t recoPdg = mcParticle->PdgCode();
2040 
2041  // Track is not a muon
2042  if ( TMath::Abs(recoPdg) != 13 ) return kRecoHadron;
2043 
2044  Int_t imother = mcParticle->GetMother();
2045 
2046  Bool_t isFirstMotherHF = kFALSE;
2047  Int_t step = 0;
2048 
2049  //printf("\n"); // REMEMBER TO CUT
2050  while ( imother >= 0 ) {
2051  TParticle* part = static_cast<AliMCParticle*>(fMCEvent->GetTrack(imother))->Particle();
2052  //printf("Particle %s -> %s\n", part->GetName(), TMCProcessName[part->GetUniqueID()]); // REMEMBER TO CUT
2053 
2054  if ( part->GetUniqueID() == kPHadronic )
2055  return kSecondaryMu;
2056 
2057  Int_t absPdg = TMath::Abs(part->GetPdgCode());
2058 
2059  step++;
2060  if ( step == 1 ) // only for first mother
2061  isFirstMotherHF = ( ( absPdg >= 400 && absPdg < 600 ) ||
2062  ( absPdg >= 4000 && absPdg < 6000 ) );
2063 
2064  if ( absPdg < 4 )
2065  return kPrimaryMu;
2066 
2067  if ( isFirstMotherHF) {
2068  if ( absPdg == 4 )
2069  return kCharmMu;
2070  else if ( absPdg == 5 )
2071  return kBeautyMu;
2072  }
2073 
2074  imother = part->GetFirstMother();
2075  }
2076 
2077  return kPrimaryMu;
2078 
2079 }
2080 
2081 //________________________________________________________________________
2083 {
2084  //
2086  //
2087  Float_t thetaDeg = ( isTheta ) ? RAtAbsEnd : TMath::ATan( RAtAbsEnd / 505. );
2088  thetaDeg *= TMath::RadToDeg();
2089  if ( thetaDeg < 2. )
2090  return 0.;
2091  else if ( thetaDeg < 3. )
2092  return 1.;
2093  else if ( thetaDeg < 10. )
2094  return 2.;
2095 
2096  return 3.;
2097 }
2098 
2099 //________________________________________________________________________
2100 void AliAnalysisTaskMuonPerformance::FillContainerInfoReco(Double_t* containerInput, AliESDMuonTrack* esdTrack,
2101  Bool_t isValid, Int_t mcID)
2102 {
2103  //
2105  //
2106 
2107  AliMCParticle* mcPart = (mcID >= 0) ? static_cast<AliMCParticle*>(fMCEvent->GetTrack(mcID)) : 0x0;
2108 
2109  if (fUseMCKinematics && mcPart) {
2110  containerInput[kVarPt] = mcPart->Pt();
2111  containerInput[kVarEta] = mcPart->Eta();
2112  containerInput[kVarPhi] = mcPart->Phi();
2113  } else {
2114  containerInput[kVarPt] = esdTrack->Pt();
2115  containerInput[kVarEta] = esdTrack->Eta();
2116  containerInput[kVarPhi] = esdTrack->Phi();
2117  }
2118  containerInput[kVarThetaZones] = GetBinThetaAbsEnd(esdTrack->GetRAtAbsorberEnd());
2119  containerInput[kVarCharge] = static_cast<Double_t>(esdTrack->Charge());
2120  containerInput[kVarHasTracker] = static_cast<Double_t>(esdTrack->ContainTrackerData() && isValid);
2121  if (esdTrack->GetMatchTrigger() == 0) containerInput[kVarTrigger] = static_cast<Double_t>(kNoMatchTrig);
2122  else if (esdTrack->GetMatchTrigger() == 1) containerInput[kVarTrigger] = static_cast<Double_t>(kAllPtTrig);
2123  else if (esdTrack->GetMatchTrigger() == 2) containerInput[kVarTrigger] = static_cast<Double_t>(kLowPtTrig);
2124  else if (esdTrack->GetMatchTrigger() == 3) containerInput[kVarTrigger] = static_cast<Double_t>(kHighPtTrig);
2125  containerInput[kVarMotherType] = static_cast<Double_t>(RecoTrackMother(mcPart));
2126 
2127 }
2128 
2129 //________________________________________________________________________
2130 void AliAnalysisTaskMuonPerformance::FillContainerInfoMC(Double_t* containerInput, AliMCParticle* mcPart)
2131 {
2132  //
2134  //
2135 
2136  containerInput[kVarPt] = mcPart->Pt();
2137  containerInput[kVarEta] = mcPart->Eta();
2138  containerInput[kVarPhi] = mcPart->Phi();
2139  containerInput[kVarThetaZones] = GetBinThetaAbsEnd(TMath::Pi()-mcPart->Theta(),kTRUE);
2140  containerInput[kVarCharge] = static_cast<Double_t>(mcPart->Charge())/3.;
2141  containerInput[kVarMotherType] = static_cast<Double_t>(RecoTrackMother(mcPart));
2142 
2143 }
2144 
2145 //________________________________________________________________________
2147 {
2148  //Fit parameters:
2149  //par[0]=Width (scale) parameter of Landau density
2150  //par[1]=Most Probable (MP, location) parameter of Landau density
2151  //par[2]=Total area (integral -inf to inf, normalization constant)
2152  //par[3]=Width (sigma) of convoluted Gaussian function
2153  //
2154  //In the Landau distribution (represented by the CERNLIB approximation),
2155  //the maximum is located at x=-0.22278298 with the location parameter=0.
2156  //This shift is corrected within this function, so that the actual
2157  //maximum is identical to the MP parameter.
2158 
2159  // Numeric constants
2160  Double_t invsq2pi = 0.3989422804014; // (2 pi)^(-1/2)
2161  Double_t mpshift = -0.22278298; // Landau maximum location
2162 
2163  // MP shift correction
2164  Double_t mpc = par[1] - mpshift * par[0];
2165 
2166  // Range of convolution integral
2167  Double_t xlow = x[0] - 5. * par[3];
2168  Double_t xupp = x[0] + 5. * par[3];
2169  Double_t step = TMath::Min(0.2*par[0],0.1*par[3]);
2170 
2171  // Convolution integral of Landau and Gaussian by sum
2172  Double_t xx = xlow + 0.5 * step;
2173  Double_t sum = 0.0;
2174  while (xx < xupp) {
2175  sum += TMath::Landau(-xx,mpc,par[0]) * TMath::Gaus(x[0],xx,par[3]);
2176  xx += step;
2177  }
2178 
2179  return (par[2] * step * sum * invsq2pi / par[3] / par[0]);
2180 
2181 }
2182 
2183 //________________________________________________________________________
2185  TGraphAsymmErrors* gMostProb, TGraphAsymmErrors* gSigma)
2186 {
2188 
2189  static TF1 *fGaus2 = 0x0;
2190  if (!fGaus2) fGaus2 = new TF1("fGaus2","gaus");
2191  static TF1* fLandauGaus = 0x0;
2192  if (!fLandauGaus) {
2193  fLandauGaus = new TF1("fLandauGaus",langaufun,h->GetYaxis()->GetBinLowEdge(1),h->GetYaxis()->GetBinLowEdge(h->GetNbinsY()+1),4);
2194  fLandauGaus->SetNpx(10000);
2195  }
2196 
2197  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2198  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2199 
2200  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2201 
2202  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2203 
2204  // first fit
2205  fGaus2->SetParameters((Double_t)tmp->GetEntries(), 0., 1.);
2206  tmp->Fit("fGaus2", "WWNQ");
2207 
2208  // rebin histo
2209  Double_t sigma = fGaus2->GetParameter(2);
2210  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*tmp->GetNbinsX(),TMath::Max(0.5*sigma/tmp->GetBinWidth(1),1.)));
2211  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2212  tmp->Rebin(rebin);
2213 
2214  // second fit
2215  Double_t mean = fGaus2->GetParameter(1);
2216  fLandauGaus->SetParameters(0.25*sigma*TMath::Sqrt(8.*log(2.)), mean, tmp->GetEntries()*tmp->GetXaxis()->GetBinWidth(1), 0.5*sigma);
2217  fLandauGaus->SetParLimits(0, 0.0025*sigma*TMath::Sqrt(8.*log(2.)), 1000.);
2218  fLandauGaus->SetParLimits(3, 0., 2.*sigma);
2219  Double_t xMin = TMath::Max(mean-50.*sigma, tmp->GetXaxis()->GetXmin());
2220  Double_t xMax = TMath::Min(mean+10.*sigma, tmp->GetXaxis()->GetXmax());
2221  if (xMin < tmp->GetXaxis()->GetXmax() && xMax > tmp->GetXaxis()->GetXmin()) fLandauGaus->SetRange(xMin, xMax);
2222  tmp->Fit("fLandauGaus","RNQ");
2223 
2224  // get fit results and fill graph
2225  Double_t fwhm = 4.*fLandauGaus->GetParameter(0);
2226  sigma = fLandauGaus->GetParameter(3);
2227  Double_t sigmaP = TMath::Sqrt(sigma*sigma + fwhm*fwhm/(8.*log(2.)));
2228  Double_t fwhmErr = fLandauGaus->GetParError(0);
2229  Double_t sigmaErr = fLandauGaus->GetParError(3);
2230  Double_t sigmaPErr = TMath::Sqrt(sigma*sigma*sigmaErr*sigmaErr + fwhm*fwhm*fwhmErr*fwhmErr/(64.*log(2.)*log(2.))) / sigmaP;
2231  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2232  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2233  h->GetXaxis()->SetRange();
2234  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2235  gMean->SetPoint(i/rebinFactorX-1, p, tmp->GetMean());
2236  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], tmp->GetMeanError(), tmp->GetMeanError());
2237  gMostProb->SetPoint(i/rebinFactorX-1, p, -fLandauGaus->GetParameter(1));
2238  gMostProb->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fLandauGaus->GetParError(1), fLandauGaus->GetParError(1));
2239  gSigma->SetPoint(i/rebinFactorX-1, p, sigmaP);
2240  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], sigmaPErr, sigmaPErr);
2241 
2242  // clean memory
2243  delete tmp;
2244  }
2245 
2246  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2247 
2248 }
2249 
2250 //________________________________________________________________________
2252  const Double_t sigma0, const char* fitting,
2253  TGraphAsymmErrors* gMean, TGraphAsymmErrors* gSigma)
2254 {
2256 
2257  static TF1* fGaus = 0x0;
2258  if (!fGaus) fGaus = new TF1("fGaus","gaus");
2259 
2260  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2261  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2262 
2263  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2264 
2265  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2266 
2267  // first fit
2268  fGaus->SetParameters(tmp->GetEntries(), mean0, sigma0);
2269  tmp->Fit("fGaus","WWNQ");
2270 
2271  // rebin histo
2272  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*tmp->GetNbinsX(),TMath::Max(0.5*fGaus->GetParameter(2)/tmp->GetBinWidth(1),1.)));
2273  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2274  tmp->Rebin(rebin);
2275 
2276  // second fit
2277  tmp->Fit("fGaus","NQ");
2278 
2279  // get fit results and fill histograms
2280  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2281  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2282  h->GetXaxis()->SetRange();
2283  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2284  gMean->SetPoint(i/rebinFactorX-1, p, fGaus->GetParameter(1));
2285  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fGaus->GetParError(1), fGaus->GetParError(1));
2286  gSigma->SetPoint(i/rebinFactorX-1, p, fGaus->GetParameter(2));
2287  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fGaus->GetParError(2), fGaus->GetParError(2));
2288 
2289  // clean memory
2290  delete tmp;
2291  }
2292 
2293  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2294 
2295 }
2296 
2297 //________________________________________________________________________
2299  TGraphAsymmErrors* gMean, TGraphAsymmErrors* gSigma)
2300 {
2302 
2303  static TF1* fPGaus = 0x0;
2304  if (!fPGaus) fPGaus = new TF1("fPGaus","x*gaus");
2305 
2306  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/fNPBins, 1);
2307  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2308 
2309  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,fNPBins)<<flush;
2310 
2311  TH1D *tmp = h->ProjectionY("tmp",i-rebinFactorX+1,i,"e");
2312 
2313  // rebin histo
2314  Int_t rebin = static_cast<Int_t>(25*(tmp->GetNbinsX()/(tmp->GetBinLowEdge(tmp->GetNbinsX()+1)-tmp->GetBinLowEdge(1))));
2315  while (tmp->GetNbinsX()%rebin!=0) rebin--;
2316  tmp->Rebin(rebin);
2317 
2318  // fit
2319  fPGaus->SetParameters(1.,0.,80.);
2320  tmp->Fit("fPGaus","NQ");
2321 
2322  // get fit results and fill histograms
2323  h->GetXaxis()->SetRange(i-rebinFactorX+1,i);
2324  Double_t p = (tmp->GetEntries() > 0) ? h->GetMean() : 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2325  h->GetXaxis()->SetRange();
2326  Double_t pErr[2] = {p-h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1), h->GetXaxis()->GetBinLowEdge(i+1)-p};
2327  gMean->SetPoint(i/rebinFactorX-1, p, fPGaus->GetParameter(1));
2328  gMean->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fPGaus->GetParError(1), fPGaus->GetParError(1));
2329  gSigma->SetPoint(i/rebinFactorX-1, p, fPGaus->GetParameter(2));
2330  gSigma->SetPointError(i/rebinFactorX-1, pErr[0], pErr[1], fPGaus->GetParError(2), fPGaus->GetParError(2));
2331 
2332  // clean memory
2333  delete tmp;
2334  }
2335 
2336  cout<<Form("\rFitting %s... %d/%d",fitting,fNPBins,fNPBins)<<endl;
2337 
2338 }
2339 
2340 //________________________________________________________________________
2342  TGraphErrors* gMean, TGraphErrors* gSigma)
2343 {
2345 
2346  static TF1* fRGaus = 0x0;
2347  Double_t mean, meanErr, sigmaErr;
2348 
2349  if (fFitResiduals) {
2350 
2351  if (!fRGaus) fRGaus = new TF1("fRGaus","gaus");
2352 
2353  // first fit
2354  Double_t xMin = h->GetXaxis()->GetXmin();
2355  Double_t xMax = h->GetXaxis()->GetXmax();
2356  fRGaus->SetRange(xMin, xMax);
2357  fRGaus->SetParameters(h->GetEntries(), 0., 0.1);
2358  fRGaus->SetParLimits(1, xMin, xMax);
2359  h->Fit("fRGaus", "WWNQ");
2360 
2361  // rebin histo
2362  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*h->GetNbinsX(),TMath::Max(0.3*fRGaus->GetParameter(2)/h->GetBinWidth(1),1.)));
2363  while (h->GetNbinsX()%rebin!=0) rebin--;
2364  h->Rebin(rebin);
2365 
2366  // second fit
2367  xMin = TMath::Max(fRGaus->GetParameter(1)-10.*fRGaus->GetParameter(2), h->GetXaxis()->GetXmin());
2368  xMax = TMath::Min(fRGaus->GetParameter(1)+10.*fRGaus->GetParameter(2), h->GetXaxis()->GetXmax());
2369  fRGaus->SetRange(xMin, xMax);
2370  fRGaus->SetParLimits(1, xMin, xMax);
2371  h->Fit("fRGaus","NQR");
2372 
2373  mean = fRGaus->GetParameter(1);
2374  meanErr = fRGaus->GetParError(1);
2375  sigma = fRGaus->GetParameter(2);
2376  sigmaErr = fRGaus->GetParError(2);
2377 
2378  } else {
2379 
2380  Zoom(h);
2381  mean = h->GetMean();
2382  meanErr = h->GetMeanError();
2383  sigma = h->GetRMS();
2384  sigmaErr = h->GetRMSError();
2385  h->GetXaxis()->SetRange(0,0);
2386 
2387  }
2388 
2389  gMean->SetPoint(i, i+1, mean);
2390  gMean->SetPointError(i, 0., meanErr);
2391 
2392  if (fCorrectForSystematics) {
2393  Double_t s = TMath::Sqrt(mean*mean + sigma*sigma);
2394  sigmaErr = (s>0.) ? TMath::Sqrt(sigma*sigma*sigmaErr*sigmaErr + mean*mean*meanErr*meanErr) / s : 0.;
2395  sigma = s;
2396  }
2397 
2398  gSigma->SetPoint(i, i+1, sigma);
2399  gSigma->SetPointError(i, 0., sigmaErr);
2400 
2401 }
2402 
2403 //________________________________________________________________________
2404 TCanvas* AliAnalysisTaskMuonPerformance::DrawVsAng(const char* name, const char* title, TH1* h1, TH2* h2)
2405 {
2407  TCanvas* c = new TCanvas(name, title);
2408  c->cd();
2409  h1->DrawCopy();
2410  TH1D *proj1 = h2->ProjectionY(Form("%s_proj_0_2",h2->GetName()),1,2);
2411  proj1->SetLineColor(2);
2412  proj1->Draw("sames");
2413  TH1D *proj2 = h2->ProjectionY(Form("%s_proj_2_3",h2->GetName()),3,3);
2414  proj2->SetLineColor(4);
2415  proj2->Draw("sames");
2416  TH1D *proj3 = h2->ProjectionY(Form("%s__proj_3_10",h2->GetName()),4,10);
2417  proj3->SetLineColor(3);
2418  proj3->Draw("sames");
2419  return c;
2420 }
2421 
2422 //________________________________________________________________________
2423 TCanvas* AliAnalysisTaskMuonPerformance::DrawVsPos(const char* name, const char* title, TH2* h1, TH2* h2, TH2* h3)
2424 {
2426  TCanvas* c = new TCanvas(name, title);
2427  c->cd();
2428  h1->SetMarkerColor(2);
2429  h1->DrawCopy();
2430  h2->SetMarkerColor(4);
2431  h2->DrawCopy("sames");
2432  h3->SetMarkerColor(3);
2433  h3->DrawCopy("sames");
2434  return c;
2435 }
2436 
2437 //________________________________________________________________________
2439  TH2* h, const Int_t nBins, const char* fitting)
2440 {
2442 
2443  static TF1 *fGaus3 = 0x0;
2444  if (!fGaus3) fGaus3 = new TF1("fGaus3","gaus");
2445  static TF1* fLandauGaus2 = 0x0;
2446  if (!fLandauGaus2) {
2447  fLandauGaus2 = new TF1("fLandauGaus2",langaufun,h->GetYaxis()->GetBinLowEdge(1),h->GetYaxis()->GetBinLowEdge(h->GetNbinsY()+1),4);
2448  fLandauGaus2->SetNpx(10000);
2449  }
2450 
2451  TLegend* l = new TLegend(0.15,0.25,0.3,0.85);
2452  TCanvas* c = new TCanvas(name, title);
2453  c->cd();
2454 
2455  h->Sumw2();
2456 
2457  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/nBins, 1);
2458  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2459 
2460  cout<<Form("\rFitting %s... %d/%d",fitting,i/rebinFactorX,nBins)<<flush;
2461 
2462  // draw projection
2463  TH1D* proj = h->ProjectionY(Form("%s_%d",h->GetName(),i/rebinFactorX),i-rebinFactorX+1,i);
2464  if (proj->GetEntries() > 0) proj->Scale(1./proj->GetEntries());
2465  proj->Draw((i==rebinFactorX)?"hist":"histsames");
2466  proj->SetLineColor(i/rebinFactorX);
2467 
2468  if (proj->GetEntries() > 10) {
2469 
2470  // first fit
2471  fGaus3->SetParameters(1., 0., 1.);
2472  proj->Fit("fGaus3", "WWNQ");
2473 
2474  // rebin histo
2475  Double_t sigma = fGaus3->GetParameter(2);
2476  Int_t rebin = static_cast<Int_t>(TMath::Min(0.1*proj->GetNbinsX(),TMath::Max(0.5*sigma/proj->GetBinWidth(1),1.)));
2477  while (proj->GetNbinsX()%rebin!=0) rebin--;
2478  proj->Rebin(rebin);
2479  proj->Scale(1./rebin);
2480 
2481  // second fit
2482  Double_t mean = fGaus3->GetParameter(1);
2483  fLandauGaus2->SetParameters(0.25*sigma*TMath::Sqrt(8.*log(2.)), mean, proj->GetXaxis()->GetBinWidth(1), 0.5*sigma);
2484  fLandauGaus2->SetParLimits(0, 0.0025*sigma*TMath::Sqrt(8.*log(2.)), 1000.);
2485  fLandauGaus2->SetParLimits(3, 0., 2.*sigma);
2486  Double_t xMin = TMath::Max(mean-50.*sigma, proj->GetXaxis()->GetXmin());
2487  Double_t xMax = TMath::Min(mean+10.*sigma, proj->GetXaxis()->GetXmax());
2488  if (xMin < proj->GetXaxis()->GetXmax() && xMax > proj->GetXaxis()->GetXmin()) fLandauGaus2->SetRange(xMin, xMax);
2489  fLandauGaus2->SetLineColor(i/rebinFactorX);
2490  proj->Fit("fLandauGaus2","RQ","sames");
2491 
2492  }
2493 
2494  // set label
2495  Double_t p = 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2496  l->AddEntry(proj,Form("%5.1f GeV",p));
2497 
2498  }
2499 
2500  cout<<Form("\rFitting %s... %d/%d",fitting,nBins,nBins)<<endl;
2501 
2502  l->Draw("same");
2503 
2504  return c;
2505 
2506 }
2507 
2508 //________________________________________________________________________
2509 TCanvas* AliAnalysisTaskMuonPerformance::DrawResPVsP(const char* name, const char* title, TH2* h, const Int_t nBins)
2510 {
2512 
2513  TLegend* l = new TLegend(0.15,0.25,0.3,0.85);
2514  TCanvas* c = new TCanvas(name, title);
2515  c->cd();
2516 
2517  h->Sumw2();
2518 
2519  Int_t rebinFactorX = TMath::Max(h->GetNbinsX()/nBins, 1);
2520  for (Int_t i = rebinFactorX; i <= h->GetNbinsX(); i+=rebinFactorX) {
2521 
2522  // draw projection
2523  TH1D* proj = h->ProjectionY(Form("%s_%d",h->GetName(),i/rebinFactorX),i-rebinFactorX+1,i);
2524  if (proj->GetEntries() > 0) proj->Scale(1./proj->GetEntries());
2525  proj->Draw((i==rebinFactorX)?"hist":"histsames");
2526  proj->SetLineColor(i/rebinFactorX);
2527  proj->SetLineWidth(2);
2528 
2529  // set label
2530  Double_t p = 0.5 * (h->GetXaxis()->GetBinLowEdge(i-rebinFactorX+1) + h->GetXaxis()->GetBinLowEdge(i+1));
2531  l->AddEntry(proj,Form("%5.1f GeV",p));
2532 
2533  }
2534 
2535  l->Draw("same");
2536 
2537  return c;
2538 }
2539 
2540 //________________________________________________________________________
2542 {
2544 
2545  Double_t maxEventsCut = fractionCut * h->GetEntries();
2546  Int_t nBins = h->GetNbinsX();
2547 
2548  // set low edge
2549  Int_t minBin;
2550  Double_t eventsCut = 0.;
2551  for (minBin = 1; minBin <= nBins; minBin++) {
2552  eventsCut += h->GetBinContent(minBin);
2553  if (eventsCut > maxEventsCut) break;
2554  }
2555 
2556  // set high edge
2557  Int_t maxBin;
2558  eventsCut = 0.;
2559  for (maxBin = nBins; maxBin >= 1; maxBin--) {
2560  eventsCut += h->GetBinContent(maxBin);
2561  if (eventsCut > maxEventsCut) break;
2562  }
2563 
2564  // set new axis range
2565  h->GetXaxis()->SetRange(--minBin, ++maxBin);
2566 }
2567 
2568 //________________________________________________________________________
2569 void AliAnalysisTaskMuonPerformance::FillEffHistos(AliCFEffGrid* efficiency, const char* suffix, TObjArray* list)
2570 {
2572 
2573  TH1* auxHisto = efficiency->Project(kVarPt);
2574  auxHisto->SetName(Form("effVsPt_%s",suffix));
2575  list->AddLast(auxHisto);
2576 
2577  auxHisto = efficiency->Project(kVarEta);
2578  auxHisto->SetName(Form("effVsEta_%s",suffix));
2579  list->AddLast(auxHisto);
2580 
2581  auxHisto = efficiency->Project(kVarPhi);
2582  auxHisto->SetName(Form("effVsPhi_%s",suffix));
2583  list->AddLast(auxHisto);
2584 
2585  auxHisto = efficiency->Project(kVarCent);
2586  auxHisto->SetName(Form("effVsCent_%s",suffix));
2587  list->AddLast(auxHisto);
2588 
2589 }
2590 
slope-Y residual at vertex versus position at absorber end for tracks with MC angle in ]3...
phi residual at vertex in 3 angular regions
TCanvas * DrawVsAng(const char *name, const char *title, TH1 *h1, TH2 *h2)
void Zoom(TH1 *h, Double_t fractionCut=0.01)
Int_t RecoTrackMother(AliMCParticle *mcParticle)
double Double_t
Definition: External.C:58
Definition: External.C:260
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 ]2,3] degrees at absorber end.
Definition: External.C:236
const char * title
Definition: MakeQAPdf.C:27
Residual of x position in first trigger chamber.
Int_t fDEIndices[1100]
index of DE in histograms refered by ID
Bool_t fUseMCKinematics
use kinematics from MC track when available
Definition: External.C:244
slope-X residual at vertex versus angle at absorber end
slope-Y residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
phi residual at vertex in 3 MC angular regions
momentum residual at vertex versus P for tracks in ]3,10[ degrees at absorber end ...
TObjArray * fSlopeAt1stClList
List of graph and canvas about slope resolution at first cluster.
centrality
P * DCA distribution versus P for tracks in ]2,3] degrees at absorber end.
TObjArray * fClusterList
List of graph and canvas about cluster resolution.
phi residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
void FillEffHistos(AliCFEffGrid *efficiency, const char *suffix, TObjArray *list)
Double_t phiMin
eta residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
void FitClusterResidual(TH1 *h, Int_t i, Double_t &sigma, TGraphErrors *gMean, TGraphErrors *gSigma)
momentum residual at first cluster versus P
mean P * MCS deviation angle versus P for tracks in ]3,10[ degrees at absorber end ...
slope-Y residual at vertex versus position at absorber end for tracks with MC angle in ]2...
P * DCA resolution versus P for tracks in ]3,10[ degrees at absorber end.
TCanvas * c
Definition: TestFitELoss.C:172
transverse momentum residual at vertex versus pT
TObjArray * fTriggerList
List of histograms for trigger resolution.
Double_t ptMin
phi residual at vertex versus position at absorber end for tracks with MC angle in ]3...
momentum residual at vertex versus angle at absorber end
Int_t nCentBins
slope-X residual at vertex in 3 angular regions
Match trigger not passing any Pt threshold (MC)
P * MCS deviation angle dispersion versus P for tracks in ]2,3] degrees at absorber end...
P * DCA distribution versus angle at absorber end.
void FitPDCAVsMom(TH2 *h, const char *fitting, TGraphAsymmErrors *gMean, TGraphAsymmErrors *gSigma)
slope-X residual at vertex versus position at absorber end in 3 MC angular regions ...
UInt_t fRequestedStationMask
mask of requested stations
momentum residual at vertex versus position at absorber end in 3 MC angular regions ...
slope-Y residual at vertex versus angle at absorber end
Double_t * sigma
const Int_t nPtBins
matched with reconstructible track and triggerable track of different ID
TCanvas * DrawResPVsP(const char *name, const char *title, TH2 *h, const Int_t nBins)
momentum residual at vertex versus angle at absorber end versus momentum
eta residual at vertex in 3 MC angular regions
slope-X residual at vertex versus position at absorber end for tracks with MC angle in ]3...
void FillContainerInfoMC(Double_t *containerInput, AliMCParticle *mcPart)
momentum residual at vertex versus P for tracks with MC angle < 2 degrees
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
TObjArray * fPhiAtVtxList
List of graph and canvas about phi resolution at vertex.
mean P * DCA versus P for tracks in ]3,10[ degrees at absorber end
float Float_t
Definition: External.C:68
not matched with either reconstructible track or triggerable track
eta residual at vertex versus angle at absorber end
Double_t fSigmaCutTrig
sigma cut to associate trigger track to triggerable track
momentum residual at vertex versus position at absorber end for tracks with MC angle in ]2...
P * MCS deviation angle distribution versus P for tracks in ]3,10[ degrees at absorber end...
eta residual at vertex in 3 angular regions
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)
P * DCA distribution versus position at absorber end for tracks with MC angle in ]3,10[ degrees.
transverse momentum residual at first cluster versus pT
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)
Double_t phiMax
TObjArray * fEfficiencyList
List of histograms for tracker/trigger efficiencies.
Definition: External.C:212
relative momentum resolution at first cluster versus P
Bool_t fMCTrigLevelFromMatchTrk
set the triggerable level from the MC matching the trk part
Bool_t fRequest2ChInSameSt45
2 fired chambers requested in the same station (4 or 5) or not
TObjArray * fPAtVtxList
List of graph and canvas about momentum resolution at vertex.
momentum residual at vertex versus position at absorber end for tracks with MC angle in ]3...
momentum residual at vertex in 3 angular regions
TObjArray * fDCAList
List of graph and canvas about DCA.
mean momentum residual at first cluster versus P
TObjArray * fTrackerList
List of histograms for tracker resolution.
P * DCA distribution versus P for tracks in ]3,10[ degrees at absorber end.
Bool_t Data(TH1F *h, Double_t *rangefit, Bool_t writefit, Double_t &sgn, Double_t &errsgn, Double_t &bkg, Double_t &errbkg, Double_t &sgnf, Double_t &errsgnf, Double_t &sigmafit, Int_t &status)
TString fDefaultStorage
location of the default OCDB storage
TObjArray * fPAt1stClList
List of graph and canvas about momentum resolution at first cluster.
most probable momentum residual at vertex versus P
relative momentum resolution at vertex versus P
eta residual at vertex versus position at absorber end for tracks with MC angle in ]3...
Double_t centMax
P * DCA distribution versus position at absorber end for tracks with MC angle in ]2,3] degrees.
slope-X residual at vertex versus position at absorber end for tracks with MC angle in ]2...
momentum residual at vertex versus position at absorber end for tracks with MC angle <= 2 degrees ...
mean P * DCA versus P for tracks in ]2,3] degrees at absorber end
slope-Y residual at vertex in 3 angular regions
eta residual at vertex versus position at absorber end in 3 MC angular regions
phi 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.
momentum residual at vertex versus P for tracks in ]2,3] degrees at absorber end
P * DCA distribution versus position at absorber end for tracks with MC angle <= 2 degrees...
Definition: External.C:220
slope-Y residual at vertex in 3 MC angular regions
slope-X residual at vertex in 3 MC angular regions
matched with reconstructible track and triggerable track of same ID
phi residual at vertex versus angle at absorber end
AliCFContainer * fCFContainer
Pointer to the CF container.
Int_t rebin
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 for tracks between 2 and 3 degrees
momentum residual at vertex in 3 MC angular regions
slope-Y residual at vertex versus position at absorber end in 3 MC angular regions ...
momentum residuals for tracks with MC angle < 2 degrees
TString fRecoParamOCDBpath
OCDB path to the recoParam file.
const char Option_t
Definition: External.C:48
Float_t GetBinThetaAbsEnd(Float_t RAtAbsEnd, Bool_t isTheta=kFALSE)
Residual of y position in first trigger chamber.
phi residual at vertex versus position at absorber end for tracks with MC angle in ]2...
const Int_t nbins
momentum residual for tracks between 3 and 10 degrees
bool Bool_t
Definition: External.C:53
eta residual at vertex versus position at absorber end for tracks with MC angle in ]2...
Double_t ptMax
Bool_t fCorrectForSystematics
add or not the systematic shifts of the residuals to the resolution
Double_t langaufun(Double_t *x, Double_t *par)
void FillContainerInfoReco(Double_t *containerInput, AliESDMuonTrack *esdTrack, Bool_t isValid, Int_t mcID)
P * MCS deviation angle dispersion 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 ...
P * DCA versus position at absorber end in 3 MC angular regions.
Definition: External.C:196
TString fAlignOCDBpath
OCDB path to the alignment file.
P * MCS deviation angle distribution versus P for tracks in ]2,3] degrees at absorber end...
mean P * MCS deviation angle versus P for tracks in ]2,3] degrees at absorber end ...
Double_t centMin
Double_t fClusterMaxRes[2]
highest chamber resolution in both directions
Bool_t GetEfficiency(AliCFEffGrid *efficiency, Double_t &calcEff, Double_t &calcEffErr)