File indexing completed on 2023-03-17 11:20:48
0001 #ifndef RecoMuon_MuonSeedValidatorHisto_H
0002 #define RecoMuon_MuonSeedValidatorHisto_H
0003
0004
0005
0006
0007
0008
0009
0010 #include "TH1F.h"
0011 #include "TH2F.h"
0012 #include "TFile.h"
0013 #include "TString.h"
0014 #include <string>
0015 #include <iostream>
0016
0017 class H2DRecHit1 {
0018 public:
0019 H2DRecHit1(std::string name_) {
0020 TString N1 = name_.c_str();
0021 name = N1;
0022
0023 heta_mu4 = new TH2F(N1 + "_heta_mu4", " Seg multiplicity vs eta_of_track", 59, -2.95, 2.95, 40, -0.25, 19.75);
0024 heta_mu3 = new TH2F(N1 + "_heta_mu3", " good Rec Seg vs eta_of_track", 59, -2.95, 2.95, 20, -0.25, 9.75);
0025 heta_mu2 = new TH2F(N1 + "_heta_mu2", " Rec Seg vs eta_of_track", 59, -2.95, 2.95, 20, -0.25, 9.75);
0026 heta_mu1 = new TH2F(N1 + "_heta_mu1", " Sim Seg vs eta_of_track", 59, -2.95, 2.95, 20, -0.25, 9.75);
0027
0028 heta_rh = new TH2F(N1 + "_heta_rh", " All rechits vs eta ", 250, 0.0, 2.5, 20, -0.25, 9.75);
0029
0030 heta_NSeed = new TH2F(N1 + "_heta_NSeed", " eta vs N of Seed", 59, -2.95, 2.95, 40, -0.25, 19.75);
0031 heta_NSta = new TH2F(N1 + "_heta_NSta", " eta vs N of Sta ", 59, -2.95, 2.95, 20, -0.25, 9.75);
0032 heta_Sim = new TH1F(N1 + "_heta_Sim", " eta of SimTrack ", 59, -2.95, 2.95);
0033 heta_simPt = new TH2F(N1 + "_heta_simPt", " eta vs sim pt ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0034 heta_simQPt = new TH2F(N1 + "_heta_simQPt", " eta vs sim Q/Pt ", 59, -2.95, 2.95, 800, -0.2, 0.2);
0035 heta_simQPa = new TH2F(N1 + "_heta_simQPa", " eta vs sim Q/Pa ", 59, -2.95, 2.95, 800, -0.2, 0.2);
0036
0037 heta_pt = new TH2F(N1 + "_heta_pt", " eta vs seed pt ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0038 heta_pa = new TH2F(N1 + "_heta_pa", " eta vs seed p ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0039 heta_bestPt = new TH2F(N1 + "_heta_bestPt", " eta vs best seed pt ", 59, -2.95, 2.95, 5000, -2500.0, 2500.0);
0040
0041 heta_pullQP = new TH2F(N1 + "_heta_pullQP", "eta vs best pull of Q/P ", 59, -2.95, 2.95, 1000, -10.0, 10.0);
0042 heta_pullQPt = new TH2F(N1 + "_heta_pullQPt", "eta vs best pull of Q/Pt ", 59, -2.95, 2.95, 1000, -10.0, 10.0);
0043 heta_QPt = new TH2F(N1 + "_heta_QPt", "eta vs best Q/Pt ", 59, -2.95, 2.95, 1000, -0.4, 0.4);
0044 heta_errQP = new TH2F(N1 + "_heta_errQP", "eta vs best Q/P error", 59, -2.95, 2.95, 500, -0.01, 0.24);
0045 heta_errQPt = new TH2F(N1 + "_heta_errQPt", "eta vs best Q/Pt error", 59, -2.95, 2.95, 500, -0.01, 0.24);
0046 heta_resolSd = new TH2F(N1 + "_heta_resolSd", "eta vs best Q/Pt resol", 59, -2.95, 2.95, 500, -1., 1.);
0047
0048 Tpt_ptLoss = new TH2F(N1 + "_Tpt_ptLoss", " trk pt vs ptLoss ratio ", 250, 0., 250., 80, 0.3, 1.1);
0049 Mpt_ptLoss = new TH2F(N1 + "_Mpt_ptLoss", " seed pt vs ptLoss ratio ", 250, 0., 250., 80, 0.3, 1.1);
0050
0051 heta_statrk = new TH1F(N1 + "_heta_statrk", " eta from sta tracks", 59, -2.95, 2.95);
0052 heta_staQPt = new TH2F(N1 + "_heta_staQPt", " eta vs sta Q/Pt ", 59, -2.95, 2.95, 800, -0.2, 0.2);
0053 heta_staQPa = new TH2F(N1 + "_heta_staQPa", " eta vs sta Q/Pa ", 59, -2.95, 2.95, 800, -0.2, 0.2);
0054 heta_staPt = new TH2F(N1 + "_heta_staPt", " eta vs sta Pt ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0055 heta_staPa = new TH2F(N1 + "_heta_staPa", " eta vs sta P ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0056 heta_simPt1 = new TH2F(N1 + "_heta_simPt1", " eta vs sim pt ", 59, -2.95, 2.95, 2500, 0.0, 2500.0);
0057 heta_resolSA = new TH2F(N1 + "_heta_resolSA", " eta vs Q/Pt resol", 59, -2.95, 2.95, 500, -1., 1.);
0058
0059 heta_errdx = new TH2F(N1 + "_heta_errdx", "eta vs dx/dz error", 59, -2.95, 2.95, 500, 0., 0.01);
0060 heta_errdy = new TH2F(N1 + "_heta_errdy", "eta vs dy/dz error", 59, -2.95, 2.95, 500, 0., 0.15);
0061 heta_errx = new TH2F(N1 + "_heta_errx", "eta vs x error", 59, -2.95, 2.95, 1000, 0., 0.1);
0062 heta_erry = new TH2F(N1 + "_heta_erry", "eta vs y error", 59, -2.95, 2.95, 1000, 0., 1.);
0063
0064 h_dh_st = new TH2F(N1 + "_h_dh_st", " dEta from vector vs station", 19, -4.75, 4.75, 400, -1.0, 1.0);
0065 h_df_st = new TH2F(N1 + "_h_df_st", " dPhi from vector vs station", 19, -4.75, 4.75, 800, -0.1, 0.1);
0066 h_dx_st = new TH2F(N1 + "_h_dx_st", " dx from position vs station", 19, -4.75, 4.75, 400, -4.0, 4.0);
0067 h_dy_st = new TH2F(N1 + "_h_dy_st", " dy from position vs station", 19, -4.75, 4.75, 400, -4.0, 4.0);
0068
0069 h_Nsta_Nseed = new TH2F(N1 + "_h_Nsta_Nseed", " N_STA vs N_Seed ", 30, -0.25, 14.75, 30, -0.25, 14.75);
0070
0071 hPt = new TH1F(N1 + "_hPt", " Pt of Tracks ", 100, 5., 205.);
0072 hPa_Pt = new TH2F(N1 + "_hPa_Pt", "P vs Pt", 50, 5., 205., 100, 0., 1000.);
0073 hPa_nSeg = new TH2F(N1 + "_hPa_nSeg", "Pt vs total CSC segments number", 40, 5., 205., 40, -0.25, 19.75);
0074
0075 hP_Showers1 = new TH2F(N1 + "_hP_Showers1", "Pt vs # of Showering segments 1", 600, 0., 1200., 50, -0.5, 49.5);
0076 hP_Showers2 = new TH2F(N1 + "_hP_Showers2", "Pt vs # of Showering segments 2", 600, 0., 1200., 50, -0.5, 49.5);
0077 hP_Showers3 = new TH2F(N1 + "_hP_Showers3", "Pt vs # of Showering segments 3", 600, 0., 1200., 50, -0.5, 49.5);
0078
0079 hP_Showers1a =
0080 new TH2F(N1 + "_hP_Showers1a", "Pt vs # of Showering segs w/o st1 1", 600, 0., 1200., 50, -0.5, 49.5);
0081 hP_Showers2a =
0082 new TH2F(N1 + "_hP_Showers2a", "Pt vs # of Showering segs w/o st1 2", 600, 0., 1200., 50, -0.5, 49.5);
0083 hP_Showers3a =
0084 new TH2F(N1 + "_hP_Showers3a", "Pt vs # of Showering segs w/o st1 3", 600, 0., 1200., 50, -0.5, 49.5);
0085
0086 hP_Showers1b = new TH2F(N1 + "_hP_Showers1b", "Pt vs # of Showering 1", 600, 0., 1200., 5, -0.5, 4.5);
0087 hP_Showers2b = new TH2F(N1 + "_hP_Showers2b", "Pt vs # of Showering 2", 600, 0., 1200., 5, -0.5, 4.5);
0088 hP_Showers3b = new TH2F(N1 + "_hP_Showers3b", "Pt vs # of Showering 3", 600, 0., 1200., 5, -0.5, 4.5);
0089
0090 hP_avShower1 = new TH2F(N1 + "_hP_avShower1", "Pt vs ave. Showering segments 1", 600, 0., 1200., 500, -0.5, 999.5);
0091 hP_avShower2 = new TH2F(N1 + "_hP_avShower2", "Pt vs ave. Showering segments 2", 600, 0., 1200., 500, -0.5, 999.5);
0092 hP_avShower3 = new TH2F(N1 + "_hP_avShower3", "Pt vs ave. Showering segments 3", 600, 0., 1200., 500, -0.5, 999.5);
0093
0094 hP_maxR1 = new TH2F(N1 + "_hP_maxR1", "Pt vs max dR of the cone 1", 600, 0., 1200., 300, -0.5, 1.);
0095 hP_maxR2 = new TH2F(N1 + "_hP_maxR2", "Pt vs max dR of the cone 2", 600, 0., 1200., 300, -0.5, 1.);
0096 hP_maxR3 = new TH2F(N1 + "_hP_maxR3", "Pt vs max dR of the cone 3", 600, 0., 1200., 300, -0.5, 1.);
0097 }
0098
0099 H2DRecHit1(TString name_, TFile *file) {
0100 name = name_;
0101 heta_mu4 = (TH2F *)file->Get(name + "_heta_mu4");
0102 heta_mu3 = (TH2F *)file->Get(name + "_heta_mu3");
0103 heta_mu2 = (TH2F *)file->Get(name + "_heta_mu2");
0104 heta_mu1 = (TH2F *)file->Get(name + "_heta_mu1");
0105
0106 heta_rh = (TH2F *)file->Get(name + "_heta_rh");
0107
0108 heta_NSeed = (TH2F *)file->Get(name + "_heta_NSeed");
0109 heta_NSta = (TH2F *)file->Get(name + "_heta_NSta");
0110 heta_Sim = (TH1F *)file->Get(name + "_heta_Sim");
0111 heta_simPt = (TH2F *)file->Get(name + "_heta_simPt");
0112 heta_simQPt = (TH2F *)file->Get(name + "_heta_simQPt");
0113 heta_simQPa = (TH2F *)file->Get(name + "_heta_simQPa");
0114
0115 heta_pt = (TH2F *)file->Get(name + "_heta_pt");
0116 heta_pa = (TH2F *)file->Get(name + "_heta_pa");
0117 heta_bestPt = (TH2F *)file->Get(name + "_heta_bestPt");
0118
0119 heta_pullQP = (TH2F *)file->Get(name + "_heta_pullQP");
0120 heta_pullQPt = (TH2F *)file->Get(name + "_heta_pullQPt");
0121 heta_QPt = (TH2F *)file->Get(name + "_heta_QPt");
0122 heta_errQP = (TH2F *)file->Get(name + "_heta_errQP");
0123 heta_errQPt = (TH2F *)file->Get(name + "_heta_errQPt");
0124 heta_resolSd = (TH2F *)file->Get(name + "_heta_resolSd");
0125 Tpt_ptLoss = (TH2F *)file->Get(name + "_Tpt_ptLoss");
0126 Mpt_ptLoss = (TH2F *)file->Get(name + "_Mpt_ptLoss");
0127
0128 heta_statrk = (TH1F *)file->Get(name + "_heta_statrk");
0129 heta_staQPt = (TH2F *)file->Get(name + "_heta_staQPt");
0130 heta_staQPa = (TH2F *)file->Get(name + "_heta_staQPa");
0131 heta_staPt = (TH2F *)file->Get(name + "_heta_staPt");
0132 heta_staPa = (TH2F *)file->Get(name + "_heta_staPa");
0133 heta_simPt1 = (TH2F *)file->Get(name + "_heta_simPt1");
0134 heta_resolSA = (TH2F *)file->Get(name + "_heta_resolSA");
0135
0136 heta_errdx = (TH2F *)file->Get(name + "_heta_errdx");
0137 heta_errdy = (TH2F *)file->Get(name + "_heta_errdy");
0138 heta_errx = (TH2F *)file->Get(name + "_heta_errx");
0139 heta_erry = (TH2F *)file->Get(name + "_heta_erry");
0140
0141 h_dh_st = (TH2F *)file->Get(name + "_h_dh_st");
0142 h_df_st = (TH2F *)file->Get(name + "_h_df_st");
0143 h_dx_st = (TH2F *)file->Get(name + "_h_dx_st");
0144 h_dy_st = (TH2F *)file->Get(name + "_h_dy_st");
0145
0146 h_Nsta_Nseed = (TH2F *)file->Get(name + "_h_Nsta_Nseed");
0147
0148 hPt = (TH1F *)file->Get(name + "_hPt");
0149 hPa_Pt = (TH2F *)file->Get(name + "_hPa_Pt");
0150 hPa_nSeg = (TH2F *)file->Get(name + "_hPa_nSeg");
0151
0152 hP_Showers1 = (TH2F *)file->Get(name + "_hP_Showers1");
0153 hP_Showers2 = (TH2F *)file->Get(name + "_hP_Showers2");
0154 hP_Showers3 = (TH2F *)file->Get(name + "_hP_Showers3");
0155
0156 hP_Showers1a = (TH2F *)file->Get(name + "_hP_Showers1a");
0157 hP_Showers2a = (TH2F *)file->Get(name + "_hP_Showers2a");
0158 hP_Showers3a = (TH2F *)file->Get(name + "_hP_Showers3a");
0159
0160 hP_Showers1b = (TH2F *)file->Get(name + "_hP_Showers1b");
0161 hP_Showers2b = (TH2F *)file->Get(name + "_hP_Showers2b");
0162 hP_Showers3b = (TH2F *)file->Get(name + "_hP_Showers3b");
0163
0164 hP_avShower1 = (TH2F *)file->Get(name + "_hP_avShower1");
0165 hP_avShower2 = (TH2F *)file->Get(name + "_hP_avShower2");
0166 hP_avShower3 = (TH2F *)file->Get(name + "_hP_avShower3");
0167
0168 hP_maxR1 = (TH2F *)file->Get(name + "_hP_maxR1");
0169 hP_maxR2 = (TH2F *)file->Get(name + "_hP_maxR2");
0170 hP_maxR3 = (TH2F *)file->Get(name + "_hP_maxR3");
0171 }
0172
0173
0174 virtual ~H2DRecHit1() {
0175 delete heta_mu4;
0176 delete heta_mu3;
0177 delete heta_mu2;
0178 delete heta_mu1;
0179
0180 delete heta_rh;
0181
0182 delete heta_NSeed;
0183 delete heta_NSta;
0184 delete heta_Sim;
0185 delete heta_simPt;
0186 delete heta_simQPt;
0187 delete heta_simQPa;
0188
0189 delete hPt;
0190 delete hPa_Pt;
0191 delete hPa_nSeg;
0192
0193 delete hP_Showers1;
0194 delete hP_Showers2;
0195 delete hP_Showers3;
0196
0197 delete hP_Showers1a;
0198 delete hP_Showers2a;
0199 delete hP_Showers3a;
0200
0201 delete hP_Showers1b;
0202 delete hP_Showers2b;
0203 delete hP_Showers3b;
0204
0205 delete hP_avShower1;
0206 delete hP_avShower2;
0207 delete hP_avShower3;
0208
0209 delete hP_maxR1;
0210 delete hP_maxR2;
0211 delete hP_maxR3;
0212
0213 delete heta_pt;
0214 delete heta_pa;
0215 delete heta_bestPt;
0216
0217 delete heta_pullQP;
0218 delete heta_pullQPt;
0219 delete heta_QPt;
0220 delete heta_errQP;
0221 delete heta_errQPt;
0222 delete heta_resolSd;
0223 delete Tpt_ptLoss;
0224 delete Mpt_ptLoss;
0225
0226 delete heta_statrk;
0227 delete heta_staQPt;
0228 delete heta_staQPa;
0229 delete heta_staPt;
0230 delete heta_staPa;
0231 delete heta_simPt1;
0232 delete heta_resolSA;
0233
0234 delete heta_errdx;
0235 delete heta_errdy;
0236 delete heta_errx;
0237 delete heta_erry;
0238
0239 delete h_dh_st;
0240 delete h_df_st;
0241 delete h_dx_st;
0242 delete h_dy_st;
0243
0244 delete h_Nsta_Nseed;
0245 }
0246
0247 void Fill1(int layer_nu, int goodlayer_nu, int seg_sum, int sim_nu, double eta_simtrk) {
0248 heta_mu4->Fill(eta_simtrk, seg_sum);
0249 heta_mu3->Fill(eta_simtrk, goodlayer_nu);
0250 heta_mu2->Fill(eta_simtrk, layer_nu);
0251 heta_mu1->Fill(eta_simtrk, sim_nu);
0252 }
0253 void Fill1a(int rh_nu, double eta_a) { heta_rh->Fill(eta_a, rh_nu); }
0254 void Fill1b(double eta_sim, int nSeed, int nSta, double simPt, double simQPt, double simQPa) {
0255 heta_NSeed->Fill(eta_sim, nSeed);
0256 heta_NSta->Fill(eta_sim, nSta);
0257 heta_Sim->Fill(eta_sim);
0258 heta_simPt->Fill(eta_sim, simPt);
0259 heta_simQPt->Fill(eta_sim, simQPt);
0260 heta_simQPa->Fill(eta_sim, simQPa);
0261 }
0262 void Fill1c(float Pt, float Pa, int nu_seg) {
0263 hPt->Fill(Pt);
0264 hPa_Pt->Fill(Pt, Pa);
0265 hPa_nSeg->Fill(Pa, nu_seg);
0266 }
0267 void Fill1d1(double P, int showers, int showersa, int showersb, double ave, double maxR) {
0268 hP_Showers1->Fill(P, showers);
0269 hP_Showers1a->Fill(P, showersa);
0270 hP_Showers1b->Fill(P, showersb);
0271 hP_avShower1->Fill(P, ave);
0272 hP_maxR1->Fill(P, maxR);
0273 }
0274 void Fill1d2(double P, int showers, int showersa, int showersb, double ave, double maxR) {
0275 hP_Showers2->Fill(P, showers);
0276 hP_Showers2a->Fill(P, showersa);
0277 hP_Showers2b->Fill(P, showersb);
0278 hP_avShower2->Fill(P, ave);
0279 hP_maxR2->Fill(P, maxR);
0280 }
0281 void Fill1d3(double P, int showers, int showersa, int showersb, double ave, double maxR) {
0282 hP_Showers3->Fill(P, showers);
0283 hP_Showers3a->Fill(P, showersa);
0284 hP_Showers3b->Fill(P, showersb);
0285 hP_avShower3->Fill(P, ave);
0286 hP_maxR3->Fill(P, maxR);
0287 }
0288 void Fill1f(double eta, double err_dx, double err_dy, double err_x, double err_y) {
0289 heta_errdx->Fill(eta, err_dx);
0290 heta_errdy->Fill(eta, err_dy);
0291 heta_errx->Fill(eta, err_x);
0292 heta_erry->Fill(eta, err_y);
0293 }
0294 void Fill1g(double seed_mT, double seed_mA, double bestSeed_mT, float eta_seed, double ptLoss, double ptTrk) {
0295 heta_pt->Fill(eta_seed, seed_mT);
0296 heta_pa->Fill(eta_seed, seed_mA);
0297 heta_bestPt->Fill(eta_seed, bestSeed_mT);
0298 Tpt_ptLoss->Fill(ptTrk, ptLoss);
0299 Mpt_ptLoss->Fill(seed_mT, ptLoss);
0300 }
0301 void Fill1i(
0302 double pull_qbp, float eta_seed, double QPt, double pull_qbpt, double errqbp, double errqbpt, double resol_qbpt) {
0303 heta_pullQP->Fill(eta_seed, pull_qbp);
0304 heta_pullQPt->Fill(eta_seed, pull_qbpt);
0305 heta_QPt->Fill(eta_seed, QPt);
0306 heta_errQP->Fill(eta_seed, errqbp);
0307 heta_errQPt->Fill(eta_seed, errqbpt);
0308 heta_resolSd->Fill(eta_seed, resol_qbpt);
0309 }
0310 void Fill1j(
0311 double eta_sta, double sta_qbp, double sta_qbpt, double sta_pt, double sta_pa, double sim_pt, double resol_qbpt) {
0312 heta_staQPt->Fill(eta_sta, sta_qbpt);
0313 heta_staQPa->Fill(eta_sta, sta_qbp);
0314 heta_staPt->Fill(eta_sta, sta_pt);
0315 heta_simPt1->Fill(eta_sta, sim_pt);
0316 heta_staPa->Fill(eta_sta, sta_pa);
0317 heta_resolSA->Fill(eta_sta, resol_qbpt);
0318 }
0319 void Fill1e(int station, double dh_v, double df_v, double dx_p, double dy_p) {
0320 h_dh_st->Fill(station, dh_v);
0321 h_df_st->Fill(station, df_v);
0322 h_dx_st->Fill(station, dx_p);
0323 h_dy_st->Fill(station, dy_p);
0324 }
0325
0326 void Fill1o(int Nsta, int Nseed) { h_Nsta_Nseed->Fill(Nseed, Nsta); }
0327
0328 void Write() {
0329 heta_mu4->Write();
0330 heta_mu3->Write();
0331 heta_mu2->Write();
0332 heta_mu1->Write();
0333
0334 heta_rh->Write();
0335
0336 heta_NSeed->Write();
0337 heta_NSta->Write();
0338 heta_Sim->Write();
0339 heta_simPt->Write();
0340 heta_simQPt->Write();
0341 heta_simQPa->Write();
0342
0343 heta_pt->Write();
0344 heta_pa->Write();
0345 heta_bestPt->Write();
0346
0347 heta_pullQP->Write();
0348 heta_pullQPt->Write();
0349 heta_QPt->Write();
0350 heta_errQP->Write();
0351 heta_errQPt->Write();
0352 heta_resolSd->Write();
0353 Tpt_ptLoss->Write();
0354 Mpt_ptLoss->Write();
0355
0356 heta_statrk->Write();
0357 heta_staQPt->Write();
0358 heta_staQPa->Write();
0359 heta_staPt->Write();
0360 heta_staPa->Write();
0361 heta_simPt1->Write();
0362 heta_resolSA->Write();
0363
0364 heta_errdx->Write();
0365 heta_errdy->Write();
0366 heta_errx->Write();
0367 heta_erry->Write();
0368
0369 h_dh_st->Write();
0370 h_df_st->Write();
0371 h_dx_st->Write();
0372 h_dy_st->Write();
0373
0374 h_Nsta_Nseed->Write();
0375
0376 hPt->Write();
0377 hPa_Pt->Write();
0378 hPa_nSeg->Write();
0379
0380 hP_Showers1->Write();
0381 hP_Showers2->Write();
0382 hP_Showers3->Write();
0383
0384 hP_Showers1a->Write();
0385 hP_Showers2a->Write();
0386 hP_Showers3a->Write();
0387
0388 hP_Showers1b->Write();
0389 hP_Showers2b->Write();
0390 hP_Showers3b->Write();
0391
0392 hP_avShower1->Write();
0393 hP_avShower2->Write();
0394 hP_avShower3->Write();
0395
0396 hP_maxR1->Write();
0397 hP_maxR2->Write();
0398 hP_maxR3->Write();
0399 }
0400
0401 TH2F *heta_mu4;
0402 TH2F *heta_mu3;
0403 TH2F *heta_mu2;
0404 TH2F *heta_mu1;
0405
0406 TH2F *heta_rh;
0407
0408 TH2F *heta_NSeed;
0409 TH2F *heta_NSta;
0410 TH1F *heta_Sim;
0411 TH2F *heta_simPt;
0412 TH2F *heta_simQPt;
0413 TH2F *heta_simQPa;
0414
0415 TH2F *heta_pt;
0416 TH2F *heta_pa;
0417 TH2F *heta_bestPt;
0418
0419 TH2F *heta_pullQP;
0420 TH2F *heta_pullQPt;
0421 TH2F *heta_QPt;
0422 TH2F *heta_errQP;
0423 TH2F *heta_errQPt;
0424 TH2F *heta_resolSd;
0425 TH2F *Tpt_ptLoss;
0426 TH2F *Mpt_ptLoss;
0427
0428 TH1F *heta_statrk;
0429 TH2F *heta_staQPt;
0430 TH2F *heta_staQPa;
0431 TH2F *heta_staPt;
0432 TH2F *heta_staPa;
0433 TH2F *heta_simPt1;
0434 TH2F *heta_resolSA;
0435
0436 TH2F *heta_errdx;
0437 TH2F *heta_errdy;
0438 TH2F *heta_errx;
0439 TH2F *heta_erry;
0440
0441 TH2F *h_dh_st;
0442 TH2F *h_df_st;
0443 TH2F *h_dx_st;
0444 TH2F *h_dy_st;
0445
0446 TH2F *h_Nsta_Nseed;
0447
0448 TH1F *hPt;
0449 TH2F *hPa_Pt;
0450 TH2F *hPa_nSeg;
0451
0452 TH2F *hP_Showers1;
0453 TH2F *hP_Showers2;
0454 TH2F *hP_Showers3;
0455
0456 TH2F *hP_Showers1a;
0457 TH2F *hP_Showers2a;
0458 TH2F *hP_Showers3a;
0459
0460 TH2F *hP_Showers1b;
0461 TH2F *hP_Showers2b;
0462 TH2F *hP_Showers3b;
0463
0464 TH2F *hP_avShower1;
0465 TH2F *hP_avShower2;
0466 TH2F *hP_avShower3;
0467
0468 TH2F *hP_maxR1;
0469 TH2F *hP_maxR2;
0470 TH2F *hP_maxR3;
0471
0472 TString name;
0473 };
0474
0475 class H2DRecHit2 {
0476 public:
0477 H2DRecHit2(std::string name_) {
0478 TString N2 = name_.c_str();
0479 name = N2;
0480
0481 heta_nSimSegs = new TH2F(N2 + "_heta_nSimSegs", "# Segs vs eta from simtrack", 59, -2.95, 2.95, 30, -0.25, 14.75);
0482 heta_nSegs = new TH2F(N2 + "_heta_nSegs", "# Segs vs eta from simtrack", 59, -2.95, 2.95, 30, -0.25, 14.75);
0483 heta_nlayers = new TH2F(N2 + "_heta_nlayers", "# layers vs eta from simtrack", 59, -2.95, 2.95, 20, -0.25, 9.75);
0484 heta_dNSegs =
0485 new TH2F(N2 + "_heta_dNSegs", "# simseg - recseg vs eta from simtrack", 59, -2.95, 2.95, 39, -9.75, 9.75);
0486 heta_dNlayers =
0487 new TH2F(N2 + "_heta_dNlayers", "# simseg - layers vs eta from simtrack", 59, -2.95, 2.95, 39, -9.75, 9.75);
0488
0489 heta_cscrh = new TH2F(N2 + "_heta_cscrh", " csc rechits vs eta ", 59, -2.95, 2.95, 40, -0.25, 19.75);
0490 heta_dtrh = new TH2F(N2 + "_heta_dtrh", " dt rechits vs eta ", 59, -2.95, 2.95, 40, -0.25, 19.75);
0491 }
0492
0493 H2DRecHit2(TString name_, TFile *file) {
0494 name = name_;
0495
0496 heta_nSimSegs = (TH2F *)file->Get(name + "_heta_nSimSegs");
0497 heta_nSegs = (TH2F *)file->Get(name + "_heta_nSegs");
0498 heta_nlayers = (TH2F *)file->Get(name + "_heta_nlayers");
0499
0500 heta_dNSegs = (TH2F *)file->Get(name + "_heta_dNSegs");
0501 heta_dNlayers = (TH2F *)file->Get(name + "_heta_dNlayers");
0502
0503 heta_cscrh = (TH2F *)file->Get(name + "_heta_cscrh");
0504 heta_dtrh = (TH2F *)file->Get(name + "_heta_dtrh");
0505 }
0506
0507
0508 virtual ~H2DRecHit2() {
0509 delete heta_nSimSegs;
0510 delete heta_nSegs;
0511 delete heta_nlayers;
0512 delete heta_dNSegs;
0513 delete heta_dNlayers;
0514
0515 delete heta_cscrh;
0516 delete heta_dtrh;
0517 }
0518
0519 void Fill2a(double eta_simTrk, int nu_layer, int nu_seg, int nu_simseg, int dNSegs, int dNlayers) {
0520 heta_nSimSegs->Fill(eta_simTrk, nu_simseg);
0521 heta_nSegs->Fill(eta_simTrk, nu_seg);
0522 heta_nlayers->Fill(eta_simTrk, nu_layer);
0523 heta_dNSegs->Fill(eta_simTrk, dNSegs);
0524 heta_dNlayers->Fill(eta_simTrk, dNlayers);
0525 }
0526 void Fill2b(double eta_c, int cscrh_nu) { heta_cscrh->Fill(eta_c, cscrh_nu); }
0527 void Fill2c(double eta_d, int dtrh_nu) { heta_dtrh->Fill(eta_d, dtrh_nu); }
0528
0529 void Write() {
0530 heta_nSimSegs->Write();
0531 heta_nSegs->Write();
0532 heta_nlayers->Write();
0533 heta_dNSegs->Write();
0534 heta_dNlayers->Write();
0535
0536 heta_cscrh->Write();
0537 heta_dtrh->Write();
0538 }
0539
0540 TH2F *heta_nSimSegs;
0541 TH2F *heta_nSegs;
0542 TH2F *heta_nlayers;
0543
0544 TH2F *heta_dNSegs;
0545 TH2F *heta_dNlayers;
0546
0547 TH2F *heta_cscrh;
0548 TH2F *heta_dtrh;
0549
0550 TString name;
0551 };
0552
0553 class H2DRecHit3 {
0554 public:
0555 H2DRecHit3(std::string name_) {
0556 TString N3 = name_.c_str();
0557 name = N3;
0558
0559
0560 hsimEta_NSeed = new TH2F(N3 + "_hsimEta_NSeed", "No sta, sim_eta vs NSeed", 59, -2.95, 2.95, 20, -0.25, 9.75);
0561
0562 hsimEta_seedEta =
0563 new TH2F(N3 + "_hsimEta_seedEta", "No sta, sim_eta vs seed_eta", 59, -2.95, 2.95, 59, -2.95, 2.95);
0564 hsimEta_seedPt = new TH2F(N3 + "_hsimEta_seedPt", "NO sta, sim_eta vs seed pt ", 59, -2.95, 2.95, 500, 0.0, 500.0);
0565 hsimEta_errPt = new TH2F(N3 + "_hsimEta_errPt", "No sta, sim_eta vs Pt error", 59, -2.95, 2.95, 500, 0.0, 500.0);
0566 hsimEta_pullQPt =
0567 new TH2F(N3 + "_hsimEta_pullQPt", "No sta, sim_eta vs q/pt pull", 59, -2.95, 2.95, 1000, -10.0, 10.0);
0568 hsimEta_NStation = new TH2F(N3 + "_hsimEta_NStation", "No sta, sim_eta vs N Seg", 59, -2.95, 2.95, 20, -0.25, 9.75);
0569 hsimEta_NSeg1 = new TH2F(N3 + "_hsimEta_NSeg1", "No sta, sim_eta vs N good Seg", 59, -2.95, 2.95, 20, -0.25, 9.75);
0570
0571 hsimEta_dh = new TH2F(N3 + "_hsimEta_dh", "No sta, sim_eta vs h_seed - h_Seg", 59, -2.95, 2.95, 400, -1., 1.);
0572 hsimEta_df = new TH2F(N3 + "_hsimEta_df", "No sta, sim_eta vs f_seed - f_Seg", 59, -2.95, 2.95, 400, -0.1, 0.1);
0573
0574 hPhi_resid = new TH1F(N3 + "_hPhi_resid", " phi residual for fail sta case ", 500, -0.25, 0.25);
0575 hEta_resid = new TH1F(N3 + "_hEta_resid", " eta residual for fail sta case ", 500, -0.25, 0.25);
0576 }
0577
0578 H2DRecHit3(TString name_, TFile *file) {
0579 name = name_;
0580
0581 hsimEta_NSeed = (TH2F *)file->Get(name + "_hsimEta_NSeed");
0582
0583 hsimEta_seedEta = (TH2F *)file->Get(name + "_hsimEta_seedEta");
0584 hsimEta_seedPt = (TH2F *)file->Get(name + "_hsimEta_seedPt");
0585 hsimEta_errPt = (TH2F *)file->Get(name + "_hsimEta_errPt");
0586 hsimEta_pullQPt = (TH2F *)file->Get(name + "_hsimEta_pullQPt");
0587 hsimEta_NStation = (TH2F *)file->Get(name + "_hsimEta_NStation");
0588 hsimEta_NSeg1 = (TH2F *)file->Get(name + "_hsimEta_NSeg1");
0589 hsimEta_dh = (TH2F *)file->Get(name + "_hsimEta_dh");
0590 hsimEta_df = (TH2F *)file->Get(name + "_hsimEta_df");
0591 hPhi_resid = (TH1F *)file->Get(name + "_hPhi_resid");
0592 hEta_resid = (TH1F *)file->Get(name + "_hEta_resid");
0593 }
0594
0595
0596 virtual ~H2DRecHit3() {
0597 delete hsimEta_NSeed;
0598 delete hsimEta_seedEta;
0599 delete hsimEta_seedPt;
0600 delete hsimEta_errPt;
0601 delete hsimEta_pullQPt;
0602 delete hsimEta_NStation;
0603 delete hsimEta_NSeg1;
0604 delete hsimEta_dh;
0605 delete hsimEta_df;
0606 delete hPhi_resid;
0607 delete hEta_resid;
0608 }
0609
0610 void Fill3a(double simEta, double seedEta, double seed_mT, double ptErr, double pullQPt) {
0611 hsimEta_seedEta->Fill(simEta, seedEta);
0612 hsimEta_seedPt->Fill(simEta, seed_mT);
0613 hsimEta_errPt->Fill(simEta, ptErr);
0614 hsimEta_pullQPt->Fill(simEta, pullQPt);
0615 }
0616 void Fill3b(double simEta, int NStation, int NGoodSeg) {
0617 hsimEta_NStation->Fill(simEta, NStation);
0618 hsimEta_NSeg1->Fill(simEta, NGoodSeg);
0619 }
0620 void Fill3c(double phi_resid, double eta_resid) {
0621 hPhi_resid->Fill(phi_resid);
0622 hEta_resid->Fill(eta_resid);
0623 }
0624 void Fill3d(double sim_eta, double dEta, double dPhi) {
0625 hsimEta_dh->Fill(sim_eta, dEta);
0626 hsimEta_df->Fill(sim_eta, dPhi);
0627 }
0628 void Fill3f(double simEta, int nu_seed) { hsimEta_NSeed->Fill(simEta, nu_seed); }
0629
0630 void Write() {
0631 hsimEta_NSeed->Write();
0632 hsimEta_seedEta->Write();
0633 hsimEta_seedPt->Write();
0634 hsimEta_errPt->Write();
0635 hsimEta_pullQPt->Write();
0636 hsimEta_NStation->Write();
0637 hsimEta_NSeg1->Write();
0638 hsimEta_dh->Write();
0639 hsimEta_df->Write();
0640 hPhi_resid->Write();
0641 hEta_resid->Write();
0642 }
0643
0644 TH2F *hsimEta_NSeed;
0645 TH2F *hsimEta_seedEta;
0646 TH2F *hsimEta_seedPt;
0647 TH2F *hsimEta_errPt;
0648 TH2F *hsimEta_pullQPt;
0649 TH2F *hsimEta_NStation;
0650 TH2F *hsimEta_NSeg1;
0651 TH2F *hsimEta_dh;
0652 TH2F *hsimEta_df;
0653 TH1F *hPhi_resid;
0654 TH1F *hEta_resid;
0655
0656 TString name;
0657 };
0658
0659 class H2DRecHit4 {
0660 public:
0661 H2DRecHit4() {
0662
0663
0664 hSeedPhi = new TH1F("hSeedPhi", " seed phi distribution ", 160, -3.1415, 3.1415);
0665 hStaPhi = new TH1F("hStaPhi", " sta phi distribution ", 160, -3.1415, 3.1415);
0666
0667 hSeedSimEta = new TH2F("hSeedSimEta", " seed - sim eta distribution ", 59, -2.95, 2.95, 59, -2.95, 2.95);
0668 hStaSimEta = new TH2F("hStaSimEta", " sta - sim eta distribution ", 59, -2.95, 2.95, 59, -2.95, 2.95);
0669
0670 h_dh_st1 = new TH2F("h_dh_st1", " dEta from vector vs station bad pt", 19, -4.75, 4.75, 400, -1.0, 1.0);
0671 h_df_st1 = new TH2F("h_df_st1", " dPhi from vector vs station bad pt", 19, -4.75, 4.75, 800, -0.1, 0.1);
0672 h_dx_st1 = new TH2F("h_dx_st1", " dx from position vs station bad pt", 19, -4.75, 4.75, 400, -4.0, 4.0);
0673 h_dy_st1 = new TH2F("h_dy_st1", " dy from position vs station bad pt", 19, -4.75, 4.75, 400, -4.0, 4.0);
0674
0675
0676 h_pt_nHits = new TH2F("h_pt_nHits", " pt vs nHits ", 104, -0.25, 51.75, 250, 0.0, 250.0);
0677 h_pt_chi2 = new TH2F("h_pt_chi2", " pt vs chi2 ", 500, 0., 50., 250, 0.0, 250.0);
0678 }
0679
0680 H2DRecHit4(TFile *file) {
0681 hSeedPhi = (TH1F *)file->Get("hSeedPhi");
0682 hStaPhi = (TH1F *)file->Get("hStaPhi");
0683 hSeedSimEta = (TH2F *)file->Get("hSeedSimEta");
0684 hStaSimEta = (TH2F *)file->Get("hStaSimEta");
0685
0686 h_dh_st1 = (TH2F *)file->Get("h_dh_st1");
0687 h_df_st1 = (TH2F *)file->Get("h_df_st1");
0688 h_dx_st1 = (TH2F *)file->Get("h_dx_st1");
0689 h_dy_st1 = (TH2F *)file->Get("h_dy_st1");
0690
0691 h_pt_nHits = (TH2F *)file->Get("h_pt_nHits");
0692 h_pt_chi2 = (TH2F *)file->Get("h_pt_chi2");
0693 }
0694
0695
0696 virtual ~H2DRecHit4() {
0697 delete hSeedPhi;
0698 delete hStaPhi;
0699
0700 delete hSeedSimEta;
0701 delete hStaSimEta;
0702
0703 delete h_dh_st1;
0704 delete h_df_st1;
0705 delete h_dx_st1;
0706 delete h_dy_st1;
0707
0708 delete h_pt_nHits;
0709 delete h_pt_chi2;
0710 }
0711
0712 void Fill4a(double sta_phi, double sta_eta, double sim_eta) {
0713 hStaPhi->Fill(sta_phi);
0714 hStaSimEta->Fill(sta_eta, sim_eta);
0715 }
0716 void Fill4b(double seed_phi, double seed_eta, double sim_eta) {
0717 hSeedPhi->Fill(seed_phi);
0718 hSeedSimEta->Fill(seed_eta, sim_eta);
0719 }
0720 void Fill4c(int station, double dh_v, double df_v, double dx_p, double dy_p) {
0721 h_dh_st1->Fill(station, dh_v);
0722 h_df_st1->Fill(station, df_v);
0723 h_dx_st1->Fill(station, dx_p);
0724 h_dy_st1->Fill(station, dy_p);
0725 }
0726 void Fill4d(double pt, int nhits, double chi2) {
0727 h_pt_nHits->Fill(nhits, pt);
0728 h_pt_chi2->Fill(chi2, pt);
0729 }
0730
0731 void Write() {
0732 hSeedPhi->Write();
0733 hStaPhi->Write();
0734
0735 hSeedSimEta->Write();
0736 hStaSimEta->Write();
0737
0738 h_dh_st1->Write();
0739 h_df_st1->Write();
0740 h_dx_st1->Write();
0741 h_dy_st1->Write();
0742
0743 h_pt_nHits->Write();
0744 h_pt_chi2->Write();
0745 }
0746
0747 TH1F *hSeedPhi;
0748 TH1F *hStaPhi;
0749
0750 TH2F *hSeedSimEta;
0751 TH2F *hStaSimEta;
0752
0753 TH2F *h_dh_st1;
0754 TH2F *h_df_st1;
0755 TH2F *h_dx_st1;
0756 TH2F *h_dy_st1;
0757
0758 TH2F *h_pt_nHits;
0759 TH2F *h_pt_chi2;
0760 };
0761
0762 class H2DRecHit5 {
0763 public:
0764 H2DRecHit5() {
0765 h_UnRelatedSeed = new TH2F("h_UnRelatedSeed", " eta of unrelated seeds ", 59, -2.95, 2.95, 2500, 0., 2500.);
0766 h_UnRelatedSta = new TH2F("h_UnRelatedSta", " eta of unrelated sta ", 59, -2.95, 2.95, 2500, 0., 2500.);
0767 h_OrphanSeed = new TH2F("h_OrphanSeed", "sim eta vs orphan seed eta", 59, -2.95, 2.95, 59, -2.95, 2.95);
0768 h_OrphanStaEta = new TH2F("h_OrphanStaEta", "sim eta vs orphanSta eta", 59, -2.95, 2.95, 59, -2.95, 2.95);
0769 h_OrphanStaPhi = new TH2F("h_OrphanStaPhi", "sim phi vs orphanSta phi", 100, -3.1416, 3.1416, 100, -3.1416, 3.1416);
0770 h_OrphanSeedPt = new TH2F("h_OrphanSeedPt", "sim eta vs orphan seed Pt", 59, -2.95, 2.95, 2500, 0., 2500);
0771 h_OrphanStaPt = new TH2F("h_OrphanStaPt", "sim eta vs orphan sta Pt", 59, -2.95, 2.95, 2500, 0., 2500);
0772 }
0773
0774 H2DRecHit5(TFile *file) {
0775 h_UnRelatedSeed = (TH2F *)file->Get("h_UnRelatedSeed");
0776 h_UnRelatedSta = (TH2F *)file->Get("h_UnRelatedSta");
0777 h_OrphanSeed = (TH2F *)file->Get("h_OrphanSeed");
0778 h_OrphanStaEta = (TH2F *)file->Get("h_OrphanStaEta");
0779 h_OrphanStaPhi = (TH2F *)file->Get("h_OrphanStaPhi");
0780
0781 h_OrphanSeedPt = (TH2F *)file->Get("h_OrphanSeedPt");
0782 h_OrphanStaPt = (TH2F *)file->Get("h_OrphanStaPt");
0783 }
0784
0785
0786 virtual ~H2DRecHit5() {
0787 delete h_UnRelatedSeed;
0788 delete h_UnRelatedSta;
0789 delete h_OrphanSeed;
0790 delete h_OrphanStaEta;
0791 delete h_OrphanStaPhi;
0792 delete h_OrphanSeedPt;
0793 delete h_OrphanStaPt;
0794 }
0795
0796 void Fill5a(double seed_eta, double seed_pt) { h_UnRelatedSeed->Fill(seed_eta, seed_pt); }
0797 void Fill5b(double seed_eta, double sim_eta, double seed_pt) {
0798 h_OrphanSeed->Fill(sim_eta, seed_eta);
0799 h_OrphanSeedPt->Fill(sim_eta, seed_pt);
0800 }
0801 void Fill5c(double sta_eta, double sta_pt) { h_UnRelatedSta->Fill(sta_eta, sta_pt); }
0802 void Fill5d(double sta_eta, double sim_eta, double sta_pt, double sta_phi, double sim_phi) {
0803 h_OrphanStaEta->Fill(sim_eta, sta_eta);
0804 h_OrphanStaPhi->Fill(sim_phi, sta_phi);
0805 h_OrphanStaPt->Fill(sim_eta, sta_pt);
0806 }
0807
0808 void Write() {
0809 h_UnRelatedSeed->Write();
0810 h_UnRelatedSta->Write();
0811 h_OrphanSeed->Write();
0812 h_OrphanStaEta->Write();
0813 h_OrphanStaPhi->Write();
0814 h_OrphanSeedPt->Write();
0815 h_OrphanStaPt->Write();
0816 }
0817
0818 TH2F *h_UnRelatedSeed;
0819 TH2F *h_UnRelatedSta;
0820 TH2F *h_OrphanSeed;
0821 TH2F *h_OrphanStaEta;
0822 TH2F *h_OrphanStaPhi;
0823 TH2F *h_OrphanSeedPt;
0824 TH2F *h_OrphanStaPt;
0825 };
0826 #endif