[u/mrichter/AliRoot.git] / EMCAL / AliEMCALQAChecker.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/*
  Checks the quality assurance. 
  By comparing with reference data

  Based on PHOS code written by
  Y. Schutz CERN July 2007

 For the moment we only implement the checking of raw data QA.
 The checked for ESD and RecPoints will be implemented later.
 

 */


// --- ROOT system ---
#include <TClass.h>
#include <TH1.h> 
#include <TF1.h> 
#include <TH1I.h> 
#include <TIterator.h> 
#include <TKey.h> 
#include <TFile.h> 
#include <TLine.h>
#include <TPaveText.h>
#include <TMath.h>

// --- Standard library ---

// --- AliRoot header files ---
#include "AliLog.h"
#include "AliQAv1.h"
#include "AliQAChecker.h"
#include "AliEMCALQAChecker.h"

ClassImp(AliEMCALQAChecker)

//__________________________________________________________________
AliEMCALQAChecker::AliEMCALQAChecker() : 
AliQACheckerBase("EMCAL","EMCAL Quality Assurance Data Maker"),
fLine(new TLine*[fknSM]),
fHref(new TLine*[fknSM]),
fText(NULL)
{
	/// ctor
  for (Int_t sm = 0 ; sm < fknSM ; sm++){
    fLine[sm] = NULL ; 
    fHref[sm] = NULL ; 
  }
}          

//__________________________________________________________________
AliEMCALQAChecker::~AliEMCALQAChecker() 
{
	/// dtor
  delete [] fLine ; 
  delete [] fHref ;
  if (fText) 
    delete fText ; 
}

//__________________________________________________________________
AliEMCALQAChecker::AliEMCALQAChecker(const AliEMCALQAChecker& qac) : 
AliQACheckerBase(qac.GetName(), qac.GetTitle()), 
fLine(new TLine*[fknSM]),
fHref(new TLine*[fknSM]),
fText(qac.fText)
{
	/// copy ctor 
  for (Int_t sm = 0 ; sm < fknSM ; sm++){
    fLine[sm] = qac.fLine[sm] ; 
    fHref[sm] = qac.fHref[sm] ; 
  }
}   
//__________________________________________________________________
AliEMCALQAChecker& AliEMCALQAChecker::operator = (const AliEMCALQAChecker &qac)
{
  fText = qac.fText;

  for (Int_t sm = 0 ; sm < fknSM ; sm++){
    fLine[sm] = qac.fLine[sm] ; 
    fHref[sm] = qac.fHref[sm] ; 
  }
  return *this ;
}

//______________________________________________________________________________
Double_t *
AliEMCALQAChecker::Check(AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/)
{
	/// Check objects in list
	
	if ( index == AliQAv1::kRAW ) 
	{
		return CheckRaws(list);
		printf ("checkers for task %d \n", index) ;		
	}
	
	if ( index == AliQAv1::kREC)
	{
		return CheckRecPoints(list);
	}
	
	if ( index == AliQAv1::kESD )
	{
		return CheckESD(list);
	}
	AliWarning(Form("Checker for task %d not implement for the moment",index));
	return NULL;
}

//______________________________________________________________________________
TH1* 
AliEMCALQAChecker::GetHisto(TObjArray* list, const char* hname, Int_t specie) const
{
	/// Get a given histo from the list
	TH1* h = static_cast<TH1*>(list->FindObject(Form("%s_%s",AliRecoParam::GetEventSpecieName(specie),hname)));
	if (!h)
	{
		AliError(Form("Did not find expected histo %s",hname));
	}
	return h;
}

//______________________________________________________________________________
Double_t 
AliEMCALQAChecker::MarkHisto(TH1& histo, Double_t value) const
{
	/// Mark histo as originator of some QA error/warning
	
	if ( value != 1.0 )
	{
		histo.SetBit(AliQAv1::GetQABit());
	}
	
	return value;
}


