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     test[specie] = 0.0 ;
 203     if ( !AliQAv1::Instance()->IsEventSpecieSet(specie)) continue ;
 204     if (list[specie]->GetEntries() == 0) test[specie] = 0. ; // nothing to check
 205     else {
 206                         //get calib histos
 207                         TH2F * hdata  = (TH2F*)list[specie]->At(k2DRatioAmp) ;
 208                         TH1F * ratio = (TH1F*)list[specie]->At(kRatioDist) ;
 209
 210                         //get L1 histos
 211                         TH2F *hL1GammaPatch = (TH2F*)list[specie]->At(kGL1);
 212                         TH2F *hL1JetPatch = (TH2F*)list[specie]->At(kJL1);
 213                         TH1I *hFrameR = (TH1I*)list[specie]->At(kSTUTRU);
 214
 215                         // ========================================================================================
 216                         //calib histo checker first:
 217
 218                         // first clean lines, text (functions)
 219                         lstF = hdata->GetListOfFunctions();
 220                         CleanListOfFunctions(lstF);
 221                         lstF = ratio->GetListOfFunctions();
 222                         CleanListOfFunctions(lstF);
 223
 224                         if(hdata->GetEntries()!=0 && ratio->GetEntries()!=0) {
 225
 226                           lstF = hdata->GetListOfFunctions();
 227
 228                                 //adding the lines to distinguish different SMs
 229                                 lstF->Add(fLineCol->Clone());
 230                                 for(Int_t iLine = 0; iLine < 5; iLine++) {
 231                                         lstF->Add(fLineRow[iLine]->Clone());
 232                                 }
 233                                 //Now adding the text to for each SM
 234                                 for(Int_t iSM = 0 ; iSM < fgknSM ; iSM++){  //number of SMs loop start
 235                                         lstF->Add(fTextSM[iSM]->Clone());
 236                                 }
 237
 238
 239                                 //now check the ratio histogram
 240                                 lstF = ratio->GetListOfFunctions();
 241
 242                                 Double_t binContent = 0. ;
 243                                 Int_t nGoodTower = 0 ;
 244                                 Double_t rv = 0. ;
 245                                 for(Int_t ix = 1; ix <= hdata->GetNbinsX(); ix++) {
 246                                         for(Int_t iy = 1; iy <= hdata->GetNbinsY(); iy++) {
 247                                                 binContent = hdata->GetBinContent(ix, iy) ;
 248                                                 if (binContent < 1.2 && binContent > 0.8) nGoodTower++ ;
 249                                         }
 250                                 }
 251                                 rv = nGoodTower/nTot ;
 252                                 //printf("%2.2f %% towers out of range [0.8, 1.2]\n", (1-rv)*100);
 253                                 if(fText){
 254                                         lstF->Add(fText->Clone()) ;
 255                                         fText->Clear() ;
 256
 257                                         fText->AddText(Form("%2.2f %% towers out of range [0.8, 1.2]", (1-rv)*100));
 258                                         if (rv < ratioThresh) {
 259                                                 test[specie] = ratioThresh;
 260                                                 // 2 lines text info for quality
 261                                                 fText->SetFillColor(2) ;
 262                                                 fText->AddText(Form("EMCAL = NOK, CALL EXPERTS!!!"));
 263                                         }
 264                                         else {
 265                                                 test[specie] = 1 - ratioThresh;
 266                                                 fText->SetFillColor(3) ;
 267                                                 fText->AddText(Form("EMCAL = OK, ENJOY..."));
 268                                         }
 269                                 }//fText
 270                         } // calib histo checking done
 271
 272                         // ========================================================================================
 273                         // now L1 checks:
 274
 275                         // first clean lines, text (functions)
 276                         lstF = hL1GammaPatch->GetListOfFunctions();
 277                         CleanListOfFunctions(lstF);
 278                         lstF = hL1JetPatch->GetListOfFunctions();
 279                         CleanListOfFunctions(lstF);
 280                         lstF = hFrameR->GetListOfFunctions();
 281                         CleanListOfFunctions(lstF);
 282
 283                         if (specie != calibSpecieId) {
 284                           //if(hL1GammaPatch->GetEntries() !=0 ) {
 285                           if(hL1GammaPatch->GetEntries() > 10) { // need some statistics for hot spot calculation
 286                                 lstF = hL1GammaPatch->GetListOfFunctions();
 287
 288                                 // Checker for L1GammaPatch
 289                                 //Double_t dL1GmeanTrig    = 1./2961.;
 290                                 //Double_t dL1GmeanTrigTRU = 1./32.;
 291                                 //Int_t sigmaG    = 100; // deviation from mean value (increased to 100)
 292                                 //Int_t sigmaGTRU = 5; // deviation from mean value for TRUs
 293                                 Double_t dL1GEntries = hL1GammaPatch->GetEntries();
 294                                 Int_t badL1G[48][64]   = {{0}} ;
 295                                 Int_t badL1GTRU[2][16] = {{0}} ;
 296                                 Int_t nBadL1G    = 0;
 297                                 Int_t nBadL1GTRU = 0;
 298                                 Double_t binContentTRU[2][16] = {{0.}};
 299                                 for(Int_t ix = 1; ix <=  hL1GammaPatch->GetNbinsX(); ix++) {
 300                                   for(Int_t iy = 1; iy <=  hL1GammaPatch->GetNbinsY(); iy++) {
 301                                     Double_t binContent = hL1GammaPatch->GetBinContent(ix, iy) ;
 302                                     if (binContent != 0) {
 303
 304                                       // fill counter for TRUs
 305                                       binContentTRU[(Int_t)((ix-1)/24)][(Int_t)((iy-1)/4)] += binContent;
 306
 307                                       // OLD METHOD (if a patch triggers > sigmaG * mean value (1/#patch positions total) says "hot spot !")
 308                                       // if ((double)binContent/(double)dL1GEntries > sigmaG*dL1GmeanTrig) {
 309                                       //        badL1G[ix-1][iy-1] += 1;
 310                                       //        nBadL1G += 1;
 311                                       // }
 312
 313                                       // NEW METHOD (if Rate > Threshold * ( (Number of towers or TRUs * Average rate) - Rate ) --> "hot spot !")
 314                                       // Thresold = how much does the noisy tower/TRU contribute to the rate
 315                                       //            1.0 --> Rate of noisy tower/TRU = Rate of all other towers/TRUs
 316                                       if (binContent/dL1GEntries > ThreshG / ( 1 + ThreshG )) {
 317                                         badL1G[ix-1][iy-1] += 1;
 318                                         nBadL1G += 1;
 319                                       }
 320                                     }
 321                                   }
 322                                 }
 323
 324                                 // check TRUs
 325                                 for(Int_t ix = 1; ix <=  2; ix++) {
 326                                   for(Int_t iy = 1; iy <=  16; iy++) {
 327                                     if(binContentTRU[ix-1][iy-1]/dL1GEntries >  ThreshG / ( 1 + ThreshG )) {
 328                                       badL1GTRU[ix-1][iy-1] += 1;
 329                                       nBadL1GTRU += 1;
 330                                     }
 331                                   }
 332                                 }
 333
 334                                 if(fTextL1[0]){
 335                                         lstF->Add(fTextL1[0]->Clone()) ;
 336                                         fTextL1[0]->Clear() ;
 337
 338                                         if (nBadL1G == 0 && nBadL1GTRU == 0 ) {
 339                                                 fTextL1[0]->SetFillColor(3) ;
 340                                                 fTextL1[0]->AddText(Form("L1 GAMMA TRIGGER = OK, ENJOY..."));
 341                                         }
 342                                         else if (nBadL1G == 0){
 343                                                 fTextL1[0]->SetFillColor(2) ;
 344                                                 fTextL1[0]->AddText(Form("HOT SPOT IN L1 GAMMA TRIGGER (TRU) = CALL EXPERT!!"));
 345
 346                                         }
 347                                         else{
 348                                           fTextL1[0]->SetFillColor(2) ;
 349                                           fTextL1[0]->AddText(Form("HOT SPOT IN L1 GAMMA TRIGGER = CALL EXPERT!!"));
 350                                           /*
 351                                             for(Int_t ix = 1; ix <=  hL1GammaPatch->GetNbinsX(); ix++) {
 352                                             for(Int_t iy = 1; iy <=  hL1GammaPatch->GetNbinsY(); iy++) {
 353                                                                 if(badL1G[ix-1][iy-1] != 0) printf("L1 Gamma patch with position x = %d, y = %d is out of range\n",ix,iy);
 354                                                                 }
 355                                                                 }
 356                                           */
 357                                         }
 358                                 }//fTextL1[0]
 359                         }// L1 gamma patch checking done
 360
 361                           //if(hL1JetPatch->GetEntries() !=0) {
 362                         if(hL1JetPatch->GetEntries() > 10) { // need some statistics for hot spot calculation
 363
 364                           lstF = hL1JetPatch->GetListOfFunctions();
 365
 366                                 // Checker for L1JetPatch
 367                                 //Double_t dL1JmeanTrig = 1/126.;
 368                                 //Int_t sigmaJ = 5; // deviation from  mean value
 369                                 Double_t dL1JEntries = hL1JetPatch->GetEntries();
 370                                 Int_t badL1J[12][16] = {{0}} ;
 371                                 Int_t nBadL1J = 0;
 372                                 for(Int_t ix = 1; ix <=  hL1JetPatch->GetNbinsX(); ix++) {
 373                                   for(Int_t iy = 1; iy <=  hL1JetPatch->GetNbinsY(); iy++) {
 374                                     Double_t binContent = hL1JetPatch->GetBinContent(ix, iy) ;
 375                                     if (binContent != 0) {
 376
 377                                       // OLD METHOD  (if a patch triggers > sigmaJ * mean value (1/#patch positions total) says "hot spot !")
 378                                       // if ((double)binContent/(double)dL1JEntries > sigmaJ*dL1JmeanTrig) {
 379                                       //        badL1J[ix-1][iy-1] += 1 ;
 380                                       //        nBadL1J += 1;
 381                                       // }
 382
 383                                       // NEW METHOD (if Rate > Threshold * ( (Number of towers or TRUs * Average rate) - Rate ) --> "hot spot !")
 384                                       // Threshold: same definitionas for Gamma
 385                                        if ((double)binContent/(double)dL1JEntries > ThreshJ / ( 1 + ThreshJ )) {
 386                                         badL1J[ix-1][iy-1] += 1 ;
 387                                         nBadL1J += 1;
 388                                        }
 389                                     }
 390                                   }
 391                                 }
 392
 393                                 if(fTextL1[1]){
 394                                         lstF->Add(fTextL1[1]->Clone()) ;
 395                                         fTextL1[1]->Clear() ;
 396
 397                                         if (nBadL1J == 0) {
 398                                                 fTextL1[1]->SetFillColor(3) ;
 399                                                 fTextL1[1]->AddText(Form("L1 JET TRIGGER = OK, ENJOY..."));
 400                                         }
 401                                         else {
 402                                                 fTextL1[1]->SetFillColor(2) ;
 403                                                 fTextL1[1]->AddText(Form("HOT SPOT IN L1 JET TRIGGER = CALL EXPERT!!"));
 404 /*
 405                                                 for(Int_t ix = 1; ix <=  hL1JetPatch->GetNbinsX(); ix++) {
 406                                                         for(Int_t iy = 1; iy <=  hL1JetPatch->GetNbinsY(); iy++) {
 407                                                                 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));
 408                                                         }
 409                                                 }
 410 */
 411
 412                                         }
 413                                 }//fTextL1[1]
 414                         } // L1 Jet patch checking done
 415       } // if (specie != calibSpecieId) ..
 416
 417                         if(hFrameR->GetEntries() !=0) {
 418                           lstF = hFrameR->GetListOfFunctions();
 419
 420                                 Int_t badLink[32] = {0};
 421                                 Int_t nBadLink = 0;
 422                                 for(Int_t ix = 1; ix <= hFrameR->GetNbinsX(); ix++) {
 423                                         Double_t binContent = hFrameR->GetBinContent(ix) ;
 424                                         if (binContent == 0) {
 425                                                 badLink[ix-1] += 1;
 426                                                 nBadLink += 1;
 427                                         }
 428                                 }
 429                                 if(fTextL1[2]){
 430                                         lstF->Add(fTextL1[2]->Clone()) ;
 431                                         fTextL1[2]->Clear() ;
 432
 433                                 if (nBadLink < badLinkThresh) {
 434                                                 fTextL1[2]->SetFillColor(3) ;
 435                                                 fTextL1[2]->AddText(Form("LINK TRU-STU = OK, ENJOY..."));
 436                                         }
 437                                         else {
 438                                                 fTextL1[2]->SetFillColor(2) ;
 439                                                 fTextL1[2]->AddText(Form("PROBLEM WITH TRU-STU LINK = CALL EXPERT!!"));
 440 /*
 441                                                 for(Int_t ix = 0; ix <= hFrameR->GetNbinsX(); ix++) {
 442                                                         if(badLink[ix] != 0) printf("STU link with TRU %d is out\n",ix);
 443                                                 }
 444 */
 445                                         }
 446                                 }//fTextL1[2]
 447                         } // Checker for link TRU-STU done
 448     }
 449
 450   } // specie
 451 }
 452 //______________________________________________________________________________
 453
 454 void AliEMCALQAChecker::Init(const AliQAv1::DETECTORINDEX_t det)
 455 {
 456         /// intialises QA and QA checker settings
 457         AliQAv1::Instance(det) ;
 458         Float_t hiValue[AliQAv1::kNBIT] ;
 459         Float_t lowValue[AliQAv1::kNBIT] ;
 460         lowValue[AliQAv1::kINFO]      = 0.0   ;
 461         hiValue[AliQAv1::kINFO]       = 0.1 ;
 462         lowValue[AliQAv1::kWARNING]   = 0.1 ;
 463   hiValue[AliQAv1::kWARNING]    = 0.5 ;
 464         lowValue[AliQAv1::kERROR]     = 0.5   ;
 465         hiValue[AliQAv1::kERROR]      = 0.8 ;
 466         lowValue[AliQAv1::kFATAL]     = 0.8   ;
 467         hiValue[AliQAv1::kFATAL]      = 1.0 ;
 468         SetHiLo(&hiValue[0], &lowValue[0]) ;
 469 }
 470
 471 //______________________________________________________________________________
 472
 473 void AliEMCALQAChecker::CleanListOfFunctions(TList *list)
 474 { // clean up
 475
 476         if (list) {
 477                 TObject *stats = list->FindObject("stats"); list->Remove(stats);
 478                 TObject *obj;
 479                 while ((obj = list->First())) { while(list->Remove(obj)) { } delete obj; }
 480                         if (stats) list->Add(stats);
 481                 }
 482         else {
 483                 AliWarning(Form("Checker : empty list of data functions; returning"));
 484                 return;
 485         }
 486
 487 }
 488
 489