added correction for events with vertex but 0 tracks
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / dNdEtaAnalysis.cxx
CommitLineData
dc740de4 1/* $Id$ */
2
75ec0f41 3#include "dNdEtaAnalysis.h"
4
ceb5d1b5 5#include <TFile.h>
45e97e28 6#include <TH3F.h>
74fd10b3 7#include <TH2F.h>
8#include <TH1F.h>
ceb5d1b5 9#include <TMath.h>
10#include <TCanvas.h>
7029240a 11#include <TCollection.h>
12#include <TIterator.h>
13#include <TList.h>
fcf2fb36 14#include <TLegend.h>
74fd10b3 15#include <TLine.h>
fcf2fb36 16
45e97e28 17#include "AlidNdEtaCorrection.h"
74fd10b3 18#include <AliCorrection.h>
19#include <AliPWG0Helper.h>
20#include <AliCorrectionMatrix2D.h>
21#include <AliCorrectionMatrix3D.h>
ceb5d1b5 22
75ec0f41 23//____________________________________________________________________
b7f4a1fd 24ClassImp(dNdEtaAnalysis)
75ec0f41 25
26//____________________________________________________________________
1afae8ff 27dNdEtaAnalysis::dNdEtaAnalysis() :
28 TNamed(),
29 fData(0),
1afae8ff 30 fPtDist(0)
31{
32 // default constructor
33
34 for (Int_t i=0; i<kVertexBinning; ++i)
35 {
3dfa46a4 36 fdNdEtaNotEventCorrected[i] = 0;
1afae8ff 37 fdNdEta[i] = 0;
38 fdNdEtaPtCutOffCorrected[i] = 0;
39 }
40}
41
42//____________________________________________________________________
770a1f1d 43dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::AnalysisMode analysisMode) :
16e24ca3 44 TNamed(name, title),
45e97e28 45 fData(0),
1afae8ff 46 fPtDist(0)
7029240a 47{
6bf0714d 48 // constructor
4dd2ad81 49
0f67a57c 50 // TODO this binning has to be the same than in AliCorrection, somehow passed?!
74fd10b3 51 Float_t binLimitsPt[] = {0.0, 0.05, 0.1, 0.125, 0.15, 0.175, 0.2, 0.225, 0.25, 0.275, 0.3, 0.325, 0.35, 0.375, 0.4, 0.425, 0.45, 0.475, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.5, 2.0, 5.0, 10.0, 100.0};
6bf0714d 52
770a1f1d 53 fData = new AliCorrection("Analysis", Form("%s Analysis", title), analysisMode);
45e97e28 54
74fd10b3 55 // do not add this hists to the directory
56 Bool_t oldStatus = TH1::AddDirectoryStatus();
57 TH1::AddDirectory(kFALSE);
1afae8ff 58
0448e811 59 fdNdEta[0] = new TH1F("dNdEta", "dN_{ch}/d#eta;#eta;dN_{ch}/d#eta", 40, -2, 2);
74fd10b3 60
3dfa46a4 61 fdNdEtaNotEventCorrected[0] = dynamic_cast<TH1F*> (fdNdEta[0]->Clone(Form("%s_noteventcorrected", fdNdEta[0]->GetName())));
74fd10b3 62 fdNdEtaPtCutOffCorrected[0] = dynamic_cast<TH1F*> (fdNdEta[0]->Clone(Form("%s_corrected", fdNdEta[0]->GetName())));
1afae8ff 63
45e97e28 64 for (Int_t i=1; i<kVertexBinning; ++i)
7029240a 65 {
74fd10b3 66 fdNdEta[i] = dynamic_cast<TH1F*> (fdNdEta[0]->Clone(Form("%s_%d", fdNdEta[0]->GetName(), i)));
3dfa46a4 67 fdNdEtaNotEventCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaNotEventCorrected[0]->Clone(Form("%s_%d", fdNdEta[0]->GetName(), i)));
74fd10b3 68 fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaPtCutOffCorrected[0]->Clone(Form("%s_%d", fdNdEtaPtCutOffCorrected[0]->GetName(), i)));
7029240a 69 }
75ec0f41 70
74fd10b3 71 fPtDist = new TH1F("Pt", "p_{T} distribution;p_{T} [GeV/c];#frac{dN}{d#eta dp_{T}} [c/GeV]", 28, binLimitsPt);
72
73 TH1::AddDirectory(oldStatus);
75ec0f41 74}
75
76//____________________________________________________________________
16e24ca3 77dNdEtaAnalysis::~dNdEtaAnalysis()
78{
79 // destructor
80
0ab29cfa 81 if (fData)
82 {
83 delete fData;
84 fData = 0;
85 }
16e24ca3 86
16e24ca3 87 for (Int_t i=0; i<kVertexBinning; ++i)
88 {
3dfa46a4 89 if (fdNdEtaNotEventCorrected[i])
90 {
91 delete fdNdEtaNotEventCorrected[i];
92 fdNdEtaNotEventCorrected[i] = 0;
93 }
0ab29cfa 94 if (fdNdEta[i])
95 {
96 delete fdNdEta[i];
97 fdNdEta[i] = 0;
98 }
99 if (fdNdEtaPtCutOffCorrected[i])
100 {
101 delete fdNdEtaPtCutOffCorrected[i];
102 fdNdEtaPtCutOffCorrected[i] = 0;
103 }
16e24ca3 104 }
1afae8ff 105
0ab29cfa 106 if (fPtDist)
107 {
108 delete fPtDist;
109 fPtDist = 0;
110 }
16e24ca3 111}
112
113//_____________________________________________________________________________
114dNdEtaAnalysis::dNdEtaAnalysis(const dNdEtaAnalysis &c) :
115 TNamed(c),
45e97e28 116 fData(0),
1afae8ff 117 fPtDist(0)
16e24ca3 118{
119 //
120 // dNdEtaAnalysis copy constructor
121 //
122
123 ((dNdEtaAnalysis &) c).Copy(*this);
124}
125
126//_____________________________________________________________________________
127dNdEtaAnalysis &dNdEtaAnalysis::operator=(const dNdEtaAnalysis &c)
128{
129 //
130 // Assignment operator
131 //
132
133 if (this != &c) ((dNdEtaAnalysis &) c).Copy(*this);
134 return *this;
135}
136
137//_____________________________________________________________________________
138void dNdEtaAnalysis::Copy(TObject &c) const
139{
140 //
141 // Copy function
142 //
143
144 dNdEtaAnalysis& target = (dNdEtaAnalysis &) c;
145
74fd10b3 146 target.fData = dynamic_cast<AliCorrection*> (fData->Clone());
16e24ca3 147
148 for (Int_t i=0; i<kVertexBinning; ++i)
1afae8ff 149 {
3dfa46a4 150 target.fdNdEtaNotEventCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaNotEventCorrected[i]->Clone());
74fd10b3 151 target.fdNdEta[i] = dynamic_cast<TH1F*> (fdNdEta[i]->Clone());
152 target.fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaPtCutOffCorrected[i]->Clone());
1afae8ff 153 }
16e24ca3 154
74fd10b3 155 target.fPtDist = dynamic_cast<TH1F*> (fPtDist->Clone());
45e97e28 156
16e24ca3 157 TNamed::Copy((TNamed &) c);
158}
159
160//____________________________________________________________________
74fd10b3 161void dNdEtaAnalysis::FillTrack(Float_t vtx, Float_t eta, Float_t pt)
6bf0714d 162{
163 // fills a track into the histograms
164
74fd10b3 165 fData->GetTrackCorrection()->FillMeas(vtx, eta, pt);
75ec0f41 166}
167
168//____________________________________________________________________
74fd10b3 169void dNdEtaAnalysis::FillEvent(Float_t vtx, Float_t n)
6bf0714d 170{
171 // fills an event into the histograms
172
74fd10b3 173 fData->GetEventCorrection()->FillMeas(vtx, n);
75ec0f41 174}
175
176//____________________________________________________________________
3dfa46a4 177void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, AlidNdEtaCorrection::CorrectionType correctionType, Int_t multCut)
fcf2fb36 178{
74fd10b3 179 //
180 // correct with the given correction values and calculate dNdEta and pT distribution
181 // the corrections that are applied can be steered by the flag correctionType
3dfa46a4 182 // the measured result is not used up to a multiplicity of multCut (the bin at multCut is the first that is used!)
74fd10b3 183 //
184
185 // TODO put tag somewhere which corrections have been applied
186
3dfa46a4 187 if (multCut > 1)
188 {
189 Printf("ERROR: A bigger multiplicity cut than 1 is not possible in the current implementation");
190 return;
191 }
192
74fd10b3 193 // set corrections to 1
194 fData->SetCorrectionToUnity();
195
196 if (correction && correctionType != AlidNdEtaCorrection::kNone)
197 {
198 TH3F* trackCorr = fData->GetTrackCorrection()->GetCorrectionHistogram();
199 TH2F* eventCorr = fData->GetEventCorrection()->GetCorrectionHistogram();
200
201 if (correctionType >= AlidNdEtaCorrection::kTrack2Particle)
202 trackCorr->Multiply(correction->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetCorrectionHistogram());
6bf0714d 203
dd367a14 204 if (correctionType >= AlidNdEtaCorrection::kVertexReco)
74fd10b3 205 {
206 trackCorr->Multiply(correction->GetVertexRecoCorrection()->GetTrackCorrection()->GetCorrectionHistogram());
207 eventCorr->Multiply(correction->GetVertexRecoCorrection()->GetEventCorrection()->GetCorrectionHistogram());
dd367a14 208 }
74fd10b3 209
210 switch (correctionType)
211 {
212 case AlidNdEtaCorrection::kINEL :
213 {
214 trackCorr->Multiply(correction->GetTriggerBiasCorrectionINEL()->GetTrackCorrection()->GetCorrectionHistogram());
215 eventCorr->Multiply(correction->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->GetCorrectionHistogram());
216 break;
217 }
218 case AlidNdEtaCorrection::kNSD :
219 {
220 trackCorr->Multiply(correction->GetTriggerBiasCorrectionNSD()->GetTrackCorrection()->GetCorrectionHistogram());
221 eventCorr->Multiply(correction->GetTriggerBiasCorrectionNSD()->GetEventCorrection()->GetCorrectionHistogram());
222 break;
223 }
224 case AlidNdEtaCorrection::kND :
225 {
226 trackCorr->Multiply(correction->GetTriggerBiasCorrectionND()->GetTrackCorrection()->GetCorrectionHistogram());
227 eventCorr->Multiply(correction->GetTriggerBiasCorrectionND()->GetEventCorrection()->GetCorrectionHistogram());
228 break;
229 }
230 default : break;
231 }
232 }
233 else
fcf2fb36 234 printf("INFO: No correction applied\n");
235
74fd10b3 236 fData->Multiply();
3dfa46a4 237 fData->PrintInfo(ptCut);
74fd10b3 238
239 TH3F* dataHist = fData->GetTrackCorrection()->GetGeneratedHistogram();
240
241 // integrate multiplicity axis out (include under/overflow bins!!!)
242 TH2F* tmp = fData->GetEventCorrection()->GetGeneratedHistogram();
3dfa46a4 243
244 // multiplicity cut correction
245 // correct for not-measured events with too small multiplicity
246 // the measured result is not used up to a multiplicity of multCut (the bin at multCut is the first that is used!)
247 TH1D* vertexHistUncorrected = tmp->ProjectionX("_px", tmp->GetYaxis()->FindBin(multCut), tmp->GetNbinsY() + 1, "e");
248 TH1D* vertexHist = (TH1D*) vertexHistUncorrected->Clone("vertexHist");
249
250 Printf("Using %d events (above a cut of %d)", (Int_t) vertexHistUncorrected->Integral(), multCut);
251 if (correction && correctionType >= AlidNdEtaCorrection::kVertexReco)
252 {
253 TH1* eventFraction = correction->GetMeasuredEventFraction(correctionType, multCut);
254 vertexHist->Divide(eventFraction);
255 Printf("Multiplicity cut correction: Corrected from %d events to %d events", (Int_t) vertexHistUncorrected->Integral(), (Int_t) vertexHist->Integral());
256 }
fcf2fb36 257
1afae8ff 258 // create pt hist
259 {
260 // reset all ranges
74fd10b3 261 dataHist->GetXaxis()->SetRange(0, 0);
262 dataHist->GetYaxis()->SetRange(0, 0);
263 dataHist->GetZaxis()->SetRange(0, 0);
1afae8ff 264
265 // vtx cut
74fd10b3 266 Int_t vertexBinBegin = dataHist->GetXaxis()->FindBin(-5);
267 Int_t vertexBinEnd = dataHist->GetXaxis()->FindBin(5);
268 dataHist->GetXaxis()->SetRange(vertexBinBegin, vertexBinEnd);
269 Float_t nEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd);
1afae8ff 270
74fd10b3 271 if (nEvents > 0)
272 {
273 // eta cut
274 dataHist->GetYaxis()->SetRange(dataHist->GetYaxis()->FindBin(-0.8), dataHist->GetYaxis()->FindBin(0.8));
275 Float_t etaWidth = 1.6;
1afae8ff 276
74fd10b3 277 TH1D* ptHist = dynamic_cast<TH1D*> (dataHist->Project3D("ze"));
1afae8ff 278
74fd10b3 279 for (Int_t i=1; i<=fPtDist->GetNbinsX(); ++i)
280 {
281 Float_t binSize = fPtDist->GetBinWidth(i);
282 fPtDist->SetBinContent(i, ptHist->GetBinContent(i) / binSize / nEvents / etaWidth);
283 fPtDist->SetBinError(i, ptHist->GetBinError(i) / binSize / nEvents / etaWidth);
284 }
1afae8ff 285
74fd10b3 286 delete ptHist;
287 }
288 else
289 printf("ERROR: nEvents is 0!\n");
1afae8ff 290 }
291
292 // reset all ranges
74fd10b3 293 dataHist->GetXaxis()->SetRange(0, 0);
294 dataHist->GetYaxis()->SetRange(0, 0);
295 dataHist->GetZaxis()->SetRange(0, 0);
1afae8ff 296
847489f7 297 // integrate over pt (with pt cut)
1afae8ff 298 Int_t ptLowBin = 1;
299 if (ptCut > 0)
74fd10b3 300 ptLowBin = dataHist->GetZaxis()->FindBin(ptCut);
1afae8ff 301
74fd10b3 302 dataHist->GetZaxis()->SetRange(ptLowBin, dataHist->GetZaxis()->GetNbins());
303 printf("range %d %d\n", ptLowBin, dataHist->GetZaxis()->GetNbins());
304 TH2D* vtxVsEta = dynamic_cast<TH2D*> (dataHist->Project3D("yx2e"));
305
306 dataHist->GetZaxis()->SetRange(0, 0);
307 vtxVsEta->GetXaxis()->SetTitle(dataHist->GetXaxis()->GetTitle());
308 vtxVsEta->GetYaxis()->SetTitle(dataHist->GetYaxis()->GetTitle());
1afae8ff 309
847489f7 310 if (vtxVsEta == 0)
311 {
312 printf("ERROR: pt integration failed\n");
313 return;
fcf2fb36 314 }
315
144ff489 316 //new TCanvas; vtxVsEta->DrawCopy();
317 //vtxVsEta->Rebin2D(1, 4);
318
0448e811 319 const Float_t vertexRange = 9.99;
b4b9cacc 320
74fd10b3 321 for (Int_t iEta=1; iEta<=vtxVsEta->GetNbinsY(); iEta++)
5af55649 322 {
7029240a 323 // do we have several histograms for different vertex positions?
b4b9cacc 324 Int_t vertexBinGlobalBegin = vertexHist->GetXaxis()->FindBin(-vertexRange);
325 Int_t vertexBinWidth = (vertexHist->GetXaxis()->FindBin(vertexRange) - vertexBinGlobalBegin + 1) / (kVertexBinning-1);
326 //printf("vertexBinGlobalBegin = %d, vertexBinWidth = %d\n", vertexBinGlobalBegin, vertexBinWidth);
7029240a 327 for (Int_t vertexPos=0; vertexPos<kVertexBinning; ++vertexPos)
328 {
2e88424e 329
b4b9cacc 330 Int_t vertexBinBegin = vertexBinGlobalBegin;
2e88424e 331 Int_t vertexBinEnd = vertexBinGlobalBegin + vertexBinWidth * (kVertexBinning-1);
fcf2fb36 332
333 // the first histogram is always for the whole vertex range
334 if (vertexPos > 0)
335 {
b4b9cacc 336 vertexBinBegin = vertexBinGlobalBegin + vertexBinWidth * (vertexPos-1);
fcf2fb36 337 vertexBinEnd = vertexBinBegin + vertexBinWidth;
338 }
7029240a 339
b4b9cacc 340 //printf("vertexBinBegin = %d, vertexBinEnd = %d\n", vertexBinBegin, vertexBinEnd);
341
3dfa46a4 342 Float_t totalEventsUncorrected = vertexHistUncorrected->Integral(vertexBinBegin, vertexBinEnd - 1);
343 if (totalEventsUncorrected == 0)
344 {
345 printf("WARNING: No events (uncorrected) for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
346 continue;
347 }
348
74fd10b3 349 Float_t totalEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd - 1);
5af55649 350 if (totalEvents == 0)
351 {
352 printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
353 continue;
7029240a 354 }
7029240a 355
5af55649 356 Float_t sum = 0;
357 Float_t sumError2 = 0;
358 for (Int_t iVtx = vertexBinBegin; iVtx < vertexBinEnd; iVtx++)
2e88424e 359 {
144ff489 360 if (vtxVsEta->GetBinContent(iVtx, iEta) != 0)
5af55649 361 {
847489f7 362 sum = sum + vtxVsEta->GetBinContent(iVtx, iEta);
2e88424e 363
144ff489 364 if (sumError2 > 10e30)
365 printf("WARNING: sum of error2 is dangerously large - be prepared for crash... ");
2e88424e 366
847489f7 367 sumError2 = sumError2 + TMath::Power(vtxVsEta->GetBinError(iVtx, iEta),2);
5af55649 368 }
369 }
7029240a 370
1afae8ff 371 Float_t ptCutOffCorrection = 1;
74fd10b3 372 if (correction && ptCut > 0)
373 ptCutOffCorrection = correction->GetMeasuredFraction(correctionType, ptCut, vtxVsEta->GetYaxis()->GetBinCenter(iEta));
1afae8ff 374
847489f7 375 if (ptCutOffCorrection <= 0)
376 {
377 printf("UNEXPECTED: ptCutOffCorrection is %f for hist %d %d %d\n", ptCutOffCorrection, vertexPos, vertexBinBegin, vertexBinEnd);
378 continue;
379 }
380
9e952c39 381 //printf("Eta: %d Vertex Range: %d %d, Event Count %f, Track Sum: %f, Track Sum corrected: %f\n", iEta, vertexBinBegin, vertexBinEnd, totalEvents, sum, sum / ptCutOffCorrection);
382
144ff489 383 Int_t bin = fdNdEta[vertexPos]->FindBin(vtxVsEta->GetYaxis()->GetBinCenter(iEta));
384 if (bin > 0 && bin < fdNdEta[vertexPos]->GetNbinsX())
385 {
3dfa46a4 386 Float_t dndeta = sum / totalEventsUncorrected;
387 Float_t error = TMath::Sqrt(sumError2) / totalEventsUncorrected;
388
389 dndeta = dndeta / fdNdEta[vertexPos]->GetBinWidth(bin);
390 error = error / fdNdEta[vertexPos]->GetBinWidth(bin);
391
392 fdNdEtaNotEventCorrected[vertexPos]->SetBinContent(bin, dndeta);
393 fdNdEtaNotEventCorrected[vertexPos]->SetBinError(bin, error);
394
395 dndeta = sum / totalEvents;
396 error = TMath::Sqrt(sumError2) / totalEvents;
7029240a 397
144ff489 398 dndeta = dndeta / fdNdEta[vertexPos]->GetBinWidth(bin);
399 error = error / fdNdEta[vertexPos]->GetBinWidth(bin);
7029240a 400
144ff489 401 fdNdEta[vertexPos]->SetBinContent(bin, dndeta);
402 fdNdEta[vertexPos]->SetBinError(bin, error);
1afae8ff 403
144ff489 404 dndeta /= ptCutOffCorrection;
405 error /= ptCutOffCorrection;
1afae8ff 406
144ff489 407 fdNdEtaPtCutOffCorrected[vertexPos]->SetBinContent(bin, dndeta);
408 fdNdEtaPtCutOffCorrected[vertexPos]->SetBinError(bin, error);
2e88424e 409
dd367a14 410 //Printf("Bin %d has dN/deta = %f", bin, dndeta);
144ff489 411 }
75ec0f41 412 }
847489f7 413 }
144ff489 414
dd367a14 415 //new TCanvas; fdNdEta[0]->DrawCopy();
75ec0f41 416}
417
75ec0f41 418//____________________________________________________________________
6bf0714d 419void dNdEtaAnalysis::SaveHistograms()
420{
421 // save the histograms to a directory with the name of this class (retrieved from TNamed)
75ec0f41 422
7029240a 423 gDirectory->mkdir(GetName());
424 gDirectory->cd(GetName());
5fbd0b17 425
1afae8ff 426 if (fData)
427 {
74fd10b3 428 fData->SaveHistograms();
1afae8ff 429 }
430 else
431 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fData is 0\n");
432
d09fb536 433 if (fPtDist)
434 fPtDist ->Write();
435 else
436 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fPtDist is 0\n");
437
7029240a 438 for (Int_t i=0; i<kVertexBinning; ++i)
1afae8ff 439 {
3dfa46a4 440 if (fdNdEtaNotEventCorrected[i])
441 fdNdEtaNotEventCorrected[i]->Write();
442 else
443 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEtaNotEventCorrected[%d] is 0\n", i);
444
1afae8ff 445 if (fdNdEta[i])
446 fdNdEta[i]->Write();
447 else
448 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEta[%d] is 0\n", i);
449
450 if (fdNdEtaPtCutOffCorrected[i])
451 fdNdEtaPtCutOffCorrected[i]->Write();
452 else
453 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEtaPtCutOffCorrected[%d] is 0\n", i);
454 }
75ec0f41 455
456 gDirectory->cd("../");
457}
458
74fd10b3 459void dNdEtaAnalysis::LoadHistograms(const Char_t* dir)
5fbd0b17 460{
6bf0714d 461 // loads the histograms from a directory with the name of this class (retrieved from TNamed)
462
74fd10b3 463 if (!dir)
464 dir = GetName();
5fbd0b17 465
74fd10b3 466 gDirectory->cd(dir);
5fbd0b17 467
74fd10b3 468 fData->LoadHistograms();
5fbd0b17 469
470 for (Int_t i=0; i<kVertexBinning; ++i)
1afae8ff 471 {
3dfa46a4 472 if (dynamic_cast<TH1F*> (gDirectory->Get(fdNdEtaNotEventCorrected[i]->GetName())))
473 fdNdEtaNotEventCorrected[i] = dynamic_cast<TH1F*> (gDirectory->Get(fdNdEtaNotEventCorrected[i]->GetName()));
474
74fd10b3 475 fdNdEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(fdNdEta[i]->GetName()));
476 fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (gDirectory->Get(fdNdEtaPtCutOffCorrected[i]->GetName()));
1afae8ff 477 }
5fbd0b17 478
74fd10b3 479 fPtDist = dynamic_cast<TH1F*> (gDirectory->Get(fPtDist->GetName()));
0ab29cfa 480
5fbd0b17 481 gDirectory->cd("../");
482}
483
ceb5d1b5 484//____________________________________________________________________
74fd10b3 485void dNdEtaAnalysis::DrawHistograms(Bool_t simple)
ceb5d1b5 486{
6bf0714d 487 // draws the histograms
74fd10b3 488
489 if (!simple)
490 {
491 if (fData)
492 fData->DrawHistograms(GetName());
6bf0714d 493
74fd10b3 494 TCanvas* canvas = new TCanvas(Form("%s_dNdEtaAnalysis", GetName()), Form("%s_dNdEtaAnalysis", GetName()), 800, 400);
495 canvas->Divide(2, 1);
ceb5d1b5 496
74fd10b3 497 canvas->cd(1);
498 if (fdNdEtaPtCutOffCorrected[0])
144ff489 499 fdNdEtaPtCutOffCorrected[0]->DrawCopy();
ceb5d1b5 500
74fd10b3 501 if (fdNdEta[0])
502 {
503 fdNdEta[0]->SetLineColor(kRed);
144ff489 504 fdNdEta[0]->DrawCopy("SAME");
74fd10b3 505 }
1afae8ff 506
74fd10b3 507 canvas->cd(2);
508 if (fPtDist)
144ff489 509 fPtDist->DrawCopy();
1afae8ff 510 }
511
fcf2fb36 512 // histograms for different vertices?
513 if (kVertexBinning > 0)
514 {
515 // doesnt work, but i dont get it, giving up...
74fd10b3 516 TCanvas* canvas2 = new TCanvas(Form("%s_dNdEtaAnalysisVtx", GetName()), Form("%s_dNdEtaAnalysisVtx", GetName()), 450, 450);
517 TCanvas* canvas3 = 0;
518 if (!simple)
519 canvas3 = new TCanvas(Form("%s_dNdEtaAnalysisVtx_noptcutoff", GetName()), Form("%s_dNdEtaAnalysisVtx_noptcutoff", GetName()), 450, 450);
520
fcf2fb36 521 //Int_t yPads = (Int_t) TMath::Ceil(((Double_t) kVertexBinning - 1) / 2);
522 //printf("%d\n", yPads);
523 //canvas2->Divide(2, yPads);
524
0448e811 525 TLegend* legend = new TLegend(0.4, 0.2, 0.6, 0.4);
fcf2fb36 526
5af55649 527 for (Int_t i=0; i<kVertexBinning; ++i)
fcf2fb36 528 {
1afae8ff 529 if (fdNdEtaPtCutOffCorrected[i])
fcf2fb36 530 {
74fd10b3 531 canvas2->cd();
532
1afae8ff 533 fdNdEtaPtCutOffCorrected[i]->SetLineColor(i+1);
144ff489 534 fdNdEtaPtCutOffCorrected[i]->DrawCopy((i == 0) ? "" : "SAME");
1afae8ff 535 legend->AddEntry(fdNdEtaPtCutOffCorrected[i], (i == 0) ? "Vtx All" : Form("Vtx Bin %d", i-1));
fcf2fb36 536 }
74fd10b3 537 if (canvas3 && fdNdEta[i])
538 {
539 canvas3->cd();
540
541 fdNdEta[i]->SetLineColor(i+1);
144ff489 542 fdNdEta[i]->DrawCopy((i == 0) ? "" : "SAME");
74fd10b3 543 }
fcf2fb36 544 }
545
74fd10b3 546 canvas2->cd();
fcf2fb36 547 legend->Draw();
74fd10b3 548 canvas2->SaveAs(Form("%s_%s.gif", canvas2->GetName(), GetName()));
549
550 if (canvas3)
551 {
552 canvas3->cd();
553 legend->Draw();
554 }
fcf2fb36 555 }
74fd10b3 556
557 if (kVertexBinning == 3)
558 {
559 TH1* clone = dynamic_cast<TH1*> (fdNdEtaPtCutOffCorrected[1]->Clone("clone"));
560 TH1* clone2 = dynamic_cast<TH1*> (fdNdEtaPtCutOffCorrected[2]->Clone("clone2"));
561
562 if (clone && clone2)
563 {
564 TCanvas* canvas4 = new TCanvas(Form("%s_dNdEtaAnalysisVtxRatios", GetName()), Form("%s_dNdEtaAnalysisVtxRatios", GetName()), 450, 450);
565
566 clone->Divide(fdNdEtaPtCutOffCorrected[0]);
567 clone->GetYaxis()->SetRangeUser(0.95, 1.05);
144ff489 568 clone->DrawCopy();
74fd10b3 569
570 clone2->Divide(fdNdEtaPtCutOffCorrected[0]);
144ff489 571 clone2->DrawCopy("SAME");
74fd10b3 572
573 TLine* line = new TLine(-1, 1, 1, 1);
574 line->Draw();
575
576 canvas4->SaveAs(Form("%s_%s.gif", canvas4->GetName(), GetName()));
577 }
578 }
7029240a 579}
580
581Long64_t dNdEtaAnalysis::Merge(TCollection* list)
582{
583 // Merges a list of dNdEtaAnalysis objects with this one.
584 // This is needed for PROOF.
585 // Returns the number of merged objects (including this)
586
587 if (!list)
588 return 0;
589
590 if (list->IsEmpty())
591 return 1;
592
593 TIterator* iter = list->MakeIterator();
594 TObject* obj;
595
596 // sub collections
3dfa46a4 597 const Int_t nCollections = 3 * kVertexBinning + 2; // 2 standalone hists, two arrays of size kVertexBinning
7029240a 598 TList* collections[nCollections];
599 for (Int_t i=0; i<nCollections; ++i)
600 collections[i] = new TList;
601
602 Int_t count = 0;
603 while ((obj = iter->Next()))
604 {
605 dNdEtaAnalysis* entry = dynamic_cast<dNdEtaAnalysis*> (obj);
606 if (entry == 0)
607 continue;
608
45e97e28 609 collections[0]->Add(entry->fData);
74fd10b3 610 collections[1]->Add(entry->fPtDist);
7029240a 611
612 for (Int_t i=0; i<kVertexBinning; ++i)
1afae8ff 613 {
74fd10b3 614 collections[2+i]->Add(entry->fdNdEta[i]);
615 collections[2+kVertexBinning+i]->Add(entry->fdNdEtaPtCutOffCorrected[i]);
3dfa46a4 616 collections[2+2*kVertexBinning+i]->Add(entry->fdNdEtaNotEventCorrected[i]);
1afae8ff 617 }
7029240a 618
619 ++count;
620 }
621
45e97e28 622 fData->Merge(collections[0]);
74fd10b3 623 fPtDist->Merge(collections[1]);
7029240a 624 for (Int_t i=0; i<kVertexBinning; ++i)
1afae8ff 625 {
74fd10b3 626 fdNdEta[i]->Merge(collections[2+i]);
3dfa46a4 627 fdNdEtaPtCutOffCorrected[i]->Merge(collections[2+kVertexBinning+i]);
628 fdNdEtaNotEventCorrected[i]->Merge(collections[2+kVertexBinning+i]);
1afae8ff 629 }
7029240a 630
631 for (Int_t i=0; i<nCollections; ++i)
632 delete collections[i];
633
634 return count+1;
ceb5d1b5 635}