Back to home page

sPhenix code displayed by LXR

 
 

    


File indexing completed on 2025-08-05 08:12:51

0001 #include <fstream>
0002 #include <iostream>
0003 #include "stdio.h"
0004 
0005 #include <TChain.h>
0006 #include <TFile.h>
0007 #include <TGraphAsymmErrors.h>
0008 #include <TGraphErrors.h>
0009 #include <TLatex.h>
0010 #include <TLegend.h>
0011 #include <TLine.h>
0012 #include <TROOT.h>
0013 #include <TMath.h>
0014 #include <TString.h>
0015 #include <TTree.h>
0016 #include <TVectorD.h>
0017 #include <TVirtualFitter.h>
0018 #include <cassert>
0019 #include <cmath>
0020 #include <vector>
0021 
0022 #include "SaveCanvas.C"
0023 #include "sPhenixStyle.C"
0024 
0025 using namespace std;
0026 
0027 Double_t v1_err(Double_t sig, Double_t v1, Double_t Res)
0028 {
0029   const Double_t Pi = 3.1415927;
0030   return 2.3 / sig * sqrt(1 - 16 * v1 * v1 / Pi / Pi) / Res;  // 2.3 - coefficient from simple two rapidity bin linear fit
0031 }
0032 
0033 TGraphErrors *GraphShiftScaling(TGraphErrors *gr_src, const double x_shift, const double err_scaling)
0034 {
0035   assert(gr_src);
0036 
0037   const int npoint = gr_src->GetN();
0038 
0039   TVectorD vx(npoint);
0040   TVectorD vy(npoint);
0041   TVectorD vex(npoint);
0042   TVectorD vey(npoint);
0043 
0044   int nfilled = 0;
0045   for (int i = 0; i < npoint; ++i)
0046   {
0047     const double &x = gr_src->GetX()[i];
0048     //        if (x<x_min or x>x_max) continue;
0049 
0050     vx[nfilled] = x + x_shift;
0051     vy[nfilled] = gr_src->GetY()[i];
0052     vex[nfilled] = gr_src->GetEX()[i];
0053     vey[nfilled] = gr_src->GetEY()[i] * err_scaling;
0054 
0055     ++nfilled;
0056   }
0057 
0058   TGraphErrors *gr = new TGraphErrors(nfilled, vx.GetMatrixArray(), vy.GetMatrixArray(),
0059                                       vex.GetMatrixArray(), vey.GetMatrixArray());
0060 
0061   gr->SetMarkerColor(gr_src->GetMarkerColor());
0062   gr->SetMarkerStyle(gr_src->GetMarkerStyle());
0063   gr->SetMarkerSize(gr_src->GetMarkerSize());
0064   gr->SetLineWidth(gr_src->GetLineWidth());
0065   gr->SetLineColor(gr_src->GetLineColor());
0066 
0067   return gr;
0068 }
0069 
0070 void makeV1_BUP2020(const Bool_t sPHENIX = 1)
0071 {
0072   //  gROOT->LoadMacro("sPhenixStyle.C");
0073   SetsPhenixStyle();
0074 
0075   const Double_t MassD = 1.86484;
0076 
0077   // STAR data points
0078   const Int_t NS = 4;
0079   const Double_t offset = 0.02;
0080   Double_t yp[NS], v1p[NS], v1pe[NS], v1pes[NS];
0081   Double_t yn[NS], v1n[NS], v1ne[NS], v1nes[NS];
0082   ifstream inData;
0083   inData.open("D0v1_STAR_final.txt");
0084   for (int i = 0; i < NS; i++)
0085   {
0086     inData >> yp[i] >> v1p[i] >> v1pe[i] >> v1pes[i];
0087     yp[i] -= offset;
0088   }
0089   inData.close();
0090   inData.open("D0barv1_STAR_final.txt");
0091   for (int i = 0; i < NS; i++)
0092   {
0093     inData >> yn[i] >> v1n[i] >> v1ne[i] >> v1nes[i];
0094     yn[i] += offset;
0095   }
0096   inData.close();
0097 
0098   //  TGraphErrors *gr_D_STAR = new TGraphErrors("D0v1_STAR_final.txt","%lg %lg %lg");
0099   TGraphErrors *gr_D_STAR = new TGraphErrors(NS, yp, v1p, 0, v1pe);
0100   gr_D_STAR->SetMarkerColor(kRed - 6);
0101   gr_D_STAR->SetMarkerStyle(kFullStar);
0102   gr_D_STAR->SetMarkerSize(2);
0103   gr_D_STAR->SetLineWidth(4);
0104   gr_D_STAR->SetLineColor(kRed - 6);
0105 
0106   //  TGraphErrors *gr_Dbar_STAR = new TGraphErrors("D0barv1_STAR_final.txt","%lg %lg %lg");
0107   TGraphErrors *gr_Dbar_STAR = new TGraphErrors(NS, yn, v1n, 0, v1ne);
0108   gr_Dbar_STAR->SetMarkerColor(kRed - 6);
0109   gr_Dbar_STAR->SetMarkerStyle(kOpenStar);
0110   gr_Dbar_STAR->SetMarkerSize(2);
0111   gr_Dbar_STAR->SetLineWidth(4);
0112   gr_Dbar_STAR->SetLineColor(kRed - 6);
0113 
0114   TFile *fin = new TFile("D0_significance.root");
0115   TGraph *gr_sig_D = (TGraph *) fin->Get("gProD0_0_80_noPid");
0116   double sig2_tot = 0.;
0117   for (int i = 0; i < gr_sig_D->GetN(); i++)
0118   {
0119     sig2_tot += TMath::Power(gr_sig_D->GetY()[i], 2.0);
0120   }
0121   double sig_tot = sqrt(sig2_tot);
0122   cout << " Total D/Dbar significance in |y|<1 = " << sig_tot << endl;
0123   double sig_tot_D = sig_tot / sqrt(2.);  // one charge sign
0124 
0125   const Double_t EPRes = 0.37;  // STAR 1st EP resolution
0126 
0127   // theory curves
0128   // Greco - D/Dbar v1 difference
0129   TGraph *gr_D_Greco = new TGraph("D0v1_RHIC_Greco.txt", "%lg %lg");
0130   TGraph *gr_Dbar_Greco = new TGraph("D0barv1_RHIC_Greco.txt", "%lg %lg");           // be careful, this calculation has a sign flip
0131   TGraph *gr_DDbar_Bozek_tmp = new TGraph("DDbarv1_eta_1712.01180.txt", "%lg %lg");  // v1 vs eta, in percentage
0132   double x[100], y[100];
0133   for (int i = 0; i < gr_DDbar_Bozek_tmp->GetN(); i++)
0134   {
0135     double eta = gr_DDbar_Bozek_tmp->GetX()[i];
0136     x[i] = eta;  // questionable!!
0137     y[i] = gr_DDbar_Bozek_tmp->GetY()[i] / 100.;
0138   }
0139   TGraph *gr_DDbar_Bozek = new TGraph(gr_DDbar_Bozek_tmp->GetN(), x, y);
0140   cout << " AAAA " << endl;
0141 
0142   // for projections v1
0143   const Int_t N_D = 8;
0144   double y_D[N_D], v1_D[N_D], v1_Dbar[N_D], v1_err_D[N_D], v1_err_Dbar[N_D];
0145 
0146   for (int i = 0; i < N_D; i++)
0147   {
0148     y_D[i] = (i + 0.5) * 2.0 / N_D - 1.0;
0149     v1_D[i] = -1.0 * gr_D_Greco->Eval(y_D[i]) + gr_DDbar_Bozek->Eval(y_D[i]);  // there is a sign flip in Greco's calculation
0150     v1_Dbar[i] = -1.0 * gr_Dbar_Greco->Eval(y_D[i]) + gr_DDbar_Bozek->Eval(y_D[i]);
0151 
0152     double sig_per_bin = sig_tot_D / sqrt((double) N_D);  // //not combine symmetric rapidity bins
0153     v1_err_D[i] = v1_err(sig_per_bin, v1_D[i], EPRes);
0154     v1_err_Dbar[i] = v1_err(sig_per_bin, v1_Dbar[i], EPRes);
0155     cout << y_D[i] << " " << v1_D[i] << "+/-" << v1_err_D[i] << endl;
0156   }
0157 
0158   TGraph *gr_D = new TGraph(N_D, y_D, v1_D);
0159   gr_D->SetMarkerSize(2);
0160   gr_D->SetMarkerColor(1);
0161   gr_D->SetMarkerStyle(20);
0162   gr_D->SetLineWidth(2.);
0163   gr_D->SetLineStyle(1);
0164   gr_D->SetLineColor(1);
0165 
0166   TGraph *gr_Dbar = new TGraph(N_D, y_D, v1_Dbar);
0167   gr_Dbar->SetMarkerSize(2);
0168   gr_Dbar->SetMarkerColor(4);
0169   gr_Dbar->SetMarkerStyle(24);
0170   gr_Dbar->SetLineWidth(2.);
0171   gr_Dbar->SetLineStyle(1);
0172   gr_Dbar->SetLineColor(4);
0173 
0174   TGraphErrors *gr_D_proj = new TGraphErrors(N_D / 2, y_D + 4, v1_D + 4, 0, v1_err_D + 4);
0175   gr_D_proj->SetMarkerSize(2);
0176   gr_D_proj->SetMarkerColor(1);
0177   gr_D_proj->SetMarkerStyle(20);
0178   gr_D_proj->SetLineWidth(4.);
0179   gr_D_proj->SetLineStyle(1);
0180   gr_D_proj->SetLineColor(1);
0181 
0182   TGraphErrors *gr_D_proj_2 = new TGraphErrors(N_D / 2, y_D, v1_D, 0, v1_err_D);
0183   gr_D_proj_2->SetMarkerSize(2);
0184   gr_D_proj_2->SetMarkerColor(1);
0185   gr_D_proj_2->SetMarkerStyle(20);
0186   gr_D_proj_2->SetLineWidth(4.);
0187   gr_D_proj_2->SetLineStyle(1);
0188   gr_D_proj_2->SetLineColor(1);
0189 
0190   TGraphErrors *gr_Dbar_proj = new TGraphErrors(N_D / 2, y_D + 4, v1_Dbar + 4, 0, v1_err_Dbar + 4);
0191   gr_Dbar_proj->SetMarkerSize(2);
0192   gr_Dbar_proj->SetMarkerColor(4);
0193   gr_Dbar_proj->SetMarkerStyle(20);
0194   gr_Dbar_proj->SetLineWidth(4.);
0195   gr_Dbar_proj->SetLineStyle(1);
0196   gr_Dbar_proj->SetLineColor(4);
0197 
0198   TGraphErrors *gr_Dbar_proj_2 = new TGraphErrors(N_D / 2, y_D, v1_Dbar, 0, v1_err_Dbar);
0199   gr_Dbar_proj_2->SetMarkerSize(2);
0200   gr_Dbar_proj_2->SetMarkerColor(4);
0201   gr_Dbar_proj_2->SetMarkerStyle(20);
0202   gr_Dbar_proj_2->SetLineWidth(4.);
0203   gr_Dbar_proj_2->SetLineStyle(1);
0204   gr_Dbar_proj_2->SetLineColor(4);
0205 
0206   ///////////////////////////////////////////////////////////////////////////
0207   //  BUP2020  lumi scaling
0208   ///////////////////////////////////////////////////////////////////////////
0209 
0210   const double refAuAuMB = 240e9;
0211   const double refAuAuXSec = 6.8252;  // b
0212 
0213   const double AuAu_Ncoll_C0_10 = 960.2;   // [DOI:?10.1103/PhysRevC.87.034911?]
0214   const double AuAu_Ncoll_C0_20 = 770.6;   // [DOI:?10.1103/PhysRevC.91.064904?]
0215   const double AuAu_Ncoll_60_70 = 29.8;    //PHYSICAL REVIEW C 87, 034911 (2013)
0216   const double AuAu_Ncoll_70_80 = 12.6;    //PHYSICAL REVIEW C 87, 034911 (2013)
0217   const double AuAu_Ncoll_C0_100 = 238.5;  // BUP2020
0218   const double pAu_Ncoll_C0_100 = 4.7;     // pb^-1 [sPH-TRG-000]
0219 
0220   const double AuAu_rec_3year = (5.7 + 15) * 1e9;       // BUP2020
0221   const double AuAu_rec_5year = AuAu_rec_3year + 30e9;  // BUP2020
0222 
0223   const double error_scale_3yr = sqrt((refAuAuMB) / (AuAu_rec_3year * refAuAuXSec));
0224   const double error_scale_5yr = sqrt((refAuAuMB) / (AuAu_rec_5year * refAuAuXSec));
0225 
0226   TGraphErrors *gr_D_proj_3yr = GraphShiftScaling(gr_D_proj, 0.05, error_scale_3yr);
0227   TGraphErrors *gr_D_proj_2_3yr = GraphShiftScaling(gr_D_proj_2, 0.05, error_scale_3yr);
0228   TGraphErrors *gr_Dbar_proj_3yr = GraphShiftScaling(gr_Dbar_proj, 0.05, error_scale_3yr);
0229   TGraphErrors *gr_Dbar_proj_2_3yr = GraphShiftScaling(gr_Dbar_proj_2, 0.05, error_scale_3yr);
0230 
0231   TGraphErrors *gr_D_proj_5yr = GraphShiftScaling(gr_D_proj, 0., error_scale_5yr);
0232   TGraphErrors *gr_D_proj_2_5yr = GraphShiftScaling(gr_D_proj_2, 0., error_scale_5yr);
0233   TGraphErrors *gr_Dbar_proj_5yr = GraphShiftScaling(gr_Dbar_proj, 0., error_scale_5yr);
0234   TGraphErrors *gr_Dbar_proj_2_5yr = GraphShiftScaling(gr_Dbar_proj_2, 0., error_scale_5yr);
0235 
0236   gr_D_proj_3yr->SetMarkerStyle(kOpenCircle);
0237   gr_D_proj_2_3yr->SetMarkerStyle(kOpenCircle);
0238   gr_Dbar_proj_3yr->SetMarkerStyle(kOpenCircle);
0239   gr_Dbar_proj_2_3yr->SetMarkerStyle(kOpenCircle);
0240 
0241   gr_D_proj_3yr->SetMarkerColor(kGray + 2);
0242   gr_D_proj_2_3yr->SetMarkerColor(kGray + 2);
0243   gr_Dbar_proj_3yr->SetMarkerColor(kBlue - 6);
0244   gr_Dbar_proj_2_3yr->SetMarkerColor(kBlue - 6);
0245 
0246   gr_D_proj_3yr->SetLineColor(kGray + 2);
0247   gr_D_proj_2_3yr->SetLineColor(kGray + 2);
0248   gr_Dbar_proj_3yr->SetLineColor(kBlue - 6);
0249   gr_Dbar_proj_2_3yr->SetLineColor(kBlue - 6);
0250 
0251   //  TCanvas *c1 = new TCanvas("v1_5y_compare", "v1_5y_compare", 0, 0, 800, 600);
0252 
0253   TCanvas *c1 = new TCanvas("D0_BUP2020_AuAu_v1_5y_compare", "D0_BUP2020_AuAu_v1_5y_compare", 1100, 800);
0254   c1->Divide(1, 1);
0255   int idx = 1;
0256   TPad *p;
0257 
0258   p = (TPad *) c1->cd(idx++);
0259   c1->Update();
0260   // gStyle->SetOptFit(0);
0261   // gStyle->SetOptStat(0);
0262   // gStyle->SetEndErrorSize(0.01);
0263   // c1->SetFillColor(10);
0264   // c1->SetBorderMode(0);
0265   // c1->SetBorderSize(2);
0266   // c1->SetFrameFillColor(0);
0267   // c1->SetFrameBorderMode(0);
0268 
0269   // // c1->SetGridx();
0270   // // c1->SetGridy();
0271   // c1->SetTickx();
0272   // c1->SetTicky();
0273 
0274   // //  c1->SetLogx();
0275   // //  c1->SetLogy();
0276 
0277   // c1->SetLeftMargin(0.13);
0278   // c1->SetBottomMargin(0.16);
0279   // c1->SetTopMargin(0.02);
0280   // c1->SetRightMargin(0.06);
0281 
0282   double x1 = -1.2;
0283   double x2 = 1.2;
0284   double y1 = -0.1;
0285   double y2 = 0.1;
0286 
0287   TH1D *d0 = new TH1D("d0", "", 1, x1, x2);
0288   d0->SetMinimum(y1);
0289   d0->SetMaximum(y2);
0290   // d0->GetXaxis()->SetNdivisions(208);
0291   // d0->GetXaxis()->CenterTitle();
0292   d0->GetXaxis()->SetTitle("Rapidity");
0293   // d0->GetXaxis()->SetTitleOffset(1.2);
0294   // d0->GetXaxis()->SetTitleSize(0.06);
0295   // d0->GetXaxis()->SetLabelOffset(0.01);
0296   // d0->GetXaxis()->SetLabelSize(0.045);
0297   // d0->GetXaxis()->SetLabelFont(42);
0298   // d0->GetXaxis()->SetTitleFont(42);
0299   // d0->GetYaxis()->SetNdivisions(505);
0300   // d0->GetYaxis()->CenterTitle();
0301   d0->GetYaxis()->SetTitle("v_{1}");
0302   // d0->GetYaxis()->SetTitleOffset(1.0);
0303   // d0->GetYaxis()->SetTitleSize(0.06);
0304   // d0->GetYaxis()->SetLabelOffset(0.005);
0305   // d0->GetYaxis()->SetLabelSize(0.045);
0306   // d0->GetYaxis()->SetLabelFont(42);
0307   // d0->GetYaxis()->SetTitleFont(42);
0308   // d0->SetLineWidth(2);
0309   d0->Draw("c");
0310 
0311   TLine *l0 = new TLine(x1, 0, x2, 0);
0312   l0->SetLineWidth(1);
0313   l0->SetLineColor(kBlack);
0314   l0->SetLineStyle(1);
0315   l0->Draw("same");
0316 
0317   gr_D_STAR->Draw("p");
0318   gr_Dbar_STAR->Draw("p");
0319 
0320   if (sPHENIX)
0321   {
0322     gr_D->Draw("c");
0323     gr_Dbar->Draw("c");
0324 
0325     gr_D_proj_2_3yr->Draw("p");
0326     gr_Dbar_proj_2_3yr->Draw("p");
0327 
0328     gr_D_proj_3yr->Draw("p");
0329     gr_Dbar_proj_3yr->Draw("p");
0330 
0331     gr_D_proj_2_5yr->Draw("p");
0332     gr_Dbar_proj_2_5yr->Draw("p");
0333 
0334     gr_D_proj_5yr->Draw("p");
0335     gr_Dbar_proj_5yr->Draw("p");
0336   }
0337 
0338   TLegend *leg = new TLegend(.3, .77, .9, .9);
0339   leg->SetFillStyle(0);
0340   leg->AddEntry("", "#it{#bf{sPHENIX}} Projection, Years 1-5", "");
0341   leg->AddEntry("", ("0-80% Au+Au, Res(#Psi_{1})=0.37"), "");
0342   //  leg->AddEntry("", Form("%.0f nb^{-1} str. O+O", OO_rec_5year/1e9), "");
0343   //  leg->AddEntry("", Form("%.0f nb^{-1} str. Ar+Ar", ArAr_rec_5year/1e9), "");
0344   leg->Draw();
0345 
0346   leg = new TLegend(.2, .18, .4, .4, "D^{0}");
0347   leg->SetFillStyle(0);
0348   leg->AddEntry(gr_D_proj_3yr, " ", "pe");
0349   leg->AddEntry(gr_D_proj_5yr, " ", "pe");
0350   leg->AddEntry(gr_D_STAR, " ", "pe");
0351   leg->Draw();
0352 
0353   leg = new TLegend(.25, .18, .5, .4, " #bar{D}^{0}");
0354   leg->SetFillStyle(0);
0355   leg->AddEntry(gr_Dbar_proj_3yr, Form("sPHENIX, Year 1-3, %.0fnb^{-1} rec. Au+Au", AuAu_rec_3year / 1e9), "pe");
0356   leg->AddEntry(gr_Dbar_proj_5yr, Form("sPHENIX, Year 1-5, %.0fnb^{-1} rec.+str. Au+Au", AuAu_rec_5year / 1e9), "pe");
0357   leg->AddEntry(gr_Dbar_STAR, "STAR, 10-80% Au+Au, PRL#bf{123}, 162301", "pe");
0358   leg->Draw();
0359 
0360   //  TLegend *leg = new TLegend(0.12, 0.2, 0.52, 0.4);
0361   //  leg->SetFillColor(10);
0362   //  leg->SetFillStyle(10);
0363   //  leg->SetLineStyle(4000);
0364   //  leg->SetLineColor(10);
0365   //  leg->SetLineWidth(0.);
0366   //  leg->SetTextFont(42);
0367   //  leg->SetTextSize(0.045);
0368   //  leg->AddEntry("", "#it{#bf{sPHENIX}} Simulation", "");
0369   //  leg->AddEntry("", "Au+Au #sqrt{s_{NN}}=200 GeV", "");
0370   //  leg->AddEntry("", "240B MB, EPRes=0.37", "");
0371   //  if (sPHENIX) leg->Draw();
0372   //
0373   //  //  TLegend *leg;
0374   //  if (sPHENIX)
0375   //    leg = new TLegend(0.62, 0.66, 0.98, 0.92);
0376   //  else
0377   //    leg = new TLegend(0.62, 0.79, 0.98, 0.92);
0378   //  leg->SetFillColor(10);
0379   //  leg->SetFillStyle(10);
0380   //  leg->SetLineStyle(4000);
0381   //  leg->SetLineColor(10);
0382   //  leg->SetLineWidth(0.);
0383   //  leg->SetTextFont(42);
0384   //  leg->SetTextSize(0.04);
0385   //  leg->AddEntry(gr_D_STAR, "D^{0} - STAR", "pl");
0386   //  leg->AddEntry(gr_Dbar_STAR, "#bar{D}^{0} - STAR", "pl");
0387   //  if (sPHENIX)
0388   //  {
0389   //    leg->AddEntry(gr_D_proj, "D^{0} - sPHENIX proj.", "pl");
0390   //    leg->AddEntry(gr_Dbar_proj, "#bar{D}^{0} - sPHENIX proj.", "pl");
0391   //  }
0392   //  // leg->AddEntry(gr_proj_D,"Uncert. of D^{0}","p");
0393   //  // leg->AddEntry(gr_proj_D_B,"Uncert. of D^{0} from B","p");
0394   //  leg->Draw();
0395 
0396   //  c1->SaveAs(Form("fig/v1_proj_240B_%d.pdf", sPHENIX));
0397   //  c1->SaveAs(Form("fig/v1_proj_240B_%d.png", sPHENIX));
0398 
0399   SaveCanvas(c1, "fig_BUP2020/" + TString(c1->GetName()), kTRUE);
0400 }