]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - STEER/AliQACheckerBase.cxx
Fix for Coverity 10007: BAD_FREE
[u/mrichter/AliRoot.git] / STEER / AliQACheckerBase.cxx
index c2e44244e0d4699d56b7e4663c98e680ec3835ae..6786975543db6a9ecb167223e0c897e827d61b10 100644 (file)
@@ -23,6 +23,7 @@
 //
 
 // --- ROOT system ---
+#include <TCanvas.h>
 #include <TClass.h>
 #include <TH1F.h> 
 #include <TH1I.h> 
 #include <TFile.h> 
 #include <TList.h>
 #include <TNtupleD.h>
+#include <TParameter.h>
+#include <TPaveText.h>
 
 // --- Standard library ---
 
 // --- AliRoot header files ---
+#include "AliCDBEntry.h"
 #include "AliLog.h"
-#include "AliQA.h"
+#include "AliQAv1.h"
 #include "AliQAChecker.h"
 #include "AliQACheckerBase.h"
 #include "AliQADataMaker.h"
+#include "AliQAManager.h"
+#include "AliDetectorRecoParam.h"
 
 ClassImp(AliQACheckerBase)
 
@@ -49,9 +55,37 @@ AliQACheckerBase::AliQACheckerBase(const char * name, const char * title) :
   TNamed(name, title), 
   fDataSubDir(0x0),
   fRefSubDir(0x0), 
-  fRefOCDBSubDir(0x0)
+  fRefOCDBSubDir(new TObjArray*[AliRecoParam::kNSpecies]), 
+  fLowTestValue(new Float_t[AliQAv1::kNBIT]),
+  fUpTestValue(new Float_t[AliQAv1::kNBIT]),
+  fImage(new TCanvas*[AliRecoParam::kNSpecies]), 
+  fPrintImage(kTRUE), 
+  fExternParamList(new TList())
 {
   // ctor
+  fLowTestValue[AliQAv1::kINFO]    =  0.5   ; 
+  fUpTestValue[AliQAv1::kINFO]     = 1.0 ; 
+  fLowTestValue[AliQAv1::kWARNING] =  0.002 ; 
+  fUpTestValue[AliQAv1::kWARNING]  = 0.5 ; 
+  fLowTestValue[AliQAv1::kERROR]   =  0.0   ; 
+  fUpTestValue[AliQAv1::kERROR]    = 0.002 ; 
+  fLowTestValue[AliQAv1::kFATAL]   = -1.0   ; 
+  fUpTestValue[AliQAv1::kFATAL]    = 0.0 ; 
+  
+  AliDebug(AliQAv1::GetQADebugLevel(), "Default setting is:") ;
+  if ( AliDebugLevel()  == AliQAv1::GetQADebugLevel() ) {
+    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", 
+                              fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO], 
+                              fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING], 
+                              fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR], 
+                              fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ; 
+    AliInfo(Form("%s", text)) ; 
+  }
+  
+  for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+    fImage[specie] = NULL ; 
+    fRefOCDBSubDir[specie] = NULL ;
+  }
 }
 
 //____________________________________________________________________________ 
@@ -59,120 +93,157 @@ AliQACheckerBase::AliQACheckerBase(const AliQACheckerBase& qac) :
   TNamed(qac.GetName(), qac.GetTitle()),
   fDataSubDir(qac.fDataSubDir), 
   fRefSubDir(qac.fRefSubDir), 
-  fRefOCDBSubDir(qac.fRefOCDBSubDir)
+  fRefOCDBSubDir(qac.fRefOCDBSubDir), 
+  fLowTestValue(new Float_t[AliQAv1::kNBIT]),
+  fUpTestValue(new Float_t[AliQAv1::kNBIT]), 
+  fImage(new TCanvas*[AliRecoParam::kNSpecies]),  
+  fPrintImage(kTRUE), 
+  fExternParamList(new TList())  
 {
   //copy ctor
-    
+  for (Int_t index = 0 ; index < AliQAv1::kNBIT ; index++) {
+    fLowTestValue[index]  = qac.fLowTestValue[index] ; 
+    fUpTestValue[index] = qac.fUpTestValue[index] ; 
+  }
+    for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+      fImage[specie] = qac.fImage[specie] ; 
+      fRefOCDBSubDir[specie] = qac.fRefOCDBSubDir[specie] ; 
+    }
+  if (qac.fExternParamList) {
+    TIter next(qac.fExternParamList) ; 
+    TParameter<double> * p ; 
+    while ( (p = (TParameter<double>*)next()) )
+      fExternParamList->Add(p) ;
+  }
 }
 
 //____________________________________________________________________________
