Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2025-08-06 08:14:31 

0001 #include "TSSAplotter.h"
0002 #include <TFile.h>
0003 #include <TDirectoryFile.h>
0004 #include <TF1.h>
0005 #include <TStyle.h>
0006 #include <TROOT.h>
0007 #include <TLatex.h>
0008 #include <TLegend.h>
0009 #include <TLegendEntry.h>
0010 #include <TLine.h>
0011 
0012 #include <string>
0013 #include <utility>
0014 #include <algorithm>
0015 #include <vector>
0016 
0017 void normalize_hist(TH1* hist, double norm=1.0) {
0018     hist->Sumw2();
0019     hist->Scale(norm/hist->Integral());
0020 }
0021 
0022 void SetupPad(TPad* pad, bool logy=false) {
0023     pad->Clear();
0024     pad->Draw();
0025     pad->cd();
0026     pad->SetFillColor(0);
0027     pad->SetBorderMode(0);
0028     pad->SetBorderSize(2);
0029     pad->SetFrameBorderMode(0);
0030     pad->SetFrameBorderMode(0);
0031     pad->SetLogy(logy);
0032     pad->SetLeftMargin(0.12);
0033 }
0034 
0035 void set_draw_style(TH1* hist, int color, int marker_style=7) {
0036     hist->SetStats(10);
0037     hist->SetMarkerStyle(marker_style);
0038     hist->SetMarkerSize(1.5);
0039     hist->SetMarkerColor(color);
0040     hist->SetLineColor(color);
0041     hist->GetXaxis()->CenterTitle(true);
0042     hist->GetXaxis()->SetLabelFont(42);
0043     hist->GetXaxis()->SetLabelSize(0.035);
0044     hist->GetXaxis()->SetTitleSize(0.05);
0045     hist->GetXaxis()->SetTitleOffset(0.8);
0046     hist->GetXaxis()->SetTitleFont(42);
0047     hist->GetYaxis()->CenterTitle(true);
0048     hist->GetYaxis()->SetLabelFont(42);
0049     hist->GetYaxis()->SetLabelSize(0.035);
0050     hist->GetYaxis()->SetTitleSize(0.05);
0051     hist->GetYaxis()->SetTitleFont(42);
0052 }
0053 
0054 void PlotSingleHist(TH1* hist, std::string outpdfname, bool logy=false, bool normalized=false) {
0055     TCanvas *c1 = new TCanvas("c1", "c1",0,50,1600,900);
0056     SetupPad(c1, logy);
0057     if (normalized) {
0058     normalize_hist(hist);
0059     hist->GetYaxis()->SetTitle("Normalized Counts");
0060     }
0061     /* else hist->GetYaxis()->SetTitle("Counts"); */
0062     set_draw_style(hist, 1, 20);
0063     hist->SetStats(0);
0064     hist->Draw("e1 x0");
0065     c1->Modified();
0066     c1->SaveAs(outpdfname.c_str());
0067     delete c1;
0068 }
0069 
0070 void PlotSingle2DHist(TH2* hist, std::string outpdfname, bool logz=false) {
0071     TCanvas *c1 = new TCanvas("c1", "c1",0,50,1600,900);
0072     SetupPad(c1, false);
0073     set_draw_style(hist, 1, 20);
0074     hist->SetStats(0);
0075     hist->Draw("colz");
0076     c1->SetLogz(logz);
0077     c1->Modified();
0078     c1->SaveAs(outpdfname.c_str());
0079     delete c1;
0080 }
0081 
0082 TSSAplotter::TSSAplotter() {}
0083 
0084 TSSAplotter::~TSSAplotter() {}
0085 
0086 void TSSAplotter::GetNMHists() {
0087     TFile* infile = new TFile(NMhistfname.c_str(), "READ");
0088     infile->cd();
0089     
0090     h_nEvents = (TH1*)infile->Get("h_nEvents");
0091     h_vtxz = (TH1*)infile->Get("h_vtxz");
0092     h_towerEta_Phi_500MeV = (TH2*)infile->Get("h_towerEta_Phi_500MeV");
0093     h_towerEta_Phi_800MeV = (TH2*)infile->Get("h_towerEta_Phi_800MeV");
0094     // cluster distributions
0095     h_nAllClusters = (TH1*)infile->Get("h_nAllClusters");
0096     h_nGoodClusters = (TH1*)infile->Get("h_nGoodClusters");
0097     h_clusterE = (TH1*)infile->Get("h_clusterE");
0098     h_clusterEta = (TH1*)infile->Get("h_clusterEta");
0099     h_clusterVtxz_Eta = (TH2*)infile->Get("h_clusterVtxz_Eta");
0100     h_clusterPhi = (TH1*)infile->Get("h_clusterPhi");
0101     h_clusterEta_Phi = (TH2*)infile->Get("h_clusterEta_Phi");
0102     h_clusterpT = (TH1*)infile->Get("h_clusterpT");
0103     h_clusterxF = (TH1*)infile->Get("h_clusterxF");
0104     h_clusterpT_xF = (TH2*)infile->Get("h_clusterpT_xF");
0105     h_clusterChi2 = (TH1*)infile->Get("h_clusterChi2");
0106     h_clusterChi2zoomed = (TH1*)infile->Get("h_clusterChi2zoomed");
0107     h_mesonClusterChi2 = (TH1*)infile->Get("h_mesonClusterChi2");
0108     h_goodClusterE = (TH1*)infile->Get("h_goodClusterE");
0109     h_goodClusterEta = (TH1*)infile->Get("h_goodClusterEta");
0110     h_goodClusterVtxz_Eta = (TH2*)infile->Get("h_goodClusterVtxz_Eta");
0111     h_goodClusterPhi = (TH1*)infile->Get("h_goodClusterPhi");
0112     h_goodClusterEta_Phi = (TH2*)infile->Get("h_goodClusterEta_Phi");
0113     h_goodClusterpT = (TH1*)infile->Get("h_goodClusterpT");
0114     h_goodClusterxF = (TH1*)infile->Get("h_goodClusterxF");
0115     h_goodClusterpT_xF = (TH2*)infile->Get("h_goodClusterpT_xF");
0116     // diphoton distributions
0117     h_nDiphotons = (TH1*)infile->Get("h_nDiphotons");
0118     h_nRecoPi0s = (TH1*)infile->Get("h_nRecoPi0s");
0119     h_nRecoEtas = (TH1*)infile->Get("h_nRecoEtas");
0120     h_diphotonE = (TH1*)infile->Get("h_diphotonE");
0121     h_diphotonMass = (TH1*)infile->Get("h_diphotonMass");
0122     h_diphotonEta = (TH1*)infile->Get("h_diphotonEta");
0123     h_diphotonVtxz_Eta = (TH2*)infile->Get("h_diphotonVtxz_Eta");
0124     h_diphotonPhi = (TH1*)infile->Get("h_diphotonPhi");
0125     h_diphotonEta_Phi = (TH2*)infile->Get("h_diphotonEta_Phi");
0126     h_diphotonpT = (TH1*)infile->Get("h_diphotonpT");
0127     h_diphotonxF = (TH1*)infile->Get("h_diphotonxF");
0128     h_diphotonpT_xF = (TH2*)infile->Get("h_diphotonpT_xF");
0129     h_diphotonAsym = (TH1*)infile->Get("h_diphotonAsym");
0130     h_diphotonDeltaR = (TH1*)infile->Get("h_diphotonDeltaR");
0131     h_diphotonAsym_DeltaR = (TH2*)infile->Get("h_diphotonAsym_DeltaR");
0132 
0133     // BHS mass
0134     double pT_upper = 20.0;
0135     double xF_upper = 0.15;
0136     std::vector<double> pTbinsMass;
0137     std::vector<double> xFbinsMass;
0138     int nbins_bhs = 20;
0139     for (int i=0; i<nbins_bhs; i++) {
0140     pTbinsMass.push_back(i*(pT_upper/nbins_bhs));
0141     xFbinsMass.push_back(2*xF_upper*i/nbins_bhs - xF_upper);
0142     }
0143     pTbinsMass.push_back(pT_upper);
0144     xFbinsMass.push_back(xF_upper);
0145     bhs_diphotonMass_pT = new BinnedHistSet(h_diphotonMass, "p_{T} (GeV)", pTbinsMass);
0146     bhs_diphotonMass_pT->GetHistsFromFile("h_diphotonMass_pT");
0147     bhs_diphotonMass_xF = new BinnedHistSet(h_diphotonMass, "x_{F}", xFbinsMass);
0148     bhs_diphotonMass_xF->GetHistsFromFile("h_diphotonMass_xF");
0149 }
0150 
0151 void TSSAplotter::GetTSSAHists() {
0152     TFile* infile = new TFile(TSSAhistfname.c_str(), "READ");
0153     infile->cd();
0154 
0155     h_nClustersMinbiasTrig = (TH1*)infile->Get("h_nClustersMinbiasTrig");
0156     h_nClustersPhotonTrig = (TH1*)infile->Get("h_nClustersPhotonTrig");
0157     h_Nevents = (TH1*)infile->Get("h_Nevents");
0158     h_Npi0s = (TH1*)infile->Get("h_Npi0s");
0159     h_Netas = (TH1*)infile->Get("h_Netas");
0160 
0161     // BHS mass
0162     /* h_diphotonMass = (TH1*)infile->Get("h_diphotonMass"); */
0163     /* bhs_diphotonMass_pT = new BinnedHistSet(h_diphotonMass, "p_{T} (GeV)", pTbinsMass); */
0164     /* bhs_diphotonMass_pT = new BinnedHistSet(h_diphotonMass, "p_{T} (GeV)", pTbins); */
0165     /* bhs_diphotonMass_pT->GetHistsFromFile("h_diphotonMass_pT"); */
0166 
0167     bhs_diphotonMass_pT_pi0binning = new BinnedHistSet(h_diphotonMass, "p_{T} (GeV)", pTbins_pi0);
0168     bhs_diphotonMass_pT_etabinning = new BinnedHistSet(h_diphotonMass, "p_{T} (GeV)", pTbins_eta);
0169     bhs_pi0_pT_pT = new BinnedHistSet(h_diphotonpT, "p_{T} (GeV)", pTbins_pi0);
0170     bhs_eta_pT_pT = new BinnedHistSet(h_diphotonpT, "p_{T} (GeV)", pTbins_eta);
0171     bhs_diphotonxF_xF = new BinnedHistSet(h_diphotonxF, "x_{F}", xFbins);
0172     bhs_diphotoneta_eta = new BinnedHistSet(h_diphotonEta, "#eta", etabins);
0173     bhs_diphotonvtxz_vtxz = new BinnedHistSet(h_vtxz, "z_{vtx} (cm)", vtxzbins);
0174     bhs_diphotonMass_pT_pi0binning->GetHistsFromFile("h_diphotonMass_pT_pi0binning");
0175     bhs_diphotonMass_pT_etabinning->GetHistsFromFile("h_diphotonMass_pT_etabinning");
0176     bhs_pi0_pT_pT->GetHistsFromFile("h_pi0_pT_pT");
0177     bhs_eta_pT_pT->GetHistsFromFile("h_eta_pT_pT");
0178     bhs_diphotonxF_xF->GetHistsFromFile("h_diphotonxF_xF");
0179     bhs_diphotoneta_eta->GetHistsFromFile("h_diphotoneta_eta");
0180     bhs_diphotonvtxz_vtxz->GetHistsFromFile("h_diphotonvtxz_vtxz");
0181 
0182     h_DiphotonMassAsym = (TH2*)infile->Get("h_DiphotonMassAsym");
0183     h_DiphotonMassdeltaR = (TH2*)infile->Get("h_DiphotonMassdeltaR");
0184     h_DiphotondeltaRAsym = (TH2*)infile->Get("h_DiphotondeltaRAsym");
0185     h_DiphotondeltaRAsym_pi0 = (TH2*)infile->Get("h_DiphotondeltaRAsym_pi0");
0186     h_DiphotondeltaRAsym_pi0_smalldR = (TH2*)infile->Get("h_DiphotondeltaRAsym_pi0_smalldR");
0187     h_DiphotondeltaRAsym_eta = (TH2*)infile->Get("h_DiphotondeltaRAsym_eta");
0188     h_DiphotondeltaRAsym_etabkgr = (TH2*)infile->Get("h_DiphotondeltaRAsym_etabkgr");
0189     h_DiphotonMassAsym_highpT = (TH2*)infile->Get("h_DiphotonMassAsym_highpT");
0190     h_DiphotonMassdeltaR_highpT = (TH2*)infile->Get("h_DiphotonMassdeltaR_highpT");
0191     h_DiphotondeltaRAsym_highpT = (TH2*)infile->Get("h_DiphotondeltaRAsym_highpT");
0192     h_DiphotondeltaRAsym_highpT_smalldR = (TH2*)infile->Get("h_DiphotondeltaRAsym_highpT_smalldR");
0193     h_DiphotondeltaRAsym_eta_highpT = (TH2*)infile->Get("h_DiphotondeltaRAsym_eta_highpT");
0194     h_DiphotondeltaRAsym_etabkgr_highpT = (TH2*)infile->Get("h_DiphotondeltaRAsym_etabkgr_highpT");
0195     // Armenteros-Podolanski plots
0196     h_armen_all = (TH2*)infile->Get("h_armen_all");
0197     h_armen_pi0 = (TH2*)infile->Get("h_armen_pi0");
0198     h_armen_pi0bkgr = (TH2*)infile->Get("h_armen_pi0bkgr");
0199     h_armen_eta = (TH2*)infile->Get("h_armen_eta");
0200     h_armen_etabkgr = (TH2*)infile->Get("h_armen_etabkgr");
0201     h_armen_eta_highpT = (TH2*)infile->Get("h_armen_eta_highpT");
0202     h_armen_etabkgr_highpT = (TH2*)infile->Get("h_armen_etabkgr_highpT");
0203     h_armen_p_L = (TH1*)infile->Get("h_armen_p_L");
0204     h_armen_p_T = (TH1*)infile->Get("h_armen_p_T");
0205     h_armen_alpha_alphaE = (TH2*)infile->Get("h_armen_alpha_alphaE");
0206     h_armen_alpha_deltaR = (TH2*)infile->Get("h_armen_alpha_deltaR");
0207     h_armen_alpha_deltaR_pi0 = (TH2*)infile->Get("h_armen_alpha_deltaR_pi0");
0208     h_armen_alpha_deltaR_eta = (TH2*)infile->Get("h_armen_alpha_deltaR_eta");
0209     h_armen_p_T_deltaR = (TH2*)infile->Get("h_armen_p_T_deltaR");
0210     h_armen_pT_p_L = (TH2*)infile->Get("h_armen_pT_p_L");
0211     h_armen_pT_p_T = (TH2*)infile->Get("h_armen_pT_p_T");
0212 
0213     // Phi hists
0214     std::string nameprefix = "h_pi0_";
0215     std::string titlewhich = "#pi^{0}";
0216     bhs_pi0_blue_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0217     bhs_pi0_blue_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0218     bhs_pi0_yellow_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0219     bhs_pi0_yellow_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0220     bhs_pi0_blue_up_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_up");
0221     bhs_pi0_blue_down_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_down");
0222     bhs_pi0_yellow_up_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_up");
0223     bhs_pi0_yellow_down_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_down");
0224     bhs_pi0_blue_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0225     bhs_pi0_blue_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0226     bhs_pi0_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0227     bhs_pi0_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0228     bhs_pi0_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_up_forward");
0229     bhs_pi0_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_down_forward");
0230     bhs_pi0_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_up_forward");
0231     bhs_pi0_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_down_forward");
0232     bhs_pi0_blue_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0233     bhs_pi0_blue_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0234     bhs_pi0_yellow_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0235     bhs_pi0_yellow_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0236     bhs_pi0_blue_up_backward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_up_backward");
0237     bhs_pi0_blue_down_backward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_blue_down_backward");
0238     bhs_pi0_yellow_up_backward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_up_backward");
0239     bhs_pi0_yellow_down_backward_phi_pT->GetHistsFromFile("h_pi0_phi_pT_yellow_down_backward");
0240     bhs_pi0_lowEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0241     bhs_pi0_lowEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0242     bhs_pi0_lowEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0243     bhs_pi0_lowEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0244     bhs_pi0_lowEta_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0_lowEta_phi_pT_blue_up_forward");
0245     bhs_pi0_lowEta_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0_lowEta_phi_pT_blue_down_forward");
0246     bhs_pi0_lowEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0_lowEta_phi_pT_yellow_up_forward");
0247     bhs_pi0_lowEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0_lowEta_phi_pT_yellow_down_forward");
0248     bhs_pi0_highEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0249     bhs_pi0_highEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0250     bhs_pi0_highEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0251     bhs_pi0_highEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0252     bhs_pi0_highEta_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0_highEta_phi_pT_blue_up_forward");
0253     bhs_pi0_highEta_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0_highEta_phi_pT_blue_down_forward");
0254     bhs_pi0_highEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0_highEta_phi_pT_yellow_up_forward");
0255     bhs_pi0_highEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0_highEta_phi_pT_yellow_down_forward");
0256 
0257     bhs_pi0_blue_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0258     bhs_pi0_blue_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0259     bhs_pi0_yellow_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0260     bhs_pi0_yellow_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0261     bhs_pi0_blue_up_phi_xF->GetHistsFromFile("h_pi0_phi_xF_blue_up");
0262     bhs_pi0_blue_down_phi_xF->GetHistsFromFile("h_pi0_phi_xF_blue_down");
0263     bhs_pi0_yellow_up_phi_xF->GetHistsFromFile("h_pi0_phi_xF_yellow_up");
0264     bhs_pi0_yellow_down_phi_xF->GetHistsFromFile("h_pi0_phi_xF_yellow_down");
0265 
0266     bhs_pi0_blue_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0267     bhs_pi0_blue_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0268     bhs_pi0_yellow_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0269     bhs_pi0_yellow_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0270     bhs_pi0_blue_up_phi_eta->GetHistsFromFile("h_pi0_phi_eta_blue_up");
0271     bhs_pi0_blue_down_phi_eta->GetHistsFromFile("h_pi0_phi_eta_blue_down");
0272     bhs_pi0_yellow_up_phi_eta->GetHistsFromFile("h_pi0_phi_eta_yellow_up");
0273     bhs_pi0_yellow_down_phi_eta->GetHistsFromFile("h_pi0_phi_eta_yellow_down");
0274 
0275     bhs_pi0_blue_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0276     bhs_pi0_blue_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0277     bhs_pi0_yellow_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0278     bhs_pi0_yellow_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", vtxzbins);
0279     bhs_pi0_blue_up_phi_vtxz->GetHistsFromFile("h_pi0_phi_vtxz_blue_up");
0280     bhs_pi0_blue_down_phi_vtxz->GetHistsFromFile("h_pi0_phi_vtxz_blue_down");
0281     bhs_pi0_yellow_up_phi_vtxz->GetHistsFromFile("h_pi0_phi_vtxz_yellow_up");
0282     bhs_pi0_yellow_down_phi_vtxz->GetHistsFromFile("h_pi0_phi_vtxz_yellow_down");
0283 
0284     nameprefix = "h_eta_";
0285     titlewhich = "#eta";
0286     bhs_eta_blue_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0287     bhs_eta_blue_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0288     bhs_eta_yellow_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0289     bhs_eta_yellow_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0290     bhs_eta_blue_up_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_up");
0291     bhs_eta_blue_down_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_down");
0292     bhs_eta_yellow_up_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_up");
0293     bhs_eta_yellow_down_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_down");
0294     bhs_eta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0295     bhs_eta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0296     bhs_eta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0297     bhs_eta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0298     bhs_eta_blue_up_forward_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_up_forward");
0299     bhs_eta_blue_down_forward_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_down_forward");
0300     bhs_eta_yellow_up_forward_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_up_forward");
0301     bhs_eta_yellow_down_forward_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_down_forward");
0302     bhs_eta_lowEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0303     bhs_eta_lowEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0304     bhs_eta_lowEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0305     bhs_eta_lowEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0306     bhs_eta_lowEta_blue_up_forward_phi_pT->GetHistsFromFile("h_eta_lowEta_phi_pT_blue_up_forward");
0307     bhs_eta_lowEta_blue_down_forward_phi_pT->GetHistsFromFile("h_eta_lowEta_phi_pT_blue_down_forward");
0308     bhs_eta_lowEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_eta_lowEta_phi_pT_yellow_up_forward");
0309     bhs_eta_lowEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_eta_lowEta_phi_pT_yellow_down_forward");
0310     bhs_eta_highEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0311     bhs_eta_highEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0312     bhs_eta_highEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0313     bhs_eta_highEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0314     bhs_eta_highEta_blue_up_forward_phi_pT->GetHistsFromFile("h_eta_highEta_phi_pT_blue_up_forward");
0315     bhs_eta_highEta_blue_down_forward_phi_pT->GetHistsFromFile("h_eta_highEta_phi_pT_blue_down_forward");
0316     bhs_eta_highEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_eta_highEta_phi_pT_yellow_up_forward");
0317     bhs_eta_highEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_eta_highEta_phi_pT_yellow_down_forward");
0318     bhs_eta_blue_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0319     bhs_eta_blue_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0320     bhs_eta_yellow_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0321     bhs_eta_yellow_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0322     bhs_eta_blue_up_backward_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_up_backward");
0323     bhs_eta_blue_down_backward_phi_pT->GetHistsFromFile("h_eta_phi_pT_blue_down_backward");
0324     bhs_eta_yellow_up_backward_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_up_backward");
0325     bhs_eta_yellow_down_backward_phi_pT->GetHistsFromFile("h_eta_phi_pT_yellow_down_backward");
0326     bhs_eta_blue_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0327     bhs_eta_blue_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0328     bhs_eta_yellow_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0329     bhs_eta_yellow_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0330     bhs_eta_blue_up_phi_xF->GetHistsFromFile("h_eta_phi_xF_blue_up");
0331     bhs_eta_blue_down_phi_xF->GetHistsFromFile("h_eta_phi_xF_blue_down");
0332     bhs_eta_yellow_up_phi_xF->GetHistsFromFile("h_eta_phi_xF_yellow_up");
0333     bhs_eta_yellow_down_phi_xF->GetHistsFromFile("h_eta_phi_xF_yellow_down");
0334     bhs_eta_blue_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0335     bhs_eta_blue_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0336     bhs_eta_yellow_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0337     bhs_eta_yellow_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0338     bhs_eta_blue_up_phi_eta->GetHistsFromFile("h_eta_phi_eta_blue_up");
0339     bhs_eta_blue_down_phi_eta->GetHistsFromFile("h_eta_phi_eta_blue_down");
0340     bhs_eta_yellow_up_phi_eta->GetHistsFromFile("h_eta_phi_eta_yellow_up");
0341     bhs_eta_yellow_down_phi_eta->GetHistsFromFile("h_eta_phi_eta_yellow_down");
0342     bhs_eta_blue_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0343     bhs_eta_blue_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0344     bhs_eta_yellow_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0345     bhs_eta_yellow_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", vtxzbins);
0346     bhs_eta_blue_up_phi_vtxz->GetHistsFromFile("h_eta_phi_vtxz_blue_up");
0347     bhs_eta_blue_down_phi_vtxz->GetHistsFromFile("h_eta_phi_vtxz_blue_down");
0348     bhs_eta_yellow_up_phi_vtxz->GetHistsFromFile("h_eta_phi_vtxz_yellow_up");
0349     bhs_eta_yellow_down_phi_vtxz->GetHistsFromFile("h_eta_phi_vtxz_yellow_down");
0350 
0351     nameprefix = "h_pi0bkgr_";
0352     titlewhich = "#pi^{0} Background";
0353     bhs_pi0bkgr_blue_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0354     bhs_pi0bkgr_blue_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0355     bhs_pi0bkgr_yellow_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0356     bhs_pi0bkgr_yellow_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0357     bhs_pi0bkgr_blue_up_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_up");
0358     bhs_pi0bkgr_blue_down_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_down");
0359     bhs_pi0bkgr_yellow_up_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_up");
0360     bhs_pi0bkgr_yellow_down_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_down");
0361     bhs_pi0bkgr_blue_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0362     bhs_pi0bkgr_blue_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0363     bhs_pi0bkgr_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0364     bhs_pi0bkgr_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0365     bhs_pi0bkgr_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_up_forward");
0366     bhs_pi0bkgr_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_down_forward");
0367     bhs_pi0bkgr_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_up_forward");
0368     bhs_pi0bkgr_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_down_forward");
0369     bhs_pi0bkgr_blue_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0370     bhs_pi0bkgr_blue_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0371     bhs_pi0bkgr_yellow_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0372     bhs_pi0bkgr_yellow_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0373     bhs_pi0bkgr_blue_up_backward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_up_backward");
0374     bhs_pi0bkgr_blue_down_backward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_blue_down_backward");
0375     bhs_pi0bkgr_yellow_up_backward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_up_backward");
0376     bhs_pi0bkgr_yellow_down_backward_phi_pT->GetHistsFromFile("h_pi0bkgr_phi_pT_yellow_down_backward");
0377     bhs_pi0bkgr_lowEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0378     bhs_pi0bkgr_lowEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0379     bhs_pi0bkgr_lowEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0380     bhs_pi0bkgr_lowEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%slowEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0381     bhs_pi0bkgr_lowEta_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_lowEta_phi_pT_blue_up_forward");
0382     bhs_pi0bkgr_lowEta_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_lowEta_phi_pT_blue_down_forward");
0383     bhs_pi0bkgr_lowEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_lowEta_phi_pT_yellow_up_forward");
0384     bhs_pi0bkgr_lowEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_lowEta_phi_pT_yellow_down_forward");
0385     bhs_pi0bkgr_highEta_blue_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0386     bhs_pi0bkgr_highEta_blue_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0387     bhs_pi0bkgr_highEta_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0388     bhs_pi0bkgr_highEta_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%shighEta_phi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_pi0);
0389     bhs_pi0bkgr_highEta_blue_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_highEta_phi_pT_blue_up_forward");
0390     bhs_pi0bkgr_highEta_blue_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_highEta_phi_pT_blue_down_forward");
0391     bhs_pi0bkgr_highEta_yellow_up_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_highEta_phi_pT_yellow_up_forward");
0392     bhs_pi0bkgr_highEta_yellow_down_forward_phi_pT->GetHistsFromFile("h_pi0bkgr_highEta_phi_pT_yellow_down_forward");
0393 
0394     bhs_pi0bkgr_blue_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0395     bhs_pi0bkgr_blue_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0396     bhs_pi0bkgr_yellow_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0397     bhs_pi0bkgr_yellow_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0398     bhs_pi0bkgr_blue_up_phi_xF->GetHistsFromFile("h_pi0bkgr_phi_xF_blue_up");
0399     bhs_pi0bkgr_blue_down_phi_xF->GetHistsFromFile("h_pi0bkgr_phi_xF_blue_down");
0400     bhs_pi0bkgr_yellow_up_phi_xF->GetHistsFromFile("h_pi0bkgr_phi_xF_yellow_up");
0401     bhs_pi0bkgr_yellow_down_phi_xF->GetHistsFromFile("h_pi0bkgr_phi_xF_yellow_down");
0402     bhs_pi0bkgr_blue_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0403     bhs_pi0bkgr_blue_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0404     bhs_pi0bkgr_yellow_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0405     bhs_pi0bkgr_yellow_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0406     bhs_pi0bkgr_blue_up_phi_eta->GetHistsFromFile("h_pi0bkgr_phi_eta_blue_up");
0407     bhs_pi0bkgr_blue_down_phi_eta->GetHistsFromFile("h_pi0bkgr_phi_eta_blue_down");
0408     bhs_pi0bkgr_yellow_up_phi_eta->GetHistsFromFile("h_pi0bkgr_phi_eta_yellow_up");
0409     bhs_pi0bkgr_yellow_down_phi_eta->GetHistsFromFile("h_pi0bkgr_phi_eta_yellow_down");
0410     bhs_pi0bkgr_blue_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0411     bhs_pi0bkgr_blue_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0412     bhs_pi0bkgr_yellow_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0413     bhs_pi0bkgr_yellow_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", vtxzbins);
0414     bhs_pi0bkgr_blue_up_phi_vtxz->GetHistsFromFile("h_pi0bkgr_phi_vtxz_blue_up");
0415     bhs_pi0bkgr_blue_down_phi_vtxz->GetHistsFromFile("h_pi0bkgr_phi_vtxz_blue_down");
0416     bhs_pi0bkgr_yellow_up_phi_vtxz->GetHistsFromFile("h_pi0bkgr_phi_vtxz_yellow_up");
0417     bhs_pi0bkgr_yellow_down_phi_vtxz->GetHistsFromFile("h_pi0bkgr_phi_vtxz_yellow_down");
0418 
0419     nameprefix = "h_etabkgr_";
0420     titlewhich = "#eta Background";
0421     bhs_etabkgr_blue_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0422     bhs_etabkgr_blue_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0423     bhs_etabkgr_yellow_up_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0424     bhs_etabkgr_yellow_down_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0425     bhs_etabkgr_blue_up_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_up");
0426     bhs_etabkgr_blue_down_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_down");
0427     bhs_etabkgr_yellow_up_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_up");
0428     bhs_etabkgr_yellow_down_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_down");
0429     bhs_etabkgr_blue_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0430     bhs_etabkgr_blue_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_forward", nameprefix.c_str()), Form("%s Blue Beam Forward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0431     bhs_etabkgr_yellow_up_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0432     bhs_etabkgr_yellow_down_forward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_forward", nameprefix.c_str()), Form("%s Yellow Beam Forward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0433     bhs_etabkgr_blue_up_forward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_up_forward");
0434     bhs_etabkgr_blue_down_forward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_down_forward");
0435     bhs_etabkgr_yellow_up_forward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_up_forward");
0436     bhs_etabkgr_yellow_down_forward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_down_forward");
0437     bhs_etabkgr_blue_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_up_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0438     bhs_etabkgr_blue_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_blue_down_backward", nameprefix.c_str()), Form("%s Blue Beam Backward-Going Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0439     bhs_etabkgr_yellow_up_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_up_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0440     bhs_etabkgr_yellow_down_backward_phi_pT = new BinnedHistSet(Form("%sphi_pT_yellow_down_backward", nameprefix.c_str()), Form("%s Yellow Beam Backward-Going Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "p_{T} (GeV)", pTbins_eta);
0441     bhs_etabkgr_blue_up_backward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_up_backward");
0442     bhs_etabkgr_blue_down_backward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_blue_down_backward");
0443     bhs_etabkgr_yellow_up_backward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_up_backward");
0444     bhs_etabkgr_yellow_down_backward_phi_pT->GetHistsFromFile("h_etabkgr_phi_pT_yellow_down_backward");
0445     bhs_etabkgr_blue_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0446     bhs_etabkgr_blue_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0447     bhs_etabkgr_yellow_up_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0448     bhs_etabkgr_yellow_down_phi_xF = new BinnedHistSet(Form("%sphi_xF_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", xFbins);
0449     bhs_etabkgr_blue_up_phi_xF->GetHistsFromFile("h_etabkgr_phi_xF_blue_up");
0450     bhs_etabkgr_blue_down_phi_xF->GetHistsFromFile("h_etabkgr_phi_xF_blue_down");
0451     bhs_etabkgr_yellow_up_phi_xF->GetHistsFromFile("h_etabkgr_phi_xF_yellow_up");
0452     bhs_etabkgr_yellow_down_phi_xF->GetHistsFromFile("h_etabkgr_phi_xF_yellow_down");
0453     bhs_etabkgr_blue_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0454     bhs_etabkgr_blue_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0455     bhs_etabkgr_yellow_up_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0456     bhs_etabkgr_yellow_down_phi_eta = new BinnedHistSet(Form("%sphi_eta_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "#eta", etabins);
0457     bhs_etabkgr_blue_up_phi_eta->GetHistsFromFile("h_etabkgr_phi_eta_blue_up");
0458     bhs_etabkgr_blue_down_phi_eta->GetHistsFromFile("h_etabkgr_phi_eta_blue_down");
0459     bhs_etabkgr_yellow_up_phi_eta->GetHistsFromFile("h_etabkgr_phi_eta_yellow_up");
0460     bhs_etabkgr_yellow_down_phi_eta->GetHistsFromFile("h_etabkgr_phi_eta_yellow_down");
0461     bhs_etabkgr_blue_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_up", nameprefix.c_str()), Form("%s Blue Beam Spin-Up #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0462     bhs_etabkgr_blue_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_blue_down", nameprefix.c_str()), Form("%s Blue Beam Spin-Down #phi^{B} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0463     bhs_etabkgr_yellow_up_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_up", nameprefix.c_str()), Form("%s Yellow Beam Spin-Up #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "z_{vtx}", vtxzbins);
0464     bhs_etabkgr_yellow_down_phi_vtxz = new BinnedHistSet(Form("%sphi_vtxz_yellow_down", nameprefix.c_str()), Form("%s Yellow Beam Spin-Down #phi^{Y} Distribution;#phi (rad);Counts", titlewhich.c_str()), nHistBins_phi, -1.0*PI, PI, "x_{F}", vtxzbins);
0465     bhs_etabkgr_blue_up_phi_vtxz->GetHistsFromFile("h_etabkgr_phi_vtxz_blue_up");
0466     bhs_etabkgr_blue_down_phi_vtxz->GetHistsFromFile("h_etabkgr_phi_vtxz_blue_down");
0467     bhs_etabkgr_yellow_up_phi_vtxz->GetHistsFromFile("h_etabkgr_phi_vtxz_yellow_up");
0468     bhs_etabkgr_yellow_down_phi_vtxz->GetHistsFromFile("h_etabkgr_phi_vtxz_yellow_down");
0469 }
0470 
0471 void TSSAplotter::PlotNMHists() {
0472     std::string outfilename = NMpdfname;
0473     std::string outfilename_start = outfilename + "(";
0474     std::string outfilename_end = outfilename + ")";
0475     TCanvas* c = new TCanvas("c", "c", 1600, 900);
0476 
0477     PlotSingleHist(h_nEvents, outfilename_start);
0478     PlotSingleHist(h_vtxz, outfilename);
0479     PlotSingle2DHist(h_towerEta_Phi_500MeV, outfilename);
0480     PlotSingle2DHist(h_towerEta_Phi_500MeV, outfilename, true);
0481     PlotSingle2DHist(h_towerEta_Phi_800MeV, outfilename);
0482     PlotSingle2DHist(h_towerEta_Phi_800MeV, outfilename, true);
0483 
0484     // clusters
0485     PlotSingleHist(h_nAllClusters, outfilename);
0486     PlotSingleHist(h_nGoodClusters, outfilename);
0487     PlotSingleHist(h_clusterE, outfilename, true);
0488     PlotSingleHist(h_clusterEta, outfilename);
0489     PlotSingle2DHist(h_clusterVtxz_Eta, outfilename);
0490     PlotSingle2DHist(h_clusterVtxz_Eta, outfilename, true);
0491     PlotSingleHist(h_clusterPhi, outfilename);
0492     PlotSingle2DHist(h_clusterEta_Phi, outfilename);
0493     PlotSingleHist(h_clusterpT, outfilename, true);
0494     PlotSingleHist(h_clusterxF, outfilename, true);
0495     PlotSingle2DHist(h_clusterpT_xF, outfilename, true);
0496     PlotSingleHist(h_clusterChi2, outfilename, true);
0497     PlotSingleHist(h_clusterChi2zoomed, outfilename, true);
0498     c->SetLogy();
0499     normalize_hist(h_clusterChi2zoomed);
0500     h_clusterChi2zoomed->Draw();
0501     normalize_hist(h_mesonClusterChi2);
0502     /* h_mesonClusterChi2->Scale(1.0/h_mesonClusterChi2->Integral()); */
0503     h_mesonClusterChi2->Draw("same");
0504     c->SaveAs(outfilename.c_str());
0505     h_clusterChi2zoomed->Scale(h_clusterChi2zoomed->GetEntries()/h_clusterChi2zoomed->Integral());
0506     h_mesonClusterChi2->Scale(h_mesonClusterChi2->GetEntries()/h_mesonClusterChi2->Integral());
0507     // good clusters
0508     PlotSingleHist(h_goodClusterE, outfilename, true);
0509     PlotSingleHist(h_goodClusterEta, outfilename);
0510     PlotSingle2DHist(h_goodClusterVtxz_Eta, outfilename);
0511     PlotSingle2DHist(h_goodClusterVtxz_Eta, outfilename, true);
0512     PlotSingleHist(h_goodClusterPhi, outfilename);
0513     PlotSingle2DHist(h_goodClusterEta_Phi, outfilename);
0514     PlotSingleHist(h_goodClusterpT, outfilename, true);
0515     PlotSingleHist(h_goodClusterxF, outfilename, true);
0516     PlotSingle2DHist(h_goodClusterpT_xF, outfilename, true);
0517 
0518     // diphotons
0519     PlotSingleHist(h_nDiphotons, outfilename);
0520     PlotSingleHist(h_nRecoPi0s, outfilename);
0521     PlotSingleHist(h_nRecoEtas, outfilename);
0522     PlotSingleHist(h_diphotonMass, outfilename);
0523     PlotSingleHist(h_diphotonE, outfilename, true);
0524     PlotSingleHist(h_diphotonEta, outfilename);
0525     PlotSingle2DHist(h_diphotonVtxz_Eta, outfilename, true);
0526     PlotSingleHist(h_diphotonPhi, outfilename);
0527     PlotSingle2DHist(h_diphotonEta_Phi, outfilename);
0528     PlotSingleHist(h_diphotonpT, outfilename, true);
0529     PlotSingleHist(h_diphotonxF, outfilename, true);
0530     PlotSingle2DHist(h_diphotonpT_xF, outfilename, true);
0531     PlotSingleHist(h_diphotonAsym, outfilename);
0532     PlotSingleHist(h_diphotonDeltaR, outfilename);
0533     PlotSingle2DHist(h_diphotonAsym_DeltaR, outfilename, true);
0534 
0535     bhs_diphotonMass_pT->PlotAllHistsWithFits(outfilename, false, "mass2");
0536     /* bhs_diphotonMass_pT->hist_vec.at(0)->Draw("x0 e1"); */
0537     bhs_diphotonMass_xF->PlotAllHistsWithFits(outfilename, false, "mass2");
0538 
0539     PlotSingleHist(h_vtxz, outfilename_end);
0540 }
0541 
0542 void TSSAplotter::PlotTSSAs() {
0543     std::string outfilename = TSSApdfname;
0544     std::string outfilename_start = outfilename + "(";
0545     std::string outfilename_end = outfilename + ")";
0546     /* TCanvas* c = new TCanvas("c", "c", 1600, 900); */
0547 
0548     PlotSingleHist(h_nClustersMinbiasTrig, outfilename_start);
0549     PlotSingleHist(h_nClustersPhotonTrig, outfilename);
0550     PlotSingleHist(h_Nevents, outfilename);
0551     PlotSingleHist(h_Npi0s, outfilename);
0552     PlotSingleHist(h_Netas, outfilename);
0553     PlotSingle2DHist(h_DiphotonMassAsym, outfilename);
0554     PlotSingle2DHist(h_DiphotonMassdeltaR, outfilename);
0555     PlotSingle2DHist(h_DiphotondeltaRAsym, outfilename);
0556     PlotSingle2DHist(h_DiphotondeltaRAsym_pi0, outfilename);
0557     PlotSingle2DHist(h_DiphotondeltaRAsym_pi0_smalldR, outfilename);
0558     PlotSingle2DHist(h_DiphotondeltaRAsym_eta, outfilename);
0559     PlotSingle2DHist(h_DiphotondeltaRAsym_etabkgr, outfilename);
0560     PlotSingle2DHist(h_DiphotonMassAsym_highpT, outfilename);
0561     PlotSingle2DHist(h_DiphotonMassdeltaR_highpT, outfilename);
0562     PlotSingle2DHist(h_DiphotondeltaRAsym_highpT, outfilename);
0563     PlotSingle2DHist(h_DiphotondeltaRAsym_highpT_smalldR, outfilename);
0564     PlotSingle2DHist(h_DiphotondeltaRAsym_eta_highpT, outfilename);
0565     PlotSingle2DHist(h_DiphotondeltaRAsym_etabkgr_highpT, outfilename);
0566     PlotSingle2DHist(h_armen_all, outfilename, true);
0567     PlotSingle2DHist(h_armen_pi0, outfilename, true);
0568     PlotSingle2DHist(h_armen_pi0bkgr, outfilename, true);
0569     PlotSingle2DHist(h_armen_eta, outfilename, true);
0570     PlotSingle2DHist(h_armen_etabkgr, outfilename, true);
0571     PlotSingle2DHist(h_armen_eta_highpT, outfilename, true);
0572     PlotSingle2DHist(h_armen_etabkgr_highpT, outfilename, true);
0573     PlotSingleHist(h_armen_p_L, outfilename, true);
0574     PlotSingleHist(h_armen_p_T, outfilename, true);
0575     PlotSingle2DHist(h_armen_alpha_alphaE, outfilename, true);
0576     PlotSingle2DHist(h_armen_alpha_deltaR, outfilename, true);
0577     PlotSingle2DHist(h_armen_alpha_deltaR_pi0, outfilename, true);
0578     PlotSingle2DHist(h_armen_alpha_deltaR_eta, outfilename, true);
0579     PlotSingle2DHist(h_armen_p_T_deltaR, outfilename, true);
0580     PlotSingle2DHist(h_armen_pT_p_L, outfilename, true);
0581     PlotSingle2DHist(h_armen_pT_p_T, outfilename, true);
0582 
0583     bhs_diphotonMass_pT_pi0binning->PlotAllHistsWithFits(outfilename, false, "mass2");
0584     bhs_diphotonMass_pT_etabinning->PlotAllHistsWithFits(outfilename, false, "mass2");
0585     PlotBackgroundFractions(outfilename);
0586     bhs_pi0_pT_pT->PlotAllHistsWithFits(outfilename, false, "");
0587     bhs_eta_pT_pT->PlotAllHistsWithFits(outfilename, false, "");
0588     bhs_diphotonxF_xF->PlotAllHistsWithFits(outfilename, false, "");
0589     bhs_diphotoneta_eta->PlotAllHistsWithFits(outfilename, false, "");
0590     bhs_diphotonvtxz_vtxz->PlotAllHistsWithFits(outfilename, false, "");
0591 
0592     PlotPi0(outfilename);
0593     PlotEta(outfilename);
0594 
0595     PlotSingleHist(h_nClustersMinbiasTrig, outfilename_end);
0596 }
0597 
0598 double TSSAplotter::RelLumiAsym(double Nup, double Ndown, double relLumi) {
0599     double asym = (Nup - relLumi*Ndown)/(Nup + relLumi*Ndown);
0600     return asym;
0601 }
0602 
0603 double TSSAplotter::RelLumiError(double Nup, double Ndown, double relLumi) {
0604     double asym = RelLumiAsym(Nup, Ndown, relLumi);
0605     double numerator = Nup - relLumi*Ndown;
0606     /* double num_err = sqrt(Nup) - relLumi*sqrt(Ndown); */
0607     double num_err = sqrt(Nup + (relLumi*relLumi*Ndown));
0608     double num_rel_err = num_err/numerator;
0609     double denominator = Nup + relLumi*Ndown;
0610     /* double den_err = sqrt(Nup) + relLumi*sqrt(Ndown); */
0611     double den_err = sqrt(Nup + (relLumi*relLumi*Ndown));
0612     double den_rel_err = den_err/denominator;
0613     
0614     double error = abs(asym)*sqrt((num_rel_err*num_rel_err) + (den_rel_err*den_rel_err));
0615     /* std::cout << Form("numerator=%f, num_err=%f, num_rel_err=%f\ndenominator=%f, den_err=%f, den_rel_err=%f\nasym=%f, asym_err=%f, asym_rel_err=%f\n", numerator, num_err, num_rel_err, denominator, den_err, den_rel_err, asym, error, (error/asym)); */
0616     return error;
0617 }
0618 
0619 double TSSAplotter::SqrtAsym(double NLup, double NLdown, double NRup, double NRdown) {
0620     double asym = (sqrt(NLup*NRdown)-sqrt(NRup*NLdown))/(sqrt(NLup*NRdown)+sqrt(NRup*NLdown));
0621     return asym;
0622 }
0623 
0624 double TSSAplotter::SqrtError(double NLup, double NLdown, double NRup, double NRdown,
0625     double NLupErr, double NLdownErr, double NRupErr, double NRdownErr) {
0626     double t1 = sqrt(NLdown*NRup*NRdown/NLup)*NLupErr;
0627     double t2 = sqrt(NLup*NRup*NRdown/NLdown)*NLdownErr;
0628     double t3 = sqrt(NLup*NLdown*NRdown/NRup)*NRupErr;
0629     double t4 = sqrt(NLup*NLdown*NRup/NRdown)*NRdownErr;
0630     double asym_err = (1/pow(sqrt(NLup*NRdown)+sqrt(NRup*NLdown),2))*sqrt(pow(t1,2)+pow(t2,2)+pow(t3,2)+pow(t4,2));
0631     return asym_err;
0632 }
0633 
0634 TGraphErrors* TSSAplotter::RelLumiGraph(TH1* phi_up, TH1* phi_down, double relLumi) {
0635     TGraphErrors* gr = new TGraphErrors();
0636     int nbins = phi_up->GetNbinsX();
0637     for (int i = 0; i < nbins; i++)
0638     {
0639     float phi = i*(2*PI/nbins) - (PI - PI/nbins);
0640     int phibin = i;
0641     
0642     double Nup = phi_up->GetBinContent(1+phibin);
0643     double Ndown = phi_down->GetBinContent(1+phibin);
0644     /* double NupErr = phi_up->GetBinError(1+phibin); */
0645     /* double NdownErr = phi_down->GetBinError(1+phibin); */
0646 
0647     double asym = RelLumiAsym(Nup, Ndown, relLumi);
0648     double err = RelLumiError(Nup, Ndown, relLumi);
0649     gr->SetPoint(i, phi, asym);
0650     gr->SetPointError(i, 0, err);
0651     }
0652     return gr;
0653 }
0654 
0655 TGraphErrors* TSSAplotter::SqrtGraph(TH1* phi_up, TH1* phi_down) {
0656     TGraphErrors* gr = new TGraphErrors();
0657     int nbins = phi_up->GetNbinsX();
0658     int npoints = (int)(nbins / 2);
0659     int offset = (int)(npoints / 2);
0660     for (int i = 0; i < npoints; i++)
0661     {
0662     float phi = i*(PI/npoints) - (PI/2 - PI/npoints/2);
0663     int phibinL = (i + offset) % nbins; // N_Left
0664     int phibinR = (phibinL + npoints + offset) % nbins; // N_Right
0665     
0666     double NLup = phi_up->GetBinContent(1+phibinL);
0667     double NLdown = phi_down->GetBinContent(1+phibinL);
0668     double NRup = phi_up->GetBinContent(1+phibinR);
0669     double NRdown = phi_down->GetBinContent(1+phibinR);
0670     double NLupErr = phi_up->GetBinError(1+phibinL);
0671     double NLdownErr = phi_down->GetBinError(1+phibinL);
0672     double NRupErr = phi_up->GetBinError(1+phibinR);
0673     double NRdownErr = phi_down->GetBinError(1+phibinR);
0674 
0675     double asym = SqrtAsym(NLup, NLdown, NRup, NRdown);
0676     double err = SqrtError(NLup, NLdown, NRup, NRdown,
0677         NLupErr, NLdownErr, NRupErr, NRdownErr);
0678     gr->SetPoint(i, phi, asym);
0679     gr->SetPointError(i, 0, err);
0680     }
0681     return gr;
0682 }
0683 
0684 void TSSAplotter::PlotOneAsym(std::string type, std::string title, TPad* pad, TGraphErrors* gr) {
0685     pad->cd();
0686     pad->Clear();
0687     double fit_lower, fit_upper;
0688     if (type == "sqrt") {
0689     fit_lower = -PI/2.0;
0690     fit_upper = PI/2.0;
0691     }
0692     else {
0693     fit_lower = -PI;
0694     fit_upper = PI;
0695     }
0696 
0697     TF1* fit = nullptr;
0698     if (type == "sqrt" || type == "rellumi") {
0699     fit = new TF1("asymfit", "[0]*sin(x - [1])", fit_lower, fit_upper);
0700     /* fit->SetParLimits(0, 0.00001, 0.4); */
0701     fit->SetParLimits(0, -0.4, 0.4);
0702     fit->SetParameter(0, 0.01);
0703     fit->SetParLimits(1, PI-0.2, PI+0.2);
0704     /* fit->SetParLimits(1, 0.0, 2*PI); */
0705     fit->SetParameter(1, 3.141592);
0706     fit->SetParName(0, "#epsilon");
0707     fit->SetParName(1, "#phi_{0}");
0708     }
0709     else {
0710     std::cout << "Invalid type option " << type << "; exiting!" << std::endl;
0711     return;
0712     }
0713 
0714     gr->SetTitle(title.c_str());
0715     gr->GetXaxis()->SetTitle("#phi");
0716     if (type == "sqrt") {
0717     gr->GetYaxis()->SetTitle("Raw Asymmetry (Square Root)");
0718     }
0719     else {
0720     gr->GetYaxis()->SetTitle("Raw Asymmetry (Rel. Lumi.)");
0721     }
0722     gStyle->SetOptFit();
0723     gStyle->SetOptStat(0);
0724     gr->Draw("ap");
0725     gr->Fit(fit, "Q");
0726     gr->GetXaxis()->SetRangeUser(fit_lower, fit_upper);
0727     gr->GetYaxis()->SetRangeUser(-0.1, 0.1);
0728 
0729     TLatex* tl = new TLatex;
0730     tl->SetTextColor(kRed);
0731     tl->DrawLatexNDC(0.20, 0.70, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit->GetParameter(0), fit->GetParError(0), fit->GetParameter(0)/fit->GetParError(0)));
0732     /* for (int i=0; i<gr->GetN(); i++) { */
0733     /* std::cout << "Point " << i << ": x=" << gr->GetPointX(i) << ", y=" << gr->GetPointY(i) << " +- " << gr->GetErrorY(i) << std::endl; */
0734     /* } */
0735     pad->Update();
0736 }
0737 
0738 void TSSAplotter::CompareOneAsym(std::string title, std::string outfilename, TGraphErrors* rl, TGraphErrors* sqrt) {
0739     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
0740     c1->cd();
0741     double fit_lower = -PI;
0742     double fit_upper = PI;
0743 
0744     TF1* fit_rl = new TF1("fit_rl", "[0]*sin(x - [1])", fit_lower, fit_upper);
0745     fit_rl->SetParLimits(0, 0.00001, 0.4);
0746     fit_rl->SetParameter(0, 0.01);
0747     fit_rl->SetParLimits(1, 0.0, 2*PI);
0748     fit_rl->SetParameter(1, 3.141592);
0749     fit_rl->SetParName(0, "#epsilon");
0750     fit_rl->SetParName(1, "#phi_{0}");
0751     fit_rl->SetLineColor(kRed);
0752 
0753     rl->SetTitle(title.c_str());
0754     rl->GetXaxis()->SetTitle("#phi");
0755     rl->GetYaxis()->SetTitle("Raw Asymmetry");
0756     rl->GetXaxis()->SetRangeUser(fit_lower, fit_upper);
0757     rl->GetYaxis()->SetRangeUser(-0.002, 0.002);
0758     /* std::cout << "RelLumi graph" << std::endl; */
0759     for (int i=0; i<rl->GetN(); i++) {
0760     std::cout << "Point " << i << ": x=" << rl->GetPointX(i) << ", y=" << rl->GetPointY(i) << " +- " << rl->GetErrorY(i) << std::endl;
0761     }
0762 
0763     int rl_marker = rl->GetMarkerStyle();
0764     int rl_mcolor = rl->GetMarkerColor();
0765     int rl_lcolor = rl->GetLineColor();
0766     rl->SetMarkerStyle(kStar);
0767     rl->SetMarkerColor(kRed);
0768     rl->SetLineColor(kRed);
0769 
0770     TF1* fit_sqrt = new TF1("fit_sqrt", "[0]*sin(x - [1])", fit_lower/2.0, fit_upper/2.0);
0771     fit_sqrt->SetParLimits(0, 0.00001, 0.4);
0772     fit_sqrt->SetParameter(0, 0.01);
0773     fit_sqrt->SetParLimits(1, 0.0, 2*PI);
0774     fit_sqrt->SetParameter(1, 3.141592);
0775     fit_sqrt->SetParName(0, "#epsilon");
0776     fit_sqrt->SetParName(1, "#phi_{0}");
0777     fit_sqrt->SetLineColor(kBlue);
0778 
0779     sqrt->SetTitle(title.c_str());
0780     sqrt->GetXaxis()->SetTitle("#phi");
0781     sqrt->GetYaxis()->SetTitle("Raw Asymmetry");
0782     sqrt->GetXaxis()->SetRangeUser(fit_lower, fit_upper);
0783     /* sqrt->GetYaxis()->SetRangeUser(-0.05, 0.05); */
0784     /* std::cout << "Sqrt graph" << std::endl; */
0785     for (int i=0; i<sqrt->GetN(); i++) {
0786     std::cout << "Point " << i << ": x=" << sqrt->GetPointX(i) << ", y=" << sqrt->GetPointY(i) << " +- " << sqrt->GetErrorY(i) << std::endl;
0787     }
0788 
0789     int sqrt_marker = sqrt->GetMarkerStyle();
0790     int sqrt_mcolor = sqrt->GetMarkerColor();
0791     int sqrt_lcolor = sqrt->GetLineColor();
0792     sqrt->SetMarkerStyle(kCircle);
0793     sqrt->SetMarkerColor(kBlue);
0794     sqrt->SetLineColor(kBlue);
0795     /* std::cout << "sqrt nPoints = " << sqrt->GetN() << std::endl; */
0796     /* std::cout << "rl nPoints = " << rl->GetN() << std::endl; */
0797 
0798     gStyle->SetOptFit(0);
0799     gStyle->SetOptStat(0);
0800     rl->Draw("ap");
0801     rl->Fit("fit_rl", "Q");
0802     sqrt->Draw("p");
0803     sqrt->Fit("fit_sqrt", "Q");
0804     TLatex* tl_rl = new TLatex;
0805     tl_rl->SetTextColor(kRed);
0806     tl_rl->DrawLatexNDC(0.20, 0.75, Form("#splitline{#epsilon = %f #pm %f}{#phi_{0} = %f #pm %f}", fit_rl->GetParameter(0), fit_rl->GetParError(0), fit_rl->GetParameter(1), fit_rl->GetParError(1)));
0807     TLatex* tl_sqrt = new TLatex;
0808     tl_sqrt->SetTextColor(kBlue);
0809     tl_sqrt->DrawLatexNDC(0.20, 0.25, Form("#splitline{#epsilon = %f #pm %f}{#phi_{0} = %f #pm %f}", fit_sqrt->GetParameter(0), fit_sqrt->GetParError(0), fit_sqrt->GetParameter(1), fit_sqrt->GetParError(1)));
0810     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
0811     leg.AddEntry(rl, "Rel. Lumi", "ep");
0812     leg.AddEntry(sqrt, "Geometric", "ep");
0813     leg.Draw();
0814 
0815     c1->SaveAs(outfilename.c_str());
0816     rl->SetMarkerStyle(rl_marker);
0817     rl->SetMarkerColor(rl_mcolor);
0818     rl->SetLineColor(rl_lcolor);
0819     sqrt->SetMarkerStyle(sqrt_marker);
0820     sqrt->SetMarkerColor(sqrt_mcolor);
0821     sqrt->SetLineColor(sqrt_lcolor);
0822 }
0823 
0824 void TSSAplotter::PlotAsymsBinned(std::string type, std::string which, BinnedHistSet* bhs_up, BinnedHistSet* bhs_down, std::string outfilename) {
0825     TCanvas* c1 = new TCanvas(which.c_str(), "", 1600, 900);
0826     c1->Divide(4, 2, 0.025, 0.025);
0827     
0828     int nbins = bhs_up->nbins;
0829     for (int i=1; i<=nbins; i++) {
0830     TPad* pad = (TPad*)c1->GetPad(i);
0831     if (!pad) return;
0832     pad->cd();
0833     gStyle->SetPadLeftMargin(0.15);
0834     gROOT->ForceStyle();
0835     c1->Update();
0836     TH1* phi_up = bhs_up->hist_vec[i];
0837     TH1* phi_down = bhs_down->hist_vec[i];
0838     TGraphErrors* gr = nullptr;
0839     if (type == "sqrt") {
0840         gr = SqrtGraph(phi_up, phi_down);
0841     }
0842     else {
0843         bool blue = (which.find("Blue") != std::string::npos);
0844         if (blue) {
0845         gr = RelLumiGraph(phi_up, phi_down, relLumiBlue);
0846         }
0847         else {
0848         gr = RelLumiGraph(phi_up, phi_down, relLumiYellow);
0849         }
0850     }
0851     std::string title_type = "";
0852     if (type == "sqrt") title_type = "Square Root Asymmetry";
0853     if (type == "rellumi") title_type = "Rel. Lumi. Asymmetry";
0854     char* range;
0855     double lower, upper;
0856     lower = bhs_up->bin_edges[i-1];
0857     upper = bhs_up->bin_edges[i];
0858     range = Form("[%.2f#leq%s<%.2f]", lower, bhs_up->binned_quantity.c_str(), upper);
0859     /* char* title = Form("%s %s %s", which.c_str(), title_type.c_str(), range); */
0860     char* title = Form("%s %s", which.c_str(), title_type.c_str());
0861     PlotOneAsym(type, std::string(title), pad, gr);
0862     TLatex* textbox = new TLatex(0, 0, "blaahhh");
0863     /* textbox->SetTextSize(textbox->GetTextSize()*1.25); */
0864     /* textbox->SetTextColor(kRed); */
0865     textbox->SetTextAlign(21);
0866     textbox->DrawLatexNDC(0.35, 0.8, range);
0867     /* std::cout << "Greg info: binned_quantity = " << bhs_up->binned_quantity << std::endl; */
0868     /* std::cout << "Greg info: gr num points = " << gr->GetN() << std::endl; */
0869     }
0870     c1->Update();
0871     c1->SaveAs(outfilename.c_str());
0872     delete c1;
0873 }
0874 
0875 void TSSAplotter::PlotAsymsBinned(std::string which, BinnedHistSet* bhs_up, BinnedHistSet* bhs_down, std::string outfilename) {
0876     TCanvas* c1 = new TCanvas(which.c_str(), "", 1600, 900);
0877     c1->Divide(4, 2, 0.025, 0.025);
0878     TCanvas* c2 = new TCanvas("c2", "", 1600, 900);
0879     
0880     int nbins = bhs_up->nbins;
0881     for (int i=1; i<=nbins; i++) {
0882     TPad* pad = (TPad*)c1->GetPad(i);
0883     if (!pad) return;
0884     pad->cd();
0885     gStyle->SetPadLeftMargin(0.15);
0886     gROOT->ForceStyle();
0887     c1->Update();
0888     TH1* phi_up = bhs_up->hist_vec[i];
0889     TH1* phi_down = bhs_down->hist_vec[i];
0890     TGraphErrors* gr_rl = nullptr;
0891     bool blue = (which.find("Blue") != std::string::npos);
0892     if (blue) {
0893         gr_rl = RelLumiGraph(phi_up, phi_down, relLumiBlue);
0894     }
0895     else {
0896         gr_rl = RelLumiGraph(phi_up, phi_down, relLumiYellow);
0897     }
0898     std::pair<double, double> asym_rl = FitGraph("rellumi", gr_rl);
0899     if (asym_rl.first == 0.0) { // needed so compiler doesn't complain about unused variables
0900     }
0901     TF1* fit_rl = gr_rl->GetFunction("asymfit");
0902     fit_rl->SetLineColor(kRed-2);
0903     TGraphErrors* gr_sqrt = SqrtGraph(phi_up, phi_down);
0904     std::pair<double, double> asym_sqrt = FitGraph("sqrt", gr_sqrt);
0905     if (asym_sqrt.first == 0.0) { // needed so compiler doesn't complain about unused variables
0906     }
0907     TF1* fit_sqrt = gr_sqrt->GetFunction("asymfit");
0908     fit_sqrt->SetLineColor(kBlue-2);
0909 
0910     std::string title_type = "Raw Asymmetry";
0911     char* range;
0912     double lower, upper;
0913     lower = bhs_up->bin_edges[i-1];
0914     upper = bhs_up->bin_edges[i];
0915     range = Form("[%.2f#leq%s<%.2f]", lower, bhs_up->binned_quantity.c_str(), upper);
0916     /* char* title = Form("%s %s %s", which.c_str(), title_type.c_str(), range); */
0917     char* title = Form("%s %s", which.c_str(), title_type.c_str());
0918     gr_rl->SetMarkerColor(kRed);
0919     gr_rl->SetMarkerStyle(kStar);
0920     gr_rl->SetLineColor(kRed);
0921     gr_rl->GetYaxis()->SetTitle("Raw Asymmetry");
0922     int npoints = gr_sqrt->GetN();
0923     double offset = (bhs_up->bin_edges[1] - bhs_up->bin_edges[0])/20.0;
0924     for (int i=0; i<npoints; i++) {
0925         gr_sqrt->SetPointX(i, gr_sqrt->GetPointX(i) + offset);
0926     }
0927     gr_sqrt->SetMarkerColor(kBlue);
0928     gr_sqrt->SetMarkerStyle(kCircle);
0929     gr_sqrt->SetLineColor(kBlue);
0930     double xmin = std::min(gr_rl->GetXaxis()->GetXmin(), gr_sqrt->GetXaxis()->GetXmin());
0931     double xmax = std::max(gr_rl->GetXaxis()->GetXmax(), gr_sqrt->GetXaxis()->GetXmax());
0932     double ymin = std::min(gr_rl->GetYaxis()->GetXmin(), gr_sqrt->GetYaxis()->GetXmin());
0933     double ymax = std::max(gr_rl->GetYaxis()->GetXmax(), gr_sqrt->GetYaxis()->GetXmax());
0934     /* TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax); */
0935     TH1F* hdummy = pad->DrawFrame(xmin, 2.0*ymin, xmax, 2.0*ymax);
0936     /* TH1F* hdummy = c1->DrawFrame(xmin, -0.1, xmax, 0.1); */
0937     hdummy->SetTitle(title);
0938     hdummy->GetXaxis()->SetTitle("#phi (rad)");
0939     hdummy->GetYaxis()->SetTitle("Raw Asymmetry");
0940     hdummy->GetXaxis()->CenterTitle(true);
0941     hdummy->GetYaxis()->CenterTitle(true);
0942     gStyle->SetOptFit(0);
0943     gStyle->SetOptStat(0);
0944     gr_rl->Draw("p");
0945     /* gr_rl->GetYaxis()->SetLimits(ymin, ymax); */
0946     /* c1->RangeAxis(xmin, ymin, xmax, ymax); */
0947     gr_sqrt->Draw("p same");
0948 
0949     TLine* line = new TLine();
0950     line->SetLineColor(kBlack);
0951     line->DrawLine(xmin, 0, xmax, 0);
0952 
0953     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
0954     leg.AddEntry(gr_rl, "Rel. Lumi", "ep");
0955     leg.AddEntry(gr_sqrt, "Geometric", "ep");
0956     leg.Draw();
0957 
0958     TLatex* textbox = new TLatex(0, 0, "blaahhh");
0959     /* textbox->SetTextSize(textbox->GetTextSize()*1.25); */
0960     /* textbox->SetTextColor(kRed); */
0961     textbox->SetTextAlign(21);
0962     textbox->DrawLatexNDC(0.5, 0.85, range);
0963     /* std::cout << "Greg info: binned_quantity = " << bhs_up->binned_quantity << std::endl; */
0964     /* std::cout << "Greg info: gr num points = " << gr->GetN() << std::endl; */
0965     TLatex* tl_rl = new TLatex;
0966     tl_rl->SetTextColor(kRed);
0967     tl_rl->DrawLatexNDC(0.30, 0.75, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_rl->GetParameter(0), fit_rl->GetParError(0), fit_rl->GetParameter(0)/fit_rl->GetParError(0)));
0968     TLatex* tl_sqrt = new TLatex;
0969     tl_sqrt->SetTextColor(kBlue);
0970     tl_sqrt->DrawLatexNDC(0.30, 0.20, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_sqrt->GetParameter(0), fit_sqrt->GetParError(0), fit_sqrt->GetParameter(0)/fit_sqrt->GetParError(0)));
0971     pad->Update();
0972 
0973     c2->cd();
0974     TH1* hdummy2 = (TH1*)hdummy->Clone("duhhh");
0975     hdummy2->Draw();
0976     gStyle->SetOptFit(0);
0977     gStyle->SetOptStat(0);
0978     gr_rl->Draw("p");
0979     gr_sqrt->Draw("p same");
0980 
0981     line->SetLineColor(kBlack);
0982     line->DrawLine(xmin, 0, xmax, 0);
0983 
0984     leg.Draw();
0985 
0986     textbox->DrawLatexNDC(0.5, 0.85, range);
0987     tl_rl->DrawLatexNDC(0.30, 0.75, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_rl->GetParameter(0), fit_rl->GetParError(0), fit_rl->GetParameter(0)/fit_rl->GetParError(0)));
0988     tl_sqrt->DrawLatexNDC(0.30, 0.20, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_sqrt->GetParameter(0), fit_sqrt->GetParError(0), fit_sqrt->GetParameter(0)/fit_sqrt->GetParError(0)));
0989     c2->Update();
0990     c2->SaveAs(outfilename.c_str());
0991     }
0992     c1->Update();
0993     c1->SaveAs(outfilename.c_str());
0994     delete c1;
0995 }
0996 
0997 void TSSAplotter::PlotAsymsBinnedBlueYellow(std::string type, std::string which, BinnedHistSet* bhs_up_blue, BinnedHistSet* bhs_down_blue, BinnedHistSet* bhs_up_yellow, BinnedHistSet* bhs_down_yellow, std::string outfilename) {
0998     TCanvas* c1 = new TCanvas(which.c_str(), "", 1600, 900);
0999     c1->Divide(4, 2, 0.025, 0.025);
1000     TCanvas* c2 = new TCanvas("c2", "", 1600, 900);
1001     
1002     int nbins = bhs_up_blue->nbins;
1003     for (int i=1; i<=nbins; i++) {
1004     TPad* pad = (TPad*)c1->GetPad(i);
1005     if (!pad) return;
1006     pad->cd();
1007     gStyle->SetPadLeftMargin(0.15);
1008     gROOT->ForceStyle();
1009     c1->Update();
1010 
1011     TH1* phi_up_blue = bhs_up_blue->hist_vec[i];
1012     TH1* phi_down_blue = bhs_down_blue->hist_vec[i];
1013     TGraphErrors* gr_blue = nullptr;
1014     if (type == "rellumi") {
1015         gr_blue = RelLumiGraph(phi_up_blue, phi_down_blue, relLumiBlue);
1016     }
1017     else {
1018         gr_blue = SqrtGraph(phi_up_blue, phi_down_blue);
1019     }
1020     std::pair<double, double> asym_blue = FitGraph(type, gr_blue);
1021     if (asym_blue.first == 0.0) { // needed so compiler doesn't complain about unused variables
1022     }
1023     TF1* fit_blue = gr_blue->GetFunction("asymfit");
1024     fit_blue->SetLineColor(kBlue-2);
1025 
1026     TH1* phi_up_yellow = bhs_up_yellow->hist_vec[i];
1027     TH1* phi_down_yellow = bhs_down_yellow->hist_vec[i];
1028     TGraphErrors* gr_yellow = nullptr;
1029     if (type == "rellumi") {
1030         gr_yellow = RelLumiGraph(phi_up_yellow, phi_down_yellow, relLumiBlue);
1031     }
1032     else {
1033         gr_yellow = SqrtGraph(phi_up_yellow, phi_down_yellow);
1034     }
1035     std::pair<double, double> asym_yellow = FitGraph(type, gr_yellow);
1036     if (asym_yellow.first == 0.0) { // needed so compiler doesn't complain about unused variables
1037     }
1038     TF1* fit_yellow = gr_yellow->GetFunction("asymfit");
1039     fit_yellow->SetLineColor(kOrange-2);
1040 
1041     std::string title_type = "Raw Asymmetry";
1042     char* range;
1043     double lower, upper;
1044     lower = bhs_up_blue->bin_edges[i-1];
1045     upper = bhs_up_blue->bin_edges[i];
1046     range = Form("[%.2f#leq%s<%.2f]", lower, bhs_up_blue->binned_quantity.c_str(), upper);
1047     /* char* title = Form("%s %s %s", which.c_str(), title_type.c_str(), range); */
1048     char* title = Form("%s %s", which.c_str(), title_type.c_str());
1049     gr_blue->SetMarkerColor(kBlue);
1050     gr_blue->SetMarkerStyle(kCircle);
1051     gr_blue->SetLineColor(kBlue);
1052     gr_blue->GetYaxis()->SetTitle("Raw Asymmetry");
1053     int npoints = gr_blue->GetN();
1054     double offset = (bhs_up_blue->bin_edges[1] - bhs_up_blue->bin_edges[0])/20.0;
1055     for (int i=0; i<npoints; i++) {
1056         gr_yellow->SetPointX(i, gr_yellow->GetPointX(i) + offset);
1057     }
1058     gr_yellow->SetMarkerColor(kOrange+1);
1059     gr_yellow->SetMarkerStyle(kStar);
1060     gr_yellow->SetLineColor(kOrange+1);
1061     double xmin = std::min(gr_blue->GetXaxis()->GetXmin(), gr_yellow->GetXaxis()->GetXmin());
1062     double xmax = std::max(gr_blue->GetXaxis()->GetXmax(), gr_yellow->GetXaxis()->GetXmax());
1063     double ymin = std::min(gr_blue->GetYaxis()->GetXmin(), gr_yellow->GetYaxis()->GetXmin());
1064     double ymax = std::max(gr_blue->GetYaxis()->GetXmax(), gr_yellow->GetYaxis()->GetXmax());
1065     /* TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax); */
1066     TH1F* hdummy = pad->DrawFrame(xmin, 2.0*ymin, xmax, 2.0*ymax);
1067     /* TH1F* hdummy = c1->DrawFrame(xmin, -0.1, xmax, 0.1); */
1068     hdummy->SetTitle(title);
1069     hdummy->GetXaxis()->SetTitle("#phi (rad)");
1070     hdummy->GetYaxis()->SetTitle("Raw Asymmetry");
1071     hdummy->GetXaxis()->CenterTitle(true);
1072     hdummy->GetYaxis()->CenterTitle(true);
1073     gStyle->SetOptFit(0);
1074     gStyle->SetOptStat(0);
1075     gr_blue->Draw("p");
1076     /* gr_blue->GetYaxis()->SetLimits(ymin, ymax); */
1077     /* c1->RangeAxis(xmin, ymin, xmax, ymax); */
1078     gr_yellow->Draw("p same");
1079 
1080     TLine* line = new TLine();
1081     line->SetLineColor(kBlack);
1082     line->DrawLine(xmin, 0, xmax, 0);
1083 
1084     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1085     leg.AddEntry(gr_blue, "Blue beam", "ep");
1086     leg.AddEntry(gr_yellow, "Yellow beam", "ep");
1087     leg.Draw();
1088 
1089     TLatex* textbox = new TLatex(0, 0, "blaahhh");
1090     /* textbox->SetTextSize(textbox->GetTextSize()*1.25); */
1091     /* textbox->SetTextColor(kRed); */
1092     textbox->SetTextAlign(21);
1093     textbox->DrawLatexNDC(0.5, 0.85, range);
1094     /* std::cout << "Greg info: binned_quantity = " << bhs_up->binned_quantity << std::endl; */
1095     /* std::cout << "Greg info: gr num points = " << gr->GetN() << std::endl; */
1096     TLatex* tl_blue = new TLatex;
1097     tl_blue->SetTextColor(kBlue);
1098     tl_blue->DrawLatexNDC(0.30, 0.75, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_blue->GetParameter(0), fit_blue->GetParError(0), fit_blue->GetParameter(0)/fit_blue->GetParError(0)));
1099     TLatex* tl_yellow = new TLatex;
1100     tl_yellow->SetTextColor(kOrange+1);
1101     tl_yellow->DrawLatexNDC(0.30, 0.20, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_yellow->GetParameter(0), fit_yellow->GetParError(0), fit_yellow->GetParameter(0)/fit_yellow->GetParError(0)));
1102     pad->Update();
1103 
1104     c2->cd();
1105     TH1* hdummy2 = (TH1*)hdummy->Clone("duhhh");
1106     hdummy2->Draw();
1107     gStyle->SetOptFit(0);
1108     gStyle->SetOptStat(0);
1109     gr_blue->Draw("p");
1110     gr_yellow->Draw("p same");
1111 
1112     line->SetLineColor(kBlack);
1113     line->DrawLine(xmin, 0, xmax, 0);
1114 
1115     leg.Draw();
1116 
1117     textbox->DrawLatexNDC(0.5, 0.85, range);
1118     tl_blue->DrawLatexNDC(0.30, 0.75, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_blue->GetParameter(0), fit_blue->GetParError(0), fit_blue->GetParameter(0)/fit_blue->GetParError(0)));
1119     tl_yellow->DrawLatexNDC(0.30, 0.20, Form("#splitline{#epsilon = %.4f #pm %.4f}{#epsilon/#delta#epsilon = %.4f}", fit_yellow->GetParameter(0), fit_yellow->GetParError(0), fit_yellow->GetParameter(0)/fit_yellow->GetParError(0)));
1120     c2->Update();
1121     c2->SaveAs(outfilename.c_str());
1122     }
1123     c1->Update();
1124     c1->SaveAs(outfilename.c_str());
1125     delete c1;
1126 }
1127 
1128 std::pair<double, double> TSSAplotter::FitGraph(std::string type, TGraphErrors* gr) {
1129     double fit_lower, fit_upper;
1130     if (type == "sqrt") {
1131     fit_lower = -PI/2.0;
1132     fit_upper = PI/2.0;
1133     }
1134     else {
1135     fit_lower = -PI;
1136     fit_upper = PI;
1137     }
1138 
1139     TF1* fit = nullptr;
1140     if (type == "sqrt" || type == "rellumi") {
1141     fit = new TF1("asymfit", "[0]*sin(x - [1])", fit_lower, fit_upper);
1142     /* fit->SetParLimits(0, 0.00001, 0.4); */
1143     fit->SetParLimits(0, -0.4, 0.4);
1144     fit->SetParameter(0, 0.001);
1145     fit->SetParLimits(1, PI-0.2, PI+0.2);
1146     /* fit->SetParLimits(1, 0.0, 2*PI); */
1147     fit->SetParameter(1, 3.141592);
1148     fit->SetParName(0, "#epsilon");
1149     fit->SetParName(1, "#phi_{0}");
1150     }
1151     else {
1152     std::cout << "Invalid type option " << type << "; exiting!" << std::endl;
1153     return std::pair<double, double>();
1154     }
1155 
1156     gr->Fit(fit, "Q");
1157     std::pair<double, double> ret;
1158     ret.first = fit->GetParameter(0);
1159     ret.second = fit->GetParError(0);
1160     return ret;
1161 }
1162 
1163 TGraphErrors* TSSAplotter::AmplitudeBinnedGraph(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up, BinnedHistSet* bhs_down) {
1164 
1165     // Check xlabel to see which mean x values to use
1166     std::vector<double> xvalues;
1167     if (xlabel == "pT_pi0") {
1168     xlabel = "p_{T} (GeV)";
1169     xvalues = pTmeans_pi0;
1170     }
1171     else if (xlabel == "pT_eta") {
1172     xlabel = "p_{T} (GeV)";
1173     xvalues = pTmeans_eta;
1174     }
1175     else if (xlabel == "xF") {
1176     xlabel = "x_{F}";
1177     xvalues = xFmeans;
1178     }
1179     else if (xlabel == "eta") {
1180     xlabel = "#eta";
1181     xvalues = etameans;
1182     }
1183     else if (xlabel == "vtxz") {
1184     xlabel = "z_{vtx}";
1185     xvalues = vtxzmeans;
1186     }
1187     else {
1188     std::cout << "In AmplitudeBinnedGraph(): unrecognized xlabel = " << xlabel << "; exiting!" << std::endl;
1189     return nullptr;
1190     }
1191 
1192     int nbins = bhs_up->nbins;
1193     /* std::vector<double> bin_centers = bhs_up->getBinCenters(bhs_up->bin_edges); */
1194     TGraphErrors* amplitudes = new TGraphErrors(nbins);
1195     amplitudes->SetTitle(title.c_str());
1196     amplitudes->GetXaxis()->SetTitle(xlabel.c_str());
1197     std::cout << "In AmplitudeBinnedGraph(): xlabel = " << xlabel << std::endl;
1198     amplitudes->GetXaxis()->CenterTitle(1);
1199     std::string ylabel;
1200     if (type == "sqrt") ylabel = "Raw Asymmetry (Square Root)";
1201     else ylabel = "Raw Asymmetry (Rel. Lumi.)";
1202     amplitudes->GetYaxis()->SetTitle(ylabel.c_str());
1203     amplitudes->GetYaxis()->CenterTitle(1);
1204 
1205     for (int i=1; i<=nbins; i++) {
1206     TH1* phi_up = bhs_up->hist_vec[i];
1207     TH1* phi_down = bhs_down->hist_vec[i];
1208     TGraphErrors* gr = nullptr;
1209     if (type == "sqrt") {
1210         gr = SqrtGraph(phi_up, phi_down);
1211     }
1212     else {
1213         bool blue = (title.find("Blue") != std::string::npos);
1214         if (blue) {
1215         gr = RelLumiGraph(phi_up, phi_down, relLumiBlue);
1216         }
1217         else {
1218         gr = RelLumiGraph(phi_up, phi_down, relLumiYellow);
1219         }
1220     }
1221     std::pair<double, double> ampl_err = FitGraph(type, gr);
1222     /* amplitudes->SetPoint(i-1, bin_centers[i-1], ampl_err.first); */
1223     /* double xrange = bhs_up->bin_edges[i] - bhs_up->bin_edges[i-1]; */
1224     /* amplitudes->SetPointError(i-1, xrange/2.0, ampl_err.second); */
1225     amplitudes->SetPoint(i-1, xvalues[i-1], ampl_err.first);
1226     amplitudes->SetPointError(i-1, 0, ampl_err.second);
1227     }
1228     return amplitudes;
1229 }
1230 
1231 void TSSAplotter::PlotAmplitudes(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up, BinnedHistSet* bhs_down, std::string outfilename) {
1232     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1233     c1->cd();
1234     TGraphErrors* gr = AmplitudeBinnedGraph(type, title, xlabel, bhs_up, bhs_down);
1235     gr->SetMarkerColor(kBlack);
1236     gr->SetMarkerStyle(kCircle);
1237     gr->SetLineColor(kBlack);
1238     gr->Draw("ap");
1239     c1->SaveAs(outfilename.c_str());
1240 }
1241 
1242 void TSSAplotter::CompareAmplitudes(std::string title, std::string xlabel, BinnedHistSet* bhs_up, BinnedHistSet* bhs_down, std::string outfilename) {
1243     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1244     c1->cd();
1245     TGraphErrors* gr_rl = AmplitudeBinnedGraph("rellumi", title, xlabel, bhs_up, bhs_down);
1246     gr_rl->SetMarkerColor(kRed);
1247     gr_rl->SetMarkerStyle(kStar);
1248     gr_rl->SetLineColor(kRed);
1249     gr_rl->GetYaxis()->SetTitle("Raw Asymmetry");
1250 
1251     TGraphErrors* gr_sqrt = AmplitudeBinnedGraph("sqrt", title, xlabel, bhs_up, bhs_down);
1252     int npoints = gr_sqrt->GetN();
1253     double offset = (bhs_up->bin_edges[1] - bhs_up->bin_edges[0])/20.0;
1254     for (int i=0; i<npoints; i++) {
1255     gr_sqrt->SetPointX(i, gr_sqrt->GetPointX(i) + offset);
1256     }
1257     gr_sqrt->SetMarkerColor(kBlue);
1258     gr_sqrt->SetMarkerStyle(kCircle);
1259     gr_sqrt->SetLineColor(kBlue);
1260 
1261     double xmin = std::min(gr_rl->GetXaxis()->GetXmin(), gr_sqrt->GetXaxis()->GetXmin());
1262     double xmax = std::max(gr_rl->GetXaxis()->GetXmax(), gr_sqrt->GetXaxis()->GetXmax());
1263     double ymin = std::min(gr_rl->GetYaxis()->GetXmin(), gr_sqrt->GetYaxis()->GetXmin());
1264     double ymax = std::max(gr_rl->GetYaxis()->GetXmax(), gr_sqrt->GetYaxis()->GetXmax());
1265     /* std::cout << "Sqrt min = " << gr_sqrt->GetYaxis()->GetXmin() << ";\tymin = " << ymin << std::endl; */
1266     TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax);
1267     hdummy->SetTitle(title.c_str());
1268     hdummy->GetXaxis()->SetTitle(xlabel.c_str());
1269     hdummy->GetYaxis()->SetTitle("Raw Asymmetry");
1270     hdummy->GetXaxis()->CenterTitle(true);
1271     hdummy->GetYaxis()->CenterTitle(true);
1272     gr_rl->Draw("p");
1273     /* gr_rl->GetYaxis()->SetLimits(ymin, ymax); */
1274     /* c1->RangeAxis(xmin, ymin, xmax, ymax); */
1275     gr_sqrt->Draw("p same");
1276 
1277     TLine* line = new TLine();
1278     line->SetLineColor(kBlack);
1279     line->DrawLine(xmin, 0, xmax, 0);
1280 
1281     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1282     leg.AddEntry(gr_rl, "Rel. Lumi", "ep");
1283     leg.AddEntry(gr_sqrt, "Geometric", "ep");
1284     leg.Draw();
1285 
1286     c1->SaveAs(outfilename.c_str());
1287 }
1288 
1289 void TSSAplotter::CompareAmplitudesEtaBinned(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up_lowEta, BinnedHistSet* bhs_down_lowEta, BinnedHistSet* bhs_up_highEta, BinnedHistSet* bhs_down_highEta, std::string outfilename) {
1290     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1291     c1->cd();
1292     TGraphErrors* gr_lowEta = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_lowEta, bhs_down_lowEta);
1293     gr_lowEta->SetMarkerColor(kRed);
1294     gr_lowEta->SetMarkerStyle(kStar);
1295     gr_lowEta->SetLineColor(kRed);
1296     gr_lowEta->GetYaxis()->SetTitle("Raw Asymmetry");
1297 
1298     TGraphErrors* gr_highEta = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_highEta, bhs_down_highEta);
1299     int npoints = gr_highEta->GetN();
1300     double offset = (bhs_up_lowEta->bin_edges[1] - bhs_up_lowEta->bin_edges[0])/20.0;
1301     for (int i=0; i<npoints; i++) {
1302     gr_highEta->SetPointX(i, gr_highEta->GetPointX(i) + offset);
1303     }
1304     gr_highEta->SetMarkerColor(kBlue);
1305     gr_highEta->SetMarkerStyle(kCircle);
1306     gr_highEta->SetLineColor(kBlue);
1307 
1308     double xmin = std::min(gr_lowEta->GetXaxis()->GetXmin(), gr_highEta->GetXaxis()->GetXmin());
1309     double xmax = std::max(gr_lowEta->GetXaxis()->GetXmax(), gr_highEta->GetXaxis()->GetXmax());
1310     double ymin = std::min(gr_lowEta->GetYaxis()->GetXmin(), gr_highEta->GetYaxis()->GetXmin());
1311     double ymax = std::max(gr_lowEta->GetYaxis()->GetXmax(), gr_highEta->GetYaxis()->GetXmax());
1312 
1313     TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax);
1314     hdummy->SetTitle(title.c_str());
1315     hdummy->GetXaxis()->SetTitle(xlabel.c_str());
1316     hdummy->GetYaxis()->SetTitle("Raw Asymmetry");
1317     hdummy->GetXaxis()->CenterTitle(true);
1318     hdummy->GetYaxis()->CenterTitle(true);
1319     gr_lowEta->Draw("p");
1320     /* gr_rl->GetYaxis()->SetLimits(ymin, ymax); */
1321     /* c1->RangeAxis(xmin, ymin, xmax, ymax); */
1322     gr_highEta->Draw("p same");
1323 
1324     TLine* line = new TLine();
1325     line->SetLineColor(kBlack);
1326     line->DrawLine(xmin, 0, xmax, 0);
1327 
1328     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1329     leg.AddEntry(gr_lowEta, "|#eta| < 0.35", "ep");
1330     leg.AddEntry(gr_highEta, "|#eta| > 0.35", "ep");
1331     leg.Draw();
1332 
1333     c1->SaveAs(outfilename.c_str());
1334 }
1335 
1336 void TSSAplotter::CompareAmplitudesBlueYellow(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up_blue, BinnedHistSet* bhs_down_blue, BinnedHistSet* bhs_up_yellow, BinnedHistSet* bhs_down_yellow, std::string outfilename) {
1337     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1338     c1->cd();
1339     TGraphErrors* gr_blue = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_blue, bhs_down_blue);
1340     gr_blue->SetMarkerColor(kBlue);
1341     gr_blue->SetMarkerStyle(kCircle);
1342     gr_blue->SetLineColor(kBlue);
1343     gr_blue->GetYaxis()->SetTitle("Raw Asymmetry");
1344 
1345     TGraphErrors* gr_yellow = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_yellow, bhs_down_yellow);
1346     int npoints = gr_yellow->GetN();
1347     double offset = (bhs_up_blue->bin_edges[1] - bhs_up_blue->bin_edges[0])/20.0;
1348     for (int i=0; i<npoints; i++) {
1349     gr_yellow->SetPointX(i, gr_yellow->GetPointX(i) + offset);
1350     }
1351     gr_yellow->SetMarkerColor(kOrange+1);
1352     gr_yellow->SetMarkerStyle(kStar);
1353     gr_yellow->SetLineColor(kOrange+1);
1354 
1355     double xmin = std::min(gr_blue->GetXaxis()->GetXmin(), gr_yellow->GetXaxis()->GetXmin());
1356     double xmax = std::max(gr_blue->GetXaxis()->GetXmax(), gr_yellow->GetXaxis()->GetXmax());
1357     double ymin = std::min(gr_blue->GetYaxis()->GetXmin(), gr_yellow->GetYaxis()->GetXmin());
1358     double ymax = std::max(gr_blue->GetYaxis()->GetXmax(), gr_yellow->GetYaxis()->GetXmax());
1359 
1360     TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax);
1361     hdummy->SetTitle(title.c_str());
1362     hdummy->GetXaxis()->SetTitle(xlabel.c_str());
1363     hdummy->GetYaxis()->SetTitle("Raw Asymmetry");
1364     hdummy->GetXaxis()->CenterTitle(true);
1365     hdummy->GetYaxis()->CenterTitle(true);
1366     gr_blue->Draw("p");
1367     /* gr_rl->GetYaxis()->SetLimits(ymin, ymax); */
1368     /* c1->RangeAxis(xmin, ymin, xmax, ymax); */
1369     gr_yellow->Draw("p same");
1370 
1371     TLine* line = new TLine();
1372     line->SetLineColor(kBlack);
1373     line->DrawLine(xmin, 0, xmax, 0);
1374 
1375     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1376     leg.AddEntry(gr_blue, "Blue beam", "ep");
1377     leg.AddEntry(gr_yellow, "Yellow beam", "ep");
1378     leg.Draw();
1379 
1380     c1->SaveAs(outfilename.c_str());
1381 }
1382 
1383 std::pair<double, double> TSSAplotter::GetCorrectedAsymErr(double A_sig, double A_sig_err, double A_bkgr, double A_bkgr_err, double r, double r_err, double pol) {
1384     double asym = (1.0/pol)*(A_sig - r*A_bkgr)/(1.0 - r);
1385     double err1 = A_sig_err*A_sig_err;
1386     double err2 = A_bkgr_err*A_bkgr_err;
1387     double err3 = ((A_sig - A_bkgr)*(A_sig - A_bkgr))/((1-r)*(1-r));
1388     double err = (1.0/pol) * (1/(1-r)) * sqrt( err1 + (r*r*err2) + (r_err*r_err*err3) );
1389     std::pair<double, double> ret;
1390     ret.first = asym;
1391     ret.second = err;
1392     return ret;
1393 }
1394 
1395 TGraphErrors* TSSAplotter::CorrectedAmplitudesGraph(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up_sig, BinnedHistSet* bhs_down_sig, BinnedHistSet* bhs_up_bkgr, BinnedHistSet* bhs_down_bkgr) {
1396     TGraphErrors* gr_sig = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_sig, bhs_down_sig);
1397     TGraphErrors* gr_bkgr = AmplitudeBinnedGraph(type, title, xlabel, bhs_up_bkgr, bhs_down_bkgr);
1398 
1399     int npoints = gr_sig->GetN();
1400     TGraphErrors* gr_corrected = new TGraphErrors(npoints);
1401     gr_corrected->SetTitle(title.c_str());
1402     gr_corrected->GetXaxis()->SetTitle(xlabel.c_str());
1403     gr_corrected->GetXaxis()->CenterTitle(1);
1404     std::string ylabel;
1405     if (type == "sqrt") ylabel = "Corrected Asymmetry (Square Root)";
1406     else ylabel = "Corrected Asymmetry (Rel. Lumi.)";
1407     gr_corrected->GetYaxis()->SetTitle(ylabel.c_str());
1408     gr_corrected->GetYaxis()->CenterTitle(1);
1409 
1410     double r, r_err, pol;
1411     bool blue = (title.find("Blue") != std::string::npos);
1412     if (blue) {
1413     pol = polBlue;
1414     }
1415     else {
1416     pol = polYellow;
1417     }
1418     bool pi0 = (title.find("pi") != std::string::npos);
1419 
1420     for (int i=0; i<npoints; i++) {
1421     if (pi0) {
1422         r = pi0_bkgr_fractions[i];
1423         r_err = pi0_bkgr_fraction_errs[i];
1424     }
1425     else {
1426         r = eta_bkgr_fractions[i];
1427         r_err =eta_bkgr_fraction_errs[i];
1428     }
1429     std::pair<double, double> asym_err = GetCorrectedAsymErr(gr_sig->GetPointY(i), gr_sig->GetErrorY(i), gr_bkgr->GetPointY(i), gr_bkgr->GetErrorY(i), r, r_err, pol);
1430     gr_corrected->SetPoint(i, gr_sig->GetPointX(i), asym_err.first);
1431     gr_corrected->SetPointError(i, gr_sig->GetErrorX(i), asym_err.second);
1432     }
1433 
1434     return gr_corrected;
1435 }
1436 
1437 void TSSAplotter::PlotCorrectedAmplitudes(std::string type, std::string title, std::string xlabel, BinnedHistSet* bhs_up_sig, BinnedHistSet* bhs_down_sig, BinnedHistSet* bhs_up_bkgr, BinnedHistSet* bhs_down_bkgr, std::string outfilename) {
1438     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1439     c1->cd();
1440     TGraphErrors* gr = CorrectedAmplitudesGraph(type, title, xlabel, bhs_up_sig, bhs_down_sig, bhs_up_bkgr, bhs_down_bkgr);
1441     gr->SetMarkerColor(kBlack);
1442     gr->SetMarkerStyle(kCircle);
1443     gr->SetLineColor(kBlack);
1444     gr->Draw("ap");
1445     TLine* line = new TLine();
1446     double xmin = gr->GetXaxis()->GetXmin();
1447     double xmax = gr->GetXaxis()->GetXmax();
1448     line->SetLineColor(kBlack);
1449     line->DrawLine(xmin, 0, xmax, 0);
1450     c1->SaveAs(outfilename.c_str());
1451 }
1452 
1453 TGraphErrors* TSSAplotter::AveragedGraphs(TGraphErrors* gr_blue, TGraphErrors* gr_yellow) {
1454     int npoints = gr_blue->GetN();
1455     TGraphErrors* gr_averaged = new TGraphErrors(npoints);
1456     for (int i=0; i<npoints; i++) {
1457     double x = gr_blue->GetPointX(i);
1458     double y_blue = gr_blue->GetPointY(i);
1459     double y_yellow = gr_yellow->GetPointY(i);
1460     double y = (y_blue + y_yellow)/2.0;
1461     double yerr_blue = gr_blue->GetErrorY(i);
1462     double yerr_yellow = gr_yellow->GetErrorY(i);
1463     double yerr = 0.5*sqrt( (yerr_blue*yerr_blue) + (yerr_yellow*yerr_yellow) );
1464     gr_averaged->SetPoint(i, x, y);
1465     gr_averaged->SetPointError(i, 0, yerr);
1466     }
1467     return gr_averaged;
1468 }
1469 
1470 void TSSAplotter::PlotAveragedAmplitudes(std::string type, std::string title, std::string xlabel, std::string outfilename, BinnedHistSet* bhs_up_sig_blue, BinnedHistSet* bhs_down_sig_blue, BinnedHistSet* bhs_up_bkgr_blue, BinnedHistSet* bhs_down_bkgr_blue, BinnedHistSet* bhs_up_sig_yellow, BinnedHistSet* bhs_down_sig_yellow, BinnedHistSet* bhs_up_bkgr_yellow, BinnedHistSet* bhs_down_bkgr_yellow) {
1471     TGraphErrors* gr_blue = CorrectedAmplitudesGraph(type, "Blue", xlabel, bhs_up_sig_blue, bhs_down_sig_blue, bhs_up_bkgr_blue, bhs_down_bkgr_blue);
1472     TGraphErrors* gr_yellow = CorrectedAmplitudesGraph(type, "Yellow", xlabel, bhs_up_sig_yellow, bhs_down_sig_yellow, bhs_up_bkgr_yellow, bhs_down_bkgr_yellow);
1473     TGraphErrors* gr_averaged = AveragedGraphs(gr_blue, gr_yellow);
1474     PlotGraph(gr_averaged, title, xlabel, "Corrected A_{N}", outfilename);
1475 }
1476 
1477 void TSSAplotter::PlotGraph(TGraphErrors* gr, std::string title, std::string xlabel, std::string ylabel, std::string outfilename, bool inset, bool ylines) {
1478     gr->SetTitle(title.c_str());
1479     gr->GetXaxis()->SetTitle(xlabel.c_str());
1480     gr->GetXaxis()->CenterTitle(true);
1481     gr->GetYaxis()->SetTitle(ylabel.c_str());
1482     gr->GetYaxis()->CenterTitle(true);
1483     gr->SetMarkerStyle(kCircle);
1484     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1485     c1->cd();
1486     gr->Draw("ap");
1487     TLine* line = new TLine();
1488     double xmin = gr->GetXaxis()->GetXmin();
1489     double xmax = gr->GetXaxis()->GetXmax();
1490     line->SetLineColor(kBlack);
1491     line->DrawLine(xmin, 0, xmax, 0);
1492     if (ylines) {
1493     TLine* line2 = new TLine();
1494     line2->SetLineStyle(kDashed);
1495     line2->DrawLine(xmin, -1.0, xmax, -1.0);
1496     line2->DrawLine(xmin, 1.0, xmax, 1.0);
1497     }
1498 
1499     if (inset) {
1500     TGraphErrors* gr_inset = new TGraphErrors(5);
1501     for (int i=0; i<5; i++) {
1502         double x = gr->GetPointX(i);
1503         double y = gr->GetPointY(i);
1504         double yerr = gr->GetErrorY(i);
1505         gr_inset->SetPoint(i, x, y);
1506         gr_inset->SetPointError(i, 0.0, yerr);
1507     }
1508     gr_inset->SetMarkerStyle(kCircle);
1509     gr_inset->SetTitle("");
1510     TPad* pad = new TPad("pad", "pad", 0.2, 0.5, 0.6, 0.7, -1, 0);
1511     pad->Draw();
1512     pad->cd();
1513     double xmin_inset = gr_inset->GetXaxis()->GetXmin();
1514     double xmax_inset = gr_inset->GetXaxis()->GetXmax();
1515     double ymin_inset = -0.01;
1516     double ymax_inset = 0.01;
1517     TH1F* h_frame = pad->DrawFrame(xmin_inset, ymin_inset, xmax_inset, ymax_inset);
1518     h_frame->GetXaxis()->SetLabelSize(0.15);
1519     /* h_frame->GetXaxis()->SetTickLength(0.5); */
1520     h_frame->GetYaxis()->SetLabelSize(0.15);
1521     h_frame->GetYaxis()->SetNdivisions(-504);
1522     gr_inset->Draw("p");
1523     line->DrawLine(xmin_inset, 0, xmax_inset, 0);
1524     pad->Modified();
1525     }
1526     c1->Modified();
1527     c1->SaveAs(outfilename.c_str());
1528 }
1529 
1530 void TSSAplotter::CompareGraphs(TGraphErrors* gr1, TGraphErrors* gr2, std::string title, std::string xlabel, std::string ylabel, std::string leglabel1, std::string leglabel2, std::string outfilename, bool inset) {
1531     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1532     c1->cd();
1533     gr1->SetMarkerColor(kRed);
1534     gr1->SetMarkerStyle(kStar);
1535     gr1->SetLineColor(kRed);
1536     /* gr1->GetYaxis()->SetTitle("Raw Asymmetry"); */
1537 
1538     int npoints = gr2->GetN();
1539     double offset = (gr2->GetPointX(1) - gr2->GetPointX(0))/20.0;
1540     for (int i=0; i<npoints; i++) {
1541     gr2->SetPointX(i, gr2->GetPointX(i) + offset);
1542     }
1543     gr2->SetMarkerColor(kBlue);
1544     gr2->SetMarkerStyle(kCircle);
1545     gr2->SetLineColor(kBlue);
1546 
1547     double xmin = std::min(gr1->GetXaxis()->GetXmin(), gr2->GetXaxis()->GetXmin());
1548     double xmax = std::max(gr1->GetXaxis()->GetXmax(), gr2->GetXaxis()->GetXmax());
1549     double ymin = std::min(gr1->GetYaxis()->GetXmin(), gr2->GetYaxis()->GetXmin());
1550     double ymax = std::max(gr1->GetYaxis()->GetXmax(), gr2->GetYaxis()->GetXmax());
1551     TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax);
1552     hdummy->SetTitle(title.c_str());
1553     hdummy->GetXaxis()->SetTitle(xlabel.c_str());
1554     hdummy->GetYaxis()->SetTitle(ylabel.c_str());
1555     hdummy->GetXaxis()->CenterTitle(true);
1556     hdummy->GetYaxis()->CenterTitle(true);
1557     gr1->Draw("p");
1558     gr2->Draw("p same");
1559 
1560     TLine* line = new TLine();
1561     line->SetLineColor(kBlack);
1562     line->DrawLine(xmin, 0, xmax, 0);
1563 
1564     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1565     leg.AddEntry(gr1, leglabel1.c_str(), "ep");
1566     leg.AddEntry(gr2, leglabel2.c_str(), "ep");
1567     leg.Draw();
1568 
1569     if (inset) {
1570     TGraphErrors* gr1_inset = new TGraphErrors(5);
1571     TGraphErrors* gr2_inset = new TGraphErrors(5);
1572     /* TGraphErrors* gr2_inset = (TGraphErrors*)gr2->Clone("gr2_inset"); */
1573     for (int i=0; i<5; i++) {
1574         double x1 = gr1->GetPointX(i);
1575         double y1 = gr1->GetPointY(i);
1576         double y1err = gr1->GetErrorY(i);
1577         gr1_inset->SetPoint(i, x1, y1);
1578         gr1_inset->SetPointError(i, 0.0, y1err);
1579         double x2 = gr2->GetPointX(i);
1580         double y2 = gr2->GetPointY(i);
1581         double x2err = gr2->GetErrorX(i);
1582         double y2err = gr2->GetErrorY(i);
1583         gr2_inset->SetPoint(i, x2, y2);
1584         gr2_inset->SetPointError(i, x2err, y2err);
1585     }
1586     gr1_inset->SetMarkerColor(kRed);
1587     gr1_inset->SetMarkerStyle(kStar);
1588     gr1_inset->SetLineColor(kRed);
1589     gr2_inset->SetMarkerColor(kBlue);
1590     gr2_inset->SetMarkerStyle(kCircle);
1591     gr2_inset->SetLineColor(kBlue);
1592 
1593     TPad* pad = new TPad("pad", "pad", 0.2, 0.5, 0.6, 0.7, -1, 0);
1594     pad->Draw();
1595     pad->cd();
1596     double xmin_inset = gr1_inset->GetXaxis()->GetXmin();
1597     double xmax_inset = gr1_inset->GetXaxis()->GetXmax();
1598     double ymin_inset = -0.01;
1599     double ymax_inset = 0.01;
1600     TH1F* h_frame = pad->DrawFrame(xmin_inset, ymin_inset, xmax_inset, ymax_inset);
1601     h_frame->GetXaxis()->SetLabelSize(0.15);
1602     /* h_frame->GetXaxis()->SetTickLength(0.5); */
1603     h_frame->GetYaxis()->SetLabelSize(0.15);
1604     h_frame->GetYaxis()->SetNdivisions(-504);
1605     gr1_inset->Draw("p");
1606     gr2_inset->Draw("p same");
1607     line->DrawLine(xmin_inset, 0, xmax_inset, 0);
1608     pad->Modified();
1609     }
1610     c1->SaveAs(outfilename.c_str());
1611     // revert gr2 x values to original
1612     for (int i=0; i<npoints; i++) {
1613     gr2->SetPointX(i, gr2->GetPointX(i) - offset);
1614     }
1615 }
1616 
1617 void TSSAplotter::CompareGraphs(TGraphErrors* gr1, TGraphAsymmErrors* gr2, std::string title, std::string xlabel, std::string ylabel, std::string leglabel1, std::string leglabel2, std::string outfilename, bool inset) {
1618     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1619     c1->cd();
1620     gr1->SetMarkerColor(kRed);
1621     gr1->SetMarkerStyle(kStar);
1622     gr1->SetLineColor(kRed);
1623     /* gr1->GetYaxis()->SetTitle("Raw Asymmetry"); */
1624 
1625     int npoints = gr2->GetN();
1626     double offset = (gr2->GetPointX(1) - gr2->GetPointX(0))/20.0;
1627     for (int i=0; i<npoints; i++) {
1628     gr2->SetPointX(i, gr2->GetPointX(i) + offset);
1629     }
1630     gr2->SetMarkerColor(kBlue);
1631     gr2->SetMarkerStyle(kCircle);
1632     gr2->SetLineColor(kBlue);
1633 
1634     double xmin = std::min(gr1->GetXaxis()->GetXmin(), gr2->GetXaxis()->GetXmin());
1635     double xmax = std::max(gr1->GetXaxis()->GetXmax(), gr2->GetXaxis()->GetXmax());
1636     double ymin = std::min(gr1->GetYaxis()->GetXmin(), gr2->GetYaxis()->GetXmin());
1637     double ymax = std::max(gr1->GetYaxis()->GetXmax(), gr2->GetYaxis()->GetXmax());
1638     ymax *= 2.0;
1639     TH1F* hdummy = c1->DrawFrame(xmin, ymin, xmax, ymax);
1640     hdummy->SetTitle(title.c_str());
1641     hdummy->GetXaxis()->SetTitle(xlabel.c_str());
1642     hdummy->GetYaxis()->SetTitle(ylabel.c_str());
1643     hdummy->GetXaxis()->CenterTitle(true);
1644     hdummy->GetYaxis()->CenterTitle(true);
1645     gr1->Draw("p");
1646     gr2->Draw("p same");
1647 
1648     TLine* line = new TLine();
1649     line->SetLineColor(kBlack);
1650     line->DrawLine(xmin, 0, xmax, 0);
1651 
1652     TLegend leg(0.75, 0.75, 0.9, 0.9, "", "nbNDC");
1653     leg.AddEntry(gr1, leglabel1.c_str(), "ep");
1654     leg.AddEntry(gr2, leglabel2.c_str(), "ep");
1655     leg.Draw();
1656 
1657     if (inset) {
1658     TGraphErrors* gr1_inset = new TGraphErrors(5);
1659     TGraphErrors* gr2_inset = new TGraphErrors(5);
1660     /* TGraphErrors* gr2_inset = (TGraphErrors*)gr2->Clone("gr2_inset"); */
1661     for (int i=0; i<5; i++) {
1662         double x = gr1->GetPointX(i);
1663         double y1 = gr1->GetPointY(i);
1664         double y1err = gr1->GetErrorY(i);
1665         gr1_inset->SetPoint(i, x, y1);
1666         gr1_inset->SetPointError(i, 0.0, y1err);
1667         double x2 = gr2->GetPointX(i);
1668         double x2err = gr2->GetErrorX(i);
1669         double y2 = gr2->GetPointY(i);
1670         double y2err = gr2->GetErrorY(i);
1671         gr2_inset->SetPoint(i, x2, y2);
1672         gr2_inset->SetPointError(i, x2err, y2err);
1673     }
1674     gr1_inset->SetMarkerColor(kRed);
1675     gr1_inset->SetMarkerStyle(kStar);
1676     gr1_inset->SetLineColor(kRed);
1677     gr2_inset->SetMarkerColor(kBlue);
1678     gr2_inset->SetMarkerStyle(kCircle);
1679     gr2_inset->SetLineColor(kBlue);
1680 
1681     TPad* pad = new TPad("pad", "pad", 0.2, 0.5, 0.6, 0.75, -1, 0);
1682     pad->Draw();
1683     pad->cd();
1684     double xmin_inset = gr1_inset->GetXaxis()->GetXmin();
1685     double xmax_inset = gr1_inset->GetXaxis()->GetXmax();
1686     double ymin_inset = -0.005;
1687     double ymax_inset = 0.005;
1688     TH1F* h_frame = pad->DrawFrame(xmin_inset, ymin_inset, xmax_inset, ymax_inset);
1689     h_frame->GetXaxis()->SetLabelSize(0.15);
1690     /* h_frame->GetXaxis()->SetTickLength(0.5); */
1691     h_frame->GetYaxis()->SetLabelSize(0.15);
1692     h_frame->GetYaxis()->SetNdivisions(-504);
1693     gr1_inset->Draw("p");
1694     gr2_inset->Draw("p same");
1695     line->DrawLine(xmin_inset, 0, xmax_inset, 0);
1696     pad->Modified();
1697     }
1698     c1->SaveAs(outfilename.c_str());
1699     // revert gr2 x values to original
1700     for (int i=0; i<npoints; i++) {
1701     gr2->SetPointX(i, gr2->GetPointX(i) - offset);
1702     }
1703 }
1704 
1705 TGraphErrors* TSSAplotter::PolCorrectedGraph(TGraphErrors* gr, double pol) {
1706     int npoints = gr->GetN();
1707     TGraphErrors* out = new TGraphErrors(npoints);
1708     for (int i=0; i<npoints; i++) {
1709     double x = gr->GetPointX(i);
1710     double y = gr->GetPointY(i);
1711     double yerr = gr->GetErrorY(i);
1712     out->SetPoint(i, x, y/pol);
1713     out->SetPointError(i, 0, yerr/pol);
1714     }
1715     return out;
1716 }
1717 
1718 void TSSAplotter::PlotBackgroundFractions(std::string outfilename) {
1719     int nbins = pi0_bkgr_fractions.size();
1720     TGraphErrors* gr_pi0_frac = new TGraphErrors(nbins - 1);
1721     TGraphErrors* gr_eta_frac = new TGraphErrors(nbins);
1722     for (int i=0; i<nbins; i++) {
1723     double x_pi0 = (pTbins_pi0[i+1] + pTbins_pi0[i])/2.0;
1724     double x_eta = (pTbins_eta[i+1] + pTbins_eta[i])/2.0;
1725     if (i < (nbins-1)) {
1726         gr_pi0_frac->SetPoint(i, x_pi0, pi0_bkgr_fractions[i]);
1727         gr_pi0_frac->SetPointError(i, 0, pi0_bkgr_fraction_errs[i]);
1728     }
1729     gr_eta_frac->SetPoint(i, x_eta, eta_bkgr_fractions[i]);
1730     gr_eta_frac->SetPointError(i, 0, eta_bkgr_fraction_errs[i]);
1731     }
1732     CompareGraphs(gr_pi0_frac, gr_eta_frac, "#pi^{0} and #eta Background Fractions", "p_{T} (GeV)", "Background Fraction", "#pi^{0}", "#eta", outfilename);
1733 }
1734 
1735 TGraphErrors* TSSAplotter::tTestGraph(std::string type, TGraphErrors* gr1, TGraphErrors* gr2) {
1736     int npoints = gr1->GetN();
1737     TGraphErrors* gr_t = new TGraphErrors(npoints);
1738     for (int i=0; i<npoints; i++) {
1739     double x = gr1->GetPointX(i);
1740     double y1 = gr1->GetPointY(i);
1741     double y2 = gr2->GetPointY(i);
1742     double y1err = gr1->GetErrorY(i);
1743     double y2err = gr2->GetErrorY(i);
1744     double num = y1 - y2;
1745     double denom = 1.0;
1746     if (type == "indep") {
1747         denom = sqrt(y1err*y1err + y2err*y2err);
1748     }
1749     else if (type == "dep") {
1750         denom = sqrt(abs(y1err*y1err - y2err*y2err));
1751     }
1752     else {
1753         std::cout << "In TSSAplotter::tTestGraph: invalid type " << type << "; returning nullptr" << std::endl;
1754         return nullptr;
1755     }
1756     double t = num/denom;
1757     gr_t->SetPoint(i, x, t);
1758     }
1759     return gr_t;
1760 }
1761 
1762 TGraphErrors* TSSAplotter::tTestZeroGraph(TGraphErrors* gr1) {
1763     int npoints = gr1->GetN();
1764     TGraphErrors* gr_t = new TGraphErrors(npoints);
1765     for (int i=0; i<npoints; i++) {
1766     double x = gr1->GetPointX(i);
1767     double y1 = gr1->GetPointY(i);
1768     double y1err = gr1->GetErrorY(i);
1769     double t = y1/y1err;
1770     gr_t->SetPoint(i, x, t);
1771     }
1772     return gr_t;
1773 }
1774 
1775 void TSSAplotter::PlottTest(TGraphErrors* gr, std::string title, std::string xlabel, std::string outfilename, bool limit_yrange) {
1776     gr->SetTitle(title.c_str());
1777     gr->GetXaxis()->SetTitle(xlabel.c_str());
1778     gr->GetXaxis()->CenterTitle(true);
1779     gr->GetYaxis()->SetTitle("t Statistic");
1780     gr->GetYaxis()->CenterTitle(true);
1781     if (limit_yrange) gr->GetYaxis()->SetRangeUser(-2.5, 2.5);
1782     gr->SetMarkerStyle(kCircle);
1783     TCanvas* c1 = new TCanvas(title.c_str(), "", 1600, 900);
1784     c1->cd();
1785     gr->Draw("ap");
1786     TLine* line = new TLine();
1787     double xmin = gr->GetXaxis()->GetXmin();
1788     double xmax = gr->GetXaxis()->GetXmax();
1789     line->SetLineColor(kBlack);
1790     line->DrawLine(xmin, 0, xmax, 0);
1791     TLine* line2 = new TLine();
1792     line2->SetLineStyle(kDashed);
1793     line2->DrawLine(xmin, -1.0, xmax, -1.0);
1794     line2->DrawLine(xmin, 1.0, xmax, 1.0);
1795     c1->SaveAs(outfilename.c_str());
1796 }
1797 
1798 void TSSAplotter::PlotPi0(std::string outfilename) {
1799     // Phi-dependent plots and fitting
1800     PlotAsymsBinned("#pi^{0} Blue Beam Forward-Going", bhs_pi0_blue_up_forward_phi_pT, bhs_pi0_blue_down_forward_phi_pT, outfilename);
1801     /* PlotAsymsBinned("#pi^{0} Blue Beam Backward-Going", bhs_pi0_blue_up_backward_phi_pT, bhs_pi0_blue_down_backward_phi_pT, outfilename); */
1802     PlotAsymsBinned("#pi^{0} Yellow Beam Forward-Going", bhs_pi0_yellow_up_forward_phi_pT, bhs_pi0_yellow_down_forward_phi_pT, outfilename);
1803     /* PlotAsymsBinned("#pi^{0} Yellow Beam Backward-Going", bhs_pi0_yellow_up_backward_phi_pT, bhs_pi0_yellow_down_backward_phi_pT, outfilename); */
1804 
1805     // Blue beam
1806     TGraphErrors* gr_rellumi_raw_forward_blue = AmplitudeBinnedGraph("rellumi", "#pi^{0} Blue Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0_blue_up_forward_phi_pT, bhs_pi0_blue_down_forward_phi_pT);
1807     TGraphErrors* gr_sqrt_raw_forward_blue = AmplitudeBinnedGraph("sqrt", "#pi^{0} Blue Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0_blue_up_forward_phi_pT, bhs_pi0_blue_down_forward_phi_pT);
1808     TGraphErrors* gr_rellumi_raw_backward_blue = AmplitudeBinnedGraph("rellumi", "#pi^{0} Blue Beam Backward-Going Asymmetry", "pT_pi0", bhs_pi0_blue_up_backward_phi_pT, bhs_pi0_blue_down_backward_phi_pT);
1809     TGraphErrors* gr_sqrt_raw_backward_blue = AmplitudeBinnedGraph("sqrt", "#pi^{0} Blue Beam Backward-Going Asymmetry", "pT_pi0", bhs_pi0_blue_up_backward_phi_pT, bhs_pi0_blue_down_backward_phi_pT);
1810     TGraphErrors* gr_t_method_blue = tTestGraph("dep", gr_rellumi_raw_forward_blue, gr_sqrt_raw_forward_blue);
1811     CompareGraphs(gr_rellumi_raw_forward_blue, gr_sqrt_raw_forward_blue, "#pi^{0} Blue Beam Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1812     PlottTest(gr_t_method_blue, "#pi^{0} Blue Beam Forward-Going Raw Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1813     CompareGraphs(gr_rellumi_raw_backward_blue, gr_sqrt_raw_backward_blue, "#pi^{0} Blue Beam Backward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1814     
1815     // Yellow beam
1816     TGraphErrors* gr_rellumi_raw_forward_yellow = AmplitudeBinnedGraph("rellumi", "#pi^{0} Yellow Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0_yellow_up_forward_phi_pT, bhs_pi0_yellow_down_forward_phi_pT);
1817     TGraphErrors* gr_sqrt_raw_forward_yellow = AmplitudeBinnedGraph("sqrt", "#pi^{0} Yellow Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0_yellow_up_forward_phi_pT, bhs_pi0_yellow_down_forward_phi_pT);
1818     TGraphErrors* gr_rellumi_raw_backward_yellow = AmplitudeBinnedGraph("rellumi", "#pi^{0} Yellow Beam Backward-Going Asymmetry", "pT_pi0", bhs_pi0_yellow_up_backward_phi_pT, bhs_pi0_yellow_down_backward_phi_pT);
1819     TGraphErrors* gr_sqrt_raw_backward_yellow = AmplitudeBinnedGraph("sqrt", "#pi^{0} Yellow Beam Backward-Going Asymmetry", "pT_pi0", bhs_pi0_yellow_up_backward_phi_pT, bhs_pi0_yellow_down_backward_phi_pT);
1820     TGraphErrors* gr_t_method_yellow = tTestGraph("dep", gr_rellumi_raw_forward_yellow, gr_sqrt_raw_forward_yellow);
1821     CompareGraphs(gr_rellumi_raw_forward_yellow, gr_sqrt_raw_forward_yellow, "#pi^{0} Yellow Beam Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1822     PlottTest(gr_t_method_yellow, "#pi^{0} Yellow Beam Forward-Going Raw Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1823     CompareGraphs(gr_rellumi_raw_backward_yellow, gr_sqrt_raw_backward_yellow, "#pi^{0} Yellow Beam Backward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1824     CompareGraphs(gr_rellumi_raw_forward_blue, gr_rellumi_raw_forward_yellow, "#pi^{0} Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry (Rel. Lumi.)", "Blue", "Yellow", outfilename);
1825     CompareGraphs(gr_sqrt_raw_forward_blue, gr_sqrt_raw_forward_yellow, "#pi^{0} Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry (Geometric)", "Blue", "Yellow", outfilename);
1826     TGraphErrors* gr_t_beam_rellumi = tTestGraph("indep", gr_rellumi_raw_forward_blue, gr_rellumi_raw_forward_yellow);
1827     PlottTest(gr_t_beam_rellumi, "#pi^{0} Forward-Going Raw Rel. Lumi. Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1828     TGraphErrors* gr_t_beam_sqrt = tTestGraph("indep", gr_sqrt_raw_forward_blue, gr_sqrt_raw_forward_yellow);
1829     PlottTest(gr_t_beam_sqrt, "#pi^{0} Forward-Going Raw Geometric Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1830     
1831     // Background
1832     TGraphErrors* gr_rellumi_raw_forward_blue_bkgr = AmplitudeBinnedGraph("rellumi", "#pi^{0} Blue Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0bkgr_blue_up_forward_phi_pT, bhs_pi0bkgr_blue_down_forward_phi_pT);
1833     TGraphErrors* gr_sqrt_raw_forward_blue_bkgr = AmplitudeBinnedGraph("sqrt", "#pi^{0} Blue Beam Forward-Going Asymmetry", "pT_pi0", bhs_pi0bkgr_blue_up_forward_phi_pT, bhs_pi0bkgr_blue_down_forward_phi_pT);
1834     TGraphErrors* gr_t_rellumi_blue_bkgr = tTestZeroGraph(gr_rellumi_raw_forward_blue_bkgr);
1835     TGraphErrors* gr_t_sqrt_blue_bkgr = tTestZeroGraph(gr_sqrt_raw_forward_blue_bkgr);
1836     PlottTest(gr_t_rellumi_blue_bkgr, "#pi^{0} Background Rel. Lumi. Asymmetry, t-Test with Zero", "p_{T} (GeV)", outfilename);
1837     PlottTest(gr_t_sqrt_blue_bkgr, "#pi^{0} Background Geometric Asymmetry, t-Test with Zero", "p_{T} (GeV)", outfilename);
1838 
1839     // Corrected asymmetries
1840     TGraphErrors* gr_corr_blue_rellumi = CorrectedAmplitudesGraph("rellumi", "#pi^{0} Blue Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_blue_up_forward_phi_pT, bhs_pi0_blue_down_forward_phi_pT, bhs_pi0bkgr_blue_up_forward_phi_pT, bhs_pi0bkgr_blue_down_forward_phi_pT);
1841     TGraphErrors* gr_corr_blue_sqrt = CorrectedAmplitudesGraph("sqrt", "#pi^{0} Blue Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_blue_up_forward_phi_pT, bhs_pi0_blue_down_forward_phi_pT, bhs_pi0bkgr_blue_up_forward_phi_pT, bhs_pi0bkgr_blue_down_forward_phi_pT);
1842     TGraphErrors* gr_corr_yellow_rellumi = CorrectedAmplitudesGraph("rellumi", "#pi^{0} Yellow Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_yellow_up_forward_phi_pT, bhs_pi0_yellow_down_forward_phi_pT, bhs_pi0bkgr_yellow_up_forward_phi_pT, bhs_pi0bkgr_yellow_down_forward_phi_pT);
1843     TGraphErrors* gr_corr_yellow_sqrt = CorrectedAmplitudesGraph("sqrt", "#pi^{0} Yellow Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_yellow_up_forward_phi_pT, bhs_pi0_yellow_down_forward_phi_pT, bhs_pi0bkgr_yellow_up_forward_phi_pT, bhs_pi0bkgr_yellow_down_forward_phi_pT);
1844     TGraphErrors* gr_t_corr_method_blue = tTestGraph("dep", gr_corr_blue_rellumi, gr_corr_blue_sqrt);
1845     TGraphErrors* gr_t_corr_method_yellow = tTestGraph("dep", gr_corr_yellow_rellumi, gr_corr_yellow_sqrt);
1846     TGraphErrors* gr_t_corr_beam_rellumi = tTestGraph("indep", gr_corr_blue_rellumi, gr_corr_yellow_rellumi);
1847     TGraphErrors* gr_t_corr_beam_sqrt = tTestGraph("dep", gr_corr_blue_sqrt, gr_corr_yellow_sqrt);
1848     CompareGraphs(gr_corr_blue_rellumi, gr_corr_blue_sqrt, "#pi^{0} Blue Beam Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename);
1849     PlotGraph(gr_t_corr_method_blue, "#pi^{0} Blue Beam Forward-Going Corrected Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", "t Statistic", outfilename, false, true);
1850     PlottTest(gr_t_corr_method_blue, "#pi^{0} Blue Beam Forward-Going Corrected Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1851     CompareGraphs(gr_corr_yellow_rellumi, gr_corr_yellow_sqrt, "#pi^{0} Yellow Beam Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename);
1852     PlottTest(gr_t_corr_method_yellow, "#pi^{0} Yellow Beam Forward-Going Corrected Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1853     CompareGraphs(gr_corr_blue_rellumi, gr_corr_yellow_rellumi, "#pi^{0} Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N} (Rel. Lumi.)", "Blue", "Yellow", outfilename);
1854     PlottTest(gr_t_corr_beam_rellumi, "#pi^{0} Forward-Going Rel. Lumi. Corrected Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1855     CompareGraphs(gr_corr_blue_sqrt, gr_corr_yellow_sqrt, "#pi^{0} Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N} (Geometric)", "Blue", "Yellow", outfilename);
1856     PlottTest(gr_t_corr_beam_sqrt, "#pi^{0} Forward-Going Geometric Corrected Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1857 
1858     // Final result
1859     TGraphErrors* gr_final_rellumi = AveragedGraphs(gr_corr_blue_rellumi, gr_corr_yellow_rellumi);
1860     TGraphErrors* gr_final_sqrt = AveragedGraphs(gr_corr_blue_sqrt, gr_corr_yellow_sqrt);
1861     TGraphErrors* gr_t_final_method = tTestGraph("dep", gr_final_rellumi, gr_final_sqrt);
1862     PlotGraph(gr_final_rellumi, "#pi^{0} Forward-Going A_{N} (Rel. Lumi.)", "p_{T} (GeV)", "A_{N}", outfilename, true);
1863     PlotGraph(gr_final_sqrt, "#pi^{0} Forward-Going A_{N} (Geometric)", "p_{T} (GeV)", "A_{N}", outfilename, true);
1864     CompareGraphs(gr_final_rellumi, gr_final_sqrt, "#pi^{0} Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename, true);
1865     PlottTest(gr_t_final_method, "#pi^{0} Forward-Going A_{N}, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1866 
1867     // PHENIX comparison
1868     TFile* f_phenix = new TFile("HEPData-ins1833997-v1-root.root", "READ");
1869     TDirectoryFile* tab1 = (TDirectoryFile*)f_phenix->Get("Table 1");
1870     TGraphAsymmErrors* gr_phenix = (TGraphAsymmErrors*)tab1->Get("Graph1D_y1");
1871     CompareGraphs(gr_final_rellumi, gr_phenix, "#pi^{0} Forward-Going Rel. Lumi. Asymmetry, PHENIX vs sPHENIX", "p_{T} (GeV)", "A_{N}", "sPHENIX", "PHENIX", outfilename, true);
1872     CompareGraphs(gr_final_sqrt, gr_phenix, "#pi^{0} Forward-Going Geometric Asymmetry, PHENIX vs sPHENIX", "p_{T} (GeV)", "A_{N}", "sPHENIX", "PHENIX", outfilename, true);
1873     // Low vs high eta
1874     TGraphErrors* gr_corr_blue_rellumi_lowEta = CorrectedAmplitudesGraph("rellumi", "#pi^{0} Blue Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_lowEta_blue_up_forward_phi_pT, bhs_pi0_lowEta_blue_down_forward_phi_pT, bhs_pi0bkgr_lowEta_blue_up_forward_phi_pT, bhs_pi0bkgr_lowEta_blue_down_forward_phi_pT);
1875     TGraphErrors* gr_corr_blue_sqrt_lowEta = CorrectedAmplitudesGraph("sqrt", "#pi^{0} Blue Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_lowEta_blue_up_forward_phi_pT, bhs_pi0_lowEta_blue_down_forward_phi_pT, bhs_pi0bkgr_lowEta_blue_up_forward_phi_pT, bhs_pi0bkgr_lowEta_blue_down_forward_phi_pT);
1876     TGraphErrors* gr_corr_yellow_rellumi_lowEta = CorrectedAmplitudesGraph("rellumi", "#pi^{0} Yellow Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_lowEta_yellow_up_forward_phi_pT, bhs_pi0_lowEta_yellow_down_forward_phi_pT, bhs_pi0bkgr_lowEta_yellow_up_forward_phi_pT, bhs_pi0bkgr_lowEta_yellow_down_forward_phi_pT);
1877     TGraphErrors* gr_corr_yellow_sqrt_lowEta = CorrectedAmplitudesGraph("sqrt", "#pi^{0} Yellow Beam Forward-Doing Asymmetry", "pT_pi0", bhs_pi0_lowEta_yellow_up_forward_phi_pT, bhs_pi0_lowEta_yellow_down_forward_phi_pT, bhs_pi0bkgr_lowEta_yellow_up_forward_phi_pT, bhs_pi0bkgr_lowEta_yellow_down_forward_phi_pT);
1878     TGraphErrors* gr_final_rellumi_lowEta = AveragedGraphs(gr_corr_blue_rellumi_lowEta, gr_corr_yellow_rellumi_lowEta);
1879     TGraphErrors* gr_final_sqrt_lowEta = AveragedGraphs(gr_corr_blue_sqrt_lowEta, gr_corr_yellow_sqrt_lowEta);
1880     CompareGraphs(gr_final_rellumi_lowEta, gr_phenix, "#pi^{0} Forward-Going Rel. Lumi. Asymmetry, PHENIX vs sPHENIX", "p_{T} (GeV)", "A_{N}", "sPHENIX |#eta|<0.35", "PHENIX", outfilename, true);
1881     CompareGraphs(gr_final_sqrt_lowEta, gr_phenix, "#pi^{0} Forward-Going Geometric Asymmetry, PHENIX vs sPHENIX", "p_{T} (GeV)", "A_{N}", "sPHENIX |#eta|<0.35", "PHENIX", outfilename, true);
1882 
1883     f_phenix->Close();
1884     delete f_phenix;
1885 }
1886 
1887 void TSSAplotter::PlotEta(std::string outfilename) {
1888     // Phi-dependent plots and fitting
1889     PlotAsymsBinned("#eta Blue Beam Forward-Going", bhs_eta_blue_up_forward_phi_pT, bhs_eta_blue_down_forward_phi_pT, outfilename);
1890     /* PlotAsymsBinned("#eta Blue Beam Backward-Going", bhs_eta_blue_up_backward_phi_pT, bhs_eta_blue_down_backward_phi_pT, outfilename); */
1891     PlotAsymsBinned("#eta Yellow Beam Forward-Going", bhs_eta_yellow_up_forward_phi_pT, bhs_eta_yellow_down_forward_phi_pT, outfilename);
1892     /* PlotAsymsBinned("#eta Yellow Beam Backward-Going", bhs_eta_yellow_up_backward_phi_pT, bhs_eta_yellow_down_backward_phi_pT, outfilename); */
1893 
1894     // Blue beam
1895     TGraphErrors* gr_rellumi_raw_forward_blue = AmplitudeBinnedGraph("rellumi", "#eta Blue Beam Forward-Going Asymmetry", "pT_eta", bhs_eta_blue_up_forward_phi_pT, bhs_eta_blue_down_forward_phi_pT);
1896     TGraphErrors* gr_sqrt_raw_forward_blue = AmplitudeBinnedGraph("sqrt", "#eta Blue Beam Forward-Going Asymmetry", "pT_eta", bhs_eta_blue_up_forward_phi_pT, bhs_eta_blue_down_forward_phi_pT);
1897     TGraphErrors* gr_rellumi_raw_backward_blue = AmplitudeBinnedGraph("rellumi", "#eta Blue Beam Backward-Going Asymmetry", "pT_eta", bhs_eta_blue_up_backward_phi_pT, bhs_eta_blue_down_backward_phi_pT);
1898     TGraphErrors* gr_sqrt_raw_backward_blue = AmplitudeBinnedGraph("sqrt", "#eta Blue Beam Backward-Going Asymmetry", "pT_eta", bhs_eta_blue_up_backward_phi_pT, bhs_eta_blue_down_backward_phi_pT);
1899     TGraphErrors* gr_t_method_blue = tTestGraph("dep", gr_rellumi_raw_forward_blue, gr_sqrt_raw_forward_blue);
1900     CompareGraphs(gr_rellumi_raw_forward_blue, gr_sqrt_raw_forward_blue, "#eta Blue Beam Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1901     PlottTest(gr_t_method_blue, "#eta Blue Beam Forward-Going Raw Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1902     CompareGraphs(gr_rellumi_raw_backward_blue, gr_sqrt_raw_backward_blue, "#eta Blue Beam Backward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1903     
1904     // Yellow beam
1905     TGraphErrors* gr_rellumi_raw_forward_yellow = AmplitudeBinnedGraph("rellumi", "#eta Yellow Beam Forward-Going Asymmetry", "pT_eta", bhs_eta_yellow_up_forward_phi_pT, bhs_eta_yellow_down_forward_phi_pT);
1906     TGraphErrors* gr_sqrt_raw_forward_yellow = AmplitudeBinnedGraph("sqrt", "#eta Yellow Beam Forward-Going Asymmetry", "pT_eta", bhs_eta_yellow_up_forward_phi_pT, bhs_eta_yellow_down_forward_phi_pT);
1907     TGraphErrors* gr_rellumi_raw_backward_yellow = AmplitudeBinnedGraph("rellumi", "#eta Yellow Beam Backward-Going Asymmetry", "pT_eta", bhs_eta_yellow_up_backward_phi_pT, bhs_eta_yellow_down_backward_phi_pT);
1908     TGraphErrors* gr_sqrt_raw_backward_yellow = AmplitudeBinnedGraph("sqrt", "#eta Yellow Beam Backward-Going Asymmetry", "pT_eta", bhs_eta_yellow_up_backward_phi_pT, bhs_eta_yellow_down_backward_phi_pT);
1909     TGraphErrors* gr_t_method_yellow = tTestGraph("dep", gr_rellumi_raw_forward_yellow, gr_sqrt_raw_forward_yellow);
1910     CompareGraphs(gr_rellumi_raw_forward_yellow, gr_sqrt_raw_forward_yellow, "#eta Yellow Beam Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1911     PlottTest(gr_t_method_yellow, "#eta Yellow Beam Forward-Going Raw Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1912     CompareGraphs(gr_rellumi_raw_backward_yellow, gr_sqrt_raw_backward_yellow, "#eta Yellow Beam Backward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry", "Rel. Lumi.", "Geometric", outfilename);
1913     CompareGraphs(gr_rellumi_raw_forward_blue, gr_rellumi_raw_forward_yellow, "#eta Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry (Rel. Lumi.)", "Blue", "Yellow", outfilename);
1914     CompareGraphs(gr_sqrt_raw_forward_blue, gr_sqrt_raw_forward_yellow, "#eta Forward-Going Asymmetry", "p_{T} (GeV)", "Raw Asymmetry (Geometric)", "Blue", "Yellow", outfilename);
1915     TGraphErrors* gr_t_beam_rellumi = tTestGraph("indep", gr_rellumi_raw_forward_blue, gr_rellumi_raw_forward_yellow);
1916     PlottTest(gr_t_beam_rellumi, "#eta Forward-Going Raw Rel. Lumi. Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1917     TGraphErrors* gr_t_beam_sqrt = tTestGraph("indep", gr_sqrt_raw_forward_blue, gr_sqrt_raw_forward_yellow);
1918     PlottTest(gr_t_beam_sqrt, "#eta Forward-Going Raw Geometric Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1919     
1920     // Background
1921     TGraphErrors* gr_rellumi_raw_forward_blue_bkgr = AmplitudeBinnedGraph("rellumi", "#eta Blue Beam Forward-Going Asymmetry", "pT_eta", bhs_etabkgr_blue_up_forward_phi_pT, bhs_etabkgr_blue_down_forward_phi_pT);
1922     TGraphErrors* gr_sqrt_raw_forward_blue_bkgr = AmplitudeBinnedGraph("sqrt", "#eta Blue Beam Forward-Going Asymmetry", "pT_eta", bhs_etabkgr_blue_up_forward_phi_pT, bhs_etabkgr_blue_down_forward_phi_pT);
1923     TGraphErrors* gr_t_rellumi_blue_bkgr = tTestZeroGraph(gr_rellumi_raw_forward_blue_bkgr);
1924     TGraphErrors* gr_t_sqrt_blue_bkgr = tTestZeroGraph(gr_sqrt_raw_forward_blue_bkgr);
1925     PlottTest(gr_t_rellumi_blue_bkgr, "#eta Background Rel. Lumi. Asymmetry, t-Test with Zero", "p_{T} (GeV)", outfilename);
1926     PlottTest(gr_t_sqrt_blue_bkgr, "#eta Background Geometric Asymmetry, t-Test with Zero", "p_{T} (GeV)", outfilename);
1927 
1928     // Corrected asymmetries
1929     TGraphErrors* gr_corr_blue_rellumi = CorrectedAmplitudesGraph("rellumi", "#eta Blue Beam Forward-Doing Asymmetry", "pT_eta", bhs_eta_blue_up_forward_phi_pT, bhs_eta_blue_down_forward_phi_pT, bhs_etabkgr_blue_up_forward_phi_pT, bhs_etabkgr_blue_down_forward_phi_pT);
1930     TGraphErrors* gr_corr_blue_sqrt = CorrectedAmplitudesGraph("sqrt", "#eta Blue Beam Forward-Doing Asymmetry", "pT_eta", bhs_eta_blue_up_forward_phi_pT, bhs_eta_blue_down_forward_phi_pT, bhs_etabkgr_blue_up_forward_phi_pT, bhs_etabkgr_blue_down_forward_phi_pT);
1931     TGraphErrors* gr_corr_yellow_rellumi = CorrectedAmplitudesGraph("rellumi", "#eta Yellow Beam Forward-Doing Asymmetry", "pT_eta", bhs_eta_yellow_up_forward_phi_pT, bhs_eta_yellow_down_forward_phi_pT, bhs_etabkgr_yellow_up_forward_phi_pT, bhs_etabkgr_yellow_down_forward_phi_pT);
1932     TGraphErrors* gr_corr_yellow_sqrt = CorrectedAmplitudesGraph("sqrt", "#eta Yellow Beam Forward-Doing Asymmetry", "pT_eta", bhs_eta_yellow_up_forward_phi_pT, bhs_eta_yellow_down_forward_phi_pT, bhs_etabkgr_yellow_up_forward_phi_pT, bhs_etabkgr_yellow_down_forward_phi_pT);
1933     TGraphErrors* gr_t_corr_method_blue = tTestGraph("dep", gr_corr_blue_rellumi, gr_corr_blue_sqrt);
1934     TGraphErrors* gr_t_corr_method_yellow = tTestGraph("dep", gr_corr_yellow_rellumi, gr_corr_yellow_sqrt);
1935     TGraphErrors* gr_t_corr_beam_rellumi = tTestGraph("indep", gr_corr_blue_rellumi, gr_corr_yellow_rellumi);
1936     TGraphErrors* gr_t_corr_beam_sqrt = tTestGraph("dep", gr_corr_blue_sqrt, gr_corr_yellow_sqrt);
1937     CompareGraphs(gr_corr_blue_rellumi, gr_corr_blue_sqrt, "#eta Blue Beam Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename);
1938     PlottTest(gr_t_corr_method_blue, "#eta Blue Beam Forward-Going Corrected Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1939     CompareGraphs(gr_corr_yellow_rellumi, gr_corr_yellow_sqrt, "#eta Yellow Beam Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename);
1940     PlottTest(gr_t_corr_method_yellow, "#eta Yellow Beam Forward-Going Corrected Asymmetry, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1941     CompareGraphs(gr_corr_blue_rellumi, gr_corr_yellow_rellumi, "#eta Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N} (Rel. Lumi.)", "Blue", "Yellow", outfilename);
1942     PlottTest(gr_t_corr_beam_rellumi, "#eta Forward-Going Rel. Lumi. Corrected Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1943     CompareGraphs(gr_corr_blue_sqrt, gr_corr_yellow_sqrt, "#eta Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N} (Geometric)", "Blue", "Yellow", outfilename);
1944     PlottTest(gr_t_corr_beam_sqrt, "#eta Forward-Going Geometric Corrected Asymmetry, t-Test Blue vs Yellow Beam", "p_{T} (GeV)", outfilename);
1945 
1946     // Final result
1947     TGraphErrors* gr_final_rellumi = AveragedGraphs(gr_corr_blue_rellumi, gr_corr_yellow_rellumi);
1948     TGraphErrors* gr_final_sqrt = AveragedGraphs(gr_corr_blue_sqrt, gr_corr_yellow_sqrt);
1949     TGraphErrors* gr_t_final_method = tTestGraph("dep", gr_final_rellumi, gr_final_sqrt);
1950     PlotGraph(gr_final_rellumi, "#eta Forward-Going A_{N} (Rel. Lumi.)", "p_{T} (GeV)", "A_{N}", outfilename, true);
1951     PlotGraph(gr_final_sqrt, "#eta Forward-Going A_{N} (Geometric)", "p_{T} (GeV)", "A_{N}", outfilename, true);
1952     CompareGraphs(gr_final_rellumi, gr_final_sqrt, "#eta Forward-Going Asymmetry", "p_{T} (GeV)", "A_{N}", "Rel. Lumi.", "Geometric", outfilename, true);
1953     PlottTest(gr_t_final_method, "#eta Forward-Going A_{N}, t-Test Rel. Lumi. vs Geometric", "p_{T} (GeV)", outfilename);
1954 }
1955 
1956 int TSSAplotter::main(std::string TSSApdf, std::string NMpdf) {
1957     TSSApdfname = TSSApdf;
1958     if (!NMpdf.empty()) {
1959     do_NMhists = true;
1960     NMpdfname = NMpdf;
1961     }
1962 
1963     if (do_NMhists) {
1964     GetNMHists();
1965     PlotNMHists();
1966     }
1967     GetTSSAHists();
1968     PlotTSSAs();
1969     return 0;
1970 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.3.7 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!