adding multiplicity correction (simple gaussian method)

adding systematic uncertainty study of sigma2vertex

diff --git a/PWG0/dNdEta/drawPlots.C b/PWG0/dNdEta/drawPlots.C
index 40b43e061a5922d523b26c7ab88308ccfa0a9da3..8902f041805f4bf6efe7139dae442a8680baac76 100644 (file)
--- a/PWG0/dNdEta/drawPlots.C
+++ b/PWG0/dNdEta/drawPlots.C
@@ -701,6 +701,117 @@ void Track2Particle3DAll()
   canvas->SaveAs("Track2Particle3DAll.eps");
 }
 
+void MultiplicityMC()
+{
+  TFile* file = TFile::Open("multiplicityMC.root");
+
+  if (!file)
+  {
+    printf("multiplicityMC.root could not be opened.\n");
+    return;
+  }
+
+  TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicityESD"));
+  TH1F* fMultiplicityMC = dynamic_cast<TH1F*> (file->Get("fMultiplicityMC"));
+  TH2F* fCorrelation = dynamic_cast<TH2F*> (file->Get("fCorrelation"));
+
+  TH1F* correction = new TH1F("MultiplicityMC_correction", "MultiplicityMC_correction;Ntracks;Npart", 76, -0.5, 75.5);
+  TH1F* correctionWidth = new TH1F("MultiplicityMC_correctionwidth", "MultiplicityMC_correctionwidth;Ntracks;Npart", 76, -0.5, 75.5);
+  //fMultiplicityMC->GetNbinsX(), fMultiplicityMC->GetXaxis()->GetXmin(), fMultiplicityMC->GetXaxis()->GetXmax());
+  for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
+  {
+    TH1D* proj = fCorrelation->ProjectionX("_px", i, i+1);
+    proj->Fit("gaus", "0");
+    correction->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(1));
+    correctionWidth->SetBinContent(i, proj->GetFunction("gaus")->GetParameter(2));
+
+    continue;
+
+    // draw for debugging
+    new TCanvas;
+    proj->DrawCopy();
+    proj->GetFunction("gaus")->DrawCopy("SAME");
+  }
+
+  TH1F* fMultiplicityESDCorrected = new TH1F("fMultiplicityESDCorrected", "fMultiplicityESDCorrected", 2010, -0.5, 200.5);
+
+  for (Int_t i=1; i<=correction->GetNbinsX(); ++i)
+  {
+    Float_t mean = correction->GetBinContent(i);
+    Float_t width = correctionWidth->GetBinContent(i);
+
+    Int_t fillBegin = fMultiplicityESDCorrected->FindBin(mean - width * 3);
+    Int_t fillEnd   = fMultiplicityESDCorrected->FindBin(mean + width * 3);
+    printf("bin %d mean %f width %f, filling from %d to %d\n", i, mean, width, fillBegin, fillEnd);
+
+    for (Int_t j=fillBegin; j <= fillEnd; ++j)
+    {
+      fMultiplicityESDCorrected->AddBinContent(j, TMath::Gaus(fMultiplicityESDCorrected->GetXaxis()->GetBinCenter(j), mean, width, kTRUE) * fMultiplicityESD->GetBinContent(i));
+    }
+  }
+
+  TH1F* fMultiplicityESDCorrectedRebinned = dynamic_cast<TH1F*> (fMultiplicityESDCorrected->Clone("fMultiplicityESDCorrectedRebinned"));
+  fMultiplicityESDCorrectedRebinned->Rebin(10);
+  fMultiplicityESDCorrectedRebinned->Scale(0.1);
+
+  TH1F* ratio2 = dynamic_cast<TH1F*> (fMultiplicityESDCorrectedRebinned->Clone("multiplicity_ratio_corrected"));
+  ratio2->Divide(fMultiplicityMC);
+
+  TCanvas* canvas = new TCanvas("MultiplicityMC", "MultiplicityMC", 1500, 1000);
+  canvas->Divide(3, 2);
+
+  canvas->cd(1);
+  fMultiplicityESD->Draw();
+  fMultiplicityMC->SetLineColor(2);
+  fMultiplicityMC->Draw("SAME");
+
+  canvas->cd(2);
+  TH1F* ratio = dynamic_cast<TH1F*> (fMultiplicityESD->Clone("multiplicity_ratio"));
+  ratio->SetTitle("ratio;Ntracks;Nreco/Ngene");
+  ratio->Divide(fMultiplicityMC);
+  ratio->Draw();
+
+  ratio2->SetLineColor(2);
+  ratio2->Draw("SAME");
+
+  canvas->cd(3);
+  fCorrelation->Draw("COLZ");
+
+  canvas->cd(4);
+  correction->Draw();
+  //correction->Fit("pol1");
+  correctionWidth->SetLineColor(2);
+  correctionWidth->Draw("SAME");
+
+  canvas->cd(5);
+  fMultiplicityESDCorrected->SetLineColor(3);
+  fMultiplicityESDCorrected->Draw();
+  fMultiplicityMC->Draw("SAME");
+  fMultiplicityESD->Draw("SAME");
+
+  canvas->cd(6);
+  fMultiplicityESDCorrectedRebinned->SetLineColor(3);
+  fMultiplicityESDCorrectedRebinned->Draw();
+  fMultiplicityMC->Draw("SAME");
+}
+
+void MultiplicityESD()
+{
+  TFile* file = TFile::Open("multiplicityESD.root");
+
+  if (!file)
+  {
+    printf("multiplicityESD.root could not be opened.\n");
+    return;
+  }
+
+  TH1F* fMultiplicityESD = dynamic_cast<TH1F*> (file->Get("fMultiplicity"));
+
+  TCanvas* canvas = new TCanvas("MultiplicityESD", "MultiplicityESD", 500, 500);
+
+  fMultiplicityESD->Draw();
+}
+
 void drawPlots(Int_t max)
 {
   gMax = max;

diff --git a/PWG0/dNdEta/drawSystematics.C b/PWG0/dNdEta/drawSystematics.C
index 7d47d7ce3a32894ccd91b9dc3493dc93db44179d..9234b82d5dd1c1c85641b62fa36bf21d37fbd5ce 100644 (file)
--- a/PWG0/dNdEta/drawSystematics.C
+++ b/PWG0/dNdEta/drawSystematics.C
@@ -361,7 +361,7 @@ TH2F* Sigma2VertexGaussianTracksHist()
   return tracks;
 }
 
