Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2023-04-04 02:11:38

0001 #include <TVector3.h>
0002 #include <TTree.h>
0003 
0004 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0005 
0006 #include "FWCore/ServiceRegistry/interface/Service.h"
0007 #include "FWCore/Framework/interface/Frameworkfwd.h"
0008 #include "FWCore/Framework/interface/stream/EDProducer.h"
0009 #include "FWCore/Framework/interface/Event.h"
0010 #include "FWCore/Framework/interface/MakerMacros.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/Utilities/interface/StreamID.h"
0013 #include "FWCore/Utilities/interface/ESGetToken.h"
0014 
0015 #include "DataFormats/Candidate/interface/Candidate.h"
0016 #include "DataFormats/PatCandidates/interface/Jet.h"
0017 #include "DataFormats/PatCandidates/interface/Muon.h"
0018 #include "DataFormats/PatCandidates/interface/Electron.h"
0019 #include "DataFormats/PatCandidates/interface/Photon.h"
0020 #include "DataFormats/PatCandidates/interface/Tau.h"
0021 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0022 
0023 #include "RecoBTag/FeatureTools/interface/TrackInfoBuilder.h"
0024 #include "RecoBTag/FeatureTools/interface/deep_helpers.h"
0025 #include "RecoBTag/FeatureTools/interface/sorting_modules.h"
0026 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
0027 
0028 #include "DataFormats/Candidate/interface/VertexCompositePtrCandidate.h"
0029 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0030 #include "DataFormats/BTauReco/interface/DeepBoostedJetTagInfo.h"
0031 
0032 using namespace btagbtvdeep;
0033 
0034 class ParticleNetFeatureEvaluator : public edm::stream::EDProducer<> {
0035 public:
0036   explicit ParticleNetFeatureEvaluator(const edm::ParameterSet &);
0037   ~ParticleNetFeatureEvaluator() override;
0038 
0039   static void fillDescriptions(edm::ConfigurationDescriptions &descriptions);
0040 
0041 private:
0042   void beginStream(edm::StreamID) override {}
0043   void produce(edm::Event &, const edm::EventSetup &) override;
0044   void endStream() override {}
0045   void fillParticleFeatures(DeepBoostedJetFeatures &fts,
0046                             const reco::Jet &jet,
0047                             const std::vector<math::XYZTLorentzVector> &tau_pfcandidates,
0048                             const pat::MuonCollection &muons,
0049                             const pat::ElectronCollection &electrons,
0050                             const pat::PhotonCollection &photons);
0051   void fillSVFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet);
0052   void fillLostTrackFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet);
0053   bool useTrackProperties(const pat::PackedCandidate *cand);
0054 
0055   const double jet_radius_;
0056   const double min_jet_pt_;
0057   const double max_jet_eta_;
0058   const double min_jet_eta_;
0059   const double min_pt_for_track_properties_;
0060   const double min_pt_for_pfcandidates_;
0061   const double min_pt_for_losttrack_;
0062   const double max_dr_for_losttrack_;
0063   const double min_pt_for_taus_;
0064   const double max_eta_for_taus_;
0065   const bool include_neutrals_;
0066 
0067   edm::EDGetTokenT<pat::MuonCollection> muon_token_;
0068   edm::EDGetTokenT<pat::ElectronCollection> electron_token_;
0069   edm::EDGetTokenT<pat::PhotonCollection> photon_token_;
0070   edm::EDGetTokenT<pat::TauCollection> tau_token_;
0071   edm::EDGetTokenT<edm::View<reco::Jet>> jet_token_;
0072   edm::EDGetTokenT<pat::PackedCandidateCollection> losttrack_token_;
0073   edm::EDGetTokenT<reco::VertexCollection> vtx_token_;
0074   edm::EDGetTokenT<reco::VertexCompositePtrCandidateCollection> sv_token_;
0075   edm::EDGetTokenT<edm::View<reco::Candidate>> pfcand_token_;
0076   edm::ESGetToken<TransientTrackBuilder, TransientTrackRecord> track_builder_token_;
0077 
0078   edm::Handle<reco::VertexCollection> vtxs_;
0079   edm::Handle<reco::VertexCompositePtrCandidateCollection> svs_;
0080   edm::Handle<edm::View<reco::Candidate>> pfcands_;
0081   edm::Handle<pat::PackedCandidateCollection> losttracks_;
0082   edm::ESHandle<TransientTrackBuilder> track_builder_;
0083 
0084   const static std::vector<std::string> particle_features_;
0085   const static std::vector<std::string> sv_features_;
0086   const static std::vector<std::string> losttrack_features_;
0087   const reco::Vertex *pv_ = nullptr;
0088 
0089   TTree *tree;
0090   unsigned int event;
0091   float jet_pt;
0092   float jet_pt_raw;
0093   float jet_eta;
0094   float jet_phi;
0095   float jet_mass;
0096   unsigned int ijet;
0097   std::vector<float> jet_pfcand_pt_log;
0098   std::vector<float> jet_pfcand_energy_log;
0099   std::vector<float> jet_pfcand_deta;
0100   std::vector<float> jet_pfcand_dphi;
0101   std::vector<float> jet_pfcand_eta;
0102   std::vector<float> jet_pfcand_charge;
0103   std::vector<float> jet_pfcand_frompv;
0104   std::vector<float> jet_pfcand_nlostinnerhits;
0105   std::vector<float> jet_pfcand_track_chi2;
0106   std::vector<float> jet_pfcand_track_qual;
0107   std::vector<float> jet_pfcand_dz;
0108   std::vector<float> jet_pfcand_dzsig;
0109   std::vector<float> jet_pfcand_dxy;
0110   std::vector<float> jet_pfcand_dxysig;
0111   std::vector<float> jet_pfcand_etarel;
0112   std::vector<float> jet_pfcand_pperp_ratio;
0113   std::vector<float> jet_pfcand_ppara_ratio;
0114   std::vector<float> jet_pfcand_trackjet_d3d;
0115   std::vector<float> jet_pfcand_trackjet_d3dsig;
0116   std::vector<float> jet_pfcand_trackjet_dist;
0117   std::vector<float> jet_pfcand_npixhits;
0118   std::vector<float> jet_pfcand_nstriphits;
0119   std::vector<float> jet_pfcand_trackjet_decayL;
0120   std::vector<float> jet_pfcand_id;
0121   std::vector<float> jet_pfcand_calofraction;
0122   std::vector<float> jet_pfcand_hcalfraction;
0123   std::vector<float> jet_pfcand_puppiw;
0124   std::vector<float> jet_pfcand_muon_id;
0125   std::vector<float> jet_pfcand_muon_isglobal;
0126   std::vector<float> jet_pfcand_muon_chi2;
0127   std::vector<float> jet_pfcand_muon_segcomp;
0128   std::vector<float> jet_pfcand_muon_nvalidhit;
0129   std::vector<float> jet_pfcand_muon_nstation;
0130   std::vector<float> jet_pfcand_electron_detaIn;
0131   std::vector<float> jet_pfcand_electron_dphiIn;
0132   std::vector<float> jet_pfcand_electron_sigIetaIeta;
0133   std::vector<float> jet_pfcand_electron_sigIphiIphi;
0134   std::vector<float> jet_pfcand_electron_r9;
0135   std::vector<float> jet_pfcand_electron_convProb;
0136   std::vector<float> jet_pfcand_photon_sigIetaIeta;
0137   std::vector<float> jet_pfcand_photon_r9;
0138   std::vector<float> jet_pfcand_photon_eVeto;
0139   std::vector<float> jet_pfcand_tau_signal;
0140   std::vector<float> jet_sv_pt_log;
0141   std::vector<float> jet_sv_mass;
0142   std::vector<float> jet_sv_deta;
0143   std::vector<float> jet_sv_dphi;
0144   std::vector<float> jet_sv_eta;
0145   std::vector<float> jet_sv_ntrack;
0146   std::vector<float> jet_sv_chi2;
0147   std::vector<float> jet_sv_dxy;
0148   std::vector<float> jet_sv_dxysig;
0149   std::vector<float> jet_sv_d3d;
0150   std::vector<float> jet_sv_d3dsig;
0151   std::vector<float> jet_losttrack_pt_log;
0152   std::vector<float> jet_losttrack_eta;
0153   std::vector<float> jet_losttrack_deta;
0154   std::vector<float> jet_losttrack_dphi;
0155   std::vector<float> jet_losttrack_charge;
0156   std::vector<float> jet_losttrack_frompv;
0157   std::vector<float> jet_losttrack_track_chi2;
0158   std::vector<float> jet_losttrack_track_qual;
0159   std::vector<float> jet_losttrack_dz;
0160   std::vector<float> jet_losttrack_dxy;
0161   std::vector<float> jet_losttrack_dzsig;
0162   std::vector<float> jet_losttrack_dxysig;
0163   std::vector<float> jet_losttrack_etarel;
0164   std::vector<float> jet_losttrack_trackjet_d3d;
0165   std::vector<float> jet_losttrack_trackjet_d3dsig;
0166   std::vector<float> jet_losttrack_trackjet_dist;
0167   std::vector<float> jet_losttrack_trackjet_decayL;
0168   std::vector<float> jet_losttrack_npixhits;
0169   std::vector<float> jet_losttrack_nstriphits;
0170 };
0171 
0172 const std::vector<std::string> ParticleNetFeatureEvaluator::particle_features_{"jet_pfcand_pt_log",
0173                                                                                "jet_pfcand_energy_log",
0174                                                                                "jet_pfcand_deta",
0175                                                                                "jet_pfcand_dphi",
0176                                                                                "jet_pfcand_eta",
0177                                                                                "jet_pfcand_charge",
0178                                                                                "jet_pfcand_frompv",
0179                                                                                "jet_pfcand_nlostinnerhits",
0180                                                                                "jet_pfcand_track_chi2",
0181                                                                                "jet_pfcand_track_qual",
0182                                                                                "jet_pfcand_dz",
0183                                                                                "jet_pfcand_dzsig",
0184                                                                                "jet_pfcand_dxy",
0185                                                                                "jet_pfcand_dxysig",
0186                                                                                "jet_pfcand_etarel",
0187                                                                                "jet_pfcand_pperp_ratio",
0188                                                                                "jet_pfcand_ppara_ratio",
0189                                                                                "jet_pfcand_trackjet_d3d",
0190                                                                                "jet_pfcand_trackjet_d3dsig",
0191                                                                                "jet_pfcand_trackjet_dist",
0192                                                                                "jet_pfcand_nhits",
0193                                                                                "jet_pfcand_npixhits",
0194                                                                                "jet_pfcand_nstriphits",
0195                                                                                "jet_pfcand_trackjet_decayL",
0196                                                                                "jet_pfcand_id",
0197                                                                                "jet_pfcand_calofraction",
0198                                                                                "jet_pfcand_hcalfraction",
0199                                                                                "jet_pfcand_puppiw",
0200                                                                                "jet_pfcand_muon_id",
0201                                                                                "jet_pfcand_muon_isglobal",
0202                                                                                "jet_pfcand_muon_segcomp",
0203                                                                                "jet_pfcand_muon_chi2",
0204                                                                                "jet_pfcand_muon_nvalidhit",
0205                                                                                "jet_pfcand_muon_nstation",
0206                                                                                "jet_pfcand_electron_detaIn",
0207                                                                                "jet_pfcand_electron_dphiIn",
0208                                                                                "jet_pfcand_electron_sigIetaIeta",
0209                                                                                "jet_pfcand_electron_sigIphiIphi",
0210                                                                                "jet_pfcand_electron_r9",
0211                                                                                "jet_pfcand_electron_convProb",
0212                                                                                "jet_pfcand_photon_sigIetaIeta",
0213                                                                                "jet_pfcand_photon_r9",
0214                                                                                "jet_pfcand_photon_eVeto",
0215                                                                                "jet_pfcand_tau_signal",
0216                                                                                "pfcand_mask"};
0217 
0218 const std::vector<std::string> ParticleNetFeatureEvaluator::sv_features_{"jet_sv_pt_log",
0219                                                                          "jet_sv_mass",
0220                                                                          "jet_sv_deta",
0221                                                                          "jet_sv_dphi",
0222                                                                          "jet_sv_eta",
0223                                                                          "jet_sv_ntrack",
0224                                                                          "jet_sv_chi2",
0225                                                                          "jet_sv_dxy",
0226                                                                          "jet_sv_dxysig",
0227                                                                          "jet_sv_d3d",
0228                                                                          "jet_sv_d3dsig",
0229                                                                          "sv_mask"};
0230 
0231 const std::vector<std::string> ParticleNetFeatureEvaluator::losttrack_features_{"jet_losttrack_pt_log",
0232                                                                                 "jet_losttrack_eta",
0233                                                                                 "jet_losttrack_deta",
0234                                                                                 "jet_losttrack_dphi",
0235                                                                                 "jet_losttrack_charge",
0236                                                                                 "jet_losttrack_frompv",
0237                                                                                 "jet_losttrack_track_chi2",
0238                                                                                 "jet_losttrack_track_qual",
0239                                                                                 "jet_losttrack_dz",
0240                                                                                 "jet_losttrack_dxy",
0241                                                                                 "jet_losttrack_dzsig",
0242                                                                                 "jet_losttrack_dxysig",
0243                                                                                 "jet_losttrack_etarel",
0244                                                                                 "jet_losttrack_trackjet_d3d",
0245                                                                                 "jet_losttrack_trackjet_d3dsig",
0246                                                                                 "jet_losttrack_trackjet_dist",
0247                                                                                 "jet_losttrack_trackjet_decayL",
0248                                                                                 "jet_losttrack_npixhits",
0249                                                                                 "jet_losttrack_nstriphits",
0250                                                                                 "lt_mask"};
0251 
0252 ParticleNetFeatureEvaluator::ParticleNetFeatureEvaluator(const edm::ParameterSet &iConfig)
0253     : jet_radius_(iConfig.getParameter<double>("jet_radius")),
0254       min_jet_pt_(iConfig.getParameter<double>("min_jet_pt")),
0255       max_jet_eta_(iConfig.getParameter<double>("max_jet_eta")),
0256       min_jet_eta_(iConfig.getParameter<double>("min_jet_eta")),
0257       min_pt_for_track_properties_(iConfig.getParameter<double>("min_pt_for_track_properties")),
0258       min_pt_for_pfcandidates_(iConfig.getParameter<double>("min_pt_for_pfcandidates")),
0259       min_pt_for_losttrack_(iConfig.getParameter<double>("min_pt_for_losttrack")),
0260       max_dr_for_losttrack_(iConfig.getParameter<double>("max_dr_for_losttrack")),
0261       min_pt_for_taus_(iConfig.getParameter<double>("min_pt_for_taus")),
0262       max_eta_for_taus_(iConfig.getParameter<double>("max_eta_for_taus")),
0263       include_neutrals_(iConfig.getParameter<bool>("include_neutrals")),
0264       muon_token_(consumes<pat::MuonCollection>(iConfig.getParameter<edm::InputTag>("muons"))),
0265       electron_token_(consumes<pat::ElectronCollection>(iConfig.getParameter<edm::InputTag>("electrons"))),
0266       photon_token_(consumes<pat::PhotonCollection>(iConfig.getParameter<edm::InputTag>("photons"))),
0267       tau_token_(consumes<pat::TauCollection>(iConfig.getParameter<edm::InputTag>("taus"))),
0268       jet_token_(consumes<edm::View<reco::Jet>>(iConfig.getParameter<edm::InputTag>("jets"))),
0269       losttrack_token_(consumes<pat::PackedCandidateCollection>(iConfig.getParameter<edm::InputTag>("losttracks"))),
0270       vtx_token_(consumes<reco::VertexCollection>(iConfig.getParameter<edm::InputTag>("vertices"))),
0271       sv_token_(consumes<reco::VertexCompositePtrCandidateCollection>(
0272           iConfig.getParameter<edm::InputTag>("secondary_vertices"))),
0273       pfcand_token_(consumes<edm::View<reco::Candidate>>(iConfig.getParameter<edm::InputTag>("pf_candidates"))),
0274       track_builder_token_(
0275           esConsumes<TransientTrackBuilder, TransientTrackRecord>(edm::ESInputTag("", "TransientTrackBuilder"))) {
0276   produces<std::vector<reco::DeepBoostedJetTagInfo>>();
0277 }
0278 
0279 ParticleNetFeatureEvaluator::~ParticleNetFeatureEvaluator() {}
0280 
0281 void ParticleNetFeatureEvaluator::fillDescriptions(edm::ConfigurationDescriptions &descriptions) {
0282   // pfDeepBoostedJetTagInfos
0283   edm::ParameterSetDescription desc;
0284   desc.add<double>("jet_radius", 0.8);
0285   desc.add<double>("min_jet_pt", 150);
0286   desc.add<double>("max_jet_eta", 99);
0287   desc.add<double>("min_jet_eta", 0.0);
0288   desc.add<double>("min_pt_for_track_properties", -1);
0289   desc.add<double>("min_pt_for_pfcandidates", -1);
0290   desc.add<double>("min_pt_for_losttrack", 1);
0291   desc.add<double>("max_dr_for_losttrack", 0.4);
0292   desc.add<double>("min_pt_for_taus", 20.);
0293   desc.add<double>("max_eta_for_taus", 2.5);
0294   desc.add<bool>("include_neutrals", true);
0295   desc.add<edm::InputTag>("vertices", edm::InputTag("offlineSlimmedPrimaryVertices"));
0296   desc.add<edm::InputTag>("secondary_vertices", edm::InputTag("slimmedSecondaryVertices"));
0297   desc.add<edm::InputTag>("pf_candidates", edm::InputTag("packedPFCandidates"));
0298   desc.add<edm::InputTag>("losttracks", edm::InputTag("lostTracks"));
0299   desc.add<edm::InputTag>("jets", edm::InputTag("slimmedJetsAK8"));
0300   desc.add<edm::InputTag>("muons", edm::InputTag("slimmedMuons"));
0301   desc.add<edm::InputTag>("taus", edm::InputTag("slimmedTaus"));
0302   desc.add<edm::InputTag>("electrons", edm::InputTag("slimmedElectrons"));
0303   desc.add<edm::InputTag>("photons", edm::InputTag("slimmedPhotons"));
0304   descriptions.add("ParticleNetFeatureEvaluator", desc);
0305 }
0306 
0307 void ParticleNetFeatureEvaluator::produce(edm::Event &iEvent, const edm::EventSetup &iSetup) {
0308   // output collection
0309   auto output_tag_infos = std::make_unique<std::vector<reco::DeepBoostedJetTagInfo>>();
0310   // Input jets
0311   auto jets = iEvent.getHandle(jet_token_);
0312   // Input muons
0313   auto muons = iEvent.getHandle(muon_token_);
0314   // Input taus
0315   auto taus = iEvent.getHandle(tau_token_);
0316   // Input electrons
0317   auto electrons = iEvent.getHandle(electron_token_);
0318   // Input photons
0319   auto photons = iEvent.getHandle(photon_token_);
0320   // Input lost tracks
0321   iEvent.getByToken(losttrack_token_, losttracks_);
0322   // Primary vertexes
0323   iEvent.getByToken(vtx_token_, vtxs_);
0324   if (vtxs_->empty()) {
0325     // produce empty TagInfos in case no primary vertex
0326     iEvent.put(std::move(output_tag_infos));
0327     return;  // exit event
0328   }
0329   // Leading vertex
0330   pv_ = &vtxs_->at(0);
0331   // Secondary vertexs
0332   iEvent.getByToken(sv_token_, svs_);
0333   // PF candidates
0334   iEvent.getByToken(pfcand_token_, pfcands_);
0335   // Track builder
0336   track_builder_ = iSetup.getHandle(track_builder_token_);
0337 
0338   // tau signal candidates
0339   std::vector<math::XYZTLorentzVector> tau_pfcandidates;
0340   for (size_t itau = 0; itau < taus->size(); itau++) {
0341     if (taus->at(itau).pt() < min_pt_for_taus_)
0342       continue;
0343     if (fabs(taus->at(itau).eta()) > max_eta_for_taus_)
0344       continue;
0345     for (unsigned ipart = 0; ipart < taus->at(itau).signalCands().size(); ipart++) {
0346       const pat::PackedCandidate *pfcand =
0347           dynamic_cast<const pat::PackedCandidate *>(taus->at(itau).signalCands()[ipart].get());
0348       tau_pfcandidates.push_back(pfcand->p4());
0349     }
0350   }
0351 
0352   // Loop over jet
0353   for (std::size_t jet_n = 0; jet_n < jets->size(); jet_n++) {
0354     const auto &jet = (*jets)[jet_n];
0355     edm::RefToBase<reco::Jet> jet_ref(jets, jet_n);
0356 
0357     // create jet features
0358     DeepBoostedJetFeatures features;
0359     for (const auto &name : particle_features_)
0360       features.add(name);
0361     for (const auto &name : sv_features_)
0362       features.add(name);
0363 
0364     // fill values only if above pt threshold and has daughters, otherwise left
0365     bool fill_vars = true;
0366     if ((jet.pt() < min_jet_pt_ and
0367          dynamic_cast<const pat::Jet *>(&jet)->correctedJet("Uncorrected").pt() < min_jet_pt_) or
0368         std::abs(jet.eta()) >= max_jet_eta_ or std::abs(jet.eta()) < min_jet_eta_)
0369       fill_vars = false;
0370     if (jet.numberOfDaughters() == 0)
0371       fill_vars = false;
0372 
0373     // fill features
0374     if (fill_vars) {
0375       fillParticleFeatures(features, jet, tau_pfcandidates, *muons, *electrons, *photons);
0376       fillSVFeatures(features, jet);
0377       fillLostTrackFeatures(features, jet);
0378       features.check_consistency(particle_features_);
0379       features.check_consistency(sv_features_);
0380       features.check_consistency(losttrack_features_);
0381     }
0382 
0383     // this should always be done even if features are not filled
0384     output_tag_infos->emplace_back(features, jet_ref);
0385   }
0386   // move output collection
0387   iEvent.put(std::move(output_tag_infos));
0388 }
0389 
0390 bool ParticleNetFeatureEvaluator::useTrackProperties(const pat::PackedCandidate *cand) {
0391   const auto *track = cand->bestTrack();
0392   return track != nullptr and track->pt() > min_pt_for_track_properties_;
0393 };
0394 
0395 void ParticleNetFeatureEvaluator::fillParticleFeatures(DeepBoostedJetFeatures &fts,
0396                                                        const reco::Jet &jet,
0397                                                        const std::vector<math::XYZTLorentzVector> &tau_pfcandidates,
0398                                                        const pat::MuonCollection &muons,
0399                                                        const pat::ElectronCollection &electrons,
0400                                                        const pat::PhotonCollection &photons) {
0401   // some jet properties
0402   math::XYZVector jet_dir = jet.momentum().Unit();
0403   TVector3 jet_direction(jet.momentum().Unit().x(), jet.momentum().Unit().y(), jet.momentum().Unit().z());
0404   GlobalVector jet_ref_track_dir(jet.px(), jet.py(), jet.pz());
0405   // vertexes
0406   reco::VertexRefProd PVRefProd(vtxs_);
0407   // track builder
0408   TrackInfoBuilder trackinfo(track_builder_);
0409 
0410   // make list of pf-candidates to be considered
0411   std::vector<const pat::PackedCandidate *> daughters;
0412   for (const auto &dau : jet.daughterPtrVector()) {
0413     // remove particles w/ extremely low puppi weights
0414     const pat::PackedCandidate *cand = dynamic_cast<const pat::PackedCandidate *>(&(*dau));
0415     if (not cand)
0416       throw edm::Exception(edm::errors::InvalidReference) << "Cannot convert to either pat::PackedCandidate";
0417     // base requirements on PF candidates
0418     if (cand->pt() < min_pt_for_pfcandidates_)
0419       continue;
0420     // charged candidate selection (for Higgs Interaction Net)
0421     if (!include_neutrals_ and (cand->charge() == 0 or cand->pt() < min_pt_for_track_properties_))
0422       continue;
0423     // filling daughters
0424     daughters.push_back(cand);
0425   }
0426 
0427   // sort by original pt (not Puppi-weighted)
0428   std::sort(daughters.begin(), daughters.end(), [](const auto &a, const auto &b) { return a->pt() > b->pt(); });
0429 
0430   // reserve space
0431   for (const auto &name : particle_features_)
0432     fts.reserve(name, daughters.size());
0433 
0434   // Build observables
0435   for (const auto &cand : daughters) {
0436     if (!include_neutrals_ and !useTrackProperties(cand))
0437       continue;
0438 
0439     // input particle is a packed PF candidate
0440     auto candP4 = cand->p4();
0441     auto candP3 = cand->momentum();
0442 
0443     // candidate track
0444     const reco::Track *track = nullptr;
0445     if (useTrackProperties(cand))
0446       track = cand->bestTrack();
0447 
0448     // reco-vertex association
0449     reco::VertexRef pv_ass = reco::VertexRef(vtxs_, 0);
0450     math::XYZPoint pv_ass_pos = pv_ass->position();
0451 
0452     TVector3 cand_direction(candP3.x(), candP3.y(), candP3.z());
0453 
0454     fts.fill("jet_pfcand_pt_log", std::isnan(std::log(candP4.pt())) ? 0 : std::log(candP4.pt()));
0455     fts.fill("jet_pfcand_energy_log", std::isnan(std::log(candP4.energy())) ? 0 : std::log(candP4.energy()));
0456     fts.fill("jet_pfcand_eta", candP4.eta());
0457     fts.fill("jet_pfcand_deta", jet_direction.Eta() - cand_direction.Eta());
0458     fts.fill("jet_pfcand_dphi", jet_direction.DeltaPhi(cand_direction));
0459     fts.fill("jet_pfcand_charge", cand->charge());
0460     fts.fill("jet_pfcand_etarel",
0461              std::isnan(reco::btau::etaRel(jet_dir, candP3)) ? 0 : reco::btau::etaRel(jet_dir, candP3));
0462     fts.fill("jet_pfcand_pperp_ratio",
0463              std::isnan(jet_direction.Perp(cand_direction) / cand_direction.Mag())
0464                  ? 0
0465                  : jet_direction.Perp(cand_direction) / cand_direction.Mag());
0466     fts.fill("jet_pfcand_ppara_ratio",
0467              std::isnan(jet_direction.Dot(cand_direction) / cand_direction.Mag())
0468                  ? 0
0469                  : jet_direction.Dot(cand_direction) / cand_direction.Mag());
0470     fts.fill("jet_pfcand_frompv", cand->fromPV());
0471     fts.fill("jet_pfcand_dz", std::isnan(cand->dz(pv_ass_pos)) ? 0 : cand->dz(pv_ass_pos));
0472     fts.fill("jet_pfcand_dxy", std::isnan(cand->dxy(pv_ass_pos)) ? 0 : cand->dxy(pv_ass_pos));
0473     fts.fill("jet_pfcand_puppiw", cand->puppiWeight());
0474     fts.fill("jet_pfcand_nlostinnerhits", cand->lostInnerHits());
0475     fts.fill("jet_pfcand_nhits", cand->numberOfHits());
0476     fts.fill("jet_pfcand_npixhits", cand->numberOfPixelHits());
0477     fts.fill("jet_pfcand_nstriphits", cand->stripLayersWithMeasurement());
0478 
0479     if (abs(cand->pdgId()) == 11 and cand->charge() != 0)
0480       fts.fill("jet_pfcand_id", 0);
0481     else if (abs(cand->pdgId()) == 13 and cand->charge() != 0)
0482       fts.fill("jet_pfcand_id", 1);
0483     else if (abs(cand->pdgId()) == 22 and cand->charge() == 0)
0484       fts.fill("jet_pfcand_id", 2);
0485     else if (abs(cand->pdgId()) != 22 and cand->charge() == 0 and abs(cand->pdgId()) != 1 and abs(cand->pdgId()) != 2)
0486       fts.fill("jet_pfcand_id", 3);
0487     else if (abs(cand->pdgId()) != 11 and abs(cand->pdgId()) != 13 and cand->charge() != 0)
0488       fts.fill("jet_pfcand_id", 4);
0489     else if (cand->charge() == 0 and abs(cand->pdgId()) == 1)
0490       fts.fill("jet_pfcand_id", 5);
0491     else if (cand->charge() == 0 and abs(cand->pdgId()) == 2)
0492       fts.fill("jet_pfcand_id", 6);
0493     else
0494       fts.fill("jet_pfcand_id", -1);
0495 
0496     fts.fill("jet_pfcand_hcalfraction", std::isnan(cand->hcalFraction()) ? 0 : cand->hcalFraction());
0497     fts.fill("jet_pfcand_calofraction", std::isnan(cand->caloFraction()) ? 0 : cand->caloFraction());
0498     fts.fill("pfcand_mask", 1);
0499 
0500     if (track) {
0501       fts.fill(
0502           "jet_pfcand_dzsig",
0503           std::isnan(fabs(cand->dz(pv_ass_pos)) / cand->dzError()) ? 0 : fabs(cand->dz(pv_ass_pos)) / cand->dzError());
0504       fts.fill("jet_pfcand_dxysig",
0505                std::isnan(fabs(cand->dxy(pv_ass_pos)) / cand->dxyError())
0506                    ? 0
0507                    : fabs(cand->dxy(pv_ass_pos)) / cand->dxyError());
0508       fts.fill("jet_pfcand_track_chi2", track->normalizedChi2());
0509       fts.fill("jet_pfcand_track_qual", track->qualityMask());
0510 
0511       reco::TransientTrack transientTrack = track_builder_->build(*track);
0512       Measurement1D meas_ip2d =
0513           IPTools::signedTransverseImpactParameter(transientTrack, jet_ref_track_dir, *pv_).second;
0514       Measurement1D meas_ip3d = IPTools::signedImpactParameter3D(transientTrack, jet_ref_track_dir, *pv_).second;
0515       Measurement1D meas_jetdist = IPTools::jetTrackDistance(transientTrack, jet_ref_track_dir, *pv_).second;
0516       Measurement1D meas_decayl = IPTools::signedDecayLength3D(transientTrack, jet_ref_track_dir, *pv_).second;
0517 
0518       fts.fill("jet_pfcand_trackjet_d3d", std::isnan(meas_ip3d.value()) ? 0 : meas_ip3d.value());
0519       fts.fill("jet_pfcand_trackjet_d3dsig",
0520                std::isnan(fabs(meas_ip3d.significance())) ? 0 : fabs(meas_ip3d.significance()));
0521       fts.fill("jet_pfcand_trackjet_dist", std::isnan(-meas_jetdist.value()) ? 0 : -meas_jetdist.value());
0522       fts.fill("jet_pfcand_trackjet_decayL", std::isnan(meas_decayl.value()) ? 0 : meas_decayl.value());
0523     } else {
0524       fts.fill("jet_pfcand_dzsig", 0);
0525       fts.fill("jet_pfcand_dxysig", 0);
0526       fts.fill("jet_pfcand_track_chi2", 0);
0527       fts.fill("jet_pfcand_track_qual", 0);
0528       fts.fill("jet_pfcand_trackjet_d3d", 0);
0529       fts.fill("jet_pfcand_trackjet_d3dsig", 0);
0530       fts.fill("jet_pfcand_trackjet_dist", 0);
0531       fts.fill("jet_pfcand_trackjet_decayL", 0);
0532     }
0533 
0534     // muons specific
0535     if (abs(cand->pdgId()) == 13) {
0536       std::vector<unsigned int> muonsToSkip;
0537       int ipos = -1;
0538       float minDR = 1000;
0539       for (size_t i = 0; i < muons.size(); i++) {
0540         if (not muons[i].isPFMuon())
0541           continue;
0542         if (std::find(muonsToSkip.begin(), muonsToSkip.end(), i) != muonsToSkip.end())
0543           continue;
0544         float dR = reco::deltaR(muons[i].p4(), candP4);
0545         if (dR < jet_radius_ and dR < minDR) {
0546           minDR = dR;
0547           ipos = i;
0548           muonsToSkip.push_back(i);
0549         }
0550       }
0551       if (ipos >= 0) {
0552         int muonId = 0;
0553         if (muons[ipos].passed(reco::Muon::CutBasedIdLoose))
0554           muonId++;
0555         if (muons[ipos].passed(reco::Muon::CutBasedIdMedium))
0556           muonId++;
0557         if (muons[ipos].passed(reco::Muon::CutBasedIdTight))
0558           muonId++;
0559         if (muons[ipos].passed(reco::Muon::CutBasedIdGlobalHighPt))
0560           muonId++;
0561         if (muons[ipos].passed(reco::Muon::CutBasedIdTrkHighPt))
0562           muonId++;
0563         fts.fill("jet_pfcand_muon_id", muonId);
0564         fts.fill("jet_pfcand_muon_isglobal", muons[ipos].isGlobalMuon());
0565         fts.fill("jet_pfcand_muon_chi2",
0566                  (muons[ipos].isGlobalMuon()) ? muons[ipos].globalTrack()->normalizedChi2() : 0);
0567         fts.fill("jet_pfcand_muon_nvalidhit",
0568                  (muons[ipos].isGlobalMuon()) ? muons[ipos].globalTrack()->hitPattern().numberOfValidMuonHits() : 0);
0569         fts.fill("jet_pfcand_muon_nstation", muons[ipos].numberOfMatchedStations());
0570         fts.fill("jet_pfcand_muon_segcomp", muon::segmentCompatibility(muons[ipos]));
0571       } else {
0572         fts.fill("jet_pfcand_muon_id", 0);
0573         fts.fill("jet_pfcand_muon_isglobal", 0);
0574         fts.fill("jet_pfcand_muon_chi2", 0);
0575         fts.fill("jet_pfcand_muon_nvalidhit", 0);
0576         fts.fill("jet_pfcand_muon_nstation", 0);
0577         fts.fill("jet_pfcand_muon_segcomp", 0);
0578       }
0579     } else {
0580       fts.fill("jet_pfcand_muon_id", 0);
0581       fts.fill("jet_pfcand_muon_isglobal", 0);
0582       fts.fill("jet_pfcand_muon_chi2", 0);
0583       fts.fill("jet_pfcand_muon_nvalidhit", 0);
0584       fts.fill("jet_pfcand_muon_nstation", 0);
0585       fts.fill("jet_pfcand_muon_segcomp", 0);
0586     }
0587 
0588     // electrons specific
0589     if (abs(cand->pdgId()) == 11) {
0590       int ipos = -1;
0591       for (size_t i = 0; i < electrons.size(); i++) {
0592         if (electrons[i].isPF()) {
0593           for (const auto &element : electrons[i].associatedPackedPFCandidates()) {
0594             if (abs(element->pdgId()) == 11 and element->p4() == candP4)
0595               ipos = i;
0596           }
0597         }
0598       }
0599       if (ipos >= 0) {
0600         fts.fill("jet_pfcand_electron_detaIn",
0601                  std::isnan(electrons[ipos].deltaEtaSuperClusterTrackAtVtx())
0602                      ? 0
0603                      : electrons[ipos].deltaEtaSuperClusterTrackAtVtx());
0604         fts.fill("jet_pfcand_electron_dphiIn",
0605                  std::isnan(electrons[ipos].deltaPhiSuperClusterTrackAtVtx())
0606                      ? 0
0607                      : electrons[ipos].deltaPhiSuperClusterTrackAtVtx());
0608         fts.fill("jet_pfcand_electron_sigIetaIeta",
0609                  std::isnan(electrons[ipos].full5x5_sigmaIetaIeta()) ? 0 : electrons[ipos].full5x5_sigmaIetaIeta());
0610         fts.fill("jet_pfcand_electron_sigIphiIphi",
0611                  std::isnan(electrons[ipos].full5x5_sigmaIphiIphi()) ? 0 : electrons[ipos].full5x5_sigmaIphiIphi());
0612         fts.fill("jet_pfcand_electron_r9", std::isnan(electrons[ipos].full5x5_r9()) ? 0 : electrons[ipos].full5x5_r9());
0613         fts.fill("jet_pfcand_electron_convProb",
0614                  std::isnan(electrons[ipos].convVtxFitProb()) ? 0 : electrons[ipos].convVtxFitProb());
0615       } else {
0616         fts.fill("jet_pfcand_electron_detaIn", 0);
0617         fts.fill("jet_pfcand_electron_dphiIn", 0);
0618         fts.fill("jet_pfcand_electron_sigIetaIeta", 0);
0619         fts.fill("jet_pfcand_electron_sigIphiIphi", 0);
0620         fts.fill("jet_pfcand_electron_r9", 0);
0621         fts.fill("jet_pfcand_electron_convProb", 0);
0622       }
0623     } else {
0624       fts.fill("jet_pfcand_electron_detaIn", 0);
0625       fts.fill("jet_pfcand_electron_dphiIn", 0);
0626       fts.fill("jet_pfcand_electron_sigIetaIeta", 0);
0627       fts.fill("jet_pfcand_electron_sigIphiIphi", 0);
0628       fts.fill("jet_pfcand_electron_r9", 0);
0629       fts.fill("jet_pfcand_electron_convProb", 0);
0630     }
0631 
0632     // photons specific
0633     if (abs(cand->pdgId()) == 22) {
0634       int ipos = -1;
0635       for (size_t i = 0; i < photons.size(); i++) {
0636         for (const auto &element : photons[i].associatedPackedPFCandidates()) {
0637           if (abs(element->pdgId()) == 22 and element->p4() == candP4)
0638             ipos = i;
0639         }
0640       }
0641       if (ipos >= 0) {
0642         fts.fill("jet_pfcand_photon_sigIetaIeta",
0643                  std::isnan(photons[ipos].full5x5_sigmaIetaIeta()) ? 0 : photons[ipos].full5x5_sigmaIetaIeta());
0644         fts.fill("jet_pfcand_photon_r9", std::isnan(photons[ipos].full5x5_r9()) ? 0 : photons[ipos].full5x5_r9());
0645         fts.fill("jet_pfcand_photon_eVeto", photons[ipos].passElectronVeto());
0646       } else {
0647         fts.fill("jet_pfcand_photon_sigIetaIeta", 0);
0648         fts.fill("jet_pfcand_photon_r9", 0);
0649         fts.fill("jet_pfcand_photon_eVeto", 0);
0650       }
0651     } else {
0652       fts.fill("jet_pfcand_photon_sigIetaIeta", 0);
0653       fts.fill("jet_pfcand_photon_r9", 0);
0654       fts.fill("jet_pfcand_photon_eVeto", 0);
0655     }
0656 
0657     // tau specific prior to any puppi weight application
0658     if (std::find(tau_pfcandidates.begin(), tau_pfcandidates.end(), cand->p4()) != tau_pfcandidates.end())
0659       fts.fill("jet_pfcand_tau_signal", 1);
0660     else
0661       fts.fill("jet_pfcand_tau_signal", 0);
0662   }
0663 }
0664 
0665 void ParticleNetFeatureEvaluator::fillSVFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet) {
0666   // secondary vertexes matching jet
0667   std::vector<const reco::VertexCompositePtrCandidate *> jetSVs;
0668   for (const auto &sv : *svs_) {
0669     if (reco::deltaR2(sv, jet) < jet_radius_ * jet_radius_) {
0670       jetSVs.push_back(&sv);
0671     }
0672   }
0673 
0674   // sort by dxy significance
0675   std::sort(jetSVs.begin(),
0676             jetSVs.end(),
0677             [&](const reco::VertexCompositePtrCandidate *sva, const reco::VertexCompositePtrCandidate *svb) {
0678               return sv_vertex_comparator(*sva, *svb, *pv_);
0679             });
0680 
0681   // reserve space
0682   for (const auto &name : sv_features_)
0683     fts.reserve(name, jetSVs.size());
0684 
0685   GlobalVector jet_global_vec(jet.px(), jet.py(), jet.pz());
0686 
0687   for (const auto *sv : jetSVs) {
0688     fts.fill("sv_mask", 1);
0689     fts.fill("jet_sv_pt_log", std::isnan(std::log(sv->pt())) ? 0 : std::log(sv->pt()));
0690     fts.fill("jet_sv_eta", sv->eta());
0691     fts.fill("jet_sv_mass", sv->mass());
0692     fts.fill("jet_sv_deta", sv->eta() - jet.eta());
0693     fts.fill("jet_sv_dphi", sv->phi() - jet.phi());
0694     fts.fill("jet_sv_ntrack", sv->numberOfDaughters());
0695     fts.fill("jet_sv_chi2", sv->vertexNormalizedChi2());
0696 
0697     reco::Vertex::CovarianceMatrix csv;
0698     sv->fillVertexCovariance(csv);
0699     reco::Vertex svtx(sv->vertex(), csv);
0700 
0701     VertexDistanceXY dxy;
0702     auto valxy = dxy.signedDistance(svtx, *pv_, jet_global_vec);
0703     fts.fill("jet_sv_dxy", std::isnan(valxy.value()) ? 0 : valxy.value());
0704     fts.fill("jet_sv_dxysig", std::isnan(fabs(valxy.significance())) ? 0 : fabs(valxy.significance()));
0705 
0706     VertexDistance3D d3d;
0707     auto val3d = d3d.signedDistance(svtx, *pv_, jet_global_vec);
0708     fts.fill("jet_sv_d3d", std::isnan(val3d.value()) ? 0 : val3d.value());
0709     fts.fill("jet_sv_d3dsig", std::isnan(fabs(val3d.significance())) ? 0 : fabs(val3d.significance()));
0710   }
0711 }
0712 
0713 void ParticleNetFeatureEvaluator::fillLostTrackFeatures(DeepBoostedJetFeatures &fts, const reco::Jet &jet) {
0714   // some jet properties
0715   TVector3 jet_direction(jet.momentum().Unit().x(), jet.momentum().Unit().y(), jet.momentum().Unit().z());
0716   GlobalVector jet_ref_track_dir(jet.px(), jet.py(), jet.pz());
0717   math::XYZVector jet_dir = jet.momentum().Unit();
0718 
0719   std::vector<pat::PackedCandidate> jet_lost_tracks;
0720   for (size_t itrk = 0; itrk < losttracks_->size(); itrk++) {
0721     if (reco::deltaR(losttracks_->at(itrk).p4(), jet.p4()) < max_dr_for_losttrack_ and
0722         losttracks_->at(itrk).pt() > min_pt_for_losttrack_) {
0723       jet_lost_tracks.push_back(losttracks_->at(itrk));
0724     }
0725   }
0726   std::sort(
0727       jet_lost_tracks.begin(), jet_lost_tracks.end(), [](const auto &a, const auto &b) { return a.pt() > b.pt(); });
0728 
0729   // reserve space
0730   for (const auto &name : losttrack_features_)
0731     fts.reserve(name, jet_lost_tracks.size());
0732 
0733   reco::VertexRef pv_ass = reco::VertexRef(vtxs_, 0);
0734   math::XYZPoint pv_ass_pos = pv_ass->position();
0735 
0736   for (auto const &ltrack : jet_lost_tracks) {
0737     fts.fill("jet_losttrack_pt_log", std::isnan(std::log(ltrack.pt())) ? 0 : std::log(ltrack.pt()));
0738     fts.fill("jet_losttrack_eta", ltrack.eta());
0739     fts.fill("jet_losttrack_charge", ltrack.charge());
0740     fts.fill("jet_losttrack_frompv", ltrack.fromPV());
0741     fts.fill("jet_losttrack_dz", std::isnan(ltrack.dz(pv_ass_pos)) ? 0 : ltrack.dz(pv_ass_pos));
0742     fts.fill("jet_losttrack_dxy", std::isnan(ltrack.dxy(pv_ass_pos)) ? 0 : ltrack.dxy(pv_ass_pos));
0743     fts.fill("jet_losttrack_npixhits", ltrack.numberOfPixelHits());
0744     fts.fill("jet_losttrack_nstriphits", ltrack.stripLayersWithMeasurement());
0745 
0746     TVector3 ltrack_momentum(ltrack.momentum().x(), ltrack.momentum().y(), ltrack.momentum().z());
0747     fts.fill("jet_losttrack_deta", jet_direction.Eta() - ltrack_momentum.Eta());
0748     fts.fill("jet_losttrack_dphi", jet_direction.DeltaPhi(ltrack_momentum));
0749     fts.fill("jet_losttrack_etarel",
0750              std::isnan(reco::btau::etaRel(jet_dir, ltrack.momentum()))
0751                  ? 0
0752                  : reco::btau::etaRel(jet_dir, ltrack.momentum()));
0753 
0754     const reco::Track *track = ltrack.bestTrack();
0755     if (track) {
0756       fts.fill("jet_losttrack_track_chi2", track->normalizedChi2());
0757       fts.fill("jet_losttrack_track_qual", track->qualityMask());
0758       fts.fill("jet_losttrack_dxysig",
0759                std::isnan(fabs(ltrack.dxy(pv_ass_pos)) / ltrack.dxyError())
0760                    ? 0
0761                    : fabs(ltrack.dxy(pv_ass_pos)) / ltrack.dxyError());
0762       fts.fill("jet_losttrack_dzsig",
0763                std::isnan(fabs(ltrack.dz(pv_ass_pos)) / ltrack.dzError())
0764                    ? 0
0765                    : fabs(ltrack.dz(pv_ass_pos)) / ltrack.dzError());
0766 
0767       reco::TransientTrack transientTrack = track_builder_->build(*track);
0768       Measurement1D meas_ip3d = IPTools::signedImpactParameter3D(transientTrack, jet_ref_track_dir, *pv_).second;
0769       Measurement1D meas_jetdist = IPTools::jetTrackDistance(transientTrack, jet_ref_track_dir, *pv_).second;
0770       Measurement1D meas_decayl = IPTools::signedDecayLength3D(transientTrack, jet_ref_track_dir, *pv_).second;
0771 
0772       fts.fill("jet_losttrack_trackjet_d3d", std::isnan(meas_ip3d.value()) ? 0 : meas_ip3d.value());
0773       fts.fill("jet_losttrack_trackjet_d3dsig",
0774                std::isnan(fabs(meas_ip3d.significance())) ? 0 : fabs(meas_ip3d.significance()));
0775       fts.fill("jet_losttrack_trackjet_dist", std::isnan(-meas_jetdist.value()) ? 0 : -meas_jetdist.value());
0776       fts.fill("jet_losttrack_trackjet_decayL", std::isnan(meas_decayl.value()) ? 0 : meas_decayl.value());
0777     } else {
0778       fts.fill("jet_losttrack_track_chi2", 0);
0779       fts.fill("jet_losttrack_track_qual", 0);
0780       fts.fill("jet_losttrack_dxysig", 0);
0781       fts.fill("jet_losttrack_dzsig", 0);
0782       fts.fill("jet_losttrack_trackjet_d3d", 0);
0783       fts.fill("jet_losttrack_trackjet_d3dsig", 0);
0784       fts.fill("jet_losttrack_trackjet_dist", 0);
0785       fts.fill("jet_losttrack_trackjet_decayL", 0);
0786     }
0787 
0788     fts.fill("lt_mask", 1);
0789   }
0790 }
0791 
0792 // define this as a plug-in
0793 DEFINE_FWK_MODULE(ParticleNetFeatureEvaluator);