]>
Commit | Line | Data |
---|---|---|
68252c5c | 1 | #if !defined(__CINT__) || defined(__MAKECINT__) |
2 | #include <TRD/AliTRDcluster.h> | |
3 | #include <AliRunLoader.h> | |
4 | #include <TCanvas.h> | |
5 | #include <TF1.h> | |
6 | #include <TH1.h> | |
7 | #include <TH2.h> | |
8 | #include <TLine.h> | |
9 | #include <TMath.h> | |
10 | #include <TObjArray.h> | |
11 | #include <TPad.h> | |
12 | #include <TROOT.h> | |
13 | #include <TStyle.h> | |
14 | #include <TSystem.h> | |
15 | #include <TTree.h> | |
16 | #endif | |
17 | ||
18 | //#define EVENT_BY_EVENT // Show amplitude event-by-event | |
19 | //#define EVENT_BY_EVENT2 // Show tracks event-by-event | |
20 | //#define PRINT_OTHER // Display various additional histograms | |
21 | ||
22 | // Number of noisy detectors | |
23 | const Int_t numNoisyDets = 10; | |
24 | ||
25 | // Holds the detector number of detectors that produce too much noise and | |
26 | // are therefore excluded from the calculation | |
27 | const Int_t exDets[numNoisyDets] = {25, 511, 515, 538, 523, 249, 294, 27, 24, 247}; | |
28 | ||
29 | Bool_t isBadCluster(Int_t det, Int_t Q, Int_t PadCol, Int_t PadRow) | |
30 | { | |
31 | if (Q < 20) return 1; | |
32 | ||
33 | // Exclude bad cols | |
34 | if (PadCol == 142 || PadCol == 136 || PadCol == 66 || PadCol == 65 || PadCol <= 2) return 1; | |
35 | //if (PadCol == 142 || PadCol == 66 || PadCol == 65 || PadCol <= 2) return 1; | |
36 | //if (PadCol == 142 || PadCol <= 2) return 1; | |
37 | ||
38 | // Exclude noisy detectors | |
39 | for (Int_t indDet = 0; indDet < numNoisyDets; indDet++) | |
40 | { | |
41 | if (det == exDets[indDet]) | |
42 | { | |
43 | // Exclude several regions of bad detectors | |
44 | if ((PadCol >= 134 && PadRow <= 5) || (PadCol == 141 && PadRow >= 10)) return 1; | |
45 | //if (PadRow == 15 || (PadCol >= 134 && PadRow <= 11) || (PadCol >= 119 && PadCol <= 123 && PadRow == 11)) return 1; | |
46 | ||
47 | break; | |
48 | } | |
49 | } | |
50 | ||
51 | // Cluster seems to be clean | |
52 | return 0; | |
53 | } | |
54 | ||
55 | ||
56 | //===================================== | |
57 | ||
58 | void AliTRDanalyseRecPoints (const char *filename="galice.root") | |
59 | { | |
60 | gROOT->SetStyle("Plain"); | |
61 | gStyle->SetPalette(1); | |
62 | ||
63 | // open run loader and load gAlice, kinematics and header | |
64 | AliRunLoader* runLoader = AliRunLoader::Open(filename); | |
65 | if (!runLoader) | |
66 | { | |
67 | printf("Error: readKine:\nGetting run loader from file \"%s\" failed", filename); | |
68 | return; | |
69 | } | |
70 | ||
71 | runLoader->LoadHeader(); | |
72 | runLoader->LoadRecPoints("TRD"); | |
73 | TObjArray *module = new TObjArray(); | |
74 | ||
75 | runLoader->CdGAFile(); | |
76 | ||
77 | ||
78 | // Flag indicating whether the signal comes from a noisy detector | |
79 | Bool_t fromNoisyDet = 0; | |
80 | ||
81 | // Counts the number of clusters per event | |
82 | Int_t clusters = 0; | |
83 | ||
84 | // Threshold for the amplitude | |
85 | const Int_t kQThreshold = 30; | |
86 | ||
87 | // Selection criterions for clusters | |
88 | const Int_t kThresholdLower = 30; | |
89 | const Int_t kThresholdUpper = 500; | |
90 | ||
91 | // Number of events passing the criterion | |
92 | Int_t numPassed = 0; | |
93 | ||
94 | // hist | |
95 | #ifdef PRINT_OTHER | |
96 | TH1D *mNCls = new TH1D("NCls", ";number of clusters", 400, 0, 400); | |
97 | TH1D *mResY = new TH1D("resY", ";sigma Y (cm)", 100, 0, 1); | |
98 | TH1D *mResZ = new TH1D("resZ", ";sigma Z (cm)", 100, 1, 3); | |
99 | ||
100 | TH1D *mNPads = new TH1D("NPads", ";number of pads", 10, -0.5, 9.5); | |
101 | TH1D *mCenter = new TH1D("center", ";center", 100, -1, 1); | |
102 | TH1D *mTime = new TH1D("timeBin", ";time bin", 25, -0.5, 24.5); | |
103 | TH1D *mPlane = new TH1D("plane", ";plane", 6, -0.5, 5.5); | |
104 | ||
105 | TH1D *mX = new TH1D("x", "; X (cm)", 100, 0, 3.5); | |
106 | TH1D *mY = new TH1D("y", "; y / chamber width (%)", 100, -1, 1); | |
107 | TH1D *mZ = new TH1D("z", "; Z (cm)", 200, -400, 400); | |
108 | ||
109 | TH2D *mPlaneY = new TH2D("planeY", ";plane;Y (cm);#Clusters", 6, -0.5, 5.5, 100, -60, 60); | |
110 | TH2D *mYX = new TH2D("yx", ";Y (cm);X (cm);#Clusters", 100, -60, 60, 120, 280, 400); | |
111 | #endif | |
112 | ||
113 | TH1D *mChrg = new TH1D("chrg", ";Charge(no pad-filter)", 200, 0, 200); | |
114 | TH1D *mChrg2 = new TH1D("chrg2", ";Charge(pads<9)", 200, 0, 200); | |
115 | TH1D *mChrg3 = new TH1D("chrg3", ";Charge(2<pads<9)", 200, 0, 200); | |
116 | ||
117 | TH2D *mQvsPads = new TH2D("QvsPads", ";Number of pads;Q;#Clusters", 15, -0.5, 14.5, 100, 0, 1000); | |
118 | ||
119 | // There are 540 detectors | |
120 | TH2D *mQvsDetector = new TH2D("QvsDet", ";Number of detector;Q;#Clusters", 540, -0.5, 539.5, 100, 0, 1000); | |
121 | ||
122 | // Supermoduleid=int(detectorNum / 30) | |
123 | TH2D *mQvsSModId = new TH2D("QvsSModId", ";Supermodule id;Q;#Clusters", 18, -0.5, 17.5, 100, 0, 1000); | |
124 | // Histograms for the already installed supermodules | |
125 | TH1D *mSMod0 = new TH1D("SMod0", ";Charge", 200, 0, 200); | |
126 | TH1D *mSMod8 = new TH1D("SMod8", ";Charge", 200, 0, 200); | |
127 | TH1D *mSMod9 = new TH1D("SMod9", ";Charge", 200, 0, 200); | |
128 | TH1D *mSMod17 = new TH1D("SMod17", ";Charge", 200, 0, 200); | |
129 | ||
130 | // Index + 2 corresponds to the number of activated pads | |
131 | TH1D *mPads[14]; | |
132 | char name[8]; | |
133 | for (Int_t ind = 0; ind < 14; ind++) | |
134 | { | |
135 | sprintf(name,"%dPads",ind + 2); | |
136 | mPads[ind] = new TH1D(name, ";Charge", 200, 0, 200); | |
137 | } | |
138 | ||
139 | // Amplitude distribution will be displayed event by event for the first 100 events | |
140 | TH1D *mChrgEBE = new TH1D("chrgEBE", ";Charge", 200, 0, 200); | |
141 | ||
142 | // Keep track of the number of clusters per event | |
143 | TH1D *mClusters = new TH1D("Clusters", ";Number of clusters per event", 300, -0.5, 299.5); | |
144 | ||
145 | // Keep track of the number of clusters per chamber with Q < kQThreshold | |
146 | char clusTitle[55]; | |
147 | sprintf(clusTitle,";Number of clusters with Q < %d vs detector", kQThreshold); | |
148 | TH1D *mClusCham = new TH1D("ClusCham", clusTitle, 540, -0.5, 539.5); | |
149 | // Same but this time noisy detectors excluded | |
150 | sprintf(clusTitle,";Number of clusters with Q < %d vs detector", kQThreshold); | |
151 | TH1D *mClusChamEx = new TH1D("ClusChamEx", clusTitle, 540, -0.5, 539.5); | |
152 | ||
153 | TH2D *mColRow = new TH2D("ColRow", ";Column;Row;#Clusters", 144, -0.5, 143.5, 18, -0.5, 17.5); | |
154 | TH2D *mColRowFiltered = new TH2D("ColRowFiltered", ";Column;Row;#Clusters", 144, -0.5, 143.5, 18, -0.5, 17.5); | |
155 | ||
156 | TH2D *mColTime = new TH2D("ColTime", ";TimeBin;Column", 18, -0.5, 17.5, 144, -0.5, 143.5); | |
157 | TH2D *mColTimeFiltered = new TH2D("ColTimeFiltered", ";TimeBin;Column;#Clusters", 18, -0.5, 17.5, 144, -0.5, 143.5); | |
158 | ||
159 | // YX histograms for the already installed supermodules | |
160 | TH2D *mTrSMod0 = new TH2D("TrSMod0", ";Col;time+layer*20;#Clusters", 144, -0.5, 143.5, 120, -0.5, 119.5); | |
161 | TH2D *mTrSMod8 = new TH2D("TrSMod8", ";Col;time+layer*20;#Clusters", 144, -0.5, 143.5, 120, -0.5, 119.5); | |
162 | TH2D *mTrSMod9 = new TH2D("TrSMod9", ";Col;time+layer*20;#Clusters", 144, -0.5, 143.5, 120, -0.5, 119.5); | |
163 | TH2D *mTrSMod17 = new TH2D("TrSMod17", ";Col;time+layer*20;#Clusters", 144, -0.5, 143.5, 120, -0.5, 119.5); | |
164 | ||
165 | ||
166 | #ifdef EVENT_BY_EVENT | |
167 | TCanvas* eventc = new TCanvas("eventc","Amplitude distribution for event 0",700,500); | |
168 | char title[40]; | |
169 | #endif | |
170 | ||
171 | #ifdef EVENT_BY_EVENT2 | |
172 | TCanvas* sTrc = new TCanvas("sTrc","Track - Supermodules",800,800); | |
173 | sTrc->SetLogz(); | |
174 | sTrc->Divide(2,2); | |
175 | char title2[40]; | |
176 | #endif | |
177 | ||
178 | ||
179 | int nEvents = runLoader->GetNumberOfEvents(); | |
180 | ||
181 | AliTRDcluster *cls = 0; | |
182 | ||
183 | for(Int_t ev = 0; ev < nEvents - 1; ev++) | |
184 | { | |
185 | clusters = 0; | |
186 | ||
187 | TTree *tree = runLoader->GetTreeR("TRD", 0); | |
188 | tree->SetBranchAddress("TRDcluster", &module); | |
189 | ||
190 | int N = tree->GetEntries(); | |
191 | ||
192 | // Check number of clusters | |
193 | for(Int_t ind = 0; ind < N; ind++) | |
194 | { | |
195 | tree->GetEntry(ind); | |
196 | Int_t m = module->GetEntries(); | |
197 | ||
198 | for (Int_t j = 0; j < m; j++) | |
199 | { | |
200 | if (cls != 0) delete cls; | |
201 | cls = (AliTRDcluster*)module->At(j); | |
202 | if (!isBadCluster(cls->GetDetector(), (Int_t)cls->GetQ(), cls->GetPadCol(), cls->GetPadRow())) clusters++; | |
203 | else | |
204 | { | |
205 | // Uncomment the following lines AND comment the next if-clause to display the histograms for the | |
206 | // EXCLUDED clusters | |
207 | ||
208 | //mQvsPads->Fill(cls->GetNPads(),cls->GetQ()); | |
209 | //mQvsDetector->Fill(cls->GetDetector(),cls->GetQ()); | |
210 | //mQvsSModId->Fill((Int_t)(cls->GetDetector() / 30.),cls->GetQ()); | |
211 | //mColTime->Fill(cls->GetLocalTimeBin(), cls->GetPadCol()); | |
212 | //mColRow->Fill(cls->GetPadCol(), cls->GetPadRow()); | |
213 | //mClusCham->Fill(cls->GetDetector()); | |
214 | //mYX->Fill(cls->GetY(), cls->GetX()); | |
215 | } | |
216 | } | |
217 | } | |
218 | ||
219 | // Process event, if number of clusters is in the correct range | |
220 | if (clusters >= kThresholdLower && clusters <= kThresholdUpper) | |
221 | { | |
222 | clusters = 0; | |
223 | numPassed++; | |
224 | ||
225 | for (Int_t i = 0; i < N; i++) | |
226 | { | |
227 | tree->GetEntry(i); | |
228 | int m = module->GetEntries(); | |
229 | ||
230 | #ifdef PRINT_OTHER | |
231 | mNCls->Fill(m); | |
232 | #endif | |
233 | ||
234 | for(Int_t j = 0; j < m; j++) | |
235 | { | |
236 | fromNoisyDet = 0; | |
237 | clusters++; | |
238 | ||
239 | if (cls != 0) delete cls; | |
240 | cls = (AliTRDcluster*)module->At(j); | |
241 | ||
242 | // Filter the clusters of the events, too | |
243 | if (cls->GetQ() < 20) continue; | |
244 | ||
245 | mColTime->Fill(cls->GetLocalTimeBin(), cls->GetPadCol()); | |
246 | ||
247 | if (cls->GetQ() <= kQThreshold) | |
248 | { | |
249 | mClusCham->Fill(cls->GetDetector()); | |
250 | } | |
251 | ||
252 | // Exclude noisy detectors | |
253 | for (Int_t indDet = 0; indDet < numNoisyDets; indDet++) | |
254 | { | |
255 | if (cls->GetDetector() == exDets[indDet]) | |
256 | { | |
257 | fromNoisyDet = 1; | |
258 | break; | |
259 | } | |
260 | } | |
261 | ||
262 | // Exclude I | |
263 | if (cls->GetPadCol() == 142 || cls->GetPadCol() == 136 || cls->GetPadCol() == 66 || cls->GetPadCol() == 65 || cls->GetPadCol() <= 2) continue; | |
264 | //if (cls->GetPadCol() == 142 || cls->GetPadCol() == 66 || cls->GetPadCol() == 65 || cls->GetPadCol() <= 2) continue; | |
265 | //if (cls->GetPadCol() == 142 || cls->GetPadCol() <= 2) continue; | |
266 | ||
267 | if (fromNoisyDet) | |
268 | { | |
269 | mColRow->Fill(cls->GetPadCol(), cls->GetPadRow()); | |
270 | ||
271 | // Exclude II | |
272 | //if (cls->GetPadRow() == 15 || (cls->GetPadCol() >= 134 && cls->GetPadRow() <= 11) || (cls->GetPadCol() >= 119 && cls->GetPadCol() <= 123 && cls->GetPadRow() == 11)) continue; | |
273 | ||
274 | if ((cls->GetPadCol() >= 134 && cls->GetPadRow() <= 5) || (cls->GetPadCol() == 141 && cls->GetPadRow() >= 10)) continue; | |
275 | ||
276 | mColRowFiltered->Fill(cls->GetPadCol(), cls->GetPadRow()); | |
277 | } | |
278 | ||
279 | mColTimeFiltered->Fill(cls->GetLocalTimeBin(), cls->GetPadCol()); | |
280 | ||
281 | if (cls->GetQ() <= kQThreshold) | |
282 | { | |
283 | mClusChamEx->Fill(cls->GetDetector()); | |
284 | } | |
285 | ||
286 | #ifdef PRINT_OTHER | |
287 | mX->Fill(cls->GetX()); | |
288 | mZ->Fill(cls->GetZ()); | |
289 | mYX->Fill(cls->GetY(), cls->GetX()); | |
290 | mResY->Fill(TMath::Sqrt(cls->GetSigmaY2())); | |
291 | mResZ->Fill(TMath::Sqrt(cls->GetSigmaZ2())); | |
292 | mNPads->Fill(cls->GetNPads()); | |
293 | mCenter->Fill(cls->GetCenter()); | |
294 | mTime->Fill(cls->GetLocalTimeBin()); | |
295 | #endif | |
296 | ||
297 | if (cls->GetNPads() > 2 && cls->GetNPads() < 9) mChrg3->Fill(cls->GetQ()); | |
298 | if (cls->GetNPads() < 9) mChrg2->Fill(cls->GetQ()); | |
299 | mChrg->Fill(cls->GetQ()); | |
300 | ||
301 | mQvsPads->Fill(cls->GetNPads(),cls->GetQ()); | |
302 | mQvsDetector->Fill(cls->GetDetector(),cls->GetQ()); | |
303 | mQvsSModId->Fill((Int_t)(cls->GetDetector() / 30.),cls->GetQ()); | |
304 | ||
305 | Int_t superModule = (Int_t)(cls->GetDetector() / 30.); | |
306 | ||
307 | // Position in supermodule = cls->GetDetector() - superModule * 30 | |
308 | // There are 6 layers per supermodule: | |
309 | Int_t layer = (cls->GetDetector() - superModule * 30) % 6; | |
310 | switch (superModule) | |
311 | { | |
312 | case 0: | |
313 | mSMod0->Fill(cls->GetQ()); | |
314 | mTrSMod0->Fill(cls->GetPadCol(), cls->GetLocalTimeBin() + layer * 20); | |
315 | break; | |
316 | case 8: | |
317 | mSMod8->Fill(cls->GetQ()); | |
318 | mTrSMod8->Fill(cls->GetPadCol(), cls->GetLocalTimeBin() + layer * 20); | |
319 | break; | |
320 | case 9: | |
321 | mSMod9->Fill(cls->GetQ()); | |
322 | mTrSMod9->Fill(cls->GetPadCol(), cls->GetLocalTimeBin() + layer * 20); | |
323 | break; | |
324 | case 17: | |
325 | mSMod17->Fill(cls->GetQ()); | |
326 | mTrSMod17->Fill(cls->GetPadCol(), cls->GetLocalTimeBin() + layer * 20); | |
327 | break; | |
328 | default: | |
329 | break; | |
330 | } | |
331 | ||
332 | // Only consider events with maximal 15 pads here | |
333 | if (cls->GetNPads() - 2 < 14 && cls->GetNPads() - 2 >= 0) | |
334 | { | |
335 | mPads[cls->GetNPads() - 2]->Fill(cls->GetQ()); | |
336 | } | |
337 | ||
338 | mChrgEBE->Fill(cls->GetQ()); | |
339 | } | |
340 | } | |
341 | ||
342 | ||
343 | #ifdef EVENT_BY_EVENT | |
344 | // Event by event display for the first 100 events | |
345 | if (numPassed < 100) | |
346 | { | |
347 | sprintf(title,"Charge distribution for event %d", ev); | |
348 | eventc->SetTitle(title); | |
349 | mChrgEBE->Draw(); | |
350 | eventc->Update(); | |
351 | gSystem->Sleep(400); | |
352 | mChrgEBE->Reset(); | |
353 | } | |
354 | #endif | |
355 | ||
356 | #ifdef EVENT_BY_EVENT2 | |
357 | // Event by event display | |
358 | sprintf(title2,"Supermodule tracks for event %d", ev); | |
359 | sTrc->SetTitle(title2); | |
360 | sTrc->cd(1); | |
361 | mTrSMod0->Draw("colz"); | |
362 | sTrc->cd(2); | |
363 | mTrSMod8->Draw("colz"); | |
364 | sTrc->cd(3); | |
365 | mTrSMod9->Draw("colz"); | |
366 | sTrc->cd(4); | |
367 | mTrSMod17->Draw("colz"); | |
368 | sTrc->Update(); | |
369 | gSystem->Sleep(400); | |
370 | mTrSMod0->Reset(); | |
371 | mTrSMod8->Reset(); | |
372 | mTrSMod9->Reset(); | |
373 | mTrSMod17->Reset(); | |
374 | #endif | |
375 | ||
376 | ||
377 | mClusters->Fill(clusters); | |
378 | } | |
379 | ||
380 | runLoader->GetNextEvent(); | |
381 | } | |
382 | ||
383 | #ifdef PRINT_OTHER | |
384 | new TCanvas(); | |
385 | gPad->SetLogy(); | |
386 | mNCls->Draw(); | |
387 | ||
388 | new TCanvas(); | |
389 | gPad->SetLogy(); | |
390 | mResY->Draw(); | |
391 | ||
392 | new TCanvas(); | |
393 | gPad->SetLogy(); | |
394 | mResZ->Draw(); | |
395 | #endif | |
396 | ||
397 | new TCanvas("chrgc", "Charge(no pad-filter)", 700, 500); | |
398 | mChrg->Draw(); | |
399 | // Try landau-fit | |
400 | TF1* lFit = new TF1("lFit", "landau", 20, 200); | |
401 | cout << "\n\nFit (no pad-filter):\n"; | |
402 | mChrg->Fit(lFit, "", "", 20, 200); | |
403 | ||
404 | new TCanvas("chrg2c", "Charge(pads<9)", 700, 500); | |
405 | mChrg2->Draw(); | |
406 | // Try landau-fit | |
407 | TF1* lFit2 = new TF1("lFit2", "landau", 20, 200); | |
408 | cout << "\n\nFit (pads<9):\n"; | |
409 | mChrg2->Fit(lFit2, "", "", 20, 200); | |
410 | ||
411 | new TCanvas("chrg3c", "Charge(2<pads<9)", 700, 500); | |
412 | mChrg3->Draw(); | |
413 | // Try landau-fit | |
414 | TF1* lFit3 = new TF1("lFit3", "landau", 20, 200); | |
415 | cout << "\n\nFit (2<pads<9):\n"; | |
416 | mChrg3->Fit(lFit3, "", "", 20, 200); | |
417 | ||
418 | #ifdef PRINT_OTHER | |
419 | new TCanvas(); | |
420 | mNPads->Draw(); | |
421 | ||
422 | new TCanvas(); | |
423 | mTime->Draw(); | |
424 | ||
425 | new TCanvas(); | |
426 | mX->Draw(); | |
427 | ||
428 | new TCanvas(); | |
429 | mY->Draw(); | |
430 | ||
431 | new TCanvas(); | |
432 | mZ->Draw(); | |
433 | ||
434 | new TCanvas(); | |
435 | mCenter->Draw(); | |
436 | ||
437 | new TCanvas(); | |
438 | mYX->Draw("colz"); | |
439 | ||
440 | new TCanvas(); | |
441 | mPlane->Draw(); | |
442 | ||
443 | new TCanvas(); | |
444 | mPlaneY->Draw("colz"); | |
445 | #endif | |
446 | ||
447 | new TCanvas(); | |
448 | gPad->SetLogz(); | |
449 | mQvsPads->Draw("colz"); | |
450 | ||
451 | new TCanvas(); | |
452 | gPad->SetLogz(); | |
453 | mQvsDetector->Draw("colz"); | |
454 | TLine* tl = 0; | |
455 | for (Int_t num = 30; num < 540; num += 30) | |
456 | { | |
457 | // Draw vertical lines every 30 bins to sketch the layer modules | |
458 | tl = new TLine(num - 0.5, 0, num - 0.5, 1000); | |
459 | tl->Draw(); | |
460 | } | |
461 | ||
462 | new TCanvas(); | |
463 | gPad->SetLogz(); | |
464 | mQvsSModId->Draw("colz"); | |
465 | ||
466 | TCanvas* sc = new TCanvas("sc","Supermodules", 800, 800); | |
467 | sc->SetLogz(); | |
468 | sc->Divide(2,2); | |
469 | sc->cd(1); | |
470 | mSMod0->Draw(); | |
471 | sc->cd(2); | |
472 | mSMod8->Draw(); | |
473 | sc->cd(3); | |
474 | mSMod9->Draw(); | |
475 | sc->cd(4); | |
476 | mSMod17->Draw(); | |
477 | ||
478 | ||
479 | TCanvas* padc = new TCanvas("padc","Pads", 1000, 1000); | |
480 | padc->Divide(4,4); | |
481 | for (Int_t ind = 0; ind < 14; ind++) | |
482 | { | |
483 | padc->cd(ind + 1); | |
484 | mPads[ind]->Draw(); | |
485 | } | |
486 | ||
487 | new TCanvas(); | |
488 | mClusters->Draw(); | |
489 | ||
490 | TCanvas* clusChamc = new TCanvas("clusChamc", "Clusters against detector with threshold", 800, 800); | |
491 | clusChamc->Divide(1,2); | |
492 | clusChamc->cd(1); | |
493 | mClusCham->Draw(); | |
494 | clusChamc->cd(2); | |
495 | mClusChamEx->Draw(); | |
496 | ||
497 | TCanvas* noisyc = new TCanvas("noisyc", "Col/Row of noisy detectors"); | |
498 | noisyc->Divide(1,2); | |
499 | noisyc->cd(1); | |
500 | gPad->SetLogz(); | |
501 | mColRow->Draw("colz"); | |
502 | noisyc->cd(2); | |
503 | gPad->SetLogz(); | |
504 | mColRowFiltered->Draw("colz"); | |
505 | ||
506 | TCanvas* ColTimec = new TCanvas("ColTimec", "Col/Time"); | |
507 | ColTimec->Divide(1,2); | |
508 | ColTimec->cd(1); | |
509 | gPad->SetLogz(); | |
510 | mColTime->Draw("colz"); | |
511 | ColTimec->cd(2); | |
512 | gPad->SetLogz(); | |
513 | mColTimeFiltered->Draw("colz"); | |
514 | ||
515 | cout << numPassed << " out of " << nEvents << " events passed the criterion (number of 'good' clusters >= "; | |
516 | cout << kThresholdLower << " && clusters <= " << kThresholdUpper << ").\n"; | |
517 | } |