Added AlidNdEtaCorrection (new procedure).
authorekman <ekman@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 13 Jun 2006 07:28:47 +0000 (07:28 +0000)
committerekman <ekman@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 13 Jun 2006 07:28:47 +0000 (07:28 +0000)
PWG0/dNdEta/AlidNdEtaCorrection.cxx [new file with mode: 0644]
PWG0/dNdEta/AlidNdEtaCorrection.h [new file with mode: 0644]
PWG0/dNdEta/makeCorrection.C

diff --git a/PWG0/dNdEta/AlidNdEtaCorrection.cxx b/PWG0/dNdEta/AlidNdEtaCorrection.cxx
new file mode 100644 (file)
index 0000000..6ee954a
--- /dev/null
@@ -0,0 +1,140 @@
+/* $Id$ */
+
+#include "AlidNdEtaCorrection.h"
+
+#include <TCanvas.h>
+
+//____________________________________________________________________
+ClassImp(AlidNdEtaCorrection)
+
+//____________________________________________________________________
+AlidNdEtaCorrection::AlidNdEtaCorrection(Char_t* name) 
+  : TNamed(name, name)
+{  
+  // constructor
+  //
+
+  fNtrackToNparticleCorrection = new CorrectionMatrix2D("nTrackToNPart", "nTrackToNPart",80,-20,20,120,-6,6);
+
+  Float_t binLimitsN[]   = {-0.5, 0.5, 1.5, 2.5, 3.5, 4.5, 5.5, 6.5, 7.5, 8.5, 9.5, 
+                           10.5, 12.5, 14.5, 16.5, 18.5, 20.5, 25.5, 30.5, 40.5, 50.5, 100.5, 300.5};
+  Float_t binLimitsVtx[] = {-20,-15,-10,-6,-3,0,3,6,10,15,20};
+  
+  fVertexRecoCorrection        = new CorrectionMatrix2D("vtxReco",       "vtxReco",10,binLimitsVtx ,22,binLimitsN);
+
+  fTriggerBiasCorrection       = new CorrectionMatrix2D("triggerBias",   "triggerBias",120,-6,6,100, 0, 10);
+
+  fNtrackToNparticleCorrection ->SetAxisTitles("vtx z [cm]", "#eta");
+  fVertexRecoCorrection        ->SetAxisTitles("vtx z [cm]", "n particles/tracks/tracklets?");
+  
+  fTriggerBiasCorrection       ->SetAxisTitles("#eta", "p_{T} [GeV/c]");
+}
+
+//____________________________________________________________________
+void
+AlidNdEtaCorrection::Finish(Int_t nEventsAll, Int_t nEventsTriggered) {  
+  //
+  // finish method
+  //
+  // divide the histograms in the CorrectionMatrix2D objects to get the corrections
+
+  
+  fNtrackToNparticleCorrection->Divide();
+
+  fVertexRecoCorrection->Divide();
+
+  fTriggerBiasCorrection->GetMeasuredHistogram()->Scale(Double_t(nEventsTriggered)/Double_t(nEventsAll));
+  fTriggerBiasCorrection->Divide();
+
+}
+
+//____________________________________________________________________
+Long64_t 
+AlidNdEtaCorrection::Merge(TCollection* list) {
+  // Merge a list of dNdEtaCorrection objects with this (needed for
+  // PROOF). 
+  // Returns the number of merged objects (including this).
+
+  if (!list)
+    return 0;
+  
+  if (list->IsEmpty())
+    return 1;
+
+  TIterator* iter = list->MakeIterator();
+  TObject* obj;
+
+  // collections of measured and generated histograms
+  TList* collectionNtrackToNparticle = new TList;
+  TList* collectionVertexReco        = new TList;
+  TList* collectionTriggerBias       = new TList;
+
+  Int_t count = 0;
+  while ((obj = iter->Next())) {
+    
+    AlidNdEtaCorrection* entry = dynamic_cast<AlidNdEtaCorrection*> (obj);
+    if (entry == 0) 
+      continue;
+
+    collectionNtrackToNparticle ->Add(entry->GetNtrackToNpraticleCorrection());
+    collectionVertexReco        ->Add(entry->GetVertexRecoCorrection());
+    collectionTriggerBias        ->Add(entry->GetTriggerBiasCorrection());
+
+    count++;
+  }
+  fNtrackToNparticleCorrection ->Merge(collectionNtrackToNparticle);
+  fVertexRecoCorrection        ->Merge(collectionVertexReco);
+  fTriggerBiasCorrection        ->Merge(collectionTriggerBias);
+  
+  delete collectionNtrackToNparticle;
+  delete collectionVertexReco;
+  delete collectionTriggerBias;
+
+  return count+1;
+}
+
+
+//____________________________________________________________________
+Bool_t
+AlidNdEtaCorrection::LoadHistograms(Char_t* fileName, Char_t* dir) {
+  //
+  // loads the histograms
+  //
+
+  fNtrackToNparticleCorrection ->LoadHistograms(fileName, dir);
+  fVertexRecoCorrection        ->LoadHistograms(fileName, dir);
+  fTriggerBiasCorrection       ->LoadHistograms(fileName, dir);
+  
+  return kTRUE;
+}
+
+
+//____________________________________________________________________
+void
+AlidNdEtaCorrection::SaveHistograms() {
+  //
+  // save the histograms
+  //
+
+  gDirectory->mkdir(fName.Data());
+  gDirectory->cd(fName.Data());
+
+  fNtrackToNparticleCorrection ->SaveHistograms();
+  fVertexRecoCorrection        ->SaveHistograms();
+  fTriggerBiasCorrection       ->SaveHistograms();
+
+  gDirectory->cd("../");
+}
+
+//____________________________________________________________________
+void AlidNdEtaCorrection::DrawHistograms()
+{
+  //
+  // call the draw histogram method of the two CorrectionMatrix2D objects
+
+  fNtrackToNparticleCorrection ->DrawHistograms();
+  fVertexRecoCorrection        ->DrawHistograms();
+  fTriggerBiasCorrection       ->DrawHistograms();
+
+  
+}
diff --git a/PWG0/dNdEta/AlidNdEtaCorrection.h b/PWG0/dNdEta/AlidNdEtaCorrection.h
new file mode 100644 (file)
index 0000000..9958d47
--- /dev/null
@@ -0,0 +1,75 @@
+/* $Id$ */
+
+#ifndef ALIDNDETACORRECTION_H
+#define ALIDNDETACORRECTION_H
+
+
+// ------------------------------------------------------
+//
+// Class to handle corrections for dN/dEta measurements
+//
+// ------------------------------------------------------
+//
+// TODO:
+// - make the ntrack to npart correction 3D
+// - add documentation
+// - add status: generate or use maps
+// - add functionality to set the bin sizes
+// 
+
+#include <TNamed.h>
+#include <TFile.h>
+
+#include <CorrectionMatrix2D.h>
+
+
+class AlidNdEtaCorrection : public TNamed
+{
+protected:  
+  
+  CorrectionMatrix2D* fNtrackToNparticleCorrection; // handles the track-to-vertex correction
+  CorrectionMatrix2D* fVertexRecoCorrection;        // handles the vertex reco (n tracks vs vtx)
+
+  CorrectionMatrix2D* fTriggerBiasCorrection;          // MB to desired sample
+
+public:
+  AlidNdEtaCorrection(Char_t* name="dndeta_correction");
+
+  void FillEvent(Float_t vtx, Float_t n)                        {fVertexRecoCorrection->FillGene(vtx, n);}
+  void FillEventWithReconstructedVertex(Float_t vtx, Float_t n) {fVertexRecoCorrection->FillMeas(vtx, n);}
+  
+  void FillParticle(Float_t vtx, Float_t eta, Float_t pt=0)                  {fNtrackToNparticleCorrection->FillGene(vtx, eta);}
+  void FillParticleWhenMeasuredTrack(Float_t vtx, Float_t eta, Float_t pt=0) {fNtrackToNparticleCorrection->FillMeas(vtx, eta);}
+  
+  void FillParticleAllEvents(Float_t eta, Float_t pt=0)          {fTriggerBiasCorrection->FillGene(eta, pt);}
+  void FillParticleWhenEventTriggered(Float_t eta, Float_t pt=0) {fTriggerBiasCorrection->FillMeas(eta, pt);}
+
+  void Finish(Int_t nEventsAll = 1, Int_t nEventsTriggered = 1);
+
+  CorrectionMatrix2D* GetNtrackToNpraticleCorrection() {return fNtrackToNparticleCorrection;}
+  CorrectionMatrix2D* GetVertexRecoCorrection()        {return fVertexRecoCorrection;}
+  CorrectionMatrix2D* GetTriggerBiasCorrection()       {return fTriggerBiasCorrection;}
+
+  virtual Long64_t Merge(TCollection* list);
+
+  void    SaveHistograms();
+  Bool_t  LoadHistograms(Char_t* fileName, Char_t* dir = "dndeta_correction");
+  Bool_t  LoadCorrection(Char_t* fileName, Char_t* dir = "dndeta_correction") 
+    {return LoadHistograms(fileName, dir);}
+  
+  void DrawHistograms();
+  
+  //  void RemoveEdges(Float_t cut=2, Int_t nBinsVtx=0, Int_t nBinsEta=0);
+  
+  Float_t GetNtracksToNpartCorrection(Float_t vtx, Float_t eta, Float_t pt) 
+    {return fNtrackToNparticleCorrection->GetCorrection(vtx, eta);}  
+  
+  Float_t GetVertexRecoCorrection(Float_t vtx, Float_t n) {return fVertexRecoCorrection->GetCorrection(vtx, n);}
+
+  Float_t GetTriggerBiasCorrection(Float_t eta, Float_t pt=0) {return fTriggerBiasCorrection->GetCorrection(eta, pt);}
+  
+
+  ClassDef(AlidNdEtaCorrection,0)
+};
+
+#endif
index c8645f3d971cd1dd976fc9ba2fa38692a731b574..b1c041c3701115e47edced4f44f8ba14fa0878f0 100644 (file)
@@ -27,7 +27,7 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
   // ########################################################
   // definition of dNdEta correction object
 