//______________________________________________________________________________
Double_t *
AliEMCALQAChecker::CheckRaws(TObjArray ** list)
{
//  Check RAW QA histograms	
//  We count the times of the response for each tower, the propability for all towers should be the same (average is given value).
//  We skip the first few cycles since the statistics is not enough, the average should be always larger than 1 at least.
//  By calculating the difference between the counts for each tower and the average, we decide whether we should recalculate 
//  the average depending the the gaus fitting on the divation distribution. 
//  During the recalutation of the average, we count how many towers are used for the calculation.
//  From the fraction of towers used, we decide whether each SM works fine or not
//  From the divation of average, we set the QA flag for the full detetcor as INFO, WARNING, ERROR or FATAL.
  
//  -- Yaxian Mao, CCNU/CERN/LPSC
					
  Float_t kThreshold = 80. ; 
  
  Double_t * test = new Double_t[AliRecoParam::kNSpecies] ;
  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 {
      TH1 * hdata = (TH1*)list[specie]->At(kTowerHG) ; 
      if(hdata->GetEntries()==0)
	continue;
      Int_t  nbin   = hdata->GetNbinsX() ;
      Int_t  nTower = nbin/4 ;  // ! number of channels for each SM
      if(fText) {
        fText->DeleteText() ; 
        fText->Clear() ;
      }
      else {
        fText = new TPaveText(0.3,0.6,0.7,0.9,"NDC") ;
      }
      fText->AddText(Form("OK if more than %2.2f %% inside aver-sigma < HG counts < aver+sigma", kThreshold));
      //TPaveText * fText[fknSM] = {0};
      Double_t rv    = 0. ;   //value to define the flag for the full detector
      for(Int_t iSM = 0 ; iSM < fknSM ; iSM++){  //number of SMs loop start
        TString projname = Form("proj_%d",iSM);
        Double_t aver  = 0. ;  //average calculated by the input histogram
        Double_t recal = 0 ;   // recalculate the average if the divation is large
        Double_t init  = 0. ;   // value to decide whether we should recalculate the average
        Double_t mean  = 0. ;  //the divation from the average, from the gaus fitting peak 
        Double_t width = 0. ;   // the sigma of the devation, from the gaus fitting
        Int_t flag     = 0 ;  //counts for channels used for recalculation of average
        Double_t ratio = 0. ;  //value to decide whether each SM works
        TH1F * proj    = NULL  ;  //a temp histogram store the difference from the average for the fitting procedure       
       
       for(Int_t iTower = iSM*nTower ; iTower<(iSM+1)*nTower ; iTower++)  //channel loop to calculate the average
          aver += hdata->GetBinContent(iTower);
        
        aver /=nTower;
        AliInfo(Form("SM: %d has average = %f\n",iSM, aver));
       Double_t ymin = hdata->GetBinContent(hdata->GetMinimumBin());
       Double_t ymax = hdata->GetBinContent(hdata->GetMaximumBin());
        proj = new TH1F(projname,projname,nbin,-aver,aver);
       fLine[iSM] = dynamic_cast<TLine*>(hdata->GetListOfFunctions()->FindObject(fLine[iSM]));
        //initialize the lines separate different SM
       if (!fLine[iSM]) {
          fLine[iSM] = new TLine((iSM+1)*nTower,ymin,(iSM+1)*nTower,ymax);
          fLine[iSM]->SetLineColor(1);
          fLine[iSM]->SetLineWidth(2);
          hdata->GetListOfFunctions()->Add(fLine[iSM]);  
          list[specie]->AddAt(hdata, kTowerHG) ; 
        }          
        else {
          fLine[iSM]->SetX1((iSM+1)*nTower) ; 
          fLine[iSM]->SetY1(ymin) ; 
          fLine[iSM]->SetX2((iSM+1)*nTower) ; 
          fLine[iSM]->SetY2(ymax) ; 
        }
        //filling the histogram with the difference and fitting with gaus function to obtain mean and sigma
       for(Int_t iTower = iSM*nTower ; iTower<(iSM+1)*nTower ; iTower++){
	 proj->Fill(hdata->GetBinContent(iTower)-aver);
       }	
       proj->Fit("gaus","","QNO");
       mean=proj->GetFunction("gaus")->GetParameter(1); // ! mean should be peaked at 0 in principal
       width=proj->GetFunction("gaus")->GetParameter(2);
       AliInfo(Form("aver = %f, mean = %f, sigma =%f\n",aver,mean, width));
        
       if(aver) init=TMath::Abs(mean/aver);  //calculate the divation from the average 
       
       //if mean or sigma is too huge or the sigma is too small, the fitting failed 
        if((aver+mean) < 1. || width/aver >2. || width < 1.e-3) 
          flag = -1 ;
       else {   //recalculate the average if the fitting didn't give the initial average  
	       aver+=mean;  //average corrected after fitting is fitting works
	     //if(aver <= 1. ) break ;
	    // if divation is too large, we conside to recalate the average by excluding channels too far from the average
         while(init>0.01 && aver >= 1.){
          if(flag) flag = 0 ;    //for each time recalculation, reset flag 
           for(Int_t iTower = iSM*nTower ; iTower < (iSM+1)*nTower ; iTower++){
             if(hdata->GetBinContent(iTower)>(aver-width) && hdata->GetBinContent(iTower)<(aver+width)){
               flag++ ;
	             recal += hdata->GetBinContent(iTower);
             }  
           }  //end of channels loop 
           
           if(flag == 0 || recal < 1.) {
             flag = -1 ;
             break ;
           }
           else {
             recal/=flag ;
	           init=(aver-recal)/(aver+recal) ; //difference between the new average and the pervious one
            aver =(aver+recal)/2 ; //recalculate the average by the last two values  
           }  
         } //out of while condition
          
         ratio=100.0*flag/nTower ;  //counting how many towers used for average recalculation
	       AliInfo(Form("SM %d has %2.2f %% chanhel works \n", iSM, ratio));
       }  // end of recalculation
        
        rv+=init ; //define the deviation from the final average
    
        //define the average line on each SM
       fHref[iSM] = dynamic_cast<TLine*>(hdata->GetListOfFunctions()->FindObject(fHref[iSM]));
       if(!fHref[iSM]) {
         fHref[iSM] = new TLine(iSM*nTower,aver,(iSM+1)*nTower,aver);
	       hdata->GetListOfFunctions()->Add(fHref[iSM]);  
	       list[specie]->AddAt(hdata, kTowerHG) ; 
       }          
       else {
         fHref[iSM]->Clear() ; 
	       fHref[iSM]->SetX1(iSM*nTower) ; 
         fHref[iSM]->SetY1(aver) ; 
         fHref[iSM]->SetX2((iSM+1)*nTower) ; 
         fHref[iSM]->SetY2(aver) ; 
       }
        
        hdata->Paint() ; 
 
        // if channels used for average recalculation smaller than the threshold,
        // then too much noise channels or channels didn't work 
        if(ratio<= kThreshold && flag >0){                             
          //fText->SetFillColor(2);
          fHref[iSM]->SetLineColor(2);
	        fHref[iSM]->SetLineWidth(2);
	        fText->SetFillColor(2);
	        fText->AddText(Form("SM %d: NOK = %2.0f %% channels OK!!!",iSM,ratio));
	        //fText[iSM]->AddText(Form("SM %d NOT OK, only %5.2f %% works!!!",iSM,flag/nTower));
        }         
        else if (flag == -1 || flag == 0 ) {  //fitting failed or recalculation didn't call
          fText->SetFillColor(2);
          fHref[iSM]->SetLineColor(2);
          fHref[iSM]->SetLineWidth(2);
          fText->SetFillColor(2);
          fText->AddText(Form("SM %d: NOK Fitting failed",iSM,ratio));
	        //fText[iSM]->AddText(Form("SM %d has %5.2f %% towers wok normally",iSM,flag/nTower));
        }  
        else {
          fText->SetFillColor(3);
	        fHref[iSM]->SetLineColor(3);
	        fHref[iSM]->SetLineWidth(2);
	        fText->AddText(Form("SM %d: OK = %2.0f %% channels OK",iSM,ratio));
	        //fText[iSM]->AddText(Form("SM %d has %5.2f %% towers wok normally",iSM,flag/nTower));
        }
        hdata->Paint() ; 
       //hdata->GetListOfFunctions()->Add(fText[iSM]);
       delete proj ; 
      }  // end of SM loop
      rv/=fknSM;
      hdata->GetListOfFunctions()->Add(fText);
      hdata->Paint() ;  
      if ( rv <=0.1 ) {
        AliInfo(Form("The deviation rv = %2.2f is small compared to average, detector works fine",rv));
        test[specie] =  0.05;
      }
      
      if ( rv <=0.5 && rv >0.1 )  {
        AliWarning(Form("The deviation rv = %2.2f is acceptable from average,  the detector works not perfect!",rv));
        test[specie] =  0.3;
      }
      
      if ( rv <=0.8 && rv >0.5 ) {
        AliError(Form("Got a large deviation of %2.2f from average, some error on the detector !!!",rv));
        test[specie] =  0.7;
      }
      
      if ( rv >0.8 ) {
        AliError(Form("Got too large deviation of %2.2f from average, detector out of control ???",rv));
        test[specie] =  0.9; 
      }
    
    } //end of checking raw
    
  }  //species loop 
  return test ;
}

