Last updates on the right version (1.44).
[u/mrichter/AliRoot.git] / MUON / AliMUONRecoDisplay.cxx
CommitLineData
2c0075de 1/**************************************************************************
2 * Copyright(c) 1998-1999, 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/*
17$Log$
18*/
19
c7ba256d 20//Authors: Mihaela Gheata, Andrei Gheata 09/10/00
2c0075de 21//////////////////////////////////////////////////////////////////////
22// //
23// AliMUONRecoDisplay //
24// //
25// This class subclasses AliDisplay and provides display of //
26// reconstructed tracks with following functionality : //
27// - front/top/side/3D display of MUON reconstructed tracks //
28// as polylines ; //
29// - context menu activated when the main pad is right-clicked //
30// The context menu contains following functions : //
31// * SetDrawHits() - switches on or off Geant hits ; //
32// * CutMomentum() - displays only tracks within Pmin - Pmax //
33// * ListTracks() - prints ID and momentum info. for all //
34// tracks within momentum range Pmin,Pmax ; //
35// * Highlight() - shows only one selected reco. track //
36// and its best matching Geant track; //
37// * UnHighlight() - self explaining; //
38// * RecoEfficiency() - compute reco. efficiency for all events//
39// from galice.root file; also fake track percentage; make //
40// plots for momentum precision //
41// * XYPlot() - make X-Y plots of reconstructed and //
42// generated tracks in all chambers //
43// //
44// Starting : generate and reconstruct events, then use the //
45// MUONrecodisplay.C macro //
46// //
47//////////////////////////////////////////////////////////////////////
c7ba256d 48
49#include <iostream.h>
50#include <AliRun.h>
51#include <TClonesArray.h>
52#include "AliMUONRecoEvent.h"
53#include "AliMUONRecoDisplay.h"
54#include <TROOT.h>
55#include <AliPoints.h>
56#include <TSlider.h>
57#include <TView.h>
58#include <TGeometry.h>
59
60ClassImp(AliMUONRecoDisplay)
61
62//-------------------------------------------------------------------
63AliMUONRecoDisplay::AliMUONRecoDisplay(Int_t nevent)
64 :AliDisplay(750)
65{
66//************ Constructor of the reco. event display**********
67 // get reconstructed event from file
68 fFile = new TFile("tree_reco.root");
69 if (!fFile) {
70 cout << "File tree_reco.root not found\n";
71 gApplication->Terminate(0);
72 }
73 fEvReco = 0;
74 fTree = (TTree *) fFile->Get("TreeRecoEvent");
75 if (!fTree) {
76 cout << "Tree of reconstructed events not found on file. Abort.\n";
77 gApplication->Terminate(0);
78 }
2c0075de 79 fEvent = nevent;
80 fEmpty = kFALSE;
c7ba256d 81 TBranch *branch = fTree->GetBranch("Event");
82 branch->SetAddress(&fEvReco);
83
84 fEvGen = new AliMUONRecoEvent();
85
86 TFile *galice_file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
87 galice_file->cd();
2c0075de 88 if (nevent > gAlice->TreeE()->GetEntries() - 1) {
89 cout << "Event number out of range !\n";
90 gApplication->Terminate(0);
91 }
c7ba256d 92 gAlice->GetEvent(nevent);
93
94
95 fRecoTracks = 0;
96 fGenTracks = 0;
97 fPolyRecoList = 0;
98 fPolyGenList = 0;
99 fHighlited = -1;
100 fMinMomentum = 0;
101 fMaxMomentum = 999;
2c0075de 102 // map this event
c7ba256d 103 MapEvent(nevent);
104}
105
106//-------------------------------------------------------------------
107AliMUONRecoDisplay::~AliMUONRecoDisplay()
108{
109// Destructor of display object
110 if (fPolyRecoList) {
111 fPolyRecoList->Delete();
112 delete fPolyRecoList;
113 }
114 if (fPolyGenList) {
115 fPolyGenList->Delete();
116 delete fPolyGenList;
117 }
118 delete fEvGen;
119}
2c0075de 120//-------------------------------------------------------------------
121Bool_t AliMUONRecoDisplay::Event(Int_t nevent)
122{
123 fEvent = nevent;
124 for (Int_t entry=0; entry<fTree->GetEntries(); entry++) {
125 fTree->GetEntry(entry);
126 if (fEvReco->GetNoEvent() == nevent) return kTRUE;
127 }
128 cout << "Event number " << nevent << " empty\n";
129 fRecoTracks = 0;
130 fPolyRecoList = 0;
131 fGenTracks = 0;
132 fPolyGenList = 0;
133 fHighlited = -1;
134 fMinMomentum = 0;
135 fMaxMomentum = 999;
136
137 return kFALSE;
138}
c7ba256d 139
140//-------------------------------------------------------------------
141void AliMUONRecoDisplay::MapEvent(Int_t nevent)
142{
143// get generated event (nevent) from galice.root and corresponding
144// reconstructed event from tree_reco.root;
145 cout << "mapping event " << nevent << endl;
2c0075de 146 fEvent = nevent;
147 fEmpty = kFALSE;
c7ba256d 148 fFile->cd();
2c0075de 149 // check if the event is not empty and make fEvReco to point to it
150 fEmpty = !Event(nevent);
c7ba256d 151 // just testing
2c0075de 152// if (!fEvReco) {
153// cout << "failed to load reco event ! Correct this.\n";
154// gApplication->Terminate(0);
155// }
156 if (fEmpty) return;
c7ba256d 157 fRecoTracks = fEvReco->TracksPtr();
158 fPolyRecoList = MakePolyLines3D(fRecoTracks);
159 // clear previous event
160 if (fEvGen) fEvGen->Clear();
161 fEvGen->SetNoEvent(nevent);
162 // get list of particles
163 TClonesArray *Particles = gAlice->Particles();
164 // connect MUON module
165 AliDetector *MUON = gAlice->GetDetector("MUON");
166 if (!MUON) {
167 cout << "MUON module not present.\n";
168 gApplication->Terminate(0);
169 }
170 // get the number of generated tracks
171 Int_t ntracks = (Int_t)gAlice->TreeH()->GetEntries();
172 // Fill the fEvGen object
173 AliMUONRecoTrack *gtrack = 0;
174 AliMUONHit *hit = 0;
175 TParticle *particle;
176 Int_t ch;
177 // loop all tracks
178 for (Int_t track=0; track<ntracks; track++) {
179 hit = (AliMUONHit *) MUON->FirstHit(track);
180 if (!hit) continue;
181 particle = (TParticle *) Particles->UncheckedAt(hit->Track());
182 if (IsReconstructible(track) && TMath::Abs(particle->GetPdgCode())==13) {
183 gtrack = fEvGen->AddEmptyTrack();
184 gtrack->SetSign(TMath::Sign((Int_t)1, -particle->GetPdgCode()));
185 // reset hits
186 for (ch=0; ch<10; ch++) gtrack->SetHitPosition(ch,0,0,0);
187 // loop all hits
188 for (AliMUONHit *muonHit=(AliMUONHit*)MUON->FirstHit(track);
189 muonHit;
190 muonHit=(AliMUONHit*)MUON->NextHit()) {
191 ch = muonHit->fChamber - 1;
192 if (ch<0 || ch>9) continue;
193 gtrack->SetHitPosition(ch, muonHit->X(), muonHit->Y(), muonHit->Z());
194 gtrack->SetMomReconstr(particle->Px(), particle->Py(), particle->Pz());
195 gtrack->SetVertexPos(particle->Vz());
196 gtrack->SetChi2r(0.0);
197 }
198 }
199 }
200 fGenTracks = fEvGen->TracksPtr();
201 fPolyGenList = MakePolyLines3D(fGenTracks);
202}
203//-------------------------------------------------------------------
204void AliMUONRecoDisplay::XYPlot()
205{
206// Plot reco. tracks hits in all chambers for current event:
207// - open blue squares : generated muons that were reconstructed accurate
208// - open cyan squares : generated muons that were not reconstructed
209// - filled green circles : reco. tracks (accurate)
210// - filled red squares : fake tracks
2c0075de 211 if (fEmpty) return;
c7ba256d 212 Double_t kMaxRadius[10];
213 kMaxRadius[0] = kMaxRadius[1] = 91.5;
214 kMaxRadius[2] = kMaxRadius[3] = 122.5;
215 kMaxRadius[4] = kMaxRadius[5] = 158.3;
216 kMaxRadius[6] = kMaxRadius[7] = 260.0;
217 kMaxRadius[8] = kMaxRadius[9] = 260.0;
218
219 TH2F *xygen_found[10];
220 TH2F *xygen_lost[10];
221 TH2F *xyreco_good[10];
222 TH2F *xyreco_fake[10];
223 Double_t x,y,r;
224 Int_t matches[500];
225 Int_t index, ch;
226
227 TPad *pad = (TPad*)gROOT->GetSelectedPad();
228 TCanvas *canvas = new TCanvas("xy", "Reconstruction efficiency");
229 canvas->Clear();
230 canvas->cd();
231 canvas->Divide(3,4);
232 canvas->cd(11);
233 gPad->Delete();
234 canvas->cd(12);
235 gPad->Delete();
236 canvas->cd(1);
237
238 // Define histograms for x-y plots
239 for (ch=0; ch<10; ch++) {
240 xygen_found[ch] = new TH2F("xygen_found","",50,-kMaxRadius[ch],kMaxRadius[ch],50,-kMaxRadius[ch],kMaxRadius[ch]);
241 xygen_lost[ch] = new TH2F("xygen_lost","",50,-kMaxRadius[ch],kMaxRadius[ch],50,-kMaxRadius[ch],kMaxRadius[ch]);
242 xyreco_good[ch] = new TH2F("xyreco_good","",50,-kMaxRadius[ch],kMaxRadius[ch],50,-kMaxRadius[ch],kMaxRadius[ch]);
243 xyreco_fake[ch] = new TH2F("xyreco_fake","",50,-kMaxRadius[ch],kMaxRadius[ch],50,-kMaxRadius[ch],kMaxRadius[ch]);
244 }
245 // find list of matching tracks
246 fPrinted = kTRUE; // no need to print
247 for (index=0; index<fRecoTracks->GetEntriesFast(); index++) {
248 matches[index] = GetBestMatch(index);
249 }
250 // and fill them
251 for (index=0; index<fGenTracks->GetEntriesFast(); index++) { // GEANT
252 Bool_t wasreconst = kFALSE;
253 for (Int_t i=0; i<fRecoTracks->GetEntriesFast(); i++) {
254 if (matches[i] == index) {
255 wasreconst = kTRUE;
256 break;
257 }
258 }
259 AliMUONRecoTrack *current = (AliMUONRecoTrack*)fGenTracks->UncheckedAt(index);
260 for (ch=0; ch<10; ch++) {
261 x = current->GetPosX(ch);
262 y = current->GetPosY(ch);
263 r = TMath::Sqrt(x*x +y*y);
264 if (r >= 10) {
265 if (wasreconst) {
266 xygen_found[ch]->Fill(x,y);
267 } else {xygen_lost[ch]->Fill(x,y);}
268 }
269 }
270 }
271 for (index=0; index<fRecoTracks->GetEntriesFast(); index++) { // reco
272 AliMUONRecoTrack *current = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(index);
273 for (ch=0; ch<10; ch++) {
274 x = current->GetPosX(ch);
275 y = current->GetPosY(ch);
276 r = TMath::Sqrt(x*x +y*y);
277 if (r >= 10) {
278 if (matches[index] >= 0) {
279 xyreco_good[ch]->Fill(x,y);
280 } else {xyreco_fake[ch]->Fill(x,y);}
281 }
282 }
283 }
284 // finally plot them
285 for (ch=0; ch<10; ch++) {
286 canvas->cd(ch+1);
287 xygen_found[ch]->SetMarkerColor(kBlue);
288 xygen_found[ch]->SetMarkerStyle(4);
289 xygen_found[ch]->SetMarkerSize(0.5);
290 xygen_found[ch]->SetStats(kFALSE);
291 xygen_found[ch]->Draw();
292 xygen_lost[ch]->SetMarkerColor(kCyan);
293 xygen_lost[ch]->SetMarkerStyle(4);
294 xygen_lost[ch]->SetMarkerSize(0.5);
295 xygen_lost[ch]->SetStats(kFALSE);
296 xygen_lost[ch]->Draw("SAME");
297 xyreco_good[ch]->SetMarkerColor(kGreen);
298 xyreco_good[ch]->SetMarkerStyle(20);
299 xyreco_good[ch]->SetMarkerSize(0.4);
300 xyreco_good[ch]->SetStats(kFALSE);
301 xyreco_good[ch]->Draw("SAME");
302 xyreco_fake[ch]->SetMarkerColor(kRed);
303 xyreco_fake[ch]->SetMarkerStyle(20);
304 xyreco_fake[ch]->SetMarkerSize(0.5);
305 xyreco_fake[ch]->SetStats(kFALSE);
306 xyreco_fake[ch]->Draw("SAME");
307 }
308 canvas->SetTitle("y vs. x for simulated and reconstructed tracks");
309 pad->cd();
310}
311
312//-------------------------------------------------------------------
313void AliMUONRecoDisplay::RecoEfficiency(Int_t first, Int_t last)
314{
315// Loop selected reconstructed events, compute efficiency as total number of
316// reconstructed tracks (accurate) over total number of generated
317// reconstructible muons. Make histogram for momentum precision and profile
318// of mean momentum precision versus total momentum (generated)
2c0075de 319 Int_t nevents = (Int_t)gAlice->TreeE()->GetEntries();
c7ba256d 320 if (last > nevents) last = nevents - 1;
321 if (first < 0) first = 0;
322 nevents = last - first + 1;
323 Int_t track;
324 Float_t efficiency;
325 Float_t generated=0, found=0, fake=0; // number of generated/found/fake tracks
326 Double_t pgen, preco, dP; // dP=preco-pgen
327 AliMUONRecoTrack *rtrack, *gtrack; // generated/reco. current tracks
328 fPrinted = kTRUE; // no need to print
329 Int_t currentEvent = gAlice->GetHeader()->GetEvent();
330
331 TH1F *pReso = new TH1F("Momentum precision", "dP = Prec - Pgen", 100, -5.0, 5.0);
332 pReso->SetXTitle("dP [GeV/c]");
333 pReso->SetYTitle("dN/dP");
334
335 TProfile *dPp = new TProfile("", "dP vs. P", 50, 0, 100, 0, 5);
336 dPp->SetXTitle("P [GeV/c]");
337 dPp->SetYTitle("<dP> [GeV/c]");
338
339 TPad *pad = (TPad*)gROOT->GetSelectedPad();
340 TCanvas *canvas = new TCanvas("c", "Reconstruction efficiency");
341 canvas->Divide(1,2);
342 canvas->cd(1);
343
344 // loop events
345 for (Int_t event=first; event<=last; event++) {
346 // get the reco. & gen. events
347 TFile *galice_file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
348 galice_file->cd();
349 gAlice->GetEvent(event);
350 MapEvent(event);
2c0075de 351 if (fEmpty) {
352 // skip empty events
353 fEmpty = kFALSE;
354 continue;
355 }
c7ba256d 356 generated += fGenTracks->GetEntriesFast();
357 // loop reco. tracks
358 for (track=0; track<fRecoTracks->GetEntriesFast(); track++) {
359 rtrack = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(track);
360 preco = rtrack->P();
361 // find correspondent generated track
362 Int_t ind = GetBestMatch(track);
363 if (ind >=0) {
364 gtrack = (AliMUONRecoTrack*)fGenTracks->UncheckedAt(ind);
365 pgen = gtrack->P();
366 found += 1;
367 dP = preco - pgen;
368 pReso->Fill(dP);
369 dPp->Fill(pgen, TMath::Abs(dP));
370 } else {
371 fake += 1;
372 }
373 }
374 }
375 efficiency = found/generated;
376 cout << "=================================================================\n";
377 cout << "|| Reconstruction efficiency and momentum precision ||\n";
378 cout << "=================================================================\n";
379 cout << "Number of events processed " << nevents << endl;
380 cout << "Total number of reconstructible muons : " << (Int_t)generated << endl;
381 cout << "RECONSTRUCTION EFFICIENCY : " << efficiency*100 << " %" << endl;
382 cout << "Fake track rate : " << 100*fake/generated << " %" << endl;
383 cout << "Momentum precision fit : \n" << endl;
384 pReso->Fit("gaus");
385 cout << "=================================================================\n";
386
387 canvas->cd(1);
388 pReso->Draw();
389 canvas->cd(2);
390 dPp->Draw("HIST");
391 pad->cd();
392 ShowNextEvent(currentEvent-last);
393}
394
395//-------------------------------------------------------------------
396void AliMUONRecoDisplay::ShowNextEvent(Int_t delta)
397{
398// overwritten from AliDisplay in order to get also mapping of next event
399 TFile *galice_file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root");
400 galice_file->cd();
401 if (delta) {
402 gAlice->Clear();
403 Int_t currentEvent = gAlice->GetHeader()->GetEvent();
404 Int_t newEvent = currentEvent + delta;
2c0075de 405 if (newEvent<0 || newEvent>(gAlice->TreeE()->GetEntries() - 1)) return;
c7ba256d 406 Int_t nparticles = gAlice->GetEvent(newEvent);
407 cout << "Event : " << newEvent << " with " << nparticles << " particles\n";
408 if (!gAlice->TreeH()) return;
409 MapEvent(newEvent);
410 fHighlited = -1;
411 }
412 LoadPoints();
413 fPad->cd();
414 Draw();
415 if (gROOT->GetListOfCanvases()->FindObject("xy")) XYPlot();
416}
417//-------------------------------------------------------------------
418Bool_t AliMUONRecoDisplay::IsReconstructible(Int_t track)
419{
420// true if at least three hits in first 2 stations, 3 in last 2 stations
421// and one in station 3
2c0075de 422 if (fEmpty) return kFALSE;
c7ba256d 423 AliDetector *MUON = gAlice->GetDetector("MUON");
424 Bool_t chHit[10];
425 Int_t ch;
426 for (ch=0; ch<10; ch++) chHit[ch] = kFALSE;
427 //loop hits
428 for (AliMUONHit *muonHit=(AliMUONHit*)MUON->FirstHit(track);
429 muonHit;
430 muonHit=(AliMUONHit*)MUON->NextHit()) {
431 ch = muonHit->fChamber - 1;
432 if (ch<0 || ch>9) continue;
433 chHit[ch] = kTRUE;
434 }
435 Int_t nhits = 0;
436 for (ch=0; ch<4; ch++) nhits += (chHit[ch])?1:0;
437 if (nhits < 3) return kFALSE;
438 nhits = 0;
439 for (ch=4; ch<6; ch++) nhits+= (chHit[ch])?1:0;
440 if (nhits < 1) return kFALSE;
441
442 for (ch=7; ch<10; ch++) nhits+= (chHit[ch])?1:0;
443 if (nhits < 3) return kFALSE;
444 return kTRUE;
445}
446
447//-------------------------------------------------------------------
448void AliMUONRecoDisplay::DrawView(Float_t theta, Float_t phi, Float_t psi)
449{
450// ovewritten from base class to change the range for MUON
451 gPad->SetCursor(kWatch);
452 gPad->SetFillColor(1);
453 gPad->Clear();
454
455 Int_t iret;
456 TView *view = new TView(1);
457 Float_t range = fRrange*fRangeSlider->GetMaximum();
458 view->SetRange(-range, -range, 0, range, range, 5*range);
459 fZoomX0[0] = -1;
460 fZoomY0[0] = -1;
461 fZoomX1[0] = 1;
462 fZoomY1[0] = 1;
463 // Display Alice geometry
464 gAlice->GetGeometry()->Draw("same");
465 //Loop on all detectors to add their products to the pad
466 DrawHits();
467 // add itself to the list (last)
468 AppendPad();
469
470 view->SetView(phi, theta, psi, iret);
471}
472//-------------------------------------------------------------------
473void AliMUONRecoDisplay::SetDrawHits(Bool_t hits)
474{
475// Turns on/off Geant hits drawing
476 fDrawHits = hits;
477 DrawView(0,-90);
478}
479
480//-------------------------------------------------------------------
481void AliMUONRecoDisplay::CutMomentum(Double_t min, Double_t max)
482{
483// Define momentum cut for displayed reconstructed tracks
484 fMinMomentum = min;
485 fMaxMomentum = max;
486 if (fHighlited >= 0) UnHighlight();
487 DrawView(0,-90);
488}
489
490//-------------------------------------------------------------------
491Int_t AliMUONRecoDisplay::GetBestMatch(Int_t indr, Float_t tolerance)
492{
493// Find the index of best Geant track matching a reconstructed track : track
494// with maximum number of compatible hits (within tolerance*sigma bending and
495// non-bending resolution) and minimum number of fake hits;
496// If no match is found within a given tolerance, the method is called recursively
497// with increasing tolerance, until tolerance = 10;
2c0075de 498 if (fEmpty) return -1;
c7ba256d 499 if (indr<0 || indr>=fRecoTracks->GetEntriesFast()) return -1;
500 AliMUONRecoTrack *rtrack = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(indr);
501 AliMUONRecoTrack *gtrack = 0;
502 Int_t bestMatch = -1;
503 Int_t maxNcompat = 0;
504 Int_t minNfakes = 10;
505 Double_t xrhit,yrhit,radius,xghit,yghit,dX,dY;
506// loop over all Geant tracks
507 for (Int_t indg=0; indg<fGenTracks->GetEntriesFast(); indg++) {
508 gtrack = (AliMUONRecoTrack*)fGenTracks->UncheckedAt(indg);
2c0075de 509 if (!gtrack) continue;
c7ba256d 510 Int_t ncompat = 0; // number of compat. hits for this track
511 Int_t nfake = 0; // number of fakes
512 // loop chambers to find compatible hits
513 for (Int_t ch=0; ch<10; ch++) {
514 xrhit = rtrack->GetPosX(ch);
515 yrhit = rtrack->GetPosY(ch);
516 radius = TMath::Sqrt(xrhit*xrhit + yrhit*yrhit);
517 if (radius<10) continue; // skip null hits
518 xghit = gtrack->GetPosX(ch);
519 yghit = gtrack->GetPosY(ch);
520 dX = TMath::Abs(xghit-xrhit);
521 dY = TMath::Abs(yghit-yrhit);
2c0075de 522 if (dX<tolerance*0.144 && dY<tolerance*0.01) {// within tol*sigma resolution
c7ba256d 523 ncompat++;
524 continue; // compatible hit
525 } else nfake++; // fake hit
526 }
527 if (ncompat && ncompat>=maxNcompat && nfake<minNfakes) { // this is best matching
528 maxNcompat = ncompat;
529 minNfakes = nfake;
530 bestMatch = indg;
531 }
532 }
533 if (bestMatch<0 && tolerance<=9.) bestMatch = GetBestMatch(indr, tolerance+=1);
534 if (!fPrinted) {
535 rtrack = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(indr);
536 Int_t sign = rtrack->GetSign();
537 cout << "Reconstructed track : " << indr << "(" << sign << ")" << endl;
538 rtrack->TrackInfo();
539 printf("Best matching Geant track within %i*sgm : %i\n", (Int_t)tolerance, bestMatch);
540 if (bestMatch >=0) {
541 gtrack = (AliMUONRecoTrack*)fGenTracks->UncheckedAt(bestMatch);
542 gtrack->TrackInfo();
543 }
544 cout << "-----------------------------------------------------------------\n";
545 }
546 fPrinted = kTRUE;
547
548 return bestMatch;
549}
550
551//-------------------------------------------------------------------
552void AliMUONRecoDisplay::Highlight(Int_t track)
553{
2c0075de 554// Highlight the specified track
555 if (fEmpty) return;
c7ba256d 556 if (fHighlited >=0) UnHighlight();
557 if (track<0 || track>fPolyRecoList->GetEntries()) return;
558 TPolyLine3D *line = (TPolyLine3D*)fPolyRecoList->UncheckedAt(track);
559 line->SetLineColor(kYellow);
560 line->SetLineWidth(1);
561 fHighlited = track;
2c0075de 562// DrawView(15,-45,135);
563 fPad->cd();
564 Draw();
c7ba256d 565}
566
567//-------------------------------------------------------------------
568void AliMUONRecoDisplay::UnHighlight()
569{
570// Unhighlight a previous highlighted track
2c0075de 571 if (fHighlited < 0 || fEmpty) return; // nothing to do
c7ba256d 572 TPolyLine3D *line = (TPolyLine3D*)fPolyRecoList->UncheckedAt(fHighlited);
573 line->SetLineColor(kRed);
574 line->SetLineWidth(1);
575 fHighlited = -1;
2c0075de 576// DrawView(0,-90);
577 fPad->cd();
578 Draw();
c7ba256d 579}
580
581//-------------------------------------------------------------------
582void AliMUONRecoDisplay::DrawHits()
583{
584// Draw hits for all ALICE detectors. Overwrites the DrawHits() method of the
585// base class for reco. track drawing
586
587 Float_t cutmin, cutmax, etamin, etamax, pmom, smin, smax, eta, theta, r;
588 const Float_t kptcutmax = 2;
589 const Float_t ketacutmax = 1.5;
590 Float_t *pxyz;
591 Int_t ntracks,track;
592 TParticle *particle;
593 TObjArray *points;
594 AliPoints *pm;
595
596 //Get cut slider
597 smax = fCutSlider->GetMaximum();
598 smin = fCutSlider->GetMinimum();
599 cutmin = kptcutmax*smin;
600 if (smax < 0.98) cutmax = kptcutmax*smax;
601 else cutmax = 100000;
602
603 //Get eta slider
604 smax = fEtaSlider->GetMaximum();
605 smin = fEtaSlider->GetMinimum();
606 etamin = ketacutmax*(2*smin-1);
607 etamax = ketacutmax*(2*smax-1);
608 if (smin < 0.02) etamin = -1000;
609 if (smax > 0.98) etamax = 1000;
610
611 TIter next(gAlice->Modules());
612 AliModule *module;
613 fHitsCuts = 0;
614 if (fDrawHits) {
615 // draw hits in all modules
616 while((module = (AliModule*)next())) {
617 if (!module->IsActive()) continue;
618 points = module->Points();
619 if (!points) continue;
620 ntracks = points->GetEntriesFast();
621 for (track=0;track<ntracks;track++) {
622 pm = (AliPoints*)points->UncheckedAt(track);
623 if (!pm) continue;
624 particle = pm->GetParticle();
625 if (!particle) continue;
626 pmom = particle->P();
627 if (pmom < cutmin) continue;
628 if (pmom > cutmax) continue;
629 // as a first approximation, take eta of first point
630 pxyz = pm->GetP();
631 r = TMath::Sqrt(pxyz[0]*pxyz[0] + pxyz[1]*pxyz[1]);
632 theta = TMath::ATan2(r,TMath::Abs(pxyz[2]));
633 if(theta) eta = -TMath::Log(TMath::Tan(0.5*theta)); else eta = 1e10;
634 if (pxyz[2] < 0) eta = -eta;
635 if (eta < etamin || eta > etamax) continue;
636 pm->Draw();
637 fHitsCuts += pm->GetN();
638 }
639 }
640 }
641 // draw reconstructed tracks
2c0075de 642 if (fEmpty) return;
c7ba256d 643 TPolyLine3D *line, *gline;
644 Int_t bestMatch;
645 Double_t px,py,pz,p;
646 AliMUONRecoTrack *rtrack;
647
648 if (fHighlited >= 0) {
649 line = (TPolyLine3D*)fPolyRecoList->UncheckedAt(fHighlited);
650 line->Draw();
651 fPrinted = kFALSE;
652 bestMatch = GetBestMatch(fHighlited);
653 if (bestMatch>=0) {
654 gline = (TPolyLine3D*)fPolyGenList->UncheckedAt(bestMatch);
655 gline->SetLineColor(kRed);
656 gline->SetLineWidth(2);
657 gline->SetLineStyle(2);
658 gline->Draw();
659 }
660 } else {
661 for (track=0; track<fPolyRecoList->GetEntries(); track++) {
662 rtrack = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(track);
663 px = rtrack->GetMomReconstr(0);
664 py = rtrack->GetMomReconstr(1);
665 pz = rtrack->GetMomReconstr(2);
666 p = rtrack->P();
667 if (p>fMinMomentum && p<fMaxMomentum) {
668 line = (TPolyLine3D*)fPolyRecoList->UncheckedAt(track);
669 line->Draw();
670 }
671 }
672 }
673}
674
675//-------------------------------------------------------------------
676void AliMUONRecoDisplay::ListTracks()
677{
678// List momentum information of all reconstructed traccks within fPmin and fPmax
679// cuts, as well as their best matching Geant tracks
2c0075de 680 if (fEmpty) return;
c7ba256d 681 cout << "================================================================\n";
682 printf("Reconstructed tracks with momentum in range : %g , %g [GeV/c]\n",
683 fMinMomentum, fMaxMomentum);
684 cout << "----------------------------------------------------------------\n";
685 AliMUONRecoTrack *rtrack;
686 Double_t p;
687 Int_t sign;
688 for (Int_t ind=0; ind<fRecoTracks->GetEntries(); ind++) {
689 rtrack = (AliMUONRecoTrack*)fRecoTracks->UncheckedAt(ind);
690 p = rtrack->P();
691 if (p>fMinMomentum && p<fMaxMomentum) {
692 fPrinted = kFALSE;
693 GetBestMatch(ind);
694 sign = rtrack->GetSign();
695 }
696 }
697 cout << "================================================================\n";
698}
699
700//-------------------------------------------------------------------
701TClonesArray* AliMUONRecoDisplay::MakePolyLines3D(TClonesArray *tracklist)
702{
703// Makes the list of polylines3D corresponding to the list of tracks
2c0075de 704 if (fEmpty) return 0;
c7ba256d 705 if (tracklist!=fRecoTracks && tracklist!=fGenTracks) return 0;
706 Bool_t reco = (tracklist==fRecoTracks)?kTRUE:kFALSE;
707 // make sure there is no other list in memory
708 if (reco) {
709 if (fPolyRecoList) {
710 fPolyRecoList->Delete();
711 delete fPolyRecoList;
712 fPolyRecoList = 0;
713 }
714 } else {
715 if (fPolyGenList) {
716 fPolyGenList->Delete();
717 delete fPolyGenList;
718 fPolyGenList = 0;
719 }
720 }
721 if (!tracklist->GetEntries()) return 0;
722
723 AliMUONRecoTrack* track = 0;
724 TClonesArray *polyLines3D = new TClonesArray("TPolyLine3D",1000);
725 TClonesArray &polylist = *polyLines3D;
726 TPolyLine3D *polyline = 0;
727
728 // loop all tracks
729 for (Int_t i=0; i<tracklist->GetEntries(); i++) {
730 track = (AliMUONRecoTrack*)tracklist->UncheckedAt(i);
731 Int_t ch = 0;
732 Double_t x,y,z,r;
733 polyline = new(polylist[i]) TPolyLine3D(2,"");
734 polyline->SetLineColor(kRed);
735 polyline->SetLineWidth(1);
736 polyline->SetNextPoint(0,0,track->GetVertexPos()); // vertex point
737 // loop chambers to fill TPolyLine3D objects
738 for (ch=0; ch<10; ch++) {
739 x = track->GetPosX(ch);
740 y = track->GetPosY(ch);
741 r = TMath::Sqrt(x*x + y*y);
742 if (r < 10) continue;
743 z = track->GetPosZ(ch);
744 polyline->SetNextPoint(x,y,z);
745 }
746 }
747 return polyLines3D;
748}
749
750//-------------------------------------------------------------------
751void AliMUONRecoDisplay::PolyLineInfo(TClonesArray *line3Dlist)
752{
753// Prints information (x, y, z coordinates) for all constructed polylines
2c0075de 754 if (fEmpty) return;
c7ba256d 755 if (line3Dlist) {
756 TPolyLine3D *polyline = 0;
757 for(Int_t trackIndex=0; trackIndex<line3Dlist->GetEntries(); trackIndex++) {
758 polyline = (TPolyLine3D*)line3Dlist->UncheckedAt(trackIndex);
759 polyline->ls();
760 Float_t *pl = polyline->GetP();
761 for (Int_t i=0; i<polyline->GetN() ;i++) {
762 printf(" x[%d]=%g, y[%d]=%g, z[%d]=%g\n",i,pl[3*i],i,pl[3*i+1],i,pl[3*i+2]);
763 }
764 }
765 }
766}
767
768