implemented analysis using 3d information (vtx_z, eta, pt)
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / AlidNdEtaCorrection.cxx
1 /* $Id$ */
2
3 #include "AlidNdEtaCorrection.h"
4
5 #include <TCanvas.h>
6 #include <TH3F.h>
7 #include <TH1D.h>
8
9 //____________________________________________________________________
10 ClassImp(AlidNdEtaCorrection)
11
12 //____________________________________________________________________
13 AlidNdEtaCorrection::AlidNdEtaCorrection(Char_t* name) 
14   : TNamed(name, name),
15   fNEvents(0),
16   fNTriggeredEvents(0)
17 {  
18   // constructor
19   //
20
21   fTrack2ParticleCorrection = new AliCorrectionMatrix3D("nTrackToNPart", "nTrackToNPart",80,-20,20,120,-6,6, 100, 0, 10);
22
23   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, 
24                             10.5, 12.5, 14.5, 16.5, 18.5, 20.5, 25.5, 30.5, 40.5, 50.5, 100.5, 300.5};
25   Float_t binLimitsVtx[] = {-20,-15,-10,-6,-3,0,3,6,10,15,20};
26   
27   fVertexRecoCorrection        = new AliCorrectionMatrix2D("vtxReco",       "vtxReco",10,binLimitsVtx ,22,binLimitsN);
28
29   fTriggerBiasCorrection       = new AliCorrectionMatrix2D("triggerBias",   "triggerBias",120,-6,6,100, 0, 10);
30
31   fTrack2ParticleCorrection ->SetAxisTitles("vtx z [cm]", "#eta", "p_{T}");
32   fVertexRecoCorrection        ->SetAxisTitles("vtx z [cm]", "n particles/tracks/tracklets?");
33
34   fTriggerBiasCorrection       ->SetAxisTitles("#eta", "p_{T} [GeV/c]");
35 }
36
37 //____________________________________________________________________
38 void
39 AlidNdEtaCorrection::Finish() {
40   //
41   // finish method
42   //
43   // divide the histograms in the AliCorrectionMatrix2D objects to get the corrections
44
45
46   fTrack2ParticleCorrection->Divide();
47
48   fVertexRecoCorrection->Divide();
49
50   fTriggerBiasCorrection->GetMeasuredHistogram()->Scale(Double_t(fNTriggeredEvents)/Double_t(fNEvents));
51   fTriggerBiasCorrection->Divide();
52
53 }
54
55 //____________________________________________________________________
56 Long64_t 
57 AlidNdEtaCorrection::Merge(TCollection* list) {
58   // Merge a list of dNdEtaCorrection objects with this (needed for
59   // PROOF). 
60   // Returns the number of merged objects (including this).
61
62   if (!list)
63     return 0;
64   
65   if (list->IsEmpty())
66     return 1;
67
68   TIterator* iter = list->MakeIterator();
69   TObject* obj;
70
71   // collections of measured and generated histograms
72   TList* collectionNtrackToNparticle = new TList;
73   TList* collectionVertexReco        = new TList;
74   TList* collectionTriggerBias       = new TList;
75
76   Int_t count = 0;
77   while ((obj = iter->Next())) {
78     
79     AlidNdEtaCorrection* entry = dynamic_cast<AlidNdEtaCorrection*> (obj);
80     if (entry == 0) 
81       continue;
82
83     collectionNtrackToNparticle ->Add(entry->GetTrack2ParticleCorrection());
84     collectionVertexReco        ->Add(entry->GetVertexRecoCorrection());
85     collectionTriggerBias        ->Add(entry->GetTriggerBiasCorrection());
86
87     count++;
88   }
89   fTrack2ParticleCorrection ->Merge(collectionNtrackToNparticle);
90   fVertexRecoCorrection        ->Merge(collectionVertexReco);
91   fTriggerBiasCorrection        ->Merge(collectionTriggerBias);
92   
93   delete collectionNtrackToNparticle;
94   delete collectionVertexReco;
95   delete collectionTriggerBias;
96
97   return count+1;
98 }
99
100
101 //____________________________________________________________________
102 Bool_t
103 AlidNdEtaCorrection::LoadHistograms(Char_t* fileName, Char_t* dir) {
104   //
105   // loads the histograms
106   //
107
108   fTrack2ParticleCorrection ->LoadHistograms(fileName, dir);
109   fVertexRecoCorrection        ->LoadHistograms(fileName, dir);
110   fTriggerBiasCorrection       ->LoadHistograms(fileName, dir);
111   
112   return kTRUE;
113 }
114
115
116 //____________________________________________________________________
117 void
118 AlidNdEtaCorrection::SaveHistograms() {
119   //
120   // save the histograms
121   //
122
123   gDirectory->mkdir(fName.Data());
124   gDirectory->cd(fName.Data());
125
126   fTrack2ParticleCorrection ->SaveHistograms();
127   fVertexRecoCorrection        ->SaveHistograms();
128   fTriggerBiasCorrection       ->SaveHistograms();
129
130   gDirectory->cd("../");
131 }
132
133 //____________________________________________________________________
134 void AlidNdEtaCorrection::DrawHistograms()
135 {
136   //
137   // call the draw histogram method of the two AliCorrectionMatrix2D objects
138
139   fTrack2ParticleCorrection ->DrawHistograms();
140   fVertexRecoCorrection        ->DrawHistograms();
141   fTriggerBiasCorrection       ->DrawHistograms();
142
143 }
144
145 //____________________________________________________________________
146 Float_t AlidNdEtaCorrection::GetMeasuredFraction(Float_t ptCutOff, Float_t eta, Bool_t debug)
147 {
148   // calculates the fraction of particles measured (some are missed due to the pt cut off)
149   // uses the generated particle histogram from fTrack2ParticleCorrection
150
151   TH3F* generated = fTrack2ParticleCorrection->GetGeneratedHistogram();
152
153   // find eta borders, if eta is negative assume -0.8 ... 0.8
154   Int_t etaBegin = 0;
155   Int_t etaEnd = 0;
156   if (eta < 0)
157   {
158     etaBegin = generated->GetYaxis()->FindBin(-0.8);
159     etaEnd = generated->GetYaxis()->FindBin(0.8);
160   }
161   else
162   {
163     etaBegin = generated->GetYaxis()->FindBin(eta);
164     etaEnd = etaBegin;
165   }
166
167   Int_t vertexBegin = generated->GetXaxis()->FindBin(-10);
168   Int_t vertexEnd = generated->GetXaxis()->FindBin(10);
169
170   TH1D* ptProj = dynamic_cast<TH1D*> (generated->ProjectionZ(Form("%s_pt", GetName()), vertexBegin, vertexEnd, etaBegin, etaEnd));
171
172   Int_t ptBin = ptProj->FindBin(ptCutOff);
173   Float_t abovePtCut = ptProj->Integral(ptBin, ptProj->GetNbinsX());
174   Float_t all = ptProj->Integral();
175
176   if (all == 0)
177     return -1;
178
179   Float_t fraction = abovePtCut / all;
180
181   if (debug)
182   {
183     new TCanvas;
184     ptProj->Draw();
185   }
186
187   return fraction;
188 }
189