Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-06-07 02:29:54

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 "DataFormats/SiPixelClusterSoA/interface/ClusteringConstants.h"
0050 
0051 #include "DataFormats/TrackSoA/interface/alpaka/TrackUtilities.h"
0052 #include "DataFormats/TrackSoA/interface/TracksHost.h"
0053 #include "DataFormats/VertexSoA/interface/ZVertexHost.h"
0054 
0055 namespace L2TauTagNNv1 {
0056   constexpr int nCellEta = 5;
0057   constexpr int nCellPhi = 5;
0058   constexpr int nVars = 31;
0059   constexpr float dR_max = 0.5;
0060   enum class NNInputs {
0061     nVertices = 0,
0062     l1Tau_pt,
0063     l1Tau_eta,
0064     l1Tau_hwIso,
0065     EcalEnergySum,
0066     EcalSize,
0067     EcalEnergyStdDev,
0068     EcalDeltaEta,
0069     EcalDeltaPhi,
0070     EcalChi2,
0071     EcalEnergySumForPositiveChi2,
0072     EcalSizeForPositiveChi2,
0073     HcalEnergySum,
0074     HcalSize,
0075     HcalEnergyStdDev,
0076     HcalDeltaEta,
0077     HcalDeltaPhi,
0078     HcalChi2,
0079     HcalEnergySumForPositiveChi2,
0080     HcalSizeForPositiveChi2,
0081     PatatrackPtSum,
0082     PatatrackSize,
0083     PatatrackSizeWithVertex,
0084     PatatrackPtSumWithVertex,
0085     PatatrackChargeSum,
0086     PatatrackDeltaEta,
0087     PatatrackDeltaPhi,
0088     PatatrackChi2OverNdof,
0089     PatatrackNdof,
0090     PatatrackDxy,
0091     PatatrackDz
0092   };
0093 
0094   const std::map<NNInputs, std::string> varNameMap = {
0095       {NNInputs::nVertices, "nVertices"},
0096       {NNInputs::l1Tau_pt, "l1Tau_pt"},
0097       {NNInputs::l1Tau_eta, "l1Tau_eta"},
0098       {NNInputs::l1Tau_hwIso, "l1Tau_hwIso"},
0099       {NNInputs::EcalEnergySum, "EcalEnergySum"},
0100       {NNInputs::EcalSize, "EcalSize"},
0101       {NNInputs::EcalEnergyStdDev, "EcalEnergyStdDev"},
0102       {NNInputs::EcalDeltaEta, "EcalDeltaEta"},
0103       {NNInputs::EcalDeltaPhi, "EcalDeltaPhi"},
0104       {NNInputs::EcalChi2, "EcalChi2"},
0105       {NNInputs::EcalEnergySumForPositiveChi2, "EcalEnergySumForPositiveChi2"},
0106       {NNInputs::EcalSizeForPositiveChi2, "EcalSizeForPositiveChi2"},
0107       {NNInputs::HcalEnergySum, "HcalEnergySum"},
0108       {NNInputs::HcalSize, "HcalSize"},
0109       {NNInputs::HcalEnergyStdDev, "HcalEnergyStdDev"},
0110       {NNInputs::HcalDeltaEta, "HcalDeltaEta"},
0111       {NNInputs::HcalDeltaPhi, "HcalDeltaPhi"},
0112       {NNInputs::HcalChi2, "HcalChi2"},
0113       {NNInputs::HcalEnergySumForPositiveChi2, "HcalEnergySumForPositiveChi2"},
0114       {NNInputs::HcalSizeForPositiveChi2, "HcalSizeForPositiveChi2"},
0115       {NNInputs::PatatrackPtSum, "PatatrackPtSum"},
0116       {NNInputs::PatatrackSize, "PatatrackSize"},
0117       {NNInputs::PatatrackSizeWithVertex, "PatatrackSizeWithVertex"},
0118       {NNInputs::PatatrackPtSumWithVertex, "PatatrackPtSumWithVertex"},
0119       {NNInputs::PatatrackChargeSum, "PatatrackChargeSum"},
0120       {NNInputs::PatatrackDeltaEta, "PatatrackDeltaEta"},
0121       {NNInputs::PatatrackDeltaPhi, "PatatrackDeltaPhi"},
0122       {NNInputs::PatatrackChi2OverNdof, "PatatrackChi2OverNdof"},
0123       {NNInputs::PatatrackNdof, "PatatrackNdof"},
0124       {NNInputs::PatatrackDxy, "PatatrackDxy"},
0125       {NNInputs::PatatrackDz, "PatatrackDz"}};
0126 }  // namespace L2TauTagNNv1
0127 namespace {
0128   inline float& getCellImpl(
0129       tensorflow::Tensor& cellGridMatrix, int tau_idx, int phi_idx, int eta_idx, L2TauTagNNv1::NNInputs NNInput_idx) {
0130     return cellGridMatrix.tensor<float, 4>()(tau_idx, phi_idx, eta_idx, static_cast<int>(NNInput_idx));
0131   }
0132 }  // namespace
0133 struct normDictElement {
0134   float mean;
0135   float std;
0136   float min;
0137   float max;
0138 };
0139 
0140 struct L2TauNNProducerAlpakaCacheData {
0141   L2TauNNProducerAlpakaCacheData() : graphDef(nullptr), session(nullptr) {}
0142   tensorflow::GraphDef* graphDef;
0143   tensorflow::Session* session;
0144   std::vector<normDictElement> normVec;
0145 };
0146 
0147 class L2TauNNProducerAlpaka : public edm::stream::EDProducer<edm::GlobalCache<L2TauNNProducerAlpakaCacheData>> {
0148 public:
0149   using TracksHost = pixelTrack::TracksHostPhase1;
0150 
0151   struct caloRecHitCollections {
0152     const HBHERecHitCollection* hbhe;
0153     const HORecHitCollection* ho;
0154     const EcalRecHitCollection* eb;
0155     const EcalRecHitCollection* ee;
0156     const CaloGeometry* geometry;
0157   };
0158 
0159   struct InputDescTau {
0160     std::string CollectionName;
0161     edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> inputToken_;
0162   };
0163 
0164   static constexpr float dR2_max = L2TauTagNNv1::dR_max * L2TauTagNNv1::dR_max;
0165   static constexpr float dEta_width = 2 * L2TauTagNNv1::dR_max / static_cast<float>(L2TauTagNNv1::nCellEta);
0166   static constexpr float dPhi_width = 2 * L2TauTagNNv1::dR_max / static_cast<float>(L2TauTagNNv1::nCellPhi);
0167 
0168   explicit L2TauNNProducerAlpaka(const edm::ParameterSet&, const L2TauNNProducerAlpakaCacheData*);
0169   static void fillDescriptions(edm::ConfigurationDescriptions&);
0170   static std::unique_ptr<L2TauNNProducerAlpakaCacheData> initializeGlobalCache(const edm::ParameterSet&);
0171   static void globalEndJob(L2TauNNProducerAlpakaCacheData*);
0172 
0173 private:
0174   void checknan(tensorflow::Tensor& tensor, int debugLevel);
0175   void standardizeTensor(tensorflow::Tensor& tensor);
0176   std::vector<float> getTauScore(const tensorflow::Tensor& cellGridMatrix);
0177   void produce(edm::Event& event, const edm::EventSetup& eventsetup) override;
0178   void fillL1TauVars(tensorflow::Tensor& cellGridMatrix, const std::vector<l1t::TauRef>& allTaus);
0179   void fillCaloRecHits(tensorflow::Tensor& cellGridMatrix,
0180                        const std::vector<l1t::TauRef>& allTaus,
0181                        const caloRecHitCollections& caloRecHits);
0182   void fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0183                       const std::vector<l1t::TauRef>& allTaus,
0184                       const TracksHost& patatracks_tsoa,
0185                       const ZVertexHost& patavtx_soa,
0186                       const reco::BeamSpot& beamspot,
0187                       const MagneticField* magfi);
0188   void selectGoodTracksAndVertices(const ZVertexHost& patavtx_soa,
0189                                    const TracksHost& patatracks_tsoa,
0190                                    std::vector<int>& trkGood,
0191                                    std::vector<int>& vtxGood);
0192   std::pair<float, float> impactParameter(int it,
0193                                           const TracksHost& patatracks_tsoa,
0194                                           float patatrackPhi,
0195                                           const reco::BeamSpot& beamspot,
0196                                           const MagneticField* magfi);
0197   template <typename VPos, typename LVec>
0198   std::tuple<float, float, int, int> getEtaPhiIndices(const VPos& position, const LVec& tau_p4);
0199   template <typename LVec>
0200   std::tuple<float, float, int, int> getEtaPhiIndices(float eta, float phi, const LVec& tau_p4);
0201 
0202 private:
0203   const int debugLevel_;
0204   const edm::EDGetTokenT<trigger::TriggerFilterObjectWithRefs> tauTriggerToken_;
0205   std::vector<InputDescTau> L1TauDesc_;
0206   const edm::EDGetTokenT<HBHERecHitCollection> hbheToken_;
0207   const edm::EDGetTokenT<HORecHitCollection> hoToken_;
0208   const edm::EDGetTokenT<EcalRecHitCollection> ebToken_;
0209   const edm::EDGetTokenT<EcalRecHitCollection> eeToken_;
0210   const edm::ESGetToken<CaloGeometry, CaloGeometryRecord> geometryToken_;
0211   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> bFieldToken_;
0212   const edm::EDGetTokenT<ZVertexHost> pataVerticesToken_;
0213   const edm::EDGetTokenT<TracksHost> pataTracksToken_;
0214   const edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;
0215   const unsigned int maxVtx_;
0216   const float fractionSumPt2_;
0217   const float minSumPt2_;
0218   const float trackPtMin_;
0219   const float trackPtMax_;
0220   const float trackChi2Max_;
0221   std::string inputTensorName_;
0222   std::string outputTensorName_;
0223   const L2TauNNProducerAlpakaCacheData* L2cacheData_;
0224 };
0225 
0226 std::unique_ptr<L2TauNNProducerAlpakaCacheData> L2TauNNProducerAlpaka::initializeGlobalCache(
0227     const edm::ParameterSet& cfg) {
0228   std::unique_ptr<L2TauNNProducerAlpakaCacheData> cacheData = std::make_unique<L2TauNNProducerAlpakaCacheData>();
0229   cacheData->normVec.reserve(L2TauTagNNv1::nVars);
0230 
0231   auto const graphPath = edm::FileInPath(cfg.getParameter<std::string>("graphPath")).fullPath();
0232 
0233   cacheData->graphDef = tensorflow::loadGraphDef(graphPath);
0234   cacheData->session = tensorflow::createSession(cacheData->graphDef);
0235 
0236   tensorflow::setLogging("2");
0237 
0238   boost::property_tree::ptree loadPtreeRoot;
0239   auto const normalizationDict = edm::FileInPath(cfg.getParameter<std::string>("normalizationDict")).fullPath();
0240   boost::property_tree::read_json(normalizationDict, loadPtreeRoot);
0241   for (const auto& [key, val] : L2TauTagNNv1::varNameMap) {
0242     boost::property_tree::ptree var = loadPtreeRoot.get_child(val);
0243     normDictElement current_element;
0244     current_element.mean = var.get_child("mean").get_value<float>();
0245     current_element.std = var.get_child("std").get_value<float>();
0246     current_element.min = var.get_child("min").get_value<float>();
0247     current_element.max = var.get_child("max").get_value<float>();
0248     cacheData->normVec.push_back(current_element);
0249   }
0250   return cacheData;
0251 }
0252 void L2TauNNProducerAlpaka::globalEndJob(L2TauNNProducerAlpakaCacheData* cacheData) {
0253   if (cacheData->graphDef != nullptr) {
0254     delete cacheData->graphDef;
0255   }
0256   tensorflow::closeSession(cacheData->session);
0257 }
0258 void L2TauNNProducerAlpaka::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0259   edm::ParameterSetDescription desc;
0260   desc.add<int>("debugLevel", 0)->setComment("set debug level for printing out info");
0261   edm::ParameterSetDescription l1TausPset;
0262   l1TausPset.add<std::string>("L1CollectionName", "DoubleTau")->setComment("Name of collections");
0263   l1TausPset.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"))
0264       ->setComment("Which trigger should the L1 Taus collection pass");
0265   edm::ParameterSet l1TausPSetDefault;
0266   l1TausPSetDefault.addParameter<std::string>("L1CollectionName", "DoubleTau");
0267   l1TausPSetDefault.addParameter<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"));
0268   desc.addVPSet("L1Taus", l1TausPset, {l1TausPSetDefault});
0269   desc.add<edm::InputTag>("hbheInput", edm::InputTag("hltHbhereco"))->setComment("HBHE recHit collection");
0270   desc.add<edm::InputTag>("hoInput", edm::InputTag("hltHoreco"))->setComment("HO recHit Collection");
0271   desc.add<edm::InputTag>("ebInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEB"))->setComment("EB recHit Collection");
0272   desc.add<edm::InputTag>("eeInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEE"))->setComment("EE recHit Collection");
0273   desc.add<edm::InputTag>("pataVertices", edm::InputTag("hltPixelVerticesSoA"))
0274       ->setComment("patatrack vertices collection");
0275   desc.add<edm::InputTag>("pataTracks", edm::InputTag("hltPixelTracksSoA"))->setComment("patatrack collection");
0276   desc.add<edm::InputTag>("BeamSpot", edm::InputTag("hltOnlineBeamSpot"))->setComment("BeamSpot Collection");
0277   desc.add<uint>("maxVtx", 100)->setComment("max output collection size (number of accepted vertices)");
0278   desc.add<double>("fractionSumPt2", 0.3)->setComment("threshold on sumPt2 fraction of the leading vertex");
0279   desc.add<double>("minSumPt2", 0.)->setComment("min sumPt2");
0280   desc.add<double>("track_pt_min", 1.0)->setComment("min track p_T");
0281   desc.add<double>("track_pt_max", 10.0)->setComment("max track p_T");
0282   desc.add<double>("track_chi2_max", 99999.)->setComment("max track chi2");
0283   desc.add<std::string>("graphPath", "RecoTauTag/TrainingFiles/data/L2TauNNTag/L2TauTag_Run3v1.pb")
0284       ->setComment("path to the saved CNN");
0285   desc.add<std::string>("normalizationDict", "RecoTauTag/TrainingFiles/data/L2TauNNTag/NormalizationDict.json")
0286       ->setComment("path to the dictionary for variable standardization");
0287   descriptions.addWithDefaultLabel(desc);
0288 }
0289 
0290 L2TauNNProducerAlpaka::L2TauNNProducerAlpaka(const edm::ParameterSet& cfg,
0291                                              const L2TauNNProducerAlpakaCacheData* 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 L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::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> L2TauNNProducerAlpaka::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> L2TauNNProducerAlpaka::getEtaPhiIndices(const VPos& position, const LVec& tau_p4) {
0434   return getEtaPhiIndices(position.eta(), position.phi(), tau_p4);
0435 }
0436 
0437 void L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::selectGoodTracksAndVertices(const ZVertexHost& patavtx_soa,
0577                                                         const TracksHost& 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> L2TauNNProducerAlpaka::impactParameter(int it,
0625                                                                const TracksHost& 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 L2TauNNProducerAlpaka::fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0656                                            const std::vector<l1t::TauRef>& allTaus,
0657                                            const TracksHost& patatracks_tsoa,
0658                                            const ZVertexHost& 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 = reco::phi(patatracks_tsoa.const_view(), it);
0687       const float patatrackEta = patatracks_tsoa.const_view()[it].eta();
0688       const float patatrackCharge = reco::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> L2TauNNProducerAlpaka::getTauScore(const tensorflow::Tensor& cellGridMatrix) {
0734   const int nTau = cellGridMatrix.shape().dim_size(0);
0735   if (nTau == 0) {
0736     return std::vector<float>();
0737   } else {
0738     std::vector<tensorflow::Tensor> pred_tensor;
0739     tensorflow::run(L2cacheData_->session, {{inputTensorName_, cellGridMatrix}}, {outputTensorName_}, &pred_tensor);
0740     std::vector<float> pred_vector(nTau);
0741     for (int tau_idx = 0; tau_idx < nTau; ++tau_idx) {
0742       pred_vector[tau_idx] = pred_tensor[0].matrix<float>()(tau_idx, 0);
0743     }
0744 
0745     return pred_vector;
0746   }
0747 }
0748 
0749 void L2TauNNProducerAlpaka::produce(edm::Event& event, const edm::EventSetup& eventsetup) {
0750   std::vector<std::vector<size_t>> TauCollectionMap(L1TauDesc_.size());
0751   l1t::TauVectorRef allTaus;
0752 
0753   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0754     l1t::TauVectorRef l1Taus;
0755     auto const& l1TriggeredTaus = event.get(L1TauDesc_[inp_idx].inputToken_);
0756     l1TriggeredTaus.getObjects(trigger::TriggerL1Tau, l1Taus);
0757     TauCollectionMap.at(inp_idx).resize(l1Taus.size());
0758 
0759     for (size_t l1_idx = 0; l1_idx < l1Taus.size(); l1_idx++) {
0760       size_t tau_idx;
0761       const auto iter = std::find(allTaus.begin(), allTaus.end(), l1Taus[l1_idx]);
0762       if (iter != allTaus.end()) {
0763         tau_idx = std::distance(allTaus.begin(), iter);
0764       } else {
0765         allTaus.push_back(l1Taus[l1_idx]);
0766         tau_idx = allTaus.size() - 1;
0767       }
0768       TauCollectionMap.at(inp_idx).at(l1_idx) = tau_idx;
0769     }
0770   }
0771   const auto ebCal = event.getHandle(ebToken_);
0772   const auto eeCal = event.getHandle(eeToken_);
0773   const auto hbhe = event.getHandle(hbheToken_);
0774   const auto ho = event.getHandle(hoToken_);
0775   auto const& patatracks_SoA = event.get(pataTracksToken_);
0776   auto const& vertices_SoA = event.get(pataVerticesToken_);
0777   const auto bsHandle = event.getHandle(beamSpotToken_);
0778 
0779   auto const fieldESH = eventsetup.getHandle(bFieldToken_);
0780   auto const geometry = eventsetup.getHandle(geometryToken_);
0781 
0782   caloRecHitCollections caloRecHits;
0783   caloRecHits.hbhe = &*hbhe;
0784   caloRecHits.ho = &*ho;
0785   caloRecHits.eb = &*ebCal;
0786   caloRecHits.ee = &*eeCal;
0787   caloRecHits.geometry = &*geometry;
0788 
0789   const int nTaus = allTaus.size();
0790   tensorflow::Tensor cellGridMatrix(tensorflow::DT_FLOAT,
0791                                     {nTaus, L2TauTagNNv1::nCellEta, L2TauTagNNv1::nCellPhi, L2TauTagNNv1::nVars});
0792   const int n_inputs = nTaus * L2TauTagNNv1::nCellEta * L2TauTagNNv1::nCellPhi * L2TauTagNNv1::nVars;
0793   for (int input_idx = 0; input_idx < n_inputs; ++input_idx) {
0794     cellGridMatrix.flat<float>()(input_idx) = 0;
0795   }
0796   fillL1TauVars(cellGridMatrix, allTaus);
0797 
0798   fillCaloRecHits(cellGridMatrix, allTaus, caloRecHits);
0799 
0800   fillPatatracks(cellGridMatrix, allTaus, patatracks_SoA, vertices_SoA, *bsHandle, fieldESH.product());
0801 
0802   standardizeTensor(cellGridMatrix);
0803 
0804   if (debugLevel_ > 0) {
0805     checknan(cellGridMatrix, debugLevel_);
0806   }
0807 
0808   std::vector<float> tau_score = getTauScore(cellGridMatrix);
0809 
0810   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0811     const size_t nTau = TauCollectionMap[inp_idx].size();
0812     auto tau_tags = std::make_unique<std::vector<float>>(nTau);
0813     for (size_t tau_pos = 0; tau_pos < nTau; ++tau_pos) {
0814       const auto tau_idx = TauCollectionMap[inp_idx][tau_pos];
0815       if (debugLevel_ > 0) {
0816         edm::LogInfo("DebugInfo") << event.id().event() << " \t " << (allTaus[tau_idx])->pt() << " \t "
0817                                   << tau_score.at(tau_idx) << std::endl;
0818       }
0819       (*tau_tags)[tau_pos] = tau_score.at(tau_idx);
0820     }
0821     event.put(std::move(tau_tags), L1TauDesc_[inp_idx].CollectionName);
0822   }
0823 }
0824 //define this as a plug-in
0825 #include "FWCore/Framework/interface/MakerMacros.h"
0826 DEFINE_FWK_MODULE(L2TauNNProducerAlpaka);