Back to home page

sPhenix code displayed by LXR
 
 

    

File indexing completed on 2026-07-16 08:11:27 

0001 #include "../util/binning.h"
0002 
0003 std::string getDate()
0004 {
0005     std::time_t t = std::time(0);   // get time now
0006     std::tm* now = std::localtime(&t);
0007 
0008     std::stringstream date;
0009     date << (now->tm_mon + 1) << '/'
0010          <<  now->tm_mday << '/'
0011          << (now->tm_year + 1900);
0012     return date.str();
0013 }
0014 
0015 /*
0016 void CutEfficiency()
0017 {
0018   string sim_inFile_Xi = "/sphenix/user/rosstom/TrackingAnalysis/Physics_Val_TF/sim/pTComp/KFParticle_2Mcascade_nopTCut.root";
0019   TFile* sim_file_Xi = new TFile(sim_inFile_Xi.c_str());
0020   TTree* sim_tree_Xi = (TTree*)sim_file_Xi->Get("DecayTree");
0021 
0022   const int nBins_pT = 5;
0023   double pTbins[nBins_pT + 1] = {0.75, 0.909375, 1.05792, 1.23212, 1.46559, 2.0};
0024 
0025   const int nBins_phi = 10;
0026   double phibins[nBins_phi + 1] = {-M_PI,-2.51177,-2.0196,-1.55623,-1.0106,-0.530045,0.0583016,1.05245,2.12079,2.59601,M_PI};
0027 
0028   const int nBins_eta = 10;
0029   double etabins[nBins_eta + 1] = {-1.1,-0.916345,-0.755196,-0.585922,-0.390321,-0.159045,0.122585,0.394623,0.637892,0.849511,1.1};
0030 
0031   std::string truthCut = "abs(Lambda0_track_1_true_ID) == 211 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3312) > 0 && abs(Lambda0_track_2_true_ID) == 2212 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3312) > 0 && abs(track_3_true_ID) == 211 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3312) > 0 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3122) == 0 && abs(Ximinus_rapidity) <= 1.0 && Lambda0_track_1_MVTX_nHits > 0 && Lambda0_track_2_MVTX_nHits > 0 && track_3_MVTX_nHits > 0";
0032   //std::string truthDataCut = "abs(Lambda0_track_1_true_ID) == 211 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3312) > 0 && abs(Lambda0_track_2_true_ID) == 2212 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3312) > 0 && abs(track_3_true_ID) == 211 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3312) > 0 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3122) == 0 && Lambda0_track_1_MVTX_nStates > 0 && Lambda0_track_2_MVTX_nStates > 0 && track_3_MVTX_nStates > 0 && Lambda0_track_1_INTT_nStates > 0 && Lambda0_track_2_INTT_nStates > 0 && track_3_INTT_nStates > 0 && Lambda0_track_1_TPC_nStates > 19 && Lambda0_track_2_TPC_nStates > 19 && track_3_TPC_nStates > 19 && abs(Lambda0_mass - 1.1157) <= 0.01 && Ximinus_decayLength_xy >= 0.15 && Lambda0_decayLength_xy >= 0.01 && (Lambda0_track_1_chi2/Lambda0_track_1_nDoF) <= 400 && (Lambda0_track_2_chi2/Lambda0_track_2_nDoF) <= 400 && (track_3_chi2/track_3_nDoF) <= 400 && track_1_track_2_DCA <= 0.5 && Ximinus_chi2/Ximinus_nDoF <= 50 && track_1_track_3_DCA <= 0.5 && track_2_track_3_DCA <= 0.5 && Lambda0_chi2/Lambda0_nDoF <= 50";
0033   std::string truthDataCut = "abs(Lambda0_track_1_true_ID) == 211 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3312) > 0 && abs(Lambda0_track_2_true_ID) == 2212 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3312) > 0 && abs(track_3_true_ID) == 211 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3312) > 0 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3122) == 0 && Lambda0_track_1_MVTX_nHits > 0 && Lambda0_track_2_MVTX_nHits > 0 && track_3_MVTX_nHits > 0 && Lambda0_track_1_INTT_nHits > 0 && Lambda0_track_2_INTT_nHits > 0 && track_3_INTT_nHits > 0 && Lambda0_track_1_TPC_nHits > 19 && Lambda0_track_2_TPC_nHits > 19 && track_3_TPC_nHits > 19 && abs(Lambda0_mass - 1.1157) <= 0.01 && Ximinus_decayLength_xy >= 0.15 && Lambda0_decayLength_xy >= 0.01 && (Lambda0_track_1_chi2/Lambda0_track_1_nDoF) <= 400 && (Lambda0_track_2_chi2/Lambda0_track_2_nDoF) <= 400 && (track_3_chi2/track_3_nDoF) <= 400 && track_1_track_2_DCA <= 0.5 && track_1_track_3_DCA <= 0.5 && track_2_track_3_DCA <= 0.5 && Ximinus_chi2/Ximinus_nDoF <= 50 && Lambda0_chi2/Lambda0_nDoF <= 50 && abs(Ximinus_rapidity) <= 1.0";
0034 
0035   TH1F* h_pT_truth = new TH1F("h_pT_truth",";#Xi^{-} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins);
0036   TH1F* h_pT_truthCuts = new TH1F("h_pT_truthCuts",";#Xi^{-} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins);
0037   
0038   TH1F* h_eta_truth = new TH1F("h_eta_truth",";#Xi^{-} #eta;Cut Efficiency",nBins_eta,etabins);
0039   TH1F* h_eta_truthCuts = new TH1F("h_eta_truthCuts",";#Xi^{-} #eta;Cut Efficiency",nBins_eta,etabins);
0040   
0041   TH1F* h_phi_truth = new TH1F("h_phi_truth",";#Xi^{-} #phi;Cut Efficiency",nBins_phi,phibins);
0042   TH1F* h_phi_truthCuts = new TH1F("h_phi_truthCuts",";#Xi^{-} #phi;Cut Efficiency",nBins_phi,phibins);
0043   
0044   TH1F* h_y_truth = new TH1F("h_y_truth",";#Xi^{-} y;Cut Efficiency",10,-1.0,1.0);
0045   TH1F* h_y_truthCuts = new TH1F("h_y_truthCuts",";#Xi^{-} y;Cut Efficiency",10,-1.0,1.0);
0046    
0047   h_pT_truth->Sumw2();
0048   h_pT_truthCuts->Sumw2();
0049   h_eta_truth->Sumw2();
0050   h_eta_truthCuts->Sumw2();
0051   h_phi_truth->Sumw2();
0052   h_phi_truthCuts->Sumw2();
0053   h_y_truth->Sumw2();
0054   h_y_truthCuts->Sumw2();
0055 
0056   sim_tree_Xi->Draw("Ximinus_pT>>+h_pT_truth",truthCut.c_str());
0057   sim_tree_Xi->Draw("Ximinus_pT>>+h_pT_truthCuts",truthDataCut.c_str());
0058   sim_tree_Xi->Draw("Ximinus_pseudorapidity>>+h_eta_truth",truthCut.c_str());
0059   sim_tree_Xi->Draw("Ximinus_pseudorapidity>>+h_eta_truthCuts",truthDataCut.c_str());
0060   sim_tree_Xi->Draw("Ximinus_phi>>+h_phi_truth",truthCut.c_str());
0061   sim_tree_Xi->Draw("Ximinus_phi>>+h_phi_truthCuts",truthDataCut.c_str());
0062   sim_tree_Xi->Draw("Ximinus_rapidity>>+h_y_truth",truthCut.c_str());
0063   sim_tree_Xi->Draw("Ximinus_rapidity>>+h_y_truthCuts",truthDataCut.c_str());
0064 
0065   TEfficiency* h_pT_eff = new TEfficiency(*h_pT_truthCuts, *h_pT_truth);
0066   TEfficiency* h_eta_eff = new TEfficiency(*h_eta_truthCuts, *h_eta_truth);
0067   TEfficiency* h_phi_eff = new TEfficiency(*h_phi_truthCuts, *h_phi_truth);
0068   TEfficiency* h_y_eff = new TEfficiency(*h_y_truthCuts, *h_y_truth);
0069   h_pT_eff->SetStatisticOption(TEfficiency::kFCP);
0070   h_pT_eff->SetMarkerStyle(20);
0071   h_pT_eff->SetMarkerSize(1.2);
0072   h_eta_eff->SetStatisticOption(TEfficiency::kFCP);
0073   h_eta_eff->SetMarkerStyle(20);
0074   h_eta_eff->SetMarkerSize(1.2);
0075   h_phi_eff->SetStatisticOption(TEfficiency::kFCP);
0076   h_phi_eff->SetMarkerStyle(20);
0077   h_phi_eff->SetMarkerSize(1.2);
0078   h_y_eff->SetStatisticOption(TEfficiency::kFCP);
0079   h_y_eff->SetMarkerStyle(20);
0080   h_y_eff->SetMarkerSize(1.2);
0081   //h_pT_eff->GetXaxis()->SetTitle("#Xi^{-} p_{T} [GeV/c]");
0082   //h_pT_eff->GetYaxis()->SetTitle("Cut Efficiency");  
0083 
0084   TCanvas *c1  = new TCanvas("myCanvas", "myCanvas",1600,1600);
0085   c1->Divide(2,2);
0086   //pT 
0087   c1->cd(1);
0088   h_pT_eff->Draw("AP");
0089   c1->Update();
0090 
0091   TGraph* g = h_pT_eff->GetPaintedGraph();
0092   if (g)
0093   {
0094     g->GetXaxis()->SetTitle("#Xi^{-} p_{T} [GeV/c]");
0095     //g->GetXaxis()->SetTitle("#Xi^{-} #eta");
0096     //g->GetXaxis()->SetTitle("#Xi^{-} #phi");
0097     //g->GetXaxis()->SetTitle("#Xi^{-} y");
0098     g->GetYaxis()->SetTitle("Cut Efficiency");
0099     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0100     g->GetXaxis()->SetLimits(pTbins[0], pTbins[nBins_pT]);
0101     //g->GetXaxis()->SetLimits(etabins[0], etabins[nBins_eta]);
0102     //g->GetXaxis()->SetLimits(phibins[0], phibins[nBins_phi]);
0103     //g->GetXaxis()->SetLimits(-1.5,1.5);
0104   }
0105   h_pT_eff->Draw("AP"); 
0106  
0107   TPaveText *pt;
0108   pt = new TPaveText(0.15,0.75,0.65,1., "NDC");
0109   pt->SetFillColor(0);
0110   pt->SetFillStyle(0);
0111   pt->SetTextFont(42);
0112   string label = "#it{#bf{sPHENIX}} Simulation";
0113   TText *pt_LaTex = pt->AddText(label.c_str());
0114   pt->SetBorderSize(0);
0115   pt->Draw();
0116   gPad->Modified();
0117 
0118   TPaveText *ptDate;
0119   ptDate = new TPaveText(0.67,0.79,0.85,0.95, "NDC");
0120   ptDate->SetFillColor(0);
0121   ptDate->SetFillStyle(0);
0122   ptDate->SetTextFont(42);
0123   std::string compilation_date = getDate();
0124   TText *pt_LaTexDate = ptDate->AddText(compilation_date.c_str());
0125   ptDate->SetBorderSize(0);
0126   ptDate->Draw();
0127   gPad->Modified();
0128  
0129   //eta
0130   c1->cd(2);
0131   h_eta_eff->Draw("AP");
0132   c1->Update();
0133 
0134   g = h_eta_eff->GetPaintedGraph();
0135   if (g)
0136   {
0137     g->GetXaxis()->SetTitle("#Xi^{-} #eta");
0138     g->GetYaxis()->SetTitle("Cut Efficiency");
0139     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0140     g->GetXaxis()->SetLimits(etabins[0], etabins[nBins_eta]);
0141   }
0142   h_eta_eff->Draw("AP"); 
0143  
0144   pt->Draw();
0145   gPad->Modified();
0146 
0147   ptDate->Draw();
0148   gPad->Modified();
0149   
0150   //phi
0151   c1->cd(3);
0152   h_phi_eff->Draw("AP");
0153   c1->Update();
0154 
0155   g = h_phi_eff->GetPaintedGraph();
0156   if (g)
0157   {
0158     g->GetXaxis()->SetTitle("#Xi^{-} #phi [rad]");
0159     g->GetYaxis()->SetTitle("Cut Efficiency");
0160     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0161     g->GetXaxis()->SetLimits(phibins[0], phibins[nBins_phi]);
0162   }
0163   h_phi_eff->Draw("AP"); 
0164  
0165   pt->Draw();
0166   gPad->Modified();
0167 
0168   ptDate->Draw();
0169   gPad->Modified();
0170   
0171   //rapidity
0172   c1->cd(4);
0173   h_y_eff->Draw("AP");
0174   c1->Update();
0175 
0176   g = h_y_eff->GetPaintedGraph();
0177   if (g)
0178   {
0179     g->GetXaxis()->SetTitle("#Xi^{-} y");
0180     g->GetYaxis()->SetTitle("Cut Efficiency");
0181     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0182     g->GetXaxis()->SetLimits(-1.0,1.0);
0183   }
0184   h_y_eff->Draw("AP"); 
0185  
0186   pt->Draw();
0187   gPad->Modified();
0188 
0189   ptDate->Draw();
0190   gPad->Modified();
0191   
0192   string extensions[] = {".png"};
0193   for (auto extension : extensions)
0194   {
0195     string output = "Xi_all_cutEff" + extension;
0196     //string output = "Xi_pT_cutEff" + extension;
0197     //string output = "Xi_eta_cutEff" + extension;
0198     //string output = "Xi_phi_cutEff" + extension;
0199     //string output = "Xi_y_cutEff" + extension;
0200     c1->SaveAs(output.c_str());
0201   }
0202 } 
0203 */
0204 /*
0205 void CutEfficiency()
0206 {
0207   string sim_inFile_KS0 = "/sphenix/user/cdean/software/analysis/LightFlavorRatios/geometric_acceptance/simulation/outputKFParticle_Kshort_reco.root";
0208   TFile* sim_file_KS0 = new TFile(sim_inFile_KS0.c_str());
0209   TTree* sim_tree_KS0 = (TTree*)sim_file_KS0->Get("DecayTree");
0210 
0211   const int nBins_pT = 5;
0212   double pTbins[nBins_pT + 1] = {0.75, 0.909375, 1.05792, 1.23212, 1.46559, 2.0};
0213 
0214   const int nBins_phi = 10;
0215   double phibins[nBins_phi + 1] = {-M_PI,-2.51177,-2.0196,-1.55623,-1.0106,-0.530045,0.0583016,1.05245,2.12079,2.59601,M_PI};
0216 
0217   const int nBins_eta = 10;
0218   double etabins[nBins_eta + 1] = {-1.1,-0.916345,-0.755196,-0.585922,-0.390321,-0.159045,0.122585,0.394623,0.637892,0.849511,1.1};
0219 
0220   std::string truthCut_KS0 = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0";
0221   std::string truthDataCut_KS0 = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && abs(track_1_IP_xy) >= 0.05 && abs(track_2_IP_xy) >= 0.05 && track_1_track_2_DCA <= 0.5 && track_1_track_2_DCA_xy <= 1.0 && K_S0_DIRA >= 0.99 && (K_S0_chi2/K_S0_nDoF) <= 20 && abs(K_S0_mass - 0.5) <= 0.1";
0222 
0223   TH1F* h_pT_truth_KS0 = new TH1F("h_pT_truth_KS0",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins);
0224   TH1F* h_pT_truthCuts_KS0 = new TH1F("h_pT_truthCuts_KS0",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins);
0225   
0226   TH1F* h_eta_truth_KS0 = new TH1F("h_eta_truth_KS0",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins);
0227   TH1F* h_eta_truthCuts_KS0 = new TH1F("h_eta_truthCuts_KS0",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins);
0228   
0229   TH1F* h_phi_truth_KS0 = new TH1F("h_phi_truth_KS0",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins);
0230   TH1F* h_phi_truthCuts_KS0 = new TH1F("h_phi_truthCuts_KS0",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins);
0231   
0232   TH1F* h_y_truth_KS0 = new TH1F("h_y_truth_KS0",";K_{S}^{0} y;Cut Efficiency",10,-1.0,1.0);
0233   TH1F* h_y_truthCuts_KS0 = new TH1F("h_y_truthCuts_KS0",";K_{S}^{0} y;Cut Efficiency",10,-1.0,1.0);
0234    
0235   h_pT_truth_KS0->Sumw2();
0236   h_pT_truthCuts_KS0->Sumw2();
0237   h_eta_truth_KS0->Sumw2();
0238   h_eta_truthCuts_KS0->Sumw2();
0239   h_phi_truth_KS0->Sumw2();
0240   h_phi_truthCuts_KS0->Sumw2();
0241   h_y_truth_KS0->Sumw2();
0242   h_y_truthCuts_KS0->Sumw2();
0243 
0244   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truth_KS0",truthCut_KS0.c_str());
0245   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truthCuts_KS0",truthDataCut_KS0.c_str());
0246   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truth_KS0",truthCut_KS0.c_str());
0247   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truthCuts_KS0",truthDataCut_KS0.c_str());
0248   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truth_KS0",truthCut_KS0.c_str());
0249   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truthCuts_KS0",truthDataCut_KS0.c_str());
0250   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truth_KS0",truthCut_KS0.c_str());
0251   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truthCuts_KS0",truthDataCut_KS0.c_str());
0252 
0253   TEfficiency* h_pT_eff_KS0 = new TEfficiency(*h_pT_truthCuts_KS0, *h_pT_truth_KS0);
0254   TEfficiency* h_eta_eff_KS0 = new TEfficiency(*h_eta_truthCuts_KS0, *h_eta_truth_KS0);
0255   TEfficiency* h_phi_eff_KS0 = new TEfficiency(*h_phi_truthCuts_KS0, *h_phi_truth_KS0);
0256   TEfficiency* h_y_eff_KS0 = new TEfficiency(*h_y_truthCuts_KS0, *h_y_truth_KS0);
0257   h_pT_eff_KS0->SetStatisticOption(TEfficiency::kFCP);
0258   h_pT_eff_KS0->SetMarkerStyle(20);
0259   h_pT_eff_KS0->SetMarkerSize(1.2);
0260   h_eta_eff_KS0->SetStatisticOption(TEfficiency::kFCP);
0261   h_eta_eff_KS0->SetMarkerStyle(20);
0262   h_eta_eff_KS0->SetMarkerSize(1.2);
0263   h_phi_eff_KS0->SetStatisticOption(TEfficiency::kFCP);
0264   h_phi_eff_KS0->SetMarkerStyle(20);
0265   h_phi_eff_KS0->SetMarkerSize(1.2);
0266   h_y_eff_KS0->SetStatisticOption(TEfficiency::kFCP);
0267   h_y_eff_KS0->SetMarkerStyle(20);
0268   h_y_eff_KS0->SetMarkerSize(1.2);
0269 
0270   TCanvas *c1  = new TCanvas("myCanvas", "myCanvas",1600,1600);
0271   c1->Divide(2,2);
0272   //pT 
0273   c1->cd(1);
0274   h_pT_eff_KS0->Draw("AP");
0275   c1->Update();
0276 
0277   TGraph* g = h_pT_eff_KS0->GetPaintedGraph();
0278   if (g)
0279   {
0280     g->GetXaxis()->SetTitle("K_{S}^{0} p_{T} [GeV/c]");
0281     g->GetYaxis()->SetTitle("Cut Efficiency");
0282     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0283     g->GetXaxis()->SetLimits(pTbins[0], pTbins[nBins_pT]);
0284   }
0285   h_pT_eff_KS0->Draw("AP"); 
0286  
0287   TPaveText *pt;
0288   pt = new TPaveText(0.15,0.75,0.65,1., "NDC");
0289   pt->SetFillColor(0);
0290   pt->SetFillStyle(0);
0291   pt->SetTextFont(42);
0292   string label = "#it{#bf{sPHENIX}} Simulation";
0293   TText *pt_LaTex = pt->AddText(label.c_str());
0294   pt->SetBorderSize(0);
0295   pt->Draw();
0296   gPad->Modified();
0297 
0298   TPaveText *ptDate;
0299   ptDate = new TPaveText(0.67,0.79,0.85,0.95, "NDC");
0300   ptDate->SetFillColor(0);
0301   ptDate->SetFillStyle(0);
0302   ptDate->SetTextFont(42);
0303   std::string compilation_date = getDate();
0304   TText *pt_LaTexDate = ptDate->AddText(compilation_date.c_str());
0305   ptDate->SetBorderSize(0);
0306   ptDate->Draw();
0307   gPad->Modified();
0308  
0309   //eta
0310   c1->cd(2);
0311   h_eta_eff_KS0->Draw("AP");
0312   c1->Update();
0313 
0314   g = h_eta_eff_KS0->GetPaintedGraph();
0315   if (g)
0316   {
0317     g->GetXaxis()->SetTitle("K_{S}^{0} #eta");
0318     g->GetYaxis()->SetTitle("Cut Efficiency");
0319     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0320     g->GetXaxis()->SetLimits(etabins[0], etabins[nBins_eta]);
0321   }
0322   h_eta_eff_KS0->Draw("AP"); 
0323  
0324   pt->Draw();
0325   gPad->Modified();
0326 
0327   ptDate->Draw();
0328   gPad->Modified();
0329   
0330   //phi
0331   c1->cd(3);
0332   h_phi_eff_KS0->Draw("AP");
0333   c1->Update();
0334 
0335   g = h_phi_eff_KS0->GetPaintedGraph();
0336   if (g)
0337   {
0338     g->GetXaxis()->SetTitle("K_{S}^{0} #phi [rad]");
0339     g->GetYaxis()->SetTitle("Cut Efficiency");
0340     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0341     g->GetXaxis()->SetLimits(phibins[0], phibins[nBins_phi]);
0342   }
0343   h_phi_eff_KS0->Draw("AP"); 
0344  
0345   pt->Draw();
0346   gPad->Modified();
0347 
0348   ptDate->Draw();
0349   gPad->Modified();
0350   
0351   //rapidity
0352   c1->cd(4);
0353   h_y_eff_KS0->Draw("AP");
0354   c1->Update();
0355 
0356   g = h_y_eff_KS0->GetPaintedGraph();
0357   if (g)
0358   {
0359     g->GetXaxis()->SetTitle("K_{S}^{0} y");
0360     g->GetYaxis()->SetTitle("Cut Efficiency");
0361     g->GetYaxis()->SetRangeUser(0.0, 1.05);
0362     g->GetXaxis()->SetLimits(-1.0,1.0);
0363   }
0364   h_y_eff_KS0->Draw("AP"); 
0365  
0366   pt->Draw();
0367   gPad->Modified();
0368 
0369   ptDate->Draw();
0370   gPad->Modified();
0371   
0372   string extensions[] = {".png"};
0373   for (auto extension : extensions)
0374   {
0375     string output = "KS0_all_cutEff" + extension;
0376     c1->SaveAs(output.c_str());
0377   }
0378 } 
0379 */
0380 void CutEfficiency()
0381 {
0382   //string sim_inFile_Xi = "/sphenix/tg/tg01/hf/mjpeters/LightFlavorResults/outputKFParticle_Lambda_reco.root";
0383   string sim_inFile_Xi = "/gpfs/mnt/gpfs02/sphenix/user/cdean/software/analysis/LightFlavorRatios/geometric_acceptance/simulation/outputKFParticle_Lambda2ppi_reco_Usman_patch.root";
0384   TFile* sim_file_Xi = new TFile(sim_inFile_Xi.c_str());
0385   TTree* sim_tree_Xi = (TTree*)sim_file_Xi->Get("DecayTree");
0386 
0387   /*
0388   const int nBins_pT = 5;
0389   double pTbins[nBins_pT + 1] = {0.75, 0.909375, 1.05792, 1.23212, 1.46559, 2.0};
0390 
0391   const int nBins_phi = 10;
0392   double phibins[nBins_phi + 1] = {-M_PI,-2.51177,-2.0196,-1.55623,-1.0106,-0.530045,0.0583016,1.05245,2.12079,2.59601,M_PI};
0393 
0394   const int nBins_eta = 10;
0395   double etabins[nBins_eta + 1] = {-1.1,-0.916345,-0.755196,-0.585922,-0.390321,-0.159045,0.122585,0.394623,0.637892,0.849511,1.1};
0396   */
0397 
0398   HistogramInfo pTbininfo = BinInfo::final_pt_bins;
0399   HistogramInfo phibininfo = BinInfo::final_phi_bins;
0400   HistogramInfo etabininfo = BinInfo::final_eta_bins;
0401   HistogramInfo ybininfo = BinInfo::final_rapidity_bins;
0402 
0403   const int nBins_pT = pTbininfo.bins.size()-1;
0404   std::vector<float> pTbins = pTbininfo.bins;
0405 
0406   const int nBins_phi = phibininfo.bins.size()-1;
0407   std::vector<float> phibins = phibininfo.bins;
0408 
0409   const int nBins_eta = etabininfo.bins.size()-1;
0410   std::vector<float> etabins = etabininfo.bins;
0411 
0412   const int nBins_rap = ybininfo.bins.size()-1;
0413   std::vector<float> rapbins = ybininfo.bins;
0414 
0415   double Lambda_PV_sigma_xy = .003553;
0416   double Lambda_track_1_sigma_xy = .002853;
0417   double Lambda_track_2_sigma_xy = .002841;
0418   double Lambda_track_1_sigma = .003507;
0419   double Lambda_track_2_sigma = .003192;
0420   double Lambda_PV_sigma = pow(9.756e-8,1./3.);
0421   double Lambda_SV_sigma = pow(2.663e-6,1./3.);
0422 
0423   double KS0_PV_sigma_xy = .003477; // sqrt(C_00 + C_11 - 2*C_01) from primary_vertex_Covariance
0424   double KS0_track_1_sigma_xy = .002584; // sqrt(C_00 + C_11 - 2*C_01) from track_1_Covariance
0425   double KS0_track_2_sigma_xy = .002585; // above but for track_2_Covariance
0426   double KS0_track_1_sigma = .002985; // track_1_IPErr
0427   double KS0_track_2_sigma = .002959; // track_2_IPErr
0428   double KS0_PV_sigma = pow(9.087e-8,1./3.); // (mother vertex volume)^1/3
0429   double KS0_SV_sigma = pow(7.531e-6,1./3.); // (primary vertex volume)^1/3
0430 
0431   std::string truthCut = "((abs(track_1_true_ID) == 211 && abs(track_2_true_ID) == 2212) || (abs(track_1_true_ID) == 2212 && abs(track_2_true_ID) == 211)) && Sum$(abs(track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(track_2_true_track_history_PDG_ID) == 3122) > 0 && abs(Lambda0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10";
0432 
0433   std::string truthDataCut = "((abs(track_1_true_ID) == 211 && abs(track_2_true_ID) == 2212) || (abs(track_1_true_ID) == 2212 && abs(track_2_true_ID) == 211)) && Sum$(abs(track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(track_2_true_track_history_PDG_ID) == 3122) > 0 && abs(Lambda0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && abs(track_1_IP_xy) >= 0.05 && abs(track_2_IP_xy) >= 0.05 && Lambda0_DIRA >= 0.99 && abs(Lambda0_mass - 1.115) <= 0.1 && Lambda0_decayLength>=0.05 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0434 
0435   std::string truthDataCut_min = "((abs(track_1_true_ID) == 211 && abs(track_2_true_ID) == 2212) || (abs(track_1_true_ID) == 2212 && abs(track_2_true_ID) == 211)) && Sum$(abs(track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(track_2_true_track_history_PDG_ID) == 3122) > 0 && abs(Lambda0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && (abs(track_1_IP_xy) - " + std::to_string(Lambda_PV_sigma_xy) + " - " + std::to_string(Lambda_track_1_sigma_xy) + ") >= 0.05 && (abs(track_2_IP_xy) - " + std::to_string(Lambda_PV_sigma_xy) + " - " + std::to_string(Lambda_track_2_sigma_xy) + ") >= 0.05 && (track_1_track_2_DCA + " + std::to_string(Lambda_track_1_sigma) + " + " + std::to_string(Lambda_track_2_sigma) + ") <= 0.5 && (track_1_track_2_DCA_xy + " + std::to_string(Lambda_track_1_sigma_xy) + " + " + std::to_string(Lambda_track_2_sigma_xy) + ") <= 1.0 && cos(acos(Lambda0_DIRA) + atan(" + std::to_string(Lambda_PV_sigma) + "/Lambda0_decayLength) + atan(" + std::to_string(Lambda_SV_sigma) + "/Lambda0_decayLength)) >= 0.99 && abs(Lambda0_mass - 1.115) <= 0.1 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0436   
0437   std::string truthDataCut_max = " ((abs(track_1_true_ID) == 211 && abs(track_2_true_ID) == 2212) || (abs(track_1_true_ID) == 2212 && abs(track_2_true_ID) == 211)) && Sum$(abs(track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(track_2_true_track_history_PDG_ID) == 3122) > 0 && abs(Lambda0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && (abs(track_1_IP_xy) + " + std::to_string(Lambda_PV_sigma_xy) + " + " + std::to_string(Lambda_track_1_sigma_xy) + ") >= 0.05 && (abs(track_2_IP_xy) + " + std::to_string(Lambda_PV_sigma_xy) + " + " + std::to_string(Lambda_track_2_sigma_xy) + ") >= 0.05 && (track_1_track_2_DCA - " + std::to_string(Lambda_track_1_sigma) + " - " + std::to_string(Lambda_track_2_sigma) + ") <= 0.5 && (track_1_track_2_DCA_xy - " + std::to_string(Lambda_track_1_sigma_xy) + " - " + std::to_string(Lambda_track_2_sigma_xy) + ") <= 1.0 && cos(acos(Lambda0_DIRA) - atan(" + std::to_string(Lambda_PV_sigma) + "/Lambda0_decayLength) - atan(" + std::to_string(Lambda_SV_sigma) + "/Lambda0_decayLength)) >= 0.99 && abs(Lambda0_mass - 1.115) <= 0.1 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0438 
0439   //std::string truthCut = "abs(Lambda0_track_1_true_ID) == 211 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3312) > 0 && abs(Lambda0_track_2_true_ID) == 2212 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3312) > 0 && abs(track_3_true_ID) == 211 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3312) > 0 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3122) == 0 && abs(Ximinus_rapidity) <= 1.0 && Lambda0_track_1_MVTX_nHits > 0 && Lambda0_track_2_MVTX_nHits > 0 && track_3_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && Lambda0_track_1_charge == track_3_charge";
0440   //std::string truthDataCut = "abs(Lambda0_track_1_true_ID) == 211 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_1_true_track_history_PDG_ID) == 3312) > 0 && abs(Lambda0_track_2_true_ID) == 2212 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3122) > 0 && Sum$(abs(Lambda0_track_2_true_track_history_PDG_ID) == 3312) > 0 && abs(track_3_true_ID) == 211 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3312) > 0 && Sum$(abs(track_3_true_track_history_PDG_ID) == 3122) == 0 && Lambda0_track_1_MVTX_nHits > 0 && Lambda0_track_2_MVTX_nHits > 0 && track_3_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && Lambda0_track_1_charge == track_3_charge && Lambda0_track_1_INTT_nHits > 0 && Lambda0_track_2_INTT_nHits > 0 && track_3_INTT_nHits > 0 && Lambda0_track_1_TPC_nHits > 19 && Lambda0_track_2_TPC_nHits > 19 && track_3_TPC_nHits > 19 && abs(Lambda0_mass - 1.1157) <= 0.01 && Ximinus_decayLength_xy >= 0.15 && Lambda0_decayLength_xy >= 0.01 && (Lambda0_track_1_chi2/Lambda0_track_1_nDoF) <= 400 && (Lambda0_track_2_chi2/Lambda0_track_2_nDoF) <= 400 && (track_3_chi2/track_3_nDoF) <= 400 && track_1_track_2_DCA <= 0.5 && track_1_track_3_DCA <= 0.5 && track_2_track_3_DCA <= 0.5 && Ximinus_chi2/Ximinus_nDoF <= 50 && Lambda0_chi2/Lambda0_nDoF <= 50 && abs(Ximinus_rapidity) <= 1.0 && Lambda0_track_1_pT >= 0.2 && Lambda0_track_2_pT >= 0.2 && track_3_pT >= 0.2";
0441 
0442   TH1F* h_pT_truth = new TH1F("h_pT_truth",";#Lambda p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0443   TH1F* h_pT_truthCuts = new TH1F("h_pT_truthCuts",";#Lambda p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0444   TH1F* h_pT_truthCuts_min = new TH1F("h_pT_truthCuts_min","#Lambda p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0445   TH1F* h_pT_truthCuts_max = new TH1F("h_pT_truthCuts_max","#Lambda p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0446   
0447   TH1F* h_eta_truth = new TH1F("h_eta_truth",";#Lambda #eta;Cut Efficiency",nBins_eta,etabins.data());
0448   TH1F* h_eta_truthCuts = new TH1F("h_eta_truthCuts",";#Lambda #eta;Cut Efficiency",nBins_eta,etabins.data());
0449   TH1F* h_eta_truthCuts_min = new TH1F("h_eta_truthCuts_min",";#Lambda #eta;Cut Efficiency",nBins_eta,etabins.data());
0450   TH1F* h_eta_truthCuts_max = new TH1F("h_eta_truthCuts_max",";#Lambda #eta;Cut Efficiency",nBins_eta,etabins.data());
0451   
0452   TH1F* h_phi_truth = new TH1F("h_phi_truth",";#Lambda #phi;Cut Efficiency",nBins_phi,phibins.data());
0453   TH1F* h_phi_truthCuts = new TH1F("h_phi_truthCuts",";#Lambda #phi;Cut Efficiency",nBins_phi,phibins.data());
0454   TH1F* h_phi_truthCuts_min = new TH1F("h_phi_truthCuts_min",";#Lambda #phi;Cut Efficiency",nBins_phi,phibins.data());
0455   TH1F* h_phi_truthCuts_max = new TH1F("h_phi_truthCuts_max",";#Lambda #phi;Cut Efficiency",nBins_phi,phibins.data());
0456   
0457   TH1F* h_y_truth = new TH1F("h_y_truth",";#Lambda y;Cut Efficiency",nBins_rap,rapbins.data());
0458   TH1F* h_y_truthCuts = new TH1F("h_y_truthCuts",";#Lambda y;Cut Efficiency",nBins_rap,rapbins.data());
0459   TH1F* h_y_truthCuts_min = new TH1F("h_y_truthCuts_min",";#Lambda y;Cut Efficiency",nBins_rap,rapbins.data());
0460   TH1F* h_y_truthCuts_max = new TH1F("h_y_truthCuts_max",";#Lambda y;Cut Efficiency",nBins_rap,rapbins.data());
0461    
0462   h_pT_truth->Sumw2();
0463   h_pT_truthCuts->Sumw2();
0464   h_pT_truthCuts_min->Sumw2();
0465   h_pT_truthCuts_max->Sumw2();
0466   h_eta_truth->Sumw2();
0467   h_eta_truthCuts->Sumw2();
0468   h_eta_truthCuts_min->Sumw2();
0469   h_eta_truthCuts_max->Sumw2();
0470   h_phi_truth->Sumw2();
0471   h_phi_truthCuts->Sumw2();
0472   h_phi_truthCuts_min->Sumw2();
0473   h_phi_truthCuts_max->Sumw2();
0474   h_y_truth->Sumw2();
0475   h_y_truthCuts->Sumw2();
0476   h_y_truthCuts_min->Sumw2();
0477   h_y_truthCuts_max->Sumw2();
0478 
0479   sim_tree_Xi->Draw("Lambda0_pT>>+h_pT_truth",truthCut.c_str());
0480   sim_tree_Xi->Draw("Lambda0_pT>>+h_pT_truthCuts",truthDataCut.c_str());
0481   sim_tree_Xi->Draw("Lambda0_pT>>+h_pT_truthCuts_min",truthDataCut_min.c_str());
0482   sim_tree_Xi->Draw("Lambda0_pT>>+h_pT_truthCuts_max",truthDataCut_max.c_str());
0483   sim_tree_Xi->Draw("Lambda0_pseudorapidity>>+h_eta_truth",truthCut.c_str());
0484   sim_tree_Xi->Draw("Lambda0_pseudorapidity>>+h_eta_truthCuts",truthDataCut.c_str());
0485   sim_tree_Xi->Draw("Lambda0_pseudorapidity>>+h_eta_truthCuts_min",truthDataCut_min.c_str());
0486   sim_tree_Xi->Draw("Lambda0_pseudorapidity>>+h_eta_truthCuts_max",truthDataCut_max.c_str());
0487   sim_tree_Xi->Draw("Lambda0_phi>>+h_phi_truth",truthCut.c_str());
0488   sim_tree_Xi->Draw("Lambda0_phi>>+h_phi_truthCuts",truthDataCut.c_str());
0489   sim_tree_Xi->Draw("Lambda0_phi>>+h_phi_truthCuts_min",truthDataCut_min.c_str());
0490   sim_tree_Xi->Draw("Lambda0_phi>>+h_phi_truthCuts_max",truthDataCut_max.c_str());
0491   sim_tree_Xi->Draw("Lambda0_rapidity>>+h_y_truth",truthCut.c_str());
0492   sim_tree_Xi->Draw("Lambda0_rapidity>>+h_y_truthCuts",truthDataCut.c_str());
0493   sim_tree_Xi->Draw("Lambda0_rapidity>>+h_y_truthCuts_min",truthDataCut_min.c_str());
0494   sim_tree_Xi->Draw("Lambda0_rapidity>>+h_y_truthCuts_max",truthDataCut_max.c_str());
0495 
0496   //string sim_inFile_KS0 = "/sphenix/user/cdean/software/analysis/LightFlavorRatios/geometric_acceptance/simulation/outputKFParticle_Kshort_reco.root";
0497   string sim_inFile_KS0 = "/gpfs/mnt/gpfs02/sphenix/user/cdean/software/analysis/LightFlavorRatios/geometric_acceptance/simulation/outputKFParticle_Kshort2pipi_reco_Usman_patch.root";
0498   TFile* sim_file_KS0 = new TFile(sim_inFile_KS0.c_str());
0499   TTree* sim_tree_KS0 = (TTree*)sim_file_KS0->Get("DecayTree");
0500 
0501   std::string truthCut_KS0 = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10";
0502   std::string truthDataCut_KS0 = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && abs(track_1_IP_xy) >= 0.05 && abs(track_2_IP_xy) >= 0.05 && K_S0_DIRA >= 0.99 && abs(K_S0_mass - 0.5) <= 0.1 && K_S0_decayLength >= 0.05 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0503   std::string truthDataCut_KS0_min = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && (abs(track_1_IP_xy) - " + std::to_string(KS0_PV_sigma_xy) + " - " + std::to_string(KS0_track_1_sigma_xy) + ") >= 0.05 && (abs(track_2_IP_xy) - " + std::to_string(KS0_PV_sigma_xy) + " - " + std::to_string(KS0_track_2_sigma_xy) + ") >= 0.05 && (track_1_track_2_DCA + " + std::to_string(KS0_track_1_sigma) + " + " + std::to_string(KS0_track_2_sigma) + ") <= 0.5 && (track_1_track_2_DCA_xy + " + std::to_string(KS0_track_1_sigma_xy) + " + " + std::to_string(KS0_track_2_sigma_xy) + ") <= 1.0 && cos(acos(K_S0_DIRA) + atan(" + std::to_string(KS0_PV_sigma) + "/K_S0_decayLength) + atan(" + std::to_string(KS0_SV_sigma) + "/K_S0_decayLength)) >= 0.99 && abs(K_S0_mass - 0.5) <= 0.1 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0504   std::string truthDataCut_KS0_max = "abs(track_1_true_ID) == 211 && Sum$(abs(track_1_true_track_history_PDG_ID) == 310) > 0 && abs(track_2_true_ID) == 211 && Sum$(abs(track_2_true_track_history_PDG_ID) == 310) > 0 && abs(K_S0_rapidity) <= 1.0 && track_1_MVTX_nHits > 0 && track_2_MVTX_nHits > 0 && abs(primary_vertex_z) < 10 && track_1_INTT_nHits > 0 && track_2_INTT_nHits > 0 && track_1_TPC_nHits > 19 && track_2_TPC_nHits > 19 && (track_1_chi2/track_1_nDoF) <= 300 && (track_2_chi2/track_2_nDoF) <= 300 && (abs(track_1_IP_xy) + " + std::to_string(KS0_PV_sigma_xy) + " + " + std::to_string(KS0_track_1_sigma_xy) + ") >= 0.05 && (abs(track_2_IP_xy) + " + std::to_string(KS0_PV_sigma_xy) + " + " + std::to_string(KS0_track_2_sigma_xy) + ") >= 0.05 && (track_1_track_2_DCA - " + std::to_string(KS0_track_1_sigma) + " - " + std::to_string(KS0_track_2_sigma) + ") <= 0.5 && (track_1_track_2_DCA_xy - " + std::to_string(KS0_track_1_sigma_xy) + " - " + std::to_string(KS0_track_2_sigma_xy) + ") <= 1.0 && cos(acos(K_S0_DIRA) - atan(" + std::to_string(KS0_PV_sigma) + "/K_S0_decayLength) - atan(" + std::to_string(KS0_SV_sigma) + "/K_S0_decayLength)) >= 0.99 && abs(K_S0_mass - 0.5) <= 0.1 && track_1_pT >= 0.2 && track_2_pT >= 0.2";
0505 
0506   TH1F* h_pT_truth_KS0 = new TH1F("h_pT_truth_KS0",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0507   TH1F* h_pT_truthCuts_KS0 = new TH1F("h_pT_truthCuts_KS0",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0508   TH1F* h_pT_truthCuts_KS0_min = new TH1F("h_pT_truthCuts_KS0_min",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0509   TH1F* h_pT_truthCuts_KS0_max = new TH1F("h_pT_truthCuts_KS0_max",";K_{S}^{0} p_{T} [GeV/c];Cut Efficiency",nBins_pT,pTbins.data());
0510   
0511   TH1F* h_eta_truth_KS0 = new TH1F("h_eta_truth_KS0",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins.data());
0512   TH1F* h_eta_truthCuts_KS0 = new TH1F("h_eta_truthCuts_KS0",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins.data());
0513   TH1F* h_eta_truthCuts_KS0_min = new TH1F("h_eta_truthCuts_KS0_min",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins.data());
0514   TH1F* h_eta_truthCuts_KS0_max = new TH1F("h_eta_truthCuts_KS0_max",";K_{S}^{0} #eta;Cut Efficiency",nBins_eta,etabins.data());
0515   
0516   TH1F* h_phi_truth_KS0 = new TH1F("h_phi_truth_KS0",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins.data());
0517   TH1F* h_phi_truthCuts_KS0 = new TH1F("h_phi_truthCuts_KS0",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins.data());
0518   TH1F* h_phi_truthCuts_KS0_min = new TH1F("h_phi_truthCuts_KS0_min",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins.data());
0519   TH1F* h_phi_truthCuts_KS0_max = new TH1F("h_phi_truthCuts_KS0_max",";K_{S}^{0} #phi;Cut Efficiency",nBins_phi,phibins.data());
0520   
0521   TH1F* h_y_truth_KS0 = new TH1F("h_y_truth_KS0",";K_{S}^{0} y;Cut Efficiency",nBins_rap,rapbins.data());
0522   TH1F* h_y_truthCuts_KS0 = new TH1F("h_y_truthCuts_KS0",";K_{S}^{0} y;Cut Efficiency",nBins_rap,rapbins.data());
0523   TH1F* h_y_truthCuts_KS0_min = new TH1F("h_y_truthCuts_KS0_min",";K_{S}^{0} y;Cut Efficiency",nBins_rap,rapbins.data());
0524   TH1F* h_y_truthCuts_KS0_max = new TH1F("h_y_truthCuts_KS0_max",";K_{S}^{0} y;Cut Efficiency",nBins_rap,rapbins.data());
0525    
0526   h_pT_truth_KS0->Sumw2();
0527   h_pT_truthCuts_KS0->Sumw2();
0528   h_pT_truthCuts_KS0_min->Sumw2();
0529   h_pT_truthCuts_KS0_max->Sumw2();
0530   h_eta_truth_KS0->Sumw2();
0531   h_eta_truthCuts_KS0->Sumw2();
0532   h_eta_truthCuts_KS0_min->Sumw2();
0533   h_eta_truthCuts_KS0_max->Sumw2();
0534   h_phi_truth_KS0->Sumw2();
0535   h_phi_truthCuts_KS0->Sumw2();
0536   h_phi_truthCuts_KS0_min->Sumw2();
0537   h_phi_truthCuts_KS0_max->Sumw2();
0538   h_y_truth_KS0->Sumw2();
0539   h_y_truthCuts_KS0->Sumw2();
0540   h_y_truthCuts_KS0_min->Sumw2();
0541   h_y_truthCuts_KS0_max->Sumw2();
0542 
0543   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truth_KS0",truthCut_KS0.c_str());
0544   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truthCuts_KS0",truthDataCut_KS0.c_str());
0545   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truthCuts_KS0_min",truthDataCut_KS0_min.c_str());
0546   sim_tree_KS0->Draw("K_S0_pT>>+h_pT_truthCuts_KS0_max",truthDataCut_KS0_max.c_str());
0547   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truth_KS0",truthCut_KS0.c_str());
0548   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truthCuts_KS0",truthDataCut_KS0.c_str());
0549   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truthCuts_KS0_min",truthDataCut_KS0_min.c_str());
0550   sim_tree_KS0->Draw("K_S0_pseudorapidity>>+h_eta_truthCuts_KS0_max",truthDataCut_KS0_max.c_str());
0551   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truth_KS0",truthCut_KS0.c_str());
0552   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truthCuts_KS0",truthDataCut_KS0.c_str());
0553   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truthCuts_KS0_min",truthDataCut_KS0_min.c_str());
0554   sim_tree_KS0->Draw("K_S0_phi>>+h_phi_truthCuts_KS0_max",truthDataCut_KS0_max.c_str());
0555   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truth_KS0",truthCut_KS0.c_str());
0556   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truthCuts_KS0",truthDataCut_KS0.c_str());
0557   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truthCuts_KS0_min",truthDataCut_KS0_min.c_str());
0558   sim_tree_KS0->Draw("K_S0_rapidity>>+h_y_truthCuts_KS0_max",truthDataCut_KS0_max.c_str());
0559 
0560   TH1F* hEff_pT = (TH1F*)h_pT_truth->Clone();
0561   TH1F* hEff_eta = (TH1F*)h_eta_truth->Clone();
0562   TH1F* hEff_phi = (TH1F*)h_phi_truth->Clone();
0563   TH1F* hEff_y = (TH1F*)h_y_truth->Clone();
0564 
0565   TH1F* hEff_pT_min = (TH1F*)h_pT_truth->Clone();
0566   TH1F* hEff_eta_min = (TH1F*)h_eta_truth->Clone();
0567   TH1F* hEff_phi_min = (TH1F*)h_phi_truth->Clone();
0568   TH1F* hEff_y_min = (TH1F*)h_y_truth->Clone();
0569 
0570   TH1F* hEff_pT_max = (TH1F*)h_pT_truth->Clone();
0571   TH1F* hEff_eta_max = (TH1F*)h_eta_truth->Clone();
0572   TH1F* hEff_phi_max = (TH1F*)h_phi_truth->Clone();
0573   TH1F* hEff_y_max = (TH1F*)h_y_truth->Clone();
0574 
0575   TH1F* hEff_pT_KS0 = (TH1F*)h_pT_truth_KS0->Clone();
0576   TH1F* hEff_eta_KS0 = (TH1F*)h_eta_truth_KS0->Clone();
0577   TH1F* hEff_phi_KS0 = (TH1F*)h_phi_truth_KS0->Clone();
0578   TH1F* hEff_y_KS0 = (TH1F*)h_y_truth_KS0->Clone();
0579 
0580   TH1F* hEff_pT_KS0_min = (TH1F*)h_pT_truth_KS0->Clone();
0581   TH1F* hEff_eta_KS0_min = (TH1F*)h_eta_truth_KS0->Clone();
0582   TH1F* hEff_phi_KS0_min = (TH1F*)h_phi_truth_KS0->Clone();
0583   TH1F* hEff_y_KS0_min = (TH1F*)h_y_truth_KS0->Clone();
0584 
0585   TH1F* hEff_pT_KS0_max = (TH1F*)h_pT_truth_KS0->Clone();
0586   TH1F* hEff_eta_KS0_max = (TH1F*)h_eta_truth_KS0->Clone();
0587   TH1F* hEff_phi_KS0_max = (TH1F*)h_phi_truth_KS0->Clone();
0588   TH1F* hEff_y_KS0_max = (TH1F*)h_y_truth_KS0->Clone();
0589 
0590   hEff_pT->SetName("eff_pT");
0591   hEff_eta->SetName("eff_eta");
0592   hEff_phi->SetName("eff_phi");
0593   hEff_y->SetName("eff_y");
0594 
0595   hEff_pT_min->SetName("eff_pT_min");
0596   hEff_eta_min->SetName("eff_eta_min");
0597   hEff_phi_min->SetName("eff_phi_min");
0598   hEff_y_min->SetName("eff_y_min");
0599 
0600   hEff_pT_max->SetName("eff_pT_max");
0601   hEff_eta_max->SetName("eff_eta_max");
0602   hEff_phi_max->SetName("eff_phi_max");
0603   hEff_y_max->SetName("eff_y_max");
0604 
0605   hEff_pT_KS0->SetName("eff_pT_KS0");
0606   hEff_eta_KS0->SetName("eff_eta_KS0");
0607   hEff_phi_KS0->SetName("eff_phi_KS0");
0608   hEff_y_KS0->SetName("eff_y_KS0");
0609 
0610   hEff_pT_KS0_min->SetName("eff_pT_KS0_min");
0611   hEff_eta_KS0_min->SetName("eff_eta_KS0_min");
0612   hEff_phi_KS0_min->SetName("eff_phi_KS0_min");
0613   hEff_y_KS0_min->SetName("eff_y_KS0_min");
0614 
0615   hEff_pT_KS0_max->SetName("eff_pT_KS0_max");
0616   hEff_eta_KS0_max->SetName("eff_eta_KS0_max");
0617   hEff_phi_KS0_max->SetName("eff_phi_KS0_max");
0618   hEff_y_KS0_max->SetName("eff_y_KS0_max");
0619 
0620   hEff_pT->Divide(h_pT_truthCuts,h_pT_truth);
0621   hEff_eta->Divide(h_eta_truthCuts,h_eta_truth);
0622   hEff_phi->Divide(h_phi_truthCuts,h_phi_truth);
0623   hEff_y->Divide(h_y_truthCuts,h_y_truth);
0624 
0625   hEff_pT_min->Divide(h_pT_truthCuts_min,h_pT_truth);
0626   hEff_eta_min->Divide(h_eta_truthCuts_min,h_eta_truth);
0627   hEff_phi_min->Divide(h_phi_truthCuts_min,h_phi_truth);
0628   hEff_y_min->Divide(h_y_truthCuts_min,h_y_truth);
0629 
0630   hEff_pT_max->Divide(h_pT_truthCuts,h_pT_truth);
0631   hEff_eta_max->Divide(h_eta_truthCuts,h_eta_truth);
0632   hEff_phi_max->Divide(h_phi_truthCuts,h_phi_truth);
0633   hEff_y_max->Divide(h_y_truthCuts,h_y_truth);
0634 
0635   hEff_pT_KS0->Divide(h_pT_truthCuts_KS0,h_pT_truth_KS0);
0636   hEff_eta_KS0->Divide(h_eta_truthCuts_KS0,h_eta_truth_KS0);
0637   hEff_phi_KS0->Divide(h_phi_truthCuts_KS0,h_phi_truth_KS0);
0638   hEff_y_KS0->Divide(h_y_truthCuts_KS0,h_y_truth_KS0);
0639 
0640   hEff_pT_KS0_min->Divide(h_pT_truthCuts_KS0_min,h_pT_truth_KS0);
0641   hEff_eta_KS0_min->Divide(h_eta_truthCuts_KS0_min,h_eta_truth_KS0);
0642   hEff_phi_KS0_min->Divide(h_phi_truthCuts_KS0_min,h_phi_truth_KS0);
0643   hEff_y_KS0_min->Divide(h_y_truthCuts_KS0_min,h_y_truth_KS0);
0644 
0645   hEff_pT_KS0_max->Divide(h_pT_truthCuts_KS0_max,h_pT_truth_KS0);
0646   hEff_eta_KS0_max->Divide(h_eta_truthCuts_KS0_max,h_eta_truth_KS0);
0647   hEff_phi_KS0_max->Divide(h_phi_truthCuts_KS0_max,h_phi_truth_KS0);
0648   hEff_y_KS0_max->Divide(h_y_truthCuts_KS0_max,h_y_truth_KS0);
0649 
0650   TH1F* hEffRatio_pT = (TH1F*)hEff_pT->Clone();
0651   TH1F* hEffRatio_eta = (TH1F*)hEff_eta->Clone();
0652   TH1F* hEffRatio_phi = (TH1F*)hEff_phi->Clone();
0653   TH1F* hEffRatio_y = (TH1F*)hEff_y->Clone();
0654 
0655   TH1F* hEffRatio_pT_min = (TH1F*)hEff_pT_min->Clone();
0656   TH1F* hEffRatio_eta_min = (TH1F*)hEff_eta_min->Clone();
0657   TH1F* hEffRatio_phi_min = (TH1F*)hEff_phi_min->Clone();
0658   TH1F* hEffRatio_y_min = (TH1F*)hEff_y_min->Clone();
0659 
0660   TH1F* hEffRatio_pT_max = (TH1F*)hEff_pT_max->Clone();
0661   TH1F* hEffRatio_eta_max = (TH1F*)hEff_eta_max->Clone();
0662   TH1F* hEffRatio_phi_max = (TH1F*)hEff_phi_max->Clone();
0663   TH1F* hEffRatio_y_max = (TH1F*)hEff_y_max->Clone();
0664 
0665   hEffRatio_pT->SetName("hEffRatio_pT");
0666   hEffRatio_eta->SetName("hEffRatio_eta");
0667   hEffRatio_phi->SetName("hEffRatio_phi");
0668   hEffRatio_y->SetName("hEffRatio_y");
0669 
0670   hEffRatio_pT_min->SetName("hEffRatio_pT_min");
0671   hEffRatio_eta_min->SetName("hEffRatio_eta_min");
0672   hEffRatio_phi_min->SetName("hEffRatio_phi_min");
0673   hEffRatio_y_min->SetName("hEffRatio_y_min");
0674 
0675   hEffRatio_pT_max->SetName("hEffRatio_pT_max");
0676   hEffRatio_eta_max->SetName("hEffRatio_eta_max");
0677   hEffRatio_phi_max->SetName("hEffRatio_phi_max");
0678   hEffRatio_y_max->SetName("hEffRatio_y_max");
0679 
0680   hEffRatio_pT->Divide(hEff_pT,hEff_pT_KS0);
0681   hEffRatio_eta->Divide(hEff_eta,hEff_eta_KS0);
0682   hEffRatio_phi->Divide(hEff_phi,hEff_phi_KS0);
0683   hEffRatio_y->Divide(hEff_y,hEff_y_KS0);
0684 
0685   hEffRatio_pT_min->Divide(hEff_pT_min,hEff_pT_KS0_min);
0686   hEffRatio_eta_min->Divide(hEff_eta_min,hEff_eta_KS0_min);
0687   hEffRatio_phi_min->Divide(hEff_phi_min,hEff_phi_KS0_min);
0688   hEffRatio_y_min->Divide(hEff_y_min,hEff_y_KS0_min);
0689 
0690   hEffRatio_pT_max->Divide(hEff_pT_max,hEff_pT_KS0_max);
0691   hEffRatio_eta_max->Divide(hEff_eta_max,hEff_eta_KS0_max);
0692   hEffRatio_phi_max->Divide(hEff_phi_max,hEff_phi_KS0_max);
0693   hEffRatio_y_max->Divide(hEff_y_max,hEff_y_KS0_max);
0694 
0695   TFile* fout = new TFile("LamdbaKsCutEfficiency_200MeV_hists.root","RECREATE");
0696 
0697   hEff_pT->Write();
0698   hEff_eta->Write();
0699   hEff_phi->Write();
0700   hEff_y->Write();
0701 
0702   hEff_pT_KS0->Write();
0703   hEff_eta_KS0->Write();
0704   hEff_phi_KS0->Write();
0705   hEff_y_KS0->Write();
0706 
0707   hEffRatio_pT->Write();
0708   hEffRatio_eta->Write();
0709   hEffRatio_phi->Write();
0710   hEffRatio_y->Write();
0711 
0712   hEff_pT_min->Write();
0713   hEff_eta_min->Write();
0714   hEff_phi_min->Write();
0715   hEff_y_min->Write();
0716 
0717   hEff_pT_KS0_min->Write();
0718   hEff_eta_KS0_min->Write();
0719   hEff_phi_KS0_min->Write();
0720   hEff_y_KS0_min->Write();
0721 
0722   hEffRatio_pT_min->Write();
0723   hEffRatio_eta_min->Write();
0724   hEffRatio_phi_min->Write();
0725   hEffRatio_y_min->Write();
0726 
0727   hEff_pT_max->Write();
0728   hEff_eta_max->Write();
0729   hEff_phi_max->Write();
0730   hEff_y_max->Write();
0731 
0732   hEff_pT_KS0_max->Write();
0733   hEff_eta_KS0_max->Write();
0734   hEff_phi_KS0_max->Write();
0735   hEff_y_KS0_max->Write();
0736 
0737   hEffRatio_pT_max->Write();
0738   hEffRatio_eta_max->Write();
0739   hEffRatio_phi_max->Write();
0740   hEffRatio_y_max->Write();
0741 
0742   TH1F* syserr_pT = new TH1F("cuteff_syserr_pT","systematic error from cut efficiency",nBins_pT,pTbins.data());
0743   TH1F* syserr_eta = new TH1F("cuteff_syserr_eta","systematic error from cut efficiency",nBins_eta,etabins.data());
0744   TH1F* syserr_phi = new TH1F("cuteff_syserr_phi","systematic error from cut efficiency",nBins_phi,phibins.data());
0745   TH1F* syserr_y = new TH1F("cuteff_syserr_y","systematic error from cut efficiency",nBins_rap,rapbins.data());
0746 
0747   for(int i=1;i<=syserr_pT->GetNbinsX();i++)
0748   {
0749     std::pair<double,double> extremes = std::minmax({hEffRatio_pT_min->GetBinContent(i),hEffRatio_pT->GetBinContent(i),hEffRatio_pT_max->GetBinContent(i)});
0750     syserr_pT->SetBinContent(i,0.5*fabs(extremes.second-extremes.first));
0751   }
0752 
0753   for(int i=1;i<=syserr_eta->GetNbinsX();i++)
0754   {
0755     std::pair<double,double> extremes = std::minmax({hEffRatio_eta_min->GetBinContent(i),hEffRatio_eta->GetBinContent(i),hEffRatio_eta_max->GetBinContent(i)});
0756     syserr_eta->SetBinContent(i,0.5*fabs(extremes.second-extremes.first));
0757   }
0758 
0759   for(int i=1;i<=syserr_phi->GetNbinsX();i++)
0760   {
0761     std::pair<double,double> extremes = std::minmax({hEffRatio_phi_min->GetBinContent(i),hEffRatio_phi->GetBinContent(i),hEffRatio_phi_max->GetBinContent(i)});
0762     syserr_phi->SetBinContent(i,0.5*fabs(extremes.second-extremes.first));
0763   }
0764 
0765   for(int i=1;i<=syserr_y->GetNbinsX();i++)
0766   {
0767     std::pair<double,double> extremes = std::minmax({hEffRatio_y_min->GetBinContent(i),hEffRatio_y->GetBinContent(i),hEffRatio_y_max->GetBinContent(i)});
0768     syserr_y->SetBinContent(i,0.5*fabs(extremes.second-extremes.first));
0769   }
0770 
0771   syserr_pT->Write();
0772   syserr_eta->Write();
0773   syserr_phi->Write();
0774   syserr_y->Write();
0775 
0776   fout->Close();
0777 
0778   // -------------------------------------------------------------------------
0779   // Build TEfficiency objects for Xi and KS0
0780   // -------------------------------------------------------------------------
0781   TEfficiency* eff_pT_Xi  = new TEfficiency(*h_pT_truthCuts,  *h_pT_truth);
0782   TEfficiency* eff_eta_Xi = new TEfficiency(*h_eta_truthCuts, *h_eta_truth);
0783   TEfficiency* eff_phi_Xi = new TEfficiency(*h_phi_truthCuts, *h_phi_truth);
0784   TEfficiency* eff_y_Xi   = new TEfficiency(*h_y_truthCuts,   *h_y_truth);
0785 
0786   TEfficiency* eff_pT_KS0  = new TEfficiency(*h_pT_truthCuts_KS0,  *h_pT_truth_KS0);
0787   TEfficiency* eff_eta_KS0 = new TEfficiency(*h_eta_truthCuts_KS0, *h_eta_truth_KS0);
0788   TEfficiency* eff_phi_KS0 = new TEfficiency(*h_phi_truthCuts_KS0, *h_phi_truth_KS0);
0789   TEfficiency* eff_y_KS0   = new TEfficiency(*h_y_truthCuts_KS0,   *h_y_truth_KS0);
0790 
0791   TEfficiency* eff_pT_Xi_min  = new TEfficiency(*h_pT_truthCuts_min,  *h_pT_truth);
0792   TEfficiency* eff_eta_Xi_min = new TEfficiency(*h_eta_truthCuts_min, *h_eta_truth);
0793   TEfficiency* eff_phi_Xi_min = new TEfficiency(*h_phi_truthCuts_min, *h_phi_truth);
0794   TEfficiency* eff_y_Xi_min   = new TEfficiency(*h_y_truthCuts_min,   *h_y_truth);
0795 
0796   TEfficiency* eff_pT_KS0_min  = new TEfficiency(*h_pT_truthCuts_KS0_min,  *h_pT_truth_KS0);
0797   TEfficiency* eff_eta_KS0_min = new TEfficiency(*h_eta_truthCuts_KS0_min, *h_eta_truth_KS0);
0798   TEfficiency* eff_phi_KS0_min = new TEfficiency(*h_phi_truthCuts_KS0_min, *h_phi_truth_KS0);
0799   TEfficiency* eff_y_KS0_min   = new TEfficiency(*h_y_truthCuts_KS0_min,   *h_y_truth_KS0);
0800 
0801   TEfficiency* eff_pT_Xi_max  = new TEfficiency(*h_pT_truthCuts_max,  *h_pT_truth);
0802   TEfficiency* eff_eta_Xi_max = new TEfficiency(*h_eta_truthCuts_max, *h_eta_truth);
0803   TEfficiency* eff_phi_Xi_max = new TEfficiency(*h_phi_truthCuts_max, *h_phi_truth);
0804   TEfficiency* eff_y_Xi_max   = new TEfficiency(*h_y_truthCuts_max,   *h_y_truth);
0805 
0806   TEfficiency* eff_pT_KS0_max  = new TEfficiency(*h_pT_truthCuts_KS0_max,  *h_pT_truth_KS0);
0807   TEfficiency* eff_eta_KS0_max = new TEfficiency(*h_eta_truthCuts_KS0_max, *h_eta_truth_KS0);
0808   TEfficiency* eff_phi_KS0_max = new TEfficiency(*h_phi_truthCuts_KS0_max, *h_phi_truth_KS0);
0809   TEfficiency* eff_y_KS0_max   = new TEfficiency(*h_y_truthCuts_KS0_max,   *h_y_truth_KS0);
0810 
0811   for (auto* e : {eff_pT_Xi, eff_eta_Xi, eff_phi_Xi, eff_y_Xi,
0812           eff_pT_KS0, eff_eta_KS0, eff_phi_KS0, eff_y_KS0,
0813                   eff_pT_Xi_min, eff_eta_Xi_min, eff_phi_Xi_min, eff_y_Xi_min,
0814                   eff_pT_KS0_min, eff_eta_KS0_min, eff_phi_KS0_min, eff_y_KS0_min,
0815                   eff_pT_Xi_max, eff_eta_Xi_max, eff_phi_Xi_max, eff_y_Xi_max,
0816                   eff_pT_KS0_max, eff_eta_KS0_max, eff_phi_KS0_max, eff_y_KS0_max})
0817       e->SetStatisticOption(TEfficiency::kFCP);
0818 
0819   auto makeEffGraph = [](TEfficiency* eff) -> TGraphAsymmErrors*
0820   {
0821       int n = eff->GetTotalHistogram()->GetNbinsX();
0822       TGraphAsymmErrors* g = new TGraphAsymmErrors(n);
0823       g->SetTitle("");
0824       for (int i = 1; i <= n; ++i)
0825       {
0826       double x    = eff->GetTotalHistogram()->GetBinCenter(i);
0827       double exlo = eff->GetTotalHistogram()->GetBinCenter(i) - eff->GetTotalHistogram()->GetBinLowEdge(i);
0828       double exhi = eff->GetTotalHistogram()->GetBinLowEdge(i+1) - eff->GetTotalHistogram()->GetBinCenter(i);
0829       double y    = eff->GetEfficiency(i);
0830       double eylo = eff->GetEfficiencyErrorLow(i);
0831       double eyhi = eff->GetEfficiencyErrorUp(i);
0832       g->SetPoint(i-1, x, y);
0833       g->SetPointError(i-1, exlo, exhi, eylo, eyhi);
0834       }
0835       return g;
0836   };
0837 
0838   auto makeRatioGraph = [](TGraphAsymmErrors* gXi, TGraphAsymmErrors* gKS0,
0839                bool XiOverKS0) -> TGraphAsymmErrors*
0840   {
0841       int n = gXi->GetN();
0842       TGraphAsymmErrors* ratio = new TGraphAsymmErrors(n);
0843       ratio->SetTitle("");
0844       for (int i = 0; i < n; ++i)
0845       {
0846       double x, yXi, yKS0;
0847       gXi ->GetPoint(i, x, yXi);
0848       double dummy;
0849       gKS0->GetPoint(i, dummy, yKS0);
0850 
0851       double exlo = gXi->GetErrorXlow(i);
0852       double exhi = gXi->GetErrorXhigh(i);
0853 
0854       // Numerator and denominator depend on which ratio we want
0855       double num    = XiOverKS0 ? yXi  : yKS0;
0856       double den    = XiOverKS0 ? yKS0 : yXi;
0857 
0858       double numElo = XiOverKS0 ? gXi->GetErrorYlow(i)  : gKS0->GetErrorYlow(i);
0859       double numEhi = XiOverKS0 ? gXi->GetErrorYhigh(i) : gKS0->GetErrorYhigh(i);
0860       double denElo = XiOverKS0 ? gKS0->GetErrorYlow(i) : gXi->GetErrorYlow(i);
0861       double denEhi = XiOverKS0 ? gKS0->GetErrorYhigh(i): gXi->GetErrorYhigh(i);
0862 
0863       if (den <= 0.)
0864       {
0865           ratio->SetPoint(i, x, 0.);
0866           ratio->SetPointError(i, exlo, exhi, 0., 0.);
0867           continue;
0868       }
0869 
0870       double r = num / den;
0871 
0872       // Error propagation for ratio f = a/b:
0873       // delta_f = f * sqrt((delta_a/a)^2 + (delta_b/b)^2)
0874       // Use asymmetric errors: propagate lo and hi separately
0875       double eylo = (num > 0.) ? r * sqrt(pow(numElo/num, 2) + pow(denEhi/den, 2)) : 0.;
0876       double eyhi = (num > 0.) ? r * sqrt(pow(numEhi/num, 2) + pow(denElo/den, 2)) : 0.;
0877 
0878       ratio->SetPoint(i, x, r);
0879       ratio->SetPointError(i, exlo, exhi, eylo, eyhi);
0880       }
0881       return ratio;
0882   };
0883 
0884   TGraphAsymmErrors* g_pT_Xi  = makeEffGraph(eff_pT_Xi);
0885   TGraphAsymmErrors* g_pT_KS0 = makeEffGraph(eff_pT_KS0);
0886   TGraphAsymmErrors* g_eta_Xi  = makeEffGraph(eff_eta_Xi);
0887   TGraphAsymmErrors* g_eta_KS0 = makeEffGraph(eff_eta_KS0);
0888   TGraphAsymmErrors* g_phi_Xi  = makeEffGraph(eff_phi_Xi);
0889   TGraphAsymmErrors* g_phi_KS0 = makeEffGraph(eff_phi_KS0);
0890   TGraphAsymmErrors* g_y_Xi    = makeEffGraph(eff_y_Xi);
0891   TGraphAsymmErrors* g_y_KS0   = makeEffGraph(eff_y_KS0);
0892 
0893   // Xi/KS0 ratios
0894   TGraphAsymmErrors* ratio_pT_XiOverKS0  = makeRatioGraph(g_pT_Xi,  g_pT_KS0,  true);
0895   TGraphAsymmErrors* ratio_eta_XiOverKS0 = makeRatioGraph(g_eta_Xi, g_eta_KS0, true);
0896   TGraphAsymmErrors* ratio_phi_XiOverKS0 = makeRatioGraph(g_phi_Xi, g_phi_KS0, true);
0897   TGraphAsymmErrors* ratio_y_XiOverKS0   = makeRatioGraph(g_y_Xi,   g_y_KS0,   true);
0898 
0899   // KS0/Xi ratios
0900   TGraphAsymmErrors* ratio_pT_KS0OverXi  = makeRatioGraph(g_pT_Xi,  g_pT_KS0,  false);
0901   TGraphAsymmErrors* ratio_eta_KS0OverXi = makeRatioGraph(g_eta_Xi, g_eta_KS0, false);
0902   TGraphAsymmErrors* ratio_phi_KS0OverXi = makeRatioGraph(g_phi_Xi, g_phi_KS0, false);
0903   TGraphAsymmErrors* ratio_y_KS0OverXi   = makeRatioGraph(g_y_Xi,   g_y_KS0,   false);
0904  
0905   TCanvas *c1  = new TCanvas("myCanvas", "myCanvas",800,800);
0906   ratio_pT_XiOverKS0->Draw("AP");
0907   ratio_pT_XiOverKS0->GetXaxis()->SetTitle("Mother p_{T} [GeV/c]");
0908   ratio_pT_XiOverKS0->GetYaxis()->SetTitle("#Lambda/K_{S}^{0} Cut Efficiency");
0909   ratio_pT_XiOverKS0->GetYaxis()->SetRangeUser(0.0, 1.05);
0910   ratio_pT_XiOverKS0->GetXaxis()->SetLimits(pTbins[0], pTbins[nBins_pT]);
0911   c1->Update();
0912   
0913   TPaveText *pt;
0914   pt = new TPaveText(0.15,0.75,0.65,1., "NDC");
0915   pt->SetFillColor(0);
0916   pt->SetFillStyle(0);
0917   pt->SetTextFont(42);
0918   string label = "#it{#bf{sPHENIX}} Simulation";
0919   TText *pt_LaTex = pt->AddText(label.c_str());
0920   pt->SetBorderSize(0);
0921   pt->Draw();
0922   gPad->Modified();
0923 
0924   TPaveText *ptDate;
0925   ptDate = new TPaveText(0.67,0.79,0.85,0.95, "NDC");
0926   ptDate->SetFillColor(0);
0927   ptDate->SetFillStyle(0);
0928   ptDate->SetTextFont(42);
0929   std::string compilation_date = getDate();
0930   TText *pt_LaTexDate = ptDate->AddText(compilation_date.c_str());
0931   ptDate->SetBorderSize(0);
0932   ptDate->Draw();
0933   gPad->Modified();
0934 
0935   string extensions[] = {"_250MeV.png","_250MeV.root"};
0936   for (auto extension : extensions)
0937   {
0938     string output = "LambdaKS0ratio_pT_cutEff" + extension;
0939     c1->SaveAs(output.c_str());
0940   }
0941   c1->Close();  
0942 
0943   TCanvas *c2  = new TCanvas("myCanvas", "myCanvas",800,800);
0944   ratio_eta_XiOverKS0->Draw("AP");
0945   ratio_eta_XiOverKS0->GetXaxis()->SetTitle("Mother #eta");
0946   ratio_eta_XiOverKS0->GetYaxis()->SetTitle("#Lambda/K_{S}^{0} Cut Efficiency");
0947   ratio_eta_XiOverKS0->GetYaxis()->SetRangeUser(0.0, 1.05);
0948   ratio_eta_XiOverKS0->GetXaxis()->SetLimits(etabins[0], etabins[nBins_eta]);
0949   c2->Update();
0950   
0951   pt->Draw();
0952   gPad->Modified();
0953 
0954   ptDate->Draw();
0955   gPad->Modified();
0956 
0957   for (auto extension : extensions)
0958   {
0959     string output = "LambdaKS0ratio_eta_cutEff" + extension;
0960     c2->SaveAs(output.c_str());
0961   }
0962   c2->Close();  
0963   
0964   TCanvas *c3  = new TCanvas("myCanvas", "myCanvas",800,800);
0965   ratio_phi_XiOverKS0->Draw("AP");
0966   ratio_phi_XiOverKS0->GetXaxis()->SetTitle("Mother #phi [rad]");
0967   ratio_phi_XiOverKS0->GetYaxis()->SetTitle("#Lambda/K_{S}^{0} Cut Efficiency");
0968   ratio_phi_XiOverKS0->GetYaxis()->SetRangeUser(0.0, 1.05);
0969   ratio_phi_XiOverKS0->GetXaxis()->SetLimits(phibins[0], phibins[nBins_phi]);
0970   c3->Update();
0971   
0972   pt->Draw();
0973   gPad->Modified();
0974 
0975   ptDate->Draw();
0976   gPad->Modified();
0977 
0978   for (auto extension : extensions)
0979   {
0980     string output = "LambdaKS0ratio_phi_cutEff" + extension;
0981     c3->SaveAs(output.c_str());
0982   }
0983   c3->Close();  
0984 
0985   TCanvas *c4  = new TCanvas("myCanvas", "myCanvas",800,800);
0986   ratio_y_XiOverKS0->Draw("AP");
0987   ratio_y_XiOverKS0->GetXaxis()->SetTitle("Mother Rapidity");
0988   ratio_y_XiOverKS0->GetYaxis()->SetTitle("#Lambda/K_{S}^{0} Cut Efficiency");
0989   ratio_y_XiOverKS0->GetYaxis()->SetRangeUser(0.0, 1.05);
0990   ratio_y_XiOverKS0->GetXaxis()->SetLimits(rapbins[0], rapbins[nBins_rap]);
0991   c4->Update();
0992   
0993   pt->Draw();
0994   gPad->Modified();
0995 
0996   ptDate->Draw();
0997   gPad->Modified();
0998 
0999   for (auto extension : extensions)
1000   {
1001     string output = "LambdaKS0ratio_y_cutEff" + extension;
1002     c4->SaveAs(output.c_str());
1003   }
1004   c4->Close();  
1005  
1006   /*
1007   c1->Divide(2,2);
1008   //pT 
1009   c1->cd(1);
1010   
1011   //ratio_pT_XiOverKS0->Draw("AP");
1012   //ratio_pT_XiOverKS0->GetXaxis()->SetTitle("Mother p_{T} [GeV/c]");
1013   //ratio_pT_XiOverKS0->GetYaxis()->SetTitle("#Xi^{-}/K_{S}^{0} Cut Efficiency");
1014   //ratio_pT_XiOverKS0->GetYaxis()->SetRangeUser(0.0, 1.05);
1015   //ratio_pT_XiOverKS0->GetXaxis()->SetLimits(pTbins[0], pTbins[nBins_pT]);
1016  
1017   ratio_pT_KS0OverXi->Draw("AP");
1018   ratio_pT_KS0OverXi->GetXaxis()->SetTitle("Mother p_{T} [GeV/c]");
1019   ratio_pT_KS0OverXi->GetYaxis()->SetTitle("K_{S}^{0}/#Xi^{-} Cut Efficiency");
1020   ratio_pT_KS0OverXi->GetYaxis()->SetRangeUser(0.0, 5.00);
1021   ratio_pT_KS0OverXi->GetXaxis()->SetLimits(pTbins[0], pTbins[nBins_pT]);
1022   c1->Update();
1023 
1024   TPaveText *pt;
1025   pt = new TPaveText(0.15,0.75,0.65,1., "NDC");
1026   pt->SetFillColor(0);
1027   pt->SetFillStyle(0);
1028   pt->SetTextFont(42);
1029   string label = "#it{#bf{sPHENIX}} Simulation";
1030   TText *pt_LaTex = pt->AddText(label.c_str());
1031   pt->SetBorderSize(0);
1032   pt->Draw();
1033   gPad->Modified();
1034 
1035   TPaveText *ptDate;
1036   ptDate = new TPaveText(0.67,0.79,0.85,0.95, "NDC");
1037   ptDate->SetFillColor(0);
1038   ptDate->SetFillStyle(0);
1039   ptDate->SetTextFont(42);
1040   std::string compilation_date = getDate();
1041   TText *pt_LaTexDate = ptDate->AddText(compilation_date.c_str());
1042   ptDate->SetBorderSize(0);
1043   ptDate->Draw();
1044   gPad->Modified();
1045  
1046   //eta
1047   c1->cd(2);
1048   ratio_eta_KS0OverXi->Draw("AP");
1049   ratio_eta_KS0OverXi->GetXaxis()->SetTitle("Mother #eta");
1050   ratio_eta_KS0OverXi->GetYaxis()->SetTitle("K_{S}^{0}/#Xi^{-} Cut Efficiency");
1051   ratio_eta_KS0OverXi->GetYaxis()->SetRangeUser(0.0, 5.0);
1052   ratio_eta_KS0OverXi->GetXaxis()->SetLimits(etabins[0], etabins[nBins_eta]);
1053   c1->Update();
1054 
1055   pt->Draw();
1056   gPad->Modified();
1057 
1058   ptDate->Draw();
1059   gPad->Modified();
1060   
1061   //phi
1062   c1->cd(3);
1063   ratio_phi_KS0OverXi->Draw("AP");
1064   ratio_phi_KS0OverXi->GetXaxis()->SetTitle("Mother #phi [rad]");
1065   ratio_phi_KS0OverXi->GetYaxis()->SetTitle("K_{S}^{0}/#Xi^{-} Cut Efficiency");
1066   ratio_phi_KS0OverXi->GetYaxis()->SetRangeUser(0.0, 5.0);
1067   ratio_phi_KS0OverXi->GetXaxis()->SetLimits(phibins[0], phibins[nBins_phi]);
1068   c1->Update();
1069 
1070   pt->Draw();
1071   gPad->Modified();
1072 
1073   ptDate->Draw();
1074   gPad->Modified();
1075   
1076   //rapidity
1077   c1->cd(4);
1078   ratio_y_KS0OverXi->Draw("AP");
1079   ratio_y_KS0OverXi->GetXaxis()->SetTitle("Mother Rapidity");
1080   ratio_y_KS0OverXi->GetYaxis()->SetTitle("K_{S}^{0}/#Xi^{-} Cut Efficiency");
1081   ratio_y_KS0OverXi->GetYaxis()->SetRangeUser(0.0, 5.0);
1082   ratio_y_KS0OverXi->GetXaxis()->SetLimits(-1.,1.);
1083   c1->Update();
1084 
1085   pt->Draw();
1086   gPad->Modified();
1087 
1088   ptDate->Draw();
1089   gPad->Modified();
1090    
1091 
1092   string extensions[] = {".png",".root"};
1093   for (auto extension : extensions)
1094   {
1095     //string output = "XiKS0ratio_all_cutEff" + extension;
1096     string output = "KS0Xiratio_all_cutEff" + extension;
1097     c1->SaveAs(output.c_str());
1098   }
1099   */
1100 }
[ 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!