-  dNdEtaCorrection* dNdEtaMap = new dNdEtaCorrection();
+  AlidNdEtaCorrection* dNdEtaMap = new AlidNdEtaCorrection();
 
   // ########################################################
   // get the data dir  
@@ -111,12 +111,15 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
     // getting number of events
     Int_t nEvents = (Int_t)runLoader->GetNumberOfEvents();
     Int_t nESDEvents = esdBranch->GetEntries();
+
+    Int_t nEventsTriggered = 0;
+    Int_t nEventsAll       = 0;
     
     if (nEvents!=nESDEvents) {
       cout << " Different number of events from runloader and esdtree!!!" << nEvents << " / " << nESDEvents << endl;
       return;
     }
-    // ########################################################
+   // ########################################################
     // loop over number of events
     cout << " looping over events..." << endl;
     for(Int_t i=0; i<nEvents; i++) {
@@ -136,15 +139,34 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
       vtx_res[1] = vtxESD->GetYRes();
       vtx_res[2] = vtxESD->GetZRes();
 
-      Bool_t goodEvent = kTRUE;
+      Bool_t vertexReconstructed = kTRUE;
 
       // the vertex should be reconstructed
       if (strcmp(vtxESD->GetName(),"default")==0) 
-       goodEvent = kFALSE;
+       vertexReconstructed = kFALSE;
 
       // the resolution should be reasonable???
       if (vtx_res[2]==0 || vtx_res[2]>0.01)
-       goodEvent = kFALSE;
+       vertexReconstructed = kFALSE;
+
+      // ########################################################
+      // get the trigger info
+      
+      Bool_t triggered = kFALSE;
+
+      // MB should be 
+      // ITS_SPD_GFO_L0  : 32
+      // VZERO_OR_LEFT   : 1
+      // VZERO_OR_RIGHT  : 2
+
+      ULong64_t triggerMask = esd->GetTriggerMask();
+
+      nEventsAll++;      
+
+      if (triggerMask&32 && ((triggerMask&1) || (triggerMask&2))) {
+       triggered = kTRUE;       
+       nEventsTriggered++;
+      }
 
       // ########################################################
       // get the MC vertex
@@ -184,21 +206,25 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
          continue;     
 
        Float_t eta = part->Eta();
+       Float_t pt  = part->Pt();
 
        if (prim) {
-         dNdEtaMap->FillParticleAllEvents(vtx[2], eta);        
+
+         dNdEtaMap->FillParticleAllEvents(eta, pt);              
+
+         if (triggered)
+           dNdEtaMap->FillParticleWhenEventTriggered(eta, pt);
          
-         if (goodEvent)
-           dNdEtaMap->FillParticleWhenGoodEvent(vtx[2], eta);  
+         if (vertexReconstructed)
+           dNdEtaMap->FillParticle(vtx[2], eta, 1.);   
        }
        
       }// end of mc particle
