Coverity stuff once more ...
[u/mrichter/AliRoot.git] / TRD / AliTRDqaBlackEvents.cxx
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  * withount fee, provided thats the abov 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 /* $Id: AliTRDqaBlackEvents.cxx 23387 2008-01-17 17:25:16Z cblume $ */
17
18 ////////////////////////////////////////////////////////////////////////////
19 //                                                                        //
20 //  QA of black events                                                    //
21 //                                                                        //
22 //  Author:                                                               //
23 //    Sylwester Radomski (radomski@physi.uni-heidelberg.de)               //
24 //                                                                        //
25 ////////////////////////////////////////////////////////////////////////////
26
27 #include "TH1D.h"
28 #include "TH2D.h"
29 #include "TH2S.h"
30 #include "TF1.h"
31 #include "TFile.h"
32 #include "TCanvas.h"
33 #include "TPad.h"
34 #include "TLatex.h"
35 #include "TStyle.h"
36 #include "TGraph.h"
37
38 #include "AliRawReader.h"
39
40 #include "AliTRDrawStreamOld.h"
41 #include "AliTRDqaBlackEvents.h"
42
43 ClassImp(AliTRDqaBlackEvents)
44
45 ///////////////////////////////////////////////////////////////////////////////////////////////////
46
47 AliTRDqaBlackEvents::AliTRDqaBlackEvents() 
48   :TObject() 
49   ,fnEvents(0)
50   ,fCreateFull(0)
51   ,fThresh(0)
52   ,fCount(0)
53   ,fRefEv(0)
54   ,fOccupancy(0)
55   ,fDetRob(0)
56   ,fTBEvent(0)
57   ,fRefHistPed(0)
58   ,fRefHistNoise(0)
59   ,fErrorHC(0)
60   ,fErrorMCM(0)
61   ,fErrorADC(0)
62   ,fErrorSMHC(0)
63   ,fErrorSMMCM(0)
64   ,fErrorSMADC(0)
65   ,fErrorGraphHC(0)
66   ,fErrorGraphMCM(0)
67   ,fErrorGraphADC(0)
68   ,fGraphMCM(0)
69   ,fMcmTracks(0)
70   ,fMapMCM(0)
71   ,fFracMCM(0)
72   ,fSMHCped(0)
73   ,fSMHCerr(0)
74   ,fNoiseTotal(0)
75   ,fPP(0)
76   ,fMinNoise(0.5)
77   ,fMaxNoise(2)
78   ,fFitType(0) 
79 {
80   //
81   // Constructor 
82   // to create the histograms call Init()
83   //
84
85   for (Int_t i = 0; i < kDET; i++) {
86     fPed[i]            = 0x0;
87     fNoise[i]          = 0x0;
88     fChPP[i]           = 0x0;
89     fNPointDist[i]     = 0x0;
90     fChPed[i]          = 0x0;
91     fChNoise[i]        = 0x0;
92     fNPoint[i]         = 0x0;
93     fData[i]           = 0x0;
94     fSignal[i]         = 0x0;
95     fnEntriesRM[i]     = 0x0;
96     fnEntriesRMDist[i] = 0x0;
97     fChPedRes[i]       = 0x0;
98     fChNoiseRes[i]     = 0x0;
99     fErrorLocHC[i]     = 0x0;
100     fErrorLocMCM[i]    = 0x0;
101     fErrorLocADC[i]    = 0x0;
102   }
103   for (Int_t i = 0; i < 3; i++) {
104     fGraphPP[i]        = 0x0;
105     fSMLink[i]         = 0x0;
106     fGrLink[i]         = 0x0;
107     fppThresh[i]       = 0;
108     fnPP[i]            = 0;
109     fnLink[i]          = 0;
110   }
111   for (Int_t i = 0; i < 2; i++) {
112     fnErrorHC[i]       = 0;
113     fnErrorMCM[i]      = 0;
114     fnErrorADC[i]      = 0;
115   }
116   for (Int_t i = 0; i < kSM; i++) {
117     fSmNoiseRms[i]     = 0x0;
118     fSmNoiseFit[i]     = 0x0;
119     fSmPP[i]           = 0x0;
120   }
121   for (Int_t i = 0; i < kSM+1; i++) {
122     fNumberADC[i]      = 0x0;
123     fnADCinSM[i]       = 0;
124   }
125   for (Int_t i = 0; i < 1000; i++) {
126     fEvNoDist[i]       = 0;
127   }
128   strncpy(fRefFileName,"",256);
129
130 }
131
132 ///////////////////////////////////////////////////////////////////////////////////////////////////
133
134 AliTRDqaBlackEvents::AliTRDqaBlackEvents(const AliTRDqaBlackEvents &qa) 
135   :TObject(qa) 
136   ,fnEvents(0)
137   ,fCreateFull(0)
138   ,fThresh(0)
139   ,fCount(0)
140   ,fRefEv(0)
141   ,fOccupancy(0)
142   ,fDetRob(0)
143   ,fTBEvent(0)
144   ,fRefHistPed(0)
145   ,fRefHistNoise(0)
146   ,fErrorHC(0)
147   ,fErrorMCM(0)
148   ,fErrorADC(0)
149   ,fErrorSMHC(0)
150   ,fErrorSMMCM(0)
151   ,fErrorSMADC(0)
152   ,fErrorGraphHC(0)
153   ,fErrorGraphMCM(0)
154   ,fErrorGraphADC(0)
155   ,fGraphMCM(0)
156   ,fMcmTracks(0)
157   ,fMapMCM(0)
158   ,fFracMCM(0)
159   ,fSMHCped(0)
160   ,fSMHCerr(0)
161   ,fNoiseTotal(0)
162   ,fPP(0)
163   ,fMinNoise(0.5)
164   ,fMaxNoise(2) 
165   ,fFitType(0)
166 {
167   //
168   // Copy constructor 
169   // to create the histograms call Init()
170   //
171   
172   for (Int_t i = 0; i < kDET; i++) {
173     fPed[i]            = 0x0;
174     fNoise[i]          = 0x0;
175     fChPP[i]           = 0x0;
176     fNPointDist[i]     = 0x0;
177     fChPed[i]          = 0x0;
178     fChNoise[i]        = 0x0;
179     fNPoint[i]         = 0x0;
180     fData[i]           = 0x0;
181     fSignal[i]         = 0x0;
182     fnEntriesRM[i]     = 0x0;
183     fnEntriesRMDist[i] = 0x0;
184     fChPedRes[i]       = 0x0;
185     fChNoiseRes[i]     = 0x0;
186     fErrorLocHC[i]     = 0x0;
187     fErrorLocMCM[i]    = 0x0;
188     fErrorLocADC[i]    = 0x0;
189   }
190   for (Int_t i = 0; i < 3; i++) {
191     fGraphPP[i]        = 0x0;
192     fSMLink[i]         = 0x0;
193     fGrLink[i]         = 0x0;
194     fppThresh[i]       = 0;
195     fnPP[i]            = 0;
196     fnLink[i]          = 0;
197   }
198   for (Int_t i = 0; i < 2; i++) {
199     fnErrorHC[i]       = 0;
200     fnErrorMCM[i]      = 0;
201     fnErrorADC[i]      = 0;
202   }
203   for (Int_t i = 0; i < kSM; i++) {
204     fSmNoiseRms[i]     = 0x0;
205     fSmNoiseFit[i]     = 0x0;
206     fSmPP[i]           = 0x0;
207   }
208   for (Int_t i = 0; i < kSM+1; i++) {
209     fNumberADC[i]      = 0x0;
210     fnADCinSM[i]       = 0;
211   }
212   for (Int_t i = 0; i < 1000; i++) {
213     fEvNoDist[i]       = 0;
214   }
215   strncpy(fRefFileName,"",256);
216
217 }
218
219 ///////////////////////////////////////////////////////////////////////////////////////////////////
220
221 void AliTRDqaBlackEvents::Init() 
222 {
223   //
224   // creates histograms 
225   // 
226
227   //TFile *file = new 
228   //Info("Init", "Statring");
229
230   fnEvents = 0;
231
232   // histograms for chambers
233   for(Int_t det=0; det<kDET; det++) {
234
235     fNPoint[det]  = new TH2D(Form("entries_%d", det), "",  16, -0.5, 15.5, 144, -0.5, 143.5);
236     //fData[det]    = new TH3F(Form("data_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5, 50, -0.5, 49.5);
237
238     // pedestal noise maps using RMS and Fit
239     fChPed[det]   = new TH2D(Form("ped_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
240     fChNoise[det] = new TH2D(Form("noise_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
241     
242     //fChPed[det]   = new TH2D(Form("ped_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
243     //fChNoise[det] = new TH2D(Form("noise_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);    
244
245     // distribution per detector
246     fPed[det]     = new TH1D(Form("pedDist_%d", det), ";pedestals (ADC counts)", 100, 5, 15);
247     fNoise[det]   = new TH1D(Form("noiseDist_%d", det), ";noise (ADC counts)", 100, 0, 5); 
248     fSignal[det]  = new TH1D(Form("signal_%d", det), ";signal (ADC counts)", 100, -0.5, 99.5);
249     fChPP[det]    = new TH1D(Form("pp_%d", det), ";pp (ADC)", 200, -0.5, 199.5);
250
251     fnEntriesRM[det] = new TH2D(Form("entriesRM_%d", det), ";ROB,MCM", 8, -0.5, 7.5, 16, -0.5, 15.5);
252     
253     // histograms after reference subtraction
254     fChPedRes[det]    = new TH2D(Form("pedRef_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
255     fChNoiseRes[det]  = new TH2D(Form("noiseRef_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
256
257
258     // error codes
259     fErrorLocMCM[det] = new TH2D(Form("errorLocMCM_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
260     fErrorLocADC[det] = new TH2D(Form("errorLocADC_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);    
261     fErrorLocHC[det]  = new TH2D(Form("errorLocHC_%d", det), "", 16, -0.5, 15.5, 144, -0.5, 143.5);
262   }
263
264   // histogram for each MCM
265   for(Int_t i=0; i < kDET * kROB * kMCM; i++)
266     fFullCounter[i] = 0;
267
268   // histograms from the whole detector
269   fOccupancy = new TH1D("occupancy", "", 20, -0.5, 19.5);
270   fDetRob    = new TH2D("DetRob", ";detector;ROB", kDET, -0.5, 539.5, 8, -0.5, 7.5);
271   fTBEvent   = new TH2D("tbEvent", ";event ID;time bin", 100, -0.5, 99.5, 30, -0.5, 29.5);
272
273   // errors statistics and location
274   fErrorHC  = new TH1D("errorHC", ";error ID;", 18, -3.5, 14.5);
275   fErrorMCM = new TH1D("errorMCM", ";error ID;", 18, -3.5, 14.5);
276   fErrorADC = new TH1D("errorADC", ";error ID;", 18, -3.5, 14.5);
277   
278   fErrorSMHC  = new TH1D("errorSM_HC", ";SM id", 18, -0.5, 17.5);
279   fErrorSMMCM = new TH1D("errorSM_MCM", ";SM id", 18, -0.5, 17.5);
280   fErrorSMADC = new TH1D("errorSM_ADC", ";SM id", 18, -0.5, 17.5);
281
282
283   fErrorGraphHC  = new TGraph();
284   fErrorGraphMCM = new TGraph();
285   fErrorGraphADC = new TGraph();
286
287   fGraphMCM = new TGraph();
288
289   for(Int_t i=0; i<3; i++) {
290     fGraphPP[i] = new TGraph();
291   }
292
293
294   fMapMCM = new TH2D("mapMCM", ";det;mcm", 540, -0.5, 539.5, kROB*kMCM, -0.5, kROB*kMCM-0.5);
295   fFracMCM = new TH1D("fracMCM", ";frequency", 100, 0, 1);
296   
297
298   fErrorGraphHC->GetHistogram()->SetTitle("fraction of events with HC error;event number");
299   fErrorGraphMCM->GetHistogram()->SetTitle("fraction of events with MCM error;event number;");
300   fErrorGraphADC->GetHistogram()->SetTitle("fraction of events with ADC error;event number;"); 
301
302
303   fSMHCped = new TH2D("smHcPed", ";super module;half chamber", 18, -0.5, 17.5, 60, -0.5, 59.5);
304   
305   // link monitor
306   const char *linkName[3] = {"smLink", "smBeaf", "smData"};
307   const char *linkGrName[3] = {"grSmLink", "grSmBeaf", "grSmData"};
308   for(Int_t i=0; i<3; i++) {
309     fSMLink[i] = new TH2D(linkName[i], ";super module;link", 18, -0.5, 17.5, 60, -0.5, 59.5);
310     fGrLink[i] = new TGraph();
311     fGrLink[i]->SetName(linkGrName[i]);
312   }
313
314   //fZSsize = new TH1D("zssizeSingle", ";threshold;nADC", 40, -0.5, 39.5);
315   
316
317   //Info("Init", "Done");
318
319   // number of ADC channels fired per SM and in total
320   for(Int_t sm=0; sm<kSM+1; sm++)
321     fNumberADC[sm] = new TGraph();
322
323   //
324   fNoiseTotal = new TH1D("noiseTotal", "noise (ADC)", 250, 0, 10);
325   fPP = new TH1D("peakPeak", "p-p (ADC)", 200, -0.5, 199.5);
326
327   for(Int_t sm=0; sm<kSM; sm++) {
328     fSmNoiseRms[sm] = new TH1D(Form("noiseRms_sm%d", sm), ";noise from RMS (ADC)", 100, 0, 10);
329     fSmNoiseFit[sm] = new TH1D(Form("noiseFit_sm%d", sm), ";noise frim Fit (ADC)", 100, 0, 10);
330     fSmPP[sm] = new TH1D(Form("peakPeak_sm%d", sm), ";peak-peak (ADC)", 200, -0.5, 199.5); 
331   }
332
333   // event number consistancy
334   for(Int_t i=0; i<1000; i++) {
335     fEvNoDist[i] = new TH1D(Form("mcmEvDist_%d", i), ";#Delta Events", 201, -100.5, 100.5);
336   }
337   fRefEv = -1;
338
339   fMcmTracks = new TObjArray();
340
341   // clean data direct
342   for(Int_t det=0; det<kDET; det++) 
343     for(Int_t row=0; row<kROW; row++)
344       for(Int_t pad=0; pad<kPAD; pad++)
345         for(Int_t ch=0; ch<kCH; ch++) {
346           fDataDirect[det][row][pad][ch] = 0;
347           fSignalDirect[det][ch] = 0;  // overdone
348         }
349 }
350
351
352 ///////////////////////////////////////////////////////////////////////////////////////////////////
353
354 void AliTRDqaBlackEvents::Reset() 
355 {
356   //
357   // Resets the histograms
358   //
359
360   for(Int_t i=0; i<kDET; i++) {
361     //fData[i]->Reset();
362     fChPed[i]->Reset();
363     fChNoise[i]->Reset();
364   }
365 }
366
367
368 ///////////////////////////////////////////////////////////////////////////////////////////////////
369
370 void AliTRDqaBlackEvents::SetRefFile(const char *filename) {
371   
372   strncpy(fRefFileName,filename,256);
373
374 }
375
376 ///////////////////////////////////////////////////////////////////////////////////////////////////
377
378 void AliTRDqaBlackEvents::ReadRefHists(Int_t det) {
379   //  
380   // Read the reference histograms
381   //
382
383   fRefHistPed = 0;
384   fRefHistNoise = 0;
385   
386   TFile *file = 0;
387   if (fRefFileName) TFile::Open(fRefFileName);
388   if (!file) return;
389
390   fRefHistPed   = (TH2D*)file->Get(Form("ped_%d",det));
391   fRefHistNoise = (TH2D*)file->Get(Form("noise_%d", det));
392
393   if (file) file->Close();
394 }
395
396 ///////////////////////////////////////////////////////////////////////////////////////////////////
397
398 void AliTRDqaBlackEvents::StartEvent()
399 {
400   //
401   // start an event
402   //
403
404   // clear the mcm data
405   for(Int_t i=0; i < kDET * kROB * kMCM; i++) {
406     if (fFullSignal[i]) fFullSignal[i]->Reset();
407     fFullCounter[i] = 0;
408   }
409
410   for(Int_t i=0; i<2; i++) {
411     fnErrorHC[i] = 0;
412     fnErrorMCM[i] = 0;
413     fnErrorADC[i] = 0;
414   }
415
416  
417   Int_t ppThresh[3] = {10, 20, 40};
418   for(Int_t i=0; i<3; i++) {
419     fppThresh[i] = ppThresh[i];
420     fnPP[i] = 0;
421     fnLink[i] = 0;
422   }
423
424   for(Int_t sm=0; sm<kSM+1; sm++) fnADCinSM[sm] = 0;
425
426   if (fRefEv > 0) fRefEv++;
427   fEvNoDist[999]->Reset();  // keep only the last event
428
429 }
430
431 ///////////////////////////////////////////////////////////////////////////////////////////////////
432
433 void AliTRDqaBlackEvents::AddBuffer(AliTRDrawStreamOld *data, AliRawReader * const reader) 
434 {
435   
436   //printf ("try to read data\n");
437   Int_t nextBuff  = data->NextBuffer();
438   //printf("done ...\n");
439
440   if (nextBuff == 0) return;
441  
442   Int_t sm = reader->GetEquipmentId() - 1024;
443   //printf("reading SM %d\n", sm);
444   AliInfo(Form("reading SM %d", sm));
445   
446   if (sm < 0 || sm > 18) return;
447
448   // lopp over stacks, links ...
449
450   for (Int_t istack = 0; istack < 5; istack++) {        
451     for (Int_t ilink = 0; ilink < 12; ilink++) {
452       
453       //printf("HC = %d %d\n", istack, ilink);
454
455       Int_t det = sm * 30 + istack * 6 + ilink/2;       
456       
457       // check if data delivered
458       if (!(data->IsLinkActiveInStack(istack, ilink))) continue;
459       fSMLink[0]->Fill(sm, istack * 12 + ilink);
460       fnLink[0]++;
461
462       // check if beaf-beaf
463       if (data->GetLinkMonitorError(istack, ilink)) {
464         fSMLink[1]->Fill(sm, istack * 12 + ilink);
465         fnLink[1]++;
466         continue;
467       }
468         
469       // fill histogram with HC header errors
470       Int_t nErrHc = 0;
471       Int_t nErrHcTot = 0;
472       
473       nErrHc = FillBits(fErrorHC, data->GetH0ErrorCode(istack, ilink), 0);
474       if (!nErrHc) fErrorHC->Fill(-3);
475       nErrHcTot += nErrHc;
476       
477       nErrHc = FillBits(fErrorHC, data->GetH1ErrorCode(istack, ilink), 2);
478       if (!nErrHc) fErrorHC->Fill(-2);
479       nErrHcTot += nErrHc;
480       
481       nErrHc = FillBits(fErrorHC, data->GetHCErrorCode(istack, ilink), 4);
482       if (!nErrHc) fErrorHC->Fill(-1);
483       nErrHcTot += nErrHc;
484       
485       // trending
486       fnErrorHC[0]++;
487       if (nErrHcTot > 0) { 
488         fnErrorHC[1]++;
489         fErrorSMHC->Fill(sm);
490       }
491       
492       // data integrity protection 
493
494       //if (data->GetHCErrorCode(istack, ilink) > 0) continue;
495       if (data->GetH0ErrorCode(istack, ilink) > 0) continue;
496       if (data->GetH1ErrorCode(istack, ilink) > 0) continue;
497       
498       fSMLink[2]->Fill(sm, istack * 12 + ilink);        
499       fnLink[2]++;
500              
501       
502       for (Int_t imcm = 0; imcm < data->GetHCMCMmax(istack, ilink); imcm++ ){
503           
504         //printf("mcm = %d %d %d\n", istack, ilink, imcm);
505         
506         // fill MCM error code
507         
508         Int_t nErrMcm = 0;
509         Int_t nErrMcmTot = 0;
510         
511         nErrMcm = FillBits(fErrorMCM, data->GetMCMhdErrorCode(istack, ilink, imcm), 0);
512         if (!nErrMcm) fErrorMCM->Fill(-3);
513         nErrMcmTot += nErrMcm;
514         
515         nErrMcm = FillBits(fErrorMCM, data->GetMCMADCMaskErrorCode(istack, ilink, imcm), 5);
516         if (!nErrMcm) fErrorMCM->Fill(-2);
517         nErrMcmTot += nErrMcm;
518         
519         nErrMcm = FillBits(fErrorMCM, data->GetMCMErrorCode(istack, ilink, imcm), 10);
520         if (!nErrMcm) fErrorMCM->Fill(-1);
521         nErrMcmTot += nErrMcm;                    
522         
523         // trending
524         fnErrorMCM[0]++;
525         if (nErrMcmTot > 0) { 
526           fnErrorMCM[1]++;
527           fErrorSMMCM->Fill(sm);
528         }
529         
530         // MCM protection
531         if ( (data->GetMCMhdErrorCode(istack,ilink,imcm)) & 2 ) continue;
532         //if ((data->GetMCMADCMaskErrorCode(istack,ilink,imcm))) continue;
533         //if ((data->GetMCMErrorCode(istack,ilink,imcm))) continue;
534         
535         Int_t mcmEvent = data->GetEventNumber(istack, ilink, imcm);
536         
537         // set the reference event number 
538         if (fRefEv < 0) {
539           fRefEv = mcmEvent;
540           printf("Reference Event Number = %d (%d %d %d)\n", fRefEv, istack, ilink, imcm);
541         }
542         
543         // fill event distribution
544         if (!(fnEvents%10)) {
545           fEvNoDist[fnEvents/10]->Fill(mcmEvent - fRefEv);
546         }
547         
548         fEvNoDist[999]->Fill(mcmEvent - fRefEv);
549         
550         Int_t mcm = data->GetMCM(istack, ilink, imcm);
551         Int_t rob = data->GetROB(istack, ilink, imcm);
552
553         // create a structure for an MCM if needed
554         Int_t mcmIndex = det * (kMCM * kROB) + rob * kMCM + mcm;
555         if (fCreateFull && !fFullSignal[mcmIndex])
556           fFullSignal[mcmIndex] = 
557             new TH2S(Form("mcm_%d_%d_%d_%d_%d", sm, istack, ilink/2, rob, mcm), 
558                      Form("mcm-%d-%d-%d-%d-%d;ADC;time bin", sm, istack, ilink/2, rob, mcm),
559                      21, -0.5, 20.5, 30, -0.5, 29.5);
560         
561         
562         //Int_t zsADC[21][40];
563         /*
564           for(Int_t ina=0; ina<21; ina++) 
565           for(Int_t th=0; th<40; th++)
566           zsADC[ina][th] = 0;
567         */
568         
569         // first loop over ADC chanels 
570         
571         for (Int_t iadc=0; iadc < data->GetADCcount(istack, ilink, imcm); iadc++) {
572           
573           //printf("ADC = %d\n", iadc);
574           
575
576           // fill ADC error bits
577           Int_t nErrAdc = FillBits(fErrorADC, data->GetADCErrorCode(), 0);
578           if (!nErrAdc) fErrorADC->Fill(-1);
579           
580           fnErrorADC[0]++;
581           if (nErrAdc > 0) {
582             fnErrorADC[1]++;
583             fErrorSMADC->Fill(sm);
584           }
585           
586           // ADC protection
587           if ((data->GetADCErrorCode(istack,ilink,imcm,iadc))) continue;
588
589           Int_t minV = 1024;
590           Int_t maxV = 0;
591             
592           Int_t *sig = data->GetSignalDirect(istack, ilink, imcm, iadc);
593
594           //Int_t adc = data->GetADCnumber(istack, ilink, imcm, iadc);
595           Int_t row = data->GetRow(istack, ilink, imcm);
596           Int_t col = data->GetCol(istack, ilink, imcm, iadc);              
597             
598           // loop over Time Bins and fill histograms
599           for(Int_t k=0; k < data->GetNumberOfTimeBins(istack, ilink); k++) { 
600                       
601             //fSignal[det]->Fill(sig[k]);
602             //fData[det]->Fill(row, col, sig[k]); // slow
603
604             if (sig[k] < kCH) {
605               fSignalDirect[det][sig[k]]++;
606               fDataDirect[det][row][col][sig[k]]++; // direct data
607             }
608
609             // peak-peak
610             minV = (minV < sig[k]) ? minV : sig[k];
611             maxV = (maxV > sig[k]) ? maxV : sig[k];
612             
613             // check for active MCMs
614             if (fCreateFull && fFullSignal[mcmIndex]) {
615               if (sig[k] > fThresh || sig[k] < 0) fFullCounter[mcmIndex]++;
616               //if (sm == 0 && istack == 0 && ilink/2 == 1 && rob == 1 && mcm == 15) fFullCounter[mcmIndex]++; // special
617               //fFullSignal[mcmIndex]->Fill(adc, k, sig[k]); // slow
618             }
619             
620             // zero suppresion tests
621             /*
622               for(Int_t th=0; th<40; th++) {
623               if (sig[k] > th) {
624               zsADC[iadc][th] = 1;
625               if (iadc > 0) zsADC[iadc-1][th] = 1;
626               if (iadc < 10) zsADC[iadc+1][th] = 1;
627               }
628               }
629             */
630             
631           } // tb
632           
633           if (maxV > 0) { 
634             fnADCinSM[sm]++;
635             fnADCinSM[kSM]++;
636           }
637           
638           Int_t adcPP = maxV - minV;
639           //if (adcPP == 100) fFullCounter[mcmIndex] += 10;
640
641           fPP->Fill(adcPP);
642           fChPP[det]->Fill(adcPP);
643           fSmPP[sm]->Fill(adcPP);
644           
645           for(Int_t i=0; i<3; i++) {
646             if ((adcPP) > fppThresh[i]) fnPP[i]++;
647           }
648           
649             
650         } // adc
651
652           // fill ZS histos
653           /*
654             for(Int_t th=0; th<40; th++) {
655             Int_t nnADC = 0;
656             for(Int_t ins=0; ins<21; ins++) 
657             nnADC += zsADC[ins][th];
658             fZSsize->Fill(th, nnADC);
659             }
660           */
661
662         // fill active MCMs
663
664         if (fCreateFull && fFullSignal[mcmIndex] && (fFullCounter[mcmIndex] > fCount)) {
665           
666           for (Int_t iadc=0; iadc < data->GetADCcount(istack, ilink, imcm); iadc++) {
667             
668             // ADC protection
669             if ((data->GetADCErrorCode(istack,ilink,imcm,iadc))) continue;        
670             
671             //Int_t row = data->GetRow(istack, ilink, imcm);
672             //Int_t col = data->GetCol(istack, ilink, imcm, iadc);                  
673             Int_t adc = data->GetADCnumber(istack, ilink, imcm, iadc);
674             
675             Int_t *sig = data->GetSignalDirect(istack, ilink, imcm, iadc);
676             
677             // loop over Time Bins and fill histograms
678             for(Int_t k=0; k < data->GetNumberOfTimeBins(istack, ilink); k++) { 
679               fFullSignal[mcmIndex]->Fill(adc, k, sig[k]); // slow
680             }
681           } // tb
682         }
683         
684       } // mcm 
685     } // link
686   } // stack  
687   
688   // printf("end of loops\n");
689
690
691
692 ///////////////////////////////////////////////////////////////////////////////////////////////////
693
694 void AliTRDqaBlackEvents::FinishEvent()
695 {
696   //
697   // Processing at the end of the current event
698   //
699
700   for(Int_t i=0; i<3; i++) {
701     fGraphPP[i]->SetPoint(fnEvents, fnEvents, fnPP[i]);
702   }
703
704   // trend of the number of links
705   for(Int_t i=0; i<3; i++) {
706     fGrLink[i]->SetPoint(fnEvents, fnEvents, fnLink[i]);
707   }
708
709   // save interesting histos
710   Int_t mcmTrackCandidate = 0;
711   for(Int_t i = 0; i < kDET * kROB * kMCM; i++) { 
712     if ((fFullCounter[i] > fCount) && fFullSignal[i] && CheckMCM(i) )  {
713       
714       fMcmTracks->AddLast(fFullSignal[i]->Clone(Form("event_%d_%s", fnEvents, fFullSignal[i]->GetName())));
715       mcmTrackCandidate++;
716       
717       Int_t mcmTrackletDet = i/(kROB * kMCM);   
718       Int_t mcmTrackletMcm = i%(kROB * kMCM);
719       fMapMCM->Fill(mcmTrackletDet, mcmTrackletMcm);
720     }
721   }
722   
723   fGraphMCM->SetPoint(fnEvents, fnEvents, mcmTrackCandidate);
724   AliInfo(Form("Number of MCM track candidates = %d\n", mcmTrackCandidate));
725   
726   
727   // update fraction of error graphs
728   Double_t err;
729   
730   err = (fnErrorHC[0] > 0)? 100.*fnErrorHC[1]/fnErrorHC[0] : -1;
731   fErrorGraphHC->SetPoint(fnEvents, fnEvents, err);
732   
733   err = (fnErrorMCM[0] > 0)? 100.*fnErrorMCM[1]/fnErrorMCM[0] : -1;
734   fErrorGraphMCM->SetPoint(fnEvents, fnEvents, err);
735   
736   err = (fnErrorADC[0] > 0)? 100.*fnErrorADC[1]/fnErrorADC[0] : -1;
737   fErrorGraphADC->SetPoint(fnEvents, fnEvents, err);
738
739   // number of fired ADC per SM
740   for(Int_t sm=0; sm<kSM+1; sm++) 
741     fNumberADC[sm]->SetPoint(fnEvents, fnEvents, fnADCinSM[sm]);
742
743   fnEvents++;
744 }
745
746 ///////////////////////////////////////////////////////////////////////////////////////////////////
747
748 void AliTRDqaBlackEvents::Process(const char *filename) 
749 {
750   //
751   // Process something
752   //
753   
754   char fn[256];
755   strncpy(fn,filename,256);
756   
757   //AliInfo(Form("FILENAME = %s (%s)\n", filename, fn));
758
759   Int_t map[kDET];
760   
761   TH1D *hist = new TH1D("fitSignal", "", 50, -0.5, 49.5);
762   TF1 *fit = new TF1("fit", "gaus(0)", 0, 20);
763   fit->SetParameters(1e3, 10, 1);
764     
765   for(Int_t det=0; det<kDET; det++) {
766
767     //AliInfo(Form("processing chamber %d\n", det));   
768
769     map[det] = 0;
770     //if (fData[det]->GetSum() < 10) continue;
771     //if (fDataDirect[det][10][10][10] < 20) continue;
772     //map[det] = 1;
773
774
775     // rewrite signal-direct
776     for(Int_t ch=0; ch<kCH; ch++) {
777       fSignal[det]->Fill(ch, fSignalDirect[det][ch]);
778     }
779
780     // read reference distributions
781     ReadRefHists(det);
782
783     //for(Int_t row=0; row<fData[det]->GetXaxis()->GetNbins(); row++) {
784     //for(Int_t pad=0; pad<fData[det]->GetYaxis()->GetNbins(); pad++) {
785         
786     for(Int_t row=0; row<kROW; row++) {
787       for(Int_t pad=0; pad<kPAD; pad++) {
788
789         // project the histogramm
790         hist->Reset();
791         //for(Int_t bb=0; bb<50; bb++) {
792         for(Int_t bb=0; bb<kCH; bb++) {
793           //Int_t dataBin = fData[det]->FindBin(row, pad, bb);
794           //Double_t v = fData[det]->GetBinContent(dataBin);
795           hist->SetBinContent(bb+1, fDataDirect[det][row][pad][bb]);
796         }
797
798         Int_t bin = fChPed[det]->FindBin(row, pad);
799
800         if (hist->GetSum() > 1) {
801           
802           map[det] = 1;
803           Double_t ped = 0, noise = 0;
804
805           if (fFitType == 0) {
806             fit->SetParameters(1e3, 10, 1);
807             hist->Fit(fit, "q0", "goff", 0, 20);
808             TF1 *f = hist->GetFunction("fit");
809             ped = TMath::Abs(f->GetParameter(1));
810             noise = TMath::Abs(f->GetParameter(2));
811             fSmNoiseFit[det/30]->Fill(noise);
812           } else {
813             ped = hist->GetMean();
814             noise = hist->GetRMS();
815             fSmNoiseRms[det/30]->Fill(noise);
816             //if (pad == 0)
817             //  AliInfo(Form("data %f %f %f\n", hist->GetSum(), ped, noise));
818           }
819
820           fChPed[det]->SetBinContent(bin, ped);
821           fChNoise[det]->SetBinContent(bin, noise);
822           fNoiseTotal->Fill(noise);
823
824           // subtract reference values
825           Double_t refped = 0;
826           Double_t refnoise = 0;
827           
828           if (fRefHistPed)   refped   = fRefHistPed->GetBinContent(bin);
829           if (fRefHistNoise) refnoise = fRefHistPed->GetBinContent(bin);
830
831           fChPedRes[det]->SetBinContent(bin, ped-refped);
832           fChNoiseRes[det]->SetBinContent(bin, noise-refnoise);
833           
834           fPed[det]->Fill(ped);
835           fNoise[det]->Fill(noise);
836
837           // fill SM-HC plot
838           Int_t sm = det / 30;
839           Int_t hc = (pad < kPAD/2) ? 2* (det % 30) : 2* (det % 30) + 1;
840           if (ped > 9. && ped < 11) fSMHCped->Fill(sm, hc, 1./1152.); // number of pads in HC
841
842         } else {
843           
844           // not enought data found 
845           fChPed[det]->SetBinContent(bin, 0);
846           fChNoise[det]->SetBinContent(bin, 0);
847           fChPedRes[det]->SetBinContent(bin, 0);
848           fChNoiseRes[det]->SetBinContent(bin, 0);
849         }
850         
851         //delete hist;
852       }
853     }
854   }
855
856
857   //AliInfo(Form("Number of events = %d\n", fnEvents));
858
859   // normalize number of entries histos
860   Int_t max = 0;
861   for(Int_t i=0; i<kDET; i++) { 
862     if (!map[i]) continue;
863     for(Int_t j=0; j<fNPoint[i]->GetXaxis()->GetNbins(); j++) {
864       for(Int_t k=0; k<fNPoint[i]->GetYaxis()->GetNbins(); k++) {
865         Int_t dataBin = fNPoint[i]->FindBin(j, k);
866         Double_t v = fNPoint[i]->GetBinContent(dataBin);
867         if (v > max) max = (Int_t)v;
868       }
869     }
870   }
871   
872   char entriesDistName[100];
873
874   for(Int_t i=0; i<kDET; i++) {
875     
876     if (!map[i]) continue;
877     
878     sprintf(entriesDistName, "entriesDist_%d", i);
879     fNPointDist[i] = new TH1D(entriesDistName, ";number of events", max+2, -0.5, max+1.5);
880     
881     for(Int_t j=0; j<fNPoint[i]->GetXaxis()->GetNbins(); j++) {
882       for(Int_t k=0; k<fNPoint[i]->GetYaxis()->GetNbins(); k++) {
883         Int_t dataBin = fNPoint[i]->FindBin(j, k);
884         Double_t v = fNPoint[i]->GetBinContent(dataBin);
885         //if (v > fnEvents) AliInfo(Form("N = %d V = %lf\n", fnEvents, v));
886         fNPointDist[i]->Fill(v); 
887       }
888     }
889     
890     fNPoint[i]->Scale(1./fnEvents);
891   }
892   
893
894   for(Int_t i=0; i<kDET; i++) {
895     fnEntriesRM[i]->SetMaximum(fnEvents * 1.5);
896   }
897
898   // save histograms
899
900   //AliInfo(Form("FILENAME 2 = %s (%d)\n", fn, fn));
901   TFile *file = new TFile(fn, "recreate");
902   for(Int_t det = 0; det < kDET; det++) {
903     if (!map[det]) continue; 
904     fChPed[det]->Write();
905     fChNoise[det]->Write();
906     fNPoint[det]->Write();
907     fNPointDist[det]->Write();
908     fPed[det]->Write();
909     fNoise[det]->Write();
910     fSignal[det]->Write();
911     fnEntriesRM[det]->Write();
912     fChPP[det]->Write();
913
914     fChPedRes[det]->Write();
915     fChNoiseRes[det]->Write();
916
917     // save error hists
918     fErrorLocMCM[det]->SetMinimum(0);
919     fErrorLocMCM[det]->SetMaximum(fnEvents);
920     fErrorLocMCM[det]->Write();
921
922     fErrorLocADC[det]->SetMinimum(0);
923     fErrorLocADC[det]->SetMaximum(fnEvents);
924     fErrorLocADC[det]->Write();
925   }
926
927   for(Int_t sm=0; sm<kSM; sm++) {
928     fSmNoiseRms[sm]->Write();
929     fSmNoiseFit[sm]->Write();
930     fSmPP[sm]->Write();
931   }
932
933
934
935   Int_t nMcm = 0;
936   for(Int_t i=0; i < kDET * kROB * kMCM; i++) {
937     if (fFullSignal[i] && fFullCounter[i] > fCount) {
938       fFullSignal[i]->Write();
939       nMcm++;
940     }
941   }
942
943   AliInfo(Form("Number of saved MCMs = %d\n", nMcm));
944   
945   fMcmTracks->Write();
946   AliInfo(Form("Number of tracks = %d\n", fMcmTracks->GetEntries()));
947   
948   // permanently problematic MCMs
949   for(Int_t det=0; det<kDET; det++) {
950     for(Int_t mcm=0; mcm<kROB*kMCM; mcm++) {
951       
952       Int_t mRob = mcm / kMCM;
953       Int_t mMcm = mcm % kMCM;
954       Int_t bin = fMapMCM->FindBin(det, mcm);
955       Double_t frac = 1. * fMapMCM->GetBinContent(bin) / fnEvents;      
956       fFracMCM->Fill(frac);
957       
958       if (frac > 0.7) {
959         AliInfo(Form("{%d, %d, %d, %f}, \n", det, mRob, mMcm, frac));
960       }      
961     }
962   }
963
964
965   fOccupancy->Write();
966   fDetRob->Write();
967   fTBEvent->Write();
968   
969   // error hists
970   fErrorHC->Write();
971   fErrorMCM->Write();
972   fErrorADC->Write();
973
974   fErrorSMHC->Write();
975   fErrorSMMCM->Write();
976   fErrorSMADC->Write();  
977   
978   // write graphs
979   fErrorGraphHC->Write("trendErrorHC");
980   fErrorGraphMCM->Write("trendErrorMCM");
981   fErrorGraphADC->Write("trendErrorADC");
982   
983   fGraphMCM->Write("trendMCM");
984   
985   for(Int_t i=0; i<3; i++) {
986     fGraphPP[i]->Write(Form("fracPP_%d", i));
987   }
988   
989   //fZSsize->Scale(1./fnEvents);
990   //fZSsize->Write();
991
992   fMapMCM->SetMaximum(fnEvents);
993   fMapMCM->Write();
994   fFracMCM->Write();
995   
996   fSMHCped->Write();
997
998   for(Int_t i=0; i<3; i++ ) {
999     fSMLink[i]->Write();
1000     fGrLink[i]->Write();
1001   }
1002
1003   for(Int_t sm=0; sm<kSM; sm++)
1004     fNumberADC[sm]->Write(Form("nADCinSM%d",sm));
1005   
1006   fNumberADC[kSM]->Write("nADCinEvent");
1007
1008   fNoiseTotal->Write();
1009   fPP->Write();
1010
1011   for(Int_t i=0; i<1000; i++) {
1012     if (fEvNoDist[i]->GetSum() > 0) fEvNoDist[i]->Write();
1013   }
1014
1015
1016   file->Close();
1017   delete file;
1018 }
1019
1020 ///////////////////////////////////////////////////////////////////////////////////////////////////
1021
1022 Int_t AliTRDqaBlackEvents::CheckMCM(Int_t index) const {
1023   //
1024   // Checks a single MCM
1025   //  
1026
1027   return 1;
1028   
1029   static Int_t data[21][3] = {
1030     {1, 0, 1}, 
1031     {242, 0, 0}, 
1032     {242, 0, 1}, 
1033     {242, 0, 2}, 
1034     {242, 0, 4}, 
1035     {242, 0, 5}, 
1036     {242, 0, 6}, 
1037     {242, 0, 8}, 
1038     {242, 0, 12}, 
1039     {251, 7, 7}, 
1040     {254, 3, 11}, 
1041     {259, 3, 14}, 
1042     {260, 1, 9}, 
1043     {260, 3, 15}, 
1044     {273, 1, 7}, 
1045     {273, 1, 15}, 
1046     {276, 5, 11}, 
1047     {280, 6, 2}, 
1048     {299, 6, 4}, 
1049     {511, 2, 9}, 
1050     {517, 7, 15}
1051   };
1052   
1053   for(Int_t i=0; i<21; i++) {
1054     Int_t wIndex = data[i][0] * kROB*kMCM + data[i][1] * kMCM + data[i][2];
1055     if (index == wIndex) return 0;
1056   }
1057
1058   return 1;
1059 }
1060
1061 ///////////////////////////////////////////////////////////////////////////////////////////////////
1062
1063
1064
1065
1066
1067
1068
1069 void AliTRDqaBlackEvents::DrawChamber(const char *filename, Int_t det, Int_t w, Int_t h) 
1070 {
1071   //
1072   // Draw raport for one chamber: 
1073   // pedestal map, noise map, distribution of pedestal and noise
1074   // 
1075   // input:
1076   // name of the file with histograms (created with Process())
1077   // detector Id (0 - 539)
1078   // 
1079
1080   // setup global style
1081   gStyle->SetPalette(1);
1082   gStyle->SetOptStat(0);
1083   gStyle->SetPadTopMargin(0.02);
1084   gStyle->SetPadBottomMargin(0.05);
1085
1086   TFile *file = new TFile(filename, "READ");
1087
1088   TCanvas *c = new TCanvas("blackEvents",Form("blackEvents %d",det), w, h);
1089   c->SetVertical(kFALSE);
1090   c->Divide(3,1, 0.01, 0.01);
1091   c->cd(3);
1092   
1093   TPad *mPad = (TPad*) gPad;
1094   mPad->Divide(1,2,0.01,0.01);
1095   
1096   c->cd(1);
1097   TH2D *h2 = (TH2D*)file->Get(Form("ped_%d",det));
1098   h2->SetMinimum(5);
1099   h2->SetMaximum(15);
1100   h2->SetTitle(";Z direction;#phi direction");
1101   h2->Draw("colz");
1102   
1103   c->cd(2);
1104   h2 = (TH2D*)file->Get(Form("noise_%d",det));
1105   h2->SetMinimum(fMinNoise);
1106   h2->SetMaximum(fMaxNoise);
1107   h2->SetTitle(";Z direction;#phi direction");
1108   h2->Draw("colz");
1109   
1110   mPad->cd(1);
1111   //gPad->SetLogy();
1112   TH1D *h1 = (TH1D*)file->Get(Form("pedDist_%d", det));
1113   h1->Draw();
1114   
1115   mPad->cd(2);
1116   gPad->SetLogy();
1117   h1 = (TH1D*)file->Get(Form("noiseDist_%d", det));
1118   h1->Draw();                    
1119   
1120   h1->Fit("gaus");
1121   TF1 *f = h1->GetFunction("gaus");
1122   const char *tt = Form("#mu = %.2f #sigma = %0.2f ", f->GetParameter(1),f->GetParameter(2));
1123   TLatex *ll = new TLatex(2, 100, tt);
1124   ll->SetTextSize(0.06);
1125   ll->Draw();
1126 }
1127
1128 ///////////////////////////////////////////////////////////////////////////////////////////////////
1129
1130 void AliTRDqaBlackEvents::DrawSm(const char *filename, Int_t sm, Int_t w, Int_t h) 
1131 {
1132   //
1133   // ????????????
1134   //
1135   
1136   gStyle->SetPalette(1);
1137   gStyle->SetOptStat(0);
1138   
1139   gStyle->SetPadTopMargin(0.02);
1140   //gStyle->SetPadBottomMargin(0.05);  
1141   //gStyle->SetPadLeftMargin(0.02);  
1142   //gStyle->SetPadRightMargin(0.02);
1143
1144   TFile *file = new TFile(filename, "READ");
1145
1146   TCanvas *c = new TCanvas("blackEventsSM",Form("blackEvents SM %d",sm), w, h);
1147   c->SetVertical(kFALSE);
1148   c->Divide(5, 6, 0.001, 0.01);
1149   
1150   for(Int_t i=0; i<30; i++) {
1151     
1152     TH2D *h2 = (TH2D*)file->Get(Form("noise_%d",i+30*sm));
1153     if (!h2) continue;
1154     h2->SetMinimum(fMinNoise);
1155     h2->SetMaximum(fMaxNoise);
1156
1157     // to be replaced by the official calculation
1158     Int_t stack = i/6;
1159     Int_t layer = i%6;
1160     Int_t index = (5-layer)*5 + stack + 1;
1161     //AliInfo(Form("%d %d %d %d\n", i, stack, layer, index));
1162     c->cd(index);
1163     gPad->SetBottomMargin(0.02);
1164     gPad->SetTopMargin(0.02);
1165
1166     h2->Draw("col");
1167   }
1168 }
1169
1170 ///////////////////////////////////////////////////////////////////////////////////////////////////
1171
1172 Int_t AliTRDqaBlackEvents::FillBits(TH1D *hist, Int_t code, Int_t offset) {
1173   //
1174   // Fill bits
1175   //
1176
1177   Int_t nb = 0;
1178   UInt_t test = 1;
1179   for(Int_t i=0; i<8; i++) {
1180     if (code & test) {
1181       hist->Fill(i+offset);
1182       nb++;
1183     }
1184     test *= 2;       
1185   }
1186   
1187   return nb;
1188 }
1189
1190 ///////////////////////////////////////////////////////////////////////////////////////////////////