1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
18 //Class to calculate the intrinsic efficiency of the detection elements of the
19 //MUON tracking chambers in function of the position in the detection element.
21 //Author: Nicolas LE BRIS - SUBATECH Nantes
22 // Modified by Matthieu LENHARDT - SUBATECH Nantes
23 // Modified by Antoine LARDEUX - SUBATECH Nantes
28 #include <THnSparse.h>
29 #include <TObjArray.h>
30 #include <TGeoGlobalMagField.h>
34 #include "AliESDEvent.h"
35 #include "AliESDMuonTrack.h"
36 #include "AliGeomManager.h"
37 #include "AliCDBManager.h"
38 #include "AliESDVZERO.h"
41 #include "AliInputEventHandler.h"
42 #include "AliAnalysisManager.h"
43 #include "AliAnalysisTaskMuonTrackingEff.h"
44 #include "AliCentrality.h"
45 #include "AliVVertex.h"
48 #include "AliMUONCDB.h"
49 #include "AliMUONESDInterface.h"
50 #include "AliMUONGeometryTransformer.h"
51 #include "AliMUONTrack.h"
52 #include "AliMUONTrackParam.h"
53 #include "AliMUONTrackExtrap.h"
54 #include "AliMUONVCluster.h"
55 #include "AliMUONConstants.h"
57 //include MUON/mapping:
58 #include "AliMpDEManager.h"
59 #include "AliMpSegmentation.h"
60 #include "AliMpVSegmentation.h"
63 ClassImp(AliAnalysisTaskMuonTrackingEff)
65 const Int_t AliAnalysisTaskMuonTrackingEff::fgkNbrOfDetectionElt[10] = {4, 4, 4, 4, 18, 18, 26, 26, 26, 26};
66 const Int_t AliAnalysisTaskMuonTrackingEff::fgkOffset = 100;
68 //________________________________________________________________________
69 AliAnalysisTaskMuonTrackingEff::AliAnalysisTaskMuonTrackingEff() :
79 fCurrentCentrality(0.),
82 fDetEltTDHistList(0x0),
83 fDetEltTTHistList(0x0),
84 fDetEltSDHistList(0x0),
85 fChamberTDHistList(0x0),
86 fChamberTTHistList(0x0),
87 fChamberSDHistList(0x0),
90 /// Default constructor
93 //________________________________________________________________________
94 AliAnalysisTaskMuonTrackingEff::AliAnalysisTaskMuonTrackingEff(TString name) :
95 AliAnalysisTaskSE(name),
104 fCurrentCentrality(100.),
107 fDetEltTDHistList(0x0),
108 fDetEltTTHistList(0x0),
109 fDetEltSDHistList(0x0),
110 fChamberTDHistList(0x0),
111 fChamberTTHistList(0x0),
112 fChamberSDHistList(0x0),
117 // Output slots 0 to 5 writes into a TClonesArray:
118 DefineOutput(1, TList::Class());
119 DefineOutput(2, TList::Class());
120 DefineOutput(3, TList::Class());
121 DefineOutput(4, TList::Class());
122 DefineOutput(5, TList::Class());
123 DefineOutput(6, TList::Class());
124 DefineOutput(7, TList::Class());
127 //________________________________________________________________________
128 AliAnalysisTaskMuonTrackingEff::~AliAnalysisTaskMuonTrackingEff()
131 if (!AliAnalysisManager::GetAnalysisManager()->IsProofMode()) {
132 delete fDetEltTDHistList;
133 delete fDetEltTTHistList;
134 delete fDetEltSDHistList;
135 delete fChamberTDHistList;
136 delete fChamberTTHistList;
137 delete fChamberSDHistList;
138 delete fExtraHistList;
143 //________________________________________________________________________
144 void AliAnalysisTaskMuonTrackingEff::NotifyRun()
146 /// Load the OCDB and the Geometry
149 if (fOCDBLoaded) return;
152 AliCDBManager* man = AliCDBManager::Instance();
153 if (man->IsDefaultStorageSet()) printf("EfficiencyTask: CDB default storage already set!\n");
154 else man->SetDefaultStorage(fOCDBpath.Data());
155 if (man->GetRun() > -1) printf("EfficiencyTask: run number already set!\n");
156 else man->SetRun(fCurrentRunNumber);
159 if (!AliGeomManager::GetGeometry()) {
160 AliGeomManager::LoadGeometry();
161 if (!AliGeomManager::GetGeometry()) return;
162 if (!AliGeomManager::ApplyAlignObjsFromCDB("MUON")) return;
164 fTransformer = new AliMUONGeometryTransformer();
165 fTransformer->LoadGeometryData();
167 // Magnetic field for track extrapolation
168 if (!TGeoGlobalMagField::Instance()->GetField()) {
169 if (!AliMUONCDB::LoadField()) return;
173 if (!AliMpSegmentation::Instance(kFALSE)) {
174 if (!AliMUONCDB::LoadMapping(kTRUE)) return;
177 // RecoParam for refitting
178 if (!AliMUONESDInterface::GetTracker()) {
179 AliMUONRecoParam* recoParam = AliMUONCDB::LoadRecoParam();
180 if (!recoParam) return;
181 AliMUONESDInterface::ResetTracker(recoParam);
187 //________________________________________________________________________
188 void AliAnalysisTaskMuonTrackingEff::UserCreateOutputObjects()
190 /// Define efficiency histograms
192 fDetEltTDHistList = new TList();
193 fDetEltTDHistList->SetOwner();
194 fDetEltTTHistList = new TList();
195 fDetEltTTHistList->SetOwner();
196 fDetEltSDHistList = new TList();
197 fDetEltSDHistList->SetOwner();
198 fChamberTDHistList = new TList();
199 fChamberTDHistList->SetOwner();
200 fChamberTTHistList = new TList();
201 fChamberTTHistList->SetOwner();
202 fChamberSDHistList = new TList();
203 fChamberSDHistList->SetOwner();
204 fExtraHistList = new TList();
205 fExtraHistList->SetOwner();
209 TString histName, histTitle;
212 Int_t nCentBins = 22;
213 Double_t centRange[2] = {-5., 105.};
215 // prepare binning for THnSparse
221 const Int_t nDims = 5;
222 Int_t nBins[nDims] = {0, nCentBins, 20, 15, 2};
223 Double_t xMin[nDims] = {0., centRange[0], 0., -4., -2.};
224 Double_t xMax[nDims] = {0., centRange[1], 20., -2.5, 2.};
226 // global index of DE in the lists
229 for (Int_t iCh = 0; iCh < 10; iCh++)
231 // histograms per chamber
232 nBins[0] = fgkNbrOfDetectionElt[iCh];
233 xMin[0] = 0.; xMax[0] = static_cast<Double_t>(fgkNbrOfDetectionElt[iCh]);
234 histTitle.Form("ChamberNbr %d", iCh+1);
235 histName.Form("TD_ChamberNbr%d", iCh+1);
236 hn = new THnSparseT<TArrayF>(histName, histTitle, nDims, nBins, xMin, xMax);
237 fChamberTDHistList->AddAt(hn, iCh);
238 histName.Form("TT_ChamberNbr%d",iCh+1);
239 hn = new THnSparseT<TArrayF>(histName, histTitle, nDims, nBins, xMin, xMax);
240 fChamberTTHistList->AddAt(hn, iCh);
241 histName.Form("SD_ChamberNbr%d", iCh+1);
242 hn = new THnSparseT<TArrayF>(histName, histTitle, nDims, nBins, xMin, xMax);
243 fChamberSDHistList->AddAt(hn, iCh);
246 for (Int_t iDE = 0; iDE < fgkNbrOfDetectionElt[iCh]; iDE++)
248 Int_t deId = FromLocalId2DetElt(iCh, iDE);
249 histTitle.Form("detEltNbr %d",deId);
252 histName.Form("TD_detEltNbr%d",deId);
253 h3 = new TH3F(histName, histTitle, 12, -10.0 , 110.0, 12, -10.0, 110.0, nCentBins, centRange[0], centRange[1]);
254 fDetEltTDHistList->AddAt(h3, iDEGlobal);
255 histName.Form("TT_detEltNbr%d",deId);
256 h3 = new TH3F(histName, histTitle, 12, -10.0 , 110.0, 12, -10.0, 110.0, nCentBins, centRange[0], centRange[1]);
257 fDetEltTTHistList->AddAt(h3, iDEGlobal);
258 histName.Form("SD_detEltNbr%d",deId);
259 h3 = new TH3F(histName, histTitle, 12, -10.0 , 110.0, 12, -10.0, 110.0, nCentBins, centRange[0], centRange[1]);
260 fDetEltSDHistList->AddAt(h3, iDEGlobal);
264 histName.Form("TD_detEltNbr%d",deId);
265 h3 = new TH3F(histName, histTitle, 28, -140.0, 140.0, 8, -40.0, 40.0, nCentBins, centRange[0], centRange[1]);
266 fDetEltTDHistList->AddAt(h3, iDEGlobal);
267 histName.Form("TT_detEltNbr%d",deId);
268 h3 = new TH3F(histName, histTitle, 28, -140.0, 140.0, 8, -40.0, 40.0, nCentBins, centRange[0], centRange[1]);
269 fDetEltTTHistList->AddAt(h3, iDEGlobal);
270 histName.Form("SD_detEltNbr%d",deId);
271 h3 = new TH3F(histName, histTitle, 28, -140.0, 140.0, 8, -40.0, 40.0, nCentBins, centRange[0], centRange[1]);
272 fDetEltSDHistList->AddAt(h3, iDEGlobal);
278 // global histograms per chamber
280 xMin[0] = 0.5; xMax[0] = 10.5;
281 hn = new THnSparseT<TArrayF>("TD_Chambers 11", "Chambers 11", nDims, nBins, xMin, xMax);
282 fChamberTDHistList->AddAt(hn, 10);
283 hn = new THnSparseT<TArrayF>("TT_Chambers 11", "Chambers 11", nDims, nBins, xMin, xMax);
284 fChamberTTHistList->AddAt(hn, 10);
285 hn = new THnSparseT<TArrayF>("SD_Chambers 11", "Chambers 11", nDims, nBins, xMin, xMax);
286 fChamberSDHistList->AddAt(hn, 10);
289 TH1F *fHistCent = new TH1F("fHistCent", "centrality distribution", nCentBins, centRange[0], centRange[1]);
290 fExtraHistList->AddAt(fHistCent,0);
291 TH1F *fHistPt = new TH1F("fHistPt", "pt distribution", 250, 0., 50.);
292 fExtraHistList->AddAt(fHistPt,1);
293 TH1F *fHistY = new TH1F("fHistY", "y distribution", 60, -4., -2.5);
294 fExtraHistList->AddAt(fHistY,2);
295 TH1F *fHistTheta = new TH1F("fHistTheta", "theta distribution", 120, 2.8, 3.2);
296 fExtraHistList->AddAt(fHistTheta,3);
297 TH1F *fHistP = new TH1F("fHistP", "momentum distribution", 250, 0., 500.);
298 fExtraHistList->AddAt(fHistP,4);
300 // post the output data at least once
301 PostData(1, fDetEltTDHistList);
302 PostData(2, fDetEltTTHistList);
303 PostData(3, fDetEltSDHistList);
304 PostData(4, fChamberTDHistList);
305 PostData(5, fChamberTTHistList);
306 PostData(6, fChamberSDHistList);
307 PostData(7, fExtraHistList);
310 //________________________________________________________________________
311 void AliAnalysisTaskMuonTrackingEff::UserExec(Option_t *)
315 // check the OCDB has been loaded properly
316 if (!fOCDBLoaded) return;
318 // get the current event
319 AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
322 // get the centrality
323 fCurrentCentrality = esd->GetCentrality()->GetCentralityPercentileUnchecked("V0M");
324 static_cast<TH1F*>(fExtraHistList->At(0))->Fill(fCurrentCentrality);
328 Int_t nTracks = (Int_t)esd->GetNumberOfMuonTracks();
329 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++)
331 AliESDMuonTrack* esdTrack = esd->GetMuonTrack(iTrack);
333 if(!esdTrack->ContainTrackerData()) continue;
335 if(fMatchTrig && !esdTrack->ContainTriggerData()) continue;
337 Double_t thetaTrackAbsEnd = TMath::ATan(esdTrack->GetRAtAbsorberEnd()/505.) * TMath::RadToDeg();
338 Double_t eta = esdTrack->Eta();
339 if(fApplyAccCut && !(thetaTrackAbsEnd >= 2. && thetaTrackAbsEnd <= 10. && eta >= -4. && eta <= -2.5)) continue;
342 const AliVVertex* primaryVertex = esd->GetPrimaryVertexSPD();
343 TVector3 trackDcaAtVz(esdTrack->GetNonBendingCoorAtDCA(), esdTrack->GetBendingCoorAtDCA(), primaryVertex->GetZ());
344 TVector3 vertex(primaryVertex->GetX(), primaryVertex->GetY(), primaryVertex->GetZ());
345 TVector3 meanDca(-0.46, -0.92, 0.); // LHC10h1
346 TVector3 dcaAtVz = trackDcaAtVz - vertex - meanDca;
347 Double_t correctedDca = dcaAtVz.Mag(); // it should also be equal to dcaAtVz.Pt().
348 Double_t pMean = 0.5 * (esdTrack->P() + esdTrack->PUncorrected());
349 Double_t cutVariable = pMean * correctedDca;
350 Double_t cutValue = (thetaTrackAbsEnd > 3.) ? 63. : 120.;
351 cutValue = TMath::Sqrt(cutValue*cutValue + 0.4*0.4*esdTrack->P()*esdTrack->P());
352 if ( cutVariable > fPDCACut*cutValue ) continue;
355 if (fChi2Cut > 0. && esdTrack->GetNormalizedChi2() > fChi2Cut) continue;
357 if (fPtCut > 0. && esdTrack->Pt() < fPtCut) continue;
359 if (fUseMCLabel && esdTrack->GetLabel() < 0) continue;
361 fCurrentTrack = esdTrack;
362 static_cast<TH1F*>(fExtraHistList->At(1))->Fill(esdTrack->Pt());
363 static_cast<TH1F*>(fExtraHistList->At(2))->Fill(esdTrack->Y());
364 static_cast<TH1F*>(fExtraHistList->At(3))->Fill(esdTrack->Theta());
365 static_cast<TH1F*>(fExtraHistList->At(4))->Fill(esdTrack->P());
367 AliMUONESDInterface::ESDToMUON(*esdTrack, track);
369 TrackParamLoop(track.GetTrackParamAtCluster());
372 // post the output data:
373 PostData(1, fDetEltTDHistList);
374 PostData(2, fDetEltTTHistList);
375 PostData(3, fDetEltSDHistList);
376 PostData(4, fChamberTDHistList);
377 PostData(5, fChamberTTHistList);
378 PostData(6, fChamberSDHistList);
379 PostData(7, fExtraHistList);
382 //________________________________________________________________________
383 void AliAnalysisTaskMuonTrackingEff::Terminate(Option_t *)
388 //________________________________________________________________________
389 void AliAnalysisTaskMuonTrackingEff::TrackParamLoop(const TObjArray* trackParams)
391 /// Loop on all the track params and fill the histos
393 Bool_t trackFilter[10];
394 memset(trackFilter, kFALSE, 10*sizeof(Bool_t));
395 Bool_t chamberResponse[10];
396 memset(chamberResponse, kFALSE, 10*sizeof(Bool_t));
398 // check if the chamber responds
399 Int_t nTrackParams = (Int_t) trackParams->GetEntriesFast();
400 for (Int_t iTrackParam = 0; iTrackParam < nTrackParams; ++iTrackParam)
402 Int_t chamberId = static_cast<AliMUONTrackParam*>(trackParams->UncheckedAt(iTrackParam))->GetClusterPtr()->GetChamberId();
403 trackFilter[chamberId] = kTRUE;
404 chamberResponse[chamberId] = kTRUE;
407 // To make sure the calculation of the efficiency of a given chamber (DE) is not biased by the tracking algorithm
408 // we must make sure the track would have been reconstructed whatever this chamber (DE) has responded or not.
409 // If the track is valid for a given chamber, the following code set trackFilter[chamberId] to kTRUE.
410 for (Int_t station = 0; station < 4; ++station)
413 Int_t ch1 = 2*station;
414 Int_t ch2 = 2*station + 1;
415 Int_t ch3 = 2*station + 2;
416 Int_t ch4 = 2*station + 3;
419 filter = trackFilter[ch1];
420 trackFilter[ch1] = trackFilter[ch2];
421 trackFilter[ch2] = filter;
425 if (chamberResponse[ch3] && chamberResponse[ch4])
427 filter = trackFilter[ch1];
428 trackFilter[ch1] = trackFilter[ch2];
429 trackFilter[ch2] = filter;
433 trackFilter[ch1] = kFALSE;
434 trackFilter[ch2] = kFALSE;
437 if (chamberResponse[ch1] && chamberResponse[ch2])
439 filter = trackFilter[ch3];
440 trackFilter[ch3] = trackFilter[ch4];
441 trackFilter[ch4] = filter;
445 trackFilter[ch3] = kFALSE;
446 trackFilter[ch4] = kFALSE;
451 // loop over track parameters
452 Int_t oldChamber = -1;
453 for (Int_t iTrackParam = 0; iTrackParam < nTrackParams; ++iTrackParam)
455 AliMUONTrackParam* trackParam = static_cast<AliMUONTrackParam*>(trackParams->UncheckedAt(iTrackParam));
456 AliMUONVCluster* cluster = trackParam->GetClusterPtr();
458 Int_t newChamber = cluster->GetChamberId();
460 Int_t detElt = cluster->GetDetElemId();
462 ///track position in the global coordinate system
463 Double_t posXG = trackParam->GetNonBendingCoor();
464 Double_t posYG = trackParam->GetBendingCoor();
465 Double_t posZG = trackParam->GetZ();
467 ///track position in the coordinate system of the DE
468 Double_t posXL, posYL, posZL;
469 fTransformer->Global2Local(detElt, posXG, posYG, posZG, posXL, posYL, posZL);
471 // fill histograms if the track is valid for this chamber
472 if(trackFilter[newChamber])
475 // fill histograms of the cluster positions on the detection element of the TRACKS DETECTED (TD)
476 FillTDHistos(newChamber, detElt, posXL, posYL);
478 // fill histograms of the cluster positions on the detection element of ALL THE TRACKS (TT)
479 FillTTHistos(newChamber, detElt, posXL, posYL);
483 FillSDHistos(newChamber, detElt, posXL, posYL);
487 // look for missing cluster(s) if any
488 if (newChamber != oldChamber)
490 if (newChamber > oldChamber + 1)
492 Int_t nbrMissChamber = newChamber - (oldChamber + 1);
494 // find the DE(s) that should have been fired and fill the corresponding histograms
495 FindAndFillMissedDetElt(trackParam, trackFilter, oldChamber+1, nbrMissChamber);
498 // in case the last chamber has not responded
499 if ( iTrackParam == nTrackParams-1 && newChamber != 9) FindAndFillMissedDetElt(trackParam, trackFilter, 9, 1);
502 oldChamber = newChamber;
506 //________________________________________________________________________
507 void AliAnalysisTaskMuonTrackingEff::FindAndFillMissedDetElt(const AliMUONTrackParam* trackParam,
508 const Bool_t* trackFilter,
509 Int_t firstMissCh, Int_t nbrMissCh)
511 /// Find which detection elements should have been hit but were missed, and fill the TT histos appropriately
513 // copy track parameters for extrapolation
514 AliMUONTrackParam extrapTrackParam(*trackParam);
516 // loop over missing chambers
517 for (Int_t iCh = 0; iCh < nbrMissCh; ++iCh)
519 Int_t chamber = firstMissCh + iCh;
521 // skip this chamber if the track is not valid for it
522 if(!trackFilter[chamber]) continue;
524 Int_t nbrOfDetElt = AliMpDEManager::GetNofDEInChamber(chamber, kTRUE);
526 Double_t pos1[6] = {0, 0, 0, 0, 0, 0};
527 Double_t pos2[6] = {0, 0, 0, 0, 0, 0};
528 Double_t posMiss[2] = {0, 0};
530 // track position at the chamber z
531 pos1[2] = AliMUONConstants::DefaultChamberZ(chamber);
532 AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam, pos1[2]);
533 pos1[0] = extrapTrackParam.GetNonBendingCoor();
534 pos1[1] = extrapTrackParam.GetBendingCoor();
536 // track position at the chamber z + dz (where dz = distance between the 2 chamber in the station)
537 pos2[2] = AliMUONConstants::DefaultChamberZ(chamber) + AliMUONConstants::DzCh();
538 AliMUONTrackExtrap::ExtrapToZ(&extrapTrackParam, pos2[2]);
539 pos2[0] = extrapTrackParam.GetNonBendingCoor();
540 pos2[1] = extrapTrackParam.GetBendingCoor();
542 // loop over all the detection element of the chamber
543 for (Int_t iDE = 0; iDE < nbrOfDetElt; iDE++)
545 Int_t deId = (chamber + 1)*fgkOffset + iDE;
547 // track positions (at chamber z and chamber z + dz) in the local coordinate system of the DE
548 fTransformer->Global2Local(deId, pos1[0], pos1[1], pos1[2], pos1[3], pos1[4], pos1[5]);
549 fTransformer->Global2Local(deId, pos2[0], pos2[1], pos2[2], pos2[3], pos2[4], pos2[5]);
551 // track position at z=0 in the local coordinate system of the DE
552 CoordinatesOfMissingCluster(pos1[3], pos1[4], pos1[5], pos2[3], pos2[4], pos2[5], posMiss[0], posMiss[1]);
554 // check if the track cross this DE and fill the corresponding histogram
555 if (CoordinatesInDetElt(deId, posMiss[0], posMiss[1])) FillTTHistos(chamber, deId, posMiss[0], posMiss[1]);
560 //________________________________________________________________________
561 void AliAnalysisTaskMuonTrackingEff::CoordinatesOfMissingCluster(Double_t x1, Double_t y1, Double_t z1,
562 Double_t x2, Double_t y2, Double_t z2,
563 Double_t& x, Double_t& y) const
565 /// Compute the coordinates of the missing cluster. They are defined by the intersection between
566 /// the straigth line joining two extrapolated points (1 and 2) and the detection element plane.
567 /// In the local coordinates, this means Z=0 in the parametric equation of the line.
568 Double_t t = - z1 / (z2 - z1);
569 x = t * (x2 - x1) + x1;
570 y = t * (y2 - y1) + y1;
573 //________________________________________________________________________
574 Bool_t AliAnalysisTaskMuonTrackingEff::CoordinatesInDetElt(Int_t DeId, Double_t x, Double_t y) const
576 /// Return kTRUE if the coordinates are in the Detection Element.
577 /// This is done by checking if a pad correspond to the (x, y) position.
578 const AliMpVSegmentation* seg1 = AliMpSegmentation::Instance()->GetMpSegmentation(DeId, AliMp::kCath0);
579 const AliMpVSegmentation* seg2 = AliMpSegmentation::Instance()->GetMpSegmentation(DeId, AliMp::kCath1);
580 if (!seg1 || !seg2) return kFALSE;
581 AliMpPad pad1 = seg1->PadByPosition(x, y, kFALSE);
582 AliMpPad pad2 = seg2->PadByPosition(x, y, kFALSE);
583 return (pad1.IsValid() && pad2.IsValid());
586 //________________________________________________________________________
587 void AliAnalysisTaskMuonTrackingEff::FillTDHistos(Int_t chamber, Int_t detElt, Double_t posXL, Double_t posYL)
589 /// Fill the histo for detected tracks
590 static_cast<TH3F*>(fDetEltTDHistList->At(FromDetElt2iDet(chamber, detElt)))->Fill(posXL, posYL, fCurrentCentrality);
591 Double_t x[5] = {0., fCurrentCentrality, fCurrentTrack->Pt(), fCurrentTrack->Y(), static_cast<Double_t>(fCurrentTrack->Charge())};
592 x[0] = static_cast<Double_t>(FromDetElt2LocalId(chamber, detElt));
593 static_cast<THnSparse*>(fChamberTDHistList->At(chamber))->Fill(x);
594 x[0] = static_cast<Double_t>(chamber+1);
595 static_cast<THnSparse*>(fChamberTDHistList->At(10))->Fill(x);
598 //________________________________________________________________________
599 void AliAnalysisTaskMuonTrackingEff::FillTTHistos(Int_t chamber, Int_t detElt, Double_t posXL, Double_t posYL)
601 /// Fill the histo for all tracks
602 static_cast<TH3F*>(fDetEltTTHistList->At(FromDetElt2iDet(chamber, detElt)))->Fill(posXL, posYL, fCurrentCentrality);
603 Double_t x[5] = {0., fCurrentCentrality, fCurrentTrack->Pt(), fCurrentTrack->Y(), static_cast<Double_t>(fCurrentTrack->Charge())};
604 x[0] = static_cast<Double_t>(FromDetElt2LocalId(chamber, detElt));
605 static_cast<THnSparse*>(fChamberTTHistList->At(chamber))->Fill(x);
606 x[0] = static_cast<Double_t>(chamber+1);
607 static_cast<THnSparse*>(fChamberTTHistList->At(10))->Fill(x);
610 //________________________________________________________________________
611 void AliAnalysisTaskMuonTrackingEff::FillSDHistos(Int_t chamber, Int_t detElt, Double_t posXL, Double_t posYL)
613 /// Fill the histo for single detected tracks
614 static_cast<TH3F*>(fDetEltSDHistList->At(FromDetElt2iDet(chamber, detElt)))->Fill(posXL, posYL, fCurrentCentrality);
615 Double_t x[5] = {0., fCurrentCentrality, fCurrentTrack->Pt(), fCurrentTrack->Y(), static_cast<Double_t>(fCurrentTrack->Charge())};
616 x[0] = static_cast<Double_t>(FromDetElt2LocalId(chamber, detElt));
617 static_cast<THnSparse*>(fChamberSDHistList->At(chamber))->Fill(x);
618 x[0] = static_cast<Double_t>(chamber+1);
619 static_cast<THnSparse*>(fChamberSDHistList->At(10))->Fill(x);
622 //________________________________________________________________________
623 Int_t AliAnalysisTaskMuonTrackingEff::FromDetElt2iDet(Int_t chamber, Int_t detElt) const
625 /// Connexion between the detection element Id and its position in the list of histograms
626 Int_t iDet = FromDetElt2LocalId(chamber, detElt);
627 for (Int_t iCh = chamber-1; iCh >=0; iCh--) iDet += fgkNbrOfDetectionElt[iCh];
631 //________________________________________________________________________
632 Int_t AliAnalysisTaskMuonTrackingEff::FromDetElt2LocalId(Int_t chamber, Int_t detElt) const
634 /// Connexion between the detection element Id and its number in the chamber
635 return detElt - fgkOffset*(chamber+1);
638 //________________________________________________________________________
639 Int_t AliAnalysisTaskMuonTrackingEff::FromLocalId2DetElt(Int_t chamber, Int_t iDet) const
641 /// Connexion between the number of the detection element in the chamber and its Id
642 return iDet + fgkOffset*(chamber+1);