-
-      if (!goodEvent)
-       continue;
-
+      
       // ########################################################
       // loop over esd tracks      
+      Int_t nGoodTracks = 0;
+
       Int_t nTracks = esd->GetNumberOfTracks();            
       for (Int_t t=0; t<nTracks; t++) {
        AliESDtrack* esdTrack = esd->GetTrack(t);      
@@ -207,15 +233,16 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
        if (!esdTrackCuts->AcceptTrack(esdTrack))
          continue;
        
-       AliTPCtrack* tpcTrack = new AliTPCtrack(*esdTrack);     
-       if (tpcTrack->GetAlpha()==0) {
-         cout << " Warning esd track alpha = 0" << endl;
-         continue;
-       }
-
-       Float_t theta = tpcTrack->Theta();
-       Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
-
+       nGoodTracks++;  
+       
+       Double_t p[3];
+       esdTrack->GetPxPyPz(p);
+       Float_t momentum = TMath::Sqrt(TMath::Power(p[0],2) + TMath::Power(p[1],2) + TMath::Power(p[2],2));
+       
+       Float_t eta = -100.;
+       if((momentum != TMath::Abs(p[2]))&&(momentum != 0))
+         eta = 0.5*TMath::Log((momentum + p[2])/(momentum - p[2]));
+       
        // using the eta of the mc particle
        Int_t label = TMath::Abs(esdTrack->GetLabel());
        if (label<0) {
@@ -224,14 +251,22 @@ makeCorrection(Char_t* dataDir, Int_t nRuns=20) {
        }
        TParticle* mcPart = particleStack->Particle(label);     
        eta = mcPart->Eta();
-
-       dNdEtaMap->FillParticleWhenMeasuredTrack(vtx[2], eta);  
-
+       Float_t pt = mcPart->Pt();
+       
+       if (vertexReconstructed)
+         dNdEtaMap->FillParticleWhenMeasuredTrack(vtx[2], eta, pt);    
+       
       } // end of track loop
+      
+      dNdEtaMap->FillEvent(vtxMC[2], nGoodTracks);
+      
+      if (vertexReconstructed)
+       dNdEtaMap->FillEventWithReconstructedVertex(vtxMC[2], nGoodTracks);
+
     } // end  of event loop
   } // end of run loop
-
-  dNdEtaMap->Finish();  
+  
+  dNdEtaMap->Finish(nEventsAll, nEventsTriggered);  
 
   TFile* fout = new TFile("correction_map.root","RECREATE");