]> git.uio.no Git - u/mrichter/AliRoot.git/blame - MUON/AliMUONClusterDrawAZ.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Doing cluster finding in local coordinates. Speeding up the
[u/mrichter/AliRoot.git] / MUON / AliMUONClusterDrawAZ.cxx
CommitLineData
1af223d7 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
d35a592d 16/* $Id$ */
17
1af223d7 18// Cluster drawing object for AZ cluster finder
19
20#include <stdlib.h>
21#include <Riostream.h>
22//#include <TROOT.h>
23#include <TCanvas.h>
24#include <TLine.h>
25//#include <TTree.h>
26#include <TH2.h>
27#include <TView.h>
28#include <TStyle.h>
29
30#include "AliMUONClusterDrawAZ.h"
31#include "AliMUONClusterFinderAZ.h"
0627f609 32#include "AliMUONGeometryModuleTransformer.h"
1af223d7 33#include "AliHeader.h"
34#include "AliRun.h"
35#include "AliMUON.h"
36//#include "AliMUONChamber.h"
37#include "AliMUONDigit.h"
38#include "AliMUONHit.h"
39#include "AliMUONRawCluster.h"
0627f609 40#include "AliMUONClusterInput.h"
1af223d7 41#include "AliMUONPixel.h"
42//#include "AliMC.h"
43#include "AliMUONLoader.h"
44#include "AliLog.h"
45
46ClassImp(AliMUONClusterDrawAZ)
47
48//_____________________________________________________________________________
49AliMUONClusterDrawAZ::AliMUONClusterDrawAZ()
50 : TObject()
51{
52// Default constructor
53 fFind = NULL; fData = NULL;
54 for (Int_t i=0; i<4; i++) fHist[i] = NULL;
55}
56
57//_____________________________________________________________________________
58AliMUONClusterDrawAZ::AliMUONClusterDrawAZ(AliMUONClusterFinderAZ *clusFinder)
59 : TObject()
60{
61// Constructor
62 fFind = clusFinder;
63 for (Int_t i=0; i<4; i++) fHist[i] = NULL;
64 fDebug = 1;
0627f609 65 fEvent = fChamber = fidDE = 0;
66 fModif = 0; //0;
1af223d7 67 Init();
68}
69
70//_____________________________________________________________________________
71AliMUONClusterDrawAZ::~AliMUONClusterDrawAZ()
72{
73 // Destructor
74}
75
76//_____________________________________________________________________________
77AliMUONClusterDrawAZ::AliMUONClusterDrawAZ(const AliMUONClusterDrawAZ& rhs)
78 : TObject(rhs)
79{
80// Protected copy constructor
81
82 AliFatal("Not implemented.");
83}
84
85
86//_____________________________________________________________________________
87AliMUONClusterDrawAZ&
88AliMUONClusterDrawAZ::operator=(const AliMUONClusterDrawAZ& rhs)
89{
90// Protected assignement operator
91
92 if (this == &rhs) return *this;
93
94 AliFatal("Not implemented.");
95 return *this;
96}
97
98//_____________________________________________________________________________
99void AliMUONClusterDrawAZ::Init()
100{
101 // Initialization
102
103 TCanvas *c1 = new TCanvas("c1","Clusters",0,0,600,700);
104 //c1->SetFillColor(10);
105 c1->Divide(1,2);
106 new TCanvas("c2","Mlem",700,0,600,350);
107
108 // Get pointer to Alice detectors
109 //AliMUON *muon = (AliMUON*) gAlice->GetModule("MUON");
110 //if (!muon) return;
111 //Loaders
112 AliRunLoader *rl = AliRunLoader::GetRunLoader();
113 AliLoader *gime = rl->GetLoader("MUONLoader");
114 fData = ((AliMUONLoader*)gime)->GetMUONData();
115
116 gime->LoadHits("READ");
117 gime->LoadRecPoints("READ");
118}
119
120//_____________________________________________________________________________
121Bool_t AliMUONClusterDrawAZ::FindEvCh(Int_t nev, Int_t ch)
122{
123 // Find requested event and chamber (skip the ones before the selected)
124
125 if (nev < fEvent) return kFALSE;
0627f609 126 else if (nev == fEvent) {
127 if (ch < fChamber) return kFALSE;
128 if (AliMUONClusterInput::Instance()->DetElemId() < fidDE) return kFALSE;
129 }
1af223d7 130 fEvent = nev;
131 fChamber = ch;
0627f609 132 fidDE = AliMUONClusterInput::Instance()->DetElemId();
1af223d7 133 return kTRUE;
134}
135
136//_____________________________________________________________________________
137void AliMUONClusterDrawAZ::DrawCluster()
138{
139 // Draw preclusters
140
141 TCanvas *c1 = (TCanvas*) gROOT->GetListOfCanvases()->FindObject("c1");
142
143 cout << " nev " << fEvent << endl;
144
145 char hName[4];
146 for (Int_t cath = 0; cath < 2; cath++) {
147 // Build histograms
0627f609 148 //if (fHist[cath*2]) {fHist[cath*2]->Delete(); fHist[cath*2] = 0;}
149 //if (fHist[cath*2+1]) {fHist[cath*2+1]->Delete(); fHist[cath*2+1] = 0;}
150 if (fHist[cath*2]) fHist[cath*2] = 0;
151 if (fHist[cath*2+1]) fHist[cath*2+1] = 0;
1af223d7 152 if (fFind->GetNPads(cath) == 0) continue; // cluster on one cathode only
153 Float_t wxMin = 999, wxMax = 0, wyMin = 999, wyMax = 0;
154 Int_t minDx = 0, maxDx = 0, minDy = 0, maxDy = 0;
155 for (Int_t i = 0; i < fFind->GetNPads(0)+fFind->GetNPads(1); i++) {
156 if (fFind->GetIJ(0,i) != cath) continue;
157 if (fFind->GetXyq(3,i) < wxMin) { wxMin = fFind->GetXyq(3,i); minDx = i; }
158 if (fFind->GetXyq(3,i) > wxMax) { wxMax = fFind->GetXyq(3,i); maxDx = i; }
159 if (fFind->GetXyq(4,i) < wyMin) { wyMin = fFind->GetXyq(4,i); minDy = i; }
160 if (fFind->GetXyq(4,i) > wyMax) { wyMax = fFind->GetXyq(4,i); maxDy = i; }
161 }
162 cout << minDx << " " << maxDx << " " << minDy << " " << maxDy << endl;
163 Int_t nx, ny, padSize;
164 Float_t xmin = 9999, xmax = -9999, ymin = 9999, ymax = -9999;
165 if (TMath::Nint(fFind->GetXyq(3,minDx)*1000) == TMath::Nint(fFind->GetXyq(3,maxDx)*1000) &&
166 TMath::Nint(fFind->GetXyq(4,minDy)*1000) == TMath::Nint(fFind->GetXyq(4,maxDy)*1000)) {
167 // the same segmentation
168 cout << " Same" << endl;
169 cout << fFind->GetXyq(3,minDx) << " " << fFind->GetXyq(3,maxDx) << " "
170 << fFind->GetXyq(4,minDy) << " " << fFind->GetXyq(4,maxDy) << endl;
171 for (Int_t i = 0; i < fFind->GetNPads(0)+fFind->GetNPads(1); i++) {
172 if (fFind->GetIJ(0,i) != cath) continue;
173 if (fFind->GetXyq(0,i) < xmin) xmin = fFind->GetXyq(0,i);
174 if (fFind->GetXyq(0,i) > xmax) xmax = fFind->GetXyq(0,i);
175 if (fFind->GetXyq(1,i) < ymin) ymin = fFind->GetXyq(1,i);
176 if (fFind->GetXyq(1,i) > ymax) ymax = fFind->GetXyq(1,i);
177 }
178 xmin -= fFind->GetXyq(3,minDx); xmax += fFind->GetXyq(3,minDx);
179 ymin -= fFind->GetXyq(4,minDy); ymax += fFind->GetXyq(4,minDy);
180 nx = TMath::Nint ((xmax-xmin)/wxMin/2);
181 ny = TMath::Nint ((ymax-ymin)/wyMin/2);
182 cout << xmin << " " << xmax << " " << nx << " " << ymin << " " << ymax << " " << ny << endl;
183 sprintf(hName,"h%d",cath*2);
184 fHist[cath*2] = new TH2D(hName,"cluster",nx,xmin,xmax,ny,ymin,ymax);
185 for (Int_t i = 0; i < fFind->GetNPads(0)+fFind->GetNPads(1); i++) {
186 if (fFind->GetIJ(0,i) != cath) continue;
187 fHist[cath*2]->Fill(fFind->GetXyq(0,i),fFind->GetXyq(1,i),fFind->GetXyq(2,i));
188 //cout << fXyq[0][i] << fXyq[1][i] << fXyq[2][i] << endl;
189 }
190 } else {
191 // different segmentation in the cluster
192 cout << " Different" << endl;
193 cout << fFind->GetXyq(3,minDx) << " " << fFind->GetXyq(3,maxDx) << " "
194 << fFind->GetXyq(4,minDy) << " " << fFind->GetXyq(4,maxDy) << endl;
195 Int_t nOK = 0;
196 Int_t indx, locMin, locMax;
197 if (TMath::Nint(fFind->GetXyq(3,minDx)*1000) != TMath::Nint(fFind->GetXyq(3,maxDx)*1000)) {
198 // different segmentation along x
199 indx = 0;
200 locMin = minDx;
201 locMax = maxDx;
202 } else {
203 // different segmentation along y
204 indx = 1;
205 locMin = minDy;
206 locMax = maxDy;
207 }
208 Int_t loc = locMin;
209 for (Int_t i = 0; i < 2; i++) {
210 // loop over different pad sizes
211 if (i > 0) loc = locMax;
212 padSize = TMath::Nint(fFind->GetXyq(indx+3,loc)*1000);
213 xmin = 9999; xmax = -9999; ymin = 9999; ymax = -9999;
214 for (Int_t j = 0; j < fFind->GetNPads(0)+fFind->GetNPads(1); j++) {
215 if (fFind->GetIJ(0,j) != cath) continue;
216 if (TMath::Nint(fFind->GetXyq(indx+3,j)*1000) != padSize) continue;
217 nOK++;
218 xmin = TMath::Min (xmin,fFind->GetXyq(0,j));
219 xmax = TMath::Max (xmax,fFind->GetXyq(0,j));
220 ymin = TMath::Min (ymin,fFind->GetXyq(1,j));
221 ymax = TMath::Max (ymax,fFind->GetXyq(1,j));
222 }
223 xmin -= fFind->GetXyq(3,loc); xmax += fFind->GetXyq(3,loc);
224 ymin -= fFind->GetXyq(4,loc); ymax += fFind->GetXyq(4,loc);
225 nx = TMath::Nint ((xmax-xmin)/fFind->GetXyq(3,loc)/2);
226 ny = TMath::Nint ((ymax-ymin)/fFind->GetXyq(4,loc)/2);
227 sprintf(hName,"h%d",cath*2+i);
228 fHist[cath*2+i] = new TH2D(hName,"cluster",nx,xmin,xmax,ny,ymin,ymax);
229 for (Int_t j = 0; j < fFind->GetNPads(0)+fFind->GetNPads(1); j++) {
230 if (fFind->GetIJ(0,j) != cath) continue;
231 if (TMath::Nint(fFind->GetXyq(indx+3,j)*1000) != padSize) continue;
232 fHist[cath*2+i]->Fill(fFind->GetXyq(0,j),fFind->GetXyq(1,j),fFind->GetXyq(2,j));
233 }
234 } // for (Int_t i=0;
235 if (nOK != fFind->GetNPads(cath)) cout << " *** Too many segmentations: nPads, nOK "
236 << fFind->GetNPads(cath) << " " << nOK << endl;
237 } // if (TMath::Nint(fFind->GetXyq(3,minDx)*1000)
238 } // for (Int_t cath = 0;
239
240 // Draw histograms and coordinates
d35a592d 241 TH2D *hFake = 0x0;
242 if (fModif) {
243 // This part works only after a modification of THistPainter::PaintLego(Option_t *)
244 // in ROOT
245 Double_t xmin = 9999, ymin = 9999, xmax = -9999, ymax = -9999, aMax = -1;
246 for (Int_t i = 0; i < 4; i++) {
247 if (!fHist[i]) continue;
248 xmin = TMath::Min (xmin, fHist[i]->GetXaxis()->GetXmin());
249 xmax = TMath::Max (xmax, fHist[i]->GetXaxis()->GetXmax());
250 ymin = TMath::Min (ymin, fHist[i]->GetYaxis()->GetXmin());
251 ymax = TMath::Max (ymax, fHist[i]->GetYaxis()->GetXmax());
252 aMax = TMath:: Max (aMax, fHist[i]->GetMaximum());
253 }
254 if (c1->FindObject("hFake")) delete c1->FindObject("hFake");
255 hFake = new TH2D ("hFake", "", 1, xmin, xmax, 1, ymin, ymax);
256 hFake->SetMaximum(aMax);
257 hFake->SetNdivisions(505,"Z");
258 hFake->SetStats(kFALSE);
259 }
260
261 TObject *stats = 0x0;
1af223d7 262 for (Int_t cath = 0; cath < 2; cath++) {
263 if (cath == 0) ModifyHistos();
264 if (fFind->GetNPads(cath) == 0) continue; // cluster on one cathode only
265 c1->cd(cath+1);
266 gPad->SetTheta(55);
267 gPad->SetPhi(30);
d35a592d 268 if (fModif) {
269 if (fHist[cath*2]) {
270 fHist[cath*2]->Draw();
271 gPad->Update();
272 stats = fHist[cath*2]->GetListOfFunctions()->FindObject("stats");
273 }
274 hFake->Draw("legoFb");
275 }
276
1af223d7 277 Double_t x, y, x0, y0, r1 = 999, r2 = 0;
278 if (fHist[cath*2+1]) {
279 //
280 x0 = fHist[cath*2]->GetXaxis()->GetXmin() - 1000*TMath::Cos(30*TMath::Pi()/180);
281 y0 = fHist[cath*2]->GetYaxis()->GetXmin() - 1000*TMath::Sin(30*TMath::Pi()/180);
282 r1 = 0;
283 Int_t ihist=cath*2;
284 for (Int_t iy = 1; iy <= fHist[ihist]->GetNbinsY(); iy++) {
285 y = fHist[ihist]->GetYaxis()->GetBinCenter(iy)
286 + fHist[ihist]->GetYaxis()->GetBinWidth(iy);
287 for (Int_t ix = 1; ix <= fHist[ihist]->GetNbinsX(); ix++) {
288 if (fHist[ihist]->GetCellContent(ix,iy) > 0.1) {
289 x = fHist[ihist]->GetXaxis()->GetBinCenter(ix)
290 + fHist[ihist]->GetXaxis()->GetBinWidth(ix);
291 r1 = TMath::Max (r1,TMath::Sqrt((x-x0)*(x-x0)+(y-y0)*(y-y0)));
292 }
293 }
294 }
295 ihist = cath*2 + 1 ;
296 for (Int_t iy = 1; iy <= fHist[ihist]->GetNbinsY(); iy++) {
297 y = fHist[ihist]->GetYaxis()->GetBinCenter(iy)
298 + fHist[ihist]->GetYaxis()->GetBinWidth(iy);
299 for (Int_t ix = 1; ix <= fHist[ihist]->GetNbinsX(); ix++) {
300 if (fHist[ihist]->GetCellContent(ix,iy) > 0.1) {
301 x = fHist[ihist]->GetXaxis()->GetBinCenter(ix)
302 + fHist[ihist]->GetXaxis()->GetBinWidth(ix);
303 r2 = TMath::Max (r2,TMath::Sqrt((x-x0)*(x-x0)+(y-y0)*(y-y0)));
304 }
305 }
306 }
307 cout << r1 << " " << r2 << endl;
308 } // if (fHist[cath*2+1])
309 if (r1 > r2) {
d35a592d 310 if (fModif) fHist[cath*2]->Draw("lego1FbSame");
311 else fHist[cath*2]->Draw("lego1Fb");
1af223d7 312 if (fHist[cath*2+1]) fHist[cath*2+1]->Draw("lego1SameAxisBbFb");
313 } else {
d35a592d 314 if (fModif) fHist[cath*2+1]->Draw("lego1FbSame");
315 else fHist[cath*2+1]->Draw("lego1Fb");
1af223d7 316 fHist[cath*2]->Draw("lego1SameAxisFbBb");
317 }
318 c1->Update();
d35a592d 319 if (fModif) stats->Draw();
1af223d7 320 } // for (Int_t cath = 0;
321
322 // Draw simulated and reconstructed hits
323 DrawHits();
324}
325
326//_____________________________________________________________________________
327void AliMUONClusterDrawAZ::DrawHits()
328{
329 // Draw simulated and reconstructed hits
330
331 TView *view[2] = { 0x0, 0x0 };
0627f609 332 Double_t p1[3]={0}, p2[3], xNDC[6], xl, yl, zl;
1af223d7 333 TLine *line[99] = {0};
334 TCanvas *c1 = (TCanvas*) gROOT->GetListOfCanvases()->FindObject("c1");
335 if (c1) {
336 c1->cd(1);
337 view[0] = c1->Pad()->GetView();
338 c1->cd(2);
339 view[1] = c1->Pad()->GetView();
340 }
341 fData->SetTreeAddress("H RC");
342
343 TH2D *hist = fHist[0] ? fHist[0] : fHist[2];
344 p2[2] = hist->GetMaximum();
345
346 // Draw simulated hits
347 cout << " *** Simulated hits *** " << endl;
0627f609 348 Int_t ntracks = 0;
349 if (fData->TreeH()) ntracks = (Int_t) fData->GetNtracks();
1af223d7 350 fnMu = 0;
351 Int_t ix, iy, iok, nLine = 0;
352 TClonesArray *hits = NULL;
353 for (Int_t i = 0; i < ntracks; i++) {
354 fData->GetTrack(i);
355 hits = fData->Hits();
356 Int_t nhits = (Int_t) hits->GetEntriesFast();
357 AliMUONHit* mHit;
358 for (Int_t ihit = 0; ihit < nhits; ihit++) {
359 mHit = (AliMUONHit*) hits->UncheckedAt(ihit);
360 if (mHit->Chamber() != fChamber+1) continue; // chamber number
0627f609 361 AliMUONClusterInput::Instance()->Segmentation2(0)->GetTransformer()->
362 Global2Local(fidDE, mHit->X(), mHit->Y(), mHit->Z(), xl, yl, zl);
363 if (TMath::Abs(zl-fFind->GetZpad()) > 1) continue; // different slat
364 p2[0] = p1[0] = xl; // x-pos of hit
365 p2[1] = p1[1] = yl; // y-pos
1af223d7 366 if (p1[0] < hist->GetXaxis()->GetXmin() ||
367 p1[0] > hist->GetXaxis()->GetXmax()) continue;
368 if (p1[1] < hist->GetYaxis()->GetXmin() ||
369 p1[1] > hist->GetYaxis()->GetXmax()) continue;
370 // Check if track comes thru pads with signal
371 iok = 0;
372 for (Int_t ihist = 0; ihist < 4; ihist++) {
373 if (!fHist[ihist]) continue;
374 ix = fHist[ihist]->GetXaxis()->FindBin(p1[0]);
375 iy = fHist[ihist]->GetYaxis()->FindBin(p1[1]);
376 if (fHist[ihist]->GetCellContent(ix,iy) > 0.5) {iok = 1; break;}
377 }
378 if (!iok) continue;
379 gStyle->SetLineColor(1);
380 if (TMath::Abs((Int_t)mHit->Particle()) == 13) {
381 gStyle->SetLineColor(4);
382 if (fnMu < 2) {
383 fxyMu[fnMu][0] = p1[0];
384 fxyMu[fnMu++][1] = p1[1];
385 }
386 }
0627f609 387 printf(" Local coord.: X=%10.4f, Y=%10.4f\n",p1[0],p1[1]);
388 printf(" Global coord.: X=%10.4f, Y=%10.4f, Z=%10.4f\n",mHit->X(),mHit->Y(),mHit->Z());
1af223d7 389 if (view[0] || view[1]) {
390 // Take into account track angles
391 p2[0] += mHit->Tlength() * TMath::Sin(mHit->Theta()/180*TMath::Pi())
392 * TMath::Cos(mHit->Phi()/180*TMath::Pi()) / 2;
393 p2[1] += mHit->Tlength() * TMath::Sin(mHit->Theta()/180*TMath::Pi())
394 * TMath::Sin(mHit->Phi()/180*TMath::Pi()) / 2;
395 for (Int_t ipad = 1; ipad < 3; ipad++) {
396 c1->cd(ipad);
d35a592d 397 if (!view[ipad-1]) continue;
1af223d7 398 view[ipad-1]->WCtoNDC(p1, &xNDC[0]);
399 view[ipad-1]->WCtoNDC(p2, &xNDC[3]);
400 //c1->DrawLine(xpad[0],xpad[1],xpad[3],xpad[4]);
401 line[nLine] = new TLine(xNDC[0],xNDC[1],xNDC[3],xNDC[4]);
402 line[nLine++]->Draw();
403 }
404 }
405 } // for (Int_t ihit = 0; ihit < nhits;
406 } // for (Int_t i = 0; i < ntracks;
407
408 // Draw reconstructed coordinates
d35a592d 409 if (fData->TreeR()) fData->GetRawClusters();
1af223d7 410 TClonesArray *rawclust = fData->RawClusters(fChamber);
411 AliMUONRawCluster *mRaw;
412 gStyle->SetLineColor(3);
413 cout << " *** Reconstructed hits *** " << endl;
414 if (rawclust) {
415 for (Int_t i = 0; i < rawclust ->GetEntriesFast(); i++) {
416 mRaw = (AliMUONRawCluster*)rawclust ->UncheckedAt(i);
0627f609 417 AliMUONClusterInput::Instance()->Segmentation2(0)->GetTransformer()->
418 Global2Local(fidDE, mRaw->GetX(0), mRaw->GetY(0), mRaw->GetZ(0), xl, yl, zl);
419 if (TMath::Abs(zl-fFind->GetZpad()) > 1) continue; // different slat
420 p2[0] = p1[0] = xl; // x-pos of hit
421 p2[1] = p1[1] = yl; // y-pos
1af223d7 422 if (p1[0] < hist->GetXaxis()->GetXmin() ||
423 p1[0] > hist->GetXaxis()->GetXmax()) continue;
424 if (p1[1] < hist->GetYaxis()->GetXmin() ||
425 p1[1] > hist->GetYaxis()->GetXmax()) continue;
426 /*
427 treeD->GetEvent(cath);
428 cout << mRaw->fMultiplicity[0] << mRaw->fMultiplicity[1] << endl;
429 for (Int_t j=0; j<mRaw->fMultiplicity[cath]; j++) {
430 Int_t digit = mRaw->fIndexMap[j][cath];
431 cout << ((AliMUONDigit*)fMuonDigits->UncheckedAt(digit))->Signal() << endl;
432 }
433 */
434 // Check if track comes thru pads with signal
435 iok = 0;
436 for (Int_t ihist = 0; ihist < 4; ihist++) {
437 if (!fHist[ihist]) continue;
438 ix = fHist[ihist]->GetXaxis()->FindBin(p1[0]);
439 iy = fHist[ihist]->GetYaxis()->FindBin(p1[1]);
440 if (fHist[ihist]->GetCellContent(ix,iy) > 0.5) {iok = 1; break;}
441 }
442 if (!iok) continue;
0627f609 443 printf(" Local coord.: X=%10.4f, Y=%10.4f\n",p1[0],p1[1]);
444 printf(" Global coord.: X=%10.4f, Y=%10.4f, Z=%10.4f\n",mRaw->GetX(0),mRaw->GetY(0),mRaw->GetZ(0));
1af223d7 445 if (view[0] || view[1]) {
446 for (Int_t ipad = 1; ipad < 3; ipad++) {
447 c1->cd(ipad);
d35a592d 448 if (!view[ipad-1]) continue;
1af223d7 449 view[ipad-1]->WCtoNDC(p1, &xNDC[0]);
450 view[ipad-1]->WCtoNDC(p2, &xNDC[3]);
451 line[nLine] = new TLine(xNDC[0],xNDC[1],xNDC[3],xNDC[4]);
452 line[nLine++]->Draw();
453 }
454 }
455 } // for (Int_t i = 0; i < rawclust ->GetEntries();
456 } // if (rawclust)
457 c1->Update();
458}
459
460//_____________________________________________________________________________
461Int_t AliMUONClusterDrawAZ::Next()
462{
463 // What to do next?
464 // File
465 FILE *lun = 0;
466 //lun = fopen("pull.dat","w");
467
468 for (Int_t i = 0; i < fnMu; i++) {
469 // Check again if muon comes thru the used pads (due to extra splitting)
470 for (Int_t j = 0; j < fFind->GetNPads(0)+fFind->GetNPads(1); j++) {
471 if (TMath::Abs(fxyMu[i][0]-fFind->GetXyq(0,j))<fFind->GetXyq(3,j) &&
472 TMath::Abs(fxyMu[i][1]-fFind->GetXyq(1,j))<fFind->GetXyq(4,j)) {
473 if (fDebug) printf("%12.3e %12.3e %12.3e %12.3e\n",fxyMu[i][2],fxyMu[i][3],fxyMu[i][4],fxyMu[i][5]);
474 if (lun) fprintf(lun,"%4d %2d %12.3e %12.3e %12.3e %12.3e\n",fEvent,fChamber,fxyMu[i][2],fxyMu[i][3],fxyMu[i][4],fxyMu[i][5]);
475 break;
476 }
477 }
478 } // for (Int_t i=0; i<fnMu;
479
480 // What's next?
481 char command[8];
482 cout << " What is next? " << endl;
483 command[0] = ' ';
484 gets(command);
0627f609 485 if (command[0] == 'n' || command[0] == 'N') { fEvent++; fChamber = fidDE = 0; } // next event
1af223d7 486 else if (command[0] == 'q' || command[0] == 'Q') { if (lun) fclose(lun); } // exit display
487 else if (command[0] == 'c' || command[0] == 'C') sscanf(command+1,"%d",&fChamber); // new chamber
0627f609 488 else if (command[0] == 'e' || command[0] == 'E') { sscanf(command+1,"%d",&fEvent); fChamber = fidDE = 0; } // new event
489 else if (command[0] == 'd' || command[0] == 'D') { sscanf(command+1,"%d",&fidDE); fChamber = fidDE / 100 - 1; } // new DetElem.
1af223d7 490 else return 1; // Next precluster
491 return 0;
492}
493
494
495//_____________________________________________________________________________
496void AliMUONClusterDrawAZ::ModifyHistos(void)
497{
498 // Modify histograms to bring them to (approximately) the same size
499 Int_t nhist = 0;
500 Float_t hlim[4][4], hbin[4][4]; // first index - xmin, xmax, ymin, ymax
501
502 Float_t binMin[4] = {999,999,999,999};
503
504 for (Int_t i = 0; i < 4; i++) {
505 if (!fHist[i]) {
506 hlim[0][i] = hlim[2][i] = 999;
507 hlim[1][i] = hlim[3][i] = -999;
508 continue;
509 }
510 hlim[0][i] = fHist[i]->GetXaxis()->GetXmin(); // xmin
511 hlim[1][i] = fHist[i]->GetXaxis()->GetXmax(); // xmax
512 hlim[2][i] = fHist[i]->GetYaxis()->GetXmin(); // ymin
513 hlim[3][i] = fHist[i]->GetYaxis()->GetXmax(); // ymax
514 hbin[0][i] = hbin[1][i] = fHist[i]->GetXaxis()->GetBinWidth(1);
515 hbin[2][i] = hbin[3][i] = fHist[i]->GetYaxis()->GetBinWidth(1);
516 binMin[0] = TMath::Min(binMin[0],hbin[0][i]);
517 binMin[2] = TMath::Min(binMin[2],hbin[2][i]);
518 nhist++;
519 }
520 binMin[1] = binMin[0];
521 binMin[3] = binMin[2];
522 cout << " Nhist: " << nhist << endl;
523
524 // Adjust histo limits for cathode with different segmentation
525 for (Int_t i = 0; i < 4; i+=2) {
526 if (!fHist[i+1]) continue;
527 Int_t imin, imax, i1 = i + 1;
528 for (Int_t lim = 0; lim < 4; lim++) {
529 while (1) {
530 if (hlim[lim][i] < hlim[lim][i1]) {
531 imin = i;
532 imax = i1;
533 } else {
534 imin = i1;
535 imax = i;
536 }
537 if (TMath::Abs(hlim[lim][imin]-hlim[lim][imax])<0.01*binMin[lim]) break;
538 if (lim == 0 || lim == 2) {
539 // find lower limit
540 hlim[lim][imax] -= hbin[lim][imax];
541 } else {
542 // find upper limit
543 hlim[lim][imin] += hbin[lim][imin];
544 }
545 } // while (1)
546 }
547 }
548
549
550 Int_t imnmx = 0, nExtra = 0;
551 for (Int_t lim = 0; lim < 4; lim++) {
552 if (lim == 0 || lim == 2) imnmx = TMath::LocMin(4,hlim[lim]); // find lower limit
553 else imnmx = TMath::LocMax(4,hlim[lim]); // find upper limit
554
555 // Adjust histogram limit
556 for (Int_t i = 0; i < 4; i++) {
557 if (!fHist[i]) continue;
558 nExtra = TMath::Nint ((hlim[lim][imnmx]-hlim[lim][i]) / hbin[lim][i]);
559 hlim[lim][i] += nExtra * hbin[lim][i];
560 }
561 }
562
563 // Rebuild histograms
564 TH2D *hist = 0;
565 Int_t nx, ny;
566 Double_t x, y, cont, cmax=0;
567 char hName[4];
568 for (Int_t ihist = 0; ihist < 4; ihist++) {
569 if (!fHist[ihist]) continue;
570 nx = TMath::Nint((hlim[1][ihist]-hlim[0][ihist])/hbin[0][ihist]);
571 ny = TMath::Nint((hlim[3][ihist]-hlim[2][ihist])/hbin[2][ihist]);
572 cout << ihist << " " << hlim[0][ihist] << " " << hlim[1][ihist] << " " << nx;
573 cout << " " << hlim[2][ihist] << " " << hlim[3][ihist] << " " << ny << endl;
574 sprintf(hName,"hh%d",ihist);
575 hist = new TH2D(hName,"hist",nx,hlim[0][ihist],hlim[1][ihist],ny,hlim[2][ihist],hlim[3][ihist]);
576 for (Int_t i=1; i<=fHist[ihist]->GetNbinsX(); i++) {
577 x = fHist[ihist]->GetXaxis()->GetBinCenter(i);
578 for (Int_t j=1; j<=fHist[ihist]->GetNbinsY(); j++) {
579 y = fHist[ihist]->GetYaxis()->GetBinCenter(j);
580 cont = fHist[ihist]->GetCellContent(i,j);
581 hist->Fill(x,y,cont);
582 }
583 }
584 cmax = TMath::Max (cmax,hist->GetMaximum());
585 sprintf(hName,"%s%d",fHist[ihist]->GetName(),ihist);
586 fHist[ihist]->Delete();
587 fHist[ihist] = new TH2D(*hist);
588 fHist[ihist]->SetName(hName);
589 fHist[ihist]->SetNdivisions(505,"Z");
590 hist->Delete();
591 }
592 if (fDebug) printf("%f \n",cmax);
593
594 for (Int_t ihist = 0; ihist < 4; ihist++) {
595 if (!fHist[ihist]) continue;
596 fHist[ihist]->SetMaximum(cmax);
597 fHist[ihist]->SetMinimum(0);
598 }
599}
600
601//_____________________________________________________________________________
602void AliMUONClusterDrawAZ::AdjustHist(Double_t *xylim, const AliMUONPixel *pixPtr)
603{
604 // Adjust histogram limits for pixel drawing
605
606 Float_t xypads[4];
607 if (fHist[0]) {
608 xypads[0] = fHist[0]->GetXaxis()->GetXmin();
609 xypads[1] = -fHist[0]->GetXaxis()->GetXmax();
610 xypads[2] = fHist[0]->GetYaxis()->GetXmin();
611 xypads[3] = -fHist[0]->GetYaxis()->GetXmax();
612 for (Int_t i = 0; i < 4; i++) {
613 while(1) {
614 if (xylim[i] < xypads[i]) break;
615 xylim[i] -= 2*pixPtr->Size(i/2);
616 }
617 }
618 }
619}
620
621//_____________________________________________________________________________
622void AliMUONClusterDrawAZ::DrawHist(const char* canvas, TH2D *hist)
623{
624 // Draw histogram in given canvas
625
626 Int_t ix = 0;
627 //((TCanvas*)gROOT->FindObject("c2"))->cd();
628 ((TCanvas*)gROOT->FindObject(canvas))->cd();
629 gPad->SetTheta(55);
630 gPad->SetPhi(30);
631 hist->Draw("lego1Fb");
632 gPad->Update();
0627f609 633 //gets((char*)&ix);
634 if (fnMu) gets((char*)&ix);
1af223d7 635}
636
637//_____________________________________________________________________________
638void AliMUONClusterDrawAZ::FillMuon(Int_t nfit, const Double_t *parOk, const Double_t *errOk)
639{
640 // Fill muon information
641
642 Int_t indx, imax;
643 Double_t cmax, rad;
644 for (Int_t i = 0; i < fnMu; i++) {
645 cmax = fxyMu[i][6];
646 for (Int_t j = 0; j < nfit; j++) {
647 indx = j<2 ? j*2 : j*2+1;
648 rad = (fxyMu[i][0]-parOk[indx])*(fxyMu[i][0]-parOk[indx]) +
649 (fxyMu[i][1]-parOk[indx+1])*(fxyMu[i][1]-parOk[indx+1]);
650 if (rad < cmax) {
651 cmax = rad;
652 imax = indx;
653 fxyMu[i][6] = cmax;
654 fxyMu[i][2] = parOk[imax] - fxyMu[i][0];
655 fxyMu[i][4] = parOk[imax+1] - fxyMu[i][1];
656 fxyMu[i][3] = errOk[imax];
657 fxyMu[i][5] = errOk[imax+1];
658 }
659 }
660 }
661}
662
d35a592d 663//_____________________________________________________________________________
664void AliMUONClusterDrawAZ::UpdateCluster(Int_t npad)
665{
666 // Update cluster after removing non-overlapped pads
667
668 Int_t cath = 0, ix = 0, iy = 0;
669 cout << " Update cluster " << endl;
670 gets((char*)&ix);
671 for (Int_t i = 0; i < npad; i++) {
672 if (TMath::Nint (fFind->GetXyq(2,i)) != -2) continue;
673 cath = fFind->GetIJ(0,i);
674 for (Int_t j = 0; j < 2; j++) {
675 Int_t ihist = cath * 2 + j;
676 if (!fHist[ihist]) continue;
677 ix = fHist[ihist]->GetXaxis()->FindBin(fFind->GetXyq(0,i));
678 iy = fHist[ihist]->GetYaxis()->FindBin(fFind->GetXyq(1,i));
679 Double_t cont = fHist[ihist]->GetCellContent(ix, iy);
680 if (cont < 0.1) continue;
681 fHist[ihist]->Fill(fFind->GetXyq(0,i), fFind->GetXyq(1,i), -cont);
682 }
683 }
684 TCanvas *c1 = (TCanvas*) gROOT->GetListOfCanvases()->FindObject("c1");
685 if (c1) {
686 c1->cd(1);
687 gPad->Modified();
688 gPad->Update();
689 c1->cd(2);
690 gPad->Modified();
691 gPad->Update();
692 }
693}