Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:27:44

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   tensorflow::setLogging("2");
0238 
0239   boost::property_tree::ptree loadPtreeRoot;
0240   auto const normalizationDict = edm::FileInPath(cfg.getParameter<std::string>("normalizationDict")).fullPath();
0241   boost::property_tree::read_json(normalizationDict, loadPtreeRoot);
0242   for (const auto& [key, val] : L2TauTagNNv1::varNameMap) {
0243     boost::property_tree::ptree var = loadPtreeRoot.get_child(val);
0244     normDictElement current_element;
0245     current_element.mean = var.get_child("mean").get_value<float>();
0246     current_element.std = var.get_child("std").get_value<float>();
0247     current_element.min = var.get_child("min").get_value<float>();
0248     current_element.max = var.get_child("max").get_value<float>();
0249     cacheData->normVec.push_back(current_element);
0250   }
0251   return cacheData;
0252 }
0253 void L2TauNNProducer::globalEndJob(L2TauNNProducerCacheData* cacheData) {
0254   if (cacheData->graphDef != nullptr) {
0255     delete cacheData->graphDef;
0256   }
0257   tensorflow::closeSession(cacheData->session);
0258 }
0259 void L2TauNNProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0260   edm::ParameterSetDescription desc;
0261   desc.add<int>("debugLevel", 0)->setComment("set debug level for printing out info");
0262   edm::ParameterSetDescription l1TausPset;
0263   l1TausPset.add<std::string>("L1CollectionName", "DoubleTau")->setComment("Name of collections");
0264   l1TausPset.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"))
0265       ->setComment("Which trigger should the L1 Taus collection pass");
0266   edm::ParameterSet l1TausPSetDefault;
0267   l1TausPSetDefault.addParameter<std::string>("L1CollectionName", "DoubleTau");
0268   l1TausPSetDefault.addParameter<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"));
0269   desc.addVPSet("L1Taus", l1TausPset, {l1TausPSetDefault});
0270   desc.add<edm::InputTag>("hbheInput", edm::InputTag("hltHbhereco"))->setComment("HBHE recHit collection");
0271   desc.add<edm::InputTag>("hoInput", edm::InputTag("hltHoreco"))->setComment("HO recHit Collection");
0272   desc.add<edm::InputTag>("ebInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEB"))->setComment("EB recHit Collection");
0273   desc.add<edm::InputTag>("eeInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEE"))->setComment("EE recHit Collection");
0274   desc.add<edm::InputTag>("pataVertices", edm::InputTag("hltPixelVerticesSoA"))
0275       ->setComment("patatrack vertices collection");
0276   desc.add<edm::InputTag>("pataTracks", edm::InputTag("hltPixelTracksSoA"))->setComment("patatrack collection");
0277   desc.add<edm::InputTag>("BeamSpot", edm::InputTag("hltOnlineBeamSpot"))->setComment("BeamSpot Collection");
0278   desc.add<uint>("maxVtx", 100)->setComment("max output collection size (number of accepted vertices)");
0279   desc.add<double>("fractionSumPt2", 0.3)->setComment("threshold on sumPt2 fraction of the leading vertex");
0280   desc.add<double>("minSumPt2", 0.)->setComment("min sumPt2");
0281   desc.add<double>("track_pt_min", 1.0)->setComment("min track p_T");
0282   desc.add<double>("track_pt_max", 10.0)->setComment("max track p_T");
0283   desc.add<double>("track_chi2_max", 99999.)->setComment("max track chi2");
0284   desc.add<std::string>("graphPath", "RecoTauTag/TrainingFiles/data/L2TauNNTag/L2TauTag_Run3v1.pb")
0285       ->setComment("path to the saved CNN");
0286   desc.add<std::string>("normalizationDict", "RecoTauTag/TrainingFiles/data/L2TauNNTag/NormalizationDict.json")
0287       ->setComment("path to the dictionary for variable standardization");
0288   descriptions.addWithDefaultLabel(desc);
0289 }
0290 
0291 L2TauNNProducer::L2TauNNProducer(const edm::ParameterSet& cfg, const L2TauNNProducerCacheData* cacheData)
0292     : debugLevel_(cfg.getParameter<int>("debugLevel")),
0293       hbheToken_(consumes<HBHERecHitCollection>(cfg.getParameter<edm::InputTag>("hbheInput"))),
0294       hoToken_(consumes<HORecHitCollection>(cfg.getParameter<edm::InputTag>("hoInput"))),
0295       ebToken_(consumes<EcalRecHitCollection>(cfg.getParameter<edm::InputTag>("ebInput"))),
0296       eeToken_(consumes<EcalRecHitCollection>(cfg.getParameter<edm::InputTag>("eeInput"))),
0297       geometryToken_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0298       bFieldToken_(esConsumes<MagneticField, IdealMagneticFieldRecord>()),
0299       pataVerticesToken_(consumes(cfg.getParameter<edm::InputTag>("pataVertices"))),
0300       pataTracksToken_(consumes(cfg.getParameter<edm::InputTag>("pataTracks"))),
0301       beamSpotToken_(consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>("BeamSpot"))),
0302       maxVtx_(cfg.getParameter<uint>("maxVtx")),
0303       fractionSumPt2_(cfg.getParameter<double>("fractionSumPt2")),
0304       minSumPt2_(cfg.getParameter<double>("minSumPt2")),
0305       trackPtMin_(cfg.getParameter<double>("track_pt_min")),
0306       trackPtMax_(cfg.getParameter<double>("track_pt_max")),
0307       trackChi2Max_(cfg.getParameter<double>("track_chi2_max")) {
0308   if (cacheData->graphDef == nullptr) {
0309     throw cms::Exception("InvalidCacheData") << "Invalid Cache Data.";
0310   }
0311   inputTensorName_ = cacheData->graphDef->node(0).name();
0312   outputTensorName_ = cacheData->graphDef->node(cacheData->graphDef->node_size() - 1).name();
0313   L2cacheData_ = cacheData;
0314   std::vector<edm::ParameterSet> L1TauCollections = cfg.getParameter<std::vector<edm::ParameterSet>>("L1Taus");
0315   L1TauDesc_.reserve(L1TauCollections.size());
0316   for (const auto& l1TauInput : L1TauCollections) {
0317     InputDescTau toInsert;
0318     toInsert.CollectionName = l1TauInput.getParameter<std::string>("L1CollectionName");
0319     toInsert.inputToken_ =
0320         consumes<trigger::TriggerFilterObjectWithRefs>(l1TauInput.getParameter<edm::InputTag>("L1TauTrigger"));
0321     L1TauDesc_.push_back(toInsert);
0322   }
0323   for (const auto& desc : L1TauDesc_)
0324     produces<std::vector<float>>(desc.CollectionName);
0325 }
0326 
0327 void L2TauNNProducer::checknan(tensorflow::Tensor& tensor, int debugLevel) {
0328   using NNInputs = L2TauTagNNv1::NNInputs;
0329   std::vector<int> tensor_shape(tensor.shape().dims());
0330   for (int d = 0; d < tensor.shape().dims(); d++) {
0331     tensor_shape.at(d) = tensor.shape().dim_size(d);
0332   }
0333   if (tensor_shape.size() != 4) {
0334     throw cms::Exception("InvalidTensor") << "Tensor shape does not have 4 dimensions!";
0335   }
0336   for (int tau_idx = 0; tau_idx < tensor_shape.at(0); tau_idx++) {
0337     for (int phi_idx = 0; phi_idx < tensor_shape.at(1); phi_idx++) {
0338       for (int eta_idx = 0; eta_idx < tensor_shape.at(2); eta_idx++) {
0339         for (int var_idx = 0; var_idx < tensor_shape.at(3); var_idx++) {
0340           auto getCell = [&](NNInputs input) -> float& {
0341             return getCellImpl(tensor, tau_idx, phi_idx, eta_idx, input);
0342           };
0343           auto nonstd_var = getCell(static_cast<NNInputs>(var_idx));
0344           if (edm::isNotFinite(nonstd_var)) {
0345             edm::LogWarning("InputVar") << "var is nan \nvar name= "
0346                                         << L2TauTagNNv1::varNameMap.at(static_cast<L2TauTagNNv1::NNInputs>(var_idx))
0347                                         << "\t var_idx = " << var_idx << "\t eta_idx = " << eta_idx
0348                                         << "\t phi_idx = " << phi_idx << "\t tau_idx = " << tau_idx;
0349             if (debugLevel > 2) {
0350               edm::LogWarning("InputVar") << "other vars in same cell \n";
0351               if (var_idx + 1 < tensor_shape.at(3))
0352                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 1))
0353                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 1));
0354               if (var_idx + 2 < tensor_shape.at(3))
0355                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 2))
0356                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 2));
0357               if (var_idx + 3 < tensor_shape.at(3))
0358                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 3))
0359                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 3));
0360               if (var_idx + 4 < tensor_shape.at(3))
0361                 edm::LogWarning("InputVar") << L2TauTagNNv1::varNameMap.at(static_cast<NNInputs>(var_idx + 4))
0362                                             << "\t = " << getCell(static_cast<NNInputs>(var_idx + 4));
0363             }
0364           }
0365         }
0366       }
0367     }
0368   }
0369 }
0370 
0371 void L2TauNNProducer::standardizeTensor(tensorflow::Tensor& tensor) {
0372   using NNInputs = L2TauTagNNv1::NNInputs;
0373   std::vector<int> tensor_shape(tensor.shape().dims());
0374   for (int d = 0; d < tensor.shape().dims(); d++) {
0375     tensor_shape.at(d) = tensor.shape().dim_size(d);
0376   }
0377   if (tensor_shape.size() != 4) {
0378     throw cms::Exception("InvalidTensor") << "Tensor shape does not have 4 dimensions!";
0379   }
0380   for (int tau_idx = 0; tau_idx < tensor_shape.at(0); tau_idx++) {
0381     for (int phi_idx = 0; phi_idx < tensor_shape.at(1); phi_idx++) {
0382       for (int eta_idx = 0; eta_idx < tensor_shape.at(2); eta_idx++) {
0383         for (int var_idx = 0; var_idx < tensor_shape.at(3); var_idx++) {
0384           auto getCell = [&](NNInputs input) -> float& {
0385             return getCellImpl(tensor, tau_idx, phi_idx, eta_idx, input);
0386           };
0387           float mean = L2cacheData_->normVec.at(var_idx).mean;
0388           float std = L2cacheData_->normVec.at(var_idx).std;
0389           float min = L2cacheData_->normVec.at(var_idx).min;
0390           float max = L2cacheData_->normVec.at(var_idx).max;
0391           float nonstd_var = getCell(static_cast<NNInputs>(var_idx));
0392           float std_var = static_cast<float>((nonstd_var - mean) / std);
0393           if (std_var > max) {
0394             std_var = static_cast<float>(max);
0395           } else if (std_var < min) {
0396             std_var = static_cast<float>(min);
0397           }
0398           getCell(static_cast<NNInputs>(var_idx)) = std_var;
0399         }
0400       }
0401     }
0402   }
0403 }
0404 
0405 void L2TauNNProducer::fillL1TauVars(tensorflow::Tensor& cellGridMatrix, const std::vector<l1t::TauRef>& allTaus) {
0406   using NNInputs = L2TauTagNNv1::NNInputs;
0407 
0408   const int nTaus = allTaus.size();
0409   for (int tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0410     for (int eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0411       for (int phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0412         auto getCell = [&](NNInputs input) -> float& {
0413           return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0414         };
0415         getCell(NNInputs::l1Tau_pt) = allTaus[tau_idx]->pt();
0416         getCell(NNInputs::l1Tau_eta) = allTaus[tau_idx]->eta();
0417         getCell(NNInputs::l1Tau_hwIso) = allTaus[tau_idx]->hwIso();
0418       }
0419     }
0420   }
0421 }
0422 
0423 template <typename LVec>
0424 std::tuple<float, float, int, int> L2TauNNProducer::getEtaPhiIndices(float eta, float phi, const LVec& tau_p4) {
0425   const float deta = eta - tau_p4.eta();
0426   const float dphi = reco::deltaPhi(phi, tau_p4.phi());
0427   const int eta_idx = static_cast<int>(floor((deta + L2TauTagNNv1::dR_max) / dEta_width));
0428   const int phi_idx = static_cast<int>(floor((dphi + L2TauTagNNv1::dR_max) / dPhi_width));
0429   return std::make_tuple(deta, dphi, eta_idx, phi_idx);
0430 }
0431 
0432 template <typename VPos, typename LVec>
0433 std::tuple<float, float, int, int> L2TauNNProducer::getEtaPhiIndices(const VPos& position, const LVec& tau_p4) {
0434   return getEtaPhiIndices(position.eta(), position.phi(), tau_p4);
0435 }
0436 
0437 void L2TauNNProducer::fillCaloRecHits(tensorflow::Tensor& cellGridMatrix,
0438                                       const std::vector<l1t::TauRef>& allTaus,
0439                                       const caloRecHitCollections& caloRecHits) {
0440   using NNInputs = L2TauTagNNv1::NNInputs;
0441 
0442   const int nTaus = allTaus.size();
0443   float deta, dphi;
0444   int eta_idx = 0;
0445   int phi_idx = 0;
0446   int tau_idx = 0;
0447 
0448   auto getCell = [&](NNInputs input) -> float& {
0449     return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0450   };
0451   for (tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0452     // calorechit_EE
0453     for (const auto& caloRecHit_ee : *caloRecHits.ee) {
0454       if (caloRecHit_ee.energy() <= 0)
0455         continue;
0456       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_ee.id())->getPosition();
0457       const float eeCalEn = caloRecHit_ee.energy();
0458       const float eeCalChi2 = caloRecHit_ee.chi2();
0459       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0460         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0461         getCell(NNInputs::EcalEnergySum) += eeCalEn;
0462         getCell(NNInputs::EcalSize) += 1.;
0463         getCell(NNInputs::EcalEnergyStdDev) += eeCalEn * eeCalEn;
0464         getCell(NNInputs::EcalDeltaEta) += deta * eeCalEn;
0465         getCell(NNInputs::EcalDeltaPhi) += dphi * eeCalEn;
0466         if (eeCalChi2 >= 0) {
0467           getCell(NNInputs::EcalChi2) += eeCalChi2 * eeCalEn;
0468           getCell(NNInputs::EcalEnergySumForPositiveChi2) += eeCalEn;
0469           getCell(NNInputs::EcalSizeForPositiveChi2) += 1.;
0470         }
0471       }
0472     }
0473 
0474     // calorechit_EB
0475     for (const auto& caloRecHit_eb : *caloRecHits.eb) {
0476       if (caloRecHit_eb.energy() <= 0)
0477         continue;
0478       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_eb.id())->getPosition();
0479       const float ebCalEn = caloRecHit_eb.energy();
0480       const float ebCalChi2 = caloRecHit_eb.chi2();
0481       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0482         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0483         getCell(NNInputs::EcalEnergySum) += ebCalEn;
0484         getCell(NNInputs::EcalSize) += 1.;
0485         getCell(NNInputs::EcalEnergyStdDev) += ebCalEn * ebCalEn;
0486         getCell(NNInputs::EcalDeltaEta) += deta * ebCalEn;
0487         getCell(NNInputs::EcalDeltaPhi) += dphi * ebCalEn;
0488         if (ebCalChi2 >= 0) {
0489           getCell(NNInputs::EcalChi2) += ebCalChi2 * ebCalEn;
0490           getCell(NNInputs::EcalEnergySumForPositiveChi2) += ebCalEn;
0491           getCell(NNInputs::EcalSizeForPositiveChi2) += 1.;
0492         }
0493       }
0494     }
0495 
0496     // calorechit_HBHE
0497     for (const auto& caloRecHit_hbhe : *caloRecHits.hbhe) {
0498       if (caloRecHit_hbhe.energy() <= 0)
0499         continue;
0500       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_hbhe.id())->getPosition();
0501       const float hbheCalEn = caloRecHit_hbhe.energy();
0502       const float hbheCalChi2 = caloRecHit_hbhe.chi2();
0503       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0504         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0505         getCell(NNInputs::HcalEnergySum) += hbheCalEn;
0506         getCell(NNInputs::HcalEnergyStdDev) += hbheCalEn * hbheCalEn;
0507         getCell(NNInputs::HcalSize) += 1.;
0508         getCell(NNInputs::HcalDeltaEta) += deta * hbheCalEn;
0509         getCell(NNInputs::HcalDeltaPhi) += dphi * hbheCalEn;
0510         if (hbheCalChi2 >= 0) {
0511           getCell(NNInputs::HcalChi2) += hbheCalChi2 * hbheCalEn;
0512           getCell(NNInputs::HcalEnergySumForPositiveChi2) += hbheCalEn;
0513           getCell(NNInputs::HcalSizeForPositiveChi2) += 1.;
0514         }
0515       }
0516     }
0517 
0518     // calorechit_HO
0519     for (const auto& caloRecHit_ho : *caloRecHits.ho) {
0520       if (caloRecHit_ho.energy() <= 0)
0521         continue;
0522       const auto& position = caloRecHits.geometry->getGeometry(caloRecHit_ho.id())->getPosition();
0523       const float hoCalEn = caloRecHit_ho.energy();
0524       if (reco::deltaR2(position, allTaus[tau_idx]->polarP4()) < dR2_max) {
0525         std::tie(deta, dphi, eta_idx, phi_idx) = getEtaPhiIndices(position, allTaus[tau_idx]->polarP4());
0526         getCell(NNInputs::HcalEnergySum) += hoCalEn;
0527         getCell(NNInputs::HcalEnergyStdDev) += hoCalEn * hoCalEn;
0528         getCell(NNInputs::HcalSize) += 1.;
0529         getCell(NNInputs::HcalDeltaEta) += deta * hoCalEn;
0530         getCell(NNInputs::HcalDeltaPhi) += dphi * hoCalEn;
0531       }
0532     }
0533 
0534     // normalize to sum and define stdDev
0535     for (eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0536       for (phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0537         /* normalize eCal vars*/
0538         if (getCell(NNInputs::EcalEnergySum) > 0.) {
0539           getCell(NNInputs::EcalDeltaEta) /= getCell(NNInputs::EcalEnergySum);
0540           getCell(NNInputs::EcalDeltaPhi) /= getCell(NNInputs::EcalEnergySum);
0541         }
0542         if (getCell(NNInputs::EcalEnergySumForPositiveChi2) > 0.) {
0543           getCell(NNInputs::EcalChi2) /= getCell(NNInputs::EcalEnergySumForPositiveChi2);
0544         }
0545         if (getCell(NNInputs::EcalSize) > 1.) {
0546           // (stdDev - (enSum*enSum)/size) / (size-1)
0547           getCell(NNInputs::EcalEnergyStdDev) =
0548               (getCell(NNInputs::EcalEnergyStdDev) -
0549                (getCell(NNInputs::EcalEnergySum) * getCell(NNInputs::EcalEnergySum)) / getCell(NNInputs::EcalSize)) /
0550               (getCell(NNInputs::EcalSize) - 1);
0551         } else {
0552           getCell(NNInputs::EcalEnergyStdDev) = 0.;
0553         }
0554         /* normalize hCal Vars */
0555         if (getCell(NNInputs::HcalEnergySum) > 0.) {
0556           getCell(NNInputs::HcalDeltaEta) /= getCell(NNInputs::HcalEnergySum);
0557           getCell(NNInputs::HcalDeltaPhi) /= getCell(NNInputs::HcalEnergySum);
0558         }
0559         if (getCell(NNInputs::HcalEnergySumForPositiveChi2) > 0.) {
0560           getCell(NNInputs::HcalChi2) /= getCell(NNInputs::HcalEnergySumForPositiveChi2);
0561         }
0562         if (getCell(NNInputs::HcalSize) > 1.) {
0563           // (stdDev - (enSum*enSum)/size) / (size-1)
0564           getCell(NNInputs::HcalEnergyStdDev) =
0565               (getCell(NNInputs::HcalEnergyStdDev) -
0566                (getCell(NNInputs::HcalEnergySum) * getCell(NNInputs::HcalEnergySum)) / getCell(NNInputs::HcalSize)) /
0567               (getCell(NNInputs::HcalSize) - 1);
0568         } else {
0569           getCell(NNInputs::HcalEnergyStdDev) = 0.;
0570         }
0571       }
0572     }
0573   }
0574 }
0575 
0576 void L2TauNNProducer::selectGoodTracksAndVertices(const ZVertexSoAHost& patavtx_soa,
0577                                                   const TrackSoAHost& patatracks_tsoa,
0578                                                   std::vector<int>& trkGood,
0579                                                   std::vector<int>& vtxGood) {
0580   using patatrackHelpers = TracksUtilities<pixelTopology::Phase1>;
0581   const auto maxTracks = patatracks_tsoa.view().metadata().size();
0582   const int nv = patavtx_soa.view().nvFinal();
0583   trkGood.clear();
0584   trkGood.reserve(maxTracks);
0585   vtxGood.clear();
0586   vtxGood.reserve(nv);
0587   auto const* quality = patatracks_tsoa.view().quality();
0588 
0589   // No need to sort either as the algorithms is just using the max (not even the location, just the max value of pt2sum).
0590   std::vector<float> pTSquaredSum(nv, 0);
0591   std::vector<int> nTrkAssociated(nv, 0);
0592 
0593   for (int32_t trk_idx = 0; trk_idx < maxTracks; ++trk_idx) {
0594     auto nHits = patatrackHelpers::nHits(patatracks_tsoa.view(), trk_idx);
0595     if (nHits == 0) {
0596       break;
0597     }
0598     int vtx_ass_to_track = patavtx_soa.view()[trk_idx].idv();
0599     if (vtx_ass_to_track >= 0 && vtx_ass_to_track < nv) {
0600       auto patatrackPt = patatracks_tsoa.view()[trk_idx].pt();
0601       ++nTrkAssociated[vtx_ass_to_track];
0602       if (patatrackPt >= trackPtMin_ && patatracks_tsoa.const_view()[trk_idx].chi2() <= trackChi2Max_) {
0603         patatrackPt = std::min(patatrackPt, trackPtMax_);
0604         pTSquaredSum[vtx_ass_to_track] += patatrackPt * patatrackPt;
0605       }
0606     }
0607     if (nHits > 0 and quality[trk_idx] >= pixelTrack::Quality::loose) {
0608       trkGood.push_back(trk_idx);
0609     }
0610   }
0611   if (nv > 0) {
0612     const auto minFOM_fromFrac = (*std::max_element(pTSquaredSum.begin(), pTSquaredSum.end())) * fractionSumPt2_;
0613     for (int j = nv - 1; j >= 0 && vtxGood.size() < maxVtx_; --j) {
0614       auto vtx_idx = patavtx_soa.view()[j].sortInd();
0615       assert(vtx_idx < nv);
0616       if (nTrkAssociated[vtx_idx] >= 2 && pTSquaredSum[vtx_idx] >= minFOM_fromFrac &&
0617           pTSquaredSum[vtx_idx] > minSumPt2_) {
0618         vtxGood.push_back(vtx_idx);
0619       }
0620     }
0621   }
0622 }
0623 
0624 std::pair<float, float> L2TauNNProducer::impactParameter(int it,
0625                                                          const TrackSoAHost& patatracks_tsoa,
0626                                                          float patatrackPhi,
0627                                                          const reco::BeamSpot& beamspot,
0628                                                          const MagneticField* magfi) {
0629   /* dxy and dz */
0630   riemannFit::Vector5d ipar, opar;
0631   riemannFit::Matrix5d icov, ocov;
0632   TracksUtilities<pixelTopology::Phase1>::copyToDense(patatracks_tsoa.view(), ipar, icov, it);
0633   riemannFit::transformToPerigeePlane(ipar, icov, opar, ocov);
0634   LocalTrajectoryParameters lpar(opar(0), opar(1), opar(2), opar(3), opar(4), 1.);
0635   float sp = std::sin(patatrackPhi);
0636   float cp = std::cos(patatrackPhi);
0637   Surface::RotationType Rotation(sp, -cp, 0, 0, 0, -1.f, cp, sp, 0);
0638   GlobalPoint BeamSpotPoint(beamspot.x0(), beamspot.y0(), beamspot.z0());
0639   Plane impPointPlane(BeamSpotPoint, Rotation);
0640   GlobalTrajectoryParameters gp(
0641       impPointPlane.toGlobal(lpar.position()), impPointPlane.toGlobal(lpar.momentum()), lpar.charge(), magfi);
0642   GlobalPoint vv = gp.position();
0643   math::XYZPoint pos(vv.x(), vv.y(), vv.z());
0644   GlobalVector pp = gp.momentum();
0645   math::XYZVector mom(pp.x(), pp.y(), pp.z());
0646   auto lambda = M_PI_2 - pp.theta();
0647   auto phi = pp.phi();
0648   float patatrackDxy = -vv.x() * std::sin(phi) + vv.y() * std::cos(phi);
0649   float patatrackDz =
0650       (vv.z() * std::cos(lambda) - (vv.x() * std::cos(phi) + vv.y() * std::sin(phi)) * std::sin(lambda)) /
0651       std::cos(lambda);
0652   return std::make_pair(patatrackDxy, patatrackDz);
0653 }
0654 
0655 void L2TauNNProducer::fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0656                                      const std::vector<l1t::TauRef>& allTaus,
0657                                      const TrackSoAHost& patatracks_tsoa,
0658                                      const ZVertexSoAHost& patavtx_soa,
0659                                      const reco::BeamSpot& beamspot,
0660                                      const MagneticField* magfi) {
0661   using NNInputs = L2TauTagNNv1::NNInputs;
0662   using patatrackHelpers = TracksUtilities<pixelTopology::Phase1>;
0663   float deta, dphi;
0664   int eta_idx = 0;
0665   int phi_idx = 0;
0666   int tau_idx = 0;
0667 
0668   auto getCell = [&](NNInputs input) -> float& {
0669     return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0670   };
0671 
0672   std::vector<int> trkGood;
0673   std::vector<int> vtxGood;
0674 
0675   selectGoodTracksAndVertices(patavtx_soa, patatracks_tsoa, trkGood, vtxGood);
0676 
0677   const int nTaus = allTaus.size();
0678   for (tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0679     const float tauEta = allTaus[tau_idx]->eta();
0680     const float tauPhi = allTaus[tau_idx]->phi();
0681 
0682     for (const auto it : trkGood) {
0683       const float patatrackPt = patatracks_tsoa.const_view()[it].pt();
0684       if (patatrackPt <= 0)
0685         continue;
0686       const float patatrackPhi = patatrackHelpers::phi(patatracks_tsoa.const_view(), it);
0687       const float patatrackEta = patatracks_tsoa.const_view()[it].eta();
0688       const float patatrackCharge = patatrackHelpers::charge(patatracks_tsoa.const_view(), it);
0689       const float patatrackChi2OverNdof = patatracks_tsoa.view()[it].chi2();
0690       const auto nHits = patatrackHelpers::nHits(patatracks_tsoa.const_view(), it);
0691       if (nHits <= 0)
0692         continue;
0693       const int patatrackNdof = 2 * std::min(6, nHits) - 5;
0694 
0695       const int vtx_idx_assTrk = patavtx_soa.view()[it].idv();
0696       if (reco::deltaR2(patatrackEta, patatrackPhi, tauEta, tauPhi) < dR2_max) {
0697         std::tie(deta, dphi, eta_idx, phi_idx) =
0698             getEtaPhiIndices(patatrackEta, patatrackPhi, allTaus[tau_idx]->polarP4());
0699         getCell(NNInputs::PatatrackPtSum) += patatrackPt;
0700         getCell(NNInputs::PatatrackSize) += 1.;
0701         getCell(NNInputs::PatatrackChargeSum) += patatrackCharge;
0702         getCell(NNInputs::PatatrackDeltaEta) += deta * patatrackPt;
0703         getCell(NNInputs::PatatrackDeltaPhi) += dphi * patatrackPt;
0704         getCell(NNInputs::PatatrackChi2OverNdof) += patatrackChi2OverNdof * patatrackPt;
0705         getCell(NNInputs::PatatrackNdof) += patatrackNdof * patatrackPt;
0706         std::pair<float, float> impactParameters = impactParameter(it, patatracks_tsoa, patatrackPhi, beamspot, magfi);
0707         getCell(NNInputs::PatatrackDxy) += impactParameters.first * patatrackPt;
0708         getCell(NNInputs::PatatrackDz) += impactParameters.second * patatrackPt;
0709         if ((std::find(vtxGood.begin(), vtxGood.end(), vtx_idx_assTrk) != vtxGood.end())) {
0710           getCell(NNInputs::PatatrackPtSumWithVertex) += patatrackPt;
0711           getCell(NNInputs::PatatrackSizeWithVertex) += 1.;
0712         }
0713       }
0714     }
0715 
0716     // normalize to sum and define stdDev
0717     for (eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0718       for (phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0719         getCell(NNInputs::nVertices) = vtxGood.size();
0720         if (getCell(NNInputs::PatatrackPtSum) > 0.) {
0721           getCell(NNInputs::PatatrackDeltaEta) /= getCell(NNInputs::PatatrackPtSum);
0722           getCell(NNInputs::PatatrackDeltaPhi) /= getCell(NNInputs::PatatrackPtSum);
0723           getCell(NNInputs::PatatrackChi2OverNdof) /= getCell(NNInputs::PatatrackPtSum);
0724           getCell(NNInputs::PatatrackNdof) /= getCell(NNInputs::PatatrackPtSum);
0725           getCell(NNInputs::PatatrackDxy) /= getCell(NNInputs::PatatrackPtSum);
0726           getCell(NNInputs::PatatrackDz) /= getCell(NNInputs::PatatrackPtSum);
0727         }
0728       }
0729     }
0730   }
0731 }
0732 
0733 std::vector<float> L2TauNNProducer::getTauScore(const tensorflow::Tensor& cellGridMatrix) {
0734   const int nTau = cellGridMatrix.shape().dim_size(0);
0735   std::vector<float> pred_vector(nTau);
0736   if (nTau > 0) {
0737     // Only run the inference if there are taus to process
0738     std::vector<tensorflow::Tensor> pred_tensor;
0739     tensorflow::run(L2cacheData_->session, {{inputTensorName_, cellGridMatrix}}, {outputTensorName_}, &pred_tensor);
0740     for (int tau_idx = 0; tau_idx < nTau; ++tau_idx) {
0741       pred_vector[tau_idx] = pred_tensor[0].matrix<float>()(tau_idx, 0);
0742     }
0743   }
0744   return pred_vector;
0745 }
0746 
0747 void L2TauNNProducer::produce(edm::Event& event, const edm::EventSetup& eventsetup) {
0748   std::vector<std::vector<size_t>> TauCollectionMap(L1TauDesc_.size());
0749   l1t::TauVectorRef allTaus;
0750 
0751   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0752     l1t::TauVectorRef l1Taus;
0753     auto const& l1TriggeredTaus = event.get(L1TauDesc_[inp_idx].inputToken_);
0754     l1TriggeredTaus.getObjects(trigger::TriggerL1Tau, l1Taus);
0755     TauCollectionMap.at(inp_idx).resize(l1Taus.size());
0756 
0757     for (size_t l1_idx = 0; l1_idx < l1Taus.size(); l1_idx++) {
0758       size_t tau_idx;
0759       const auto iter = std::find(allTaus.begin(), allTaus.end(), l1Taus[l1_idx]);
0760       if (iter != allTaus.end()) {
0761         tau_idx = std::distance(allTaus.begin(), iter);
0762       } else {
0763         allTaus.push_back(l1Taus[l1_idx]);
0764         tau_idx = allTaus.size() - 1;
0765       }
0766       TauCollectionMap.at(inp_idx).at(l1_idx) = tau_idx;
0767     }
0768   }
0769   const auto ebCal = event.getHandle(ebToken_);
0770   const auto eeCal = event.getHandle(eeToken_);
0771   const auto hbhe = event.getHandle(hbheToken_);
0772   const auto ho = event.getHandle(hoToken_);
0773   auto const& patatracks_SoA = event.get(pataTracksToken_);
0774   auto const& vertices_SoA = event.get(pataVerticesToken_);
0775   const auto bsHandle = event.getHandle(beamSpotToken_);
0776 
0777   auto const fieldESH = eventsetup.getHandle(bFieldToken_);
0778   auto const geometry = eventsetup.getHandle(geometryToken_);
0779 
0780   caloRecHitCollections caloRecHits;
0781   caloRecHits.hbhe = &*hbhe;
0782   caloRecHits.ho = &*ho;
0783   caloRecHits.eb = &*ebCal;
0784   caloRecHits.ee = &*eeCal;
0785   caloRecHits.geometry = &*geometry;
0786 
0787   const int nTaus = allTaus.size();
0788   tensorflow::Tensor cellGridMatrix(tensorflow::DT_FLOAT,
0789                                     {nTaus, L2TauTagNNv1::nCellEta, L2TauTagNNv1::nCellPhi, L2TauTagNNv1::nVars});
0790   const int n_inputs = nTaus * L2TauTagNNv1::nCellEta * L2TauTagNNv1::nCellPhi * L2TauTagNNv1::nVars;
0791   for (int input_idx = 0; input_idx < n_inputs; ++input_idx) {
0792     cellGridMatrix.flat<float>()(input_idx) = 0;
0793   }
0794   fillL1TauVars(cellGridMatrix, allTaus);
0795 
0796   fillCaloRecHits(cellGridMatrix, allTaus, caloRecHits);
0797 
0798   fillPatatracks(cellGridMatrix, allTaus, patatracks_SoA, vertices_SoA, *bsHandle, fieldESH.product());
0799 
0800   standardizeTensor(cellGridMatrix);
0801 
0802   if (debugLevel_ > 0) {
0803     checknan(cellGridMatrix, debugLevel_);
0804   }
0805 
0806   std::vector<float> tau_score = getTauScore(cellGridMatrix);
0807 
0808   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0809     const size_t nTau = TauCollectionMap[inp_idx].size();
0810     auto tau_tags = std::make_unique<std::vector<float>>(nTau);
0811     for (size_t tau_pos = 0; tau_pos < nTau; ++tau_pos) {
0812       const auto tau_idx = TauCollectionMap[inp_idx][tau_pos];
0813       if (debugLevel_ > 0) {
0814         edm::LogInfo("DebugInfo") << event.id().event() << " \t " << (allTaus[tau_idx])->pt() << " \t "
0815                                   << tau_score.at(tau_idx) << std::endl;
0816       }
0817       (*tau_tags)[tau_pos] = tau_score.at(tau_idx);
0818     }
0819     event.put(std::move(tau_tags), L1TauDesc_[inp_idx].CollectionName);
0820   }
0821 }
0822 //define this as a plug-in
0823 #include "FWCore/Framework/interface/MakerMacros.h"
0824 DEFINE_FWK_MODULE(L2TauNNProducer);