File indexing completed on 2023-02-04 03:48:10
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #include "RecoTauTag/RecoTau/interface/DeepTauBase.h"
0011 #include "RecoTauTag/RecoTau/interface/DeepTauScaling.h"
0012 #include "FWCore/Utilities/interface/isFinite.h"
0013 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0014 #include "DataFormats/TauReco/interface/PFTauTransverseImpactParameterAssociation.h"
0015
0016 #include <fstream>
0017 #include "oneapi/tbb/concurrent_unordered_set.h"
0018
0019 namespace deep_tau {
0020 constexpr int NumberOfOutputs = 4;
0021 }
0022
0023 namespace {
0024
0025 namespace dnn_inputs_v2 {
0026 constexpr int number_of_inner_cell = 11;
0027 constexpr int number_of_outer_cell = 21;
0028 constexpr int number_of_conv_features = 64;
0029 namespace TauBlockInputs {
0030 enum vars {
0031 rho = 0,
0032 tau_pt,
0033 tau_eta,
0034 tau_phi,
0035 tau_mass,
0036 tau_E_over_pt,
0037 tau_charge,
0038 tau_n_charged_prongs,
0039 tau_n_neutral_prongs,
0040 chargedIsoPtSum,
0041 chargedIsoPtSumdR03_over_dR05,
0042 footprintCorrection,
0043 neutralIsoPtSum,
0044 neutralIsoPtSumWeight_over_neutralIsoPtSum,
0045 neutralIsoPtSumWeightdR03_over_neutralIsoPtSum,
0046 neutralIsoPtSumdR03_over_dR05,
0047 photonPtSumOutsideSignalCone,
0048 puCorrPtSum,
0049 tau_dxy_pca_x,
0050 tau_dxy_pca_y,
0051 tau_dxy_pca_z,
0052 tau_dxy_valid,
0053 tau_dxy,
0054 tau_dxy_sig,
0055 tau_ip3d_valid,
0056 tau_ip3d,
0057 tau_ip3d_sig,
0058 tau_dz,
0059 tau_dz_sig_valid,
0060 tau_dz_sig,
0061 tau_flightLength_x,
0062 tau_flightLength_y,
0063 tau_flightLength_z,
0064 tau_flightLength_sig,
0065 tau_pt_weighted_deta_strip,
0066 tau_pt_weighted_dphi_strip,
0067 tau_pt_weighted_dr_signal,
0068 tau_pt_weighted_dr_iso,
0069 tau_leadingTrackNormChi2,
0070 tau_e_ratio_valid,
0071 tau_e_ratio,
0072 tau_gj_angle_diff_valid,
0073 tau_gj_angle_diff,
0074 tau_n_photons,
0075 tau_emFraction,
0076 tau_inside_ecal_crack,
0077 leadChargedCand_etaAtEcalEntrance_minus_tau_eta,
0078 NumberOfInputs
0079 };
0080 std::vector<int> varsToDrop = {
0081 tau_phi, tau_dxy_pca_x, tau_dxy_pca_y, tau_dxy_pca_z};
0082 }
0083
0084 namespace EgammaBlockInputs {
0085 enum vars {
0086 rho = 0,
0087 tau_pt,
0088 tau_eta,
0089 tau_inside_ecal_crack,
0090 pfCand_ele_valid,
0091 pfCand_ele_rel_pt,
0092 pfCand_ele_deta,
0093 pfCand_ele_dphi,
0094 pfCand_ele_pvAssociationQuality,
0095 pfCand_ele_puppiWeight,
0096 pfCand_ele_charge,
0097 pfCand_ele_lostInnerHits,
0098 pfCand_ele_numberOfPixelHits,
0099 pfCand_ele_vertex_dx,
0100 pfCand_ele_vertex_dy,
0101 pfCand_ele_vertex_dz,
0102 pfCand_ele_vertex_dx_tauFL,
0103 pfCand_ele_vertex_dy_tauFL,
0104 pfCand_ele_vertex_dz_tauFL,
0105 pfCand_ele_hasTrackDetails,
0106 pfCand_ele_dxy,
0107 pfCand_ele_dxy_sig,
0108 pfCand_ele_dz,
0109 pfCand_ele_dz_sig,
0110 pfCand_ele_track_chi2_ndof,
0111 pfCand_ele_track_ndof,
0112 ele_valid,
0113 ele_rel_pt,
0114 ele_deta,
0115 ele_dphi,
0116 ele_cc_valid,
0117 ele_cc_ele_rel_energy,
0118 ele_cc_gamma_rel_energy,
0119 ele_cc_n_gamma,
0120 ele_rel_trackMomentumAtVtx,
0121 ele_rel_trackMomentumAtCalo,
0122 ele_rel_trackMomentumOut,
0123 ele_rel_trackMomentumAtEleClus,
0124 ele_rel_trackMomentumAtVtxWithConstraint,
0125 ele_rel_ecalEnergy,
0126 ele_ecalEnergy_sig,
0127 ele_eSuperClusterOverP,
0128 ele_eSeedClusterOverP,
0129 ele_eSeedClusterOverPout,
0130 ele_eEleClusterOverPout,
0131 ele_deltaEtaSuperClusterTrackAtVtx,
0132 ele_deltaEtaSeedClusterTrackAtCalo,
0133 ele_deltaEtaEleClusterTrackAtCalo,
0134 ele_deltaPhiEleClusterTrackAtCalo,
0135 ele_deltaPhiSuperClusterTrackAtVtx,
0136 ele_deltaPhiSeedClusterTrackAtCalo,
0137 ele_mvaInput_earlyBrem,
0138 ele_mvaInput_lateBrem,
0139 ele_mvaInput_sigmaEtaEta,
0140 ele_mvaInput_hadEnergy,
0141 ele_mvaInput_deltaEta,
0142 ele_gsfTrack_normalizedChi2,
0143 ele_gsfTrack_numberOfValidHits,
0144 ele_rel_gsfTrack_pt,
0145 ele_gsfTrack_pt_sig,
0146 ele_has_closestCtfTrack,
0147 ele_closestCtfTrack_normalizedChi2,
0148 ele_closestCtfTrack_numberOfValidHits,
0149 pfCand_gamma_valid,
0150 pfCand_gamma_rel_pt,
0151 pfCand_gamma_deta,
0152 pfCand_gamma_dphi,
0153 pfCand_gamma_pvAssociationQuality,
0154 pfCand_gamma_fromPV,
0155 pfCand_gamma_puppiWeight,
0156 pfCand_gamma_puppiWeightNoLep,
0157 pfCand_gamma_lostInnerHits,
0158 pfCand_gamma_numberOfPixelHits,
0159 pfCand_gamma_vertex_dx,
0160 pfCand_gamma_vertex_dy,
0161 pfCand_gamma_vertex_dz,
0162 pfCand_gamma_vertex_dx_tauFL,
0163 pfCand_gamma_vertex_dy_tauFL,
0164 pfCand_gamma_vertex_dz_tauFL,
0165 pfCand_gamma_hasTrackDetails,
0166 pfCand_gamma_dxy,
0167 pfCand_gamma_dxy_sig,
0168 pfCand_gamma_dz,
0169 pfCand_gamma_dz_sig,
0170 pfCand_gamma_track_chi2_ndof,
0171 pfCand_gamma_track_ndof,
0172 NumberOfInputs
0173 };
0174 }
0175
0176 namespace MuonBlockInputs {
0177 enum vars {
0178 rho = 0,
0179 tau_pt,
0180 tau_eta,
0181 tau_inside_ecal_crack,
0182 pfCand_muon_valid,
0183 pfCand_muon_rel_pt,
0184 pfCand_muon_deta,
0185 pfCand_muon_dphi,
0186 pfCand_muon_pvAssociationQuality,
0187 pfCand_muon_fromPV,
0188 pfCand_muon_puppiWeight,
0189 pfCand_muon_charge,
0190 pfCand_muon_lostInnerHits,
0191 pfCand_muon_numberOfPixelHits,
0192 pfCand_muon_vertex_dx,
0193 pfCand_muon_vertex_dy,
0194 pfCand_muon_vertex_dz,
0195 pfCand_muon_vertex_dx_tauFL,
0196 pfCand_muon_vertex_dy_tauFL,
0197 pfCand_muon_vertex_dz_tauFL,
0198 pfCand_muon_hasTrackDetails,
0199 pfCand_muon_dxy,
0200 pfCand_muon_dxy_sig,
0201 pfCand_muon_dz,
0202 pfCand_muon_dz_sig,
0203 pfCand_muon_track_chi2_ndof,
0204 pfCand_muon_track_ndof,
0205 muon_valid,
0206 muon_rel_pt,
0207 muon_deta,
0208 muon_dphi,
0209 muon_dxy,
0210 muon_dxy_sig,
0211 muon_normalizedChi2_valid,
0212 muon_normalizedChi2,
0213 muon_numberOfValidHits,
0214 muon_segmentCompatibility,
0215 muon_caloCompatibility,
0216 muon_pfEcalEnergy_valid,
0217 muon_rel_pfEcalEnergy,
0218 muon_n_matches_DT_1,
0219 muon_n_matches_DT_2,
0220 muon_n_matches_DT_3,
0221 muon_n_matches_DT_4,
0222 muon_n_matches_CSC_1,
0223 muon_n_matches_CSC_2,
0224 muon_n_matches_CSC_3,
0225 muon_n_matches_CSC_4,
0226 muon_n_matches_RPC_1,
0227 muon_n_matches_RPC_2,
0228 muon_n_matches_RPC_3,
0229 muon_n_matches_RPC_4,
0230 muon_n_hits_DT_1,
0231 muon_n_hits_DT_2,
0232 muon_n_hits_DT_3,
0233 muon_n_hits_DT_4,
0234 muon_n_hits_CSC_1,
0235 muon_n_hits_CSC_2,
0236 muon_n_hits_CSC_3,
0237 muon_n_hits_CSC_4,
0238 muon_n_hits_RPC_1,
0239 muon_n_hits_RPC_2,
0240 muon_n_hits_RPC_3,
0241 muon_n_hits_RPC_4,
0242 NumberOfInputs
0243 };
0244 }
0245
0246 namespace HadronBlockInputs {
0247 enum vars {
0248 rho = 0,
0249 tau_pt,
0250 tau_eta,
0251 tau_inside_ecal_crack,
0252 pfCand_chHad_valid,
0253 pfCand_chHad_rel_pt,
0254 pfCand_chHad_deta,
0255 pfCand_chHad_dphi,
0256 pfCand_chHad_leadChargedHadrCand,
0257 pfCand_chHad_pvAssociationQuality,
0258 pfCand_chHad_fromPV,
0259 pfCand_chHad_puppiWeight,
0260 pfCand_chHad_puppiWeightNoLep,
0261 pfCand_chHad_charge,
0262 pfCand_chHad_lostInnerHits,
0263 pfCand_chHad_numberOfPixelHits,
0264 pfCand_chHad_vertex_dx,
0265 pfCand_chHad_vertex_dy,
0266 pfCand_chHad_vertex_dz,
0267 pfCand_chHad_vertex_dx_tauFL,
0268 pfCand_chHad_vertex_dy_tauFL,
0269 pfCand_chHad_vertex_dz_tauFL,
0270 pfCand_chHad_hasTrackDetails,
0271 pfCand_chHad_dxy,
0272 pfCand_chHad_dxy_sig,
0273 pfCand_chHad_dz,
0274 pfCand_chHad_dz_sig,
0275 pfCand_chHad_track_chi2_ndof,
0276 pfCand_chHad_track_ndof,
0277 pfCand_chHad_hcalFraction,
0278 pfCand_chHad_rawCaloFraction,
0279 pfCand_nHad_valid,
0280 pfCand_nHad_rel_pt,
0281 pfCand_nHad_deta,
0282 pfCand_nHad_dphi,
0283 pfCand_nHad_puppiWeight,
0284 pfCand_nHad_puppiWeightNoLep,
0285 pfCand_nHad_hcalFraction,
0286 NumberOfInputs
0287 };
0288 }
0289 }
0290
0291 float getTauID(const pat::Tau& tau, const std::string& tauID, float default_value = -999., bool assert_input = true) {
0292 static tbb::concurrent_unordered_set<std::string> isFirstWarning;
0293 if (tau.isTauIDAvailable(tauID)) {
0294 return tau.tauID(tauID);
0295 } else {
0296 if (assert_input) {
0297 throw cms::Exception("DeepTauId")
0298 << "Exception in <getTauID>: No tauID '" << tauID << "' available in pat::Tau given as function argument.";
0299 }
0300 if (isFirstWarning.insert(tauID).second) {
0301 edm::LogWarning("DeepTauID") << "Warning in <getTauID>: No tauID '" << tauID
0302 << "' available in pat::Tau given as function argument."
0303 << " Using default_value = " << default_value << " instead." << std::endl;
0304 }
0305 return default_value;
0306 }
0307 }
0308
0309 struct TauFunc {
0310 const reco::TauDiscriminatorContainer* basicTauDiscriminatorCollection;
0311 const reco::TauDiscriminatorContainer* basicTauDiscriminatordR03Collection;
0312 const edm::AssociationVector<reco::PFTauRefProd, std::vector<reco::PFTauTransverseImpactParameterRef>>*
0313 pfTauTransverseImpactParameters;
0314
0315 using BasicDiscr = deep_tau::DeepTauBase::BasicDiscriminator;
0316 std::map<BasicDiscr, size_t> indexMap;
0317 std::map<BasicDiscr, size_t> indexMapdR03;
0318
0319 const float getChargedIsoPtSum(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0320 return (*basicTauDiscriminatorCollection)[tau_ref].rawValues.at(indexMap.at(BasicDiscr::ChargedIsoPtSum));
0321 }
0322 const float getChargedIsoPtSum(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0323 return getTauID(tau, "chargedIsoPtSum");
0324 }
0325 const float getChargedIsoPtSumdR03(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0326 return (*basicTauDiscriminatordR03Collection)[tau_ref].rawValues.at(indexMapdR03.at(BasicDiscr::ChargedIsoPtSum));
0327 }
0328 const float getChargedIsoPtSumdR03(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0329 return getTauID(tau, "chargedIsoPtSumdR03");
0330 }
0331 const float getFootprintCorrectiondR03(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0332 return (*basicTauDiscriminatordR03Collection)[tau_ref].rawValues.at(
0333 indexMapdR03.at(BasicDiscr::FootprintCorrection));
0334 }
0335 const float getFootprintCorrectiondR03(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0336 return getTauID(tau, "footprintCorrectiondR03");
0337 }
0338 const float getFootprintCorrection(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0339 return getTauID(tau, "footprintCorrection");
0340 }
0341 const float getNeutralIsoPtSum(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0342 return (*basicTauDiscriminatorCollection)[tau_ref].rawValues.at(indexMap.at(BasicDiscr::NeutralIsoPtSum));
0343 }
0344 const float getNeutralIsoPtSum(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0345 return getTauID(tau, "neutralIsoPtSum");
0346 }
0347 const float getNeutralIsoPtSumdR03(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0348 return (*basicTauDiscriminatordR03Collection)[tau_ref].rawValues.at(indexMapdR03.at(BasicDiscr::NeutralIsoPtSum));
0349 }
0350 const float getNeutralIsoPtSumdR03(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0351 return getTauID(tau, "neutralIsoPtSumdR03");
0352 }
0353 const float getNeutralIsoPtSumWeight(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0354 return (*basicTauDiscriminatorCollection)[tau_ref].rawValues.at(indexMap.at(BasicDiscr::NeutralIsoPtSumWeight));
0355 }
0356 const float getNeutralIsoPtSumWeight(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0357 return getTauID(tau, "neutralIsoPtSumWeight");
0358 }
0359 const float getNeutralIsoPtSumdR03Weight(const reco::PFTau& tau,
0360 const edm::RefToBase<reco::BaseTau> tau_ref) const {
0361 return (*basicTauDiscriminatordR03Collection)[tau_ref].rawValues.at(
0362 indexMapdR03.at(BasicDiscr::NeutralIsoPtSumWeight));
0363 }
0364 const float getNeutralIsoPtSumdR03Weight(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0365 return getTauID(tau, "neutralIsoPtSumWeightdR03");
0366 }
0367 const float getPhotonPtSumOutsideSignalCone(const reco::PFTau& tau,
0368 const edm::RefToBase<reco::BaseTau> tau_ref) const {
0369 return (*basicTauDiscriminatorCollection)[tau_ref].rawValues.at(
0370 indexMap.at(BasicDiscr::PhotonPtSumOutsideSignalCone));
0371 }
0372 const float getPhotonPtSumOutsideSignalCone(const pat::Tau& tau,
0373 const edm::RefToBase<reco::BaseTau> tau_ref) const {
0374 return getTauID(tau, "photonPtSumOutsideSignalCone");
0375 }
0376 const float getPhotonPtSumOutsideSignalConedR03(const reco::PFTau& tau,
0377 const edm::RefToBase<reco::BaseTau> tau_ref) const {
0378 return (*basicTauDiscriminatordR03Collection)[tau_ref].rawValues.at(
0379 indexMapdR03.at(BasicDiscr::PhotonPtSumOutsideSignalCone));
0380 }
0381 const float getPhotonPtSumOutsideSignalConedR03(const pat::Tau& tau,
0382 const edm::RefToBase<reco::BaseTau> tau_ref) const {
0383 return getTauID(tau, "photonPtSumOutsideSignalConedR03");
0384 }
0385 const float getPuCorrPtSum(const reco::PFTau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0386 return (*basicTauDiscriminatorCollection)[tau_ref].rawValues.at(indexMap.at(BasicDiscr::PUcorrPtSum));
0387 }
0388 const float getPuCorrPtSum(const pat::Tau& tau, const edm::RefToBase<reco::BaseTau> tau_ref) const {
0389 return getTauID(tau, "puCorrPtSum");
0390 }
0391
0392 auto getdxyPCA(const reco::PFTau& tau, const size_t tau_index) const {
0393 return pfTauTransverseImpactParameters->value(tau_index)->dxy_PCA();
0394 }
0395 auto getdxyPCA(const pat::Tau& tau, const size_t tau_index) const { return tau.dxy_PCA(); }
0396 auto getdxy(const reco::PFTau& tau, const size_t tau_index) const {
0397 return pfTauTransverseImpactParameters->value(tau_index)->dxy();
0398 }
0399 auto getdxy(const pat::Tau& tau, const size_t tau_index) const { return tau.dxy(); }
0400 auto getdxyError(const reco::PFTau& tau, const size_t tau_index) const {
0401 return pfTauTransverseImpactParameters->value(tau_index)->dxy_error();
0402 }
0403 auto getdxyError(const pat::Tau& tau, const size_t tau_index) const { return tau.dxy_error(); }
0404 auto getdxySig(const reco::PFTau& tau, const size_t tau_index) const {
0405 return pfTauTransverseImpactParameters->value(tau_index)->dxy_Sig();
0406 }
0407 auto getdxySig(const pat::Tau& tau, const size_t tau_index) const { return tau.dxy_Sig(); }
0408 auto getip3d(const reco::PFTau& tau, const size_t tau_index) const {
0409 return pfTauTransverseImpactParameters->value(tau_index)->ip3d();
0410 }
0411 auto getip3d(const pat::Tau& tau, const size_t tau_index) const { return tau.ip3d(); }
0412 auto getip3dError(const reco::PFTau& tau, const size_t tau_index) const {
0413 return pfTauTransverseImpactParameters->value(tau_index)->ip3d_error();
0414 }
0415 auto getip3dError(const pat::Tau& tau, const size_t tau_index) const { return tau.ip3d_error(); }
0416 auto getip3dSig(const reco::PFTau& tau, const size_t tau_index) const {
0417 return pfTauTransverseImpactParameters->value(tau_index)->ip3d_Sig();
0418 }
0419 auto getip3dSig(const pat::Tau& tau, const size_t tau_index) const { return tau.ip3d_Sig(); }
0420 auto getHasSecondaryVertex(const reco::PFTau& tau, const size_t tau_index) const {
0421 return pfTauTransverseImpactParameters->value(tau_index)->hasSecondaryVertex();
0422 }
0423 auto getHasSecondaryVertex(const pat::Tau& tau, const size_t tau_index) const { return tau.hasSecondaryVertex(); }
0424 auto getFlightLength(const reco::PFTau& tau, const size_t tau_index) const {
0425 return pfTauTransverseImpactParameters->value(tau_index)->flightLength();
0426 }
0427 auto getFlightLength(const pat::Tau& tau, const size_t tau_index) const { return tau.flightLength(); }
0428 auto getFlightLengthSig(const reco::PFTau& tau, const size_t tau_index) const {
0429 return pfTauTransverseImpactParameters->value(tau_index)->flightLengthSig();
0430 }
0431 auto getFlightLengthSig(const pat::Tau& tau, const size_t tau_index) const { return tau.flightLengthSig(); }
0432
0433 auto getLeadingTrackNormChi2(const reco::PFTau& tau) { return reco::tau::lead_track_chi2(tau); }
0434 auto getLeadingTrackNormChi2(const pat::Tau& tau) { return tau.leadingTrackNormChi2(); }
0435 auto getEmFraction(const pat::Tau& tau) { return tau.emFraction_MVA(); }
0436 auto getEmFraction(const reco::PFTau& tau) { return tau.emFraction(); }
0437 auto getEtaAtEcalEntrance(const pat::Tau& tau) { return tau.etaAtEcalEntranceLeadChargedCand(); }
0438 auto getEtaAtEcalEntrance(const reco::PFTau& tau) {
0439 return tau.leadPFChargedHadrCand()->positionAtECALEntrance().eta();
0440 }
0441 auto getEcalEnergyLeadingChargedHadr(const reco::PFTau& tau) { return tau.leadPFChargedHadrCand()->ecalEnergy(); }
0442 auto getEcalEnergyLeadingChargedHadr(const pat::Tau& tau) { return tau.ecalEnergyLeadChargedHadrCand(); }
0443 auto getHcalEnergyLeadingChargedHadr(const reco::PFTau& tau) { return tau.leadPFChargedHadrCand()->hcalEnergy(); }
0444 auto getHcalEnergyLeadingChargedHadr(const pat::Tau& tau) { return tau.hcalEnergyLeadChargedHadrCand(); }
0445
0446 template <typename PreDiscrType>
0447 bool passPrediscriminants(const PreDiscrType prediscriminants,
0448 const size_t andPrediscriminants,
0449 const edm::RefToBase<reco::BaseTau> tau_ref) {
0450 bool passesPrediscriminants = (andPrediscriminants ? 1 : 0);
0451
0452 size_t nPrediscriminants = prediscriminants.size();
0453 for (size_t iDisc = 0; iDisc < nPrediscriminants; ++iDisc) {
0454
0455 double discResult = (*prediscriminants[iDisc].handle)[tau_ref];
0456 uint8_t thisPasses = (discResult > prediscriminants[iDisc].cut) ? 1 : 0;
0457
0458
0459
0460
0461
0462
0463
0464
0465
0466
0467
0468
0469
0470
0471
0472
0473 if (thisPasses ^ andPrediscriminants)
0474 {
0475 passesPrediscriminants = (andPrediscriminants ? 0 : 1);
0476 break;
0477 }
0478 }
0479 return passesPrediscriminants;
0480 }
0481 };
0482
0483 namespace candFunc {
0484 auto getTauDz(const reco::PFCandidate& cand) { return cand.bestTrack()->dz(); }
0485 auto getTauDz(const pat::PackedCandidate& cand) { return cand.dz(); }
0486 auto getTauDZSigValid(const reco::PFCandidate& cand) {
0487 return cand.bestTrack() != nullptr && std::isnormal(cand.bestTrack()->dz()) && std::isnormal(cand.dzError()) &&
0488 cand.dzError() > 0;
0489 }
0490 auto getTauDZSigValid(const pat::PackedCandidate& cand) {
0491 return cand.hasTrackDetails() && std::isnormal(cand.dz()) && std::isnormal(cand.dzError()) && cand.dzError() > 0;
0492 }
0493 auto getTauDxy(const reco::PFCandidate& cand) { return cand.bestTrack()->dxy(); }
0494 auto getTauDxy(const pat::PackedCandidate& cand) { return cand.dxy(); }
0495 auto getPvAssocationQuality(const reco::PFCandidate& cand) { return 0.7013f; }
0496 auto getPvAssocationQuality(const pat::PackedCandidate& cand) { return cand.pvAssociationQuality(); }
0497 auto getPuppiWeight(const reco::PFCandidate& cand, const float aod_value) { return aod_value; }
0498 auto getPuppiWeight(const pat::PackedCandidate& cand, const float aod_value) { return cand.puppiWeight(); }
0499 auto getPuppiWeightNoLep(const reco::PFCandidate& cand, const float aod_value) { return aod_value; }
0500 auto getPuppiWeightNoLep(const pat::PackedCandidate& cand, const float aod_value) {
0501 return cand.puppiWeightNoLep();
0502 }
0503 auto getLostInnerHits(const reco::PFCandidate& cand, float default_value) {
0504 return cand.bestTrack() != nullptr
0505 ? cand.bestTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
0506 : default_value;
0507 }
0508 auto getLostInnerHits(const pat::PackedCandidate& cand, float default_value) { return cand.lostInnerHits(); }
0509 auto getNumberOfPixelHits(const reco::PFCandidate& cand, float default_value) {
0510 return cand.bestTrack() != nullptr
0511 ? cand.bestTrack()->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
0512 : default_value;
0513 }
0514 auto getNumberOfPixelHits(const pat::PackedCandidate& cand, float default_value) {
0515 return cand.numberOfPixelHits();
0516 }
0517 auto getHasTrackDetails(const reco::PFCandidate& cand) { return cand.bestTrack() != nullptr; }
0518 auto getHasTrackDetails(const pat::PackedCandidate& cand) { return cand.hasTrackDetails(); }
0519 auto getPseudoTrack(const reco::PFCandidate& cand) { return *cand.bestTrack(); }
0520 auto getPseudoTrack(const pat::PackedCandidate& cand) { return cand.pseudoTrack(); }
0521 auto getFromPV(const reco::PFCandidate& cand) { return 0.9994f; }
0522 auto getFromPV(const pat::PackedCandidate& cand) { return cand.fromPV(); }
0523 auto getHCalFraction(const reco::PFCandidate& cand, bool disable_hcalFraction_workaround) {
0524 return cand.rawHcalEnergy() / (cand.rawHcalEnergy() + cand.rawEcalEnergy());
0525 }
0526 auto getHCalFraction(const pat::PackedCandidate& cand, bool disable_hcalFraction_workaround) {
0527 float hcal_fraction = 0.;
0528 if (disable_hcalFraction_workaround) {
0529
0530
0531 hcal_fraction = cand.hcalFraction();
0532 } else {
0533
0534 if (cand.pdgId() == 1 || cand.pdgId() == 130) {
0535 hcal_fraction = cand.hcalFraction();
0536 } else if (cand.isIsolatedChargedHadron()) {
0537 hcal_fraction = cand.rawHcalFraction();
0538 }
0539 }
0540 return hcal_fraction;
0541 }
0542 auto getRawCaloFraction(const reco::PFCandidate& cand) {
0543 return (cand.rawEcalEnergy() + cand.rawHcalEnergy()) / cand.energy();
0544 }
0545 auto getRawCaloFraction(const pat::PackedCandidate& cand) { return cand.rawCaloFraction(); }
0546 };
0547
0548 template <typename LVector1, typename LVector2>
0549 float dEta(const LVector1& p4, const LVector2& tau_p4) {
0550 return static_cast<float>(p4.eta() - tau_p4.eta());
0551 }
0552
0553 template <typename LVector1, typename LVector2>
0554 float dPhi(const LVector1& p4_1, const LVector2& p4_2) {
0555 return static_cast<float>(reco::deltaPhi(p4_2.phi(), p4_1.phi()));
0556 }
0557
0558 struct MuonHitMatchV2 {
0559 static constexpr size_t n_muon_stations = 4;
0560 static constexpr int first_station_id = 1;
0561 static constexpr int last_station_id = first_station_id + n_muon_stations - 1;
0562 using CountArray = std::array<unsigned, n_muon_stations>;
0563 using CountMap = std::map<int, CountArray>;
0564
0565 const std::vector<int>& consideredSubdets() {
0566 static const std::vector<int> subdets = {MuonSubdetId::DT, MuonSubdetId::CSC, MuonSubdetId::RPC};
0567 return subdets;
0568 }
0569
0570 const std::string& subdetName(int subdet) {
0571 static const std::map<int, std::string> subdet_names = {
0572 {MuonSubdetId::DT, "DT"}, {MuonSubdetId::CSC, "CSC"}, {MuonSubdetId::RPC, "RPC"}};
0573 if (!subdet_names.count(subdet))
0574 throw cms::Exception("MuonHitMatch") << "Subdet name for subdet id " << subdet << " not found.";
0575 return subdet_names.at(subdet);
0576 }
0577
0578 size_t getStationIndex(int station, bool throw_exception) const {
0579 if (station < first_station_id || station > last_station_id) {
0580 if (throw_exception)
0581 throw cms::Exception("MuonHitMatch") << "Station id is out of range";
0582 return std::numeric_limits<size_t>::max();
0583 }
0584 return static_cast<size_t>(station - 1);
0585 }
0586
0587 MuonHitMatchV2(const pat::Muon& muon) {
0588 for (int subdet : consideredSubdets()) {
0589 n_matches[subdet].fill(0);
0590 n_hits[subdet].fill(0);
0591 }
0592
0593 countMatches(muon, n_matches);
0594 countHits(muon, n_hits);
0595 }
0596
0597 void countMatches(const pat::Muon& muon, CountMap& n_matches) {
0598 for (const auto& segment : muon.matches()) {
0599 if (segment.segmentMatches.empty() && segment.rpcMatches.empty())
0600 continue;
0601 if (n_matches.count(segment.detector())) {
0602 const size_t station_index = getStationIndex(segment.station(), true);
0603 ++n_matches.at(segment.detector()).at(station_index);
0604 }
0605 }
0606 }
0607
0608 void countHits(const pat::Muon& muon, CountMap& n_hits) {
0609 if (muon.outerTrack().isNonnull()) {
0610 const auto& hit_pattern = muon.outerTrack()->hitPattern();
0611 for (int hit_index = 0; hit_index < hit_pattern.numberOfAllHits(reco::HitPattern::TRACK_HITS); ++hit_index) {
0612 auto hit_id = hit_pattern.getHitPattern(reco::HitPattern::TRACK_HITS, hit_index);
0613 if (hit_id == 0)
0614 break;
0615 if (hit_pattern.muonHitFilter(hit_id) && (hit_pattern.getHitType(hit_id) == TrackingRecHit::valid ||
0616 hit_pattern.getHitType(hit_id) == TrackingRecHit::bad)) {
0617 const size_t station_index = getStationIndex(hit_pattern.getMuonStation(hit_id), false);
0618 if (station_index < n_muon_stations) {
0619 CountArray* muon_n_hits = nullptr;
0620 if (hit_pattern.muonDTHitFilter(hit_id))
0621 muon_n_hits = &n_hits.at(MuonSubdetId::DT);
0622 else if (hit_pattern.muonCSCHitFilter(hit_id))
0623 muon_n_hits = &n_hits.at(MuonSubdetId::CSC);
0624 else if (hit_pattern.muonRPCHitFilter(hit_id))
0625 muon_n_hits = &n_hits.at(MuonSubdetId::RPC);
0626
0627 if (muon_n_hits)
0628 ++muon_n_hits->at(station_index);
0629 }
0630 }
0631 }
0632 }
0633 }
0634
0635 unsigned nMatches(int subdet, int station) const {
0636 if (!n_matches.count(subdet))
0637 throw cms::Exception("MuonHitMatch") << "Subdet " << subdet << " not found.";
0638 const size_t station_index = getStationIndex(station, true);
0639 return n_matches.at(subdet).at(station_index);
0640 }
0641
0642 unsigned nHits(int subdet, int station) const {
0643 if (!n_hits.count(subdet))
0644 throw cms::Exception("MuonHitMatch") << "Subdet " << subdet << " not found.";
0645 const size_t station_index = getStationIndex(station, true);
0646 return n_hits.at(subdet).at(station_index);
0647 }
0648
0649 unsigned countMuonStationsWithMatches(int first_station, int last_station) const {
0650 static const std::map<int, std::vector<bool>> masks = {
0651 {MuonSubdetId::DT, {false, false, false, false}},
0652 {MuonSubdetId::CSC, {true, false, false, false}},
0653 {MuonSubdetId::RPC, {false, false, false, false}},
0654 };
0655 const size_t first_station_index = getStationIndex(first_station, true);
0656 const size_t last_station_index = getStationIndex(last_station, true);
0657 unsigned cnt = 0;
0658 for (size_t n = first_station_index; n <= last_station_index; ++n) {
0659 for (const auto& match : n_matches) {
0660 if (!masks.at(match.first).at(n) && match.second.at(n) > 0)
0661 ++cnt;
0662 }
0663 }
0664 return cnt;
0665 }
0666
0667 unsigned countMuonStationsWithHits(int first_station, int last_station) const {
0668 static const std::map<int, std::vector<bool>> masks = {
0669 {MuonSubdetId::DT, {false, false, false, false}},
0670 {MuonSubdetId::CSC, {false, false, false, false}},
0671 {MuonSubdetId::RPC, {false, false, false, false}},
0672 };
0673
0674 const size_t first_station_index = getStationIndex(first_station, true);
0675 const size_t last_station_index = getStationIndex(last_station, true);
0676 unsigned cnt = 0;
0677 for (size_t n = first_station_index; n <= last_station_index; ++n) {
0678 for (const auto& hit : n_hits) {
0679 if (!masks.at(hit.first).at(n) && hit.second.at(n) > 0)
0680 ++cnt;
0681 }
0682 }
0683 return cnt;
0684 }
0685
0686 private:
0687 CountMap n_matches, n_hits;
0688 };
0689
0690 enum class CellObjectType {
0691 PfCand_electron,
0692 PfCand_muon,
0693 PfCand_chargedHadron,
0694 PfCand_neutralHadron,
0695 PfCand_gamma,
0696 Electron,
0697 Muon,
0698 Other
0699 };
0700
0701 template <typename Object>
0702 CellObjectType GetCellObjectType(const Object&);
0703 template <>
0704 CellObjectType GetCellObjectType(const pat::Electron&) {
0705 return CellObjectType::Electron;
0706 }
0707 template <>
0708 CellObjectType GetCellObjectType(const pat::Muon&) {
0709 return CellObjectType::Muon;
0710 }
0711
0712 template <>
0713 CellObjectType GetCellObjectType(reco::Candidate const& cand) {
0714 static const std::map<int, CellObjectType> obj_types = {{11, CellObjectType::PfCand_electron},
0715 {13, CellObjectType::PfCand_muon},
0716 {22, CellObjectType::PfCand_gamma},
0717 {130, CellObjectType::PfCand_neutralHadron},
0718 {211, CellObjectType::PfCand_chargedHadron}};
0719
0720 auto iter = obj_types.find(std::abs(cand.pdgId()));
0721 if (iter == obj_types.end())
0722 return CellObjectType::Other;
0723 return iter->second;
0724 }
0725
0726 using Cell = std::map<CellObjectType, size_t>;
0727 struct CellIndex {
0728 int eta, phi;
0729
0730 bool operator<(const CellIndex& other) const {
0731 if (eta != other.eta)
0732 return eta < other.eta;
0733 return phi < other.phi;
0734 }
0735 };
0736
0737 class CellGrid {
0738 public:
0739 using Map = std::map<CellIndex, Cell>;
0740 using const_iterator = Map::const_iterator;
0741
0742 CellGrid(unsigned n_cells_eta,
0743 unsigned n_cells_phi,
0744 double cell_size_eta,
0745 double cell_size_phi,
0746 bool disable_CellIndex_workaround)
0747 : nCellsEta(n_cells_eta),
0748 nCellsPhi(n_cells_phi),
0749 nTotal(nCellsEta * nCellsPhi),
0750 cellSizeEta(cell_size_eta),
0751 cellSizePhi(cell_size_phi),
0752 disable_CellIndex_workaround_(disable_CellIndex_workaround) {
0753 if (nCellsEta % 2 != 1 || nCellsEta < 1)
0754 throw cms::Exception("DeepTauId") << "Invalid number of eta cells.";
0755 if (nCellsPhi % 2 != 1 || nCellsPhi < 1)
0756 throw cms::Exception("DeepTauId") << "Invalid number of phi cells.";
0757 if (cellSizeEta <= 0 || cellSizePhi <= 0)
0758 throw cms::Exception("DeepTauId") << "Invalid cell size.";
0759 }
0760
0761 int maxEtaIndex() const { return static_cast<int>((nCellsEta - 1) / 2); }
0762 int maxPhiIndex() const { return static_cast<int>((nCellsPhi - 1) / 2); }
0763 double maxDeltaEta() const { return cellSizeEta * (0.5 + maxEtaIndex()); }
0764 double maxDeltaPhi() const { return cellSizePhi * (0.5 + maxPhiIndex()); }
0765 int getEtaTensorIndex(const CellIndex& cellIndex) const { return cellIndex.eta + maxEtaIndex(); }
0766 int getPhiTensorIndex(const CellIndex& cellIndex) const { return cellIndex.phi + maxPhiIndex(); }
0767
0768 bool tryGetCellIndex(double deltaEta, double deltaPhi, CellIndex& cellIndex) const {
0769 const auto getCellIndex = [this](double x, double maxX, double size, int& index) {
0770 const double absX = std::abs(x);
0771 if (absX > maxX)
0772 return false;
0773 double absIndex;
0774 if (disable_CellIndex_workaround_) {
0775
0776
0777 absIndex = std::floor(absX / size + 0.5);
0778 } else {
0779
0780 absIndex = std::floor(std::abs(absX / size - 0.5));
0781 }
0782 index = static_cast<int>(std::copysign(absIndex, x));
0783 return true;
0784 };
0785
0786 return getCellIndex(deltaEta, maxDeltaEta(), cellSizeEta, cellIndex.eta) &&
0787 getCellIndex(deltaPhi, maxDeltaPhi(), cellSizePhi, cellIndex.phi);
0788 }
0789
0790 size_t num_valid_cells() const { return cells.size(); }
0791 Cell& operator[](const CellIndex& cellIndex) { return cells[cellIndex]; }
0792 const Cell& at(const CellIndex& cellIndex) const { return cells.at(cellIndex); }
0793 size_t count(const CellIndex& cellIndex) const { return cells.count(cellIndex); }
0794 const_iterator find(const CellIndex& cellIndex) const { return cells.find(cellIndex); }
0795 const_iterator begin() const { return cells.begin(); }
0796 const_iterator end() const { return cells.end(); }
0797
0798 public:
0799 const unsigned nCellsEta, nCellsPhi, nTotal;
0800 const double cellSizeEta, cellSizePhi;
0801
0802 private:
0803 std::map<CellIndex, Cell> cells;
0804 const bool disable_CellIndex_workaround_;
0805 };
0806 }
0807
0808 using bd = deep_tau::DeepTauBase::BasicDiscriminator;
0809 const std::map<bd, std::string> deep_tau::DeepTauBase::stringFromDiscriminator_{
0810 {bd::ChargedIsoPtSum, "ChargedIsoPtSum"},
0811 {bd::NeutralIsoPtSum, "NeutralIsoPtSum"},
0812 {bd::NeutralIsoPtSumWeight, "NeutralIsoPtSumWeight"},
0813 {bd::FootprintCorrection, "TauFootprintCorrection"},
0814 {bd::PhotonPtSumOutsideSignalCone, "PhotonPtSumOutsideSignalCone"},
0815 {bd::PUcorrPtSum, "PUcorrPtSum"}};
0816 const std::vector<bd> deep_tau::DeepTauBase::requiredBasicDiscriminators_ = {bd::ChargedIsoPtSum,
0817 bd::NeutralIsoPtSum,
0818 bd::NeutralIsoPtSumWeight,
0819 bd::PhotonPtSumOutsideSignalCone,
0820 bd::PUcorrPtSum};
0821 const std::vector<bd> deep_tau::DeepTauBase::requiredBasicDiscriminatorsdR03_ = {bd::ChargedIsoPtSum,
0822 bd::NeutralIsoPtSum,
0823 bd::NeutralIsoPtSumWeight,
0824 bd::PhotonPtSumOutsideSignalCone,
0825 bd::FootprintCorrection};
0826
0827 class DeepTauId : public deep_tau::DeepTauBase {
0828 public:
0829 static constexpr float default_value = -999.;
0830
0831 static const OutputCollection& GetOutputs() {
0832 static constexpr size_t e_index = 0, mu_index = 1, tau_index = 2, jet_index = 3;
0833 static const OutputCollection outputs_ = {
0834 {"VSe", Output({tau_index}, {e_index, tau_index})},
0835 {"VSmu", Output({tau_index}, {mu_index, tau_index})},
0836 {"VSjet", Output({tau_index}, {jet_index, tau_index})},
0837 };
0838 return outputs_;
0839 }
0840
0841 const std::map<BasicDiscriminator, size_t> matchDiscriminatorIndices(
0842 edm::Event& event,
0843 edm::EDGetTokenT<reco::TauDiscriminatorContainer> discriminatorContainerToken,
0844 std::vector<BasicDiscriminator> requiredDiscr) {
0845 std::map<std::string, size_t> discrIndexMapStr;
0846 auto const aHandle = event.getHandle(discriminatorContainerToken);
0847 auto const aProv = aHandle.provenance();
0848 if (aProv == nullptr)
0849 aHandle.whyFailed()->raise();
0850 const auto& psetsFromProvenance = edm::parameterSet(aProv->stable(), event.processHistory());
0851 auto const idlist = psetsFromProvenance.getParameter<std::vector<edm::ParameterSet>>("IDdefinitions");
0852 for (size_t j = 0; j < idlist.size(); ++j) {
0853 std::string idname = idlist[j].getParameter<std::string>("IDname");
0854 if (discrIndexMapStr.count(idname)) {
0855 throw cms::Exception("DeepTauId")
0856 << "basic discriminator " << idname << " appears more than once in the input.";
0857 }
0858 discrIndexMapStr[idname] = j;
0859 }
0860
0861
0862 std::map<BasicDiscriminator, size_t> discrIndexMap;
0863 for (size_t i = 0; i < requiredDiscr.size(); i++) {
0864 if (discrIndexMapStr.find(stringFromDiscriminator_.at(requiredDiscr[i])) == discrIndexMapStr.end())
0865 throw cms::Exception("DeepTauId") << "Basic Discriminator " << stringFromDiscriminator_.at(requiredDiscr[i])
0866 << " was not provided in the config file.";
0867 else
0868 discrIndexMap[requiredDiscr[i]] = discrIndexMapStr[stringFromDiscriminator_.at(requiredDiscr[i])];
0869 }
0870 return discrIndexMap;
0871 }
0872
0873 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0874 edm::ParameterSetDescription desc;
0875 desc.add<edm::InputTag>("electrons", edm::InputTag("slimmedElectrons"));
0876 desc.add<edm::InputTag>("muons", edm::InputTag("slimmedMuons"));
0877 desc.add<edm::InputTag>("taus", edm::InputTag("slimmedTaus"));
0878 desc.add<edm::InputTag>("pfcands", edm::InputTag("packedPFCandidates"));
0879 desc.add<edm::InputTag>("vertices", edm::InputTag("offlineSlimmedPrimaryVertices"));
0880 desc.add<edm::InputTag>("rho", edm::InputTag("fixedGridRhoAll"));
0881 desc.add<std::vector<std::string>>("graph_file",
0882 {"RecoTauTag/TrainingFiles/data/DeepTauId/deepTau_2017v2p6_e6.pb"});
0883 desc.add<bool>("mem_mapped", false);
0884 desc.add<unsigned>("year", 2017);
0885 desc.add<unsigned>("version", 2);
0886 desc.add<unsigned>("sub_version", 1);
0887 desc.add<int>("debug_level", 0);
0888 desc.add<bool>("disable_dxy_pca", false);
0889 desc.add<bool>("disable_hcalFraction_workaround", false);
0890 desc.add<bool>("disable_CellIndex_workaround", false);
0891 desc.add<bool>("save_inputs", false);
0892 desc.add<bool>("is_online", false);
0893
0894 desc.add<std::vector<std::string>>("VSeWP", {"-1."});
0895 desc.add<std::vector<std::string>>("VSmuWP", {"-1."});
0896 desc.add<std::vector<std::string>>("VSjetWP", {"-1."});
0897
0898 desc.addUntracked<edm::InputTag>("basicTauDiscriminators", edm::InputTag("basicTauDiscriminators"));
0899 desc.addUntracked<edm::InputTag>("basicTauDiscriminatorsdR03", edm::InputTag("basicTauDiscriminatorsdR03"));
0900 desc.add<edm::InputTag>("pfTauTransverseImpactParameters", edm::InputTag("hpsPFTauTransverseImpactParameters"));
0901
0902 {
0903 edm::ParameterSetDescription pset_Prediscriminants;
0904 pset_Prediscriminants.add<std::string>("BooleanOperator", "and");
0905 {
0906 edm::ParameterSetDescription psd1;
0907 psd1.add<double>("cut");
0908 psd1.add<edm::InputTag>("Producer");
0909 pset_Prediscriminants.addOptional<edm::ParameterSetDescription>("decayMode", psd1);
0910 }
0911 desc.add<edm::ParameterSetDescription>("Prediscriminants", pset_Prediscriminants);
0912 }
0913
0914 descriptions.add("DeepTau", desc);
0915 }
0916
0917 public:
0918 explicit DeepTauId(const edm::ParameterSet& cfg, const deep_tau::DeepTauCache* cache)
0919 : DeepTauBase(cfg, GetOutputs(), cache),
0920 electrons_token_(consumes<std::vector<pat::Electron>>(cfg.getParameter<edm::InputTag>("electrons"))),
0921 muons_token_(consumes<std::vector<pat::Muon>>(cfg.getParameter<edm::InputTag>("muons"))),
0922 rho_token_(consumes<double>(cfg.getParameter<edm::InputTag>("rho"))),
0923 basicTauDiscriminators_inputToken_(consumes<reco::TauDiscriminatorContainer>(
0924 cfg.getUntrackedParameter<edm::InputTag>("basicTauDiscriminators"))),
0925 basicTauDiscriminatorsdR03_inputToken_(consumes<reco::TauDiscriminatorContainer>(
0926 cfg.getUntrackedParameter<edm::InputTag>("basicTauDiscriminatorsdR03"))),
0927 pfTauTransverseImpactParameters_token_(
0928 consumes<edm::AssociationVector<reco::PFTauRefProd, std::vector<reco::PFTauTransverseImpactParameterRef>>>(
0929 cfg.getParameter<edm::InputTag>("pfTauTransverseImpactParameters"))),
0930 year_(cfg.getParameter<unsigned>("year")),
0931 version_(cfg.getParameter<unsigned>("version")),
0932 sub_version_(cfg.getParameter<unsigned>("sub_version")),
0933 debug_level(cfg.getParameter<int>("debug_level")),
0934 disable_dxy_pca_(cfg.getParameter<bool>("disable_dxy_pca")),
0935 disable_hcalFraction_workaround_(cfg.getParameter<bool>("disable_hcalFraction_workaround")),
0936 disable_CellIndex_workaround_(cfg.getParameter<bool>("disable_CellIndex_workaround")),
0937 save_inputs_(cfg.getParameter<bool>("save_inputs")),
0938 json_file_(nullptr),
0939 file_counter_(0) {
0940 if (version_ == 2) {
0941 using namespace dnn_inputs_v2;
0942 namespace sc = deep_tau::Scaling;
0943 tauInputs_indices_.resize(TauBlockInputs::NumberOfInputs);
0944 std::iota(std::begin(tauInputs_indices_), std::end(tauInputs_indices_), 0);
0945
0946 if (sub_version_ == 1) {
0947 tauBlockTensor_ = std::make_unique<tensorflow::Tensor>(
0948 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, TauBlockInputs::NumberOfInputs});
0949 scalingParamsMap_ = &sc::scalingParamsMap_v2p1;
0950 } else if (sub_version_ == 5) {
0951 std::sort(TauBlockInputs::varsToDrop.begin(), TauBlockInputs::varsToDrop.end());
0952 for (auto v : TauBlockInputs::varsToDrop) {
0953 tauInputs_indices_.at(v) = -1;
0954 for (std::size_t i = v + 1; i < TauBlockInputs::NumberOfInputs; ++i)
0955 tauInputs_indices_.at(i) -= 1;
0956 }
0957 tauBlockTensor_ = std::make_unique<tensorflow::Tensor>(
0958 tensorflow::DT_FLOAT,
0959 tensorflow::TensorShape{1,
0960 static_cast<int>(TauBlockInputs::NumberOfInputs) -
0961 static_cast<int>(TauBlockInputs::varsToDrop.size())});
0962 if (year_ == 2026) {
0963 scalingParamsMap_ = &sc::scalingParamsMap_PhaseIIv2p5;
0964 } else {
0965 scalingParamsMap_ = &sc::scalingParamsMap_v2p5;
0966 }
0967 } else
0968 throw cms::Exception("DeepTauId") << "subversion " << sub_version_ << " is not supported.";
0969
0970 std::map<std::vector<bool>, std::vector<sc::FeatureT>> GridFeatureTypes_map = {
0971 {{false}, {sc::FeatureT::TauFlat, sc::FeatureT::GridGlobal}},
0972 {{false, true},
0973 {sc::FeatureT::PfCand_electron,
0974 sc::FeatureT::PfCand_muon,
0975 sc::FeatureT::PfCand_chHad,
0976 sc::FeatureT::PfCand_nHad,
0977 sc::FeatureT::PfCand_gamma,
0978 sc::FeatureT::Electron,
0979 sc::FeatureT::Muon}}};
0980
0981
0982 for (const auto& p : GridFeatureTypes_map) {
0983 for (auto is_inner : p.first) {
0984 for (auto featureType : p.second) {
0985 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(featureType, is_inner));
0986 if (!(sp.mean_.size() == sp.std_.size() && sp.mean_.size() == sp.lim_min_.size() &&
0987 sp.mean_.size() == sp.lim_max_.size()))
0988 throw cms::Exception("DeepTauId") << "sizes of scaling parameter vectors do not match between each other";
0989 }
0990 }
0991 }
0992
0993 for (size_t n = 0; n < 2; ++n) {
0994 const bool is_inner = n == 0;
0995 const auto n_cells = is_inner ? number_of_inner_cell : number_of_outer_cell;
0996 eGammaTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
0997 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, 1, 1, EgammaBlockInputs::NumberOfInputs});
0998 muonTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
0999 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, 1, 1, MuonBlockInputs::NumberOfInputs});
1000 hadronsTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1001 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, 1, 1, HadronBlockInputs::NumberOfInputs});
1002 convTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1003 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, n_cells, n_cells, number_of_conv_features});
1004 zeroOutputTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1005 tensorflow::DT_FLOAT, tensorflow::TensorShape{1, 1, 1, number_of_conv_features});
1006
1007 eGammaTensor_[is_inner]->flat<float>().setZero();
1008 muonTensor_[is_inner]->flat<float>().setZero();
1009 hadronsTensor_[is_inner]->flat<float>().setZero();
1010
1011 setCellConvFeatures(*zeroOutputTensor_[is_inner], getPartialPredictions(is_inner), 0, 0, 0);
1012 }
1013 } else {
1014 throw cms::Exception("DeepTauId") << "version " << version_ << " is not supported.";
1015 }
1016 }
1017
1018 static std::unique_ptr<deep_tau::DeepTauCache> initializeGlobalCache(const edm::ParameterSet& cfg) {
1019 return DeepTauBase::initializeGlobalCache(cfg);
1020 }
1021
1022 static void globalEndJob(const deep_tau::DeepTauCache* cache_) { return DeepTauBase::globalEndJob(cache_); }
1023
1024 private:
1025 static constexpr float pi = M_PI;
1026
1027 template <typename T>
1028 static float getValue(T value) {
1029 return std::isnormal(value) ? static_cast<float>(value) : 0.f;
1030 }
1031
1032 template <typename T>
1033 static float getValueLinear(T value, float min_value, float max_value, bool positive) {
1034 const float fixed_value = getValue(value);
1035 const float clamped_value = std::clamp(fixed_value, min_value, max_value);
1036 float transformed_value = (clamped_value - min_value) / (max_value - min_value);
1037 if (!positive)
1038 transformed_value = transformed_value * 2 - 1;
1039 return transformed_value;
1040 }
1041
1042 template <typename T>
1043 static float getValueNorm(T value, float mean, float sigma, float n_sigmas_max = 5) {
1044 const float fixed_value = getValue(value);
1045 const float norm_value = (fixed_value - mean) / sigma;
1046 return std::clamp(norm_value, -n_sigmas_max, n_sigmas_max);
1047 }
1048
1049 static bool isAbove(double value, double min) { return std::isnormal(value) && value > min; }
1050
1051 static bool calculateElectronClusterVarsV2(const pat::Electron& ele,
1052 float& cc_ele_energy,
1053 float& cc_gamma_energy,
1054 int& cc_n_gamma) {
1055 cc_ele_energy = cc_gamma_energy = 0;
1056 cc_n_gamma = 0;
1057 const auto& superCluster = ele.superCluster();
1058 if (superCluster.isNonnull() && superCluster.isAvailable() && superCluster->clusters().isNonnull() &&
1059 superCluster->clusters().isAvailable()) {
1060 for (auto iter = superCluster->clustersBegin(); iter != superCluster->clustersEnd(); ++iter) {
1061 const float energy = static_cast<float>((*iter)->energy());
1062 if (iter == superCluster->clustersBegin())
1063 cc_ele_energy += energy;
1064 else {
1065 cc_gamma_energy += energy;
1066 ++cc_n_gamma;
1067 }
1068 }
1069 return true;
1070 } else
1071 return false;
1072 }
1073
1074 inline void checkInputs(const tensorflow::Tensor& inputs,
1075 const std::string& block_name,
1076 int n_inputs,
1077 const CellGrid* grid = nullptr) const {
1078 if (debug_level >= 1) {
1079 std::cout << "<checkInputs>: block_name = " << block_name << std::endl;
1080 if (block_name == "input_tau") {
1081 for (int input_index = 0; input_index < n_inputs; ++input_index) {
1082 float input = inputs.matrix<float>()(0, input_index);
1083 if (edm::isNotFinite(input)) {
1084 throw cms::Exception("DeepTauId")
1085 << "in the " << block_name
1086 << ", input is not finite, i.e. infinite or NaN, for input_index = " << input_index;
1087 }
1088 if (debug_level >= 2) {
1089 std::cout << block_name << "[var = " << input_index << "] = " << std::setprecision(5) << std::fixed << input
1090 << std::endl;
1091 }
1092 }
1093 } else {
1094 assert(grid);
1095 int n_eta, n_phi;
1096 if (block_name.find("input_inner") != std::string::npos) {
1097 n_eta = 5;
1098 n_phi = 5;
1099 } else if (block_name.find("input_outer") != std::string::npos) {
1100 n_eta = 10;
1101 n_phi = 10;
1102 } else
1103 assert(0);
1104 int eta_phi_index = 0;
1105 for (int eta = -n_eta; eta <= n_eta; ++eta) {
1106 for (int phi = -n_phi; phi <= n_phi; ++phi) {
1107 const CellIndex cell_index{eta, phi};
1108 const auto cell_iter = grid->find(cell_index);
1109 if (cell_iter != grid->end()) {
1110 for (int input_index = 0; input_index < n_inputs; ++input_index) {
1111 float input = inputs.tensor<float, 4>()(eta_phi_index, 0, 0, input_index);
1112 if (edm::isNotFinite(input)) {
1113 throw cms::Exception("DeepTauId")
1114 << "in the " << block_name << ", input is not finite, i.e. infinite or NaN, for eta = " << eta
1115 << ", phi = " << phi << ", input_index = " << input_index;
1116 }
1117 if (debug_level >= 2) {
1118 std::cout << block_name << "[eta = " << eta << "][phi = " << phi << "][var = " << input_index
1119 << "] = " << std::setprecision(5) << std::fixed << input << std::endl;
1120 }
1121 }
1122 eta_phi_index += 1;
1123 }
1124 }
1125 }
1126 }
1127 }
1128 }
1129
1130 inline void saveInputs(const tensorflow::Tensor& inputs,
1131 const std::string& block_name,
1132 int n_inputs,
1133 const CellGrid* grid = nullptr) {
1134 if (debug_level >= 1) {
1135 std::cout << "<saveInputs>: block_name = " << block_name << std::endl;
1136 }
1137 if (!is_first_block_)
1138 (*json_file_) << ", ";
1139 (*json_file_) << "\"" << block_name << "\": [";
1140 if (block_name == "input_tau") {
1141 for (int input_index = 0; input_index < n_inputs; ++input_index) {
1142 float input = inputs.matrix<float>()(0, input_index);
1143 if (input_index != 0)
1144 (*json_file_) << ", ";
1145 (*json_file_) << input;
1146 }
1147 } else {
1148 assert(grid);
1149 int n_eta, n_phi;
1150 if (block_name.find("input_inner") != std::string::npos) {
1151 n_eta = 5;
1152 n_phi = 5;
1153 } else if (block_name.find("input_outer") != std::string::npos) {
1154 n_eta = 10;
1155 n_phi = 10;
1156 } else
1157 assert(0);
1158 int eta_phi_index = 0;
1159 for (int eta = -n_eta; eta <= n_eta; ++eta) {
1160 if (eta != -n_eta)
1161 (*json_file_) << ", ";
1162 (*json_file_) << "[";
1163 for (int phi = -n_phi; phi <= n_phi; ++phi) {
1164 if (phi != -n_phi)
1165 (*json_file_) << ", ";
1166 (*json_file_) << "[";
1167 const CellIndex cell_index{eta, phi};
1168 const auto cell_iter = grid->find(cell_index);
1169 for (int input_index = 0; input_index < n_inputs; ++input_index) {
1170 float input = 0.;
1171 if (cell_iter != grid->end()) {
1172 input = inputs.tensor<float, 4>()(eta_phi_index, 0, 0, input_index);
1173 }
1174 if (input_index != 0)
1175 (*json_file_) << ", ";
1176 (*json_file_) << input;
1177 }
1178 if (cell_iter != grid->end()) {
1179 eta_phi_index += 1;
1180 }
1181 (*json_file_) << "]";
1182 }
1183 (*json_file_) << "]";
1184 }
1185 }
1186 (*json_file_) << "]";
1187 is_first_block_ = false;
1188 }
1189
1190 private:
1191 tensorflow::Tensor getPredictions(edm::Event& event, edm::Handle<TauCollection> taus) override {
1192
1193 const std::vector<pat::Electron> electron_collection_default;
1194 const std::vector<pat::Muon> muon_collection_default;
1195 const reco::TauDiscriminatorContainer basicTauDiscriminators_default;
1196 const reco::TauDiscriminatorContainer basicTauDiscriminatorsdR03_default;
1197 const edm::AssociationVector<reco::PFTauRefProd, std::vector<reco::PFTauTransverseImpactParameterRef>>
1198 pfTauTransverseImpactParameters_default;
1199
1200 const std::vector<pat::Electron>* electron_collection;
1201 const std::vector<pat::Muon>* muon_collection;
1202 const reco::TauDiscriminatorContainer* basicTauDiscriminators;
1203 const reco::TauDiscriminatorContainer* basicTauDiscriminatorsdR03;
1204 const edm::AssociationVector<reco::PFTauRefProd, std::vector<reco::PFTauTransverseImpactParameterRef>>*
1205 pfTauTransverseImpactParameters;
1206
1207 if (!is_online_) {
1208 electron_collection = &event.get(electrons_token_);
1209 muon_collection = &event.get(muons_token_);
1210 pfTauTransverseImpactParameters = &pfTauTransverseImpactParameters_default;
1211 basicTauDiscriminators = &basicTauDiscriminators_default;
1212 basicTauDiscriminatorsdR03 = &basicTauDiscriminatorsdR03_default;
1213 } else {
1214 electron_collection = &electron_collection_default;
1215 muon_collection = &muon_collection_default;
1216 pfTauTransverseImpactParameters = &event.get(pfTauTransverseImpactParameters_token_);
1217 basicTauDiscriminators = &event.get(basicTauDiscriminators_inputToken_);
1218 basicTauDiscriminatorsdR03 = &event.get(basicTauDiscriminatorsdR03_inputToken_);
1219
1220
1221 if (!discrIndicesMapped_) {
1222 basicDiscrIndexMap_ =
1223 matchDiscriminatorIndices(event, basicTauDiscriminators_inputToken_, requiredBasicDiscriminators_);
1224 basicDiscrdR03IndexMap_ =
1225 matchDiscriminatorIndices(event, basicTauDiscriminatorsdR03_inputToken_, requiredBasicDiscriminatorsdR03_);
1226 discrIndicesMapped_ = true;
1227 }
1228 }
1229
1230 TauFunc tauIDs = {basicTauDiscriminators,
1231 basicTauDiscriminatorsdR03,
1232 pfTauTransverseImpactParameters,
1233 basicDiscrIndexMap_,
1234 basicDiscrdR03IndexMap_};
1235
1236 edm::Handle<edm::View<reco::Candidate>> pfCands;
1237 event.getByToken(pfcandToken_, pfCands);
1238
1239 edm::Handle<reco::VertexCollection> vertices;
1240 event.getByToken(vtxToken_, vertices);
1241
1242 edm::Handle<double> rho;
1243 event.getByToken(rho_token_, rho);
1244
1245 auto const& eventnr = event.id().event();
1246
1247 tensorflow::Tensor predictions(tensorflow::DT_FLOAT, {static_cast<int>(taus->size()), deep_tau::NumberOfOutputs});
1248
1249 for (size_t tau_index = 0; tau_index < taus->size(); ++tau_index) {
1250 const edm::RefToBase<reco::BaseTau> tauRef = taus->refAt(tau_index);
1251
1252 std::vector<tensorflow::Tensor> pred_vector;
1253
1254 bool passesPrediscriminants;
1255 if (is_online_) {
1256 passesPrediscriminants = tauIDs.passPrediscriminants<std::vector<TauDiscInfo<reco::PFTauDiscriminator>>>(
1257 recoPrediscriminants_, andPrediscriminants_, tauRef);
1258 } else {
1259 passesPrediscriminants = tauIDs.passPrediscriminants<std::vector<TauDiscInfo<pat::PATTauDiscriminator>>>(
1260 patPrediscriminants_, andPrediscriminants_, tauRef);
1261 }
1262
1263 if (passesPrediscriminants) {
1264 if (version_ == 2) {
1265 if (is_online_) {
1266 getPredictionsV2<reco::PFCandidate, reco::PFTau>(taus->at(tau_index),
1267 tau_index,
1268 tauRef,
1269 electron_collection,
1270 muon_collection,
1271 *pfCands,
1272 vertices->at(0),
1273 *rho,
1274 eventnr,
1275 pred_vector,
1276 tauIDs);
1277 } else
1278 getPredictionsV2<pat::PackedCandidate, pat::Tau>(taus->at(tau_index),
1279 tau_index,
1280 tauRef,
1281 electron_collection,
1282 muon_collection,
1283 *pfCands,
1284 vertices->at(0),
1285 *rho,
1286 eventnr,
1287 pred_vector,
1288 tauIDs);
1289 } else {
1290 throw cms::Exception("DeepTauId") << "version " << version_ << " is not supported.";
1291 }
1292
1293 for (int k = 0; k < deep_tau::NumberOfOutputs; ++k) {
1294 const float pred = pred_vector[0].flat<float>()(k);
1295 if (!(pred >= 0 && pred <= 1))
1296 throw cms::Exception("DeepTauId")
1297 << "invalid prediction = " << pred << " for tau_index = " << tau_index << ", pred_index = " << k;
1298 predictions.matrix<float>()(tau_index, k) = pred;
1299 }
1300 } else {
1301
1302
1303
1304
1305
1306 for (int k = 0; k < deep_tau::NumberOfOutputs; ++k) {
1307 predictions.matrix<float>()(tau_index, k) = (k == 2) ? -1.f : 2.f;
1308 }
1309 }
1310 }
1311 return predictions;
1312 }
1313
1314 template <typename CandidateCastType, typename TauCastType>
1315 void getPredictionsV2(TauCollection::const_reference& tau,
1316 const size_t tau_index,
1317 const edm::RefToBase<reco::BaseTau> tau_ref,
1318 const std::vector<pat::Electron>* electrons,
1319 const std::vector<pat::Muon>* muons,
1320 const edm::View<reco::Candidate>& pfCands,
1321 const reco::Vertex& pv,
1322 double rho,
1323 const edm::EventNumber_t& eventnr,
1324 std::vector<tensorflow::Tensor>& pred_vector,
1325 TauFunc tau_funcs) {
1326 using namespace dnn_inputs_v2;
1327 if (debug_level >= 2) {
1328 std::cout << "<DeepTauId::getPredictionsV2 (moduleLabel = " << moduleDescription().moduleLabel()
1329 << ")>:" << std::endl;
1330 std::cout << " tau: pT = " << tau.pt() << ", eta = " << tau.eta() << ", phi = " << tau.phi()
1331 << ", eventnr = " << eventnr << std::endl;
1332 }
1333 CellGrid inner_grid(number_of_inner_cell, number_of_inner_cell, 0.02, 0.02, disable_CellIndex_workaround_);
1334 CellGrid outer_grid(number_of_outer_cell, number_of_outer_cell, 0.05, 0.05, disable_CellIndex_workaround_);
1335 fillGrids(dynamic_cast<const TauCastType&>(tau), *electrons, inner_grid, outer_grid);
1336 fillGrids(dynamic_cast<const TauCastType&>(tau), *muons, inner_grid, outer_grid);
1337 fillGrids(dynamic_cast<const TauCastType&>(tau), pfCands, inner_grid, outer_grid);
1338
1339 createTauBlockInputs<CandidateCastType>(
1340 dynamic_cast<const TauCastType&>(tau), tau_index, tau_ref, pv, rho, tau_funcs);
1341 checkInputs(*tauBlockTensor_, "input_tau", static_cast<int>(tauBlockTensor_->shape().dim_size(1)));
1342 createConvFeatures<CandidateCastType>(dynamic_cast<const TauCastType&>(tau),
1343 tau_index,
1344 tau_ref,
1345 pv,
1346 rho,
1347 electrons,
1348 muons,
1349 pfCands,
1350 inner_grid,
1351 tau_funcs,
1352 true);
1353 checkInputs(*eGammaTensor_[true], "input_inner_egamma", EgammaBlockInputs::NumberOfInputs, &inner_grid);
1354 checkInputs(*muonTensor_[true], "input_inner_muon", MuonBlockInputs::NumberOfInputs, &inner_grid);
1355 checkInputs(*hadronsTensor_[true], "input_inner_hadrons", HadronBlockInputs::NumberOfInputs, &inner_grid);
1356 createConvFeatures<CandidateCastType>(dynamic_cast<const TauCastType&>(tau),
1357 tau_index,
1358 tau_ref,
1359 pv,
1360 rho,
1361 electrons,
1362 muons,
1363 pfCands,
1364 outer_grid,
1365 tau_funcs,
1366 false);
1367 checkInputs(*eGammaTensor_[false], "input_outer_egamma", EgammaBlockInputs::NumberOfInputs, &outer_grid);
1368 checkInputs(*muonTensor_[false], "input_outer_muon", MuonBlockInputs::NumberOfInputs, &outer_grid);
1369 checkInputs(*hadronsTensor_[false], "input_outer_hadrons", HadronBlockInputs::NumberOfInputs, &outer_grid);
1370
1371 if (save_inputs_) {
1372 std::string json_file_name = "DeepTauId_" + std::to_string(eventnr) + "_" + std::to_string(tau_index) + ".json";
1373 json_file_ = new std::ofstream(json_file_name.data());
1374 is_first_block_ = true;
1375 (*json_file_) << "{";
1376 saveInputs(*tauBlockTensor_, "input_tau", static_cast<int>(tauBlockTensor_->shape().dim_size(1)));
1377 saveInputs(
1378 *eGammaTensor_[true], "input_inner_egamma", dnn_inputs_v2::EgammaBlockInputs::NumberOfInputs, &inner_grid);
1379 saveInputs(*muonTensor_[true], "input_inner_muon", dnn_inputs_v2::MuonBlockInputs::NumberOfInputs, &inner_grid);
1380 saveInputs(
1381 *hadronsTensor_[true], "input_inner_hadrons", dnn_inputs_v2::HadronBlockInputs::NumberOfInputs, &inner_grid);
1382 saveInputs(
1383 *eGammaTensor_[false], "input_outer_egamma", dnn_inputs_v2::EgammaBlockInputs::NumberOfInputs, &outer_grid);
1384 saveInputs(*muonTensor_[false], "input_outer_muon", dnn_inputs_v2::MuonBlockInputs::NumberOfInputs, &outer_grid);
1385 saveInputs(
1386 *hadronsTensor_[false], "input_outer_hadrons", dnn_inputs_v2::HadronBlockInputs::NumberOfInputs, &outer_grid);
1387 (*json_file_) << "}";
1388 delete json_file_;
1389 ++file_counter_;
1390 }
1391
1392 tensorflow::run(&(cache_->getSession("core")),
1393 {{"input_tau", *tauBlockTensor_},
1394 {"input_inner", *convTensor_.at(true)},
1395 {"input_outer", *convTensor_.at(false)}},
1396 {"main_output/Softmax"},
1397 &pred_vector);
1398 if (debug_level >= 1) {
1399 std::cout << "output = { ";
1400 for (int idx = 0; idx < deep_tau::NumberOfOutputs; ++idx) {
1401 if (idx > 0)
1402 std::cout << ", ";
1403 std::string label;
1404 if (idx == 0)
1405 label = "e";
1406 else if (idx == 1)
1407 label = "mu";
1408 else if (idx == 2)
1409 label = "tau";
1410 else if (idx == 3)
1411 label = "jet";
1412 else
1413 assert(0);
1414 std::cout << label << " = " << pred_vector[0].flat<float>()(idx);
1415 }
1416 std::cout << " }" << std::endl;
1417 }
1418 }
1419
1420 template <typename Collection, typename TauCastType>
1421 void fillGrids(const TauCastType& tau, const Collection& objects, CellGrid& inner_grid, CellGrid& outer_grid) {
1422 static constexpr double outer_dR2 = 0.25;
1423 const double inner_radius = getInnerSignalConeRadius(tau.polarP4().pt());
1424 const double inner_dR2 = std::pow(inner_radius, 2);
1425
1426 const auto addObject = [&](size_t n, double deta, double dphi, CellGrid& grid) {
1427 const auto& obj = objects.at(n);
1428 const CellObjectType obj_type = GetCellObjectType(obj);
1429 if (obj_type == CellObjectType::Other)
1430 return;
1431 CellIndex cell_index;
1432 if (grid.tryGetCellIndex(deta, dphi, cell_index)) {
1433 Cell& cell = grid[cell_index];
1434 auto iter = cell.find(obj_type);
1435 if (iter != cell.end()) {
1436 const auto& prev_obj = objects.at(iter->second);
1437 if (obj.polarP4().pt() > prev_obj.polarP4().pt())
1438 iter->second = n;
1439 } else {
1440 cell[obj_type] = n;
1441 }
1442 }
1443 };
1444
1445 for (size_t n = 0; n < objects.size(); ++n) {
1446 const auto& obj = objects.at(n);
1447 const double deta = obj.polarP4().eta() - tau.polarP4().eta();
1448 const double dphi = reco::deltaPhi(obj.polarP4().phi(), tau.polarP4().phi());
1449 const double dR2 = std::pow(deta, 2) + std::pow(dphi, 2);
1450 if (dR2 < inner_dR2)
1451 addObject(n, deta, dphi, inner_grid);
1452 if (dR2 < outer_dR2)
1453 addObject(n, deta, dphi, outer_grid);
1454 }
1455 }
1456
1457 tensorflow::Tensor getPartialPredictions(bool is_inner) {
1458 std::vector<tensorflow::Tensor> pred_vector;
1459 if (is_inner) {
1460 tensorflow::run(&(cache_->getSession("inner")),
1461 {
1462 {"input_inner_egamma", *eGammaTensor_.at(is_inner)},
1463 {"input_inner_muon", *muonTensor_.at(is_inner)},
1464 {"input_inner_hadrons", *hadronsTensor_.at(is_inner)},
1465 },
1466 {"inner_all_dropout_4/Identity"},
1467 &pred_vector);
1468 } else {
1469 tensorflow::run(&(cache_->getSession("outer")),
1470 {
1471 {"input_outer_egamma", *eGammaTensor_.at(is_inner)},
1472 {"input_outer_muon", *muonTensor_.at(is_inner)},
1473 {"input_outer_hadrons", *hadronsTensor_.at(is_inner)},
1474 },
1475 {"outer_all_dropout_4/Identity"},
1476 &pred_vector);
1477 }
1478 return pred_vector.at(0);
1479 }
1480
1481 template <typename CandidateCastType, typename TauCastType>
1482 void createConvFeatures(const TauCastType& tau,
1483 const size_t tau_index,
1484 const edm::RefToBase<reco::BaseTau> tau_ref,
1485 const reco::Vertex& pv,
1486 double rho,
1487 const std::vector<pat::Electron>* electrons,
1488 const std::vector<pat::Muon>* muons,
1489 const edm::View<reco::Candidate>& pfCands,
1490 const CellGrid& grid,
1491 TauFunc tau_funcs,
1492 bool is_inner) {
1493 if (debug_level >= 2) {
1494 std::cout << "<DeepTauId::createConvFeatures (is_inner = " << is_inner << ")>:" << std::endl;
1495 }
1496 tensorflow::Tensor& convTensor = *convTensor_.at(is_inner);
1497 eGammaTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1498 tensorflow::DT_FLOAT,
1499 tensorflow::TensorShape{
1500 (long long int)grid.num_valid_cells(), 1, 1, dnn_inputs_v2::EgammaBlockInputs::NumberOfInputs});
1501 muonTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1502 tensorflow::DT_FLOAT,
1503 tensorflow::TensorShape{
1504 (long long int)grid.num_valid_cells(), 1, 1, dnn_inputs_v2::MuonBlockInputs::NumberOfInputs});
1505 hadronsTensor_[is_inner] = std::make_unique<tensorflow::Tensor>(
1506 tensorflow::DT_FLOAT,
1507 tensorflow::TensorShape{
1508 (long long int)grid.num_valid_cells(), 1, 1, dnn_inputs_v2::HadronBlockInputs::NumberOfInputs});
1509
1510 eGammaTensor_[is_inner]->flat<float>().setZero();
1511 muonTensor_[is_inner]->flat<float>().setZero();
1512 hadronsTensor_[is_inner]->flat<float>().setZero();
1513
1514 unsigned idx = 0;
1515 for (int eta = -grid.maxEtaIndex(); eta <= grid.maxEtaIndex(); ++eta) {
1516 for (int phi = -grid.maxPhiIndex(); phi <= grid.maxPhiIndex(); ++phi) {
1517 if (debug_level >= 2) {
1518 std::cout << "processing ( eta = " << eta << ", phi = " << phi << " )" << std::endl;
1519 }
1520 const CellIndex cell_index{eta, phi};
1521 const auto cell_iter = grid.find(cell_index);
1522 if (cell_iter != grid.end()) {
1523 if (debug_level >= 2) {
1524 std::cout << " creating inputs for ( eta = " << eta << ", phi = " << phi << " ): idx = " << idx
1525 << std::endl;
1526 }
1527 const Cell& cell = cell_iter->second;
1528 createEgammaBlockInputs<CandidateCastType>(
1529 idx, tau, tau_index, tau_ref, pv, rho, electrons, pfCands, cell, tau_funcs, is_inner);
1530 createMuonBlockInputs<CandidateCastType>(
1531 idx, tau, tau_index, tau_ref, pv, rho, muons, pfCands, cell, tau_funcs, is_inner);
1532 createHadronsBlockInputs<CandidateCastType>(
1533 idx, tau, tau_index, tau_ref, pv, rho, pfCands, cell, tau_funcs, is_inner);
1534 idx += 1;
1535 } else {
1536 if (debug_level >= 2) {
1537 std::cout << " skipping creation of inputs, because ( eta = " << eta << ", phi = " << phi
1538 << " ) is not in the grid !!" << std::endl;
1539 }
1540 }
1541 }
1542 }
1543
1544 const auto predTensor = getPartialPredictions(is_inner);
1545 idx = 0;
1546 for (int eta = -grid.maxEtaIndex(); eta <= grid.maxEtaIndex(); ++eta) {
1547 for (int phi = -grid.maxPhiIndex(); phi <= grid.maxPhiIndex(); ++phi) {
1548 const CellIndex cell_index{eta, phi};
1549 const int eta_index = grid.getEtaTensorIndex(cell_index);
1550 const int phi_index = grid.getPhiTensorIndex(cell_index);
1551
1552 const auto cell_iter = grid.find(cell_index);
1553 if (cell_iter != grid.end()) {
1554 setCellConvFeatures(convTensor, predTensor, idx, eta_index, phi_index);
1555 idx += 1;
1556 } else {
1557 setCellConvFeatures(convTensor, *zeroOutputTensor_[is_inner], 0, eta_index, phi_index);
1558 }
1559 }
1560 }
1561 }
1562
1563 void setCellConvFeatures(tensorflow::Tensor& convTensor,
1564 const tensorflow::Tensor& features,
1565 unsigned batch_idx,
1566 int eta_index,
1567 int phi_index) {
1568 for (int n = 0; n < dnn_inputs_v2::number_of_conv_features; ++n) {
1569 convTensor.tensor<float, 4>()(0, eta_index, phi_index, n) = features.tensor<float, 4>()(batch_idx, 0, 0, n);
1570 }
1571 }
1572
1573 template <typename CandidateCastType, typename TauCastType>
1574 void createTauBlockInputs(const TauCastType& tau,
1575 const size_t& tau_index,
1576 const edm::RefToBase<reco::BaseTau> tau_ref,
1577 const reco::Vertex& pv,
1578 double rho,
1579 TauFunc tau_funcs) {
1580 namespace dnn = dnn_inputs_v2::TauBlockInputs;
1581 namespace sc = deep_tau::Scaling;
1582 sc::FeatureT ft = sc::FeatureT::TauFlat;
1583 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft, false));
1584
1585 tensorflow::Tensor& inputs = *tauBlockTensor_;
1586 inputs.flat<float>().setZero();
1587
1588 const auto& get = [&](int var_index) -> float& {
1589 return inputs.matrix<float>()(0, tauInputs_indices_.at(var_index));
1590 };
1591
1592 auto leadChargedHadrCand = dynamic_cast<const CandidateCastType*>(tau.leadChargedHadrCand().get());
1593
1594 get(dnn::rho) = sp.scale(rho, tauInputs_indices_[dnn::rho]);
1595 get(dnn::tau_pt) = sp.scale(tau.polarP4().pt(), tauInputs_indices_[dnn::tau_pt]);
1596 get(dnn::tau_eta) = sp.scale(tau.polarP4().eta(), tauInputs_indices_[dnn::tau_eta]);
1597 if (sub_version_ == 1) {
1598 get(dnn::tau_phi) = getValueLinear(tau.polarP4().phi(), -pi, pi, false);
1599 }
1600 get(dnn::tau_mass) = sp.scale(tau.polarP4().mass(), tauInputs_indices_[dnn::tau_mass]);
1601 get(dnn::tau_E_over_pt) = sp.scale(tau.p4().energy() / tau.p4().pt(), tauInputs_indices_[dnn::tau_E_over_pt]);
1602 get(dnn::tau_charge) = sp.scale(tau.charge(), tauInputs_indices_[dnn::tau_charge]);
1603 get(dnn::tau_n_charged_prongs) = sp.scale(tau.decayMode() / 5 + 1, tauInputs_indices_[dnn::tau_n_charged_prongs]);
1604 get(dnn::tau_n_neutral_prongs) = sp.scale(tau.decayMode() % 5, tauInputs_indices_[dnn::tau_n_neutral_prongs]);
1605 get(dnn::chargedIsoPtSum) =
1606 sp.scale(tau_funcs.getChargedIsoPtSum(tau, tau_ref), tauInputs_indices_[dnn::chargedIsoPtSum]);
1607 get(dnn::chargedIsoPtSumdR03_over_dR05) =
1608 sp.scale(tau_funcs.getChargedIsoPtSumdR03(tau, tau_ref) / tau_funcs.getChargedIsoPtSum(tau, tau_ref),
1609 tauInputs_indices_[dnn::chargedIsoPtSumdR03_over_dR05]);
1610 if (sub_version_ == 1)
1611 get(dnn::footprintCorrection) =
1612 sp.scale(tau_funcs.getFootprintCorrectiondR03(tau, tau_ref), tauInputs_indices_[dnn::footprintCorrection]);
1613 else if (sub_version_ == 5)
1614 get(dnn::footprintCorrection) =
1615 sp.scale(tau_funcs.getFootprintCorrection(tau, tau_ref), tauInputs_indices_[dnn::footprintCorrection]);
1616
1617 get(dnn::neutralIsoPtSum) =
1618 sp.scale(tau_funcs.getNeutralIsoPtSum(tau, tau_ref), tauInputs_indices_[dnn::neutralIsoPtSum]);
1619 get(dnn::neutralIsoPtSumWeight_over_neutralIsoPtSum) =
1620 sp.scale(tau_funcs.getNeutralIsoPtSumWeight(tau, tau_ref) / tau_funcs.getNeutralIsoPtSum(tau, tau_ref),
1621 tauInputs_indices_[dnn::neutralIsoPtSumWeight_over_neutralIsoPtSum]);
1622 get(dnn::neutralIsoPtSumWeightdR03_over_neutralIsoPtSum) =
1623 sp.scale(tau_funcs.getNeutralIsoPtSumdR03Weight(tau, tau_ref) / tau_funcs.getNeutralIsoPtSum(tau, tau_ref),
1624 tauInputs_indices_[dnn::neutralIsoPtSumWeightdR03_over_neutralIsoPtSum]);
1625 get(dnn::neutralIsoPtSumdR03_over_dR05) =
1626 sp.scale(tau_funcs.getNeutralIsoPtSumdR03(tau, tau_ref) / tau_funcs.getNeutralIsoPtSum(tau, tau_ref),
1627 tauInputs_indices_[dnn::neutralIsoPtSumdR03_over_dR05]);
1628 get(dnn::photonPtSumOutsideSignalCone) = sp.scale(tau_funcs.getPhotonPtSumOutsideSignalCone(tau, tau_ref),
1629 tauInputs_indices_[dnn::photonPtSumOutsideSignalCone]);
1630 get(dnn::puCorrPtSum) = sp.scale(tau_funcs.getPuCorrPtSum(tau, tau_ref), tauInputs_indices_[dnn::puCorrPtSum]);
1631
1632
1633
1634
1635 if (sub_version_ == 1) {
1636 if (!disable_dxy_pca_) {
1637 auto const pca = tau_funcs.getdxyPCA(tau, tau_index);
1638 get(dnn::tau_dxy_pca_x) = sp.scale(pca.x(), tauInputs_indices_[dnn::tau_dxy_pca_x]);
1639 get(dnn::tau_dxy_pca_y) = sp.scale(pca.y(), tauInputs_indices_[dnn::tau_dxy_pca_y]);
1640 get(dnn::tau_dxy_pca_z) = sp.scale(pca.z(), tauInputs_indices_[dnn::tau_dxy_pca_z]);
1641 } else {
1642 get(dnn::tau_dxy_pca_x) = 0;
1643 get(dnn::tau_dxy_pca_y) = 0;
1644 get(dnn::tau_dxy_pca_z) = 0;
1645 }
1646 }
1647
1648 const bool tau_dxy_valid =
1649 isAbove(tau_funcs.getdxy(tau, tau_index), -10) && isAbove(tau_funcs.getdxyError(tau, tau_index), 0);
1650 if (tau_dxy_valid) {
1651 get(dnn::tau_dxy_valid) = sp.scale(tau_dxy_valid, tauInputs_indices_[dnn::tau_dxy_valid]);
1652 get(dnn::tau_dxy) = sp.scale(tau_funcs.getdxy(tau, tau_index), tauInputs_indices_[dnn::tau_dxy]);
1653 get(dnn::tau_dxy_sig) =
1654 sp.scale(std::abs(tau_funcs.getdxy(tau, tau_index)) / tau_funcs.getdxyError(tau, tau_index),
1655 tauInputs_indices_[dnn::tau_dxy_sig]);
1656 }
1657 const bool tau_ip3d_valid =
1658 isAbove(tau_funcs.getip3d(tau, tau_index), -10) && isAbove(tau_funcs.getip3dError(tau, tau_index), 0);
1659 if (tau_ip3d_valid) {
1660 get(dnn::tau_ip3d_valid) = sp.scale(tau_ip3d_valid, tauInputs_indices_[dnn::tau_ip3d_valid]);
1661 get(dnn::tau_ip3d) = sp.scale(tau_funcs.getip3d(tau, tau_index), tauInputs_indices_[dnn::tau_ip3d]);
1662 get(dnn::tau_ip3d_sig) =
1663 sp.scale(std::abs(tau_funcs.getip3d(tau, tau_index)) / tau_funcs.getip3dError(tau, tau_index),
1664 tauInputs_indices_[dnn::tau_ip3d_sig]);
1665 }
1666 if (leadChargedHadrCand) {
1667 const bool hasTrackDetails = candFunc::getHasTrackDetails(*leadChargedHadrCand);
1668 const float tau_dz = (is_online_ && !hasTrackDetails) ? 0 : candFunc::getTauDz(*leadChargedHadrCand);
1669 get(dnn::tau_dz) = sp.scale(tau_dz, tauInputs_indices_[dnn::tau_dz]);
1670 get(dnn::tau_dz_sig_valid) =
1671 sp.scale(candFunc::getTauDZSigValid(*leadChargedHadrCand), tauInputs_indices_[dnn::tau_dz_sig_valid]);
1672 const double dzError = hasTrackDetails ? leadChargedHadrCand->dzError() : -999.;
1673 get(dnn::tau_dz_sig) = sp.scale(std::abs(tau_dz) / dzError, tauInputs_indices_[dnn::tau_dz_sig]);
1674 }
1675 get(dnn::tau_flightLength_x) =
1676 sp.scale(tau_funcs.getFlightLength(tau, tau_index).x(), tauInputs_indices_[dnn::tau_flightLength_x]);
1677 get(dnn::tau_flightLength_y) =
1678 sp.scale(tau_funcs.getFlightLength(tau, tau_index).y(), tauInputs_indices_[dnn::tau_flightLength_y]);
1679 get(dnn::tau_flightLength_z) =
1680 sp.scale(tau_funcs.getFlightLength(tau, tau_index).z(), tauInputs_indices_[dnn::tau_flightLength_z]);
1681 if (sub_version_ == 1)
1682 get(dnn::tau_flightLength_sig) = 0.55756444;
1683 else if (sub_version_ == 5)
1684 get(dnn::tau_flightLength_sig) =
1685 sp.scale(tau_funcs.getFlightLengthSig(tau, tau_index), tauInputs_indices_[dnn::tau_flightLength_sig]);
1686
1687 get(dnn::tau_pt_weighted_deta_strip) = sp.scale(reco::tau::pt_weighted_deta_strip(tau, tau.decayMode()),
1688 tauInputs_indices_[dnn::tau_pt_weighted_deta_strip]);
1689
1690 get(dnn::tau_pt_weighted_dphi_strip) = sp.scale(reco::tau::pt_weighted_dphi_strip(tau, tau.decayMode()),
1691 tauInputs_indices_[dnn::tau_pt_weighted_dphi_strip]);
1692 get(dnn::tau_pt_weighted_dr_signal) = sp.scale(reco::tau::pt_weighted_dr_signal(tau, tau.decayMode()),
1693 tauInputs_indices_[dnn::tau_pt_weighted_dr_signal]);
1694 get(dnn::tau_pt_weighted_dr_iso) =
1695 sp.scale(reco::tau::pt_weighted_dr_iso(tau, tau.decayMode()), tauInputs_indices_[dnn::tau_pt_weighted_dr_iso]);
1696 get(dnn::tau_leadingTrackNormChi2) =
1697 sp.scale(tau_funcs.getLeadingTrackNormChi2(tau), tauInputs_indices_[dnn::tau_leadingTrackNormChi2]);
1698 const auto eratio = reco::tau::eratio(tau);
1699 const bool tau_e_ratio_valid = std::isnormal(eratio) && eratio > 0.f;
1700 get(dnn::tau_e_ratio_valid) = sp.scale(tau_e_ratio_valid, tauInputs_indices_[dnn::tau_e_ratio_valid]);
1701 get(dnn::tau_e_ratio) = tau_e_ratio_valid ? sp.scale(eratio, tauInputs_indices_[dnn::tau_e_ratio]) : 0.f;
1702 const double gj_angle_diff = calculateGottfriedJacksonAngleDifference(tau, tau_index, tau_funcs);
1703 const bool tau_gj_angle_diff_valid = (std::isnormal(gj_angle_diff) || gj_angle_diff == 0) && gj_angle_diff >= 0;
1704 get(dnn::tau_gj_angle_diff_valid) =
1705 sp.scale(tau_gj_angle_diff_valid, tauInputs_indices_[dnn::tau_gj_angle_diff_valid]);
1706 get(dnn::tau_gj_angle_diff) =
1707 tau_gj_angle_diff_valid ? sp.scale(gj_angle_diff, tauInputs_indices_[dnn::tau_gj_angle_diff]) : 0;
1708 get(dnn::tau_n_photons) = sp.scale(reco::tau::n_photons_total(tau), tauInputs_indices_[dnn::tau_n_photons]);
1709 get(dnn::tau_emFraction) = sp.scale(tau_funcs.getEmFraction(tau), tauInputs_indices_[dnn::tau_emFraction]);
1710
1711 get(dnn::tau_inside_ecal_crack) =
1712 sp.scale(isInEcalCrack(tau.p4().eta()), tauInputs_indices_[dnn::tau_inside_ecal_crack]);
1713 get(dnn::leadChargedCand_etaAtEcalEntrance_minus_tau_eta) =
1714 sp.scale(tau_funcs.getEtaAtEcalEntrance(tau) - tau.p4().eta(),
1715 tauInputs_indices_[dnn::leadChargedCand_etaAtEcalEntrance_minus_tau_eta]);
1716 }
1717
1718 template <typename CandidateCastType, typename TauCastType>
1719 void createEgammaBlockInputs(unsigned idx,
1720 const TauCastType& tau,
1721 const size_t tau_index,
1722 const edm::RefToBase<reco::BaseTau> tau_ref,
1723 const reco::Vertex& pv,
1724 double rho,
1725 const std::vector<pat::Electron>* electrons,
1726 const edm::View<reco::Candidate>& pfCands,
1727 const Cell& cell_map,
1728 TauFunc tau_funcs,
1729 bool is_inner) {
1730 namespace dnn = dnn_inputs_v2::EgammaBlockInputs;
1731 namespace sc = deep_tau::Scaling;
1732 sc::FeatureT ft_global = sc::FeatureT::GridGlobal;
1733 sc::FeatureT ft_PFe = sc::FeatureT::PfCand_electron;
1734 sc::FeatureT ft_PFg = sc::FeatureT::PfCand_gamma;
1735 sc::FeatureT ft_e = sc::FeatureT::Electron;
1736
1737
1738 int PFe_index_offset = scalingParamsMap_->at(std::make_pair(ft_global, false)).mean_.size();
1739 int e_index_offset = PFe_index_offset + scalingParamsMap_->at(std::make_pair(ft_PFe, false)).mean_.size();
1740 int PFg_index_offset = e_index_offset + scalingParamsMap_->at(std::make_pair(ft_e, false)).mean_.size();
1741
1742
1743 int fill_index_offset_e = 0;
1744 int fill_index_offset_PFg = 0;
1745 if (sub_version_ == 5) {
1746 fill_index_offset_e =
1747 scalingParamsMap_->at(std::make_pair(ft_PFg, false)).mean_.size();
1748 fill_index_offset_PFg =
1749 -scalingParamsMap_->at(std::make_pair(ft_e, false)).mean_.size();
1750 }
1751
1752 tensorflow::Tensor& inputs = *eGammaTensor_.at(is_inner);
1753
1754 const auto& get = [&](int var_index) -> float& { return inputs.tensor<float, 4>()(idx, 0, 0, var_index); };
1755
1756 const bool valid_index_pf_ele = cell_map.count(CellObjectType::PfCand_electron);
1757 const bool valid_index_pf_gamma = cell_map.count(CellObjectType::PfCand_gamma);
1758 const bool valid_index_ele = cell_map.count(CellObjectType::Electron);
1759
1760 if (!cell_map.empty()) {
1761 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_global, false));
1762 get(dnn::rho) = sp.scale(rho, dnn::rho);
1763 get(dnn::tau_pt) = sp.scale(tau.polarP4().pt(), dnn::tau_pt);
1764 get(dnn::tau_eta) = sp.scale(tau.polarP4().eta(), dnn::tau_eta);
1765 get(dnn::tau_inside_ecal_crack) = sp.scale(isInEcalCrack(tau.polarP4().eta()), dnn::tau_inside_ecal_crack);
1766 }
1767 if (valid_index_pf_ele) {
1768 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_PFe, is_inner));
1769 size_t index_pf_ele = cell_map.at(CellObjectType::PfCand_electron);
1770 const auto& ele_cand = dynamic_cast<const CandidateCastType&>(pfCands.at(index_pf_ele));
1771
1772 get(dnn::pfCand_ele_valid) = sp.scale(valid_index_pf_ele, dnn::pfCand_ele_valid - PFe_index_offset);
1773 get(dnn::pfCand_ele_rel_pt) =
1774 sp.scale(ele_cand.polarP4().pt() / tau.polarP4().pt(), dnn::pfCand_ele_rel_pt - PFe_index_offset);
1775 get(dnn::pfCand_ele_deta) =
1776 sp.scale(ele_cand.polarP4().eta() - tau.polarP4().eta(), dnn::pfCand_ele_deta - PFe_index_offset);
1777 get(dnn::pfCand_ele_dphi) =
1778 sp.scale(dPhi(tau.polarP4(), ele_cand.polarP4()), dnn::pfCand_ele_dphi - PFe_index_offset);
1779 get(dnn::pfCand_ele_pvAssociationQuality) = sp.scale<int>(
1780 candFunc::getPvAssocationQuality(ele_cand), dnn::pfCand_ele_pvAssociationQuality - PFe_index_offset);
1781 get(dnn::pfCand_ele_puppiWeight) = is_inner ? sp.scale(candFunc::getPuppiWeight(ele_cand, 0.9906834f),
1782 dnn::pfCand_ele_puppiWeight - PFe_index_offset)
1783 : sp.scale(candFunc::getPuppiWeight(ele_cand, 0.9669586f),
1784 dnn::pfCand_ele_puppiWeight - PFe_index_offset);
1785 get(dnn::pfCand_ele_charge) = sp.scale(ele_cand.charge(), dnn::pfCand_ele_charge - PFe_index_offset);
1786 get(dnn::pfCand_ele_lostInnerHits) =
1787 sp.scale<int>(candFunc::getLostInnerHits(ele_cand, 0), dnn::pfCand_ele_lostInnerHits - PFe_index_offset);
1788 get(dnn::pfCand_ele_numberOfPixelHits) =
1789 sp.scale(candFunc::getNumberOfPixelHits(ele_cand, 0), dnn::pfCand_ele_numberOfPixelHits - PFe_index_offset);
1790 get(dnn::pfCand_ele_vertex_dx) =
1791 sp.scale(ele_cand.vertex().x() - pv.position().x(), dnn::pfCand_ele_vertex_dx - PFe_index_offset);
1792 get(dnn::pfCand_ele_vertex_dy) =
1793 sp.scale(ele_cand.vertex().y() - pv.position().y(), dnn::pfCand_ele_vertex_dy - PFe_index_offset);
1794 get(dnn::pfCand_ele_vertex_dz) =
1795 sp.scale(ele_cand.vertex().z() - pv.position().z(), dnn::pfCand_ele_vertex_dz - PFe_index_offset);
1796 get(dnn::pfCand_ele_vertex_dx_tauFL) =
1797 sp.scale(ele_cand.vertex().x() - pv.position().x() - tau_funcs.getFlightLength(tau, tau_index).x(),
1798 dnn::pfCand_ele_vertex_dx_tauFL - PFe_index_offset);
1799 get(dnn::pfCand_ele_vertex_dy_tauFL) =
1800 sp.scale(ele_cand.vertex().y() - pv.position().y() - tau_funcs.getFlightLength(tau, tau_index).y(),
1801 dnn::pfCand_ele_vertex_dy_tauFL - PFe_index_offset);
1802 get(dnn::pfCand_ele_vertex_dz_tauFL) =
1803 sp.scale(ele_cand.vertex().z() - pv.position().z() - tau_funcs.getFlightLength(tau, tau_index).z(),
1804 dnn::pfCand_ele_vertex_dz_tauFL - PFe_index_offset);
1805
1806 const bool hasTrackDetails = candFunc::getHasTrackDetails(ele_cand);
1807 if (hasTrackDetails) {
1808 get(dnn::pfCand_ele_hasTrackDetails) =
1809 sp.scale(hasTrackDetails, dnn::pfCand_ele_hasTrackDetails - PFe_index_offset);
1810 get(dnn::pfCand_ele_dxy) = sp.scale(candFunc::getTauDxy(ele_cand), dnn::pfCand_ele_dxy - PFe_index_offset);
1811 get(dnn::pfCand_ele_dxy_sig) = sp.scale(std::abs(candFunc::getTauDxy(ele_cand)) / ele_cand.dxyError(),
1812 dnn::pfCand_ele_dxy_sig - PFe_index_offset);
1813 get(dnn::pfCand_ele_dz) = sp.scale(candFunc::getTauDz(ele_cand), dnn::pfCand_ele_dz - PFe_index_offset);
1814 get(dnn::pfCand_ele_dz_sig) = sp.scale(std::abs(candFunc::getTauDz(ele_cand)) / ele_cand.dzError(),
1815 dnn::pfCand_ele_dz_sig - PFe_index_offset);
1816 get(dnn::pfCand_ele_track_chi2_ndof) =
1817 candFunc::getPseudoTrack(ele_cand).ndof() > 0
1818 ? sp.scale(candFunc::getPseudoTrack(ele_cand).chi2() / candFunc::getPseudoTrack(ele_cand).ndof(),
1819 dnn::pfCand_ele_track_chi2_ndof - PFe_index_offset)
1820 : 0;
1821 get(dnn::pfCand_ele_track_ndof) =
1822 candFunc::getPseudoTrack(ele_cand).ndof() > 0
1823 ? sp.scale(candFunc::getPseudoTrack(ele_cand).ndof(), dnn::pfCand_ele_track_ndof - PFe_index_offset)
1824 : 0;
1825 }
1826 }
1827 if (valid_index_pf_gamma) {
1828 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_PFg, is_inner));
1829 size_t index_pf_gamma = cell_map.at(CellObjectType::PfCand_gamma);
1830 const auto& gamma_cand = dynamic_cast<const CandidateCastType&>(pfCands.at(index_pf_gamma));
1831
1832 get(dnn::pfCand_gamma_valid + fill_index_offset_PFg) =
1833 sp.scale(valid_index_pf_gamma, dnn::pfCand_gamma_valid - PFg_index_offset);
1834 get(dnn::pfCand_gamma_rel_pt + fill_index_offset_PFg) =
1835 sp.scale(gamma_cand.polarP4().pt() / tau.polarP4().pt(), dnn::pfCand_gamma_rel_pt - PFg_index_offset);
1836 get(dnn::pfCand_gamma_deta + fill_index_offset_PFg) =
1837 sp.scale(gamma_cand.polarP4().eta() - tau.polarP4().eta(), dnn::pfCand_gamma_deta - PFg_index_offset);
1838 get(dnn::pfCand_gamma_dphi + fill_index_offset_PFg) =
1839 sp.scale(dPhi(tau.polarP4(), gamma_cand.polarP4()), dnn::pfCand_gamma_dphi - PFg_index_offset);
1840 get(dnn::pfCand_gamma_pvAssociationQuality + fill_index_offset_PFg) = sp.scale<int>(
1841 candFunc::getPvAssocationQuality(gamma_cand), dnn::pfCand_gamma_pvAssociationQuality - PFg_index_offset);
1842 get(dnn::pfCand_gamma_fromPV + fill_index_offset_PFg) =
1843 sp.scale<int>(candFunc::getFromPV(gamma_cand), dnn::pfCand_gamma_fromPV - PFg_index_offset);
1844 get(dnn::pfCand_gamma_puppiWeight + fill_index_offset_PFg) =
1845 is_inner ? sp.scale(candFunc::getPuppiWeight(gamma_cand, 0.9084110f),
1846 dnn::pfCand_gamma_puppiWeight - PFg_index_offset)
1847 : sp.scale(candFunc::getPuppiWeight(gamma_cand, 0.4211567f),
1848 dnn::pfCand_gamma_puppiWeight - PFg_index_offset);
1849 get(dnn::pfCand_gamma_puppiWeightNoLep + fill_index_offset_PFg) =
1850 is_inner ? sp.scale(candFunc::getPuppiWeightNoLep(gamma_cand, 0.8857716f),
1851 dnn::pfCand_gamma_puppiWeightNoLep - PFg_index_offset)
1852 : sp.scale(candFunc::getPuppiWeightNoLep(gamma_cand, 0.3822604f),
1853 dnn::pfCand_gamma_puppiWeightNoLep - PFg_index_offset);
1854 get(dnn::pfCand_gamma_lostInnerHits + fill_index_offset_PFg) =
1855 sp.scale<int>(candFunc::getLostInnerHits(gamma_cand, 0), dnn::pfCand_gamma_lostInnerHits - PFg_index_offset);
1856 get(dnn::pfCand_gamma_numberOfPixelHits + fill_index_offset_PFg) = sp.scale(
1857 candFunc::getNumberOfPixelHits(gamma_cand, 0), dnn::pfCand_gamma_numberOfPixelHits - PFg_index_offset);
1858 get(dnn::pfCand_gamma_vertex_dx + fill_index_offset_PFg) =
1859 sp.scale(gamma_cand.vertex().x() - pv.position().x(), dnn::pfCand_gamma_vertex_dx - PFg_index_offset);
1860 get(dnn::pfCand_gamma_vertex_dy + fill_index_offset_PFg) =
1861 sp.scale(gamma_cand.vertex().y() - pv.position().y(), dnn::pfCand_gamma_vertex_dy - PFg_index_offset);
1862 get(dnn::pfCand_gamma_vertex_dz + fill_index_offset_PFg) =
1863 sp.scale(gamma_cand.vertex().z() - pv.position().z(), dnn::pfCand_gamma_vertex_dz - PFg_index_offset);
1864 get(dnn::pfCand_gamma_vertex_dx_tauFL + fill_index_offset_PFg) =
1865 sp.scale(gamma_cand.vertex().x() - pv.position().x() - tau_funcs.getFlightLength(tau, tau_index).x(),
1866 dnn::pfCand_gamma_vertex_dx_tauFL - PFg_index_offset);
1867 get(dnn::pfCand_gamma_vertex_dy_tauFL + fill_index_offset_PFg) =
1868 sp.scale(gamma_cand.vertex().y() - pv.position().y() - tau_funcs.getFlightLength(tau, tau_index).y(),
1869 dnn::pfCand_gamma_vertex_dy_tauFL - PFg_index_offset);
1870 get(dnn::pfCand_gamma_vertex_dz_tauFL + fill_index_offset_PFg) =
1871 sp.scale(gamma_cand.vertex().z() - pv.position().z() - tau_funcs.getFlightLength(tau, tau_index).z(),
1872 dnn::pfCand_gamma_vertex_dz_tauFL - PFg_index_offset);
1873 const bool hasTrackDetails = candFunc::getHasTrackDetails(gamma_cand);
1874 if (hasTrackDetails) {
1875 get(dnn::pfCand_gamma_hasTrackDetails + fill_index_offset_PFg) =
1876 sp.scale(hasTrackDetails, dnn::pfCand_gamma_hasTrackDetails - PFg_index_offset);
1877 get(dnn::pfCand_gamma_dxy + fill_index_offset_PFg) =
1878 sp.scale(candFunc::getTauDxy(gamma_cand), dnn::pfCand_gamma_dxy - PFg_index_offset);
1879 get(dnn::pfCand_gamma_dxy_sig + fill_index_offset_PFg) =
1880 sp.scale(std::abs(candFunc::getTauDxy(gamma_cand)) / gamma_cand.dxyError(),
1881 dnn::pfCand_gamma_dxy_sig - PFg_index_offset);
1882 get(dnn::pfCand_gamma_dz + fill_index_offset_PFg) =
1883 sp.scale(candFunc::getTauDz(gamma_cand), dnn::pfCand_gamma_dz - PFg_index_offset);
1884 get(dnn::pfCand_gamma_dz_sig + fill_index_offset_PFg) =
1885 sp.scale(std::abs(candFunc::getTauDz(gamma_cand)) / gamma_cand.dzError(),
1886 dnn::pfCand_gamma_dz_sig - PFg_index_offset);
1887 get(dnn::pfCand_gamma_track_chi2_ndof + fill_index_offset_PFg) =
1888 candFunc::getPseudoTrack(gamma_cand).ndof() > 0
1889 ? sp.scale(candFunc::getPseudoTrack(gamma_cand).chi2() / candFunc::getPseudoTrack(gamma_cand).ndof(),
1890 dnn::pfCand_gamma_track_chi2_ndof - PFg_index_offset)
1891 : 0;
1892 get(dnn::pfCand_gamma_track_ndof + fill_index_offset_PFg) =
1893 candFunc::getPseudoTrack(gamma_cand).ndof() > 0
1894 ? sp.scale(candFunc::getPseudoTrack(gamma_cand).ndof(), dnn::pfCand_gamma_track_ndof - PFg_index_offset)
1895 : 0;
1896 }
1897 }
1898 if (valid_index_ele) {
1899 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_e, is_inner));
1900 size_t index_ele = cell_map.at(CellObjectType::Electron);
1901 const auto& ele = electrons->at(index_ele);
1902
1903 get(dnn::ele_valid + fill_index_offset_e) = sp.scale(valid_index_ele, dnn::ele_valid - e_index_offset);
1904 get(dnn::ele_rel_pt + fill_index_offset_e) =
1905 sp.scale(ele.polarP4().pt() / tau.polarP4().pt(), dnn::ele_rel_pt - e_index_offset);
1906 get(dnn::ele_deta + fill_index_offset_e) =
1907 sp.scale(ele.polarP4().eta() - tau.polarP4().eta(), dnn::ele_deta - e_index_offset);
1908 get(dnn::ele_dphi + fill_index_offset_e) =
1909 sp.scale(dPhi(tau.polarP4(), ele.polarP4()), dnn::ele_dphi - e_index_offset);
1910
1911 float cc_ele_energy, cc_gamma_energy;
1912 int cc_n_gamma;
1913 const bool cc_valid = calculateElectronClusterVarsV2(ele, cc_ele_energy, cc_gamma_energy, cc_n_gamma);
1914 if (cc_valid) {
1915 get(dnn::ele_cc_valid + fill_index_offset_e) = sp.scale(cc_valid, dnn::ele_cc_valid - e_index_offset);
1916 get(dnn::ele_cc_ele_rel_energy + fill_index_offset_e) =
1917 sp.scale(cc_ele_energy / ele.polarP4().pt(), dnn::ele_cc_ele_rel_energy - e_index_offset);
1918 get(dnn::ele_cc_gamma_rel_energy + fill_index_offset_e) =
1919 sp.scale(cc_gamma_energy / cc_ele_energy, dnn::ele_cc_gamma_rel_energy - e_index_offset);
1920 get(dnn::ele_cc_n_gamma + fill_index_offset_e) = sp.scale(cc_n_gamma, dnn::ele_cc_n_gamma - e_index_offset);
1921 }
1922 get(dnn::ele_rel_trackMomentumAtVtx + fill_index_offset_e) =
1923 sp.scale(ele.trackMomentumAtVtx().R() / ele.polarP4().pt(), dnn::ele_rel_trackMomentumAtVtx - e_index_offset);
1924 get(dnn::ele_rel_trackMomentumAtCalo + fill_index_offset_e) = sp.scale(
1925 ele.trackMomentumAtCalo().R() / ele.polarP4().pt(), dnn::ele_rel_trackMomentumAtCalo - e_index_offset);
1926 get(dnn::ele_rel_trackMomentumOut + fill_index_offset_e) =
1927 sp.scale(ele.trackMomentumOut().R() / ele.polarP4().pt(), dnn::ele_rel_trackMomentumOut - e_index_offset);
1928 get(dnn::ele_rel_trackMomentumAtEleClus + fill_index_offset_e) = sp.scale(
1929 ele.trackMomentumAtEleClus().R() / ele.polarP4().pt(), dnn::ele_rel_trackMomentumAtEleClus - e_index_offset);
1930 get(dnn::ele_rel_trackMomentumAtVtxWithConstraint + fill_index_offset_e) =
1931 sp.scale(ele.trackMomentumAtVtxWithConstraint().R() / ele.polarP4().pt(),
1932 dnn::ele_rel_trackMomentumAtVtxWithConstraint - e_index_offset);
1933 get(dnn::ele_rel_ecalEnergy + fill_index_offset_e) =
1934 sp.scale(ele.ecalEnergy() / ele.polarP4().pt(), dnn::ele_rel_ecalEnergy - e_index_offset);
1935 get(dnn::ele_ecalEnergy_sig + fill_index_offset_e) =
1936 sp.scale(ele.ecalEnergy() / ele.ecalEnergyError(), dnn::ele_ecalEnergy_sig - e_index_offset);
1937 get(dnn::ele_eSuperClusterOverP + fill_index_offset_e) =
1938 sp.scale(ele.eSuperClusterOverP(), dnn::ele_eSuperClusterOverP - e_index_offset);
1939 get(dnn::ele_eSeedClusterOverP + fill_index_offset_e) =
1940 sp.scale(ele.eSeedClusterOverP(), dnn::ele_eSeedClusterOverP - e_index_offset);
1941 get(dnn::ele_eSeedClusterOverPout + fill_index_offset_e) =
1942 sp.scale(ele.eSeedClusterOverPout(), dnn::ele_eSeedClusterOverPout - e_index_offset);
1943 get(dnn::ele_eEleClusterOverPout + fill_index_offset_e) =
1944 sp.scale(ele.eEleClusterOverPout(), dnn::ele_eEleClusterOverPout - e_index_offset);
1945 get(dnn::ele_deltaEtaSuperClusterTrackAtVtx + fill_index_offset_e) =
1946 sp.scale(ele.deltaEtaSuperClusterTrackAtVtx(), dnn::ele_deltaEtaSuperClusterTrackAtVtx - e_index_offset);
1947 get(dnn::ele_deltaEtaSeedClusterTrackAtCalo + fill_index_offset_e) =
1948 sp.scale(ele.deltaEtaSeedClusterTrackAtCalo(), dnn::ele_deltaEtaSeedClusterTrackAtCalo - e_index_offset);
1949 get(dnn::ele_deltaEtaEleClusterTrackAtCalo + fill_index_offset_e) =
1950 sp.scale(ele.deltaEtaEleClusterTrackAtCalo(), dnn::ele_deltaEtaEleClusterTrackAtCalo - e_index_offset);
1951 get(dnn::ele_deltaPhiEleClusterTrackAtCalo + fill_index_offset_e) =
1952 sp.scale(ele.deltaPhiEleClusterTrackAtCalo(), dnn::ele_deltaPhiEleClusterTrackAtCalo - e_index_offset);
1953 get(dnn::ele_deltaPhiSuperClusterTrackAtVtx + fill_index_offset_e) =
1954 sp.scale(ele.deltaPhiSuperClusterTrackAtVtx(), dnn::ele_deltaPhiSuperClusterTrackAtVtx - e_index_offset);
1955 get(dnn::ele_deltaPhiSeedClusterTrackAtCalo + fill_index_offset_e) =
1956 sp.scale(ele.deltaPhiSeedClusterTrackAtCalo(), dnn::ele_deltaPhiSeedClusterTrackAtCalo - e_index_offset);
1957 const bool mva_valid =
1958 (ele.mvaInput().earlyBrem > -2) ||
1959 (year_ !=
1960 2026);
1961 if (mva_valid) {
1962 get(dnn::ele_mvaInput_earlyBrem + fill_index_offset_e) =
1963 sp.scale(ele.mvaInput().earlyBrem, dnn::ele_mvaInput_earlyBrem - e_index_offset);
1964 get(dnn::ele_mvaInput_lateBrem + fill_index_offset_e) =
1965 sp.scale(ele.mvaInput().lateBrem, dnn::ele_mvaInput_lateBrem - e_index_offset);
1966 get(dnn::ele_mvaInput_sigmaEtaEta + fill_index_offset_e) =
1967 sp.scale(ele.mvaInput().sigmaEtaEta, dnn::ele_mvaInput_sigmaEtaEta - e_index_offset);
1968 get(dnn::ele_mvaInput_hadEnergy + fill_index_offset_e) =
1969 sp.scale(ele.mvaInput().hadEnergy, dnn::ele_mvaInput_hadEnergy - e_index_offset);
1970 get(dnn::ele_mvaInput_deltaEta + fill_index_offset_e) =
1971 sp.scale(ele.mvaInput().deltaEta, dnn::ele_mvaInput_deltaEta - e_index_offset);
1972 }
1973 const auto& gsfTrack = ele.gsfTrack();
1974 if (gsfTrack.isNonnull()) {
1975 get(dnn::ele_gsfTrack_normalizedChi2 + fill_index_offset_e) =
1976 sp.scale(gsfTrack->normalizedChi2(), dnn::ele_gsfTrack_normalizedChi2 - e_index_offset);
1977 get(dnn::ele_gsfTrack_numberOfValidHits + fill_index_offset_e) =
1978 sp.scale(gsfTrack->numberOfValidHits(), dnn::ele_gsfTrack_numberOfValidHits - e_index_offset);
1979 get(dnn::ele_rel_gsfTrack_pt + fill_index_offset_e) =
1980 sp.scale(gsfTrack->pt() / ele.polarP4().pt(), dnn::ele_rel_gsfTrack_pt - e_index_offset);
1981 get(dnn::ele_gsfTrack_pt_sig + fill_index_offset_e) =
1982 sp.scale(gsfTrack->pt() / gsfTrack->ptError(), dnn::ele_gsfTrack_pt_sig - e_index_offset);
1983 }
1984 const auto& closestCtfTrack = ele.closestCtfTrackRef();
1985 const bool has_closestCtfTrack = closestCtfTrack.isNonnull();
1986 if (has_closestCtfTrack) {
1987 get(dnn::ele_has_closestCtfTrack + fill_index_offset_e) =
1988 sp.scale(has_closestCtfTrack, dnn::ele_has_closestCtfTrack - e_index_offset);
1989 get(dnn::ele_closestCtfTrack_normalizedChi2 + fill_index_offset_e) =
1990 sp.scale(closestCtfTrack->normalizedChi2(), dnn::ele_closestCtfTrack_normalizedChi2 - e_index_offset);
1991 get(dnn::ele_closestCtfTrack_numberOfValidHits + fill_index_offset_e) =
1992 sp.scale(closestCtfTrack->numberOfValidHits(), dnn::ele_closestCtfTrack_numberOfValidHits - e_index_offset);
1993 }
1994 }
1995 }
1996
1997 template <typename CandidateCastType, typename TauCastType>
1998 void createMuonBlockInputs(unsigned idx,
1999 const TauCastType& tau,
2000 const size_t tau_index,
2001 const edm::RefToBase<reco::BaseTau> tau_ref,
2002 const reco::Vertex& pv,
2003 double rho,
2004 const std::vector<pat::Muon>* muons,
2005 const edm::View<reco::Candidate>& pfCands,
2006 const Cell& cell_map,
2007 TauFunc tau_funcs,
2008 bool is_inner) {
2009 namespace dnn = dnn_inputs_v2::MuonBlockInputs;
2010 namespace sc = deep_tau::Scaling;
2011 sc::FeatureT ft_global = sc::FeatureT::GridGlobal;
2012 sc::FeatureT ft_PFmu = sc::FeatureT::PfCand_muon;
2013 sc::FeatureT ft_mu = sc::FeatureT::Muon;
2014
2015
2016 int PFmu_index_offset = scalingParamsMap_->at(std::make_pair(ft_global, false)).mean_.size();
2017 int mu_index_offset = PFmu_index_offset + scalingParamsMap_->at(std::make_pair(ft_PFmu, false)).mean_.size();
2018
2019 tensorflow::Tensor& inputs = *muonTensor_.at(is_inner);
2020
2021 const auto& get = [&](int var_index) -> float& { return inputs.tensor<float, 4>()(idx, 0, 0, var_index); };
2022
2023 const bool valid_index_pf_muon = cell_map.count(CellObjectType::PfCand_muon);
2024 const bool valid_index_muon = cell_map.count(CellObjectType::Muon);
2025
2026 if (!cell_map.empty()) {
2027 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_global, false));
2028 get(dnn::rho) = sp.scale(rho, dnn::rho);
2029 get(dnn::tau_pt) = sp.scale(tau.polarP4().pt(), dnn::tau_pt);
2030 get(dnn::tau_eta) = sp.scale(tau.polarP4().eta(), dnn::tau_eta);
2031 get(dnn::tau_inside_ecal_crack) = sp.scale(isInEcalCrack(tau.polarP4().eta()), dnn::tau_inside_ecal_crack);
2032 }
2033 if (valid_index_pf_muon) {
2034 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_PFmu, is_inner));
2035 size_t index_pf_muon = cell_map.at(CellObjectType::PfCand_muon);
2036 const auto& muon_cand = dynamic_cast<const CandidateCastType&>(pfCands.at(index_pf_muon));
2037
2038 get(dnn::pfCand_muon_valid) = sp.scale(valid_index_pf_muon, dnn::pfCand_muon_valid - PFmu_index_offset);
2039 get(dnn::pfCand_muon_rel_pt) =
2040 sp.scale(muon_cand.polarP4().pt() / tau.polarP4().pt(), dnn::pfCand_muon_rel_pt - PFmu_index_offset);
2041 get(dnn::pfCand_muon_deta) =
2042 sp.scale(muon_cand.polarP4().eta() - tau.polarP4().eta(), dnn::pfCand_muon_deta - PFmu_index_offset);
2043 get(dnn::pfCand_muon_dphi) =
2044 sp.scale(dPhi(tau.polarP4(), muon_cand.polarP4()), dnn::pfCand_muon_dphi - PFmu_index_offset);
2045 get(dnn::pfCand_muon_pvAssociationQuality) = sp.scale<int>(
2046 candFunc::getPvAssocationQuality(muon_cand), dnn::pfCand_muon_pvAssociationQuality - PFmu_index_offset);
2047 get(dnn::pfCand_muon_fromPV) =
2048 sp.scale<int>(candFunc::getFromPV(muon_cand), dnn::pfCand_muon_fromPV - PFmu_index_offset);
2049 get(dnn::pfCand_muon_puppiWeight) = is_inner ? sp.scale(candFunc::getPuppiWeight(muon_cand, 0.9786588f),
2050 dnn::pfCand_muon_puppiWeight - PFmu_index_offset)
2051 : sp.scale(candFunc::getPuppiWeight(muon_cand, 0.8132477f),
2052 dnn::pfCand_muon_puppiWeight - PFmu_index_offset);
2053 get(dnn::pfCand_muon_charge) = sp.scale(muon_cand.charge(), dnn::pfCand_muon_charge - PFmu_index_offset);
2054 get(dnn::pfCand_muon_lostInnerHits) =
2055 sp.scale<int>(candFunc::getLostInnerHits(muon_cand, 0), dnn::pfCand_muon_lostInnerHits - PFmu_index_offset);
2056 get(dnn::pfCand_muon_numberOfPixelHits) = sp.scale(candFunc::getNumberOfPixelHits(muon_cand, 0),
2057 dnn::pfCand_muon_numberOfPixelHits - PFmu_index_offset);
2058 get(dnn::pfCand_muon_vertex_dx) =
2059 sp.scale(muon_cand.vertex().x() - pv.position().x(), dnn::pfCand_muon_vertex_dx - PFmu_index_offset);
2060 get(dnn::pfCand_muon_vertex_dy) =
2061 sp.scale(muon_cand.vertex().y() - pv.position().y(), dnn::pfCand_muon_vertex_dy - PFmu_index_offset);
2062 get(dnn::pfCand_muon_vertex_dz) =
2063 sp.scale(muon_cand.vertex().z() - pv.position().z(), dnn::pfCand_muon_vertex_dz - PFmu_index_offset);
2064 get(dnn::pfCand_muon_vertex_dx_tauFL) =
2065 sp.scale(muon_cand.vertex().x() - pv.position().x() - tau_funcs.getFlightLength(tau, tau_index).x(),
2066 dnn::pfCand_muon_vertex_dx_tauFL - PFmu_index_offset);
2067 get(dnn::pfCand_muon_vertex_dy_tauFL) =
2068 sp.scale(muon_cand.vertex().y() - pv.position().y() - tau_funcs.getFlightLength(tau, tau_index).y(),
2069 dnn::pfCand_muon_vertex_dy_tauFL - PFmu_index_offset);
2070 get(dnn::pfCand_muon_vertex_dz_tauFL) =
2071 sp.scale(muon_cand.vertex().z() - pv.position().z() - tau_funcs.getFlightLength(tau, tau_index).z(),
2072 dnn::pfCand_muon_vertex_dz_tauFL - PFmu_index_offset);
2073
2074 const bool hasTrackDetails = candFunc::getHasTrackDetails(muon_cand);
2075 if (hasTrackDetails) {
2076 get(dnn::pfCand_muon_hasTrackDetails) =
2077 sp.scale(hasTrackDetails, dnn::pfCand_muon_hasTrackDetails - PFmu_index_offset);
2078 get(dnn::pfCand_muon_dxy) = sp.scale(candFunc::getTauDxy(muon_cand), dnn::pfCand_muon_dxy - PFmu_index_offset);
2079 get(dnn::pfCand_muon_dxy_sig) = sp.scale(std::abs(candFunc::getTauDxy(muon_cand)) / muon_cand.dxyError(),
2080 dnn::pfCand_muon_dxy_sig - PFmu_index_offset);
2081 get(dnn::pfCand_muon_dz) = sp.scale(candFunc::getTauDz(muon_cand), dnn::pfCand_muon_dz - PFmu_index_offset);
2082 get(dnn::pfCand_muon_dz_sig) = sp.scale(std::abs(candFunc::getTauDz(muon_cand)) / muon_cand.dzError(),
2083 dnn::pfCand_muon_dz_sig - PFmu_index_offset);
2084 get(dnn::pfCand_muon_track_chi2_ndof) =
2085 candFunc::getPseudoTrack(muon_cand).ndof() > 0
2086 ? sp.scale(candFunc::getPseudoTrack(muon_cand).chi2() / candFunc::getPseudoTrack(muon_cand).ndof(),
2087 dnn::pfCand_muon_track_chi2_ndof - PFmu_index_offset)
2088 : 0;
2089 get(dnn::pfCand_muon_track_ndof) =
2090 candFunc::getPseudoTrack(muon_cand).ndof() > 0
2091 ? sp.scale(candFunc::getPseudoTrack(muon_cand).ndof(), dnn::pfCand_muon_track_ndof - PFmu_index_offset)
2092 : 0;
2093 }
2094 }
2095 if (valid_index_muon) {
2096 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_mu, is_inner));
2097 size_t index_muon = cell_map.at(CellObjectType::Muon);
2098 const auto& muon = muons->at(index_muon);
2099
2100 get(dnn::muon_valid) = sp.scale(valid_index_muon, dnn::muon_valid - mu_index_offset);
2101 get(dnn::muon_rel_pt) = sp.scale(muon.polarP4().pt() / tau.polarP4().pt(), dnn::muon_rel_pt - mu_index_offset);
2102 get(dnn::muon_deta) = sp.scale(muon.polarP4().eta() - tau.polarP4().eta(), dnn::muon_deta - mu_index_offset);
2103 get(dnn::muon_dphi) = sp.scale(dPhi(tau.polarP4(), muon.polarP4()), dnn::muon_dphi - mu_index_offset);
2104 get(dnn::muon_dxy) = sp.scale(muon.dB(pat::Muon::PV2D), dnn::muon_dxy - mu_index_offset);
2105 get(dnn::muon_dxy_sig) =
2106 sp.scale(std::abs(muon.dB(pat::Muon::PV2D)) / muon.edB(pat::Muon::PV2D), dnn::muon_dxy_sig - mu_index_offset);
2107
2108 const bool normalizedChi2_valid = muon.globalTrack().isNonnull() && muon.normChi2() >= 0;
2109 if (normalizedChi2_valid) {
2110 get(dnn::muon_normalizedChi2_valid) =
2111 sp.scale(normalizedChi2_valid, dnn::muon_normalizedChi2_valid - mu_index_offset);
2112 get(dnn::muon_normalizedChi2) = sp.scale(muon.normChi2(), dnn::muon_normalizedChi2 - mu_index_offset);
2113 if (muon.innerTrack().isNonnull())
2114 get(dnn::muon_numberOfValidHits) =
2115 sp.scale(muon.numberOfValidHits(), dnn::muon_numberOfValidHits - mu_index_offset);
2116 }
2117 get(dnn::muon_segmentCompatibility) =
2118 sp.scale(muon.segmentCompatibility(), dnn::muon_segmentCompatibility - mu_index_offset);
2119 get(dnn::muon_caloCompatibility) =
2120 sp.scale(muon.caloCompatibility(), dnn::muon_caloCompatibility - mu_index_offset);
2121
2122 const bool pfEcalEnergy_valid = muon.pfEcalEnergy() >= 0;
2123 if (pfEcalEnergy_valid) {
2124 get(dnn::muon_pfEcalEnergy_valid) =
2125 sp.scale(pfEcalEnergy_valid, dnn::muon_pfEcalEnergy_valid - mu_index_offset);
2126 get(dnn::muon_rel_pfEcalEnergy) =
2127 sp.scale(muon.pfEcalEnergy() / muon.polarP4().pt(), dnn::muon_rel_pfEcalEnergy - mu_index_offset);
2128 }
2129
2130 MuonHitMatchV2 hit_match(muon);
2131 static const std::map<int, std::pair<int, int>> muonMatchHitVars = {
2132 {MuonSubdetId::DT, {dnn::muon_n_matches_DT_1, dnn::muon_n_hits_DT_1}},
2133 {MuonSubdetId::CSC, {dnn::muon_n_matches_CSC_1, dnn::muon_n_hits_CSC_1}},
2134 {MuonSubdetId::RPC, {dnn::muon_n_matches_RPC_1, dnn::muon_n_hits_RPC_1}}};
2135
2136 for (int subdet : hit_match.MuonHitMatchV2::consideredSubdets()) {
2137 const auto& matchHitVar = muonMatchHitVars.at(subdet);
2138 for (int station = MuonHitMatchV2::first_station_id; station <= MuonHitMatchV2::last_station_id; ++station) {
2139 const unsigned n_matches = hit_match.nMatches(subdet, station);
2140 const unsigned n_hits = hit_match.nHits(subdet, station);
2141 get(matchHitVar.first + station - 1) = sp.scale(n_matches, matchHitVar.first + station - 1 - mu_index_offset);
2142 get(matchHitVar.second + station - 1) = sp.scale(n_hits, matchHitVar.second + station - 1 - mu_index_offset);
2143 }
2144 }
2145 }
2146 }
2147
2148 template <typename CandidateCastType, typename TauCastType>
2149 void createHadronsBlockInputs(unsigned idx,
2150 const TauCastType& tau,
2151 const size_t tau_index,
2152 const edm::RefToBase<reco::BaseTau> tau_ref,
2153 const reco::Vertex& pv,
2154 double rho,
2155 const edm::View<reco::Candidate>& pfCands,
2156 const Cell& cell_map,
2157 TauFunc tau_funcs,
2158 bool is_inner) {
2159 namespace dnn = dnn_inputs_v2::HadronBlockInputs;
2160 namespace sc = deep_tau::Scaling;
2161 sc::FeatureT ft_global = sc::FeatureT::GridGlobal;
2162 sc::FeatureT ft_PFchH = sc::FeatureT::PfCand_chHad;
2163 sc::FeatureT ft_PFnH = sc::FeatureT::PfCand_nHad;
2164
2165
2166 int PFchH_index_offset = scalingParamsMap_->at(std::make_pair(ft_global, false)).mean_.size();
2167 int PFnH_index_offset = PFchH_index_offset + scalingParamsMap_->at(std::make_pair(ft_PFchH, false)).mean_.size();
2168
2169 tensorflow::Tensor& inputs = *hadronsTensor_.at(is_inner);
2170
2171 const auto& get = [&](int var_index) -> float& { return inputs.tensor<float, 4>()(idx, 0, 0, var_index); };
2172
2173 const bool valid_chH = cell_map.count(CellObjectType::PfCand_chargedHadron);
2174 const bool valid_nH = cell_map.count(CellObjectType::PfCand_neutralHadron);
2175
2176 if (!cell_map.empty()) {
2177 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_global, false));
2178 get(dnn::rho) = sp.scale(rho, dnn::rho);
2179 get(dnn::tau_pt) = sp.scale(tau.polarP4().pt(), dnn::tau_pt);
2180 get(dnn::tau_eta) = sp.scale(tau.polarP4().eta(), dnn::tau_eta);
2181 get(dnn::tau_inside_ecal_crack) = sp.scale(isInEcalCrack(tau.polarP4().eta()), dnn::tau_inside_ecal_crack);
2182 }
2183 if (valid_chH) {
2184 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_PFchH, is_inner));
2185 size_t index_chH = cell_map.at(CellObjectType::PfCand_chargedHadron);
2186 const auto& chH_cand = dynamic_cast<const CandidateCastType&>(pfCands.at(index_chH));
2187
2188 get(dnn::pfCand_chHad_valid) = sp.scale(valid_chH, dnn::pfCand_chHad_valid - PFchH_index_offset);
2189 get(dnn::pfCand_chHad_rel_pt) =
2190 sp.scale(chH_cand.polarP4().pt() / tau.polarP4().pt(), dnn::pfCand_chHad_rel_pt - PFchH_index_offset);
2191 get(dnn::pfCand_chHad_deta) =
2192 sp.scale(chH_cand.polarP4().eta() - tau.polarP4().eta(), dnn::pfCand_chHad_deta - PFchH_index_offset);
2193 get(dnn::pfCand_chHad_dphi) =
2194 sp.scale(dPhi(tau.polarP4(), chH_cand.polarP4()), dnn::pfCand_chHad_dphi - PFchH_index_offset);
2195 get(dnn::pfCand_chHad_leadChargedHadrCand) =
2196 sp.scale(&chH_cand == dynamic_cast<const CandidateCastType*>(tau.leadChargedHadrCand().get()),
2197 dnn::pfCand_chHad_leadChargedHadrCand - PFchH_index_offset);
2198 get(dnn::pfCand_chHad_pvAssociationQuality) = sp.scale<int>(
2199 candFunc::getPvAssocationQuality(chH_cand), dnn::pfCand_chHad_pvAssociationQuality - PFchH_index_offset);
2200 get(dnn::pfCand_chHad_fromPV) =
2201 sp.scale<int>(candFunc::getFromPV(chH_cand), dnn::pfCand_chHad_fromPV - PFchH_index_offset);
2202 const float default_chH_pw_inner = 0.7614090f;
2203 const float default_chH_pw_outer = 0.1974930f;
2204 get(dnn::pfCand_chHad_puppiWeight) = is_inner ? sp.scale(candFunc::getPuppiWeight(chH_cand, default_chH_pw_inner),
2205 dnn::pfCand_chHad_puppiWeight - PFchH_index_offset)
2206 : sp.scale(candFunc::getPuppiWeight(chH_cand, default_chH_pw_outer),
2207 dnn::pfCand_chHad_puppiWeight - PFchH_index_offset);
2208 get(dnn::pfCand_chHad_puppiWeightNoLep) =
2209 is_inner ? sp.scale(candFunc::getPuppiWeightNoLep(chH_cand, default_chH_pw_inner),
2210 dnn::pfCand_chHad_puppiWeightNoLep - PFchH_index_offset)
2211 : sp.scale(candFunc::getPuppiWeightNoLep(chH_cand, default_chH_pw_outer),
2212 dnn::pfCand_chHad_puppiWeightNoLep - PFchH_index_offset);
2213 get(dnn::pfCand_chHad_charge) = sp.scale(chH_cand.charge(), dnn::pfCand_chHad_charge - PFchH_index_offset);
2214 get(dnn::pfCand_chHad_lostInnerHits) =
2215 sp.scale<int>(candFunc::getLostInnerHits(chH_cand, 0), dnn::pfCand_chHad_lostInnerHits - PFchH_index_offset);
2216 get(dnn::pfCand_chHad_numberOfPixelHits) = sp.scale(candFunc::getNumberOfPixelHits(chH_cand, 0),
2217 dnn::pfCand_chHad_numberOfPixelHits - PFchH_index_offset);
2218 get(dnn::pfCand_chHad_vertex_dx) =
2219 sp.scale(chH_cand.vertex().x() - pv.position().x(), dnn::pfCand_chHad_vertex_dx - PFchH_index_offset);
2220 get(dnn::pfCand_chHad_vertex_dy) =
2221 sp.scale(chH_cand.vertex().y() - pv.position().y(), dnn::pfCand_chHad_vertex_dy - PFchH_index_offset);
2222 get(dnn::pfCand_chHad_vertex_dz) =
2223 sp.scale(chH_cand.vertex().z() - pv.position().z(), dnn::pfCand_chHad_vertex_dz - PFchH_index_offset);
2224 get(dnn::pfCand_chHad_vertex_dx_tauFL) =
2225 sp.scale(chH_cand.vertex().x() - pv.position().x() - tau_funcs.getFlightLength(tau, tau_index).x(),
2226 dnn::pfCand_chHad_vertex_dx_tauFL - PFchH_index_offset);
2227 get(dnn::pfCand_chHad_vertex_dy_tauFL) =
2228 sp.scale(chH_cand.vertex().y() - pv.position().y() - tau_funcs.getFlightLength(tau, tau_index).y(),
2229 dnn::pfCand_chHad_vertex_dy_tauFL - PFchH_index_offset);
2230 get(dnn::pfCand_chHad_vertex_dz_tauFL) =
2231 sp.scale(chH_cand.vertex().z() - pv.position().z() - tau_funcs.getFlightLength(tau, tau_index).z(),
2232 dnn::pfCand_chHad_vertex_dz_tauFL - PFchH_index_offset);
2233
2234 const bool hasTrackDetails = candFunc::getHasTrackDetails(chH_cand);
2235 if (hasTrackDetails) {
2236 get(dnn::pfCand_chHad_hasTrackDetails) =
2237 sp.scale(hasTrackDetails, dnn::pfCand_chHad_hasTrackDetails - PFchH_index_offset);
2238 get(dnn::pfCand_chHad_dxy) =
2239 sp.scale(candFunc::getTauDxy(chH_cand), dnn::pfCand_chHad_dxy - PFchH_index_offset);
2240 get(dnn::pfCand_chHad_dxy_sig) = sp.scale(std::abs(candFunc::getTauDxy(chH_cand)) / chH_cand.dxyError(),
2241 dnn::pfCand_chHad_dxy_sig - PFchH_index_offset);
2242 get(dnn::pfCand_chHad_dz) = sp.scale(candFunc::getTauDz(chH_cand), dnn::pfCand_chHad_dz - PFchH_index_offset);
2243 get(dnn::pfCand_chHad_dz_sig) = sp.scale(std::abs(candFunc::getTauDz(chH_cand)) / chH_cand.dzError(),
2244 dnn::pfCand_chHad_dz_sig - PFchH_index_offset);
2245 get(dnn::pfCand_chHad_track_chi2_ndof) =
2246 candFunc::getPseudoTrack(chH_cand).ndof() > 0
2247 ? sp.scale(candFunc::getPseudoTrack(chH_cand).chi2() / candFunc::getPseudoTrack(chH_cand).ndof(),
2248 dnn::pfCand_chHad_track_chi2_ndof - PFchH_index_offset)
2249 : 0;
2250 get(dnn::pfCand_chHad_track_ndof) =
2251 candFunc::getPseudoTrack(chH_cand).ndof() > 0
2252 ? sp.scale(candFunc::getPseudoTrack(chH_cand).ndof(), dnn::pfCand_chHad_track_ndof - PFchH_index_offset)
2253 : 0;
2254 }
2255 float hcal_fraction = candFunc::getHCalFraction(chH_cand, disable_hcalFraction_workaround_);
2256 get(dnn::pfCand_chHad_hcalFraction) =
2257 sp.scale(hcal_fraction, dnn::pfCand_chHad_hcalFraction - PFchH_index_offset);
2258 get(dnn::pfCand_chHad_rawCaloFraction) =
2259 sp.scale(candFunc::getRawCaloFraction(chH_cand), dnn::pfCand_chHad_rawCaloFraction - PFchH_index_offset);
2260 }
2261 if (valid_nH) {
2262 const sc::ScalingParams& sp = scalingParamsMap_->at(std::make_pair(ft_PFnH, is_inner));
2263 size_t index_nH = cell_map.at(CellObjectType::PfCand_neutralHadron);
2264 const auto& nH_cand = dynamic_cast<const CandidateCastType&>(pfCands.at(index_nH));
2265
2266 get(dnn::pfCand_nHad_valid) = sp.scale(valid_nH, dnn::pfCand_nHad_valid - PFnH_index_offset);
2267 get(dnn::pfCand_nHad_rel_pt) =
2268 sp.scale(nH_cand.polarP4().pt() / tau.polarP4().pt(), dnn::pfCand_nHad_rel_pt - PFnH_index_offset);
2269 get(dnn::pfCand_nHad_deta) =
2270 sp.scale(nH_cand.polarP4().eta() - tau.polarP4().eta(), dnn::pfCand_nHad_deta - PFnH_index_offset);
2271 get(dnn::pfCand_nHad_dphi) =
2272 sp.scale(dPhi(tau.polarP4(), nH_cand.polarP4()), dnn::pfCand_nHad_dphi - PFnH_index_offset);
2273 get(dnn::pfCand_nHad_puppiWeight) = is_inner ? sp.scale(candFunc::getPuppiWeight(nH_cand, 0.9798355f),
2274 dnn::pfCand_nHad_puppiWeight - PFnH_index_offset)
2275 : sp.scale(candFunc::getPuppiWeight(nH_cand, 0.7813260f),
2276 dnn::pfCand_nHad_puppiWeight - PFnH_index_offset);
2277 get(dnn::pfCand_nHad_puppiWeightNoLep) = is_inner
2278 ? sp.scale(candFunc::getPuppiWeightNoLep(nH_cand, 0.9046796f),
2279 dnn::pfCand_nHad_puppiWeightNoLep - PFnH_index_offset)
2280 : sp.scale(candFunc::getPuppiWeightNoLep(nH_cand, 0.6554860f),
2281 dnn::pfCand_nHad_puppiWeightNoLep - PFnH_index_offset);
2282 float hcal_fraction = candFunc::getHCalFraction(nH_cand, disable_hcalFraction_workaround_);
2283 get(dnn::pfCand_nHad_hcalFraction) = sp.scale(hcal_fraction, dnn::pfCand_nHad_hcalFraction - PFnH_index_offset);
2284 }
2285 }
2286
2287 static void calculateElectronClusterVars(const pat::Electron* ele, float& elecEe, float& elecEgamma) {
2288 if (ele) {
2289 elecEe = elecEgamma = 0;
2290 auto superCluster = ele->superCluster();
2291 if (superCluster.isNonnull() && superCluster.isAvailable() && superCluster->clusters().isNonnull() &&
2292 superCluster->clusters().isAvailable()) {
2293 for (auto iter = superCluster->clustersBegin(); iter != superCluster->clustersEnd(); ++iter) {
2294 const double energy = (*iter)->energy();
2295 if (iter == superCluster->clustersBegin())
2296 elecEe += energy;
2297 else
2298 elecEgamma += energy;
2299 }
2300 }
2301 } else {
2302 elecEe = elecEgamma = default_value;
2303 }
2304 }
2305
2306 template <typename CandidateCollection, typename TauCastType>
2307 static void processSignalPFComponents(const TauCastType& tau,
2308 const CandidateCollection& candidates,
2309 LorentzVectorXYZ& p4_inner,
2310 LorentzVectorXYZ& p4_outer,
2311 float& pt_inner,
2312 float& dEta_inner,
2313 float& dPhi_inner,
2314 float& m_inner,
2315 float& pt_outer,
2316 float& dEta_outer,
2317 float& dPhi_outer,
2318 float& m_outer,
2319 float& n_inner,
2320 float& n_outer) {
2321 p4_inner = LorentzVectorXYZ(0, 0, 0, 0);
2322 p4_outer = LorentzVectorXYZ(0, 0, 0, 0);
2323 n_inner = 0;
2324 n_outer = 0;
2325
2326 const double innerSigCone_radius = getInnerSignalConeRadius(tau.pt());
2327 for (const auto& cand : candidates) {
2328 const double dR = reco::deltaR(cand->p4(), tau.leadChargedHadrCand()->p4());
2329 const bool isInside_innerSigCone = dR < innerSigCone_radius;
2330 if (isInside_innerSigCone) {
2331 p4_inner += cand->p4();
2332 ++n_inner;
2333 } else {
2334 p4_outer += cand->p4();
2335 ++n_outer;
2336 }
2337 }
2338
2339 pt_inner = n_inner != 0 ? p4_inner.Pt() : default_value;
2340 dEta_inner = n_inner != 0 ? dEta(p4_inner, tau.p4()) : default_value;
2341 dPhi_inner = n_inner != 0 ? dPhi(p4_inner, tau.p4()) : default_value;
2342 m_inner = n_inner != 0 ? p4_inner.mass() : default_value;
2343
2344 pt_outer = n_outer != 0 ? p4_outer.Pt() : default_value;
2345 dEta_outer = n_outer != 0 ? dEta(p4_outer, tau.p4()) : default_value;
2346 dPhi_outer = n_outer != 0 ? dPhi(p4_outer, tau.p4()) : default_value;
2347 m_outer = n_outer != 0 ? p4_outer.mass() : default_value;
2348 }
2349
2350 template <typename CandidateCollection, typename TauCastType>
2351 static void processIsolationPFComponents(const TauCastType& tau,
2352 const CandidateCollection& candidates,
2353 LorentzVectorXYZ& p4,
2354 float& pt,
2355 float& d_eta,
2356 float& d_phi,
2357 float& m,
2358 float& n) {
2359 p4 = LorentzVectorXYZ(0, 0, 0, 0);
2360 n = 0;
2361
2362 for (const auto& cand : candidates) {
2363 p4 += cand->p4();
2364 ++n;
2365 }
2366
2367 pt = n != 0 ? p4.Pt() : default_value;
2368 d_eta = n != 0 ? dEta(p4, tau.p4()) : default_value;
2369 d_phi = n != 0 ? dPhi(p4, tau.p4()) : default_value;
2370 m = n != 0 ? p4.mass() : default_value;
2371 }
2372
2373 static double getInnerSignalConeRadius(double pt) {
2374 static constexpr double min_pt = 30., min_radius = 0.05, cone_opening_coef = 3.;
2375
2376 return std::max(cone_opening_coef / std::max(pt, min_pt), min_radius);
2377 }
2378
2379
2380 template <typename TauCastType>
2381 static bool calculateGottfriedJacksonAngleDifference(const TauCastType& tau,
2382 const size_t tau_index,
2383 double& gj_diff,
2384 TauFunc tau_funcs) {
2385 if (tau_funcs.getHasSecondaryVertex(tau, tau_index)) {
2386 static constexpr double mTau = 1.77682;
2387 const double mAOne = tau.p4().M();
2388 const double pAOneMag = tau.p();
2389 const double argumentThetaGJmax = (std::pow(mTau, 2) - std::pow(mAOne, 2)) / (2 * mTau * pAOneMag);
2390 const double argumentThetaGJmeasured = tau.p4().Vect().Dot(tau_funcs.getFlightLength(tau, tau_index)) /
2391 (pAOneMag * tau_funcs.getFlightLength(tau, tau_index).R());
2392 if (std::abs(argumentThetaGJmax) <= 1. && std::abs(argumentThetaGJmeasured) <= 1.) {
2393 double thetaGJmax = std::asin(argumentThetaGJmax);
2394 double thetaGJmeasured = std::acos(argumentThetaGJmeasured);
2395 gj_diff = thetaGJmeasured - thetaGJmax;
2396 return true;
2397 }
2398 }
2399 return false;
2400 }
2401
2402 template <typename TauCastType>
2403 static float calculateGottfriedJacksonAngleDifference(const TauCastType& tau,
2404 const size_t tau_index,
2405 TauFunc tau_funcs) {
2406 double gj_diff;
2407 if (calculateGottfriedJacksonAngleDifference(tau, tau_index, gj_diff, tau_funcs))
2408 return static_cast<float>(gj_diff);
2409 return default_value;
2410 }
2411
2412 static bool isInEcalCrack(double eta) {
2413 const double abs_eta = std::abs(eta);
2414 return abs_eta > 1.46 && abs_eta < 1.558;
2415 }
2416
2417 template <typename TauCastType>
2418 static const pat::Electron* findMatchedElectron(const TauCastType& tau,
2419 const std::vector<pat::Electron>* electrons,
2420 double deltaR) {
2421 const double dR2 = deltaR * deltaR;
2422 const pat::Electron* matched_ele = nullptr;
2423 for (const auto& ele : *electrons) {
2424 if (reco::deltaR2(tau.p4(), ele.p4()) < dR2 && (!matched_ele || matched_ele->pt() < ele.pt())) {
2425 matched_ele = &ele;
2426 }
2427 }
2428 return matched_ele;
2429 }
2430
2431 private:
2432 edm::EDGetTokenT<std::vector<pat::Electron>> electrons_token_;
2433 edm::EDGetTokenT<std::vector<pat::Muon>> muons_token_;
2434 edm::EDGetTokenT<double> rho_token_;
2435 edm::EDGetTokenT<reco::TauDiscriminatorContainer> basicTauDiscriminators_inputToken_;
2436 edm::EDGetTokenT<reco::TauDiscriminatorContainer> basicTauDiscriminatorsdR03_inputToken_;
2437 edm::EDGetTokenT<edm::AssociationVector<reco::PFTauRefProd, std::vector<reco::PFTauTransverseImpactParameterRef>>>
2438 pfTauTransverseImpactParameters_token_;
2439 std::string input_layer_, output_layer_;
2440 const unsigned year_;
2441 const unsigned version_;
2442 const unsigned sub_version_;
2443 const int debug_level;
2444 const bool disable_dxy_pca_;
2445 const bool disable_hcalFraction_workaround_;
2446 const bool disable_CellIndex_workaround_;
2447 std::unique_ptr<tensorflow::Tensor> tauBlockTensor_;
2448 std::array<std::unique_ptr<tensorflow::Tensor>, 2> eGammaTensor_, muonTensor_, hadronsTensor_, convTensor_,
2449 zeroOutputTensor_;
2450 const std::map<std::pair<deep_tau::Scaling::FeatureT, bool>, deep_tau::Scaling::ScalingParams>* scalingParamsMap_;
2451 const bool save_inputs_;
2452 std::ofstream* json_file_;
2453 bool is_first_block_;
2454 int file_counter_;
2455 std::vector<int> tauInputs_indices_;
2456
2457
2458 bool discrIndicesMapped_ = false;
2459 std::map<BasicDiscriminator, size_t> basicDiscrIndexMap_;
2460 std::map<BasicDiscriminator, size_t> basicDiscrdR03IndexMap_;
2461 };
2462
2463 #include "FWCore/Framework/interface/MakerMacros.h"
2464 DEFINE_FWK_MODULE(DeepTauId);