-AliQACheckerBase& AliQACheckerBase::operator = (const AliQACheckerBase& qadm )
+AliQACheckerBase& AliQACheckerBase::operator = (const AliQACheckerBase& qac )
 {
   // Equal operator.
   this->~AliQACheckerBase();
-  new(this) AliQACheckerBase(qadm);
+  new(this) AliQACheckerBase(qac);
   return *this;
 }
 
+//____________________________________________________________________________ 
+AliQACheckerBase::~AliQACheckerBase()
+{
+  delete [] fLowTestValue ; 
+  delete [] fUpTestValue ; 
+  for (Int_t esIndex = 0 ; esIndex < AliRecoParam::kNSpecies ; esIndex++) {
+    if ( fImage[esIndex] ) 
+      delete fImage[esIndex] ;
+    if ( fRefOCDBSubDir[esIndex] ) 
+      delete fRefOCDBSubDir[esIndex] ; 
+  }
+  delete[] fImage ; 
+  delete[] fRefOCDBSubDir ; 
+  AliQAv1::GetQAResultFile()->Close() ; 
+  if (fExternParamList) {
+    fExternParamList->Clear() ; 
+    delete fExternParamList ; 
+  }
+}
+
 //____________________________________________________________________________
-const Double_t AliQACheckerBase::Check(AliQA::ALITASK_t /*index*/
+void AliQACheckerBase::Check(Double_t * test, AliQAv1::ALITASK_t index, const AliDetectorRecoParam * recoParam
 {
   // Performs a basic checking
   // Compares all the histograms stored in the directory
   // With reference histograms either in a file of in OCDB  
 
-       Double_t test = 0.0  ;
-       Int_t count = 0 ; 
-
-       if (!fDataSubDir)  
-               test = 1. ; // nothing to check
-       else 
-               if (!fRefSubDir && !fRefOCDBSubDir)
-                       test = -1 ; // no reference data
-               else {
-                       TList * keyList = fDataSubDir->GetListOfKeys() ; 
-                       TIter next(keyList) ; 
-                       TKey * key ;
-                       count = 0 ; 
-                       while ( (key = static_cast<TKey *>(next())) ) {
-                               TObject * odata = fRefSubDir->Get(key->GetName()) ; 
-                               if ( odata->IsA()->InheritsFrom("TH1") ) {
-                                       TH1 * hdata = static_cast<TH1*>(odata) ;
-                                       TH1 * href = NULL ; 
-                                       if (fRefSubDir) 
-                                               href  = static_cast<TH1*>(fRefSubDir->Get(key->GetName())) ;
-                                       else if (fRefOCDBSubDir) {
-                                               href  = static_cast<TH1*>(fRefOCDBSubDir->FindObject(key->GetName())) ;
-                                       }
-                                       if (!href) 
-                                               test = -1 ; // no reference data ; 
-                                       else {
-                                               Double_t rv =  DiffK(hdata, href) ;
-                                               AliInfo(Form("%s ->Test = %f", hdata->GetName(), rv)) ; 
-                                               test += rv ; 
-                                               count++ ; 
-                                       }
-                               } else
-                                       AliError(Form("%s Is a Classname that cannot be processed", key->GetClassName())) ;
-                       }
-               } 
-       
-       if (count != 0) 
-               test /= count ;
-   
-       return test ;
+  TObjArray ** list = new TObjArray *[AliRecoParam::kNSpecies] ; 
+  Int_t specie ;
+  for (specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+    list[specie] =  new TObjArray(AliQAv1::GetMaxQAObj()) ; 
+    if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) 
+      continue ; 
+    if (fDataSubDir) {
+      TList * keyList = fDataSubDir->GetListOfKeys() ; 
+      TIter next(keyList) ; 
+      TKey * key ;
+      while ( (key = static_cast<TKey *>(next())) ) {
+        TDirectory * specieDir = fDataSubDir->GetDirectory(key->GetName()) ; 
+        TList * keykeyList = specieDir->GetListOfKeys() ; 
+        TIter next2(keykeyList) ; 
+        TKey * keykey ;
+        while ( (keykey = static_cast<TKey *>(next2())) ) {
+          TObject * odata = specieDir->Get(keykey->GetName()) ; 
+          if ( odata->IsA()->InheritsFrom("TH1") ) {
+            TH1 * hdata = static_cast<TH1*>(odata) ;
+            list[specie]->Add(hdata) ; 
+          } else if (!odata->IsA()->InheritsFrom("TDirectory")) // skip the expert directory
+            AliError(Form("%s Is a Classname that cannot be processed", key->GetClassName())) ;
+        }
+      }
+    }
+  }
+  Check(test, index, list, recoParam) ;
+  
+  delete[] list ; 
+    
 }  
 
 //____________________________________________________________________________
-const Double_t AliQACheckerBase::Check(AliQA::ALITASK_t /*index*/, TObjArray * list
+void AliQACheckerBase::Check(Double_t * test, AliQAv1::ALITASK_t task, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/
 {
   // Performs a basic checking
   // Compares all the histograms in the list
 
-       Double_t test = 0.0  ;
-       Int_t count = 0 ; 
-
-       if (list->GetEntries() == 0)  
-               test = 1. ; // nothing to check
-       else {
-               if (!fRefSubDir)
-                       test = -1 ; // no reference data
-               else {
-                       TIter next(list) ; 
-                       TH1 * hdata ;
-                       count = 0 ; 
-                       while ( (hdata = dynamic_cast<TH1 *>(next())) ) {
-                               if ( hdata) { 
-                                       TH1 * href = NULL ; 
-                                       if (fRefSubDir) 
-                                               href  = static_cast<TH1*>(fRefSubDir->Get(hdata->GetName())) ;
-                                       else if (fRefOCDBSubDir)
-                                               href  = static_cast<TH1*>(fRefOCDBSubDir->FindObject(hdata->GetName())) ;
-                                       if (!href) 
-                                               test = -1 ; // no reference data ; 
-                                       else {
-                                               Double_t rv =  DiffK(hdata, href) ;
-                                               AliInfo(Form("%s ->Test = %f", hdata->GetName(), rv)) ; 
-                                               test += rv ; 
-                                               count++ ; 
-                                       }
-                               } 
-                               else
-                                       AliError("Data type cannot be processed") ;
-                       }
-               }
-       }
-       if (count != 0) 
-               test /= count ;
-       return test ;
+       Int_t count[AliRecoParam::kNSpecies]   = { 0 }; 
+
+  GetRefSubDir(GetName(), AliQAv1::GetTaskName(task), fRefSubDir, fRefOCDBSubDir) ;
+ // SetRefandData(refDir, refOCDBDir) ; 
+  
+  for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+    test[specie] = 1.0 ; 
+    if ( !AliQAv1::Instance()->IsEventSpecieSet(specie)) 
+      continue ; 
+    if (list[specie]->GetEntries() == 0)  
+      test[specie] = 0. ; // nothing to check
+    else {
+      if (!fRefSubDir && !fRefOCDBSubDir)
+        test[specie] = -1 ; // no reference data
+      else {
+        TIter next(list[specie]) ; 
+        TH1 * hdata ;
+        count[specie] = 0 ; 
+        while ( (hdata = static_cast<TH1 *>(next())) ) {
+          if ( hdata->IsA()->InheritsFrom("TH1") ) {
+            if ( hdata->TestBit(AliQAv1::GetExpertBit()) )  // does not perform the test for expert data
+              continue ; 
+            TH1 * href = NULL ; 
+            if (fRefSubDir) 
+              href  = static_cast<TH1*>(fRefSubDir->Get(hdata->GetName())) ;
+            else if (fRefOCDBSubDir[specie])
+              href  = static_cast<TH1*>(fRefOCDBSubDir[specie]->FindObject(hdata->GetName())) ;
+            if (!href) 
+              test[specie] = -1 ; // no reference data ; 
+            else {
+              Double_t rv =  DiffK(hdata, href) ;
+              AliDebug(AliQAv1::GetQADebugLevel(), Form("%s ->Test = %f", hdata->GetName(), rv)) ; 
+              test[specie] += rv ; 
+              count[specie]++ ;
+            }
+          } else
+            AliError("Data type cannot be processed") ;
+          if (count[specie] != 0) 
+            test[specie] /= count[specie] ;
+        }
+      }
+    }
+  }
 }  
 
+
 //____________________________________________________________________________ 
-const Double_t AliQACheckerBase::DiffC(const TH1 * href, const TH1 * hin) const
+Double_t AliQACheckerBase::DiffC(const TH1 * href, const TH1 * hin) const
 {
   // compares two histograms using the Chi2 test
   if ( hin->Integral() == 0 ) {
-    AliWarning(Form("Spectrum %s is empty", hin->GetName())) ; 
+    AliDebug(AliQAv1::GetQADebugLevel(), Form("Spectrum %s is empty", hin->GetName())) ; 
     return 0. ;
   }
     
@@ -180,47 +251,98 @@ const Double_t AliQACheckerBase::DiffC(const TH1 * href, const TH1 * hin) const
 }
 
 //____________________________________________________________________________ 
-const Double_t AliQACheckerBase::DiffK(const TH1 * href, const TH1 * hin) const
+Double_t AliQACheckerBase::DiffK(const TH1 * href, const TH1 * hin) const
 {
   // compares two histograms using the Kolmogorov test
-  if ( hin->Integral() == 0 ) {
-    AliWarning(Form("Spectrum %s is empty", hin->GetName())) ; 
+  if ( hin->Integral() == 0 || href->Integral() == 0) {
+    AliDebug(AliQAv1::GetQADebugLevel(), Form("Spectrum %s or its reference is empty", hin->GetName())) ; 
     return 0. ;
   }
     
   return hin->KolmogorovTest(href) ;  
 }
 
-//____________________________________________________________________________ 
-void AliQACheckerBase::Init(const AliQA::DETECTORINDEX_t det)
-{
-  AliQA::Instance(det) ; 
+  //_____________________________________________________________________________
+void AliQACheckerBase::GetRefSubDir(const char * det, const char * task, TDirectory *& dirFile, TObjArray **& dirOCDB)     
+{ 
+    // Opens and returns the file with the reference data 
+  dirFile = NULL ; 
+  TString refStorage(AliQAv1::GetQARefStorage()) ;
+  if (!refStorage.Contains(AliQAv1::GetLabLocalOCDB()) && !refStorage.Contains(AliQAv1::GetLabAliEnOCDB())) {
+    AliError(Form("%s is not a valid location for reference data", refStorage.Data())) ; 
+    return ; 
+  } else {
+    AliQAManager* manQA = AliQAManager::QAManager(AliQAv1::GetTaskIndex(task)) ;
+      //    dirOCDB = new TObjArray*[AliRecoParam::kNSpecies] ;        
+    for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+      dirOCDB[specie] = NULL ; 
+      if ( !AliQAv1::Instance()->IsEventSpecieSet(specie) ) 
+        continue ; 
+      AliQAv1::SetQARefDataDirName(specie) ;
+      if ( ! manQA->GetLock() ) { 
+        manQA->SetDefaultStorage(AliQAv1::GetQARefStorage()) ; 
+        manQA->SetSpecificStorage("*", AliQAv1::GetQARefStorage()) ;
+        manQA->SetRun(AliCDBManager::Instance()->GetRun()) ; 
+        manQA->SetLock() ; 
+      }
+      char * detOCDBDir = Form("%s/%s/%s", det, AliQAv1::GetRefOCDBDirName(), AliQAv1::GetRefDataDirName()) ; 
+      AliCDBEntry * entry = manQA->Get(detOCDBDir, manQA->GetRun()) ;
+      if (entry) {
+        TList * listDetQAD =static_cast<TList *>(entry->GetObject()) ;
+        if ( strcmp(listDetQAD->ClassName(), "TList") != 0 ) {
+          AliError(Form("Expected a Tlist and found a %s for detector %s", listDetQAD->ClassName(), det)) ; 
+          listDetQAD = NULL ; 
+          continue ; 
+        } 
+        if ( listDetQAD ) {
+          TIter next(listDetQAD) ;
+          while ( (TObjArray*)next() ) 
+            dirOCDB[specie] = static_cast<TObjArray *>(listDetQAD->FindObject(Form("%s/%s", task, AliRecoParam::GetEventSpecieName(specie)))) ;             
+        }
+      }
+    }
+  }
 }
 
 //____________________________________________________________________________
-void AliQACheckerBase::Run(AliQA::ALITASK_t index, TObject * obj) 
+void AliQACheckerBase::PrintExternParam() 
+{
+    // Print the list of external parameter list
+  TIter next(fExternParamList) ; 
+  TParameter<double> *pp ; 
+  TString printit("\n") ;
+  while( (pp = (TParameter<double>*)next()) )
+    printit += Form("%s = %f\n", pp->GetName(), pp->GetVal());  
+  AliInfo(Form("%s", printit.Data())) ;
+}
+  
+//____________________________________________________________________________
+void AliQACheckerBase::Run(AliQAv1::ALITASK_t index, AliDetectorRecoParam * recoParam) 
 { 
-       AliDebug(1, Form("Processing %s", AliQA::GetAliTaskName(index))) ; 
+       AliDebug(AliQAv1::GetQADebugLevel(), Form("Processing %s", AliQAv1::GetAliTaskName(index))) ; 
   
-       Double_t rv = -1 ;      
-  if ( !obj ) {
-    rv = Check(index) ;
-  } else { 
-    TString className(obj->IsA()->GetName()) ; 
-    if (className.Contains("TObjArray")) {
-      rv = Check(index, static_cast<TObjArray *>(obj)) ;
-    } else if (className.Contains("TNtupleD")) { 
-      rv = Check(index, static_cast<TNtupleD *>(obj)) ;
-    }  else {
-      AliError(Form("%s class not implemented", className.Data())) ; 
-    }
-  }
+       Double_t * rv = new Double_t[AliRecoParam::kNSpecies] ;
+  Check(rv, index, recoParam) ;
+       SetQA(index, rv) ;      
+       
+  AliDebug(AliQAv1::GetQADebugLevel(), Form("Test result of %s", AliQAv1::GetAliTaskName(index))) ;
+       
+  delete [] rv ; 
+  Finish() ; 
+}
+
+//____________________________________________________________________________
+void AliQACheckerBase::Run(AliQAv1::ALITASK_t index, TObjArray ** list, AliDetectorRecoParam * recoParam) 
+{ 
+       AliDebug(AliQAv1::GetQADebugLevel(), Form("Processing %s", AliQAv1::GetAliTaskName(index))) ; 
   
+       Double_t * rv = new Double_t[AliRecoParam::kNSpecies] ;
+  Check(rv, index, list, recoParam) ;
        SetQA(index, rv) ;      
        
-  AliDebug(1, Form("Test result of %s", AliQA::GetAliTaskName(index))) ;
+  AliDebug(AliQAv1::GetQADebugLevel(), Form("Test result of %s", AliQAv1::GetAliTaskName(index))) ;
        
+  delete [] rv ; 
   Finish() ; 
 }
 
@@ -228,26 +350,130 @@ void AliQACheckerBase::Run(AliQA::ALITASK_t index, TObject * obj)
 void AliQACheckerBase::Finish() const 
 {
        // wrap up and save QA in proper file
-       AliQA * qa = AliQA::Instance() ; 
-       qa->Show() ;
-       AliQA::GetQAResultFile()->cd() ; 
-       qa->Write(qa->GetName(), kWriteDelete) ;   
-       AliQA::GetQAResultFile()->Close() ; 
+       AliQAv1::GetQAResultFile() ; 
+       AliQAv1 * qa = AliQAv1::Instance() ; 
+       qa->Write(AliQAv1::GetQAName(), kWriteDelete) ;   
+}
+
+//____________________________________________________________________________ 
+void AliQACheckerBase::MakeImage( TObjArray ** list, AliQAv1::TASKINDEX_t task, AliQAv1::MODE_t mode) 
+{
+  // makes the QA image for sim and rec
+  Int_t nImages = 0 ;
+  for (Int_t esIndex = 0 ; esIndex < AliRecoParam::kNSpecies ; esIndex++) {
+    if (! AliQAv1::Instance(AliQAv1::GetDetIndex(GetName()))->IsEventSpecieSet(AliRecoParam::ConvertIndex(esIndex)) ) 
+      continue ;
+    TIter next(list[esIndex]) ;  
+    TH1 * hdata = NULL ; 
+    while ( (hdata=static_cast<TH1 *>(next())) ) {
+      TString cln(hdata->ClassName()) ; 
+      if ( ! cln.Contains("TH") )
+        continue ; 
+      if ( hdata->TestBit(AliQAv1::GetImageBit()) )
+        nImages++; 
+    }
+    break ; 
+  }
+  if ( nImages == 0 ) {
+    AliDebug(AliQAv1::GetQADebugLevel(), Form("No histogram will be plotted for %s %s\n", GetName(), AliQAv1::GetTaskName(task).Data())) ;  
+  } else {
+    AliDebug(AliQAv1::GetQADebugLevel(), Form("%d histograms will be plotted for %s %s\n", nImages, GetName(), AliQAv1::GetTaskName(task).Data())) ;  
+    for (Int_t esIndex = 0 ; esIndex < AliRecoParam::kNSpecies ; esIndex++) {
+      if (! AliQAv1::Instance(AliQAv1::GetDetIndex(GetName()))->IsEventSpecieSet(AliRecoParam::ConvertIndex(esIndex)) || list[esIndex]->GetEntries() == 0) 
+        continue ;
+      const Char_t * title = Form("QA_%s_%s_%s", GetName(), AliQAv1::GetTaskName(task).Data(), AliRecoParam::GetEventSpecieName(esIndex)) ; 
+      if ( !fImage[esIndex] ) {
+        fImage[esIndex] = new TCanvas(title, title) ;
+      }
+      fImage[esIndex]->Clear() ; 
+      fImage[esIndex]->SetTitle(title) ; 
+      fImage[esIndex]->cd() ; 
+      TPaveText someText(0.015, 0.015, 0.98, 0.98) ;
+      someText.AddText(title) ;
+      someText.Draw() ; 
+      fImage[esIndex]->Print(Form("%s%s%d.%s", AliQAv1::GetImageFileName(), AliQAv1::GetModeName(mode), AliQAChecker::Instance()->GetRunNumber(), AliQAv1::GetImageFileFormat()), "ps") ; 
+      fImage[esIndex]->Clear() ; 
+      Int_t nx = TMath::Nint(TMath::Sqrt(nImages));
+      Int_t ny = nx  ; 
+      if (nx < TMath::Sqrt(nImages))
+        ny++ ; 
+      
+      fImage[esIndex]->Divide(nx, ny) ; 
+      TIter nexthist(list[esIndex]) ; 
+      TH1* hist = NULL ;
+      Int_t npad = 1 ; 
+      fImage[esIndex]->cd(npad) ; 
+      while ( (hist=static_cast<TH1*>(nexthist())) ) {
+        TString cln(hist->ClassName()) ; 
+        if ( ! cln.Contains("TH") )
+          continue ; 
+        if(hist->TestBit(AliQAv1::GetImageBit())) {
+          TString opts = hist->GetDrawOption();
+          if (opts.Contains("logy",TString::kIgnoreCase)) {
+            gPad->SetLogy();
+            opts.ReplaceAll("logy", "");
+          }
+          if (opts.Contains("logx", TString::kIgnoreCase)) {
+            gPad->SetLogx();
+            opts.ReplaceAll("logx", "");
+          }
+          hist->DrawCopy() ; 
+          fImage[esIndex]->cd(++npad) ; 
+        }
+      }
+      fImage[esIndex]->Print(Form("%s%s%d.%s", AliQAv1::GetImageFileName(), AliQAv1::GetModeName(mode), AliQAChecker::Instance()->GetRunNumber(), AliQAv1::GetImageFileFormat()), "ps") ; 
+    }
+  }  
 }
 
 //____________________________________________________________________________
-void AliQACheckerBase::SetQA(AliQA::ALITASK_t index, const Double_t value) const
+void AliQACheckerBase::SetHiLo(Float_t * hiValue, Float_t * lowValue) 
+{
+  AliDebug(AliQAv1::GetQADebugLevel(), "Previous setting was:") ;
+  if ( AliDebugLevel() == AliQAv1::GetQADebugLevel() ) {
+    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", 
+                              fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO], 
+                              fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING], 
+                              fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR], 
+                              fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ; 
+    AliInfo(Form("%s", text)) ; 
+  }
+  
+  for (Int_t index = 0 ; index < AliQAv1::kNBIT ; index++) {
+    fLowTestValue[index]  = lowValue[index] ; 
+    fUpTestValue[index]   = hiValue[index] ; 
+  }
+  AliDebug(AliQAv1::GetQADebugLevel(), "Current setting is:") ;
+  if ( AliDebugLevel()  == AliQAv1::GetQADebugLevel() ) {
+    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", 
+                              fLowTestValue[AliQAv1::kINFO], fUpTestValue[AliQAv1::kINFO], 
+                              fLowTestValue[AliQAv1::kWARNING], fUpTestValue[AliQAv1::kWARNING], 
+                              fLowTestValue[AliQAv1::kERROR], fUpTestValue[AliQAv1::kERROR], 
+                              fLowTestValue[AliQAv1::kFATAL], fUpTestValue[AliQAv1::kFATAL]) ;     AliInfo(Form("%s", text)) ; 
+  }
+}
+
+//____________________________________________________________________________
+void AliQACheckerBase::SetQA(AliQAv1::ALITASK_t index, Double_t * value) const
 {
        // sets the QA according the return value of the Check
 
-       AliQA * qa = AliQA::Instance(index) ; 
+  AliQAv1 * qa = AliQAv1::Instance(index) ;
 
-       if ( value <= 0.) 
-               qa->Set(AliQA::kFATAL) ; 
-       else if ( value > 0 && value <= 0.0002 )
-               qa->Set(AliQA::kERROR) ; 
-       else if ( value > 0.0002 && value <= 0.5 )
-               qa->Set(AliQA::kWARNING) ;
-       else if ( value > 0.5 && value < 1 ) 
-               qa->Set(AliQA::kINFO) ;                 
+  for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+    if (! qa->IsEventSpecieSet(AliRecoParam::ConvertIndex(specie)))
+      continue ;
+    if (  value == NULL ) { // No checker is implemented, set all QA to Fatal
+      qa->Set(AliQAv1::kFATAL, specie) ; 
+    } else {
+      if ( value[specie] >= fLowTestValue[AliQAv1::kFATAL] && value[specie] < fUpTestValue[AliQAv1::kFATAL] ) 
+        qa->Set(AliQAv1::kFATAL, AliRecoParam::ConvertIndex(specie)) ; 
+      else if ( value[specie] > fLowTestValue[AliQAv1::kERROR] && value[specie] <= fUpTestValue[AliQAv1::kERROR]  )
+        qa->Set(AliQAv1::kERROR, AliRecoParam::ConvertIndex(specie)) ; 
+      else if ( value[specie] > fLowTestValue[AliQAv1::kWARNING] && value[specie] <= fUpTestValue[AliQAv1::kWARNING]  )
+        qa->Set(AliQAv1::kWARNING, AliRecoParam::ConvertIndex(specie)) ;
+      else if ( value[specie] > fLowTestValue[AliQAv1::kINFO] && value[specie] <= fUpTestValue[AliQAv1::kINFO] ) 
+        qa->Set(AliQAv1::kINFO, AliRecoParam::ConvertIndex(specie)) ;  
+    }
+  }
 }