//______________________________________________________________________________
void AliEMCALQAChecker::Init(const AliQAv1::DETECTORINDEX_t det) 
{
	/// intialises QA and QA checker settings
	AliQAv1::Instance(det) ; 
	Float_t hiValue[AliQAv1::kNBIT] ; 
	Float_t lowValue[AliQAv1::kNBIT] ;
	lowValue[AliQAv1::kINFO]      = 0.0   ; 
	hiValue[AliQAv1::kINFO]       = 0.1 ; 
	lowValue[AliQAv1::kWARNING]   = 0.1 ; 
  hiValue[AliQAv1::kWARNING]    = 0.5 ; 
	lowValue[AliQAv1::kERROR]     = 0.5   ; 
	hiValue[AliQAv1::kERROR]      = 0.8 ; 
	lowValue[AliQAv1::kFATAL]     = 0.8   ; 
	hiValue[AliQAv1::kFATAL]      = 1.0 ; 
	SetHiLo(&hiValue[0], &lowValue[0]) ; 
}
Commit	Line	Data
94594e5d	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	/*
	17	Checks the quality assurance.
	18	By comparing with reference data
	19
	20	Based on PHOS code written by
	21	Y. Schutz CERN July 2007
9e47432c	22
a2655076	23	For the moment we only implement the checking of raw data QA.
a2655076	24	The checked for ESD and RecPoints will be implemented later.
9e47432c	25
9e47432c	26
9e47432c	27	*/
9e47432c	28
94594e5d	29
	30	// --- ROOT system ---
	31	#include <TClass.h>
