1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
20 // Checks the quality assurance.
21 // By comparing with reference data
25 // --- ROOT system ---
31 #include <TIterator.h>
35 #include <TParameter.h>
36 #include <TPaveText.h>
37 // --- Standard library ---
39 // --- AliRoot header files ---
42 #include "AliQAChecker.h"
43 #include "AliHMPIDQAChecker.h"
44 #include "AliCDBEntry.h"
45 #include "AliQAManager.h"
46 #include "AliQAThresholds.h"
48 ClassImp(AliHMPIDQAChecker)
49 //_________________________________________________________________
50 AliHMPIDQAChecker::AliHMPIDQAChecker() :
51 AliQACheckerBase("HMPID","HMPID Quality Assurance Data Checker"),
55 fHmpQaThr_NumberOfExcludedDDL(0),
56 fHmpQaThr_DataSizeLowerThreshold(900),
57 fHmpQaThr_DataSizeUpperThreshold(1500),
58 fHmpQaThr_PadOccupancyLowerThreshold(0.005),
59 fHmpQaThr_PadOccupancyUpperThreshold(0.8),
60 fHmpQaThr_SectorGainLossWarningThreshold(3),
61 fHmpQaThr_SectorGainLossErrorThreshold(6),
62 fHmpQaThr_MissingPadFractionWarningThreshold(0.3),
63 fHmpQaThr_MissingPadFractionErrorThreshold(0.5),
69 //ctor, fetches the reference data from OCDB
70 char * detOCDBDir = Form("HMPID/%s/%s", AliQAv1::GetRefOCDBDirName(), AliQAv1::GetRefDataDirName()) ;
71 AliCDBEntry * QARefRec = AliQAManager::QAManager()->Get(detOCDBDir);
73 fQARefRec = dynamic_cast<TObjArray*> (QARefRec->GetObject()) ;
75 if (fQARefRec->GetEntries())
76 fNoReference = kFALSE ;
78 AliInfo("QA reference data NOT retrieved for Reconstruction check. No HMPIDChecker!");
81 InitOnlineThresholds();
85 //_________________________________________________________________
86 AliHMPIDQAChecker::~AliHMPIDQAChecker()
88 if(fQARefRec) { fQARefRec->Delete() ; delete fQARefRec ; }
90 //_________________________________________________________________
91 void AliHMPIDQAChecker::Check(Double_t * check, AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/)
94 // Main check function: Depending on the TASK, different checks are applied
95 // At the moment: check for empty histograms and checks for RecPoints
99 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
101 //printf("+++++++++++++++++++++ specie %d name: %s \n",specie,AliRecoParam::GetEventSpecieName(specie));
102 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) continue ;
103 // checking for empy histograms
104 if(CheckEntries(list[specie]) == 0) {
105 AliWarning("histograms are empty");
106 check[specie] = 0.4;//-> Corresponds to kWARNING see AliQACheckerBase::Run
109 check[specie] = AliQAv1::kINFO ;
113 if(index == AliQAv1::kSIM) check[specie] = CheckSim(list[specie], fQARefRec);
115 // checking rec points
116 if(index == AliQAv1::kREC) check[specie] = CheckRec(list[specie], fQARefRec);
119 if(index == AliQAv1::kRAW) { check[specie] = CheckRaw(specie,list[specie]); }
126 //_________________________________________________________________
127 Double_t AliHMPIDQAChecker::CheckEntries(TObjArray * list) const
130 // check on the QA histograms on the input list:
133 Double_t test = 0.0 ;
136 if (list->GetEntries() == 0){
137 test = 1. ; // nothing to check
143 while ( (hdata = dynamic_cast<TH1 *>(next())) ) {
146 if(hdata->GetEntries()>0)rv=1;
151 AliError("Data type cannot be processed") ;
157 AliWarning("Histograms are booked for THIS specific Task, but they are all empty: setting flag to kWARNING");
158 test = 0.; //upper limit value to set kWARNING flag for a task
166 //_________________________________________________________________
167 Double_t AliHMPIDQAChecker::CheckSim(TObjArray *listsim, TObjArray *listref) const
170 // check on the HMPID RecPoints by using expo fit and Kolmogorov Test:
173 Float_t checkresponse = 0;
175 Float_t counter = 0 ;
176 TIter next(listsim) ;
178 while ( (histo = dynamic_cast<TH1 *>(next())) ) {
179 //PH The histogram should have at least 10 bins with at least 5 entries
180 Int_t nbinsabove = 0;
181 for (Int_t ibin=histo->FindBin(1); ibin<=histo->FindBin(50); ibin++) {
182 if (histo->GetBinContent(ibin)>5) nbinsabove++;
185 if( nbinsabove < 10 ) counter++;
187 TString h = histo->GetTitle();
188 if(h.Contains("Zoom")){
189 histo->Fit("expo","LQ0","",5,50);
190 if(histo->GetFunction("expo")->GetParameter(1) !=0 ) if(TMath::Abs((-1./(histo->GetFunction("expo"))->GetParameter(1)) - 35 ) > 5) counter++;
192 if(h.Contains("size MIP")) if(TMath::Abs(histo->GetMean()-5) > 2) counter++;
193 if(h.Contains("size Phots")) if(TMath::Abs(histo->GetMean()-2) > 2) counter++;
194 if(h.Contains("distribution")) if(histo->KolmogorovTest((TH1F *)listref->At(0))<0.8) counter++;
195 AliDebug(AliQAv1::GetQADebugLevel(),Form(" Kolm. test : %f ",histo->KolmogorovTest((TH1F *)listref->At(0))));
198 Float_t response = counter/(7.+7.+42.+42.); // 7.+7.+42 +42 = N checked histograms (-> To be replaced by listsim->GetEntries())
200 if(response < 0.1) checkresponse = 0.7; // <10% of the check histograms show a failing check -> Corresponds to kINFO see AliQACheckerBase::Run
201 else if(response < 0.5) checkresponse = 0.4; // 50% of the check histograms show a failing check -> Corresponds to kWARNING see AliQACheckerBase::Run
202 else checkresponse = 0.001; // > 50% of the check histograms show a failing check -> Corresponds to kERROR see AliQACheckerBase::Run
203 return checkresponse;
206 //___________________________________________________________________________________________________
207 Double_t AliHMPIDQAChecker::CheckRec(TObjArray *listrec, TObjArray *listref) const
210 // check on the HMPID RecPoints by using expo fit and Kolmogorov Test:
213 Float_t checkresponse = 0;
215 Float_t counter = 0 ;
216 TIter next(listrec) ;
218 while ( (histo = dynamic_cast<TH1 *>(next())) ) {
219 //PH The histogram should have at least 10 bins with at least 5 entries
220 Int_t nbinsabove = 0;
221 for (Int_t ibin=histo->FindBin(1); ibin<=histo->FindBin(50); ibin++) {
222 if (histo->GetBinContent(ibin)>5) nbinsabove++;
225 if( nbinsabove < 10 ) counter++;
227 TString h = histo->GetTitle();
228 if(h.Contains("Zoom")){
229 histo->Fit("expo","LQ0","",5,50);
230 if(histo->GetFunction("expo")->GetParameter(1) !=0 ) if(TMath::Abs((-1./(histo->GetFunction("expo"))->GetParameter(1)) - 35 ) > 5) counter++;
232 if(h.Contains("size MIP")) if(TMath::Abs(histo->GetMean()-5) > 2) counter++;
233 if(h.Contains("size Phots")) if(TMath::Abs(histo->GetMean()-2) > 2) counter++;
234 if(h.Contains("distribution")) if(histo->KolmogorovTest((TH1F *)listref->At(0))<0.8) counter++;
235 AliDebug(AliQAv1::GetQADebugLevel(),Form(" Kolm. test : %f ",histo->KolmogorovTest((TH1F *)listref->At(0))));
238 Float_t response = counter/(7.+7.+42.+42.); // 7.+7.+42 +42 = N checked histograms (-> To be replaced by listrec->GetEntries())
240 if(response < 0.1) checkresponse = 0.7; // <10% of the check histograms show a failing check -> Corresponds to kINFO see AliQACheckerBase::Run
241 else if(response < 0.5) checkresponse = 0.4; // 50% of the check histograms show a failing check -> Corresponds to kWARNING see AliQACheckerBase::Run
242 else checkresponse = 0.001; // > 50% of the check histograms show a failing check -> Corresponds to kERROR see AliQACheckerBase::Run
243 return checkresponse;
245 //___________________________________________________________________________________________________
246 Double_t AliHMPIDQAChecker::CheckRaw(Int_t specie, TObjArray* list)
249 // Check the raw data for offline / online using default or updated thresholds from AMORE
250 // As of now (06/07/2012) the quality flag of all histos in raw will be est to the result of teh CheckRaw and displayed
251 // in AMORE. But we can pu undividual labels.
253 Int_t raqQualFlag = AliQAv1::kNULLBit;
255 Int_t hmpQaFlags[4]={-1}; //init for the 4 shifter histos
257 TString histname ="null";
258 TPaveText text(0.65,0.8,0.9,0.99,"NDC");
259 TPaveText text1(0.65,0.8,0.9,0.99,"NDC");
262 Int_t entries = list->GetEntriesFast();
264 if ( entries == 0 ) {
265 AliWarning(Form("HMPID QA Checker RAWS: no object to analyse! Exiting..."));
266 return AliQAv1::kFATAL;
269 TLine* lineDdlDataSizeMin = new TLine(1536,fHmpQaThr_DataSizeUpperThreshold,1548,fHmpQaThr_DataSizeUpperThreshold);
270 TLine* lineDdlDataSizeMax = new TLine(1536,fHmpQaThr_DataSizeLowerThreshold,1548,fHmpQaThr_DataSizeLowerThreshold);
272 TLine* linePadOccupMin = new TLine(1536,fHmpQaThr_PadOccupancyUpperThreshold,1548,fHmpQaThr_PadOccupancyUpperThreshold);
273 TLine* linePadOccupMax = new TLine(1536,fHmpQaThr_PadOccupancyLowerThreshold,1548,fHmpQaThr_PadOccupancyLowerThreshold);
276 Int_t badDdlCnt = 0, badOccCnt = 0;
278 //___ check data size per ddl
279 if( list->FindObject(Form("%s_hHmpDdlDataSize",AliRecoParam::GetEventSpecieName(specie)) ))
281 TH1* h1 = dynamic_cast<TH1*>( list->FindObject(Form("%s_hHmpDdlDataSize",AliRecoParam::GetEventSpecieName(specie))));
283 if( h1->Integral() > 1 ) {
284 // no entres -> fatal
285 if( h1->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;}
289 // clean up the text, stat, lines, ...
290 if( h1->GetListOfFunctions() ) h1->GetListOfFunctions()->Clear();
292 for ( Int_t iddl = 1 ; iddl <= 14; iddl++)
294 if( h1->GetBinContent(iddl) < fHmpQaThr_DataSizeLowerThreshold || h1->GetBinContent(iddl) > fHmpQaThr_DataSizeUpperThreshold ) badDdlCnt++;
296 //___ check if one or more DDLs are excluded
299 badDdlCnt -= fHmpQaThr_NumberOfExcludedDDL;
303 //___ check how many are bad
304 if ( badDdlCnt == 0 )
306 hmpQaFlags[0] = AliQAv1::kINFO;
308 text1.AddText(Form("OK (%d)",fIsOnlineThr));
309 text1.SetFillColor(kGreen);
310 h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone());
311 lineDdlDataSizeMin->SetLineColor(kGreen);
312 lineDdlDataSizeMax->SetLineColor(kGreen);
313 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone());
314 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone());
316 else if ( badDdlCnt == 1 )
318 hmpQaFlags[0] = AliQAv1::kWARNING;
320 text1.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr));
321 text1.SetFillColor(kOrange);
322 h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone());
323 lineDdlDataSizeMin->SetLineColor(kOrange);
324 lineDdlDataSizeMax->SetLineColor(kOrange);
325 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone());
326 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone());
328 else if ( badDdlCnt >= 2 )
330 hmpQaFlags[0] = AliQAv1::kERROR;
332 text1.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr));
333 text1.SetFillColor(kRed);
334 h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone());
335 lineDdlDataSizeMin->SetLineColor(kRed);
336 lineDdlDataSizeMax->SetLineColor(kRed);
337 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone());
338 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone());
342 hmpQaFlags[0] = AliQAv1::kFATAL;
344 text1.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr));
345 text1.SetFillColor(kRed);
346 h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone());
347 lineDdlDataSizeMin->SetLineColor(kRed);
348 lineDdlDataSizeMax->SetLineColor(kRed);
349 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone());
350 h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone());
353 }//the histo is filled
354 }//___hHmpDdlDataSize
358 if( list->FindObject(Form("%s_fHmpPadOcc",AliRecoParam::GetEventSpecieName(specie)) ))
362 TH1* h1 = dynamic_cast<TH1*>( list->FindObject(Form("%s_fHmpPadOcc",AliRecoParam::GetEventSpecieName(specie))));
364 if( h1->Integral() > 1 ) {
365 // no entres -> fatal
366 if( h1->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;}
368 // clean up the text, stat, lines, ...
369 if( h1->GetListOfFunctions() ) h1->GetListOfFunctions()->Clear();
371 for ( Int_t iddl = 1 ; iddl <= 14; iddl++)
373 if( h1->GetBinContent(iddl) < fHmpQaThr_PadOccupancyLowerThreshold || h1->GetBinContent(iddl) > fHmpQaThr_PadOccupancyUpperThreshold ) badOccCnt++;
376 badOccCnt -= fHmpQaThr_NumberOfExcludedDDL;
378 //___ check how many are bad
379 if ( badOccCnt == 0 )
381 hmpQaFlags[1] = AliQAv1::kINFO;
383 text.AddText(Form("OK (%d)",fIsOnlineThr));
384 text.SetFillColor(kGreen);
385 h1->GetListOfFunctions()->Add((TPaveText*)text.Clone());
386 linePadOccupMin->SetLineColor(kGreen);
387 linePadOccupMax->SetLineColor(kGreen);
388 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone());
389 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone());
391 else if ( badOccCnt == 1 )
393 hmpQaFlags[1] = AliQAv1::kWARNING;
395 text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr));
396 text.SetFillColor(kOrange);
397 h1->GetListOfFunctions()->Add((TPaveText*)text.Clone());
398 linePadOccupMin->SetLineColor(kGreen);
399 linePadOccupMax->SetLineColor(kGreen);
400 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone());
401 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone());
403 else if ( badOccCnt == 2 )
405 hmpQaFlags[1] = AliQAv1::kERROR;
407 text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr));
408 text.SetFillColor(kRed);
409 h1->GetListOfFunctions()->Add((TPaveText*)text.Clone());
410 linePadOccupMin->SetLineColor(kGreen);
411 linePadOccupMax->SetLineColor(kGreen);
412 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone());
413 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone());
417 hmpQaFlags[1] = AliQAv1::kFATAL;
419 text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr));
420 text.SetFillColor(kRed);
421 h1->GetListOfFunctions()->Add((TPaveText*)text.Clone());
422 linePadOccupMin->SetLineColor(kGreen);
423 linePadOccupMax->SetLineColor(kGreen);
424 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone());
425 h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone());
431 Int_t sumPadMapChROR[7]={0};
432 Int_t sumPadMapChROL[7]={0};
433 Int_t bigMapFlag = AliQAv1::kINFO;
434 Int_t errCntBigMap = 0;
436 if( list->FindObject(Form("%s_hHmpBigMap",AliRecoParam::GetEventSpecieName(specie)) ))
438 TH2* h2 = dynamic_cast<TH2*>( list->FindObject(Form("%s_hHmpBigMap",AliRecoParam::GetEventSpecieName(specie))));
440 if( h2->Integral() > 1 ) {
441 // no entres -> fatal
442 if( h2->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;}
444 // clean up the text, stat, lines, ...
445 if( h2->GetListOfFunctions() ) h2->GetListOfFunctions()->Clear();
447 //calculate missing pad fraction
448 for(Int_t ich = 0; ich < 7; ich++)
450 for(Int_t iy=1+ich*144;iy<=144+ich*144;iy++) {
451 for(Int_t ix=1;ix<=80;ix++) if(h2->GetBinContent(ix,iy) > 0) sumPadMapChROL[ich]++;
452 for(Int_t ix=81;ix<=160;ix++) if(h2->GetBinContent(ix,iy) > 0) sumPadMapChROR[ich]++;
455 //check the calculated missing pad fraction
456 for(Int_t ich = 0; ich < 7; ich++)
459 bigMapFlag = AliQAv1::kINFO;
460 if( (1-sumPadMapChROL[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionWarningThreshold ||
461 (1-sumPadMapChROR[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionWarningThreshold )
463 bigMapFlag = AliQAv1::kWARNING;
465 if( (1-sumPadMapChROL[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionErrorThreshold ||
466 (1-sumPadMapChROR[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionErrorThreshold ) {
467 bigMapFlag = AliQAv1::kERROR;
472 if( errCntBigMap > 0 ) bigMapFlag = AliQAv1::kERROR;
475 if ( bigMapFlag == AliQAv1::kINFO )
477 hmpQaFlags[2] = AliQAv1::kINFO;
479 text.AddText(Form("OK (%d)",fIsOnlineThr));
480 text.SetFillColor(kGreen);
481 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
483 else if ( bigMapFlag == AliQAv1::kWARNING )
485 hmpQaFlags[2] = AliQAv1::kWARNING;
487 text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr));
488 text.SetFillColor(kOrange);
489 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
491 else if ( bigMapFlag == AliQAv1::kERROR )
493 hmpQaFlags[2] = AliQAv1::kERROR;
495 text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr));
496 text.SetFillColor(kRed);
497 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
501 hmpQaFlags[2] = AliQAv1::kFATAL;
503 text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr));
504 text.SetFillColor(kRed);
505 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
512 Int_t numSectorsMissing = 0, numSectorsGainLoss = 0;
513 Double_t hSumQPerSector[42]={0};
515 if( list->FindObject(Form("%s_fHmpHvSectorQ",AliRecoParam::GetEventSpecieName(specie))))
518 TH2* h2 = dynamic_cast<TH2*>( list->FindObject(Form("%s_fHmpHvSectorQ",AliRecoParam::GetEventSpecieName(specie))));
520 if(h2->Integral() > 0 ) {
521 // no entres -> fatal
522 if( h2->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;}
524 // clean up the text, stat, lines, ...
525 if( h2->GetListOfFunctions() ) h2->GetListOfFunctions()->Clear();
528 for(Int_t isec = 1 ; isec <= 42; isec++)
530 for(Int_t ibiny=100;ibiny<410;ibiny++) {hSumQPerSector[isec-1] += h2->GetBinContent(ibiny,isec); }
531 if(hSumQPerSector[isec-1] < 0.001) {numSectorsGainLoss++;} // there is no photon and mip peak, gain loss
532 if(h2->GetBinContent(1,isec) < 0.01 ) {numSectorsMissing++; } //practically there is no charge , the sector is missing
534 Int_t sectorErrors = numSectorsGainLoss+numSectorsMissing;
537 if ( sectorErrors <= 3)
539 hmpQaFlags[3] = AliQAv1::kINFO;
541 text.AddText(Form("OK (%d)",fIsOnlineThr));
542 text.SetFillColor(kGreen);
543 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
545 else if ( sectorErrors > fHmpQaThr_SectorGainLossWarningThreshold)
547 hmpQaFlags[3] = AliQAv1::kWARNING;
549 text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr));
550 if(numSectorsMissing > 0 ) text.AddText(Form("MISSING SECTOR?"));
551 text.SetFillColor(kOrange);
552 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
554 else if ( sectorErrors > fHmpQaThr_SectorGainLossErrorThreshold)
556 hmpQaFlags[3] = AliQAv1::kERROR;
558 text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr));
559 if(numSectorsMissing > 0 ) text.AddText(Form("MISSING SECTOR?"));
560 text.SetFillColor(kRed);
561 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
565 hmpQaFlags[3] = AliQAv1::kFATAL;
567 text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr));
568 text.SetFillColor(kRed);
569 h2->GetListOfFunctions()->Add((TPaveText*)text.Clone());
577 lineDdlDataSizeMin->Delete();
578 lineDdlDataSizeMax->Delete();
579 linePadOccupMin->Delete();
580 linePadOccupMax->Delete();
583 switch ( TMath::MaxElement(4,hmpQaFlags))
589 case AliQAv1::kWARNING:
593 case AliQAv1::kERROR:
597 case AliQAv1::kFATAL:
602 dflag = AliQAv1::kNULLBit;
612 //___________________________________________________________________________________________________
613 void AliHMPIDQAChecker::InitOnlineThresholds()
616 // Init the online thresholds from GRP generated by AMORE
619 AliCDBManager* man = AliCDBManager::Instance();
620 if(!man) { fIsOnlineThr = kFALSE; return; }
622 // man->SetRun(AliQAChecker::Instance()->GetRunNumber()); //this should be the normal code
623 // man->SetSpecificStorage("GRP/Calib/QAThresholds","local://$ALICE_ROOT/OCDB/");// need for development machine testing
625 if(!man->Get("GRP/Calib/QAThresholds")) { fIsOnlineThr = kFALSE; return; }
627 AliCDBEntry* entry = man->Get("GRP/Calib/QAThresholds");
628 if(!entry) { fIsOnlineThr = kFALSE; return; }
630 TObjArray* branch = (TObjArray*) entry->GetObject();
631 if(!branch ) { fIsOnlineThr = kFALSE; return; }
633 AliQAThresholds* thresholds = (AliQAThresholds*) branch->FindObject("HMP");
634 if(!thresholds) { fIsOnlineThr = kFALSE; return; }
636 fIsOnlineThr = kTRUE;
639 TString parName = "zero";
640 while ( thresholds->GetThreshold(teCnt))
643 parName = (thresholds->GetThreshold(teCnt))->GetName();
646 if( parName.Contains("HmpNumberOfExcludedDDLthreshold") )
648 TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt);
649 fHmpQaThr_NumberOfExcludedDDL = myParam->GetVal();
654 if( parName.Contains("HmpDataSizeLowerThreshold") )
656 TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt);
657 fHmpQaThr_DataSizeLowerThreshold = myParam->GetVal();
661 if( parName.Contains("HmpDataSizeUpperThreshold") )
663 TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt);
664 fHmpQaThr_DataSizeUpperThreshold = myParam->GetVal();
669 if( parName.Contains("HmpPadOccupancyLowerThreshold") )
671 TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt);
672 fHmpQaThr_PadOccupancyLowerThreshold = myParam->GetVal();
676 if( parName.Contains("HmpPadOccupancyUpperThreshold") )
678 TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt);
679 fHmpQaThr_PadOccupancyUpperThreshold = myParam->GetVal();
684 if( parName.Contains("HmpSectorGainLossWarningThreshold") )
686 TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt);
687 fHmpQaThr_SectorGainLossWarningThreshold = myParam->GetVal();
691 if( parName.Contains("HmpSectorGainLossErrorThreshold") )
693 TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt);
694 fHmpQaThr_SectorGainLossErrorThreshold = myParam->GetVal();
698 if( parName.Contains("HmpMissingPadFractionWarningThreshold") )
700 TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt);
701 fHmpQaThr_MissingPadFractionWarningThreshold = myParam->GetVal();
705 if( parName.Contains("HmpMissingPadFractionErrorThreshold") )
707 TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt);
708 fHmpQaThr_MissingPadFractionErrorThreshold = myParam->GetVal();
716 // PrintThresholds();
718 //___________________________________________________________________________________________________
719 void AliHMPIDQAChecker::PrintThresholds()
721 Printf("--- Printing thresholds ---");
722 Printf(Form("--- Default or online: %d ---",fIsOnlineThr));
723 Printf(Form("--- fHmpQaThr_NumberOfExcludedDDL %d ---",fHmpQaThr_NumberOfExcludedDDL));
724 Printf(Form("--- fHmpQaThr_DataSizeLowerThreshold %d ---",fHmpQaThr_DataSizeLowerThreshold));
725 Printf(Form("--- fHmpQaThr_DataSizeUpperThreshold %d ---",fHmpQaThr_DataSizeUpperThreshold));
726 Printf(Form("--- fHmpQaThr_PadOccupancyLowerThreshold %f ---",fHmpQaThr_PadOccupancyLowerThreshold));
727 Printf(Form("--- fHmpQaThr_PadOccupancyUpperThreshold %f ---",fHmpQaThr_PadOccupancyUpperThreshold));
728 Printf(Form("--- fHmpQaThr_SectorGainLossWarningThreshold %d ---",fHmpQaThr_SectorGainLossWarningThreshold));
729 Printf(Form("--- fHmpQaThr_SectorGainLossErrorThreshold %d ---",fHmpQaThr_SectorGainLossErrorThreshold));
730 Printf(Form("--- fHmpQaThr_MissingPadFractionWarningThreshold %f ---",fHmpQaThr_MissingPadFractionWarningThreshold));
731 Printf(Form("--- fHmpQaThr_MissingPadFractionErrorThreshold %f ---",fHmpQaThr_MissingPadFractionErrorThreshold));
732 Printf("--- Printing thresholds done ---");
738 //___________________________________________________________________________________________________