-void Sigma2VertexGaussian()
+TH1F* Sigma2VertexGaussian()
 {
   TH2F* tracks = Sigma2VertexGaussianTracksHist();
 
@@ -377,7 +377,7 @@ void Sigma2VertexGaussian()
   ratio->SetMarkerStyle(21);
 
   canvas->cd(2);
-  ratio->Draw("P");
+  ratio->DrawCopy("P");
 
   TH1F* ratio2 = new TH1F("Sigma2Vertex_ratio2", "Sigma2Vertex_ratio2;nSigma;% included 3 sigma / % included n sigma", 10, 0.25, 5.25);
   Double_t sigma3 = Sigma2VertexCount(tracks, 3);
@@ -386,12 +386,14 @@ void Sigma2VertexGaussian()
   ratio2->SetMarkerStyle(21);
 
   canvas->cd(3);
-  ratio2->Draw("P");
+  ratio2->DrawCopy("P");
 
   canvas->SaveAs("Sigma2Vertex.eps");
+
+  return ratio2;
 }
 
-void Sigma2VertexSimulation()
+TH1F* Sigma2VertexSimulation()
 {
   TFile* file = TFile::Open("systematics.root");
 
@@ -402,7 +404,7 @@ void Sigma2VertexSimulation()
     return;
   }
 
-  TH1F* ratio = new TH1F("sigmavertex_ratio", "sigmavertex_ratio;Nsigma;% included 3 sigma / % included n sigma", sigmavertex->GetNbinsX(), sigmavertex->GetXaxis()->GetXmin(), sigmavertex->GetXaxis()->GetXmax());
+  TH1F* ratio = new TH1F("sigmavertexsimulation_ratio", "sigmavertexsimulation_ratio;Nsigma;% included 3 sigma / % included n sigma", sigmavertex->GetNbinsX(), sigmavertex->GetXaxis()->GetXmin(), sigmavertex->GetXaxis()->GetXmax());
 
   for (Int_t i=1; i<=sigmavertex->GetNbinsX(); ++i)
     ratio->SetBinContent(i, sigmavertex->GetBinContent(sigmavertex->GetXaxis()->FindBin(3)) / sigmavertex->GetBinContent(i));
@@ -416,9 +418,37 @@ void Sigma2VertexSimulation()
 
   canvas->cd(2);
   ratio->SetMarkerStyle(21);
-  ratio->Draw("P");
+  ratio->DrawCopy("P");
+
+  return ratio;
 }
 
+void Sigma2VertexCompare()
+{
+  TH1F* ratio1 = Sigma2VertexGaussian();
+  TH1F* ratio2 = Sigma2VertexSimulation();
+
+  ratio1->SetStats(kFALSE);
+  ratio2->SetStats(kFALSE);
+
+  ratio1->SetMarkerStyle(0);
+  ratio2->SetMarkerStyle(0);
+
+  TLegend* legend = new TLegend(0.647177,0.775424,0.961694,0.966102);
+  legend->AddEntry(ratio1, "Gaussian");
+  legend->AddEntry(ratio2, "Simulation");
+
+  ratio1->GetXaxis()->SetTitleOffset(1.5);
+
+  TCanvas* canvas = new TCanvas("Sigma2VertexCompare", "Sigma2VertexCompare", 500, 500);
+  InitPad();
+
+  ratio1->Draw();
+  ratio2->SetLineColor(kRed);
+  ratio2->Draw("SAME");
+
+  legend->Draw();
+}
 
 void drawSystematics()
 {