]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGPP/TRD/AliTRDonlineTrackletQA.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGPP / TRD / AliTRDonlineTrackletQA.cxx
CommitLineData
6419bebb 1// AliTRDonlineTrackletQA implements the standard QA for the TRD
2// on-line tracklets. It reads the tracklets provided by
3// AliTRDonlineTrackletFilter and analyses them. The QA can be run
4// on tracklets from simulation or from real data. In addition it is
5// possible to compare the tracklets from real data to re-simulation
6// based on the ADC data.
7
b220e595 8#include "TFile.h"
9#include "TTree.h"
10#include "TChain.h"
11#include "TH1F.h"
12#include "TH2F.h"
13#include "TProfile.h"
14#include "TCanvas.h"
6911c7f9 15#include "Math/Factory.h"
b220e595 16
17#include "AliAnalysisManager.h"
18#include "AliESDEvent.h"
19#include "AliAODEvent.h"
20#include "AliMCEvent.h"
21#include "AliESDInputHandler.h"
22#include "AliAODHandler.h"
23#include "AliMCEventHandler.h"
24#include "AliLog.h"
25#include "AliESDTrdTrack.h"
6911c7f9 26#include "AliESDTrdTracklet.h"
b220e595 27#include "AliTRDtrackletMCM.h"
6911c7f9 28#include "AliTRDtrackletWord.h"
29#include "AliTRDtrapConfig.h"
30#include "AliTRDmcmSim.h"
31#include "AliTRDcalibDB.h"
32
b220e595 33#include "AliVParticle.h"
6911c7f9 34#include "AliMCParticle.h"
b220e595 35#include "AliTrackReference.h"
36#include "AliMagF.h"
37#include "TGeoGlobalMagField.h"
b220e595 38
39#include "AliTRDonlineTrackletQA.h"
40
41ClassImp(AliTRDonlineTrackletQA)
42
43AliTRDonlineTrackletQA::AliTRDonlineTrackletQA(const char *name) :
44 AliAnalysisTask(name, ""),
45 fESD(0x0),
46 fInputHandler(0x0),
47 fInputEvent(0x0),
48 fOutputAOD(0x0),
49 fMCEvent(0x0),
50 fTrackletsRaw(0x0),
6911c7f9 51 fTrackletsSim(0x0),
b220e595 52 fOutputList(0x0),
53 fHistYpos(0x0),
54 fHistYres(0x0),
55 fHistYresDy(0x0),
56 fHistdY(0x0),
57 fHistdYres(0x0),
58 fHistYresESD(0x0),
59 fHistdYresESD(0x0),
60 fHistCanddY(0x0),
61 fHistFounddY(0x0),
62 fHistTrklPerRef(0x0),
63 fHistdYdYref(0x0),
64 fHistYposRaw(0x0),
65 fHistdYRaw(0x0),
6911c7f9 66 fHistAlphaRaw(0x0),
67 fHistYdYRaw(0x0),
b220e595 68 fHistZrow(0x0),
69 fHistZrowRaw(0x0),
70 fHistPid(0x0),
71 fHistPidRaw(0x0),
6911c7f9 72 fHistPidDiff(0x0),
b220e595 73 fHistYdiff(0x0),
74 fHistdYdiff(0x0),
75 fHistdYdYraw(0x0),
6911c7f9 76 fHistFitYres(0x0),
77 fHistFitDyresEven(0x0),
78 fHistFitDyresOdd(0x0),
79 fHistNoMatchSim(0x0),
80 fHistNoMatchRaw(0x0),
81 fHistResY(0x0),
82 fHistResZ(0x0),
b220e595 83 fTreeTracklets(0x0),
54d34aac 84 fY(0.),
85 fDY(0.),
86 fYdiff(0.),
87 fDYdiff(0.),
88 fQ0(0),
89 fQ1(0),
90 fNHits(0),
b220e595 91 fMinPt(1.),
92 fGeo(new AliTRDgeometry),
93 fNevent(0),
6911c7f9 94 fTrackletTree(0x0)
b220e595 95{
96 // ctor
97
98 DefineInput(0, TChain::Class());
99 DefineInput(1, TTree::Class());
100
6911c7f9 101 DefineOutput(0, TTree::Class());
b220e595 102 DefineOutput(1, TList::Class());
103}
104
105AliTRDonlineTrackletQA::~AliTRDonlineTrackletQA()
106{
107 // dtor
108
109 delete fGeo;
110}
111
54d34aac 112void AliTRDonlineTrackletQA::ConnectInputData(const Option_t * /* option */)
b220e595 113{
54d34aac 114 // connect input data
115
b220e595 116 fInputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
117 if (fInputHandler)
118 fInputEvent = fInputHandler->GetEvent();
119
120 AliMCEventHandler *mcH = (AliMCEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler();
121 if (mcH)
122 fMCEvent = mcH->MCEvent();
123}
124
125void AliTRDonlineTrackletQA::CreateOutputObjects()
126{
54d34aac 127 // create output objects
128
6911c7f9 129 OpenFile(1);
130
b220e595 131 fOutputList = new TList();
6911c7f9 132 fOutputList->SetOwner(kTRUE);
b220e595 133
6911c7f9 134 fHistYpos = new TH1F("ypos",
135 "Tracklet (sim) y-position;y (cm);count",
136 8192/32, -4096*160e-4, 4095*160e-4);
137 fHistYposRaw = new TH1F("ypos_raw",
138 "Tracklet (raw) y-position;y (cm);count",
b220e595 139 8192/32, -4096*160e-4, 4095*160e-4);
6911c7f9 140 fHistYres = new TH1F("yres",
141 "Tracklet (sim) #Deltay;y_{tracklet}-y_{MC} (cm);count",
b220e595 142 8192/32, -4096/32*160e-4, 4095/32*160e-4);
6911c7f9 143 fHistYresDy = new TH2F("yresdy",
144 "Tracklet (sim) #Deltay;y_{tracklet}-y_{MC} (cm);deflection (bin)",
b220e595 145 8192/32, -4096/32*160e-4, 4095/32*160e-4,
146 128, -64.5, 63.5);
6911c7f9 147 fHistYresESD = new TH1F("yresesd",
148 "Tracklet #Deltay;y (cm);count",
b220e595 149 100, -10, 10);
6911c7f9 150 fHistYdiff = new TH1F("ydiff",
151 "Tracklet #Deltay (sim - raw);y_{sim}-y_{raw} (160 #mum);count",
152 200, -100, 100);
b220e595 153 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
6911c7f9 154 fHistYlocal[iLayer] = new TH2F(Form("ylocal_%i", iLayer),
b220e595 155 Form("Tracklet local y, layer %i;y_{MC} (pad width);y_{trkl} (pad width)", iLayer),
156 100, -1, 1, 100, -1, 1);
157 }
158
6911c7f9 159 fHistdY = new TH1F("dy",
160 "deflection (sim);dy (140 #mum)",
b220e595 161 128, -64.5, 63.5);
6911c7f9 162 fHistdYRaw = new TH1F("dy_raw",
163 "deflection (raw);dy (140 #mum)",
b220e595 164 128, -64.5, 63.5);
6911c7f9 165 fHistAlphaRaw = new TH1F("fHistAlphaRaw",
166 "angle w.r.t. to straight line",
167 256, -128.5, 127.5);
168 fHistdYres = new TH1F("dyres",
169 "deflection residual;dy (cm)",
b220e595 170 128, -1., 1.);
6911c7f9 171 fHistdYresESD = new TH1F("dyresesd",
172 "deflection residual;dy (cm)",
b220e595 173 128, -1., 1.);
6911c7f9 174 fHistCanddY = new TH1F("dycand",
175 "deflection;dy (140 #mum)",
b220e595 176 128, -64.5, 63.5);
6911c7f9 177 fHistFounddY = new TH1F("dyfound",
178 "deflection;dy (140 #mum)",
b220e595 179 128, -64.5, 63.5);
6911c7f9 180 fHistdYdiff = new TH1F("dydiff",
181 "deflection #Deltady;dy_{sim}-dy_{raw} (140 #mum)",
182 100, -2., 2.);
183 fHistdYdYraw = new TH2F("dydyraw",
b220e595 184 "deflection from sim. vs raw;dy_{sim} (140 #mum);dy_{raw} (140 #mum)",
185 128, -64.5, 63.5, 128, -64.5, 63.5);
186
6911c7f9 187 fHistTrklPerRef = new TH1F("trklperref",
b220e595 188 "No. of tracklets per track reference;no. of tracklets",
189 10, -0.5, 9.5);
190
6911c7f9 191 fHistdYdYref = new TH2F("dydyref",
b220e595 192 "deflection vs. deflection from track reference;dy_{ref} (140 #mum);dy (140 #mum)",
193 128, -64.5, 63.5, 128, -64.5, 63.5);
194
6911c7f9 195 fHistZrow = new TH1F("zrow",
196 "z-position;pad row",
b220e595 197 16, -0.5, 15.5);
6911c7f9 198 fHistZrowRaw = new TH1F("zrow-raw",
199 "z-position;pad row",
b220e595 200 16, -0.5, 15.5);
201
6911c7f9 202 fHistPid = new TH1F("pid",
203 "pid",
b220e595 204 256, -0.5, 255.5);
6911c7f9 205 fHistPidRaw = new TH1F("pid-raw",
206 "pid",
b220e595 207 256, -0.5, 255.5);
6911c7f9 208 fHistPidDiff = new TH1F("piddiff",
209 "piddiff",
210 256, -127.5, 128.5);
b220e595 211
6911c7f9 212 fHistYdYRaw = new TH2F("ydyraw",
213 "y vs dy (raw tracklets);y (cm);dy (140 #mum)",
214 8192/32, -4096*160e-4, 4095*160e-4,
b220e595 215 128, -64.5, 63.5);
216
6911c7f9 217 fHistFitYres = new TH1F("fityres",
218 "Tracklet #Deltay;y_{tracklet}-y_{track} (cm);count",
219 8192/32, -4096/32*160e-4, 4095/32*160e-4);
220
221 fHistFitDyresEven= new TH1F("fitdyreseven",
222 "Tracklet #Deltady;dy_{tracklet}-dy_{track} (cm);count",
223 100, -2., 2.);
224 fHistFitDyresOdd = new TH1F("fitdyresodd",
225 "Tracklet #Deltady;dy_{tracklet}-dy_{track} (cm);count",
226 100, -2., 2.);
227
228 fHistNoMatchSim = new TH2F("nomatchsim",
229 "Unmatched tracklets from Simulation",
230 8192/32, -4096*160e-4, 4095*160e-4,
231 540, -0.5, 539.5);
232 fHistNoMatchRaw = new TH2F("nomatchraw",
233 "Unmatched tracklets from raw data",
234 8192/32, -4096*160e-4, 4095*160e-4,
235 540, -0.5, 539.5);
236
237 fHistResY = new TH1F("resy",
238 "Residuals to GTU track",
239 100, -10., 10.);
240 fHistResZ = new TH1F("resz",
241 "Residuals to GTU track",
242 100, -20., 20.);
243
b220e595 244 fTreeTracklets = new TTree("trkl", "trkl");
245 fTreeTracklets->Branch("y", &fY);
246 fTreeTracklets->Branch("dy", &fDY);
247 fTreeTracklets->Branch("ydiff", &fYdiff);
248 fTreeTracklets->Branch("dydiff", &fDYdiff);
249 fTreeTracklets->Branch("q0", &fQ0);
250 fTreeTracklets->Branch("q1", &fQ1);
251 fTreeTracklets->Branch("nhits", &fNHits);
252
253 fOutputList->Add(fHistYpos);
254 fOutputList->Add(fHistdY);
255 fOutputList->Add(fHistZrow);
256 fOutputList->Add(fHistPid);
257
258 fOutputList->Add(fHistYres);
259 fOutputList->Add(fHistYresDy);
260 fOutputList->Add(fHistCanddY);
261 fOutputList->Add(fHistFounddY);
262 fOutputList->Add(fHistTrklPerRef);
263 fOutputList->Add(fHistdYres);
264 fOutputList->Add(fHistYresESD);
265 fOutputList->Add(fHistdYresESD);
266 fOutputList->Add(fHistdYdYref);
267
6911c7f9 268 for (Int_t iLayer = 0; iLayer < 6; iLayer++)
b220e595 269 fOutputList->Add(fHistYlocal[iLayer]);
270
271 fOutputList->Add(fHistYposRaw);
272 fOutputList->Add(fHistdYRaw);
6911c7f9 273 fOutputList->Add(fHistAlphaRaw);
b220e595 274 fOutputList->Add(fHistZrowRaw);
275 fOutputList->Add(fHistPidRaw);
276 fOutputList->Add(fHistYdYRaw);
277
278 fOutputList->Add(fHistYdiff);
279 fOutputList->Add(fHistdYdiff);
6911c7f9 280 fOutputList->Add(fHistPidDiff);
b220e595 281 fOutputList->Add(fHistdYdYraw);
282
6911c7f9 283 fOutputList->Add(fHistFitYres);
284 fOutputList->Add(fHistFitDyresEven);
285 fOutputList->Add(fHistFitDyresOdd);
286
287 fOutputList->Add(fHistNoMatchSim);
288 fOutputList->Add(fHistNoMatchRaw);
289
290 fOutputList->Add(fHistResY);
291 fOutputList->Add(fHistResZ);
292
b220e595 293 fOutputList->Add(fTreeTracklets);
6911c7f9 294
295 PostData(1, fOutputList);
b220e595 296}
297
54d34aac 298void AliTRDonlineTrackletQA::Exec(const Option_t * /* option */)
b220e595 299{
300 // execute this for each event
301
302 // connect input data
303 fTrackletTree = (TTree*) GetInputData(1);
6911c7f9 304 if (fTrackletTree) {
305 fTrackletTree->SetBranchAddress("tracklets_sim", &fTrackletsSim);
306 fTrackletTree->SetBranchAddress("tracklets_raw", &fTrackletsRaw);
307 fTrackletTree->GetEntry(fTrackletTree->GetEntriesFast()-1);
308 }
b220e595 309
6911c7f9 310 fESD = dynamic_cast<AliESDEvent*> (fInputEvent);
311
312 TList trackletsSim[540];
313 TList trackletsRaw[540];
b220e595 314
315 // prepare raw tracklets for comparison
316 Int_t detRaw;
6911c7f9 317 Int_t robRaw;
b220e595 318 Int_t mcmRaw;
6911c7f9 319 Int_t yRaw;
b220e595 320 Int_t dyRaw;
321 TTree trklRaw("raw tracklets", "raw tracklets");
322 trklRaw.Branch("det", &detRaw);
323 trklRaw.Branch("rob", &robRaw);
324 trklRaw.Branch("mcm", &mcmRaw);
325 trklRaw.Branch("y", &yRaw);
326 trklRaw.Branch("dy", &dyRaw);
327 trklRaw.SetDirectory(0x0);
328 // prepare simulated tracklets for comparison
329 Int_t detSim;
6911c7f9 330 Int_t robSim;
b220e595 331 Int_t mcmSim;
6911c7f9 332 Int_t ySim;
b220e595 333 Int_t dySim;
334 TTree trklSim("sim tracklets", "sim tracklets");
335 trklSim.Branch("det", &detSim);
336 trklSim.Branch("rob", &robSim);
337 trklSim.Branch("mcm", &mcmSim);
338 trklSim.Branch("y", &ySim);
339 trklSim.Branch("dy", &dySim);
340 trklSim.SetDirectory(0x0);
341
6911c7f9 342 // ----- read ESD tracklets -----
343 if (fESD) {
344 Int_t nTracklets = fESD->GetNumberOfTrdTracklets();
345 Int_t nTrackletsSim = 0;
346 for (Int_t iTracklet = 0; iTracklet < nTracklets; iTracklet++) {
347 AliESDTrdTracklet *tracklet = fESD->GetTrdTracklet(iTracklet);
348
349 if (tracklet) {
350 // check for tracklet from simulation (label >= -1)
351 if (tracklet->GetLabel() >= -1) {
352 fHistYpos->Fill(tracklet->GetLocalY());
353 fHistdY->Fill(tracklet->GetBinDy());
354 fHistZrow->Fill(tracklet->GetBinZ());
355 fHistPid->Fill(tracklet->GetPID());
356
357 detSim = tracklet->GetDetector();
358 robSim = tracklet->GetROB();
359 mcmSim = tracklet->GetMCM();
360 ySim = tracklet->GetBinY();
361 dySim = tracklet->GetBinDy();
362 trklSim.Fill();
363
364 AliDebug(1, Form("trkl sim 0x%08x in %4i (ESD)", tracklet->GetTrackletWord(), tracklet->GetHCId()));
365 ++nTrackletsSim;
366
367 if (tracklet->GetDetector() > -1) {
368 trackletsSim[tracklet->GetDetector()].Add(tracklet);
369 }
370
371 PlotMC(tracklet);
372 PlotESD(tracklet);
373 }
374 // otherwise raw tracklet
375 else {
376 if (tracklet->GetDetector() > -1)
377 trackletsRaw[tracklet->GetDetector()].Add(tracklet);
378
379 fHistYposRaw->Fill(tracklet->GetLocalY());
380 fHistdYRaw->Fill(tracklet->GetBinDy());
381 Float_t alpha = tracklet->GetBinDy() - 3./140.e-4*tracklet->GetLocalY()/GetX(tracklet);
382 fHistAlphaRaw->Fill(alpha);
383 fHistZrowRaw->Fill(tracklet->GetBinZ());
384 fHistPidRaw->Fill(tracklet->GetPID());
385 fHistYdYRaw->Fill(tracklet->GetLocalY(), tracklet->GetBinDy());
386
387 detRaw = tracklet->GetDetector();
388 robRaw = tracklet->GetROB();
389 mcmRaw = tracklet->GetMCM();
390 yRaw = tracklet->GetBinY();
391 dyRaw = tracklet->GetBinDy();
392 trklRaw.Fill();
393 }
394 }
395 }
396 AliDebug(1, Form("no. of simulated tracklets in ESDs: %i", nTrackletsSim));
397 }
398
399 // ----- simulated tracklet from TRD.Tracklets.root -----
b220e595 400 if (fTrackletsSim) {
6911c7f9 401 AliDebug(1, Form("no. of simulated tracklets in TRD.Tracklets.root: %i", fTrackletsSim->GetEntries()));
b220e595 402 for (Int_t iTracklet = 0; iTracklet < fTrackletsSim->GetEntries(); iTracklet++) {
6911c7f9 403 AliTRDtrackletMCM *trkl = (AliTRDtrackletMCM*) (*fTrackletsSim)[iTracklet];
404 AliDebug(1, Form("trkl sim 0x%08x in %4i (TRD)", trkl->GetTrackletWord(), trkl->GetHCId()));
b220e595 405 }
406 }
407
6911c7f9 408 // ----- match simulated and raw tracklets (detector-wise) -----
409 for (Int_t iDetector = 0; iDetector < 540; iDetector++) {
410 for (Int_t iTracklet = 0; iTracklet < trackletsRaw[iDetector].GetEntries(); iTracklet++) {
411 AliESDTrdTracklet *trackletRaw = (AliESDTrdTracklet*) trackletsRaw[iDetector].At(iTracklet);
412
413 // search for the best match
414 AliESDTrdTracklet *trackletMatch = 0x0;
415 Int_t dyMax = 20;
416 Int_t yMax = 100;
417 for (Int_t iTrackletSim = 0; iTrackletSim < trackletsSim[iDetector].GetEntries(); iTrackletSim++) {
418 AliESDTrdTracklet *trackletSim = (AliESDTrdTracklet*) trackletsSim[iDetector].At(iTrackletSim);
419 AliDebug(10, Form("comparing 0x%08x with 0x%08x, det: %i %i",
420 trackletRaw->GetTrackletWord(), trackletSim->GetTrackletWord(),
421 trackletRaw->GetDetector(), trackletSim->GetDetector() ));
422 // require equal z
423 if (trackletRaw->GetBinZ() != trackletSim->GetBinZ())
424 continue;
425 // match on y
426 if (TMath::Abs(trackletRaw->GetBinY() - trackletSim->GetBinY()) > yMax)
427 continue;
b220e595 428
6911c7f9 429 trackletMatch = trackletSim;
430 yMax = TMath::Abs(trackletRaw->GetBinY() - trackletSim->GetBinY());
431 dyMax = TMath::Abs(trackletRaw->GetBinDy() - trackletSim->GetBinDy());
432 }
433 if (trackletMatch) {
434 fHistYdiff->Fill(trackletRaw->GetBinY() - trackletMatch->GetBinY());
435 fHistdYdiff->Fill(3.*trackletRaw->GetDyDx() - 3.*trackletMatch->GetDyDx());
436 fHistdYdYraw->Fill(trackletMatch->GetBinDy(), trackletRaw->GetBinDy());
437 fHistPidDiff->Fill(trackletRaw->GetPID() - trackletMatch->GetPID());
438 trackletsSim[iDetector].Remove(trackletMatch);
439 }
440 else {
441 // store unmatched raw tracklets
442 fHistNoMatchRaw->Fill(trackletRaw->GetLocalY(), trackletRaw->GetDetector());
443 // AliInfo(Form("unmatched raw tracklet: 0x%08x", trackletRaw->GetTrackletWord()));
444 }
445 }
446 // store the unmatched sim tracklets
447 for (Int_t iTracklet = 0; iTracklet < trackletsSim[iDetector].GetEntries(); iTracklet++) {
448 AliESDTrdTracklet *tracklet = (AliESDTrdTracklet*) trackletsSim[iDetector].At(iTracklet);
449 fHistNoMatchSim->Fill(tracklet->GetLocalY(), tracklet->GetDetector());
450 // AliInfo(Form("unmatched simulated tracklet: 0x%08x", tracklet->GetTrackletWord()));
b220e595 451 }
452 }
453
b220e595 454 // ----- MC tracks and track references -----
455 // determine tracklet efficiency
456 if (fMCEvent) {
457 Int_t nTracksMC = fMCEvent->GetNumberOfTracks();
6911c7f9 458 // AliInfo(Form("no. of MC tracks: %i", nTracksMC));
b220e595 459 for (Int_t iTrack = 0; iTrack < nTracksMC; iTrack++) {
460 // we want primaries
461 if (!fMCEvent->IsPhysicalPrimary(iTrack))
462 continue;
463
464 AliMCParticle *mcpart = (AliMCParticle*) fMCEvent->GetTrack(iTrack);
465
466 // don't look at tracks with too low pt
467 if (TMath::Abs(mcpart->Pt()) < fMinPt)
468 continue;
469
470 // look for two track references in a chamber
471 Int_t nTrackRefs = mcpart->GetNumberOfTrackReferences();
472 AliTrackReference *tr[2] = { 0x0 };
473 Int_t nRef = 0;
474
475 for (Int_t iTrackRef = 0; iTrackRef < nTrackRefs; iTrackRef++) {
476 AliTrackReference *trackRef = mcpart->GetTrackReference(iTrackRef);
477 if (trackRef->DetectorId() != AliTrackReference::kTRD)
478 continue;
479 if (trackRef->Pt() < fMinPt)
480 continue;
481 Int_t label = trackRef->Label();
482 if (label < 0)
483 continue;
484
485 // first track reference, remember it
486 if (nRef == 0) {
487 tr[nRef] = trackRef;
488 nRef++;
489 continue;
490 }
491 else {
492 // next one is too far away, remember it but forget the previous one
493 if (TMath::Abs(trackRef->LocalX() - tr[0]->LocalX()) > 5.) {
494 tr[0] = trackRef;
495 nRef = 1;
496 continue;
497 }
498 // too close to previous track reference
499 // we don't want it
500 else if (TMath::Abs(trackRef->LocalX() - tr[0]->LocalX()) < .5) {
501 continue;
502 }
503 // then it must be ok, so we take it
504 else {
505 tr[1] = trackRef;
506 nRef++;
507 }
508 }
509
510 // calculation deflection from track references
511 Float_t deflLength = 3. * (tr[1]->LocalY() - tr[0]->LocalY()) / (tr[1]->LocalX() - tr[0]->LocalX());
512 // if it is too large we reject it
513 if (deflLength < 1.) {
514 fHistCanddY->Fill(deflLength/140e-4);
515 }
516 else {
517 nRef = 0;
518 continue;
519 }
6911c7f9 520
b220e595 521 // now search for tracklets belonging to this track reference
e068d4e7 522 Int_t nTrackletsPerRef(0);
523 Int_t defl(0);
6911c7f9 524 TIter nextTracklet(trackletsSim);
525 while (AliESDTrdTracklet *trkl = (AliESDTrdTracklet*) nextTracklet()) {
b220e595 526 // they must have the same label
6911c7f9 527 if (trkl->GetLabel() != label)
b220e595 528 continue;
6911c7f9 529
b220e595 530 // and be close enough in radial position
6911c7f9 531 if (TMath::Abs(trackRef->LocalX() - GetX(trkl)) > 5.)
b220e595 532 continue;
6911c7f9 533
b220e595 534 // if they are close in position we accept it
6911c7f9 535 if ((TMath::Abs(trackRef->LocalY() - trkl->GetLocalY()) < 5.) &&
536 (TMath::Abs(trackRef->Z() - GetZ(trkl)) < 5.)) {
537 defl = trkl->GetDyDx();
b220e595 538 nTrackletsPerRef++;
539 }
540 }
541 fHistTrklPerRef->Fill(nTrackletsPerRef);
542 if (nTrackletsPerRef == 0) {
543 AliInfo(Form("Track ref without assigned tracklet: x=%4.2f, y=%4.2f, z=%4.2f, pt=%4.2f (%i)",
544 trackRef->X(), trackRef->Y(), trackRef->Z(), trackRef->Pt(), trackRef->Label()));
545 }
546 if (nTrackletsPerRef == 1) {
547 fHistdYdYref->Fill(deflLength/140.e-4, defl);
548 fHistFounddY->Fill(deflLength/140.e-4);
549 }
550 nRef = 0;
551 }
552 }
553 }
554
555 // ----- ESD tracks -----
556 if (fESD) {
557 for (Int_t iTrack = 0; iTrack < fESD->GetNumberOfTracks(); iTrack++) {
558 AliESDtrack *esdTrack = fESD->GetTrack(iTrack);
54d34aac 559 AliDebug(1, Form("ESD track pt: %7.2f", esdTrack->Pt()));
b220e595 560 }
561 }
562
6911c7f9 563 // ----- ESD TRD tracks -----
564 if (fESD) {
565 // AliTRDtrapConfig *trapcfg = AliTRDcalibDB::Instance()->GetTrapConfig();
566 AliDebug(1, Form("no. of TRD tracks: %i", fESD->GetNumberOfTrdTracks()));
567 for (Int_t iTrack = 0; iTrack < fESD->GetNumberOfTrdTracks(); iTrack++) {
568 AliESDTrdTrack *trdTrack = fESD->GetTrdTrack(iTrack);
569 AliInfo(Form("TRD track pt: %7.2f", trdTrack->Pt()));
570 // AliTRDtrackOnline track;
571 AliESDTrdTracklet tracklets[6];
572 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
573 AliESDTrdTracklet *trkl = trdTrack->GetTracklet(iLayer);
574 if (trkl) {
575 tracklets[iLayer].SetTrackletWord(trkl->GetTrackletWord());
576 tracklets[iLayer].SetHCId(trkl->GetHCId());
577 // track.AddTracklet(&(tracklets[iLayer])); //???
578 }
579 }
580 // track.AddParametrization(new AliTRDtrackParametrizationCurved());
581 // ROOT::Math::Minimizer *minim = ROOT::Math::Factory::CreateMinimizer("Minuit");
582 // track.Fit(minim);
583 // const TList &paramList = track.GetParametrizations();
584 // TIter param(&paramList, kFALSE);
585 // while (AliTRDtrackParametrization *trackParam = (AliTRDtrackParametrization*) param()) {
586 // trackParam->Print();
587 // AliTRDtrackParametrizationCurved *curvedParam = dynamic_cast<AliTRDtrackParametrizationCurved*> (trackParam);
588 // if (curvedParam) {
589 // // calculate residuals
590 // for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
591 // AliESDTrdTracklet *trkl = trdTrack->GetTracklet(iLayer);
592 // if (trkl) {
593 // AliTRDtrackPosition pos = curvedParam->ExtrapolateToLayer(iLayer);
594 // fHistFitYres->Fill(pos.GetY() - trkl->GetLocalY());
595 // Int_t rob = (trkl->GetBinZ()/4) + (trkl->GetHCId() % 0x1 ? 1 : 0);
596 // Int_t mcm = trkl->GetBinZ() & 0x3;
597 // Float_t deflcorr = ((Int_t) trapcfg->GetDmemUnsigned(AliTRDmcmSim::fgkDmemAddrDeflCorr, trkl->GetDetector(), rob, mcm))
598 // * (0.635 + 0.03 * iLayer)/(256.0) / 32.;
599 // Float_t scale = 1.;
600 // Float_t defl = (3.*trkl->GetDyDx() - deflcorr) * scale + deflcorr;
601 // if (iLayer & 0x1) {
602 // fHistFitDyresOdd->Fill(pos.GetdY() - defl);
603 // }
604 // else {
605 // fHistFitDyresEven->Fill(pos.GetdY() - defl);
606 // }
607 // }
608 // }
609 // }
610 // }
611 }
612 }
613
614 PostData(1, fOutputList);
b220e595 615}
616
617void AliTRDonlineTrackletQA::LocalInit()
618{
619
620}
621
54d34aac 622void AliTRDonlineTrackletQA::Terminate(const Option_t * /* option */)
b220e595 623{
54d34aac 624 // upon terminate
625
b220e595 626 fOutputList = dynamic_cast<TList*> (GetOutputData(1));
627
628 if (!fOutputList) {
629 AliError("No output objects found!");
630 return;
631 }
632
633// fHistYpos = dynamic_cast<TH1F*> (fOutputList->FindObject("ypos"));
634// if (fHistYpos) {
635// TCanvas *c = new TCanvas;
636// c->cd();
637// gPad->SetGridx();
638// gPad->SetGridy();
639// fHistYpos->DrawCopy();
640// }
641}
642
6911c7f9 643void AliTRDonlineTrackletQA::PlotMC(AliESDTrdTracklet *trkl)
b220e595 644{
645 // compare the given tracklet to the MC information,
646 // i.e. track references
6911c7f9 647
b220e595 648 Int_t label = trkl->GetLabel();
649
650 // check for valid label
651 if (label < 0 ) {
6911c7f9 652 AliDebug(1, "MC tracklet has no label");
b220e595 653 return;
654 }
655 if (label >= fMCEvent->GetNumberOfTracks()) {
6911c7f9 656 AliError("MC tracklet has invalid label");
b220e595 657 return;
658 }
659
660 // some cuts on the tracklet
661// if (!fMCEvent->IsPhysicalPrimary(label))
662// return;
663// if (TMath::Abs(trkl->GetdYdX()) > 0.05)
664// return;
665// if (trkl->GetDetector() % 6 != 5)
666// return;
6911c7f9 667
b220e595 668 // get MC particle for this tracklet
669 AliMCParticle *track = (AliMCParticle*) fMCEvent->GetTrack(label);
670
671 // don't look at tracks with too low pt
672 if (TMath::Abs(track->Pt() < fMinPt))
673 return;
674
675 // serach track references corresponding to the current tracklet
676 AliTrackReference *tr[2] = { 0x0 };
677 Int_t nTrackRefs = 0;
678
679 for (Int_t iTrackRef = 0; iTrackRef < track->GetNumberOfTrackReferences(); iTrackRef++) {
680 AliTrackReference *trackRef = track->GetTrackReference(iTrackRef);
681 if (trackRef->DetectorId() != AliTrackReference::kTRD)
682 continue;
683 if (trackRef->Pt() < fMinPt)
684 continue;
6911c7f9 685
686 if (TMath::Abs(GetX(trkl) - trackRef->LocalX()) > 5.)
b220e595 687 continue;
688
6911c7f9 689 tr[nTrackRefs++] = trackRef;
b220e595 690
691 if (nTrackRefs == 2)
692 break;
693 }
694
6911c7f9 695 // if there were exactly 2 track references
b220e595 696 // (otherwise something is strange and we want to look at clean cases)
697 // compare tracklet to them
698 if (nTrackRefs == 2) {
699 // sanity check to be in the same sector
700 if ( TMath::Abs((tr[0]->Alpha()*180./TMath::Pi()-10.)/20. - (trkl->GetDetector()/30)) > .1) {
701 AliError("Track reference in different sector");
702 }
703 // require minimal distance in X and maximum deflection in Y
704 // for the track references
6911c7f9 705 else if ((tr[1]->LocalX() - tr[0]->LocalX()) > 0.1 && TMath::Abs(tr[1]->LocalY() - tr[0]->LocalY()) < 1.) {
b220e595 706 // calculate slope from track references
707 // and check whether it's in the allowed range
708 Float_t slope = 3. * (tr[1]->LocalY() - tr[0]->LocalY()) / (tr[1]->LocalX() - tr[0]->LocalX());
709 if (TMath::Abs(slope) < 64*140e-4) {
6911c7f9 710 AliDebug(1,Form("x1: %f, x0: %f, y1: %f, y0:%f",
b220e595 711 tr[1]->LocalX(), tr[0]->LocalX(), tr[1]->LocalY(), tr[0]->LocalY() ));
712 // calculate y-position scaled to radial position of the tracklet
713 // and consider the tilting angle of the pads
714 // since the tracklets are affected by it
6911c7f9 715 Float_t yMC = (tr[1]->LocalY() + (-0.5+GetX(trkl) - tr[1]->LocalX()) * (tr[1]->LocalY() - tr[0]->LocalY()) / (tr[1]->LocalX() - tr[0]->LocalX()));
716 Float_t yMCtilt = yMC + (TMath::Tan(TMath::Power(-1, (trkl->GetDetector() % 6))*2.*TMath::Pi()/180.) * (tr[1]->Z() - GetZ(trkl)));
717 if (TMath::Abs(trkl->GetLocalY() - yMCtilt) > 10.) {
718 AliError(Form("Deviation too large for tracklet: 0x%08x in det. %i at x = %f, y = %f, z = %f, alpha = %f",
719 trkl->GetTrackletWord(), trkl->GetDetector(), GetX(trkl), trkl->GetLocalY(), GetZ(trkl), tr[0]->Alpha()));
b220e595 720 }
6911c7f9 721 fHistYres->Fill(trkl->GetLocalY() - yMCtilt);
722 fHistYresDy->Fill(trkl->GetLocalY() - yMCtilt, trkl->GetDyDx());
b220e595 723 // what about tilt correction here ???
6911c7f9 724 fHistdYres->Fill(3. * trkl->GetDyDx() -
b220e595 725 3. * (tr[1]->LocalY() - tr[0]->LocalY()) / (tr[1]->LocalX() - tr[0]->LocalX()));
726 // plot position deviation in pad-coordinates
727 // to study the influence of the position LUT
728 Float_t padWidth = fGeo->GetPadPlane(trkl->GetDetector())->GetWidthIPad();
6911c7f9 729 Float_t yMClocal = yMCtilt/padWidth - floor(yMCtilt/padWidth) - padWidth/2.;
b220e595 730 Int_t layer = trkl->GetDetector() % 6;
6911c7f9 731 fHistYlocal[layer]->Fill(yMClocal,
732 trkl->GetLocalY()/padWidth - floor(trkl->GetLocalY()/padWidth) - padWidth/2. - yMClocal);
b220e595 733 // and fill everything to the tree
6911c7f9 734 fYdiff = trkl->GetLocalY() - yMCtilt;
735 fDYdiff = 3. * trkl->GetDyDx() -
b220e595 736 3. * (tr[1]->LocalY() - tr[0]->LocalY()) / (tr[1]->LocalX() - tr[0]->LocalX());
6911c7f9 737 fY = trkl->GetLocalY();
738 fDY = trkl->GetDyDx();
b220e595 739 fTreeTracklets->Fill();
740 // output tracklets with large deviation
741 if (TMath::Abs(fYdiff) > 0.5) {
6911c7f9 742 AliWarning(Form("tracklet: y=%4.2f, dy=%4.2f, ydiff=%4.2f, dydiff=%4.2f, q0=%5d, q1=%5d, nhits=%2d, label=%i",
743 trkl->GetLocalY(), trkl->GetDyDx(), fYdiff, fDYdiff, fQ0, fQ1, fNHits, label));
b220e595 744 }
745 }
746 }
747 }
748}
749
750
6911c7f9 751void AliTRDonlineTrackletQA::PlotESD(AliESDTrdTracklet *trkl)
b220e595 752{
54d34aac 753 // plot comparison to ESD
b220e595 754
6911c7f9 755 Float_t xTrkl = GetX(trkl);
756 Float_t yTrkl = trkl->GetLocalY();
757 Float_t zTrkl = GetZ(trkl);
b220e595 758
759 Float_t alpha = (trkl->GetDetector() / 30) * 20. + 10.;
760 alpha *= TMath::Pi() / 180.;
761
762 AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*> (fInputEvent);
763 if (!esdEvent)
764 return;
765
6911c7f9 766 Float_t mag = 0.5; //((AliMagF*) TGeoGlobalMagField::Instance()->GetField())->SolenoidField();
b220e595 767
768 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) {
769 AliESDtrack *track = esdEvent->GetTrack(iTrack);
770
771 if (!track->GetOuterParam())
772 continue;
773
774 AliExternalTrackParam *param = new AliExternalTrackParam(*(track->GetOuterParam()));
775
6911c7f9 776 AliDebug(10, Form("track %i at x = %f, y = %f",
b220e595 777 iTrack, param->GetX(), param->GetY()));
778 param->Propagate(alpha, xTrkl, mag);
6911c7f9 779 AliDebug(10, Form("after propagating track %i at x = %f, y = %f",
b220e595 780 iTrack, param->GetX(), param->GetY()));
781
782 if ((TMath::Abs(xTrkl - param->GetX()) < 10.) &&
783 (TMath::Abs(yTrkl - param->GetY()) < 5.) &&
784 (TMath::Abs(zTrkl - param->GetZ()) < 10.)) {
6911c7f9 785 AliDebug(2, Form("match of tracklet-track: %i <-> %i",
786 trkl->GetLabel(), track->GetLabel()));
b220e595 787 AliDebug(5, Form("tracklet position: det: %3i x = %f, y = %f, z = %f, alpha = %f",
6911c7f9 788 trkl->GetDetector(), GetX(trkl), trkl->GetLocalY(), GetZ(trkl), alpha));
789 AliDebug(5, Form("after propagating track %i at x = %f, y = %f, z = %f",
b220e595 790 iTrack, param->GetX(), param->GetY(), param->GetZ()));
791
792 fHistYresESD->Fill(yTrkl - param->GetY());
793 }
794
795 delete param;
796 }
797
798}
799
6419bebb 800Int_t AliTRDonlineTrackletQA::GetTrackletsForMC(Int_t /* label */, Int_t /*idx*/ []) const
b220e595 801{
54d34aac 802 // get tracklets for MC label
803 // not implemented
804
b220e595 805 return 0;
806}