//____________________________________________________________________
dNdEtaAnalysis::dNdEtaAnalysis(Char_t* name, Char_t* title) :
-TNamed(name, title)
+ TNamed(name, title),
+ fEtaVsVtx(0),
+ fEtaVsVtxUncorrected(0),
+ fVtx(0)
{
- hEtaVsVtx = new TH2F(Form("%s_eta_vs_vtx", name),"",80,-20,20,120,-6,6);
- hEtaVsVtx->SetXTitle("vtx z [cm]");
- hEtaVsVtx->SetYTitle("#eta");
+ // constructor
- hEtaVsVtxUncorrected = dynamic_cast<TH2F*> (hEtaVsVtx->Clone(Form("%s_eta_vs_vtx_uncorrected", name)));
- hEtaVsVtxCheck = dynamic_cast<TH2F*> (hEtaVsVtx->Clone(Form("%s_eta_vs_vtx_check", name)));
- hVtx = hEtaVsVtx->ProjectionX(Form("%s_vtx", name));
+ fEtaVsVtx = new TH2F(Form("%s_eta_vs_vtx", name),"",80,-20,20,120,-6,6);
+ fEtaVsVtx->SetXTitle("vtx z [cm]");
+ fEtaVsVtx->SetYTitle("#eta");
+
+ fEtaVsVtxUncorrected = dynamic_cast<TH2F*> (fEtaVsVtx->Clone(Form("%s_eta_vs_vtx_uncorrected", name)));
+ fVtx = fEtaVsVtx->ProjectionX(Form("%s_vtx", name));
for (Int_t i=0; i<kVertexBinning; ++i)
{
- hdNdEta[i] = hEtaVsVtx->ProjectionY(Form("%s_dNdEta_%d", name, i));
- hdNdEta[i]->SetYTitle("dN/d#eta");
+ fdNdEta[i] = fEtaVsVtx->ProjectionY(Form("%s_dNdEta_%d", name, i));
+ fdNdEta[i]->SetYTitle("dN/d#eta");
}
- hEtaVsVtx->Sumw2();
- hVtx->Sumw2();
+ fEtaVsVtx->Sumw2();
+ fVtx->Sumw2();
}
//____________________________________________________________________
-void
-dNdEtaAnalysis::FillTrack(Float_t vtx, Float_t eta, Float_t c) {
- hEtaVsVtxUncorrected->Fill(vtx,eta);
- hEtaVsVtxCheck->Fill(vtx, eta, c);
+dNdEtaAnalysis::~dNdEtaAnalysis()
+{
+ // destructor
+
+ delete fEtaVsVtx;
+ fEtaVsVtx = 0;
+
+ delete fEtaVsVtxUncorrected;
+ fEtaVsVtxUncorrected = 0;
+
+ delete fVtx;
+ fVtx = 0;
+
+ for (Int_t i=0; i<kVertexBinning; ++i)
+ {
+ delete fdNdEta[i];
+ fdNdEta[i] = 0;
+ }
+}
+
+//_____________________________________________________________________________
+dNdEtaAnalysis::dNdEtaAnalysis(const dNdEtaAnalysis &c) :
+ TNamed(c),
+ fEtaVsVtx(0),
+ fEtaVsVtxUncorrected(0),
+ fVtx(0)
+{
+ //
+ // dNdEtaAnalysis copy constructor
+ //
+
+ ((dNdEtaAnalysis &) c).Copy(*this);
+}
+
+//_____________________________________________________________________________
+dNdEtaAnalysis &dNdEtaAnalysis::operator=(const dNdEtaAnalysis &c)
+{
+ //
+ // Assignment operator
+ //
+
+ if (this != &c) ((dNdEtaAnalysis &) c).Copy(*this);
+ return *this;
+}
+
+//_____________________________________________________________________________
+void dNdEtaAnalysis::Copy(TObject &c) const
+{
+ //
+ // Copy function
+ //
+
+ dNdEtaAnalysis& target = (dNdEtaAnalysis &) c;
+
+ target.fEtaVsVtx = dynamic_cast<TH2F*> (fEtaVsVtx->Clone());
+ target.fEtaVsVtxUncorrected = dynamic_cast<TH2F*> (fEtaVsVtxUncorrected->Clone());
+ target.fVtx = dynamic_cast<TH1D*> (fVtx->Clone());
+
+ for (Int_t i=0; i<kVertexBinning; ++i)
+ target.fdNdEta[i] = dynamic_cast<TH1D*> (fdNdEta[i]->Clone());
+
+ TNamed::Copy((TNamed &) c);
}
//____________________________________________________________________
-void
-dNdEtaAnalysis::FillEvent(Float_t vtx) {
- hVtx->Fill(vtx);
+void dNdEtaAnalysis::FillTrack(Float_t vtx, Float_t eta)
+{
+ // fills a track into the histograms
+
+ fEtaVsVtxUncorrected->Fill(vtx,eta);
+}
+
+//____________________________________________________________________
+void dNdEtaAnalysis::FillEvent(Float_t vtx)
+{
+ // fills an event into the histograms
+
+ fVtx->Fill(vtx);
}
//____________________________________________________________________
void dNdEtaAnalysis::Finish(dNdEtaCorrection* correction)
{
// correct with correction values if available
+
// TODO what do we do with the error?
if (!correction)
printf("INFO: No correction applied\n");
// this can be replaced by TH2F::Divide if we agree that the binning will be always the same
- for (Int_t iVtx=0; iVtx<=hEtaVsVtxUncorrected->GetNbinsX(); iVtx++)
+ for (Int_t iVtx=0; iVtx<=fEtaVsVtxUncorrected->GetNbinsX(); iVtx++)
{
- for (Int_t iEta=0; iEta<=hEtaVsVtxUncorrected->GetNbinsY(); iEta++)
+ for (Int_t iEta=0; iEta<=fEtaVsVtxUncorrected->GetNbinsY(); iEta++)
{
Float_t correctionValue = 1;
if (correction)
- correctionValue = correction->GetCorrection(hEtaVsVtxUncorrected->GetXaxis()->GetBinCenter(iVtx), hEtaVsVtxUncorrected->GetYaxis()->GetBinCenter(iEta));
+ correctionValue = correction->GetCorrection(fEtaVsVtxUncorrected->GetXaxis()->GetBinCenter(iVtx), fEtaVsVtxUncorrected->GetYaxis()->GetBinCenter(iEta));
- Float_t value = hEtaVsVtxUncorrected->GetBinContent(iVtx, iEta);
- Float_t error = hEtaVsVtxUncorrected->GetBinError(iVtx, iEta);
+ Float_t value = fEtaVsVtxUncorrected->GetBinContent(iVtx, iEta);
+ Float_t error = fEtaVsVtxUncorrected->GetBinError(iVtx, iEta);
Float_t correctedValue = value * correctionValue;
Float_t correctedError = error * correctionValue;
if (correctedValue != 0)
{
- hEtaVsVtx->SetBinContent(iVtx, iEta, correctedValue);
- hEtaVsVtx->SetBinError(iVtx, iEta, correctedError);
+ fEtaVsVtx->SetBinContent(iVtx, iEta, correctedValue);
+ fEtaVsVtx->SetBinError(iVtx, iEta, correctedError);
}
}
}
- // normalize with n events (per vtx)
- for (Int_t iVtx=0; iVtx<=hVtx->GetNbinsX(); iVtx++) {
- Float_t nEvents = hVtx->GetBinContent(iVtx);
- Float_t nEventsError = hVtx->GetBinError(iVtx);
-
- if (nEvents==0) continue;
-
- for (Int_t iEta=0; iEta<=hEtaVsVtx->GetNbinsY(); iEta++) {
- Float_t value = hEtaVsVtx->GetBinContent(iVtx, iEta) / nEvents;
- if (value==0) continue;
- Float_t error = hEtaVsVtx->GetBinError(iVtx, iEta)/nEvents;
- error = TMath::Sqrt(TMath::Power(hEtaVsVtx->GetBinError(iVtx, iEta)/
- hEtaVsVtx->GetBinContent(iVtx, iEta),2) +
- TMath::Power(nEventsError/nEvents,2));
- hEtaVsVtx->SetBinContent(iVtx, iEta, value);
- hEtaVsVtx->SetBinError(iVtx, iEta, error);
- }
-
- //debug
- for (Int_t iEta=0; iEta<=hEtaVsVtxCheck->GetNbinsY(); iEta++) {
- Float_t value = hEtaVsVtxCheck->GetBinContent(iVtx, iEta) / nEvents;
- if (value==0) continue;
- Float_t error = hEtaVsVtxCheck->GetBinError(iVtx, iEta)/nEvents;
- error = TMath::Sqrt(TMath::Power(hEtaVsVtxCheck->GetBinError(iVtx, iEta)/
- hEtaVsVtxCheck->GetBinContent(iVtx, iEta),2) +
- TMath::Power(nEventsError/nEvents,2));
- 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<=fEtaVsVtx->GetNbinsY(); iEta++)
+ {
// do we have several histograms for different vertex positions?
- Int_t vertexBinWidth = hVtx->GetNbinsX() / (kVertexBinning-1);
+ Int_t vertexBinWidth = fVtx->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 = fVtx->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 (hEtaVsVtx->GetBinContent(iVtx, iEta)==0) continue;
-
- Float_t w = 1/TMath::Power(hEtaVsVtx->GetBinError(iVtx, iEta),2);
- sumXw = sumXw + hEtaVsVtx->GetBinContent(iVtx, iEta)*w;
- sumW = sumW + w;
-
- sum = sum + hEtaVsVtx->GetBinContent(iVtx, iEta);
- sumError2 = sumError2 + TMath::Power(hEtaVsVtx->GetBinError(iVtx, iEta),2);
- nMeasurements++;
+ Float_t totalEvents = fVtx->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 (fEtaVsVtx->GetBinContent(iVtx, iEta) != 0)
+ {
+ sum = sum + fEtaVsVtx->GetBinContent(iVtx, iEta);
+ sumError2 = sumError2 + TMath::Power(fEtaVsVtx->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/fdNdEta[vertexPos]->GetBinWidth(iEta);
+ error = error/fdNdEta[vertexPos]->GetBinWidth(iEta);
- hdNdEta[vertexPos]->SetBinContent(iEta, dndeta);
- hdNdEta[vertexPos]->SetBinError(iEta, error);
- }
+ fdNdEta[vertexPos]->SetBinContent(iEta, dndeta);
+ fdNdEta[vertexPos]->SetBinError(iEta, error);
}
}
}
//____________________________________________________________________
-void
-dNdEtaAnalysis::SaveHistograms() {
+void dNdEtaAnalysis::SaveHistograms()
+{
+ // save the histograms to a directory with the name of this class (retrieved from TNamed)
gDirectory->mkdir(GetName());
gDirectory->cd(GetName());
-
- hEtaVsVtx ->Write();
- hEtaVsVtxUncorrected->Write();
- hVtx ->Write();
+
+ fEtaVsVtx ->Write();
+ fEtaVsVtxUncorrected->Write();
+ fVtx ->Write();
for (Int_t i=0; i<kVertexBinning; ++i)
- hdNdEta[i] ->Write();
+ fdNdEta[i] ->Write();
+
+ gDirectory->cd("../");
+}
+
+void dNdEtaAnalysis::LoadHistograms()
+{
+ // loads the histograms from a directory with the name of this class (retrieved from TNamed)
+
+ gDirectory->cd(GetName());
+
+ fEtaVsVtx = dynamic_cast<TH2F*> (gDirectory->Get(fEtaVsVtx->GetName()));
+ fEtaVsVtxUncorrected = dynamic_cast<TH2F*> (gDirectory->Get(fEtaVsVtxUncorrected->GetName()));
+
+ fVtx = dynamic_cast<TH1D*> (gDirectory->Get(fVtx->GetName()));
+
+ for (Int_t i=0; i<kVertexBinning; ++i)
+ fdNdEta[i] = dynamic_cast<TH1D*> (gDirectory->Get(fdNdEta[i]->GetName()));
gDirectory->cd("../");
}
//____________________________________________________________________
void dNdEtaAnalysis::DrawHistograms()
{
+ // draws the histograms
+
TCanvas* canvas = new TCanvas("dNdEtaAnalysis", "dNdEtaAnalysis", 800, 800);
canvas->Divide(2, 2);
canvas->cd(1);
- if (hEtaVsVtx)
- hEtaVsVtx->Draw("COLZ");
+ if (fEtaVsVtx)
+ fEtaVsVtx->Draw("COLZ");
canvas->cd(2);
- if (hEtaVsVtxCheck)
- hEtaVsVtxCheck->Draw("COLZ");
+ if (fEtaVsVtxUncorrected)
+ fEtaVsVtxUncorrected->Draw("COLZ");
canvas->cd(3);
- if (hEtaVsVtxUncorrected)
- hEtaVsVtxUncorrected->Draw("COLZ");
-
- /*canvas->cd(3);
- if (hVtx)
- hVtx->Draw();*/
+ if (fVtx)
+ fVtx->Draw();
canvas->cd(4);
- if (hdNdEta[0])
- hdNdEta[0]->Draw();
+ if (fdNdEta[0])
+ fdNdEta[0]->Draw();
// histograms for different vertices?
if (kVertexBinning > 0)
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])
+ if (fdNdEta[i])
{
- hdNdEta[i]->SetLineColor(i);
- hdNdEta[i]->Draw((i == 1) ? "" : "SAME");
- legend->AddEntry(hdNdEta[i], Form("Vtx Bin %d", i-1));
+ fdNdEta[i]->SetLineColor(i+1);
+ fdNdEta[i]->Draw((i == 0) ? "" : "SAME");
+ legend->AddEntry(fdNdEta[i], (i == 0) ? "Vtx All" : Form("Vtx Bin %d", i-1));
}
}
if (entry == 0)
continue;
- collections[0]->Add(entry->hEtaVsVtx);
- collections[1]->Add(entry->hEtaVsVtxUncorrected);
- collections[2]->Add(entry->hVtx);
+ collections[0]->Add(entry->fEtaVsVtx);
+ collections[1]->Add(entry->fEtaVsVtxUncorrected);
+ collections[2]->Add(entry->fVtx);
for (Int_t i=0; i<kVertexBinning; ++i)
- collections[3+i]->Add(entry->hdNdEta[i]);
+ collections[3+i]->Add(entry->fdNdEta[i]);
++count;
}
- hEtaVsVtx->Merge(collections[0]);
- hEtaVsVtxUncorrected->Merge(collections[1]);
- hVtx->Merge(collections[2]);
+ fEtaVsVtx->Merge(collections[0]);
+ fEtaVsVtxUncorrected->Merge(collections[1]);
+ fVtx->Merge(collections[2]);
for (Int_t i=0; i<kVertexBinning; ++i)
- hdNdEta[i]->Merge(collections[3+i]);
+ fdNdEta[i]->Merge(collections[3+i]);
for (Int_t i=0; i<nCollections; ++i)
delete collections[i];