STEER/CDB/AliBaseCalibViewer.cxx

   1 ///////////////////////////////////////////////////////////////////////////////
   2 //                                                                           //
   3 //  Base class for the AliTPCCalibViewer and AliTRDCalibViewer               //
   4 //  used for the calibration monitor                                         //
   5 //                                                                           //
   6 //  Authors:     Marian Ivanov (Marian.Ivanov@cern.ch)                       //
   7 //               Jens Wiechula (Jens.Wiechula@cern.ch)                       //
   8 //               Ionut Arsene  (iarsene@cern.ch)                             //
   9 //                                                                           //
  10 ///////////////////////////////////////////////////////////////////////////////
  11
  12
  13 #include <iostream>
  14 #include <fstream>
  15 #include <TString.h>
  16 #include <TRandom.h>
  17 #include <TLegend.h>
  18 #include <TLine.h>
  19 //#include <TCanvas.h>
  20 #include <TROOT.h>
  21 #include <TStyle.h>
  22 #include <TH1.h>
  23 #include <TH1F.h>
  24 #include <TMath.h>
  25 #include <THashTable.h>
  26 #include <TObjString.h>
  27 #include <TLinearFitter.h>
  28 #include <TFile.h>
  29 #include <TKey.h>
  30 #include <TGraph.h>
  31 #include <TDirectory.h>
  32 #include <TFriendElement.h>
  33
  34 #include "TTreeStream.h"
  35 #include "AliBaseCalibViewer.h"
  36
  37 ClassImp(AliBaseCalibViewer)
  38
  39 AliBaseCalibViewer::AliBaseCalibViewer()
  40 :TObject(),
  41   fTree(0),
  42   fFile(0),
  43   fListOfObjectsToBeDeleted(0),
  44   fTreeMustBeDeleted(0),
  45   fAbbreviation(0),
  46   fAppendString(0)
  47 {
  48   //
  49   // Default constructor
  50   //
  51 }
  52
  53 //_____________________________________________________________________________
  54 AliBaseCalibViewer::AliBaseCalibViewer(const AliBaseCalibViewer &c)
  55 :TObject(c),
  56   fTree(0),
  57   fFile(0),
  58   fListOfObjectsToBeDeleted(0),
  59   fTreeMustBeDeleted(0),
  60   fAbbreviation(0),
  61   fAppendString(0)
  62 {
  63   //
  64   // dummy AliBaseCalibViewer copy constructor
  65   // not yet working!!!
  66   //
  67   fTree = c.fTree;
  68   fTreeMustBeDeleted = c.fTreeMustBeDeleted;
  69   fListOfObjectsToBeDeleted = c.fListOfObjectsToBeDeleted;
  70   fAbbreviation = c.fAbbreviation;
  71   fAppendString = c.fAppendString;
  72 }
  73
  74 //_____________________________________________________________________________
  75 AliBaseCalibViewer::AliBaseCalibViewer(TTree* tree)
  76 :TObject(),
  77   fTree(0),
  78   fFile(0),
  79   fListOfObjectsToBeDeleted(0),
  80   fTreeMustBeDeleted(0),
  81   fAbbreviation(0),
  82   fAppendString(0)
  83 {
  84   //
  85   // Constructor that initializes the calibration viewer
  86   //
  87   fTree = tree;
  88   fTreeMustBeDeleted = kFALSE;
  89   fListOfObjectsToBeDeleted = new TObjArray();
  90   fAbbreviation = "~";
  91   fAppendString = ".fElements";
  92 }
  93
  94 //_____________________________________________________________________________
  95 AliBaseCalibViewer::AliBaseCalibViewer(const Char_t* fileName, const Char_t* treeName)
  96 :TObject(),
  97   fTree(0),
  98   fFile(0),
  99   fListOfObjectsToBeDeleted(0),
 100   fTreeMustBeDeleted(0),
 101   fAbbreviation(0),
 102   fAppendString(0)
 103
 104 {
 105   //
 106   // Constructor to initialize the calibration viewer
 107   // the file 'fileName' contains the tree 'treeName'
 108   //
 109   fFile = new TFile(fileName, "read");
 110   fTree = (TTree*) fFile->Get(treeName);
 111   fTreeMustBeDeleted = kTRUE;
 112   fListOfObjectsToBeDeleted = new TObjArray();
 113   fAbbreviation = "~";
 114   fAppendString = ".fElements";
 115 }
 116
 117 //____________________________________________________________________________
 118 AliBaseCalibViewer & AliBaseCalibViewer::operator =(const AliBaseCalibViewer & param)
 119 {
 120   //
 121   // assignment operator - dummy
 122   // not yet working!!!
 123   //
 124   if(&param == this) return *this;
 125   TObject::operator=(param);
 126
 127
 128   fTree = param.fTree;
 129   fTreeMustBeDeleted = param.fTreeMustBeDeleted;
 130   fListOfObjectsToBeDeleted = param.fListOfObjectsToBeDeleted;
 131   fAbbreviation = param.fAbbreviation;
 132   fAppendString = param.fAppendString;
 133   return (*this);
 134 }
 135
 136 //_____________________________________________________________________________
 137 AliBaseCalibViewer::~AliBaseCalibViewer()
 138 {
 139   //
 140   // AliBaseCalibViewer destructor
 141   // all objects will be deleted, the file will be closed, the pictures will disappear
 142   //
 143   if (fTree && fTreeMustBeDeleted) {
 144     fTree->SetCacheSize(0);
 145     fTree->Delete();
 146   }
 147   if (fFile) {
 148     fFile->Close();
 149     fFile = 0;
 150   }
 151
 152   for (Int_t i = fListOfObjectsToBeDeleted->GetEntriesFast()-1; i >= 0; i--) {
 153     delete fListOfObjectsToBeDeleted->At(i);
 154   }
 155   delete fListOfObjectsToBeDeleted;
 156 }
 157
 158 //_____________________________________________________________________________
 159 void AliBaseCalibViewer::Delete(Option_t* /*option*/) {
 160   //
 161   // Should be called from AliBaseCalibViewerGUI class only.
 162   // If you use Delete() do not call the destructor.
 163   // All objects (except those contained in fListOfObjectsToBeDeleted) will be deleted, the file will be closed.
 164   //
 165
 166   if (fTree && fTreeMustBeDeleted) {
 167     fTree->SetCacheSize(0);
 168     fTree->Delete();
 169   }
 170   delete fFile;
 171   delete fListOfObjectsToBeDeleted;
 172 }
 173
 174 //_____________________________________________________________________________
 175 void AliBaseCalibViewer::FormatHistoLabels(TH1 *histo) const {
 176   //
 177   // formats title and axis labels of histo
 178   // removes '.fElements'
 179   //
 180   if (!histo) return;
 181   TString replaceString(fAppendString.Data());
 182   TString *str = new TString(histo->GetTitle());
 183   str->ReplaceAll(replaceString, "");
 184   histo->SetTitle(str->Data());
 185   delete str;
 186   if (histo->GetXaxis()) {
 187     str = new TString(histo->GetXaxis()->GetTitle());
 188     str->ReplaceAll(replaceString, "");
 189     histo->GetXaxis()->SetTitle(str->Data());
 190     delete str;
 191   }
 192   if (histo->GetYaxis()) {
 193     str = new TString(histo->GetYaxis()->GetTitle());
 194     str->ReplaceAll(replaceString, "");
 195     histo->GetYaxis()->SetTitle(str->Data());
 196     delete str;
 197   }
 198   if (histo->GetZaxis()) {
 199     str = new TString(histo->GetZaxis()->GetTitle());
 200     str->ReplaceAll(replaceString, "");
 201     histo->GetZaxis()->SetTitle(str->Data());
 202     delete str;
 203   }
 204 }
 205
 206 //_____________________________________________________________________________
 207 TFriendElement* AliBaseCalibViewer::AddReferenceTree(const Char_t* filename, const Char_t* treename, const Char_t* refname){
 208   //
 209   // add a reference tree to the current tree
 210   // by default the treename is 'tree' and the reference treename is 'R'
 211   //
 212   TFile *file = new TFile(filename);
 213   fListOfObjectsToBeDeleted->Add(file);
 214   TTree * tree = (TTree*)file->Get(treename);
 215   return AddFriend(tree, refname);
 216 }
 217
 218 //_____________________________________________________________________________
 219 TString* AliBaseCalibViewer::Fit(const Char_t* drawCommand, const Char_t* formula, const Char_t* cuts,
 220     Double_t & chi2, TVectorD &fitParam, TMatrixD &covMatrix){
 221   //
 222   // fit an arbitrary function, specified by formula into the data, specified by drawCommand and cuts
 223   // returns chi2, fitParam and covMatrix
 224   // returns TString with fitted formula
 225   //
 226
 227   TString formulaStr(formula);
 228   TString drawStr(drawCommand);
 229   TString cutStr(cuts);
 230
 231   // abbreviations:
 232   drawStr.ReplaceAll(fAbbreviation, fAppendString);
 233   cutStr.ReplaceAll(fAbbreviation, fAppendString);
 234   formulaStr.ReplaceAll(fAbbreviation, fAppendString);
 235
 236   formulaStr.ReplaceAll("++", fAbbreviation);
 237   TObjArray* formulaTokens = formulaStr.Tokenize(fAbbreviation.Data());
 238   Int_t dim = formulaTokens->GetEntriesFast();
 239
 240   fitParam.ResizeTo(dim);
 241   covMatrix.ResizeTo(dim,dim);
 242
 243   TLinearFitter* fitter = new TLinearFitter(dim+1, Form("hyp%d",dim));
 244   fitter->StoreData(kTRUE);
 245   fitter->ClearPoints();
 246
 247   Int_t entries = Draw(drawStr.Data(), cutStr.Data(), "goff");
 248   if (entries == -1) {
 249     delete fitter;
 250     delete formulaTokens;
 251     return new TString("An ERROR has occured during fitting!");
 252   }
 253   Double_t **values = new Double_t*[dim+1] ;
 254
 255   for (Int_t i = 0; i < dim + 1; i++){
 256     Int_t centries = 0;
 257     if (i < dim) centries = fTree->Draw(((TObjString*)formulaTokens->At(i))->GetName(), cutStr.Data(), "goff");
 258     else  centries = fTree->Draw(drawStr.Data(), cutStr.Data(), "goff");
 259
 260     if (entries != centries) {
 261       delete fitter;
 262       delete [] values;
 263       delete formulaTokens;
 264       return new TString("An ERROR has occured during fitting!");
 265     }
 266     values[i] = new Double_t[entries];
 267     memcpy(values[i],  fTree->GetV1(), entries*sizeof(Double_t));
 268   }
 269
 270   // add points to the fitter
 271   for (Int_t i = 0; i < entries; i++){
 272     Double_t x[1000];
 273     for (Int_t j=0; j<dim;j++) x[j]=values[j][i];
 274     fitter->AddPoint(x, values[dim][i], 1);
 275   }
 276
 277   fitter->Eval();
 278   fitter->GetParameters(fitParam);
 279   fitter->GetCovarianceMatrix(covMatrix);
 280   chi2 = fitter->GetChisquare();
 281
 282   TString *preturnFormula = new TString(Form("( %e+",fitParam[0])), &returnFormula = *preturnFormula;
 283
 284   for (Int_t iparam = 0; iparam < dim; iparam++) {
 285     returnFormula.Append(Form("%s*(%e)",((TObjString*)formulaTokens->At(iparam))->GetName(),fitParam[iparam+1]));
 286     if (iparam < dim-1) returnFormula.Append("+");
 287   }
 288   returnFormula.Append(" )");
 289   delete formulaTokens;
 290   delete fitter;
 291   for (Int_t i = 0; i < dim + 1; i++) delete [] values[i];
 292   delete[] values;
 293   return preturnFormula;
 294 }
 295
 296 //_____________________________________________________________________________
 297 Double_t AliBaseCalibViewer::GetLTM(Int_t n, Double_t *array, Double_t *sigma, Double_t fraction){
 298   //
 299   //  returns the LTM and sigma
 300   //
 301   Double_t *ddata = new Double_t[n];
 302   Double_t mean = 0, lsigma = 0;
 303   UInt_t nPoints = 0;
 304   for (UInt_t i = 0; i < (UInt_t)n; i++) {
 305     ddata[nPoints]= array[nPoints];
 306     nPoints++;
 307   }
 308   Int_t hh = TMath::Min(TMath::Nint(fraction * nPoints), Int_t(n));
 309   AliMathBase::EvaluateUni(nPoints, ddata, mean, lsigma, hh);
 310   if (sigma) *sigma = lsigma;
 311   delete [] ddata;
 312   return mean;
 313 }
 314
 315 //_____________________________________________________________________________
 316 Int_t AliBaseCalibViewer::GetBin(Float_t value, Int_t nbins, Double_t binLow, Double_t binUp){
 317   // Returns the 'bin' for 'value'
 318   // The interval between 'binLow' and 'binUp' is divided into 'nbins' equidistant bins
 319   // avoid index out of bounds error: 'if (bin < binLow) bin = binLow' and vice versa
 320   /* Begin_Latex
 321      GetBin(value) = #frac{nbins - 1}{binUp - binLow} #upoint (value - binLow) +1
 322      End_Latex
 323      */
 324
 325   Int_t bin =  TMath::Nint( (Float_t)(value - binLow) / (Float_t)(binUp - binLow) * (nbins-1) ) + 1;
 326   // avoid index out of bounds:
 327   if (value < binLow) bin = 0;
 328   if (value > binUp)  bin = nbins + 1;
 329   return bin;
 330
 331 }
 332
 333 //_____________________________________________________________________________
 334 TH1F* AliBaseCalibViewer::SigmaCut(TH1F *histogram, Float_t mean, Float_t sigma, Float_t sigmaMax,
 335     Float_t sigmaStep, Bool_t pm) {
 336   //
 337   // Creates a cumulative histogram Begin_Latex S(t, #mu, #sigma) End_Latex, where you can see, how much of the data are inside sigma-intervals around the mean value
 338   // The data of the distribution Begin_Latex f(x, #mu, #sigma) End_Latex are given in 'histogram'
 339   // 'mean' and 'sigma' are Begin_Latex #mu End_Latex and  Begin_Latex #sigma End_Latex of the distribution in 'histogram', to be specified by the user
 340   // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma, Begin_Latex t #sigma End_Latex)
 341   // sigmaStep: the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
 342   // pm: Decide weather Begin_Latex t > 0 End_Latex (first case) or Begin_Latex t End_Latex arbitrary (secound case)
 343   // The actual work is done on the array.
 344   /* Begin_Latex
 345      f(x, #mu, #sigma)     #Rightarrow       S(t, #mu, #sigma) = (#int_{#mu}^{#mu + t #sigma} f(x, #mu, #sigma) dx + #int_{#mu}^{#mu - t #sigma} f(x, #mu, #sigma) dx) / (#int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx),    for  t > 0
 346      or
 347      f(x, #mu, #sigma)     #Rightarrow       S(t, #mu, #sigma) = #int_{#mu}^{#mu + t #sigma} f(x, #mu, #sigma) dx / #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx
 348      End_Latex
 349      Begin_Macro(source)
 350      {
 351      Float_t mean = 0;
 352      Float_t sigma = 1.5;
 353      Float_t sigmaMax = 4;
 354      gROOT->SetStyle("Plain");
 355      TH1F *distribution = new TH1F("Distribution1", "Distribution f(x, #mu, #sigma)", 1000,-5,5);
 356      TRandom rand(23);
 357      for (Int_t i = 0; i <50000;i++) distribution->Fill(rand.Gaus(mean, sigma));
 358      Float_t *ar = distribution->GetArray();
 359
 360      TCanvas* macro_example_canvas = new TCanvas("cAliBaseCalibViewer", "", 350, 350);
 361      macro_example_canvas->Divide(0,3);
 362      TVirtualPad *pad1 = macro_example_canvas->cd(1);
 363      pad1->SetGridy();
 364      pad1->SetGridx();
 365      distribution->Draw();
 366      TVirtualPad *pad2 = macro_example_canvas->cd(2);
 367      pad2->SetGridy();
 368      pad2->SetGridx();
 369
 370      TH1F *shist = AliTPCCalibViewer::SigmaCut(distribution, mean, sigma, sigmaMax);
 371      shist->SetNameTitle("Cumulative","Cumulative S(t, #mu, #sigma)");
 372      shist->Draw();
 373      TVirtualPad *pad3 = macro_example_canvas->cd(3);
 374      pad3->SetGridy();
 375      pad3->SetGridx();
 376      TH1F *shistPM = AliTPCCalibViewer::SigmaCut(distribution, mean, sigma, sigmaMax, -1, kTRUE);
 377      shistPM->Draw();
 378      return macro_example_canvas;
 379      }
 380      End_Macro
 381      */
 382
 383   Float_t *array = histogram->GetArray();
 384   Int_t    nbins = histogram->GetXaxis()->GetNbins();
 385   Float_t binLow = histogram->GetXaxis()->GetXmin();
 386   Float_t binUp  = histogram->GetXaxis()->GetXmax();
 387   return AliBaseCalibViewer::SigmaCut(nbins, array, mean, sigma, nbins, binLow, binUp, sigmaMax, sigmaStep, pm);
 388 }
 389
 390 //_____________________________________________________________________________
 391 TH1F* AliBaseCalibViewer::SigmaCut(Int_t n, Float_t *array, Float_t mean, Float_t sigma, Int_t nbins, Float_t binLow, Float_t binUp, Float_t sigmaMax, Float_t sigmaStep, Bool_t pm){
 392   //
 393   // Creates a histogram Begin_Latex S(t, #mu, #sigma) End_Latex, where you can see, how much of the data are inside sigma-intervals around the mean value
 394   // The data of the distribution Begin_Latex f(x, #mu, #sigma) End_Latex are given in 'array', 'n' specifies the length of the array
 395   // 'mean' and 'sigma' are Begin_Latex #mu End_Latex and  Begin_Latex #sigma End_Latex of the distribution in 'array', to be specified by the user
 396   // 'nbins': number of bins, 'binLow': first bin, 'binUp': last bin
 397   // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma, Begin_Latex t #sigma End_Latex)
 398   // sigmaStep: the binsize of the generated histogram
 399   // Here the actual work is done.
 400
 401   if (TMath::Abs(sigma) < 1.e-10) return 0;
 402   Float_t binWidth = (binUp-binLow)/(nbins - 1);
 403   if (sigmaStep <= 0) sigmaStep = binWidth;
 404   Int_t kbins = (Int_t)(sigmaMax * sigma / sigmaStep) + 1; // + 1  due to overflow bin in histograms
 405   if (pm) kbins = 2 * (Int_t)(sigmaMax * sigma / sigmaStep) + 1;
 406   Float_t kbinLow = !pm ? 0 : -sigmaMax;
 407   Float_t kbinUp  = sigmaMax;
 408   TH1F *hist = new TH1F("sigmaCutHisto","Cumulative; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
 409   hist->SetDirectory(0);
 410   hist->Reset();
 411
 412   // calculate normalization
 413   Double_t normalization = 0;
 414   for (Int_t i = 0; i <= n; i++) {
 415     normalization += array[i];
 416   }
 417
 418   // given units: units from given histogram
 419   // sigma units: in units of sigma
 420   // iDelta: integrate in interval (mean +- iDelta), given units
 421   // x:      ofset from mean for integration, given units
 422   // hist:   needs
 423
 424   // fill histogram
 425   for (Float_t iDelta = 0; iDelta <= sigmaMax * sigma; iDelta += sigmaStep) {
 426     // integrate array
 427     Double_t valueP = array[GetBin(mean, nbins, binLow, binUp)];
 428     Double_t valueM = array[GetBin(mean-binWidth, nbins, binLow, binUp)];
 429     // add bin of mean value only once to the histogram
 430     for (Float_t x = binWidth; x <= iDelta; x += binWidth) {
 431       valueP += (mean + x <= binUp)  ? array[GetBin(mean + x, nbins, binLow, binUp)] : 0;
 432       valueM += (mean-binWidth - x >= binLow) ? array[GetBin(mean-binWidth - x, nbins, binLow, binUp)] : 0;
 433     }
 434
 435     if (valueP / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail  +++ \n", valueP, normalization);
 436     if (valueP / normalization > 100) return hist;
 437     if (valueM / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail  +++ \n", valueM, normalization);
 438     if (valueM / normalization > 100) return hist;
 439     valueP = (valueP / normalization);
 440     valueM = (valueM / normalization);
 441     if (pm) {
 442       Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
 443       hist->SetBinContent(bin, valueP);
 444       bin = GetBin(-iDelta/sigma, kbins, kbinLow, kbinUp);
 445       hist->SetBinContent(bin, valueM);
 446     }
 447     else { // if (!pm)
 448       Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
 449       hist->SetBinContent(bin, valueP + valueM);
 450     }
 451   }
 452   if (!pm) hist->SetMaximum(1.2);
 453   return hist;
 454 }
 455
 456 //_____________________________________________________________________________
 457 TH1F* AliBaseCalibViewer::SigmaCut(Int_t /*n*/, Double_t */*array*/, Double_t /*mean*/, Double_t /*sigma*/,
 458     Int_t /*nbins*/, Double_t */*xbins*/, Double_t /*sigmaMax*/){
 459   //
 460   // SigmaCut for variable binsize
 461   // NOT YET IMPLEMENTED !!!
 462   //
 463   printf("SigmaCut with variable binsize, Not yet implemented\n");
 464
 465   return 0;
 466 }
 467
 468 //_____________________________________________________________________________
 469 Int_t  AliBaseCalibViewer::DrawHisto1D(const Char_t* drawCommand, const Char_t* sector, const Char_t* cuts,
 470     const Char_t *sigmas, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM) const
 471 {
 472   //
 473   // Easy drawing of data, in principle the same as EasyDraw1D
 474   // Difference: A line for the mean / median / LTM is drawn
 475   // in 'sigmas' you can specify in which distance to the mean/median/LTM you want to see a line in sigma-units, separated by ';'
 476   // example: sigmas = "2; 4; 6;"  at Begin_Latex 2 #sigma End_Latex, Begin_Latex 4 #sigma End_Latex and Begin_Latex 6 #sigma End_Latex  a line is drawn.
 477   // "plotMean", "plotMedian" and "plotLTM": what kind of lines do you want to see?
 478   //
 479   Int_t oldOptStat = gStyle->GetOptStat();
 480   gStyle->SetOptStat(0000000);
 481   Double_t ltmFraction = 0.8;
 482
 483   TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
 484   TVectorF nsigma(sigmasTokens->GetEntriesFast());
 485   for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
 486     TString str(((TObjString*)sigmasTokens->At(i))->GetString());
 487     Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
 488     nsigma[i] = sig;
 489   }
 490   delete sigmasTokens;
 491   //
 492   TString drawStr(drawCommand);
 493   Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
 494   if (dangerousToDraw) {
 495     Warning("DrawHisto1D", "The draw string must not contain ':' or '>>'.");
 496     return -1;
 497   }
 498   drawStr += " >> tempHist";
 499   Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts);
 500   TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
 501   // FIXME is this histogram deleted automatically?
 502   Double_t *values = fTree->GetV1();  // value is the array containing 'entries' numbers
 503
 504   Double_t mean = TMath::Mean(entries, values);
 505   Double_t median = TMath::Median(entries, values);
 506   Double_t sigma = TMath::RMS(entries, values);
 507   Double_t maxY = htemp->GetMaximum();
 508
 509   TLegend * legend = new TLegend(.7,.7, .99, .99, "Statistical information");
 510   //fListOfObjectsToBeDeleted->Add(legend);
 511
 512   if (plotMean) {
 513     // draw Mean
 514     TLine* line = new TLine(mean, 0, mean, maxY);
 515     //fListOfObjectsToBeDeleted->Add(line);
 516     line->SetLineColor(kRed);
 517     line->SetLineWidth(2);
 518     line->SetLineStyle(1);
 519     line->Draw();
 520     legend->AddEntry(line, Form("Mean: %f", mean), "l");
 521     // draw sigma lines
 522     for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
 523       TLine* linePlusSigma = new TLine(mean + nsigma[i] * sigma, 0, mean + nsigma[i] * sigma, maxY);
 524       //fListOfObjectsToBeDeleted->Add(linePlusSigma);
 525       linePlusSigma->SetLineColor(kRed);
 526       linePlusSigma->SetLineStyle(2 + i);
 527       linePlusSigma->Draw();
 528       TLine* lineMinusSigma = new TLine(mean - nsigma[i] * sigma, 0, mean - nsigma[i] * sigma, maxY);
 529       //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
 530       lineMinusSigma->SetLineColor(kRed);
 531       lineMinusSigma->SetLineStyle(2 + i);
 532       lineMinusSigma->Draw();
 533       legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma)), "l");
 534     }
 535   }
 536   if (plotMedian) {
 537     // draw median
 538     TLine* line = new TLine(median, 0, median, maxY);
 539     //fListOfObjectsToBeDeleted->Add(line);
 540     line->SetLineColor(kBlue);
 541     line->SetLineWidth(2);
 542     line->SetLineStyle(1);
 543     line->Draw();
 544     legend->AddEntry(line, Form("Median: %f", median), "l");
 545     // draw sigma lines
 546     for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
 547       TLine* linePlusSigma = new TLine(median + nsigma[i] * sigma, 0, median + nsigma[i]*sigma, maxY);
 548       //fListOfObjectsToBeDeleted->Add(linePlusSigma);
 549       linePlusSigma->SetLineColor(kBlue);
 550       linePlusSigma->SetLineStyle(2 + i);
 551       linePlusSigma->Draw();
 552       TLine* lineMinusSigma = new TLine(median - nsigma[i] * sigma, 0, median - nsigma[i]*sigma, maxY);
 553       //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
 554       lineMinusSigma->SetLineColor(kBlue);
 555       lineMinusSigma->SetLineStyle(2 + i);
 556       lineMinusSigma->Draw();
 557       legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f",(Int_t)(nsigma[i]), (Float_t)(nsigma[i] * sigma)), "l");
 558     }
 559   }
 560   if (plotLTM) {
 561     // draw LTM
 562     Double_t ltmRms = 0;
 563     Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
 564     TLine* line = new TLine(ltm, 0, ltm, maxY);
 565     //fListOfObjectsToBeDeleted->Add(line);
 566     line->SetLineColor(kGreen+2);
 567     line->SetLineWidth(2);
 568     line->SetLineStyle(1);
 569     line->Draw();
 570     legend->AddEntry(line, Form("LTM: %f", ltm), "l");
 571     // draw sigma lines
 572     for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
 573       TLine* linePlusSigma = new TLine(ltm + nsigma[i] * ltmRms, 0, ltm + nsigma[i] * ltmRms, maxY);
 574       //fListOfObjectsToBeDeleted->Add(linePlusSigma);
 575       linePlusSigma->SetLineColor(kGreen+2);
 576       linePlusSigma->SetLineStyle(2+i);
 577       linePlusSigma->Draw();
 578
 579       TLine* lineMinusSigma = new TLine(ltm - nsigma[i] * ltmRms, 0, ltm - nsigma[i] * ltmRms, maxY);
 580       //fListOfObjectsToBeDeleted->Add(lineMinusSigma);
 581       lineMinusSigma->SetLineColor(kGreen+2);
 582       lineMinusSigma->SetLineStyle(2+i);
 583       lineMinusSigma->Draw();
 584       legend->AddEntry(lineMinusSigma, Form("%i #sigma = %f", (Int_t)(nsigma[i]), (Float_t)(nsigma[i] * ltmRms)), "l");
 585     }
 586   }
 587   if (!plotMean && !plotMedian && !plotLTM) return -1;
 588   legend->Draw();
 589   gStyle->SetOptStat(oldOptStat);
 590   return 1;
 591 }
 592
 593 //_____________________________________________________________________________
 594 Int_t AliBaseCalibViewer::SigmaCut(const Char_t* drawCommand, const Char_t* sector, const Char_t* cuts,
 595     Float_t sigmaMax, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM, Bool_t pm,
 596     const Char_t *sigmas, Float_t sigmaStep) const {
 597   //
 598   // Creates a histogram, where you can see, how much of the data are inside sigma-intervals
 599   // around the mean/median/LTM
 600   // with drawCommand, sector and cuts you specify your input data, see EasyDraw
 601   // sigmaMax: up to which sigma around the mean/median/LTM the histogram is generated (in units of sigma)
 602   // sigmaStep: the binsize of the generated histogram
 603   // plotMean/plotMedian/plotLTM: specifies where to put the center
 604   //
 605
 606   Double_t ltmFraction = 0.8;
 607
 608   TString drawStr(drawCommand);
 609   Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
 610   if (dangerousToDraw) {
 611     Warning("SigmaCut", "The draw string must not contain ':' or '>>'.");
 612     return -1;
 613   }
 614   drawStr += " >> tempHist";
 615
 616   Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
 617   TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
 618   // FIXME is this histogram deleted automatically?
 619   Double_t *values = fTree->GetV1();  // value is the array containing 'entries' numbers
 620
 621   Double_t mean = TMath::Mean(entries, values);
 622   Double_t median = TMath::Median(entries, values);
 623   Double_t sigma = TMath::RMS(entries, values);
 624
 625   TLegend * legend = new TLegend(.7,.7, .99, .99, "Cumulative");
 626   //fListOfObjectsToBeDeleted->Add(legend);
 627   TH1F *cutHistoMean = 0;
 628   TH1F *cutHistoMedian = 0;
 629   TH1F *cutHistoLTM = 0;
 630
 631   TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
 632   TVectorF nsigma(sigmasTokens->GetEntriesFast());
 633   for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
 634     TString str(((TObjString*)sigmasTokens->At(i))->GetString());
 635     Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
 636     nsigma[i] = sig;
 637   }
 638   delete sigmasTokens;
 639   //
 640   if (plotMean) {
 641     cutHistoMean = SigmaCut(htemp, mean, sigma, sigmaMax, sigmaStep, pm);
 642     if (cutHistoMean) {
 643       //fListOfObjectsToBeDeleted->Add(cutHistoMean);
 644       cutHistoMean->SetLineColor(kRed);
 645       legend->AddEntry(cutHistoMean, "Mean", "l");
 646       cutHistoMean->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 647       cutHistoMean->Draw();
 648       DrawLines(cutHistoMean, nsigma, legend, kRed, pm);
 649     } // if (cutHistoMean)
 650
 651   }
 652   if (plotMedian) {
 653     cutHistoMedian = SigmaCut(htemp, median, sigma, sigmaMax, sigmaStep, pm);
 654     if (cutHistoMedian) {
 655       //fListOfObjectsToBeDeleted->Add(cutHistoMedian);
 656       cutHistoMedian->SetLineColor(kBlue);
 657       legend->AddEntry(cutHistoMedian, "Median", "l");
 658       cutHistoMedian->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 659       if (plotMean && cutHistoMean) cutHistoMedian->Draw("same");
 660       else cutHistoMedian->Draw();
 661       DrawLines(cutHistoMedian, nsigma, legend, kBlue, pm);
 662     }  // if (cutHistoMedian)
 663   }
 664   if (plotLTM) {
 665     Double_t ltmRms = 0;
 666     Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
 667     cutHistoLTM = SigmaCut(htemp, ltm, ltmRms, sigmaMax, sigmaStep, pm);
 668     if (cutHistoLTM) {
 669       //fListOfObjectsToBeDeleted->Add(cutHistoLTM);
 670       cutHistoLTM->SetLineColor(kGreen+2);
 671       legend->AddEntry(cutHistoLTM, "LTM", "l");
 672       cutHistoLTM->SetTitle(Form("%s, cumulative; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 673       if ((plotMean && cutHistoMean) || (plotMedian && cutHistoMedian)) cutHistoLTM->Draw("same");
 674       else cutHistoLTM->Draw();
 675       DrawLines(cutHistoLTM, nsigma, legend, kGreen+2, pm);
 676     }
 677   }
 678   if (!plotMean && !plotMedian && !plotLTM) return -1;
 679   legend->Draw();
 680   return 1;
 681 }
 682
 683 //_____________________________________________________________________________
 684 Int_t AliBaseCalibViewer::Integrate(const Char_t* drawCommand, const Char_t* sector, const Char_t* cuts,
 685     Float_t /*sigmaMax*/, Bool_t plotMean, Bool_t plotMedian, Bool_t plotLTM,
 686     const Char_t *sigmas, Float_t /*sigmaStep*/) const {
 687   //
 688   // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
 689   // "mean" and "sigma" are Begin_Latex #mu End_Latex and  Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
 690   // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
 691   // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
 692   // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
 693   // The actual work is done on the array.
 694   /* Begin_Latex
 695      f(x, #mu, #sigma)     #Rightarrow       S(t, #mu, #sigma) = #int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx / #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx
 696      End_Latex
 697      */
 698
 699   Double_t ltmFraction = 0.8;
 700
 701   TString drawStr(drawCommand);
 702   Bool_t dangerousToDraw = drawStr.Contains(":") || drawStr.Contains(">>");
 703   if (dangerousToDraw) {
 704     Warning("Integrate", "The draw string must not contain ':' or '>>'.");
 705     return -1;
 706   }
 707   drawStr += " >> tempHist";
 708
 709   Int_t entries = EasyDraw1D(drawStr.Data(), sector, cuts, "goff");
 710   TH1F *htemp = (TH1F*)gDirectory->Get("tempHist");
 711   TGraph *integralGraphMean   = 0;
 712   TGraph *integralGraphMedian = 0;
 713   TGraph *integralGraphLTM    = 0;
 714   Double_t *values = fTree->GetV1();  // value is the array containing 'entries' numbers
 715   Int_t    *index  = new Int_t[entries];
 716   Float_t  *xarray = new Float_t[entries];
 717   Float_t  *yarray = new Float_t[entries];
 718   TMath::Sort(entries, values, index, kFALSE);
 719
 720   Double_t mean = TMath::Mean(entries, values);
 721   Double_t median = TMath::Median(entries, values);
 722   Double_t sigma = TMath::RMS(entries, values);
 723
 724   // parse sigmas string
 725   TObjArray *sigmasTokens = TString(sigmas).Tokenize(";");
 726   TVectorF nsigma(sigmasTokens->GetEntriesFast());
 727   for (Int_t i = 0; i < sigmasTokens->GetEntriesFast(); i++) {
 728     TString str(((TObjString*)sigmasTokens->At(i))->GetString());
 729     Double_t sig = (str.IsFloat()) ? str.Atof() : 0;
 730     nsigma[i] = sig;
 731   }
 732   delete sigmasTokens;
 733   TLegend * legend = new TLegend(.7,.7, .99, .99, "Integrated histogram");
 734   //fListOfObjectsToBeDeleted->Add(legend);
 735
 736   if (plotMean) {
 737     for (Int_t i = 0; i < entries; i++) {
 738       xarray[i] = (values[index[i]] - mean) / sigma;
 739       yarray[i] = float(i) / float(entries);
 740     }
 741     integralGraphMean = new TGraph(entries, xarray, yarray);
 742     if (integralGraphMean) {
 743       //fListOfObjectsToBeDeleted->Add(integralGraphMean);
 744       integralGraphMean->SetLineColor(kRed);
 745       legend->AddEntry(integralGraphMean, "Mean", "l");
 746       integralGraphMean->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 747       integralGraphMean->Draw("alu");
 748       DrawLines(integralGraphMean, nsigma, legend, kRed, kTRUE);
 749     }
 750   }
 751   if (plotMedian) {
 752     for (Int_t i = 0; i < entries; i++) {
 753       xarray[i] = (values[index[i]] - median) / sigma;
 754       yarray[i] = float(i) / float(entries);
 755     }
 756     integralGraphMedian = new TGraph(entries, xarray, yarray);
 757     if (integralGraphMedian) {
 758       //fListOfObjectsToBeDeleted->Add(integralGraphMedian);
 759       integralGraphMedian->SetLineColor(kBlue);
 760       legend->AddEntry(integralGraphMedian, "Median", "l");
 761       integralGraphMedian->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 762       if (plotMean && integralGraphMean) integralGraphMedian->Draw("samelu");
 763       else integralGraphMedian->Draw("alu");
 764       DrawLines(integralGraphMedian, nsigma, legend, kBlue, kTRUE);
 765     }
 766   }
 767   if (plotLTM) {
 768     Double_t ltmRms = 0;
 769     Double_t ltm = GetLTM(entries, values, &ltmRms, ltmFraction);
 770     for (Int_t i = 0; i < entries; i++) {
 771       xarray[i] = (values[index[i]] - ltm) / ltmRms;
 772       yarray[i] = float(i) / float(entries);
 773     }
 774     integralGraphLTM = new TGraph(entries, xarray, yarray);
 775     if (integralGraphLTM) {
 776       //fListOfObjectsToBeDeleted->Add(integralGraphLTM);
 777       integralGraphLTM->SetLineColor(kGreen+2);
 778       legend->AddEntry(integralGraphLTM, "LTM", "l");
 779       integralGraphLTM->SetTitle(Form("%s, integrated; Multiples of #sigma; Fraction of included data", htemp->GetTitle()));
 780       if ((plotMean && integralGraphMean) || (plotMedian && integralGraphMedian)) integralGraphLTM->Draw("samelu");
 781       else integralGraphLTM->Draw("alu");
 782       DrawLines(integralGraphLTM, nsigma, legend, kGreen+2, kTRUE);
 783     }
 784   }
 785   delete [] index;
 786   delete [] xarray;
 787   delete [] yarray;
 788   if (!plotMean && !plotMedian && !plotLTM) return -1;
 789   legend->Draw();
 790   return entries;
 791 }
 792
 793 //_____________________________________________________________________________
 794 TH1F* AliBaseCalibViewer::Integrate(TH1F *histogram, Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
 795   //
 796   // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
 797   // "mean" and "sigma" are Begin_Latex #mu End_Latex and  Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
 798   // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
 799   // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
 800   // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
 801   // The actual work is done on the array.
 802   /* Begin_Latex
 803      f(x, #mu, #sigma)     #Rightarrow       S(t, #mu, #sigma) = #int_{-#infty}^{#mu + t #sigma} f(x, #mu, #sigma) dx / #int_{-#infty}^{+#infty} f(x, #mu, #sigma) dx
 804      End_Latex
 805      Begin_Macro(source)
 806      {
 807      Float_t mean = 0;
 808      Float_t sigma = 1.5;
 809      Float_t sigmaMax = 4;
 810      gROOT->SetStyle("Plain");
 811      TH1F *distribution = new TH1F("Distribution2", "Distribution f(x, #mu, #sigma)", 1000,-5,5);
 812      TRandom rand(23);
 813      for (Int_t i = 0; i <50000;i++) distribution->Fill(rand.Gaus(mean, sigma));
 814      Float_t *ar = distribution->GetArray();
 815
 816      TCanvas* macro_example_canvas = new TCanvas("macro_example_canvas_Integrate", "", 350, 350);
 817      macro_example_canvas->Divide(0,2);
 818      TVirtualPad *pad1 = macro_example_canvas->cd(1);
 819      pad1->SetGridy();
 820      pad1->SetGridx();
 821      distribution->Draw();
 822      TVirtualPad *pad2 = macro_example_canvas->cd(2);
 823      pad2->SetGridy();
 824      pad2->SetGridx();
 825      TH1F *shist = AliTPCCalibViewer::Integrate(distribution, mean, sigma, sigmaMax);
 826      shist->SetNameTitle("Cumulative","Cumulative S(t, #mu, #sigma)");
 827      shist->Draw();
 828
 829      return macro_example_canvas_Integrate;
 830      }
 831      End_Macro
 832      */
 833
 834
 835   Float_t *array = histogram->GetArray();
 836   Int_t    nbins = histogram->GetXaxis()->GetNbins();
 837   Float_t binLow = histogram->GetXaxis()->GetXmin();
 838   Float_t binUp  = histogram->GetXaxis()->GetXmax();
 839   return Integrate(nbins, array, nbins, binLow, binUp, mean, sigma, sigmaMax, sigmaStep);
 840 }
 841
 842 //_____________________________________________________________________________
 843 TH1F* AliBaseCalibViewer::Integrate(Int_t n, Float_t *array, Int_t nbins, Float_t binLow, Float_t binUp,
 844     Float_t mean, Float_t sigma, Float_t sigmaMax, Float_t sigmaStep){
 845   // Creates an integrated histogram Begin_Latex S(t, #mu, #sigma) End_Latex, out of the input distribution distribution Begin_Latex f(x, #mu, #sigma) End_Latex, given in "histogram"
 846   // "mean" and "sigma" are Begin_Latex #mu End_Latex and  Begin_Latex #sigma End_Latex of the distribution in "histogram", to be specified by the user
 847   // sigmaMax: up to which sigma around the mean/median/LTM you want to integrate
 848   // if "igma == 0" and "sigmaMax == 0" the whole histogram is integrated
 849   // "sigmaStep": the binsize of the generated histogram, -1 means, that the maximal reasonable stepsize is used
 850   // Here the actual work is done.
 851
 852   Bool_t givenUnits = kTRUE;
 853   if (TMath::Abs(sigma) < 1.e-10 && TMath::Abs(sigmaMax) < 1.e-10) givenUnits = kFALSE;
 854   if (givenUnits) {
 855     sigma = 1;
 856     sigmaMax = (binUp - binLow) / 2.;
 857   }
 858
 859   Float_t binWidth = (binUp-binLow)/(nbins - 1);
 860   if (sigmaStep <= 0) sigmaStep = binWidth;
 861   Int_t kbins =  (Int_t)(sigmaMax * sigma / sigmaStep) + 1;  // + 1  due to overflow bin in histograms
 862   Float_t kbinLow = givenUnits ? binLow : -sigmaMax;
 863   Float_t kbinUp  = givenUnits ? binUp  : sigmaMax;
 864   TH1F *hist = 0;
 865   if (givenUnits)  hist = new TH1F("integratedHisto","Integrated Histogram; Given x; Fraction of included data", kbins, kbinLow, kbinUp);
 866   if (!givenUnits) hist = new TH1F("integratedHisto","Integrated Histogram; Multiples of #sigma; Fraction of included data", kbins, kbinLow, kbinUp);
 867   hist->SetDirectory(0);
 868   hist->Reset();
 869
 870   // calculate normalization
 871   //  printf("calculating normalization, integrating from bin 1 to %i \n", n);
 872   Double_t normalization = 0;
 873   for (Int_t i = 1; i <= n; i++) {
 874     normalization += array[i];
 875   }
 876   //  printf("normalization: %f \n", normalization);
 877
 878   // given units: units from given histogram
 879   // sigma units: in units of sigma
 880   // iDelta: integrate in interval (mean +- iDelta), given units
 881   // x:      ofset from mean for integration, given units
 882   // hist:   needs
 883
 884   // fill histogram
 885   for (Float_t iDelta = mean - sigmaMax * sigma; iDelta <= mean + sigmaMax * sigma; iDelta += sigmaStep) {
 886     // integrate array
 887     Double_t value = 0;
 888     for (Float_t x = mean - sigmaMax * sigma; x <= iDelta; x += binWidth) {
 889       value += (x <= binUp && x >= binLow)  ? array[GetBin(x, nbins, binLow, binUp)] : 0;
 890     }
 891     if (value / normalization > 100) printf("+++ Error, value to big: %f, normalization with %f will fail  +++ \n", value, normalization);
 892     if (value / normalization > 100) return hist;
 893     Int_t bin = GetBin(iDelta/sigma, kbins, kbinLow, kbinUp);
 894     //  printf("first integration bin: %i, last integration bin: %i \n", GetBin(mean - sigmaMax * sigma, nbins, binLow, binUp), GetBin(iDelta, nbins, binLow, binUp));
 895     //  printf("value: %f, normalization: %f, normalized value: %f, iDelta: %f, Bin: %i \n", value, normalization, value/normalization, iDelta, bin);
 896     value = (value / normalization);
 897     hist->SetBinContent(bin, value);
 898   }
 899   return hist;
 900 }
 901
 902 //_____________________________________________________________________________
 903 void AliBaseCalibViewer::DrawLines(TH1F *histogram, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
 904   //
 905   // Private function for SigmaCut(...) and Integrate(...)
 906   // Draws lines into the given histogram, specified by "nsigma", the lines are addeed to the legend
 907   //
 908
 909   // start to draw the lines, loop over requested sigmas
 910   for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
 911     if (!pm) {
 912       Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
 913       TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
 914       //fListOfObjectsToBeDeleted->Add(lineUp);
 915       lineUp->SetLineColor(color);
 916       lineUp->SetLineStyle(2 + i);
 917       lineUp->Draw();
 918       TLine* lineLeft = new TLine(nsigma[i], histogram->GetBinContent(bin), 0, histogram->GetBinContent(bin));
 919       //fListOfObjectsToBeDeleted->Add(lineLeft);
 920       lineLeft->SetLineColor(color);
 921       lineLeft->SetLineStyle(2 + i);
 922       lineLeft->Draw();
 923       legend->AddEntry(lineLeft, Form("Fraction(%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
 924     }
 925     else { // if (pm)
 926       Int_t bin = histogram->GetXaxis()->FindBin(nsigma[i]);
 927       TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], histogram->GetBinContent(bin));
 928       //fListOfObjectsToBeDeleted->Add(lineUp1);
 929       lineUp1->SetLineColor(color);
 930       lineUp1->SetLineStyle(2 + i);
 931       lineUp1->Draw();
 932       TLine* lineLeft1 = new TLine(nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
 933       //fListOfObjectsToBeDeleted->Add(lineLeft1);
 934       lineLeft1->SetLineColor(color);
 935       lineLeft1->SetLineStyle(2 + i);
 936       lineLeft1->Draw();
 937       legend->AddEntry(lineLeft1, Form("Fraction(+%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
 938       bin = histogram->GetXaxis()->FindBin(-nsigma[i]);
 939       TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], histogram->GetBinContent(bin));
 940       //fListOfObjectsToBeDeleted->Add(lineUp2);
 941       lineUp2->SetLineColor(color);
 942       lineUp2->SetLineStyle(2 + i);
 943       lineUp2->Draw();
 944       TLine* lineLeft2 = new TLine(-nsigma[i], histogram->GetBinContent(bin), histogram->GetBinLowEdge(0)+histogram->GetBinWidth(0), histogram->GetBinContent(bin));
 945       //fListOfObjectsToBeDeleted->Add(lineLeft2);
 946       lineLeft2->SetLineColor(color);
 947       lineLeft2->SetLineStyle(2 + i);
 948       lineLeft2->Draw();
 949       legend->AddEntry(lineLeft2, Form("Fraction(-%f #sigma) = %f",nsigma[i], histogram->GetBinContent(bin)), "l");
 950     }
 951   }  // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
 952 }
 953
 954 //_____________________________________________________________________________
 955 void AliBaseCalibViewer::DrawLines(TGraph *graph, TVectorF nsigma, TLegend *legend, Int_t color, Bool_t pm) const {
 956   //
 957   // Private function for SigmaCut(...) and Integrate(...)
 958   // Draws lines into the given histogram, specified by "nsigma", the lines are addeed to the legend
 959   //
 960
 961   // start to draw the lines, loop over requested sigmas
 962   for (Int_t i = 0; i < nsigma.GetNoElements(); i++) {
 963     if (!pm) {
 964       TLine* lineUp = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
 965       //fListOfObjectsToBeDeleted->Add(lineUp);
 966       lineUp->SetLineColor(color);
 967       lineUp->SetLineStyle(2 + i);
 968       lineUp->Draw();
 969       TLine* lineLeft = new TLine(nsigma[i], graph->Eval(nsigma[i]), 0, graph->Eval(nsigma[i]));
 970       //fListOfObjectsToBeDeleted->Add(lineLeft);
 971       lineLeft->SetLineColor(color);
 972       lineLeft->SetLineStyle(2 + i);
 973       lineLeft->Draw();
 974       legend->AddEntry(lineLeft, Form("Fraction(%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i])), "l");
 975     }
 976     else { // if (pm)
 977       TLine* lineUp1 = new TLine(nsigma[i], 0, nsigma[i], graph->Eval(nsigma[i]));
 978       //fListOfObjectsToBeDeleted->Add(lineUp1);
 979       lineUp1->SetLineColor(color);
 980       lineUp1->SetLineStyle(2 + i);
 981       lineUp1->Draw();
 982       TLine* lineLeft1 = new TLine(nsigma[i], graph->Eval(nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(nsigma[i]));
 983       //fListOfObjectsToBeDeleted->Add(lineLeft1);
 984       lineLeft1->SetLineColor(color);
 985       lineLeft1->SetLineStyle(2 + i);
 986       lineLeft1->Draw();
 987       legend->AddEntry(lineLeft1, Form("Fraction(+%f #sigma) = %f",nsigma[i], graph->Eval(nsigma[i])), "l");
 988       TLine* lineUp2 = new TLine(-nsigma[i], 0, -nsigma[i], graph->Eval(-nsigma[i]));
 989       //fListOfObjectsToBeDeleted->Add(lineUp2);
 990       lineUp2->SetLineColor(color);
 991       lineUp2->SetLineStyle(2 + i);
 992       lineUp2->Draw();
 993       TLine* lineLeft2 = new TLine(-nsigma[i], graph->Eval(-nsigma[i]), graph->GetHistogram()->GetXaxis()->GetBinLowEdge(0), graph->Eval(-nsigma[i]));
 994       //fListOfObjectsToBeDeleted->Add(lineLeft2);
 995       lineLeft2->SetLineColor(color);
 996       lineLeft2->SetLineStyle(2 + i);
 997       lineLeft2->Draw();
 998       legend->AddEntry(lineLeft2, Form("Fraction(-%f #sigma) = %f",nsigma[i], graph->Eval(-nsigma[i])), "l");
 999     }
1000   }  // for (Int_t i = 0; i < nsigma.GetNoElements(); i++)
1001 }