3 #include "dNdEtaAnalysis.h"
11 #include <TCollection.h>
12 #include <TIterator.h>
16 #include <TParameter.h>
18 #include "AlidNdEtaCorrection.h"
19 #include <AliCorrection.h>
20 #include <AliPWG0Helper.h>
21 #include <AliCorrectionMatrix2D.h>
22 #include <AliCorrectionMatrix3D.h>
24 //____________________________________________________________________
25 ClassImp(dNdEtaAnalysis)
27 //____________________________________________________________________
28 dNdEtaAnalysis::dNdEtaAnalysis() :
33 fAnalysisMode(AliPWG0Helper::kInvalid),
36 // default constructor
38 for (Int_t i=0; i<kVertexBinning; ++i)
41 fdNdEtaPtCutOffCorrected[i] = 0;
45 //____________________________________________________________________
46 dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title, AliPWG0Helper::AnalysisMode analysisMode) :
51 fAnalysisMode(analysisMode),
56 fData = new AliCorrection("Analysis", Form("%s Analysis", title), analysisMode);
58 // do not add this hists to the directory
59 Bool_t oldStatus = TH1::AddDirectoryStatus();
60 TH1::AddDirectory(kFALSE);
62 fMult = new TH1F("TriggeredMultiplicity", "Triggered Events;raw multiplicity;entries", 1000, -0.5, 999.5);
64 TH1* histForBinning = fData->GetTrackCorrection()->GetGeneratedHistogram();
65 fdNdEta[0] = new TH1F("dNdEta", "dN_{ch}/d#eta;#eta;dN_{ch}/d#eta", histForBinning->GetNbinsY(), histForBinning->GetYaxis()->GetXbins()->GetArray());
67 fdNdEtaPtCutOffCorrected[0] = dynamic_cast<TH1F*> (fdNdEta[0]->Clone(Form("%s_corrected", fdNdEta[0]->GetName())));
69 for (Int_t i=1; i<kVertexBinning; ++i)
71 fdNdEta[i] = dynamic_cast<TH1F*> (fdNdEta[0]->Clone(Form("%s_%d", fdNdEta[0]->GetName(), i)));
72 fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaPtCutOffCorrected[0]->Clone(Form("%s_%d", fdNdEtaPtCutOffCorrected[0]->GetName(), i)));
75 fPtDist = new TH1F("Pt", "p_{T} distribution;p_{T} [GeV/c];#frac{dN}{d#eta dp_{T}} [c/GeV]", histForBinning->GetNbinsZ(), histForBinning->GetZaxis()->GetXbins()->GetArray());
77 TH1::AddDirectory(oldStatus);
80 //____________________________________________________________________
81 dNdEtaAnalysis::~dNdEtaAnalysis()
97 for (Int_t i=0; i<kVertexBinning; ++i)
104 if (fdNdEtaPtCutOffCorrected[i])
106 delete fdNdEtaPtCutOffCorrected[i];
107 fdNdEtaPtCutOffCorrected[i] = 0;
118 //_____________________________________________________________________________
119 dNdEtaAnalysis::dNdEtaAnalysis(const dNdEtaAnalysis &c) :
124 fAnalysisMode(AliPWG0Helper::kInvalid),
128 // dNdEtaAnalysis copy constructor
131 ((dNdEtaAnalysis &) c).Copy(*this);
134 //_____________________________________________________________________________
135 dNdEtaAnalysis &dNdEtaAnalysis::operator=(const dNdEtaAnalysis &c)
138 // Assignment operator
141 if (this != &c) ((dNdEtaAnalysis &) c).Copy(*this);
145 //_____________________________________________________________________________
146 void dNdEtaAnalysis::Copy(TObject &c) const
152 dNdEtaAnalysis& target = (dNdEtaAnalysis &) c;
154 target.fData = dynamic_cast<AliCorrection*> (fData->Clone());
155 target.fMult = dynamic_cast<TH1F*> (fMult->Clone());
157 for (Int_t i=0; i<kVertexBinning; ++i)
159 target.fdNdEta[i] = dynamic_cast<TH1F*> (fdNdEta[i]->Clone());
160 target.fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (fdNdEtaPtCutOffCorrected[i]->Clone());
163 target.fPtDist = dynamic_cast<TH1F*> (fPtDist->Clone());
165 target.fAnalysisMode = fAnalysisMode;
168 TNamed::Copy((TNamed &) c);
171 //____________________________________________________________________
172 void dNdEtaAnalysis::FillTrack(Float_t vtx, Float_t eta, Float_t pt)
174 // fills a track into the histograms
176 fData->GetTrackCorrection()->FillMeas(vtx, eta, pt);
179 //____________________________________________________________________
180 void dNdEtaAnalysis::FillEvent(Float_t vtx, Float_t n)
182 // fills an event into the histograms
184 fData->GetEventCorrection()->FillMeas(vtx, n);
187 //____________________________________________________________________
188 void dNdEtaAnalysis::FillTriggeredEvent(Float_t n)
190 // fills a triggered event into the histograms
195 //____________________________________________________________________
196 void dNdEtaAnalysis::Finish(AlidNdEtaCorrection* correction, Float_t ptCut, AlidNdEtaCorrection::CorrectionType correctionType, const char* tag)
199 // correct with the given correction values and calculate dNdEta and pT distribution
200 // the corrections that are applied can be steered by the flag correctionType
201 // the measured result is not used up to a multiplicity of multCut (the bin at multCut is the first that is used!)
204 fTag.Form("Correcting dN/deta spectrum >>> %s <<<. Correction type: %d, pt cut: %.2f.", tag, (Int_t) correctionType, ptCut);
205 Printf("\n\n%s", fTag.Data());
207 // set corrections to 1
208 fData->SetCorrectionToUnity();
210 if (correction && correctionType != AlidNdEtaCorrection::kNone)
212 TH3F* trackCorr = fData->GetTrackCorrection()->GetCorrectionHistogram();
213 TH2F* eventCorr = fData->GetEventCorrection()->GetCorrectionHistogram();
215 if (correctionType >= AlidNdEtaCorrection::kTrack2Particle)
216 trackCorr->Multiply(correction->GetTrack2ParticleCorrection()->GetTrackCorrection()->GetCorrectionHistogram());
218 if (correctionType >= AlidNdEtaCorrection::kVertexReco)
220 trackCorr->Multiply(correction->GetVertexRecoCorrection()->GetTrackCorrection()->GetCorrectionHistogram());
221 eventCorr->Multiply(correction->GetVertexRecoCorrection()->GetEventCorrection()->GetCorrectionHistogram());
223 // set bin with multiplicity 0 to unity (correction has no meaning in this bin)
224 for (Int_t i=0; i<=eventCorr->GetNbinsX()+1; i++)
225 eventCorr->SetBinContent(i, 1, 1);
228 switch (correctionType)
230 case AlidNdEtaCorrection::kINEL :
232 trackCorr->Multiply(correction->GetTriggerBiasCorrectionINEL()->GetTrackCorrection()->GetCorrectionHistogram());
233 eventCorr->Multiply(correction->GetTriggerBiasCorrectionINEL()->GetEventCorrection()->GetCorrectionHistogram());
236 case AlidNdEtaCorrection::kNSD :
238 trackCorr->Multiply(correction->GetTriggerBiasCorrectionNSD()->GetTrackCorrection()->GetCorrectionHistogram());
239 eventCorr->Multiply(correction->GetTriggerBiasCorrectionNSD()->GetEventCorrection()->GetCorrectionHistogram());
242 case AlidNdEtaCorrection::kND :
244 trackCorr->Multiply(correction->GetTriggerBiasCorrectionND()->GetTrackCorrection()->GetCorrectionHistogram());
245 eventCorr->Multiply(correction->GetTriggerBiasCorrectionND()->GetEventCorrection()->GetCorrectionHistogram());
252 printf("INFO: No correction applied\n");
254 TH2F* rawMeasured = (TH2F*) fData->GetEventCorrection()->GetMeasuredHistogram()->Clone("rawMeasured");
258 if (correctionType >= AlidNdEtaCorrection::kVertexReco)
260 // There are no events with vertex that have 0 multiplicity, therefore
261 // populate bin with 0 multiplicity with the following idea:
262 // alpha = triggered events with vertex at a given vertex position / all triggered events with vertex
263 // triggered events without vertex and 0 multiplicity at a given vertex position = alpha * all triggered events with 0 multiplicity
264 // afterwards we still correct for the trigger efficiency
266 //TH2* measuredEvents = fData->GetEventCorrection()->GetMeasuredHistogram();
267 TH2* correctedEvents = fData->GetEventCorrection()->GetGeneratedHistogram();
269 TH2* eTrig = correction->GetVertexRecoCorrection()->GetEventCorrection()->GetGeneratedHistogram();
270 TH2* eTrigVtx = correction->GetVertexRecoCorrection()->GetEventCorrection()->GetMeasuredHistogram();
271 TH1* eTrigVtx_projx = eTrigVtx->ProjectionX("eTrigVtx_projx", 2, rawMeasured->GetNbinsY()+1);
273 //new TCanvas; correctedEvents->DrawCopy("TEXT");
275 // start above 0 mult. bin with integration
276 TH1* vertexDist = rawMeasured->ProjectionX("vertexdist_measured", 2, rawMeasured->GetNbinsY()+1);
277 //new TCanvas; vertexDist->DrawCopy();
279 Int_t allEventsWithVertex = (Int_t) vertexDist->Integral(0, vertexDist->GetNbinsX()+1); // include under/overflow!
280 Int_t triggeredEventsWith0Mult = (Int_t) fMult->GetBinContent(1);
282 Printf("%d triggered events with 0 mult. -- %d triggered events with vertex", triggeredEventsWith0Mult, allEventsWithVertex);
284 for (Int_t i = 1; i <= rawMeasured->GetNbinsX(); i++)
286 Double_t alpha = (Double_t) vertexDist->GetBinContent(i) / allEventsWithVertex;
287 Double_t events = alpha * triggeredEventsWith0Mult;
289 if (eTrigVtx_projx->GetBinContent(i) == 0)
292 Double_t fZ = eTrigVtx_projx->Integral(0, eTrigVtx_projx->GetNbinsX()+1) / eTrigVtx_projx->GetBinContent(i) *
293 eTrig->GetBinContent(i, 1) / eTrig->Integral(0, eTrig->GetNbinsX()+1, 1, 1);
297 // multiply with trigger correction if set above
298 events *= fData->GetEventCorrection()->GetCorrectionHistogram()->GetBinContent(i, 1);
300 Printf("Bin %d, alpha is %.2f, fZ is %.3f, number of events with 0 mult.: %.2f", i, alpha * 100., fZ, events);
302 correctedEvents->SetBinContent(i, 1, events);
305 //new TCanvas; correctedEvents->DrawCopy("TEXT");
308 fData->PrintInfo(ptCut);
310 TH3F* dataHist = fData->GetTrackCorrection()->GetGeneratedHistogram();
312 // integrate multiplicity axis out (include under/overflow bins!!!)
313 TH2F* tmp = fData->GetEventCorrection()->GetGeneratedHistogram();
315 TH1D* vertexHist = (TH1D*) tmp->ProjectionX("_px", 0, tmp->GetNbinsY() + 1, "e");
320 dataHist->GetXaxis()->SetRange(0, 0);
321 dataHist->GetYaxis()->SetRange(0, 0);
322 dataHist->GetZaxis()->SetRange(0, 0);
325 Int_t vertexBinBegin = dataHist->GetXaxis()->FindBin(-5);
326 Int_t vertexBinEnd = dataHist->GetXaxis()->FindBin(5);
327 dataHist->GetXaxis()->SetRange(vertexBinBegin, vertexBinEnd);
328 Float_t nEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd);
333 dataHist->GetYaxis()->SetRange(dataHist->GetYaxis()->FindBin(-0.8), dataHist->GetYaxis()->FindBin(0.8));
334 Float_t etaWidth = 1.6;
336 TH1D* ptHist = dynamic_cast<TH1D*> (dataHist->Project3D("ze"));
338 for (Int_t i=1; i<=fPtDist->GetNbinsX(); ++i)
340 Float_t binSize = fPtDist->GetBinWidth(i);
341 fPtDist->SetBinContent(i, ptHist->GetBinContent(i) / binSize / nEvents / etaWidth);
342 fPtDist->SetBinError(i, ptHist->GetBinError(i) / binSize / nEvents / etaWidth);
348 printf("ERROR: nEvents is 0!\n");
352 dataHist->GetXaxis()->SetRange(0, 0);
353 dataHist->GetYaxis()->SetRange(0, 0);
354 dataHist->GetZaxis()->SetRange(0, 0);
356 // integrate over pt (with pt cut)
359 ptLowBin = dataHist->GetZaxis()->FindBin(ptCut);
361 dataHist->GetZaxis()->SetRange(ptLowBin, dataHist->GetZaxis()->GetNbins()+1);
362 printf("pt range %d %d\n", ptLowBin, dataHist->GetZaxis()->GetNbins()+1);
363 TH2D* vtxVsEta = dynamic_cast<TH2D*> (dataHist->Project3D("yx2e"));
365 dataHist->GetZaxis()->SetRange(0, 0);
366 vtxVsEta->GetXaxis()->SetTitle(dataHist->GetXaxis()->GetTitle());
367 vtxVsEta->GetYaxis()->SetTitle(dataHist->GetYaxis()->GetTitle());
371 printf("ERROR: pt integration failed\n");
375 // clear result histograms
376 for (Int_t vertexPos=0; vertexPos<kVertexBinning; ++vertexPos)
378 fdNdEta[vertexPos]->Reset();
379 fdNdEtaPtCutOffCorrected[vertexPos]->Reset();
382 const Float_t vertexRangeBegin[kVertexBinning] = { -9.99, -9.99, 0.01 };
383 const Float_t vertexRangeEnd[kVertexBinning] = { 9.99, -0.01, 9.99 };
385 for (Int_t iEta=1; iEta<=vtxVsEta->GetNbinsY(); iEta++)
387 // loop over vertex ranges
388 for (Int_t vertexPos=0; vertexPos<kVertexBinning; ++vertexPos)
390 Int_t vertexBinBegin = vertexHist->GetXaxis()->FindBin(vertexRangeBegin[vertexPos]);
391 Int_t vertexBinEnd = vertexHist->GetXaxis()->FindBin(vertexRangeEnd[vertexPos]);
393 const Int_t *binBegin = 0;
394 // adjust acceptance range for SPD
395 if (fAnalysisMode == AliPWG0Helper::kSPD)
397 //const Int_t binBegin[30] = { 18, 16, 15, 14, 13, 13, 12, 11, 9, 7, 6, 5, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 3, 3, 3, 2, 2, 2, 2, 1 };
398 const Int_t binBeginSPD[30] = { -1, -1, -1, -1, 16, 14, 12, 10, 9, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 2, 2, 2, 2, 2, -1, -1, -1, -1}; // limit in correction map is 5
399 //const Int_t binBegin[30] = { -1, -1, -1, 17, 15, 14, 12, 10, 8, 7, 6, 5, 4, 4, 4, 4, 3, 3, 3, 3, 3, 2, 2, 2, 2, 2, 2, -1, -1, -1}; // limit in correction map is 10
400 binBegin = binBeginSPD;
403 else if (fAnalysisMode == AliPWG0Helper::kTPC)
405 // const Int_t binBegin[30] = { -1, -1, -1, -1, 16, 14, 12, 10, 9, 7, 6, 5, 5, 4, 4, 4, 4, 3, 3, 3, 3, 2, 2, 2, 2, 2, -1, -1, -1, -1}; // limit in correction map is 5
406 const Int_t binBeginTPC[30] = {-1, -1, -1, -1, -1, -1, -1, -1, 9, 3, 3, 3, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, -1, -1, -1, -1, -1, -1, -1, -1};
407 binBegin = binBeginTPC;
409 else if (fAnalysisMode == AliPWG0Helper::kTPCITS)
414 Int_t vtxBegin = binBegin[iEta - 1];
415 Int_t vtxEnd = 18 + 1 - binBegin[30 - iEta];
417 // eta range not accessible
421 //Float_t eta = vtxVsEta->GetYaxis()->GetBinCenter(iEta);
422 //Printf("%d %d | %d %d", vtxBegin, vertexHist->GetXaxis()->FindBin(GetVtxMin(eta)), vtxEnd, vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta)));
423 //vtxBegin = vertexHist->GetXaxis()->FindBin(GetVtxMin(eta));
424 //vtxEnd = vertexHist->GetXaxis()->FindBin(-GetVtxMin(-eta));
425 //printf("Eta bin: %d (%f) Vertex range: %d; before: %d %d ", iEta, eta, vertexPos, vertexBinBegin, vertexBinEnd);
426 vertexBinBegin = TMath::Max(vertexBinBegin, vtxBegin);
427 vertexBinEnd = TMath::Min(vertexBinEnd, vtxEnd);
428 //Printf("after: %d %d", vertexBinBegin, vertexBinEnd);
430 // no data for this bin
431 if (vertexBinBegin > vertexBinEnd)
433 //Printf("Bin empty. Skipped");
437 Float_t totalEvents = vertexHist->Integral(vertexBinBegin, vertexBinEnd);
438 if (totalEvents == 0)
440 printf("WARNING: No events for hist %d %d %d\n", vertexPos, vertexBinBegin, vertexBinEnd);
445 Float_t sumError2 = 0;
446 for (Int_t iVtx = vertexBinBegin; iVtx <= vertexBinEnd; iVtx++)
448 if (vtxVsEta->GetBinContent(iVtx, iEta) != 0)
450 sum = sum + vtxVsEta->GetBinContent(iVtx, iEta);
452 if (sumError2 > 10e30)
453 Printf("WARNING: sum of error2 is dangerously large - be prepared for crash... ");
455 sumError2 = sumError2 + TMath::Power(vtxVsEta->GetBinError(iVtx, iEta),2);
459 Float_t ptCutOffCorrection = 1;
460 if (correction && ptCut > 0)
461 ptCutOffCorrection = correction->GetMeasuredFraction(correctionType, ptCut, vtxVsEta->GetYaxis()->GetBinCenter(iEta));
463 if (ptCutOffCorrection <= 0)
465 printf("UNEXPECTED: ptCutOffCorrection is %f for hist %d %d %d\n", ptCutOffCorrection, vertexPos, vertexBinBegin, vertexBinEnd);
469 //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);
471 Int_t bin = fdNdEta[vertexPos]->FindBin(vtxVsEta->GetYaxis()->GetBinCenter(iEta));
472 if (bin > 0 && bin <= fdNdEta[vertexPos]->GetNbinsX())
474 Float_t dndeta = sum / totalEvents;
475 Float_t error = TMath::Sqrt(sumError2) / totalEvents;
477 dndeta = dndeta / fdNdEta[vertexPos]->GetBinWidth(bin);
478 error = error / fdNdEta[vertexPos]->GetBinWidth(bin);
480 fdNdEta[vertexPos]->SetBinContent(bin, dndeta);
481 fdNdEta[vertexPos]->SetBinError(bin, error);
483 dndeta /= ptCutOffCorrection;
484 error /= ptCutOffCorrection;
486 fdNdEtaPtCutOffCorrected[vertexPos]->SetBinContent(bin, dndeta);
487 fdNdEtaPtCutOffCorrected[vertexPos]->SetBinError(bin, error);
489 //Printf("Bin %d has dN/deta = %f +- %f; %f +- %f", bin, dndeta, error, fdNdEta[vertexPos]->GetBinContent(bin), fdNdEta[vertexPos]->GetBinError(bin));
495 //____________________________________________________________________
496 Float_t dNdEtaAnalysis::GetVtxMin(Float_t eta)
498 // helper function for the SPD acceptance
499 // the function returns the beginning of the acceptance window in z vertex position as function of eta
500 // to get the maximum: -GetVtxMin(-eta)
502 Float_t a[2] = { -15, 0 };
503 Float_t b[2] = { 0, -1.4 };
504 Float_t c[2] = { 15, -2.2 };
507 meanAB[0] = (b[0] + a[0]) / 2;
508 meanAB[1] = (b[1] + a[1]) / 2;
511 meanBC[0] = (c[0] + b[0]) / 2;
512 meanBC[1] = (c[1] + b[1]) / 2;
514 Float_t mAB = (b[1] - a[1]) / (b[0] - a[0]);
515 Float_t mBC = (c[1] - b[1]) / (c[0] - b[0]);
517 Float_t bAB = meanAB[1] - mAB * meanAB[0];
518 Float_t bBC = meanBC[1] - mBC * meanBC[0];
521 return 1.0 / mAB * eta - bAB / mAB;
523 return 1.0 / mBC * eta - bBC / mBC;
526 //____________________________________________________________________
527 void dNdEtaAnalysis::SaveHistograms()
529 // save the histograms to a directory with the name of this class (retrieved from TNamed)
531 gDirectory->mkdir(GetName());
532 gDirectory->cd(GetName());
536 fData->SaveHistograms();
539 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fData is 0\n");
546 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fMult is 0\n");
551 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fPtDist is 0\n");
553 for (Int_t i=0; i<kVertexBinning; ++i)
558 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEta[%d] is 0\n", i);
560 if (fdNdEtaPtCutOffCorrected[i])
561 fdNdEtaPtCutOffCorrected[i]->Write();
563 printf("dNdEtaAnalysis::SaveHistograms: UNEXPECTED: fdNdEtaPtCutOffCorrected[%d] is 0\n", i);
566 TNamed named("fTag", fTag.Data());
569 TParameter<Int_t> param("fAnalysisMode", fAnalysisMode);
572 gDirectory->cd("../");
575 void dNdEtaAnalysis::LoadHistograms(const Char_t* dir)
577 // loads the histograms from a directory with the name of this class (retrieved from TNamed)
584 fData->LoadHistograms();
585 fMult = dynamic_cast<TH1F*> (gDirectory->Get(fMult->GetName()));
587 for (Int_t i=0; i<kVertexBinning; ++i)
589 fdNdEta[i] = dynamic_cast<TH1F*> (gDirectory->Get(fdNdEta[i]->GetName()));
590 fdNdEtaPtCutOffCorrected[i] = dynamic_cast<TH1F*> (gDirectory->Get(fdNdEtaPtCutOffCorrected[i]->GetName()));
593 fPtDist = dynamic_cast<TH1F*> (gDirectory->Get(fPtDist->GetName()));
595 if (dynamic_cast<TNamed*> (gFile->Get("fTag")))
596 fTag = (dynamic_cast<TNamed*> (gFile->Get("fTag")))->GetTitle();
598 if (dynamic_cast<TParameter<Int_t>*> (gFile->Get("fAnalysisMode")))
599 fAnalysisMode = (AliPWG0Helper::AnalysisMode) (dynamic_cast<TParameter<Int_t>*> (gFile->Get("fAnalysisMode")))->GetVal();
601 gDirectory->cd("../");
604 //____________________________________________________________________
605 void dNdEtaAnalysis::DrawHistograms(Bool_t simple)
607 // draws the histograms
612 fData->DrawHistograms(GetName());
614 TCanvas* canvas = new TCanvas(Form("%s_dNdEtaAnalysis", GetName()), Form("%s_dNdEtaAnalysis", GetName()), 800, 400);
615 canvas->Divide(2, 1);
618 if (fdNdEtaPtCutOffCorrected[0])
619 fdNdEtaPtCutOffCorrected[0]->DrawCopy();
623 fdNdEta[0]->SetLineColor(kRed);
624 fdNdEta[0]->DrawCopy("SAME");
632 // histograms for different vertices?
633 if (kVertexBinning > 0)
635 // doesnt work, but i dont get it, giving up...
636 TCanvas* canvas2 = new TCanvas(Form("%s_dNdEtaAnalysisVtx", GetName()), Form("%s_dNdEtaAnalysisVtx", GetName()), 450, 450);
637 TCanvas* canvas3 = 0;
639 canvas3 = new TCanvas(Form("%s_dNdEtaAnalysisVtx_noptcutoff", GetName()), Form("%s_dNdEtaAnalysisVtx_noptcutoff", GetName()), 450, 450);
641 //Int_t yPads = (Int_t) TMath::Ceil(((Double_t) kVertexBinning - 1) / 2);
642 //printf("%d\n", yPads);
643 //canvas2->Divide(2, yPads);
645 TLegend* legend = new TLegend(0.4, 0.2, 0.6, 0.4);
647 for (Int_t i=0; i<kVertexBinning; ++i)
649 if (fdNdEtaPtCutOffCorrected[i])
653 fdNdEtaPtCutOffCorrected[i]->SetLineColor(i+1);
654 fdNdEtaPtCutOffCorrected[i]->DrawCopy((i == 0) ? "" : "SAME");
655 legend->AddEntry(fdNdEtaPtCutOffCorrected[i], (i == 0) ? "Vtx All" : Form("Vtx Bin %d", i-1));
657 if (canvas3 && fdNdEta[i])
661 fdNdEta[i]->SetLineColor(i+1);
662 fdNdEta[i]->DrawCopy((i == 0) ? "" : "SAME");
668 canvas2->SaveAs(Form("%s_%s.gif", canvas2->GetName(), GetName()));
677 if (kVertexBinning == 3)
679 TH1* clone = dynamic_cast<TH1*> (fdNdEtaPtCutOffCorrected[1]->Clone("clone"));
680 TH1* clone2 = dynamic_cast<TH1*> (fdNdEtaPtCutOffCorrected[2]->Clone("clone2"));
684 TCanvas* canvas4 = new TCanvas(Form("%s_dNdEtaAnalysisVtxRatios", GetName()), Form("%s_dNdEtaAnalysisVtxRatios", GetName()), 450, 450);
686 clone->Divide(fdNdEtaPtCutOffCorrected[0]);
687 clone->GetYaxis()->SetRangeUser(0.95, 1.05);
690 clone2->Divide(fdNdEtaPtCutOffCorrected[0]);
691 clone2->DrawCopy("SAME");
693 TLine* line = new TLine(-1, 1, 1, 1);
696 canvas4->SaveAs(Form("%s_%s.gif", canvas4->GetName(), GetName()));
701 Long64_t dNdEtaAnalysis::Merge(TCollection* list)
703 // Merges a list of dNdEtaAnalysis objects with this one.
704 // This is needed for PROOF.
705 // Returns the number of merged objects (including this)
713 TIterator* iter = list->MakeIterator();
717 const Int_t nCollections = 2 * kVertexBinning + 3; // 3 standalone hists, 3 arrays of size kVertexBinning
718 TList* collections[nCollections];
719 for (Int_t i=0; i<nCollections; ++i)
720 collections[i] = new TList;
723 while ((obj = iter->Next()))
725 dNdEtaAnalysis* entry = dynamic_cast<dNdEtaAnalysis*> (obj);
729 collections[0]->Add(entry->fData);
730 collections[1]->Add(entry->fMult);
731 collections[2]->Add(entry->fPtDist);
733 for (Int_t i=0; i<kVertexBinning; ++i)
735 collections[3+i]->Add(entry->fdNdEta[i]);
736 collections[3+kVertexBinning+i]->Add(entry->fdNdEtaPtCutOffCorrected[i]);
742 fData->Merge(collections[0]);
743 fMult->Merge(collections[1]);
744 fPtDist->Merge(collections[2]);
746 for (Int_t i=0; i<kVertexBinning; ++i)
748 fdNdEta[i]->Merge(collections[3+i]);
749 fdNdEtaPtCutOffCorrected[i]->Merge(collections[3+kVertexBinning+i]);
752 for (Int_t i=0; i<nCollections; ++i)
753 delete collections[i];