TRD/AliTRDqaBlackEvents.cxx

   1 /**************************************************************************\r
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   3  *                                                                        *\r
   4  * Author: The ALICE Off-line Project.                                    *\r
   5  * Contributors are mentioned in the code where appropriate.              *\r
   6  *                                                                        *\r
   7  * Permission to use, copy, modify and distribute this software and its   *\r
   8  * documentation strictly for non-commercial purposes is hereby granted   *\r
   9  * withount fee, provided that the abov copyright notice appears in all   *\r
  10  * copies and that both the copyright notice and this permission notice   *\r
  11  * appear in the supporting documentation. The authors make no claims     *\r
  12  * about the suitability of this software for any purpose. It is          *\r
  13  * provided "as is" without express or implied warranty.                  *\r
  14  **************************************************************************/\r
  15 \r
  16 /* $Id: AliTRDqaBlackEvents.cxx 23387 2008-01-17 17:25:16Z cblume $ */\r
  17 \r
  18 ////////////////////////////////////////////////////////////////////////////\r
  19 //                                                                        //\r
  20 //  QA of black events                                                    //\r
  21 //                                                                        //\r
  22 //  Author:                                                               //\r
  23 //    Sylwester Radomski (radomski@physi.uni-heidelberg.de)               //\r
  24 //                                                                        //\r
  25 ////////////////////////////////////////////////////////////////////////////\r
  26 \r
  27 #include "TH1D.h"\r
  28 #include "TH2D.h"\r
  29 #include "TH2S.h"\r
  30 #include "TH3F.h"\r
  31 #include "TF1.h"\r
  32 #include "TFile.h"\r
  33 #include "TCanvas.h"\r
  34 #include "TPad.h"\r
  35 #include "TLatex.h"\r
  36 #include "TStyle.h"\r
  37 \r
  38 #include "AliTRDgeometry.h"\r
  39 #include "AliTRDrawStreamTB.h"\r
  40 #include "AliTRDqaBlackEvents.h"\r
  41 \r
  42 ClassImp(AliTRDqaBlackEvents)\r
  43 \r
  44 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
  45 \r
  46 AliTRDqaBlackEvents::AliTRDqaBlackEvents() \r
  47   :TObject() \r
  48   ,fnEvents(0)\r
  49   ,fCreateFull(0)\r
  50   ,fThresh(0)\r
  51   ,fCount(0)\r
  52   ,fOccupancy(0)\r
  53   ,fDetRob(0)\r
  54   ,fTBEvent(0)\r
  55   ,fFitType(0)\r
  56   ,fMinNoise(0.5)\r
  57   ,fMaxNoise(2) \r
  58 {\r
  59   //\r
  60   // Constructor \r
  61   // to create the histograms call Init()\r
  62   //\r
  63 }\r
  64 \r
  65 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
  66 \r
  67 AliTRDqaBlackEvents::AliTRDqaBlackEvents(const AliTRDqaBlackEvents &qa) \r
  68   :TObject(qa) \r
  69   ,fnEvents(0)\r
  70   ,fCreateFull(0)\r
  71   ,fThresh(0)\r
  72   ,fCount(0)\r
  73   ,fOccupancy(0)\r
  74   ,fDetRob(0)\r
  75   ,fTBEvent(0)\r
  76   ,fFitType(0)\r
  77   ,fMinNoise(0.5)\r
  78   ,fMaxNoise(2) \r
  79 {\r
  80   //\r
  81   // Copy constructor \r
  82   // to create the histograms call Init()\r
  83   //\r
  84 }\r
  85 \r
  86 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
  87 \r
  88 void AliTRDqaBlackEvents::Init() \r
  89 {\r
  90   //\r
  91   // creates histograms \r
  92   // \r
  93 \r
  94   //TFile *file = new \r
  95   //Info("Init", "Statring");\r
  96 \r
  97   fnEvents = 0;\r
  98 \r
  99   // histograms for chambers\r
 100   for(Int_t i=0; i<kDET; i++) {\r
 101     fNPoint[i]  = new TH2D(Form("entries_%d", i), "",  16, -0.5, 15.5, 144, -0.5, 143.5);\r
 102     fData[i]    = new TH3F(Form("data_%d", i), "", 16, -0.5, 15.5, 144, -0.5, 143.5, 50, -0.5, 49.5);\r
 103     fChPed[i]   = new TH2D(Form("ped_%d", i), "", 16, -0.5, 15.5, 144, -0.5, 143.5);\r
 104     fChNoise[i] = new TH2D(Form("noise_%d", i), "", 16, -0.5, 15.5, 144, -0.5, 143.5);\r
 105     fPed[i]     = new TH1D(Form("pedDist_%d", i), ";pedestals (ADC counts)", 100, 5, 15);\r
 106 \r
 107     fNoise[i]   = new TH1D(Form("noiseDist_%d", i), ";noise (ADC counts)", 100, 0, 5); \r
 108     fSignal[i]  = new TH1D(Form("signal_%d", i), ";signal (ADC counts)", 100, -0.5, 99.5);\r
 109 \r
 110     fnEntriesRM[i] = new TH2D(Form("entriesRM_%d", i), ";ROB,MCM", 8, -0.5, 7.5, 16, -0.5, 15.5);\r
 111   }\r
 112 \r
 113   // histogram for each MCM\r
 114   for(Int_t i=0; i < kDET * kROB * kMCM; i++)\r
 115     fFullCounter[i] = 0;\r
 116 \r
 117   // histograms from the whole detector\r
 118   fOccupancy = new TH1D("occupancy", "", 20, -0.5, 19.5);\r
 119   fDetRob    = new TH2D("DetRob", ";detector;ROB", kDET, -0.5, 539.5, 8, -0.5, 7.5);\r
 120   fTBEvent   = new TH2D("tbEvent", ";event ID;time bin", 100, -0.5, 99.5, 30, -0.5, 29.5);\r
 121 \r
 122   //Info("Init", "Done");\r
 123 }\r
 124 \r
 125 \r
 126 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
 127 \r
 128 void AliTRDqaBlackEvents::Reset() \r
 129 {\r
 130   //\r
 131   // Resets the histograms\r
 132   //\r
 133 \r
 134   for(Int_t i=0; i<kDET; i++) {\r
 135     fData[i]->Reset();\r
 136     fChPed[i]->Reset();\r
 137     fChNoise[i]->Reset();\r
 138   }\r
 139 }\r
 140 \r
 141 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
 142 \r
 143 Int_t AliTRDqaBlackEvents::AddEvent(AliTRDrawStreamTB *data) \r
 144 {\r
 145   //\r
 146   // Add an event\r
 147   //\r
 148 \r
 149   // structure to keep track if particular chanel is used\r
 150   Char_t isUsed[kDET][kCOL][kPAD]; \r
 151   for(Int_t i=0; i<kDET; i++)\r
 152     for(Int_t j=0; j<kCOL; j++)\r
 153       for(Int_t k=0; k<kPAD; k++)\r
 154         isUsed[i][j][k] = 0;\r
 155   \r
 156   Int_t nb = 0;\r
 157   Int_t rob_last = -1;\r
 158   Int_t mcm_last = -1;\r
 159   Int_t sm_01 = -1;\r
 160 \r
 161   while (data->Next()) {\r
 162 \r
 163     Int_t det = data->GetDet();\r
 164 \r
 165     Int_t row = data->GetRow();\r
 166     Int_t col = data->GetCol();\r
 167 \r
 168     Int_t rob = data->GetROB();\r
 169     Int_t mcm = data->GetMCM();\r
 170     Int_t adc = data->GetADC();\r
 171 \r
 172     Int_t *sig = data->GetSignals();\r
 173     nb++;\r
 174 \r
 175     // ugly hook\r
 176     if (det == 0 && row == 0 && col == 0) {\r
 177       if (sm_01 > -1) printf("\t\t\t!!! !!! second data set !!! !!!\n");\r
 178       sm_01++; \r
 179     }\r
 180 \r
 181     det += sm_01 * 30;   /// ugly\r
 182 \r
 183     // memory coruption protection\r
 184     if (det<0 || det>=kDET) continue;\r
 185     \r
 186     // check the ROBs\r
 187     fDetRob->Fill(det, rob, 1./(kMCM*18));\r
 188 \r
 189     isUsed[det][row][col]++;\r
 190 \r
 191     // check if mcm signal is continuus\r
 192     if ((rob_last != rob) || (mcm_last != mcm)) {\r
 193       rob_last = rob;\r
 194       mcm_last = mcm;\r
 195       fnEntriesRM[det]->Fill(rob,mcm);\r
 196     }\r
 197     \r
 198     // number of entries for each channels\r
 199     fNPoint[det]->Fill(row, col);\r
 200     \r
 201 \r
 202     // create a structure for an MCM if needed\r
 203     Int_t mcmIndex = det * (kMCM * kROB) + rob * kMCM + mcm;\r
 204     if (fCreateFull && !fFullSignal[mcmIndex])\r
 205       fFullSignal[mcmIndex] = new TH2S(Form("mcm_%d_%d_%d", det, rob, mcm), \r
 206                                        Form("mcm-%d-%d-%d;ADC;time bin", det, rob,mcm),\r
 207                                        21, -0.5, 20.5, 30, -0.5, 29.5);\r
 208     \r
 209 \r
 210     // loop over Time Bins and fill histograms\r
 211     for(Int_t k=0; k<kTB; k++) { /// to be corrected\r
 212 \r
 213       //if (sig[k] < 1) \r
 214       //printf("det = %d rob = %d mcm = %d adc = %d k = %d S = %d\n", det, rob, mcm, adc, k, sig[k]);\r
 215       \r
 216       fSignal[det]->Fill(sig[k]);\r
 217       fData[det]->Fill(row, col, sig[k]);\r
 218       \r
 219       // check if data strange enought\r
 220       if (fCreateFull && fFullSignal[mcmIndex]) {\r
 221         if (sig[k] > fThresh || sig[k] < 1) fFullCounter[mcmIndex]++;\r
 222         fFullSignal[mcmIndex]->Fill(adc, k, sig[k]);\r
 223       }\r
 224       \r
 225       // noisy chamber\r
 226       if (det == 29) {\r
 227         fTBEvent->Fill(fnEvents, k, sig[k]);\r
 228       }\r
 229 \r
 230     }\r
 231   }\r
 232   \r
 233   // is the dead-alive status changing during the run\r
 234   for(Int_t i=0; i<kDET; i++) {\r
 235     for(Int_t j=0; j<kCOL; j++)\r
 236       for(Int_t k=0; k<kPAD; k++)\r
 237         fOccupancy->Fill(isUsed[i][j][k]);\r
 238   }\r
 239 \r
 240   fnEvents++;\r
 241   return nb;\r
 242 }\r
 243 \r
 244 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
 245 \r
 246 void AliTRDqaBlackEvents::Process(const char *filename) \r
 247 {\r
 248   //\r
 249   // Process something\r
 250   //\r
 251 \r
 252   Int_t map[kDET];\r
 253   \r
 254   TH1D *hist = new TH1D("fitSignal", "", 50, -0.5, 49.5);\r
 255   TF1 *fit = new TF1("fit", "gaus(0)", 0, 20);\r
 256   fit->SetParameters(1e3, 10, 1);\r
 257     \r
 258   for(Int_t i=0; i<kDET; i++) {\r
 259     \r
 260     map[i] = 0;\r
 261     if (fData[i]->GetSum() < 10) continue;\r
 262     map[i] = 1;\r
 263 \r
 264     Info("process", "processing chamber %d", i);\r
 265 \r
 266     for(Int_t j=0; j<fData[i]->GetXaxis()->GetNbins(); j++) {\r
 267       for(Int_t k=0; k<fData[i]->GetYaxis()->GetNbins(); k++) {\r
 268         \r
 269         // project the histogramm\r
 270         hist->Reset();\r
 271         for(Int_t bb=0; bb<50; bb++) {\r
 272           Int_t dataBin = fData[i]->FindBin(j, k, bb);\r
 273           Double_t v = fData[i]->GetBinContent(dataBin);\r
 274           hist->SetBinContent(bb+1, v);\r
 275         }\r
 276 \r
 277         //TH1D *hist = fData[i]->ProjectionZ(Form("pad_%_%d_%d", i, j, k), j+1, j+1, k+1, k+1);\r
 278         \r
 279         Int_t bin = fChPed[i]->FindBin(j, k);\r
 280 \r
 281         if (hist->GetSum() > 1) {\r
 282           \r
 283           Double_t ped = 0, noise = 0;\r
 284 \r
 285           if (fFitType == 0) {\r
 286             fit->SetParameters(1e3, 10, 1);\r
 287             hist->Fit(fit, "q0", "goff", 0, 20);\r
 288             TF1 *f = hist->GetFunction("fit");\r
 289             ped = TMath::Abs(f->GetParameter(1));\r
 290             noise = TMath::Abs(f->GetParameter(2));\r
 291           } else {\r
 292             ped = hist->GetMean();\r
 293             noise = hist->GetRMS();\r
 294           }\r
 295 \r
 296           fChPed[i]->SetBinContent(bin, ped);\r
 297           fChNoise[i]->SetBinContent(bin, noise);\r
 298           \r
 299           fPed[i]->Fill(ped);\r
 300           fNoise[i]->Fill(noise);\r
 301 \r
 302         } else {\r
 303           fChPed[i]->SetBinContent(bin, 0);\r
 304           fChNoise[i]->SetBinContent(bin, 0);\r
 305         }\r
 306         \r
 307         //delete hist;\r
 308       }\r
 309     }\r
 310   }\r
 311 \r
 312   Info("Process", "Number of events = %d", fnEvents);\r
 313 \r
 314   // normalize number of entries histos\r
 315   Int_t max = 0;\r
 316   for(Int_t i=0; i<kDET; i++) { \r
 317     if (!map[i]) continue;\r
 318     for(Int_t j=0; j<fNPoint[i]->GetXaxis()->GetNbins(); j++) {\r
 319       for(Int_t k=0; k<fNPoint[i]->GetYaxis()->GetNbins(); k++) {\r
 320         Int_t dataBin = fNPoint[i]->FindBin(j, k);\r
 321         Double_t v = fNPoint[i]->GetBinContent(dataBin);\r
 322         if (v > max) max = (Int_t)v;\r
 323       }\r
 324     }\r
 325   }\r
 326   \r
 327   char entriesDistName[100];\r
 328   \r
 329   for(Int_t i=0; i<kDET; i++) {\r
 330     \r
 331     if (!map[i]) continue;\r
 332     \r
 333     sprintf(entriesDistName, "entriesDist_%d", i);\r
 334     fNPointDist[i] = new TH1D(entriesDistName, ";number of events", max+2, -0.5, max+1.5);\r
 335     \r
 336     for(Int_t j=0; j<fNPoint[i]->GetXaxis()->GetNbins(); j++) {\r
 337       for(Int_t k=0; k<fNPoint[i]->GetYaxis()->GetNbins(); k++) {\r
 338         Int_t dataBin = fNPoint[i]->FindBin(j, k);\r
 339         Double_t v = fNPoint[i]->GetBinContent(dataBin);\r
 340         //if (v > fnEvents) printf("N = %d V = %lf\n", fnEvents, v);\r
 341         fNPointDist[i]->Fill(v); \r
 342       }\r
 343     }\r
 344     \r
 345     fNPoint[i]->Scale(1./fnEvents);\r
 346   }\r
 347   \r
 348 \r
 349   for(Int_t i=0; i<kDET; i++) {\r
 350     fnEntriesRM[i]->SetMaximum(fnEvents * 1.5);\r
 351   }\r
 352 \r
 353   // save histograms\r
 354 \r
 355   TFile *file = new TFile(filename, "recreate");\r
 356   for(Int_t i=0; i<kDET; i++) {\r
 357     if (!map[i]) continue; \r
 358     fChPed[i]->Write();\r
 359     fChNoise[i]->Write();\r
 360     fNPoint[i]->Write();\r
 361     fNPointDist[i]->Write();\r
 362     fPed[i]->Write();\r
 363     fNoise[i]->Write();\r
 364     fSignal[i]->Write();\r
 365     fnEntriesRM[i]->Write();\r
 366   }\r
 367 \r
 368   Int_t nMcm = 0;\r
 369   for(Int_t i=0; i < kDET * kROB * kMCM; i++) {\r
 370     if (fFullSignal[i] && fFullCounter[i] > fCount) {\r
 371       fFullSignal[i]->Write();\r
 372       nMcm++;\r
 373     }\r
 374   }\r
 375 \r
 376   printf("Number of saved MCMs = %d\n", nMcm);\r
 377 \r
 378   fOccupancy->Write();\r
 379   fDetRob->Write();\r
 380   fTBEvent->Write();\r
 381   file->Close();\r
 382   delete file;\r
 383 }\r
 384 \r
 385 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
 386 \r
 387 void AliTRDqaBlackEvents::DrawChamber(const char *filename, Int_t det, Int_t w, Int_t h) \r
 388 {\r
 389   //\r
 390   // Draw raport for one chamber: \r
 391   // pedestal map, noise map, distribution of pedestal and noise\r
 392   // \r
 393   // input:\r
 394   // name of the file with histograms (created with Process())\r
 395   // detector Id (0 - 539)\r
 396   // \r
 397 \r
 398   // setup global style\r
 399   gStyle->SetPalette(1);\r
 400   gStyle->SetOptStat(0);\r
 401   gStyle->SetPadTopMargin(0.02);\r
 402   gStyle->SetPadBottomMargin(0.05);\r
 403 \r
 404   TFile *file = new TFile(filename, "READ");\r
 405 \r
 406   TCanvas *c = new TCanvas("blackEvents",Form("blackEvents %d",det), w, h);\r
 407   c->SetVertical(kFALSE);\r
 408   c->Divide(3,1, 0.01, 0.01);\r
 409   c->cd(3);\r
 410   \r
 411   TPad *mPad = (TPad*) gPad;\r
 412   mPad->Divide(1,2,0.01,0.01);\r
 413   \r
 414   c->cd(1);\r
 415   TH2D *h2 = (TH2D*)file->Get(Form("ped_%d",det));\r
 416   h2->SetMinimum(5);\r
 417   h2->SetMaximum(15);\r
 418   h2->SetTitle(";Z direction;#phi direction");\r
 419   h2->Draw("colz");\r
 420   \r
 421   c->cd(2);\r
 422   h2 = (TH2D*)file->Get(Form("noise_%d",det));\r
 423   h2->SetMinimum(fMinNoise);\r
 424   h2->SetMaximum(fMaxNoise);\r
 425   h2->SetTitle(";Z direction;#phi direction");\r
 426   h2->Draw("colz");\r
 427   \r
 428   mPad->cd(1);\r
 429   //gPad->SetLogy();\r
 430   TH1D *h1 = (TH1D*)file->Get(Form("pedDist_%d", det));\r
 431   h1->Draw();\r
 432   \r
 433   mPad->cd(2);\r
 434   gPad->SetLogy();\r
 435   h1 = (TH1D*)file->Get(Form("noiseDist_%d", det));\r
 436   h1->Draw();                    \r
 437   \r
 438   h1->Fit("gaus");\r
 439   TF1 *f = h1->GetFunction("gaus");\r
 440   const char *tt = Form("#mu = %.2f #sigma = %0.2f ", f->GetParameter(1),f->GetParameter(2));\r
 441   TLatex *ll = new TLatex(2, 100, tt);\r
 442   ll->SetTextSize(0.06);\r
 443   ll->Draw();\r
 444 }\r
 445 \r
 446 ///////////////////////////////////////////////////////////////////////////////////////////////////\r
 447 \r
 448 void AliTRDqaBlackEvents::DrawSm(const char *filename, Int_t sm, Int_t w, Int_t h) \r
 449 {\r
 450   //\r
 451   // ????????????\r
 452   //\r
 453   \r
 454   gStyle->SetPalette(1);\r
 455   gStyle->SetOptStat(0);\r
 456   \r
 457   gStyle->SetPadTopMargin(0.02);\r
 458   //gStyle->SetPadBottomMargin(0.05);  \r
 459   //gStyle->SetPadLeftMargin(0.02);  \r
 460   //gStyle->SetPadRightMargin(0.02);\r
 461 \r
 462   TFile *file = new TFile(filename, "READ");\r
 463 \r
 464   TCanvas *c = new TCanvas("blackEventsSM",Form("blackEvents SM %d",sm), w, h);\r
 465   c->SetVertical(kFALSE);\r
 466   c->Divide(5, 6, 0.001, 0.01);\r
 467   \r
 468   for(Int_t i=0; i<30; i++) {\r
 469     \r
 470     TH2D *h2 = (TH2D*)file->Get(Form("noise_%d",i+30*sm));\r
 471     if (!h2) continue;\r
 472     h2->SetMinimum(fMinNoise);\r
 473     h2->SetMaximum(fMaxNoise);\r
 474 \r
 475     // to be replaced by the official calculation\r
 476     Int_t stack = i/6;\r
 477     Int_t layer = i%6;\r
 478     Int_t index = (5-layer)*5 + stack + 1;\r
 479     //printf("%d %d %d %d\n", i, stack, layer, index);\r
 480     c->cd(index);\r
 481     gPad->SetBottomMargin(0.02);\r
 482     gPad->SetTopMargin(0.02);\r
 483 \r
 484     h2->Draw("col");\r
 485   }\r
 486 }\r
 487 \r
 488 ///////////////////////////////////////////////////////////////////////////////////////////////////\r