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