]>
Commit | Line | Data |
---|---|---|
e4e6ae61 | 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 | /* $Id$ */ | |
18 | ||
19 | //... | |
20 | // Checks the quality assurance. | |
21 | // By comparing with reference data | |
22 | // Skeleton for HMPID | |
23 | //... | |
24 | ||
25 | // --- ROOT system --- | |
26 | #include <TClass.h> | |
27 | #include <TH1F.h> | |
74db01e8 | 28 | #include <TH1I.h> |
29 | #include <TH2.h> | |
fdd7f404 | 30 | #include <TF1.h> |
e4e6ae61 | 31 | #include <TIterator.h> |
32 | #include <TKey.h> | |
74db01e8 | 33 | #include <TFile.h> |
34 | #include <TLine.h> | |
35 | #include <TParameter.h> | |
36 | #include <TPaveText.h> | |
e4e6ae61 | 37 | // --- Standard library --- |
38 | ||
39 | // --- AliRoot header files --- | |
40 | #include "AliLog.h" | |
4e25ac79 | 41 | #include "AliQAv1.h" |
b990acb3 | 42 | #include "AliQAChecker.h" |
43 | #include "AliHMPIDQAChecker.h" | |
fdd7f404 | 44 | #include "AliCDBEntry.h" |
fc7e0df2 | 45 | #include "AliQAManager.h" |
74db01e8 | 46 | #include "AliQAThresholds.h" |
e4e6ae61 | 47 | |
b990acb3 | 48 | ClassImp(AliHMPIDQAChecker) |
34fbdb5b | 49 | //_________________________________________________________________ |
50 | AliHMPIDQAChecker::AliHMPIDQAChecker() : | |
51 | AliQACheckerBase("HMPID","HMPID Quality Assurance Data Checker"), | |
09b0a300 | 52 | fNoReference(kTRUE), |
74db01e8 | 53 | fQARefRec(NULL), |
54 | ||
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), | |
64 | fIsOnlineThr(kFALSE) | |
65 | ||
66 | ||
67 | ||
34fbdb5b | 68 | { |
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); | |
72 | if(QARefRec) { | |
73 | fQARefRec = dynamic_cast<TObjArray*> (QARefRec->GetObject()) ; | |
74 | if (fQARefRec) | |
75 | if (fQARefRec->GetEntries()) | |
76 | fNoReference = kFALSE ; | |
77 | if (fNoReference) | |
78 | AliInfo("QA reference data NOT retrieved for Reconstruction check. No HMPIDChecker!"); | |
79 | } | |
303d90aa | 80 | |
34fbdb5b | 81 | } |
82 | ||
34fbdb5b | 83 | //_________________________________________________________________ |
84 | AliHMPIDQAChecker::~AliHMPIDQAChecker() | |
85 | { | |
71edfefd | 86 | if(fQARefRec) { fQARefRec->Delete() ; delete fQARefRec ; } |
34fbdb5b | 87 | } |
fdd7f404 | 88 | //_________________________________________________________________ |
86017bd8 | 89 | void AliHMPIDQAChecker::Check(Double_t * check, AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/) |
fdd7f404 | 90 | { |
91 | // | |
92 | // Main check function: Depending on the TASK, different checks are applied | |
93 | // At the moment: check for empty histograms and checks for RecPoints | |
94 | ||
303d90aa | 95 | InitOnlineThresholds(); |
96 | ||
34fbdb5b | 97 | if(fNoReference) |
fdd7f404 | 98 | |
57acd2d2 | 99 | for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) { |
74db01e8 | 100 | check[specie] = 1.0; |
101 | //printf("+++++++++++++++++++++ specie %d name: %s \n",specie,AliRecoParam::GetEventSpecieName(specie)); | |
102 | if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) continue ; | |
2ac899f2 | 103 | // checking for empy histograms |
57acd2d2 | 104 | if(CheckEntries(list[specie]) == 0) { |
105 | AliWarning("histograms are empty"); | |
106 | check[specie] = 0.4;//-> Corresponds to kWARNING see AliQACheckerBase::Run | |
107 | } | |
108 | ||
74db01e8 | 109 | check[specie] = AliQAv1::kINFO ; |
110 | ||
111 | ||
2ac899f2 | 112 | //check sim |
34fbdb5b | 113 | if(index == AliQAv1::kSIM) check[specie] = CheckSim(list[specie], fQARefRec); |
2ac899f2 | 114 | |
57acd2d2 | 115 | // checking rec points |
34fbdb5b | 116 | if(index == AliQAv1::kREC) check[specie] = CheckRec(list[specie], fQARefRec); |
117 | ||
74db01e8 | 118 | //checking raw data |
119 | if(index == AliQAv1::kRAW) { check[specie] = CheckRaw(specie,list[specie]); } | |
120 | ||
121 | ||
2ac899f2 | 122 | } // species loop |
fdd7f404 | 123 | } |
124 | //_________________________________________________________________ | |
125 | Double_t AliHMPIDQAChecker::CheckEntries(TObjArray * list) const | |
126 | { | |
127 | // | |
128 | // check on the QA histograms on the input list: | |
129 | // | |
130 | ||
131 | Double_t test = 0.0 ; | |
132 | Int_t count = 0 ; | |
133 | ||
134 | if (list->GetEntries() == 0){ | |
135 | test = 1. ; // nothing to check | |
136 | } | |
137 | else { | |
138 | TIter next(list) ; | |
139 | TH1 * hdata ; | |
140 | count = 0 ; | |
141 | while ( (hdata = dynamic_cast<TH1 *>(next())) ) { | |
142 | if (hdata) { | |
143 | Double_t rv = 0.; | |
fdd7f404 | 144 | if(hdata->GetEntries()>0)rv=1; |
145 | count++ ; | |
146 | test += rv ; | |
147 | } | |
148 | else{ | |
149 | AliError("Data type cannot be processed") ; | |
150 | } | |
151 | ||
152 | } | |
153 | if (count != 0) { | |
154 | if (test==0) { | |
155 | AliWarning("Histograms are booked for THIS specific Task, but they are all empty: setting flag to kWARNING"); | |
156 | test = 0.; //upper limit value to set kWARNING flag for a task | |
157 | } | |
158 | else test = 1 ; | |
159 | } | |
160 | } | |
161 | ||
162 | return test ; | |
163 | } | |
164 | //_________________________________________________________________ | |
2ac899f2 | 165 | Double_t AliHMPIDQAChecker::CheckSim(TObjArray *listsim, TObjArray *listref) const |
166 | { | |
167 | // | |
168 | // check on the HMPID RecPoints by using expo fit and Kolmogorov Test: | |
169 | // | |
170 | ||
171 | Float_t checkresponse = 0; | |
172 | ||
173 | Float_t counter = 0 ; | |
174 | TIter next(listsim) ; | |
175 | TH1* histo; | |
176 | while ( (histo = dynamic_cast<TH1 *>(next())) ) { | |
5089b327 | 177 | //PH The histogram should have at least 10 bins with at least 5 entries |
2ac899f2 | 178 | Int_t nbinsabove = 0; |
5089b327 | 179 | for (Int_t ibin=histo->FindBin(1); ibin<=histo->FindBin(50); ibin++) { |
180 | if (histo->GetBinContent(ibin)>5) nbinsabove++; | |
2ac899f2 | 181 | } |
182 | ||
5089b327 | 183 | if( nbinsabove < 10 ) counter++; |
2ac899f2 | 184 | else { |
185 | TString h = histo->GetTitle(); | |
186 | if(h.Contains("Zoom")){ | |
5089b327 | 187 | histo->Fit("expo","LQ0","",5,50); |
2ac899f2 | 188 | if(histo->GetFunction("expo")->GetParameter(1) !=0 ) if(TMath::Abs((-1./(histo->GetFunction("expo"))->GetParameter(1)) - 35 ) > 5) counter++; |
189 | } | |
190 | if(h.Contains("size MIP")) if(TMath::Abs(histo->GetMean()-5) > 2) counter++; | |
191 | if(h.Contains("size Phots")) if(TMath::Abs(histo->GetMean()-2) > 2) counter++; | |
192 | if(h.Contains("distribution")) if(histo->KolmogorovTest((TH1F *)listref->At(0))<0.8) counter++; | |
193 | AliDebug(AliQAv1::GetQADebugLevel(),Form(" Kolm. test : %f ",histo->KolmogorovTest((TH1F *)listref->At(0)))); | |
194 | } | |
195 | } | |
196 | Float_t response = counter/(7.+7.+42.+42.); // 7.+7.+42 +42 = N checked histograms (-> To be replaced by listsim->GetEntries()) | |
197 | ||
198 | if(response < 0.1) checkresponse = 0.7; // <10% of the check histograms show a failing check -> Corresponds to kINFO see AliQACheckerBase::Run | |
199 | else if(response < 0.5) checkresponse = 0.4; // 50% of the check histograms show a failing check -> Corresponds to kWARNING see AliQACheckerBase::Run | |
200 | else checkresponse = 0.001; // > 50% of the check histograms show a failing check -> Corresponds to kERROR see AliQACheckerBase::Run | |
201 | return checkresponse; | |
202 | } | |
fdd7f404 | 203 | |
2ac899f2 | 204 | //___________________________________________________________________________________________________ |
205 | Double_t AliHMPIDQAChecker::CheckRec(TObjArray *listrec, TObjArray *listref) const | |
fdd7f404 | 206 | { |
0f46609e | 207 | // |
fbe0560f | 208 | // check on the HMPID RecPoints by using expo fit and Kolmogorov Test: |
209 | // | |
fdd7f404 | 210 | |
211 | Float_t checkresponse = 0; | |
212 | ||
213 | Float_t counter = 0 ; | |
214 | TIter next(listrec) ; | |
215 | TH1* histo; | |
216 | while ( (histo = dynamic_cast<TH1 *>(next())) ) { | |
5089b327 | 217 | //PH The histogram should have at least 10 bins with at least 5 entries |
3762c8fb | 218 | Int_t nbinsabove = 0; |
5089b327 | 219 | for (Int_t ibin=histo->FindBin(1); ibin<=histo->FindBin(50); ibin++) { |
220 | if (histo->GetBinContent(ibin)>5) nbinsabove++; | |
3762c8fb | 221 | } |
222 | ||
5089b327 | 223 | if( nbinsabove < 10 ) counter++; |
fdd7f404 | 224 | else { |
225 | TString h = histo->GetTitle(); | |
226 | if(h.Contains("Zoom")){ | |
5089b327 | 227 | histo->Fit("expo","LQ0","",5,50); |
fdd7f404 | 228 | if(histo->GetFunction("expo")->GetParameter(1) !=0 ) if(TMath::Abs((-1./(histo->GetFunction("expo"))->GetParameter(1)) - 35 ) > 5) counter++; |
229 | } | |
230 | if(h.Contains("size MIP")) if(TMath::Abs(histo->GetMean()-5) > 2) counter++; | |
231 | if(h.Contains("size Phots")) if(TMath::Abs(histo->GetMean()-2) > 2) counter++; | |
232 | if(h.Contains("distribution")) if(histo->KolmogorovTest((TH1F *)listref->At(0))<0.8) counter++; | |
5379c4a3 | 233 | AliDebug(AliQAv1::GetQADebugLevel(),Form(" Kolm. test : %f ",histo->KolmogorovTest((TH1F *)listref->At(0)))); |
fdd7f404 | 234 | } |
235 | } | |
236 | Float_t response = counter/(7.+7.+42.+42.); // 7.+7.+42 +42 = N checked histograms (-> To be replaced by listrec->GetEntries()) | |
237 | ||
238 | if(response < 0.1) checkresponse = 0.7; // <10% of the check histograms show a failing check -> Corresponds to kINFO see AliQACheckerBase::Run | |
239 | else if(response < 0.5) checkresponse = 0.4; // 50% of the check histograms show a failing check -> Corresponds to kWARNING see AliQACheckerBase::Run | |
240 | else checkresponse = 0.001; // > 50% of the check histograms show a failing check -> Corresponds to kERROR see AliQACheckerBase::Run | |
241 | return checkresponse; | |
242 | } | |
74db01e8 | 243 | //___________________________________________________________________________________________________ |
244 | Double_t AliHMPIDQAChecker::CheckRaw(Int_t specie, TObjArray* list) | |
245 | { | |
246 | // | |
247 | // Check the raw data for offline / online using default or updated thresholds from AMORE | |
248 | // 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 | |
249 | // in AMORE. But we can pu undividual labels. | |
250 | // | |
303d90aa | 251 | |
23874fc6 | 252 | //Int_t raqQualFlag = AliQAv1::kNULLBit; |
74db01e8 | 253 | |
254 | Int_t hmpQaFlags[4]={-1}; //init for the 4 shifter histos | |
255 | ||
256 | TString histname ="null"; | |
257 | TPaveText text(0.65,0.8,0.9,0.99,"NDC"); | |
258 | TPaveText text1(0.65,0.8,0.9,0.99,"NDC"); | |
259 | ||
260 | ||
261 | Int_t entries = list->GetEntriesFast(); | |
262 | ||
263 | if ( entries == 0 ) { | |
264 | AliWarning(Form("HMPID QA Checker RAWS: no object to analyse! Exiting...")); | |
265 | return AliQAv1::kFATAL; | |
266 | } | |
267 | ||
268 | TLine* lineDdlDataSizeMin = new TLine(1536,fHmpQaThr_DataSizeUpperThreshold,1548,fHmpQaThr_DataSizeUpperThreshold); | |
269 | TLine* lineDdlDataSizeMax = new TLine(1536,fHmpQaThr_DataSizeLowerThreshold,1548,fHmpQaThr_DataSizeLowerThreshold); | |
270 | ||
271 | TLine* linePadOccupMin = new TLine(1536,fHmpQaThr_PadOccupancyUpperThreshold,1548,fHmpQaThr_PadOccupancyUpperThreshold); | |
272 | TLine* linePadOccupMax = new TLine(1536,fHmpQaThr_PadOccupancyLowerThreshold,1548,fHmpQaThr_PadOccupancyLowerThreshold); | |
273 | ||
274 | ||
275 | Int_t badDdlCnt = 0, badOccCnt = 0; | |
276 | ||
277 | //___ check data size per ddl | |
278 | if( list->FindObject(Form("%s_hHmpDdlDataSize",AliRecoParam::GetEventSpecieName(specie)) )) | |
279 | { | |
280 | TH1* h1 = dynamic_cast<TH1*>( list->FindObject(Form("%s_hHmpDdlDataSize",AliRecoParam::GetEventSpecieName(specie)))); | |
281 | if(h1) { | |
282 | if( h1->Integral() > 1 ) { | |
283 | // no entres -> fatal | |
23874fc6 | 284 | //if( h1->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;} |
74db01e8 | 285 | h1->SetStats(0); |
286 | ||
287 | ||
288 | // clean up the text, stat, lines, ... | |
289 | if( h1->GetListOfFunctions() ) h1->GetListOfFunctions()->Clear(); | |
290 | ||
291 | for ( Int_t iddl = 1 ; iddl <= 14; iddl++) | |
292 | { | |
293 | if( h1->GetBinContent(iddl) < fHmpQaThr_DataSizeLowerThreshold || h1->GetBinContent(iddl) > fHmpQaThr_DataSizeUpperThreshold ) badDdlCnt++; | |
294 | } | |
295 | //___ check if one or more DDLs are excluded | |
296 | ||
297 | ||
298 | badDdlCnt -= fHmpQaThr_NumberOfExcludedDDL; | |
299 | ||
300 | ||
301 | ||
302 | //___ check how many are bad | |
303 | if ( badDdlCnt == 0 ) | |
304 | { | |
305 | hmpQaFlags[0] = AliQAv1::kINFO; | |
306 | text1.Clear(); | |
307 | text1.AddText(Form("OK (%d)",fIsOnlineThr)); | |
308 | text1.SetFillColor(kGreen); | |
309 | h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone()); | |
310 | lineDdlDataSizeMin->SetLineColor(kGreen); | |
311 | lineDdlDataSizeMax->SetLineColor(kGreen); | |
312 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone()); | |
313 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone()); | |
314 | } | |
315 | else if ( badDdlCnt == 1 ) | |
316 | { | |
317 | hmpQaFlags[0] = AliQAv1::kWARNING; | |
318 | text1.Clear(); | |
319 | text1.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr)); | |
320 | text1.SetFillColor(kOrange); | |
321 | h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone()); | |
322 | lineDdlDataSizeMin->SetLineColor(kOrange); | |
323 | lineDdlDataSizeMax->SetLineColor(kOrange); | |
324 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone()); | |
325 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone()); | |
326 | } | |
327 | else if ( badDdlCnt >= 2 ) | |
328 | { | |
329 | hmpQaFlags[0] = AliQAv1::kERROR; | |
330 | text1.Clear(); | |
331 | text1.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr)); | |
332 | text1.SetFillColor(kRed); | |
333 | h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone()); | |
334 | lineDdlDataSizeMin->SetLineColor(kRed); | |
335 | lineDdlDataSizeMax->SetLineColor(kRed); | |
336 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone()); | |
337 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone()); | |
338 | } | |
339 | else | |
340 | { | |
341 | hmpQaFlags[0] = AliQAv1::kFATAL; | |
342 | text1.Clear(); | |
343 | text1.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr)); | |
344 | text1.SetFillColor(kRed); | |
345 | h1->GetListOfFunctions()->Add((TPaveText*)text1.Clone()); | |
346 | lineDdlDataSizeMin->SetLineColor(kRed); | |
347 | lineDdlDataSizeMax->SetLineColor(kRed); | |
348 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMin->Clone()); | |
349 | h1->GetListOfFunctions()->Add((TLine*)lineDdlDataSizeMax->Clone()); | |
350 | } | |
351 | } | |
352 | }//the histo is filled | |
353 | }//___hHmpDdlDataSize | |
354 | ||
355 | ||
356 | ||
357 | if( list->FindObject(Form("%s_fHmpPadOcc",AliRecoParam::GetEventSpecieName(specie)) )) | |
358 | { | |
359 | ||
360 | ||
361 | TH1* h1 = dynamic_cast<TH1*>( list->FindObject(Form("%s_fHmpPadOcc",AliRecoParam::GetEventSpecieName(specie)))); | |
362 | if(h1) { | |
363 | if( h1->Integral() > 1 ) { | |
364 | // no entres -> fatal | |
23874fc6 | 365 | //if( h1->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;} |
74db01e8 | 366 | h1->SetStats(0); |
367 | // clean up the text, stat, lines, ... | |
368 | if( h1->GetListOfFunctions() ) h1->GetListOfFunctions()->Clear(); | |
369 | ||
370 | for ( Int_t iddl = 1 ; iddl <= 14; iddl++) | |
371 | { | |
372 | if( h1->GetBinContent(iddl) < fHmpQaThr_PadOccupancyLowerThreshold || h1->GetBinContent(iddl) > fHmpQaThr_PadOccupancyUpperThreshold ) badOccCnt++; | |
373 | } | |
374 | ||
375 | badOccCnt -= fHmpQaThr_NumberOfExcludedDDL; | |
376 | ||
377 | //___ check how many are bad | |
378 | if ( badOccCnt == 0 ) | |
379 | { | |
380 | hmpQaFlags[1] = AliQAv1::kINFO; | |
381 | text.Clear(); | |
382 | text.AddText(Form("OK (%d)",fIsOnlineThr)); | |
383 | text.SetFillColor(kGreen); | |
384 | h1->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
385 | linePadOccupMin->SetLineColor(kGreen); | |
386 | linePadOccupMax->SetLineColor(kGreen); | |
387 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone()); | |
388 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone()); | |
389 | } | |
390 | else if ( badOccCnt == 1 ) | |
391 | { | |
392 | hmpQaFlags[1] = AliQAv1::kWARNING; | |
393 | text.Clear(); | |
394 | text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr)); | |
395 | text.SetFillColor(kOrange); | |
396 | h1->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
397 | linePadOccupMin->SetLineColor(kGreen); | |
398 | linePadOccupMax->SetLineColor(kGreen); | |
399 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone()); | |
400 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone()); | |
401 | } | |
402 | else if ( badOccCnt == 2 ) | |
403 | { | |
404 | hmpQaFlags[1] = AliQAv1::kERROR; | |
405 | text.Clear(); | |
406 | text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr)); | |
407 | text.SetFillColor(kRed); | |
408 | h1->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
409 | linePadOccupMin->SetLineColor(kGreen); | |
410 | linePadOccupMax->SetLineColor(kGreen); | |
411 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone()); | |
412 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone()); | |
413 | } | |
414 | else | |
415 | { | |
416 | hmpQaFlags[1] = AliQAv1::kFATAL; | |
417 | text.Clear(); | |
418 | text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr)); | |
419 | text.SetFillColor(kRed); | |
420 | h1->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
421 | linePadOccupMin->SetLineColor(kGreen); | |
422 | linePadOccupMax->SetLineColor(kGreen); | |
423 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMin->Clone()); | |
424 | h1->GetListOfFunctions()->Add((TLine*)linePadOccupMax->Clone()); | |
425 | } | |
426 | } | |
427 | } | |
428 | }//___HmpPadOcc | |
429 | ||
430 | Int_t sumPadMapChROR[7]={0}; | |
431 | Int_t sumPadMapChROL[7]={0}; | |
432 | Int_t bigMapFlag = AliQAv1::kINFO; | |
433 | Int_t errCntBigMap = 0; | |
434 | ||
435 | if( list->FindObject(Form("%s_hHmpBigMap",AliRecoParam::GetEventSpecieName(specie)) )) | |
436 | { | |
437 | TH2* h2 = dynamic_cast<TH2*>( list->FindObject(Form("%s_hHmpBigMap",AliRecoParam::GetEventSpecieName(specie)))); | |
438 | if(h2) { | |
439 | if( h2->Integral() > 1 ) { | |
440 | // no entres -> fatal | |
23874fc6 | 441 | // if( h2->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;} |
74db01e8 | 442 | h2->SetStats(0); |
443 | // clean up the text, stat, lines, ... | |
444 | if( h2->GetListOfFunctions() ) h2->GetListOfFunctions()->Clear(); | |
445 | ||
446 | //calculate missing pad fraction | |
447 | for(Int_t ich = 0; ich < 7; ich++) | |
448 | { | |
449 | for(Int_t iy=1+ich*144;iy<=144+ich*144;iy++) { | |
450 | for(Int_t ix=1;ix<=80;ix++) if(h2->GetBinContent(ix,iy) > 0) sumPadMapChROL[ich]++; | |
451 | for(Int_t ix=81;ix<=160;ix++) if(h2->GetBinContent(ix,iy) > 0) sumPadMapChROR[ich]++; | |
452 | } | |
453 | }//ch loop | |
454 | //check the calculated missing pad fraction | |
455 | for(Int_t ich = 0; ich < 7; ich++) | |
456 | { | |
457 | ||
458 | bigMapFlag = AliQAv1::kINFO; | |
459 | if( (1-sumPadMapChROL[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionWarningThreshold || | |
460 | (1-sumPadMapChROR[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionWarningThreshold ) | |
461 | { | |
462 | bigMapFlag = AliQAv1::kWARNING; | |
463 | } | |
464 | if( (1-sumPadMapChROL[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionErrorThreshold || | |
465 | (1-sumPadMapChROR[ich]/1.0/11520) > fHmpQaThr_MissingPadFractionErrorThreshold ) { | |
466 | bigMapFlag = AliQAv1::kERROR; | |
467 | errCntBigMap++; | |
468 | } | |
469 | ||
470 | }//ch loop | |
471 | if( errCntBigMap > 0 ) bigMapFlag = AliQAv1::kERROR; | |
472 | ||
473 | //update labels | |
474 | if ( bigMapFlag == AliQAv1::kINFO ) | |
475 | { | |
476 | hmpQaFlags[2] = AliQAv1::kINFO; | |
477 | text.Clear(); | |
478 | text.AddText(Form("OK (%d)",fIsOnlineThr)); | |
479 | text.SetFillColor(kGreen); | |
480 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
481 | } | |
482 | else if ( bigMapFlag == AliQAv1::kWARNING ) | |
483 | { | |
484 | hmpQaFlags[2] = AliQAv1::kWARNING; | |
485 | text.Clear(); | |
486 | text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr)); | |
487 | text.SetFillColor(kOrange); | |
488 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
489 | } | |
490 | else if ( bigMapFlag == AliQAv1::kERROR ) | |
491 | { | |
492 | hmpQaFlags[2] = AliQAv1::kERROR; | |
493 | text.Clear(); | |
494 | text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr)); | |
495 | text.SetFillColor(kRed); | |
496 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
497 | } | |
498 | else | |
499 | { | |
500 | hmpQaFlags[2] = AliQAv1::kFATAL; | |
501 | text.Clear(); | |
502 | text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr)); | |
503 | text.SetFillColor(kRed); | |
504 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
505 | } | |
506 | } | |
507 | } | |
508 | }//___HmpBigMap | |
509 | ||
510 | ||
511 | Int_t numSectorsMissing = 0, numSectorsGainLoss = 0; | |
512 | Double_t hSumQPerSector[42]={0}; | |
513 | ||
514 | if( list->FindObject(Form("%s_fHmpHvSectorQ",AliRecoParam::GetEventSpecieName(specie)))) | |
515 | { | |
516 | ||
517 | TH2* h2 = dynamic_cast<TH2*>( list->FindObject(Form("%s_fHmpHvSectorQ",AliRecoParam::GetEventSpecieName(specie)))); | |
518 | if(h2) { | |
519 | if(h2->Integral() > 0 ) { | |
520 | // no entres -> fatal | |
23874fc6 | 521 | // if( h2->GetEntries() == 0) {raqQualFlag = AliQAv1::kFATAL;} |
74db01e8 | 522 | h2->SetStats(0); |
523 | // clean up the text, stat, lines, ... | |
524 | if( h2->GetListOfFunctions() ) h2->GetListOfFunctions()->Clear(); | |
525 | ||
526 | //___ check sectors | |
527 | for(Int_t isec = 1 ; isec <= 42; isec++) | |
528 | { | |
529 | for(Int_t ibiny=100;ibiny<410;ibiny++) {hSumQPerSector[isec-1] += h2->GetBinContent(ibiny,isec); } | |
530 | if(hSumQPerSector[isec-1] < 0.001) {numSectorsGainLoss++;} // there is no photon and mip peak, gain loss | |
531 | if(h2->GetBinContent(1,isec) < 0.01 ) {numSectorsMissing++; } //practically there is no charge , the sector is missing | |
532 | } | |
533 | Int_t sectorErrors = numSectorsGainLoss+numSectorsMissing; | |
534 | ||
535 | ||
536 | if ( sectorErrors <= 3) | |
537 | { | |
538 | hmpQaFlags[3] = AliQAv1::kINFO; | |
539 | text.Clear(); | |
540 | text.AddText(Form("OK (%d)",fIsOnlineThr)); | |
541 | text.SetFillColor(kGreen); | |
542 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
543 | } | |
544 | else if ( sectorErrors > fHmpQaThr_SectorGainLossWarningThreshold) | |
545 | { | |
546 | hmpQaFlags[3] = AliQAv1::kWARNING; | |
547 | text.Clear(); | |
548 | text.AddText(Form("WARNING CHECK TWIKI (%d)",fIsOnlineThr)); | |
549 | if(numSectorsMissing > 0 ) text.AddText(Form("MISSING SECTOR?")); | |
550 | text.SetFillColor(kOrange); | |
551 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
552 | } | |
553 | else if ( sectorErrors > fHmpQaThr_SectorGainLossErrorThreshold) | |
554 | { | |
555 | hmpQaFlags[3] = AliQAv1::kERROR; | |
556 | text.Clear(); | |
557 | text.AddText(Form("ERROR CALL ONCALL (%d)",fIsOnlineThr)); | |
558 | if(numSectorsMissing > 0 ) text.AddText(Form("MISSING SECTOR?")); | |
559 | text.SetFillColor(kRed); | |
560 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
561 | } | |
562 | else | |
563 | { | |
564 | hmpQaFlags[3] = AliQAv1::kFATAL; | |
565 | text.Clear(); | |
566 | text.AddText(Form("FATAL CALL ONCALL (%d)",fIsOnlineThr)); | |
567 | text.SetFillColor(kRed); | |
568 | h2->GetListOfFunctions()->Add((TPaveText*)text.Clone()); | |
569 | } | |
570 | } | |
571 | } | |
572 | } | |
573 | ||
574 | ||
575 | //del lines, ... | |
576 | lineDdlDataSizeMin->Delete(); | |
577 | lineDdlDataSizeMax->Delete(); | |
578 | linePadOccupMin->Delete(); | |
579 | linePadOccupMax->Delete(); | |
580 | ||
581 | Double_t dflag = -1; | |
582 | switch ( TMath::MaxElement(4,hmpQaFlags)) | |
583 | { | |
584 | case AliQAv1::kINFO: | |
585 | dflag = 1.0; | |
586 | ||
587 | break; | |
588 | case AliQAv1::kWARNING: | |
589 | dflag = 0.75; | |
590 | ||
591 | break; | |
592 | case AliQAv1::kERROR: | |
593 | dflag = 0.25; | |
594 | ||
595 | break; | |
596 | case AliQAv1::kFATAL: | |
597 | dflag = -1.0; | |
598 | ||
599 | break; | |
600 | default: | |
601 | dflag = AliQAv1::kNULLBit; | |
602 | ||
603 | break; | |
604 | } | |
605 | ||
606 | return dflag; | |
607 | ||
608 | ||
609 | ||
610 | } | |
611 | //___________________________________________________________________________________________________ | |
612 | void AliHMPIDQAChecker::InitOnlineThresholds() | |
613 | { | |
614 | // | |
615 | // Init the online thresholds from GRP generated by AMORE | |
616 | // | |
303d90aa | 617 | |
74db01e8 | 618 | AliCDBManager* man = AliCDBManager::Instance(); |
619 | if(!man) { fIsOnlineThr = kFALSE; return; } | |
e683ec90 | 620 | |
74db01e8 | 621 | |
e683ec90 | 622 | if(!man->Get("GRP/Calib/QAThresholds")) { fIsOnlineThr = kFALSE; return; } |
74db01e8 | 623 | |
624 | AliCDBEntry* entry = man->Get("GRP/Calib/QAThresholds"); | |
625 | if(!entry) { fIsOnlineThr = kFALSE; return; } | |
626 | ||
627 | TObjArray* branch = (TObjArray*) entry->GetObject(); | |
628 | if(!branch ) { fIsOnlineThr = kFALSE; return; } | |
629 | ||
630 | AliQAThresholds* thresholds = (AliQAThresholds*) branch->FindObject("HMP"); | |
631 | if(!thresholds) { fIsOnlineThr = kFALSE; return; } | |
632 | else | |
633 | fIsOnlineThr = kTRUE; | |
634 | ||
303d90aa | 635 | |
74db01e8 | 636 | Int_t teCnt = 0; |
637 | TString parName = "zero"; | |
638 | while ( thresholds->GetThreshold(teCnt)) | |
639 | { | |
303d90aa | 640 | if(!((thresholds->GetThreshold(teCnt))->GetName())) return; |
641 | ||
74db01e8 | 642 | parName = (thresholds->GetThreshold(teCnt))->GetName(); |
303d90aa | 643 | |
74db01e8 | 644 | if( parName.Contains("HmpNumberOfExcludedDDLthreshold") ) |
645 | { | |
646 | TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt); | |
647 | fHmpQaThr_NumberOfExcludedDDL = myParam->GetVal(); | |
74db01e8 | 648 | } |
649 | ||
650 | ||
651 | if( parName.Contains("HmpDataSizeLowerThreshold") ) | |
652 | { | |
653 | TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt); | |
654 | fHmpQaThr_DataSizeLowerThreshold = myParam->GetVal(); | |
74db01e8 | 655 | } |
656 | ||
657 | if( parName.Contains("HmpDataSizeUpperThreshold") ) | |
658 | { | |
659 | TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt); | |
660 | fHmpQaThr_DataSizeUpperThreshold = myParam->GetVal(); | |
74db01e8 | 661 | } |
303d90aa | 662 | |
74db01e8 | 663 | if( parName.Contains("HmpPadOccupancyLowerThreshold") ) |
664 | { | |
665 | TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt); | |
666 | fHmpQaThr_PadOccupancyLowerThreshold = myParam->GetVal(); | |
74db01e8 | 667 | } |
668 | ||
669 | if( parName.Contains("HmpPadOccupancyUpperThreshold") ) | |
670 | { | |
671 | TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt); | |
672 | fHmpQaThr_PadOccupancyUpperThreshold = myParam->GetVal(); | |
74db01e8 | 673 | } |
303d90aa | 674 | |
74db01e8 | 675 | if( parName.Contains("HmpSectorGainLossWarningThreshold") ) |
676 | { | |
677 | TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt); | |
678 | fHmpQaThr_SectorGainLossWarningThreshold = myParam->GetVal(); | |
74db01e8 | 679 | } |
680 | ||
681 | if( parName.Contains("HmpSectorGainLossErrorThreshold") ) | |
682 | { | |
683 | TParameter<int>* myParam = (TParameter<int>*) thresholds->GetThreshold(teCnt); | |
684 | fHmpQaThr_SectorGainLossErrorThreshold = myParam->GetVal(); | |
74db01e8 | 685 | } |
686 | ||
687 | if( parName.Contains("HmpMissingPadFractionWarningThreshold") ) | |
688 | { | |
689 | TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt); | |
690 | fHmpQaThr_MissingPadFractionWarningThreshold = myParam->GetVal(); | |
74db01e8 | 691 | } |
692 | ||
693 | if( parName.Contains("HmpMissingPadFractionErrorThreshold") ) | |
694 | { | |
695 | TParameter<float>* myParam = (TParameter<float>*) thresholds->GetThreshold(teCnt); | |
696 | fHmpQaThr_MissingPadFractionErrorThreshold = myParam->GetVal(); | |
74db01e8 | 697 | } |
698 | ||
699 | teCnt++; | |
700 | }//while | |
701 | ||
702 | ||
703 | // PrintThresholds(); | |
704 | } | |
705 | //___________________________________________________________________________________________________ | |
706 | void AliHMPIDQAChecker::PrintThresholds() | |
707 | { | |
708 | Printf("--- Printing thresholds ---"); | |
94fb00f0 | 709 | Printf("--- Default or online: %i ---",fIsOnlineThr); |
710 | Printf("--- fHmpQaThr_NumberOfExcludedDDL: %i ---",fHmpQaThr_NumberOfExcludedDDL); | |
711 | Printf("--- fHmpQaThr_DataSizeLowerThreshold: %i ---",fHmpQaThr_DataSizeLowerThreshold); | |
712 | Printf("--- fHmpQaThr_DataSizeUpperThreshold: %i ---",fHmpQaThr_DataSizeUpperThreshold); | |
713 | Printf("--- fHmpQaThr_PadOccupancyLowerThreshold: %f ---",fHmpQaThr_PadOccupancyLowerThreshold); | |
714 | Printf("--- fHmpQaThr_PadOccupancyUpperThreshold: %f ---",fHmpQaThr_PadOccupancyUpperThreshold); | |
715 | Printf("--- fHmpQaThr_SectorGainLossWarningThreshold: %i ---",fHmpQaThr_SectorGainLossWarningThreshold); | |
716 | Printf("--- fHmpQaThr_SectorGainLossErrorThreshold: %i ---",fHmpQaThr_SectorGainLossErrorThreshold); | |
717 | Printf("--- fHmpQaThr_MissingPadFractionWarningThreshold: %f ---",fHmpQaThr_MissingPadFractionWarningThreshold); | |
718 | Printf("--- fHmpQaThr_MissingPadFractionErrorThreshold: %f ---",fHmpQaThr_MissingPadFractionErrorThreshold); | |
74db01e8 | 719 | Printf("--- Printing thresholds done ---"); |
720 | ||
721 | ||
722 | ||
723 | ||
724 | } | |
725 | //___________________________________________________________________________________________________ | |
726 | ||
e4e6ae61 | 727 |