EMCAL/AliEMCALQAChecker.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  * without fee, provided that the above 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 /*
  17   Checks the quality assurance.
  18   By comparing with reference data
  19
  20   Based on PHOS code written by
  21   Y. Schutz CERN July 2007
  22
  23  For the moment we only implement the checking of raw data QA.
  24  The checked for ESD and RecPoints will be implemented later.
  25
  26
  27  */
  28
  29
  30 // --- ROOT system ---
  31 #include <TClass.h>
  32 #include <TH1.h>
  33 #include <TF1.h>
  34 #include <TH1I.h>
  35 #include <TH2F.h>
  36 #include <TIterator.h>
  37 #include <TKey.h>
  38 #include <TFile.h>
  39 #include <TLine.h>
  40 #include <TText.h>
  41 #include <TPaveText.h>
  42 #include <TMath.h>
  43
  44 // --- Standard library ---
  45
  46 // --- AliRoot header files ---
  47 #include "AliLog.h"
  48 #include "AliCDBManager.h"
  49 #include "AliCDBEntry.h"
  50 #include "AliQAv1.h"
  51 #include "AliQAChecker.h"
  52 #include "AliQAThresholds.h"
  53 #include "AliEMCALQAChecker.h"
  54
  55 ClassImp(AliEMCALQAChecker)
  56
  57 //__________________________________________________________________
  58 AliEMCALQAChecker::AliEMCALQAChecker() :
  59 AliQACheckerBase("EMCAL","EMCAL Quality Assurance Data Maker"),
  60 fTextSM(new TText*[fgknSM]),
  61 fLineCol(new TLine(47.5,-0.5,47.5,119.5)),
  62 fText(new TPaveText(0.2,0.7,0.8,0.9,"NDC"))
  63 {
  64   // ctor
  65   fLineCol->SetLineColor(1);
  66   fLineCol->SetLineWidth(2);
  67
  68   fTextSM[0]= new TText(20, 12, "SM A0");
  69   fTextSM[1]= new TText(20, 36, "SM A1");
  70   fTextSM[2]= new TText(20, 60, "SM A2");
  71   fTextSM[3]= new TText(20, 84, "SM A3");
  72   fTextSM[4]= new TText(20, 108,"SM A4");
  73   fTextSM[5]= new TText(20, 132,"SM A5");
  74
  75   fTextSM[6]= new TText(64, 12, "SM C0");
  76   fTextSM[7]= new TText(64, 36, "SM C1");
  77   fTextSM[8]= new TText(64, 60, "SM C2");
  78   fTextSM[9]= new TText(64, 84, "SM C3");
  79   fTextSM[10]= new TText(64, 108,"SM C4");
  80   fTextSM[11]= new TText(64, 132,"SM C5");
  81
  82         for(int i = 0; i < 5; i++) {
  83                 fLineRow[i] = new TLine(-0.5,23.5+(24*i),95.5,23.5+(24*i));
  84                 fLineRow[i]->SetLineColor(1);
  85                 fLineRow[i]->SetLineWidth(2);
  86         }
  87
  88         for(int i = 0; i < 3; i++) {
  89                 fTextL1[i] = new TPaveText(0.2,0.8,0.8,0.9,"NDC");
  90         }
  91
  92
  93 }
  94
  95 //__________________________________________________________________
  96 AliEMCALQAChecker::~AliEMCALQAChecker()
  97 {
  98         /// dtor
  99   delete [] fTextSM ;
 100   delete fLineCol ;
 101   for (Int_t i=0; i<5; ++i) delete fLineRow[i] ;
 102   delete fText  ;
 103 }
 104
 105 //______________________________________________________________________________
 106 void AliEMCALQAChecker::Check(Double_t * test, AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/)
 107 {
 108         /// Check objects in list
 109
 110         if ( index == AliQAv1::kRAW )
 111         {
 112     CheckRaws(test, list);
 113     //printf ("checkers for task %d \n", index) ;
 114         }
 115
 116         if ( index == AliQAv1::kREC)
 117         {
 118     CheckRecPoints(test, list);
 119         }
 120
 121         if ( index == AliQAv1::kESD )
 122         {
 123     CheckESD(test, list);
 124         }
 125         //AliWarning(Form("Checker for task %d not implement for the moment",index));
 126 }
 127
 128 //______________________________________________________________________________
 129 TH1*
 130 AliEMCALQAChecker::GetHisto(TObjArray* list, const char* hname, Int_t specie) const
 131 {
 132         /// Get a given histo from the list
 133         TH1* h = static_cast<TH1*>(list->FindObject(Form("%s_%s",AliRecoParam::GetEventSpecieName(specie),hname)));
 134         if (!h)
 135         {
 136                 AliError(Form("Did not find expected histo %s",hname));
 137         }
 138         return h;
 139 }
 140
 141 //______________________________________________________________________________
 142 Double_t
 143 AliEMCALQAChecker::MarkHisto(TH1& histo, Double_t value) const
 144 {
 145         /// Mark histo as originator of some QA error/warning
 146
 147         if ( value != 1.0 )
 148         {
 149                 histo.SetBit(AliQAv1::GetQABit());
 150         }
 151
 152         return value;
 153 }
 154
 155
 156 //______________________________________________________________________________
 157 void AliEMCALQAChecker::CheckRaws(Double_t * test, TObjArray ** list)
 158 {
 159   //  Check RAW QA histograms
 160   //  -- Yaxian Mao, CCNU/CERN/LPSC
 161   //adding new checking method: 25/04/2010, Yaxian Mao
 162   //Comparing the amplitude from current run to the reference run, if the ratio in the range [0.8, .12], count as a good tower.
 163   //If more than 90% towers are good, EMCAL works fine, otherwise experts should be contacted.
 164
 165
 166   // Retrieve the thresholds
 167   Float_t ratioThresh = 0.9;   // threshold for calibration ratio = good towers/all towers (default 0.9)
 168   Float_t ThreshG     = 0.5;   // threshold for L1 Gamma triggers (default 0.5)
 169   Float_t ThreshJ     = 0.5;   // threshold for L1 Jet triggers (default 0.5)
 170   Int_t badLinkThresh = 1;     // threshold for bad links (default 1)
 171
 172   AliCDBManager* man = AliCDBManager::Instance();
 173   if(man){
 174     AliCDBEntry* entry = man->Get("GRP/Calib/QAThresholds");
 175     if(entry){
 176       TObjArray* branch = (TObjArray*) entry->GetObject();
 177       if(branch){
 178         AliQAThresholds* thresholds = (AliQAThresholds*) branch->FindObject("EMC");
 179         if(thresholds){
 180           TParameter<float>* ParamR  = (TParameter<float>*) thresholds->GetThreshold(0);
 181           TParameter<float>* ParamG  = (TParameter<float>*) thresholds->GetThreshold(1);
 182           TParameter<float>* ParamJ  = (TParameter<float>*) thresholds->GetThreshold(2);
 183           TParameter<int>* ParamL  = (TParameter<int>*) thresholds->GetThreshold(3);
 184           if(ParamR)
 185             ratioThresh = ParamR->GetVal();
 186           if(ParamG)
 187             ThreshG = ParamG->GetVal();
 188           if(ParamJ)
 189             ThreshJ = ParamJ->GetVal();
 190           if(ParamL)
 191             badLinkThresh = ParamL->GetVal();
 192         }
 193       }
 194     }
 195   }
 196
 197   Int_t nTowersPerSM = 24*48; // number of towers in a SuperModule; 24x48
 198   Double_t nTot = fgknSM * nTowersPerSM ;
 199   TList *lstF = 0;
 200   Int_t calibSpecieId = (Int_t)TMath::Log2( AliRecoParam::kCalib );
 201   for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
 202
 203     test[specie] = 0.0 ;
 204     if ( !AliQAv1::Instance()->IsEventSpecieSet(specie)) continue ;
 205     if (list[specie]->GetEntries() == 0) test[specie] = 0. ; // nothing to check
 206     else {
 207       //get calib histos
 208       TH2F * hdata  = (TH2F*)list[specie]->At(k2DRatioAmp) ;
 209       TH1F * ratio = (TH1F*)list[specie]->At(kRatioDist) ;
 210
 211       //get L1 histos
 212       TH2F *hL1GammaPatch = (TH2F*)list[specie]->At(kGL1);
 213       TH2F *hL1JetPatch = (TH2F*)list[specie]->At(kJL1);
 214       TH1I *hFrameR = (TH1I*)list[specie]->At(kSTUTRU);
 215
 216
 217       // =======================================================================================
 218       //calib histo checker first:
 219       if( hdata && ratio ){
 220
 221         // first clean lines, text (functions)
 222         lstF = hdata->GetListOfFunctions();
 223         CleanListOfFunctions(lstF);
 224         lstF = ratio->GetListOfFunctions();
 225         CleanListOfFunctions(lstF);
 226
 227         if(hdata->GetEntries()!=0 && ratio->GetEntries()!=0) {
 228
 229           lstF = hdata->GetListOfFunctions();
 230
 231           //adding the lines to distinguish different SMs
 232           lstF->Add(fLineCol->Clone());
 233           for(Int_t iLine = 0; iLine < 5; iLine++) {
 234             lstF->Add(fLineRow[iLine]->Clone());
 235           }
 236           //Now adding the text to for each SM
 237           for(Int_t iSM = 0 ; iSM < fgknSM ; iSM++){  //number of SMs loop start
 238             lstF->Add(fTextSM[iSM]->Clone());
 239           }
 240
 241
 242           //now check the ratio histogram
 243           lstF = ratio->GetListOfFunctions();
 244
 245           Double_t binContent = 0. ;
 246           Int_t nGoodTower = 0 ;
 247           Double_t rv = 0. ;
 248           for(Int_t ix = 1; ix <= hdata->GetNbinsX(); ix++) {
 249             for(Int_t iy = 1; iy <= hdata->GetNbinsY(); iy++) {
 250               binContent = hdata->GetBinContent(ix, iy) ;
 251                                 if (binContent < 1.2 && binContent > 0.8) nGoodTower++ ;
 252             }
 253           }
 254           rv = nGoodTower/nTot ;
 255           //printf("%2.2f %% towers out of range [0.8, 1.2]\n", (1-rv)*100);
 256           if(fText){
 257             lstF->Add(fText->Clone()) ;
 258             fText->Clear() ;
 259
 260             fText->AddText(Form("%2.2f %% towers out of range [0.8, 1.2]", (1-rv)*100));
 261             if (rv < ratioThresh) {
 262               test[specie] = ratioThresh;
 263               // 2 lines text info for quality
 264               fText->SetFillColor(2) ;
 265               fText->AddText(Form("EMCAL = NOK, CALL EXPERTS!!!"));
 266             }
 267             else {
 268               test[specie] = 1 - ratioThresh;
 269               fText->SetFillColor(3) ;
 270               fText->AddText(Form("EMCAL = OK, ENJOY..."));
 271             }
 272           }//fText
 273         }//calib histo checking done
 274       }//histograms NOT NULL
 275
 276       // ========================================================================================
 277       // now L1 checks:
 278
 279       if( hL1GammaPatch ){
 280
 281         // first clean lines, text (functions)
 282         lstF = hL1GammaPatch->GetListOfFunctions();
 283         CleanListOfFunctions(lstF);
 284
 285         if (specie != calibSpecieId) {
 286           //if(hL1GammaPatch->GetEntries() !=0 ) {
 287           if(hL1GammaPatch->GetEntries() > 10) { // need some statistics for hot spot calculation
 288             lstF = hL1GammaPatch->GetListOfFunctions();
 289
 290             // Checker for L1GammaPatch
 291             //Double_t dL1GmeanTrig    = 1./2961.;
 292             //Double_t dL1GmeanTrigTRU = 1./32.;
 293             //Int_t sigmaG    = 100; // deviation from mean value (increased to 100)
 294             //Int_t sigmaGTRU = 5; // deviation from mean value for TRUs
 295             Double_t dL1GEntries = hL1GammaPatch->GetEntries();
 296             Int_t badL1G[48][64]   = {{0}} ;
 297             Int_t badL1GTRU[2][16] = {{0}} ;
 298             Int_t nBadL1G    = 0;
 299             Int_t nBadL1GTRU = 0;
 300             Double_t binContentTRU[2][16] = {{0.}};
 301             for(Int_t ix = 1; ix <=  hL1GammaPatch->GetNbinsX(); ix++) {
 302               for(Int_t iy = 1; iy <=  hL1GammaPatch->GetNbinsY(); iy++) {
 303                 Double_t binContent = hL1GammaPatch->GetBinContent(ix, iy) ;
 304                 if (binContent != 0) {
 305
 306                   // fill counter for TRUs
 307                   binContentTRU[(Int_t)((ix-1)/24)][(Int_t)((iy-1)/4)] += binContent;
 308
 309                   // OLD METHOD (if a patch triggers > sigmaG * mean value (1/#patch positions total) says "hot spot !")
 310                   // if ((double)binContent/(double)dL1GEntries > sigmaG*dL1GmeanTrig) {
 311                   //    badL1G[ix-1][iy-1] += 1;
 312                   //    nBadL1G += 1;
 313                   // }
 314
 315                   // NEW METHOD (if Rate > Threshold * ( (Number of towers or TRUs * Average rate) - Rate ) --> "hot spot !")
 316                   // Thresold = how much does the noisy tower/TRU contribute to the rate
 317                   //            1.0 --> Rate of noisy tower/TRU = Rate of all other towers/TRUs
 318                   if (binContent/dL1GEntries > ThreshG / ( 1 + ThreshG )) {
 319                     badL1G[ix-1][iy-1] += 1;
 320                     nBadL1G += 1;
 321                   }
 322                 }
 323               }
 324             }
 325
 326             // check TRUs
 327             for(Int_t ix = 1; ix <=  2; ix++) {
 328               for(Int_t iy = 1; iy <=  16; iy++) {
 329                 if(binContentTRU[ix-1][iy-1]/dL1GEntries >  ThreshG / ( 1 + ThreshG )) {
 330                   badL1GTRU[ix-1][iy-1] += 1;
 331                   nBadL1GTRU += 1;
 332                 }
 333               }
 334             }
 335
 336             if(fTextL1[0]){
 337               lstF->Add(fTextL1[0]->Clone()) ;
 338               fTextL1[0]->Clear() ;
 339
 340               if (nBadL1G == 0 && nBadL1GTRU == 0 ) {
 341                 fTextL1[0]->SetFillColor(3) ;
 342                 fTextL1[0]->AddText(Form("L1 GAMMA TRIGGER = OK, ENJOY..."));
 343               }
 344               else if (nBadL1G == 0){
 345                 fTextL1[0]->SetFillColor(2) ;
 346                 fTextL1[0]->AddText(Form("HOT SPOT IN L1 GAMMA TRIGGER (TRU) = CALL EXPERT!!"));
 347
 348               }
 349               else{
 350                 fTextL1[0]->SetFillColor(2) ;
 351                 fTextL1[0]->AddText(Form("HOT SPOT IN L1 GAMMA TRIGGER = CALL EXPERT!!"));
 352                 /*
 353                   for(Int_t ix = 1; ix <=  hL1GammaPatch->GetNbinsX(); ix++) {
 354                   for(Int_t iy = 1; iy <=  hL1GammaPatch->GetNbinsY(); iy++) {
 355                   if(badL1G[ix-1][iy-1] != 0) printf("L1 Gamma patch with position x = %d, y = %d is out of range\n",ix,iy);
 356                   }
 357                   }
 358                 */
 359               }
 360             }//fTextL1[0]
 361           }// L1 gamma patch checking done
 362         }// if (specie != calibSpecieId) ..
 363       }//hL1GammaPatch NOT NULL
 364
 365       if( hL1JetPatch ){
 366
 367         lstF = hL1JetPatch->GetListOfFunctions();
 368         CleanListOfFunctions(lstF);
 369
 370         if (specie != calibSpecieId) {
 371
 372           //if(hL1JetPatch->GetEntries() !=0) {
 373           if(hL1JetPatch->GetEntries() > 10) { // need some statistics for hot spot calculation
 374
 375             lstF = hL1JetPatch->GetListOfFunctions();
 376
 377             // Checker for L1JetPatch
 378             //Double_t dL1JmeanTrig = 1/126.;
 379             //Int_t sigmaJ = 5; // deviation from  mean value
 380             Double_t dL1JEntries = hL1JetPatch->GetEntries();
 381             Int_t badL1J[12][16] = {{0}} ;
 382             Int_t nBadL1J = 0;
 383             for(Int_t ix = 1; ix <=  hL1JetPatch->GetNbinsX(); ix++) {
 384               for(Int_t iy = 1; iy <=  hL1JetPatch->GetNbinsY(); iy++) {
 385                 Double_t binContent = hL1JetPatch->GetBinContent(ix, iy) ;
 386                 if (binContent != 0) {
 387
 388                   // OLD METHOD  (if a patch triggers > sigmaJ * mean value (1/#patch positions total) says "hot spot !")
 389                   // if ((double)binContent/(double)dL1JEntries > sigmaJ*dL1JmeanTrig) {
 390                   //    badL1J[ix-1][iy-1] += 1 ;
 391                   //    nBadL1J += 1;
 392                   // }
 393
 394                   // NEW METHOD (if Rate > Threshold * ( (Number of towers or TRUs * Average rate) - Rate ) --> "hot spot !")
 395                   // Threshold: same definitionas for Gamma
 396                   if ((double)binContent/(double)dL1JEntries > ThreshJ / ( 1 + ThreshJ )) {
 397                     badL1J[ix-1][iy-1] += 1 ;
 398                     nBadL1J += 1;
 399                   }
 400                 }
 401               }
 402             }
 403
 404             if(fTextL1[1]){
 405               lstF->Add(fTextL1[1]->Clone()) ;
 406               fTextL1[1]->Clear() ;
 407
 408               if (nBadL1J == 0) {
 409                 fTextL1[1]->SetFillColor(3) ;
 410                 fTextL1[1]->AddText(Form("L1 JET TRIGGER = OK, ENJOY..."));
 411               }
 412               else {
 413                 fTextL1[1]->SetFillColor(2) ;
 414                 fTextL1[1]->AddText(Form("HOT SPOT IN L1 JET TRIGGER = CALL EXPERT!!"));
 415                 /*
 416                   for(Int_t ix = 1; ix <=  hL1JetPatch->GetNbinsX(); ix++) {
 417                   for(Int_t iy = 1; iy <=  hL1JetPatch->GetNbinsY(); iy++) {
 418                   if(badL1J[ix-1][iy-1] != 0) printf("L1 Jet patch with position x = %d, y = %d is out of range\n",(4*ix-4),(4*iy-4));
 419                   }
 420                   }
 421                 */
 422
 423               }
 424             }//fTextL1[1]
 425           } // L1 Jet patch checking done
 426         } // if (specie != calibSpecieId) ..
 427       }// hL1JetPatch NOT NULL
 428
 429       if(hFrameR){
 430
 431         lstF = hFrameR->GetListOfFunctions();
 432         CleanListOfFunctions(lstF);
 433
 434         if(hFrameR->GetEntries() !=0) {
 435           lstF = hFrameR->GetListOfFunctions();
 436
 437           Int_t badLink[32] = {0};
 438           Int_t nBadLink = 0;
 439           for(Int_t ix = 1; ix <= hFrameR->GetNbinsX(); ix++) {
 440             Double_t binContent = hFrameR->GetBinContent(ix) ;
 441             if (binContent == 0) {
 442               badLink[ix-1] += 1;
 443               nBadLink += 1;
 444             }
 445           }
 446           if(fTextL1[2]){
 447             lstF->Add(fTextL1[2]->Clone()) ;
 448             fTextL1[2]->Clear() ;
 449
 450             if (nBadLink < badLinkThresh) {
 451               fTextL1[2]->SetFillColor(3) ;
 452               fTextL1[2]->AddText(Form("LINK TRU-STU = OK, ENJOY..."));
 453             }
 454             else {
 455               fTextL1[2]->SetFillColor(2) ;
 456               fTextL1[2]->AddText(Form("PROBLEM WITH TRU-STU LINK = CALL EXPERT!!"));
 457               /*
 458                 for(Int_t ix = 0; ix <= hFrameR->GetNbinsX(); ix++) {
 459                 if(badLink[ix] != 0) printf("STU link with TRU %d is out\n",ix);
 460                 }
 461               */
 462             }
 463           }//fTextL1[2]
 464         } // Checker for link TRU-STU done
 465       }//hFrameR NOT NULL
 466     } // species processed
 467   } // specie
 468 }
 469 //______________________________________________________________________________
 470
 471 void AliEMCALQAChecker::Init(const AliQAv1::DETECTORINDEX_t det)
 472 {
 473         /// intialises QA and QA checker settings
 474         AliQAv1::Instance(det) ;
 475         Float_t hiValue[AliQAv1::kNBIT] ;
 476         Float_t lowValue[AliQAv1::kNBIT] ;
 477         lowValue[AliQAv1::kINFO]      = 0.0   ;
 478         hiValue[AliQAv1::kINFO]       = 0.1 ;
 479         lowValue[AliQAv1::kWARNING]   = 0.1 ;
 480   hiValue[AliQAv1::kWARNING]    = 0.5 ;
 481         lowValue[AliQAv1::kERROR]     = 0.5   ;
 482         hiValue[AliQAv1::kERROR]      = 0.8 ;
 483         lowValue[AliQAv1::kFATAL]     = 0.8   ;
 484         hiValue[AliQAv1::kFATAL]      = 1.0 ;
 485         SetHiLo(&hiValue[0], &lowValue[0]) ;
 486 }
 487
 488 //______________________________________________________________________________
 489
 490 void AliEMCALQAChecker::CleanListOfFunctions(TList *list)
 491 { // clean up
 492
 493         if (list) {
 494                 TObject *stats = list->FindObject("stats"); list->Remove(stats);
 495                 TObject *obj;
 496                 while ((obj = list->First())) { while(list->Remove(obj)) { } delete obj; }
 497                         if (stats) list->Add(stats);
 498                 }
 499         else {
 500                 AliWarning(Form("Checker : empty list of data functions; returning"));
 501                 return;
 502         }
 503
 504 }
 505
 506