AliRoot Core  v5-06-15 (45dab64)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliTPCCalibViewer.cxx
Go to the documentation of this file.
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 
35 
36 #include <fstream>
37 #include <iostream>
38 
39 #include <TFile.h>
40 #include <TFriendElement.h>
41 #include <TGraph.h>
42 #include <TKey.h>
43 #include <TPad.h>
44 //#include <TCanvas.h>
45 #include <TH1.h>
46 #include <TH1F.h>
47 #include <TLegend.h>
48 #include <TLine.h>
49 #include <TMath.h>
50 #include <TObjString.h>
51 //#include <TROOT.h>
52 #include <TRandom.h>
53 #include <TString.h>
54 #include <TStyle.h>
55 #include <TTreeStream.h>
56 
57 #include "AliTPCCalibCE.h"
58 #include "AliMathBase.h"
59 #include "AliTPCCalPad.h"
60 #include "AliTPCCalROC.h"
61 #include "AliTPCCalibPedestal.h"
62 #include "AliTPCCalibPulser.h"
63 
64 //
65 // AliRoot includes
66 //
67 #include "AliTPCCalibViewer.h"
68 
69 using std::ifstream;
73 
74 
76  :TObject(),
77  fTree(0),
78  fFile(0),
79  fListOfObjectsToBeDeleted(0),
80  fTreeMustBeDeleted(0),
81  fAbbreviation(0),
82  fAppendString(0)
83 {
84  //
85  // Default constructor
86  //
87 
88 }
89 
90 //_____________________________________________________________________________
92  :TObject(c),
93  fTree(0),
94  fFile(0),
95  fListOfObjectsToBeDeleted(0),
96  fTreeMustBeDeleted(0),
97  fAbbreviation(0),
98  fAppendString(0)
99 {
102 
103  fTree = c.fTree;
105  //fFile = new TFile(*(c.fFile));
109 }
110 
111 //_____________________________________________________________________________
113  :TObject(),
114  fTree(0),
115  fFile(0),
116  fListOfObjectsToBeDeleted(0),
117  fTreeMustBeDeleted(0),
118  fAbbreviation(0),
119  fAppendString(0)
120 {
122 
123  fTree = tree;
124  fTreeMustBeDeleted = kFALSE;
126  fAbbreviation = "~";
127  fAppendString = ".fElements";
128 }
129 
130 //_____________________________________________________________________________
131 AliTPCCalibViewer::AliTPCCalibViewer(const char* fileName, const char* treeName)
132  :TObject(),
133  fTree(0),
134  fFile(0),
135  fListOfObjectsToBeDeleted(0),
136  fTreeMustBeDeleted(0),
137  fAbbreviation(0),
138  fAppendString(0)
139 
140 {
143 
144  fFile = new TFile(fileName, "read");
145  fTree = (TTree*) fFile->Get(treeName);
146  fTreeMustBeDeleted = kTRUE;
148  fAbbreviation = "~";
149  fAppendString = ".fElements";
150 }
151 
152 //____________________________________________________________________________
154 {
157 
158  if (this == &param) return (*this);
159 
160  fTree = param.fTree;
162  //fFile = new TFile(*(param.fFile));
166  return (*this);
167 }
168 
169 //_____________________________________________________________________________
171 {
174 
175  if (fTree && fTreeMustBeDeleted) {
176  fTree->SetCacheSize(0);
177  fTree->Delete();
178  //delete fTree;
179  }
180  if (fFile) {
181  fFile->Close();
182  fFile = 0;
183  }
184 
185  for (Int_t i = fListOfObjectsToBeDeleted->GetEntriesFast()-1; i >= 0; i--) {
186  //cout << "Index " << i << " trying to delete the following object: " << fListOfObjectsToBeDeleted->At(i)->GetName() << "..."<< endl;
187  delete fListOfObjectsToBeDeleted->At(i);
188  }
190 }
191 
192 //_____________________________________________________________________________
193 void AliTPCCalibViewer::Delete(Option_t* /*option*/) {
197 
198  if (fTree && fTreeMustBeDeleted) {
199  fTree->SetCacheSize(0);
200  fTree->Delete();
201  }
202  delete fFile;
204 }
205 
206 
207 const char* AliTPCCalibViewer::AddAbbreviations(const Char_t *c, Bool_t printDrawCommand){
230 
231  TString str(c);
232  TString removeString = "!#"; // very unpropable combination of chars
233  TString replaceString = "";
234  TString searchString = "";
235  TString normString = "";
236  TObjArray *listOfVariables = GetListOfVariables();
237  listOfVariables->Add(new TObjString("channel"));
238  listOfVariables->Add(new TObjString("gx"));
239  listOfVariables->Add(new TObjString("gy"));
240  listOfVariables->Add(new TObjString("lx"));
241  listOfVariables->Add(new TObjString("ly"));
242  listOfVariables->Add(new TObjString("pad"));
243  listOfVariables->Add(new TObjString("row"));
244  listOfVariables->Add(new TObjString("rpad"));
245  listOfVariables->Add(new TObjString("sector"));
246  TObjArray *listOfNormalizationVariables = GetListOfNormalizationVariables();
247  Int_t nVariables = listOfVariables->GetEntriesFast();
248  Int_t nNorm = listOfNormalizationVariables->GetEntriesFast();
249 
250  Int_t *varLengths = new Int_t[nVariables];
251  for (Int_t i = 0; i < nVariables; i++) {
252  varLengths[i] = ((TObjString*)listOfVariables->At(i))->String().Length();
253  }
254  Int_t *normLengths = new Int_t[nNorm];
255  for (Int_t i = 0; i < nNorm; i++) {
256  normLengths[i] = ((TObjString*)listOfNormalizationVariables->At(i))->String().Length();
257  // printf("normLengths[%i] (%s) = %i \n", i,((TObjString*)listOfNormalizationVariables->At(i))->String().Data(), normLengths[i]);
258  }
259  Int_t *varSort = new Int_t[nVariables];
260  TMath::Sort(nVariables, varLengths, varSort, kTRUE);
261  Int_t *normSort = new Int_t[nNorm];
262  TMath::Sort(nNorm, normLengths, normSort, kTRUE);
263  // for (Int_t i = 0; i<nNorm; i++) printf("normLengths: %i\n", normLengths[normSort[i]]);
264  // for (Int_t i = 0; i<nVariables; i++) printf("varLengths: %i\n", varLengths[varSort[i]]);
265 
266  for (Int_t ivar = 0; ivar < nVariables; ivar++) {
267  // ***** save correct tokens *****
268  // first get the next variable:
269  searchString = ((TObjString*)listOfVariables->At(varSort[ivar]))->String();
270  // printf("searchString: %s ++++++++++++++\n", searchString.Data());
271  // form replaceString:
272  replaceString = "";
273  for (Int_t i = 0; i < searchString.Length(); i++) {
274  replaceString.Append(searchString[i]);
275  if (i == 0) replaceString.Append(removeString);
276  }
277  // go through normalization:
278  // printf("go through normalization\n");
279  for (Int_t inorm = 0; inorm < nNorm; inorm++) {
280  // printf(" inorm=%i, nNorm=%i, normSort[inorm]=%i \n", inorm, nNorm, normSort[inorm]);
281  normString = ((TObjString*)listOfNormalizationVariables->At(normSort[inorm]))->String();
282  // printf(" walking in normalization, i=%i, normString=%s \n", inorm, normString.Data());
283  str.ReplaceAll(searchString + normString, replaceString + normString);
284  // like: str.ReplaceAll("CEQmean_Mean", "C!EQmean_Mean");
285  }
286  str.ReplaceAll(searchString + fAbbreviation, replaceString + fAbbreviation);
287  // like: str.ReplaceAll("CEQmean~", "C!EQmean~");
288  str.ReplaceAll(searchString + fAppendString, replaceString + fAppendString);
289  // like: str.ReplaceAll("CEQmean.fElements", "C!EQmean.fElements");
290 
291  // ***** add missing extensions *****
292  str.ReplaceAll(searchString, replaceString + fAbbreviation);
293  // like: str.ReplaceAll("CEQmean", "C!EQmean~");
294  }
295 
296  // ***** undo saving *****
297  str.ReplaceAll(removeString, "");
298 
299  if (printDrawCommand) std::cout << "The string looks now like: " << str.Data() << std::endl;
300  delete [] varLengths;
301  delete [] normLengths;
302  delete [] varSort;
303  delete [] normSort;
304  return str.Data();
305 }
306 
307 
308 
309 
310 //_____________________________________________________________________________
311 Int_t AliTPCCalibViewer::EasyDraw(const char* drawCommand, const char* sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
320 
321  TString drawStr(drawCommand);
322  TString sectorStr(sector);
323  sectorStr.ToUpper();
324  TString cutStr("");
325  //TString drawOptionsStr("profcolz ");
326  Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
327  if (dangerousToDraw) {
328  Warning("EasyDraw", "The draw string must not contain ':' or '>>'. Using only first variable for drawing!");
329 // return -1;
330 // drawStr.Resize(drawStr.First(">"));
331  drawStr.Resize(drawStr.First(":"));
332  }
333 
334  TString drawOptionsStr("");
335  TRandom rnd(0);
336  Int_t rndNumber = rnd.Integer(10000);
337 
338  if (drawOptions && strcmp(drawOptions, "") != 0)
339  drawOptionsStr += drawOptions;
340  else
341  drawOptionsStr += "profcolz";
342 
343  if (sectorStr == "A") {
344  drawStr += Form(":gy%s:gx%s>>prof", fAppendString.Data(), fAppendString.Data());
345  drawStr += rndNumber;
346  drawStr += "(330,-250,250,330,-250,250)";
347  cutStr += "(sector/18)%2==0 ";
348  }
349  else if (sectorStr == "C") {
350  drawStr += Form(":gy%s:gx%s>>prof", fAppendString.Data(), fAppendString.Data());
351  drawStr += rndNumber;
352  drawStr += "(330,-250,250,330,-250,250)";
353  cutStr += "(sector/18)%2==1 ";
354  }
355  else if (sectorStr == "ALL") {
356  drawStr += Form(":gy%s:gx%s>>prof", fAppendString.Data(), fAppendString.Data());
357  drawStr += rndNumber;
358  drawStr += "(330,-250,250,330,-250,250)";
359  }
360  else if (sectorStr.Contains("S")) {
361  drawStr += Form(":rpad%s:row%s+(sector>35)*63>>prof", fAppendString.Data(), fAppendString.Data());
362  drawStr += rndNumber;
363  drawStr += "(159,0,159,140,-70,70)";
364  TString sec=sectorStr;
365  sec.Remove(0,1);
366  cutStr += "sector%36=="+sec+" ";
367  }
368  else if (sectorStr.IsDigit()) {
369  Int_t isec = sectorStr.Atoi();
370  drawStr += Form(":rpad%s:row%s>>prof", fAppendString.Data(), fAppendString.Data());
371  drawStr += rndNumber;
372  if (isec < 36 && isec >= 0)
373  drawStr += "(63,0,63,108,-54,54)";
374  else if (isec < 72 && isec >= 36)
375  drawStr += "(96,0,96,140,-70,70)";
376  else {
377  Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
378  return -1;
379  }
380  cutStr += "(sector==";
381  cutStr += isec;
382  cutStr += ") ";
383  }
384 
385  if (cuts && cuts[0] != 0) {
386  if (cutStr.Length() != 0) cutStr += "&& ";
387  cutStr += "(";
388  cutStr += cuts;
389  cutStr += ")";
390  }
391  drawStr.ReplaceAll(fAbbreviation, fAppendString);
392  cutStr.ReplaceAll(fAbbreviation, fAppendString);
393  if (writeDrawCommand) std::cout << "fTree->Draw(\"" << drawStr << "\", \"" << cutStr << "\", \"" << drawOptionsStr << "\");" << std::endl;
394  Int_t returnValue = fTree->Draw(drawStr.Data(), cutStr.Data(), drawOptionsStr.Data());
395  TString profName("prof");
396  profName += rndNumber;
397  TObject *obj = gDirectory->Get(profName.Data());
398  if (obj && obj->InheritsFrom("TH1")) FormatHistoLabels((TH1*)obj);
399  return returnValue;
400 }
401 
402 
403 Int_t AliTPCCalibViewer::EasyDraw(const char* drawCommand, Int_t sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
410 
411  if (sector >= 0 && sector < 72) {
412  return EasyDraw(drawCommand, Form("%i", sector), cuts, drawOptions, writeDrawCommand);
413  }
414  Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
415  return -1;
416 }
417 
418 
419 //_____________________________________________________________________________
420 Int_t AliTPCCalibViewer::EasyDraw1D(const char* drawCommand, const char* sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
429 
430  TString drawStr(drawCommand);
431  TString sectorStr(sector);
432  TString drawOptionsStr(drawOptions);
433  sectorStr.ToUpper();
434  TString cutStr("");
435 
436  if (sectorStr == "A")
437  cutStr += "(sector/18)%2==0 ";
438  else if (sectorStr == "C")
439  cutStr += "(sector/18)%2==1 ";
440  else if (sectorStr.IsDigit()) {
441  Int_t isec = sectorStr.Atoi();
442  if (isec < 0 || isec > 71) {
443  Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
444  return -1;
445  }
446  cutStr += "(sector==";
447  cutStr += isec;
448  cutStr += ") ";
449  }
450  else if (sectorStr.Contains("S")) {
451  TString sec=sectorStr;
452  sec.Remove(0,1);
453  cutStr += "sector%36=="+sec+" ";
454  }
455 
456  if (cuts && cuts[0] != 0) {
457  if (cutStr.Length() != 0) cutStr += "&& ";
458  cutStr += "(";
459  cutStr += cuts;
460  cutStr += ")";
461  }
462 
463  drawStr.ReplaceAll(fAbbreviation, fAppendString);
464  cutStr.ReplaceAll(fAbbreviation, fAppendString);
465  if (writeDrawCommand) std::cout << "fTree->Draw(\"" << drawStr << "\", \"" << cutStr << "\", \"" << drawOptionsStr << "\");" << std::endl;
466  Int_t returnValue = fTree->Draw(drawStr.Data(), cutStr.Data(), drawOptionsStr.Data());
467  if (returnValue == -1) return -1;
468 
469  TObject *obj = (gPad) ? gPad->GetPrimitive("htemp") : 0;
470  if (!obj) obj = (TH1F*)gDirectory->Get("htemp");
471  if (!obj) obj = gPad->GetPrimitive("tempHist");
472  if (!obj) obj = (TH1F*)gDirectory->Get("tempHist");
473  if (!obj) obj = gPad->GetPrimitive("Graph");
474  if (!obj) obj = (TH1F*)gDirectory->Get("Graph");
475  if (obj && obj->InheritsFrom("TH1")) FormatHistoLabels((TH1*)obj);
476  return returnValue;
477 }
478 
479 
480 Int_t AliTPCCalibViewer::EasyDraw1D(const char* drawCommand, Int_t sector, const char* cuts, const char* drawOptions, Bool_t writeDrawCommand) const {
487 
488  if (sector >= 0 && sector < 72) {
489  return EasyDraw1D(drawCommand, Form("%i",sector), cuts, drawOptions, writeDrawCommand);
490  }
491  Error("EasyDraw","The TPC contains only sectors between 0 and 71.");
492  return -1;
493 }
494 
495 
496 void AliTPCCalibViewer::FormatHistoLabels(TH1 *histo) const {
499 
500  if (!histo) return;
501  TString replaceString(fAppendString.Data());
502  TString *str = new TString(histo->GetTitle());
503  str->ReplaceAll(replaceString, "");
504  histo->SetTitle(str->Data());
505  delete str;
506  if (histo->GetXaxis()) {
507  str = new TString(histo->GetXaxis()->GetTitle());
508  str->ReplaceAll(replaceString, "");
509  histo->GetXaxis()->SetTitle(str->Data());
510  delete str;
511  }
512  if (histo->GetYaxis()) {
513  str = new TString(histo->GetYaxis()->GetTitle());
514  str->ReplaceAll(replaceString, "");
515  histo->GetYaxis()->SetTitle(str->Data());
516  delete str;
517  }
518  if (histo->GetZaxis()) {
519  str = new TString(histo->GetZaxis()->GetTitle());
520  str->ReplaceAll(replaceString, "");
521  histo->GetZaxis()->SetTitle(str->Data());
522  delete str;
523  }
524 }
525 
526 
527 Int_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 {
533 
534  if (sector >= 0 && sector < 72) {
535  return DrawHisto1D(drawCommand, Form("%i", sector), cuts, sigmas, plotMean, plotMedian, plotLTM);
536  }
537  Error("DrawHisto1D","The TPC contains only sectors between 0 and 71.");
538  return -1;
539 }
540 
541 
542 Int_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 {
548 
549  Int_t oldOptStat = gStyle->GetOptStat();
550  gStyle->SetOptStat(0000000);
551  Double_t ltmFraction = 0.8;
552 
553  TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
554  TVectorF nsigma(sigmasTokens->GetEntriesFast());
555  for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
556  TString str(((TObjString*)sigmasTokens->At(i))->GetString());
557  Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
558  nsigma[i] = sig;
559  }
560  delete sigmasTokens;
561 
562  TString drawStr(drawCommand);
563  Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
564  if (dangerousToDraw) {
565  Warning("DrawHisto1D", "The draw string must not contain ':' or '>>'.");
566  return -1;
567  }
568  drawStr += " >> tempHist";
569  Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts);
570  TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
571  // FIXME is this histogram deleted automatically?
572  Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
573 
574  Double_t mean = TMath::Mean(entries, values);
575  Double_t median = TMath::Median(entries, values);
576  Double_t sigma = TMath::RMS(entries, values);
577  Double_t maxY = htemp->GetMaximum();
578 
579  TLegend * legend = new TLegend(.7,.7, .99, .99, "Statistical information");
580  //fListOfObjectsToBeDeleted->Add(legend);
581 
582  if (plotMean) {
583  // draw Mean
584  TLine* line = new TLine(mean, 0, mean, maxY);
585  //fListOfObjectsToBeDeleted->Add(line);
586  line->SetLineColor(kRed);
587  line->SetLineWidth(2);
588  line->SetLineStyle(1);
589  line->Draw();
590  legend->AddEntry(line, Form("Mean: %f", mean), "l");
591  // draw sigma lines
592  for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
593  TLine* linePlusSigma = new TLine(mean + nsigma[i] * sigma, 0, mean + nsigma[i] * sigma, maxY);
594  //fListOfObjectsToBeDeleted->Add(linePlusSigma);
595  linePlusSigma->SetLineColor(kRed);
596  linePlusSigma->SetLineStyle(2 + i);
597  linePlusSigma->Draw();
598  TLine* lineMinusSigma = new TLine(mean - nsigma[i] * sigma, 0, mean - nsigma[i] * sigma, maxY);
599  //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
600  lineMinusSigma->SetLineColor(kRed);
601  lineMinusSigma->SetLineStyle(2 + i);
602  lineMinusSigma->Draw();
603  legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma)), "l");
604  }
605  }
606  if (plotMedian) {
607  // draw median
608  TLine* line = new TLine(median, 0, median, maxY);
609  //fListOfObjectsToBeDeleted->Add(line);
610  line->SetLineColor(kBlue);
611  line->SetLineWidth(2);
612  line->SetLineStyle(1);
613  line->Draw();
614  legend->AddEntry(line, Form("Median: %f", median), "l");
615  // draw sigma lines
616  for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
617  TLine* linePlusSigma = new TLine(median + nsigma[i] * sigma, 0, median + nsigma[i]*sigma, maxY);
618  //fListOfObjectsToBeDeleted->Add(linePlusSigma);
619  linePlusSigma->SetLineColor(kBlue);
620  linePlusSigma->SetLineStyle(2 + i);
621  linePlusSigma->Draw();
622  TLine* lineMinusSigma = new TLine(median - nsigma[i] * sigma, 0, median - nsigma[i]*sigma, maxY);
623  //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
624  lineMinusSigma->SetLineColor(kBlue);
625  lineMinusSigma->SetLineStyle(2 + i);
626  lineMinusSigma->Draw();
627  legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma)), "l");
628  }
629  }
630  if (plotLTM) {
631  // draw LTM
632  Double_t ltmRms = 0;
633  Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
634  TLine* line = new TLine(ltm, 0, ltm, maxY);
635  //fListOfObjectsToBeDeleted->Add(line);
636  line->SetLineColor(kGreen+2);
637  line->SetLineWidth(2);
638  line->SetLineStyle(1);
639  line->Draw();
640  legend->AddEntry(line, Form("LTM: %f", ltm), "l");
641  // draw sigma lines
642  for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
643  TLine* linePlusSigma = new TLine(ltm + nsigma[i] * ltmRms, 0, ltm + nsigma[i] * ltmRms, maxY);
644  //fListOfObjectsToBeDeleted->Add(linePlusSigma);
645  linePlusSigma->SetLineColor(kGreen+2);
646  linePlusSigma->SetLineStyle(2+i);
647  linePlusSigma->Draw();
648 
649  TLine* lineMinusSigma = new TLine(ltm - nsigma[i] * ltmRms, 0, ltm - nsigma[i] * ltmRms, maxY);
650  //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
651  lineMinusSigma->SetLineColor(kGreen+2);
652  lineMinusSigma->SetLineStyle(2+i);
653  lineMinusSigma->Draw();
654  legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f", (Int_t)(nsigma[i]), (Float_t)(nsigma[i] * ltmRms)), "l");
655  }
656  }
657  if (!plotMean && !plotMedian && !plotLTM) return -1;
658  legend->Draw();
659  gStyle->SetOptStat(oldOptStat);
660  return 1;
661 }
662 
663 
664 Int_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 {
681 
682  if (sector >= 0 && sector < 72) {
683  return SigmaCut(drawCommand, Form("%i", sector), cuts, sigmaMax, plotMean, plotMedian, plotLTM, pm, sigmas, sigmaStep);
684  }
685  Error("SigmaCut","The TPC contains only sectors between 0 and 71.");
686  return -1;
687 }
688 
689 
690 Int_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 {
697 
698  Double_t ltmFraction = 0.8;
699 
700  TString drawStr(drawCommand);
701  Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
702  if (dangerousToDraw) {
703  Warning("SigmaCut", "The draw string must not contain ':' or '>>'.");
704  return -1;
705  }
706  drawStr += " >> tempHist";
707 
708  Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
709  TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
710  // FIXME is this histogram deleted automatically?
711  Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
712 
713  Double_t mean = TMath::Mean(entries, values);
714  Double_t median = TMath::Median(entries, values);
715  Double_t sigma = TMath::RMS(entries, values);
716 
717  TLegend * legend = new TLegend(.7,.7, .99, .99, "Cumulative");
718  //fListOfObjectsToBeDeleted->Add(legend);
719  TH1F *cutHistoMean = 0;
720  TH1F *cutHistoMedian = 0;
721  TH1F *cutHistoLTM = 0;
722 
723  TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
724  TVectorF nsigma(sigmasTokens->GetEntriesFast());
725  for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
726  TString str(((TObjString*)sigmasTokens->At(i))->GetString());
727  Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
728  nsigma[i] = sig;
729  }
730  delete sigmasTokens;
731 
732  if (plotMean) {
733  cutHistoMean = AliTPCCalibViewer::SigmaCut(htemp, mean, sigma, sigmaMax, sigmaStep, pm);
734  if (cutHistoMean) {
735  //fListOfObjectsToBeDeleted->Add(cutHistoMean);
736  cutHistoMean->SetLineColor(kRed);
737  legend->AddEntry(cutHistoMean, "Mean", "l");
738  cutHistoMean->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
739  cutHistoMean->Draw();
740  DrawLines(cutHistoMean, nsigma, legend, kRed, pm);
741  } // if (cutHistoMean)
742 
743  }
744  if (plotMedian) {
745  cutHistoMedian = AliTPCCalibViewer::SigmaCut(htemp, median, sigma, sigmaMax, sigmaStep, pm);
746  if (cutHistoMedian) {
747  //fListOfObjectsToBeDeleted->Add(cutHistoMedian);
748  cutHistoMedian->SetLineColor(kBlue);
749  legend->AddEntry(cutHistoMedian, "Median", "l");
750  cutHistoMedian->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
751  if (plotMean && cutHistoMean) cutHistoMedian->Draw("same");
752  else cutHistoMedian->Draw();
753  DrawLines(cutHistoMedian, nsigma, legend, kBlue, pm);
754  } // if (cutHistoMedian)
755  }
756  if (plotLTM) {
757  Double_t ltmRms = 0;
758  Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
759  cutHistoLTM = AliTPCCalibViewer::SigmaCut(htemp, ltm, ltmRms, sigmaMax, sigmaStep, pm);
760  if (cutHistoLTM) {
761  //fListOfObjectsToBeDeleted->Add(cutHistoLTM);
762  cutHistoLTM->SetLineColor(kGreen+2);
763  legend->AddEntry(cutHistoLTM, "LTM", "l");
764  cutHistoLTM->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
765  if ((plotMean && cutHistoMean) || (plotMedian && cutHistoMedian)) cutHistoLTM->Draw("same");
766  else cutHistoLTM->Draw();
767  DrawLines(cutHistoLTM, nsigma, legend, kGreen+2, pm);
768  }
769  }
770  if (!plotMean && !plotMedian && !plotLTM) return -1;
771  legend->Draw();
772  return 1;
773 }
774 
775 
776 Int_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 {
783 
784  // Double_t ltmFraction = 0.8; //unused
785 
786  TString drawStr(drawCommand);
787  drawStr += " >> tempHist";
788 
789  Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
790  TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
791  TGraph *cutGraphMean = 0;
792  // TGraph *cutGraphMedian = 0;
793  // TGraph *cutGraphLTM = 0;
794  Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
795  Int_t *index = new Int_t[entries];
796  Float_t *xarray = new Float_t[entries];
797  Float_t *yarray = new Float_t[entries];
798  TMath::Sort(entries, values, index, kFALSE);
799 
800  Double_t mean = TMath::Mean(entries, values);
801  // Double_t median = TMath::Median(entries, values);
802  Double_t sigma = TMath::RMS(entries, values);
803 
804  TLegend * legend = new TLegend(.7,.7, .99, .99, "Cumulative");
805  //fListOfObjectsToBeDeleted->Add(legend);
806 
807  // parse sigmas string
808  TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
809  TVectorF nsigma(sigmasTokens->GetEntriesFast());
810  for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
811  TString str(((TObjString*)sigmasTokens->At(i))->GetString());
812  Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
813  nsigma[i] = sig;
814  }
815  delete sigmasTokens;
816 
817  if (plotMean) {
818  for (Int_t i = 0; i < entries; i++) {
819  xarray[i] = TMath::Abs(values[index[i]] - mean) / sigma;
820  yarray[i] = float(i) / float(entries);
821  }
822  cutGraphMean = new TGraph(entries, xarray, yarray);
823  if (cutGraphMean) {
824  //fListOfObjectsToBeDeleted->Add(cutGraphMean);
825  cutGraphMean->SetLineColor(kRed);
826  legend->AddEntry(cutGraphMean, "Mean", "l");
827  cutGraphMean->SetTitle(Form("%s, Cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
828  cutGraphMean->Draw("alu");
829  DrawLines(cutGraphMean, nsigma, legend, kRed, kTRUE);
830  }
831  }
832  delete [] index;
833  delete [] xarray;
834  delete [] yarray;
835  /*
836  if (plotMedian) {
837  cutHistoMedian = AliTPCCalibViewer::SigmaCut(htemp, median, sigma, sigmaMax, sigmaStep, pm);
838  if (cutHistoMedian) {
839  fListOfObjectsToBeDeleted->Add(cutHistoMedian);
840  cutHistoMedian->SetLineColor(kBlue);
841  legend->AddEntry(cutHistoMedian, "Median", "l");
842  cutHistoMedian->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
843  if (plotMean && cutHistoMean) cutHistoMedian->Draw("same");
844  else cutHistoMedian->Draw();
845  DrawLines(cutHistoMedian, nsigma, legend, kBlue, pm);
846  } // if (cutHistoMedian)
847  }
848  if (plotLTM) {
849  Double_t ltmRms = 0;
850  Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
851  cutHistoLTM = AliTPCCalibViewer::SigmaCut(htemp, ltm, ltmRms, sigmaMax, sigmaStep, pm);
852  if (cutHistoLTM) {
853  fListOfObjectsToBeDeleted->Add(cutHistoLTM);
854  cutHistoLTM->SetLineColor(kGreen+2);
855  legend->AddEntry(cutHistoLTM, "LTM", "l");
856  cutHistoLTM->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
857  if (plotMean && cutHistoMean || plotMedian && cutHistoMedian) cutHistoLTM->Draw("same");
858  else cutHistoLTM->Draw();
859  DrawLines(cutHistoLTM, nsigma, legend, kGreen+2, pm);
860  }
861  }*/
862  if (!plotMean && !plotMedian && !plotLTM) return -1;
863  legend->Draw();
864  return 1;
865 }
866 
867 
868 Int_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 {
879 
880  if (sector >= 0 && sector < 72) {
881  return Integrate(drawCommand, Form("%i", sector), cuts, sigmaMax, plotMean, plotMedian, plotLTM, sigmas, sigmaStep);
882  }
883  Error("Integrate","The TPC contains only sectors between 0 and 71.");
884  return -1;
885 
886 }
887 
888 
889 Int_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 {
900 
901 
902  Double_t ltmFraction = 0.8;
903 
904  TString drawStr(drawCommand);
905  drawStr += " >> tempHist";
906 
907  Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
908  TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
909  // FIXME is this histogram deleted automatically?
910  Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
911 
912  Double_t mean = TMath::Mean(entries, values);
913  Double_t median = TMath::Median(entries, values);
914  Double_t sigma = TMath::RMS(entries, values);
915 
916  TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
917  TVectorF nsigma(sigmasTokens->GetEntriesFast());
918  for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
919  TString str(((TObjString*)sigmasTokens->At(i))->GetString());
920  Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
921  nsigma[i] = sig;
922  }
923  delete sigmasTokens;
924 
925  TLegend * legend = new TLegend(.7,.7, .99, .99, "Integrated histogram");
926  //fListOfObjectsToBeDeleted->Add(legend);
927  TH1F *integralHistoMean = 0;
928  TH1F *integralHistoMedian = 0;
929  TH1F *integralHistoLTM = 0;
930 
931  if (plotMean) {
932  integralHistoMean = AliTPCCalibViewer::Integrate(htemp, mean, sigma, sigmaMax, sigmaStep);
933  if (integralHistoMean) {
934  //fListOfObjectsToBeDeleted->Add(integralHistoMean);
935  integralHistoMean->SetLineColor(kRed);
936  legend->AddEntry(integralHistoMean, "Mean", "l");
937  integralHistoMean->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
938  integralHistoMean->Draw();
939  DrawLines(integralHistoMean, nsigma, legend, kRed, kTRUE);
940  }
941  }
942  if (plotMedian) {
943  integralHistoMedian = AliTPCCalibViewer::Integrate(htemp, median, sigma, sigmaMax, sigmaStep);
944  if (integralHistoMedian) {
945  //fListOfObjectsToBeDeleted->Add(integralHistoMedian);
946  integralHistoMedian->SetLineColor(kBlue);
947  legend->AddEntry(integralHistoMedian, "Median", "l");
948  integralHistoMedian->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
949  if (plotMean && integralHistoMean) integralHistoMedian->Draw("same");
950  else integralHistoMedian->Draw();
951  DrawLines(integralHistoMedian, nsigma, legend, kBlue, kTRUE);
952  }
953  }
954  if (plotLTM) {
955  Double_t ltmRms = 0;
956  Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
957  integralHistoLTM = AliTPCCalibViewer::Integrate(htemp, ltm, ltmRms, sigmaMax, sigmaStep);
958  if (integralHistoLTM) {
959  //fListOfObjectsToBeDeleted->Add(integralHistoLTM);
960  integralHistoLTM->SetLineColor(kGreen+2);
961  legend->AddEntry(integralHistoLTM, "LTM", "l");
962  integralHistoLTM->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
963  if ((plotMean && integralHistoMean) || (plotMedian && integralHistoMedian)) integralHistoLTM->Draw("same");
964  else integralHistoLTM->Draw();
965  DrawLines(integralHistoLTM, nsigma, legend, kGreen+2, kTRUE);
966  }
967  }
968  if (!plotMean && !plotMedian && !plotLTM) return -1;
969  legend->Draw();
970  return 1;
971 }
972 
973 
974 Int_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 {
985 
986  Double_t ltmFraction = 0.8;
987 
988  TString drawStr(drawCommand);
989  Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
990  if (dangerousToDraw) {
991  Warning("Integrate", "The draw string must not contain ':' or '>>'.");
992  return -1;
993  }
994  drawStr += " >> tempHist";
995 
996  Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
997  TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
998  TGraph *integralGraphMean = 0;
999  TGraph *integralGraphMedian = 0;
1000  TGraph *integralGraphLTM = 0;
1001  Double_t *values = fTree->GetV1(); // value is the array containing 'entries' numbers
1002  Int_t *index = new Int_t[entries];
1003  Float_t *xarray = new Float_t[entries];
1004  Float_t *yarray = new Float_t[entries];
1005  TMath::Sort(entries, values, index, kFALSE);
1006 
1007  Double_t mean = TMath::Mean(entries, values);
1008  Double_t median = TMath::Median(entries, values);
1009  Double_t sigma = TMath::RMS(entries, values);
1010 
1011  // parse sigmas string
1012  TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
1013  TVectorF nsigma(sigmasTokens->GetEntriesFast());
1014  for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
1015  TString str(((TObjString*)sigmasTokens->At(i))->GetString());
1016  Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
1017  nsigma[i] = sig;
1018  }
1019  delete sigmasTokens;
1020 
1021  TLegend * legend = new TLegend(.7,.7, .99, .99, "Integrated histogram");
1022  //fListOfObjectsToBeDeleted->Add(legend);
1023 
1024  if (plotMean) {
1025  for (Int_t i = 0; i < entries; i++) {
1026  xarray[i] = (values[index[i]] - mean) / sigma;
1027  yarray[i] = float(i) / float(entries);
1028  }
1029  integralGraphMean = new TGraph(entries, xarray, yarray);
1030  if (integralGraphMean) {
1031  //fListOfObjectsToBeDeleted->Add(integralGraphMean);
1032  integralGraphMean->SetLineColor(kRed);
1033  legend->AddEntry(integralGraphMean, "Mean", "l");
1034  integralGraphMean->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
1035  integralGraphMean->Draw("alu");
1036  DrawLines(integralGraphMean, nsigma, legend, kRed, kTRUE);
1037  }
1038  }
1039  if (plotMedian) {
1040  for (Int_t i = 0; i < entries; i++) {
1041  xarray[i] = (values[index[i]] - median) / sigma;
1042  yarray[i] = float(i) / float(entries);
1043  }
1044  integralGraphMedian = new TGraph(entries, xarray, yarray);
1045  if (integralGraphMedian) {
1046  //fListOfObjectsToBeDeleted->Add(integralGraphMedian);
1047  integralGraphMedian->SetLineColor(kBlue);
1048  legend->AddEntry(integralGraphMedian, "Median", "l");
1049  integralGraphMedian->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
1050  if (plotMean && integralGraphMean) integralGraphMedian->Draw("samelu");
1051  else integralGraphMedian->Draw("alu");
1052  DrawLines(integralGraphMedian, nsigma, legend, kBlue, kTRUE);
1053  }
1054  }
1055  if (plotLTM) {
1056  Double_t ltmRms = 0;
1057  Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
1058  for (Int_t i = 0; i < entries; i++) {
1059  xarray[i] = (values[index[i]] - ltm) / ltmRms;
1060  yarray[i] = float(i) / float(entries);
1061  }
1062  integralGraphLTM = new TGraph(entries, xarray, yarray);
1063  if (integralGraphLTM) {
1064  //fListOfObjectsToBeDeleted->Add(integralGraphLTM);
1065  integralGraphLTM->SetLineColor(kGreen+2);
1066  legend->AddEntry(integralGraphLTM, "LTM", "l");
1067  integralGraphLTM->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
1068  if ((plotMean && integralGraphMean) || (plotMedian && integralGraphMedian)) integralGraphLTM->Draw("samelu");
1069  else integralGraphLTM->Draw("alu");
1070  DrawLines(integralGraphLTM, nsigma, legend, kGreen+2, kTRUE);
1071  }
1072  }
1073  delete [] index;
1074  delete [] xarray;
1075  delete [] yarray;
1076 
1077  if (!plotMean && !plotMedian && !plotLTM) return -1;
1078  legend->Draw();
1079  return entries;
1080 }
1081 
1082 
1083 void AliTPCCalibViewer::DrawLines(TH1F *histogram, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
1086 
1087  // start to draw the lines, loop over requested sigmas
1088  for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
1089  if (!pm) {
1090  Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
1091  TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
1092  //fListOfObjectsToBeDeleted->Add(lineUp);
1093  lineUp->SetLineColor(color);
1094  lineUp->SetLineStyle(2 + i);
1095  lineUp->Draw();
1096  TLine* lineLeft = new TLine(nsigma[i], histogram->GetBinContent(bin), 0, histogram->GetBinContent(bin));
1097  //fListOfObjectsToBeDeleted->Add(lineLeft);
1098  lineLeft->SetLineColor(color);
1099  lineLeft->SetLineStyle(2 + i);
1100  lineLeft->Draw();
1101  legend->AddEntry(lineLeft, Form("Fraction(%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
1102  }
1103  else { // if (pm)
1104  Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
1105  TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
1106  //fListOfObjectsToBeDeleted->Add(lineUp1);
1107  lineUp1->SetLineColor(color);
1108  lineUp1->SetLineStyle(2 + i);
1109  lineUp1->Draw();
1110  TLine* lineLeft1 = new TLine(nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
1111  //fListOfObjectsToBeDeleted->Add(lineLeft1);
1112  lineLeft1->SetLineColor(color);
1113  lineLeft1->SetLineStyle(2 + i);
1114  lineLeft1->Draw();
1115  legend->AddEntry(lineLeft1, Form("Fraction(+%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
1116  bin = histogram->GetXaxis()->FindBin(-nsigma[i]);
1117  TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], histogram->GetBinContent(bin));
1118  //fListOfObjectsToBeDeleted->Add(lineUp2);
1119  lineUp2->SetLineColor(color);
1120  lineUp2->SetLineStyle(2 + i);
1121  lineUp2->Draw();
1122  TLine* lineLeft2 = new TLine(-nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
1123  //fListOfObjectsToBeDeleted->Add(lineLeft2);
1124  lineLeft2->SetLineColor(color);
1125  lineLeft2->SetLineStyle(2 + i);
1126  lineLeft2->Draw();
1127  legend->AddEntry(lineLeft2, Form("Fraction(-%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
1128  }
1129  } // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
1130 }
1131 
1132 
1133 void AliTPCCalibViewer::DrawLines(TGraph *graph, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
1136 
1137  // start to draw the lines, loop over requested sigmas
1138  for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
1139  if (!pm) {
1140  TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
1141  //fListOfObjectsToBeDeleted->Add(lineUp);
1142  lineUp->SetLineColor(color);
1143  lineUp->SetLineStyle(2 + i);
1144  lineUp->Draw();
1145  TLine* lineLeft = new TLine(nsigma[i], graph->Eval(nsigma[i]), 0, graph->Eval(nsigma[i]));
1146  //fListOfObjectsToBeDeleted->Add(lineLeft);
1147  lineLeft->SetLineColor(color);
1148  lineLeft->SetLineStyle(2 + i);
1149  lineLeft->Draw();
1150  legend->AddEntry(lineLeft, Form("Fraction(%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i])), "l");
1151  }
1152  else { // if (pm)
1153  TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
1154  //fListOfObjectsToBeDeleted->Add(lineUp1);
1155  lineUp1->SetLineColor(color);
1156  lineUp1->SetLineStyle(2 + i);
1157  lineUp1->Draw();
1158  TLine* lineLeft1 = new TLine(nsigma[i], graph->Eval(nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(nsigma[i]));
1159  //fListOfObjectsToBeDeleted->Add(lineLeft1);
1160  lineLeft1->SetLineColor(color);
1161  lineLeft1->SetLineStyle(2 + i);
1162  lineLeft1->Draw();
1163  legend->AddEntry(lineLeft1, Form("Fraction(+%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i])), "l");
1164  TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], graph->Eval(-nsigma[i]));
1165  //fListOfObjectsToBeDeleted->Add(lineUp2);
1166  lineUp2->SetLineColor(color);
1167  lineUp2->SetLineStyle(2 + i);
1168  lineUp2->Draw();
1169  TLine* lineLeft2 = new TLine(-nsigma[i], graph->Eval(-nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(-nsigma[i]));
1170  //fListOfObjectsToBeDeleted->Add(lineLeft2);
1171  lineLeft2->SetLineColor(color);
1172  lineLeft2->SetLineStyle(2 + i);
1173  lineLeft2->Draw();
1174  legend->AddEntry(lineLeft2, Form("Fraction(-%f #sigma) = %f",nsigma[i], graph->Eval(-nsigma[i])), "l");
1175  }
1176  } // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
1177 }
1178 
1179 
1180 
1181 
1182 
1184 // Array tools //
1186 
1187 
1188 Int_t AliTPCCalibViewer::GetBin(Float_t value, Int_t nbins, Double_t binLow, Double_t binUp){
1196 
1197  Int_t bin = TMath::Nint( (Float_t)(value - binLow) / (Float_t)(binUp - binLow) * (nbins-1) ) + 1;
1198  // avoid index out of bounds:
1199  if (value < binLow) bin = 0;
1200  if (value > binUp) bin = nbins + 1;
1201  return bin;
1202 
1203 }
1204 
1205 
1206 Double_t AliTPCCalibViewer::GetLTM(Int_t n, const Double_t *const array, Double_t *const sigma, Double_t fraction){
1208 
1209  Double_t *ddata = new Double_t[n];
1210  Double_t mean = 0, lsigma = 0;
1211  UInt_t nPoints = 0;
1212  for (UInt_t i = 0; i < (UInt_t)n; i++) {
1213  ddata[nPoints]= array[nPoints];
1214  nPoints++;
1215  }
1216  Int_t hh = TMath::Min(TMath::Nint(fraction * nPoints), Int_t(n));
1217  AliMathBase::EvaluateUni(nPoints, ddata, mean, lsigma, hh);
1218  if (sigma) *sigma = lsigma;
1219  delete [] ddata;
1220  return mean;
1221 }
1222 
1223 
1224 TH1F* AliTPCCalibViewer::SigmaCut(TH1F *const histogram, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep, Bool_t pm) {
1239 
1240  Float_t *array = histogram->GetArray();
1241  Int_t nbins = histogram->GetXaxis()->GetNbins();
1242  Float_t binLow = histogram->GetXaxis()->GetXmin();
1243  Float_t binUp = histogram->GetXaxis()->GetXmax();
1244  return AliTPCCalibViewer::SigmaCut(nbins, array, mean, sigma, nbins, binLow, binUp, sigmaMax, sigmaStep, pm);
1245 }
1246 
1247 
1248 TH1F* AliTPCCalibViewer::SigmaCut(Int_t n, const 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){
1256 
1257  if (sigma == 0) return 0;
1258  Float_t binWidth = (binUp-binLow)/(nbins - 1);
1259  if (sigmaStep <= 0) sigmaStep = binWidth;
1260  Int_t kbins = (Int_t)(sigmaMax * sigma / sigmaStep) + 1; // + 1 due to overflow bin in histograms
1261  if (pm) kbins = 2 * (Int_t)(sigmaMax * sigma / sigmaStep) + 1;
1262  Float_t kbinLow = !pm ? 0 : -sigmaMax;
1263  Float_t kbinUp = sigmaMax;
1264  TH1F *hist = new TH1F("sigmaCutHisto","Cumulative; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
1265  hist->SetDirectory(0);
1266  hist->Reset();
1267 
1268  // calculate normalization
1269  Double_t normalization = 0;
1270  for (Int_t i = 0; i <= n; i++) {
1271  normalization += array[i];
1272  }
1273 
1274  // given units: units from given histogram
1275  // sigma units: in units of sigma
1276  // iDelta: integrate in interval (mean +- iDelta), given units
1277  // x: ofset from mean for integration, given units
1278  // hist: needs
1279 
1280 // printf("nbins: %i, binLow: %f, binUp: %f \n", nbins, binLow, binUp);
1281  // fill histogram
1282  for (Float_t iDelta = 0; iDelta <= sigmaMax * sigma; iDelta += sigmaStep) {
1283  // integrate array
1284  Double_t valueP = array[GetBin(mean, nbins, binLow, binUp)];
1285  Double_t valueM = array[GetBin(mean-binWidth, nbins, binLow, binUp)];
1286  // add bin of mean value only once to the histogram
1287 // printf("++ adding bins: ");
1288  for (Float_t x = binWidth; x <= iDelta; x += binWidth) {
1289  valueP += (mean + x <= binUp) ? array[GetBin(mean + x, nbins, binLow, binUp)] : 0;
1290  valueM += (mean-binWidth - x >= binLow) ? array[GetBin(mean-binWidth - x, nbins, binLow, binUp)] : 0;
1291 // printf("%i, ", GetBin(mean + x, nbins, binLow, binUp));
1292  }
1293 // printf("\n");
1294  if (valueP / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", valueP, normalization);
1295  if (valueP / normalization > 100) return hist;
1296  if (valueM / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", valueM, normalization);
1297  if (valueM / normalization > 100) return hist;
1298  valueP = (valueP / normalization);
1299  valueM = (valueM / normalization);
1300  if (pm) {
1301  Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1302  hist->SetBinContent(bin, valueP);
1303  bin = GetBin(-iDelta/sigma, kbins, kbinLow, kbinUp);
1304  hist->SetBinContent(bin, valueM);
1305  }
1306  else { // if (!pm)
1307  Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1308  hist->SetBinContent(bin, valueP + valueM);
1309 // printf(" first integration bin: %i, last integration bin in + direction: %i \n", GetBin(mean+binWidth, nbins, binLow, binUp), GetBin(iDelta, nbins, binLow, binUp));
1310 // printf(" first integration bin: %i, last integration bin in - direction: %i \n", GetBin(mean+binWidth, nbins, binLow, binUp), GetBin(-iDelta, nbins, binLow, binUp));
1311 // printf(" value: %f, normalization: %f, iDelta: %f, Bin: %i \n", valueP+valueM, normalization, iDelta, bin);
1312  }
1313  }
1314  //hist->SetMaximum(0.7);
1315  if (!pm) hist->SetMaximum(1.2);
1316  return hist;
1317 }
1318 
1319 
1320 TH1F* AliTPCCalibViewer::SigmaCut(Int_t /*n*/, const Double_t */*array*/, Double_t /*mean*/, Double_t /*sigma*/, Int_t /*nbins*/, const Double_t */*xbins*/, Double_t /*sigmaMax*/){
1323 
1324  printf("SigmaCut with variable binsize, Not yet implemented\n");
1325 
1326  return 0;
1327 }
1328 
1329 
1330 TH1F* AliTPCCalibViewer::Integrate(TH1F *const histogram, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
1341 
1342 
1343  Float_t *array = histogram->GetArray();
1344  Int_t nbins = histogram->GetXaxis()->GetNbins();
1345  Float_t binLow = histogram->GetXaxis()->GetXmin();
1346  Float_t binUp = histogram->GetXaxis()->GetXmax();
1347  return AliTPCCalibViewer::Integrate(nbins, array, nbins, binLow, binUp, mean, sigma, sigmaMax, sigmaStep);
1348 }
1349 
1350 
1351 TH1F* AliTPCCalibViewer::Integrate(Int_t n, const Float_t *const array, Int_t nbins, Float_t binLow, Float_t binUp, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
1358 
1359  Bool_t givenUnits = kTRUE;
1360  if (sigma != 0 && sigmaMax != 0) givenUnits = kFALSE;
1361  if (givenUnits) {
1362  sigma = 1;
1363  sigmaMax = (binUp - binLow) / 2.;
1364  }
1365 
1366  Float_t binWidth = (binUp-binLow)/(nbins - 1);
1367  if (sigmaStep <= 0) sigmaStep = binWidth;
1368  Int_t kbins = (Int_t)(sigmaMax * sigma / sigmaStep) + 1; // + 1 due to overflow bin in histograms
1369  Float_t kbinLow = givenUnits ? binLow : -sigmaMax;
1370  Float_t kbinUp = givenUnits ? binUp : sigmaMax;
1371  TH1F *hist = 0;
1372  if (givenUnits) hist = new TH1F("integratedHisto","Integrated Histogram; Given x; Fraction of included data", kbins, kbinLow, kbinUp);
1373  if (!givenUnits) hist = new TH1F("integratedHisto","Integrated Histogram; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
1374  hist->SetDirectory(0);
1375  hist->Reset();
1376 
1377  // calculate normalization
1378  // printf("calculating normalization, integrating from bin 1 to %i \n", n);
1379  Double_t normalization = 0;
1380  for (Int_t i = 1; i <= n; i++) {
1381  normalization += array[i];
1382  }
1383  // printf("normalization: %f \n", normalization);
1384 
1385  // given units: units from given histogram
1386  // sigma units: in units of sigma
1387  // iDelta: integrate in interval (mean +- iDelta), given units
1388  // x: ofset from mean for integration, given units
1389  // hist: needs
1390 
1391  // fill histogram
1392  for (Float_t iDelta = mean - sigmaMax * sigma; iDelta <= mean + sigmaMax * sigma; iDelta += sigmaStep) {
1393  // integrate array
1394  Double_t value = 0;
1395  for (Float_t x = mean - sigmaMax * sigma; x <= iDelta; x += binWidth) {
1396  value += (x <= binUp && x >= binLow) ? array[GetBin(x, nbins, binLow, binUp)] : 0;
1397  }
1398  if (value / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail +++ \n", value, normalization);
1399  if (value / normalization > 100) return hist;
1400  Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
1401  // printf("first integration bin: %i, last integration bin: %i \n", GetBin(mean - sigmaMax * sigma, nbins, binLow, binUp), GetBin(iDelta, nbins, binLow, binUp));
1402  // printf("value: %f, normalization: %f, normalized value: %f, iDelta: %f, Bin: %i \n", value, normalization, value/normalization, iDelta, bin);
1403  value = (value / normalization);
1404  hist->SetBinContent(bin, value);
1405  }
1406  return hist;
1407 }
1408 
1409 
1410 
1411 
1412 
1414 // end of Array tools //
1416 
1417 
1418 
1419 //_____________________________________________________________________________
1420 AliTPCCalPad* AliTPCCalibViewer::GetCalPadOld(const char* desiredData, const char* cuts, const char* calPadName) const {
1425 
1426  TString drawStr(desiredData);
1427  drawStr.Append(":channel");
1428  drawStr.Append(fAbbreviation);
1429  AliTPCCalPad * createdCalPad = new AliTPCCalPad(calPadName, calPadName);
1430  Int_t entries = 0;
1431  for (Int_t sec = 0; sec < 72; sec++) {
1432  AliTPCCalROC * roc = createdCalPad->GetCalROC(sec);
1433  entries = EasyDraw1D(drawStr.Data(), (Int_t)sec, cuts, "goff");
1434  if (entries == -1) return 0;
1435  const Double_t *pchannel = fTree->GetV2();
1436  const Double_t *pvalue = fTree->GetV1();
1437  for (Int_t i = 0; i < entries; i++)
1438  roc->SetValue((UInt_t)(pchannel[i]), (Float_t)(pvalue[i]));
1439  }
1440  return createdCalPad;
1441 }
1442 
1443 
1444 //_____________________________________________________________________________
1445 AliTPCCalPad* AliTPCCalibViewer::GetCalPad(const char* desiredData, const char* cuts, const char* calPadName) const {
1450 
1451  TString drawStr(desiredData);
1452  drawStr.Append(":channel.fElements:sector");
1453  AliTPCCalPad * createdCalPad = new AliTPCCalPad(calPadName, calPadName);
1454  //
1455  Int_t entries = fTree->Draw(drawStr, cuts,"goff");
1456  const Double_t *pvalue = fTree->GetV1();
1457  const Double_t *pchannel = fTree->GetV2();
1458  const Double_t *psector = fTree->GetV3();
1459 
1460  for (Int_t ientry=0; ientry<entries; ientry++){
1461  Int_t sector= TMath::Nint(psector[ientry]);
1462  AliTPCCalROC * roc = createdCalPad->GetCalROC(sector);
1463  if (roc) roc->SetValue((UInt_t)(pchannel[ientry]), (Float_t)(pvalue[ientry]));
1464  }
1465 
1466  // for (Int_t sec = 0; sec < 72; sec++) {
1467 // AliTPCCalROC * roc = createdCalPad->GetCalROC(sec);
1468 // entries = EasyDraw1D(drawStr.Data(), (Int_t)sec, cuts, "goff");
1469 // if (entries == -1) return 0;
1470 // for (Int_t i = 0; i < entries; i++)
1471 // roc->SetValue((UInt_t)(pchannel[i]), (Float_t)(pvalue[i]));
1472 // }
1473  return createdCalPad;
1474 }
1475 
1476 //_____________________________________________________________________________
1477 AliTPCCalROC* AliTPCCalibViewer::GetCalROC(const char* desiredData, UInt_t sector, const char* cuts) const {
1481 
1482  TString drawStr(desiredData);
1483  drawStr.Append(":channel");
1484  drawStr.Append(fAbbreviation);
1485  Int_t entries = EasyDraw1D(drawStr.Data(), (Int_t)sector, cuts, "goff");
1486  if (entries == -1) return 0;
1487  AliTPCCalROC * createdROC = new AliTPCCalROC(sector);
1488  for (Int_t i = 0; i < entries; i++)
1489  createdROC->SetValue((UInt_t)(fTree->GetV2()[i]), fTree->GetV1()[i]);
1490  return createdROC;
1491 }
1492 
1493 
1497 
1498  TObjArray* arr = new TObjArray();
1499  TObjString* str = 0;
1500  if (!fTree) return 0;
1501  Int_t nentries = fTree->GetListOfBranches()->GetEntries();
1502  for (Int_t i = 0; i < nentries; i++) {
1503  str = new TObjString(fTree->GetListOfBranches()->At(i)->GetName());
1504  str->String().ReplaceAll("_Median", "");
1505  str->String().ReplaceAll("_Mean", "");
1506  str->String().ReplaceAll("_RMS", "");
1507  str->String().ReplaceAll("_LTM", "");
1508  str->String().ReplaceAll("_OutlierCutted", "");
1509  str->String().ReplaceAll(".", "");
1510  if (!arr->FindObject(str) &&
1511  !(str->String() == "channel" || str->String() == "gx" || str->String() == "gy" ||
1512  str->String() == "lx" || str->String() == "ly" || str->String() == "pad" ||
1513  str->String() == "row" || str->String() == "rpad" || str->String() == "sector" ))
1514  arr->Add(str);
1515  }
1516 
1517  // loop over all friends (if there are some) and add them to the list
1518  if (fTree->GetListOfFriends()) {
1519  for (Int_t ifriend = 0; ifriend < fTree->GetListOfFriends()->GetEntries(); ifriend++){
1520  // printf("iterating through friendlist, currently at %i\n", ifriend);
1521  // printf("working with %s\n", fTree->GetListOfFriends()->At(ifriend)->ClassName());
1522  if (TString(fTree->GetListOfFriends()->At(ifriend)->ClassName()) != "TFriendElement") continue; // no friendElement found
1523  TFriendElement *friendElement = (TFriendElement*)fTree->GetListOfFriends()->At(ifriend);
1524  if (friendElement->GetTree() == 0) continue; // no tree found in friendElement
1525  // printf("friend found \n");
1526  for (Int_t i = 0; i < friendElement->GetTree()->GetListOfBranches()->GetEntries(); i++) {
1527  // printf("iterating through friendelement entries, currently at %i\n", i);
1528  str = new TObjString(friendElement->GetTree()->GetListOfBranches()->At(i)->GetName());
1529  str->String().ReplaceAll("_Median", "");
1530  str->String().ReplaceAll("_Mean", "");
1531  str->String().ReplaceAll("_RMS", "");
1532  str->String().ReplaceAll("_LTM", "");
1533  str->String().ReplaceAll("_OutlierCutted", "");
1534  str->String().ReplaceAll(".", "");
1535  if (!(str->String() == "channel" || str->String() == "gx" || str->String() == "gy" ||
1536  str->String() == "lx" || str->String() == "ly" || str->String() == "pad" ||
1537  str->String() == "row" || str->String() == "rpad" || str->String() == "sector" )){
1538  // insert "<friendName>." at the beginning: (<friendName> is per default "R")
1539  str->String().Insert(0, ".");
1540  str->String().Insert(0, friendElement->GetName());
1541  if (!arr->FindObject(str)) arr->Add(str);
1542  // printf("added string %s \n", str->String().Data());
1543  }
1544  }
1545  }
1546  } // if (fTree->GetListOfFriends())
1547 
1548  arr->Sort();
1549 // ((TFriendElement*)gui->GetViewer()->GetTree()->GetListOfFriends()->At(0))->GetTree()->GetListOfBranches()->At(0)->GetName()
1550 // ((TFriendElement*)gui->GetViewer()->GetTree()->GetListOfFriends()->At(0))->GetTree()->GetListOfBranches()
1551 
1552 
1553  if (printList) {
1554  TIterator* iter = arr->MakeIterator();
1555  iter->Reset();
1556  TObjString* currentStr = 0;
1557  while ( (currentStr = (TObjString*)(iter->Next())) ) {
1558  std::cout << currentStr->GetString().Data() << std::endl;
1559  }
1560  delete iter;
1561  }
1562  return arr;
1563 }
1564 
1565 
1569 
1570  TObjArray* arr = new TObjArray();
1571  arr->Add(new TObjString("_Mean"));
1572  arr->Add(new TObjString("_Mean_OutlierCutted"));
1573  arr->Add(new TObjString("_Median"));
1574  arr->Add(new TObjString("_Median_OutlierCutted"));
1575  arr->Add(new TObjString("_LTM"));
1576  arr->Add(new TObjString("_LTM_OutlierCutted"));
1577  arr->Add(new TObjString(Form("LFitIntern_4_8%s", fAppendString.Data())));
1578  arr->Add(new TObjString(Form("GFitIntern_Lin%s", fAppendString.Data())));
1579  arr->Add(new TObjString(Form("GFitIntern_Par%s", fAppendString.Data())));
1580  arr->Add(new TObjString("FitLinLocal"));
1581  arr->Add(new TObjString("FitLinGlobal"));
1582  arr->Add(new TObjString("FitParLocal"));
1583  arr->Add(new TObjString("FitParGlobal"));
1584 
1585  if (printList) {
1586  TIterator* iter = arr->MakeIterator();
1587  iter->Reset();
1588  TObjString* currentStr = 0;
1589  while ((currentStr = (TObjString*)(iter->Next()))) {
1590  std::cout << currentStr->GetString().Data() << std::endl;
1591  }
1592  delete iter;
1593  }
1594  return arr;
1595 }
1596 
1597 
1598 TFriendElement* AliTPCCalibViewer::AddReferenceTree(const char* filename, const char* treename, const char* refname){
1601 
1602  TFile *file = new TFile(filename);
1603  fListOfObjectsToBeDeleted->Add(file);
1604  TTree * tree = (TTree*)file->Get(treename);
1605  return AddFriend(tree, refname);
1606 }
1607 
1608 
1612 
1613  TObjArray *listOfCalPads = GetListOfVariables();
1614  TObjArray *calPadsArray = new TObjArray();
1615  Int_t numberOfCalPads = listOfCalPads->GetEntries();
1616  for (Int_t i = 0; i < numberOfCalPads; i++) {
1617  std::cout << "Creating calPad " << (i+1) << " of " << numberOfCalPads << "\r" << std::flush;
1618  char* calPadName = (char*)((TObjString*)(listOfCalPads->At(i)))->GetString().Data();
1619  TString drawCommand = ((TObjString*)(listOfCalPads->At(i)))->GetString();
1620  drawCommand.Append(fAbbreviation.Data());
1621  AliTPCCalPad* calPad = GetCalPad(drawCommand.Data(), "", calPadName);
1622  calPadsArray->Add(calPad);
1623  }
1624  std::cout << std::endl;
1625  listOfCalPads->Delete();
1626  delete listOfCalPads;
1627  return calPadsArray;
1628 }
1629 
1630 
1631 TString* AliTPCCalibViewer::Fit(const char* drawCommand, const char* formula, const char* cuts, Double_t & chi2, TVectorD &fitParam, TMatrixD &covMatrix){
1635 
1636  TString formulaStr(formula);
1637  TString drawStr(drawCommand);
1638  TString cutStr(cuts);
1639 
1640  // abbreviations:
1641  drawStr.ReplaceAll(fAbbreviation, fAppendString);
1642  cutStr.ReplaceAll(fAbbreviation, fAppendString);
1643  formulaStr.ReplaceAll(fAbbreviation, fAppendString);
1644 
1645  formulaStr.ReplaceAll("++", fAbbreviation);
1646  TObjArray* formulaTokens = formulaStr.Tokenize(fAbbreviation.Data());
1647  Int_t dim = formulaTokens->GetEntriesFast();
1648 
1649  fitParam.ResizeTo(dim);
1650  covMatrix.ResizeTo(dim,dim);
1651 
1652  TLinearFitter* fitter = new TLinearFitter(dim+1, Form("hyp%d",dim));
1653  fitter->StoreData(kTRUE);
1654  fitter->ClearPoints();
1655 
1656  Int_t entries = Draw(drawStr.Data(), cutStr.Data(), "goff");
1657  if (entries == -1) {
1658  delete formulaTokens;
1659  return new TString("An ERROR has occured during fitting!");
1660  }
1661  Double_t **values = new Double_t*[dim+1] ;
1662 
1663  for (Int_t i = 0; i < dim + 1; i++){
1664  Int_t centries = 0;
1665  if (i < dim) centries = fTree->Draw(((TObjString*)formulaTokens->At(i))->GetName(), cutStr.Data(), "goff");
1666  else centries = fTree->Draw(drawStr.Data(), cutStr.Data(), "goff");
1667 
1668  if (entries != centries) {
1669  delete [] values;
1670  delete formulaTokens;
1671  return new TString("An ERROR has occured during fitting!");
1672  }
1673  values[i] = new Double_t[entries];
1674  memcpy(values[i], fTree->GetV1(), entries*sizeof(Double_t));
1675  }
1676 
1677  // add points to the fitter
1678  for (Int_t i = 0; i < entries; i++){
1679  Double_t x[1000];
1680  for (Int_t j=0; j<dim;j++) x[j]=values[j][i];
1681  fitter->AddPoint(x, values[dim][i], 1);
1682  }
1683 
1684  fitter->Eval();
1685  fitter->GetParameters(fitParam);
1686  fitter->GetCovarianceMatrix(covMatrix);
1687  chi2 = fitter->GetChisquare();
1688 // chi2 = chi2;
1689 
1690  TString *preturnFormula = new TString(Form("( %e+",fitParam[0])), &returnFormula = *preturnFormula;
1691 
1692  for (Int_t iparam = 0; iparam < dim; iparam++) {
1693  returnFormula.Append(Form("%s*(%e)",((TObjString*)formulaTokens->At(iparam))->GetName(),fitParam[iparam+1]));
1694  if (iparam < dim-1) returnFormula.Append("+");
1695  }
1696  returnFormula.Append(" )");
1697  delete formulaTokens;
1698  delete fitter;
1699  for (Int_t i = 0; i < dim + 1; i++){
1700  delete [] values[i];
1701  }
1702  delete [] values;
1703  return preturnFormula;
1704 }
1705 
1706 
1707 void AliTPCCalibViewer::MakeTreeWithObjects(const char *fileName, const TObjArray *const array, const char * mapFileName) {
1712 
1713  AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
1714 
1715  TObjArray* mapIROCs = 0;
1716  TObjArray* mapOROCs = 0;
1717  TVectorF *mapIROCArray = 0;
1718  TVectorF *mapOROCArray = 0;
1719  Int_t mapEntries = 0;
1720  TString* mapNames = 0;
1721 
1722  if (mapFileName) {
1723  TFile mapFile(mapFileName, "read");
1724 
1725  TList* listOfROCs = mapFile.GetListOfKeys();
1726  mapEntries = listOfROCs->GetEntries()/2;
1727  mapIROCs = new TObjArray(mapEntries*2);
1728  mapOROCs = new TObjArray(mapEntries*2);
1729  mapIROCArray = new TVectorF[mapEntries];
1730  mapOROCArray = new TVectorF[mapEntries];
1731 
1732  mapNames = new TString[mapEntries];
1733  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1734  TString rocName(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
1735  rocName.Remove(rocName.Length()-4, 4);
1736  mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "IROC").Data()), ivalue);
1737  mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "OROC").Data()), ivalue);
1738  mapNames[ivalue].Append(rocName);
1739  }
1740 
1741  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1742  mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
1743  mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
1744 
1745  for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
1746  (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
1747  for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
1748  (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
1749  }
1750 
1751  } // if (mapFileName)
1752 
1753  TTreeSRedirector cstream(fileName);
1754  Int_t arrayEntries = array->GetEntries();
1755 
1756  // Read names of AliTPCCalPads and save them in names[]
1757  TString* names = new TString[arrayEntries];
1758  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1759  names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
1760 
1761  for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
1762 
1763  TVectorF *vectorArray = new TVectorF[arrayEntries];
1764  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1765  vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1766 
1767 
1768  //
1769  // fill vectors of variable per pad
1770  //
1771  TVectorF *posArray = new TVectorF[8];
1772  for (Int_t ivalue = 0; ivalue < 8; ivalue++)
1773  posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1774 
1775  Float_t posG[3] = {0};
1776  Float_t posL[3] = {0};
1777  Int_t ichannel = 0;
1778  for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
1779  for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
1780  tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
1781  tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
1782  posArray[0][ichannel] = irow;
1783  posArray[1][ichannel] = ipad;
1784  posArray[2][ichannel] = posL[0];
1785  posArray[3][ichannel] = posL[1];
1786  posArray[4][ichannel] = posG[0];
1787  posArray[5][ichannel] = posG[1];
1788  posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
1789  posArray[7][ichannel] = ichannel;
1790 
1791  // loop over array containing AliTPCCalPads
1792  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1793  AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1794  AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1795  if (calROC)
1796  (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
1797  else
1798  (vectorArray[ivalue])[ichannel] = 0;
1799  }
1800  ichannel++;
1801  }
1802  }
1803  AliTPCCalROC dummyROC(0);
1804  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1805  AliTPCCalROC *roc = ((AliTPCCalPad*)array->At(ivalue))->GetCalROC(isector);
1806  if (!roc) roc = &dummyROC;
1807  cstream << "calPads" <<
1808  (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
1809  cstream << "calPads" <<
1810  (Char_t*)((names[ivalue] + "Pad.=").Data()) << roc;
1811  }
1812 
1813  if (mapFileName) {
1814  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1815  if (isector < 36)
1816  cstream << "calPads" <<
1817  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
1818  else
1819  cstream << "calPads" <<
1820  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
1821  }
1822  }
1823 
1824  cstream << "calPads" <<
1825  "sector=" << isector;
1826 
1827  cstream << "calPads" <<
1828  "row.=" << &posArray[0] <<
1829  "pad.=" << &posArray[1] <<
1830  "lx.=" << &posArray[2] <<
1831  "ly.=" << &posArray[3] <<
1832  "gx.=" << &posArray[4] <<
1833  "gy.=" << &posArray[5] <<
1834  "rpad.=" << &posArray[6] <<
1835  "channel.=" << &posArray[7];
1836 
1837  cstream << "calPads" <<
1838  "\n";
1839 
1840  delete[] posArray;
1841  delete[] vectorArray;
1842  } //for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++)
1843 
1844  delete[] names;
1845  if (mapFileName) {
1846  delete mapIROCs;
1847  delete mapOROCs;
1848  delete[] mapIROCArray;
1849  delete[] mapOROCArray;
1850  delete[] mapNames;
1851  }
1852 }
1853 
1854 
1855 void AliTPCCalibViewer::MakeTree(const char * fileName, TObjArray * array, const char * mapFileName, AliTPCCalPad *const outlierPad, Float_t ltmFraction) {
1866 
1867  AliTPCROC* tpcROCinstance = AliTPCROC::Instance();
1868 
1869  TObjArray* mapIROCs = 0;
1870  TObjArray* mapOROCs = 0;
1871  TVectorF *mapIROCArray = 0;
1872  TVectorF *mapOROCArray = 0;
1873  Int_t mapEntries = 0;
1874  TString* mapNames = 0;
1875 
1876  if (mapFileName) {
1877  TFile mapFile(mapFileName, "read");
1878 
1879  TList* listOfROCs = mapFile.GetListOfKeys();
1880  mapEntries = listOfROCs->GetEntries()/2;
1881  mapIROCs = new TObjArray(mapEntries*2);
1882  mapOROCs = new TObjArray(mapEntries*2);
1883  mapIROCArray = new TVectorF[mapEntries];
1884  mapOROCArray = new TVectorF[mapEntries];
1885 
1886  mapNames = new TString[mapEntries];
1887  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1888  TString rocName(((TKey*)(listOfROCs->At(ivalue*2)))->GetName());
1889  rocName.Remove(rocName.Length()-4, 4);
1890  mapIROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "IROC").Data()), ivalue);
1891  mapOROCs->AddAt((AliTPCCalROC*)mapFile.Get((rocName + "OROC").Data()), ivalue);
1892  mapNames[ivalue].Append(rocName);
1893  }
1894 
1895  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
1896  mapIROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(0));
1897  mapOROCArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(36));
1898 
1899  for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(0); ichannel++)
1900  (mapIROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapIROCs->At(ivalue)))->GetValue(ichannel);
1901  for (UInt_t ichannel = 0; ichannel < tpcROCinstance->GetNChannels(36); ichannel++)
1902  (mapOROCArray[ivalue])[ichannel] = ((AliTPCCalROC*)(mapOROCs->At(ivalue)))->GetValue(ichannel);
1903  }
1904 
1905  } // if (mapFileName)
1906 
1907  TTreeSRedirector cstream(fileName,"recreate");
1908  Int_t arrayEntries = 0;
1909  if (array) arrayEntries = array->GetEntries();
1910 
1911  TString* names = new TString[arrayEntries];
1912  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++)
1913  names[ivalue].Append(((AliTPCCalPad*)array->At(ivalue))->GetName());
1914 
1915  for (UInt_t isector = 0; isector < tpcROCinstance->GetNSectors(); isector++) {
1916  //
1917  // get statistic for given sector
1918  //
1919  TVectorF median(arrayEntries);
1920  TVectorF mean(arrayEntries);
1921  TVectorF rms(arrayEntries);
1922  TVectorF ltm(arrayEntries);
1923  TVectorF ltmrms(arrayEntries);
1924  TVectorF medianWithOut(arrayEntries);
1925  TVectorF meanWithOut(arrayEntries);
1926  TVectorF rmsWithOut(arrayEntries);
1927  TVectorF ltmWithOut(arrayEntries);
1928  TVectorF ltmrmsWithOut(arrayEntries);
1929 
1930  TVectorF *vectorArray = new TVectorF[arrayEntries];
1931  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++){
1932  vectorArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1933  vectorArray[ivalue].SetUniqueID(array->UncheckedAt(ivalue)->GetUniqueID());
1934 // printf("UniqueID: %d\n",vectorArray[ivalue].GetUniqueID());
1935  }
1936 
1937  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1938  AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1939  AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1940  AliTPCCalROC* outlierROC = 0;
1941  if (outlierPad) outlierROC = outlierPad->GetCalROC(isector);
1942  if (calROC) {
1943  median[ivalue] = calROC->GetMedian();
1944  mean[ivalue] = calROC->GetMean();
1945  rms[ivalue] = calROC->GetRMS();
1946  Double_t ltmrmsValue = 0;
1947  ltm[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction);
1948  ltmrms[ivalue] = ltmrmsValue;
1949  if (outlierROC) {
1950  medianWithOut[ivalue] = calROC->GetMedian(outlierROC);
1951  meanWithOut[ivalue] = calROC->GetMean(outlierROC);
1952  rmsWithOut[ivalue] = calROC->GetRMS(outlierROC);
1953  ltmrmsValue = 0;
1954  ltmWithOut[ivalue] = calROC->GetLTM(&ltmrmsValue, ltmFraction, outlierROC);
1955  ltmrmsWithOut[ivalue] = ltmrmsValue;
1956  }
1957  }
1958  else {
1959  median[ivalue] = 0.;
1960  mean[ivalue] = 0.;
1961  rms[ivalue] = 0.;
1962  ltm[ivalue] = 0.;
1963  ltmrms[ivalue] = 0.;
1964  medianWithOut[ivalue] = 0.;
1965  meanWithOut[ivalue] = 0.;
1966  rmsWithOut[ivalue] = 0.;
1967  ltmWithOut[ivalue] = 0.;
1968  ltmrmsWithOut[ivalue] = 0.;
1969  }
1970  }
1971 
1972  //
1973  // fill vectors of variable per pad
1974  //
1975  TVectorF *posArray = new TVectorF[8];
1976  for (Int_t ivalue = 0; ivalue < 8; ivalue++)
1977  posArray[ivalue].ResizeTo(tpcROCinstance->GetNChannels(isector));
1978 
1979  Float_t posG[3] = {0};
1980  Float_t posL[3] = {0};
1981  Int_t ichannel = 0;
1982  for (UInt_t irow = 0; irow < tpcROCinstance->GetNRows(isector); irow++) {
1983  for (UInt_t ipad = 0; ipad < tpcROCinstance->GetNPads(isector, irow); ipad++) {
1984  tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
1985  tpcROCinstance->GetPositionGlobal(isector, irow, ipad, posG);
1986  posArray[0][ichannel] = irow;
1987  posArray[1][ichannel] = ipad;
1988  posArray[2][ichannel] = posL[0];
1989  posArray[3][ichannel] = posL[1];
1990  posArray[4][ichannel] = posG[0];
1991  posArray[5][ichannel] = posG[1];
1992  posArray[6][ichannel] = (Int_t)(ipad - (Double_t)(tpcROCinstance->GetNPads(isector, irow))/2);
1993  posArray[7][ichannel] = ichannel;
1994 
1995  // loop over array containing AliTPCCalPads
1996  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
1997  AliTPCCalPad* calPad = (AliTPCCalPad*) array->At(ivalue);
1998  AliTPCCalROC* calROC = calPad->GetCalROC(isector);
1999  if (calROC)
2000  (vectorArray[ivalue])[ichannel] = calROC->GetValue(irow, ipad);
2001  else
2002  (vectorArray[ivalue])[ichannel] = 0;
2003  }
2004  ichannel++;
2005  }
2006  }
2007 
2008  cstream << "calPads" <<
2009  "sector=" << isector;
2010 
2011  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
2012  if (outlierPad==0) {
2013  cstream << "calPads" <<
2014  (Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
2015  (Char_t*)((names[ivalue] + "_Mean=").Data()) << mean[ivalue] <<
2016  (Char_t*)((names[ivalue] + "_RMS=").Data()) << rms[ivalue] <<
2017  (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltm[ivalue] <<
2018  (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrms[ivalue];
2019  }
2020  if (outlierPad) {
2021  cstream << "calPads" <<
2022  (Char_t*)((names[ivalue] + "_Median=").Data()) << medianWithOut[ivalue] <<
2023  (Char_t*)((names[ivalue] + "_Mean").Data()) << meanWithOut[ivalue] <<
2024  (Char_t*)((names[ivalue] + "_RMS=").Data()) << rmsWithOut[ivalue] <<
2025  (Char_t*)((names[ivalue] + "_LTM=").Data()) << ltmWithOut[ivalue] <<
2026  (Char_t*)((names[ivalue] + "_RMS_LTM=").Data()) << ltmrmsWithOut[ivalue];
2027  }
2028  //timestamp and run, if given in title
2029  /* TString title(((AliTPCCalPad*) array->At(ivalue))->GetTitle());
2030  TObjArray *arrtitle=title.Tokenize(",");
2031  Int_t run=-1;
2032  UInt_t time=0;
2033  TIter next(arrtitle);
2034  TObject *o=0;
2035  while ( (o=next()) ){
2036  TString &entry=((TObjString*)o)->GetString();
2037  entry.Remove(TString::kBoth,' ');
2038  entry.Remove(TString::kBoth,'\t');
2039  if (entry.BeginsWith("Run:")) {
2040  run=entry(4,entry.Length()).Atoi();
2041  } else if (entry.BeginsWith("Time:")) {
2042  time=entry(6,entry.Length()).Atoi();
2043  }
2044  }
2045  delete arrtitle;*/
2046 
2047  }
2048 
2049  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
2050  cstream << "calPads" <<
2051  (Char_t*)((names[ivalue] + ".=").Data()) << &vectorArray[ivalue];
2052  }
2053 
2054  if (mapFileName) {
2055  for (Int_t ivalue = 0; ivalue < mapEntries; ivalue++) {
2056  if (isector < 36)
2057  cstream << "calPads" <<
2058  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapIROCArray[ivalue];
2059  else
2060  cstream << "calPads" <<
2061  (Char_t*)((mapNames[ivalue] + ".=").Data()) << &mapOROCArray[ivalue];
2062  }
2063  }
2064 
2065  cstream << "calPads" <<
2066  "row.=" << &posArray[0] <<
2067  "pad.=" << &posArray[1] <<
2068  "lx.=" << &posArray[2] <<
2069  "ly.=" << &posArray[3] <<
2070  "gx.=" << &posArray[4] <<
2071  "gy.=" << &posArray[5] <<
2072  "rpad.=" << &posArray[6] <<
2073  "channel.=" << &posArray[7];
2074 
2075  cstream << "calPads" <<
2076  "\n";
2077 
2078  delete[] posArray;
2079  delete[] vectorArray;
2080  }
2081  TTree * tree = (cstream << "calPads").GetTree();
2082  MakeCalPadAliases(tree);
2083  for (Int_t ivalue = 0; ivalue < arrayEntries; ivalue++) {
2084  //(Char_t*)((names[ivalue] + "_Median=").Data()) << median[ivalue] <<
2085  tree->SetAlias(TString::Format("%s_LTMRatio",names[ivalue].Data()).Data(), TString::Format("%s.fElements/%s_LTM",names[ivalue].Data(),names[ivalue].Data()).Data());
2086  tree->SetAlias(TString::Format("%s_MeanRatio",names[ivalue].Data()).Data(), TString::Format("%s.fElements/%s_Mean",names[ivalue].Data(),names[ivalue].Data()).Data());
2087  tree->SetAlias(TString::Format("%s_MedianRatio",names[ivalue].Data()).Data(), TString::Format("%s.fElements/%s_Median",names[ivalue].Data(),names[ivalue].Data()).Data());
2088  tree->SetAlias(names[ivalue].Data(), TString::Format("%s.fElements",names[ivalue].Data()).Data());
2089  }
2090  delete[] names;
2091  if (mapFileName) {
2092  delete mapIROCs;
2093  delete mapOROCs;
2094  delete[] mapIROCArray;
2095  delete[] mapOROCArray;
2096  delete[] mapNames;
2097  }
2098 }
2101 
2102  tree->SetAlias("Pad0","MapPadLength.fElements==7.5"); // pad types
2103  tree->SetAlias("Pad1","MapPadLength.fElements==10.0");
2104  tree->SetAlias("Pad2","MapPadLength.fElements==15.0");
2105  tree->SetAlias("ly","(ly.fElements+((sector<36)*(2*((sector%36)<18)-1)*0.2)+((sector>=36)*(2*((sector%36)<18)-1)*0.3))"); // correction for bug in tpcROCinstance->GetPositionLocal(isector, irow, ipad, posL);
2106  tree->SetAlias("dRowNorm0","(row.fElements/32-1)"); // normalized distance to the center of the chamber in lx (-1,1)
2107  tree->SetAlias("dRowNorm1","(row.fElements/32-1)"); //
2108  tree->SetAlias("dRowNorm2","((row.fElements-63)/16-1)"); //
2109  tree->SetAlias("dRowNorm","(row.fElements<63)*(row.fElements/32-1)+(row.fElements>63)*((row.fElements-63)/16-1)");
2110  tree->SetAlias("dPhiNorm","(ly/lx.fElements)/(pi/18)"); // normalized distance to the center of the chamber in phi (-1,1)
2111  //
2112  tree->SetAlias("fsector","(sector%36)+(0.5+9*atan2(ly,lx.fElements)/pi)"); // float sector number
2113  tree->SetAlias("posEdge","((pi/18.)-abs(atan(ly/lx.fElements)))*lx.fElements"); //distance to the edge
2114  tree->SetAlias("lphi","atan(ly/lx.fElements)"); //local phi angle
2115 }
2116 
2117 void AliTPCCalibViewer::MakeTree(const char *outPutFileName, const Char_t *inputFileName, AliTPCCalPad *outlierPad, Float_t ltmFraction, const char *mapFileName ){
2136 
2137  TObjArray objArray;
2138  CreateObjectList(inputFileName, &objArray);
2139  MakeTree(outPutFileName, &objArray, mapFileName, outlierPad, ltmFraction);
2140 }
2141 
2142 
2143 void AliTPCCalibViewer::CreateObjectList(const Char_t *filename, TObjArray *calibObjects){
2161 
2162  if ( calibObjects == 0x0 ) return;
2163  ifstream in;
2164  in.open(filename);
2165  if ( !in.is_open() ){
2166  fprintf(stderr,"Error: cannot open list file '%s'", filename);
2167  return;
2168  }
2169 
2170  AliTPCCalPad *calPad=0x0;
2171 
2172  TString sFile;
2173  sFile.ReadFile(in);
2174  in.close();
2175 
2176  TObjArray *arrFileLine = sFile.Tokenize("\n");
2177  TIter nextLine(arrFileLine);
2178 
2179  TObjString *sObjLine = 0x0;
2180  while ( (sObjLine = (TObjString*)nextLine()) ){
2181  TString sLine(sObjLine->GetString());
2182 
2183  TObjArray *arrCol = sLine.Tokenize("\t");
2184  Int_t nCols = arrCol->GetEntriesFast();
2185 
2186  TObjString *sObjType = (TObjString*)(arrCol->At(0));
2187  TObjString *sObjFileName = (TObjString*)(arrCol->At(1));
2188  TObjString *sObjName = 0x0;
2189 
2190  if ( !sObjType || !sObjFileName ) continue;
2191  TString sType(sObjType->GetString());
2192  TString sFileName(sObjFileName->GetString());
2193  printf("Type %s, opening %s \n", sType.Data(), sFileName.Data());
2194  TFile *fIn = TFile::Open(sFileName);
2195  if ( !fIn ){
2196  fprintf(stderr,"File not found: '%s'", sFileName.Data());
2197  continue;
2198  }
2199 
2200  if ( sType == "CE" ){ // next three colums are the names for CETmean, CEQmean and CETrms
2201  AliTPCCalibCE *ce = (AliTPCCalibCE*)fIn->Get("AliTPCCalibCE");
2202  calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadT0());
2203  if (nCols > 2) { // check, if the name is provided
2204  sObjName = (TObjString*)(arrCol->At(2));
2205  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2206  }
2207  else calPad->SetNameTitle("CETmean","CETmean");
2208  calibObjects->Add(calPad);
2209 
2210  calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadQ());
2211  if (nCols > 3) { // check, if the name is provided
2212  sObjName = (TObjString*)(arrCol->At(3));
2213  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2214  }
2215  else calPad->SetNameTitle("CEQmean","CEQmean");
2216  calibObjects->Add(calPad);
2217 
2218  calPad = new AliTPCCalPad((TObjArray*)ce->GetCalPadRMS());
2219  if (nCols > 4) { // check, if the name is provided
2220  sObjName = (TObjString*)(arrCol->At(4));
2221  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2222  }
2223  else calPad->SetNameTitle("CETrms","CETrms");
2224  calibObjects->Add(calPad);
2225 
2226  } else if ( sType == "Pulser") {
2227  AliTPCCalibPulser *sig = (AliTPCCalibPulser*)fIn->Get("AliTPCCalibPulser");
2228 
2229  calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadT0());
2230  if (nCols > 2) {
2231  sObjName = (TObjString*)(arrCol->At(2));
2232  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2233  }
2234  else calPad->SetNameTitle("PulserTmean","PulserTmean");
2235  calibObjects->Add(calPad);
2236 
2237  calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadQ());
2238  if (nCols > 3) {
2239  sObjName = (TObjString*)(arrCol->At(3));
2240  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2241  }
2242  else calPad->SetNameTitle("PulserQmean","PulserQmean");
2243  calibObjects->Add(calPad);
2244 
2245  calPad = new AliTPCCalPad((TObjArray*)sig->GetCalPadRMS());
2246  if (nCols > 4) {
2247  sObjName = (TObjString*)(arrCol->At(4));
2248  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2249  }
2250  else calPad->SetNameTitle("PulserTrms","PulserTrms");
2251  calibObjects->Add(calPad);
2252 
2253  } else if ( sType == "Pedestals") {
2254  AliTPCCalibPedestal *ped = (AliTPCCalibPedestal*)fIn->Get("AliTPCCalibPedestal");
2255 
2256  calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadPedestal());
2257  if (nCols > 2) {
2258  sObjName = (TObjString*)(arrCol->At(2));
2259  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2260  }
2261  else calPad->SetNameTitle("Pedestals","Pedestals");
2262  calibObjects->Add(calPad);
2263 
2264  calPad = new AliTPCCalPad((TObjArray*)ped->GetCalPadRMS());
2265  if (nCols > 3) {
2266  sObjName = (TObjString*)(arrCol->At(3));
2267  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2268  }
2269  else calPad->SetNameTitle("Noise","Noise");
2270  calibObjects->Add(calPad);
2271 
2272  } else if ( sType == "CalPad") {
2273  if (nCols > 2) sObjName = (TObjString*)(arrCol->At(2));
2274  else continue;
2275  calPad = new AliTPCCalPad(*(AliTPCCalPad*)fIn->Get(sObjName->GetString().Data()));
2276  if (nCols > 3) {
2277  sObjName = (TObjString*)(arrCol->At(3));
2278  calPad->SetNameTitle(sObjName->GetString().Data(), sObjName->GetString().Data());
2279  }
2280  calibObjects->Add(calPad);
2281  } else {
2282  fprintf(stderr,"Undefined Type: '%s'",sType.Data());
2283  }
2284  delete fIn;
2285  }
2286  delete arrFileLine;
2287 }
TString fAppendString
'.fElements', stored in a TStrig
static Int_t GetBin(Float_t value, Int_t nbins, Double_t binLow, Double_t binUp)
printf("Chi2/npoints = %f\n", TMath::Sqrt(chi2/npoints))
UInt_t GetNPads(UInt_t sector, UInt_t row) const
Definition: AliTPCROC.h:29
const TObjArray * GetCalPadRMS() const
TTree * fTree
Definition: MakeTreeStat.C:55
AliTPCCalROC * GetCalROC(Int_t sector) const
Definition: AliTPCCalPad.h:42
const TObjArray * GetCalPadRMS() const
Definition: AliTPCCalibCE.h:57
TStyle * gStyle
static void MakeCalPadAliases(TTree *tree)
#define TObjArray
static Double_t GetLTM(Int_t n, const Double_t *const array, Double_t *const sigma=0, Double_t fraction=0.9)
void DrawLines(TH1F *cutHistoMean, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const
UInt_t GetNRows(UInt_t sector) const
Definition: AliTPCROC.h:27
Float_t GetValue(UInt_t row, UInt_t pad) const
Definition: AliTPCCalROC.h:33
Implementation of the TPC pedestal calibration.
const TObjArray * GetCalPadT0() const
Definition: AliTPCCalibCE.h:54
Double_t GetRMS(AliTPCCalROC *const outlierROC=0) const
void FormatHistoLabels(TH1 *histo) const
TFriendElement * AddFriend(const char *treename, const char *filename)
UInt_t GetNChannels(UInt_t sector) const
Definition: AliTPCROC.h:28
TTree * fTree
tree containing visualization data (e.g. written by AliTPCCalPad::MakeTree(...)
Class for viewing/visualizing TPC calibration data.
Int_t Integrate(const char *drawCommand, const char *sector, const char *cuts=0, Float_t sigmaMax=5, Bool_t plotMean=kTRUE, Bool_t plotMedian=kTRUE, Bool_t plotLTM=kTRUE, const char *sigmas="", Float_t sigmaStep=-1) const
void GetPositionLocal(UInt_t sector, UInt_t row, UInt_t pad, Float_t *pos)
Definition: AliTPCROC.cxx:414
TObjArray * GetListOfNormalizationVariables(Bool_t printList=kFALSE) const
TTree * tree
const char * AddAbbreviations(const Char_t *c, Bool_t printDrawCommand=kFALSE)
virtual void Draw(Option_t *opt="")
Int_t SigmaCutNew(const char *drawCommand, const char *sector, const char *cuts=0, Float_t sigmaMax=5, Bool_t plotMean=kTRUE, Bool_t plotMedian=kTRUE, Bool_t plotLTM=kTRUE, Bool_t pm=kFALSE, const char *sigmas="", Float_t sigmaStep=-1) const
Int_t SigmaCut(const char *drawCommand, Int_t sector, const char *cuts=0, Float_t sigmaMax=5, Bool_t plotMean=kTRUE, Bool_t plotMedian=kTRUE, Bool_t plotLTM=kTRUE, Bool_t pm=kFALSE, const char *sigmas="", Float_t sigmaStep=-1) const
const TObjArray * GetCalPadRMS() const
TObjArray * GetListOfVariables(Bool_t printList=kFALSE)
TFile * fFile
file that contains a calPads tree (e.g. written by AliTPCCalPad::MakeTree(...)
TObjArray * array
Definition: AnalyzeLaser.C:12
ClassImp(TPCGenInfo)
Definition: AliTPCCmpNG.C:254
Int_t EasyDraw(const char *drawCommand, const char *sector, const char *cuts=0, const char *drawOptions=0, Bool_t writeDrawCommand=kFALSE) const
Implementation of the TPC pulser calibration.
TString * Fit(const char *drawCommand, const char *formula, const char *cuts, Double_t &chi2, TVectorD &fitParam, TMatrixD &covMatrix)
Double_t chi2
Definition: AnalyzeLaser.C:7
static void CreateObjectList(const Char_t *filename, TObjArray *calibObjects)
void GetPositionGlobal(UInt_t sector, UInt_t row, UInt_t pad, Float_t *pos)
Definition: AliTPCROC.cxx:432
const TObjArray * GetCalPadT0() const
Geometry class for a single ROC.
Definition: AliTPCROC.h:14
void SetValue(UInt_t row, UInt_t pad, Float_t vd)
Definition: AliTPCCalROC.h:35
void Add(Float_t c1)
TFriendElement * AddReferenceTree(const char *filename, const char *treename="calPads", const char *refname="R")
Bool_t fTreeMustBeDeleted
decides weather the tree must be deleted in destructor or not
Double_t GetMedian(AliTPCCalROC *const outlierROC=0) const
Double_t GetMean(AliTPCCalROC *const outlierROC=0) const
TPC calibration base class for one ROC.
Definition: AliTPCCalROC.h:19
AliTPCCalROC * GetCalROC(const char *desiredData, UInt_t sector, const char *cuts="") const
TObject * htemp
Definition: PlotSys.C:37
const TObjArray * GetCalPadPedestal() const
AliTPCCalPad * GetCalPadOld(const char *desiredData, const char *cuts="", const char *calPadName="NoName") const
static void MakeTree(const char *fileName, TObjArray *array, const char *mapFileName=0, AliTPCCalPad *const outlierPad=0, Float_t ltmFraction=0.9)
AliTPCCalibViewer & operator=(const AliTPCCalibViewer &param)
virtual void Delete(Option_t *option="")
TObjArray * GetArrayOfCalPads()
AliTPCCalPad * GetCalPad(const char *desiredData, const char *cuts="", const char *calPadName="NoName") const
const TObjArray * GetCalPadQ() const
Definition: AliTPCCalibCE.h:56
Int_t EasyDraw1D(const char *drawCommand, const char *sector, const char *cuts=0, const char *drawOptions=0, Bool_t writeDrawCommand=kFALSE) const
static AliTPCROC * Instance()
Definition: AliTPCROC.cxx:34
Implementation of the TPC Central Electrode calibration.
Definition: AliTPCCalibCE.h:29
UInt_t GetNSectors() const
Definition: AliTPCROC.h:26
Double_t GetLTM(Double_t *const sigma=0, Double_t fraction=0.9, AliTPCCalROC *const outlierROC=0)
TTree * GetTree() const
static void MakeTreeWithObjects(const char *fileName, const TObjArray *const array, const char *mapFileName=0)
const TObjArray * GetCalPadQ() const
Int_t IntegrateOld(const char *drawCommand, const char *sector, const char *cuts=0, Float_t sigmaMax=5, Bool_t plotMean=kTRUE, Bool_t plotMedian=kTRUE, Bool_t plotLTM=kTRUE, const char *sigmas="", Float_t sigmaStep=-1) const
Int_t DrawHisto1D(const char *drawCommand, Int_t sector, const char *cuts=0, const char *sigmas="2;4;6", Bool_t plotMean=kTRUE, Bool_t plotMedian=kTRUE, Bool_t plotLTM=kTRUE) const
TObjArray * fListOfObjectsToBeDeleted
Objects, that will be deleted when the destructor ist called.
TString fAbbreviation
the abreviation for '.fElements'