ROOT logo
/* $Id$ */

// ------------------------------------------------------
//
// Class to handle 3d-corrections.
//
// ------------------------------------------------------
//

#include <TDirectory.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TH3F.h>
#include <TString.h>
#include <TMath.h>

#include <AliLog.h>

#include "AliCorrectionMatrix3D.h"
#include "AliCorrectionMatrix2D.h"
#include "AliPWG0Helper.h"

//____________________________________________________________________
ClassImp(AliCorrectionMatrix3D)

//____________________________________________________________________
AliCorrectionMatrix3D::AliCorrectionMatrix3D() :
  AliCorrectionMatrix()
{
  // default constructor
}

//____________________________________________________________________
AliCorrectionMatrix3D::AliCorrectionMatrix3D(const AliCorrectionMatrix3D& c)
  : AliCorrectionMatrix(c)
{
  // copy constructor
  ((AliCorrectionMatrix3D &)c).Copy(*this);
}

//____________________________________________________________________
AliCorrectionMatrix3D &AliCorrectionMatrix3D::operator=(const AliCorrectionMatrix3D &c)
{
  // assigment operator

  if (this != &c)
    ((AliCorrectionMatrix3D &) c).Copy(*this);

  return *this;
}

//____________________________________________________________________
AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
              Int_t nBinX, Float_t Xmin, Float_t Xmax,
              Int_t nBinY, Float_t Ymin, Float_t Ymax,
              Int_t nBinZ, Float_t Zmin, Float_t Zmax)
  : AliCorrectionMatrix(name, title)
{
  //
  // constructor
  //

  Float_t* binLimitsX = new Float_t[nBinX+1];
  for (Int_t i=0; i<=nBinX; ++i)
    binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;

  Float_t* binLimitsY = new Float_t[nBinY+1];
  for (Int_t i=0; i<=nBinY; ++i)
    binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;

  Float_t* binLimitsZ = new Float_t[nBinZ+1];
  for (Int_t i=0; i<=nBinZ; ++i)
    binLimitsZ[i] = Zmin + (Zmax - Zmin) / nBinZ * i;

  CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);

  delete[] binLimitsX;
  delete[] binLimitsY;
  delete[] binLimitsZ;
}

AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
      Int_t nBinX, Float_t Xmin, Float_t Xmax,
      Int_t nBinY, Float_t Ymin, Float_t Ymax,
      Int_t nBinZ, const Float_t* zbins)
  : AliCorrectionMatrix(name, title)
{
  // constructor with variable bin sizes

  Float_t* binLimitsX = new Float_t[nBinX+1];
  for (Int_t i=0; i<=nBinX; ++i)
    binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;

  Float_t* binLimitsY = new Float_t[nBinY+1];
  for (Int_t i=0; i<=nBinY; ++i)
    binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;

  CreateHists(nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);

  delete[] binLimitsX;
  delete[] binLimitsY;
}

AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title, TH3* hBinning)
  : AliCorrectionMatrix(name, title)
{
  // constructor with variable bin sizes (uses binning of hBinning)

  // do not add this hists to the directory
  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);

  fhMeas = (TH3F*)hBinning->Clone("measured");
  fhGene = (TH3F*)hBinning->Clone("generated");
  fhCorr = (TH3F*)hBinning->Clone("correction");

  fhMeas->SetTitle(Form("%s measured", GetTitle()));
  fhGene->SetTitle(Form("%s generated", GetTitle()));
  fhCorr->SetTitle(Form("%s correction", GetTitle()));

  fhMeas->Reset();
  fhGene->Reset();
  fhCorr->Reset();

  TH1::AddDirectory(oldStatus);

  fhMeas->Sumw2();
  fhGene->Sumw2();
  fhCorr->Sumw2();
}

