File indexing completed on 2024-10-17 22:59:09
0001 #include "TTreeValidation.h"
0002 #include "Event.h"
0003 #include "RecoTracker/MkFitCore/interface/Config.h"
0004 #include "RecoTracker/MkFitCore/standalone/ConfigStandalone.h"
0005 #include "RecoTracker/MkFitCore/interface/IterationConfig.h"
0006
0007 #ifdef WITH_ROOT
0008
0009 namespace mkfit {
0010
0011 TTreeValidation::TTreeValidation(std::string fileName, const TrackerInfo* trk_info) {
0012 std::lock_guard<std::mutex> locker(glock_);
0013 gROOT->ProcessLine("#include <vector>");
0014
0015 ntotallayers_fit_ = trk_info->n_layers();
0016
0017
0018
0019
0020
0021
0022 TDirectory::TContext contextEraser;
0023 f_ = std::unique_ptr<TFile>(TFile::Open(fileName.c_str(), "recreate"));
0024
0025 if (Config::sim_val_for_cmssw || Config::sim_val) {
0026 TTreeValidation::initializeEfficiencyTree();
0027 TTreeValidation::initializeFakeRateTree();
0028 }
0029 if (Config::cmssw_val) {
0030 TTreeValidation::initializeCMSSWEfficiencyTree();
0031 TTreeValidation::initializeCMSSWFakeRateTree();
0032 }
0033 if (Config::fit_val) {
0034 for (int i = 0; i < nfvs_; ++i)
0035 fvs_[i].resize(ntotallayers_fit_);
0036 TTreeValidation::initializeFitTree();
0037 }
0038 TTreeValidation::initializeConfigTree();
0039 }
0040
0041 void TTreeValidation::initializeEfficiencyTree() {
0042
0043 efftree_ = std::make_unique<TTree>("efftree", "efftree");
0044 efftree_->SetDirectory(0);
0045
0046 efftree_->Branch("evtID", &evtID_eff_);
0047 efftree_->Branch("mcID", &mcID_eff_);
0048
0049 efftree_->Branch("nHits_mc", &nHits_mc_eff_);
0050 efftree_->Branch("nLayers_mc", &nLayers_mc_eff_);
0051 efftree_->Branch("lastlyr_mc", &lastlyr_mc_eff_);
0052
0053 efftree_->Branch("seedID_seed", &seedID_seed_eff_);
0054 efftree_->Branch("seedID_build", &seedID_build_eff_);
0055 efftree_->Branch("seedID_fit", &seedID_fit_eff_);
0056
0057 efftree_->Branch("x_mc_gen", &x_mc_gen_eff_);
0058 efftree_->Branch("y_mc_gen", &y_mc_gen_eff_);
0059 efftree_->Branch("z_mc_gen", &z_mc_gen_eff_);
0060
0061 efftree_->Branch("pt_mc_gen", &pt_mc_gen_eff_);
0062 efftree_->Branch("phi_mc_gen", &phi_mc_gen_eff_);
0063 efftree_->Branch("eta_mc_gen", &eta_mc_gen_eff_);
0064
0065 efftree_->Branch("mcmask_seed", &mcmask_seed_eff_);
0066 efftree_->Branch("mcmask_build", &mcmask_build_eff_);
0067 efftree_->Branch("mcmask_fit", &mcmask_fit_eff_);
0068
0069 efftree_->Branch("mcTSmask_seed", &mcTSmask_seed_eff_);
0070 efftree_->Branch("mcTSmask_build", &mcTSmask_build_eff_);
0071 efftree_->Branch("mcTSmask_fit", &mcTSmask_fit_eff_);
0072
0073 efftree_->Branch("xhit_seed", &xhit_seed_eff_);
0074 efftree_->Branch("xhit_build", &xhit_build_eff_);
0075 efftree_->Branch("xhit_fit", &xhit_fit_eff_);
0076
0077 efftree_->Branch("yhit_seed", &yhit_seed_eff_);
0078 efftree_->Branch("yhit_build", &yhit_build_eff_);
0079 efftree_->Branch("yhit_fit", &yhit_fit_eff_);
0080
0081 efftree_->Branch("zhit_seed", &zhit_seed_eff_);
0082 efftree_->Branch("zhit_build", &zhit_build_eff_);
0083 efftree_->Branch("zhit_fit", &zhit_fit_eff_);
0084
0085 efftree_->Branch("pt_mc_seed", &pt_mc_seed_eff_);
0086 efftree_->Branch("pt_seed", &pt_seed_eff_);
0087 efftree_->Branch("ept_seed", &ept_seed_eff_);
0088 efftree_->Branch("pt_mc_build", &pt_mc_build_eff_);
0089 efftree_->Branch("pt_build", &pt_build_eff_);
0090 efftree_->Branch("ept_build", &ept_build_eff_);
0091 efftree_->Branch("pt_mc_fit", &pt_mc_fit_eff_);
0092 efftree_->Branch("pt_fit", &pt_fit_eff_);
0093 efftree_->Branch("ept_fit", &ept_fit_eff_);
0094
0095 efftree_->Branch("phi_mc_seed", &phi_mc_seed_eff_);
0096 efftree_->Branch("phi_seed", &phi_seed_eff_);
0097 efftree_->Branch("ephi_seed", &ephi_seed_eff_);
0098 efftree_->Branch("phi_mc_build", &phi_mc_build_eff_);
0099 efftree_->Branch("phi_build", &phi_build_eff_);
0100 efftree_->Branch("ephi_build", &ephi_build_eff_);
0101 efftree_->Branch("phi_mc_fit", &phi_mc_fit_eff_);
0102 efftree_->Branch("phi_fit", &phi_fit_eff_);
0103 efftree_->Branch("ephi_fit", &ephi_fit_eff_);
0104
0105 efftree_->Branch("eta_mc_seed", &eta_mc_seed_eff_);
0106 efftree_->Branch("eta_seed", &eta_seed_eff_);
0107 efftree_->Branch("eeta_seed", &eeta_seed_eff_);
0108 efftree_->Branch("eta_mc_build", &eta_mc_build_eff_);
0109 efftree_->Branch("eta_build", &eta_build_eff_);
0110 efftree_->Branch("eeta_build", &eeta_build_eff_);
0111 efftree_->Branch("eta_mc_fit", &eta_mc_fit_eff_);
0112 efftree_->Branch("eta_fit", &eta_fit_eff_);
0113 efftree_->Branch("eeta_fit", &eeta_fit_eff_);
0114
0115 efftree_->Branch("nHits_seed", &nHits_seed_eff_);
0116 efftree_->Branch("nHits_build", &nHits_build_eff_);
0117 efftree_->Branch("nHits_fit", &nHits_fit_eff_);
0118
0119 efftree_->Branch("nLayers_seed", &nLayers_seed_eff_);
0120 efftree_->Branch("nLayers_build", &nLayers_build_eff_);
0121 efftree_->Branch("nLayers_fit", &nLayers_fit_eff_);
0122
0123 efftree_->Branch("nHitsMatched_seed", &nHitsMatched_seed_eff_);
0124 efftree_->Branch("nHitsMatched_build", &nHitsMatched_build_eff_);
0125 efftree_->Branch("nHitsMatched_fit", &nHitsMatched_fit_eff_);
0126
0127 efftree_->Branch("fracHitsMatched_seed", &fracHitsMatched_seed_eff_);
0128 efftree_->Branch("fracHitsMatched_build", &fracHitsMatched_build_eff_);
0129 efftree_->Branch("fracHitsMatched_fit", &fracHitsMatched_fit_eff_);
0130
0131 efftree_->Branch("lastlyr_seed", &lastlyr_seed_eff_);
0132 efftree_->Branch("lastlyr_build", &lastlyr_build_eff_);
0133 efftree_->Branch("lastlyr_fit", &lastlyr_fit_eff_);
0134
0135 efftree_->Branch("dphi_seed", &dphi_seed_eff_);
0136 efftree_->Branch("dphi_build", &dphi_build_eff_);
0137 efftree_->Branch("dphi_fit", &dphi_fit_eff_);
0138
0139 efftree_->Branch("hitchi2_seed", &hitchi2_seed_eff_);
0140 efftree_->Branch("hitchi2_build", &hitchi2_build_eff_);
0141 efftree_->Branch("hitchi2_fit", &hitchi2_fit_eff_);
0142
0143 efftree_->Branch("score_seed", &score_seed_eff_);
0144 efftree_->Branch("score_build", &score_build_eff_);
0145 efftree_->Branch("score_fit", &score_fit_eff_);
0146
0147 efftree_->Branch("helixchi2_seed", &helixchi2_seed_eff_);
0148 efftree_->Branch("helixchi2_build", &helixchi2_build_eff_);
0149 efftree_->Branch("helixchi2_fit", &helixchi2_fit_eff_);
0150
0151 efftree_->Branch("duplmask_seed", &duplmask_seed_eff_);
0152 efftree_->Branch("duplmask_build", &duplmask_build_eff_);
0153 efftree_->Branch("duplmask_fit", &duplmask_fit_eff_);
0154
0155 efftree_->Branch("nTkMatches_seed", &nTkMatches_seed_eff_);
0156 efftree_->Branch("nTkMatches_build", &nTkMatches_build_eff_);
0157 efftree_->Branch("nTkMatches_fit", &nTkMatches_fit_eff_);
0158
0159 efftree_->Branch("itermask_seed", &itermask_seed_eff_);
0160 efftree_->Branch("itermask_build", &itermask_build_eff_);
0161 efftree_->Branch("itermask_fit", &itermask_fit_eff_);
0162 efftree_->Branch("iterduplmask_seed", &iterduplmask_seed_eff_);
0163 efftree_->Branch("iterduplmask_build", &iterduplmask_build_eff_);
0164 efftree_->Branch("iterduplmask_fit", &iterduplmask_fit_eff_);
0165 efftree_->Branch("algo_seed", &algo_seed_eff_);
0166
0167 if (Config::keepHitInfo) {
0168 efftree_->Branch("hitlyrs_mc", &hitlyrs_mc_eff_);
0169 efftree_->Branch("hitlyrs_seed", &hitlyrs_seed_eff_);
0170 efftree_->Branch("hitlyrs_build", &hitlyrs_build_eff_);
0171 efftree_->Branch("hitlyrs_fit", &hitlyrs_fit_eff_);
0172
0173 efftree_->Branch("hitidxs_mc", &hitidxs_mc_eff_);
0174 efftree_->Branch("hitidxs_seed", &hitidxs_seed_eff_);
0175 efftree_->Branch("hitidxs_build", &hitidxs_build_eff_);
0176 efftree_->Branch("hitidxs_fit", &hitidxs_fit_eff_);
0177
0178 efftree_->Branch("hitmcTkIDs_mc", &hitmcTkIDs_mc_eff_);
0179 efftree_->Branch("hitmcTkIDs_seed", &hitmcTkIDs_seed_eff_);
0180 efftree_->Branch("hitmcTkIDs_build", &hitmcTkIDs_build_eff_);
0181 efftree_->Branch("hitmcTkIDs_fit", &hitmcTkIDs_fit_eff_);
0182
0183 efftree_->Branch("hitxs_mc", &hitxs_mc_eff_);
0184 efftree_->Branch("hitxs_seed", &hitxs_seed_eff_);
0185 efftree_->Branch("hitxs_build", &hitxs_build_eff_);
0186 efftree_->Branch("hitxs_fit", &hitxs_fit_eff_);
0187
0188 efftree_->Branch("hitys_mc", &hitys_mc_eff_);
0189 efftree_->Branch("hitys_seed", &hitys_seed_eff_);
0190 efftree_->Branch("hitys_build", &hitys_build_eff_);
0191 efftree_->Branch("hitys_fit", &hitys_fit_eff_);
0192
0193 efftree_->Branch("hitzs_mc", &hitzs_mc_eff_);
0194 efftree_->Branch("hitzs_seed", &hitzs_seed_eff_);
0195 efftree_->Branch("hitzs_build", &hitzs_build_eff_);
0196 efftree_->Branch("hitzs_fit", &hitzs_fit_eff_);
0197 }
0198 }
0199
0200 void TTreeValidation::initializeFakeRateTree() {
0201
0202 frtree_ = std::make_unique<TTree>("frtree", "frtree");
0203 frtree_->SetDirectory(0);
0204
0205 frtree_->Branch("evtID", &evtID_FR_);
0206 frtree_->Branch("seedID", &seedID_FR_);
0207
0208 frtree_->Branch("seedmask_seed", &seedmask_seed_FR_);
0209 frtree_->Branch("seedmask_build", &seedmask_build_FR_);
0210 frtree_->Branch("seedmask_fit", &seedmask_fit_FR_);
0211
0212 frtree_->Branch("xhit_seed", &xhit_seed_FR_);
0213 frtree_->Branch("xhit_build", &xhit_build_FR_);
0214 frtree_->Branch("xhit_fit", &xhit_fit_FR_);
0215
0216 frtree_->Branch("yhit_seed", &yhit_seed_FR_);
0217 frtree_->Branch("yhit_build", &yhit_build_FR_);
0218 frtree_->Branch("yhit_fit", &yhit_fit_FR_);
0219
0220 frtree_->Branch("zhit_seed", &zhit_seed_FR_);
0221 frtree_->Branch("zhit_build", &zhit_build_FR_);
0222 frtree_->Branch("zhit_fit", &zhit_fit_FR_);
0223
0224 frtree_->Branch("pt_seed", &pt_seed_FR_);
0225 frtree_->Branch("ept_seed", &ept_seed_FR_);
0226 frtree_->Branch("pt_build", &pt_build_FR_);
0227 frtree_->Branch("ept_build", &ept_build_FR_);
0228 frtree_->Branch("pt_fit", &pt_fit_FR_);
0229 frtree_->Branch("ept_fit", &ept_fit_FR_);
0230
0231 frtree_->Branch("phi_seed", &phi_seed_FR_);
0232 frtree_->Branch("ephi_seed", &ephi_seed_FR_);
0233 frtree_->Branch("phi_build", &phi_build_FR_);
0234 frtree_->Branch("ephi_build", &ephi_build_FR_);
0235 frtree_->Branch("phi_fit", &phi_fit_FR_);
0236 frtree_->Branch("ephi_fit", &ephi_fit_FR_);
0237
0238 frtree_->Branch("eta_seed", &eta_seed_FR_);
0239 frtree_->Branch("eeta_seed", &eeta_seed_FR_);
0240 frtree_->Branch("eta_build", &eta_build_FR_);
0241 frtree_->Branch("eeta_build", &eeta_build_FR_);
0242 frtree_->Branch("eta_fit", &eta_fit_FR_);
0243 frtree_->Branch("eeta_fit", &eeta_fit_FR_);
0244
0245 frtree_->Branch("nHits_seed", &nHits_seed_FR_);
0246 frtree_->Branch("nHits_build", &nHits_build_FR_);
0247 frtree_->Branch("nHits_fit", &nHits_fit_FR_);
0248
0249 frtree_->Branch("nLayers_seed", &nLayers_seed_FR_);
0250 frtree_->Branch("nLayers_build", &nLayers_build_FR_);
0251 frtree_->Branch("nLayers_fit", &nLayers_fit_FR_);
0252
0253 frtree_->Branch("nHitsMatched_seed", &nHitsMatched_seed_FR_);
0254 frtree_->Branch("nHitsMatched_build", &nHitsMatched_build_FR_);
0255 frtree_->Branch("nHitsMatched_fit", &nHitsMatched_fit_FR_);
0256
0257 frtree_->Branch("fracHitsMatched_seed", &fracHitsMatched_seed_FR_);
0258 frtree_->Branch("fracHitsMatched_build", &fracHitsMatched_build_FR_);
0259 frtree_->Branch("fracHitsMatched_fit", &fracHitsMatched_fit_FR_);
0260
0261 frtree_->Branch("lastlyr_seed", &lastlyr_seed_FR_);
0262 frtree_->Branch("lastlyr_build", &lastlyr_build_FR_);
0263 frtree_->Branch("lastlyr_fit", &lastlyr_fit_FR_);
0264
0265 frtree_->Branch("dphi_seed", &dphi_seed_FR_);
0266 frtree_->Branch("dphi_build", &dphi_build_FR_);
0267 frtree_->Branch("dphi_fit", &dphi_fit_FR_);
0268
0269 frtree_->Branch("hitchi2_seed", &hitchi2_seed_FR_);
0270 frtree_->Branch("hitchi2_build", &hitchi2_build_FR_);
0271 frtree_->Branch("hitchi2_fit", &hitchi2_fit_FR_);
0272
0273 frtree_->Branch("score_seed", &score_seed_FR_);
0274 frtree_->Branch("score_build", &score_build_FR_);
0275 frtree_->Branch("score_fit", &score_fit_FR_);
0276
0277
0278 frtree_->Branch("mcID_seed", &mcID_seed_FR_);
0279 frtree_->Branch("mcID_build", &mcID_build_FR_);
0280 frtree_->Branch("mcID_fit", &mcID_fit_FR_);
0281
0282 frtree_->Branch("mcmask_seed", &mcmask_seed_FR_);
0283 frtree_->Branch("mcmask_build", &mcmask_build_FR_);
0284 frtree_->Branch("mcmask_fit", &mcmask_fit_FR_);
0285
0286 frtree_->Branch("mcTSmask_seed", &mcTSmask_seed_FR_);
0287 frtree_->Branch("mcTSmask_build", &mcTSmask_build_FR_);
0288 frtree_->Branch("mcTSmask_fit", &mcTSmask_fit_FR_);
0289
0290 frtree_->Branch("pt_mc_seed", &pt_mc_seed_FR_);
0291 frtree_->Branch("pt_mc_build", &pt_mc_build_FR_);
0292 frtree_->Branch("pt_mc_fit", &pt_mc_fit_FR_);
0293
0294 frtree_->Branch("phi_mc_seed", &phi_mc_seed_FR_);
0295 frtree_->Branch("phi_mc_build", &phi_mc_build_FR_);
0296 frtree_->Branch("phi_mc_fit", &phi_mc_fit_FR_);
0297
0298 frtree_->Branch("eta_mc_seed", &eta_mc_seed_FR_);
0299 frtree_->Branch("eta_mc_build", &eta_mc_build_FR_);
0300 frtree_->Branch("eta_mc_fit", &eta_mc_fit_FR_);
0301
0302 frtree_->Branch("nHits_mc_seed", &nHits_mc_seed_FR_);
0303 frtree_->Branch("nHits_mc_build", &nHits_mc_build_FR_);
0304 frtree_->Branch("nHits_mc_fit", &nHits_mc_fit_FR_);
0305
0306 frtree_->Branch("nLayers_mc_seed", &nLayers_mc_seed_FR_);
0307 frtree_->Branch("nLayers_mc_build", &nLayers_mc_build_FR_);
0308 frtree_->Branch("nLayers_mc_fit", &nLayers_mc_fit_FR_);
0309
0310 frtree_->Branch("lastlyr_mc_seed", &lastlyr_mc_seed_FR_);
0311 frtree_->Branch("lastlyr_mc_build", &lastlyr_mc_build_FR_);
0312 frtree_->Branch("lastlyr_mc_fit", &lastlyr_mc_fit_FR_);
0313
0314 frtree_->Branch("helixchi2_seed", &helixchi2_seed_FR_);
0315 frtree_->Branch("helixchi2_build", &helixchi2_build_FR_);
0316 frtree_->Branch("helixchi2_fit", &helixchi2_fit_FR_);
0317
0318 frtree_->Branch("duplmask_seed", &duplmask_seed_FR_);
0319 frtree_->Branch("duplmask_build", &duplmask_build_FR_);
0320 frtree_->Branch("duplmask_fit", &duplmask_fit_FR_);
0321
0322 frtree_->Branch("iTkMatches_seed", &iTkMatches_seed_FR_);
0323 frtree_->Branch("iTkMatches_build", &iTkMatches_build_FR_);
0324 frtree_->Branch("iTkMatches_fit", &iTkMatches_fit_FR_);
0325
0326 frtree_->Branch("algorithm", &algorithm_FR_);
0327
0328 if (Config::keepHitInfo) {
0329 frtree_->Branch("hitlyrs_seed", &hitlyrs_seed_FR_);
0330 frtree_->Branch("hitlyrs_mc_seed", &hitlyrs_mc_seed_FR_);
0331 frtree_->Branch("hitlyrs_build", &hitlyrs_build_FR_);
0332 frtree_->Branch("hitlyrs_mc_build", &hitlyrs_mc_build_FR_);
0333 frtree_->Branch("hitlyrs_fit", &hitlyrs_fit_FR_);
0334 frtree_->Branch("hitlyrs_mc_fit", &hitlyrs_mc_fit_FR_);
0335
0336 frtree_->Branch("hitidxs_seed", &hitidxs_seed_FR_);
0337 frtree_->Branch("hitidxs_mc_seed", &hitidxs_mc_seed_FR_);
0338 frtree_->Branch("hitidxs_build", &hitidxs_build_FR_);
0339 frtree_->Branch("hitidxs_mc_build", &hitidxs_mc_build_FR_);
0340 frtree_->Branch("hitidxs_fit", &hitidxs_fit_FR_);
0341 frtree_->Branch("hitidxs_mc_fit", &hitidxs_mc_fit_FR_);
0342
0343 frtree_->Branch("hitmcTkIDs_seed", &hitmcTkIDs_seed_FR_);
0344 frtree_->Branch("hitmcTkIDs_mc_seed", &hitmcTkIDs_mc_seed_FR_);
0345 frtree_->Branch("hitmcTkIDs_build", &hitmcTkIDs_build_FR_);
0346 frtree_->Branch("hitmcTkIDs_mc_build", &hitmcTkIDs_mc_build_FR_);
0347 frtree_->Branch("hitmcTkIDs_fit", &hitmcTkIDs_fit_FR_);
0348 frtree_->Branch("hitmcTkIDs_mc_fit", &hitmcTkIDs_mc_fit_FR_);
0349
0350 frtree_->Branch("hitxs_seed", &hitxs_seed_FR_);
0351 frtree_->Branch("hitxs_mc_seed", &hitxs_mc_seed_FR_);
0352 frtree_->Branch("hitxs_build", &hitxs_build_FR_);
0353 frtree_->Branch("hitxs_mc_build", &hitxs_mc_build_FR_);
0354 frtree_->Branch("hitxs_fit", &hitxs_fit_FR_);
0355 frtree_->Branch("hitxs_mc_fit", &hitxs_mc_fit_FR_);
0356
0357 frtree_->Branch("hitys_seed", &hitys_seed_FR_);
0358 frtree_->Branch("hitys_mc_seed", &hitys_mc_seed_FR_);
0359 frtree_->Branch("hitys_build", &hitys_build_FR_);
0360 frtree_->Branch("hitys_mc_build", &hitys_mc_build_FR_);
0361 frtree_->Branch("hitys_fit", &hitys_fit_FR_);
0362 frtree_->Branch("hitys_mc_fit", &hitys_mc_fit_FR_);
0363
0364 frtree_->Branch("hitzs_seed", &hitzs_seed_FR_);
0365 frtree_->Branch("hitzs_mc_seed", &hitzs_mc_seed_FR_);
0366 frtree_->Branch("hitzs_build", &hitzs_build_FR_);
0367 frtree_->Branch("hitzs_mc_build", &hitzs_mc_build_FR_);
0368 frtree_->Branch("hitzs_fit", &hitzs_fit_FR_);
0369 frtree_->Branch("hitzs_mc_fit", &hitzs_mc_fit_FR_);
0370 }
0371 }
0372
0373 void TTreeValidation::initializeConfigTree() {
0374
0375 configtree_ = std::make_unique<TTree>("configtree", "configtree");
0376 configtree_->SetDirectory(0);
0377
0378 configtree_->Branch("Ntracks", &Ntracks_);
0379 configtree_->Branch("Nevents", &Nevents_);
0380
0381 configtree_->Branch("nLayers", &nLayers_);
0382
0383 configtree_->Branch("nlayers_per_seed", &nlayers_per_seed_);
0384 configtree_->Branch("maxCand", &maxCand_);
0385 configtree_->Branch("chi2Cut_min", &chi2Cut_min_);
0386 configtree_->Branch("nSigma", &nSigma_);
0387 configtree_->Branch("minDPhi", &minDPhi_);
0388 configtree_->Branch("maxDPhi", &maxDPhi_);
0389 configtree_->Branch("minDEta", &minDEta_);
0390 configtree_->Branch("maxDEta", &maxDEta_);
0391
0392 configtree_->Branch("beamspotX", &beamspotX_);
0393 configtree_->Branch("beamspotY", &beamspotY_);
0394 configtree_->Branch("beamspotZ", &beamspotZ_);
0395
0396 configtree_->Branch("minSimPt", &minSimPt_);
0397 configtree_->Branch("maxSimPt", &maxSimPt_);
0398
0399 configtree_->Branch("hitposerrXY", &hitposerrXY_);
0400 configtree_->Branch("hitposerrZ", &hitposerrZ_);
0401 configtree_->Branch("hitposerrR", &hitposerrR_);
0402
0403 configtree_->Branch("varXY", &varXY_);
0404 configtree_->Branch("varZ", &varZ_);
0405
0406 configtree_->Branch("ptinverr049", &ptinverr049_);
0407 configtree_->Branch("phierr049", &phierr049_);
0408 configtree_->Branch("thetaerr049", &thetaerr049_);
0409 configtree_->Branch("ptinverr012", &ptinverr012_);
0410 configtree_->Branch("phierr012", &phierr012_);
0411 configtree_->Branch("thetaerr012", &thetaerr012_);
0412 }
0413
0414 void TTreeValidation::initializeCMSSWEfficiencyTree() {
0415
0416 cmsswefftree_ = std::make_unique<TTree>("cmsswefftree", "cmsswefftree");
0417 cmsswefftree_->SetDirectory(0);
0418
0419 cmsswefftree_->Branch("evtID", &evtID_ceff_);
0420 cmsswefftree_->Branch("cmsswID", &cmsswID_ceff_);
0421 cmsswefftree_->Branch("seedID_cmssw", &seedID_cmssw_ceff_);
0422
0423
0424 cmsswefftree_->Branch("x_cmssw", &x_cmssw_ceff_);
0425 cmsswefftree_->Branch("y_cmssw", &y_cmssw_ceff_);
0426 cmsswefftree_->Branch("z_cmssw", &z_cmssw_ceff_);
0427
0428 cmsswefftree_->Branch("pt_cmssw", &pt_cmssw_ceff_);
0429 cmsswefftree_->Branch("phi_cmssw", &phi_cmssw_ceff_);
0430 cmsswefftree_->Branch("eta_cmssw", &eta_cmssw_ceff_);
0431
0432 cmsswefftree_->Branch("nHits_cmssw", &nHits_cmssw_ceff_);
0433 cmsswefftree_->Branch("nLayers_cmssw", &nLayers_cmssw_ceff_);
0434 cmsswefftree_->Branch("lastlyr_cmssw", &lastlyr_cmssw_ceff_);
0435
0436
0437 cmsswefftree_->Branch("cmsswmask_build", &cmsswmask_build_ceff_);
0438 cmsswefftree_->Branch("seedID_build", &seedID_build_ceff_);
0439 cmsswefftree_->Branch("mcTrackID_build", &mcTrackID_build_ceff_);
0440
0441 cmsswefftree_->Branch("pt_build", &pt_build_ceff_);
0442 cmsswefftree_->Branch("ept_build", &ept_build_ceff_);
0443 cmsswefftree_->Branch("phi_build", &phi_build_ceff_);
0444 cmsswefftree_->Branch("ephi_build", &ephi_build_ceff_);
0445 cmsswefftree_->Branch("eta_build", &eta_build_ceff_);
0446 cmsswefftree_->Branch("eeta_build", &eeta_build_ceff_);
0447
0448 cmsswefftree_->Branch("x_mc_build", &x_mc_build_ceff_);
0449 cmsswefftree_->Branch("y_mc_build", &y_mc_build_ceff_);
0450 cmsswefftree_->Branch("z_mc_build", &z_mc_build_ceff_);
0451 cmsswefftree_->Branch("pt_mc_build", &pt_mc_build_ceff_);
0452 cmsswefftree_->Branch("phi_mc_build", &phi_mc_build_ceff_);
0453 cmsswefftree_->Branch("eta_mc_build", &eta_mc_build_ceff_);
0454
0455 cmsswefftree_->Branch("nHits_build", &nHits_build_ceff_);
0456 cmsswefftree_->Branch("nLayers_build", &nLayers_build_ceff_);
0457 cmsswefftree_->Branch("nHitsMatched_build", &nHitsMatched_build_ceff_);
0458 cmsswefftree_->Branch("fracHitsMatched_build", &fracHitsMatched_build_ceff_);
0459 cmsswefftree_->Branch("lastlyr_build", &lastlyr_build_ceff_);
0460
0461 cmsswefftree_->Branch("xhit_build", &xhit_build_ceff_);
0462 cmsswefftree_->Branch("yhit_build", &yhit_build_ceff_);
0463 cmsswefftree_->Branch("zhit_build", &zhit_build_ceff_);
0464
0465 cmsswefftree_->Branch("hitchi2_build", &hitchi2_build_ceff_);
0466 cmsswefftree_->Branch("helixchi2_build", &helixchi2_build_ceff_);
0467 cmsswefftree_->Branch("score_build", &score_build_ceff_);
0468 cmsswefftree_->Branch("dphi_build", &dphi_build_ceff_);
0469
0470 cmsswefftree_->Branch("duplmask_build", &duplmask_build_ceff_);
0471 cmsswefftree_->Branch("nTkMatches_build", &nTkMatches_build_ceff_);
0472
0473 cmsswefftree_->Branch("itermask_build", &itermask_build_ceff_);
0474 cmsswefftree_->Branch("iterduplmask_build", &iterduplmask_build_ceff_);
0475
0476
0477 cmsswefftree_->Branch("cmsswmask_fit", &cmsswmask_fit_ceff_);
0478 cmsswefftree_->Branch("seedID_fit", &seedID_fit_ceff_);
0479 cmsswefftree_->Branch("mcTrackID_fit", &mcTrackID_fit_ceff_);
0480
0481 cmsswefftree_->Branch("pt_fit", &pt_fit_ceff_);
0482 cmsswefftree_->Branch("ept_fit", &ept_fit_ceff_);
0483 cmsswefftree_->Branch("phi_fit", &phi_fit_ceff_);
0484 cmsswefftree_->Branch("ephi_fit", &ephi_fit_ceff_);
0485 cmsswefftree_->Branch("eta_fit", &eta_fit_ceff_);
0486 cmsswefftree_->Branch("eeta_fit", &eeta_fit_ceff_);
0487
0488 cmsswefftree_->Branch("x_mc_fit", &x_mc_fit_ceff_);
0489 cmsswefftree_->Branch("y_mc_fit", &y_mc_fit_ceff_);
0490 cmsswefftree_->Branch("z_mc_fit", &z_mc_fit_ceff_);
0491 cmsswefftree_->Branch("pt_mc_fit", &pt_mc_fit_ceff_);
0492 cmsswefftree_->Branch("phi_mc_fit", &phi_mc_fit_ceff_);
0493 cmsswefftree_->Branch("eta_mc_fit", &eta_mc_fit_ceff_);
0494
0495 cmsswefftree_->Branch("nHits_fit", &nHits_fit_ceff_);
0496 cmsswefftree_->Branch("nLayers_fit", &nLayers_fit_ceff_);
0497 cmsswefftree_->Branch("nHitsMatched_fit", &nHitsMatched_fit_ceff_);
0498 cmsswefftree_->Branch("fracHitsMatched_fit", &fracHitsMatched_fit_ceff_);
0499 cmsswefftree_->Branch("lastlyr_fit", &lastlyr_fit_ceff_);
0500
0501 cmsswefftree_->Branch("xhit_fit", &xhit_fit_ceff_);
0502 cmsswefftree_->Branch("yhit_fit", &yhit_fit_ceff_);
0503 cmsswefftree_->Branch("zhit_fit", &zhit_fit_ceff_);
0504
0505 cmsswefftree_->Branch("hitchi2_fit", &hitchi2_fit_ceff_);
0506 cmsswefftree_->Branch("helixchi2_fit", &helixchi2_fit_ceff_);
0507 cmsswefftree_->Branch("score_fit", &score_fit_ceff_);
0508 cmsswefftree_->Branch("dphi_fit", &dphi_fit_ceff_);
0509
0510 cmsswefftree_->Branch("duplmask_fit", &duplmask_fit_ceff_);
0511 cmsswefftree_->Branch("nTkMatches_fit", &nTkMatches_fit_ceff_);
0512
0513 cmsswefftree_->Branch("itermask_fit", &itermask_fit_ceff_);
0514 cmsswefftree_->Branch("iterduplmask_fit", &iterduplmask_fit_ceff_);
0515
0516 cmsswefftree_->Branch("algo_seed", &algo_seed_ceff_);
0517
0518 if (Config::keepHitInfo) {
0519 cmsswefftree_->Branch("hitlyrs_cmssw", &hitlyrs_cmssw_ceff_);
0520 cmsswefftree_->Branch("hitlyrs_build", &hitlyrs_build_ceff_);
0521 cmsswefftree_->Branch("hitlyrs_mc_build", &hitlyrs_mc_build_ceff_);
0522 cmsswefftree_->Branch("hitlyrs_fit", &hitlyrs_fit_ceff_);
0523 cmsswefftree_->Branch("hitlyrs_mc_fit", &hitlyrs_mc_fit_ceff_);
0524
0525 cmsswefftree_->Branch("hitidxs_cmssw", &hitidxs_cmssw_ceff_);
0526 cmsswefftree_->Branch("hitidxs_build", &hitidxs_build_ceff_);
0527 cmsswefftree_->Branch("hitidxs_mc_build", &hitidxs_mc_build_ceff_);
0528 cmsswefftree_->Branch("hitidxs_fit", &hitidxs_fit_ceff_);
0529 cmsswefftree_->Branch("hitidxs_mc_fit", &hitidxs_mc_fit_ceff_);
0530 }
0531 }
0532
0533 void TTreeValidation::initializeCMSSWFakeRateTree() {
0534
0535 cmsswfrtree_ = std::make_unique<TTree>("cmsswfrtree", "cmsswfrtree");
0536 cmsswfrtree_->SetDirectory(0);
0537
0538 cmsswfrtree_->Branch("evtID", &evtID_cFR_);
0539 cmsswfrtree_->Branch("seedID", &seedID_cFR_);
0540 cmsswfrtree_->Branch("mcTrackID", &mcTrackID_cFR_);
0541
0542
0543 cmsswfrtree_->Branch("x_mc", &x_mc_cFR_);
0544 cmsswfrtree_->Branch("y_mc", &y_mc_cFR_);
0545 cmsswfrtree_->Branch("z_mc", &z_mc_cFR_);
0546 cmsswfrtree_->Branch("pt_mc", &pt_mc_cFR_);
0547 cmsswfrtree_->Branch("phi_mc", &phi_mc_cFR_);
0548 cmsswfrtree_->Branch("eta_mc", &eta_mc_cFR_);
0549
0550
0551 cmsswfrtree_->Branch("cmsswID_build", &cmsswID_build_cFR_);
0552 cmsswfrtree_->Branch("cmsswmask_build", &cmsswmask_build_cFR_);
0553
0554 cmsswfrtree_->Branch("pt_build", &pt_build_cFR_);
0555 cmsswfrtree_->Branch("ept_build", &ept_build_cFR_);
0556 cmsswfrtree_->Branch("phi_build", &phi_build_cFR_);
0557 cmsswfrtree_->Branch("ephi_build", &ephi_build_cFR_);
0558 cmsswfrtree_->Branch("eta_build", &eta_build_cFR_);
0559 cmsswfrtree_->Branch("eeta_build", &eeta_build_cFR_);
0560
0561 cmsswfrtree_->Branch("nHits_build", &nHits_build_cFR_);
0562 cmsswfrtree_->Branch("nLayers_build", &nLayers_build_cFR_);
0563 cmsswfrtree_->Branch("nHitsMatched_build", &nHitsMatched_build_cFR_);
0564 cmsswfrtree_->Branch("fracHitsMatched_build", &fracHitsMatched_build_cFR_);
0565 cmsswfrtree_->Branch("lastlyr_build", &lastlyr_build_cFR_);
0566
0567 cmsswfrtree_->Branch("xhit_build", &xhit_build_cFR_);
0568 cmsswfrtree_->Branch("yhit_build", &yhit_build_cFR_);
0569 cmsswfrtree_->Branch("zhit_build", &zhit_build_cFR_);
0570
0571 cmsswfrtree_->Branch("hitchi2_build", &hitchi2_build_cFR_);
0572 cmsswfrtree_->Branch("helixchi2_build", &helixchi2_build_cFR_);
0573 cmsswfrtree_->Branch("score_build", &score_build_cFR_);
0574 cmsswfrtree_->Branch("dphi_build", &dphi_build_cFR_);
0575
0576 cmsswfrtree_->Branch("duplmask_build", &duplmask_build_cFR_);
0577 cmsswfrtree_->Branch("iTkMatches_build", &iTkMatches_build_cFR_);
0578
0579 cmsswfrtree_->Branch("seedID_cmssw_build", &seedID_cmssw_build_cFR_);
0580
0581 cmsswfrtree_->Branch("x_cmssw_build", &x_cmssw_build_cFR_);
0582 cmsswfrtree_->Branch("y_cmssw_build", &y_cmssw_build_cFR_);
0583 cmsswfrtree_->Branch("z_cmssw_build", &z_cmssw_build_cFR_);
0584
0585 cmsswfrtree_->Branch("pt_cmssw_build", &pt_cmssw_build_cFR_);
0586 cmsswfrtree_->Branch("phi_cmssw_build", &phi_cmssw_build_cFR_);
0587 cmsswfrtree_->Branch("eta_cmssw_build", &eta_cmssw_build_cFR_);
0588
0589 cmsswfrtree_->Branch("nHits_cmssw_build", &nHits_cmssw_build_cFR_);
0590 cmsswfrtree_->Branch("nLayers_cmssw_build", &nLayers_cmssw_build_cFR_);
0591 cmsswfrtree_->Branch("lastlyr_cmssw_build", &lastlyr_cmssw_build_cFR_);
0592
0593
0594 cmsswfrtree_->Branch("cmsswID_fit", &cmsswID_fit_cFR_);
0595 cmsswfrtree_->Branch("cmsswmask_fit", &cmsswmask_fit_cFR_);
0596
0597 cmsswfrtree_->Branch("pt_fit", &pt_fit_cFR_);
0598 cmsswfrtree_->Branch("ept_fit", &ept_fit_cFR_);
0599 cmsswfrtree_->Branch("phi_fit", &phi_fit_cFR_);
0600 cmsswfrtree_->Branch("ephi_fit", &ephi_fit_cFR_);
0601 cmsswfrtree_->Branch("eta_fit", &eta_fit_cFR_);
0602 cmsswfrtree_->Branch("eeta_fit", &eeta_fit_cFR_);
0603
0604 cmsswfrtree_->Branch("nHits_fit", &nHits_fit_cFR_);
0605 cmsswfrtree_->Branch("nLayers_fit", &nLayers_fit_cFR_);
0606 cmsswfrtree_->Branch("nHitsMatched_fit", &nHitsMatched_fit_cFR_);
0607 cmsswfrtree_->Branch("fracHitsMatched_fit", &fracHitsMatched_fit_cFR_);
0608 cmsswfrtree_->Branch("lastlyr_fit", &lastlyr_fit_cFR_);
0609
0610 cmsswfrtree_->Branch("xhit_fit", &xhit_fit_cFR_);
0611 cmsswfrtree_->Branch("yhit_fit", &yhit_fit_cFR_);
0612 cmsswfrtree_->Branch("zhit_fit", &zhit_fit_cFR_);
0613
0614 cmsswfrtree_->Branch("hitchi2_fit", &hitchi2_fit_cFR_);
0615 cmsswfrtree_->Branch("helixchi2_fit", &helixchi2_fit_cFR_);
0616 cmsswfrtree_->Branch("score_fit", &score_fit_cFR_);
0617 cmsswfrtree_->Branch("dphi_fit", &dphi_fit_cFR_);
0618
0619 cmsswfrtree_->Branch("duplmask_fit", &duplmask_fit_cFR_);
0620 cmsswfrtree_->Branch("iTkMatches_fit", &iTkMatches_fit_cFR_);
0621
0622 cmsswfrtree_->Branch("seedID_cmssw_fit", &seedID_cmssw_fit_cFR_);
0623
0624 cmsswfrtree_->Branch("x_cmssw_fit", &x_cmssw_fit_cFR_);
0625 cmsswfrtree_->Branch("y_cmssw_fit", &y_cmssw_fit_cFR_);
0626 cmsswfrtree_->Branch("z_cmssw_fit", &z_cmssw_fit_cFR_);
0627
0628 cmsswfrtree_->Branch("pt_cmssw_fit", &pt_cmssw_fit_cFR_);
0629 cmsswfrtree_->Branch("phi_cmssw_fit", &phi_cmssw_fit_cFR_);
0630 cmsswfrtree_->Branch("eta_cmssw_fit", &eta_cmssw_fit_cFR_);
0631
0632 cmsswfrtree_->Branch("nHits_cmssw_fit", &nHits_cmssw_fit_cFR_);
0633 cmsswfrtree_->Branch("nLayers_cmssw_fit", &nLayers_cmssw_fit_cFR_);
0634 cmsswfrtree_->Branch("lastlyr_cmssw_fit", &lastlyr_cmssw_fit_cFR_);
0635
0636 cmsswfrtree_->Branch("algorithm", &algorithm_cFR_);
0637
0638 if (Config::keepHitInfo) {
0639 cmsswfrtree_->Branch("hitlyrs_mc", &hitlyrs_mc_cFR_);
0640 cmsswfrtree_->Branch("hitlyrs_build", &hitlyrs_build_cFR_);
0641 cmsswfrtree_->Branch("hitlyrs_cmssw_build", &hitlyrs_cmssw_build_cFR_);
0642 cmsswfrtree_->Branch("hitlyrs_fit", &hitlyrs_fit_cFR_);
0643 cmsswfrtree_->Branch("hitlyrs_cmssw_fit", &hitlyrs_cmssw_fit_cFR_);
0644
0645 cmsswfrtree_->Branch("hitidxs_mc", &hitidxs_mc_cFR_);
0646 cmsswfrtree_->Branch("hitidxs_build", &hitidxs_build_cFR_);
0647 cmsswfrtree_->Branch("hitidxs_cmssw_build", &hitidxs_cmssw_build_cFR_);
0648 cmsswfrtree_->Branch("hitidxs_fit", &hitidxs_fit_cFR_);
0649 cmsswfrtree_->Branch("hitidxs_cmssw_fit", &hitidxs_cmssw_fit_cFR_);
0650 }
0651 }
0652
0653 void TTreeValidation::initializeFitTree() {
0654 fittree_ = std::make_unique<TTree>("fittree", "fittree");
0655 fittree_->SetDirectory(0);
0656
0657 fittree_->Branch("ntotallayers", &ntotallayers_fit_, "ntotallayers_fit_/I");
0658 fittree_->Branch("tkid", &tkid_fit_, "tkid_fit_/I");
0659 fittree_->Branch("evtid", &evtid_fit_, "evtid_fit_/I");
0660
0661 fittree_->Branch("z_prop", &z_prop_fit_, "z_prop_fit_[ntotallayers_fit_]/F");
0662 fittree_->Branch("ez_prop", &ez_prop_fit_, "ez_prop_fit_[ntotallayers_fit_]/F");
0663 fittree_->Branch("z_hit", &z_hit_fit_, "z_hit_fit_[ntotallayers_fit_]/F");
0664 fittree_->Branch("ez_hit", &ez_hit_fit_, "ez_hit_fit_[ntotallayers_fit_]/F");
0665 fittree_->Branch("z_sim", &z_sim_fit_, "z_sim_fit_[ntotallayers_fit_]/F");
0666 fittree_->Branch("ez_sim", &ez_sim_fit_, "ez_sim_fit_[ntotallayers_fit_]/F");
0667
0668 fittree_->Branch("pphi_prop", &pphi_prop_fit_, "pphi_prop_fit_[ntotallayers_fit_]/F");
0669 fittree_->Branch("epphi_prop", &epphi_prop_fit_, "epphi_prop_fit_[ntotallayers_fit_]/F");
0670 fittree_->Branch("pphi_hit", &pphi_hit_fit_, "pphi_hit_fit_[ntotallayers_fit_]/F");
0671 fittree_->Branch("epphi_hit", &epphi_hit_fit_, "epphi_hit_fit_[ntotallayers_fit_]/F");
0672 fittree_->Branch("pphi_sim", &pphi_sim_fit_, "pphi_sim_fit_[ntotallayers_fit_]/F");
0673 fittree_->Branch("epphi_sim", &epphi_sim_fit_, "epphi_sim_fit_[ntotallayers_fit_]/F");
0674
0675 fittree_->Branch("pt_up", &pt_up_fit_, "pt_up_fit_[ntotallayers_fit_]/F");
0676 fittree_->Branch("ept_up", &ept_up_fit_, "ept_up_fit_[ntotallayers_fit_]/F");
0677 fittree_->Branch("pt_sim", &pt_sim_fit_, "pt_sim_fit_[ntotallayers_fit_]/F");
0678 fittree_->Branch("ept_sim", &ept_sim_fit_, "ept_sim_fit_[ntotallayers_fit_]/F");
0679
0680 fittree_->Branch("mphi_up", &mphi_up_fit_, "mphi_up_fit_[ntotallayers_fit_]/F");
0681 fittree_->Branch("emphi_up", &emphi_up_fit_, "emphi_up_fit_[ntotallayers_fit_]/F");
0682 fittree_->Branch("mphi_sim", &mphi_sim_fit_, "mphi_sim_fit_[ntotallayers_fit_]/F");
0683 fittree_->Branch("emphi_sim", &emphi_sim_fit_, "emphi_sim_fit_[ntotallayers_fit_]/F");
0684
0685 fittree_->Branch("meta_up", &meta_up_fit_, "meta_up_fit_[ntotallayers_fit_]/F");
0686 fittree_->Branch("emeta_up", &emeta_up_fit_, "emeta_up_fit_[ntotallayers_fit_]/F");
0687 fittree_->Branch("meta_sim", &meta_sim_fit_, "meta_sim_fit_[ntotallayers_fit_]/F");
0688 fittree_->Branch("emeta_sim", &emeta_sim_fit_, "emeta_sim_fit_[ntotallayers_fit_]/F");
0689 }
0690
0691 void TTreeValidation::alignTracks(TrackVec& evt_tracks, TrackExtraVec& evt_extras, bool alignExtra) {
0692 std::lock_guard<std::mutex> locker(glock_);
0693
0694
0695 if (alignExtra) {
0696 TrackExtraVec trackExtra_tmp(evt_tracks.size());
0697
0698
0699 for (int itrack = 0; itrack < (int)evt_tracks.size(); itrack++) {
0700 trackExtra_tmp[itrack] = evt_extras[evt_tracks[itrack].label()];
0701 }
0702
0703
0704 evt_extras = trackExtra_tmp;
0705 }
0706
0707
0708 for (int itrack = 0; itrack < (int)evt_tracks.size(); itrack++) {
0709 evt_tracks[itrack].setLabel(itrack);
0710 }
0711 }
0712
0713 void TTreeValidation::collectFitInfo(const FitVal& tmpfitval, int tkid, int layer) {
0714 std::lock_guard<std::mutex> locker(glock_);
0715
0716 fitValTkMapMap_[tkid][layer] = tmpfitval;
0717 }
0718
0719 void TTreeValidation::resetValidationMaps() {
0720 std::lock_guard<std::mutex> locker(glock_);
0721
0722 fitValTkMapMap_.clear();
0723
0724
0725 simToSeedMap_.clear();
0726 simToBuildMap_.clear();
0727 simToFitMap_.clear();
0728
0729
0730 seedToBuildMap_.clear();
0731 seedToFitMap_.clear();
0732
0733
0734 cmsswToBuildMap_.clear();
0735 cmsswToFitMap_.clear();
0736
0737
0738 seedToCmsswMap_.clear();
0739 cmsswToSeedMap_.clear();
0740
0741
0742 buildToCmsswMap_.clear();
0743
0744
0745 buildToFitMap_.clear();
0746 fitToBuildMap_.clear();
0747
0748
0749 candToSeedMapDumbCMSSW_.clear();
0750 fitToSeedMapDumbCMSSW_.clear();
0751 }
0752
0753 void TTreeValidation::setTrackExtras(Event& ev) {
0754 std::lock_guard<std::mutex> locker(glock_);
0755
0756 const auto& layerhits = ev.layerHits_;
0757
0758 if (Config::sim_val_for_cmssw || Config::sim_val) {
0759 const auto& simhits = ev.simHitsInfo_;
0760 const auto& simtracks = ev.simTracks_;
0761 const auto& seedtracks = ev.seedTracks_;
0762 auto& seedextras = ev.seedTracksExtra_;
0763 const auto& buildtracks = ev.candidateTracks_;
0764 auto& buildextras = ev.candidateTracksExtra_;
0765 const auto& fittracks = ev.fitTracks_;
0766 auto& fitextras = ev.fitTracksExtra_;
0767
0768
0769 for (int itrack = 0; itrack < (int)seedtracks.size(); itrack++) {
0770 const auto& track = seedtracks[itrack];
0771 auto& extra = seedextras[itrack];
0772
0773 extra.findMatchingSeedHits(track, track, layerhits);
0774 extra.setMCTrackIDInfo(
0775 track,
0776 layerhits,
0777 simhits,
0778 simtracks,
0779 true,
0780 (Config::seedInput == simSeeds));
0781 }
0782
0783
0784 for (int itrack = 0; itrack < (int)buildtracks.size(); itrack++) {
0785 const auto& track = buildtracks[itrack];
0786 auto& extra = buildextras[itrack];
0787
0788 if (Config::sim_val) {
0789 extra.findMatchingSeedHits(track, seedtracks[track.label()], layerhits);
0790 } else if (Config::sim_val_for_cmssw) {
0791 extra.findMatchingSeedHits(track, seedtracks[candToSeedMapDumbCMSSW_[track.label()]], layerhits);
0792 }
0793
0794 extra.setMCTrackIDInfo(track, layerhits, simhits, simtracks, false, (Config::seedInput == simSeeds));
0795 }
0796
0797
0798 for (int itrack = 0; itrack < (int)fittracks.size(); itrack++) {
0799 const auto& track = fittracks[itrack];
0800 auto& extra = fitextras[itrack];
0801
0802 if (Config::sim_val) {
0803 extra.findMatchingSeedHits(track, seedtracks[track.label()], layerhits);
0804 } else if (Config::sim_val_for_cmssw) {
0805 extra.findMatchingSeedHits(track, seedtracks[fitToSeedMapDumbCMSSW_[track.label()]], layerhits);
0806 }
0807
0808 extra.setMCTrackIDInfo(track, layerhits, simhits, simtracks, false, (Config::seedInput == simSeeds));
0809 }
0810 }
0811
0812 if (Config::cmssw_val) {
0813
0814 storeSeedAndMCID(ev);
0815
0816 const auto& cmsswtracks = ev.cmsswTracks_;
0817 const auto& cmsswextras = ev.cmsswTracksExtra_;
0818 const auto& seedtracks = ev.seedTracks_;
0819 const auto& buildtracks = ev.candidateTracks_;
0820 auto& buildextras = ev.candidateTracksExtra_;
0821 const auto& fittracks = ev.fitTracks_;
0822 auto& fitextras = ev.fitTracksExtra_;
0823
0824
0825 RedTrackVec reducedCMSSW;
0826 LayIdxIDVecMapMap cmsswHitIDMap;
0827 setupCMSSWMatching(ev, reducedCMSSW, cmsswHitIDMap);
0828
0829
0830 for (int itrack = 0; itrack < (int)buildtracks.size(); itrack++) {
0831 const auto& track = buildtracks[itrack];
0832 auto& extra = buildextras[itrack];
0833
0834
0835 extra.findMatchingSeedHits(track,
0836 seedtracks[track.label()],
0837 layerhits);
0838
0839 if (Config::cmsswMatchingFW == trkParamBased) {
0840 extra.setCMSSWTrackIDInfoByTrkParams(track, layerhits, cmsswtracks, reducedCMSSW, true);
0841 } else if (Config::cmsswMatchingFW == hitBased) {
0842 extra.setCMSSWTrackIDInfoByHits(track,
0843 cmsswHitIDMap,
0844 cmsswtracks,
0845 cmsswextras,
0846 reducedCMSSW,
0847 -1);
0848 } else if (Config::cmsswMatchingFW == labelBased)
0849 {
0850 extra.setCMSSWTrackIDInfoByHits(track,
0851 cmsswHitIDMap,
0852 cmsswtracks,
0853 cmsswextras,
0854 reducedCMSSW,
0855 reducedCMSSW[cmsswtracks[buildToCmsswMap_[track.label()]].label()].label());
0856 } else {
0857 std::cerr << "Specified CMSSW validation, but using an incorrect matching option! Exiting..." << std::endl;
0858 exit(1);
0859 }
0860 }
0861
0862
0863 for (int itrack = 0; itrack < (int)fittracks.size(); itrack++) {
0864 const auto& track = fittracks[itrack];
0865 auto& extra = fitextras[itrack];
0866
0867
0868 extra.findMatchingSeedHits(track,
0869 seedtracks[track.label()],
0870 layerhits);
0871
0872 if (Config::cmsswMatchingBK == trkParamBased) {
0873 extra.setCMSSWTrackIDInfoByTrkParams(track, layerhits, cmsswtracks, reducedCMSSW, true);
0874 } else if (Config::cmsswMatchingBK == hitBased) {
0875 extra.setCMSSWTrackIDInfoByHits(track,
0876 cmsswHitIDMap,
0877 cmsswtracks,
0878 cmsswextras,
0879 reducedCMSSW,
0880 -1);
0881 } else if (Config::cmsswMatchingBK == labelBased)
0882 {
0883 extra.setCMSSWTrackIDInfoByHits(
0884 track,
0885 cmsswHitIDMap,
0886 cmsswtracks,
0887 cmsswextras,
0888 reducedCMSSW,
0889 reducedCMSSW[cmsswtracks[buildToCmsswMap_[fitToBuildMap_[track.label()]]].label()].label());
0890 } else {
0891 std::cerr << "Specified CMSSW validation, but using an incorrect matching option! Exiting..." << std::endl;
0892 exit(1);
0893 }
0894 }
0895 }
0896 }
0897
0898 void TTreeValidation::makeSimTkToRecoTksMaps(Event& ev) {
0899 std::lock_guard<std::mutex> locker(glock_);
0900
0901 TTreeValidation::mapRefTkToRecoTks(ev.seedTracks_, ev.seedTracksExtra_, simToSeedMap_);
0902 TTreeValidation::mapRefTkToRecoTks(ev.candidateTracks_, ev.candidateTracksExtra_, simToBuildMap_);
0903 TTreeValidation::mapRefTkToRecoTks(ev.fitTracks_, ev.fitTracksExtra_, simToFitMap_);
0904 }
0905
0906 void TTreeValidation::mapRefTkToRecoTks(const TrackVec& evt_tracks,
0907 TrackExtraVec& evt_extras,
0908 TkIDToTkIDVecMap& refTkMap) {
0909 for (auto itrack = 0; itrack < (int)evt_tracks.size(); ++itrack) {
0910 auto&& track(evt_tracks[itrack]);
0911 auto&& extra(evt_extras[itrack]);
0912 if (Config::sim_val_for_cmssw || Config::sim_val) {
0913 if (extra.mcTrackID() >= 0)
0914 {
0915 refTkMap[extra.mcTrackID()].push_back(
0916 track.label());
0917 }
0918 }
0919 if (Config::cmssw_val) {
0920 if (extra.cmsswTrackID() >= 0)
0921 {
0922 refTkMap[extra.cmsswTrackID()].push_back(
0923 track.label());
0924 }
0925 }
0926 }
0927
0928 for (auto&& refTkMatches : refTkMap) {
0929 if (refTkMatches.second.size() < 2)
0930 {
0931 auto& extra(evt_extras[refTkMatches.second[0]]);
0932 extra.setDuplicateInfo(0, bool(false));
0933 } else
0934 {
0935
0936
0937 TrackVec tmpMatches;
0938 for (auto&& label :
0939 refTkMatches.second)
0940 {
0941 tmpMatches.emplace_back(evt_tracks[label]);
0942 }
0943
0944 std::sort(tmpMatches.begin(), tmpMatches.end(), sortByScoreCand);
0945 for (auto itrack = 0; itrack < (int)tmpMatches.size();
0946 itrack++)
0947 {
0948 refTkMatches.second[itrack] = tmpMatches[itrack].label();
0949 }
0950
0951 int duplicateID = 0;
0952 for (auto&& label : refTkMatches.second)
0953 {
0954 auto& extra(evt_extras[label]);
0955 extra.setDuplicateInfo(duplicateID, bool(true));
0956 duplicateID++;
0957 }
0958 }
0959 }
0960 }
0961
0962 void TTreeValidation::makeSeedTkToRecoTkMaps(Event& ev) {
0963 std::lock_guard<std::mutex> locker(glock_);
0964
0965 TTreeValidation::mapSeedTkToRecoTk(ev.candidateTracks_, ev.candidateTracksExtra_, seedToBuildMap_);
0966 TTreeValidation::mapSeedTkToRecoTk(ev.fitTracks_, ev.fitTracksExtra_, seedToFitMap_);
0967 }
0968
0969 void TTreeValidation::mapSeedTkToRecoTk(const TrackVec& evt_tracks,
0970 const TrackExtraVec& evt_extras,
0971 TkIDToTkIDMap& seedTkMap) {
0972 for (auto&& track : evt_tracks) {
0973 seedTkMap[evt_extras[track.label()].seedID()] = track.label();
0974 }
0975 }
0976
0977 void TTreeValidation::makeRecoTkToRecoTkMaps(Event& ev) {
0978 std::lock_guard<std::mutex> locker(glock_);
0979 TTreeValidation::makeRecoTkToRecoTkMap(
0980 buildToFitMap_, ev.candidateTracks_, ev.candidateTracksExtra_, ev.fitTracks_, ev.fitTracksExtra_);
0981 TTreeValidation::makeRecoTkToRecoTkMap(
0982 fitToBuildMap_, ev.fitTracks_, ev.fitTracksExtra_, ev.candidateTracks_, ev.candidateTracksExtra_);
0983 }
0984
0985 void TTreeValidation::makeRecoTkToRecoTkMap(TkIDToTkIDMap& refToPairMap,
0986 const TrackVec& reftracks,
0987 const TrackExtraVec& refextras,
0988 const TrackVec& pairtracks,
0989 const TrackExtraVec& pairextras) {
0990
0991
0992 for (auto&& reftrack : reftracks) {
0993 const auto& refextra = refextras[reftrack.label()];
0994 for (auto&& pairtrack : pairtracks) {
0995 const auto& pairextra = pairextras[pairtrack.label()];
0996 if (refextra.seedID() == pairextra.seedID()) {
0997 refToPairMap[reftrack.label()] = pairtrack.label();
0998 break;
0999 }
1000 }
1001 }
1002 }
1003
1004 void TTreeValidation::makeCMSSWTkToRecoTksMaps(Event& ev) {
1005 std::lock_guard<std::mutex> locker(glock_);
1006
1007 TTreeValidation::mapRefTkToRecoTks(ev.candidateTracks_, ev.candidateTracksExtra_, cmsswToBuildMap_);
1008 TTreeValidation::mapRefTkToRecoTks(ev.fitTracks_, ev.fitTracksExtra_, cmsswToFitMap_);
1009 }
1010
1011 void TTreeValidation::makeSeedTkToCMSSWTkMap(Event& ev) {
1012 const auto& seedtracks = ev.seedTracks_;
1013 const auto& cmsswtracks = ev.cmsswTracks_;
1014 for (int itrack = 0; itrack < (int)seedtracks.size(); itrack++) {
1015 for (auto&& cmsswtrack : cmsswtracks) {
1016 if (cmsswtrack.label() == itrack) {
1017 seedToCmsswMap_[seedtracks[itrack].label()] = cmsswtrack.label();
1018 break;
1019 }
1020 }
1021 }
1022 }
1023
1024 void TTreeValidation::makeCMSSWTkToSeedTkMap(Event& ev) {
1025 const auto& seedtracks = ev.seedTracks_;
1026
1027 for (const auto& seedToCmsswPair : seedToCmsswMap_) {
1028 const auto seedlabel =
1029 seedToCmsswPair
1030 .first;
1031 const auto cmsswlabel = seedToCmsswPair.second;
1032
1033 for (int itrack = 0; itrack < (int)seedtracks.size(); itrack++) {
1034 const auto& seedtrack = seedtracks[itrack];
1035 if (seedtrack.label() == seedlabel) {
1036 cmsswToSeedMap_[cmsswlabel] = itrack;
1037 break;
1038 }
1039 }
1040 }
1041 }
1042
1043 void TTreeValidation::makeRecoTkToSeedTkMapsDumbCMSSW(Event& ev) {
1044 std::lock_guard<std::mutex> locker(glock_);
1045
1046 TTreeValidation::makeRecoTkToSeedTkMapDumbCMSSW(
1047 ev.candidateTracksExtra_, ev.seedTracksExtra_, candToSeedMapDumbCMSSW_);
1048 TTreeValidation::makeRecoTkToSeedTkMapDumbCMSSW(ev.fitTracksExtra_, ev.seedTracksExtra_, fitToSeedMapDumbCMSSW_);
1049 }
1050
1051 void TTreeValidation::makeRecoTkToSeedTkMapDumbCMSSW(const TrackExtraVec& recoextras,
1052 const TrackExtraVec& seedextras,
1053 TkIDToTkIDMap& recoToSeedMap) {
1054 for (int itrack = 0; itrack < (int)recoextras.size(); itrack++) {
1055 const auto reco_seedID = recoextras[itrack].seedID();
1056 for (int jtrack = 0; jtrack < (int)seedextras.size(); jtrack++) {
1057 const auto seed_seedID = seedextras[jtrack].seedID();
1058 if (reco_seedID == seed_seedID) {
1059 recoToSeedMap[itrack] = jtrack;
1060 break;
1061 }
1062 }
1063 }
1064 }
1065
1066 void TTreeValidation::setTrackScoresDumbCMSSW(Event& ev) {
1067 auto& seedtracks = ev.seedTracks_;
1068 auto& candtracks = ev.candidateTracks_;
1069 auto& fittracks = ev.fitTracks_;
1070 auto score_calc = IterationConfig::get_track_scorer("default");
1071
1072
1073 for (auto& seedtrack : seedtracks) {
1074 seedtrack.setScore(getScoreCand(score_calc, seedtrack));
1075 }
1076
1077
1078 for (const auto& candToSeedPair : candToSeedMapDumbCMSSW_) {
1079 auto& candtrack = candtracks[candToSeedPair.first];
1080
1081 candtrack.setScore(getScoreCand(score_calc, candtrack));
1082 }
1083 for (const auto& fitToSeedPair : fitToSeedMapDumbCMSSW_) {
1084 auto& fittrack = fittracks[fitToSeedPair.first];
1085
1086 fittrack.setScore(getScoreCand(score_calc, fittrack));
1087 }
1088 }
1089
1090 void TTreeValidation::storeSeedAndMCID(Event& ev) {
1091 const auto& buildtracks = ev.candidateTracks_;
1092 auto& buildextras = ev.candidateTracksExtra_;
1093
1094 const auto& fittracks = ev.fitTracks_;
1095 auto& fitextras = ev.fitTracksExtra_;
1096
1097 const auto& cmsswtracks = ev.cmsswTracks_;
1098 auto& cmsswextras = ev.cmsswTracksExtra_;
1099
1100
1101 int newlabel = -1;
1102 for (int itrack = 0; itrack < (int)buildtracks.size(); itrack++) {
1103 auto& extra = buildextras[itrack];
1104 const int seedID = extra.seedID();
1105
1106 extra.setmcTrackID(seedID);
1107
1108 if (seedToCmsswMap_.count(seedID)) {
1109 extra.setseedID(seedToCmsswMap_[seedID]);
1110 if (Config::cmsswMatchingFW == labelBased || Config::cmsswMatchingBK == labelBased) {
1111 for (int ctrack = 0; ctrack < (int)cmsswextras.size(); ctrack++) {
1112 if (cmsswextras[ctrack].seedID() == extra.seedID()) {
1113 buildToCmsswMap_[itrack] = cmsswtracks[ctrack].label();
1114 break;
1115 }
1116 }
1117 }
1118 } else {
1119 extra.setseedID(--newlabel);
1120 }
1121 }
1122
1123
1124 for (int itrack = 0; itrack < (int)fittracks.size(); itrack++) {
1125 auto& extra = fitextras[itrack];
1126
1127 extra.setmcTrackID(buildextras[fitToBuildMap_[itrack]].mcTrackID());
1128 extra.setseedID(buildextras[fitToBuildMap_[itrack]].seedID());
1129 }
1130 }
1131
1132 void TTreeValidation::setupCMSSWMatching(const Event& ev,
1133 RedTrackVec& reducedCMSSW,
1134 LayIdxIDVecMapMap& cmsswHitIDMap) {
1135
1136 const auto& layerhits = ev.layerHits_;
1137 const auto& cmsswtracks = ev.cmsswTracks_;
1138 auto& cmsswextras = ev.cmsswTracksExtra_;
1139 const auto& seedtracks = ev.seedTracks_;
1140
1141
1142 reducedCMSSW.resize(cmsswtracks.size());
1143
1144 for (int itrack = 0; itrack < (int)cmsswtracks.size(); itrack++) {
1145
1146 auto& cmsswextra = cmsswextras[itrack];
1147 const auto& cmsswtrack = cmsswtracks[itrack];
1148 const auto& seedtrack = seedtracks[cmsswToSeedMap_[cmsswtrack.label()]];
1149
1150
1151 cmsswextra.findMatchingSeedHits(cmsswtrack, seedtrack, layerhits);
1152
1153
1154 const auto seedID = cmsswextra.seedID();
1155 const auto& params = cmsswtrack.parameters();
1156 SVector2 tmpv(params[3], params[5]);
1157
1158 HitLayerMap tmpmap;
1159 for (int ihit = 0; ihit < cmsswtrack.nTotalHits(); ihit++) {
1160 const int lyr = cmsswtrack.getHitLyr(ihit);
1161 const int idx = cmsswtrack.getHitIdx(ihit);
1162
1163
1164 if (cmsswextra.isSeedHit(lyr, idx))
1165 continue;
1166
1167 if (lyr >= 0 && idx >= 0) {
1168 tmpmap[lyr].push_back(idx);
1169 cmsswHitIDMap[lyr][idx].push_back(cmsswtrack.label());
1170 }
1171 }
1172
1173
1174 reducedCMSSW[itrack] = ReducedTrack(cmsswtrack.label(), seedID, tmpv, cmsswtrack.momPhi(), tmpmap);
1175 }
1176 }
1177
1178 int TTreeValidation::getLastFoundHit(const int trackMCHitID, const int mcTrackID, const Event& ev) {
1179 int mcHitID = -1;
1180 if (ev.simHitsInfo_[trackMCHitID].mcTrackID() == mcTrackID) {
1181 mcHitID = trackMCHitID;
1182 } else {
1183 mcHitID = ev.simTracks_[mcTrackID].getMCHitIDFromLayer(ev.layerHits_, ev.simHitsInfo_[trackMCHitID].layer());
1184 }
1185 return mcHitID;
1186 }
1187
1188 int TTreeValidation::getMaskAssignment(const int refID) {
1189
1190 auto refmask = -99;
1191
1192 if (refID >= 0)
1193 {
1194 refmask = 1;
1195 } else if (refID == -10) {
1196 refmask = -2;
1197 } else {
1198 if (Config::inclusiveShorts)
1199 {
1200 if (refID == -1 || refID == -5 || refID == -8 || refID == -9) {
1201 refmask = 0;
1202 } else if (refID == -2) {
1203 refmask = 2;
1204 } else
1205 {
1206 refmask = -1;
1207 }
1208 } else
1209 {
1210 if (refID == -1 || refID == -9) {
1211 refmask = 0;
1212 } else if (Config::mtvLikeValidation && refID == -4) {
1213 refmask = 2;
1214 } else
1215 {
1216 refmask = -1;
1217 }
1218 }
1219 }
1220
1221 return refmask;
1222 }
1223
1224 void TTreeValidation::resetFitBranches() {
1225 for (int ilayer = 0; ilayer < ntotallayers_fit_; ++ilayer) {
1226 z_prop_fit_[ilayer] = -1000.f;
1227 ez_prop_fit_[ilayer] = -1000.f;
1228 z_hit_fit_[ilayer] = -1000.f;
1229 ez_hit_fit_[ilayer] = -1000.f;
1230 z_sim_fit_[ilayer] = -1000.f;
1231 ez_sim_fit_[ilayer] = -1000.f;
1232
1233 pphi_prop_fit_[ilayer] = -1000.f;
1234 epphi_prop_fit_[ilayer] = -1000.f;
1235 pphi_hit_fit_[ilayer] = -1000.f;
1236 epphi_hit_fit_[ilayer] = -1000.f;
1237 pphi_sim_fit_[ilayer] = -1000.f;
1238 epphi_sim_fit_[ilayer] = -1000.f;
1239
1240 pt_up_fit_[ilayer] = -1000.f;
1241 ept_up_fit_[ilayer] = -1000.f;
1242 pt_sim_fit_[ilayer] = -1000.f;
1243 ept_sim_fit_[ilayer] = -1000.f;
1244
1245 mphi_up_fit_[ilayer] = -1000.f;
1246 emphi_up_fit_[ilayer] = -1000.f;
1247 mphi_sim_fit_[ilayer] = -1000.f;
1248 emphi_sim_fit_[ilayer] = -1000.f;
1249
1250 meta_up_fit_[ilayer] = -1000.f;
1251 emeta_up_fit_[ilayer] = -1000.f;
1252 meta_sim_fit_[ilayer] = -1000.f;
1253 emeta_sim_fit_[ilayer] = -1000.f;
1254 }
1255 }
1256
1257 void TTreeValidation::fillFitTree(const Event& ev) {
1258 std::lock_guard<std::mutex> locker(glock_);
1259
1260 evtid_fit_ = ev.evtID();
1261 const auto& simtracks = ev.simTracks_;
1262 const auto& layerhits = ev.layerHits_;
1263 const auto& simtrackstates = ev.simTrackStates_;
1264
1265 for (auto&& fitvalmapmap : fitValTkMapMap_) {
1266 TTreeValidation::resetFitBranches();
1267
1268 tkid_fit_ = fitvalmapmap.first;
1269
1270 const auto& simtrack = simtracks[tkid_fit_];
1271 const auto& fitvalmap = fitvalmapmap.second;
1272 for (int ilayer = 0; ilayer < ntotallayers_fit_; ++ilayer) {
1273 if (fitvalmap.count(ilayer)) {
1274 const auto& hit = layerhits[ilayer][simtrack.getHitIdx(ilayer)];
1275 const auto& initTS = simtrackstates.at(hit.mcHitID());
1276 const auto& fitval = fitvalmap.at(ilayer);
1277
1278 z_hit_fit_[ilayer] = hit.z();
1279 ez_hit_fit_[ilayer] = std::sqrt(hit.ezz());
1280 z_sim_fit_[ilayer] = initTS.z();
1281 ez_sim_fit_[ilayer] = initTS.ezz();
1282 z_prop_fit_[ilayer] = fitval.ppz;
1283 ez_prop_fit_[ilayer] = fitval.eppz;
1284
1285 pphi_hit_fit_[ilayer] = hit.phi();
1286 epphi_hit_fit_[ilayer] = std::sqrt(hit.ephi());
1287 pphi_sim_fit_[ilayer] = initTS.posPhi();
1288 epphi_sim_fit_[ilayer] = initTS.eposPhi();
1289 pphi_prop_fit_[ilayer] = fitval.ppphi;
1290 epphi_prop_fit_[ilayer] = fitval.eppphi;
1291
1292 pt_up_fit_[ilayer] = fitval.upt;
1293 ept_up_fit_[ilayer] = fitval.eupt;
1294 pt_sim_fit_[ilayer] = initTS.pT();
1295 ept_sim_fit_[ilayer] = initTS.epT();
1296
1297 mphi_up_fit_[ilayer] = fitval.umphi;
1298 emphi_up_fit_[ilayer] = fitval.eumphi;
1299 mphi_sim_fit_[ilayer] = initTS.momPhi();
1300 emphi_sim_fit_[ilayer] = initTS.emomPhi();
1301
1302 meta_up_fit_[ilayer] = fitval.umeta;
1303 emeta_up_fit_[ilayer] = fitval.eumeta;
1304 meta_sim_fit_[ilayer] = initTS.momEta();
1305 emeta_sim_fit_[ilayer] = initTS.emomEta();
1306 }
1307 }
1308 fittree_->Fill();
1309 }
1310 }
1311
1312 void TTreeValidation::fillFullHitInfo(const Event& ev,
1313 const Track& track,
1314 std::vector<int>& lyrs,
1315 std::vector<int>& idxs,
1316 std::vector<int>& mcTkIDs,
1317 std::vector<float>& xs,
1318 std::vector<float>& ys,
1319 std::vector<float>& zs) {
1320
1321 const auto& layerHits = ev.layerHits_;
1322 const auto& simHitsInfo = ev.simHitsInfo_;
1323
1324
1325 const auto nTotalHits = track.nTotalHits();
1326 lyrs.resize(nTotalHits);
1327 idxs.resize(nTotalHits);
1328 mcTkIDs.resize(nTotalHits, -99);
1329 xs.resize(nTotalHits, -9999.f);
1330 ys.resize(nTotalHits, -9999.f);
1331 zs.resize(nTotalHits, -9999.f);
1332
1333
1334 for (auto ihit = 0; ihit < nTotalHits; ihit++) {
1335 const auto lyr = track.getHitLyr(ihit);
1336 const auto idx = track.getHitIdx(ihit);
1337
1338 lyrs[ihit] = lyr;
1339 idxs[ihit] = idx;
1340
1341 if (lyr < 0)
1342 continue;
1343 if (idx < 0)
1344 continue;
1345
1346 const auto& hit = layerHits[lyr][idx];
1347 mcTkIDs[ihit] = hit.mcTrackID(simHitsInfo);
1348 xs[ihit] = hit.x();
1349 ys[ihit] = hit.y();
1350 zs[ihit] = hit.z();
1351 }
1352 }
1353
1354 void TTreeValidation::fillMinHitInfo(const Track& track, std::vector<int>& lyrs, std::vector<int>& idxs) {
1355 for (int ihit = 0; ihit < track.nTotalHits(); ihit++) {
1356 lyrs.emplace_back(track.getHitLyr(ihit));
1357 idxs.emplace_back(track.getHitIdx(ihit));
1358 }
1359 }
1360
1361 void TTreeValidation::fillEfficiencyTree(const Event& ev) {
1362 std::lock_guard<std::mutex> locker(glock_);
1363
1364 const auto ievt = ev.evtID();
1365 const auto& evt_sim_tracks = ev.simTracks_;
1366 const auto& evt_seed_tracks = ev.seedTracks_;
1367 const auto& evt_seed_extras = ev.seedTracksExtra_;
1368 const auto& evt_build_tracks = ev.candidateTracks_;
1369 const auto& evt_build_extras = ev.candidateTracksExtra_;
1370 const auto& evt_fit_tracks = ev.fitTracks_;
1371 const auto& evt_fit_extras = ev.fitTracksExtra_;
1372 const auto& evt_layer_hits = ev.layerHits_;
1373 const auto& evt_sim_trackstates = ev.simTrackStates_;
1374
1375 unsigned int count = 0;
1376 for (const auto& simtrack : evt_sim_tracks) {
1377
1378 if (Config::keepHitInfo) {
1379 hitlyrs_mc_eff_.clear();
1380 hitlyrs_seed_eff_.clear();
1381 hitlyrs_build_eff_.clear();
1382 hitlyrs_fit_eff_.clear();
1383
1384 hitidxs_mc_eff_.clear();
1385 hitidxs_seed_eff_.clear();
1386 hitidxs_build_eff_.clear();
1387 hitidxs_fit_eff_.clear();
1388
1389 hitmcTkIDs_mc_eff_.clear();
1390 hitmcTkIDs_seed_eff_.clear();
1391 hitmcTkIDs_build_eff_.clear();
1392 hitmcTkIDs_fit_eff_.clear();
1393
1394 hitxs_mc_eff_.clear();
1395 hitxs_seed_eff_.clear();
1396 hitxs_build_eff_.clear();
1397 hitxs_fit_eff_.clear();
1398
1399 hitys_mc_eff_.clear();
1400 hitys_seed_eff_.clear();
1401 hitys_build_eff_.clear();
1402 hitys_fit_eff_.clear();
1403
1404 hitzs_mc_eff_.clear();
1405 hitzs_seed_eff_.clear();
1406 hitzs_build_eff_.clear();
1407 hitzs_fit_eff_.clear();
1408 }
1409
1410 evtID_eff_ = ievt;
1411 mcID_eff_ = simtrack.label();
1412
1413
1414 x_mc_gen_eff_ = simtrack.x();
1415 y_mc_gen_eff_ = simtrack.y();
1416 z_mc_gen_eff_ = simtrack.z();
1417
1418 pt_mc_gen_eff_ = simtrack.pT();
1419 phi_mc_gen_eff_ = simtrack.momPhi();
1420 eta_mc_gen_eff_ = simtrack.momEta();
1421 nHits_mc_eff_ = simtrack.nFoundHits();
1422 nLayers_mc_eff_ = simtrack.nUniqueLayers();
1423 lastlyr_mc_eff_ = simtrack.getLastFoundHitLyr();
1424
1425 itermask_seed_eff_ = 0;
1426 itermask_build_eff_ = 0;
1427 itermask_fit_eff_ = 0;
1428 iterduplmask_seed_eff_ = 0;
1429 iterduplmask_build_eff_ = 0;
1430 iterduplmask_fit_eff_ = 0;
1431 algo_seed_eff_ = 0;
1432
1433 if (Config::mtvRequireSeeds) {
1434 for (auto aa : ev.simTracksExtra_[count].seedAlgos()) {
1435 algo_seed_eff_ = (algo_seed_eff_ | (1 << aa));
1436 }
1437 }
1438 count++;
1439
1440
1441 if (Config::keepHitInfo)
1442 TTreeValidation::fillFullHitInfo(ev,
1443 simtrack,
1444 hitlyrs_mc_eff_,
1445 hitidxs_mc_eff_,
1446 hitmcTkIDs_mc_eff_,
1447 hitxs_mc_eff_,
1448 hitys_mc_eff_,
1449 hitzs_mc_eff_);
1450
1451
1452 if (simToSeedMap_.count(mcID_eff_) &&
1453 simtrack
1454 .isFindable())
1455 {
1456 for (unsigned int ii = 0; ii < simToSeedMap_[mcID_eff_].size(); ii++) {
1457 const int theAlgo = evt_seed_tracks[simToSeedMap_[mcID_eff_][ii]].algoint();
1458 if ((itermask_seed_eff_ >> theAlgo) & 1)
1459 iterduplmask_seed_eff_ = (iterduplmask_seed_eff_ | (1 << theAlgo));
1460 itermask_seed_eff_ = (itermask_seed_eff_ | (1 << theAlgo));
1461 }
1462 const auto& seedtrack =
1463 evt_seed_tracks[simToSeedMap_[mcID_eff_][0]];
1464 const auto& seedextra = evt_seed_extras[seedtrack.label()];
1465 mcmask_seed_eff_ = 1;
1466
1467 seedID_seed_eff_ = seedextra.seedID();
1468
1469
1470 const int mcHitID =
1471 TTreeValidation::getLastFoundHit(seedtrack.getLastFoundMCHitID(evt_layer_hits), mcID_eff_, ev);
1472 if (mcHitID >= 0 && Config::readSimTrackStates) {
1473 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
1474
1475 pt_mc_seed_eff_ = initLayTS.pT();
1476 phi_mc_seed_eff_ = initLayTS.momPhi();
1477 eta_mc_seed_eff_ = initLayTS.momEta();
1478 helixchi2_seed_eff_ = computeHelixChi2(initLayTS.parameters, seedtrack.parameters(), seedtrack.errors());
1479
1480 mcTSmask_seed_eff_ = 1;
1481 } else if (Config::tryToSaveSimInfo)
1482 {
1483
1484 pt_mc_seed_eff_ = pt_mc_gen_eff_;
1485 phi_mc_seed_eff_ = phi_mc_gen_eff_;
1486 eta_mc_seed_eff_ = eta_mc_gen_eff_;
1487 helixchi2_seed_eff_ = computeHelixChi2(simtrack.parameters(), seedtrack.parameters(), seedtrack.errors());
1488
1489 mcTSmask_seed_eff_ = 0;
1490 } else {
1491 pt_mc_seed_eff_ = -101;
1492 phi_mc_seed_eff_ = -101;
1493 eta_mc_seed_eff_ = -101;
1494 helixchi2_seed_eff_ = -101;
1495
1496 mcTSmask_seed_eff_ = -2;
1497 }
1498
1499
1500 const Hit& lasthit = evt_layer_hits[seedtrack.getLastFoundHitLyr()][seedtrack.getLastFoundHitIdx()];
1501 xhit_seed_eff_ = lasthit.x();
1502 yhit_seed_eff_ = lasthit.y();
1503 zhit_seed_eff_ = lasthit.z();
1504
1505 pt_seed_eff_ = seedtrack.pT();
1506 ept_seed_eff_ = seedtrack.epT();
1507 phi_seed_eff_ = seedtrack.momPhi();
1508 ephi_seed_eff_ = seedtrack.emomPhi();
1509 eta_seed_eff_ = seedtrack.momEta();
1510 eeta_seed_eff_ = seedtrack.emomEta();
1511
1512
1513 nHits_seed_eff_ = seedtrack.nFoundHits();
1514 nLayers_seed_eff_ = seedtrack.nUniqueLayers();
1515 nHitsMatched_seed_eff_ = seedextra.nHitsMatched();
1516 fracHitsMatched_seed_eff_ = seedextra.fracHitsMatched();
1517 lastlyr_seed_eff_ = seedtrack.getLastFoundHitLyr();
1518
1519
1520 dphi_seed_eff_ = seedextra.dPhi();
1521
1522
1523 hitchi2_seed_eff_ = seedtrack.chi2();
1524 score_seed_eff_ = seedtrack.score();
1525
1526 duplmask_seed_eff_ = seedextra.isDuplicate();
1527 nTkMatches_seed_eff_ = simToSeedMap_[mcID_eff_].size();
1528
1529
1530 if (Config::keepHitInfo)
1531 TTreeValidation::fillFullHitInfo(ev,
1532 seedtrack,
1533 hitlyrs_seed_eff_,
1534 hitidxs_seed_eff_,
1535 hitmcTkIDs_seed_eff_,
1536 hitxs_seed_eff_,
1537 hitys_seed_eff_,
1538 hitzs_seed_eff_);
1539 } else
1540 {
1541 mcmask_seed_eff_ = (simtrack.isFindable() ? 0 : -1);
1542
1543 seedID_seed_eff_ = -99;
1544
1545 pt_mc_seed_eff_ = -99;
1546 phi_mc_seed_eff_ = -99;
1547 eta_mc_seed_eff_ = -99;
1548 helixchi2_seed_eff_ = -99;
1549
1550 mcTSmask_seed_eff_ = -1;
1551
1552 xhit_seed_eff_ = -2000;
1553 yhit_seed_eff_ = -2000;
1554 zhit_seed_eff_ = -2000;
1555
1556 pt_seed_eff_ = -99;
1557 ept_seed_eff_ = -99;
1558 phi_seed_eff_ = -99;
1559 ephi_seed_eff_ = -99;
1560 eta_seed_eff_ = -99;
1561 eeta_seed_eff_ = -99;
1562
1563 nHits_seed_eff_ = -99;
1564 nLayers_seed_eff_ = -99;
1565 nHitsMatched_seed_eff_ = -99;
1566 fracHitsMatched_seed_eff_ = -99;
1567 lastlyr_seed_eff_ = -99;
1568
1569 dphi_seed_eff_ = -99;
1570
1571 hitchi2_seed_eff_ = -99;
1572 score_seed_eff_ = -17000;
1573
1574 duplmask_seed_eff_ = -1;
1575 nTkMatches_seed_eff_ = -99;
1576 }
1577
1578
1579 if (simToBuildMap_.count(mcID_eff_) &&
1580 simtrack
1581 .isFindable())
1582 {
1583 for (unsigned int ii = 0; ii < simToBuildMap_[mcID_eff_].size(); ii++) {
1584 const int theAlgo = evt_build_tracks[simToBuildMap_[mcID_eff_][ii]].algoint();
1585 if ((itermask_build_eff_ >> theAlgo) & 1)
1586 iterduplmask_build_eff_ = (iterduplmask_build_eff_ | (1 << theAlgo));
1587 itermask_build_eff_ = (itermask_build_eff_ | (1 << theAlgo));
1588 }
1589 const auto& buildtrack =
1590 evt_build_tracks[simToBuildMap_[mcID_eff_][0]];
1591 const auto& buildextra =
1592 evt_build_extras[buildtrack.label()];
1593 mcmask_build_eff_ = 1;
1594
1595 seedID_build_eff_ = buildextra.seedID();
1596
1597
1598 const int mcHitID =
1599 TTreeValidation::getLastFoundHit(buildtrack.getLastFoundMCHitID(evt_layer_hits), mcID_eff_, ev);
1600 if (mcHitID >= 0 && Config::readSimTrackStates) {
1601 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
1602
1603 pt_mc_build_eff_ = initLayTS.pT();
1604 phi_mc_build_eff_ = initLayTS.momPhi();
1605 eta_mc_build_eff_ = initLayTS.momEta();
1606 helixchi2_build_eff_ = computeHelixChi2(initLayTS.parameters, buildtrack.parameters(), buildtrack.errors());
1607
1608 mcTSmask_build_eff_ = 1;
1609 } else if (Config::tryToSaveSimInfo)
1610 {
1611
1612 pt_mc_build_eff_ = pt_mc_gen_eff_;
1613 phi_mc_build_eff_ = phi_mc_gen_eff_;
1614 eta_mc_build_eff_ = eta_mc_gen_eff_;
1615 helixchi2_build_eff_ = computeHelixChi2(simtrack.parameters(), buildtrack.parameters(), buildtrack.errors());
1616
1617 mcTSmask_build_eff_ = 0;
1618 } else {
1619 pt_mc_build_eff_ = -101;
1620 phi_mc_build_eff_ = -101;
1621 eta_mc_build_eff_ = -101;
1622 helixchi2_build_eff_ = -101;
1623
1624 mcTSmask_build_eff_ = -2;
1625 }
1626
1627
1628 const Hit& lasthit = evt_layer_hits[buildtrack.getLastFoundHitLyr()][buildtrack.getLastFoundHitIdx()];
1629 xhit_build_eff_ = lasthit.x();
1630 yhit_build_eff_ = lasthit.y();
1631 zhit_build_eff_ = lasthit.z();
1632
1633 pt_build_eff_ = buildtrack.pT();
1634 ept_build_eff_ = buildtrack.epT();
1635 phi_build_eff_ = buildtrack.momPhi();
1636 ephi_build_eff_ = buildtrack.emomPhi();
1637 eta_build_eff_ = buildtrack.momEta();
1638 eeta_build_eff_ = buildtrack.emomEta();
1639
1640 nHits_build_eff_ = buildtrack.nFoundHits();
1641 nLayers_build_eff_ = buildtrack.nUniqueLayers();
1642 nHitsMatched_build_eff_ = buildextra.nHitsMatched();
1643 fracHitsMatched_build_eff_ = buildextra.fracHitsMatched();
1644 lastlyr_build_eff_ = buildtrack.getLastFoundHitLyr();
1645
1646
1647 dphi_build_eff_ = buildextra.dPhi();
1648
1649
1650 hitchi2_build_eff_ = buildtrack.chi2();
1651 score_build_eff_ = buildtrack.score();
1652
1653 duplmask_build_eff_ = buildextra.isDuplicate();
1654 nTkMatches_build_eff_ = simToBuildMap_[mcID_eff_].size();
1655
1656
1657 if (Config::keepHitInfo)
1658 TTreeValidation::fillFullHitInfo(ev,
1659 buildtrack,
1660 hitlyrs_build_eff_,
1661 hitidxs_build_eff_,
1662 hitmcTkIDs_build_eff_,
1663 hitxs_build_eff_,
1664 hitys_build_eff_,
1665 hitzs_build_eff_);
1666 } else
1667 {
1668 mcmask_build_eff_ = (simtrack.isFindable() ? 0 : -1);
1669
1670 seedID_build_eff_ = -99;
1671
1672 pt_mc_build_eff_ = -99;
1673 phi_mc_build_eff_ = -99;
1674 eta_mc_build_eff_ = -99;
1675 helixchi2_build_eff_ = -99;
1676
1677 mcTSmask_build_eff_ = -1;
1678
1679 xhit_build_eff_ = -2000;
1680 yhit_build_eff_ = -2000;
1681 zhit_build_eff_ = -2000;
1682
1683 pt_build_eff_ = -99;
1684 ept_build_eff_ = -99;
1685 phi_build_eff_ = -99;
1686 ephi_build_eff_ = -99;
1687 eta_build_eff_ = -99;
1688 eeta_build_eff_ = -99;
1689
1690 nHits_build_eff_ = -99;
1691 nLayers_build_eff_ = -99;
1692 nHitsMatched_build_eff_ = -99;
1693 fracHitsMatched_build_eff_ = -99;
1694 lastlyr_build_eff_ = -99;
1695
1696 dphi_build_eff_ = -99;
1697
1698 hitchi2_build_eff_ = -99;
1699 score_build_eff_ = -17000;
1700
1701 duplmask_build_eff_ = -1;
1702 nTkMatches_build_eff_ = -99;
1703 }
1704
1705
1706 if (simToFitMap_.count(mcID_eff_) &&
1707 simtrack
1708 .isFindable())
1709 {
1710 for (unsigned int ii = 0; ii < simToFitMap_[mcID_eff_].size(); ii++) {
1711 const int theAlgo = evt_fit_tracks[simToFitMap_[mcID_eff_][ii]].algoint();
1712 if ((itermask_fit_eff_ >> theAlgo) & 1)
1713 iterduplmask_fit_eff_ = (iterduplmask_fit_eff_ | (1 << theAlgo));
1714 itermask_fit_eff_ = (itermask_fit_eff_ | (1 << theAlgo));
1715 }
1716 const auto& fittrack =
1717 evt_fit_tracks[simToFitMap_[mcID_eff_][0]];
1718 const auto& fitextra = evt_fit_extras[fittrack.label()];
1719 mcmask_fit_eff_ = 1;
1720
1721 seedID_fit_eff_ = fitextra.seedID();
1722
1723
1724 const int mcHitID =
1725 TTreeValidation::getLastFoundHit(fittrack.getLastFoundMCHitID(evt_layer_hits), mcID_eff_, ev);
1726 if (mcHitID >= 0 && Config::readSimTrackStates) {
1727 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
1728
1729 pt_mc_fit_eff_ = initLayTS.pT();
1730 phi_mc_fit_eff_ = initLayTS.momPhi();
1731 eta_mc_fit_eff_ = initLayTS.momEta();
1732 helixchi2_fit_eff_ = computeHelixChi2(initLayTS.parameters, fittrack.parameters(), fittrack.errors());
1733
1734 mcTSmask_fit_eff_ = 1;
1735 } else if (Config::tryToSaveSimInfo)
1736 {
1737
1738 pt_mc_fit_eff_ = pt_mc_gen_eff_;
1739 phi_mc_fit_eff_ = phi_mc_gen_eff_;
1740 eta_mc_fit_eff_ = eta_mc_gen_eff_;
1741 helixchi2_fit_eff_ = computeHelixChi2(simtrack.parameters(), fittrack.parameters(), fittrack.errors());
1742
1743 mcTSmask_fit_eff_ = 0;
1744 } else {
1745 pt_mc_fit_eff_ = -101;
1746 phi_mc_fit_eff_ = -101;
1747 eta_mc_fit_eff_ = -101;
1748 helixchi2_fit_eff_ = -101;
1749
1750 mcTSmask_fit_eff_ = -2;
1751 }
1752
1753
1754 const Hit& lasthit = evt_layer_hits[fittrack.getLastFoundHitLyr()][fittrack.getLastFoundHitIdx()];
1755 xhit_fit_eff_ = lasthit.x();
1756 yhit_fit_eff_ = lasthit.y();
1757 zhit_fit_eff_ = lasthit.z();
1758
1759 pt_fit_eff_ = fittrack.pT();
1760 ept_fit_eff_ = fittrack.epT();
1761 phi_fit_eff_ = fittrack.momPhi();
1762 ephi_fit_eff_ = fittrack.emomPhi();
1763 eta_fit_eff_ = fittrack.momEta();
1764 eeta_fit_eff_ = fittrack.emomEta();
1765
1766
1767 nHits_fit_eff_ = fittrack.nFoundHits();
1768 nLayers_fit_eff_ = fittrack.nUniqueLayers();
1769 nHitsMatched_fit_eff_ = fitextra.nHitsMatched();
1770 fracHitsMatched_fit_eff_ = fitextra.fracHitsMatched();
1771 lastlyr_fit_eff_ = fittrack.getLastFoundHitLyr();
1772
1773
1774 dphi_fit_eff_ = fitextra.dPhi();
1775
1776
1777 hitchi2_fit_eff_ = fittrack.chi2();
1778 score_fit_eff_ = fittrack.score();
1779
1780 duplmask_fit_eff_ = fitextra.isDuplicate();
1781 nTkMatches_fit_eff_ = simToFitMap_[mcID_eff_].size();
1782
1783
1784 if (Config::keepHitInfo)
1785 TTreeValidation::fillFullHitInfo(ev,
1786 fittrack,
1787 hitlyrs_fit_eff_,
1788 hitidxs_fit_eff_,
1789 hitmcTkIDs_fit_eff_,
1790 hitxs_fit_eff_,
1791 hitys_fit_eff_,
1792 hitzs_fit_eff_);
1793 } else
1794 {
1795 mcmask_fit_eff_ = (simtrack.isFindable() ? 0 : -1);
1796
1797 seedID_fit_eff_ = -99;
1798
1799 pt_mc_fit_eff_ = -99;
1800 phi_mc_fit_eff_ = -99;
1801 eta_mc_fit_eff_ = -99;
1802 helixchi2_fit_eff_ = -99;
1803
1804 mcTSmask_fit_eff_ = -1;
1805
1806 xhit_fit_eff_ = -2000;
1807 yhit_fit_eff_ = -2000;
1808 zhit_fit_eff_ = -2000;
1809
1810 pt_fit_eff_ = -99;
1811 ept_fit_eff_ = -99;
1812 phi_fit_eff_ = -99;
1813 ephi_fit_eff_ = -99;
1814 eta_fit_eff_ = -99;
1815 eeta_fit_eff_ = -99;
1816
1817 nHits_fit_eff_ = -99;
1818 nLayers_fit_eff_ = -99;
1819 nHitsMatched_fit_eff_ = -99;
1820 fracHitsMatched_fit_eff_ = -99;
1821 lastlyr_fit_eff_ = -99;
1822
1823 dphi_fit_eff_ = -99;
1824
1825 hitchi2_fit_eff_ = -99;
1826 score_fit_eff_ = -17000;
1827
1828 duplmask_fit_eff_ = -1;
1829 nTkMatches_fit_eff_ = -99;
1830 }
1831
1832 efftree_->Fill();
1833 }
1834 }
1835
1836 void TTreeValidation::fillFakeRateTree(const Event& ev) {
1837 std::lock_guard<std::mutex> locker(glock_);
1838
1839 const auto ievt = ev.evtID();
1840 const auto& evt_sim_tracks =
1841 ev.simTracks_;
1842 const auto& evt_seed_tracks = ev.seedTracks_;
1843 const auto& evt_seed_extras = ev.seedTracksExtra_;
1844 const auto& evt_build_tracks = ev.candidateTracks_;
1845 const auto& evt_build_extras = ev.candidateTracksExtra_;
1846 const auto& evt_fit_tracks = ev.fitTracks_;
1847 const auto& evt_fit_extras = ev.fitTracksExtra_;
1848 const auto& evt_layer_hits = ev.layerHits_;
1849 const auto& evt_sim_trackstates = ev.simTrackStates_;
1850
1851 for (const auto& seedtrack : evt_seed_tracks) {
1852 if (Config::keepHitInfo) {
1853 hitlyrs_seed_FR_.clear();
1854 hitlyrs_mc_seed_FR_.clear();
1855 hitlyrs_build_FR_.clear();
1856 hitlyrs_mc_build_FR_.clear();
1857 hitlyrs_fit_FR_.clear();
1858 hitlyrs_mc_fit_FR_.clear();
1859
1860 hitidxs_seed_FR_.clear();
1861 hitidxs_mc_seed_FR_.clear();
1862 hitidxs_build_FR_.clear();
1863 hitidxs_mc_build_FR_.clear();
1864 hitidxs_fit_FR_.clear();
1865 hitidxs_mc_fit_FR_.clear();
1866
1867 hitmcTkIDs_seed_FR_.clear();
1868 hitmcTkIDs_mc_seed_FR_.clear();
1869 hitmcTkIDs_build_FR_.clear();
1870 hitmcTkIDs_mc_build_FR_.clear();
1871 hitmcTkIDs_fit_FR_.clear();
1872 hitmcTkIDs_mc_fit_FR_.clear();
1873
1874 hitxs_seed_FR_.clear();
1875 hitxs_mc_seed_FR_.clear();
1876 hitxs_build_FR_.clear();
1877 hitxs_mc_build_FR_.clear();
1878 hitxs_fit_FR_.clear();
1879 hitxs_mc_fit_FR_.clear();
1880
1881 hitys_seed_FR_.clear();
1882 hitys_mc_seed_FR_.clear();
1883 hitys_build_FR_.clear();
1884 hitys_mc_build_FR_.clear();
1885 hitys_fit_FR_.clear();
1886 hitys_mc_fit_FR_.clear();
1887
1888 hitzs_seed_FR_.clear();
1889 hitzs_mc_seed_FR_.clear();
1890 hitzs_build_FR_.clear();
1891 hitzs_mc_build_FR_.clear();
1892 hitzs_fit_FR_.clear();
1893 hitzs_mc_fit_FR_.clear();
1894 }
1895
1896 evtID_FR_ = ievt;
1897
1898
1899 const auto& seedextra = evt_seed_extras[seedtrack.label()];
1900 seedID_FR_ = seedextra.seedID();
1901 seedmask_seed_FR_ =
1902 1;
1903
1904
1905
1906 xhit_seed_FR_ = 0;
1907 yhit_seed_FR_ = 0;
1908 zhit_seed_FR_ = 0;
1909
1910 pt_seed_FR_ = seedtrack.pT();
1911 ept_seed_FR_ = seedtrack.epT();
1912 phi_seed_FR_ = seedtrack.momPhi();
1913 ephi_seed_FR_ = seedtrack.emomPhi();
1914 eta_seed_FR_ = seedtrack.momEta();
1915 eeta_seed_FR_ = seedtrack.emomEta();
1916
1917 nHits_seed_FR_ = seedtrack.nFoundHits();
1918 nLayers_seed_FR_ = seedtrack.nUniqueLayers();
1919 nHitsMatched_seed_FR_ = seedextra.nHitsMatched();
1920 fracHitsMatched_seed_FR_ = seedextra.fracHitsMatched();
1921 lastlyr_seed_FR_ = seedtrack.getLastFoundHitLyr();
1922
1923 algorithm_FR_ = seedtrack.algoint();
1924
1925
1926 dphi_seed_FR_ = seedextra.dPhi();
1927
1928
1929 hitchi2_seed_FR_ = seedtrack.chi2();
1930 score_seed_FR_ = seedtrack.score();
1931
1932 if (Config::keepHitInfo)
1933 TTreeValidation::fillFullHitInfo(ev,
1934 seedtrack,
1935 hitlyrs_seed_FR_,
1936 hitidxs_seed_FR_,
1937 hitmcTkIDs_seed_FR_,
1938 hitxs_seed_FR_,
1939 hitys_seed_FR_,
1940 hitzs_seed_FR_);
1941
1942
1943 mcID_seed_FR_ = seedextra.mcTrackID();
1944 mcmask_seed_FR_ = TTreeValidation::getMaskAssignment(mcID_seed_FR_);
1945
1946 if (mcmask_seed_FR_ == 1)
1947 {
1948 const auto& simtrack = evt_sim_tracks[mcID_seed_FR_];
1949
1950 const int mcHitID =
1951 TTreeValidation::getLastFoundHit(seedtrack.getLastFoundMCHitID(evt_layer_hits), mcID_seed_FR_, ev);
1952 if (mcHitID >= 0 && Config::readSimTrackStates) {
1953 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
1954 pt_mc_seed_FR_ = initLayTS.pT();
1955 phi_mc_seed_FR_ = initLayTS.momPhi();
1956 eta_mc_seed_FR_ = initLayTS.momEta();
1957 helixchi2_seed_FR_ = computeHelixChi2(initLayTS.parameters, seedtrack.parameters(), seedtrack.errors());
1958
1959 mcTSmask_seed_FR_ = 1;
1960 } else if (Config::tryToSaveSimInfo) {
1961 pt_mc_seed_FR_ = simtrack.pT();
1962 phi_mc_seed_FR_ = simtrack.momPhi();
1963 eta_mc_seed_FR_ = simtrack.momEta();
1964 helixchi2_seed_FR_ = computeHelixChi2(simtrack.parameters(), seedtrack.parameters(), seedtrack.errors());
1965
1966 mcTSmask_seed_FR_ = 0;
1967 } else {
1968 pt_mc_seed_FR_ = -101;
1969 phi_mc_seed_FR_ = -101;
1970 eta_mc_seed_FR_ = -101;
1971 helixchi2_seed_FR_ = -101;
1972
1973 mcTSmask_seed_FR_ = -2;
1974 }
1975
1976 nHits_mc_seed_FR_ = simtrack.nFoundHits();
1977 nLayers_mc_seed_FR_ = simtrack.nUniqueLayers();
1978 lastlyr_mc_seed_FR_ = simtrack.getLastFoundHitLyr();
1979
1980 duplmask_seed_FR_ = seedextra.isDuplicate();
1981 iTkMatches_seed_FR_ =
1982 seedextra
1983 .duplicateID();
1984
1985 if (Config::keepHitInfo)
1986 TTreeValidation::fillFullHitInfo(ev,
1987 simtrack,
1988 hitlyrs_mc_seed_FR_,
1989 hitidxs_mc_seed_FR_,
1990 hitmcTkIDs_mc_seed_FR_,
1991 hitxs_mc_seed_FR_,
1992 hitys_mc_seed_FR_,
1993 hitzs_mc_seed_FR_);
1994 } else {
1995
1996 pt_mc_seed_FR_ = -99;
1997 phi_mc_seed_FR_ = -99;
1998 eta_mc_seed_FR_ = -99;
1999 helixchi2_seed_FR_ = -99;
2000
2001 mcTSmask_seed_FR_ = -1;
2002
2003 nHits_mc_seed_FR_ = -99;
2004 nLayers_mc_seed_FR_ = -99;
2005 lastlyr_mc_seed_FR_ = -99;
2006
2007 duplmask_seed_FR_ = -1;
2008 iTkMatches_seed_FR_ = -99;
2009 }
2010
2011
2012
2013
2014 if (seedToBuildMap_.count(seedID_FR_)) {
2015 seedmask_build_FR_ = 1;
2016
2017 const auto& buildtrack = evt_build_tracks[seedToBuildMap_[seedID_FR_]];
2018 const auto& buildextra = evt_build_extras[buildtrack.label()];
2019
2020
2021 const Hit& lasthit = evt_layer_hits[buildtrack.getLastFoundHitLyr()][buildtrack.getLastFoundHitIdx()];
2022 xhit_build_FR_ = lasthit.x();
2023 yhit_build_FR_ = lasthit.y();
2024 zhit_build_FR_ = lasthit.z();
2025
2026 pt_build_FR_ = buildtrack.pT();
2027 ept_build_FR_ = buildtrack.epT();
2028 phi_build_FR_ = buildtrack.momPhi();
2029 ephi_build_FR_ = buildtrack.emomPhi();
2030 eta_build_FR_ = buildtrack.momEta();
2031 eeta_build_FR_ = buildtrack.emomEta();
2032
2033 nHits_build_FR_ = buildtrack.nFoundHits();
2034 nLayers_build_FR_ = buildtrack.nUniqueLayers();
2035 nHitsMatched_build_FR_ = buildextra.nHitsMatched();
2036 fracHitsMatched_build_FR_ = buildextra.fracHitsMatched();
2037 lastlyr_build_FR_ = buildtrack.getLastFoundHitLyr();
2038
2039
2040 dphi_build_FR_ = buildextra.dPhi();
2041
2042
2043 hitchi2_build_FR_ = buildtrack.chi2();
2044 score_build_FR_ = buildtrack.score();
2045
2046 if (Config::keepHitInfo)
2047 TTreeValidation::fillFullHitInfo(ev,
2048 buildtrack,
2049 hitlyrs_build_FR_,
2050 hitidxs_build_FR_,
2051 hitmcTkIDs_build_FR_,
2052 hitxs_build_FR_,
2053 hitys_build_FR_,
2054 hitzs_build_FR_);
2055
2056
2057 mcID_build_FR_ = buildextra.mcTrackID();
2058 mcmask_build_FR_ = TTreeValidation::getMaskAssignment(mcID_build_FR_);
2059
2060 if (mcmask_build_FR_ == 1)
2061 {
2062 const auto& simtrack = evt_sim_tracks[mcID_build_FR_];
2063
2064 const int mcHitID =
2065 TTreeValidation::getLastFoundHit(buildtrack.getLastFoundMCHitID(evt_layer_hits), mcID_build_FR_, ev);
2066 if (mcHitID >= 0 && Config::readSimTrackStates) {
2067 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
2068 pt_mc_build_FR_ = initLayTS.pT();
2069 phi_mc_build_FR_ = initLayTS.momPhi();
2070 eta_mc_build_FR_ = initLayTS.momEta();
2071 helixchi2_build_FR_ = computeHelixChi2(initLayTS.parameters, buildtrack.parameters(), buildtrack.errors());
2072
2073 mcTSmask_build_FR_ = 1;
2074 } else if (Config::tryToSaveSimInfo) {
2075 pt_mc_build_FR_ = simtrack.pT();
2076 phi_mc_build_FR_ = simtrack.momPhi();
2077 eta_mc_build_FR_ = simtrack.momEta();
2078 helixchi2_build_FR_ = computeHelixChi2(simtrack.parameters(), buildtrack.parameters(), buildtrack.errors());
2079
2080 mcTSmask_build_FR_ = 0;
2081 } else {
2082 pt_mc_build_FR_ = -101;
2083 phi_mc_build_FR_ = -101;
2084 eta_mc_build_FR_ = -101;
2085 helixchi2_build_FR_ = -101;
2086
2087 mcTSmask_build_FR_ = -2;
2088 }
2089
2090 nHits_mc_build_FR_ = simtrack.nFoundHits();
2091 nLayers_mc_build_FR_ = simtrack.nUniqueLayers();
2092 lastlyr_mc_build_FR_ = simtrack.getLastFoundHitLyr();
2093
2094 duplmask_build_FR_ = buildextra.isDuplicate();
2095 iTkMatches_build_FR_ =
2096 buildextra
2097 .duplicateID();
2098
2099 if (Config::keepHitInfo)
2100 TTreeValidation::fillFullHitInfo(ev,
2101 simtrack,
2102 hitlyrs_mc_build_FR_,
2103 hitidxs_mc_build_FR_,
2104 hitmcTkIDs_mc_build_FR_,
2105 hitxs_mc_build_FR_,
2106 hitys_mc_build_FR_,
2107 hitzs_mc_build_FR_);
2108 } else
2109 {
2110
2111 pt_mc_build_FR_ = -99;
2112 phi_mc_build_FR_ = -99;
2113 eta_mc_build_FR_ = -99;
2114 helixchi2_build_FR_ = -99;
2115
2116 mcTSmask_build_FR_ = -1;
2117
2118 nHits_mc_build_FR_ = -99;
2119 nLayers_mc_build_FR_ = -99;
2120 lastlyr_mc_build_FR_ = -99;
2121
2122 duplmask_build_FR_ = -1;
2123 iTkMatches_build_FR_ = -99;
2124 }
2125 }
2126
2127 else
2128 {
2129 seedmask_build_FR_ = 0;
2130
2131
2132 xhit_build_FR_ = -3000;
2133 yhit_build_FR_ = -3000;
2134 zhit_build_FR_ = -3000;
2135
2136
2137 pt_build_FR_ = -100;
2138 ept_build_FR_ = -100;
2139 phi_build_FR_ = -100;
2140 ephi_build_FR_ = -100;
2141 eta_build_FR_ = -100;
2142 eeta_build_FR_ = -100;
2143
2144 nHits_build_FR_ = -100;
2145 nLayers_build_FR_ = -100;
2146 nHitsMatched_build_FR_ = -100;
2147 fracHitsMatched_build_FR_ = -100;
2148 lastlyr_build_FR_ = -100;
2149
2150 dphi_build_FR_ = -100;
2151
2152 hitchi2_build_FR_ = -100;
2153 score_build_FR_ = -5001;
2154
2155
2156 mcmask_build_FR_ = -2;
2157 mcID_build_FR_ = -100;
2158
2159 pt_mc_build_FR_ = -100;
2160 phi_mc_build_FR_ = -100;
2161 eta_mc_build_FR_ = -100;
2162 helixchi2_build_FR_ = -100;
2163
2164 mcTSmask_build_FR_ = -3;
2165
2166 nHits_mc_build_FR_ = -100;
2167 nLayers_mc_build_FR_ = -100;
2168 lastlyr_mc_build_FR_ = -100;
2169
2170 duplmask_build_FR_ = -2;
2171 iTkMatches_build_FR_ = -100;
2172 }
2173
2174
2175 if (seedToFitMap_.count(seedID_FR_)) {
2176 seedmask_fit_FR_ = 1;
2177
2178 const auto& fittrack = evt_fit_tracks[seedToFitMap_[seedID_FR_]];
2179 const auto& fitextra = evt_fit_extras[fittrack.label()];
2180
2181
2182 const Hit& lasthit = evt_layer_hits[fittrack.getLastFoundHitLyr()][fittrack.getLastFoundHitIdx()];
2183 xhit_fit_FR_ = lasthit.x();
2184 yhit_fit_FR_ = lasthit.y();
2185 zhit_fit_FR_ = lasthit.z();
2186
2187 pt_fit_FR_ = fittrack.pT();
2188 ept_fit_FR_ = fittrack.epT();
2189 phi_fit_FR_ = fittrack.momPhi();
2190 ephi_fit_FR_ = fittrack.emomPhi();
2191 eta_fit_FR_ = fittrack.momEta();
2192 eeta_fit_FR_ = fittrack.emomEta();
2193
2194 nHits_fit_FR_ = fittrack.nFoundHits();
2195 nLayers_fit_FR_ = fittrack.nUniqueLayers();
2196 nHitsMatched_fit_FR_ = fitextra.nHitsMatched();
2197 fracHitsMatched_fit_FR_ = fitextra.fracHitsMatched();
2198 lastlyr_fit_FR_ = fittrack.getLastFoundHitLyr();
2199
2200
2201 dphi_fit_FR_ = fitextra.dPhi();
2202
2203
2204 hitchi2_fit_FR_ = fittrack.chi2();
2205 score_fit_FR_ = fittrack.score();
2206
2207 if (Config::keepHitInfo)
2208 TTreeValidation::fillFullHitInfo(ev,
2209 fittrack,
2210 hitlyrs_fit_FR_,
2211 hitidxs_fit_FR_,
2212 hitmcTkIDs_fit_FR_,
2213 hitxs_fit_FR_,
2214 hitys_fit_FR_,
2215 hitzs_fit_FR_);
2216
2217
2218 mcID_fit_FR_ = fitextra.mcTrackID();
2219 mcmask_fit_FR_ = TTreeValidation::getMaskAssignment(mcID_fit_FR_);
2220
2221 if (mcmask_fit_FR_ == 1)
2222 {
2223 const auto& simtrack = evt_sim_tracks[mcID_fit_FR_];
2224
2225 const int mcHitID = TTreeValidation::getLastFoundHit(
2226 fittrack.getLastFoundMCHitID(evt_layer_hits), mcID_fit_FR_, ev);
2227 if (mcHitID >= 0 && Config::readSimTrackStates) {
2228 const TrackState& initLayTS = evt_sim_trackstates[mcHitID];
2229 pt_mc_fit_FR_ = initLayTS.pT();
2230 phi_mc_fit_FR_ = initLayTS.momPhi();
2231 eta_mc_fit_FR_ = initLayTS.momEta();
2232 helixchi2_fit_FR_ = computeHelixChi2(initLayTS.parameters, fittrack.parameters(), fittrack.errors());
2233
2234 mcTSmask_fit_FR_ = 1;
2235 } else if (Config::tryToSaveSimInfo) {
2236 pt_mc_fit_FR_ = simtrack.pT();
2237 phi_mc_fit_FR_ = simtrack.momPhi();
2238 eta_mc_fit_FR_ = simtrack.momEta();
2239 helixchi2_fit_FR_ = computeHelixChi2(simtrack.parameters(), fittrack.parameters(), fittrack.errors());
2240
2241 mcTSmask_fit_FR_ = 0;
2242 } else {
2243 pt_mc_fit_FR_ = -101;
2244 phi_mc_fit_FR_ = -101;
2245 eta_mc_fit_FR_ = -101;
2246 helixchi2_fit_FR_ = -101;
2247
2248 mcTSmask_fit_FR_ = -2;
2249 }
2250
2251 nHits_mc_fit_FR_ = simtrack.nFoundHits();
2252 nLayers_mc_fit_FR_ = simtrack.nUniqueLayers();
2253 lastlyr_mc_fit_FR_ = simtrack.getLastFoundHitLyr();
2254
2255 duplmask_fit_FR_ = fitextra.isDuplicate();
2256 iTkMatches_fit_FR_ =
2257 fitextra
2258 .duplicateID();
2259
2260 if (Config::keepHitInfo)
2261 TTreeValidation::fillFullHitInfo(ev,
2262 simtrack,
2263 hitlyrs_mc_fit_FR_,
2264 hitidxs_mc_fit_FR_,
2265 hitmcTkIDs_mc_fit_FR_,
2266 hitxs_mc_fit_FR_,
2267 hitys_mc_fit_FR_,
2268 hitzs_mc_fit_FR_);
2269 } else
2270 {
2271
2272 pt_mc_fit_FR_ = -99;
2273 phi_mc_fit_FR_ = -99;
2274 eta_mc_fit_FR_ = -99;
2275 helixchi2_fit_FR_ = -99;
2276
2277 mcTSmask_fit_FR_ = -1;
2278
2279 nHits_mc_fit_FR_ = -99;
2280 nLayers_mc_fit_FR_ = -99;
2281 lastlyr_mc_fit_FR_ = -99;
2282
2283 duplmask_fit_FR_ = -1;
2284 iTkMatches_fit_FR_ = -99;
2285 }
2286 }
2287
2288 else
2289 {
2290 seedmask_fit_FR_ = 0;
2291
2292
2293 xhit_fit_FR_ = -3000;
2294 yhit_fit_FR_ = -3000;
2295 zhit_fit_FR_ = -3000;
2296
2297
2298 pt_fit_FR_ = -100;
2299 ept_fit_FR_ = -100;
2300 phi_fit_FR_ = -100;
2301 ephi_fit_FR_ = -100;
2302 eta_fit_FR_ = -100;
2303 eeta_fit_FR_ = -100;
2304
2305 nHits_fit_FR_ = -100;
2306 nLayers_fit_FR_ = -100;
2307 nHitsMatched_fit_FR_ = -100;
2308 fracHitsMatched_fit_FR_ = -100;
2309 lastlyr_fit_FR_ = -100;
2310
2311 dphi_fit_FR_ = -100;
2312
2313 hitchi2_fit_FR_ = -100;
2314 score_fit_FR_ = -5001;
2315
2316
2317 mcmask_fit_FR_ = -2;
2318 mcID_fit_FR_ = -100;
2319
2320 pt_mc_fit_FR_ = -100;
2321 phi_mc_fit_FR_ = -100;
2322 eta_mc_fit_FR_ = -100;
2323 helixchi2_fit_FR_ = -100;
2324
2325 mcTSmask_fit_FR_ = -3;
2326
2327 nHits_mc_fit_FR_ = -100;
2328 nLayers_mc_fit_FR_ = -100;
2329 lastlyr_mc_fit_FR_ = -100;
2330
2331 duplmask_fit_FR_ = -2;
2332 iTkMatches_fit_FR_ = -100;
2333 }
2334
2335 frtree_->Fill();
2336 }
2337 }
2338
2339 void TTreeValidation::fillConfigTree() {
2340 std::lock_guard<std::mutex> locker(glock_);
2341
2342 Ntracks_ = Config::nTracks;
2343 Nevents_ = Config::nEvents;
2344
2345 nLayers_ = Config::nLayers;
2346
2347 nlayers_per_seed_ = Config::ItrInfo[0].m_params.nlayers_per_seed;
2348 maxCand_ = Config::ItrInfo[0].m_params.maxCandsPerSeed;
2349 chi2Cut_min_ = Config::ItrInfo[0].m_params.chi2Cut_min;
2350 nSigma_ = Config::nSigma;
2351 minDPhi_ = Config::minDPhi;
2352 maxDPhi_ = Config::maxDPhi;
2353 minDEta_ = Config::minDEta;
2354 maxDEta_ = Config::maxDEta;
2355
2356 beamspotX_ = Config::beamspotX;
2357 beamspotY_ = Config::beamspotY;
2358 beamspotZ_ = Config::beamspotZ;
2359
2360 minSimPt_ = Config::minSimPt;
2361 maxSimPt_ = Config::maxSimPt;
2362
2363 hitposerrXY_ = Config::hitposerrXY;
2364 hitposerrZ_ = Config::hitposerrZ;
2365 hitposerrR_ = Config::hitposerrR;
2366 varXY_ = Config::varXY;
2367 varZ_ = Config::varZ;
2368
2369 ptinverr049_ = Config::ptinverr049;
2370 phierr049_ = Config::phierr049;
2371 thetaerr049_ = Config::thetaerr049;
2372 ptinverr012_ = Config::ptinverr012;
2373 phierr012_ = Config::phierr012;
2374 thetaerr012_ = Config::thetaerr012;
2375
2376 configtree_->Fill();
2377 }
2378
2379 void TTreeValidation::fillCMSSWEfficiencyTree(const Event& ev) {
2380 std::lock_guard<std::mutex> locker(glock_);
2381
2382 const auto ievt = ev.evtID();
2383 const auto& evt_sim_tracks = ev.simTracks_;
2384 const auto& evt_cmssw_tracks = ev.cmsswTracks_;
2385 const auto& evt_cmssw_extras = ev.cmsswTracksExtra_;
2386 const auto& evt_build_tracks = ev.candidateTracks_;
2387 const auto& evt_build_extras = ev.candidateTracksExtra_;
2388 const auto& evt_fit_tracks = ev.fitTracks_;
2389 const auto& evt_fit_extras = ev.fitTracksExtra_;
2390 const auto& evt_layer_hits = ev.layerHits_;
2391
2392 for (const auto& cmsswtrack : evt_cmssw_tracks) {
2393
2394 if (Config::keepHitInfo) {
2395 hitlyrs_cmssw_ceff_.clear();
2396 hitlyrs_build_ceff_.clear();
2397 hitlyrs_mc_build_ceff_.clear();
2398 hitlyrs_fit_ceff_.clear();
2399 hitlyrs_mc_fit_ceff_.clear();
2400
2401 hitidxs_cmssw_ceff_.clear();
2402 hitidxs_build_ceff_.clear();
2403 hitidxs_mc_build_ceff_.clear();
2404 hitidxs_fit_ceff_.clear();
2405 hitidxs_mc_fit_ceff_.clear();
2406 }
2407
2408 const auto& cmsswextra = evt_cmssw_extras[cmsswtrack.label()];
2409
2410 evtID_ceff_ = ievt;
2411 cmsswID_ceff_ = cmsswtrack.label();
2412 seedID_cmssw_ceff_ = cmsswextra.seedID();
2413
2414
2415 x_cmssw_ceff_ = cmsswtrack.x();
2416 y_cmssw_ceff_ = cmsswtrack.y();
2417 z_cmssw_ceff_ = cmsswtrack.z();
2418
2419 pt_cmssw_ceff_ = cmsswtrack.pT();
2420 phi_cmssw_ceff_ = cmsswtrack.momPhi();
2421 eta_cmssw_ceff_ = cmsswtrack.momEta();
2422
2423 nHits_cmssw_ceff_ = cmsswtrack.nFoundHits();
2424 nLayers_cmssw_ceff_ = cmsswtrack.nUniqueLayers();
2425 lastlyr_cmssw_ceff_ = cmsswtrack.getLastFoundHitLyr();
2426
2427 itermask_build_ceff_ = 0;
2428 itermask_fit_ceff_ = 0;
2429 iterduplmask_build_ceff_ = 0;
2430 iterduplmask_fit_ceff_ = 0;
2431 algo_seed_ceff_ = 0;
2432
2433 for (auto aa : cmsswextra.seedAlgos())
2434 algo_seed_ceff_ = (algo_seed_ceff_ | (1 << aa));
2435
2436 if (Config::keepHitInfo)
2437 TTreeValidation::fillMinHitInfo(cmsswtrack, hitlyrs_cmssw_ceff_, hitidxs_cmssw_ceff_);
2438
2439
2440 if (cmsswToBuildMap_.count(cmsswID_ceff_) &&
2441 cmsswtrack
2442 .isFindable())
2443 {
2444 for (unsigned int ii = 0; ii < cmsswToBuildMap_[cmsswID_ceff_].size(); ii++) {
2445 const int theAlgo = evt_build_tracks[cmsswToBuildMap_[cmsswID_ceff_][ii]].algoint();
2446 if ((itermask_build_ceff_ >> theAlgo) & 1)
2447 iterduplmask_build_ceff_ = (iterduplmask_build_ceff_ | (1 << theAlgo));
2448 itermask_build_ceff_ = (itermask_build_ceff_ | (1 << theAlgo));
2449 }
2450
2451 const auto& buildtrack =
2452 evt_build_tracks[cmsswToBuildMap_[cmsswID_ceff_][0]];
2453 const auto& buildextra =
2454 evt_build_extras[buildtrack.label()];
2455 cmsswmask_build_ceff_ = 1;
2456
2457 seedID_build_ceff_ = buildextra.seedID();
2458 mcTrackID_build_ceff_ = buildextra.mcTrackID();
2459
2460
2461 pt_build_ceff_ = buildtrack.pT();
2462 ept_build_ceff_ = buildtrack.epT();
2463 phi_build_ceff_ = buildtrack.momPhi();
2464 ephi_build_ceff_ = buildtrack.emomPhi();
2465 eta_build_ceff_ = buildtrack.momEta();
2466 eeta_build_ceff_ = buildtrack.emomEta();
2467
2468
2469 if (mcTrackID_build_ceff_ >= 0) {
2470 const auto& simtrack = evt_sim_tracks[mcTrackID_build_ceff_];
2471 x_mc_build_ceff_ = simtrack.x();
2472 y_mc_build_ceff_ = simtrack.y();
2473 z_mc_build_ceff_ = simtrack.z();
2474 pt_mc_build_ceff_ = simtrack.pT();
2475 phi_mc_build_ceff_ = simtrack.momPhi();
2476 eta_mc_build_ceff_ = simtrack.momEta();
2477
2478 if (Config::keepHitInfo)
2479 TTreeValidation::fillMinHitInfo(simtrack, hitlyrs_mc_build_ceff_, hitidxs_mc_build_ceff_);
2480 } else {
2481 x_mc_build_ceff_ = -1000;
2482 y_mc_build_ceff_ = -1000;
2483 z_mc_build_ceff_ = -1000;
2484 pt_mc_build_ceff_ = -99;
2485 phi_mc_build_ceff_ = -99;
2486 eta_mc_build_ceff_ = -99;
2487 }
2488
2489
2490 nHits_build_ceff_ = buildtrack.nFoundHits();
2491 nLayers_build_ceff_ = buildtrack.nUniqueLayers();
2492 nHitsMatched_build_ceff_ = buildextra.nHitsMatched();
2493 fracHitsMatched_build_ceff_ = buildextra.fracHitsMatched();
2494 lastlyr_build_ceff_ = buildtrack.getLastFoundHitLyr();
2495
2496
2497 const Hit& lasthit = evt_layer_hits[buildtrack.getLastFoundHitLyr()][buildtrack.getLastFoundHitIdx()];
2498 xhit_build_ceff_ = lasthit.x();
2499 yhit_build_ceff_ = lasthit.y();
2500 zhit_build_ceff_ = lasthit.z();
2501
2502
2503 hitchi2_build_ceff_ = buildtrack.chi2();
2504 helixchi2_build_ceff_ = buildextra.helixChi2();
2505 score_build_ceff_ = buildtrack.score();
2506
2507
2508 dphi_build_ceff_ = buildextra.dPhi();
2509
2510
2511 duplmask_build_ceff_ = buildextra.isDuplicate();
2512 nTkMatches_build_ceff_ = cmsswToBuildMap_[cmsswID_ceff_].size();
2513
2514 if (Config::keepHitInfo)
2515 TTreeValidation::fillMinHitInfo(buildtrack, hitlyrs_build_ceff_, hitidxs_build_ceff_);
2516 } else
2517 {
2518 cmsswmask_build_ceff_ = (cmsswtrack.isFindable() ? 0 : -1);
2519
2520 seedID_build_ceff_ = -99;
2521 mcTrackID_build_ceff_ = -99;
2522
2523 pt_build_ceff_ = -99;
2524 ept_build_ceff_ = -99;
2525 phi_build_ceff_ = -99;
2526 ephi_build_ceff_ = -99;
2527 eta_build_ceff_ = -99;
2528 eeta_build_ceff_ = -99;
2529
2530 x_mc_build_ceff_ = -2000;
2531 y_mc_build_ceff_ = -2000;
2532 z_mc_build_ceff_ = -2000;
2533 pt_mc_build_ceff_ = -99;
2534 phi_mc_build_ceff_ = -99;
2535 eta_mc_build_ceff_ = -99;
2536
2537 nHits_build_ceff_ = -99;
2538 nLayers_build_ceff_ = -99;
2539 nHitsMatched_build_ceff_ = -99;
2540 fracHitsMatched_build_ceff_ = -99;
2541 lastlyr_build_ceff_ = -99;
2542
2543 xhit_build_ceff_ = -2000;
2544 yhit_build_ceff_ = -2000;
2545 zhit_build_ceff_ = -2000;
2546
2547 hitchi2_build_ceff_ = -99;
2548 helixchi2_build_ceff_ = -99;
2549 score_build_ceff_ = -17000;
2550
2551 dphi_build_ceff_ = -99;
2552
2553 duplmask_build_ceff_ = -1;
2554 nTkMatches_build_ceff_ = -99;
2555 }
2556
2557
2558 if (cmsswToFitMap_.count(cmsswID_ceff_) &&
2559 cmsswtrack
2560 .isFindable())
2561 {
2562 for (unsigned int ii = 0; ii < cmsswToFitMap_[cmsswID_ceff_].size(); ii++) {
2563 const int theAlgo = evt_build_tracks[cmsswToFitMap_[cmsswID_ceff_][ii]].algoint();
2564 if ((itermask_fit_ceff_ >> theAlgo) & 1)
2565 iterduplmask_fit_ceff_ = (iterduplmask_fit_ceff_ | (1 << theAlgo));
2566 itermask_fit_ceff_ = (itermask_fit_ceff_ | (1 << theAlgo));
2567 }
2568
2569 const auto& fittrack =
2570 evt_fit_tracks[cmsswToFitMap_[cmsswID_ceff_][0]];
2571 const auto& fitextra = evt_fit_extras[fittrack.label()];
2572 cmsswmask_fit_ceff_ = 1;
2573
2574 seedID_fit_ceff_ = fitextra.seedID();
2575 mcTrackID_fit_ceff_ = fitextra.mcTrackID();
2576
2577
2578 pt_fit_ceff_ = fittrack.pT();
2579 ept_fit_ceff_ = fittrack.epT();
2580 phi_fit_ceff_ = fittrack.momPhi();
2581 ephi_fit_ceff_ = fittrack.emomPhi();
2582 eta_fit_ceff_ = fittrack.momEta();
2583 eeta_fit_ceff_ = fittrack.emomEta();
2584
2585
2586 if (mcTrackID_fit_ceff_ >= 0) {
2587 const auto& simtrack = evt_sim_tracks[mcTrackID_fit_ceff_];
2588 x_mc_fit_ceff_ = simtrack.x();
2589 y_mc_fit_ceff_ = simtrack.y();
2590 z_mc_fit_ceff_ = simtrack.z();
2591 pt_mc_fit_ceff_ = simtrack.pT();
2592 phi_mc_fit_ceff_ = simtrack.momPhi();
2593 eta_mc_fit_ceff_ = simtrack.momEta();
2594
2595 if (Config::keepHitInfo)
2596 TTreeValidation::fillMinHitInfo(simtrack, hitlyrs_mc_fit_ceff_, hitidxs_mc_fit_ceff_);
2597 } else {
2598 x_mc_fit_ceff_ = -1000;
2599 y_mc_fit_ceff_ = -1000;
2600 z_mc_fit_ceff_ = -1000;
2601 pt_mc_fit_ceff_ = -99;
2602 phi_mc_fit_ceff_ = -99;
2603 eta_mc_fit_ceff_ = -99;
2604 }
2605
2606
2607 nHits_fit_ceff_ = fittrack.nFoundHits();
2608 nLayers_fit_ceff_ = fittrack.nUniqueLayers();
2609 nHitsMatched_fit_ceff_ = fitextra.nHitsMatched();
2610 fracHitsMatched_fit_ceff_ = fitextra.fracHitsMatched();
2611 lastlyr_fit_ceff_ = fittrack.getLastFoundHitLyr();
2612
2613
2614 const Hit& lasthit = evt_layer_hits[fittrack.getLastFoundHitLyr()][fittrack.getLastFoundHitIdx()];
2615 xhit_fit_ceff_ = lasthit.x();
2616 yhit_fit_ceff_ = lasthit.y();
2617 zhit_fit_ceff_ = lasthit.z();
2618
2619
2620 hitchi2_fit_ceff_ = fittrack.chi2();
2621 helixchi2_fit_ceff_ = fitextra.helixChi2();
2622 score_fit_ceff_ = fittrack.score();
2623
2624
2625 dphi_fit_ceff_ = fitextra.dPhi();
2626
2627
2628 duplmask_fit_ceff_ = fitextra.isDuplicate();
2629 nTkMatches_fit_ceff_ = cmsswToFitMap_[cmsswID_ceff_].size();
2630
2631 if (Config::keepHitInfo)
2632 TTreeValidation::fillMinHitInfo(fittrack, hitlyrs_fit_ceff_, hitidxs_fit_ceff_);
2633 } else
2634 {
2635 cmsswmask_fit_ceff_ = (cmsswtrack.isFindable() ? 0 : -1);
2636
2637 seedID_fit_ceff_ = -99;
2638 mcTrackID_fit_ceff_ = -99;
2639
2640 pt_fit_ceff_ = -99;
2641 ept_fit_ceff_ = -99;
2642 phi_fit_ceff_ = -99;
2643 ephi_fit_ceff_ = -99;
2644 eta_fit_ceff_ = -99;
2645 eeta_fit_ceff_ = -99;
2646
2647 x_mc_fit_ceff_ = -2000;
2648 y_mc_fit_ceff_ = -2000;
2649 z_mc_fit_ceff_ = -2000;
2650 pt_mc_fit_ceff_ = -99;
2651 phi_mc_fit_ceff_ = -99;
2652 eta_mc_fit_ceff_ = -99;
2653
2654 nHits_fit_ceff_ = -99;
2655 nLayers_fit_ceff_ = -99;
2656 nHitsMatched_fit_ceff_ = -99;
2657 fracHitsMatched_fit_ceff_ = -99;
2658 lastlyr_fit_ceff_ = -99;
2659
2660 xhit_fit_ceff_ = -2000;
2661 yhit_fit_ceff_ = -2000;
2662 zhit_fit_ceff_ = -2000;
2663
2664 hitchi2_fit_ceff_ = -99;
2665 helixchi2_fit_ceff_ = -99;
2666 score_fit_ceff_ = -17000;
2667
2668 dphi_fit_ceff_ = -99;
2669
2670 duplmask_fit_ceff_ = -1;
2671 nTkMatches_fit_ceff_ = -99;
2672 }
2673
2674 cmsswefftree_->Fill();
2675 }
2676 }
2677
2678 void TTreeValidation::fillCMSSWFakeRateTree(const Event& ev) {
2679 std::lock_guard<std::mutex> locker(glock_);
2680
2681 auto ievt = ev.evtID();
2682 const auto& evt_sim_tracks = ev.simTracks_;
2683 const auto& evt_cmssw_tracks = ev.cmsswTracks_;
2684 const auto& evt_cmssw_extras = ev.cmsswTracksExtra_;
2685 const auto& evt_build_tracks = ev.candidateTracks_;
2686 const auto& evt_build_extras = ev.candidateTracksExtra_;
2687 const auto& evt_fit_tracks = ev.fitTracks_;
2688 const auto& evt_fit_extras = ev.fitTracksExtra_;
2689 const auto& evt_layer_hits = ev.layerHits_;
2690
2691 for (const auto& buildtrack : evt_build_tracks) {
2692 if (Config::keepHitInfo) {
2693 hitlyrs_mc_cFR_.clear();
2694 hitlyrs_build_cFR_.clear();
2695 hitlyrs_cmssw_build_cFR_.clear();
2696 hitlyrs_fit_cFR_.clear();
2697 hitlyrs_cmssw_fit_cFR_.clear();
2698
2699 hitidxs_mc_cFR_.clear();
2700 hitidxs_build_cFR_.clear();
2701 hitidxs_cmssw_build_cFR_.clear();
2702 hitidxs_fit_cFR_.clear();
2703 hitidxs_cmssw_fit_cFR_.clear();
2704 }
2705
2706 algorithm_cFR_ = buildtrack.algoint();
2707
2708 const auto& buildextra = evt_build_extras[buildtrack.label()];
2709
2710
2711 evtID_cFR_ = ievt;
2712 seedID_cFR_ = buildextra.seedID();
2713 mcTrackID_cFR_ = buildextra.mcTrackID();
2714
2715
2716 pt_build_cFR_ = buildtrack.pT();
2717 ept_build_cFR_ = buildtrack.epT();
2718 phi_build_cFR_ = buildtrack.momPhi();
2719 ephi_build_cFR_ = buildtrack.emomPhi();
2720 eta_build_cFR_ = buildtrack.momEta();
2721 eeta_build_cFR_ = buildtrack.emomEta();
2722
2723
2724 if (mcTrackID_cFR_ >= 0) {
2725 const auto& simtrack = evt_sim_tracks[mcTrackID_cFR_];
2726 x_mc_cFR_ = simtrack.x();
2727 y_mc_cFR_ = simtrack.y();
2728 z_mc_cFR_ = simtrack.z();
2729 pt_mc_cFR_ = simtrack.pT();
2730 phi_mc_cFR_ = simtrack.momPhi();
2731 eta_mc_cFR_ = simtrack.momEta();
2732
2733 if (Config::keepHitInfo)
2734 TTreeValidation::fillMinHitInfo(simtrack, hitlyrs_mc_cFR_, hitidxs_mc_cFR_);
2735 } else {
2736 x_mc_cFR_ = -1000;
2737 y_mc_cFR_ = -1000;
2738 z_mc_cFR_ = -1000;
2739 pt_mc_cFR_ = -99;
2740 phi_mc_cFR_ = -99;
2741 eta_mc_cFR_ = -99;
2742 }
2743
2744
2745 nHits_build_cFR_ = buildtrack.nFoundHits();
2746 nLayers_build_cFR_ = buildtrack.nUniqueLayers();
2747 nHitsMatched_build_cFR_ = buildextra.nHitsMatched();
2748 fracHitsMatched_build_cFR_ = buildextra.fracHitsMatched();
2749 lastlyr_build_cFR_ = buildtrack.getLastFoundHitLyr();
2750
2751
2752 const Hit& lasthit = evt_layer_hits[buildtrack.getLastFoundHitLyr()][buildtrack.getLastFoundHitIdx()];
2753 xhit_build_cFR_ = lasthit.x();
2754 yhit_build_cFR_ = lasthit.y();
2755 zhit_build_cFR_ = lasthit.z();
2756
2757
2758 hitchi2_build_cFR_ = buildtrack.chi2();
2759 helixchi2_build_cFR_ = buildextra.helixChi2();
2760 score_build_cFR_ = buildtrack.score();
2761
2762
2763 dphi_build_cFR_ = buildextra.dPhi();
2764
2765 if (Config::keepHitInfo)
2766 TTreeValidation::fillMinHitInfo(buildtrack, hitlyrs_build_cFR_, hitidxs_build_cFR_);
2767
2768
2769 cmsswID_build_cFR_ = buildextra.cmsswTrackID();
2770 cmsswmask_build_cFR_ = TTreeValidation::getMaskAssignment(cmsswID_build_cFR_);
2771
2772 if (cmsswmask_build_cFR_ == 1)
2773 {
2774 const auto& cmsswtrack = evt_cmssw_tracks[cmsswID_build_cFR_];
2775 const auto& cmsswextra = evt_cmssw_extras[cmsswtrack.label()];
2776
2777 seedID_cmssw_build_cFR_ = cmsswextra.seedID();
2778
2779 x_cmssw_build_cFR_ = cmsswtrack.x();
2780 y_cmssw_build_cFR_ = cmsswtrack.y();
2781 z_cmssw_build_cFR_ = cmsswtrack.z();
2782
2783 pt_cmssw_build_cFR_ = cmsswtrack.pT();
2784 phi_cmssw_build_cFR_ = cmsswtrack.momPhi();
2785 eta_cmssw_build_cFR_ = cmsswtrack.momEta();
2786
2787 nHits_cmssw_build_cFR_ = cmsswtrack.nFoundHits();
2788 nLayers_cmssw_build_cFR_ = cmsswtrack.nUniqueLayers();
2789 lastlyr_cmssw_build_cFR_ = cmsswtrack.getLastFoundHitLyr();
2790
2791
2792 duplmask_build_cFR_ = buildextra.isDuplicate();
2793 iTkMatches_build_cFR_ = buildextra.duplicateID();
2794
2795 if (Config::keepHitInfo)
2796 TTreeValidation::fillMinHitInfo(cmsswtrack, hitlyrs_cmssw_build_cFR_, hitidxs_cmssw_build_cFR_);
2797 } else
2798 {
2799 seedID_cmssw_build_cFR_ = -99;
2800
2801 x_cmssw_build_cFR_ = -2000;
2802 y_cmssw_build_cFR_ = -2000;
2803 z_cmssw_build_cFR_ = -2000;
2804
2805 pt_cmssw_build_cFR_ = -99;
2806 phi_cmssw_build_cFR_ = -99;
2807 eta_cmssw_build_cFR_ = -99;
2808
2809 nHits_cmssw_build_cFR_ = -99;
2810 nLayers_cmssw_build_cFR_ = -99;
2811 lastlyr_cmssw_build_cFR_ = -99;
2812
2813 duplmask_build_cFR_ = -1;
2814 iTkMatches_build_cFR_ = -99;
2815 }
2816
2817
2818 if (buildToFitMap_.count(buildtrack.label())) {
2819 const auto& fittrack = evt_fit_tracks[buildToFitMap_[buildtrack.label()]];
2820 const auto& fitextra = evt_fit_extras[fittrack.label()];
2821
2822
2823 pt_fit_cFR_ = fittrack.pT();
2824 ept_fit_cFR_ = fittrack.epT();
2825 phi_fit_cFR_ = fittrack.momPhi();
2826 ephi_fit_cFR_ = fittrack.emomPhi();
2827 eta_fit_cFR_ = fittrack.momEta();
2828 eeta_fit_cFR_ = fittrack.emomEta();
2829
2830
2831 nHits_fit_cFR_ = fittrack.nFoundHits();
2832 nLayers_fit_cFR_ = fittrack.nUniqueLayers();
2833 nHitsMatched_fit_cFR_ = fitextra.nHitsMatched();
2834 fracHitsMatched_fit_cFR_ = fitextra.fracHitsMatched();
2835 lastlyr_fit_cFR_ = fittrack.getLastFoundHitLyr();
2836
2837
2838 const Hit& lasthit = evt_layer_hits[fittrack.getLastFoundHitLyr()][fittrack.getLastFoundHitIdx()];
2839 xhit_fit_cFR_ = lasthit.x();
2840 yhit_fit_cFR_ = lasthit.y();
2841 zhit_fit_cFR_ = lasthit.z();
2842
2843
2844 hitchi2_fit_cFR_ = fittrack.chi2();
2845 helixchi2_fit_cFR_ = fitextra.helixChi2();
2846 score_fit_cFR_ = fittrack.score();
2847
2848
2849 dphi_fit_cFR_ = fitextra.dPhi();
2850
2851 if (Config::keepHitInfo)
2852 TTreeValidation::fillMinHitInfo(buildtrack, hitlyrs_fit_cFR_, hitidxs_fit_cFR_);
2853
2854
2855 cmsswID_fit_cFR_ = fitextra.cmsswTrackID();
2856 cmsswmask_fit_cFR_ = TTreeValidation::getMaskAssignment(cmsswID_fit_cFR_);
2857
2858 if (cmsswmask_fit_cFR_ == 1)
2859 {
2860 const auto& cmsswtrack = evt_cmssw_tracks[cmsswID_fit_cFR_];
2861 const auto& cmsswextra = evt_cmssw_extras[cmsswtrack.label()];
2862
2863 seedID_cmssw_fit_cFR_ = cmsswextra.seedID();
2864
2865 x_cmssw_fit_cFR_ = cmsswtrack.x();
2866 y_cmssw_fit_cFR_ = cmsswtrack.y();
2867 z_cmssw_fit_cFR_ = cmsswtrack.z();
2868
2869 pt_cmssw_fit_cFR_ = cmsswtrack.pT();
2870 phi_cmssw_fit_cFR_ = cmsswtrack.momPhi();
2871 eta_cmssw_fit_cFR_ = cmsswtrack.momEta();
2872
2873 nHits_cmssw_fit_cFR_ = cmsswtrack.nFoundHits();
2874 nLayers_cmssw_fit_cFR_ = cmsswtrack.nUniqueLayers();
2875 lastlyr_cmssw_fit_cFR_ = cmsswtrack.getLastFoundHitLyr();
2876
2877
2878 duplmask_fit_cFR_ = fitextra.isDuplicate();
2879 iTkMatches_fit_cFR_ = fitextra.duplicateID();
2880
2881 if (Config::keepHitInfo)
2882 TTreeValidation::fillMinHitInfo(fittrack, hitlyrs_cmssw_fit_cFR_, hitidxs_cmssw_fit_cFR_);
2883 } else
2884 {
2885 seedID_cmssw_fit_cFR_ = -99;
2886
2887 x_cmssw_fit_cFR_ = -2000;
2888 y_cmssw_fit_cFR_ = -2000;
2889 z_cmssw_fit_cFR_ = -2000;
2890
2891 pt_cmssw_fit_cFR_ = -99;
2892 phi_cmssw_fit_cFR_ = -99;
2893 eta_cmssw_fit_cFR_ = -99;
2894
2895 nHits_cmssw_fit_cFR_ = -99;
2896 nLayers_cmssw_fit_cFR_ = -99;
2897 lastlyr_cmssw_fit_cFR_ = -99;
2898
2899 duplmask_fit_cFR_ = -1;
2900 iTkMatches_fit_cFR_ = -99;
2901 }
2902 } else
2903 {
2904 pt_fit_cFR_ = -100;
2905 ept_fit_cFR_ = -100;
2906 phi_fit_cFR_ = -100;
2907 ephi_fit_cFR_ = -100;
2908 eta_fit_cFR_ = -100;
2909 eeta_fit_cFR_ = -100;
2910
2911 nHits_fit_cFR_ = -100;
2912 nLayers_fit_cFR_ = -100;
2913 nHitsMatched_fit_cFR_ = -100;
2914 fracHitsMatched_fit_cFR_ = -100;
2915 lastlyr_fit_cFR_ = -100;
2916
2917 xhit_fit_cFR_ = -3000;
2918 yhit_fit_cFR_ = -3000;
2919 zhit_fit_cFR_ = -3000;
2920
2921 hitchi2_fit_cFR_ = -100;
2922 helixchi2_fit_cFR_ = -100;
2923 score_fit_cFR_ = -5001;
2924 dphi_fit_cFR_ = -100;
2925
2926 cmsswID_fit_cFR_ = -100;
2927 cmsswmask_fit_cFR_ = -2;
2928
2929 seedID_cmssw_fit_cFR_ = -100;
2930
2931 x_cmssw_fit_cFR_ = -3000;
2932 y_cmssw_fit_cFR_ = -3000;
2933 z_cmssw_fit_cFR_ = -3000;
2934
2935 pt_cmssw_fit_cFR_ = -100;
2936 phi_cmssw_fit_cFR_ = -100;
2937 eta_cmssw_fit_cFR_ = -100;
2938
2939 nHits_cmssw_fit_cFR_ = -100;
2940 nLayers_cmssw_fit_cFR_ = -100;
2941 lastlyr_cmssw_fit_cFR_ = -100;
2942
2943 duplmask_fit_cFR_ = -2;
2944 iTkMatches_fit_cFR_ = -100;
2945 }
2946
2947 cmsswfrtree_->Fill();
2948 }
2949 }
2950
2951 void TTreeValidation::saveTTrees() {
2952 std::lock_guard<std::mutex> locker(glock_);
2953 f_->cd();
2954
2955 if (Config::sim_val_for_cmssw || Config::sim_val) {
2956 efftree_->SetDirectory(f_.get());
2957 efftree_->Write();
2958
2959 frtree_->SetDirectory(f_.get());
2960 frtree_->Write();
2961 }
2962 if (Config::cmssw_val) {
2963 cmsswefftree_->SetDirectory(f_.get());
2964 cmsswefftree_->Write();
2965
2966 cmsswfrtree_->SetDirectory(f_.get());
2967 cmsswfrtree_->Write();
2968 }
2969 if (Config::fit_val) {
2970 fittree_->SetDirectory(f_.get());
2971 fittree_->Write();
2972 }
2973
2974 configtree_->SetDirectory(f_.get());
2975 configtree_->Write();
2976 }
2977
2978 }
2979 #endif