1 /**************************************************************************
2 * Copyright(c) 2004, 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 **************************************************************************/
15 //__________________________________________________________________
25 //__________________________________________________________________
27 // --- ROOT system ---
32 #include <TIterator.h>
37 #include <TPaveText.h>
41 // --- AliRoot header files ---
44 #include "AliQAChecker.h"
45 #include "AliFMDQAChecker.h"
46 #include "AliRecoParam.h"
48 ClassImp(AliFMDQAChecker)
50 ; // This is for Emacs! - do not delete
53 //__________________________________________________________________
55 AliFMDQAChecker::CheckOne(AliQAv1::ALITASK_t what,
56 AliRecoParam::EventSpecie_t specie,
59 if(what == AliQAv1::kESD) return CheckESD(specie, hist);
60 if(what == AliQAv1::kRAW) return CheckRaw(specie, hist);
61 if(what == AliQAv1::kSIM) return CheckSim(specie, hist);
62 if(what == AliQAv1::kREC) return CheckRec(specie, hist);
65 //__________________________________________________________________
67 AliFMDQAChecker::CheckESD(AliRecoParam::EventSpecie_t /* specie*/,
70 return (hist->GetMean() > 0 ? 1 : 0);
72 //__________________________________________________________________
74 AliFMDQAChecker::CheckRaw(AliRecoParam::EventSpecie_t /* specie*/,
77 return (hist->GetMean() > 0 ? 1 : 0);
79 //__________________________________________________________________
81 AliFMDQAChecker::CheckSim(AliRecoParam::EventSpecie_t /* specie*/,
84 return (hist->GetMean() > 0 ? 1 : 0);
86 //__________________________________________________________________
88 AliFMDQAChecker::CheckRec(AliRecoParam::EventSpecie_t /* specie*/,
91 return (hist->GetMean() > 0 ? 1 : 0);
94 //__________________________________________________________________
95 void AliFMDQAChecker::Check(Double_t* rv,
96 AliQAv1::ALITASK_t what,
98 const AliDetectorRecoParam* /*t*/)
101 // Member function called to do the actual checking
104 // rv Array of return values.
105 // what What to check
106 // list Array of arrays of histograms. There's one arrat for
108 // t Reconstruction parameters - not used.
111 // Double_t* rv = new Double_t[AliRecoParam::kNSpecies] ;
112 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
116 if (!AliQAv1::Instance()->IsEventSpecieSet(specie) )
119 if(!list[specie]) continue;
122 Int_t nHist = list[specie]->GetEntriesFast();
123 for(Int_t i= 0; i< nHist; i++) {
125 if (!(hist = static_cast<TH1*>(list[specie]->At(i)))) continue;
127 rv[specie] += CheckOne(what, AliRecoParam::ConvertIndex(specie), hist);
129 // if (count != 0) rv[specie] /= count;
135 Int_t CheckForLog(TAxis* axis,
140 TString t(axis->GetTitle());
141 if (!t.Contains("[log]", TString::kIgnoreCase)) return 0;
142 t.ReplaceAll("[log]", "");
144 case 1: pad->SetLogx(); ret |= 0x1; break;
145 case 2: pad->SetLogy(); ret |= 0x2; break;
146 case 3: pad->SetLogz(); ret |= 0x4; break;
151 void RestoreLog(TAxis* axis, Bool_t log)
154 TString t(axis->GetTitle());
161 void FindMinMax(TH1* h, Double_t& min, Double_t& max)
165 for (Int_t i = 1; i <= h->GetNbinsX(); i++) {
166 Double_t c = h->GetBinContent(i);
167 if (c < 1e-8) continue;
168 tmin = TMath::Min(tmin, c);
169 tmax = TMath::Max(tmax, c);
175 //____________________________________________________________________________
177 AliFMDQAChecker::MakeImage(TObjArray** list,
178 AliQAv1::TASKINDEX_t task,
179 AliQAv1::MODE_t mode)
181 // makes the QA image for sim and rec
184 // task What to check
185 // list Array of arrays of histograms. There's one array for
187 // t Reconstruction parameters - not used.
191 Double_t min = 10000;
193 // Loop over all species
194 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
195 AliRecoParam::EventSpecie_t spe = AliRecoParam::ConvertIndex(specie);
196 if (!AliQAv1::Instance(AliQAv1::GetDetIndex(GetName()))
197 ->IsEventSpecieSet(spe))
200 // If nothing is defined for this specie, go on.
201 if(!list[specie] || list[specie]->GetEntriesFast() == 0) continue;
203 // Loop over the histograms and figure out how many histograms we
204 // have and the min/max
206 Int_t nHist = list[specie]->GetEntriesFast();
207 for(Int_t i= 0; i< nHist; i++) {
208 hist = static_cast<TH1F*>(list[specie]->At(i));
209 if (hist && hist->TestBit(AliQAv1::GetImageBit())) {
211 TString name(hist->GetName());
212 if (name.Contains("readouterrors", TString::kIgnoreCase)) continue;
214 // Double_t hMax = hist->GetMaximum();
215 // hist->GetBinContent(hist->GetMaximumBin());
216 // Double_t hMin = hist->GetMinimum();
217 // hist->GetBinContent(hist->GetMinimumBin());
219 FindMinMax(hist, hMin, hMax);
220 max = TMath::Max(max, hMax);
221 min = TMath::Min(min, hMin);
226 min = TMath::Max(1e-6, min);
227 max = TMath::Max(1.0, max);
229 // IF no images, go on.
231 AliDebug(AliQAv1::GetQADebugLevel(),
232 Form("No histogram will be plotted for %s %s\n", GetName(),
233 AliQAv1::GetTaskName(task).Data()));
237 AliDebug(AliQAv1::GetQADebugLevel(),
238 Form("%d histograms will be plotted for %s %s\n",
239 nImages, GetName(), AliQAv1::GetTaskName(task).Data()));
240 gStyle->SetOptStat(0);
242 // Again loop over species and draw a canvas
243 for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
244 if (!AliQAv1::Instance(AliQAv1::GetDetIndex(GetName()))
245 ->IsEventSpecieSet(specie)) continue;
247 // if Nothing here, go on
248 if(!list[specie] || list[specie]->GetEntries() <= 0 ||
249 nImages <= 0) continue;
252 const Char_t * title = Form("QA_%s_%s_%s", GetName(),
253 AliQAv1::GetTaskName(task).Data(),
254 AliRecoParam::GetEventSpecieName(specie));
255 if (!fImage[specie]) fImage[specie] = new TCanvas(title, title) ;
256 fImage[specie]->Clear() ;
257 fImage[specie]->SetTitle(title) ;
258 fImage[specie]->cd() ;
260 // Put something in the canvas - even if empty
261 TPaveText someText(0.015, 0.015, 0.98, 0.98) ;
262 someText.AddText(title) ;
263 someText.SetFillColor(0);
264 someText.SetFillStyle(0);
265 someText.SetBorderSize(0);
266 someText.SetTextColor(kRed+1);
268 TString outName(Form("%s%s%d.%s", AliQAv1::GetImageFileName(),
269 AliQAv1::GetModeName(mode),
270 AliQAChecker::Instance()->GetRunNumber(),
271 AliQAv1::GetImageFileFormat()));
272 fImage[specie]->Print(outName, "ps") ;
275 // Now set some parameters on the canvas
276 fImage[specie]->Clear();
277 fImage[specie]->SetTopMargin(0.10);
278 fImage[specie]->SetBottomMargin(0.15);
279 fImage[specie]->SetLeftMargin(0.15);
280 fImage[specie]->SetRightMargin(0.05);
283 const char* topT = Form("Mode: %s, Task: %s, Specie: %s, Run: %d",
284 AliQAv1::GetModeName(mode),
285 AliQAv1::GetTaskName(task).Data(),
286 AliRecoParam::GetEventSpecieName(specie),
287 AliQAChecker::Instance()->GetRunNumber());
288 TLatex* topText = new TLatex(.5, .99, topT);
289 topText->SetTextAlign(23);
290 topText->SetTextSize(.038);
291 topText->SetTextFont(42);
292 topText->SetTextColor(kBlue+3);
297 Int_t nx = int(nImages + .5) / 2;
299 fImage[specie]->Divide(nx, ny, 0, 0);
302 // Loop over histograms
304 Int_t nHist = list[specie]->GetEntriesFast();
306 for (Int_t i = 0; i < nHist; i++) {
307 hist = static_cast<TH1*>(list[specie]->At(i));
308 if (!hist || !hist->TestBit(AliQAv1::GetImageBit())) continue;
311 TVirtualPad* pad = fImage[specie]->cd(++j);
312 pad->SetRightMargin(0.01);
314 // Check for log scale
316 logOpts |= CheckForLog(hist->GetXaxis(), pad, 1);
317 logOpts |= CheckForLog(hist->GetYaxis(), pad, 2);
318 logOpts |= CheckForLog(hist->GetZaxis(), pad, 3);
320 // Figure out special cases
322 TString name(hist->GetName());
323 if (name.Contains("readouterrors", TString::kIgnoreCase)) {
324 pad->SetRightMargin(0.15);
325 pad->SetBottomMargin(0.10);
326 pad->SetTopMargin(0.02);
332 hist->SetMinimum(min);
333 hist->SetMaximum(max);
340 if (name.Contains("readouterrors", TString::kIgnoreCase)) {
341 for (Int_t kk = 1; kk <= 3; kk++) {
342 TH1* proj = static_cast<TH2*>(hist)->ProjectionY("",kk,kk);
343 Double_t m = proj->GetMean();
344 TLatex* l = new TLatex(kk, 30, Form("Mean: %f", m));
346 l->SetTextColor(m > 10 ? kRed+1 : m > .7 ? kOrange+2 :kGreen+2);
351 gStyle->SetOptTitle(0);
352 TPad* insert = new TPad("insert", "Zoom",
353 .4,.4, .99, .95, 0, 0, 0);
354 insert->SetTopMargin(0.01);
355 insert->SetRightMargin(0.01);
356 insert->SetFillColor(0);
357 insert->SetBorderSize(1);
358 insert->SetBorderMode(0);
361 if (logOpts & 0x1) insert->SetLogx();
362 if (logOpts & 0x2) insert->SetLogy();
363 if (logOpts & 0x4) insert->SetLogz();
364 hist->GetXaxis()->SetRange(1, hist->GetNbinsX()/8);
365 TH1* copy = hist->DrawCopy(opt);
366 copy->GetXaxis()->SetNdivisions(408, false);
367 // Restore full range
368 hist->GetXaxis()->SetRange(0, 0);
369 gStyle->SetOptTitle(1);
372 // Possibly restore the log options
373 RestoreLog(hist->GetXaxis(), logOpts & 0x1);
374 RestoreLog(hist->GetYaxis(), logOpts & 0x2);
375 RestoreLog(hist->GetZaxis(), logOpts & 0x4);
377 // Print to a post-script file
378 fImage[specie]->Print(outName, "ps");
382 //__________________________________________________________________