PWG0/multPbPb/AliAnalysisMultPbTrackHistoManager.cxx

   1 #include "AliAnalysisMultPbTrackHistoManager.h"
   2 #include "AliLog.h"
   3 #include "TH1.h"
   4 #include "TH3D.h"
   5 #include "TH1I.h"
   6 #include "TROOT.h"
   7 #include <iostream>
   8 using namespace std;
   9 ClassImp(AliAnalysisMultPbTrackHistoManager)
  10
  11 const char * AliAnalysisMultPbTrackHistoManager::kStatStepNames[]     = { "All Events", "After centrality selection", "With Vertex" };
  12 const char * AliAnalysisMultPbTrackHistoManager::kHistoPtEtaVzNames[] = { "hGenPtEtaVz", "hRecPtEtaVz", "hRecPtEtaVzPrim",
  13                                                                           "hRecPtEtaVzSecWeak", "hRecPtEtaVzSecMaterial", "hRecPtEtaVzFake"};
  14 const char * AliAnalysisMultPbTrackHistoManager::kHistoDCANames[]     = { "hGenDCA", "hRecDCA", "hRecDCAPrim", "hRecDCASecWeak","hRecDCASecMaterial", "hRecDCAFake"};
  15
  16
  17
  18 AliAnalysisMultPbTrackHistoManager::AliAnalysisMultPbTrackHistoManager() : AliHistoListWrapper(), fHNameSuffix(""){
  19   // standard ctor
  20
  21 }
  22
  23 AliAnalysisMultPbTrackHistoManager::AliAnalysisMultPbTrackHistoManager(const char * name, const char * title): AliHistoListWrapper(name,title), fHNameSuffix("")  {
  24   //named ctor
  25 };
  26
  27 AliAnalysisMultPbTrackHistoManager::AliAnalysisMultPbTrackHistoManager(const AliAnalysisMultPbTrackHistoManager& obj) : AliHistoListWrapper (obj) {
  28   // copy ctor
  29   AliError("Not Implemented");
  30 };
  31
  32 AliAnalysisMultPbTrackHistoManager::~AliAnalysisMultPbTrackHistoManager() {
  33   // dtor
  34
  35 }
  36
  37 TH3D * AliAnalysisMultPbTrackHistoManager::GetHistoPtEtaVz(Histo_t id) {
  38   // Returns a 3D histo of Pt/eta/vtx. It it does not exist, books it.
  39
  40   TH3D * h = (TH3D*) GetHisto(kHistoPtEtaVzNames[id]);
  41   if (!h) {
  42     h = BookHistoPtEtaVz(kHistoPtEtaVzNames[id], Form("Pt Eta Vz distribution (%s)",kHistoPtEtaVzNames[id]));
  43   }
  44
  45   return h;
  46
  47 }
  48
  49 TH1D * AliAnalysisMultPbTrackHistoManager::GetHistoDCA(Histo_t id) {
  50   // Returns a 3D histo of Pt/eta/vtx. It it does not exist, books it.
  51
  52   TH1D * h = (TH1D*) GetHisto(kHistoDCANames[id]);
  53   if (!h) {
  54     h = BookHistoDCA(kHistoDCANames[id], Form("Pt Eta Vz distribution (%s)",kHistoDCANames[id]));
  55   }
  56
  57   return h;
  58
  59 }
  60
  61 TH1D * AliAnalysisMultPbTrackHistoManager::GetHistoPt (Histo_t id,
  62                                                        Float_t minEta, Float_t maxEta,
  63                                                        Float_t minVz, Float_t maxVz,
  64                                                        Bool_t scaleWidth) {
  65   // Returns a projection of the 3D histo pt/eta/vz.
  66   // WARNING: since that is a histo, the requested range will be discretized to the binning.
  67   // Always avoids under (over) flows
  68   // If scaleWidth = kTRUE, the projection is scaled for the bin width (default)
  69
  70   TH3D * h3D = GetHistoPtEtaVz(id);
  71
  72   // Get range in terms of bin numners.  If the float range is
  73   // less than -11111 take the range from the first to the last bin (i.e. no
  74   // under/over-flows)
  75   Int_t min1 = minEta  < -11111 ? 1 : h3D ->GetYaxis()->FindBin(minEta);
  76   Int_t min2  = minVz  < -11111 ? 1 : h3D ->GetZaxis()->FindBin(minVz) ;
  77
  78   Int_t max1 = maxEta  < -11111 ? h3D->GetNbinsY() : h3D ->GetYaxis()->FindBin(maxEta-0.00001);
  79   Int_t max2  = maxVz  < -11111 ? h3D->GetNbinsZ() : h3D ->GetZaxis()->FindBin(maxVz -0.00001);
  80
  81
  82   TString hname = h3D->GetName();
  83   hname = hname +  "_pt_" + long (min1)  + "_" + long(max1) + "_" + long (min2)  + "_" + long(max2);
  84
  85
  86   if (gROOT->FindObjectAny(hname.Data())){
  87     AliError(Form("An object called %s already exists",hname.Data()));
  88   }
  89
  90   TH1D * h = h3D->ProjectionX(hname.Data(), min1, max1, min2, max2, "E");
  91   if(scaleWidth) h ->Scale(1.,"width");
  92
  93   return h;
  94
  95 }
  96
  97 TH1D * AliAnalysisMultPbTrackHistoManager::GetHistoVz (Histo_t id,
  98                                                        Float_t minPt, Float_t maxPt,
  99                                                        Float_t minEta, Float_t maxEta,
 100                                                        Bool_t scaleWidth) {
 101   // Returns a projection of the 3D histo pt/eta/vz.
 102   // WARNING: since that is a histo, the requested range will be discretized to the binning.
 103   // Always avoids under (over) flows
 104   // If scaleWidth = kTRUE, the projection is scaled for the bin width (default)
 105
 106   TH3D * h3D = GetHistoPtEtaVz(id);
 107
 108   // Get range in terms of bin numners.  If the float range is
 109   // less than -11111 take the range from the first to the last bin (i.e. no
 110   // under/over-flows)
 111   Int_t min1  = minPt  < -11111 ? 1 : h3D ->GetXaxis()->FindBin(minPt) ;
 112   Int_t min2  = minEta < -11111 ? 1 : h3D ->GetYaxis()->FindBin(minEta);
 113
 114   Int_t max1  = maxPt  < -11111 ? h3D->GetNbinsX() : h3D ->GetXaxis()->FindBin(maxPt -0.00001);
 115   Int_t max2  = maxEta < -11111 ? h3D->GetNbinsY() : h3D ->GetYaxis()->FindBin(maxEta-0.00001);
 116
 117
 118   TString hname = h3D->GetName();
 119   hname = hname +  "_Vz_" + long (min1)  + "_" + long(max1) + "_" + long (min2)  + "_" + long(max2);
 120
 121   if (gROOT->FindObjectAny(hname.Data())){
 122     AliError(Form("An object called %s already exists",hname.Data()));
 123   }
 124
 125   TH1D * h = h3D->ProjectionZ(hname.Data(), min1, max1, min2, max2, "E");
 126   if(scaleWidth) h ->Scale(1.,"width");
 127   return h;
 128
 129
 130 }
 131
 132 TH1D * AliAnalysisMultPbTrackHistoManager::GetHistoEta (Histo_t id,
 133                                                         Float_t minPt, Float_t maxPt,
 134                                                         Float_t minVz, Float_t maxVz,
 135                                                         Bool_t scaleWidth) {
 136   // Returns a projection of the 3D histo pt/eta/vz.
 137   // WARNING: since that is a histo, the requested range will be discretized to the binning.
 138   // Always avoids under (over) flows
 139   // If scaleWidth = kTRUE, the projection is scaled for the bin width (default)
 140
 141   TH3D * h3D = GetHistoPtEtaVz(id);
 142
 143   // Get range in terms of bin numners.  If the float range is
 144   // less than -11111 take the range from the first to the last bin (i.e. no
 145   // under/over-flows)
 146   Int_t min1 = minPt < -11111 ? 1 : h3D ->GetXaxis()->FindBin(minPt) ;
 147   Int_t min2 = minVz < -11111 ? 1 : h3D ->GetYaxis()->FindBin(minVz);
 148
 149   Int_t max1 = maxPt < -11111 ? h3D->GetNbinsX() : h3D ->GetXaxis()->FindBin(maxPt -0.00001);
 150   Int_t max2 = maxVz < -11111 ? h3D->GetNbinsY() : h3D ->GetYaxis()->FindBin(maxVz-0.00001);
 151
 152   TString hname = h3D->GetName();
 153   hname = hname +  "_Eta_" + long (min1)  + "_" + long(max1) + "_" + long (min2)  + "_" + long(max2);
 154
 155   if (gROOT->FindObjectAny(hname.Data())){
 156     AliError(Form("An object called %s already exists",hname.Data()));
 157   }
 158
 159   TH1D * h = h3D->ProjectionY(hname.Data(), min1, max1, min2, max2, "E");
 160   if(scaleWidth) h ->Scale(1.,"width");
 161   return h;
 162 }
 163
 164
 165 TH1I * AliAnalysisMultPbTrackHistoManager::GetHistoStats() {
 166   // Returns histogram with event statistiscs (processed events at each step)
 167
 168   TH1I * h =  (TH1I*) GetHisto("hStats");
 169   if (!h) h = BookHistoStats();
 170   return h;
 171
 172 }
 173
 174
 175
 176 TH1 * AliAnalysisMultPbTrackHistoManager::GetHisto(const char * name) {
 177   //Search list for histo
 178   // TODO: keep track of histo id rather than searching by name?
 179   return (TH1*) fList->FindObject(TString(name)+fHNameSuffix);
 180
 181 }
 182
 183 void AliAnalysisMultPbTrackHistoManager::ScaleHistos(Double_t nev, Option_t * option) {
 184   // Scales all histos in the list for nev
 185   // option can be used to pass further options to TH1::Scale
 186   TH1 * h = 0;
 187   TIter iter = fList->MakeIterator();
 188   while ((h = (TH1*) iter.Next())) {
 189     if (!h->InheritsFrom("TH1")) {
 190       AliFatal (Form("%s does not inherits from TH1, cannot scale",h->GetName()));
 191     }
 192     cout << "Scaling " << h->GetName() << " " << nev << endl;
 193
 194     h->Scale(1./nev,option);
 195   }
 196
 197 }
 198
 199 TH3D * AliAnalysisMultPbTrackHistoManager::BookHistoPtEtaVz(const char * name, const char * title) {
 200   // Books a 3D histo of Pt/eta/vtx
 201   // TODO: make the binning settable, variable binning?
 202
 203   Bool_t oldStatus = TH1::AddDirectoryStatus();
 204   TH1::AddDirectory(kFALSE);
 205
 206   TString hname = name;
 207   hname+=fHNameSuffix;
 208
 209   AliInfo(Form("Booking %s",hname.Data()));
 210
 211
 212   TH3D * h = new TH3D (hname,title,
 213                        50,0,10, // Pt
 214                        20,-1, 1,   // Eta
 215                        10,-10,10 // Vz
 216                        );
 217
 218   h->SetYTitle("#eta");
 219   h->SetXTitle("p_{T}");
 220   h->SetZTitle("V_{z} (cm)");
 221   h->Sumw2();
 222
 223   fList->Add(h);
 224
 225   TH1::AddDirectory(oldStatus);
 226   return h;
 227 }
 228
 229 TH1D * AliAnalysisMultPbTrackHistoManager::BookHistoDCA(const char * name, const char * title) {
 230   // Books a DCA histo
 231
 232   Bool_t oldStatus = TH1::AddDirectoryStatus();
 233   TH1::AddDirectory(kFALSE);
 234
 235   TString hname = name;
 236   hname+=fHNameSuffix;
 237
 238   AliInfo(Form("Booking %s",hname.Data()));
 239
 240
 241   TH1D * h = new TH1D (hname,title, 200,0,200);
 242
 243   h->SetXTitle("#Delta DCA");
 244   h->Sumw2();
 245
 246   fList->Add(h);
 247
 248   TH1::AddDirectory(oldStatus);
 249   return h;
 250 }
 251
 252 TH1I * AliAnalysisMultPbTrackHistoManager::BookHistoStats() {
 253   // Books histogram with event statistiscs (processed events at each step)
 254
 255   Bool_t oldStatus = TH1::AddDirectoryStatus();
 256   TH1::AddDirectory(kFALSE);
 257
 258   AliInfo(Form("Booking Stat histo"));
 259
 260   TH1I * h = new TH1I (TString("hStats")+fHNameSuffix, "Number of processed events", kNStatBins, -0.5, kNStatBins-0.5);
 261   for(Int_t istatbin = 0; istatbin < kNStatBins; istatbin++){
 262     h->GetXaxis()->SetBinLabel(istatbin+1,kStatStepNames[istatbin]);
 263   }
 264   TH1::AddDirectory(oldStatus);
 265   fList->Add(h);
 266   return h;
 267 }
 268
 269
 270
 271
 272