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