//____________________________________________________________________
void AliCorrectionMatrix3D::CreateHists(Int_t nBinX, const Float_t* binLimitsX,
      Int_t nBinY, const Float_t* binLimitsY,
      Int_t nBinZ, const Float_t* binLimitsZ)
{
  // create the histograms

  // do not add this hists to the directory
  Bool_t oldStatus = TH1::AddDirectoryStatus();
  TH1::AddDirectory(kFALSE);

  fhMeas  = new TH3F("measured",   Form("%s measured",GetTitle()),   nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
  fhGene  = new TH3F("generated",  Form("%s generated",GetTitle()),  nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);
  fhCorr  = new TH3F("correction", Form("%s correction",GetTitle()), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, binLimitsZ);

  fhMeas->Sumw2();
  fhGene->Sumw2();
  fhCorr->Sumw2();

  TH1::AddDirectory(oldStatus);
}


//____________________________________________________________________
AliCorrectionMatrix3D::~AliCorrectionMatrix3D()
{
  //
  // destructor
  //

  // histograms already deleted in base class
}

//____________________________________________________________________
TH3* AliCorrectionMatrix3D::GetGeneratedHistogram()
{
  // return generated histogram casted to correct type
  return dynamic_cast<TH3F*> (fhGene);
}

//____________________________________________________________________
TH3* AliCorrectionMatrix3D::GetMeasuredHistogram()
{
  // return measured histogram casted to correct type
  return dynamic_cast<TH3F*> (fhMeas);
}

//____________________________________________________________________
TH3* AliCorrectionMatrix3D::GetCorrectionHistogram()
{
  // return correction histogram casted to correct type
  return dynamic_cast<TH3F*> (fhCorr);
}

//____________________________________________________________________
AliCorrectionMatrix2D* AliCorrectionMatrix3D::Get2DCorrection(Option_t* opt, Float_t aMin, Float_t aMax)
{
  // returns a 2D projection of this correction

  TString option = opt;

  // unzoom
  fhMeas->GetXaxis()->SetRange(0, 0);
  fhMeas->GetYaxis()->SetRange(0, 0);
  fhMeas->GetZaxis()->SetRange(0, 0);

  fhGene->GetXaxis()->SetRange(0, 0);
  fhGene->GetYaxis()->SetRange(0, 0);
  fhGene->GetZaxis()->SetRange(0, 0);

  if (aMin<aMax) {
    if (option.Contains("xy") || option.Contains("yx")) {
      Int_t bMin = fhMeas->GetZaxis()->FindBin(aMin);
      Int_t bMax = fhMeas->GetZaxis()->FindBin(aMax);
      fhGene->GetZaxis()->SetRange(bMin, bMax);
      fhMeas->GetZaxis()->SetRange(bMin, bMax);
    }
    else if (option.Contains("xz") || option.Contains("zx")) {
      Int_t bMin = fhMeas->GetYaxis()->FindBin(aMin);
      Int_t bMax = fhMeas->GetYaxis()->FindBin(aMax);
      fhGene->GetYaxis()->SetRange(bMin, bMax);
      fhMeas->GetYaxis()->SetRange(bMin, bMax);
    }
    else if (option.Contains("yz") || option.Contains("zy")) {
      Int_t bMin = fhMeas->GetXaxis()->FindBin(aMin);
      Int_t bMax = fhMeas->GetXaxis()->FindBin(aMax);
      fhGene->GetXaxis()->SetRange(bMin, bMax);
      fhMeas->GetXaxis()->SetRange(bMin, bMax);
    }
    else {
      AliDebug(AliLog::kWarning, Form("WARNING: unknown projection option %s", opt));
      return 0;
    }
  }

  AliCorrectionMatrix2D* corr2D = new AliCorrectionMatrix2D(Form("%s_%s",GetName(),opt),Form("%s projection %s",GetName(),opt),100,0,100,100,0,100);

  TH2F* meas = (TH2F*) ((TH3F*)fhMeas)->Project3D(option)->Clone(Form("%s_meas", corr2D->GetName()));
  TH2F* gene = (TH2F*) ((TH3F*)fhGene)->Project3D(option)->Clone(Form("%s_gene", corr2D->GetName()));
  
  // set errors
  for (Int_t x = 0; x<=meas->GetNbinsX()+1; x++)
    for (Int_t y = 0; y<=meas->GetNbinsY()+1; y++)
    {
      gene->SetBinError(x, y, TMath::Sqrt(gene->GetBinContent(x, y)));
      meas->SetBinError(x, y, TMath::Sqrt(meas->GetBinContent(x, y)));
    }

  TH2F* corr = (TH2F*)gene->Clone(Form("%s_corr", corr2D->GetName()));
  corr->Reset();

  corr2D->SetGeneratedHistogram(gene);
  corr2D->SetMeasuredHistogram(meas);
  corr2D->SetCorrectionHistogram(corr);

  corr2D->Divide();

  // unzoom
  fhMeas->GetXaxis()->SetRange(0, 0);
  fhMeas->GetYaxis()->SetRange(0, 0);
  fhMeas->GetZaxis()->SetRange(0, 0);

  fhGene->GetXaxis()->SetRange(0, 0);
  fhGene->GetYaxis()->SetRange(0, 0);
  fhGene->GetZaxis()->SetRange(0, 0);

  return corr2D;
}

//____________________________________________________________________
TH1* AliCorrectionMatrix3D::Get1DCorrectionHistogram(Option_t* opt, Float_t aMin1, Float_t aMax1, Float_t aMin2, Float_t aMax2)
{
  // returns a 1D projection of this correction
  
  AliCorrectionMatrix2D* corr2D = 0;
  if (strcmp(opt,"x")==0) {
    corr2D = Get2DCorrection("yx",aMin2,aMax2);
    return corr2D->Get1DCorrectionHistogram("x",aMin1,aMax1);
  }
  if (strcmp(opt,"y")==0) {
    corr2D = Get2DCorrection("xy",aMin2,aMax2);
    return corr2D->Get1DCorrectionHistogram("x",aMin1,aMax1);
  }  
  if (strcmp(opt,"z")==0) {
    corr2D = Get2DCorrection("yz",aMin1,aMax1);
    return corr2D->Get1DCorrectionHistogram("x",aMin2,aMax2);
  }  
  AliDebug(AliLog::kWarning, Form("WARNING: unknown projection option %s (should be x,y or z)", opt));
  
  return 0;
}

//____________________________________________________________________
void AliCorrectionMatrix3D::FillMeas(Float_t ax, Float_t ay, Float_t az, Double_t weight)
{
  // add value to measured histogram
  GetMeasuredHistogram()->Fill(ax, ay, az, weight);
}

//____________________________________________________________________
void AliCorrectionMatrix3D::FillGene(Float_t ax, Float_t ay, Float_t az, Double_t weight)
{
  // add value to generated histogram
  GetGeneratedHistogram()->Fill(ax, ay, az, weight);
}

//____________________________________________________________________
Float_t AliCorrectionMatrix3D::GetCorrection(Float_t ax, Float_t ay, Float_t az) const
{
  // returns a value of the correction map
  return fhCorr->GetBinContent(fhCorr->FindBin(ax, ay, az));
}

//____________________________________________________________________
//void AliCorrectionMatrix3D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge, Int_t nBinsZedge)
void AliCorrectionMatrix3D::RemoveEdges(Float_t, Int_t, Int_t, Int_t)
{
  // so what do we do here...
}

//____________________________________________________________________
void AliCorrectionMatrix3D::SaveHistograms()
{
  //
  // saves the histograms
  //

  AliCorrectionMatrix::SaveHistograms();

  if (GetGeneratedHistogram() && GetMeasuredHistogram())
  {
    gDirectory->cd(GetName());
    
    AliPWG0Helper::CreateDividedProjections(GetGeneratedHistogram(), GetMeasuredHistogram(), 0, kFALSE, kTRUE);
    
    gDirectory->cd("..");
  }
}

//____________________________________________________________________
Int_t AliCorrectionMatrix3D::CheckEmptyBins(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax, Bool_t quiet)
{
  //
  // counts the number of empty Bins in a given region
  //

  TH3* hist = GetGeneratedHistogram();
  if (!hist)
    return -1;

  Int_t emptyBins = 0;
  for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
    for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
      for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
        if (hist->GetBinContent(x, y, z) == 0)
        {
          if (!quiet)
            printf("Empty bin in %s at vtx = %f, eta = %f, pt = %f\n", GetName(), hist->GetXaxis()->GetBinCenter(x), hist->GetYaxis()->GetBinCenter(y), hist->GetZaxis()->GetBinCenter(z));
          ++emptyBins;
        }

  return emptyBins;
}

//____________________________________________________________________
TH1* AliCorrectionMatrix3D::PlotBinErrors(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax)
{
  //
  // makes a 1d plots of the relative bin errors
  //

  TH3* hist = GetCorrectionHistogram();
  if (!hist)
    return 0;

  TH1F* target = new TH1F("relerrors", "relerrors", 100, 0, 10);

  for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
    for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
      for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
        if (hist->GetBinContent(x, y, z) != 0)
          target->Fill(100 * hist->GetBinError(x, y, z) / hist->GetBinContent(x, y, z));

  return target;
}

 AliCorrectionMatrix3D.cxx:1
 AliCorrectionMatrix3D.cxx:2
 AliCorrectionMatrix3D.cxx:3
 AliCorrectionMatrix3D.cxx:4
 AliCorrectionMatrix3D.cxx:5
 AliCorrectionMatrix3D.cxx:6
 AliCorrectionMatrix3D.cxx:7
 AliCorrectionMatrix3D.cxx:8
 AliCorrectionMatrix3D.cxx:9
 AliCorrectionMatrix3D.cxx:10
 AliCorrectionMatrix3D.cxx:11
 AliCorrectionMatrix3D.cxx:12
 AliCorrectionMatrix3D.cxx:13
 AliCorrectionMatrix3D.cxx:14
 AliCorrectionMatrix3D.cxx:15
 AliCorrectionMatrix3D.cxx:16
 AliCorrectionMatrix3D.cxx:17
 AliCorrectionMatrix3D.cxx:18
 AliCorrectionMatrix3D.cxx:19
 AliCorrectionMatrix3D.cxx:20
 AliCorrectionMatrix3D.cxx:21
 AliCorrectionMatrix3D.cxx:22
 AliCorrectionMatrix3D.cxx:23
 AliCorrectionMatrix3D.cxx:24
 AliCorrectionMatrix3D.cxx:25
 AliCorrectionMatrix3D.cxx:26
 AliCorrectionMatrix3D.cxx:27
 AliCorrectionMatrix3D.cxx:28
 AliCorrectionMatrix3D.cxx:29
 AliCorrectionMatrix3D.cxx:30
 AliCorrectionMatrix3D.cxx:31
 AliCorrectionMatrix3D.cxx:32
 AliCorrectionMatrix3D.cxx:33
 AliCorrectionMatrix3D.cxx:34
 AliCorrectionMatrix3D.cxx:35
 AliCorrectionMatrix3D.cxx:36
 AliCorrectionMatrix3D.cxx:37
 AliCorrectionMatrix3D.cxx:38
 AliCorrectionMatrix3D.cxx:39
 AliCorrectionMatrix3D.cxx:40
 AliCorrectionMatrix3D.cxx:41
 AliCorrectionMatrix3D.cxx:42
 AliCorrectionMatrix3D.cxx:43
 AliCorrectionMatrix3D.cxx:44
 AliCorrectionMatrix3D.cxx:45
 AliCorrectionMatrix3D.cxx:46
 AliCorrectionMatrix3D.cxx:47
 AliCorrectionMatrix3D.cxx:48
 AliCorrectionMatrix3D.cxx:49
 AliCorrectionMatrix3D.cxx:50
 AliCorrectionMatrix3D.cxx:51
 AliCorrectionMatrix3D.cxx:52
 AliCorrectionMatrix3D.cxx:53
 AliCorrectionMatrix3D.cxx:54
 AliCorrectionMatrix3D.cxx:55
 AliCorrectionMatrix3D.cxx:56
 AliCorrectionMatrix3D.cxx:57
 AliCorrectionMatrix3D.cxx:58
 AliCorrectionMatrix3D.cxx:59
 AliCorrectionMatrix3D.cxx:60
 AliCorrectionMatrix3D.cxx:61
 AliCorrectionMatrix3D.cxx:62
 AliCorrectionMatrix3D.cxx:63
 AliCorrectionMatrix3D.cxx:64
 AliCorrectionMatrix3D.cxx:65
 AliCorrectionMatrix3D.cxx:66
 AliCorrectionMatrix3D.cxx:67
 AliCorrectionMatrix3D.cxx:68
 AliCorrectionMatrix3D.cxx:69
 AliCorrectionMatrix3D.cxx:70
 AliCorrectionMatrix3D.cxx:71
 AliCorrectionMatrix3D.cxx:72
 AliCorrectionMatrix3D.cxx:73
 AliCorrectionMatrix3D.cxx:74
 AliCorrectionMatrix3D.cxx:75
 AliCorrectionMatrix3D.cxx:76
 AliCorrectionMatrix3D.cxx:77
 AliCorrectionMatrix3D.cxx:78
 AliCorrectionMatrix3D.cxx:79
 AliCorrectionMatrix3D.cxx:80
 AliCorrectionMatrix3D.cxx:81
 AliCorrectionMatrix3D.cxx:82
 AliCorrectionMatrix3D.cxx:83
 AliCorrectionMatrix3D.cxx:84
 AliCorrectionMatrix3D.cxx:85
 AliCorrectionMatrix3D.cxx:86
 AliCorrectionMatrix3D.cxx:87
 AliCorrectionMatrix3D.cxx:88
 AliCorrectionMatrix3D.cxx:89
 AliCorrectionMatrix3D.cxx:90
 AliCorrectionMatrix3D.cxx:91
 AliCorrectionMatrix3D.cxx:92
 AliCorrectionMatrix3D.cxx:93
 AliCorrectionMatrix3D.cxx:94
 AliCorrectionMatrix3D.cxx:95
 AliCorrectionMatrix3D.cxx:96
 AliCorrectionMatrix3D.cxx:97
 AliCorrectionMatrix3D.cxx:98
 AliCorrectionMatrix3D.cxx:99
 AliCorrectionMatrix3D.cxx:100
 AliCorrectionMatrix3D.cxx:101
 AliCorrectionMatrix3D.cxx:102
 AliCorrectionMatrix3D.cxx:103
 AliCorrectionMatrix3D.cxx:104
 AliCorrectionMatrix3D.cxx:105
 AliCorrectionMatrix3D.cxx:106
 AliCorrectionMatrix3D.cxx:107
 AliCorrectionMatrix3D.cxx:108
 AliCorrectionMatrix3D.cxx:109
 AliCorrectionMatrix3D.cxx:110
 AliCorrectionMatrix3D.cxx:111
 AliCorrectionMatrix3D.cxx:112
 AliCorrectionMatrix3D.cxx:113
 AliCorrectionMatrix3D.cxx:114
 AliCorrectionMatrix3D.cxx:115
 AliCorrectionMatrix3D.cxx:116
 AliCorrectionMatrix3D.cxx:117
 AliCorrectionMatrix3D.cxx:118
 AliCorrectionMatrix3D.cxx:119
 AliCorrectionMatrix3D.cxx:120
 AliCorrectionMatrix3D.cxx:121
 AliCorrectionMatrix3D.cxx:122
 AliCorrectionMatrix3D.cxx:123
 AliCorrectionMatrix3D.cxx:124
 AliCorrectionMatrix3D.cxx:125
 AliCorrectionMatrix3D.cxx:126
 AliCorrectionMatrix3D.cxx:127
 AliCorrectionMatrix3D.cxx:128
 AliCorrectionMatrix3D.cxx:129
 AliCorrectionMatrix3D.cxx:130
 AliCorrectionMatrix3D.cxx:131
 AliCorrectionMatrix3D.cxx:132
 AliCorrectionMatrix3D.cxx:133
 AliCorrectionMatrix3D.cxx:134
 AliCorrectionMatrix3D.cxx:135
 AliCorrectionMatrix3D.cxx:136
 AliCorrectionMatrix3D.cxx:137
 AliCorrectionMatrix3D.cxx:138
 AliCorrectionMatrix3D.cxx:139
 AliCorrectionMatrix3D.cxx:140
 AliCorrectionMatrix3D.cxx:141
 AliCorrectionMatrix3D.cxx:142
 AliCorrectionMatrix3D.cxx:143
 AliCorrectionMatrix3D.cxx:144
 AliCorrectionMatrix3D.cxx:145
 AliCorrectionMatrix3D.cxx:146
 AliCorrectionMatrix3D.cxx:147
 AliCorrectionMatrix3D.cxx:148
 AliCorrectionMatrix3D.cxx:149
 AliCorrectionMatrix3D.cxx:150
 AliCorrectionMatrix3D.cxx:151
 AliCorrectionMatrix3D.cxx:152
 AliCorrectionMatrix3D.cxx:153
 AliCorrectionMatrix3D.cxx:154
 AliCorrectionMatrix3D.cxx:155
 AliCorrectionMatrix3D.cxx:156
 AliCorrectionMatrix3D.cxx:157
 AliCorrectionMatrix3D.cxx:158
 AliCorrectionMatrix3D.cxx:159
 AliCorrectionMatrix3D.cxx:160
 AliCorrectionMatrix3D.cxx:161
 AliCorrectionMatrix3D.cxx:162
 AliCorrectionMatrix3D.cxx:163
 AliCorrectionMatrix3D.cxx:164
 AliCorrectionMatrix3D.cxx:165
 AliCorrectionMatrix3D.cxx:166
 AliCorrectionMatrix3D.cxx:167
 AliCorrectionMatrix3D.cxx:168
 AliCorrectionMatrix3D.cxx:169
 AliCorrectionMatrix3D.cxx:170
 AliCorrectionMatrix3D.cxx:171
 AliCorrectionMatrix3D.cxx:172
 AliCorrectionMatrix3D.cxx:173
 AliCorrectionMatrix3D.cxx:174
 AliCorrectionMatrix3D.cxx:175
 AliCorrectionMatrix3D.cxx:176
 AliCorrectionMatrix3D.cxx:177
 AliCorrectionMatrix3D.cxx:178
 AliCorrectionMatrix3D.cxx:179
 AliCorrectionMatrix3D.cxx:180
 AliCorrectionMatrix3D.cxx:181
 AliCorrectionMatrix3D.cxx:182
 AliCorrectionMatrix3D.cxx:183
 AliCorrectionMatrix3D.cxx:184
 AliCorrectionMatrix3D.cxx:185
 AliCorrectionMatrix3D.cxx:186
 AliCorrectionMatrix3D.cxx:187
 AliCorrectionMatrix3D.cxx:188
 AliCorrectionMatrix3D.cxx:189
 AliCorrectionMatrix3D.cxx:190
 AliCorrectionMatrix3D.cxx:191
 AliCorrectionMatrix3D.cxx:192
 AliCorrectionMatrix3D.cxx:193
 AliCorrectionMatrix3D.cxx:194
 AliCorrectionMatrix3D.cxx:195
 AliCorrectionMatrix3D.cxx:196
 AliCorrectionMatrix3D.cxx:197
 AliCorrectionMatrix3D.cxx:198
 AliCorrectionMatrix3D.cxx:199
 AliCorrectionMatrix3D.cxx:200
 AliCorrectionMatrix3D.cxx:201
 AliCorrectionMatrix3D.cxx:202
 AliCorrectionMatrix3D.cxx:203
 AliCorrectionMatrix3D.cxx:204
 AliCorrectionMatrix3D.cxx:205
 AliCorrectionMatrix3D.cxx:206
 AliCorrectionMatrix3D.cxx:207
 AliCorrectionMatrix3D.cxx:208
 AliCorrectionMatrix3D.cxx:209
 AliCorrectionMatrix3D.cxx:210
 AliCorrectionMatrix3D.cxx:211
 AliCorrectionMatrix3D.cxx:212
 AliCorrectionMatrix3D.cxx:213
 AliCorrectionMatrix3D.cxx:214
 AliCorrectionMatrix3D.cxx:215
 AliCorrectionMatrix3D.cxx:216
 AliCorrectionMatrix3D.cxx:217
 AliCorrectionMatrix3D.cxx:218
 AliCorrectionMatrix3D.cxx:219
 AliCorrectionMatrix3D.cxx:220
 AliCorrectionMatrix3D.cxx:221
 AliCorrectionMatrix3D.cxx:222
 AliCorrectionMatrix3D.cxx:223
 AliCorrectionMatrix3D.cxx:224
 AliCorrectionMatrix3D.cxx:225
 AliCorrectionMatrix3D.cxx:226
 AliCorrectionMatrix3D.cxx:227
 AliCorrectionMatrix3D.cxx:228
 AliCorrectionMatrix3D.cxx:229
 AliCorrectionMatrix3D.cxx:230
 AliCorrectionMatrix3D.cxx:231
 AliCorrectionMatrix3D.cxx:232
 AliCorrectionMatrix3D.cxx:233
 AliCorrectionMatrix3D.cxx:234
 AliCorrectionMatrix3D.cxx:235
 AliCorrectionMatrix3D.cxx:236
 AliCorrectionMatrix3D.cxx:237
 AliCorrectionMatrix3D.cxx:238
 AliCorrectionMatrix3D.cxx:239
 AliCorrectionMatrix3D.cxx:240
 AliCorrectionMatrix3D.cxx:241
 AliCorrectionMatrix3D.cxx:242
 AliCorrectionMatrix3D.cxx:243
 AliCorrectionMatrix3D.cxx:244
 AliCorrectionMatrix3D.cxx:245
 AliCorrectionMatrix3D.cxx:246
 AliCorrectionMatrix3D.cxx:247
 AliCorrectionMatrix3D.cxx:248
 AliCorrectionMatrix3D.cxx:249
 AliCorrectionMatrix3D.cxx:250
 AliCorrectionMatrix3D.cxx:251
 AliCorrectionMatrix3D.cxx:252
 AliCorrectionMatrix3D.cxx:253
 AliCorrectionMatrix3D.cxx:254
 AliCorrectionMatrix3D.cxx:255
 AliCorrectionMatrix3D.cxx:256
 AliCorrectionMatrix3D.cxx:257
 AliCorrectionMatrix3D.cxx:258
 AliCorrectionMatrix3D.cxx:259
 AliCorrectionMatrix3D.cxx:260
 AliCorrectionMatrix3D.cxx:261
 AliCorrectionMatrix3D.cxx:262
 AliCorrectionMatrix3D.cxx:263
 AliCorrectionMatrix3D.cxx:264
 AliCorrectionMatrix3D.cxx:265
 AliCorrectionMatrix3D.cxx:266
 AliCorrectionMatrix3D.cxx:267
 AliCorrectionMatrix3D.cxx:268
 AliCorrectionMatrix3D.cxx:269
 AliCorrectionMatrix3D.cxx:270
 AliCorrectionMatrix3D.cxx:271
 AliCorrectionMatrix3D.cxx:272
 AliCorrectionMatrix3D.cxx:273
 AliCorrectionMatrix3D.cxx:274
 AliCorrectionMatrix3D.cxx:275
 AliCorrectionMatrix3D.cxx:276
 AliCorrectionMatrix3D.cxx:277
 AliCorrectionMatrix3D.cxx:278
 AliCorrectionMatrix3D.cxx:279
 AliCorrectionMatrix3D.cxx:280
 AliCorrectionMatrix3D.cxx:281
 AliCorrectionMatrix3D.cxx:282
 AliCorrectionMatrix3D.cxx:283
 AliCorrectionMatrix3D.cxx:284
 AliCorrectionMatrix3D.cxx:285
 AliCorrectionMatrix3D.cxx:286
 AliCorrectionMatrix3D.cxx:287
 AliCorrectionMatrix3D.cxx:288
 AliCorrectionMatrix3D.cxx:289
 AliCorrectionMatrix3D.cxx:290
 AliCorrectionMatrix3D.cxx:291
 AliCorrectionMatrix3D.cxx:292
 AliCorrectionMatrix3D.cxx:293
 AliCorrectionMatrix3D.cxx:294
 AliCorrectionMatrix3D.cxx:295
 AliCorrectionMatrix3D.cxx:296
 AliCorrectionMatrix3D.cxx:297
 AliCorrectionMatrix3D.cxx:298
 AliCorrectionMatrix3D.cxx:299
 AliCorrectionMatrix3D.cxx:300
 AliCorrectionMatrix3D.cxx:301
 AliCorrectionMatrix3D.cxx:302
 AliCorrectionMatrix3D.cxx:303
 AliCorrectionMatrix3D.cxx:304
 AliCorrectionMatrix3D.cxx:305
 AliCorrectionMatrix3D.cxx:306
 AliCorrectionMatrix3D.cxx:307
 AliCorrectionMatrix3D.cxx:308
 AliCorrectionMatrix3D.cxx:309
 AliCorrectionMatrix3D.cxx:310
 AliCorrectionMatrix3D.cxx:311
 AliCorrectionMatrix3D.cxx:312
 AliCorrectionMatrix3D.cxx:313
 AliCorrectionMatrix3D.cxx:314
 AliCorrectionMatrix3D.cxx:315
 AliCorrectionMatrix3D.cxx:316
 AliCorrectionMatrix3D.cxx:317
 AliCorrectionMatrix3D.cxx:318
 AliCorrectionMatrix3D.cxx:319
 AliCorrectionMatrix3D.cxx:320
 AliCorrectionMatrix3D.cxx:321
 AliCorrectionMatrix3D.cxx:322
 AliCorrectionMatrix3D.cxx:323
 AliCorrectionMatrix3D.cxx:324
 AliCorrectionMatrix3D.cxx:325
 AliCorrectionMatrix3D.cxx:326
 AliCorrectionMatrix3D.cxx:327
 AliCorrectionMatrix3D.cxx:328
 AliCorrectionMatrix3D.cxx:329
 AliCorrectionMatrix3D.cxx:330
 AliCorrectionMatrix3D.cxx:331
 AliCorrectionMatrix3D.cxx:332
 AliCorrectionMatrix3D.cxx:333
 AliCorrectionMatrix3D.cxx:334
 AliCorrectionMatrix3D.cxx:335
 AliCorrectionMatrix3D.cxx:336
 AliCorrectionMatrix3D.cxx:337
 AliCorrectionMatrix3D.cxx:338
 AliCorrectionMatrix3D.cxx:339
 AliCorrectionMatrix3D.cxx:340
 AliCorrectionMatrix3D.cxx:341
 AliCorrectionMatrix3D.cxx:342
 AliCorrectionMatrix3D.cxx:343
 AliCorrectionMatrix3D.cxx:344
 AliCorrectionMatrix3D.cxx:345
 AliCorrectionMatrix3D.cxx:346
 AliCorrectionMatrix3D.cxx:347
 AliCorrectionMatrix3D.cxx:348
 AliCorrectionMatrix3D.cxx:349
 AliCorrectionMatrix3D.cxx:350
 AliCorrectionMatrix3D.cxx:351
 AliCorrectionMatrix3D.cxx:352
 AliCorrectionMatrix3D.cxx:353
 AliCorrectionMatrix3D.cxx:354
 AliCorrectionMatrix3D.cxx:355
 AliCorrectionMatrix3D.cxx:356
 AliCorrectionMatrix3D.cxx:357
 AliCorrectionMatrix3D.cxx:358
 AliCorrectionMatrix3D.cxx:359
 AliCorrectionMatrix3D.cxx:360
 AliCorrectionMatrix3D.cxx:361
 AliCorrectionMatrix3D.cxx:362
 AliCorrectionMatrix3D.cxx:363
 AliCorrectionMatrix3D.cxx:364
 AliCorrectionMatrix3D.cxx:365
 AliCorrectionMatrix3D.cxx:366
 AliCorrectionMatrix3D.cxx:367
 AliCorrectionMatrix3D.cxx:368
 AliCorrectionMatrix3D.cxx:369
 AliCorrectionMatrix3D.cxx:370
 AliCorrectionMatrix3D.cxx:371
 AliCorrectionMatrix3D.cxx:372
 AliCorrectionMatrix3D.cxx:373
 AliCorrectionMatrix3D.cxx:374
 AliCorrectionMatrix3D.cxx:375
 AliCorrectionMatrix3D.cxx:376
 AliCorrectionMatrix3D.cxx:377