Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-05-07 01:50:03

0001 #include "performance.h"
0002 
0003 enum { pT5 = 7, pT3 = 5, T5 = 4, pLS = 8 };
0004 
0005 //__________________________________________________________________________________________________________________________________________________________________________
0006 int main(int argc, char** argv) {
0007   // Parse arguments
0008   parseArguments(argc, argv);
0009 
0010   // Initialize input and output root files
0011   initializeInputsAndOutputs();
0012 
0013   // Set of pdgids
0014   std::vector<int> pdgids = {0, 11, 211, 13, 321};
0015 
0016   // Set of charges
0017   std::vector<int> charges = {0, 1, -1};
0018 
0019   // Set of extra selections for efficiency plots
0020   std::vector<TString> selnames = {
0021       "base",    // default baseline that is more inline with MTV
0022       "loweta",  // When the eta cut is restricted to 2.4
0023       "xtr",     // When the eta cut is restricted to transition regions
0024       "vtr"      // When the eta cut is vetoing transition regions
0025   };
0026   std::vector<std::function<bool(unsigned int)>> sels = {
0027       [&](unsigned int isim) { return 1.; },
0028       [&](unsigned int isim) { return std::abs(lstEff.getVF("sim_eta").at(isim)) < 2.4; },
0029       [&](unsigned int isim) {
0030         return std::abs(lstEff.getVF("sim_eta").at(isim)) > 1.1 and std::abs(lstEff.getVF("sim_eta").at(isim)) < 1.7;
0031       },
0032       [&](unsigned int isim) {
0033         float sim_eta = lstEff.getVF("sim_eta").at(isim);
0034         return (std::abs(sim_eta) < 1.1 or std::abs(sim_eta) > 1.7) and std::abs(sim_eta) < 2.4;
0035       }};
0036   pdgids.insert(pdgids.end(), ana.pdgids.begin(), ana.pdgids.end());
0037 
0038   std::vector<SimTrackSetDefinition> list_effSetDef;
0039 
0040   // creating a set of efficiency plots for each pdgids being considered
0041   for (auto& pdgid : pdgids) {
0042     for (auto& charge : charges) {
0043       for (unsigned int isel = 0; isel < sels.size(); ++isel) {
0044         list_effSetDef.push_back(
0045             SimTrackSetDefinition(/* name  */
0046                                   TString("TC_") + selnames[isel],
0047                                   /* pdgid */ pdgid,
0048                                   /* q     */ charge,
0049                                   /* pass  */
0050                                   [&](unsigned int isim) { return lstEff.getVI("sim_TC_matched").at(isim) > 0; },
0051                                   /* sel   */ sels[isel]));
0052         list_effSetDef.push_back(SimTrackSetDefinition(
0053             /* name  */
0054             TString("pT5_") + selnames[isel],
0055             /* pdgid */ pdgid,
0056             /* q     */ charge,
0057             /* pass  */
0058             [&](unsigned int isim) { return lstEff.getVI("sim_TC_matched_mask").at(isim) & (1 << pT5); },
0059             /* sel   */ sels[isel]));
0060         list_effSetDef.push_back(SimTrackSetDefinition(
0061             /* name  */
0062             TString("pT3_") + selnames[isel],
0063             /* pdgid */ pdgid,
0064             /* q     */ charge,
0065             /* pass  */
0066             [&](unsigned int isim) { return lstEff.getVI("sim_TC_matched_mask").at(isim) & (1 << pT3); },
0067             /* sel   */ sels[isel]));
0068         list_effSetDef.push_back(SimTrackSetDefinition(
0069             /* name  */
0070             TString("T5_") + selnames[isel],
0071             /* pdgid */ pdgid,
0072             /* q     */ charge,
0073             /* pass  */
0074             [&](unsigned int isim) { return lstEff.getVI("sim_TC_matched_mask").at(isim) & (1 << T5); },
0075             /* sel   */ sels[isel]));
0076         list_effSetDef.push_back(SimTrackSetDefinition(
0077             /* name  */
0078             TString("pLS_") + selnames[isel],
0079             /* pdgid */ pdgid,
0080             /* q     */ charge,
0081             /* pass  */
0082             [&](unsigned int isim) { return lstEff.getVI("sim_TC_matched_mask").at(isim) & (1 << pLS); },
0083             /* sel   */ sels[isel]));
0084 
0085         if (ana.do_lower_level) {
0086           //lower objects - name will have pT5_lower_, T5_lower_, pT3_lower_
0087           list_effSetDef.push_back(SimTrackSetDefinition(
0088               /* name  */
0089               TString("pT5_lower_") + selnames[isel],
0090               /* pdgid */ pdgid,
0091               /* q     */ charge,
0092               /* pass  */ [&](unsigned int isim) { return lstEff.getVI("sim_pT5_matched").at(isim) > 0; },
0093               /* sel   */ sels[isel]));
0094           list_effSetDef.push_back(
0095               SimTrackSetDefinition(/* name  */
0096                                     TString("T5_lower_") + selnames[isel],
0097                                     /* pdgid */ pdgid,
0098                                     /* q     */ charge,
0099                                     /* pass  */
0100                                     [&](unsigned int isim) { return lstEff.getVI("sim_T5_matched").at(isim) > 0; },
0101                                     /* sel   */ sels[isel]));
0102           list_effSetDef.push_back(SimTrackSetDefinition(
0103               /* name  */
0104               TString("pT3_lower_") + selnames[isel],
0105               /* pdgid */ pdgid,
0106               /* q     */ charge,
0107               /* pass  */ [&](unsigned int isim) { return lstEff.getVI("sim_pT3_matched").at(isim) > 0; },
0108               /* sel   */ sels[isel]));
0109         }
0110       }
0111     }
0112   }
0113 
0114   bookEfficiencySets(list_effSetDef);
0115 
0116   // creating a set of fake rate plots
0117   std::vector<RecoTrackSetDefinition> list_FRSetDef;
0118   list_FRSetDef.push_back(
0119       RecoTrackSetDefinition(/* name  */
0120                              "TC",
0121                              /* pass  */
0122                              [&](unsigned int itc) { return lstEff.getVI("tc_isFake").at(itc) > 0; },
0123                              /* sel   */ [&](unsigned int itc) { return 1; },
0124                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0125                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0126                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0127                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0128   list_FRSetDef.push_back(
0129       RecoTrackSetDefinition(/* name  */
0130                              "pT5",
0131                              /* pass  */
0132                              [&](unsigned int itc) { return lstEff.getVI("tc_isFake").at(itc) > 0; },
0133                              /* sel   */
0134                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pT5; },
0135                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0136                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0137                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0138                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0139   list_FRSetDef.push_back(
0140       RecoTrackSetDefinition(/* name  */
0141                              "pT3",
0142                              /* pass  */
0143                              [&](unsigned int itc) { return lstEff.getVI("tc_isFake").at(itc) > 0; },
0144                              /* sel   */
0145                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pT3; },
0146                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0147                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0148                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0149                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0150   list_FRSetDef.push_back(
0151       RecoTrackSetDefinition(/* name  */
0152                              "T5",
0153                              /* pass  */
0154                              [&](unsigned int itc) { return lstEff.getVI("tc_isFake").at(itc) > 0; },
0155                              /* sel   */
0156                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == T5; },
0157                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0158                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0159                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0160                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0161   list_FRSetDef.push_back(
0162       RecoTrackSetDefinition(/* name  */
0163                              "pLS",
0164                              /* pass  */
0165                              [&](unsigned int itc) { return lstEff.getVI("tc_isFake").at(itc) > 0; },
0166                              /* sel   */
0167                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pLS; },
0168                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0169                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0170                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0171                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0172 
0173   if (ana.do_lower_level) {
0174     list_FRSetDef.push_back(RecoTrackSetDefinition(
0175         /* name  */
0176         "pT5_lower",
0177         /* pass  */ [&](unsigned int ipT5) { return lstEff.getVI("pT5_isFake").at(ipT5) > 0; },
0178         /* sel   */ [&](unsigned int ipT5) { return 1; },
0179         /* pt    */ [&]() { return lstEff.getVF("pT5_pt"); },
0180         /* eta   */ [&]() { return lstEff.getVF("pT5_eta"); },
0181         /* phi   */ [&]() { return lstEff.getVF("pT5_phi"); },
0182         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("pT5_pt").size(), 1); }));
0183     list_FRSetDef.push_back(RecoTrackSetDefinition(
0184         /* name  */
0185         "T5_lower",
0186         /* pass  */ [&](unsigned int it5) { return lstEff.getVI("t5_isFake").at(it5) > 0; },
0187         /* sel   */ [&](unsigned int it5) { return 1; },
0188         /* pt    */ [&]() { return lstEff.getVF("t5_pt"); },
0189         /* eta   */ [&]() { return lstEff.getVF("t5_eta"); },
0190         /* phi   */ [&]() { return lstEff.getVF("t5_phi"); },
0191         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("t5_pt").size(), 1); }));
0192     list_FRSetDef.push_back(RecoTrackSetDefinition(
0193         /* name  */
0194         "pT3_lower",
0195         /* pass  */ [&](unsigned int ipT3) { return lstEff.getVI("pT3_isFake").at(ipT3) > 0; },
0196         /* sel   */ [&](unsigned int ipT3) { return 1; },
0197         /* pt    */ [&]() { return lstEff.getVF("pT3_pt"); },
0198         /* eta   */ [&]() { return lstEff.getVF("pT3_eta"); },
0199         /* phi   */ [&]() { return lstEff.getVF("pT3_phi"); },
0200         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("pT3_pt").size(), 1); }));
0201   }
0202 
0203   bookFakeRateSets(list_FRSetDef);
0204 
0205   // creating a set of duplicate rate plots
0206   std::vector<RecoTrackSetDefinition> list_DRSetDef;
0207   list_DRSetDef.push_back(
0208       RecoTrackSetDefinition(/* name  */
0209                              "TC",
0210                              /* pass  */
0211                              [&](unsigned int itc) { return lstEff.getVI("tc_isDuplicate").at(itc) > 0; },
0212                              /* sel   */ [&](unsigned int itc) { return 1; },
0213                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0214                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0215                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0216                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0217   list_DRSetDef.push_back(
0218       RecoTrackSetDefinition(/* name  */
0219                              "pT5",
0220                              /* pass  */
0221                              [&](unsigned int itc) { return lstEff.getVI("tc_isDuplicate").at(itc) > 0; },
0222                              /* sel   */
0223                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pT5; },
0224                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0225                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0226                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0227                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0228   list_DRSetDef.push_back(
0229       RecoTrackSetDefinition(/* name  */
0230                              "pT3",
0231                              /* pass  */
0232                              [&](unsigned int itc) { return lstEff.getVI("tc_isDuplicate").at(itc) > 0; },
0233                              /* sel   */
0234                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pT3; },
0235                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0236                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0237                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0238                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0239   list_DRSetDef.push_back(
0240       RecoTrackSetDefinition(/* name  */
0241                              "T5",
0242                              /* pass  */
0243                              [&](unsigned int itc) { return lstEff.getVI("tc_isDuplicate").at(itc) > 0; },
0244                              /* sel   */
0245                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == T5; },
0246                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0247                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0248                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0249                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0250   list_DRSetDef.push_back(
0251       RecoTrackSetDefinition(/* name  */
0252                              "pLS",
0253                              /* pass  */
0254                              [&](unsigned int itc) { return lstEff.getVI("tc_isDuplicate").at(itc) > 0; },
0255                              /* sel   */
0256                              [&](unsigned int itc) { return lstEff.getVI("tc_type").at(itc) == pLS; },
0257                              /* pt    */ [&]() { return lstEff.getVF("tc_pt"); },
0258                              /* eta   */ [&]() { return lstEff.getVF("tc_eta"); },
0259                              /* phi   */ [&]() { return lstEff.getVF("tc_phi"); },
0260                              /* type  */ [&]() { return lstEff.getVI("tc_type"); }));
0261 
0262   if (ana.do_lower_level) {
0263     list_DRSetDef.push_back(RecoTrackSetDefinition(
0264         /* name  */
0265         "pT5_lower",
0266         /* pass  */ [&](unsigned int ipT5) { return lstEff.getVI("pT5_isDuplicate").at(ipT5) > 0; },
0267         /* sel   */ [&](unsigned int ipT5) { return 1; },
0268         /* pt    */ [&]() { return lstEff.getVF("pT5_pt"); },
0269         /* eta   */ [&]() { return lstEff.getVF("pT5_eta"); },
0270         /* phi   */ [&]() { return lstEff.getVF("pT5_phi"); },
0271         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("pT5_pt").size(), 1); }));
0272     list_DRSetDef.push_back(RecoTrackSetDefinition(
0273         /* name  */
0274         "T5_lower",
0275         /* pass  */ [&](unsigned int it5) { return lstEff.getVI("t5_isDuplicate").at(it5) > 0; },
0276         /* sel   */ [&](unsigned int it5) { return 1; },
0277         /* pt    */ [&]() { return lstEff.getVF("t5_pt"); },
0278         /* eta   */ [&]() { return lstEff.getVF("t5_eta"); },
0279         /* phi   */ [&]() { return lstEff.getVF("t5_phi"); },
0280         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("t5_pt").size(), 1); }));
0281     list_DRSetDef.push_back(RecoTrackSetDefinition(
0282         /* name  */
0283         "pT3_lower",
0284         /* pass  */ [&](unsigned int ipT3) { return lstEff.getVI("pT3_isDuplicate").at(ipT3) > 0; },
0285         /* sel   */ [&](unsigned int ipT3) { return 1; },
0286         /* pt    */ [&]() { return lstEff.getVF("pT3_pt"); },
0287         /* eta   */ [&]() { return lstEff.getVF("pT3_eta"); },
0288         /* phi   */ [&]() { return lstEff.getVF("pT3_phi"); },
0289         /* type  */ [&]() { return std::vector<int>(lstEff.getVF("pT3_pt").size(), 1); }));
0290   }
0291 
0292   bookDuplicateRateSets(list_DRSetDef);
0293 
0294   // Book Histograms
0295   ana.cutflow.bookHistograms(ana.histograms);  // if just want to book everywhere
0296 
0297   int nevts = 0;
0298 
0299   // Looping input file
0300   while (ana.looper.nextEvent()) {
0301     // If splitting jobs are requested then determine whether to process the event or not based on remainder
0302     if (ana.job_index != -1 and ana.nsplit_jobs != -1) {
0303       if (ana.looper.getNEventsProcessed() % ana.nsplit_jobs != (unsigned int)ana.job_index)
0304         continue;
0305     }
0306 
0307     // Reset all temporary variables necessary for histogramming
0308     ana.tx.clear();
0309 
0310     // Compute all temporary variables and pack the vector quantities that will get filled to the histograms
0311     fillEfficiencySets(list_effSetDef);
0312     fillFakeRateSets(list_FRSetDef);
0313     fillDuplicateRateSets(list_DRSetDef);
0314 
0315     // Reset all temporary variables necessary for histogramming
0316     ana.cutflow.fill();
0317 
0318     // Counting number of events processed
0319     nevts++;
0320   }
0321 
0322   // Write number of events processed
0323   ana.output_tfile->cd();
0324   TH1F* h_nevts = new TH1F("nevts", "nevts", 1, 0, 1);
0325   h_nevts->SetBinContent(1, nevts);
0326   h_nevts->Write();
0327 
0328   // Writing output file
0329   ana.cutflow.saveOutput();
0330 
0331   // The below can be sometimes crucial
0332   delete ana.output_tfile;
0333 
0334   // Done
0335   return 0;
0336 }
0337 
0338 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0339 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0340 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0341 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0342 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0343 
0344 //__________________________________________________________________________________________________________________________________________________________________________
0345 void bookEfficiencySets(std::vector<SimTrackSetDefinition>& effsets) {
0346   for (auto& effset : effsets)
0347     bookEfficiencySet(effset);
0348 }
0349 
0350 //__________________________________________________________________________________________________________________________________________________________________________
0351 void bookEfficiencySet(SimTrackSetDefinition& effset) {
0352   TString category_name = TString::Format("%s_%d_%d", effset.set_name.Data(), effset.pdgid, effset.q);
0353 
0354   // Denominator tracks' quantities
0355   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_pt");
0356   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_eta");
0357   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_dxy");
0358   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_vxy");
0359   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_dz");
0360   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_denom_phi");
0361 
0362   // Numerator tracks' quantities
0363   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_pt");
0364   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_eta");
0365   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_dxy");
0366   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_vxy");
0367   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_dz");
0368   ana.tx.createBranch<std::vector<float>>(category_name + "_ef_numer_phi");
0369 
0370   // Inefficiencies
0371   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_pt");
0372   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_eta");
0373   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_dxy");
0374   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_vxy");
0375   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_dz");
0376   ana.tx.createBranch<std::vector<float>>(category_name + "_ie_numer_phi");
0377 
0378   ana.histograms.addVecHistogram(category_name + "_ef_denom_pt", getPtBounds(0), [&, category_name]() {
0379     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_pt");
0380   });
0381   ana.histograms.addVecHistogram(category_name + "_ef_denom_ptlow", getPtBounds(4), [&, category_name]() {
0382     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_pt");
0383   });
0384   ana.histograms.addVecHistogram(category_name + "_ef_denom_ptmtv", getPtBounds(9), [&, category_name]() {
0385     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_pt");
0386   });
0387   ana.histograms.addVecHistogram(category_name + "_ef_denom_ptflatbin", 180, 0., 100, [&, category_name]() {
0388     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_pt");
0389   });
0390   ana.histograms.addVecHistogram(category_name + "_ef_denom_eta", 180, -4.5, 4.5, [&, category_name]() {
0391     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_eta");
0392   });
0393   ana.histograms.addVecHistogram(category_name + "_ef_denom_dxy", 180, -30., 30., [&, category_name]() {
0394     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_dxy");
0395   });
0396   ana.histograms.addVecHistogram(category_name + "_ef_denom_vxy", 180, -30., 30., [&, category_name]() {
0397     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_vxy");
0398   });
0399   ana.histograms.addVecHistogram(category_name + "_ef_denom_dz", 180, -30., 30., [&, category_name]() {
0400     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_dz");
0401   });
0402   ana.histograms.addVecHistogram(category_name + "_ef_denom_phi", 180, -M_PI, M_PI, [&, category_name]() {
0403     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_denom_phi");
0404   });
0405 
0406   ana.histograms.addVecHistogram(category_name + "_ef_numer_pt", getPtBounds(0), [&, category_name]() {
0407     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_pt");
0408   });
0409   ana.histograms.addVecHistogram(category_name + "_ef_numer_ptlow", getPtBounds(4), [&, category_name]() {
0410     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_pt");
0411   });
0412   ana.histograms.addVecHistogram(category_name + "_ef_numer_ptmtv", getPtBounds(9), [&, category_name]() {
0413     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_pt");
0414   });
0415   ana.histograms.addVecHistogram(category_name + "_ef_numer_ptflatbin", 180, 0., 100, [&, category_name]() {
0416     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_pt");
0417   });
0418   ana.histograms.addVecHistogram(category_name + "_ef_numer_eta", 180, -4.5, 4.5, [&, category_name]() {
0419     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_eta");
0420   });
0421   ana.histograms.addVecHistogram(category_name + "_ef_numer_dxy", 180, -30., 30., [&, category_name]() {
0422     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_dxy");
0423   });
0424   ana.histograms.addVecHistogram(category_name + "_ef_numer_vxy", 180, -30., 30., [&, category_name]() {
0425     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_vxy");
0426   });
0427   ana.histograms.addVecHistogram(category_name + "_ef_numer_dz", 180, -30., 30., [&, category_name]() {
0428     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_dz");
0429   });
0430   ana.histograms.addVecHistogram(category_name + "_ef_numer_phi", 180, -M_PI, M_PI, [&, category_name]() {
0431     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ef_numer_phi");
0432   });
0433 
0434   ana.histograms.addVecHistogram(category_name + "_ie_numer_pt", getPtBounds(0), [&, category_name]() {
0435     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_pt");
0436   });
0437   ana.histograms.addVecHistogram(category_name + "_ie_numer_ptlow", getPtBounds(4), [&, category_name]() {
0438     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_pt");
0439   });
0440   ana.histograms.addVecHistogram(category_name + "_ie_numer_ptmtv", getPtBounds(9), [&, category_name]() {
0441     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_pt");
0442   });
0443   ana.histograms.addVecHistogram(category_name + "_ie_numer_ptflatbin", 180, 0., 100, [&, category_name]() {
0444     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_pt");
0445   });
0446   ana.histograms.addVecHistogram(category_name + "_ie_numer_eta", 180, -4.5, 4.5, [&, category_name]() {
0447     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_eta");
0448   });
0449   ana.histograms.addVecHistogram(category_name + "_ie_numer_dxy", 180, -30., 30., [&, category_name]() {
0450     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_dxy");
0451   });
0452   ana.histograms.addVecHistogram(category_name + "_ie_numer_vxy", 180, -30., 30., [&, category_name]() {
0453     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_vxy");
0454   });
0455   ana.histograms.addVecHistogram(category_name + "_ie_numer_dz", 180, -30., 30., [&, category_name]() {
0456     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_dz");
0457   });
0458   ana.histograms.addVecHistogram(category_name + "_ie_numer_phi", 180, -M_PI, M_PI, [&, category_name]() {
0459     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_ie_numer_phi");
0460   });
0461 }
0462 
0463 //__________________________________________________________________________________________________________________________________________________________________________
0464 void bookDuplicateRateSets(std::vector<RecoTrackSetDefinition>& DRsets) {
0465   for (auto& DRset : DRsets) {
0466     bookDuplicateRateSet(DRset);
0467   }
0468 }
0469 
0470 //__________________________________________________________________________________________________________________________________________________________________________
0471 void bookDuplicateRateSet(RecoTrackSetDefinition& DRset) {
0472   TString category_name = DRset.set_name;
0473 
0474   // Denominator tracks' quantities
0475   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_denom_pt");
0476   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_denom_eta");
0477   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_denom_phi");
0478 
0479   // Numerator tracks' quantities
0480   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_numer_pt");
0481   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_numer_eta");
0482   ana.tx.createBranch<std::vector<float>>(category_name + "_dr_numer_phi");
0483 
0484   // Histogram utility object that is used to define the histograms
0485   ana.histograms.addVecHistogram(category_name + "_dr_denom_pt", getPtBounds(0), [&, category_name]() {
0486     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_denom_pt");
0487   });
0488   ana.histograms.addVecHistogram(category_name + "_dr_denom_ptlow", getPtBounds(4), [&, category_name]() {
0489     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_denom_pt");
0490   });
0491   ana.histograms.addVecHistogram(category_name + "_dr_denom_ptmtv", getPtBounds(9), [&, category_name]() {
0492     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_denom_pt");
0493   });
0494   ana.histograms.addVecHistogram(category_name + "_dr_denom_eta", 180, -4.5, 4.5, [&, category_name]() {
0495     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_denom_eta");
0496   });
0497   ana.histograms.addVecHistogram(category_name + "_dr_denom_phi", 180, -M_PI, M_PI, [&, category_name]() {
0498     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_denom_phi");
0499   });
0500   ana.histograms.addVecHistogram(category_name + "_dr_numer_pt", getPtBounds(0), [&, category_name]() {
0501     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_numer_pt");
0502   });
0503   ana.histograms.addVecHistogram(category_name + "_dr_numer_ptlow", getPtBounds(4), [&, category_name]() {
0504     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_numer_pt");
0505   });
0506   ana.histograms.addVecHistogram(category_name + "_dr_numer_ptmtv", getPtBounds(9), [&, category_name]() {
0507     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_numer_pt");
0508   });
0509   ana.histograms.addVecHistogram(category_name + "_dr_numer_eta", 180, -4.5, 4.5, [&, category_name]() {
0510     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_numer_eta");
0511   });
0512   ana.histograms.addVecHistogram(category_name + "_dr_numer_phi", 180, -M_PI, M_PI, [&, category_name]() {
0513     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_dr_numer_phi");
0514   });
0515 }
0516 
0517 //__________________________________________________________________________________________________________________________________________________________________________
0518 void bookFakeRateSets(std::vector<RecoTrackSetDefinition>& FRsets) {
0519   for (auto& FRset : FRsets) {
0520     bookFakeRateSet(FRset);
0521   }
0522 }
0523 
0524 //__________________________________________________________________________________________________________________________________________________________________________
0525 void bookFakeRateSet(RecoTrackSetDefinition& FRset) {
0526   TString category_name = FRset.set_name;
0527 
0528   // Denominator tracks' quantities
0529   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_denom_pt");
0530   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_denom_eta");
0531   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_denom_phi");
0532 
0533   // Numerator tracks' quantities
0534   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_numer_pt");
0535   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_numer_eta");
0536   ana.tx.createBranch<std::vector<float>>(category_name + "_fr_numer_phi");
0537 
0538   // Histogram utility object that is used to define the histograms
0539   ana.histograms.addVecHistogram(category_name + "_fr_denom_pt", getPtBounds(0), [&, category_name]() {
0540     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_denom_pt");
0541   });
0542   ana.histograms.addVecHistogram(category_name + "_fr_denom_ptlow", getPtBounds(4), [&, category_name]() {
0543     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_denom_pt");
0544   });
0545   ana.histograms.addVecHistogram(category_name + "_fr_denom_ptmtv", getPtBounds(9), [&, category_name]() {
0546     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_denom_pt");
0547   });
0548   ana.histograms.addVecHistogram(category_name + "_fr_denom_eta", 180, -4.5, 4.5, [&, category_name]() {
0549     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_denom_eta");
0550   });
0551   ana.histograms.addVecHistogram(category_name + "_fr_numer_phi", 180, -M_PI, M_PI, [&, category_name]() {
0552     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_numer_phi");
0553   });
0554   ana.histograms.addVecHistogram(category_name + "_fr_numer_pt", getPtBounds(0), [&, category_name]() {
0555     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_numer_pt");
0556   });
0557   ana.histograms.addVecHistogram(category_name + "_fr_numer_ptlow", getPtBounds(4), [&, category_name]() {
0558     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_numer_pt");
0559   });
0560   ana.histograms.addVecHistogram(category_name + "_fr_numer_ptmtv", getPtBounds(9), [&, category_name]() {
0561     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_numer_pt");
0562   });
0563   ana.histograms.addVecHistogram(category_name + "_fr_numer_eta", 180, -4.5, 4.5, [&, category_name]() {
0564     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_numer_eta");
0565   });
0566   ana.histograms.addVecHistogram(category_name + "_fr_denom_phi", 180, -M_PI, M_PI, [&, category_name]() {
0567     return ana.tx.getBranchLazy<std::vector<float>>(category_name + "_fr_denom_phi");
0568   });
0569 }
0570 
0571 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0572 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0573 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0574 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0575 // ---------------------------------------------------------=============================================-------------------------------------------------------------------
0576 
0577 //__________________________________________________________________________________________________________________________________________________________________________
0578 void fillEfficiencySets(std::vector<SimTrackSetDefinition>& effsets) {
0579   std::vector<float> const& pt = lstEff.getVF("sim_pt");
0580   std::vector<float> const& eta = lstEff.getVF("sim_eta");
0581   std::vector<float> const& dz = lstEff.getVF("sim_pca_dz");
0582   std::vector<float> const& dxy = lstEff.getVF("sim_pca_dxy");
0583   std::vector<float> const& phi = lstEff.getVF("sim_phi");
0584   std::vector<int> const& pdgidtrk = lstEff.getVI("sim_pdgId");
0585   std::vector<int> const& q = lstEff.getVI("sim_q");
0586   std::vector<float> const& vtx_x = lstEff.getVF("sim_vx");
0587   std::vector<float> const& vtx_y = lstEff.getVF("sim_vy");
0588   std::vector<float> const& vtx_z = lstEff.getVF("sim_vz");
0589 
0590   for (auto& effset : effsets) {
0591     for (unsigned int isimtrk = 0; isimtrk < lstEff.getVF("sim_pt").size(); ++isimtrk) {
0592       fillEfficiencySet(isimtrk,
0593                         effset,
0594                         pt.at(isimtrk),
0595                         eta.at(isimtrk),
0596                         dz.at(isimtrk),
0597                         dxy.at(isimtrk),
0598                         phi.at(isimtrk),
0599                         pdgidtrk.at(isimtrk),
0600                         q.at(isimtrk),
0601                         vtx_x.at(isimtrk),
0602                         vtx_y.at(isimtrk),
0603                         vtx_z.at(isimtrk));
0604     }
0605   }
0606 }
0607 
0608 //__________________________________________________________________________________________________________________________________________________________________________
0609 void fillEfficiencySet(int isimtrk,
0610                        SimTrackSetDefinition& effset,
0611                        float pt,
0612                        float eta,
0613                        float dz,
0614                        float dxy,
0615                        float phi,
0616                        int pdgidtrk,
0617                        int q,
0618                        float vtx_x,
0619                        float vtx_y,
0620                        float vtx_z) {
0621   //=========================================================
0622   // NOTE: The following is not applied as the LSTNtuple no longer writes this.
0623   // const int &bunch = lstEff.getVI("sim_bunchCrossing")[isimtrk];
0624   // const int &event = lstEff.getVI("sim_event")[isimtrk];
0625   // if (bunch != 0)
0626   //     return;
0627   // if (event != 0)
0628   //     return;
0629   //=========================================================
0630 
0631   const float vtx_perp = sqrt(vtx_x * vtx_x + vtx_y * vtx_y);
0632   bool pass = effset.pass(isimtrk);
0633   bool sel = effset.sel(isimtrk);
0634 
0635   if (effset.pdgid != 0) {
0636     if (std::abs(pdgidtrk) != std::abs(effset.pdgid))
0637       return;
0638   }
0639 
0640   if (effset.q != 0) {
0641     if (q != effset.q)
0642       return;
0643   }
0644 
0645   if (effset.pdgid == 0 and q == 0)
0646     return;
0647 
0648   if (not sel)
0649     return;
0650 
0651   TString category_name = TString::Format("%s_%d_%d", effset.set_name.Data(), effset.pdgid, effset.q);
0652 
0653   // https://github.com/cms-sw/cmssw/blob/7cbdb18ec6d11d5fd17ca66c1153f0f4e869b6b0/SimTracker/Common/python/trackingParticleSelector_cfi.py
0654   // https://github.com/cms-sw/cmssw/blob/7cbdb18ec6d11d5fd17ca66c1153f0f4e869b6b0/SimTracker/Common/interface/TrackingParticleSelector.h#L122-L124
0655   const float vtx_z_thresh = 30;
0656   const float vtx_perp_thresh = 2.5;
0657 
0658   // N minus eta cut
0659   if (pt > ana.pt_cut and std::abs(vtx_z) < vtx_z_thresh and std::abs(vtx_perp) < vtx_perp_thresh) {
0660     // vs. eta plot
0661     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_eta", eta);
0662     if (pass)
0663       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_eta", eta);
0664     else
0665       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_eta", eta);
0666   }
0667 
0668   // N minus pt cut
0669   if (std::abs(eta) < ana.eta_cut and std::abs(vtx_z) < vtx_z_thresh and std::abs(vtx_perp) < vtx_perp_thresh) {
0670     // vs. pt plot
0671     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_pt", pt);
0672     if (pass)
0673       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_pt", pt);
0674     else
0675       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_pt", pt);
0676   }
0677 
0678   // N minus dxy cut
0679   if (std::abs(eta) < ana.eta_cut and pt > ana.pt_cut and std::abs(vtx_z) < vtx_z_thresh) {
0680     // vs. dxy plot
0681     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_dxy", dxy);
0682     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_vxy", vtx_perp);
0683     if (pass) {
0684       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_dxy", dxy);
0685       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_vxy", vtx_perp);
0686     } else {
0687       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_dxy", dxy);
0688       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_vxy", vtx_perp);
0689     }
0690   }
0691 
0692   // N minus dz cut
0693   if (std::abs(eta) < ana.eta_cut and pt > ana.pt_cut and std::abs(vtx_perp) < vtx_perp_thresh) {
0694     // vs. dz plot
0695     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_dz", dz);
0696     if (pass)
0697       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_dz", dz);
0698     else
0699       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_dz", dz);
0700   }
0701 
0702   // All phase-space cuts
0703   if (std::abs(eta) < ana.eta_cut and pt > ana.pt_cut and std::abs(vtx_z) < vtx_z_thresh and
0704       std::abs(vtx_perp) < vtx_perp_thresh) {
0705     // vs. Phi plot
0706     ana.tx.pushbackToBranch<float>(category_name + "_ef_denom_phi", phi);
0707     if (pass)
0708       ana.tx.pushbackToBranch<float>(category_name + "_ef_numer_phi", phi);
0709     else
0710       ana.tx.pushbackToBranch<float>(category_name + "_ie_numer_phi", phi);
0711   }
0712 }
0713 
0714 //__________________________________________________________________________________________________________________________________________________________________________
0715 void fillFakeRateSets(std::vector<RecoTrackSetDefinition>& FRsets) {
0716   for (auto& FRset : FRsets) {
0717     std::vector<float> const& pt = FRset.pt();
0718     std::vector<float> const& eta = FRset.eta();
0719     std::vector<float> const& phi = FRset.phi();
0720     for (unsigned int itc = 0; itc < FRset.pt().size(); ++itc) {
0721       fillFakeRateSet(itc, FRset, pt.at(itc), eta.at(itc), phi.at(itc));
0722     }
0723   }
0724 }
0725 
0726 //__________________________________________________________________________________________________________________________________________________________________________
0727 void fillFakeRateSet(int itc, RecoTrackSetDefinition& FRset, float pt, float eta, float phi) {
0728   TString category_name = FRset.set_name;
0729   bool pass = FRset.pass(itc);
0730   bool sel = FRset.sel(itc);
0731 
0732   if (not sel)
0733     return;
0734 
0735   if (pt > ana.pt_cut) {
0736     ana.tx.pushbackToBranch<float>(category_name + "_fr_denom_eta", eta);
0737     if (pass)
0738       ana.tx.pushbackToBranch<float>(category_name + "_fr_numer_eta", eta);
0739   }
0740   if (std::abs(eta) < ana.eta_cut) {
0741     ana.tx.pushbackToBranch<float>(category_name + "_fr_denom_pt", pt);
0742     if (pass)
0743       ana.tx.pushbackToBranch<float>(category_name + "_fr_numer_pt", pt);
0744   }
0745   if (std::abs(eta) < ana.eta_cut and pt > ana.pt_cut) {
0746     ana.tx.pushbackToBranch<float>(category_name + "_fr_denom_phi", phi);
0747     if (pass)
0748       ana.tx.pushbackToBranch<float>(category_name + "_fr_numer_phi", phi);
0749   }
0750 }
0751 
0752 //__________________________________________________________________________________________________________________________________________________________________________
0753 void fillDuplicateRateSets(std::vector<RecoTrackSetDefinition>& DRsets) {
0754   for (auto& DRset : DRsets) {
0755     std::vector<float> const& pt = DRset.pt();
0756     std::vector<float> const& eta = DRset.eta();
0757     std::vector<float> const& phi = DRset.phi();
0758     for (unsigned int itc = 0; itc < DRset.pt().size(); ++itc) {
0759       fillDuplicateRateSet(itc, DRset, pt.at(itc), eta.at(itc), phi.at(itc));
0760     }
0761   }
0762 }
0763 
0764 //__________________________________________________________________________________________________________________________________________________________________________
0765 void fillDuplicateRateSet(int itc, RecoTrackSetDefinition& DRset, float pt, float eta, float phi) {
0766   TString category_name = DRset.set_name;
0767   bool pass = DRset.pass(itc);
0768   bool sel = DRset.sel(itc);
0769 
0770   if (not sel)
0771     return;
0772 
0773   if (pt > ana.pt_cut) {
0774     ana.tx.pushbackToBranch<float>(category_name + "_dr_denom_eta", eta);
0775     if (pass)
0776       ana.tx.pushbackToBranch<float>(category_name + "_dr_numer_eta", eta);
0777   }
0778   if (std::abs(eta) < ana.eta_cut) {
0779     ana.tx.pushbackToBranch<float>(category_name + "_dr_denom_pt", pt);
0780     if (pass)
0781       ana.tx.pushbackToBranch<float>(category_name + "_dr_numer_pt", pt);
0782   }
0783   if (std::abs(eta) < ana.eta_cut and pt > ana.pt_cut) {
0784     ana.tx.pushbackToBranch<float>(category_name + "_dr_denom_phi", phi);
0785     if (pass)
0786       ana.tx.pushbackToBranch<float>(category_name + "_dr_numer_phi", phi);
0787   }
0788 }