o) new way of calculating the final dndeta

o) debugged the dndeta histograms which show only a selection of the vertex range

diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx
index b486000878b749ec03cd9aaff26f6439856b2176..5b511787501ac68939270564584378e140d8345a 100644 (file)
--- a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx
+++ b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.cxx
@@ -8,6 +8,7 @@
 #include <TParticle.h>
 #include <TParticlePDG.h>
 #include <TVector3.h>
+#include <TH1F.h>
 #include <TH3F.h>
 #include <TTree.h>
 
@@ -22,7 +23,9 @@ ClassImp(AlidNdEtaAnalysisMCSelector)
 
 AlidNdEtaAnalysisMCSelector::AlidNdEtaAnalysisMCSelector() :
   AlidNdEtaAnalysisSelector(),
-  fVertex(0)
+  fVertex(0),
+  fPartEta(0),
+  fEvents(0)
 {
   //
   // Constructor. Initialization of pointers
@@ -42,7 +45,9 @@ void AlidNdEtaAnalysisMCSelector::Init(TTree *tree)
 
   tree->SetBranchStatus("ESD", 0);
 
-  fVertex = new TH3F("vertex", "vertex", 50, -50, 50, 50, -50, 50, 50, -50, 50);
+  fVertex = new TH3F("vertex_check", "vertex_check", 50, -50, 50, 50, -50, 50, 50, -50, 50);
+  fPartEta = new TH1F("dndeta_check", "dndeta_check", 120, -6, 6);
+  fPartEta->Sumw2();
 }
 
 Bool_t AlidNdEtaAnalysisMCSelector::Process(Long64_t entry)
@@ -102,9 +107,13 @@ Bool_t AlidNdEtaAnalysisMCSelector::Process(Long64_t entry)
 
     fdNdEtaAnalysis->FillTrack(vtxMC[2], particle->Eta(), 1);
     fVertex->Fill(particle->Vx(), particle->Vy(), particle->Vz());
+
+    fPartEta->Fill(particle->Eta());
   }
   fdNdEtaAnalysis->FillEvent(vtxMC[2]);
 
+  ++fEvents;
+
   return kTRUE;
 }
 
@@ -112,6 +121,15 @@ void AlidNdEtaAnalysisMCSelector::Terminate()
 {
   AlidNdEtaAnalysisSelector::Terminate();
 
-  new TCanvas;
+  fPartEta->Scale(1.0/fEvents);
+  fPartEta->Scale(1.0/fPartEta->GetBinWidth(1));
+
+  TCanvas* canvas = new TCanvas("control", "control", 900, 450);
+  canvas->Divide(2, 1);
+
+  canvas->cd(1);
   fVertex->Draw();
+
+  canvas->cd(2);
+  fPartEta->Draw();
 }

diff --git a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h
index 4f54eb8d77add1547ee502b2a7dbfe5bff1c1588..3aae19ec7917b787bfde779ba611e8cd2ffc014a 100644 (file)
--- a/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h
+++ b/PWG0/dNdEta/AlidNdEtaAnalysisMCSelector.h
@@ -6,6 +6,7 @@
 #include "AlidNdEtaAnalysisSelector.h"
 
 class TH3F;
+class TH1F;
 
 class AlidNdEtaAnalysisMCSelector : public AlidNdEtaAnalysisSelector {
   public:
@@ -19,9 +20,11 @@ class AlidNdEtaAnalysisMCSelector : public AlidNdEtaAnalysisSelector {
  protected:
 
  private:
-  TH3F* fVertex;
+    TH3F* fVertex;  //! vertex of counted particles
+    TH1F* fPartEta; //! counted particles as function of eta
+    Int_t fEvents;  //! number of processed events
 
-  ClassDef(AlidNdEtaAnalysisMCSelector, 0);
+    ClassDef(AlidNdEtaAnalysisMCSelector, 0);
 };
 
 #endif

diff --git a/PWG0/dNdEta/dNdEtaAnalysis.cxx b/PWG0/dNdEta/dNdEtaAnalysis.cxx
index ecb6f7bf67d97a949731072eeedee875978c25b8..38d878de605e25f485c8f786ac91817af5401641 100644 (file)
--- a/PWG0/dNdEta/dNdEtaAnalysis.cxx
+++ b/PWG0/dNdEta/dNdEtaAnalysis.cxx
@@ -83,7 +83,7 @@ void dNdEtaAnalysis::Finish(dNdEtaCorrection* correction)
   }
 
   // normalize with n events (per vtx)