9e47432c	32	#include <TH1.h>
9e47432c	33	#include <TF1.h>
94594e5d	34	#include <TH1I.h>
	35	#include <TIterator.h>
	36	#include <TKey.h>
	37	#include <TFile.h>
9e47432c	38	#include <TLine.h>
	39	#include <TPaveText.h>
	40	#include <TMath.h>
94594e5d	41
	42	// --- Standard library ---
	43
	44	// --- AliRoot header files ---
	45	#include "AliLog.h"
4e25ac79	46	#include "AliQAv1.h"
94594e5d	47	#include "AliQAChecker.h"
	48	#include "AliEMCALQAChecker.h"
	49
	50	ClassImp(AliEMCALQAChecker)
	51
9e47432c	52	//__________________________________________________________________
	53	AliEMCALQAChecker::AliEMCALQAChecker() :
	54	AliQACheckerBase("EMCAL","EMCAL Quality Assurance Data Maker"),
	55	fLine(new TLine*[fknSM]),
	56	fHref(new TLine*[fknSM]),
	57	fText(NULL)
	58	{
	59	/// ctor
	60	for (Int_t sm = 0 ; sm < fknSM ; sm++){
	61	fLine[sm] = NULL ;
	62	fHref[sm] = NULL ;
	63	}
	64	}
	65
	66	//__________________________________________________________________
	67	AliEMCALQAChecker::~AliEMCALQAChecker()
	68	{
	69	/// dtor
	70	delete [] fLine ;
	71	delete [] fHref ;
	72	if (fText)
	73	delete fText ;
	74	}
	75
	76	//__________________________________________________________________
	77	AliEMCALQAChecker::AliEMCALQAChecker(const AliEMCALQAChecker& qac) :
	78	AliQACheckerBase(qac.GetName(), qac.GetTitle()),
	79	fLine(new TLine*[fknSM]),
	80	fHref(new TLine*[fknSM]),
	81	fText(qac.fText)
	82	{
	83	/// copy ctor
	84	for (Int_t sm = 0 ; sm < fknSM ; sm++){
	85	fLine[sm] = qac.fLine[sm] ;
	86	fHref[sm] = qac.fHref[sm] ;
	87	}
	88	}
	89	//__________________________________________________________________
	90	AliEMCALQAChecker& AliEMCALQAChecker::operator = (const AliEMCALQAChecker &qac)
	91	{
	92	fText = qac.fText;
	93
	94	for (Int_t sm = 0 ; sm < fknSM ; sm++){
	95	fLine[sm] = qac.fLine[sm] ;
	96	fHref[sm] = qac.fHref[sm] ;
	97	}
	98	return *this ;
	99	}
	100
	101	//______________________________________________________________________________
	102	Double_t *
