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 **************************************************************************/
16 // Checks the quality assurance.
17 // By comparing with reference data
19 //---------------------------------------------
20 //checkig without reference data:
21 //for RAW QA all histograms should have approximatly the same
22 //number of entries as RefPoint
23 //for Rec Points checks
24 // - amplitude measured by 2 methos
25 // - online and offline T0 measurements
26 // for ESD quality of reconstruction ( and measurements):
27 // RMS of vertex and T0 less than 75ps
29 // Alla.Maevskaya@cern.ch
32 // --- ROOT system ---
33 #include <Riostream.h>
37 #include <TFitResultPtr.h>
39 #include <TIterator.h>
44 #include <TPaveText.h>
46 // --- Standard library ---
48 // --- AliRoot header files ---
51 #include "AliQAChecker.h"
52 #include "AliCDBEntry.h"
53 #include "AliQAManager.h"
54 #include "AliT0QAChecker.h"
56 ClassImp(AliT0QAChecker)
57 //____________________________________________________________________________
58 AliT0QAChecker::AliT0QAChecker() :
59 AliQACheckerBase("T0","T0 Quality Assurance Checker")
62 // Standard constructor
66 //____________________________________________________________________________
67 AliT0QAChecker::AliT0QAChecker(const AliT0QAChecker& qac):
68 AliQACheckerBase(qac.GetName(), qac.GetTitle())
71 AliError("Copy should not be used with this class\n");
73 //____________________________________________________________________________
74 AliT0QAChecker& AliT0QAChecker::operator=(const AliT0QAChecker& qac){
75 // assignment operator
76 this->~AliT0QAChecker();
77 new(this)AliT0QAChecker(qac);
82 //____________________________________________________________________________
83 AliT0QAChecker::~AliT0QAChecker(){
88 //__________________________________________________________________
89 void AliT0QAChecker::Check(Double_t * test, AliQAv1::ALITASK_t index, TObjArray ** list, const AliDetectorRecoParam * /*recoParam*/)
91 // Super-basic check on the QA histograms on the input list:
92 // look whether they are empty!
94 char * detOCDBDir = Form("T0/%s/%s", AliQAv1::GetRefOCDBDirName(), AliQAv1::GetRefDataDirName()) ;
96 AliCDBEntry *QARefRec = AliQAManager::QAManager()->Get(detOCDBDir);
99 AliInfo("QA reference data NOT retrieved for Reconstruction check. No T0 reference distribution");
103 for(Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++){
104 test[specie] = 1.0; //FK// initiate qa flag for the whole set of histograms as good
108 for(Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
109 if(!(AliQAv1::Instance()->IsEventSpecieSet(specie) && list[specie]) || list[specie]->GetEntries() == 0) {
113 if(index == AliQAv1::kRAW){
115 if(AliRecoParam::ConvertIndex(specie) == AliRecoParam::kCalib){// if (index == AliQAv1::kRAW )
116 //check laser data efficiencies
117 Double_t qaFlag = CheckLaser(list[specie]);
118 if(qaFlag < test[specie]) test[specie] = qaFlag;
121 if(AliRecoParam::ConvertIndex(specie) == AliRecoParam::kCalib ||
122 //AliRecoParam::ConvertIndex(specie) == AliRecoParam::kHighMult ||
123 AliRecoParam::ConvertIndex(specie) == AliRecoParam::kLowMult){
124 //AliRecoParam::ConvertIndex(specie) == AliRecoParam::kDefault ||
127 Double_t qaFlag = CheckBCID(list[specie]);
128 if(qaFlag < test[specie]) test[specie] = qaFlag;
131 if(//AliRecoParam::ConvertIndex(specie) == AliRecoParam::kHighMult ||
132 AliRecoParam::ConvertIndex(specie) == AliRecoParam::kLowMult){
133 //AliRecoParam::ConvertIndex(specie) == AliRecoParam::kDefault ||
135 Double_t qaFlag = CheckRaw(list[specie]);
136 if(qaFlag < test[specie]) test[specie] = qaFlag;
140 if(index == AliQAv1::kESD && AliRecoParam::Convert(specie) != AliRecoParam::kCalib){
141 test[specie] = CheckESD(list[specie]);
146 //--------------------------------------------------------------------------
147 Double_t AliT0QAChecker::CheckLaser(TObjArray *listrec) const {
153 //thresholds for warning and error on efficiencies
154 Float_t thrWarning = 0.5; //FK// warning level
155 Float_t thrError = 0.2; //FK// error level
157 const int kNumberOfHistos = 3;
158 Int_t consecutiveHistoNumber[kNumberOfHistos] = { 207, 208, 209}; //Checked histos fhCDFeff, hEffLED, hEffQTC
159 Int_t qualityFlag[kNumberOfHistos]; //quality flag for a given histogram
161 for(Int_t ir=0; ir<kNumberOfHistos; ir++){
162 qualityFlag[ir] = kT0Info; //init quality flag for a given histogram
164 hdata = (TH1*) listrec->UncheckedAt(consecutiveHistoNumber[ir]);
166 fhRawEff[nEffHistos] = hdata;
171 TLine linelowyellow(0, thrWarning, 24, thrWarning);
172 linelowyellow.SetLineColor(5);
173 linelowyellow.SetLineStyle(3);
174 linelowyellow.SetLineWidth(4);
175 TLine linelowred(0, thrError, 24, thrError);
176 linelowred.SetLineColor(2);
177 linelowred.SetLineStyle(3);
178 linelowred.SetLineWidth(4);
180 Bool_t bEffHistosNotEmpty = kFALSE; //check if all histograms have some counts
182 for(Int_t ih = 0; ih < nEffHistos; ih++){
184 EraseOldMessages((TH1*) fhRawEff[ih]);// clean objects added at previous checks
186 fhRawEff[ih]->SetLineWidth(2);
187 fhRawEff[ih]->SetMaximum(2.);
188 fhRawEff[ih]->SetMinimum(0.);
189 fhRawEff[ih]->GetListOfFunctions()->Add((TLine*)linelowyellow.Clone());
190 fhRawEff[ih]->GetListOfFunctions()->Add((TLine*)linelowred.Clone());
192 if(fhRawEff[ih]->Integral()>0) bEffHistosNotEmpty = kTRUE; //this histo does have some counts in it
194 Int_t nbins= fhRawEff[ih]->GetNbinsX();
195 for(Int_t ib=1; ib<=nbins; ib++){ //loop over bins and check if the efficiency is above level
197 Float_t chcont = fhRawEff[ih]->GetBinContent(ib);
198 if(chcont < thrWarning && qualityFlag[ih] > kT0Error ) qualityFlag[ih] = kT0Warning;//Warning level
199 if(chcont < thrError) qualityFlag[ih] = kT0Error;//Error level
202 if(qualityFlag[ih] == kT0Info ){
203 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 efficiency %s is good", fhRawEff[ih]->GetName() ));
204 }else if(qualityFlag[ih] == kT0Warning){
205 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 efficiency %s is not so good", fhRawEff[ih]->GetName() ));
206 }else if(qualityFlag[ih] == kT0Error){
207 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 efficiency %s is not good", fhRawEff[ih]->GetName() ));
212 int lowestQualityFlag = (int) kT0Info;
213 for(Int_t ih = 0; ih < nEffHistos; ih++){
215 if(!bEffHistosNotEmpty){ //all laser efficiency plots are empty
216 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
217 text.AddText(Form("1) T0 is in BEAMTUNIG: empty plots are ok"));
218 text.AddText(Form("2) T0 is in READY: check calibriation trigger"));
219 text.AddText(Form("if also physics data are empty report"));
220 text.AddText(Form("readout problem to the T0 on-call expert"));
221 fhRawEff[ih]->GetListOfFunctions()->Add((TPaveText*)text.Clone());
224 if( qualityFlag[ih] <lowestQualityFlag ) lowestQualityFlag = qualityFlag[ih];
227 return ConvertQualityFlagToDouble(lowestQualityFlag);
229 //--------------------------------------------------------------------------
230 Double_t AliT0QAChecker::CheckBCID(TObjArray *listrec) const {
232 Int_t qualityFlagBCID = kT0Info; //init quality flag for a given histogram;
234 TH2F *hBCID = (TH2F*) listrec->UncheckedAt(224); //BCID versus TRM BCID
237 // clean objects added at previous checks
238 EraseOldMessages((TH1*)hBCID);
240 if(hBCID->Integral()>0){
241 //BCID does have some counts in it
243 Int_t nbinsX = hBCID->GetNbinsX();
244 Int_t nbinsY = hBCID->GetNbinsY();
245 double entriesOnDiagonal = 0; //count diagonal and off diagonal entries
246 double entriesOffDiagonal = 0;
248 for(Int_t ix=1; ix<=nbinsX; ix++){
249 for(Int_t iy=1; iy<=nbinsY; iy++){
250 if(TMath::Abs(ix-iy)<6) entriesOnDiagonal += hBCID->GetBinContent(ix,iy); //On Diagonal
251 else entriesOffDiagonal += hBCID->GetBinContent(ix,iy); //Off Diagonal
254 if(entriesOnDiagonal<1 || entriesOffDiagonal>0){
255 qualityFlagBCID = kT0Error; //no entries on diagonal
256 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 %s is not diagonal", hBCID->GetName() ));
258 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
259 text.AddText(Form("Check if entries are on a diagonal."));
260 text.AddText(Form("Report readout problem to the T0 on-call expert"));
261 hBCID->GetListOfFunctions()->Add((TPaveText*)text.Clone());
265 qualityFlagBCID = kT0Error;
266 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 : %s has NO entries", hBCID->GetName() ));
268 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
269 text.AddText(Form("NO ENTRIES!!!"));
270 text.AddText(Form("If T0 is READY report"));
271 text.AddText(Form("readout problem to the T0 on-call expert"));
272 hBCID->GetListOfFunctions()->Add((TPaveText*)text.Clone());
276 int lowestQualityFlag = (int) qualityFlagBCID;
278 return ConvertQualityFlagToDouble(lowestQualityFlag);
282 //--------------------------------------------------------------------------
283 Double_t AliT0QAChecker::CheckRaw(TObjArray *listrec) const {
286 Int_t qualityFlagTrigger = kT0Info; //init quality flag for a given histogram;
288 TH1F *hTrigger = (TH1F*) listrec->UncheckedAt(169);//hRawTrigger
291 // clean objects added at previous checks
292 EraseOldMessages((TH1*) hTrigger);
294 if(hTrigger->Integral()>0){
295 //trigger plot does have some counts in it
296 //are Mean, ORA and ORC not empty?
297 if( hTrigger->GetBinContent(1)<0.001 || hTrigger->GetBinContent(3)<0.001 || hTrigger->GetBinContent(4)<0.001){
298 qualityFlagTrigger = kT0Error; //no entries on diagonal
299 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0: too little ORA and ORC in %s", hTrigger->GetName() ));
301 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
302 text.AddText(Form("Check ORA and ORC"));
303 text.AddText(Form("Report problem to the T0 on-call expert"));
304 hTrigger->GetListOfFunctions()->Add((TPaveText*)text.Clone());
307 }else{ //Trigger histo empty
309 qualityFlagTrigger = kT0Error;
310 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 histogram %s has NO entries", hTrigger->GetName() ));
312 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
313 text.AddText(Form("NO ENTRIES!!!"));
314 text.AddText(Form("If T0 is READY report"));
315 text.AddText(Form("readout problem to the T0 on-call expert"));
316 hTrigger->GetListOfFunctions()->Add((TPaveText*)text.Clone());
320 //--------- timing plots ---------
321 Int_t qualityFlagMeanBest = kT0Info; //init quality flag for a given histogram;
322 TH1F *hMeanBest = (TH1F*) listrec->UncheckedAt(223); //fhMeanBest time
323 EraseOldMessages((TH1*) hMeanBest);
325 if(hMeanBest->Integral()<1){
326 qualityFlagMeanBest = kT0Error;
327 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 histogram %s has NO entries", hMeanBest->GetName() ));
329 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
330 text.AddText(Form("NO ENTRIES!!!"));
331 text.AddText(Form("If T0 is READY and beam is on report"));
332 text.AddText(Form("the problem to the T0 on-call expert"));
333 hMeanBest->GetListOfFunctions()->Add((TPaveText*)text.Clone());
339 //----------- vetext plots ---------
341 Int_t qualityFlagVetrexFirst = kT0Info; //init quality flag for a given histogram;
342 Int_t qualityFlagVetrexBest = kT0Info; //init quality flag for a given histogram;
344 TH1F *hVertexFirst = (TH1F*) listrec->UncheckedAt(225);//fhVertex1st
345 TH1F *hVertexBest = (TH1F*) listrec->UncheckedAt(226);//fhVertexBest
346 TH1F *hOrCminOrATvdcOn = (TH1F*) listrec->UncheckedAt(217);//vertex with TVDC on
348 // clean objects added at previous checks
349 EraseOldMessages((TH1*) hVertexFirst);
350 EraseOldMessages((TH1*) hVertexBest);
351 float errorLevelDifference = 200;
352 float warningLevelDifference = 100;
354 if(hVertexFirst->Integral()>0 && hOrCminOrATvdcOn->Integral()>0){
355 float meanVertexFirst = hVertexFirst->GetMean();
356 float meanVertexTVDC = hOrCminOrATvdcOn->GetMean();
357 float diff = TMath::Abs(meanVertexFirst - meanVertexTVDC);
358 if(diff > errorLevelDifference){
359 qualityFlagVetrexFirst = kT0Error; //init quality flag for a given histogram;
361 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
362 text.AddText(Form("large diff. between TVDC vertex and"));
363 text.AddText(Form("first vertex. Alert the T0 on-call expert"));
364 hVertexFirst->GetListOfFunctions()->Add((TPaveText*)text.Clone());
366 }else if(diff > warningLevelDifference){
367 qualityFlagVetrexFirst = kT0Warning; //init quality flag for a given histogram;
369 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
370 text.AddText(Form("large diff. between TVDC vertex and"));
371 text.AddText(Form("first vertex. Inform the T0 on-call expert"));
372 hVertexFirst->GetListOfFunctions()->Add((TPaveText*)text.Clone());
376 qualityFlagVetrexFirst = kT0Error;
377 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 histogram %s has NO entries", hVertexFirst->GetName() ));
379 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
380 text.AddText(Form("NO ENTRIES!!!"));
381 text.AddText(Form("If T0 is READY and beam is on report"));
382 text.AddText(Form("the problem to the T0 on-call expert"));
383 hVertexFirst->GetListOfFunctions()->Add((TPaveText*)text.Clone());
386 if(hVertexBest->Integral()>0 && hOrCminOrATvdcOn->Integral()>0){
387 float meanVertexBest = hVertexBest->GetMean();
388 float meanVertexTVDC = hOrCminOrATvdcOn->GetMean();
389 float diff = TMath::Abs(meanVertexBest - meanVertexTVDC);
390 if(diff > errorLevelDifference){
391 qualityFlagVetrexBest = kT0Error; //init quality flag for a given histogram;
393 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
394 text.AddText(Form("Large. diff. between TVDC vertex and"));
395 text.AddText(Form("best vertex. Alert the T0 on-call expert"));
396 hVertexBest->GetListOfFunctions()->Add((TPaveText*)text.Clone());
398 }else if(diff > warningLevelDifference){
399 qualityFlagVetrexBest = kT0Warning; //init quality flag for a given histogram;
401 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
402 text.AddText(Form("Large. diff. between TVDC vertex and"));
403 text.AddText(Form("best vertex. Infrom the T0 on-call expert"));
404 hVertexBest->GetListOfFunctions()->Add((TPaveText*)text.Clone());
407 qualityFlagVetrexBest = kT0Error;
408 AliDebug(AliQAv1::GetQADebugLevel(), Form("T0 histogram %s has NO entries", hVertexBest->GetName() ));
410 TPaveText text(0.20,0.50,0.99,0.99,"NDC");
411 text.AddText(Form("NO ENTRIES!!!"));
412 text.AddText(Form("If T0 is READY and beam is on report"));
413 text.AddText(Form("the problem to the T0 on-call expert"));
414 hVertexBest->GetListOfFunctions()->Add((TPaveText*)text.Clone());
419 int lowestQualityFlag = (int) qualityFlagTrigger;
420 if(qualityFlagVetrexBest<lowestQualityFlag) lowestQualityFlag = qualityFlagVetrexBest;
421 if(qualityFlagVetrexFirst<lowestQualityFlag) lowestQualityFlag = qualityFlagVetrexFirst;
422 if(qualityFlagMeanBest <lowestQualityFlag) lowestQualityFlag = qualityFlagMeanBest;
425 return ConvertQualityFlagToDouble(lowestQualityFlag);
429 //--------------------------------------------------------------------------
430 Double_t AliT0QAChecker::CheckESD(TObjArray *listrec ) const
435 fhESD = (TH1*) listrec->UncheckedAt(2);
437 AliDebug(AliQAv1::GetQADebugLevel(), Form("count %s ", fhESD->GetName()) );
438 TF1 *f1 = new TF1("f1","gaus",-1,1);
439 fhESD->Fit("f1","R","Q", -1,1);
441 f1->GetParameters(&par[0]);
443 TPaveText text(0.30,0.50,0.99,0.99,"NDC");
445 text.AddText(Form("T0 RUN %d ",AliCDBManager::Instance()->GetRun()));
447 AliDebug(AliQAv1::GetQADebugLevel(), Form("numentries %d mean %f #sigma %f", (int)fhESD->GetEntries(),par[1], par[2]));
450 if (par[2] > 0.07 && par[2] < 1.) {
452 text.AddText(Form("not good resolution :\n %f ns\n", par[2] ));
453 text.SetFillColor(5);
454 printf("T0 detector resolution is not good enouph: %f ns\n",par[2] );
456 if(TMath::Abs(par[1])>0.05) {
458 text.AddText(Form(" Check clock shift on %f ns", par[1]));
459 text.SetFillColor(5);
461 if (par[2] > 1. || TMath::Abs(par[1])>0.1) {
463 text.AddText(Form(" Bad resolution:\n mean %f ns sigma %f ns", par[1], par[2]));
464 text.SetFillColor(2);
465 { // RS Clean previous additions
466 TList* lstF = fhESD->GetListOfFunctions();
468 TObject *stats = lstF->FindObject("stats");
471 while ((obj = lstF->First())) {
472 while(lstF->Remove(obj)) { }
475 if (stats) lstF->Add(stats);
478 fhESD->GetListOfFunctions()->Add((TPaveText*)text.Clone());
479 AliDebug(AliQAv1::GetQADebugLevel(),
480 Form("Please, check calibration: shift= %f resolution %f test=%f\n",
481 par[1], par[2], checkr) ) ;
486 AliDebug(AliQAv1::GetQADebugLevel(),
487 Form("No ESD QA histogram found, nothing to check"));
496 //--------------------------------------------------------------------------
497 void AliT0QAChecker::EraseOldMessages(TH1* h) const
499 //erase the old captions
500 TList* lstF = h->GetListOfFunctions();
502 TObject *stats = lstF->FindObject("stats");
505 while ((obj = lstF->First())) {
506 while(lstF->Remove(obj)) { }
509 if (stats) lstF->Add(stats);
512 //--------------------------------------------------------------------------
513 Double_t AliT0QAChecker::ConvertQualityFlagToDouble(int qualityFlag) const
515 //covert quality flag to double
518 switch ( qualityFlag ){
522 checkr = 0.75; break;
524 checkr = 0.25; break;
526 checkr = -1.0; break;
528 AliError("Invalid ecc value. FIXME !");
529 checkr = 0.25; break;