PWG0/dNdEta/dNdEtaCorrection.cxx

   1 /* $Id$ */
   2
   3 #include "dNdEtaCorrection.h"
   4
   5 #include <TCanvas.h>
   6
   7 //____________________________________________________________________
   8 ClassImp(dNdEtaCorrection)
   9
  10 //____________________________________________________________________
  11 dNdEtaCorrection::dNdEtaCorrection(Char_t* name) {
  12
  13   fName = TString(name);
  14
  15   hEtaVsVtx_meas  = new TH2F("etaVsVtx_meas", "etaVsVtx_meas" ,80,-20,20,120,-6,6);
  16   hEtaVsVtx_gene  = new TH2F("etaVsVtx_gene", "etaVsVtx_gene" ,80,-20,20,120,-6,6);
  17   hEtaVsVtx_corr  = new TH2F("etaVsVtx_corr", "etaVsVtx_corr", 80,-20,20,120,-6,6);
  18
  19   hEtaVsVtx_ratio = new TH2F("etaVsVtx_ratio","etaVsVtx_ratio",80,-20,20,120,-6,6);
  20
  21   hEtaVsVtx_meas ->SetXTitle("vtx z [cm]");  hEtaVsVtx_meas ->SetYTitle("#eta");
  22   hEtaVsVtx_gene ->SetXTitle("vtx z [cm]");  hEtaVsVtx_gene ->SetYTitle("#eta");
  23   hEtaVsVtx_corr ->SetXTitle("vtx z [cm]");  hEtaVsVtx_corr ->SetYTitle("#eta");
  24   hEtaVsVtx_ratio->SetXTitle("vtx z [cm]");  hEtaVsVtx_ratio->SetYTitle("#eta");
  25
  26 }
  27
  28 //____________________________________________________________________
  29 void
  30 dNdEtaCorrection::Finish() {
  31
  32   hEtaVsVtx_ratio->Divide(hEtaVsVtx_meas, hEtaVsVtx_gene, 1,1,"B");
  33   hEtaVsVtx_corr->Divide(hEtaVsVtx_gene, hEtaVsVtx_meas, 1,1,"B");
  34
  35   Int_t nBinsVtx = hEtaVsVtx_corr->GetNbinsX();
  36   Int_t nBinsEta = hEtaVsVtx_corr->GetNbinsY();
  37
  38   TH2F* tmp = (TH2F*)hEtaVsVtx_corr->Clone("tmp");
  39
  40   // cut at 0.2
  41   for (Int_t bx=0; bx<=nBinsVtx; bx++) {
  42     for (Int_t by=0; by<=nBinsEta; by++) {
  43       if (tmp->GetBinContent(bx,by)<0.2) {
  44         hEtaVsVtx_corr->SetBinContent(bx,by,0);
  45         hEtaVsVtx_corr->SetBinError(bx,by,0);
  46
  47         tmp->SetBinContent(bx,by,0);
  48       }
  49       else
  50         tmp->SetBinContent(bx,by,1);
  51     }
  52   }
  53
  54 }
  55
  56 //____________________________________________________________________
  57 void
  58 dNdEtaCorrection::RemoveEdges(Float_t cut, Int_t nBinsVtx, Int_t nBinsEta) {
  59
  60   // remove edges of correction histogram by removing
  61   // - bins with content bigger than cut
  62   // - bins next to bins with zero bin content
  63
  64   Int_t nBinsX = hEtaVsVtx_corr->GetNbinsX();
  65   Int_t nBinsY = hEtaVsVtx_corr->GetNbinsY();
  66
  67   // set bin content to zero for bins with content bigger than cut
  68   for (Int_t bx=0; bx<=nBinsX; bx++) {
  69     for (Int_t by=0; by<=nBinsY; by++) {
  70       if (hEtaVsVtx_corr->GetBinContent(bx,by)>cut) {
  71         hEtaVsVtx_corr->SetBinContent(bx,by,0);
  72         hEtaVsVtx_corr->SetBinError(bx,by,0);
  73       }
  74     }
  75   }
  76
  77   // set bin content to zero for bins next to bins with zero
  78   TH2F* tmp = (TH2F*)hEtaVsVtx_corr->Clone("tmp");
  79   tmp->Reset();
  80
  81   Bool_t done = kFALSE;
  82   Int_t nBinsVtxCount = 0;
  83   Int_t nBinsEtaCount = 0;
  84   while (!done) {
  85     if (nBinsVtxCount<nBinsVtx)
  86       for (Int_t bx=0; bx<=nBinsX; bx++) {
  87         for (Int_t by=0; by<=nBinsY; by++) {
  88           if ((hEtaVsVtx_corr->GetBinContent(bx+1,by)==0)||
  89               (hEtaVsVtx_corr->GetBinContent(bx-1,by)==0))
  90             tmp->SetBinContent(bx,by,1);
  91
  92         }
  93       }
  94     if (nBinsEtaCount<nBinsEta)
  95       for (Int_t bx=0; bx<=nBinsX; bx++) {
  96         for (Int_t by=0; by<=nBinsY; by++) {
  97           if ((hEtaVsVtx_corr->GetBinContent(bx,by+1)==0)||
  98               (hEtaVsVtx_corr->GetBinContent(bx,by-1)==0))
  99             tmp->SetBinContent(bx,by,1);
 100         }
 101       }
 102     for (Int_t bx=0; bx<=nBinsX; bx++) {
 103       for (Int_t by=0; by<=nBinsY; by++) {
 104         if (tmp->GetBinContent(bx,by)==1) {
 105           hEtaVsVtx_corr->SetBinContent(bx,by,0);
 106           hEtaVsVtx_corr->SetBinError(bx,by,0);
 107         }
 108       }
 109     }
 110     nBinsVtxCount++;
 111     nBinsEtaCount++;
 112     if ((nBinsVtxCount>=nBinsVtx)&&(nBinsEtaCount>=nBinsEta)) done=kTRUE;
 113   }
 114   tmp->Delete();
 115
 116 }
 117
 118 //____________________________________________________________________
 119 Bool_t
 120 dNdEtaCorrection::LoadHistograms(Char_t* fileName, Char_t* dir) {
 121
 122   TFile* fin = TFile::Open(fileName);
 123
 124   // add test of file
 125   // return kFALSE
 126
 127   hEtaVsVtx_meas  = (TH2F*)fin->Get(Form("%s/etaVsVtx_meas",dir));
 128   hEtaVsVtx_gene  = (TH2F*)fin->Get(Form("%s/etaVsVtx_gene",dir));
 129   hEtaVsVtx_corr  = (TH2F*)fin->Get(Form("%s/etaVsVtx_corr",dir));
 130
 131   hEtaVsVtx_ratio = (TH2F*)fin->Get(Form("%s/etaVsVtx_ratio",dir));
 132
 133   return kTRUE;
 134 }
 135
 136
 137 //____________________________________________________________________
 138 void
 139 dNdEtaCorrection::SaveHistograms() {
 140
 141   gDirectory->mkdir(fName.Data());
 142   gDirectory->cd(fName.Data());
 143
 144   hEtaVsVtx_meas ->Write();
 145   hEtaVsVtx_gene ->Write();
 146
 147   if (hEtaVsVtx_corr)
 148     hEtaVsVtx_corr->Write();
 149
 150   if (hEtaVsVtx_ratio)
 151     hEtaVsVtx_ratio->Write();
 152
 153
 154   gDirectory->cd("../");
 155 }
 156
 157 //____________________________________________________________________
 158 void dNdEtaCorrection::DrawHistograms()
 159 {
 160   TCanvas* canvas = new TCanvas("dNdEtaCorrection", "dNdEtaCorrection", 800, 800);
 161   canvas->Divide(2, 2);
 162
 163   canvas->cd(1);
 164   if (hEtaVsVtx_meas)
 165     hEtaVsVtx_meas->Draw("COLZ");
 166
 167   canvas->cd(2);
 168   if (hEtaVsVtx_gene)
 169     hEtaVsVtx_gene->Draw("COLZ");
 170
 171   canvas->cd(3);
 172   if (hEtaVsVtx_ratio)
 173     hEtaVsVtx_ratio->Draw("COLZ");
 174
 175   canvas->cd(4);
 176   if (hEtaVsVtx_corr)
 177     hEtaVsVtx_corr->Draw("COLZ");
 178 }