-  for (Int_t iVtx=0; iVtx<=hVtx->GetNbinsX(); iVtx++) {
+/*  for (Int_t iVtx=0; iVtx<=hVtx->GetNbinsX(); iVtx++) {
     Float_t nEvents      = hVtx->GetBinContent(iVtx);
     Float_t nEventsError = hVtx->GetBinError(iVtx);
 
@@ -111,60 +111,52 @@ void dNdEtaAnalysis::Finish(dNdEtaCorrection* correction)
       hEtaVsVtxCheck->SetBinContent(iVtx, iEta, value);
       hEtaVsVtxCheck->SetBinError(iVtx,   iEta, error);
     }
-  }
+  }*/
 
   // then take the wieghted average for each eta
   // is this the right way to do it???
-  for (Int_t iEta=0; iEta<=hEtaVsVtx->GetNbinsY(); iEta++) {
-    Float_t sum           = 0;
-    Float_t sumError2     = 0;
-    Int_t   nMeasurements = 0;    
-
-    Float_t sumXw = 0;
-    Float_t sumW  = 0;
-    
+  for (Int_t iEta=0; iEta<=hEtaVsVtx->GetNbinsY(); iEta++)
+  {
     // do we have several histograms for different vertex positions?
     Int_t vertexBinWidth = hVtx->GetNbinsX() / (kVertexBinning-1);
     for (Int_t vertexPos=0; vertexPos<kVertexBinning; ++vertexPos)
     {
-      Int_t vertexBinBegin = 0;
-      Int_t vertexBinEnd = hVtx->GetNbinsX();
+      Int_t vertexBinBegin = 1;
+      Int_t vertexBinEnd = hVtx->GetNbinsX() + 1;
 
       // the first histogram is always for the whole vertex range
       if (vertexPos > 0)
       {
-        vertexBinBegin = vertexBinWidth * (vertexPos-1);
+        vertexBinBegin = 1 + vertexBinWidth * (vertexPos-1);
         vertexBinEnd = vertexBinBegin + vertexBinWidth;
       }
 
-      for (Int_t iVtx=vertexBinBegin; iVtx<=vertexBinEnd; iVtx++) {
-        if (hVtx->GetBinContent(iVtx)==0)             continue;
-        if (hEtaVsVtxCheck->GetBinContent(iVtx, iEta)==0) continue;
-
-        Float_t w = 1/TMath::Power(hEtaVsVtx->GetBinError(iVtx, iEta),2);
-        sumXw = sumXw + hEtaVsVtxCheck->GetBinContent(iVtx, iEta)*w;
-        sumW  = sumW + w;
-
-        sum = sum + hEtaVsVtxCheck->GetBinContent(iVtx, iEta);
-        sumError2 = sumError2 + TMath::Power(hEtaVsVtxCheck->GetBinError(iVtx, iEta),2);      
-        nMeasurements++;
+      Float_t totalEvents = hVtx->Integral(vertexBinBegin, vertexBinEnd - 1);
+      if (totalEvents == 0)
+      {
+        printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
+        continue;
       }
-      Float_t dndeta = 0;
-      Float_t error  = 0;
 
-      if (nMeasurements!=0) {
-        dndeta = sum/Float_t(nMeasurements);
-        error  = TMath::Sqrt(sumError2)/Float_t(nMeasurements);
+      Float_t sum = 0;
+      Float_t sumError2 = 0;
+      for (Int_t iVtx = vertexBinBegin; iVtx < vertexBinEnd; iVtx++)
+      {
+        if (hEtaVsVtx->GetBinContent(iVtx, iEta) != 0)
+        {
+          sum = sum + hEtaVsVtx->GetBinContent(iVtx, iEta);
+          sumError2 = sumError2 + TMath::Power(hEtaVsVtx->GetBinError(iVtx, iEta),2);
+        }
+      }
 
-        dndeta = sumXw/sumW;
-        error  = 1/TMath::Sqrt(sumW);
+      Float_t dndeta = sum / totalEvents;
+      Float_t error  = TMath::Sqrt(sumError2) / totalEvents;
 
-        dndeta = dndeta/hdNdEta[vertexPos]->GetBinWidth(iEta);
-        error  = error/hdNdEta[vertexPos]->GetBinWidth(iEta);
+      dndeta = dndeta/hdNdEta[vertexPos]->GetBinWidth(iEta);
+      error  = error/hdNdEta[vertexPos]->GetBinWidth(iEta);
 
-        hdNdEta[vertexPos]->SetBinContent(iEta, dndeta);
-        hdNdEta[vertexPos]->SetBinError(iEta, error);
-      }
+      hdNdEta[vertexPos]->SetBinContent(iEta, dndeta);
+      hdNdEta[vertexPos]->SetBinError(iEta, error);
     }
   }
 }
@@ -227,15 +219,15 @@ void dNdEtaAnalysis::DrawHistograms()
 
     TLegend* legend = new TLegend(0.7, 0.7, 0.9, 0.9);
 
-    for (Int_t i=1; i<kVertexBinning; ++i)
+    for (Int_t i=0; i<kVertexBinning; ++i)
     {
       //canvas2->cd(i-1);
       //printf("%d\n", i);
       if (hdNdEta[i])
       {
-        hdNdEta[i]->SetLineColor(i);
-        hdNdEta[i]->Draw((i == 1) ? "" : "SAME");
-        legend->AddEntry(hdNdEta[i], Form("Vtx Bin %d", i-1));
+        hdNdEta[i]->SetLineColor(i+1);
+        hdNdEta[i]->Draw((i == 0) ? "" : "SAME");
+        legend->AddEntry(hdNdEta[i], (i == 0) ? "Vtx All" : Form("Vtx Bin %d", i-1));
       }
     }
 

diff --git a/PWG0/dNdEta/dNdEtaAnalysis.h b/PWG0/dNdEta/dNdEtaAnalysis.h
index 754f7978206fcddd4bc1ed1d1e40b78e46751e5e..ef90a1358fa9a5bbdac909258c1b94ab6d9a25d4 100644 (file)
--- a/PWG0/dNdEta/dNdEtaAnalysis.h
+++ b/PWG0/dNdEta/dNdEtaAnalysis.h
@@ -28,7 +28,7 @@ class dNdEtaCorrection;
 class dNdEtaAnalysis : public TNamed
 {
 public:
-  enum { kVertexBinning = 1+6 }; // the first is for the whole vertex range, the others divide the vertex range
+  enum { kVertexBinning = 1+4 }; // the first is for the whole vertex range, the others divide the vertex range
 
   dNdEtaAnalysis(Char_t* name, Char_t* title);
 

diff --git a/PWG0/dNdEta/drawVertexRecEff.C b/PWG0/dNdEta/drawVertexRecEff.C
new file mode 100644 (file)
index 0000000..0b5ae24
--- /dev/null
+++ b/PWG0/dNdEta/drawVertexRecEff.C
@@ -0,0 +1,57 @@
+/* $Id$ */
+
+// draws the result of AlidNdEtaVertexRecEffSelector
+
+void drawVertexRecEff()
+{
+  TFile* fout = TFile::Open("vertexRecEff.root");
+
+  TH1F* dNGen = dynamic_cast<TH1F*> fout->Get("dNGen");
+  TH1F* dNRec = dynamic_cast<TH1F*> fout->Get("dNRec");
+
+  TH1F* vtxGen = dynamic_cast<TH1F*> fout->Get("VtxGen");
+  TH1F* vtxRec = dynamic_cast<TH1F*> fout->Get("VtxRec");
+
+  // calculate ratios
+  TH1F* dNRatiodN = dynamic_cast<TH1F*> (dNRec->Clone("dNRatiodN"));
+  dNRatiodN->SetTitle("Ratio");
+  dNRatiodN->Divide(dNGen);
+
+  //vtxGen->Rebin(4);
+  //vtxRec->Rebin(4);
+
+  // calculate correction ratio number
+  Float_t sumGen = 0;
+  Float_t sumRec = 0;
+  for (Int_t i=1; i<=dNGen->GetNbinsX(); ++i)
+  {
+    sumGen += dNGen->GetBinCenter(i) * dNGen->GetBinContent(i);
+    sumRec += dNRec->GetBinCenter(i) * dNRec->GetBinContent(i);
+  }
+  Float_t ratio = sumRec / dNRec->Integral(1, dNGen->GetNbinsX()) * dNGen->Integral(1, dNGen->GetNbinsX()) / sumGen;
+
+  cout << "Ratio: " << ratio << endl;
+
+  TH1F* dNRatioVtx = dynamic_cast<TH1F*> (vtxRec->Clone("dNRatioVtx"));
+  dNRatioVtx->SetTitle("Ratio");
+  dNRatioVtx->Divide(vtxGen);
+
+  TCanvas* canvas = new TCanvas("dN", "dN", 1000, 1000);
+  canvas->Divide(2, 2);
+
+  canvas->cd(1);
+  dNGen->Draw();
+  dNRec->SetLineColor(kRed);
+  dNRec->Draw("SAME");
+
+  canvas->cd(2);
+  dNRatiodN->Draw();
+
+  canvas->cd(3);
+  vtxGen->Draw();
+  vtxRec->SetLineColor(kRed);
+  vtxRec->Draw("SAME");
+
+  canvas->cd(4);
+  dNRatioVtx->Draw();
+}