CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Public Attributes
CTPPSProtonReconstructionSimulationValidator::PlotGroup Struct Reference

Public Member Functions

 PlotGroup ()
 
void write () const
 

Static Public Member Functions

static TGraphErrors profileToRMSGraph (TProfile *p, const char *name="")
 

Public Attributes

std::unique_ptr< TH1D > h_de_t
 
std::unique_ptr< TH1D > h_de_th_x
 
std::unique_ptr< TH1D > h_de_th_y
 
std::unique_ptr< TH1D > h_de_vtx_y
 
std::unique_ptr< TH1D > h_de_xi
 
std::unique_ptr< TH2D > h_xi_reco_vs_xi_simu
 
std::unique_ptr< TProfile > p_de_t_vs_t_simu
 
std::unique_ptr< TProfile > p_de_t_vs_xi_simu
 
std::unique_ptr< TProfile > p_de_th_x_vs_xi_simu
 
std::unique_ptr< TProfile > p_de_th_y_vs_xi_simu
 
std::unique_ptr< TProfile > p_de_vtx_y_vs_xi_simu
 
std::unique_ptr< TProfile > p_de_xi_vs_xi_simu
 

Detailed Description

Definition at line 57 of file CTPPSProtonReconstructionSimulationValidator.cc.

Constructor & Destructor Documentation

CTPPSProtonReconstructionSimulationValidator::PlotGroup::PlotGroup ( )
inline

Definition at line 75 of file CTPPSProtonReconstructionSimulationValidator.cc.

75  :
76  h_de_xi(new TH1D("", ";#xi_{reco} - #xi_{simu}", 100, 0., 0.)),
77  p_de_xi_vs_xi_simu(new TProfile("", ";#xi_{simu};#xi_{reco} - #xi_{simu}", 19, 0.015, 0.205)),
78  h_xi_reco_vs_xi_simu(new TH2D("", ";#xi_{simu};#xi_{reco}", 100, 0., 0.30, 100, 0., 0.30)),
79 
80  h_de_th_x(new TH1D("", ";#theta_{x,reco} - #theta_{x,simu}", 100, 0., 0.)),
81  p_de_th_x_vs_xi_simu(new TProfile("", ";#xi_{simu};#theta_{x,reco} - #theta_{x,simu}", 19, 0.015, 0.205)),
82 
83  h_de_th_y(new TH1D("", ";#theta_{y,reco} - #theta_{y,simu}", 100, 0., 0.)),
84  p_de_th_y_vs_xi_simu(new TProfile("", ";#xi_{simu};#theta_{y,reco} - #theta_{y,simu}", 19, 0.015, 0.205)),
85 
86  h_de_vtx_y(new TH1D("", ";vtx_{y,reco} - vtx_{y,simu} (mm)", 100, 0., 0.)),
87  p_de_vtx_y_vs_xi_simu(new TProfile("", ";#xi_{simu};vtx_{y,reco} - vtx_{y,simu} (mm)", 19, 0.015, 0.205)),
88 
89  h_de_t(new TH1D("", ";t_{reco} - t_{simu}", 100, -1., +1.)),
90  p_de_t_vs_xi_simu(new TProfile("", ";xi_{simu};t_{reco} - t_{simu}", 19, 0.015, 0.205)),
91  p_de_t_vs_t_simu(new TProfile("", ";t_{simu};t_{reco} - t_{simu}", 20, 0., 5.))
92  {}

Member Function Documentation

static TGraphErrors CTPPSProtonReconstructionSimulationValidator::PlotGroup::profileToRMSGraph ( TProfile *  p,
const char *  name = "" 
)
inlinestatic

Definition at line 94 of file CTPPSProtonReconstructionSimulationValidator.cc.

References EnergyCorrector::c, N, dataset::name, mathSSE::sqrt(), and w.

Referenced by write().

