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