486788fc	103	AliEMCALQAChecker::Check(AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /recoParam/)
9e47432c	104	{
	105	/// Check objects in list
	106
	107	if ( index == AliQAv1::kRAW )
	108	{
	109	return CheckRaws(list);
	110	printf ("checkers for task %d \n", index) ;
	111	}
	112
	113	if ( index == AliQAv1::kREC)
	114	{
	115	return CheckRecPoints(list);
	116	}
	117
	118	if ( index == AliQAv1::kESD )
	119	{
	120	return CheckESD(list);
	121	}
	122	AliWarning(Form("Checker for task %d not implement for the moment",index));
	123	return NULL;
	124	}
	125
	126	//______________________________________________________________________________
	127	TH1*
	128	AliEMCALQAChecker::GetHisto(TObjArray* list, const char* hname, Int_t specie) const
	129	{
	130	/// Get a given histo from the list
	131	TH1* h = static_cast<TH1*>(list->FindObject(Form("%s_%s",AliRecoParam::GetEventSpecieName(specie),hname)));
	132	if (!h)
	133	{
	134	AliError(Form("Did not find expected histo %s",hname));
	135	}
	136	return h;
	137	}
	138
	139	//______________________________________________________________________________
	140	Double_t
	141	AliEMCALQAChecker::MarkHisto(TH1& histo, Double_t value) const
	142	{
	143	/// Mark histo as originator of some QA error/warning
	144
	145	if ( value != 1.0 )
	146	{
	147	histo.SetBit(AliQAv1::GetQABit());
	148	}
	149
	150	return value;
	151	}
	152
	153
	154	//______________________________________________________________________________
	155	Double_t *
	156	AliEMCALQAChecker::CheckRaws(TObjArray ** list)
	157	{
a2655076	158	// Check RAW QA histograms
	159	// We count the times of the response for each tower, the propability for all towers should be the same (average is given value).
	160	// We skip the first few cycles since the statistics is not enough, the average should be always larger than 1 at least.
	161	// By calculating the difference between the counts for each tower and the average, we decide whether we should recalculate
	162	// the average depending the the gaus fitting on the divation distribution.
	163	// During the recalutation of the average, we count how many towers are used for the calculation.
	164	// From the fraction of towers used, we decide whether each SM works fine or not
	165	// From the divation of average, we set the QA flag for the full detetcor as INFO, WARNING, ERROR or FATAL.
6a754398	166
a2655076	167	// -- Yaxian Mao, CCNU/CERN/LPSC
a2655076	168
9e47432c	169	Float_t kThreshold = 80. ;
a2655076	170
	171	Double_t * test = new Double_t[AliRecoParam::kNSpecies] ;
	172	for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
	173	test[specie] = 1.0 ;
	174	if ( !AliQAv1::Instance()->IsEventSpecieSet(specie))
	175	continue ;
	176	if (list[specie]->GetEntries() == 0)
	177	test[specie] = 0. ; // nothing to check
	178	else {
	179	TH1 * hdata = (TH1*)list[specie]->At(kTowerHG) ;
	180	if(hdata->GetEntries()==0)
	181	continue;
9e47432c	182	Int_t nbin = hdata->GetNbinsX() ;
a2655076	183	Int_t nTower = nbin/4 ; // ! number of channels for each SM
9e47432c	184	if(fText) {
	185	fText->DeleteText() ;
	186	fText->Clear() ;
	187	}
	188	else {
	189	fText = new TPaveText(0.3,0.6,0.7,0.9,"NDC") ;
	190	}
	191	fText->AddText(Form("OK if more than %2.2f %% inside aver-sigma < HG counts < aver+sigma", kThreshold));
a2655076	192	//TPaveText * fText[fknSM] = {0};
	193	Double_t rv = 0. ; //value to define the flag for the full detector
	194	for(Int_t iSM = 0 ; iSM < fknSM ; iSM++){ //number of SMs loop start
9e47432c	195	TString projname = Form("proj_%d",iSM);
a2655076	196	Double_t aver = 0. ; //average calculated by the input histogram
	197	Double_t recal = 0 ; // recalculate the average if the divation is large
	198	Double_t init = 0. ; // value to decide whether we should recalculate the average
	199	Double_t mean = 0. ; //the divation from the average, from the gaus fitting peak
	200	Double_t width = 0. ; // the sigma of the devation, from the gaus fitting
	201	Int_t flag = 0 ; //counts for channels used for recalculation of average
	202	Double_t ratio = 0. ; //value to decide whether each SM works
	203	TH1F * proj = NULL ; //a temp histogram store the difference from the average for the fitting procedure
	204
	205	for(Int_t iTower = iSMnTower ; iTower<(iSM+1)nTower ; iTower++) //channel loop to calculate the average
	206	aver += hdata->GetBinContent(iTower);
9e47432c	207
a2655076	208	aver /=nTower;
	209	AliInfo(Form("SM: %d has average = %f\n",iSM, aver));
	210	Double_t ymin = hdata->GetBinContent(hdata->GetMinimumBin());
	211	Double_t ymax = hdata->GetBinContent(hdata->GetMaximumBin());
	212	proj = new TH1F(projname,projname,nbin,-aver,aver);
	213	fLine[iSM] = dynamic_cast<TLine*>(hdata->GetListOfFunctions()->FindObject(fLine[iSM]));
	214	//initialize the lines separate different SM
	215	if (!fLine[iSM]) {
9e47432c	216	fLine[iSM] = new TLine((iSM+1)nTower,ymin,(iSM+1)nTower,ymax);
	217	fLine[iSM]->SetLineColor(1);
	218	fLine[iSM]->SetLineWidth(2);
	219	hdata->GetListOfFunctions()->Add(fLine[iSM]);
	220	list[specie]->AddAt(hdata, kTowerHG) ;
	221	}
	222	else {
	223	fLine[iSM]->SetX1((iSM+1)*nTower) ;
	224	fLine[iSM]->SetY1(ymin) ;
	225	fLine[iSM]->SetX2((iSM+1)*nTower) ;
	226	fLine[iSM]->SetY2(ymax) ;
	227	}
a2655076	228	//filling the histogram with the difference and fitting with gaus function to obtain mean and sigma
	229	for(Int_t iTower = iSMnTower ; iTower<(iSM+1)nTower ; iTower++){
	230	proj->Fill(hdata->GetBinContent(iTower)-aver);
	231	}
	232	proj->Fit("gaus","","QNO");
	233	mean=proj->GetFunction("gaus")->GetParameter(1); // ! mean should be peaked at 0 in principal
	234	width=proj->GetFunction("gaus")->GetParameter(2);
	235	AliInfo(Form("aver = %f, mean = %f, sigma =%f\n",aver,mean, width));
6a754398	236
6a754398	237	if(aver) init=TMath::Abs(mean/aver); //calculate the divation from the average
a2655076	238
a2655076	239	//if mean or sigma is too huge or the sigma is too small, the fitting failed
6a754398	240	if((aver+mean) < 1. \|\| width/aver >2. \|\| width < 1.e-3)
6a754398	241	flag = -1 ;
a2655076	242	else { //recalculate the average if the fitting didn't give the initial average
6a754398	243	aver+=mean; //average corrected after fitting is fitting works
	244	//if(aver <= 1. ) break ;
	245	// if divation is too large, we conside to recalate the average by excluding channels too far from the average
	246	while(init>0.01 && aver >= 1.){
	247	if(flag) flag = 0 ; //for each time recalculation, reset flag
	248	for(Int_t iTower = iSMnTower ; iTower < (iSM+1)nTower ; iTower++){
	249	if(hdata->GetBinContent(iTower)>(aver-width) && hdata->GetBinContent(iTower)<(aver+width)){
	250	flag++ ;
	251	recal += hdata->GetBinContent(iTower);
	252	}
	253	} //end of channels loop
	254
	255	if(flag == 0 \|\| recal < 1.) {
	256	flag = -1 ;
	257	break ;
	258	}
	259	else {
	260	recal/=flag ;
	261	init=(aver-recal)/(aver+recal) ; //difference between the new average and the pervious one
	262	aver =(aver+recal)/2 ; //recalculate the average by the last two values
	263	}
	264	} //out of while condition
a2655076	265
6a754398	266	ratio=100.0*flag/nTower ; //counting how many towers used for average recalculation
6a754398	267	AliInfo(Form("SM %d has %2.2f %% chanhel works \n", iSM, ratio));
a2655076	268	} // end of recalculation
6a754398	269
	270	rv+=init ; //define the deviation from the final average
	271
	272	//define the average line on each SM
a2655076	273	fHref[iSM] = dynamic_cast<TLine*>(hdata->GetListOfFunctions()->FindObject(fHref[iSM]));
a2655076	274	if(!fHref[iSM]) {
6a754398	275	fHref[iSM] = new TLine(iSMnTower,aver,(iSM+1)nTower,aver);
	276	hdata->GetListOfFunctions()->Add(fHref[iSM]);
	277	list[specie]->AddAt(hdata, kTowerHG) ;
a2655076	278	}
a2655076	279	else {
6a754398	280	fHref[iSM]->Clear() ;
	281	fHref[iSM]->SetX1(iSM*nTower) ;
	282	fHref[iSM]->SetY1(aver) ;
	283	fHref[iSM]->SetX2((iSM+1)*nTower) ;
	284	fHref[iSM]->SetY2(aver) ;
a2655076	285	}
6a754398	286
	287	hdata->Paint() ;
	288
	289	// if channels used for average recalculation smaller than the threshold,
	290	// then too much noise channels or channels didn't work
	291	if(ratio<= kThreshold && flag >0){
	292	//fText->SetFillColor(2);
	293	fHref[iSM]->SetLineColor(2);
	294	fHref[iSM]->SetLineWidth(2);
	295	fText->SetFillColor(2);
	296	fText->AddText(Form("SM %d: NOK = %2.0f %% channels OK!!!",iSM,ratio));
	297	//fText[iSM]->AddText(Form("SM %d NOT OK, only %5.2f %% works!!!",iSM,flag/nTower));
	298	}
	299	else if (flag == -1 \|\| flag == 0 ) { //fitting failed or recalculation didn't call
	300	fText->SetFillColor(2);
	301	fHref[iSM]->SetLineColor(2);
	302	fHref[iSM]->SetLineWidth(2);
	303	fText->SetFillColor(2);
	304	fText->AddText(Form("SM %d: NOK Fitting failed",iSM,ratio));
	305	//fText[iSM]->AddText(Form("SM %d has %5.2f %% towers wok normally",iSM,flag/nTower));
	306	}
	307	else {
	308	fText->SetFillColor(3);
	309	fHref[iSM]->SetLineColor(3);
	310	fHref[iSM]->SetLineWidth(2);
	311	fText->AddText(Form("SM %d: OK = %2.0f %% channels OK",iSM,ratio));
	312	//fText[iSM]->AddText(Form("SM %d has %5.2f %% towers wok normally",iSM,flag/nTower));
	313	}
	314	hdata->Paint() ;
a2655076	315	//hdata->GetListOfFunctions()->Add(fText[iSM]);
	316	delete proj ;
	317	} // end of SM loop
	318	rv/=fknSM;
	319	hdata->GetListOfFunctions()->Add(fText);
	320	hdata->Paint() ;
6a754398	321	if ( rv <=0.1 ) {
	322	AliInfo(Form("The deviation rv = %2.2f is small compared to average, detector works fine",rv));
	323	test[specie] = 0.05;
	324	}
a2655076	325
6a754398	326	if ( rv <=0.5 && rv >0.1 ) {
	327	AliWarning(Form("The deviation rv = %2.2f is acceptable from average, the detector works not perfect!",rv));
	328	test[specie] = 0.3;
	329	}
a2655076	330
6a754398	331	if ( rv <=0.8 && rv >0.5 ) {
8976d7bd	332	AliError(Form("Got a large deviation of %2.2f from average, some error on the detector !!!",rv));
6a754398	333	test[specie] = 0.7;
	334	}
	335
	336	if ( rv >0.8 ) {
8976d7bd	337	AliError(Form("Got too large deviation of %2.2f from average, detector out of control ???",rv));
6a754398	338	test[specie] = 0.9;
	339	}
	340
a2655076	341	} //end of checking raw
	342
	343	} //species loop
	344	return test ;
9e47432c	345	}
	346
	347	//______________________________________________________________________________
	348	void AliEMCALQAChecker::Init(const AliQAv1::DETECTORINDEX_t det)
	349	{
	350	/// intialises QA and QA checker settings
	351	AliQAv1::Instance(det) ;
	352	Float_t hiValue[AliQAv1::kNBIT] ;
	353	Float_t lowValue[AliQAv1::kNBIT] ;
	354	lowValue[AliQAv1::kINFO] = 0.0 ;
	355	hiValue[AliQAv1::kINFO] = 0.1 ;
9e47432c	356	lowValue[AliQAv1::kWARNING] = 0.1 ;
6a754398	357	hiValue[AliQAv1::kWARNING] = 0.5 ;
9e47432c	358	lowValue[AliQAv1::kERROR] = 0.5 ;
	359	hiValue[AliQAv1::kERROR] = 0.8 ;
	360	lowValue[AliQAv1::kFATAL] = 0.8 ;
	361	hiValue[AliQAv1::kFATAL] = 1.0 ;
	362	SetHiLo(&hiValue[0], &lowValue[0]) ;
	363	}
94594e5d	364