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