Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-09-24 22:51:38 

0001 /*
0002  * \class L2TauTagProducer
0003  *
0004  * L2Tau identification using Convolutional NN.
0005  *
0006  * \author Valeria D'Amante, Università di Siena and INFN Pisa
0007  *         Konstantin Androsov, EPFL and ETHZ
0008 */
0009 #include <memory>
0010 #include <boost/property_tree/json_parser.hpp>
0011 #include <cmath>
0012 #include "FWCore/Framework/interface/stream/EDProducer.h"
0013 #include "FWCore/Framework/interface/ESHandle.h"
0014 #include "FWCore/Framework/interface/Event.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/Framework/interface/Frameworkfwd.h"
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 #include "DataFormats/Math/interface/deltaR.h"
0019 #include "DataFormats/Common/interface/Handle.h"
0020 #include "FWCore/Utilities/interface/InputTag.h"
0021 #include "FWCore/Utilities/interface/isFinite.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0024 #include "PhysicsTools/TensorFlow/interface/TensorFlow.h"
0025 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0026 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0027 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0028 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0029 #include "DataFormats/CaloRecHit/interface/CaloRecHit.h"
0030 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
0031 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0032 #include "DataFormats/EcalDetId/interface/EcalDetIdCollections.h"
0033 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0034 #include "DataFormats/HcalRecHit/interface/HBHERecHit.h"
0035 #include "DataFormats/HcalRecHit/interface/HcalRecHitDefs.h"
0036 #include "DataFormats/HcalRecHit/interface/HFRecHit.h"
0037 #include "DataFormats/HcalRecHit/interface/HORecHit.h"
0038 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0039 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0040 #include "TrackingTools/TrajectoryParametrization/interface/CurvilinearTrajectoryError.h"
0041 #include "RecoTracker/PixelTrackFitting/interface/FitUtils.h"
0042 #include "TrackingTools/TrajectoryParametrization/interface/GlobalTrajectoryParameters.h"
0043 #include "DataFormats/TrackReco/interface/HitPattern.h"
0044 #include "TrackingTools/AnalyticalJacobians/interface/JacobianLocalToCurvilinear.h"
0045 #include "DataFormats/TrajectoryState/interface/LocalTrajectoryParameters.h"
0046 #include "DataFormats/GeometrySurface/interface/Plane.h"
0047 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0048 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0049 #include "CUDADataFormats/SiPixelCluster/interface/gpuClusteringConstants.h"
0050 
0051 #include "CUDADataFormats/Track/interface/PixelTrackUtilities.h"
0052 #include "CUDADataFormats/Track/interface/TrackSoAHeterogeneousHost.h"
0053 #include "CUDADataFormats/Vertex/interface/ZVertexUtilities.h"
0054 #include "CUDADataFormats/Vertex/interface/ZVertexSoAHeterogeneousHost.h"
0055 
0056 namespace L2TauTagNNv1 {
0057   constexpr int nCellEta = 5;
0058   constexpr int nCellPhi = 5;
0059   constexpr int nVars = 31;
0060   constexpr float dR_max = 0.5;
0061   enum class NNInputs {
0062     nVertices = 0,
0063     l1Tau_pt,
0064     l1Tau_eta,
0065     l1Tau_hwIso,
0066     EcalEnergySum,
0067     EcalSize,
0068     EcalEnergyStdDev,
0069     EcalDeltaEta,
0070     EcalDeltaPhi,
0071     EcalChi2,
0072     EcalEnergySumForPositiveChi2,
0073     EcalSizeForPositiveChi2,
0074     HcalEnergySum,
0075     HcalSize,
0076     HcalEnergyStdDev,
0077     HcalDeltaEta,
0078     HcalDeltaPhi,
0079     HcalChi2,
0080     HcalEnergySumForPositiveChi2,
0081     HcalSizeForPositiveChi2,
0082     PatatrackPtSum,
0083     PatatrackSize,
0084     PatatrackSizeWithVertex,
0085     PatatrackPtSumWithVertex,
0086     PatatrackChargeSum,
0087     PatatrackDeltaEta,
0088     PatatrackDeltaPhi,
0089     PatatrackChi2OverNdof,
0090     PatatrackNdof,
0091     PatatrackDxy,
0092     PatatrackDz
0093   };
0094 
0095   const std::map<NNInputs, std::string> varNameMap = {
0096       {NNInputs::nVertices, "nVertices"},
0097       {NNInputs::l1Tau_pt, "l1Tau_pt"},
0098       {NNInputs::l1Tau_eta, "l1Tau_eta"},
0099       {NNInputs::l1Tau_hwIso, "l1Tau_hwIso"},
0100       {NNInputs::EcalEnergySum, "EcalEnergySum"},
0101       {NNInputs::EcalSize, "EcalSize"},
0102       {NNInputs::EcalEnergyStdDev, "EcalEnergyStdDev"},
0103       {NNInputs::EcalDeltaEta, "EcalDeltaEta"},
0104       {NNInputs::EcalDeltaPhi, "EcalDeltaPhi"},
0105       {NNInputs::EcalChi2, "EcalChi2"},
0106       {NNInputs::EcalEnergySumForPositiveChi2, "EcalEnergySumForPositiveChi2"},
0107       {NNInputs::EcalSizeForPositiveChi2, "EcalSizeForPositiveChi2"},
0108       {NNInputs::HcalEnergySum, "HcalEnergySum"},
0109       {NNInputs::HcalSize, "HcalSize"},
0110       {NNInputs::HcalEnergyStdDev, "HcalEnergyStdDev"},
0111       {NNInputs::HcalDeltaEta, "HcalDeltaEta"},
0112       {NNInputs::HcalDeltaPhi, "HcalDeltaPhi"},
0113       {NNInputs::HcalChi2, "HcalChi2"},
0114       {NNInputs::HcalEnergySumForPositiveChi2, "HcalEnergySumForPositiveChi2"},
0115       {NNInputs::HcalSizeForPositiveChi2, "HcalSizeForPositiveChi2"},
0116       {NNInputs::PatatrackPtSum, "PatatrackPtSum"},
0117       {NNInputs::PatatrackSize, "PatatrackSize"},
0118       {NNInputs::PatatrackSizeWithVertex, "PatatrackSizeWithVertex"},
0119       {NNInputs::PatatrackPtSumWithVertex, "PatatrackPtSumWithVertex"},
0120       {NNInputs::PatatrackChargeSum, "PatatrackChargeSum"},
0121       {NNInputs::PatatrackDeltaEta, "PatatrackDeltaEta"},
0122       {NNInputs::PatatrackDeltaPhi, "PatatrackDeltaPhi"},
0123       {NNInputs::PatatrackChi2OverNdof, "PatatrackChi2OverNdof"},
0124       {NNInputs::PatatrackNdof, "PatatrackNdof"},
0125       {NNInputs::PatatrackDxy, "PatatrackDxy"},
0126       {NNInputs::PatatrackDz, "PatatrackDz"}};
0127 }  // namespace L2TauTagNNv1
0128 namespace {
0129   inline float& getCellImpl(
0130       tensorflow::Tensor& cellGridMatrix, int tau_idx, int phi_idx, int eta_idx, L2TauTagNNv1::NNInputs NNInput_idx) {
0131     return cellGridMatrix.tensor<float, 4>()(tau_idx, phi_idx, eta_idx, static_cast<int>(NNInput_idx));
0132   }
0133 }  // namespace
0134 struct normDictElement {
0135   float mean;
0136   float std;
0137   float min;
0138   float max;
0139 };
0140 
0141 struct L2TauNNProducerCacheData {
0142   L2TauNNProducerCacheData() : graphDef(nullptr), session(nullptr) {}
0143   tensorflow::GraphDef* graphDef;
0144   tensorflow::Session* session;
0145   std::vector<normDictElement> normVec;
0146 };
0147 
0148 class L2TauNNProducer : public edm::stream::EDProducer<edm::GlobalCache<L2TauNNProducerCacheData>> {
0149 public:
0150   using TrackSoAHost = pixelTrack::TrackSoAHostPhase1;
0151 
0152   struct caloRecHitCollections {
0153     const HBHERecHitCollection* hbhe;
0154     const HORecHitCollection* ho;
0155     const EcalRecHitCollection* eb;
0156     const EcalRecHitCollection* ee;
0157     const CaloGeometry* geometry;
0158   };
0159 
0160   struct InputDescTau {
0161     std::string CollectionName;
0162     edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> inputToken_;
0163   };
0164 
0165   static constexpr float dR2_max = L2TauTagNNv1::dR_max * L2TauTagNNv1::dR_max;
0166   static constexpr float dEta_width = 2 * L2TauTagNNv1::dR_max / static_cast<float>(L2TauTagNNv1::nCellEta);
0167   static constexpr float dPhi_width = 2 * L2TauTagNNv1::dR_max / static_cast<float>(L2TauTagNNv1::nCellPhi);
0168 
0169   explicit L2TauNNProducer(const edm::ParameterSet&, const L2TauNNProducerCacheData*);
0170   static void fillDescriptions(edm::ConfigurationDescriptions&);
0171   static std::unique_ptr<L2TauNNProducerCacheData> initializeGlobalCache(const edm::ParameterSet&);
0172   static void globalEndJob(L2TauNNProducerCacheData*);
0173 
0174 private:
0175   void checknan(tensorflow::Tensor& tensor, int debugLevel);
0176   void standardizeTensor(tensorflow::Tensor& tensor);
0177   std::vector<float> getTauScore(const tensorflow::Tensor& cellGridMatrix);
0178   void produce(edm::Event& event, const edm::EventSetup& eventsetup) override;
0179   void fillL1TauVars(tensorflow::Tensor& cellGridMatrix, const std::vector<l1t::TauRef>& allTaus);
0180   void fillCaloRecHits(tensorflow::Tensor& cellGridMatrix,
0181                        const std::vector<l1t::TauRef>& allTaus,
0182                        const caloRecHitCollections& caloRecHits);
0183   void fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0184                       const std::vector<l1t::TauRef>& allTaus,
0185                       const TrackSoAHost& patatracks_tsoa,
0186                       const ZVertexSoAHost& patavtx_soa,
0187                       const reco::BeamSpot& beamspot,
0188                       const MagneticField* magfi);
0189   void selectGoodTracksAndVertices(const ZVertexSoAHost& patavtx_soa,
0190                                    const TrackSoAHost& patatracks_tsoa,
0191                                    std::vector<int>& trkGood,
0192                                    std::vector<int>& vtxGood);
0193 
0194   std::pair<float, float> impactParameter(int it,
0195                                           const TrackSoAHost& patatracks_tsoa,
0196                                           float patatrackPhi,
0197                                           const reco::BeamSpot& beamspot,
0198                                           const MagneticField* magfi);
0199   template <typename VPos, typename LVec>
0200   std::tuple<float, float, int, int> getEtaPhiIndices(const VPos& position, const LVec& tau_p4);
0201   template <typename LVec>
0202   std::tuple<float, float, int, int> getEtaPhiIndices(float eta, float phi, const LVec& tau_p4);
0203 
0204 private:
0205   const int debugLevel_;
0206   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tauTriggerToken_;
0207   std::vector<InputDescTau> L1TauDesc_;
0208   const edm::EDGetTokenT<HBHERecHitCollection> hbheToken_;
0209   const edm::EDGetTokenT<HORecHitCollection> hoToken_;
0210   const edm::EDGetTokenT<EcalRecHitCollection> ebToken_;
0211   const edm::EDGetTokenT<EcalRecHitCollection> eeToken_;
0212   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometryToken_;
0213   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> bFieldToken_;
0214   const edm::EDGetTokenT<ZVertexSoAHost> pataVerticesToken_;
0215   const edm::EDGetTokenT<TrackSoAHost> pataTracksToken_;
0216   const edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;
0217   const unsigned int maxVtx_;
0218   const float fractionSumPt2_;
0219   const float minSumPt2_;
0220   const float trackPtMin_;
0221   const float trackPtMax_;
0222   const float trackChi2Max_;
0223   std::string inputTensorName_;
0224   std::string outputTensorName_;
0225   const L2TauNNProducerCacheData* L2cacheData_;
0226 };
0227 
0228 std::unique_ptr<L2TauNNProducerCacheData> L2TauNNProducer::initializeGlobalCache(const edm::ParameterSet& cfg) {
0229   std::unique_ptr<L2TauNNProducerCacheData> cacheData = std::make_unique<L2TauNNProducerCacheData>();
0230   cacheData->normVec.reserve(L2TauTagNNv1::nVars);
0231 
0232   auto const graphPath = edm::FileInPath(cfg.getParameter<std::string>("graphPath")).fullPath();
0233 
0234   cacheData->graphDef = tensorflow::loadGraphDef(graphPath);
0235   cacheData->session = tensorflow::createSession(cacheData->graphDef);
0236 
0237   boost::property_tree::ptree loadPtreeRoot;
0238   auto const normalizationDict = edm::FileInPath(cfg.getParameter<std::string>("normalizationDict")).fullPath();
0239   boost::property_tree::read_json(normalizationDict, loadPtreeRoot);
0240   for (const auto& [key, val] : L2TauTagNNv1::varNameMap) {
0241     boost::property_tree::ptree var = loadPtreeRoot.get_child(val);
0242     normDictElement current_element;
0243     current_element.mean = var.get_child("mean").get_value<float>();
0244     current_element.std = var.get_child("std").get_value<float>();
0245     current_element.min = var.get_child("min").get_value<float>();
0246     current_element.max = var.get_child("max").get_value<float>();
0247     cacheData->normVec.push_back(current_element);
0248   }
0249   return cacheData;
0250 }
0251 void L2TauNNProducer::globalEndJob(L2TauNNProducerCacheData* cacheData) {
0252   if (cacheData->graphDef != nullptr) {
0253     delete cacheData->graphDef;
0254   }
0255   tensorflow::closeSession(cacheData->session);
0256 }
0257 void L2TauNNProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0258   edm::ParameterSetDescription desc;
0259   desc.add<int>("debugLevel", 0)->setComment("set debug level for printing out info");
0260   edm::ParameterSetDescription l1TausPset;
0261   l1TausPset.add<std::string>("L1CollectionName", "DoubleTau")->setComment("Name of collections");
0262   l1TausPset.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"))
0263       ->setComment("Which trigger should the L1 Taus collection pass");
0264   edm::ParameterSet l1TausPSetDefault;
0265   l1TausPSetDefault.addParameter<std::string>("L1CollectionName", "DoubleTau");
0266   l1TausPSetDefault.addParameter<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"));
0267   desc.addVPSet("L1Taus", l1TausPset, {l1TausPSetDefault});
0268   desc.add<edm::InputTag>("hbheInput", edm::InputTag("hltHbhereco"))->setComment("HBHE recHit collection");
0269   desc.add<edm::InputTag>("hoInput", edm::InputTag("hltHoreco"))->setComment("HO recHit Collection");
0270   desc.add<edm::InputTag>("ebInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEB"))->setComment("EB recHit Collection");
0271   desc.add<edm::InputTag>("eeInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEE"))->setComment("EE recHit Collection");
0272   desc.add<edm::InputTag>("pataVertices", edm::InputTag("hltPixelVerticesSoA"))
0273       ->setComment("patatrack vertices collection");
0274   desc.add<edm::InputTag>("pataTracks", edm::InputTag("hltPixelTracksSoA"))->setComment("patatrack collection");
0275   desc.add<edm::InputTag>("BeamSpot", edm::InputTag("hltOnlineBeamSpot"))->setComment("BeamSpot Collection");
0276   desc.add<uint>("maxVtx", 100)->setComment("max output collection size (number of accepted vertices)");
0277   desc.add<double>("fractionSumPt2", 0.3)->setComment("threshold on sumPt2 fraction of the leading vertex");
0278   desc.add<double>("minSumPt2", 0.)->setComment("min sumPt2");
0279   desc.add<double>("track_pt_min", 1.0)->setComment("min track p_T");
0280   desc.add<double>("track_pt_max", 10.0)->setComment("max track p_T");
0281   desc.add<double>("track_chi2_max", 99999.)->setComment("max track chi2");
0282   desc.add<std::string>("graphPath", "RecoTauTag/TrainingFiles/data/L2TauNNTag/L2TauTag_Run3v1.pb")
0283       ->setComment("path to the saved CNN");
0284   desc.add<std::string>("normalizationDict", "RecoTauTag/TrainingFiles/data/L2TauNNTag/NormalizationDict.json")
0285       ->setComment("path to the dictionary for variable standardization");
0286   descriptions.addWithDefaultLabel(desc);
0287 }
0288 
0289 L2TauNNProducer::L2TauNNProducer(const edm::ParameterSet& cfg, const L2TauNNProducerCacheData* cacheData)
0290     : debugLevel_(cfg.getParameter<int>("debugLevel")),
0291       hbheToken_(consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheInput"))),
0292       hoToken_(consumes<HORecHitCollection>(cfg.getParameter<edm::InputTag>("hoInput"))),
0293       ebToken_(consumes<EcalRecHitCollection>(cfg.getParameter<edm::InputTag>("ebInput"))),
0294       eeToken_(consumes<EcalRecHitCollection>(cfg.getParameter<edm::InputTag>("eeInput"))),
0295       geometryToken_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0296       bFieldToken_(esConsumes<MagneticField, IdealMagneticFieldRecord>()),
0297       pataVerticesToken_(consumes(cfg.getParameter<edm::InputTag>("pataVertices"))),
0298       pataTracksToken_(consumes(cfg.getParameter<edm::InputTag>("pataTracks"))),
0299       beamSpotToken_(consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>("BeamSpot"))),
0300       maxVtx_(cfg.getParameter<uint>("maxVtx")),
0301       fractionSumPt2_(cfg.getParameter<double>("fractionSumPt2")),
0302       minSumPt2_(cfg.getParameter<double>("minSumPt2")),
0303       trackPtMin_(cfg.getParameter<double>("track_pt_min")),
0304       trackPtMax_(cfg.getParameter<double>("track_pt_max")),
0305       trackChi2Max_(cfg.getParameter<double>("track_chi2_max")) {
0306   if (cacheData->graphDef == nullptr) {
0307     throw cms::Exception("InvalidCacheData") << "Invalid Cache Data.";
0308   }
0309   inputTensorName_ = cacheData->graphDef->node(0).name();
0310   outputTensorName_ = cacheData->graphDef->node(cacheData->graphDef->node_size() - 1).name();
0311   L2cacheData_ = cacheData;
0312   std::vector<edm::ParameterSet> L1TauCollections = cfg.getParameter<std::vector<edm::ParameterSet>>("L1Taus");
0313   L1TauDesc_.reserve(L1TauCollections.size());
0314   for (const auto& l1TauInput : L1TauCollections) {
0315     InputDescTau toInsert;
0316     toInsert.CollectionName = l1TauInput.getParameter<std::string>("L1CollectionName");
0317     toInsert.inputToken_ =
0318         consumes<trigger::TriggerFilterObjectWithRefs>(l1TauInput.getParameter<edm::InputTag>("L1TauTrigger"));
0319     L1TauDesc_.push_back(toInsert);
0320   }
0321   for (const auto& desc : L1TauDesc_)
0322     produces<std::vector<float>>(desc.CollectionName);
0323 }
0324 
0325 void L2TauNNProducer::checknan(tensorflow::Tensor& tensor, int debugLevel) {
0326   using NNInputs = L2TauTagNNv1::NNInputs;
0327   std::vector<int> tensor_shape(tensor.shape().dims());
0328   for (int d = 0; d < tensor.shape().dims(); d++) {
0329     tensor_shape.at(d) = tensor.shape().dim_size(d);
0330   }
0331   if (tensor_shape.size() != 4) {
0332     throw cms::Exception("InvalidTensor") << "Tensor shape does not have 4 dimensions!";
0333   }
0334   for (int tau_idx = 0; tau_idx < tensor_shape.at(0); tau_idx++) {
0335     for (int phi_idx = 0; phi_idx < tensor_shape.at(1); phi_idx++) {
0336       for (int eta_idx = 0; eta_idx < tensor_shape.at(2); eta_idx++) {
0337         for (int var_idx = 0; var_idx < tensor_shape.at(3); var_idx++) {
0338           auto getCell = [&](NNInputs input) -> float& {
0339             return getCellImpl(tensor, tau_idx, phi_idx, eta_idx, input);
0340           };
0341           auto nonstd_var = getCell(static_cast<NNInputs>(var_idx));
0342           if (edm::isNotFinite(nonstd_var)) {
0343             edm::LogWarning("InputVar") << "var is nan \nvar name= "
0344                                         << L2TauTagNNv1::varNameMap.at(static_cast<L2TauTagNNv1::NNInputs>(var_idx))
0345                                         << "\t var_idx = " << var_idx << "\t eta_idx = " << eta_idx
0346                                         << "\t phi_idx = " << phi_idx << "\t tau_idx = " << tau_idx;
0347             if (debugLevel > 2) {
0348               edm::LogWarning("InputVar") << "other vars in same cell \n";
0349               if (var_idx + 1 < tensor_shape.at(3))
0350                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 1))
0351                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 1));
0352               if (var_idx + 2 < tensor_shape.at(3))
0353                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 2))
0354                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 2));
0355               if (var_idx + 3 < tensor_shape.at(3))
0356                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 3))
0357                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 3));
0358               if (var_idx + 4 < tensor_shape.at(3))
0359                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 4))
0360                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 4));
0361             }
0362           }
0363         }
0364       }
0365     }
0366   }
0367 }
0368 
0369 void L2TauNNProducer::standardizeTensor(tensorflow::Tensor& tensor) {
0370   using NNInputs = L2TauTagNNv1::NNInputs;
0371   std::vector<int> tensor_shape(tensor.shape().dims());
0372   for (int d = 0; d < tensor.shape().dims(); d++) {
0373     tensor_shape.at(d) = tensor.shape().dim_size(d);
0374   }
0375   if (tensor_shape.size() != 4) {
0376     throw cms::Exception("InvalidTensor") << "Tensor shape does not have 4 dimensions!";
0377   }
0378   for (int tau_idx = 0; tau_idx < tensor_shape.at(0); tau_idx++) {
0379     for (int phi_idx = 0; phi_idx < tensor_shape.at(1); phi_idx++) {
0380       for (int eta_idx = 0; eta_idx < tensor_shape.at(2); eta_idx++) {
0381         for (int var_idx = 0; var_idx < tensor_shape.at(3); var_idx++) {
0382           auto getCell = [&](NNInputs input) -> float& {
0383             return getCellImpl(tensor, tau_idx, phi_idx, eta_idx, input);
0384           };
0385           float mean = L2cacheData_->normVec.at(var_idx).mean;
0386           float std = L2cacheData_->normVec.at(var_idx).std;
0387           float min = L2cacheData_->normVec.at(var_idx).min;
0388           float max = L2cacheData_->normVec.at(var_idx).max;
0389           float nonstd_var = getCell(static_cast<NNInputs>(var_idx));
0390           float std_var = static_cast<float>((nonstd_var - mean) / std);
0391           if (std_var > max) {
0392             std_var = static_cast<float>(max);
0393           } else if (std_var < min) {
0394             std_var = static_cast<float>(min);
0395           }
0396           getCell(static_cast<NNInputs>(var_idx)) = std_var;
0397         }
0398       }
0399     }
0400   }
0401 }
0402 
0403 void L2TauNNProducer::fillL1TauVars(tensorflow::Tensor& cellGridMatrix, const std::vector<l1t::TauRef>& allTaus) {
0404   using NNInputs = L2TauTagNNv1::NNInputs;
0405 
0406   const int nTaus = allTaus.size();
0407   for (int tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0408     for (int eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0409       for (int phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0410         auto getCell = [&](NNInputs input) -> float& {
0411           return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0412         };
0413         getCell(NNInputs::l1Tau_pt) = allTaus[tau_idx]->pt();
0414         getCell(NNInputs::l1Tau_eta) = allTaus[tau_idx]->eta();
0415         getCell(NNInputs::l1Tau_hwIso) = allTaus[tau_idx]->hwIso();
0416       }
0417     }
0418   }
0419 }
0420 
0421 template <typename LVec>
0422 std::tuple<float, float, int, int> L2TauNNProducer::getEtaPhiIndices(float eta, float phi, const LVec& tau_p4) {
0423   const float deta = eta - tau_p4.eta();
0424   const float dphi = reco::deltaPhi(phi, tau_p4.phi());
0425   const int eta_idx = static_cast<int>(floor((deta + L2TauTagNNv1::dR_max) / dEta_width));
0426   const int phi_idx = static_cast<int>(floor((dphi + L2TauTagNNv1::dR_max) / dPhi_width));
0427   return std::make_tuple(deta, dphi, eta_idx, phi_idx);
0428 }
0429 
0430 template <typename VPos, typename LVec>
0431 std::tuple<float, float, int, int> L2TauNNProducer::getEtaPhiIndices(const VPos& position, const LVec& tau_p4) {
0432   return getEtaPhiIndices(position.eta(), position.phi(), tau_p4);
0433 }
0434 
0435 void L2TauNNProducer::fillCaloRecHits(tensorflow::Tensor& cellGridMatrix,
0436                                       const std::vector<l1t::TauRef>& allTaus,
0437                                       const caloRecHitCollections& caloRecHits) {
0438   using NNInputs = L2TauTagNNv1::NNInputs;
0439 
0440   const int nTaus = allTaus.size();
0441   float deta, dphi;
0442   int eta_idx = 0;
0443   int phi_idx = 0;
0444   int tau_idx = 0;
0445 
0446   auto getCell = [&](NNInputs input) -> float& {
0447     return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0448   };
0449   for (tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0450     // calorechit_EE
0451     for (const auto& caloRecHit_ee : *caloRecHits.ee) {
0452       if (caloRecHit_ee.energy() <= 0)
0453         continue;
0454       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_ee.id())->getPosition();
0455       const float eeCalEn = caloRecHit_ee.energy();
0456       const float eeCalChi2 = caloRecHit_ee.chi2();
0457       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0458         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0459         getCell(NNInputs::EcalEnergySum) += eeCalEn;
0460         getCell(NNInputs::EcalSize) += 1.;
0461         getCell(NNInputs::EcalEnergyStdDev) += eeCalEn * eeCalEn;
0462         getCell(NNInputs::EcalDeltaEta) += deta * eeCalEn;
0463         getCell(NNInputs::EcalDeltaPhi) += dphi * eeCalEn;
0464         if (eeCalChi2 >= 0) {
0465           getCell(NNInputs::EcalChi2) += eeCalChi2 * eeCalEn;
0466           getCell(NNInputs::EcalEnergySumForPositiveChi2) += eeCalEn;
0467           getCell(NNInputs::EcalSizeForPositiveChi2) += 1.;
0468         }
0469       }
0470     }
0471 
0472     // calorechit_EB
0473     for (const auto& caloRecHit_eb : *caloRecHits.eb) {
0474       if (caloRecHit_eb.energy() <= 0)
0475         continue;
0476       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_eb.id())->getPosition();
0477       const float ebCalEn = caloRecHit_eb.energy();
0478       const float ebCalChi2 = caloRecHit_eb.chi2();
0479       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0480         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0481         getCell(NNInputs::EcalEnergySum) += ebCalEn;
0482         getCell(NNInputs::EcalSize) += 1.;
0483         getCell(NNInputs::EcalEnergyStdDev) += ebCalEn * ebCalEn;
0484         getCell(NNInputs::EcalDeltaEta) += deta * ebCalEn;
0485         getCell(NNInputs::EcalDeltaPhi) += dphi * ebCalEn;
0486         if (ebCalChi2 >= 0) {
0487           getCell(NNInputs::EcalChi2) += ebCalChi2 * ebCalEn;
0488           getCell(NNInputs::EcalEnergySumForPositiveChi2) += ebCalEn;
0489           getCell(NNInputs::EcalSizeForPositiveChi2) += 1.;
0490         }
0491       }
0492     }
0493 
0494     // calorechit_HBHE
0495     for (const auto& caloRecHit_hbhe : *caloRecHits.hbhe) {
0496       if (caloRecHit_hbhe.energy() <= 0)
0497         continue;
0498       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_hbhe.id())->getPosition();
0499       const float hbheCalEn = caloRecHit_hbhe.energy();
0500       const float hbheCalChi2 = caloRecHit_hbhe.chi2();
0501       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0502         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0503         getCell(NNInputs::HcalEnergySum) += hbheCalEn;
0504         getCell(NNInputs::HcalEnergyStdDev) += hbheCalEn * hbheCalEn;
0505         getCell(NNInputs::HcalSize) += 1.;
0506         getCell(NNInputs::HcalDeltaEta) += deta * hbheCalEn;
0507         getCell(NNInputs::HcalDeltaPhi) += dphi * hbheCalEn;
0508         if (hbheCalChi2 >= 0) {
0509           getCell(NNInputs::HcalChi2) += hbheCalChi2 * hbheCalEn;
0510           getCell(NNInputs::HcalEnergySumForPositiveChi2) += hbheCalEn;
0511           getCell(NNInputs::HcalSizeForPositiveChi2) += 1.;
0512         }
0513       }
0514     }
0515 
0516     // calorechit_HO
0517     for (const auto& caloRecHit_ho : *caloRecHits.ho) {
0518       if (caloRecHit_ho.energy() <= 0)
0519         continue;
0520       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_ho.id())->getPosition();
0521       const float hoCalEn = caloRecHit_ho.energy();
0522       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0523         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0524         getCell(NNInputs::HcalEnergySum) += hoCalEn;
0525         getCell(NNInputs::HcalEnergyStdDev) += hoCalEn * hoCalEn;
0526         getCell(NNInputs::HcalSize) += 1.;
0527         getCell(NNInputs::HcalDeltaEta) += deta * hoCalEn;
0528         getCell(NNInputs::HcalDeltaPhi) += dphi * hoCalEn;
0529       }
0530     }
0531 
0532     // normalize to sum and define stdDev
0533     for (eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0534       for (phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0535         /* normalize eCal vars*/
0536         if (getCell(NNInputs::EcalEnergySum) > 0.) {
0537           getCell(NNInputs::EcalDeltaEta) /= getCell(NNInputs::EcalEnergySum);
0538           getCell(NNInputs::EcalDeltaPhi) /= getCell(NNInputs::EcalEnergySum);
0539         }
0540         if (getCell(NNInputs::EcalEnergySumForPositiveChi2) > 0.) {
0541           getCell(NNInputs::EcalChi2) /= getCell(NNInputs::EcalEnergySumForPositiveChi2);
0542         }
0543         if (getCell(NNInputs::EcalSize) > 1.) {
0544           // (stdDev - (enSum*enSum)/size) / (size-1)
0545           getCell(NNInputs::EcalEnergyStdDev) =
0546               (getCell(NNInputs::EcalEnergyStdDev) -
0547                (getCell(NNInputs::EcalEnergySum) * getCell(NNInputs::EcalEnergySum)) / getCell(NNInputs::EcalSize)) /
0548               (getCell(NNInputs::EcalSize) - 1);
0549         } else {
0550           getCell(NNInputs::EcalEnergyStdDev) = 0.;
0551         }
0552         /* normalize hCal Vars */
0553         if (getCell(NNInputs::HcalEnergySum) > 0.) {
0554           getCell(NNInputs::HcalDeltaEta) /= getCell(NNInputs::HcalEnergySum);
0555           getCell(NNInputs::HcalDeltaPhi) /= getCell(NNInputs::HcalEnergySum);
0556         }
0557         if (getCell(NNInputs::HcalEnergySumForPositiveChi2) > 0.) {
0558           getCell(NNInputs::HcalChi2) /= getCell(NNInputs::HcalEnergySumForPositiveChi2);
0559         }
0560         if (getCell(NNInputs::HcalSize) > 1.) {
0561           // (stdDev - (enSum*enSum)/size) / (size-1)
0562           getCell(NNInputs::HcalEnergyStdDev) =
0563               (getCell(NNInputs::HcalEnergyStdDev) -
0564                (getCell(NNInputs::HcalEnergySum) * getCell(NNInputs::HcalEnergySum)) / getCell(NNInputs::HcalSize)) /
0565               (getCell(NNInputs::HcalSize) - 1);
0566         } else {
0567           getCell(NNInputs::HcalEnergyStdDev) = 0.;
0568         }
0569       }
0570     }
0571   }
0572 }
0573 
0574 void L2TauNNProducer::selectGoodTracksAndVertices(const ZVertexSoAHost& patavtx_soa,
0575                                                   const TrackSoAHost& patatracks_tsoa,
0576                                                   std::vector<int>& trkGood,
0577                                                   std::vector<int>& vtxGood) {
0578   using patatrackHelpers = TracksUtilities<pixelTopology::Phase1>;
0579   const auto maxTracks = patatracks_tsoa.view().metadata().size();
0580   const int nv = patavtx_soa.view().nvFinal();
0581   trkGood.clear();
0582   trkGood.reserve(maxTracks);
0583   vtxGood.clear();
0584   vtxGood.reserve(nv);
0585   auto const* quality = patatracks_tsoa.view().quality();
0586 
0587   // No need to sort either as the algorithms is just using the max (not even the location, just the max value of pt2sum).
0588   std::vector<float> pTSquaredSum(nv, 0);
0589   std::vector<int> nTrkAssociated(nv, 0);
0590 
0591   for (int32_t trk_idx = 0; trk_idx < maxTracks; ++trk_idx) {
0592     auto nHits = patatrackHelpers::nHits(patatracks_tsoa.view(), trk_idx);
0593     if (nHits == 0) {
0594       break;
0595     }
0596     int vtx_ass_to_track = patavtx_soa.view()[trk_idx].idv();
0597     if (vtx_ass_to_track >= 0 && vtx_ass_to_track < nv) {
0598       auto patatrackPt = patatracks_tsoa.view()[trk_idx].pt();
0599       ++nTrkAssociated[vtx_ass_to_track];
0600       if (patatrackPt >= trackPtMin_ && patatracks_tsoa.const_view()[trk_idx].chi2() <= trackChi2Max_) {
0601         patatrackPt = std::min(patatrackPt, trackPtMax_);
0602         pTSquaredSum[vtx_ass_to_track] += patatrackPt * patatrackPt;
0603       }
0604     }
0605     if (nHits > 0 and quality[trk_idx] >= pixelTrack::Quality::loose) {
0606       trkGood.push_back(trk_idx);
0607     }
0608   }
0609   if (nv > 0) {
0610     const auto minFOM_fromFrac = (*std::max_element(pTSquaredSum.begin(), pTSquaredSum.end())) * fractionSumPt2_;
0611     for (int j = nv - 1; j >= 0 && vtxGood.size() < maxVtx_; --j) {
0612       auto vtx_idx = patavtx_soa.view()[j].sortInd();
0613       assert(vtx_idx < nv);
0614       if (nTrkAssociated[vtx_idx] >= 2 && pTSquaredSum[vtx_idx] >= minFOM_fromFrac &&
0615           pTSquaredSum[vtx_idx] > minSumPt2_) {
0616         vtxGood.push_back(vtx_idx);
0617       }
0618     }
0619   }
0620 }
0621 
0622 std::pair<float, float> L2TauNNProducer::impactParameter(int it,
0623                                                          const TrackSoAHost& patatracks_tsoa,
0624                                                          float patatrackPhi,
0625                                                          const reco::BeamSpot& beamspot,
0626                                                          const MagneticField* magfi) {
0627   /* dxy and dz */
0628   riemannFit::Vector5d ipar, opar;
0629   riemannFit::Matrix5d icov, ocov;
0630   TracksUtilities<pixelTopology::Phase1>::copyToDense(patatracks_tsoa.view(), ipar, icov, it);
0631   riemannFit::transformToPerigeePlane(ipar, icov, opar, ocov);
0632   LocalTrajectoryParameters lpar(opar(0), opar(1), opar(2), opar(3), opar(4), 1.);
0633   float sp = std::sin(patatrackPhi);
0634   float cp = std::cos(patatrackPhi);
0635   Surface::RotationType Rotation(sp, -cp, 0, 0, 0, -1.f, cp, sp, 0);
0636   GlobalPoint BeamSpotPoint(beamspot.x0(), beamspot.y0(), beamspot.z0());
0637   Plane impPointPlane(BeamSpotPoint, Rotation);
0638   GlobalTrajectoryParameters gp(
0639       impPointPlane.toGlobal(lpar.position()), impPointPlane.toGlobal(lpar.momentum()), lpar.charge(), magfi);
0640   GlobalPoint vv = gp.position();
0641   math::XYZPoint pos(vv.x(), vv.y(), vv.z());
0642   GlobalVector pp = gp.momentum();
0643   math::XYZVector mom(pp.x(), pp.y(), pp.z());
0644   auto lambda = M_PI_2 - pp.theta();
0645   auto phi = pp.phi();
0646   float patatrackDxy = -vv.x() * std::sin(phi) + vv.y() * std::cos(phi);
0647   float patatrackDz =
0648       (vv.z() * std::cos(lambda) - (vv.x() * std::cos(phi) + vv.y() * std::sin(phi)) * std::sin(lambda)) /
0649       std::cos(lambda);
0650   return std::make_pair(patatrackDxy, patatrackDz);
0651 }
0652 
0653 void L2TauNNProducer::fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0654                                      const std::vector<l1t::TauRef>& allTaus,
0655                                      const TrackSoAHost& patatracks_tsoa,
0656                                      const ZVertexSoAHost& patavtx_soa,
0657                                      const reco::BeamSpot& beamspot,
0658                                      const MagneticField* magfi) {
0659   using NNInputs = L2TauTagNNv1::NNInputs;
0660   using patatrackHelpers = TracksUtilities<pixelTopology::Phase1>;
0661   float deta, dphi;
0662   int eta_idx = 0;
0663   int phi_idx = 0;
0664   int tau_idx = 0;
0665 
0666   auto getCell = [&](NNInputs input) -> float& {
0667     return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0668   };
0669 
0670   std::vector<int> trkGood;
0671   std::vector<int> vtxGood;
0672 
0673   selectGoodTracksAndVertices(patavtx_soa, patatracks_tsoa, trkGood, vtxGood);
0674 
0675   const int nTaus = allTaus.size();
0676   for (tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0677     const float tauEta = allTaus[tau_idx]->eta();
0678     const float tauPhi = allTaus[tau_idx]->phi();
0679 
0680     for (const auto it : trkGood) {
0681       const float patatrackPt = patatracks_tsoa.const_view()[it].pt();
0682       if (patatrackPt <= 0)
0683         continue;
0684       const float patatrackPhi = patatrackHelpers::phi(patatracks_tsoa.const_view(), it);
0685       const float patatrackEta = patatracks_tsoa.const_view()[it].eta();
0686       const float patatrackCharge = patatrackHelpers::charge(patatracks_tsoa.const_view(), it);
0687       const float patatrackChi2OverNdof = patatracks_tsoa.view()[it].chi2();
0688       const auto nHits = patatrackHelpers::nHits(patatracks_tsoa.const_view(), it);
0689       if (nHits <= 0)
0690         continue;
0691       const int patatrackNdof = 2 * std::min(6, nHits) - 5;
0692 
0693       const int vtx_idx_assTrk = patavtx_soa.view()[it].idv();
0694       if (reco::deltaR2(patatrackEta, patatrackPhi, tauEta, tauPhi) < dR2_max) {
0695         std::tie(deta, dphi, eta_idx, phi_idx) =
0696             getEtaPhiIndices(patatrackEta, patatrackPhi, allTaus[tau_idx]->polarP4());
0697         getCell(NNInputs::PatatrackPtSum) += patatrackPt;
0698         getCell(NNInputs::PatatrackSize) += 1.;
0699         getCell(NNInputs::PatatrackChargeSum) += patatrackCharge;
0700         getCell(NNInputs::PatatrackDeltaEta) += deta * patatrackPt;
0701         getCell(NNInputs::PatatrackDeltaPhi) += dphi * patatrackPt;
0702         getCell(NNInputs::PatatrackChi2OverNdof) += patatrackChi2OverNdof * patatrackPt;
0703         getCell(NNInputs::PatatrackNdof) += patatrackNdof * patatrackPt;
0704         std::pair<float, float> impactParameters = impactParameter(it, patatracks_tsoa, patatrackPhi, beamspot, magfi);
0705         getCell(NNInputs::PatatrackDxy) += impactParameters.first * patatrackPt;
0706         getCell(NNInputs::PatatrackDz) += impactParameters.second * patatrackPt;
0707         if ((std::find(vtxGood.begin(), vtxGood.end(), vtx_idx_assTrk) != vtxGood.end())) {
0708           getCell(NNInputs::PatatrackPtSumWithVertex) += patatrackPt;
0709           getCell(NNInputs::PatatrackSizeWithVertex) += 1.;
0710         }
0711       }
0712     }
0713 
0714     // normalize to sum and define stdDev
0715     for (eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0716       for (phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0717         getCell(NNInputs::nVertices) = vtxGood.size();
0718         if (getCell(NNInputs::PatatrackPtSum) > 0.) {
0719           getCell(NNInputs::PatatrackDeltaEta) /= getCell(NNInputs::PatatrackPtSum);
0720           getCell(NNInputs::PatatrackDeltaPhi) /= getCell(NNInputs::PatatrackPtSum);
0721           getCell(NNInputs::PatatrackChi2OverNdof) /= getCell(NNInputs::PatatrackPtSum);
0722           getCell(NNInputs::PatatrackNdof) /= getCell(NNInputs::PatatrackPtSum);
0723           getCell(NNInputs::PatatrackDxy) /= getCell(NNInputs::PatatrackPtSum);
0724           getCell(NNInputs::PatatrackDz) /= getCell(NNInputs::PatatrackPtSum);
0725         }
0726       }
0727     }
0728   }
0729 }
0730 
0731 std::vector<float> L2TauNNProducer::getTauScore(const tensorflow::Tensor& cellGridMatrix) {
0732   const int nTau = cellGridMatrix.shape().dim_size(0);
0733   std::vector<float> pred_vector(nTau);
0734   if (nTau > 0) {
0735     // Only run the inference if there are taus to process
0736     std::vector<tensorflow::Tensor> pred_tensor;
0737     tensorflow::run(L2cacheData_->session, {{inputTensorName_, cellGridMatrix}}, {outputTensorName_}, &pred_tensor);
0738     for (int tau_idx = 0; tau_idx < nTau; ++tau_idx) {
0739       pred_vector[tau_idx] = pred_tensor[0].matrix<float>()(tau_idx, 0);
0740     }
0741   }
0742   return pred_vector;
0743 }
0744 
0745 void L2TauNNProducer::produce(edm::Event& event, const edm::EventSetup& eventsetup) {
0746   std::vector<std::vector<size_t>> TauCollectionMap(L1TauDesc_.size());
0747   l1t::TauVectorRef allTaus;
0748 
0749   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0750     l1t::TauVectorRef l1Taus;
0751     auto const& l1TriggeredTaus = event.get(L1TauDesc_[inp_idx].inputToken_);
0752     l1TriggeredTaus.getObjects(trigger::TriggerL1Tau, l1Taus);
0753     TauCollectionMap.at(inp_idx).resize(l1Taus.size());
0754 
0755     for (size_t l1_idx = 0; l1_idx < l1Taus.size(); l1_idx++) {
0756       size_t tau_idx;
0757       const auto iter = std::find(allTaus.begin(), allTaus.end(), l1Taus[l1_idx]);
0758       if (iter != allTaus.end()) {
0759         tau_idx = std::distance(allTaus.begin(), iter);
0760       } else {
0761         allTaus.push_back(l1Taus[l1_idx]);
0762         tau_idx = allTaus.size() - 1;
0763       }
0764       TauCollectionMap.at(inp_idx).at(l1_idx) = tau_idx;
0765     }
0766   }
0767   const auto ebCal = event.getHandle(ebToken_);
0768   const auto eeCal = event.getHandle(eeToken_);
0769   const auto hbhe = event.getHandle(hbheToken_);
0770   const auto ho = event.getHandle(hoToken_);
0771   auto const& patatracks_SoA = event.get(pataTracksToken_);
0772   auto const& vertices_SoA = event.get(pataVerticesToken_);
0773   const auto bsHandle = event.getHandle(beamSpotToken_);
0774 
0775   auto const fieldESH = eventsetup.getHandle(bFieldToken_);
0776   auto const geometry = eventsetup.getHandle(geometryToken_);
0777 
0778   caloRecHitCollections caloRecHits;
0779   caloRecHits.hbhe = &*hbhe;
0780   caloRecHits.ho = &*ho;
0781   caloRecHits.eb = &*ebCal;
0782   caloRecHits.ee = &*eeCal;
0783   caloRecHits.geometry = &*geometry;
0784 
0785   const int nTaus = allTaus.size();
0786   tensorflow::Tensor cellGridMatrix(tensorflow::DT_FLOAT,
0787                                     {nTaus, L2TauTagNNv1::nCellEta, L2TauTagNNv1::nCellPhi, L2TauTagNNv1::nVars});
0788   const int n_inputs = nTaus * L2TauTagNNv1::nCellEta * L2TauTagNNv1::nCellPhi * L2TauTagNNv1::nVars;
0789   for (int input_idx = 0; input_idx < n_inputs; ++input_idx) {
0790     cellGridMatrix.flat<float>()(input_idx) = 0;
0791   }
0792   fillL1TauVars(cellGridMatrix, allTaus);
0793 
0794   fillCaloRecHits(cellGridMatrix, allTaus, caloRecHits);
0795 
0796   fillPatatracks(cellGridMatrix, allTaus, patatracks_SoA, vertices_SoA, *bsHandle, fieldESH.product());
0797 
0798   standardizeTensor(cellGridMatrix);
0799 
0800   if (debugLevel_ > 0) {
0801     checknan(cellGridMatrix, debugLevel_);
0802   }
0803 
0804   std::vector<float> tau_score = getTauScore(cellGridMatrix);
0805 
0806   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0807     const size_t nTau = TauCollectionMap[inp_idx].size();
0808     auto tau_tags = std::make_unique<std::vector<float>>(nTau);
0809     for (size_t tau_pos = 0; tau_pos < nTau; ++tau_pos) {
0810       const auto tau_idx = TauCollectionMap[inp_idx][tau_pos];
0811       if (debugLevel_ > 0) {
0812         edm::LogInfo("DebugInfo") << event.id().event() << " \t " << (allTaus[tau_idx])->pt() << " \t "
0813                                   << tau_score.at(tau_idx) << std::endl;
0814       }
0815       (*tau_tags)[tau_pos] = tau_score.at(tau_idx);
0816     }
0817     event.put(std::move(tau_tags), L1TauDesc_[inp_idx].CollectionName);
0818   }
0819 }
0820 //define this as a plug-in
0821 #include "FWCore/Framework/interface/MakerMacros.h"
0822 DEFINE_FWK_MODULE(L2TauNNProducer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!