Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:24:25 

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