PWG0/AliCorrectionMatrix3D.cxx

   1 /* $Id$ */
   2
   3 // ------------------------------------------------------
   4 //
   5 // Class to handle 3d-corrections.
   6 //
   7 // ------------------------------------------------------
   8 //
   9
  10 #include <TH3F.h>
  11 #include <TH1F.h>
  12
  13 #include <AliLog.h>
  14
  15 #include "AliCorrectionMatrix3D.h"
  16 #include "AliPWG0Helper.h"
  17
  18 //____________________________________________________________________
  19 ClassImp(AliCorrectionMatrix3D)
  20
  21 //____________________________________________________________________
  22 AliCorrectionMatrix3D::AliCorrectionMatrix3D() :
  23   AliCorrectionMatrix()
  24 {
  25   // default constructor
  26 }
  27
  28 //____________________________________________________________________
  29 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const AliCorrectionMatrix3D& c)
  30   : AliCorrectionMatrix(c)
  31 {
  32   // copy constructor
  33   ((AliCorrectionMatrix3D &)c).Copy(*this);
  34 }
  35
  36 //____________________________________________________________________
  37 AliCorrectionMatrix3D &AliCorrectionMatrix3D::operator=(const AliCorrectionMatrix3D &c)
  38 {
  39   // assigment operator
  40
  41   if (this != &c)
  42     ((AliCorrectionMatrix3D &) c).Copy(*this);
  43
  44   return *this;
  45 }
  46
  47 //____________________________________________________________________
  48 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
  49               Int_t nBinX, Float_t Xmin, Float_t Xmax,
  50               Int_t nBinY, Float_t Ymin, Float_t Ymax,
  51               Int_t nBinZ, Float_t Zmin, Float_t Zmax)
  52   : AliCorrectionMatrix(name, title)
  53 {
  54   //
  55   // constructor
  56   //
  57
  58   // do not add this hists to the directory
  59   Bool_t oldStatus = TH1::AddDirectoryStatus();
  60   TH1::AddDirectory(kFALSE);
  61
  62   fhMeas  = new TH3F("measured",   Form("%s measured",title),   nBinX, Xmin, Xmax, nBinY, Ymin, Ymax, nBinZ, Zmin, Zmax);
  63   fhGene  = new TH3F("generated",  Form("%s generated",title),  nBinX, Xmin, Xmax, nBinY, Ymin, Ymax, nBinZ, Zmin, Zmax);
  64   fhCorr  = new TH3F("correction", Form("%s correction",title), nBinX, Xmin, Xmax, nBinY, Ymin, Ymax, nBinZ, Zmin, Zmax);
  65
  66   TH1::AddDirectory(oldStatus);
  67
  68   fhMeas->Sumw2();
  69   fhGene->Sumw2();
  70   fhCorr->Sumw2();
  71 }
  72
  73 AliCorrectionMatrix3D::AliCorrectionMatrix3D(const Char_t* name, const Char_t* title,
  74       Int_t nBinX, Float_t Xmin, Float_t Xmax,
  75       Int_t nBinY, Float_t Ymin, Float_t Ymax,
  76       Int_t nBinZ, const Float_t* zbins)
  77   : AliCorrectionMatrix(name, title)
  78 {
  79   // constructor with variable bin sizes
  80
  81   Float_t* binLimitsX = new Float_t[nBinX+1];
  82   for (Int_t i=0; i<=nBinX; ++i)
  83     binLimitsX[i] = Xmin + (Xmax - Xmin) / nBinX * i;
  84
  85   Float_t* binLimitsY = new Float_t[nBinY+1];
  86   for (Int_t i=0; i<=nBinY; ++i)
  87     binLimitsY[i] = Ymin + (Ymax - Ymin) / nBinY * i;
  88
  89   // do not add this hists to the directory
  90   Bool_t oldStatus = TH1::AddDirectoryStatus();
  91   TH1::AddDirectory(kFALSE);
  92
  93   fhMeas  = new TH3F("measured",   Form("%s measured",title),   nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);
  94   fhGene  = new TH3F("generated",  Form("%s generated",title),  nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);
  95   fhCorr  = new TH3F("correction", Form("%s correction",title), nBinX, binLimitsX, nBinY, binLimitsY, nBinZ, zbins);
  96
  97   TH1::AddDirectory(oldStatus);
  98
  99   delete[] binLimitsX;
 100   delete[] binLimitsY;
 101
 102   fhMeas->Sumw2();
 103   fhGene->Sumw2();
 104   fhCorr->Sumw2();
 105 }
 106
 107 //____________________________________________________________________
 108 AliCorrectionMatrix3D::~AliCorrectionMatrix3D()
 109 {
 110   //
 111   // destructor
 112   //
 113
 114   // histograms already deleted in base class
 115 }
 116
 117 //____________________________________________________________________
 118 TH3F* AliCorrectionMatrix3D::GetGeneratedHistogram()
 119 {
 120   // return generated histogram casted to correct type
 121   return dynamic_cast<TH3F*> (fhGene);
 122 }
 123
 124 //____________________________________________________________________
 125 TH3F* AliCorrectionMatrix3D::GetMeasuredHistogram()
 126 {
 127   // return measured histogram casted to correct type
 128   return dynamic_cast<TH3F*> (fhMeas);
 129 }
 130
 131 //____________________________________________________________________
 132 TH3F* AliCorrectionMatrix3D::GetCorrectionHistogram()
 133 {
 134   // return correction histogram casted to correct type
 135   return dynamic_cast<TH3F*> (fhCorr);
 136 }
 137
 138 //____________________________________________________________________
 139 void AliCorrectionMatrix3D::FillMeas(Float_t ax, Float_t ay, Float_t az)
 140 {
 141   // add value to measured histogram
 142   GetMeasuredHistogram()->Fill(ax, ay, az);
 143 }
 144
 145 //____________________________________________________________________
 146 void AliCorrectionMatrix3D::FillGene(Float_t ax, Float_t ay, Float_t az)
 147 {
 148   // add value to generated histogram
 149   GetGeneratedHistogram()->Fill(ax, ay, az);
 150 }
 151
 152 //____________________________________________________________________
 153 Float_t AliCorrectionMatrix3D::GetCorrection(Float_t ax, Float_t ay, Float_t az) const
 154 {
 155   // returns a value of the correction map
 156   return fhCorr->GetBinContent(fhCorr->FindBin(ax, ay, az));
 157 }
 158
 159 //____________________________________________________________________
 160 //void AliCorrectionMatrix3D::RemoveEdges(Float_t cut, Int_t nBinsXedge, Int_t nBinsYedge, Int_t nBinsZedge)
 161 void AliCorrectionMatrix3D::RemoveEdges(Float_t, Int_t, Int_t, Int_t)
 162 {
 163   // so what do we do here...
 164 }
 165
 166 //____________________________________________________________________
 167 void AliCorrectionMatrix3D::SaveHistograms()
 168 {
 169   //
 170   // saves the histograms
 171   //
 172
 173   AliCorrectionMatrix::SaveHistograms();
 174
 175   if (GetGeneratedHistogram() && GetMeasuredHistogram())
 176   {
 177     gDirectory->cd(GetName());
 178
 179     AliPWG0Helper::CreateDividedProjections(GetGeneratedHistogram(), GetMeasuredHistogram(), 0, kFALSE, kTRUE);
 180
 181     gDirectory->cd("..");
 182   }
 183 }
 184
 185 //____________________________________________________________________
 186 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)
 187 {
 188   //
 189   // counts the number of empty Bins in a given region
 190   //
 191
 192   TH3F* hist = GetGeneratedHistogram();
 193   if (!hist)
 194     return -1;
 195
 196   Int_t emptyBins = 0;
 197   for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
 198     for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
 199       for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
 200         if (hist->GetBinContent(x, y, z) == 0)
 201         {
 202           if (!quiet)
 203             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));
 204           ++emptyBins;
 205         }
 206
 207   return emptyBins;
 208 }
 209
 210 //____________________________________________________________________
 211 TH1F* AliCorrectionMatrix3D::PlotBinErrors(Float_t xmin, Float_t xmax, Float_t ymin, Float_t ymax, Float_t zmin, Float_t zmax)
 212 {
 213   //
 214   // makes a 1d plots of the relative bin errors
 215   //
 216
 217   TH3F* hist = GetCorrectionHistogram();
 218   if (!hist)
 219     return 0;
 220
 221   TH1F* target = new TH1F("relerrors", "relerrors", 100, 0, 10);
 222
 223   for (Int_t x=hist->GetXaxis()->FindBin(xmin); x<=hist->GetXaxis()->FindBin(xmax); ++x)
 224     for (Int_t y=hist->GetYaxis()->FindBin(ymin); y<=hist->GetYaxis()->FindBin(ymax); ++y)
 225       for (Int_t z=hist->GetZaxis()->FindBin(zmin); z<=hist->GetZaxis()->FindBin(zmax); ++z)
 226         if (hist->GetBinContent(x, y, z) != 0)
 227           target->Fill(100 * hist->GetBinError(x, y, z) / hist->GetBinContent(x, y, z));
 228
 229   return target;
 230 }
 231