analysis_code/final_asymmetries/CreateVtxPlot.C

0001 #include "/sphenix/u/virgilemahaut/style/sPhenixStyle_Greg.C"
0002
0003 #include <TCanvas.h>
0004 #include <TGraphErrors.h>
0005 #include <TLegend.h>
0006 #include <TAxis.h>
0007 #include <TStyle.h>
0008
0009 #include <fstream>
0010 #include <sstream>
0011 #include <string>
0012 #include <vector>
0013 #include <algorithm>
0014 #include <iostream>
0015
0016 double GetAbsMaxInRange(TGraphErrors *graph, double xMin, double xMax)
0017 {
0018   if (!graph ) return -1; // Error handling
0019
0020   double maxVal = -1;
0021   int n = graph->GetN();
0022
0023   for (int i = 0; i < n; i++)
0024   {
0025     double x, y, ey, yVal;
0026     graph->GetPoint(i, x, y);
0027     ey = graph->GetErrorY(i);
0028     yVal = std::max(std::abs(y-ey), std::abs(y+ey));
0029     if (x >= xMin && x <= xMax)
0030     {
0031       if (yVal > maxVal)
0032       {
0033         maxVal = yVal;
0034       }
0035     }
0036   }
0037   return maxVal;
0038 }
0039
0040
0041 void parse_file(const std::string& filename,
0042                 std::vector<double>& x,
0043                 std::vector<double>& y,
0044                 std::vector<double>& estat,
0045                 std::vector<double>& esyst,
0046                 std::vector<double>& exsyst,
0047                 const double width = 0.15
0048                 )
0049 {
0050   std::ifstream fin(filename);
0051   if (!fin.is_open()) {
0052     std::cerr << "Cannot open file: " << filename << std::endl;
0053     return;
0054   }
0055
0056   std::string line;
0057   while (std::getline(fin, line)) {
0058     // Skip empty lines
0059     if (line.empty()) continue;
0060
0061     // Replace '&' with spaces so stringstream can parse numbers
0062     std::replace(line.begin(), line.end(), '&', ' ');
0063
0064     std::stringstream ss(line);
0065     double xv, yv, statv, systv;
0066     if (!(ss >> xv >> yv >> statv >> systv)) {
0067       std::cerr << "Skipping malformed line: " << line << std::endl;
0068       continue;
0069     }
0070
0071     x.push_back(xv);
0072     y.push_back(yv);
0073     estat.push_back(statv);
0074     esyst.push_back(systv);
0075
0076     // Constant half-width for the systematic boxes in x.
0077     // Choose something visually reasonable for your spacing.
0078     exsyst.push_back(width);
0079   }
0080 }
0081
0082 void CreateVtxPlot()
0083 {
0084   std::string plots_folder_xf = "xf_asymmetries_low_vtx/SYSTEMATIC/";
0085
0086   gSystem->Exec(("mkdir -p " + plots_folder_xf).c_str());
0087
0088
0089   const std::string inputfolder = "asym_values/";
0090   const int nParticles = 2;
0091   const std::string particles[nParticles] = {"pi0", "eta"};
0092   const int nDirections = 2;
0093   const std::string directions[nDirections] = {"forward", "backward"};
0094   for (int iP = 0; iP < nParticles; iP++) {
0095     {
0096       std::stringstream inputfilename_low;
0097       inputfilename_low << inputfolder << "asym_summary_" << particles[iP] << "_xf_low_vtx" << ".txt";
0098       std::stringstream inputfilename_high;
0099       inputfilename_high << inputfolder << "asym_summary_" << particles[iP] << "_xf_high_vtx" << ".txt";
0100       std::stringstream canvas_name;
0101       canvas_name << "canvas_asym_" << particles[iP] << "_xf_vs_vertex";
0102       std::vector<double> x_low, y_low, estat_low, esyst_low, exsyst_low;
0103       std::vector<double> x_high, y_high, estat_high, esyst_high, exsyst_high;
0104
0105       parse_file(inputfilename_low.str(), x_low, y_low, estat_low, esyst_low, exsyst_low, 0.002);
0106       parse_file(inputfilename_high.str(), x_high, y_high, estat_high, esyst_high, exsyst_high, 0.002);
0107
0108       const int n = x_low.size();
0109       auto gSyst_low = new TGraphErrors(n);
0110       auto gStat_low = new TGraphErrors(n);
0111       auto gSyst_high = new TGraphErrors(n);
0112       auto gStat_high = new TGraphErrors(n);
0113
0114       for (int i = 0; i < n; ++i) {
0115         gSyst_low->SetPoint(i, x_low[i], y_low[i]);
0116         gSyst_low->SetPointError(i, exsyst_low[i], esyst_low[i]);
0117
0118         gStat_low->SetPoint(i, x_low[i], y_low[i]);
0119         gStat_low->SetPointError(i, 0.0, estat_low[i]);
0120
0121         gSyst_high->SetPoint(i, x_high[i], y_high[i]);
0122         gSyst_high->SetPointError(i, exsyst_high[i], esyst_high[i]);
0123
0124         gStat_high->SetPoint(i, x_high[i], y_high[i]);
0125         gStat_high->SetPointError(i, 0.0, estat_high[i]);
0126       }
0127       SetsPhenixStyle();
0128
0129       TCanvas *canvas = new TCanvas(canvas_name.str().c_str(), "", 1600, 900);
0130       gStat_low->SetTitle(";x_{F};A_{N}");
0131
0132       double y_bound = 0;
0133       y_bound = GetAbsMaxInRange(gStat_low, -0.3, 0.3);
0134
0135       int color = (iP == 0 ? kRed : kBlue);
0136
0137       gStat_low->SetLineWidth(2);
0138       gStat_low->SetLineColor(color);
0139       gStat_low->SetMarkerColor(color);
0140       gStat_low->SetMarkerStyle(kFullCircle);
0141       gStat_low->SetMarkerSize(2.2);
0142       gStat_low->SetMinimum(-2 * y_bound);
0143       gStat_low->SetMaximum(2 * y_bound);
0144       gStat_low->Draw("AP E1");
0145
0146       gStat_high->SetLineWidth(2);
0147       gStat_high->SetLineColor(kViolet);
0148       gStat_high->SetMarkerColor(kViolet);
0149       gStat_high->SetMarkerStyle(kFullCircle);
0150       gStat_high->SetMarkerSize(2.2);
0151       gStat_high->SetMinimum(-2 * y_bound);
0152       gStat_high->SetMaximum(2 * y_bound);
0153       gStat_high->Draw("P E1 SAME");
0154
0155        // Style systematic boxes (transparent fill)
0156       if (iP == 0) {
0157         gSyst_low->SetFillColorAlpha(kPink-9, 0.35); // alpha controls transparency
0158         gSyst_low->SetLineColor(kPink-9);
0159         gSyst_high->SetFillColorAlpha(kViolet-9, 0.35); // alpha controls transparency
0160         gSyst_high->SetLineColor(kViolet-9);
0161       } else {
0162         gSyst_low->SetFillColorAlpha(kAzure-4, 0.35); // alpha controls transparency
0163         gSyst_low->SetLineColor(kAzure-4);
0164       }
0165       gSyst_low->SetMarkerStyle(0);
0166       gSyst_low->Draw("P E2 SAME");
0167       gSyst_high->SetMarkerStyle(0);
0168       gSyst_high->Draw("P E2 SAME");
0169
0170       TLegend *legend = new TLegend(0.25, 0.25, 0.5, 0.4);
0171       legend->SetTextSize(0.05);
0172       legend->AddEntry(gStat_low, "|z_{vtx}| < 30 cm");
0173       legend->AddEntry(gStat_high, "|z_{vtx}| > 30 cm");
0174       legend->Draw();
0175
0176       gPad->Modified();
0177       gPad->Update();
0178       double min_x = gPad->GetUxmin();
0179       double max_x = gPad->GetUxmax();
0180
0181       // Add "sPHENIX internal"
0182       TPad *p = new TPad("p","p",0.,0.,1.,1.); p->SetFillStyle(0); p->Draw(); p->cd();
0183       TBox *whiteBox = new TBox(0.17, 0.72, 0.46, 0.90);
0184       whiteBox->Draw();
0185       canvas->cd();
0186       whiteBox->SetFillColorAlpha(kWhite, 1);
0187       std::stringstream stream;
0188       stream.str("");
0189       TLatex latex;
0190       latex.SetNDC();
0191       latex.SetTextColor(kBlack);
0192       latex.DrawLatex(0.22, 0.85, "#font[72]{sPHENIX} Internal");
0193       latex.DrawLatex(0.22, 0.75, "p^{#uparrow}+p #sqrt{s} = 200 GeV");
0194       latex.SetTextSize(0.03);
0195       latex.DrawLatex(0.19, 0.68, "7% polarization scale uncertainty not shown");
0196       latex.SetTextSize(0.05);
0197
0198       if (iP == 0) {
0199         latex.DrawLatex(0.5, 0.85, "p^{#uparrow}+p #rightarrow #pi^{0} X");
0200       } else {
0201         latex.DrawLatex(0.5, 0.85, "p^{#uparrow}+p #rightarrow #eta X");
0202       }
0203
0204       TLine *tline = new TLine();
0205       tline->SetLineWidth(2);
0206       tline->SetLineColor(kBlack);
0207       tline->SetLineStyle(kDashed);
0208       tline->DrawLine(min_x, 0, max_x, 0);
0209
0210       canvas->Update();
0211       canvas->Draw();
0212       canvas->SaveAs((plots_folder_xf +  "/" + canvas_name.str() + ".png").c_str());
0213       canvas->SaveAs((plots_folder_xf +  "/" + canvas_name.str() + ".pdf").c_str());
0214     }
0215   }
0216   gSystem->Exit(0);
0217 }