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