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