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 // Base class for detectors quality assurance checkers
21 // Compares Data made by QADataMakers with reference data
22 // Y. Schutz CERN August 2007
25 // --- ROOT system ---
30 #include <TIterator.h>
35 #include <TParameter.h>
36 #include <TPaveText.h>
38 // --- Standard library ---
40 // --- AliRoot header files ---
41 #include "AliCDBEntry.h"
44 #include "AliQAChecker.h"
45 #include "AliQACheckerBase.h"
46 #include "AliQADataMaker.h"
47 #include "AliQAManager.h"
48 #include "AliDetectorRecoParam.h"
50 ClassImp(AliQACheckerBase)
53 //____________________________________________________________________________
54 AliQACheckerBase::AliQACheckerBase(const char * name, const char * title) :
58 fRefOCDBSubDir(new TObjArray*[AliRecoParam::kNSpecies]),
59 fLowTestValue(new Float_t[AliQAv1::kNBIT]),
60 fUpTestValue(new Float_t[AliQAv1::kNBIT]),
61 fImage(new TCanvas*[AliRecoParam::kNSpecies]),
63 fExternParamList(new TList())
66 fLowTestValue[AliQAv1::kINFO] = 0.5 ;
67 fUpTestValue[AliQAv1::kINFO] = 1.0 ;
68 fLowTestValue[AliQAv1::kWARNING] = 0.002 ;
69 fUpTestValue[AliQAv1::kWARNING] = 0.5 ;
70 fLowTestValue[AliQAv1::kERROR] = 0.0 ;
71 fUpTestValue[AliQAv1::kERROR] = 0.002 ;
72 fLowTestValue[AliQAv1::kFATAL] = -1.0 ;
73 fUpTestValue[AliQAv1::kFATAL] = 0.0 ;
75 AliDebug(AliQAv1::GetQADebugLevel(), "Default setting is:") ;
76 if ( AliDebugLevel() == AliQAv1::GetQADebugLevel() ) {
77 const Char_t * text= Form(" INFO -> %1.5f < value < %1.5f WARNING -> %1.5f < value <= %1.5f \n ERROR -> %1.5f < value <= %1.5f \n FATAL -> %1.5f <= value < %1.5f \n",
78 fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO],
79 fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING],
80 fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR],
81 fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ;
82 AliInfo(Form("%s", text)) ;
85 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
86 fImage[specie] = NULL ;
87 fRefOCDBSubDir[specie] = NULL ;
91 //____________________________________________________________________________
92 AliQACheckerBase::~AliQACheckerBase()
94 delete [] fLowTestValue ;
95 delete [] fUpTestValue ;
98 delete[] fRefOCDBSubDir ;
99 AliQAv1::GetQAResultFile()->Close() ;
100 if (fExternParamList) {
101 fExternParamList->Clear() ;
102 delete fExternParamList ;
106 //____________________________________________________________________________
107 void AliQACheckerBase::PrivateCheck(Double_t * test, AliQAv1::ALITASK_t index, const AliDetectorRecoParam * recoParam)
109 // Performs a basic checking
110 // Compares all the histograms stored in the directory
111 // With reference histograms either in a file of in OCDB
113 TObjArray ** list = new TObjArray *[AliRecoParam::kNSpecies] ;
115 for (specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
116 list[specie] = new TObjArray(AliQAv1::GetMaxQAObj()) ;
117 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) )
120 TList * keyList = fDataSubDir->GetListOfKeys() ;
121 TIter next(keyList) ;
123 while ( (key = static_cast<TKey *>(next())) ) {
124 TDirectory * specieDir = fDataSubDir->GetDirectory(key->GetName()) ;
125 TList * keykeyList = specieDir->GetListOfKeys() ;
126 TIter next2(keykeyList) ;
128 while ( (keykey = static_cast<TKey *>(next2())) ) {
129 TObject * odata = specieDir->Get(keykey->GetName()) ;
130 if ( odata->IsA()->InheritsFrom("TH1") ) {
131 TH1 * hdata = static_cast<TH1*>(odata) ;
132 list[specie]->Add(hdata) ;
133 } else if (!odata->IsA()->InheritsFrom("TDirectory")) // skip the expert directory
134 AliError(Form("%s Is a Classname that cannot be processed", key->GetClassName())) ;
140 Check(test, index, list, recoParam) ;
146 //____________________________________________________________________________
147 void AliQACheckerBase::Check(Double_t * test, AliQAv1::ALITASK_t task, TObjArray ** list, const AliDetectorRecoParam * /* recoParam */)
149 // Performs a basic checking
150 // Compares all the histograms in the list
152 Int_t count[AliRecoParam::kNSpecies] = { 0 };
154 GetRefSubDir(GetName(), AliQAv1::GetTaskName(task), fRefSubDir, fRefOCDBSubDir) ;
155 // SetRefandData(refDir, refOCDBDir) ;
157 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
159 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie)) continue ;
160 if (list[specie]->GetEntries() == 0)
161 test[specie] = 0. ; // nothing to check
163 if (!fRefSubDir && !fRefOCDBSubDir)
164 test[specie] = -1 ; // no reference data
166 TIter next(list[specie]) ;
169 while ( (hdata = static_cast<TH1 *>(next())) ) {
170 if ( hdata->IsA()->InheritsFrom("TH1") ) {
171 if ( hdata->TestBit(AliQAv1::GetExpertBit()) ) // does not perform the test for expert data
174 // First try to find the ref histo with exact name (with possible trigger clon ending)
175 TString hname = hdata->GetName();
177 if (fRefSubDir) href = static_cast<TH1*>(fRefSubDir->Get(hname.Data())) ;
178 else if (fRefOCDBSubDir[specie]) href = static_cast<TH1*>(fRefOCDBSubDir[specie]->FindObject(hname.Data()));
180 if (!href && hdata->TestBit(AliQAv1::GetClonedBit())) { // try to find the histo for the base name (w/o trigger ending
181 int ind = hname.Index(AliQADataMaker::GetTriggerPrefix());
184 if (fRefSubDir) href = static_cast<TH1*>(fRefSubDir->Get(hname.Data())) ;
185 else if (fRefOCDBSubDir[specie]) href = static_cast<TH1*>(fRefOCDBSubDir[specie]->FindObject(hname.Data()));
190 test[specie] = -1 ; // no reference data ;
192 Double_t rv = DiffK(hdata, href) ;
193 AliDebug(AliQAv1::GetQADebugLevel(), Form("%s ->Test = %f", hdata->GetName(), rv)) ;
198 AliError("Data type cannot be processed") ;
199 if (count[specie] != 0)
200 test[specie] /= count[specie] ;
208 //____________________________________________________________________________
209 void AliQACheckerBase::DeleteImages()
212 for (Int_t esIndex = 0 ; esIndex < AliRecoParam::kNSpecies ; esIndex++) {
213 if ( fImage[esIndex] ) {delete fImage[esIndex]; fImage[esIndex] = 0;}
214 if ( fRefOCDBSubDir[esIndex] ) {delete fRefOCDBSubDir[esIndex]; fRefOCDBSubDir[esIndex] = 0;}
218 //____________________________________________________________________________
219 Double_t AliQACheckerBase::DiffC(const TH1 * href, const TH1 * hin) const
221 // compares two histograms using the Chi2 test
222 if ( hin->Integral() == 0 ) {
223 AliDebug(AliQAv1::GetQADebugLevel(), Form("Spectrum %s is empty", hin->GetName())) ;
227 return hin->Chi2Test(href) ;
230 //____________________________________________________________________________
231 Double_t AliQACheckerBase::DiffK(const TH1 * href, const TH1 * hin) const
233 // compares two histograms using the Kolmogorov test
234 if ( hin->Integral() == 0 || href->Integral() == 0) {
235 AliDebug(AliQAv1::GetQADebugLevel(), Form("Spectrum %s or its reference is empty", hin->GetName())) ;
239 return hin->KolmogorovTest(href) ;
242 //_____________________________________________________________________________
243 void AliQACheckerBase::GetRefSubDir(const char * det, const char * task, TDirectory *& dirFile, TObjArray **& dirOCDB)
245 // Opens and returns the file with the reference data
247 TString refStorage(AliQAv1::GetQARefStorage()) ;
248 if (!refStorage.Contains(AliQAv1::GetLabLocalOCDB()) && !refStorage.Contains(AliQAv1::GetLabAliEnOCDB())) {
249 AliError(Form("%s is not a valid location for reference data", refStorage.Data())) ;
252 AliQAManager* manQA = AliQAManager::QAManager(AliQAv1::GetTaskIndex(task)) ;
253 // dirOCDB = new TObjArray*[AliRecoParam::kNSpecies] ;
254 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
255 dirOCDB[specie] = NULL ;
256 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) )
258 AliQAv1::SetQARefDataDirName(specie) ;
259 if ( ! manQA->GetLock() ) {
260 manQA->SetDefaultStorage(AliQAv1::GetQARefStorage()) ;
261 manQA->SetSpecificStorage("*", AliQAv1::GetQARefStorage()) ;
262 manQA->SetRun(AliCDBManager::Instance()->GetRun()) ;
265 char * detOCDBDir = Form("%s/%s/%s", det, AliQAv1::GetRefOCDBDirName(), AliQAv1::GetRefDataDirName()) ;
266 AliCDBEntry * entry = manQA->Get(detOCDBDir, manQA->GetRun()) ;
268 TList * listDetQAD =static_cast<TList *>(entry->GetObject()) ;
269 if ( listDetQAD && strcmp(listDetQAD->ClassName(), "TList") != 0 ) {
270 AliError(Form("Expected a Tlist and found a %s for detector %s", listDetQAD->ClassName(), det)) ;
275 TIter next(listDetQAD) ;
276 while ( (TObjArray*)next() )
277 dirOCDB[specie] = static_cast<TObjArray *>(listDetQAD->FindObject(Form("%s/%s", task, AliRecoParam::GetEventSpecieName(specie)))) ;
284 //____________________________________________________________________________
285 void AliQACheckerBase::PrintExternParam()
287 // Print the list of external parameter list
288 TIter next(fExternParamList) ;
289 TParameter<double> *pp ;
290 TString printit("\n") ;
291 while( (pp = (TParameter<double>*)next()) )
292 printit += Form("%s = %f\n", pp->GetName(), pp->GetVal());
293 AliInfo(Form("%s", printit.Data())) ;
296 //____________________________________________________________________________
297 void AliQACheckerBase::Run(AliQAv1::ALITASK_t index, const AliDetectorRecoParam * recoParam)
299 //Run the checker for all kind of species
300 AliDebug(AliQAv1::GetQADebugLevel(), Form("Processing %s", AliQAv1::GetAliTaskName(index))) ;
302 Double_t * rv = new Double_t[AliRecoParam::kNSpecies] ;
303 for (int i=AliRecoParam::kNSpecies;i--;) rv[i] = 0.0;
304 PrivateCheck(rv, index, recoParam) ;
307 AliDebug(AliQAv1::GetQADebugLevel(), Form("Test result of %s", AliQAv1::GetAliTaskName(index))) ;
313 //____________________________________________________________________________
314 void AliQACheckerBase::Run(AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * recoParam)
316 // RS: perform check for all trigger classes in loop
317 Double_t * rv = new Double_t[AliRecoParam::kNSpecies] ;
319 TObjArray ** listTrig = new TObjArray *[AliRecoParam::kNSpecies];
321 for (int itc=-1;itc<AliQADataMaker::GetNTrigClasses();itc++) {
323 // RS: fetch the histograms for each specie and for given trigger
324 //AliInfo(Form("Processing %s for trigger: %s", AliQAv1::GetAliTaskName(index),AliQADataMaker::GetTrigClassName(itc)));
326 for (int specie=0;specie<AliRecoParam::kNSpecies;specie++) {
327 listTrig[specie] = 0;
328 if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) || !list[specie]) continue;
329 listTrig[specie] = new TObjArray( list[specie]->GetSize() ); // destination for clones of this trigger
330 AliQADataMaker::GetDataOfTrigClass(list[specie],itc, listTrig[specie]);
332 AliDebug(AliQAv1::GetQADebugLevel(), Form("Processing %s for trigger: %s", AliQAv1::GetAliTaskName(index),AliQADataMaker::GetTrigClassName(itc)));
333 Check(rv, index, listTrig, recoParam) ;
335 AliDebug(AliQAv1::GetQADebugLevel(), Form("Test result of %s for trigger: %s", AliQAv1::GetAliTaskName(index),AliQADataMaker::GetTrigClassName(itc)));
337 for (int specie=0;specie<AliRecoParam::kNSpecies;specie++) if (listTrig[specie]) delete listTrig[specie]; // clean temporary container
344 //____________________________________________________________________________
345 void AliQACheckerBase::Finish() const
347 // wrap up and save QA in proper file
348 AliQAv1::GetQAResultFile() ;
349 AliQAv1 * qa = AliQAv1::Instance() ;
350 qa->Write(AliQAv1::GetQAName(), kWriteDelete) ;
353 //____________________________________________________________________________
354 void AliQACheckerBase::MakeImage( TObjArray ** list, AliQAv1::TASKINDEX_t task, AliQAv1::MODE_t mode)
356 // makes the QA image for sim and rec
357 TObjArray tmpArr; // array to store flat version of original array (which may contain clones)
359 for (Int_t esIndex = 0; esIndex < AliRecoParam::kNSpecies; esIndex++) {
360 if (! AliQAv1::Instance(AliQAv1::GetDetIndex(GetName()))->IsEventSpecieSet(AliRecoParam::ConvertIndex(esIndex)) || list[esIndex]->GetEntries() == 0) continue;
362 TIter next(list[esIndex]);
363 TObject* hdata = NULL;
365 while ( (hdata=(next())) ) { // count histos and transfere to flat array
366 if (hdata->InheritsFrom(TH1::Class()) && hdata->TestBit(AliQAv1::GetImageBit()) ) { // histo, not cloned
368 tmpArr.AddLast(hdata);
371 if (!hdata->TestBit(AliQAv1::GetClonedBit())) continue; // not an array of clones, unknown object
372 TIter nextCl((TObjArray*)hdata); // array of histo clones
374 while ((hcl=nextCl())) if (hcl->InheritsFrom(TH1::Class()) && hcl->TestBit(AliQAv1::GetImageBit())) {tmpArr.AddLast(hcl); nImages++;}
377 if ( nImages == 0 ) {
378 AliDebug(AliQAv1::GetQADebugLevel(), Form("No histogram will be plotted for %s %s %s\n", GetName(), AliQAv1::GetTaskName(task).Data(), AliRecoParam::GetEventSpecieName(esIndex)));
381 AliDebug(AliQAv1::GetQADebugLevel(), Form("%d histograms will be plotted for %s %s %s\n", nImages, GetName(), AliQAv1::GetTaskName(task).Data(),AliRecoParam::GetEventSpecieName(esIndex)));
383 const Char_t * title = Form("QA_%s_%s_%s", GetName(), AliQAv1::GetTaskName(task).Data(), AliRecoParam::GetEventSpecieName(esIndex));
385 if ( !fImage[esIndex] ) fImage[esIndex] = new TCanvas(title, title);
387 fImage[esIndex]->Clear();
388 fImage[esIndex]->SetTitle(title);
389 fImage[esIndex]->cd();
390 TPaveText someText(0.015, 0.015, 0.98, 0.98);
391 someText.AddText(title);
393 fImage[esIndex]->Print(Form("%s%s%d.%s", AliQAv1::GetImageFileName(), AliQAv1::GetModeName(mode), AliQAChecker::Instance()->GetRunNumber(), AliQAv1::GetImageFileFormat()));
394 fImage[esIndex]->Clear();
395 Int_t nx = TMath::Nint(TMath::Sqrt(nImages));
397 if (nx < TMath::Sqrt(nImages)) ny++;
399 fImage[esIndex]->Divide(nx, ny);
400 TIter nexthist(&tmpArr);
402 fImage[esIndex]->cd(npad);
404 while ( (histo=(TH1*)nexthist()) ) { // tmpArr is guaranteed to contain only plottable histos, no checks needed
405 TString opts = histo->GetDrawOption();
406 if (opts.Contains("logy",TString::kIgnoreCase)) gPad->SetLogy();
407 if (opts.Contains("logx",TString::kIgnoreCase)) gPad->SetLogx();
409 fImage[esIndex]->cd(++npad);
411 fImage[esIndex]->Print(Form("%s%s%d.%s", AliQAv1::GetImageFileName(), AliQAv1::GetModeName(mode), AliQAChecker::Instance()->GetRunNumber(), AliQAv1::GetImageFileFormat()), "ps");
415 //____________________________________________________________________________
416 void AliQACheckerBase::SetHiLo(Float_t * hiValue, Float_t * lowValue)
418 AliDebug(AliQAv1::GetQADebugLevel(), "Previous setting was:") ;
419 if ( AliDebugLevel() == AliQAv1::GetQADebugLevel() ) {
420 const Char_t * text= Form(" INFO -> %1.5f < value < %1.5f WARNING -> %1.5f < value <= %1.5f \n ERROR -> %1.5f < value <= %1.5f \n FATAL -> %1.5f <= value < %1.5f \n",
421 fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO],
422 fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING],
423 fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR],
424 fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ;
425 AliInfo(Form("%s", text)) ;
428 for (Int_t index = 0 ; index < AliQAv1::kNBIT ; index++) {
429 fLowTestValue[index] = lowValue[index] ;
430 fUpTestValue[index] = hiValue[index] ;
432 AliDebug(AliQAv1::GetQADebugLevel(), "Current setting is:") ;
433 if ( AliDebugLevel() == AliQAv1::GetQADebugLevel() ) {
434 const Char_t * text= Form(" INFO -> %1.5f < value < %1.5f WARNING -> %1.5f < value <= %1.5f \n ERROR -> %1.5f < value <= %1.5f \n FATAL -> %1.5f <= value < %1.5f \n",
435 fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO],
436 fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING],
437 fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR],
438 fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ; AliInfo(Form("%s", text)) ;
442 //____________________________________________________________________________
443 void AliQACheckerBase::SetQA(AliQAv1::ALITASK_t index, Double_t * value) const
445 // sets the QA according the return value of the Check
447 AliQAv1 * qa = AliQAv1::Instance(index) ;
449 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
450 if (! qa->IsEventSpecieSet(AliRecoParam::ConvertIndex(specie)))
452 if ( value == NULL ) { // No checker is implemented, set all QA to Fatal
453 qa->Set(AliQAv1::kFATAL, specie) ;
455 if ( value[specie] >= fLowTestValue[AliQAv1::kFATAL] && value[specie] < fUpTestValue[AliQAv1::kFATAL] )
456 qa->Set(AliQAv1::kFATAL, AliRecoParam::ConvertIndex(specie)) ;
457 else if ( value[specie] > fLowTestValue[AliQAv1::kERROR] && value[specie] <= fUpTestValue[AliQAv1::kERROR] )
458 qa->Set(AliQAv1::kERROR, AliRecoParam::ConvertIndex(specie)) ;
459 else if ( value[specie] > fLowTestValue[AliQAv1::kWARNING] && value[specie] <= fUpTestValue[AliQAv1::kWARNING] )
460 qa->Set(AliQAv1::kWARNING, AliRecoParam::ConvertIndex(specie)) ;
461 else if ( value[specie] > fLowTestValue[AliQAv1::kINFO] && value[specie] <= fUpTestValue[AliQAv1::kINFO] )
462 qa->Set(AliQAv1::kINFO, AliRecoParam::ConvertIndex(specie)) ;