Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-04-22 06:27:41 

0001 #include <TVector3.h>
0002 
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004 #include "FWCore/Framework/interface/stream/EDProducer.h"
0005 
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/MakerMacros.h"
0008 
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/Utilities/interface/StreamID.h"
0011 #include "FWCore/Utilities/interface/ESGetToken.h"
0012 
0013 #include "DataFormats/Candidate/interface/Candidate.h"
0014 #include "DataFormats/PatCandidates/interface/Jet.h"
0015 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0016 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0017 
0018 #include "RecoBTag/FeatureTools/interface/TrackInfoBuilder.h"
0019 #include "RecoBTag/FeatureTools/interface/deep_helpers.h"
0020 #include "RecoBTag/FeatureTools/interface/sorting_modules.h"
0021 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0022 
0023 #include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h"
0024 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0025 
0026 #include "DataFormats/BTauReco/interface/DeepBoostedJetTagInfo.h"
0027 #include "DataFormats/Common/interface/AssociationMap.h"
0028 
0029 using namespace btagbtvdeep;
0030 
0031 class DeepBoostedJetTagInfoProducer : public edm::stream::EDProducer<> {
0032 public:
0033   explicit DeepBoostedJetTagInfoProducer(const edm::ParameterSet &);
0034   ~DeepBoostedJetTagInfoProducer() override;
0035 
0036   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0037 
0038 private:
0039   typedef std::vector<reco::DeepBoostedJetTagInfo> DeepBoostedJetTagInfoCollection;
0040   typedef reco::VertexCompositePtrCandidateCollection SVCollection;
0041   typedef reco::VertexCollection VertexCollection;
0042   typedef edm::View<reco::Candidate> CandidateView;
0043   typedef edm::AssociationMap<edm::OneToOne<reco::JetView, reco::JetView>> JetMatchMap;
0044 
0045   void beginStream(edm::StreamID) override {}
0046   void produce(edm::Event &, const edm::EventSetup &) override;
0047   void endStream() override {}
0048 
0049   void fillParticleFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &unsubJet, const reco::Jet &jet);
0050   void fillSVFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet);
0051   void fillParticleFeaturesHLT(DeepBoostedJetFeatures &fts, const reco::Jet &jet, const reco::VertexRefProd &PVRefProd);
0052   void fillSVFeaturesHLT(DeepBoostedJetFeatures &fts, const reco::Jet &jet);
0053 
0054   float puppiWgt(const reco::CandidatePtr &cand);
0055   bool useTrackProperties(const reco::PFCandidate *reco_cand);
0056 
0057   const double jet_radius_;
0058   const double min_jet_pt_;
0059   const double max_jet_eta_;
0060   const double min_pt_for_track_properties_;
0061   const double min_pt_for_pfcandidates_;
0062   const bool use_puppiP4_;
0063   const double min_puppi_wgt_;
0064   const bool include_neutrals_;
0065   const bool sort_by_sip2dsig_;
0066   const bool flip_ip_sign_;
0067   const double max_sip3dsig_;
0068   const bool use_hlt_features_;
0069 
0070   edm::EDGetTokenT<edm::View<reco::Jet>> jet_token_;
0071   edm::EDGetTokenT<JetMatchMap> unsubjet_map_token_;
0072   edm::EDGetTokenT<VertexCollection> vtx_token_;
0073   edm::EDGetTokenT<SVCollection> sv_token_;
0074   edm::EDGetTokenT<CandidateView> pfcand_token_;
0075 
0076   bool use_puppi_value_map_;
0077   bool use_pvasq_value_map_;
0078   bool is_packed_pf_candidate_collection_;
0079   bool use_unsubjet_map_;
0080 
0081   edm::EDGetTokenT<edm::ValueMap<float>> puppi_value_map_token_;
0082   edm::EDGetTokenT<edm::ValueMap<int>> pvasq_value_map_token_;
0083   edm::EDGetTokenT<edm::Association<VertexCollection>> pvas_token_;
0084   edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> track_builder_token_;
0085 
0086   edm::Handle<VertexCollection> vtxs_;
0087   edm::Handle<SVCollection> svs_;
0088   edm::Handle<CandidateView> pfcands_;
0089   edm::ESHandle<TransientTrackBuilder> track_builder_;
0090   edm::Handle<edm::ValueMap<float>> puppi_value_map_;
0091   edm::Handle<edm::ValueMap<int>> pvasq_value_map_;
0092   edm::Handle<edm::Association<VertexCollection>> pvas_;
0093 
0094   const static std::vector<std::string> particle_features_;
0095   const static std::vector<std::string> sv_features_;
0096   const static std::vector<std::string> particle_features_hlt_;
0097   const static std::vector<std::string> sv_features_hlt_;
0098   const reco::Vertex *pv_ = nullptr;
0099   const static float min_track_pt_property_;
0100   const static int min_valid_pixel_hits_;
0101   const int covarianceVersion_;
0102   const std::vector<int> covariancePackingSchemas_;
0103 
0104   std::map<reco::CandidatePtr::key_type, float> puppi_wgt_cache;
0105 
0106   const static std::vector<std::string> particle_features_scouting_;
0107   const bool use_scouting_features_;
0108   edm::EDGetTokenT<edm::ValueMap<float>> normchi2_value_map_token_;
0109   edm::EDGetTokenT<edm::ValueMap<float>> dz_value_map_token_;
0110   edm::EDGetTokenT<edm::ValueMap<float>> dxy_value_map_token_;
0111   edm::EDGetTokenT<edm::ValueMap<float>> dzsig_value_map_token_;
0112   edm::EDGetTokenT<edm::ValueMap<float>> dxysig_value_map_token_;
0113   edm::EDGetTokenT<edm::ValueMap<int>> lostInnerHits_value_map_token_;
0114   edm::EDGetTokenT<edm::ValueMap<int>> quality_value_map_token_;
0115   edm::EDGetTokenT<edm::ValueMap<float>> trkPt_value_map_token_;
0116   edm::EDGetTokenT<edm::ValueMap<float>> trkEta_value_map_token_;
0117   edm::EDGetTokenT<edm::ValueMap<float>> trkPhi_value_map_token_;
0118 
0119   edm::Handle<edm::ValueMap<float>> normchi2_value_map_;
0120   edm::Handle<edm::ValueMap<float>> dz_value_map_;
0121   edm::Handle<edm::ValueMap<float>> dxy_value_map_;
0122   edm::Handle<edm::ValueMap<float>> dzsig_value_map_;
0123   edm::Handle<edm::ValueMap<float>> dxysig_value_map_;
0124   edm::Handle<edm::ValueMap<int>> lostInnerHits_value_map_;
0125   edm::Handle<edm::ValueMap<int>> quality_value_map_;
0126   edm::Handle<edm::ValueMap<float>> trkPt_value_map_;
0127   edm::Handle<edm::ValueMap<float>> trkEta_value_map_;
0128   edm::Handle<edm::ValueMap<float>> trkPhi_value_map_;
0129 };
0130 
0131 const std::vector<std::string> DeepBoostedJetTagInfoProducer::particle_features_{
0132     "pfcand_puppiw",        "pfcand_hcalFrac",       "pfcand_VTX_ass",      "pfcand_lostInnerHits",
0133     "pfcand_quality",       "pfcand_charge",         "pfcand_isEl",         "pfcand_isMu",
0134     "pfcand_isChargedHad",  "pfcand_isGamma",        "pfcand_isNeutralHad", "pfcand_phirel",
0135     "pfcand_etarel",        "pfcand_deltaR",         "pfcand_abseta",       "pfcand_ptrel_log",
0136     "pfcand_erel_log",      "pfcand_pt_log",         "pfcand_drminsv",      "pfcand_drsubjet1",
0137     "pfcand_drsubjet2",     "pfcand_normchi2",       "pfcand_dz",           "pfcand_dzsig",
0138     "pfcand_dxy",           "pfcand_dxysig",         "pfcand_dptdpt",       "pfcand_detadeta",
0139     "pfcand_dphidphi",      "pfcand_dxydxy",         "pfcand_dzdz",         "pfcand_dxydz",
0140     "pfcand_dphidxy",       "pfcand_dlambdadz",      "pfcand_btagEtaRel",   "pfcand_btagPtRatio",
0141     "pfcand_btagPParRatio", "pfcand_btagSip2dVal",   "pfcand_btagSip2dSig", "pfcand_btagSip3dVal",
0142     "pfcand_btagSip3dSig",  "pfcand_btagJetDistVal", "pfcand_mask",         "pfcand_pt_log_nopuppi",
0143     "pfcand_e_log_nopuppi", "pfcand_ptrel",          "pfcand_erel"};
0144 
0145 const std::vector<std::string> DeepBoostedJetTagInfoProducer::particle_features_hlt_{"jet_pfcand_pt_log",
0146                                                                                      "jet_pfcand_energy_log",
0147                                                                                      "jet_pfcand_deta",
0148                                                                                      "jet_pfcand_dphi",
0149                                                                                      "jet_pfcand_eta",
0150                                                                                      "jet_pfcand_charge",
0151                                                                                      "jet_pfcand_frompv",
0152                                                                                      "jet_pfcand_nlostinnerhits",
0153                                                                                      "jet_pfcand_track_chi2",
0154                                                                                      "jet_pfcand_track_qual",
0155                                                                                      "jet_pfcand_dz",
0156                                                                                      "jet_pfcand_dzsig",
0157                                                                                      "jet_pfcand_dxy",
0158                                                                                      "jet_pfcand_dxysig",
0159                                                                                      "jet_pfcand_etarel",
0160                                                                                      "jet_pfcand_pperp_ratio",
0161                                                                                      "jet_pfcand_ppara_ratio",
0162                                                                                      "jet_pfcand_trackjet_d3d",
0163                                                                                      "jet_pfcand_trackjet_d3dsig",
0164                                                                                      "jet_pfcand_trackjet_dist",
0165                                                                                      "jet_pfcand_nhits",
0166                                                                                      "jet_pfcand_npixhits",
0167                                                                                      "jet_pfcand_nstriphits",
0168                                                                                      "jet_pfcand_trackjet_decayL",
0169                                                                                      "jet_pfcand_puppiw",
0170                                                                                      "pfcand_mask"};
0171 
0172 const std::vector<std::string> DeepBoostedJetTagInfoProducer::particle_features_scouting_{"pfcand_px",
0173                                                                                           "pfcand_py",
0174                                                                                           "pfcand_pz",
0175                                                                                           "pfcand_energy",
0176                                                                                           "pfcand_quality",
0177                                                                                           "pfcand_charge",
0178                                                                                           "pfcand_isEl",
0179                                                                                           "pfcand_isMu",
0180                                                                                           "pfcand_isChargedHad",
0181                                                                                           "pfcand_isGamma",
0182                                                                                           "pfcand_isNeutralHad",
0183                                                                                           "pfcand_phirel",
0184                                                                                           "pfcand_etarel",
0185                                                                                           "pfcand_deltaR",
0186                                                                                           "pfcand_abseta",
0187                                                                                           "pfcand_ptrel_log",
0188                                                                                           "pfcand_erel_log",
0189                                                                                           "pfcand_pt_log",
0190                                                                                           "pfcand_normchi2",
0191                                                                                           "pfcand_dz",
0192                                                                                           "pfcand_dxy",
0193                                                                                           "pfcand_dxysig",
0194                                                                                           "pfcand_btagEtaRel",
0195                                                                                           "pfcand_btagPtRatio",
0196                                                                                           "pfcand_btagPParRatio",
0197                                                                                           "pfcand_mask",
0198                                                                                           "pfcand_pt_log_nopuppi",
0199                                                                                           "pfcand_dzsig",
0200                                                                                           "pfcand_e_log_nopuppi",
0201                                                                                           "pfcand_ptrel",
0202                                                                                           "pfcand_erel",
0203                                                                                           "pfcand_lostInnerHits"};
0204 
0205 const std::vector<std::string> DeepBoostedJetTagInfoProducer::sv_features_{"sv_mask",
0206                                                                            "sv_ptrel",
0207                                                                            "sv_erel",
0208                                                                            "sv_phirel",
0209                                                                            "sv_etarel",
0210                                                                            "sv_deltaR",
0211                                                                            "sv_abseta",
0212                                                                            "sv_mass",
0213                                                                            "sv_ptrel_log",
0214                                                                            "sv_erel_log",
0215                                                                            "sv_pt_log",
0216                                                                            "sv_pt",
0217                                                                            "sv_ntracks",
0218                                                                            "sv_normchi2",
0219                                                                            "sv_dxy",
0220                                                                            "sv_dxysig",
0221                                                                            "sv_d3d",
0222                                                                            "sv_d3dsig",
0223                                                                            "sv_costhetasvpv"};
0224 
0225 const std::vector<std::string> DeepBoostedJetTagInfoProducer::sv_features_hlt_{"jet_sv_pt_log",
0226                                                                                "jet_sv_mass",
0227                                                                                "jet_sv_deta",
0228                                                                                "jet_sv_dphi",
0229                                                                                "jet_sv_eta",
0230                                                                                "jet_sv_ntrack",
0231                                                                                "jet_sv_chi2",
0232                                                                                "jet_sv_dxy",
0233                                                                                "jet_sv_dxysig",
0234                                                                                "jet_sv_d3d",
0235                                                                                "jet_sv_d3dsig",
0236                                                                                "sv_mask"};
0237 
0238 const float DeepBoostedJetTagInfoProducer::min_track_pt_property_ = 0.5;
0239 const int DeepBoostedJetTagInfoProducer::min_valid_pixel_hits_ = 0;
0240 
0241 DeepBoostedJetTagInfoProducer::DeepBoostedJetTagInfoProducer(const edm::ParameterSet &iConfig)
0242     : jet_radius_(iConfig.getParameter<double>("jet_radius")),
0243       min_jet_pt_(iConfig.getParameter<double>("min_jet_pt")),
0244       max_jet_eta_(iConfig.getParameter<double>("max_jet_eta")),
0245       min_pt_for_track_properties_(iConfig.getParameter<double>("min_pt_for_track_properties")),
0246       min_pt_for_pfcandidates_(iConfig.getParameter<double>("min_pt_for_pfcandidates")),
0247       use_puppiP4_(iConfig.getParameter<bool>("use_puppiP4")),
0248       min_puppi_wgt_(iConfig.getParameter<double>("min_puppi_wgt")),
0249       include_neutrals_(iConfig.getParameter<bool>("include_neutrals")),
0250       sort_by_sip2dsig_(iConfig.getParameter<bool>("sort_by_sip2dsig")),
0251       flip_ip_sign_(iConfig.getParameter<bool>("flip_ip_sign")),
0252       max_sip3dsig_(iConfig.getParameter<double>("sip3dSigMax")),
0253       use_hlt_features_(iConfig.getParameter<bool>("use_hlt_features")),
0254       jet_token_(consumes<edm::View<reco::Jet>>(iConfig.getParameter<edm::InputTag>("jets"))),
0255       vtx_token_(consumes<VertexCollection>(iConfig.getParameter<edm::InputTag>("vertices"))),
0256       sv_token_(consumes<SVCollection>(iConfig.getParameter<edm::InputTag>("secondary_vertices"))),
0257       pfcand_token_(consumes<CandidateView>(iConfig.getParameter<edm::InputTag>("pf_candidates"))),
0258       use_puppi_value_map_(false),
0259       use_pvasq_value_map_(false),
0260       use_unsubjet_map_(false),
0261       track_builder_token_(
0262           esConsumes<TransientTrackBuilder, TransientTrackRecord>(edm::ESInputTag("", "TransientTrackBuilder"))),
0263       covarianceVersion_(iConfig.getParameter<int>("covarianceVersion")),
0264       covariancePackingSchemas_(iConfig.getParameter<std::vector<int>>("covariancePackingSchemas")),
0265       use_scouting_features_(iConfig.getParameter<bool>("use_scouting_features")) {
0266   const auto &puppi_value_map_tag = iConfig.getParameter<edm::InputTag>("puppi_value_map");
0267   if (!puppi_value_map_tag.label().empty()) {
0268     puppi_value_map_token_ = consumes<edm::ValueMap<float>>(puppi_value_map_tag);
0269     use_puppi_value_map_ = true;
0270   } else if (use_puppiP4_) {
0271     throw edm::Exception(edm::errors::Configuration,
0272                          "puppi_value_map is not set but use_puppiP4 is set to True. Must also set puppi_value_map.");
0273   }
0274 
0275   const auto &pvas_tag = iConfig.getParameter<edm::InputTag>("vertex_associator");
0276   if (!pvas_tag.label().empty()) {
0277     pvasq_value_map_token_ = consumes<edm::ValueMap<int>>(pvas_tag);
0278     pvas_token_ = consumes<edm::Association<VertexCollection>>(pvas_tag);
0279     use_pvasq_value_map_ = true;
0280   }
0281 
0282   const auto &normchi2_value_map_tag = iConfig.getParameter<edm::InputTag>("normchi2_value_map");
0283   if (!normchi2_value_map_tag.label().empty()) {
0284     normchi2_value_map_token_ = consumes<edm::ValueMap<float>>(normchi2_value_map_tag);
0285   }
0286 
0287   const auto &dz_value_map_tag = iConfig.getParameter<edm::InputTag>("dz_value_map");
0288   if (!dz_value_map_tag.label().empty()) {
0289     dz_value_map_token_ = consumes<edm::ValueMap<float>>(dz_value_map_tag);
0290   }
0291 
0292   const auto &dxy_value_map_tag = iConfig.getParameter<edm::InputTag>("dxy_value_map");
0293   if (!dxy_value_map_tag.label().empty()) {
0294     dxy_value_map_token_ = consumes<edm::ValueMap<float>>(dxy_value_map_tag);
0295   }
0296 
0297   const auto &dzsig_value_map_tag = iConfig.getParameter<edm::InputTag>("dzsig_value_map");
0298   if (!dzsig_value_map_tag.label().empty()) {
0299     dzsig_value_map_token_ = consumes<edm::ValueMap<float>>(dzsig_value_map_tag);
0300   }
0301 
0302   const auto &dxysig_value_map_tag = iConfig.getParameter<edm::InputTag>("dxysig_value_map");
0303   if (!dxysig_value_map_tag.label().empty()) {
0304     dxysig_value_map_token_ = consumes<edm::ValueMap<float>>(dxysig_value_map_tag);
0305   }
0306 
0307   const auto &lostInnerHits_value_map_tag = iConfig.getParameter<edm::InputTag>("lostInnerHits_value_map");
0308   if (!lostInnerHits_value_map_tag.label().empty()) {
0309     lostInnerHits_value_map_token_ = consumes<edm::ValueMap<int>>(lostInnerHits_value_map_tag);
0310   }
0311 
0312   const auto &quality_value_map_tag = iConfig.getParameter<edm::InputTag>("quality_value_map");
0313   if (!quality_value_map_tag.label().empty()) {
0314     quality_value_map_token_ = consumes<edm::ValueMap<int>>(quality_value_map_tag);
0315   }
0316 
0317   const auto &trkPt_value_map_tag = iConfig.getParameter<edm::InputTag>("trkPt_value_map");
0318   if (!trkPt_value_map_tag.label().empty()) {
0319     trkPt_value_map_token_ = consumes<edm::ValueMap<float>>(trkPt_value_map_tag);
0320   }
0321 
0322   const auto &trkEta_value_map_tag = iConfig.getParameter<edm::InputTag>("trkEta_value_map");
0323   if (!trkEta_value_map_tag.label().empty()) {
0324     trkEta_value_map_token_ = consumes<edm::ValueMap<float>>(trkEta_value_map_tag);
0325   }
0326 
0327   const auto &trkPhi_value_map_tag = iConfig.getParameter<edm::InputTag>("trkPhi_value_map");
0328   if (!trkPhi_value_map_tag.label().empty()) {
0329     trkPhi_value_map_token_ = consumes<edm::ValueMap<float>>(trkPhi_value_map_tag);
0330   }
0331 
0332   const auto &unsubjet_map_tag = iConfig.getParameter<edm::InputTag>("unsubjet_map");
0333   if (!unsubjet_map_tag.label().empty()) {
0334     unsubjet_map_token_ = consumes<JetMatchMap>(unsubjet_map_tag);
0335     use_unsubjet_map_ = true;
0336   }
0337 
0338   produces<DeepBoostedJetTagInfoCollection>();
0339 }
0340 
0341 DeepBoostedJetTagInfoProducer::~DeepBoostedJetTagInfoProducer() {}
0342 
0343 void DeepBoostedJetTagInfoProducer::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0344   // pfDeepBoostedJetTagInfos
0345   edm::ParameterSetDescription desc;
0346   desc.add<double>("jet_radius", 0.8);
0347   desc.add<double>("min_jet_pt", 150);
0348   desc.add<double>("max_jet_eta", 99);
0349   desc.add<double>("min_pt_for_track_properties", -1);
0350   desc.add<double>("min_pt_for_pfcandidates", -1);
0351   desc.add<bool>("use_puppiP4", true);
0352   desc.add<bool>("include_neutrals", true);
0353   desc.add<bool>("sort_by_sip2dsig", false);
0354   desc.add<double>("min_puppi_wgt", 0.01);
0355   desc.add<bool>("flip_ip_sign", false);
0356   desc.add<double>("sip3dSigMax", -1);
0357   desc.add<bool>("use_hlt_features", false);
0358   desc.add<edm::InputTag>("vertices", edm::InputTag("offlinePrimaryVertices"));
0359   desc.add<edm::InputTag>("secondary_vertices", edm::InputTag("inclusiveCandidateSecondaryVertices"));
0360   desc.add<edm::InputTag>("pf_candidates", edm::InputTag("particleFlow"));
0361   desc.add<edm::InputTag>("jets", edm::InputTag("ak8PFJetsPuppi"));
0362   desc.add<edm::InputTag>("unsubjet_map", {});
0363   desc.add<edm::InputTag>("puppi_value_map", edm::InputTag("puppi"));
0364   desc.add<edm::InputTag>("vertex_associator", edm::InputTag("primaryVertexAssociation", "original"));
0365   desc.add<bool>("use_scouting_features", false);
0366   desc.add<edm::InputTag>("normchi2_value_map", edm::InputTag(""));
0367   desc.add<edm::InputTag>("dz_value_map", edm::InputTag(""));
0368   desc.add<edm::InputTag>("dxy_value_map", edm::InputTag(""));
0369   desc.add<edm::InputTag>("dzsig_value_map", edm::InputTag(""));
0370   desc.add<edm::InputTag>("dxysig_value_map", edm::InputTag(""));
0371   desc.add<edm::InputTag>("lostInnerHits_value_map", edm::InputTag(""));
0372   desc.add<edm::InputTag>("quality_value_map", edm::InputTag(""));
0373   desc.add<edm::InputTag>("trkPt_value_map", edm::InputTag(""));
0374   desc.add<edm::InputTag>("trkEta_value_map", edm::InputTag(""));
0375   desc.add<edm::InputTag>("trkPhi_value_map", edm::InputTag(""));
0376   desc.add<int>("covarianceVersion", 0)->setComment("so far: 0 is Phase0, 1 is Phase1");
0377   desc.add<std::vector<int>>("covariancePackingSchemas", {8, 264, 520, 776, 0});
0378   descriptions.add("pfDeepBoostedJetTagInfos", desc);
0379 }
0380 
0381 void DeepBoostedJetTagInfoProducer::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0382   // output collection
0383   auto output_tag_infos = std::make_unique<DeepBoostedJetTagInfoCollection>();
0384   // Input jets
0385   auto jets = iEvent.getHandle(jet_token_);
0386   auto unsubjet_map = use_unsubjet_map_ ? iEvent.getHandle(unsubjet_map_token_) : edm::Handle<JetMatchMap>();
0387   // Primary vertexes
0388   if (!use_scouting_features_) {
0389     iEvent.getByToken(vtx_token_, vtxs_);
0390     if (vtxs_->empty()) {
0391       // produce empty TagInfos in case no primary vertex
0392       iEvent.put(std::move(output_tag_infos));
0393       return;  // exit event
0394     }
0395     // Leading vertex
0396     pv_ = &vtxs_->at(0);
0397     // Secondary vertexs
0398     iEvent.getByToken(sv_token_, svs_);
0399     // Track builder
0400     track_builder_ = iSetup.getHandle(track_builder_token_);
0401   }
0402   // PF candidates
0403   iEvent.getByToken(pfcand_token_, pfcands_);
0404   is_packed_pf_candidate_collection_ = false;
0405   if (dynamic_cast<const pat::PackedCandidateCollection *>(pfcands_.product()))
0406     is_packed_pf_candidate_collection_ = true;
0407   // Puppi weight value map
0408   if (use_puppi_value_map_)
0409     iEvent.getByToken(puppi_value_map_token_, puppi_value_map_);
0410   // Vertex associator map
0411   if (use_pvasq_value_map_) {
0412     iEvent.getByToken(pvasq_value_map_token_, pvasq_value_map_);
0413     iEvent.getByToken(pvas_token_, pvas_);
0414   }
0415   if (use_scouting_features_) {
0416     iEvent.getByToken(normchi2_value_map_token_, normchi2_value_map_);
0417     iEvent.getByToken(dz_value_map_token_, dz_value_map_);
0418     iEvent.getByToken(dxy_value_map_token_, dxy_value_map_);
0419     iEvent.getByToken(dzsig_value_map_token_, dzsig_value_map_);
0420     iEvent.getByToken(dxysig_value_map_token_, dxysig_value_map_);
0421     iEvent.getByToken(lostInnerHits_value_map_token_, lostInnerHits_value_map_);
0422     iEvent.getByToken(quality_value_map_token_, quality_value_map_);
0423     iEvent.getByToken(trkPt_value_map_token_, trkPt_value_map_);
0424     iEvent.getByToken(trkEta_value_map_token_, trkEta_value_map_);
0425     iEvent.getByToken(trkPhi_value_map_token_, trkPhi_value_map_);
0426   }
0427 
0428   // Loop over jet
0429   for (std::size_t jet_n = 0; jet_n < jets->size(); jet_n++) {
0430     const auto &jet = (*jets)[jet_n];
0431     edm::RefToBase<reco::Jet> jet_ref(jets, jet_n);
0432     const auto &unsubJet =
0433         (use_unsubjet_map_ && (*unsubjet_map)[jet_ref].isNonnull()) ? *(*unsubjet_map)[jet_ref] : jet;
0434 
0435     // create jet features
0436     DeepBoostedJetFeatures features;
0437     if (use_hlt_features_) {
0438       // declare all the feature variables (init as empty vector)
0439       for (const auto &name : particle_features_hlt_) {
0440         features.add(name);
0441       }
0442       for (const auto &name : sv_features_hlt_) {
0443         features.add(name);
0444       }
0445     } else if (use_scouting_features_) {
0446       for (const auto &name : particle_features_scouting_) {
0447         features.add(name);
0448       }
0449     } else {
0450       for (const auto &name : particle_features_) {
0451         features.add(name);
0452       }
0453       for (const auto &name : sv_features_) {
0454         features.add(name);
0455       }
0456     }
0457 
0458     // fill values only if above pt threshold and has daughters, otherwise left
0459     bool fill_vars = true;
0460     if (jet.pt() < min_jet_pt_ or std::abs(jet.eta()) > max_jet_eta_) {
0461       fill_vars = false;
0462     }
0463     if (unsubJet.numberOfDaughters() == 0 and !use_scouting_features_) {
0464       fill_vars = false;
0465     }
0466 
0467     // fill features
0468     if (fill_vars) {
0469       fillParticleFeatures(features, unsubJet, jet);
0470       if (!use_scouting_features_) {
0471         fillSVFeatures(features, jet);
0472       }
0473       if (use_hlt_features_) {
0474         features.check_consistency(particle_features_hlt_);
0475         features.check_consistency(sv_features_hlt_);
0476       } else if (use_scouting_features_) {
0477         features.check_consistency(particle_features_scouting_);
0478       } else {
0479         features.check_consistency(particle_features_);
0480         features.check_consistency(sv_features_);
0481       }
0482     }
0483     // this should always be done even if features are not filled
0484     output_tag_infos->emplace_back(features, jet_ref);
0485   }
0486   // move output collection
0487   iEvent.put(std::move(output_tag_infos));
0488 }
0489 
0490 float DeepBoostedJetTagInfoProducer::puppiWgt(const reco::CandidatePtr &cand) {
0491   const auto *pack_cand = dynamic_cast<const pat::PackedCandidate *>(&(*cand));
0492   const auto *reco_cand = dynamic_cast<const reco::PFCandidate *>(&(*cand));
0493 
0494   //
0495   // Access puppi weight from ValueMap.
0496   //
0497   float wgt = 1.;  // Set to fallback value
0498 
0499   if (pack_cand) {
0500     if (use_puppi_value_map_)
0501       wgt = (*puppi_value_map_)[cand];
0502   } else if (reco_cand) {
0503     if (use_puppi_value_map_)
0504       wgt = (*puppi_value_map_)[cand];
0505   } else
0506     throw edm::Exception(edm::errors::InvalidReference)
0507         << "Cannot convert to either pat::PackedCandidate or reco::PFCandidate";
0508   puppi_wgt_cache[cand.key()] = wgt;
0509   return wgt;
0510 }
0511 
0512 bool DeepBoostedJetTagInfoProducer::useTrackProperties(const reco::PFCandidate *reco_cand) {
0513   const auto *track = reco_cand->bestTrack();
0514   return track != nullptr and track->pt() > min_pt_for_track_properties_;
0515 };
0516 
0517 void DeepBoostedJetTagInfoProducer::fillParticleFeatures(DeepBoostedJetFeatures &fts,
0518                                                          const reco::Jet &unsubJet,
0519                                                          const reco::Jet &jet) {
0520   // some jet properties
0521   math::XYZVector jet_dir = jet.momentum().Unit();
0522   TVector3 jet_direction(jet.momentum().Unit().x(), jet.momentum().Unit().y(), jet.momentum().Unit().z());
0523   GlobalVector jet_ref_track_dir(jet.px(), jet.py(), jet.pz());
0524   const float etasign = jet.eta() > 0 ? 1 : -1;
0525   std::unique_ptr<reco::VertexRefProd> PVRefProd;
0526   std::unique_ptr<TrackInfoBuilder> trackinfo;
0527   if (not use_scouting_features_) {
0528     // vertexes
0529     PVRefProd = std::make_unique<reco::VertexRefProd>(vtxs_);
0530     // track builder
0531     trackinfo = std::make_unique<TrackInfoBuilder>(track_builder_);
0532   }
0533 
0534   // make list of pf-candidates to be considered
0535   std::vector<reco::CandidatePtr> daughters;
0536   for (const auto &dau : unsubJet.daughterPtrVector()) {
0537     // remove particles w/ extremely low puppi weights
0538     // [Note] use jet daughters here to get the puppiWgt correctly
0539     if ((puppiWgt(dau)) < min_puppi_wgt_)
0540       continue;
0541     // from here: get the original reco/packed candidate not scaled by the puppi weight
0542     auto cand = pfcands_->ptrAt(dau.key());
0543     // base requirements on PF candidates
0544     if (use_hlt_features_ and cand->pt() < min_pt_for_pfcandidates_)
0545       continue;
0546     // charged candidate selection (for Higgs Interaction Net)
0547     if (!include_neutrals_ and (cand->charge() == 0 or cand->pt() < min_pt_for_track_properties_))
0548       continue;
0549     // only when computing the nagative tagger: remove charged candidates with high sip3d
0550     if (flip_ip_sign_ and cand->charge()) {
0551       (*trackinfo).buildTrackInfo(&(*cand), jet_dir, jet_ref_track_dir, *pv_);
0552       if ((*trackinfo).getTrackSip3dSig() > max_sip3dsig_)
0553         continue;
0554     }
0555     // filling daughters
0556     daughters.push_back(cand);
0557   }
0558 
0559   // Sorting of PF-candidates
0560   std::vector<btagbtvdeep::SortingClass<reco::CandidatePtr>> c_sorted;
0561   if (sort_by_sip2dsig_) {
0562     // sort charged pf candidates by 2d impact parameter significance
0563     for (const auto &cand : daughters) {
0564       (*trackinfo).buildTrackInfo(&(*cand), jet_dir, jet_ref_track_dir, *pv_);
0565       c_sorted.emplace_back(cand,
0566                             (*trackinfo).getTrackSip2dSig(),
0567                             -btagbtvdeep::mindrsvpfcand(*svs_, &(*cand), jet_radius_),
0568                             cand->pt() / jet.pt());
0569     }
0570     std::sort(c_sorted.begin(), c_sorted.end(), btagbtvdeep::SortingClass<reco::CandidatePtr>::compareByABCInv);
0571     for (unsigned int i = 0; i < c_sorted.size(); i++) {
0572       const auto &c = c_sorted.at(i);
0573       const auto &cand = c.get();
0574       daughters.at(i) = cand;
0575     }
0576   } else {
0577     // sort by Puppi-weighted pt
0578     if (use_puppiP4_)
0579       std::sort(daughters.begin(), daughters.end(), [&](const reco::CandidatePtr &a, const reco::CandidatePtr &b) {
0580         return puppi_wgt_cache.at(a.key()) * a->pt() > puppi_wgt_cache.at(b.key()) * b->pt();
0581       });
0582     // sort by original pt (not Puppi-weighted)
0583     else
0584       std::sort(daughters.begin(), daughters.end(), [](const auto &a, const auto &b) { return a->pt() > b->pt(); });
0585   }
0586 
0587   // reserve space
0588   if (use_hlt_features_) {
0589     for (const auto &name : particle_features_hlt_)
0590       fts.reserve(name, daughters.size());
0591   } else if (use_scouting_features_) {
0592     for (const auto &name : particle_features_scouting_)
0593       fts.reserve(name, daughters.size());
0594   } else {
0595     for (const auto &name : particle_features_)
0596       fts.reserve(name, daughters.size());
0597   }
0598 
0599   // build white list of candidates i.e. particles and tracks belonging to SV. Needed when input particles are not packed PF candidates
0600   std::vector<unsigned int> whiteListSV;
0601   std::vector<reco::TrackRef> whiteListTk;
0602   if (!use_scouting_features_) {
0603     if (not is_packed_pf_candidate_collection_) {
0604       for (size_t isv = 0; isv < svs_->size(); isv++) {
0605         for (size_t icand = 0; icand < svs_->at(isv).numberOfSourceCandidatePtrs(); icand++) {
0606           const edm::Ptr<reco::Candidate> &cand = svs_->at(isv).sourceCandidatePtr(icand);
0607           if (cand.id() == pfcands_.id())
0608             whiteListSV.push_back(cand.key());
0609         }
0610         for (auto cand = svs_->at(isv).begin(); cand != svs_->at(isv).end(); cand++) {
0611           const reco::RecoChargedCandidate *chCand = dynamic_cast<const reco::RecoChargedCandidate *>(&(*cand));
0612           if (chCand != nullptr) {
0613             whiteListTk.push_back(chCand->track());
0614           }
0615         }
0616       }
0617     }
0618   }
0619 
0620   // Build observables
0621   size_t icand = 0;
0622   for (const auto &cand : daughters) {
0623     const auto *packed_cand = dynamic_cast<const pat::PackedCandidate *>(&(*cand));
0624     const auto *reco_cand = dynamic_cast<const reco::PFCandidate *>(&(*cand));
0625 
0626     if (not packed_cand and not reco_cand)
0627       throw edm::Exception(edm::errors::InvalidReference)
0628           << "Cannot convert to either reco::PFCandidate or pat::PackedCandidate";
0629 
0630     if (!include_neutrals_ and
0631         ((packed_cand and !packed_cand->hasTrackDetails()) or (reco_cand and !useTrackProperties(reco_cand)))) {
0632       icand++;
0633       continue;
0634     }
0635 
0636     const float ip_sign = flip_ip_sign_ ? -1 : 1;
0637 
0638     // input particle is a packed PF candidate
0639     auto candP4 = use_puppiP4_ ? puppi_wgt_cache.at(cand.key()) * cand->p4() : cand->p4();
0640     auto candP3 = use_puppiP4_ ? puppi_wgt_cache.at(cand.key()) * cand->momentum() : cand->momentum();
0641 
0642     if (use_scouting_features_) {
0643       fts.fill("pfcand_px", candP4.px());
0644       fts.fill("pfcand_py", candP4.py());
0645       fts.fill("pfcand_pz", candP4.pz());
0646       fts.fill("pfcand_energy", candP4.energy());
0647       fts.fill("pfcand_charge", reco_cand->charge());
0648       fts.fill("pfcand_isEl", std::abs(reco_cand->pdgId()) == 11);
0649       fts.fill("pfcand_isMu", std::abs(reco_cand->pdgId()) == 13);
0650       fts.fill("pfcand_isChargedHad", std::abs(reco_cand->pdgId()) == 211);
0651       fts.fill("pfcand_isGamma", std::abs(reco_cand->pdgId()) == 22);
0652       fts.fill("pfcand_isNeutralHad", std::abs(reco_cand->pdgId()) == 130);
0653       fts.fill("pfcand_phirel", reco::deltaPhi(candP4, jet));
0654       fts.fill("pfcand_etarel", etasign * (candP4.eta() - jet.eta()));
0655       fts.fill("pfcand_deltaR", reco::deltaR(candP4, jet));
0656       fts.fill("pfcand_abseta", std::abs(candP4.eta()));
0657       fts.fill("pfcand_ptrel_log", std::log(candP4.pt() / jet.pt()));
0658       fts.fill("pfcand_ptrel", candP4.pt() / jet.pt());
0659       fts.fill("pfcand_erel_log", std::log(candP4.energy() / jet.energy()));
0660       fts.fill("pfcand_erel", candP4.energy() / jet.energy());
0661       fts.fill("pfcand_pt_log", std::log(candP4.pt()));
0662       fts.fill("pfcand_mask", 1);
0663       fts.fill("pfcand_pt_log_nopuppi", std::log(cand->pt()));
0664       fts.fill("pfcand_e_log_nopuppi", std::log(cand->energy()));
0665 
0666       // if normchi2 is 999.f then the pfcand has no "best track" and therefore no track variables
0667       if ((*normchi2_value_map_)[cand] > 900) {
0668         fts.fill("pfcand_normchi2", 0);
0669         fts.fill("pfcand_lostInnerHits", 0);
0670         fts.fill("pfcand_quality", 0);
0671         fts.fill("pfcand_dz", 0);
0672         fts.fill("pfcand_dzsig", 0);
0673         fts.fill("pfcand_dxy", 0);
0674         fts.fill("pfcand_dxysig", 0);
0675         fts.fill("pfcand_btagEtaRel", 0);
0676         fts.fill("pfcand_btagPtRatio", 0);
0677         fts.fill("pfcand_btagPParRatio", 0);
0678       } else {
0679         fts.fill("pfcand_normchi2", (*normchi2_value_map_)[cand]);
0680         fts.fill("pfcand_lostInnerHits", (*lostInnerHits_value_map_)[cand]);
0681         fts.fill("pfcand_quality", (*quality_value_map_)[cand]);
0682         fts.fill("pfcand_dz", (*dz_value_map_)[cand]);
0683         fts.fill("pfcand_dzsig", (*dzsig_value_map_)[cand]);
0684         fts.fill("pfcand_dxy", (*dxy_value_map_)[cand]);
0685         fts.fill("pfcand_dxysig", (*dxysig_value_map_)[cand]);
0686         float trk_px = (*trkPt_value_map_)[cand] * std::cos((*trkPhi_value_map_)[cand]);
0687         float trk_py = (*trkPt_value_map_)[cand] * std::sin((*trkPhi_value_map_)[cand]);
0688         float trk_pz = (*trkPt_value_map_)[cand] * std::sinh((*trkEta_value_map_)[cand]);
0689         math::XYZVector track_mom(trk_px, trk_py, trk_pz);
0690         TVector3 track_direction(trk_px, trk_py, trk_pz);
0691         double track_mag = sqrt(trk_px * trk_px + trk_py * trk_py + trk_pz * trk_pz);
0692         fts.fill("pfcand_btagEtaRel", reco::btau::etaRel(jet_dir, track_mom));
0693         fts.fill("pfcand_btagPtRatio", track_direction.Perp(jet_direction) / track_mag);
0694         fts.fill("pfcand_btagPParRatio", jet_dir.Dot(track_mom) / track_mag);
0695       }
0696     } else {
0697       // candidate track
0698       const reco::Track *track = nullptr;
0699       if (packed_cand)
0700         track = packed_cand->bestTrack();
0701       else if (reco_cand and useTrackProperties(reco_cand))
0702         track = reco_cand->bestTrack();
0703 
0704       // reco-vertex association
0705       int pv_ass_quality = 0;
0706       reco::VertexRef pv_ass;
0707       float vtx_ass = 0;
0708       if (reco_cand) {
0709         if (use_pvasq_value_map_) {
0710           pv_ass_quality = (*pvasq_value_map_)[cand];
0711           pv_ass = (*pvas_)[cand];
0712           vtx_ass = vtx_ass_from_pfcand(*reco_cand, pv_ass_quality, pv_ass);
0713         } else
0714           throw edm::Exception(edm::errors::InvalidReference) << "Vertex association missing";
0715       }
0716 
0717       // Building offline features
0718       if (not use_hlt_features_) {
0719         // in case of packed candidate
0720         if (packed_cand) {
0721           float hcal_fraction = 0.;
0722           if (packed_cand->pdgId() == 1 or packed_cand->pdgId() == 130)
0723             hcal_fraction = packed_cand->hcalFraction();
0724           else if (packed_cand->isIsolatedChargedHadron())
0725             hcal_fraction = packed_cand->rawHcalFraction();
0726 
0727           fts.fill("pfcand_hcalFrac", hcal_fraction);
0728           fts.fill("pfcand_VTX_ass", packed_cand->pvAssociationQuality());
0729           fts.fill("pfcand_lostInnerHits", packed_cand->lostInnerHits());
0730           fts.fill("pfcand_quality", track ? track->qualityMask() : 0);
0731           fts.fill("pfcand_charge", packed_cand->charge());
0732           fts.fill("pfcand_isEl", std::abs(packed_cand->pdgId()) == 11);
0733           fts.fill("pfcand_isMu", std::abs(packed_cand->pdgId()) == 13);
0734           fts.fill("pfcand_isChargedHad", std::abs(packed_cand->pdgId()) == 211);
0735           fts.fill("pfcand_isGamma", std::abs(packed_cand->pdgId()) == 22);
0736           fts.fill("pfcand_isNeutralHad", std::abs(packed_cand->pdgId()) == 130);
0737           fts.fill("pfcand_dz", ip_sign * packed_cand->dz());
0738           fts.fill("pfcand_dxy", ip_sign * packed_cand->dxy());
0739           fts.fill("pfcand_dzsig", track ? ip_sign * packed_cand->dz() / packed_cand->dzError() : 0);
0740           fts.fill("pfcand_dxysig", track ? ip_sign * packed_cand->dxy() / packed_cand->dxyError() : 0);
0741 
0742         }
0743         // in the case of reco candidate
0744         else if (reco_cand) {
0745           fts.fill("pfcand_hcalFrac", reco_cand->hcalEnergy() / (reco_cand->ecalEnergy() + reco_cand->hcalEnergy()));
0746           fts.fill("pfcand_VTX_ass", vtx_ass);
0747           fts.fill("pfcand_lostInnerHits", useTrackProperties(reco_cand) ? lost_inner_hits_from_pfcand(*reco_cand) : 0);
0748           fts.fill("pfcand_quality", useTrackProperties(reco_cand) ? quality_from_pfcand(*reco_cand) : 0);
0749           fts.fill("pfcand_charge", reco_cand->charge());
0750           fts.fill("pfcand_isEl", std::abs(reco_cand->pdgId()) == 11);
0751           fts.fill("pfcand_isMu", std::abs(reco_cand->pdgId()) == 13);
0752           fts.fill("pfcand_isChargedHad", std::abs(reco_cand->pdgId()) == 211);
0753           fts.fill("pfcand_isGamma", std::abs(reco_cand->pdgId()) == 22);
0754           fts.fill("pfcand_isNeutralHad", std::abs(reco_cand->pdgId()) == 130);
0755           fts.fill("pfcand_dz", track ? ip_sign * track->dz(pv_->position()) : 0);
0756           fts.fill("pfcand_dzsig", track ? ip_sign * track->dz(pv_->position()) / track->dzError() : 0);
0757           fts.fill("pfcand_dxy", track ? ip_sign * track->dxy(pv_->position()) : 0);
0758           fts.fill("pfcand_dxysig", track ? ip_sign * track->dxy(pv_->position()) / track->dxyError() : 0);
0759         }
0760 
0761         // generic candidate observables
0762         fts.fill("pfcand_puppiw", puppi_wgt_cache.at(cand.key()));
0763         fts.fill("pfcand_phirel", reco::deltaPhi(candP4, jet));
0764         fts.fill("pfcand_etarel", etasign * (candP4.eta() - jet.eta()));
0765         fts.fill("pfcand_deltaR", reco::deltaR(candP4, jet));
0766         fts.fill("pfcand_abseta", std::abs(candP4.eta()));
0767 
0768         fts.fill("pfcand_ptrel_log", std::log(candP4.pt() / jet.pt()));
0769         fts.fill("pfcand_ptrel", candP4.pt() / jet.pt());
0770         fts.fill("pfcand_erel_log", std::log(candP4.energy() / jet.energy()));
0771         fts.fill("pfcand_erel", candP4.energy() / jet.energy());
0772         fts.fill("pfcand_pt_log", std::log(candP4.pt()));
0773 
0774         fts.fill("pfcand_mask", 1);
0775         fts.fill("pfcand_pt_log_nopuppi", std::log(cand->pt()));
0776         fts.fill("pfcand_e_log_nopuppi", std::log(cand->energy()));
0777 
0778         float drminpfcandsv = btagbtvdeep::mindrsvpfcand(*svs_, &(*cand), std::numeric_limits<float>::infinity());
0779         fts.fill("pfcand_drminsv", drminpfcandsv);
0780 
0781         if (track) {
0782           auto cov = [&](unsigned i, unsigned j) { return track->covariance(i, j); };
0783           fts.fill("pfcand_dptdpt", cov(0, 0));
0784           fts.fill("pfcand_detadeta", cov(1, 1));
0785           fts.fill("pfcand_dphidphi", cov(2, 2));
0786           fts.fill("pfcand_dxydxy", cov(3, 3));
0787           fts.fill("pfcand_dzdz", cov(4, 4));
0788           fts.fill("pfcand_dxydz", cov(3, 4));
0789           fts.fill("pfcand_dphidxy", cov(2, 3));
0790           fts.fill("pfcand_dlambdadz", cov(1, 4));
0791 
0792           fts.fill("pfcand_normchi2", std::floor(track->normalizedChi2()));
0793 
0794           (*trackinfo).buildTrackInfo(&(*cand), jet_dir, jet_ref_track_dir, *pv_);
0795           fts.fill("pfcand_btagEtaRel", (*trackinfo).getTrackEtaRel());
0796           fts.fill("pfcand_btagPtRatio", (*trackinfo).getTrackPtRatio());
0797           fts.fill("pfcand_btagPParRatio", (*trackinfo).getTrackPParRatio());
0798           fts.fill("pfcand_btagSip2dVal", ip_sign * (*trackinfo).getTrackSip2dVal());
0799           fts.fill("pfcand_btagSip2dSig", ip_sign * (*trackinfo).getTrackSip2dSig());
0800           fts.fill("pfcand_btagSip3dVal", ip_sign * (*trackinfo).getTrackSip3dVal());
0801           fts.fill("pfcand_btagSip3dSig", ip_sign * (*trackinfo).getTrackSip3dSig());
0802           fts.fill("pfcand_btagJetDistVal", (*trackinfo).getTrackJetDistVal());
0803         } else {
0804           fts.fill("pfcand_normchi2", 999);
0805           fts.fill("pfcand_dptdpt", 0);
0806           fts.fill("pfcand_detadeta", 0);
0807           fts.fill("pfcand_dphidphi", 0);
0808           fts.fill("pfcand_dxydxy", 0);
0809           fts.fill("pfcand_dzdz", 0);
0810           fts.fill("pfcand_dxydz", 0);
0811           fts.fill("pfcand_dphidxy", 0);
0812           fts.fill("pfcand_dlambdadz", 0);
0813           fts.fill("pfcand_btagEtaRel", 0);
0814           fts.fill("pfcand_btagPtRatio", 0);
0815           fts.fill("pfcand_btagPParRatio", 0);
0816           fts.fill("pfcand_btagSip2dVal", 0);
0817           fts.fill("pfcand_btagSip2dSig", 0);
0818           fts.fill("pfcand_btagSip3dVal", 0);
0819           fts.fill("pfcand_btagSip3dSig", 0);
0820           fts.fill("pfcand_btagJetDistVal", 0);
0821         }
0822 
0823         // subjets only if the incomming jets is a PAT one
0824         const auto *patJet = dynamic_cast<const pat::Jet *>(&jet);
0825         if (patJet and patJet->nSubjetCollections() > 0) {
0826           auto subjets = patJet->subjets();
0827           std::sort(subjets.begin(), subjets.end(), [](const edm::Ptr<pat::Jet> &p1, const edm::Ptr<pat::Jet> &p2) {
0828             return p1->pt() > p2->pt();
0829           });
0830           fts.fill("pfcand_drsubjet1", !subjets.empty() ? reco::deltaR(*cand, *subjets.at(0)) : -1);
0831           fts.fill("pfcand_drsubjet2", subjets.size() > 1 ? reco::deltaR(*cand, *subjets.at(1)) : -1);
0832         } else {
0833           fts.fill("pfcand_drsubjet1", -1);
0834           fts.fill("pfcand_drsubjet2", -1);
0835         }
0836       }
0837       // using HLT features
0838       else {
0839         pat::PackedCandidate candidate;
0840         math::XYZPoint pv_ass_pos;
0841         // In case input is a packed candidate (evaluate HLT network on offline)
0842         if (packed_cand) {
0843           candidate = *packed_cand;
0844           pv_ass = reco::VertexRef(vtxs_, 0);
0845           pv_ass_pos = pv_ass->position();
0846         }
0847         // In case input is a reco::PFCandidate
0848         else if (reco_cand) {
0849           // follow what is done in PhysicsTools/PatAlgos/plugins/PATPackedCandidateProducer.cc to minimize differences between HLT and RECO choices of input observables
0850           // no reference to vertex, take closest in dz or position 0
0851           if (not pv_ass.isNonnull()) {
0852             if (track) {
0853               float z_dist = 99999;
0854               int pv_pos = -1;
0855               for (size_t iv = 0; iv < vtxs_->size(); iv++) {
0856                 float dz = std::abs(track->dz(((*vtxs_)[iv]).position()));
0857                 if (dz < z_dist) {
0858                   z_dist = dz;
0859                   pv_pos = iv;
0860                 }
0861               }
0862               pv_ass = reco::VertexRef(vtxs_, pv_pos);
0863             } else
0864               pv_ass = reco::VertexRef(vtxs_, 0);
0865           }
0866           pv_ass_pos = pv_ass->position();
0867 
0868           // create a transient packed candidate
0869           if (track) {
0870             candidate = pat::PackedCandidate(cand->polarP4(),
0871                                              track->referencePoint(),
0872                                              track->pt(),
0873                                              track->eta(),
0874                                              track->phi(),
0875                                              cand->pdgId(),
0876                                              (*PVRefProd),
0877                                              pv_ass.key());
0878             candidate.setAssociationQuality(pat::PackedCandidate::PVAssociationQuality(
0879                 btagbtvdeep::vtx_ass_from_pfcand(*reco_cand, pv_ass_quality, pv_ass)));
0880             candidate.setCovarianceVersion(covarianceVersion_);
0881             pat::PackedCandidate::LostInnerHits lostHits = pat::PackedCandidate::noLostInnerHits;
0882             int nlost = track->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
0883             if (nlost == 0) {
0884               if (track->hitPattern().hasValidHitInPixelLayer(PixelSubdetector::SubDetector::PixelBarrel, 1))
0885                 lostHits = pat::PackedCandidate::validHitInFirstPixelBarrelLayer;
0886             } else
0887               lostHits = (nlost == 1 ? pat::PackedCandidate::oneLostInnerHit : pat::PackedCandidate::moreLostInnerHits);
0888             candidate.setLostInnerHits(lostHits);
0889             candidate.setTrkAlgo(static_cast<uint8_t>(track->algo()), static_cast<uint8_t>(track->originalAlgo()));
0890 
0891             if (useTrackProperties(reco_cand) or
0892                 std::find(whiteListSV.begin(), whiteListSV.end(), icand) != whiteListSV.end() or
0893                 std::find(whiteListTk.begin(), whiteListTk.end(), reco_cand->trackRef()) != whiteListTk.end()) {
0894               candidate.setFirstHit(track->hitPattern().getHitPattern(reco::HitPattern::TRACK_HITS, 0));
0895               if (abs(cand->pdgId()) == 22)
0896                 candidate.setTrackProperties(*track, covariancePackingSchemas_[4], covarianceVersion_);
0897               else {
0898                 if (track->hitPattern().numberOfValidPixelHits() > min_valid_pixel_hits_)
0899                   candidate.setTrackProperties(*track, covariancePackingSchemas_[0], covarianceVersion_);
0900                 else
0901                   candidate.setTrackProperties(*track, covariancePackingSchemas_[1], covarianceVersion_);
0902               }
0903             } else {
0904               if (candidate.pt() > min_track_pt_property_) {
0905                 if (track->hitPattern().numberOfValidPixelHits() > 0)
0906                   candidate.setTrackProperties(*track, covariancePackingSchemas_[2], covarianceVersion_);
0907                 else
0908                   candidate.setTrackProperties(*track, covariancePackingSchemas_[3], covarianceVersion_);
0909               }
0910             }
0911             candidate.setTrackHighPurity(reco_cand->trackRef().isNonnull() and
0912                                          reco_cand->trackRef()->quality(reco::Track::highPurity));
0913           } else {
0914             candidate = pat::PackedCandidate(cand->polarP4(),
0915                                              pv_ass_pos,
0916                                              cand->pt(),
0917                                              cand->eta(),
0918                                              cand->phi(),
0919                                              cand->pdgId(),
0920                                              (*PVRefProd),
0921                                              pv_ass.key());
0922             candidate.setAssociationQuality(
0923                 pat::PackedCandidate::PVAssociationQuality(pat::PackedCandidate::UsedInFitTight));
0924           }
0925           /// override track
0926           track = candidate.bestTrack();
0927         }
0928 
0929         TVector3 cand_direction(candP3.x(), candP3.y(), candP3.z());
0930 
0931         fts.fill("jet_pfcand_pt_log", std::log(candP4.pt()));
0932         fts.fill("jet_pfcand_energy_log", std::log(candP4.energy()));
0933         fts.fill("jet_pfcand_eta", candP4.eta());
0934         fts.fill("jet_pfcand_deta", jet_direction.Eta() - cand_direction.Eta());
0935         fts.fill("jet_pfcand_dphi", jet_direction.DeltaPhi(cand_direction));
0936         fts.fill("jet_pfcand_charge", cand->charge());
0937         fts.fill("jet_pfcand_etarel", reco::btau::etaRel(jet_dir, candP3));
0938         fts.fill("jet_pfcand_pperp_ratio", jet_direction.Perp(cand_direction) / cand_direction.Mag());
0939         fts.fill("jet_pfcand_ppara_ratio", jet_direction.Dot(cand_direction) / cand_direction.Mag());
0940         fts.fill("jet_pfcand_frompv", candidate.fromPV());
0941         fts.fill("jet_pfcand_dz", candidate.dz(pv_ass_pos));
0942         fts.fill("jet_pfcand_dxy", candidate.dxy(pv_ass_pos));
0943         fts.fill("jet_pfcand_puppiw", puppi_wgt_cache.at(cand.key()));
0944         fts.fill("jet_pfcand_nlostinnerhits", candidate.lostInnerHits());
0945         fts.fill("jet_pfcand_nhits", candidate.numberOfHits());
0946         fts.fill("jet_pfcand_npixhits", candidate.numberOfPixelHits());
0947         fts.fill("jet_pfcand_nstriphits", candidate.stripLayersWithMeasurement());
0948         fts.fill("pfcand_mask", 1);
0949 
0950         if (track) {
0951           fts.fill("jet_pfcand_dzsig", fabs(candidate.dz(pv_ass_pos)) / candidate.dzError());
0952           fts.fill("jet_pfcand_dxysig", fabs(candidate.dxy(pv_ass_pos)) / candidate.dxyError());
0953           fts.fill("jet_pfcand_track_chi2", track->normalizedChi2());
0954           fts.fill("jet_pfcand_track_qual", track->qualityMask());
0955 
0956           reco::TransientTrack transientTrack = track_builder_->build(*track);
0957           Measurement1D meas_ip2d =
0958               IPTools::signedTransverseImpactParameter(transientTrack, jet_ref_track_dir, *pv_).second;
0959           Measurement1D meas_ip3d = IPTools::signedImpactParameter3D(transientTrack, jet_ref_track_dir, *pv_).second;
0960           Measurement1D meas_jetdist = IPTools::jetTrackDistance(transientTrack, jet_ref_track_dir, *pv_).second;
0961           Measurement1D meas_decayl = IPTools::signedDecayLength3D(transientTrack, jet_ref_track_dir, *pv_).second;
0962 
0963           fts.fill("jet_pfcand_trackjet_d3d", meas_ip3d.value());
0964           fts.fill("jet_pfcand_trackjet_d3dsig", fabs(meas_ip3d.significance()));
0965           fts.fill("jet_pfcand_trackjet_dist", -meas_jetdist.value());
0966           fts.fill("jet_pfcand_trackjet_decayL", meas_decayl.value());
0967         } else {
0968           fts.fill("jet_pfcand_dzsig", 0);
0969           fts.fill("jet_pfcand_dxysig", 0);
0970           fts.fill("jet_pfcand_track_chi2", 0);
0971           fts.fill("jet_pfcand_track_qual", 0);
0972           fts.fill("jet_pfcand_trackjet_d3d", 0);
0973           fts.fill("jet_pfcand_trackjet_d3dsig", 0);
0974           fts.fill("jet_pfcand_trackjet_dist", 0);
0975           fts.fill("jet_pfcand_trackjet_decayL", 0);
0976         }
0977       }
0978     }
0979     icand++;
0980   }
0981 }
0982 
0983 void DeepBoostedJetTagInfoProducer::fillSVFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet) {
0984   // secondary vertexes matching jet
0985   std::vector<const reco::VertexCompositePtrCandidate *> jetSVs;
0986   for (const auto &sv : *svs_) {
0987     if (reco::deltaR2(sv, jet) < jet_radius_ * jet_radius_) {
0988       jetSVs.push_back(&sv);
0989     }
0990   }
0991   // sort by dxy significance
0992   std::sort(jetSVs.begin(),
0993             jetSVs.end(),
0994             [&](const reco::VertexCompositePtrCandidate *sva, const reco::VertexCompositePtrCandidate *svb) {
0995               return sv_vertex_comparator(*sva, *svb, *pv_);
0996             });
0997 
0998   // reserve space
0999   if (not use_hlt_features_) {
1000     for (const auto &name : sv_features_) {
1001       fts.reserve(name, jetSVs.size());
1002     }
1003   } else {
1004     for (const auto &name : sv_features_hlt_) {
1005       fts.reserve(name, jetSVs.size());
1006     }
1007   }
1008 
1009   const float etasign = jet.eta() > 0 ? 1 : -1;
1010 
1011   GlobalVector jet_global_vec(jet.px(), jet.py(), jet.pz());
1012 
1013   for (const auto *sv : jetSVs) {
1014     // features for reco
1015     if (not use_hlt_features_) {
1016       fts.fill("sv_mask", 1);
1017       fts.fill("sv_phirel", reco::deltaPhi(*sv, jet));
1018       fts.fill("sv_etarel", etasign * (sv->eta() - jet.eta()));
1019       fts.fill("sv_deltaR", reco::deltaR(*sv, jet));
1020       fts.fill("sv_abseta", std::abs(sv->eta()));
1021       fts.fill("sv_mass", sv->mass());
1022 
1023       fts.fill("sv_ptrel_log", std::log(sv->pt() / jet.pt()));
1024       fts.fill("sv_ptrel", sv->pt() / jet.pt());
1025       fts.fill("sv_erel_log", std::log(sv->energy() / jet.energy()));
1026       fts.fill("sv_erel", sv->energy() / jet.energy());
1027       fts.fill("sv_pt_log", std::log(sv->pt()));
1028       fts.fill("sv_pt", sv->pt());
1029 
1030       fts.fill("sv_ntracks", sv->numberOfDaughters());
1031       fts.fill("sv_normchi2", sv->vertexNormalizedChi2());
1032       const auto &dxy = vertexDxy(*sv, *pv_);
1033       fts.fill("sv_dxy", dxy.value());
1034       fts.fill("sv_dxysig", dxy.significance());
1035       const auto &d3d = vertexD3d(*sv, *pv_);
1036       fts.fill("sv_d3d", d3d.value());
1037       fts.fill("sv_d3dsig", d3d.significance());
1038 
1039       fts.fill("sv_costhetasvpv", (flip_ip_sign_ ? -1.f : 1.f) * vertexDdotP(*sv, *pv_));
1040     } else {
1041       fts.fill("sv_mask", 1);
1042       fts.fill("jet_sv_pt_log", log(sv->pt()));
1043       fts.fill("jet_sv_eta", sv->eta());
1044       fts.fill("jet_sv_mass", sv->mass());
1045       fts.fill("jet_sv_deta", sv->eta() - jet.eta());
1046       fts.fill("jet_sv_dphi", sv->phi() - jet.phi());
1047       fts.fill("jet_sv_ntrack", sv->numberOfDaughters());
1048       fts.fill("jet_sv_chi2", sv->vertexNormalizedChi2());
1049 
1050       reco::Vertex::CovarianceMatrix csv;
1051       sv->fillVertexCovariance(csv);
1052       reco::Vertex svtx(sv->vertex(), csv);
1053 
1054       VertexDistanceXY dxy;
1055       auto valxy = dxy.signedDistance(svtx, *pv_, jet_global_vec);
1056       fts.fill("jet_sv_dxy", valxy.value());
1057       fts.fill("jet_sv_dxysig", fabs(valxy.significance()));
1058 
1059       VertexDistance3D d3d;
1060       auto val3d = d3d.signedDistance(svtx, *pv_, jet_global_vec);
1061       fts.fill("jet_sv_d3d", val3d.value());
1062       fts.fill("jet_sv_d3dsig", fabs(val3d.significance()));
1063     }
1064   }
1065 }
1066 
1067 // define this as a plug-in
1068 DEFINE_FWK_MODULE(DeepBoostedJetTagInfoProducer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!