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
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
0309 auto output_tag_infos = std::make_unique<std::vector<reco::DeepBoostedJetTagInfo>>();
0310
0311 auto jets = iEvent.getHandle(jet_token_);
0312
0313 auto muons = iEvent.getHandle(muon_token_);
0314
0315 auto taus = iEvent.getHandle(tau_token_);
0316
0317 auto electrons = iEvent.getHandle(electron_token_);
0318
0319 auto photons = iEvent.getHandle(photon_token_);
0320
0321 iEvent.getByToken(losttrack_token_, losttracks_);
0322
0323 iEvent.getByToken(vtx_token_, vtxs_);
0324 if (vtxs_->empty()) {
0325
0326 iEvent.put(std::move(output_tag_infos));
0327 return;
0328 }
0329
0330 pv_ = &vtxs_->at(0);
0331
0332 iEvent.getByToken(sv_token_, svs_);
0333
0334 iEvent.getByToken(pfcand_token_, pfcands_);
0335
0336 track_builder_ = iSetup.getHandle(track_builder_token_);
0337
0338
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
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
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
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
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
0384 output_tag_infos->emplace_back(features, jet_ref);
0385 }
0386
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
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
0406 reco::VertexRefProd PVRefProd(vtxs_);
0407
0408 TrackInfoBuilder trackinfo(track_builder_);
0409
0410
0411 std::vector<const pat::PackedCandidate *> daughters;
0412 for (const auto &dau : jet.daughterPtrVector()) {
0413
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
0418 if (cand->pt() < min_pt_for_pfcandidates_)
0419 continue;
0420
0421 if (!include_neutrals_ and (cand->charge() == 0 or cand->pt() < min_pt_for_track_properties_))
0422 continue;
0423
0424 daughters.push_back(cand);
0425 }
0426
0427
0428 std::sort(daughters.begin(), daughters.end(), [](const auto &a, const auto &b) { return a->pt() > b->pt(); });
0429
0430
0431 for (const auto &name : particle_features_)
0432 fts.reserve(name, daughters.size());
0433
0434
0435 for (const auto &cand : daughters) {
0436 if (!include_neutrals_ and !useTrackProperties(cand))
0437 continue;
0438
0439
0440 auto candP4 = cand->p4();
0441 auto candP3 = cand->momentum();
0442
0443
0444 const reco::Track *track = nullptr;
0445 if (useTrackProperties(cand))
0446 track = cand->bestTrack();
0447
0448
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
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
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
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
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
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
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
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
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
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 <rack : 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
0793 DEFINE_FWK_MODULE(ParticleNetFeatureEvaluator);