95  {
96  TGraphErrors gr_err;
97  gr_err.SetName(name);
98 
99  for (int bi = 1; bi <= p->GetNbinsX(); ++bi) {
100  double c = p->GetBinCenter(bi);
101  double w = p->GetBinWidth(bi);
102 
103  double N = p->GetBinEntries(bi);
104  double Sy = p->GetBinContent(bi) * N;
105  double Syy = p->GetSumw2()->At(bi);
106 
107  double si_sq = Syy/N - Sy*Sy/N/N;
108  double si = (si_sq >= 0.) ? sqrt(si_sq) : 0.;
109  double si_unc_sq = si_sq / 2. / N; // Gaussian approximation
110  double si_unc = (si_unc_sq >= 0.) ? sqrt(si_unc_sq) : 0.;
111 
112  int idx = gr_err.GetN();
113  gr_err.SetPoint(idx, c, si);
114  gr_err.SetPointError(idx, w/2., si_unc);
115  }
116 
117  return gr_err;
118  }
const double w
Definition: UKUtility.cc:23
T sqrt(T t)
Definition: SSEVec.h:18
#define N
Definition: blowfish.cc:9
void CTPPSProtonReconstructionSimulationValidator::PlotGroup::write ( ) const
inline

Definition at line 120 of file CTPPSProtonReconstructionSimulationValidator.cc.

References profileToRMSGraph().

121  {
122  h_xi_reco_vs_xi_simu->Write("h_xi_reco_vs_xi_simu");
123  h_de_xi->Write("h_de_xi");
124  p_de_xi_vs_xi_simu->Write("p_de_xi_vs_xi_simu");
125  profileToRMSGraph(p_de_xi_vs_xi_simu.get(), "g_rms_de_xi_vs_xi_simu").Write();
126 
127  h_de_th_x->Write("h_de_th_x");
128  p_de_th_x_vs_xi_simu->Write("p_de_th_x_vs_xi_simu");
129  profileToRMSGraph(p_de_th_x_vs_xi_simu.get(), "g_rms_de_th_x_vs_xi_simu").Write();
130 
131  h_de_th_y->Write("h_de_th_y");
132  p_de_th_y_vs_xi_simu->Write("p_de_th_y_vs_xi_simu");
133  profileToRMSGraph(p_de_th_y_vs_xi_simu.get(), "g_rms_de_th_y_vs_xi_simu").Write();
134 
135  h_de_vtx_y->Write("h_de_vtx_y");
136  p_de_vtx_y_vs_xi_simu->Write("p_de_vtx_y_vs_xi_simu");
137  profileToRMSGraph(p_de_vtx_y_vs_xi_simu.get(), "g_rms_de_vtx_y_vs_xi_simu").Write();
138 
139  h_de_t->Write("h_de_t");
140  p_de_t_vs_xi_simu->Write("p_de_t_vs_xi_simu");
141  profileToRMSGraph(p_de_t_vs_xi_simu.get(), "g_rms_de_t_vs_xi_simu").Write();
142  p_de_t_vs_t_simu->Write("p_de_t_vs_t_simu");
143  profileToRMSGraph(p_de_t_vs_t_simu.get(), "g_rms_de_t_vs_t_simu").Write();
144  }
static TGraphErrors profileToRMSGraph(TProfile *p, const char *name="")

Member Data Documentation

std::unique_ptr<TH1D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_de_t
std::unique_ptr<TH1D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_de_th_x
std::unique_ptr<TH1D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_de_th_y
std::unique_ptr<TH1D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_de_vtx_y
std::unique_ptr<TH1D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_de_xi
std::unique_ptr<TH2D> CTPPSProtonReconstructionSimulationValidator::PlotGroup::h_xi_reco_vs_xi_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_t_vs_t_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_t_vs_xi_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_th_x_vs_xi_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_th_y_vs_xi_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_vtx_y_vs_xi_simu
std::unique_ptr<TProfile> CTPPSProtonReconstructionSimulationValidator::PlotGroup::p_de_xi_vs_xi_simu