]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TPC/AliTPCCalibViewer.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge Jens and Christian changes
[u/mrichter/AliRoot.git] / TPC / AliTPCCalibViewer.cxx
CommitLineData
39bcd65d 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///////////////////////////////////////////////////////////////////////////////
18// //
19// Class for viewing/visualizing TPC calibration data //
20// base on TTree functionality for visualization //
72d0ab7e 21// //
22// Create a list of AliTPCCalPads, arrange them in an TObjArray. //
23// Pass this TObjArray to MakeTree and create the calibration Tree //
24// While craating this tree some statistical information are calculated //
25// Open the viewer with this Tree: AliTPCCalibViewer v("CalibTree.root") //
26// Have fun! //
27// EasyDraw("CETmean~-CETmean_mean", "A", "(CETmean~-CETmean_mean)>0") //
28// //
29// If you like to click, we recommand you the //
30// AliTPCCalibViewerGUI //
31// //
32// THE DOCUMENTATION IS STILL NOT COMPLETED !!!! //
33// //
39bcd65d 34///////////////////////////////////////////////////////////////////////////////
35
36//
37// ROOT includes
38//
39#include <iostream>
4f3934a1 40#include <fstream>
39bcd65d 41#include <TString.h>
42#include <TRandom.h>
43#include <TLegend.h>
44#include <TLine.h>
45#include <TCanvas.h>
46#include <TROOT.h>
47#include <TStyle.h>
72d0ab7e 48#include <TH1.h>
39bcd65d 49#include <TH1F.h>
50#include <THashTable.h>
51#include <TObjString.h>
52#include "TTreeStream.h"
53#include "TFile.h"
54#include "TKey.h"
4f3934a1 55#include "TGraph.h"
56#include "AliTPCCalibPulser.h"
57#include "AliTPCCalibPedestal.h"
58#include "AliTPCCalibCE.h"
59// #include "TObjArray.h"
60// #include "TObjString.h"
61// #include "TString.h"
62// #include "AliTPCCalPad.h"
39bcd65d 63
64
65//
66// AliRoot includes
67//
68#include "AliTPCCalibViewer.h"
69
70ClassImp(AliTPCCalibViewer)
71
72AliTPCCalibViewer::AliTPCCalibViewer()
73 :TObject(),
74 fTree(0),
75 fFile(0),
a6d2bd0c 76 fListOfObjectsToBeDeleted(0),
77 fTreeMustBeDeleted(0)
39bcd65d 78{
79 //
80 // Default constructor
81 //
82
83}
84
85//_____________________________________________________________________________
86AliTPCCalibViewer::AliTPCCalibViewer(const AliTPCCalibViewer &c)
87 :TObject(c),
88 fTree(0),
89 fFile(0),
a6d2bd0c 90 fListOfObjectsToBeDeleted(0),
91 fTreeMustBeDeleted(0)
39bcd65d 92{
93 //
94 // dummy AliTPCCalibViewer copy constructor
95 // not yet working!!!
96 //
97 fTree = c.fTree;
a6d2bd0c 98 fTreeMustBeDeleted = c.fTreeMustBeDeleted;
39bcd65d 99 //fFile = new TFile(*(c.fFile));
100 fListOfObjectsToBeDeleted = c.fListOfObjectsToBeDeleted;
101}
102
103//_____________________________________________________________________________
104AliTPCCalibViewer::AliTPCCalibViewer(TTree* tree)
105 :TObject(),
106 fTree(0),
107 fFile(0),
a6d2bd0c 108 fListOfObjectsToBeDeleted(0),
109 fTreeMustBeDeleted(0)
39bcd65d 110{
111 //
112 // Constructor that initializes the calibration viewer
113 //
114 fTree = tree;
a6d2bd0c 115 fTreeMustBeDeleted = kFALSE;
39bcd65d 116 fListOfObjectsToBeDeleted = new TObjArray();
117}
118
119//_____________________________________________________________________________
120AliTPCCalibViewer::AliTPCCalibViewer(char* fileName, char* treeName)
121 :TObject(),
122 fTree(0),
123 fFile(0),
a6d2bd0c 124 fListOfObjectsToBeDeleted(0),
125 fTreeMustBeDeleted(0)
39bcd65d 126{
127 //
128 // Constructor to initialize the calibration viewer
129 // the file 'fileName' contains the tree 'treeName'
130 //
131 fFile = new TFile(fileName, "read");
132 fTree = (TTree*) fFile->Get(treeName);
a6d2bd0c 133 fTreeMustBeDeleted = kTRUE;
39bcd65d 134 fListOfObjectsToBeDeleted = new TObjArray();
135}
136
137//____________________________________________________________________________
138AliTPCCalibViewer & AliTPCCalibViewer::operator =(const AliTPCCalibViewer & param)
139{
140 //
141 // assignment operator - dummy
142 // not yet working!!!
143 //
144 fTree = param.fTree;
a6d2bd0c 145 fTreeMustBeDeleted = param.fTreeMustBeDeleted;
39bcd65d 146 //fFile = new TFile(*(param.fFile));
147 fListOfObjectsToBeDeleted = param.fListOfObjectsToBeDeleted;
148 return (*this);
149}
150
151//_____________________________________________________________________________
152AliTPCCalibViewer::~AliTPCCalibViewer()
153{
154 //
155 // AliTPCCalibViewer destructor
a6d2bd0c 156 // all objects will be deleted, the file will be closed, the pictures will disappear
39bcd65d 157 //
a6d2bd0c 158 if (fTree && fTreeMustBeDeleted) {
159 fTree->SetCacheSize(0);
160 fTree->Delete();
161 //delete fTree;
162 }
39bcd65d 163 if (fFile) {
164 fFile->Close();
165 fFile = 0;
166 }
167
168 for (Int_t i = fListOfObjectsToBeDeleted->GetEntriesFast()-1; i >= 0; i--) {
169 //cout << "Index " << i << " trying to delete the following object: " << fListOfObjectsToBeDeleted->At(i)->GetName() << "..."<< endl;
170 delete fListOfObjectsToBeDeleted->At(i);
171 }
172 delete fListOfObjectsToBeDeleted;
173}
174
a6d2bd0c 175//_____________________________________________________________________________
176void AliTPCCalibViewer::Delete(Option_t* option) {
177 //
178 // Should be called from AliTPCCalibViewerGUI class only.
179 // If you use Delete() do not call the destructor.
180 // All objects (except those contained in fListOfObjectsToBeDeleted) will be deleted, the file will be closed.
181 //
182
183 option = option; // to avoid warnings on compiling
184 if (fTree && fTreeMustBeDeleted) {
185 fTree->SetCacheSize(0);
186 fTree->Delete();
187 }
188 if (fFile)
189 delete fFile;
190 delete fListOfObjectsToBeDeleted;
191}
192
39bcd65d 193//_____________________________________________________________________________
194Int_t AliTPCCalibViewer::EasyDraw(const char* drawCommand, const char* sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
195 //
196 // easy drawing of data, use '~' for abbreviation of '.fElements'
197 // example: EasyDraw("CETmean~-CETmean_mean", "A", "(CETmean~-CETmean_mean)>0")
198 // sector: sector-number - only the specified sector will be drwawn
199 // 'A'/'C' or 'a'/'c' - side A/C will be drawn
200 // 'ALL' - whole TPC will be drawn, projected on one side
201 // cuts: specifies cuts
202 // drawOptions: draw options like 'same'
203 // writeDrawCommand: write the command, that is passed to TTree::Draw
204 //
72d0ab7e 205
39bcd65d 206 TString drawStr(drawCommand);
207 TString sectorStr(sector);
208 sectorStr.ToUpper();
209 TString cutStr("");
72d0ab7e 210 //TString drawOptionsStr("profcolz ");
211 TString drawOptionsStr("");
39bcd65d 212 TRandom rnd(0);
213 Int_t rndNumber = rnd.Integer(10000);
72d0ab7e 214
215 if (drawOptions && strcmp(drawOptions, "") != 0)
39bcd65d 216 drawOptionsStr += drawOptions;
72d0ab7e 217 else
218 drawOptionsStr += "profcolz";
39bcd65d 219
220 if (sectorStr == "A") {
221 drawStr += ":gy.fElements:gx.fElements>>prof";
222 drawStr += rndNumber;
223 drawStr += "(330,-250,250,330,-250,250)";
224 cutStr += "(sector/18)%2==0 ";
225 }
226 else if (sectorStr == "C") {
227 drawStr += ":gy.fElements:gx.fElements>>prof";
228 drawStr += rndNumber;
229 drawStr += "(330,-250,250,330,-250,250)";
230 cutStr += "(sector/18)%2==1 ";
231 }
232 else if (sectorStr == "ALL") {
233 drawStr += ":gy.fElements:gx.fElements>>prof";
234 drawStr += rndNumber;
235 drawStr += "(330,-250,250,330,-250,250)";
236 }
237 else if (sectorStr.IsDigit()) {
238 Int_t isec = sectorStr.Atoi();
239 drawStr += ":rpad.fElements:row.fElements>>prof";
240 drawStr += rndNumber;
241 if (isec < 36 && isec >= 0)
242 drawStr += "(63,0,63,108,-54,54)";
243 else if (isec < 72 && isec >= 36)
244 drawStr += "(96,0,96,140,-70,70)";
245 else {
246 Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
247 return -1;
248 }
249 cutStr += "(sector==";
250 cutStr += isec;
251 cutStr += ") ";
252 }
253
254 if (cuts && cuts[0] != 0) {
255 if (cutStr.Length() != 0) cutStr += "&& ";
256 cutStr += "(";
257 cutStr += cuts;
258 cutStr += ")";
259 }
260 drawStr.ReplaceAll("~", ".fElements");
261 cutStr.ReplaceAll("~", ".fElements");
262 if (writeDrawCommand) std::cout << "fTree->Draw(\"" << drawStr << "\", \"" << cutStr << "\", \"" << drawOptionsStr << "\");" << std::endl;
263 return fTree->Draw(drawStr.Data(), cutStr.Data(), drawOptionsStr.Data());
264}
265
72d0ab7e 266
39bcd65d 267Int_t AliTPCCalibViewer::EasyDraw(const char* drawCommand, Int_t sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
268 //
269 // easy drawing of data, use '~' for abbreviation of '.fElements'
270 // example: EasyDraw("CETmean~-CETmean_mean", 34, "(CETmean~-CETmean_mean)>0")
271 // sector: sector-number - only the specified sector will be drwawn
272 // cuts: specifies cuts
273 // drawOptions: draw options like 'same'
274 // writeDrawCommand: write the command, that is passed to TTree::Draw
275 //
276 if (sector >= 0 && sector < 72) {
277 char sectorChr[3];
278 sprintf(sectorChr, "%i", sector);
279 return EasyDraw(drawCommand, sectorChr, cuts, drawOptions, writeDrawCommand);
280 }
281 Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
282 return -1;
283}
284
72d0ab7e 285
39bcd65d 286//_____________________________________________________________________________
287Int_t AliTPCCalibViewer::EasyDraw1D(const char* drawCommand, const char* sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
288 //
289 // easy drawing of data, use '~' for abbreviation of '.fElements'
290 // example: EasyDraw("CETmean~-CETmean_mean", "A", "(CETmean~-CETmean_mean)>0")
291 // sector: sector-number - the specified sector will be drwawn
292 // 'A'/'C' or 'a'/'c' - side A/C will be drawn
293 // 'ALL' - whole TPC will be drawn, projected on one side
294 // cuts: specifies cuts
295 // drawOptions: draw options like 'same'
296 // writeDrawCommand: write the command, that is passed to TTree::Draw
297 //
298
299 TString drawStr(drawCommand);
300 TString sectorStr(sector);
301 TString drawOptionsStr(drawOptions);
302 sectorStr.ToUpper();
303 TString cutStr("");
304
305 if (sectorStr == "A")
306 cutStr += "(sector/18)%2==0 ";
307 else if (sectorStr == "C")
308 cutStr += "(sector/18)%2==1 ";
309 else if (sectorStr.IsDigit()) {
310 Int_t isec = sectorStr.Atoi();
311 if (isec < 0 || isec > 71) {
312 Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
313 return -1;
314 }
315 cutStr += "(sector==";
316 cutStr += isec;
317 cutStr += ") ";
318 }
319
320 if (cuts && cuts[0] != 0) {
321 if (cutStr.Length() != 0) cutStr += "&& ";
322 cutStr += "(";
323 cutStr += cuts;
324 cutStr += ")";
325 }
326
327 drawStr.ReplaceAll("~", ".fElements");
328 cutStr.ReplaceAll("~", ".fElements");
329 if (writeDrawCommand) std::cout << "fTree->Draw(\"" << drawStr << "\", \"" << cutStr << "\", \"" << drawOptionsStr << "\");" << std::endl;
330 return fTree->Draw(drawStr.Data(), cutStr.Data(), drawOptionsStr.Data());
331}
332
72d0ab7e 333
39bcd65d 334Int_t AliTPCCalibViewer::EasyDraw1D(const char* drawCommand, Int_t sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
335 //
336 // easy drawing of data, use '~' for abbreviation of '.fElements'
337 // example: EasyDraw("CETmean~-CETmean_mean", 34, "(CETmean~-CETmean_mean)>0")
338 // sector: sector-number - the specified sector will be drwawn
339 // cuts: specifies cuts
340 // drawOptions: draw options like 'same'
341 // writeDrawCommand: write the command, that is passed to TTree::Draw
342 //
343
344 if (sector >= 0 && sector < 72) {
345 char sectorChr[3];
346 sprintf(sectorChr, "%i", sector);
347 return EasyDraw1D(drawCommand, sectorChr, cuts, drawOptions, writeDrawCommand);
348 }
349 Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
350 return -1;
351}
352
39bcd65d 353
72d0ab7e 354Int_t AliTPCCalibViewer::DrawHisto1D(const char* drawCommand, Int_t sector, const char* cuts, const char *sigmas, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM) const {
355 //
356 // Easy drawing of data, in principle the same as EasyDraw1D
357 // Difference: A line for the mean / median / LTM is drawn
358 // in 'sigmas' you can specify in which distance to the mean/median/LTM you want to see a line in sigma-units, separated by ';'
359 // example: sigmas = "2; 4; 6;" at Begin_Latex 2 #sigma End_Latex, Begin_Latex 4 #sigma End_Latex and Begin_Latex 6 #sigma End_Latex a line is drawn.
360 // "plotMean", "plotMedian" and "plotLTM": what kind of lines do you want to see?
361 //
362 if (sector >= 0 && sector < 72) {
363 char sectorChr[3];
364 sprintf(sectorChr, "%i", sector);
365 return DrawHisto1D(drawCommand, sectorChr, cuts, sigmas, plotMean, plotMedian, plotLTM);
366 }
367 Error("DrawHisto1D","The TPC contains only sectors between 0 and 71.");
368 return -1;
369}
39bcd65d 370
72d0ab7e 371
372Int_t AliTPCCalibViewer::DrawHisto1D(const char* drawCommand, const char* sector, const char* cuts, const char *sigmas, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM) const {
373 //
374 // Easy drawing of data, in principle the same as EasyDraw1D
375 // Difference: A line for the mean / median / LTM is drawn
376 // in 'sigmas' you can specify in which distance to the mean/median/LTM you want to see a line in sigma-units, separated by ';'
377 // example: sigmas = "2; 4; 6;" at Begin_Latex 2 #sigma End_Latex, Begin_Latex 4 #sigma End_Latex and Begin_Latex 6 #sigma End_Latex a line is drawn.
378 // "plotMean", "plotMedian" and "plotLTM": what kind of lines do you want to see?
379 //
39bcd65d 380 Int_t oldOptStat = gStyle->GetOptStat();
381 gStyle->SetOptStat(0000000);
72d0ab7e 382 Double_t ltmFraction = 0.8;
383
384 TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
385 TVectorF nsigma(sigmasTokens->GetEntriesFast());
386 for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
387 TString str(((TObjString*)sigmasTokens->At(i))->GetString());
388 Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
389 nsigma[i] = sig;
39bcd65d 390 }
391
72d0ab7e 392 TString drawStr(drawCommand);
393 drawStr += " >> tempHist";
394 Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts);
395 TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
396 // FIXME is this histogram deleted automatically?
397 Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
398
399 Double_t mean = TMath::Mean(entries, values);
400 Double_t median = TMath::Median(entries, values);
401 Double_t sigma = TMath::RMS(entries, values);
402 Double_t maxY = htemp->GetMaximum();
403
39bcd65d 404 char c[500];
72d0ab7e 405 TLegend * legend = new TLegend(.7,.7, .99, .99, "Statistical information");
406// sprintf(c, "%s, sector: %i", type, sector);
39bcd65d 407 fListOfObjectsToBeDeleted->Add(legend);
408
39bcd65d 409 if (plotMean) {
72d0ab7e 410 // draw Mean
411 TLine* line = new TLine(mean, 0, mean, maxY);
39bcd65d 412 fListOfObjectsToBeDeleted->Add(line);
413 line->SetLineColor(kRed);
414 line->SetLineWidth(2);
415 line->SetLineStyle(1);
416 line->Draw();
72d0ab7e 417 sprintf(c, "Mean: %f", mean);
39bcd65d 418 legend->AddEntry(line, c, "l");
72d0ab7e 419 // draw sigma lines
39bcd65d 420 for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
72d0ab7e 421 TLine* linePlusSigma = new TLine(mean + nsigma[i] * sigma, 0, mean + nsigma[i] * sigma, maxY);
39bcd65d 422 fListOfObjectsToBeDeleted->Add(linePlusSigma);
423 linePlusSigma->SetLineColor(kRed);
72d0ab7e 424 linePlusSigma->SetLineStyle(2 + i);
39bcd65d 425 linePlusSigma->Draw();
72d0ab7e 426 TLine* lineMinusSigma = new TLine(mean - nsigma[i] * sigma, 0, mean - nsigma[i] * sigma, maxY);
39bcd65d 427 fListOfObjectsToBeDeleted->Add(lineMinusSigma);
428 lineMinusSigma->SetLineColor(kRed);
72d0ab7e 429 lineMinusSigma->SetLineStyle(2 + i);
39bcd65d 430 lineMinusSigma->Draw();
72d0ab7e 431 sprintf(c, "%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma));
39bcd65d 432 legend->AddEntry(lineMinusSigma, c, "l");
433 }
434 }
39bcd65d 435 if (plotMedian) {
72d0ab7e 436 // draw median
437 TLine* line = new TLine(median, 0, median, maxY);
39bcd65d 438 fListOfObjectsToBeDeleted->Add(line);
439 line->SetLineColor(kBlue);
440 line->SetLineWidth(2);
441 line->SetLineStyle(1);
442 line->Draw();
72d0ab7e 443 sprintf(c, "Median: %f", median);
39bcd65d 444 legend->AddEntry(line, c, "l");
72d0ab7e 445 // draw sigma lines
39bcd65d 446 for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
72d0ab7e 447 TLine* linePlusSigma = new TLine(median + nsigma[i] * sigma, 0, median + nsigma[i]*sigma, maxY);
39bcd65d 448 fListOfObjectsToBeDeleted->Add(linePlusSigma);
449 linePlusSigma->SetLineColor(kBlue);
72d0ab7e 450 linePlusSigma->SetLineStyle(2 + i);
39bcd65d 451 linePlusSigma->Draw();
72d0ab7e 452 TLine* lineMinusSigma = new TLine(median - nsigma[i] * sigma, 0, median - nsigma[i]*sigma, maxY);
39bcd65d 453 fListOfObjectsToBeDeleted->Add(lineMinusSigma);
454 lineMinusSigma->SetLineColor(kBlue);
72d0ab7e 455 lineMinusSigma->SetLineStyle(2 + i);
39bcd65d 456 lineMinusSigma->Draw();
72d0ab7e 457 sprintf(c, "%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma));
39bcd65d 458 legend->AddEntry(lineMinusSigma, c, "l");
459 }
460 }
39bcd65d 461 if (plotLTM) {
72d0ab7e 462 // draw LTM
463 Double_t ltmRms = 0;
464 Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
465 TLine* line = new TLine(ltm, 0, ltm, maxY);
39bcd65d 466 fListOfObjectsToBeDeleted->Add(line);
467 line->SetLineColor(kGreen+2);
468 line->SetLineWidth(2);
469 line->SetLineStyle(1);
470 line->Draw();
72d0ab7e 471 sprintf(c, "LTM: %f", ltm);
39bcd65d 472 legend->AddEntry(line, c, "l");
72d0ab7e 473 // draw sigma lines
39bcd65d 474 for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
72d0ab7e 475 TLine* linePlusSigma = new TLine(ltm + nsigma[i] * ltmRms, 0, ltm + nsigma[i] * ltmRms, maxY);
39bcd65d 476 fListOfObjectsToBeDeleted->Add(linePlusSigma);
477 linePlusSigma->SetLineColor(kGreen+2);
478 linePlusSigma->SetLineStyle(2+i);
479 linePlusSigma->Draw();
480
72d0ab7e 481 TLine* lineMinusSigma = new TLine(ltm - nsigma[i] * ltmRms, 0, ltm - nsigma[i] * ltmRms, maxY);
39bcd65d 482 fListOfObjectsToBeDeleted->Add(lineMinusSigma);
483 lineMinusSigma->SetLineColor(kGreen+2);
484 lineMinusSigma->SetLineStyle(2+i);
485 lineMinusSigma->Draw();
72d0ab7e 486 sprintf(c, "%i #sigma = %f", (Int_t)(nsigma[i]), (Float_t)(nsigma[i] * ltmRms));
39bcd65d 487 legend->AddEntry(lineMinusSigma, c, "l");
488 }
489 }
72d0ab7e 490 if (!plotMean && !plotMedian && !plotLTM) return -1;
39bcd65d 491 legend->Draw();
492 gStyle->SetOptStat(oldOptStat);
72d0ab7e 493 return 1;
39bcd65d 494}
495
39bcd65d 496
72d0ab7e 497Int_t AliTPCCalibViewer::SigmaCut(const char* drawCommand, Int_t sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, Bool_t pm, const char *sigmas, Float_t sigmaStep) const {
498 //
499 // Creates a histogram Begin_Latex S(t, #mu, #sigma) End_Latex, where you can see, how much of the data are inside sigma-intervals around the mean value
500 // The data of the distribution Begin_Latex f(x, #mu, #sigma) End_Latex are given in 'array', 'n' specifies the length of the array
501 // 'mean' and 'sigma' are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in 'array', to be specified by the user
502 // 'nbins': number of bins, 'binLow': first bin, 'binUp': last bin
503 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma, Begin_Latex t #sigma End_Latex)
504 // sigmaStep: the binsize of the generated histogram
505 // Begin_Latex
506 // f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{#mu}^{#mu + t #sigma} f(x, #mu, #sigma) dx + #int_{#mu}^{#mu - t #sigma} f(x, #mu, #sigma) dx }{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
507 // End_Latex
508 //
509 //
510 // Creates a histogram, where you can see, how much of the data are inside sigma-intervals
511 // around the mean/median/LTM
512 // with drawCommand, sector and cuts you specify your input data, see EasyDraw
513 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma)
514 // sigmaStep: the binsize of the generated histogram
515 // plotMean/plotMedian/plotLTM: specifies where to put the center
516 //
517 if (sector >= 0 && sector < 72) {
518 char sectorChr[3];
519 sprintf(sectorChr, "%i", sector);
520 return SigmaCut(drawCommand, sectorChr, cuts, sigmaMax, plotMean, plotMedian, plotLTM, pm, sigmas, sigmaStep);
521 }
522 Error("SigmaCut","The TPC contains only sectors between 0 and 71.");
523 return -1;
524}
525
526
527Int_t AliTPCCalibViewer::SigmaCut(const char* drawCommand, const char* sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, Bool_t pm, const char *sigmas, Float_t sigmaStep) const {
528 //
529 // Creates a histogram, where you can see, how much of the data are inside sigma-intervals
530 // around the mean/median/LTM
531 // with drawCommand, sector and cuts you specify your input data, see EasyDraw
532 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma)
533 // sigmaStep: the binsize of the generated histogram
534 // plotMean/plotMedian/plotLTM: specifies where to put the center
535 //
536
537 Double_t ltmFraction = 0.8;
538
539 TString drawStr(drawCommand);
540 drawStr += " >> tempHist";
541
542 Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
543 TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
544 // FIXME is this histogram deleted automatically?
545 Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
546
547 Double_t mean = TMath::Mean(entries, values);
548 Double_t median = TMath::Median(entries, values);
549 Double_t sigma = TMath::RMS(entries, values);
550
551 TLegend * legend = new TLegend(.7,.7, .99, .99, "Cumulative");
552 fListOfObjectsToBeDeleted->Add(legend);
553 TH1F *cutHistoMean = 0;
554 TH1F *cutHistoMedian = 0;
555 TH1F *cutHistoLTM = 0;
556
557 TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
558 TVectorF nsigma(sigmasTokens->GetEntriesFast());
559 for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
560 TString str(((TObjString*)sigmasTokens->At(i))->GetString());
561 Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
562 nsigma[i] = sig;
563 }
564
565 if (plotMean) {
566 cutHistoMean = AliTPCCalibViewer::SigmaCut(htemp, mean, sigma, sigmaMax, sigmaStep, pm);
567 if (cutHistoMean) {
568 fListOfObjectsToBeDeleted->Add(cutHistoMean);
569 cutHistoMean->SetLineColor(kRed);
570 legend->AddEntry(cutHistoMean, "Mean", "l");
571 cutHistoMean->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
572 cutHistoMean->Draw();
573 DrawLines(cutHistoMean, nsigma, legend, kRed, pm);
574 } // if (cutHistoMean)
575
576 }
577 if (plotMedian) {
578 cutHistoMedian = AliTPCCalibViewer::SigmaCut(htemp, median, sigma, sigmaMax, sigmaStep, pm);
579 if (cutHistoMedian) {
580 fListOfObjectsToBeDeleted->Add(cutHistoMedian);
581 cutHistoMedian->SetLineColor(kBlue);
582 legend->AddEntry(cutHistoMedian, "Median", "l");
583 cutHistoMedian->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
584 if (plotMean && cutHistoMean) cutHistoMedian->Draw("same");
585 else cutHistoMedian->Draw();
586 DrawLines(cutHistoMedian, nsigma, legend, kBlue, pm);
587 } // if (cutHistoMedian)
588 }
589 if (plotLTM) {
590 Double_t ltmRms = 0;
591 Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
592 cutHistoLTM = AliTPCCalibViewer::SigmaCut(htemp, ltm, ltmRms, sigmaMax, sigmaStep, pm);
593 if (cutHistoLTM) {
594 fListOfObjectsToBeDeleted->Add(cutHistoLTM);
595 cutHistoLTM->SetLineColor(kGreen+2);
596 legend->AddEntry(cutHistoLTM, "LTM", "l");
597 cutHistoLTM->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
598 if (plotMean && cutHistoMean || plotMedian && cutHistoMedian) cutHistoLTM->Draw("same");
599 else cutHistoLTM->Draw();
600 DrawLines(cutHistoLTM, nsigma, legend, kGreen+2, pm);
601 }
602 }
603 if (!plotMean && !plotMedian && !plotLTM) return -1;
604 legend->Draw();
605 return 1;
606}
607
4f3934a1 608Int_t AliTPCCalibViewer::SigmaCutNew(const char* drawCommand, const char* sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, Bool_t pm, const char *sigmas, Float_t sigmaStep) const {
609 //
610 // Creates a histogram, where you can see, how much of the data are inside sigma-intervals
611 // around the mean/median/LTM
612 // with drawCommand, sector and cuts you specify your input data, see EasyDraw
613 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma)
614 // sigmaStep: the binsize of the generated histogram
615 // plotMean/plotMedian/plotLTM: specifies where to put the center
616 //
617
618 // Double_t ltmFraction = 0.8; //unused
619 // avoid compiler warnings:
620 sigmaMax = sigmaMax;
621 pm = pm;
622 sigmaStep = sigmaStep;
623
624 TString drawStr(drawCommand);
625 drawStr += " >> tempHist";
626
627 Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
628 TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
629 TGraph *cutGraphMean = 0;
630 // TGraph *cutGraphMedian = 0;
631 // TGraph *cutGraphLTM = 0;
632 Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
633 Int_t *index = new Int_t[entries];
634 Float_t *xarray = new Float_t[entries];
635 Float_t *yarray = new Float_t[entries];
636 TMath::Sort(entries, values, index, kFALSE);
637
638 Double_t mean = TMath::Mean(entries, values);
639 // Double_t median = TMath::Median(entries, values);
640 Double_t sigma = TMath::RMS(entries, values);
641
642 TLegend * legend = new TLegend(.7,.7, .99, .99, "Cumulative");
643 fListOfObjectsToBeDeleted->Add(legend);
644
645 // parse sigmas string
646 TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
647 TVectorF nsigma(sigmasTokens->GetEntriesFast());
648 for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
649 TString str(((TObjString*)sigmasTokens->At(i))->GetString());
650 Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
651 nsigma[i] = sig;
652 }
653
654 if (plotMean) {
655 for (Int_t i = 0; i < entries; i++) {
656 xarray[i] = TMath::Abs(values[index[i]] - mean) / sigma;
657 yarray[i] = float(i) / float(entries);
658 }
659 cutGraphMean = new TGraph(entries, xarray, yarray);
660 if (cutGraphMean) {
661 fListOfObjectsToBeDeleted->Add(cutGraphMean);
662 cutGraphMean->SetLineColor(kRed);
663 legend->AddEntry(cutGraphMean, "Mean", "l");
664 cutGraphMean->SetTitle(Form("%s, Cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
665 cutGraphMean->Draw("alu");
666 DrawLines(cutGraphMean, nsigma, legend, kRed, kTRUE);
667 }
668 }
669 /*
670 if (plotMedian) {
671 cutHistoMedian = AliTPCCalibViewer::SigmaCut(htemp, median, sigma, sigmaMax, sigmaStep, pm);
672 if (cutHistoMedian) {
673 fListOfObjectsToBeDeleted->Add(cutHistoMedian);
674 cutHistoMedian->SetLineColor(kBlue);
675 legend->AddEntry(cutHistoMedian, "Median", "l");
676 cutHistoMedian->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
677 if (plotMean && cutHistoMean) cutHistoMedian->Draw("same");
678 else cutHistoMedian->Draw();
679 DrawLines(cutHistoMedian, nsigma, legend, kBlue, pm);
680 } // if (cutHistoMedian)
681 }
682 if (plotLTM) {
683 Double_t ltmRms = 0;
684 Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
685 cutHistoLTM = AliTPCCalibViewer::SigmaCut(htemp, ltm, ltmRms, sigmaMax, sigmaStep, pm);
686 if (cutHistoLTM) {
687 fListOfObjectsToBeDeleted->Add(cutHistoLTM);
688 cutHistoLTM->SetLineColor(kGreen+2);
689 legend->AddEntry(cutHistoLTM, "LTM", "l");
690 cutHistoLTM->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
691 if (plotMean && cutHistoMean || plotMedian && cutHistoMedian) cutHistoLTM->Draw("same");
692 else cutHistoLTM->Draw();
693 DrawLines(cutHistoLTM, nsigma, legend, kGreen+2, pm);
694 }
695 }*/
696 if (!plotMean && !plotMedian && !plotLTM) return -1;
697 legend->Draw();
698 return 1;
699}
700
72d0ab7e 701
702
703
704Int_t AliTPCCalibViewer::Integrate(const char* drawCommand, Int_t sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, const char *sigmas, Float_t sigmaStep) const {
705 //
706 // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
707 // "mean" and "sigma" are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
708 // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
709 // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
710 // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
711 // The actual work is done on the array.
712 /* Begin_Latex
713 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx}{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
714 End_Latex
715 */
716 if (sector >= 0 && sector < 72) {
717 char sectorChr[3];
718 sprintf(sectorChr, "%i", sector);
719 return Integrate(drawCommand, sectorChr, cuts, sigmaMax, plotMean, plotMedian, plotLTM, sigmas, sigmaStep);
720 }
721 Error("Integrate","The TPC contains only sectors between 0 and 71.");
722 return -1;
723
724}
725
726
4f3934a1 727Int_t AliTPCCalibViewer::IntegrateOld(const char* drawCommand, const char* sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, const char *sigmas, Float_t sigmaStep) const {
72d0ab7e 728 //
729 // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
730 // "mean" and "sigma" are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
731 // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
732 // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
733 // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
734 // The actual work is done on the array.
735 /* Begin_Latex
736 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx}{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
737 End_Latex
738 */
739
740 Double_t ltmFraction = 0.8;
741
742 TString drawStr(drawCommand);
743 drawStr += " >> tempHist";
744
745 Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
746 TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
747 // FIXME is this histogram deleted automatically?
748 Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
749
750 Double_t mean = TMath::Mean(entries, values);
751 Double_t median = TMath::Median(entries, values);
752 Double_t sigma = TMath::RMS(entries, values);
753
754 TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
755 TVectorF nsigma(sigmasTokens->GetEntriesFast());
756 for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
757 TString str(((TObjString*)sigmasTokens->At(i))->GetString());
758 Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
759 nsigma[i] = sig;
760 }
761
762 TLegend * legend = new TLegend(.7,.7, .99, .99, "Integrated histogram");
763 fListOfObjectsToBeDeleted->Add(legend);
764 TH1F *integralHistoMean = 0;
765 TH1F *integralHistoMedian = 0;
766 TH1F *integralHistoLTM = 0;
767
768 if (plotMean) {
769 integralHistoMean = AliTPCCalibViewer::Integrate(htemp, mean, sigma, sigmaMax, sigmaStep);
770 if (integralHistoMean) {
771 fListOfObjectsToBeDeleted->Add(integralHistoMean);
772 integralHistoMean->SetLineColor(kRed);
773 legend->AddEntry(integralHistoMean, "Mean", "l");
774 integralHistoMean->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
775 integralHistoMean->Draw();
776 DrawLines(integralHistoMean, nsigma, legend, kRed, kTRUE);
39bcd65d 777 }
72d0ab7e 778 }
779 if (plotMedian) {
780 integralHistoMedian = AliTPCCalibViewer::Integrate(htemp, median, sigma, sigmaMax, sigmaStep);
781 if (integralHistoMedian) {
782 fListOfObjectsToBeDeleted->Add(integralHistoMedian);
783 integralHistoMedian->SetLineColor(kBlue);
784 legend->AddEntry(integralHistoMedian, "Median", "l");
785 integralHistoMedian->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
786 if (plotMean && integralHistoMean) integralHistoMedian->Draw("same");
787 else integralHistoMedian->Draw();
788 DrawLines(integralHistoMedian, nsigma, legend, kBlue, kTRUE);
39bcd65d 789 }
72d0ab7e 790 }
791 if (plotLTM) {
792 Double_t ltmRms = 0;
793 Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
794 integralHistoLTM = AliTPCCalibViewer::Integrate(htemp, ltm, ltmRms, sigmaMax, sigmaStep);
795 if (integralHistoLTM) {
796 fListOfObjectsToBeDeleted->Add(integralHistoLTM);
797 integralHistoLTM->SetLineColor(kGreen+2);
798 legend->AddEntry(integralHistoLTM, "LTM", "l");
799 integralHistoLTM->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
800 if (plotMean && integralHistoMean || plotMedian && integralHistoMedian) integralHistoLTM->Draw("same");
801 else integralHistoLTM->Draw();
802 DrawLines(integralHistoLTM, nsigma, legend, kGreen+2, kTRUE);
39bcd65d 803 }
804 }
72d0ab7e 805 if (!plotMean && !plotMedian && !plotLTM) return -1;
806 legend->Draw();
807 return 1;
808}
809
810
4f3934a1 811Int_t AliTPCCalibViewer::Integrate(const char* drawCommand, const char* sector, const char* cuts, Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, const char *sigmas, Float_t sigmaStep) const {
812 //
813 // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
814 // "mean" and "sigma" are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
815 // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
816 // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
817 // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
818 // The actual work is done on the array.
819 /* Begin_Latex
820 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx}{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
821 End_Latex
822 */
823
824 Double_t ltmFraction = 0.8;
825 // avoid compiler warnings:
826 sigmaMax = sigmaMax;
827 sigmaStep = sigmaStep;
828
829 TString drawStr(drawCommand);
830 drawStr += " >> tempHist";
831
832 Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
833 TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
834 TGraph *integralGraphMean = 0;
835 TGraph *integralGraphMedian = 0;
836 TGraph *integralGraphLTM = 0;
837 Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
838 Int_t *index = new Int_t[entries];
839 Float_t *xarray = new Float_t[entries];
840 Float_t *yarray = new Float_t[entries];
841 TMath::Sort(entries, values, index, kFALSE);
842
843 Double_t mean = TMath::Mean(entries, values);
844 Double_t median = TMath::Median(entries, values);
845 Double_t sigma = TMath::RMS(entries, values);
846
847 // parse sigmas string
848 TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
849 TVectorF nsigma(sigmasTokens->GetEntriesFast());
850 for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
851 TString str(((TObjString*)sigmasTokens->At(i))->GetString());
852 Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
853 nsigma[i] = sig;
854 }
855
856 TLegend * legend = new TLegend(.7,.7, .99, .99, "Integrated histogram");
857 fListOfObjectsToBeDeleted->Add(legend);
858
859 if (plotMean) {
860 for (Int_t i = 0; i < entries; i++) {
861 xarray[i] = (values[index[i]] - mean) / sigma;
862 yarray[i] = float(i) / float(entries);
863 }
864 integralGraphMean = new TGraph(entries, xarray, yarray);
865 if (integralGraphMean) {
866 fListOfObjectsToBeDeleted->Add(integralGraphMean);
867 integralGraphMean->SetLineColor(kRed);
868 legend->AddEntry(integralGraphMean, "Mean", "l");
869 integralGraphMean->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
870 integralGraphMean->Draw("alu");
871 DrawLines(integralGraphMean, nsigma, legend, kRed, kTRUE);
872 }
873 }
874 if (plotMedian) {
875 for (Int_t i = 0; i < entries; i++) {
876 xarray[i] = (values[index[i]] - median) / sigma;
877 yarray[i] = float(i) / float(entries);
878 }
879 integralGraphMedian = new TGraph(entries, xarray, yarray);
880 if (integralGraphMedian) {
881 fListOfObjectsToBeDeleted->Add(integralGraphMedian);
882 integralGraphMedian->SetLineColor(kBlue);
883 legend->AddEntry(integralGraphMedian, "Median", "l");
884 integralGraphMedian->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
885 if (plotMean && integralGraphMean) integralGraphMedian->Draw("samelu");
886 else integralGraphMedian->Draw("alu");
887 DrawLines(integralGraphMedian, nsigma, legend, kBlue, kTRUE);
888 }
889 }
890 if (plotLTM) {
891 Double_t ltmRms = 0;
892 Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
893 for (Int_t i = 0; i < entries; i++) {
894 xarray[i] = (values[index[i]] - ltm) / ltmRms;
895 yarray[i] = float(i) / float(entries);
896 }
897 integralGraphLTM = new TGraph(entries, xarray, yarray);
898 if (integralGraphLTM) {
899 fListOfObjectsToBeDeleted->Add(integralGraphLTM);
900 integralGraphLTM->SetLineColor(kGreen+2);
901 legend->AddEntry(integralGraphLTM, "LTM", "l");
902 integralGraphLTM->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
903 if (plotMean && integralGraphMean || plotMedian && integralGraphMedian) integralGraphLTM->Draw("samelu");
904 else integralGraphLTM->Draw("alu");
905 DrawLines(integralGraphLTM, nsigma, legend, kGreen+2, kTRUE);
906 }
907 }
908 if (!plotMean && !plotMedian && !plotLTM) return -1;
909 legend->Draw();
910 return entries;
911}
912
913
72d0ab7e 914void AliTPCCalibViewer::DrawLines(TH1F *histogram, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
915 //
916 // Private function for SigmaCut(...) and Integrate(...)
917 // Draws lines into the given histogram, specified by "nsigma", the lines are addeed to the legend
918 //
39bcd65d 919
72d0ab7e 920 // start to draw the lines, loop over requested sigmas
39bcd65d 921 char c[500];
72d0ab7e 922 for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
923 if (!pm) {
924 Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
925 TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
926 fListOfObjectsToBeDeleted->Add(lineUp);
927 lineUp->SetLineColor(color);
928 lineUp->SetLineStyle(2 + i);
929 lineUp->Draw();
930 TLine* lineLeft = new TLine(nsigma[i], histogram->GetBinContent(bin), 0, histogram->GetBinContent(bin));
931 fListOfObjectsToBeDeleted->Add(lineLeft);
932 lineLeft->SetLineColor(color);
933 lineLeft->SetLineStyle(2 + i);
934 lineLeft->Draw();
935 sprintf(c, "Fraction(%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin));
936 legend->AddEntry(lineLeft, c, "l");
937 }
938 else { // if (pm)
939 Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
940 TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
941 fListOfObjectsToBeDeleted->Add(lineUp1);
942 lineUp1->SetLineColor(color);
943 lineUp1->SetLineStyle(2 + i);
944 lineUp1->Draw();
945 TLine* lineLeft1 = new TLine(nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
946 fListOfObjectsToBeDeleted->Add(lineLeft1);
947 lineLeft1->SetLineColor(color);
948 lineLeft1->SetLineStyle(2 + i);
949 lineLeft1->Draw();
950 sprintf(c, "Fraction(+%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin));
951 legend->AddEntry(lineLeft1, c, "l");
952 bin = histogram->GetXaxis()->FindBin(-nsigma[i]);
953 TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], histogram->GetBinContent(bin));
954 fListOfObjectsToBeDeleted->Add(lineUp2);
955 lineUp2->SetLineColor(color);
956 lineUp2->SetLineStyle(2 + i);
957 lineUp2->Draw();
958 TLine* lineLeft2 = new TLine(-nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
959 fListOfObjectsToBeDeleted->Add(lineLeft2);
960 lineLeft2->SetLineColor(color);
961 lineLeft2->SetLineStyle(2 + i);
962 lineLeft2->Draw();
963 sprintf(c, "Fraction(-%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin));
964 legend->AddEntry(lineLeft2, c, "l");
965 }
966 } // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
967}
968
969
4f3934a1 970void AliTPCCalibViewer::DrawLines(TGraph *graph, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
971 //
972 // Private function for SigmaCut(...) and Integrate(...)
973 // Draws lines into the given histogram, specified by "nsigma", the lines are addeed to the legend
974 //
975
976 // start to draw the lines, loop over requested sigmas
977 char c[500];
978 for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
979 if (!pm) {
980 TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
981 fListOfObjectsToBeDeleted->Add(lineUp);
982 lineUp->SetLineColor(color);
983 lineUp->SetLineStyle(2 + i);
984 lineUp->Draw();
985 TLine* lineLeft = new TLine(nsigma[i], graph->Eval(nsigma[i]), 0, graph->Eval(nsigma[i]));
986 fListOfObjectsToBeDeleted->Add(lineLeft);
987 lineLeft->SetLineColor(color);
988 lineLeft->SetLineStyle(2 + i);
989 lineLeft->Draw();
990 sprintf(c, "Fraction(%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i]));
991 legend->AddEntry(lineLeft, c, "l");
992 }
993 else { // if (pm)
994 TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
995 fListOfObjectsToBeDeleted->Add(lineUp1);
996 lineUp1->SetLineColor(color);
997 lineUp1->SetLineStyle(2 + i);
998 lineUp1->Draw();
999 TLine* lineLeft1 = new TLine(nsigma[i], graph->Eval(nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(nsigma[i]));
1000 fListOfObjectsToBeDeleted->Add(lineLeft1);
1001 lineLeft1->SetLineColor(color);
1002 lineLeft1->SetLineStyle(2 + i);
1003 lineLeft1->Draw();
1004 sprintf(c, "Fraction(+%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i]));
1005 legend->AddEntry(lineLeft1, c, "l");
1006 TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], graph->Eval(-nsigma[i]));
1007 fListOfObjectsToBeDeleted->Add(lineUp2);
1008 lineUp2->SetLineColor(color);
1009 lineUp2->SetLineStyle(2 + i);
1010 lineUp2->Draw();
1011 TLine* lineLeft2 = new TLine(-nsigma[i], graph->Eval(-nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(-nsigma[i]));
1012 fListOfObjectsToBeDeleted->Add(lineLeft2);
1013 lineLeft2->SetLineColor(color);
1014 lineLeft2->SetLineStyle(2 + i);
1015 lineLeft2->Draw();
1016 sprintf(c, "Fraction(-%f #sigma) = %f",nsigma[i], graph->Eval(-nsigma[i]));
1017 legend->AddEntry(lineLeft2, c, "l");
1018 }
1019 } // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
1020}
1021
72d0ab7e 1022
1023
1024
1025
1026/////////////////
1027// Array tools //
1028/////////////////
1029
1030
1031Int_t AliTPCCalibViewer::GetBin(Float_t value, Int_t nbins, Double_t binLow, Double_t binUp){
1032 // Returns the 'bin' for 'value'
1033 // The interval between 'binLow' and 'binUp' is divided into 'nbins' equidistant bins
1034 // avoid index out of bounds error: 'if (bin < binLow) bin = binLow' and vice versa
1035 /* Begin_Latex
1036 GetBin(value) = #frac{nbins - 1}{binUp - binLow} #upoint (value - binLow) +1
1037 End_Latex
1038 */
1039
1040 Int_t bin = TMath::Nint( (Float_t)(value - binLow) / (Float_t)(binUp - binLow) * (nbins-1) ) + 1;
1041 // avoid index out of bounds:
1042 if (value < binLow) bin = 0;
1043 if (value > binUp) bin = nbins + 1;
1044 return bin;
39bcd65d 1045
72d0ab7e 1046}
39bcd65d 1047
72d0ab7e 1048
1049Double_t AliTPCCalibViewer::GetLTM(Int_t n, Double_t *array, Double_t *sigma, Double_t fraction){
1050 //
1051 // returns the LTM and sigma
1052 //
1053 Double_t *ddata = new Double_t[n];
1054 Double_t mean = 0, lsigma = 0;
1055 UInt_t nPoints = 0;
1056 for (UInt_t i = 0; i < (UInt_t)n; i++) {
1057 ddata[nPoints]= array[nPoints];
1058 nPoints++;
1059 }
1060 Int_t hh = TMath::Min(TMath::Nint(fraction * nPoints), Int_t(n));
1061 AliMathBase::EvaluateUni(nPoints, ddata, mean, lsigma, hh);
1062 if (sigma) *sigma = lsigma;
1063 delete [] ddata;
1064 return mean;
39bcd65d 1065}
1066
1067
72d0ab7e 1068TH1F* AliTPCCalibViewer::SigmaCut(TH1F *histogram, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep, Bool_t pm) {
1069 //
1070 // Creates a cumulative histogram Begin_Latex S(t, #mu, #sigma) End_Latex, where you can see, how much of the data are inside sigma-intervals around the mean value
1071 // The data of the distribution Begin_Latex f(x, #mu, #sigma) End_Latex are given in 'histogram'
1072 // 'mean' and 'sigma' are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in 'histogram', to be specified by the user
1073 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma, Begin_Latex t #sigma End_Latex)
1074 // sigmaStep: the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
1075 // pm: Decide weather Begin_Latex t > 0 End_Latex (first case) or Begin_Latex t End_Latex arbitrary (secound case)
1076 // The actual work is done on the array.
1077 /* Begin_Latex
1078 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{#mu}^{#mu + t #sigma} f(x, #mu, #sigma) dx + #int_{#mu}^{#mu - t #sigma} f(x, #mu, #sigma) dx }{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx } , for t > 0
1079 or
1080 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{#mu}^{#mu + t #sigma} f(x, #mu, #sigma) dx}{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
1081 End_Latex
1082 begin_macro(source)
1083 {
1084 Float_t mean = 0;
1085 Float_t sigma = 1.5;
1086 Float_t sigmaMax = 4;
1087 gROOT->SetStyle("Plain");
1088 TH1F *distribution = new TH1F("Distribution", "Distribution f(x, #mu, #sigma)", 1000,-5,5);
1089 TRandom rand(23);
1090 for (Int_t i = 0; i <50000;i++) distribution->Fill(rand.Gaus(mean, sigma));
1091 Float_t *ar = distribution->GetArray();
1092
1093 TCanvas* macro_example_canvas = new TCanvas("macro_example_canvas", "", 350, 350);
1094 macro_example_canvas->Divide(0,3);
1095 TVirtualPad *pad1 = macro_example_canvas->cd(1);
1096 pad1->SetGridy();
1097 pad1->SetGridx();
1098 distribution->Draw();
1099 TVirtualPad *pad2 = macro_example_canvas->cd(2);
1100 pad2->SetGridy();
1101 pad2->SetGridx();
1102
1103 TH1F *shist = AliTPCCalibViewer::SigmaCut(distribution, mean, sigma, sigmaMax);
1104 shist->SetNameTitle("Cumulative","Cumulative S(t, #mu, #sigma)");
1105 shist->Draw();
1106 TVirtualPad *pad3 = macro_example_canvas->cd(3);
1107 pad3->SetGridy();
1108 pad3->SetGridx();
1109 TH1F *shistPM = AliTPCCalibViewer::SigmaCut(distribution, mean, sigma, sigmaMax, -1, kTRUE);
1110 shistPM->Draw();
1111 return macro_example_canvas;
1112 }
1113 end_macro
1114 */
1115
1116 Float_t *array = histogram->GetArray();
1117 Int_t nbins = histogram->GetXaxis()->GetNbins();
1118 Float_t binLow = histogram->GetXaxis()->GetXmin();
1119 Float_t binUp = histogram->GetXaxis()->GetXmax();
1120 return AliTPCCalibViewer::SigmaCut(nbins, array, mean, sigma, nbins, binLow, binUp, sigmaMax, sigmaStep, pm);
1121}
1122
1123
1124
1125TH1F* AliTPCCalibViewer::SigmaCut(Int_t n, Float_t *array, Float_t mean, Float_t sigma, Int_t nbins, Float_t binLow, Float_t binUp, Float_t sigmaMax, Float_t sigmaStep, Bool_t pm){
1126 //
1127 // Creates a histogram Begin_Latex S(t, #mu, #sigma) End_Latex, where you can see, how much of the data are inside sigma-intervals around the mean value
1128 // The data of the distribution Begin_Latex f(x, #mu, #sigma) End_Latex are given in 'array', 'n' specifies the length of the array
1129 // 'mean' and 'sigma' are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in 'array', to be specified by the user
1130 // 'nbins': number of bins, 'binLow': first bin, 'binUp': last bin
1131 // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma, Begin_Latex t #sigma End_Latex)
1132 // sigmaStep: the binsize of the generated histogram
1133 // Here the actual work is done.
1134
1135 if (sigma == 0) return 0;
1136 Float_t binWidth = (binUp-binLow)/(nbins - 1);
1137 if (sigmaStep <= 0) sigmaStep = binWidth;
1138 Int_t kbins = (Int_t)(sigmaMax * sigma / sigmaStep) + 1; // + 1 due to overflow bin in histograms
1139 if (pm) kbins = 2 * (Int_t)(sigmaMax * sigma / sigmaStep) + 1;
1140 Float_t kbinLow = !pm ? 0 : -sigmaMax;
1141 Float_t kbinUp = sigmaMax;
1142 TH1F *hist = new TH1F("sigmaCutHisto","Cumulative; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
1143 hist->SetDirectory(0);
1144 hist->Reset();
1145
1146 // calculate normalization
1147 Double_t normalization = 0;
1148 for (Int_t i = 0; i <= n; i++) {
1149 normalization += array[i];
1150 }
1151
1152 // given units: units from given histogram
1153 // sigma units: in units of sigma
1154 // iDelta: integrate in interval (mean +- iDelta), given units
1155 // x: ofset from mean for integration, given units
1156 // hist: needs
1157
1158// printf("nbins: %i, binLow: %f, binUp: %f \n", nbins, binLow, binUp);
1159 // fill histogram
1160 for (Float_t iDelta = 0; iDelta <= sigmaMax * sigma; iDelta += sigmaStep) {
1161 // integrate array
1162 Double_t valueP = array[GetBin(mean, nbins, binLow, binUp)];
1163 Double_t valueM = array[GetBin(mean-binWidth, nbins, binLow, binUp)];
1164 // add bin of mean value only once to the histogram
1165// printf("++ adding bins: ");
1166 for (Float_t x = binWidth; x <= iDelta; x += binWidth) {
1167 valueP += (mean + x <= binUp) ? array[GetBin(mean + x, nbins, binLow, binUp)] : 0;
1168 valueM += (mean-binWidth - x >= binLow) ? array[GetBin(mean-binWidth - x, nbins, binLow, binUp)] : 0;
1169// printf("%i, ", GetBin(mean + x, nbins, binLow, binUp));
1170 }
1171// printf("\n");
1172 if (valueP / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", valueP, normalization);
1173 if (valueP / normalization > 100) return hist;
1174 if (valueM / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", valueM, normalization);
1175 if (valueM / normalization > 100) return hist;
1176 valueP = (valueP / normalization);
1177 valueM = (valueM / normalization);
1178 if (pm) {
1179 Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1180 hist->SetBinContent(bin, valueP);
1181 bin = GetBin(-iDelta/sigma, kbins, kbinLow, kbinUp);
1182 hist->SetBinContent(bin, valueM);
1183 }
1184 else { // if (!pm)
1185 Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1186 hist->SetBinContent(bin, valueP + valueM);
1187// printf(" first integration bin: %i, last integration bin in + direction: %i \n", GetBin(mean+binWidth, nbins, binLow, binUp), GetBin(iDelta, nbins, binLow, binUp));
1188// printf(" first integration bin: %i, last integration bin in - direction: %i \n", GetBin(mean+binWidth, nbins, binLow, binUp), GetBin(-iDelta, nbins, binLow, binUp));
1189// printf(" value: %f, normalization: %f, iDelta: %f, Bin: %i \n", valueP+valueM, normalization, iDelta, bin);
1190 }
1191 }
1192 //hist->SetMaximum(0.7);
1193 if (!pm) hist->SetMaximum(1.2);
1194 return hist;
1195}
1196
1197
1198TH1F* AliTPCCalibViewer::SigmaCut(Int_t n, Double_t *array, Double_t mean, Double_t sigma, Int_t nbins, Double_t *xbins, Double_t sigmaMax){
1199 //
1200 // SigmaCut for variable binsize
1201 // NOT YET IMPLEMENTED !!!
1202 //
1203 printf("SigmaCut with variable binsize, Not yet implemented\n");
1204 // avoid compiler warnings:
1205 n=n;
1206 mean=mean;
1207 sigma=sigma;
1208 nbins=nbins;
1209 sigmaMax=sigmaMax;
1210 array=array;
1211 xbins=xbins;
1212
1213 return 0;
1214}
1215
1216
1217TH1F* AliTPCCalibViewer::Integrate(TH1F *histogram, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
1218 //
1219 // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
1220 // "mean" and "sigma" are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
1221 // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
1222 // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
1223 // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
1224 // The actual work is done on the array.
1225 /* Begin_Latex
1226 f(x, #mu, #sigma) #Rightarrow S(t, #mu, #sigma) = #frac{#int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx}{ #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx }
1227 End_Latex
1228 begin_macro(source)
1229 {
1230 Float_t mean = 0;
1231 Float_t sigma = 1.5;
1232 Float_t sigmaMax = 4;
1233 gROOT->SetStyle("Plain");
1234 TH1F *distribution = new TH1F("Distribution", "Distribution f(x, #mu, #sigma)", 1000,-5,5);
1235 TRandom rand(23);
1236 for (Int_t i = 0; i <50000;i++) distribution->Fill(rand.Gaus(mean, sigma));
1237 Float_t *ar = distribution->GetArray();
1238
1239 TCanvas* macro_example_canvas = new TCanvas("macro_example_canvas", "", 350, 350);
1240 macro_example_canvas->Divide(0,2);
1241 TVirtualPad *pad1 = macro_example_canvas->cd(1);
1242 pad1->SetGridy();
1243 pad1->SetGridx();
1244 distribution->Draw();
1245 TVirtualPad *pad2 = macro_example_canvas->cd(2);
1246 pad2->SetGridy();
1247 pad2->SetGridx();
1248 TH1F *shist = AliTPCCalibViewer::Integrate(distribution, mean, sigma, sigmaMax);
1249 shist->SetNameTitle("Cumulative","Cumulative S(t, #mu, #sigma)");
1250 shist->Draw();
1251
1252 return macro_example_canvas;
1253 }
1254 end_macro
1255 */
1256
1257
1258 Float_t *array = histogram->GetArray();
1259 Int_t nbins = histogram->GetXaxis()->GetNbins();
1260 Float_t binLow = histogram->GetXaxis()->GetXmin();
1261 Float_t binUp = histogram->GetXaxis()->GetXmax();
1262 return AliTPCCalibViewer::Integrate(nbins, array, nbins, binLow, binUp, mean, sigma, sigmaMax, sigmaStep);
1263}
1264
1265
1266TH1F* AliTPCCalibViewer::Integrate(Int_t n, Float_t *array, Int_t nbins, Float_t binLow, Float_t binUp, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
1267 // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
1268 // "mean" and "sigma" are Begin_Latex #mu End_Latex and Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
1269 // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
1270 // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
1271 // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
1272 // Here the actual work is done.
1273
1274 Bool_t givenUnits = kTRUE;
1275 if (sigma != 0 && sigmaMax != 0) givenUnits = kFALSE;
1276 if (givenUnits) {
1277 sigma = 1;
1278 sigmaMax = (binUp - binLow) / 2.;
1279 }
1280
1281 Float_t binWidth = (binUp-binLow)/(nbins - 1);
1282 if (sigmaStep <= 0) sigmaStep = binWidth;
1283 Int_t kbins = (Int_t)(sigmaMax * sigma / sigmaStep) + 1; // + 1 due to overflow bin in histograms
1284 Float_t kbinLow = givenUnits ? binLow : -sigmaMax;
1285 Float_t kbinUp = givenUnits ? binUp : sigmaMax;
1286 TH1F *hist = 0;
1287 if (givenUnits) hist = new TH1F("integratedHisto","Integrated Histogram; Given x; Fraction of included data", kbins, kbinLow, kbinUp);
1288 if (!givenUnits) hist = new TH1F("integratedHisto","Integrated Histogram; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
1289 hist->SetDirectory(0);
1290 hist->Reset();
1291
1292 // calculate normalization
1293 // printf("calculating normalization, integrating from bin 1 to %i \n", n);
1294 Double_t normalization = 0;
1295 for (Int_t i = 1; i <= n; i++) {
1296 normalization += array[i];
1297 }
1298 // printf("normalization: %f \n", normalization);
1299
1300 // given units: units from given histogram
1301 // sigma units: in units of sigma
1302 // iDelta: integrate in interval (mean +- iDelta), given units
1303 // x: ofset from mean for integration, given units
1304 // hist: needs
1305
1306 // fill histogram
1307 for (Float_t iDelta = mean - sigmaMax * sigma; iDelta <= mean + sigmaMax * sigma; iDelta += sigmaStep) {
1308 // integrate array
1309 Double_t value = 0;
1310 for (Float_t x = mean - sigmaMax * sigma; x <= iDelta; x += binWidth) {
1311 value += (x <= binUp && x >= binLow) ? array[GetBin(x, nbins, binLow, binUp)] : 0;
1312 }
1313 if (value / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", value, normalization);
1314 if (value / normalization > 100) return hist;
1315 Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1316 // printf("first integration bin: %i, last integration bin: %i \n", GetBin(mean - sigmaMax * sigma, nbins, binLow, binUp), GetBin(iDelta, nbins, binLow, binUp));
1317 // printf("value: %f, normalization: %f, normalized value: %f, iDelta: %f, Bin: %i \n", value, normalization, value/normalization, iDelta, bin);
1318 value = (value / normalization);
1319 hist->SetBinContent(bin, value);
1320 }
1321 return hist;
1322}
1323
1324
1325
1326
1327
1328////////////////////////
1329// end of Array tools //
1330////////////////////////
1331
1332
1333
39bcd65d 1334//_____________________________________________________________________________
1335AliTPCCalPad* AliTPCCalibViewer::GetCalPad(const char* desiredData, char* cuts, char* calPadName) const {
1336 //
1337 // creates a AliTPCCalPad out of the 'desiredData'
1338 // the functionality of EasyDraw1D is used
1339 // calPadName specifies the name of the created AliTPCCalPad
1340 // - this takes a while -
1341 //
1342 TString drawStr(desiredData);
1343 drawStr.Append(":channel~");
1344 AliTPCCalPad * createdCalPad = new AliTPCCalPad(calPadName, calPadName);
1345 Int_t entries = 0;
1346 for (Int_t sec = 0; sec < 72; sec++) {
1347 entries = EasyDraw1D(drawStr.Data(), (Int_t)sec, cuts, "goff");
a6d2bd0c 1348 if (entries == -1) return 0;
39bcd65d 1349 for (Int_t i = 0; i < entries; i++)
1350 createdCalPad->GetCalROC(sec)->SetValue((UInt_t)(fTree->GetV2()[i]), (Float_t)(fTree->GetV1()[i]));
1351 }
1352 return createdCalPad;
1353}
1354
1355//_____________________________________________________________________________
1356AliTPCCalROC* AliTPCCalibViewer::GetCalROC(const char* desiredData, UInt_t sector, char* cuts) const {
1357 //
1358 // creates a AliTPCCalROC out of the desiredData
1359 // the functionality of EasyDraw1D is used
1360 // sector specifies the sector of the created AliTPCCalROC
1361 //
1362 TString drawStr(desiredData);
1363 drawStr.Append(":channel~");
1364 Int_t entries = EasyDraw1D(drawStr.Data(), (Int_t)sector, cuts, "goff");
a6d2bd0c 1365 if (entries == -1) return 0;
39bcd65d 1366 AliTPCCalROC * createdROC = new AliTPCCalROC(sector);
1367 for (Int_t i = 0; i < entries; i++)
1368 createdROC->SetValue((UInt_t)(fTree->GetV2()[i]), fTree->GetV1()[i]);
1369 return createdROC;
1370}
1371
1372
1373TObjArray* AliTPCCalibViewer::GetListOfVariables(Bool_t printList) {
1374 //
1375 // scan the tree - produces a list of available variables in the tree
1376 // printList: print the list to the screen, after the scan is done
1377 //
1378 TObjArray* arr = new TObjArray();
1379 TObjString* str = 0;
1380 Int_t nentries = fTree->GetListOfBranches()->GetEntries();
1381 for (Int_t i = 0; i < nentries; i++) {
1382 str = new TObjString(fTree->GetListOfBranches()->At(i)->GetName());
1383 str->String().ReplaceAll("_Median", "");
1384 str->String().ReplaceAll("_Mean", "");
1385 str->String().ReplaceAll("_RMS", "");
1386 str->String().ReplaceAll("_LTM", "");
1387 str->String().ReplaceAll("_OutlierCutted", "");
1388 str->String().ReplaceAll(".", "");
1389 if (!arr->FindObject(str) &&
1390 !(str->String() == "channel" || str->String() == "gx" || str->String() == "gy" ||
1391 str->String() == "lx" || str->String() == "ly" || str->String() == "pad" ||
1392 str->String() == "row" || str->String() == "rpad" || str->String() == "sector" ))
1393 arr->Add(str);
1394 }
1395 arr->Sort();
1396
1397 if (printList) {
1398 TIterator* iter = arr->MakeIterator();
1399 iter->Reset();
1400 TObjString* currentStr = 0;
1401 while ( (currentStr = (TObjString*)(iter->Next())) ) {
1402 std::cout << currentStr->GetString().Data() << std::endl;
1403 }
1404 delete iter;
1405 }
1406 return arr;
1407}
1408
1409
72d0ab7e 1410TObjArray* AliTPCCalibViewer::GetListOfNormalizationVariables(Bool_t printList) const{
39bcd65d 1411 //
1412 // produces a list of available variables for normalization in the tree
1413 // printList: print the list to the screen, after the scan is done
1414 //
1415 TObjArray* arr = new TObjArray();
1416 arr->Add(new TObjString("_Mean"));
1417 arr->Add(new TObjString("_Mean_OutlierCutted"));
1418 arr->Add(new TObjString("_Median"));
1419 arr->Add(new TObjString("_Median_OutlierCutted"));
1420 arr->Add(new TObjString("_LTM"));
1421 arr->Add(new TObjString("_LTM_OutlierCutted"));
1422 arr->Add(new TObjString("LFitIntern_4_8.fElements"));
1423 arr->Add(new TObjString("GFitIntern_Lin.fElements"));
1424 arr->Add(new TObjString("GFitIntern_Par.fElements"));
a6d2bd0c 1425 arr->Add(new TObjString("FitLinLocal"));
1426 arr->Add(new TObjString("FitLinGlobal"));
72d0ab7e 1427 arr->Add(new TObjString("FitParLocal"));
1428 arr->Add(new TObjString("FitParGlobal"));
39bcd65d 1429
1430 if (printList) {
1431 TIterator* iter = arr->MakeIterator();
1432 iter->Reset();
1433 TObjString* currentStr = 0;
1434 while ((currentStr = (TObjString*)(iter->Next()))) {
1435 std::cout << currentStr->GetString().Data() << std::endl;
1436 }
1437 delete iter;
1438 }
1439 return arr;
1440}
1441
1442
1443TFriendElement* AliTPCCalibViewer::AddReferenceTree(const char* filename, const char* treename, const char* refname){
1444 //
1445 // add a reference tree to the current tree
1446 // by default the treename is 'calPads' and the reference treename is 'R'
1447 //
1448 TFile *file = new TFile(filename);
1449 fListOfObjectsToBeDeleted->Add(file);
1450 TTree * tree = (TTree*)file->Get(treename);
1451 return AddFriend(tree, refname);
1452}
1453
1454
1455TObjArray* AliTPCCalibViewer::GetArrayOfCalPads(){
1456 //
1457 // Returns a TObjArray with all AliTPCCalPads that are stored in the tree
1458 // - this takes a while -
1459 //
1460 TObjArray *listOfCalPads = GetListOfVariables();
1461 TObjArray *calPadsArray = new TObjArray();
1462 Int_t numberOfCalPads = listOfCalPads->GetEntries();
1463 for (Int_t i = 0; i < numberOfCalPads; i++) {
1464 std::cout << "Creating calPad " << (i+1) << " of " << numberOfCalPads << "\r" << std::flush;
1465 char* calPadName = (char*)((TObjString*)(listOfCalPads->At(i)))->GetString().Data();
1466 TString drawCommand = ((TObjString*)(listOfCalPads->At(i)))->GetString();
1467 drawCommand.Append("~");
1468 AliTPCCalPad* calPad = GetCalPad(drawCommand.Data(), "", calPadName);
1469 calPadsArray->Add(calPad);
1470 }
1471 std::cout << std::endl;
1472 listOfCalPads->Delete();
1473 delete listOfCalPads;
1474 return calPadsArray;
1475}
1476
1477
a6d2bd0c 1478TString* AliTPCCalibViewer::Fit(const char* drawCommand, const char* formula, const char* cuts, Double_t & chi2, TVectorD &fitParam, TMatrixD &covMatrix){
1479 //
1480 // fit an arbitrary function, specified by formula into the data, specified by drawCommand and cuts
1481 // returns chi2, fitParam and covMatrix
1482 // returns TString with fitted formula
1483 //
1484
1485 TString formulaStr(formula);
1486 TString drawStr(drawCommand);
1487 TString cutStr(cuts);
1488
1489 // abbreviations:
1490 drawStr.ReplaceAll("~",".fElements");
1491 cutStr.ReplaceAll("~",".fElements");
1492 formulaStr.ReplaceAll("~", ".fElements");
1493
1494 formulaStr.ReplaceAll("++", "~");
1495 TObjArray* formulaTokens = formulaStr.Tokenize("~");
1496 Int_t dim = formulaTokens->GetEntriesFast();
1497
1498 fitParam.ResizeTo(dim);
1499 covMatrix.ResizeTo(dim,dim);
1500
1501 TLinearFitter* fitter = new TLinearFitter(dim+1, Form("hyp%d",dim));
1502 fitter->StoreData(kTRUE);
1503 fitter->ClearPoints();
1504
1505 Int_t entries = Draw(drawStr.Data(), cutStr.Data(), "goff");
1506 if (entries == -1) return new TString("An ERROR has occured during fitting!");
1507 Double_t **values = new Double_t*[dim+1] ;
1508
1509 for (Int_t i = 0; i < dim + 1; i++){
1510 Int_t centries = 0;
1511 if (i < dim) centries = fTree->Draw(((TObjString*)formulaTokens->At(i))->GetName(), cutStr.Data(), "goff");
1512 else centries = fTree->Draw(drawStr.Data(), cutStr.Data(), "goff");
1513
1514 if (entries != centries) return new TString("An ERROR has occured during fitting!");
1515 values[i] = new Double_t[entries];
1516 memcpy(values[i], fTree->GetV1(), entries*sizeof(Double_t));
1517 }
1518
1519 // add points to the fitter
1520 for (Int_t i = 0; i < entries; i++){
1521 Double_t x[1000];
1522 for (Int_t j=0; j<dim;j++) x[j]=values[j][i];
1523 fitter->AddPoint(x, values[dim][i], 1);
1524 }
1525
1526 fitter->Eval();
1527 fitter->GetParameters(fitParam);
1528 fitter->GetCovarianceMatrix(covMatrix);
1529 chi2 = fitter->GetChisquare();
1530 chi2 = chi2;
1531
1532 TString *preturnFormula = new TString(Form("( %f+",fitParam[0])), &returnFormula = *preturnFormula;
1533
1534 for (Int_t iparam = 0; iparam < dim; iparam++) {
1535 returnFormula.Append(Form("%s*(%f)",((TObjString*)formulaTokens->At(iparam))->GetName(),fitParam[iparam+1]));
1536 if (iparam < dim-1) returnFormula.Append("+");
1537 }
1538 returnFormula.Append(" )");
1539 delete formulaTokens;
1540 delete fitter;
1541 delete[] values;
1542 return preturnFormula;
1543}
1544
1545
39bcd65d 1546void AliTPCCalibViewer::MakeTreeWithObjects(const char * fileName, TObjArray * array, const char * mapFileName) {
1547 //
1548 // Write tree with all available information
1549 // im mapFileName is speciefied, the Map information are also written to the tree
1550 // AliTPCCalPad-Objects are written directly to the tree, so that they can be accessd later on
1551 // (does not work!!!)
1552 //
1553 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
1554
1555 TObjArray* mapIROCs = 0;
1556 TObjArray* mapOROCs = 0;
1557 TVectorF *mapIROCArray = 0;
1558 TVectorF *mapOROCArray = 0;
1559 Int_t mapEntries = 0;
1560 TString* mapNames = 0;
1561
1562 if (mapFileName) {
1563 TFile mapFile(mapFileName, "read");
1564
1565 TList* listOfROCs = mapFile.GetListOfKeys();
1566 mapEntries = listOfROCs->GetEntries()/2;
1567 mapIROCs = new TObjArray(mapEntries*2);
1568 mapOROCs = new TObjArray(mapEntries*2);
1569 mapIROCArray = new TVectorF[mapEntries];
1570 mapOROCArray = new TVectorF[mapEntries];
1571
1572 mapNames = new TString[mapEntries];
1573 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
72d0ab7e 1574 TString rocName(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
1575 rocName.Remove(rocName.Length()-4, 4);
1576 mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "IROC").Data()), ivalue);
1577 mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "OROC").Data()), ivalue);
1578 mapNames[ivalue].Append(rocName);
39bcd65d 1579 }
1580
1581 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1582 mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
1583 mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
1584
1585 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
1586 (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
1587 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
1588 (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
1589 }
1590
1591 } // if (mapFileName)
1592
1593 TTreeSRedirector cstream(fileName);
1594 Int_t arrayEntries = array->GetEntries();
1595
1596 // Read names of AliTPCCalPads and save them in names[]
1597 TString* names = new TString[arrayEntries];
1598 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1599 names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
1600
1601 for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
1602
1603 TVectorF *vectorArray = new TVectorF[arrayEntries];
1604 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1605 vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1606
1607
1608 //
1609 // fill vectors of variable per pad
1610 //
1611 TVectorF *posArray = new TVectorF[8];
1612 for (Int_t ivalue = 0; ivalue < 8; ivalue++)
1613 posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1614
1615 Float_t posG[3] = {0};
1616 Float_t posL[3] = {0};
1617 Int_t ichannel = 0;
1618 for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
1619 for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
1620 tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
1621 tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
1622 posArray[0][ichannel] = irow;
1623 posArray[1][ichannel] = ipad;
1624 posArray[2][ichannel] = posL[0];
1625 posArray[3][ichannel] = posL[1];
1626 posArray[4][ichannel] = posG[0];
1627 posArray[5][ichannel] = posG[1];
1628 posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
1629 posArray[7][ichannel] = ichannel;
1630
1631 // loop over array containing AliTPCCalPads
1632 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1633 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1634 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1635 if (calROC)
1636 (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
1637 else
1638 (vectorArray[ivalue])[ichannel] = 0;
1639 }
1640 ichannel++;
1641 }
1642 }
1643 AliTPCCalROC dummyROC(0);
1644 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1645 AliTPCCalROC *roc = ((AliTPCCalPad*)array->At(ivalue))->GetCalROC(isector);
1646 if (!roc) roc = &dummyROC;
1647 cstream << "calPads" <<
1648 (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
1649 cstream << "calPads" <<
1650 (Char_t*)((names[ivalue] + "Pad.=").Data()) << roc;
1651 }
1652
1653 if (mapFileName) {
1654 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1655 if (isector < 36)
1656 cstream << "calPads" <<
1657 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
1658 else
1659 cstream << "calPads" <<
1660 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
1661 }
1662 }
1663
1664 cstream << "calPads" <<
1665 "sector=" << isector;
1666
1667 cstream << "calPads" <<
1668 "row.=" << &posArray[0] <<
1669 "pad.=" << &posArray[1] <<
1670 "lx.=" << &posArray[2] <<
1671 "ly.=" << &posArray[3] <<
1672 "gx.=" << &posArray[4] <<
1673 "gy.=" << &posArray[5] <<
1674 "rpad.=" << &posArray[6] <<
1675 "channel.=" << &posArray[7];
1676
1677 cstream << "calPads" <<
1678 "\n";
1679
1680 delete[] posArray;
1681 delete[] vectorArray;
1682 } //for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++)
1683
1684 delete[] names;
1685 if (mapFileName) {
1686 delete mapIROCs;
1687 delete mapOROCs;
1688 delete[] mapIROCArray;
1689 delete[] mapOROCArray;
1690 delete[] mapNames;
1691 }
1692}
1693
72d0ab7e 1694
39bcd65d 1695void AliTPCCalibViewer::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad* outlierPad, Float_t ltmFraction) {
1696 //
1697 // Write a tree with all available information
72d0ab7e 1698 // if mapFileName is speciefied, the Map information are also written to the tree
39bcd65d 1699 // pads specified in outlierPad are not used for calculating statistics
72d0ab7e 1700 // The following statistical information on the basis of a ROC are calculated:
1701 // "_Median", "_Mean", "_LTM", "_RMS_LTM"
1702 // "_Median_OutlierCutted", "_Mean_OutlierCutted", "_RMS_OutlierCutted", "_LTM_OutlierCutted", "_RMS_LTM_OutlierCutted"
1703 // The following position variables are available:
1704 // "row", "pad", "lx", "ly", "gx", "gy", "rpad", "channel"
1705 //
1706 // The tree out of this function is the basis for the AliTPCCalibViewer and the AliTPCCalibViewerGUI.
1707
39bcd65d 1708 AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
1709
1710 TObjArray* mapIROCs = 0;
1711 TObjArray* mapOROCs = 0;
1712 TVectorF *mapIROCArray = 0;
1713 TVectorF *mapOROCArray = 0;
1714 Int_t mapEntries = 0;
1715 TString* mapNames = 0;
1716
1717 if (mapFileName) {
1718 TFile mapFile(mapFileName, "read");
1719
1720 TList* listOfROCs = mapFile.GetListOfKeys();
1721 mapEntries = listOfROCs->GetEntries()/2;
1722 mapIROCs = new TObjArray(mapEntries*2);
1723 mapOROCs = new TObjArray(mapEntries*2);
1724 mapIROCArray = new TVectorF[mapEntries];
1725 mapOROCArray = new TVectorF[mapEntries];
1726
1727 mapNames = new TString[mapEntries];
1728 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
72d0ab7e 1729 TString rocName(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
1730 rocName.Remove(rocName.Length()-4, 4);
1731 mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "IROC").Data()), ivalue);
1732 mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "OROC").Data()), ivalue);
1733 mapNames[ivalue].Append(rocName);
39bcd65d 1734 }
1735
1736 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1737 mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
1738 mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
1739
1740 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
1741 (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
1742 for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
1743 (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
1744 }
1745
1746 } // if (mapFileName)
1747
1748 TTreeSRedirector cstream(fileName);
72d0ab7e 1749 Int_t arrayEntries = 0;
1750 if (array) arrayEntries = array->GetEntries();
39bcd65d 1751
1752 TString* names = new TString[arrayEntries];
1753 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1754 names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
1755
1756 for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
1757 //
1758 // get statistic for given sector
1759 //
1760 TVectorF median(arrayEntries);
1761 TVectorF mean(arrayEntries);
1762 TVectorF rms(arrayEntries);
1763 TVectorF ltm(arrayEntries);
1764 TVectorF ltmrms(arrayEntries);
1765 TVectorF medianWithOut(arrayEntries);
1766 TVectorF meanWithOut(arrayEntries);
1767 TVectorF rmsWithOut(arrayEntries);
1768 TVectorF ltmWithOut(arrayEntries);
1769 TVectorF ltmrmsWithOut(arrayEntries);
1770
1771 TVectorF *vectorArray = new TVectorF[arrayEntries];
1772 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1773 vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1774
1775 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1776 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1777 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1778 AliTPCCalROC* outlierROC = 0;
1779 if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
1780 if (calROC) {
1781 median[ivalue] = calROC->GetMedian();
1782 mean[ivalue] = calROC->GetMean();
1783 rms[ivalue] = calROC->GetRMS();
1784 Double_t ltmrmsValue = 0;
1785 ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
1786 ltmrms[ivalue] = ltmrmsValue;
1787 if (outlierROC) {
1788 medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
1789 meanWithOut[ivalue] = calROC->GetMean(outlierROC);
1790 rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
1791 ltmrmsValue = 0;
1792 ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
1793 ltmrmsWithOut[ivalue] = ltmrmsValue;
1794 }
1795 }
1796 else {
1797 median[ivalue] = 0.;
1798 mean[ivalue] = 0.;
1799 rms[ivalue] = 0.;
1800 ltm[ivalue] = 0.;
1801 ltmrms[ivalue] = 0.;
1802 medianWithOut[ivalue] = 0.;
1803 meanWithOut[ivalue] = 0.;
1804 rmsWithOut[ivalue] = 0.;
1805 ltmWithOut[ivalue] = 0.;
1806 ltmrmsWithOut[ivalue] = 0.;
1807 }
1808 }
1809
1810 //
1811 // fill vectors of variable per pad
1812 //
1813 TVectorF *posArray = new TVectorF[8];
1814 for (Int_t ivalue = 0; ivalue < 8; ivalue++)
1815 posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1816
1817 Float_t posG[3] = {0};
1818 Float_t posL[3] = {0};
1819 Int_t ichannel = 0;
1820 for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
1821 for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
1822 tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
1823 tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
1824 posArray[0][ichannel] = irow;
1825 posArray[1][ichannel] = ipad;
1826 posArray[2][ichannel] = posL[0];
1827 posArray[3][ichannel] = posL[1];
1828 posArray[4][ichannel] = posG[0];
1829 posArray[5][ichannel] = posG[1];
1830 posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
1831 posArray[7][ichannel] = ichannel;
1832
1833 // loop over array containing AliTPCCalPads
1834 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1835 AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1836 AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1837 if (calROC)
1838 (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
1839 else
1840 (vectorArray[ivalue])[ichannel] = 0;
1841 }
1842 ichannel++;
1843 }
1844 }
1845
1846 cstream << "calPads" <<
1847 "sector=" << isector;
1848
1849 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1850 cstream << "calPads" <<
1851 (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
1852 (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
1853 (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
1854 (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
1855 (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
1856 if (outlierPad) {
1857 cstream << "calPads" <<
1858 (Char_t*)((names[ivalue] + "_Median_OutlierCutted=").Data()) << medianWithOut[ivalue] <<
1859 (Char_t*)((names[ivalue] + "_Mean_OutlierCutted=").Data()) << meanWithOut[ivalue] <<
1860 (Char_t*)((names[ivalue] + "_RMS_OutlierCutted=").Data()) << rmsWithOut[ivalue] <<
1861 (Char_t*)((names[ivalue] + "_LTM_OutlierCutted=").Data()) << ltmWithOut[ivalue] <<
1862 (Char_t*)((names[ivalue] + "_RMS_LTM_OutlierCutted=").Data()) << ltmrmsWithOut[ivalue];
1863 }
1864 }
1865
1866 for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1867 cstream << "calPads" <<
1868 (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
1869 }
1870
1871 if (mapFileName) {
1872 for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1873 if (isector < 36)
1874 cstream << "calPads" <<
1875 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
1876 else
1877 cstream << "calPads" <<
1878 (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
1879 }
1880 }
1881
1882 cstream << "calPads" <<
1883 "row.=" << &posArray[0] <<
1884 "pad.=" << &posArray[1] <<
1885 "lx.=" << &posArray[2] <<
1886 "ly.=" << &posArray[3] <<
1887 "gx.=" << &posArray[4] <<
1888 "gy.=" << &posArray[5] <<
1889 "rpad.=" << &posArray[6] <<
1890 "channel.=" << &posArray[7];
1891
1892 cstream << "calPads" <<
1893 "\n";
1894
1895 delete[] posArray;
1896 delete[] vectorArray;
1897 }
1898
1899
1900 delete[] names;
1901 if (mapFileName) {
1902 delete mapIROCs;
1903 delete mapOROCs;
1904 delete[] mapIROCArray;
1905 delete[] mapOROCArray;
1906 delete[] mapNames;
1907 }
1908}
1909
4f3934a1 1910
1911void AliTPCCalibViewer::MakeTree(const char *outPutFileName, const Char_t *inputFileName, AliTPCCalPad *outlierPad, Float_t ltmFraction, const char *mapFileName ){
1912 //
1913 // Function to create a calibration Tree with all available information.
1914 // See also documentation to MakeTree
1915 // the file "inputFileName" must be in the following format (see also CreateObjectList):
1916 // (each colum separated by tabs, "dependingOnType" can have 2 or 3 colums)
1917 //
1918 // type path dependingOnType
1919 //
1920 // type == "CE":
1921 // dependingOnType = CETmean CEQmean CETrms
1922 //
1923 // type == "Pulser":
1924 // dependingOnType = PulserTmean PulsterQmean PulserTrms
1925 //
1926 // type == "Pedestals":
1927 // dependingOnType = Pedestals Noise
1928 //
1929 // type == "CalPad":
1930 // dependingOnType = NameInFile NameToWriteToFile
1931 //
1932 //
1933 TObjArray objArray;
1934 CreateObjectList(inputFileName, &objArray);
1935 MakeTree(outPutFileName, &objArray, mapFileName, outlierPad, ltmFraction);
1936}
1937
1938void AliTPCCalibViewer::CreateObjectList(const Char_t *filename, TObjArray *calibObjects){
1939 //
1940 // Function to create a TObjArray out of a given file
1941 // the file must be in the following format:
1942 // (each colum separated by tabs, "dependingOnType" can have 2 or 3 colums)
1943 //
1944 //
1945 // type path dependingOnType
1946 //
1947 // type == "CE":
1948 // dependingOnType = CETmean CEQmean CETrms
1949 //
1950 // type == "Pulser":
1951 // dependingOnType = PulserTmean PulsterQmean PulserTrms
1952 //
1953 // type == "Pedestals":
1954 // dependingOnType = Pedestals Noise
1955 //
1956 // type == "CalPad":
1957 // dependingOnType = NameInFile NameToWriteToFile
1958 //
1959 //
1960 //
1961 if ( calibObjects == 0x0 ) return;
1962 ifstream in;
1963 in.open(filename);
1964 if ( !in.is_open() ){
1965 fprintf(stderr,"Error: cannot open list file '%s'", filename);
1966 return;
1967 }
1968
1969 AliTPCCalPad *calPad=0x0;
1970
1971 TString sFile;
1972 sFile.ReadFile(in);
1973 in.close();
1974
1975 TObjArray *arrFileLine = sFile.Tokenize("\n");
1976 TIter nextLine(arrFileLine);
1977
1978 TObjString *sObjLine = 0x0;
1979 while ( (sObjLine = (TObjString*)nextLine()) ){
1980 TString sLine(sObjLine->GetString());
1981
1982 TObjArray *arrCol = sLine.Tokenize("\t");
1983 Int_t nCols = arrCol->GetEntriesFast();
1984
1985 TObjString *sObjType = (TObjString*)(arrCol->At(0));
1986 TObjString *sObjFileName = (TObjString*)(arrCol->At(1));
1987 TObjString *sObjName = 0x0;
1988
1989 if ( !sObjType || !sObjFileName ) continue;
1990 TString sType(sObjType->GetString());
1991 TString sFileName(sObjFileName->GetString());
1992 printf("Type %s, opening %s \n", sType.Data(), sFileName.Data());
1993 TFile *fIn = TFile::Open(sFileName);
1994 if ( !fIn ){
1995 fprintf(stderr,"File not found: '%s'", sFileName.Data());
1996 continue;
1997 }
1998
1999 if ( sType == "CE" ){ // next three colums are the names for CETmean, CEQmean and CETrms
2000 AliTPCCalibCE *ce = (AliTPCCalibCE*)fIn->Get("AliTPCCalibCE");
2001 calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadT0());
2002 if (nCols > 2) { // check, if the name is provided
2003 sObjName = (TObjString*)(arrCol->At(2));
2004 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2005 }
2006 else calPad->SetNameTitle("CETmean","CETmean");
2007 calibObjects->Add(calPad);
2008
2009 calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadQ());
2010 if (nCols > 3) { // check, if the name is provided
2011 sObjName = (TObjString*)(arrCol->At(3));
2012 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2013 }
2014 else calPad->SetNameTitle("CEQmean","CEQmean");
2015 calibObjects->Add(calPad);
2016
2017 calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadRMS());
2018 if (nCols > 4) { // check, if the name is provided
2019 sObjName = (TObjString*)(arrCol->At(4));
2020 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2021 }
2022 else calPad->SetNameTitle("CETrms","CETrms");
2023 calibObjects->Add(calPad);
2024
2025 } else if ( sType == "Pulser") {
2026 AliTPCCalibPulser *sig = (AliTPCCalibPulser*)fIn->Get("AliTPCCalibPulser");
2027
2028 calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadT0());
2029 if (nCols > 2) {
2030 sObjName = (TObjString*)(arrCol->At(2));
2031 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2032 }
2033 else calPad->SetNameTitle("PulserTmean","PulserTmean");
2034 calibObjects->Add(calPad);
2035
2036 calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadQ());
2037 if (nCols > 3) {
2038 sObjName = (TObjString*)(arrCol->At(3));
2039 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2040 }
2041 else calPad->SetNameTitle("PulserQmean","PulserQmean");
2042 calibObjects->Add(calPad);
2043
2044 calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadRMS());
2045 if (nCols > 4) {
2046 sObjName = (TObjString*)(arrCol->At(4));
2047 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2048 }
2049 else calPad->SetNameTitle("PulserTrms","PulserTrms");
2050 calibObjects->Add(calPad);
2051
2052 } else if ( sType == "Pedestals") {
2053 AliTPCCalibPedestal *ped = (AliTPCCalibPedestal*)fIn->Get("AliTPCCalibPedestal");
2054
2055 calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadPedestal());
2056 if (nCols > 2) {
2057 sObjName = (TObjString*)(arrCol->At(2));
2058 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2059 }
2060 else calPad->SetNameTitle("Pedestals","Pedestals");
2061 calibObjects->Add(calPad);
2062
2063 calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadRMS());
2064 if (nCols > 3) {
2065 sObjName = (TObjString*)(arrCol->At(3));
2066 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2067 }
2068 else calPad->SetNameTitle("Noise","Noise");
2069 calibObjects->Add(calPad);
2070
2071 } else if ( sType == "CalPad") {
2072 if (nCols > 2) sObjName = (TObjString*)(arrCol->At(2));
2073 else continue;
2074 calPad = new AliTPCCalPad(*(AliTPCCalPad*)fIn->Get(sObjName->GetString().Data()));
2075 if (nCols > 3) {
2076 sObjName = (TObjString*)(arrCol->At(3));
2077 calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2078 }
2079 calibObjects->Add(calPad);
2080 } else {
2081 fprintf(stderr,"Undefined Type: '%s'",sType.Data());
2082 }
2083 delete fIn;
2084 }
2085}
2086
2087
2088