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