Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-07-03 00:42:39 

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/TracksSoA.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 = reco::TracksHost;
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   boost::property_tree::ptree loadPtreeRoot;
0237   auto const normalizationDict = edm::FileInPath(cfg.getParameter<std::string>("normalizationDict")).fullPath();
0238   boost::property_tree::read_json(normalizationDict, loadPtreeRoot);
0239   for (const auto& [key, val] : L2TauTagNNv1::varNameMap) {
0240     boost::property_tree::ptree var = loadPtreeRoot.get_child(val);
0241     normDictElement current_element;
0242     current_element.mean = var.get_child("mean").get_value<float>();
0243     current_element.std = var.get_child("std").get_value<float>();
0244     current_element.min = var.get_child("min").get_value<float>();
0245     current_element.max = var.get_child("max").get_value<float>();
0246     cacheData->normVec.push_back(current_element);
0247   }
0248   return cacheData;
0249 }
0250 void L2TauNNProducerAlpaka::globalEndJob(L2TauNNProducerAlpakaCacheData* cacheData) {
0251   if (cacheData->graphDef != nullptr) {
0252     delete cacheData->graphDef;
0253   }
0254   tensorflow::closeSession(cacheData->session);
0255 }
0256 void L2TauNNProducerAlpaka::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0257   edm::ParameterSetDescription desc;
0258   desc.add<int>("debugLevel", 0)->setComment("set debug level for printing out info");
0259   edm::ParameterSetDescription l1TausPset;
0260   l1TausPset.add<std::string>("L1CollectionName", "DoubleTau")->setComment("Name of collections");
0261   l1TausPset.add<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"))
0262       ->setComment("Which trigger should the L1 Taus collection pass");
0263   edm::ParameterSet l1TausPSetDefault;
0264   l1TausPSetDefault.addParameter<std::string>("L1CollectionName", "DoubleTau");
0265   l1TausPSetDefault.addParameter<edm::InputTag>("L1TauTrigger", edm::InputTag("hltL1sDoubleTauBigOR"));
0266   desc.addVPSet("L1Taus", l1TausPset, {l1TausPSetDefault});
0267   desc.add<edm::InputTag>("hbheInput", edm::InputTag("hltHbhereco"))->setComment("HBHE recHit collection");
0268   desc.add<edm::InputTag>("hoInput", edm::InputTag("hltHoreco"))->setComment("HO recHit Collection");
0269   desc.add<edm::InputTag>("ebInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEB"))->setComment("EB recHit Collection");
0270   desc.add<edm::InputTag>("eeInput", edm::InputTag("hltEcalRecHit:EcalRecHitsEE"))->setComment("EE recHit Collection");
0271   desc.add<edm::InputTag>("pataVertices", edm::InputTag("hltPixelVerticesSoA"))
0272       ->setComment("patatrack vertices collection");
0273   desc.add<edm::InputTag>("pataTracks", edm::InputTag("hltPixelTracksSoA"))->setComment("patatrack collection");
0274   desc.add<edm::InputTag>("BeamSpot", edm::InputTag("hltOnlineBeamSpot"))->setComment("BeamSpot Collection");
0275   desc.add<uint>("maxVtx", 100)->setComment("max output collection size (number of accepted vertices)");
0276   desc.add<double>("fractionSumPt2", 0.3)->setComment("threshold on sumPt2 fraction of the leading vertex");
0277   desc.add<double>("minSumPt2", 0.)->setComment("min sumPt2");
0278   desc.add<double>("track_pt_min", 1.0)->setComment("min track p_T");
0279   desc.add<double>("track_pt_max", 10.0)->setComment("max track p_T");
0280   desc.add<double>("track_chi2_max", 99999.)->setComment("max track chi2");
0281   desc.add<std::string>("graphPath", "RecoTauTag/TrainingFiles/data/L2TauNNTag/L2TauTag_Run3v1.pb")
0282       ->setComment("path to the saved CNN");
0283   desc.add<std::string>("normalizationDict", "RecoTauTag/TrainingFiles/data/L2TauNNTag/NormalizationDict.json")
0284       ->setComment("path to the dictionary for variable standardization");
0285   descriptions.addWithDefaultLabel(desc);
0286 }
0287 
0288 L2TauNNProducerAlpaka::L2TauNNProducerAlpaka(const edm::ParameterSet& cfg,
0289                                              const L2TauNNProducerAlpakaCacheData* 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 L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::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> L2TauNNProducerAlpaka::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> L2TauNNProducerAlpaka::getEtaPhiIndices(const VPos& position, const LVec& tau_p4) {
0432   return getEtaPhiIndices(position.eta(), position.phi(), tau_p4);
0433 }
0434 
0435 void L2TauNNProducerAlpaka::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 L2TauNNProducerAlpaka::selectGoodTracksAndVertices(const ZVertexHost& patavtx_soa,
0575                                                         const TracksHost& patatracks_tsoa,
0576                                                         std::vector<int>& trkGood,
0577                                                         std::vector<int>& vtxGood) {
0578   const auto maxTracks = patatracks_tsoa.view().metadata().size();
0579   const int nv = patavtx_soa.view().nvFinal();
0580   trkGood.clear();
0581   trkGood.reserve(maxTracks);
0582   vtxGood.clear();
0583   vtxGood.reserve(nv);
0584   auto const* quality = patatracks_tsoa.view().quality();
0585 
0586   // No need to sort either as the algorithms is just using the max (not even the location, just the max value of pt2sum).
0587   std::vector<float> pTSquaredSum(nv, 0);
0588   std::vector<int> nTrkAssociated(nv, 0);
0589 
0590   for (int32_t trk_idx = 0; trk_idx < maxTracks; ++trk_idx) {
0591     auto n_hits = nHits(patatracks_tsoa.view(), trk_idx);
0592     if (n_hits == 0) {
0593       break;
0594     }
0595     int vtx_ass_to_track = patavtx_soa.view<reco::ZVertexTracksSoA>()[trk_idx].idv();
0596     if (vtx_ass_to_track >= 0 && vtx_ass_to_track < nv) {
0597       auto patatrackPt = patatracks_tsoa.view()[trk_idx].pt();
0598       ++nTrkAssociated[vtx_ass_to_track];
0599       if (patatrackPt >= trackPtMin_ && patatracks_tsoa.const_view()[trk_idx].chi2() <= trackChi2Max_) {
0600         patatrackPt = std::min(patatrackPt, trackPtMax_);
0601         pTSquaredSum[vtx_ass_to_track] += patatrackPt * patatrackPt;
0602       }
0603     }
0604     if (n_hits > 0 and quality[trk_idx] >= pixelTrack::Quality::loose) {
0605       trkGood.push_back(trk_idx);
0606     }
0607   }
0608   if (nv > 0) {
0609     const auto minFOM_fromFrac = (*std::max_element(pTSquaredSum.begin(), pTSquaredSum.end())) * fractionSumPt2_;
0610     for (int j = nv - 1; j >= 0 && vtxGood.size() < maxVtx_; --j) {
0611       auto vtx_idx = patavtx_soa.view()[j].sortInd();
0612       assert(vtx_idx < nv);
0613       if (nTrkAssociated[vtx_idx] >= 2 && pTSquaredSum[vtx_idx] >= minFOM_fromFrac &&
0614           pTSquaredSum[vtx_idx] > minSumPt2_) {
0615         vtxGood.push_back(vtx_idx);
0616       }
0617     }
0618   }
0619 }
0620 
0621 std::pair<float, float> L2TauNNProducerAlpaka::impactParameter(int it,
0622                                                                const TracksHost& patatracks_tsoa,
0623                                                                float patatrackPhi,
0624                                                                const reco::BeamSpot& beamspot,
0625                                                                const MagneticField* magfi) {
0626   /* dxy and dz */
0627   riemannFit::Vector5d ipar, opar;
0628   riemannFit::Matrix5d icov, ocov;
0629   copyToDense(patatracks_tsoa.view(), ipar, icov, it);
0630   riemannFit::transformToPerigeePlane(ipar, icov, opar, ocov);
0631   LocalTrajectoryParameters lpar(opar(0), opar(1), opar(2), opar(3), opar(4), 1.);
0632   float sp = std::sin(patatrackPhi);
0633   float cp = std::cos(patatrackPhi);
0634   Surface::RotationType Rotation(sp, -cp, 0, 0, 0, -1.f, cp, sp, 0);
0635   GlobalPoint BeamSpotPoint(beamspot.x0(), beamspot.y0(), beamspot.z0());
0636   Plane impPointPlane(BeamSpotPoint, Rotation);
0637   GlobalTrajectoryParameters gp(
0638       impPointPlane.toGlobal(lpar.position()), impPointPlane.toGlobal(lpar.momentum()), lpar.charge(), magfi);
0639   GlobalPoint vv = gp.position();
0640   math::XYZPoint pos(vv.x(), vv.y(), vv.z());
0641   GlobalVector pp = gp.momentum();
0642   math::XYZVector mom(pp.x(), pp.y(), pp.z());
0643   auto lambda = M_PI_2 - pp.theta();
0644   auto phi = pp.phi();
0645   float patatrackDxy = -vv.x() * std::sin(phi) + vv.y() * std::cos(phi);
0646   float patatrackDz =
0647       (vv.z() * std::cos(lambda) - (vv.x() * std::cos(phi) + vv.y() * std::sin(phi)) * std::sin(lambda)) /
0648       std::cos(lambda);
0649   return std::make_pair(patatrackDxy, patatrackDz);
0650 }
0651 
0652 void L2TauNNProducerAlpaka::fillPatatracks(tensorflow::Tensor& cellGridMatrix,
0653                                            const std::vector<l1t::TauRef>& allTaus,
0654                                            const TracksHost& patatracks_tsoa,
0655                                            const ZVertexHost& patavtx_soa,
0656                                            const reco::BeamSpot& beamspot,
0657                                            const MagneticField* magfi) {
0658   using NNInputs = L2TauTagNNv1::NNInputs;
0659 
0660   float deta, dphi;
0661   int eta_idx = 0;
0662   int phi_idx = 0;
0663   int tau_idx = 0;
0664 
0665   auto getCell = [&](NNInputs input) -> float& {
0666     return getCellImpl(cellGridMatrix, tau_idx, phi_idx, eta_idx, input);
0667   };
0668 
0669   std::vector<int> trkGood;
0670   std::vector<int> vtxGood;
0671 
0672   selectGoodTracksAndVertices(patavtx_soa, patatracks_tsoa, trkGood, vtxGood);
0673 
0674   const int nTaus = allTaus.size();
0675   for (tau_idx = 0; tau_idx < nTaus; tau_idx++) {
0676     const float tauEta = allTaus[tau_idx]->eta();
0677     const float tauPhi = allTaus[tau_idx]->phi();
0678 
0679     for (const auto it : trkGood) {
0680       const float patatrackPt = patatracks_tsoa.const_view()[it].pt();
0681       if (patatrackPt <= 0)
0682         continue;
0683       const float patatrackPhi = reco::phi(patatracks_tsoa.const_view(), it);
0684       const float patatrackEta = patatracks_tsoa.const_view()[it].eta();
0685       const float patatrackCharge = reco::charge(patatracks_tsoa.const_view(), it);
0686       const float patatrackChi2OverNdof = patatracks_tsoa.view()[it].chi2();
0687       const auto n_hits = nHits(patatracks_tsoa.const_view(), it);
0688       if (n_hits <= 0)
0689         continue;
0690       const int patatrackNdof = 2 * std::min(6, n_hits) - 5;
0691 
0692       const int vtx_idx_assTrk = patavtx_soa.view<reco::ZVertexTracksSoA>()[it].idv();
0693       if (reco::deltaR2(patatrackEta, patatrackPhi, tauEta, tauPhi) < dR2_max) {
0694         std::tie(deta, dphi, eta_idx, phi_idx) =
0695             getEtaPhiIndices(patatrackEta, patatrackPhi, allTaus[tau_idx]->polarP4());
0696         getCell(NNInputs::PatatrackPtSum) += patatrackPt;
0697         getCell(NNInputs::PatatrackSize) += 1.;
0698         getCell(NNInputs::PatatrackChargeSum) += patatrackCharge;
0699         getCell(NNInputs::PatatrackDeltaEta) += deta * patatrackPt;
0700         getCell(NNInputs::PatatrackDeltaPhi) += dphi * patatrackPt;
0701         getCell(NNInputs::PatatrackChi2OverNdof) += patatrackChi2OverNdof * patatrackPt;
0702         getCell(NNInputs::PatatrackNdof) += patatrackNdof * patatrackPt;
0703         std::pair<float, float> impactParameters = impactParameter(it, patatracks_tsoa, patatrackPhi, beamspot, magfi);
0704         getCell(NNInputs::PatatrackDxy) += impactParameters.first * patatrackPt;
0705         getCell(NNInputs::PatatrackDz) += impactParameters.second * patatrackPt;
0706         if ((std::find(vtxGood.begin(), vtxGood.end(), vtx_idx_assTrk) != vtxGood.end())) {
0707           getCell(NNInputs::PatatrackPtSumWithVertex) += patatrackPt;
0708           getCell(NNInputs::PatatrackSizeWithVertex) += 1.;
0709         }
0710       }
0711     }
0712 
0713     // normalize to sum and define stdDev
0714     for (eta_idx = 0; eta_idx < L2TauTagNNv1::nCellEta; eta_idx++) {
0715       for (phi_idx = 0; phi_idx < L2TauTagNNv1::nCellPhi; phi_idx++) {
0716         getCell(NNInputs::nVertices) = vtxGood.size();
0717         if (getCell(NNInputs::PatatrackPtSum) > 0.) {
0718           getCell(NNInputs::PatatrackDeltaEta) /= getCell(NNInputs::PatatrackPtSum);
0719           getCell(NNInputs::PatatrackDeltaPhi) /= getCell(NNInputs::PatatrackPtSum);
0720           getCell(NNInputs::PatatrackChi2OverNdof) /= getCell(NNInputs::PatatrackPtSum);
0721           getCell(NNInputs::PatatrackNdof) /= getCell(NNInputs::PatatrackPtSum);
0722           getCell(NNInputs::PatatrackDxy) /= getCell(NNInputs::PatatrackPtSum);
0723           getCell(NNInputs::PatatrackDz) /= getCell(NNInputs::PatatrackPtSum);
0724         }
0725       }
0726     }
0727   }
0728 }
0729 
0730 std::vector<float> L2TauNNProducerAlpaka::getTauScore(const tensorflow::Tensor& cellGridMatrix) {
0731   const int nTau = cellGridMatrix.shape().dim_size(0);
0732   if (nTau == 0) {
0733     return std::vector<float>();
0734   } else {
0735     std::vector<tensorflow::Tensor> pred_tensor;
0736     tensorflow::run(L2cacheData_->session, {{inputTensorName_, cellGridMatrix}}, {outputTensorName_}, &pred_tensor);
0737     std::vector<float> pred_vector(nTau);
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 
0746 void L2TauNNProducerAlpaka::produce(edm::Event& event, const edm::EventSetup& eventsetup) {
0747   std::vector<std::vector<size_t>> TauCollectionMap(L1TauDesc_.size());
0748   l1t::TauVectorRef allTaus;
0749 
0750   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0751     l1t::TauVectorRef l1Taus;
0752     auto const& l1TriggeredTaus = event.get(L1TauDesc_[inp_idx].inputToken_);
0753     l1TriggeredTaus.getObjects(trigger::TriggerL1Tau, l1Taus);
0754     TauCollectionMap.at(inp_idx).resize(l1Taus.size());
0755 
0756     for (size_t l1_idx = 0; l1_idx < l1Taus.size(); l1_idx++) {
0757       size_t tau_idx;
0758       const auto iter = std::find(allTaus.begin(), allTaus.end(), l1Taus[l1_idx]);
0759       if (iter != allTaus.end()) {
0760         tau_idx = std::distance(allTaus.begin(), iter);
0761       } else {
0762         allTaus.push_back(l1Taus[l1_idx]);
0763         tau_idx = allTaus.size() - 1;
0764       }
0765       TauCollectionMap.at(inp_idx).at(l1_idx) = tau_idx;
0766     }
0767   }
0768   const auto ebCal = event.getHandle(ebToken_);
0769   const auto eeCal = event.getHandle(eeToken_);
0770   const auto hbhe = event.getHandle(hbheToken_);
0771   const auto ho = event.getHandle(hoToken_);
0772   auto const& patatracks_SoA = event.get(pataTracksToken_);
0773   auto const& vertices_SoA = event.get(pataVerticesToken_);
0774   const auto bsHandle = event.getHandle(beamSpotToken_);
0775 
0776   auto const fieldESH = eventsetup.getHandle(bFieldToken_);
0777   auto const geometry = eventsetup.getHandle(geometryToken_);
0778 
0779   caloRecHitCollections caloRecHits;
0780   caloRecHits.hbhe = &*hbhe;
0781   caloRecHits.ho = &*ho;
0782   caloRecHits.eb = &*ebCal;
0783   caloRecHits.ee = &*eeCal;
0784   caloRecHits.geometry = &*geometry;
0785 
0786   const int nTaus = allTaus.size();
0787   tensorflow::Tensor cellGridMatrix(tensorflow::DT_FLOAT,
0788                                     {nTaus, L2TauTagNNv1::nCellEta, L2TauTagNNv1::nCellPhi, L2TauTagNNv1::nVars});
0789   const int n_inputs = nTaus * L2TauTagNNv1::nCellEta * L2TauTagNNv1::nCellPhi * L2TauTagNNv1::nVars;
0790   for (int input_idx = 0; input_idx < n_inputs; ++input_idx) {
0791     cellGridMatrix.flat<float>()(input_idx) = 0;
0792   }
0793   fillL1TauVars(cellGridMatrix, allTaus);
0794 
0795   fillCaloRecHits(cellGridMatrix, allTaus, caloRecHits);
0796 
0797   fillPatatracks(cellGridMatrix, allTaus, patatracks_SoA, vertices_SoA, *bsHandle, fieldESH.product());
0798 
0799   standardizeTensor(cellGridMatrix);
0800 
0801   if (debugLevel_ > 0) {
0802     checknan(cellGridMatrix, debugLevel_);
0803   }
0804 
0805   std::vector<float> tau_score = getTauScore(cellGridMatrix);
0806 
0807   for (size_t inp_idx = 0; inp_idx < L1TauDesc_.size(); inp_idx++) {
0808     const size_t nTau = TauCollectionMap[inp_idx].size();
0809     auto tau_tags = std::make_unique<std::vector<float>>(nTau);
0810     for (size_t tau_pos = 0; tau_pos < nTau; ++tau_pos) {
0811       const auto tau_idx = TauCollectionMap[inp_idx][tau_pos];
0812       if (debugLevel_ > 0) {
0813         edm::LogInfo("DebugInfo") << event.id().event() << " \t " << (allTaus[tau_idx])->pt() << " \t "
0814                                   << tau_score.at(tau_idx) << std::endl;
0815       }
0816       (*tau_tags)[tau_pos] = tau_score.at(tau_idx);
0817     }
0818     event.put(std::move(tau_tags), L1TauDesc_[inp_idx].CollectionName);
0819   }
0820 }
0821 //define this as a plug-in
0822 #include "FWCore/Framework/interface/MakerMacros.h"
0823 DEFINE_FWK_MODULE(L2TauNNProducerAlpaka);
[ 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!