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 ///////////////////////////////////////////////////////////////////////////////
18 // Class to process a tree and create alarms based on thresholds //
19 // origin: jens wiechula: jens.wiechula@cern.ch //
21 ///////////////////////////////////////////////////////////////////////////////
24 #include <TObjArray.h>
26 #include <TObjString.h>
36 #include <TIterator.h>
41 #include <TStopwatch.h>
45 #include "AliTPCCalibQAChecker.h"
49 AliTPCCalibQAChecker::AliTPCCalibQAChecker() :
50 TNamed("AliTPCCalibQAChecker","AliTPCCalibQAChecker"),
56 fIterSubCheckers(0x0),
58 fArrAlarmDescriptions(0x0),
71 ResetAlarmThresholds();
73 //_________________________________________________________________________
74 AliTPCCalibQAChecker::AliTPCCalibQAChecker(const char* name, const char *title) :
81 fIterSubCheckers(0x0),
83 fArrAlarmDescriptions(0x0),
96 ResetAlarmThresholds();
98 //_________________________________________________________________________
99 AliTPCCalibQAChecker::~AliTPCCalibQAChecker()
104 if (fHistRep) delete fHistRep;
105 if (fIterSubCheckers) delete fIterSubCheckers;
106 if (fArrAlarmDescriptions) delete fArrAlarmDescriptions;
108 //_________________________________________________________________________
109 void AliTPCCalibQAChecker::AddSubChecker(AliTPCCalibQAChecker *alarm)
112 // add a sub checker to this checker
115 if (!fArrSubCheckers) {
116 fArrSubCheckers=new TObjArray;
117 fArrSubCheckers->SetOwner();
119 fArrSubCheckers->Add(alarm);
121 //_________________________________________________________________________
122 void AliTPCCalibQAChecker::Process()
125 // Process the alarm thresholds, decide the alarm level, create the representation histogram
128 //reset quality level
133 //decide which type of checker to use
134 if (fArrSubCheckers && fArrSubCheckers->GetEntries()>0) ProcessSub();
135 else if (fTreePtr && *fTreePtr) ProcessTree();
136 else if (fHistPtr && *fHistPtr) ProcessHist();
137 else if (fGraphPtr && *fGraphPtr) ProcessGraph();
138 else if (fNumberPtr) ProcessNumber();
140 AliInfo(Form("Processing Time (%s): %fs",GetName(),s.RealTime()));
142 //_________________________________________________________________________
143 void AliTPCCalibQAChecker::ProcessSub()
146 // sub checker type checker
148 QualityFlag_t quality=kINFO;
149 if (fArrSubCheckers && fArrSubCheckers->GetEntries()>0){
150 TIter next(fArrSubCheckers);
152 while ( (o=next()) ) {
153 AliTPCCalibQAChecker *al=(AliTPCCalibQAChecker*)o;
155 QualityFlag_t subQuality=al->GetQuality();
156 if (subQuality>quality) quality=subQuality;
159 fQualityLevel=quality;
161 //_________________________________________________________________________
162 void AliTPCCalibQAChecker::ProcessTree()
165 // process tree type checker
168 //Create Representation Histogram
169 CreateRepresentationHist();
171 // if (!fTree) return;
172 //chek for the quality
188 //_________________________________________________________________________
189 void AliTPCCalibQAChecker::ProcessHist()
192 // process histogram type checker
195 if (!(fHistPtr && *fHistPtr)) return;
208 //create representation histogram if we are not in tree mode
211 fHistRep=(*fHistPtr)->Clone();
212 TH1* h=(TH1*)fHistRep;
214 h->SetNameTitle(Form("h%s",GetName()),GetTitle());
216 TH1* h=(TH1*)fHistRep;
219 if (!fHistRep->InheritsFrom(TH2::Class())) AddQualityLines(h);
222 //_________________________________________________________________________
223 void AliTPCCalibQAChecker::ProcessGraph()
226 // process graph type checker
228 if (!(fGraphPtr && *fGraphPtr)) return;
229 Int_t npoints=(*fGraphPtr)->GetN();
230 fQualityLevel=GetQuality(npoints,(*fGraphPtr)->GetY());
232 //_________________________________________________________________________
233 void AliTPCCalibQAChecker::ProcessNumber()
236 // process number type checker
238 if (!fNumberPtr) return;
239 fQualityLevel=GetQuality(*fNumberPtr);
240 //create representation histogram
242 TH1* h=new TH1F(Form("h%s",GetName()),GetTitle(),1,0,1);
243 h->GetXaxis()->SetBinLabel(1,"Value");
248 TH1 *h=(TH1*)fHistRep;
249 TMarker *marker=(TMarker*)h->GetListOfFunctions()->FindObject("TMarker");
251 marker=new TMarker(.5,0,20);
252 h->GetListOfFunctions()->Add(marker);
254 marker->SetMarkerColor(GetQualityColor());
255 marker->SetY(*fNumberPtr);
257 //_________________________________________________________________________
258 void AliTPCCalibQAChecker::ProcessEntries()
261 // Processing function which analyses the number of affected rows of a tree draw
263 TString draw=fStrDraw;
264 if (draw.IsNull()) return;
266 TString cuts=fStrCuts;
268 TString opt=fStrDrawOpt;
271 //draw and get the histogram
272 Int_t res=(*fTreePtr)->Draw(draw.Data(),cuts.Data(),opt.Data());
273 fQualityLevel=GetQuality(res);
275 //_________________________________________________________________________
276 void AliTPCCalibQAChecker::ProcessMean()
279 // Processing function which analyses the mean of the resulting histogram
283 Double_t value=h->GetMean();
284 if (fAlarmType==kNentries) value=h->GetEntries();
285 fQualityLevel=GetQuality(value);
287 //_________________________________________________________________________
288 void AliTPCCalibQAChecker::ProcessBin()
291 // Process a histogram bin by bin and check for thresholds
294 //bin quality counters
295 Int_t nquality[kNQualityFlags];
296 for (Int_t iquality=(Int_t)kINFO; iquality<kNQualityFlags; ++iquality) nquality[iquality]=0;
300 Int_t nbinsX=h->GetNbinsX();
301 Int_t nbinsY=h->GetNbinsY();
302 Int_t nbinsZ=h->GetNbinsZ();
303 Int_t nbinsTotal=nbinsX*nbinsY*nbinsZ;
306 for (Int_t ibinZ=1;ibinZ<nbinsZ+1;++ibinZ){
307 for (Int_t ibinY=1;ibinY<nbinsY+1;++ibinY){
308 for (Int_t ibinX=1;ibinX<nbinsX+1;++ibinX){
309 Double_t value = (*fHistPtr)->GetBinContent(ibinX, ibinY, ibinZ);
310 QualityFlag_t quality=GetQuality(value);
316 //loop over Quality levels and set quality
317 for (Int_t iquality=(Int_t)kINFO; iquality<kNQualityFlags; ++iquality){
318 if (fAlarmType==kBinAny){
319 if (nquality[iquality]) fQualityLevel=(QualityFlag_t)iquality;
320 } else if (fAlarmType==kBinAll){
321 if (nquality[iquality]==nbinsTotal) fQualityLevel=(QualityFlag_t)iquality;
325 //_________________________________________________________________________
326 void AliTPCCalibQAChecker::CreateRepresentationHist()
329 // Create the representation histogram which will be shown in the draw function
331 ResetRepresentationHist();
333 TString draw=fStrDrawRep;
336 fStrDrawRepOpt=fStrDrawOpt;
338 draw.ReplaceAll("%alarm%",fStrDraw.Data());
340 if (draw.IsNull()) return;
342 TString cuts=fStrCuts;
344 TString opt=fStrDrawRepOpt;
347 Int_t res=(*fTreePtr)->Draw(draw.Data(),cuts.Data(),opt.Data());
348 TH1 *hist=(*fTreePtr)->GetHistogram();
350 AliError(Form("Could not create representation histogram of alarm '%s'",GetName()));
353 fHistRep=hist->Clone();
354 TH1 *h=(TH1*)fHistRep;
356 h->SetNameTitle(Form("h%s",GetName()),GetTitle());
358 //_________________________________________________________________________
359 void AliTPCCalibQAChecker::CreateAlarmHist()
362 // create alarm histogram from the tree
365 TString draw=fStrDraw;
366 if (draw.IsNull()) return;
368 TString cuts=fStrCuts;
370 TString opt=fStrDrawOpt;
373 //draw and get the histogram
374 Int_t res=(*fTreePtr)->Draw(draw.Data(),cuts.Data(),opt.Data());
375 fHist=(*fTreePtr)->GetHistogram();
376 if (res<0 || !fHist){
377 AliError(Form("Could not create alarm histogram of alarm '%s'",GetName()));
380 fHist->SetDirectory(0);
383 //_________________________________________________________________________
384 void AliTPCCalibQAChecker::ResetAlarmHist()
387 // delete the alarm histogram and reset the pointer
390 if (*fHistPtr) delete *fHistPtr;
394 //_________________________________________________________________________
395 void AliTPCCalibQAChecker::Draw(Option_t *option)
398 // object draw function
399 // by default the pad backgound color is set to the quality level color
400 // use 'nobc' to change this
403 //case of a node with subcheckers and no specific representation histogram
404 if (HasSubCheckers()&&!fHistRep) {DrawSubNodes(option); return;}
405 if (fHistRep) {DrawRepresentationHist(option); return;}
407 //_________________________________________________________________________
408 void AliTPCCalibQAChecker::Print(Option_t *option) const
411 // print the quality status. If we have sub checkers print recursively
413 TString sOpt(option);
414 cout << sOpt << GetName() << ": " << GetQualityName() << endl;
415 if (fArrSubCheckers && fArrSubCheckers->GetEntries()>0){
416 sOpt.ReplaceAll("+-"," ");
418 TIter next(fArrSubCheckers);
420 while ( (o=next()) ) o->Print(sOpt.Data());
423 //_________________________________________________________________________
424 void AliTPCCalibQAChecker::SetAlarmThreshold(const Double_t min, const Double_t max, const QualityFlag_t quality)
427 //set the alarm thresholds for a specific quality level
429 if ((Int_t)quality<(Int_t)kINFO||(Int_t)quality>=kNQualityFlags) return;
430 fThresMin[quality]=min;
431 fThresMax[quality]=max;
433 //_________________________________________________________________________
434 void AliTPCCalibQAChecker::ResetAlarmThreshold(const QualityFlag_t quality)
437 //set the alarm thresholds for a specific quality level
439 if ((Int_t)quality<(Int_t)kINFO||(Int_t)quality>=kNQualityFlags) return;
440 fThresMin[quality]=0;
441 fThresMax[quality]=0;
443 //_________________________________________________________________________
444 void AliTPCCalibQAChecker::ResetAlarmThresholds()
447 //reset all the alarm thresholds
449 for (Int_t i=0;i<kNQualityFlags;++i){
454 //_________________________________________________________________________
455 void AliTPCCalibQAChecker::SetQualityDescription(const char* text, const QualityFlag_t quality)
458 // set an description for the quality level
459 // %min and %max will be replaced by the min and max values of the alarm, when the quality
460 // description is queried (see QualityDescription)
463 if (quality<kINFO||quality>kFATAL) return;
464 if (! fArrAlarmDescriptions ) fArrAlarmDescriptions=new TObjArray(kNQualityFlags);
465 TObjString *s=(TObjString*)fArrAlarmDescriptions->At(quality);
466 if (!s) fArrAlarmDescriptions->AddAt(s=new TObjString,quality);
470 //_________________________________________________________________________
471 const AliTPCCalibQAChecker* AliTPCCalibQAChecker::GetSubChecker(const char* name, Bool_t recursive) const
477 if (sname==GetName()) return this;
478 if (!fArrSubCheckers || !fArrSubCheckers->GetEntries()) return 0x0;
479 const AliTPCCalibQAChecker *al=0x0;
481 TIter next(fArrSubCheckers);
483 while ( (o=next()) ){
484 AliTPCCalibQAChecker *sal=(AliTPCCalibQAChecker*)o;
485 al=sal->GetSubChecker(name);
489 al=dynamic_cast<AliTPCCalibQAChecker*>(fArrSubCheckers->FindObject(name));
493 //_________________________________________________________________________
494 Int_t AliTPCCalibQAChecker::GetNumberOfSubCheckers(Bool_t recursive) const
497 // get the number of sub checkers
498 // if recursive get total number of non subchecker type sub checkers
502 if (!fArrSubCheckers) return 1;
503 if (!fArrSubCheckers->GetEntries()) return 0;
504 TIter next(fArrSubCheckers);
506 while ( (o=next()) ){
507 AliTPCCalibQAChecker *al=(AliTPCCalibQAChecker*)o;
508 nsub+=al->GetNumberOfSubCheckers();
511 if (fArrSubCheckers) nsub=fArrSubCheckers->GetEntries();
515 //_________________________________________________________________________
516 AliTPCCalibQAChecker* AliTPCCalibQAChecker::NextSubChecker()
519 // loop over sub checkers
520 // if recursive, recursively return the pointers of non subchecker type sub checkers
522 if (!fArrSubCheckers && !fArrSubCheckers->GetEntries()) return 0;
523 if (!fIterSubCheckers) fIterSubCheckers=fArrSubCheckers->MakeIterator();
524 AliTPCCalibQAChecker *al=(AliTPCCalibQAChecker*)fIterSubCheckers->Next();
526 delete fIterSubCheckers;
527 fIterSubCheckers=0x0;
529 // if (recursive && al->GetNumberOfSubCheckers(kFALSE)) al=al->NextSubChecker();
532 //_________________________________________________________________________
533 const char* AliTPCCalibQAChecker::QualityName(const AliTPCCalibQAChecker::QualityFlag_t quality)
536 // get quality name for quality
555 //_________________________________________________________________________
556 Color_t AliTPCCalibQAChecker::QualityColor(const AliTPCCalibQAChecker::QualityFlag_t quality)
559 // get quality color for quality
561 Color_t info = kSpring-8;
562 Color_t warning = kOrange;
563 Color_t error = kRed;
564 Color_t fatal = kRed+2;
565 Color_t none = kWhite;
586 //_________________________________________________________________________
587 const char* AliTPCCalibQAChecker::QualityDescription(const QualityFlag_t quality) const
590 // return description for quality
592 if (!fArrAlarmDescriptions || !fArrAlarmDescriptions->At(quality)) return "";
593 TString s(fArrAlarmDescriptions->At(quality)->GetName());
595 min+=fThresMin[quality];
596 max+=fThresMax[quality];
597 s.ReplaceAll("%min",min);
598 s.ReplaceAll("%max",max);
601 //_________________________________________________________________________
602 Int_t AliTPCCalibQAChecker::DrawInPad(TVirtualPad *pad, Int_t sub)
608 if (fArrSubCheckers){
609 if (fArrSubCheckers->GetEntries()>0){
610 TIter next(fArrSubCheckers);
612 while ( (o=next()) ) {
613 AliTPCCalibQAChecker *al=(AliTPCCalibQAChecker*)o;
614 sub=al->DrawInPad(pad,sub);
624 //_________________________________________________________________________
625 void AliTPCCalibQAChecker::DrawSubNodes(Option_t */*option*/)
628 // divide the current pad in sub pads and draw all sub nodes
630 if (!gPad) new TCanvas;
631 TPad *mother=(TPad*)gPad;
633 //find number of sub pads needed
634 Int_t nPads = GetNumberOfSubCheckers();
635 Int_t nCols = (Int_t)TMath::Ceil( TMath::Sqrt(nPads) );
636 Int_t nRows = (Int_t)TMath::Ceil( (Double_t)nPads/(Double_t)nCols );
637 gPad->Divide(nCols,nRows);
643 while ( (pad=(TPad*)mother->GetPad(ipad)) ){
644 TFrame* frame=(TFrame*)pad->GetPrimitive("TFrame");
645 if (frame) frame->SetFillColor(kWhite);
646 else {printf("no frame\n");}
652 //_________________________________________________________________________
653 void AliTPCCalibQAChecker::DrawRepresentationHist(Option_t *option)
656 // Draw the representation histogram
658 if (!fHistRep) return;
659 if (!gPad) new TCanvas;
660 Bool_t withBackColor=kTRUE;
665 if (opt.Contains("nobc")) withBackColor=kFALSE;
666 opt.ReplaceAll("nobc","");
668 if (opt.IsNull()) opt=fStrDrawRepOpt;
671 opt.ReplaceAll("prof","");
673 fHistRep->Draw(opt.Data());
675 gPad->SetFillColor(GetQualityColor());
680 //_________________________________________________________________________
681 void AliTPCCalibQAChecker::AddQualityLines(TH1 *hist)
684 // add lines indicating the quality level to hist
687 Double_t xmin=hist->GetXaxis()->GetXmin();
688 Double_t xmax=hist->GetXaxis()->GetXmax();
692 for (Int_t i=(Int_t)kINFO; i<kNQualityFlags; ++i){
693 if (fThresMin[i]>=fThresMax[i]) continue;
696 TLine *lmin=new TLine(xmin,min,xmax,min);
697 lmin->SetLineWidth(2);
698 lmin->SetLineColor(QualityColor((QualityFlag_t)i));
699 TLine *lmax=new TLine(xmin,max,xmax,max);
700 lmax->SetLineWidth(2);
701 lmax->SetLineColor(QualityColor((QualityFlag_t)i));
702 hist->GetListOfFunctions()->Add(lmin);
703 hist->GetListOfFunctions()->Add(lmax);
705 hist->SetAxisRange(min,max,"Y");