diff --git a/DataFormats/L1TrackTrigger/interface/TTTrack.h b/DataFormats/L1TrackTrigger/interface/TTTrack.h index 8428d694ac779..19ceec4c53453 100644 --- a/DataFormats/L1TrackTrigger/interface/TTTrack.h +++ b/DataFormats/L1TrackTrigger/interface/TTTrack.h @@ -443,7 +443,6 @@ void TTTrack::setTrackWordBits() { } unsigned int valid = true; - unsigned int mvaQuality = 0; unsigned int mvaOther = 0; // missing conversion of global phi to difference from sector center phi @@ -457,7 +456,7 @@ void TTTrack::setTrackWordBits() { 0, theStubPtConsistency_, theHitPattern_, - mvaQuality, + theTrkMVA1_, mvaOther, thePhiSector_); } else { @@ -469,7 +468,7 @@ void TTTrack::setTrackWordBits() { chi2ZRed(), chi2BendRed(), theHitPattern_, - mvaQuality, + theTrkMVA1_, mvaOther, thePhiSector_); } diff --git a/DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h b/DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h index e4fce81139fe6..c235c7f545ba3 100644 --- a/DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h +++ b/DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h @@ -103,6 +103,8 @@ class TTTrack_TrackWord { {0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 8.0, 10.0, 20.0, 50.0}}; static constexpr std::array bendChi2Bins = { {0.0, 0.75, 1.0, 1.5, 2.25, 3.5, 5.0, 20.0}}; + static constexpr std::array mvaQualityBins = { + {0.0, 0.125, 0.25, 0.375, 0.5, 0.625, 0.750, 0.875}}; // Sector constants static constexpr unsigned int nSectors = 9; @@ -141,7 +143,7 @@ class TTTrack_TrackWord { double chi2RZ, double bendChi2, unsigned int hitPattern, - unsigned int mvaQuality, + double mvaQuality, unsigned int mvaOther, unsigned int sector); TTTrack_TrackWord(unsigned int valid, @@ -224,7 +226,7 @@ class TTTrack_TrackWord { double getBendChi2() const { return bendChi2Bins[getBendChi2Bits()]; } unsigned int getHitPattern() const { return getHitPatternBits(); } unsigned int getNStubs() const { return countSetBits(getHitPatternBits()); } - unsigned int getMVAQuality() const { return getMVAQualityBits(); } + double getMVAQuality() const { return mvaQualityBins[getMVAQualityBits()]; } unsigned int getMVAOther() const { return getMVAOtherBits(); } // ----------member functions (setters) ------------ @@ -236,7 +238,7 @@ class TTTrack_TrackWord { double chi2RZ, double bendChi2, unsigned int hitPattern, - unsigned int mvaQuality, + double mvaQuality, unsigned int mvaOther, unsigned int sector); diff --git a/DataFormats/L1TrackTrigger/src/TTTrack_TrackWord.cc b/DataFormats/L1TrackTrigger/src/TTTrack_TrackWord.cc index 0b1420d2ca09f..f4428e0cf8d55 100644 --- a/DataFormats/L1TrackTrigger/src/TTTrack_TrackWord.cc +++ b/DataFormats/L1TrackTrigger/src/TTTrack_TrackWord.cc @@ -35,7 +35,7 @@ TTTrack_TrackWord::TTTrack_TrackWord(unsigned int valid, double chi2RZ, double bendChi2, unsigned int hitPattern, - unsigned int mvaQuality, + double mvaQuality, unsigned int mvaOther, unsigned int sector) { setTrackWord(valid, momentum, POCA, rInv, chi2RPhi, chi2RZ, bendChi2, hitPattern, mvaQuality, mvaOther, sector); @@ -65,7 +65,7 @@ void TTTrack_TrackWord::setTrackWord(unsigned int valid, double chi2RZ, double bendChi2, unsigned int hitPattern, - unsigned int mvaQuality, + double mvaQuality, unsigned int mvaOther, unsigned int sector) { // first, derive quantities to be packed @@ -85,7 +85,7 @@ void TTTrack_TrackWord::setTrackWord(unsigned int valid, chi2rz_t chi2RZ_ = getBin(chi2RZ, chi2RZBins); bendChi2_t bendChi2_ = getBin(bendChi2, bendChi2Bins); hit_t hitPattern_ = hitPattern; - qualityMVA_t mvaQuality_ = mvaQuality; + qualityMVA_t mvaQuality_ = getBin(mvaQuality, mvaQualityBins); otherMVA_t mvaOther_ = mvaOther; // pack the track word diff --git a/L1Trigger/L1TTrackMatch/BuildFile.xml b/L1Trigger/L1TTrackMatch/BuildFile.xml index c8d624409711b..2caaa0df5472c 100644 --- a/L1Trigger/L1TTrackMatch/BuildFile.xml +++ b/L1Trigger/L1TTrackMatch/BuildFile.xml @@ -3,6 +3,7 @@ + diff --git a/L1Trigger/L1TTrackMatch/plugins/L1TrackVertexAssociationProducer.cc b/L1Trigger/L1TTrackMatch/plugins/L1TrackVertexAssociationProducer.cc index 36438dc052f5b..34d7e0970abeb 100644 --- a/L1Trigger/L1TTrackMatch/plugins/L1TrackVertexAssociationProducer.cc +++ b/L1Trigger/L1TTrackMatch/plugins/L1TrackVertexAssociationProducer.cc @@ -12,8 +12,8 @@ std::vector - Each floating point TTTrack inside this collection inherits from a bit-accurate TTTrack_TrackWord, used for emulation purposes. Outputs: - std::vector - A collection of TTTracks selected from cuts on the TTTrack properties - std::vector - A collection of TTTracks selected from cuts on the TTTrack_TrackWord properties + std::vector - A collection of TTTrack Refs selected from cuts on the TTTrack properties + std::vector - A collection of TTTrack Refs selected from cuts on the TTTrack_TrackWord properties */ // // Original Author: Alexx Perloff @@ -62,7 +62,9 @@ #include "FWCore/Utilities/interface/EDMException.h" #include "FWCore/Utilities/interface/StreamID.h" #include "Geometry/Records/interface/TrackerTopologyRcd.h" +#include "PhysicsTools/TensorFlow/interface/TensorFlow.h" #include "L1Trigger/DemonstratorTools/interface/codecs/tracks.h" +#include "FWCore/ParameterSet/interface/FileInPath.h" // // class declaration @@ -149,6 +151,76 @@ class L1TrackVertexAssociationProducer : public edm::global::EDProducer<> { std::vector deltaZMaxEtaBounds_; std::vector deltaZMax_; }; + + struct NNTrackWordSelector { + NNTrackWordSelector(tensorflow::Session* AssociationSesh, + const double AssociationThreshold, + const std::vector& AssociationNetworkZ0binning, + const std::vector& AssociationNetworkEtaBounds, + const std::vector& AssociationNetworkZ0ResBins) + : AssociationSesh_(AssociationSesh), + AssociationThreshold_(AssociationThreshold), + z0_binning_(AssociationNetworkZ0binning), + eta_bins_(AssociationNetworkEtaBounds), + res_bins_(AssociationNetworkZ0ResBins) {} + + bool operator()(const TTTrackType& t, const l1t::VertexWord& v) const { + tensorflow::Tensor inputAssoc(tensorflow::DT_FLOAT, {1, 4}); + std::vector outputAssoc; + + TTTrack_TrackWord::tanl_t etaEmulationBits = t.getTanlWord(); + ap_fixed<16, 3> etaEmulation; + etaEmulation.V = (etaEmulationBits.range()); + + auto lower = std::lower_bound(eta_bins_.begin(), eta_bins_.end(), etaEmulation.to_double()); + + int resbin = std::distance(eta_bins_.begin(), lower); + float binWidth = z0_binning_[2]; + // Calculate integer dZ from track z0 and vertex z0 (use floating point version and convert internally allowing use of both emulator and simulator vertex and track) + float dZ = + abs(floor(((t.getZ0() + z0_binning_[1]) / (binWidth))) - floor(((v.z0() + z0_binning_[1]) / (binWidth)))); + + // The following constants <14, 9>, <22, 9> are defined by the quantisation of the Neural Network + ap_uint<14> ptEmulationBits = t.getTrackWord()(TTTrack_TrackWord::TrackBitLocations::kRinvMSB - 1, + TTTrack_TrackWord::TrackBitLocations::kRinvLSB); + ap_ufixed<14, 9> ptEmulation; + ptEmulation.V = (ptEmulationBits.range()); + + ap_ufixed<22, 9> ptEmulation_rescale; + ptEmulation_rescale = ptEmulation.to_double(); + + ap_ufixed<22, 9> resBinEmulation_rescale; + resBinEmulation_rescale = res_bins_[resbin]; + + ap_ufixed<22, 9> MVAEmulation_rescale; + MVAEmulation_rescale = t.getMVAQualityBits(); + + ap_ufixed<22, 9> dZEmulation_rescale; + dZEmulation_rescale = dZ; + + inputAssoc.tensor()(0, 0) = ptEmulation_rescale.to_double(); + inputAssoc.tensor()(0, 1) = MVAEmulation_rescale.to_double(); + inputAssoc.tensor()(0, 2) = resBinEmulation_rescale.to_double() / 16.0; + inputAssoc.tensor()(0, 3) = dZEmulation_rescale.to_double(); + + // Run Association Network: + tensorflow::run(AssociationSesh_, {{"assoc:0", inputAssoc}}, {"Identity:0"}, &outputAssoc); + + double NNOutput = (double)outputAssoc[0].tensor()(0, 0); + + double NNOutput_exp = 1.0 / (1.0 + exp(-1.0 * (NNOutput))); + + return NNOutput_exp >= AssociationThreshold_; + } + + private: + tensorflow::Session* AssociationSesh_; + double AssociationThreshold_; + std::vector z0_binning_; + std::vector eta_bins_; + std::vector res_bins_; + }; + struct TTTrackWordLinkLimitSelector { TTTrackWordLinkLimitSelector(const unsigned int fwNTrackSetsTVA) : fwNTrackSetsTVA_(fwNTrackSetsTVA) { //create a counter for all 18 GTT input links, 2 per phiSector of the TrackFindingProcessors @@ -209,6 +281,15 @@ class L1TrackVertexAssociationProducer : public edm::global::EDProducer<> { const double useDisplacedTracksDeltaZOverride_; // corresponds to N_TRACK_SETS_TVA in LibHLS https://gitlab.cern.ch/GTT/LibHLS/-/blob/master/DataFormats/Track/interface/TrackConstants.h const unsigned int fwNTrackSetsTVA_; + + //NNVtx: + edm::FileInPath associationGraphPath_; + const double associationThreshold_; + bool useAssociationNetwork_; + tensorflow::GraphDef* associationGraph_; + tensorflow::Session* associationSesh_; + std::vector associationNetworkZ0binning_, associationNetworkEtaBounds_, associationNetworkZ0ResBins_; + int debug_; }; @@ -241,7 +322,17 @@ L1TrackVertexAssociationProducer::L1TrackVertexAssociationProducer(const edm::Pa deltaZMax_(cutSet_.getParameter>("deltaZMax")), useDisplacedTracksDeltaZOverride_(iConfig.getParameter("useDisplacedTracksDeltaZOverride")), fwNTrackSetsTVA_(iConfig.getParameter("fwNTrackSetsTVA")), + associationThreshold_(iConfig.getParameter("associationThreshold")), + useAssociationNetwork_(iConfig.getParameter("useAssociationNetwork")), + associationNetworkZ0binning_(iConfig.getParameter>("associationNetworkZ0binning")), + associationNetworkEtaBounds_(iConfig.getParameter>("associationNetworkEtaBounds")), + associationNetworkZ0ResBins_(iConfig.getParameter>("associationNetworkZ0ResBins")), debug_(iConfig.getParameter("debug")) { + if (useAssociationNetwork_) { + associationGraphPath_ = iConfig.getParameter("associationGraph"); + associationGraph_ = tensorflow::loadGraphDef(associationGraphPath_.fullPath()); + associationSesh_ = tensorflow::createSession(associationGraph_); + } // Confirm the the configuration makes sense if (!processSimulatedTracks_ && !processEmulatedTracks_) { throw cms::Exception("You must process at least one of the track collections (simulated or emulated)."); @@ -408,6 +499,12 @@ void L1TrackVertexAssociationProducer::produce(edm::StreamID, edm::Event& iEvent TTTrackDeltaZMaxSelector deltaZSel(deltaZMaxEtaBounds_, deltaZMax_); TTTrackWordDeltaZMaxSelector deltaZSelEmu(deltaZMaxEtaBounds_, deltaZMax_); + NNTrackWordSelector TTTrackNetworkSelector(associationSesh_, + associationThreshold_, + associationNetworkZ0binning_, + associationNetworkEtaBounds_, + associationNetworkZ0ResBins_); + iEvent.getByToken(l1TracksToken_, l1TracksHandle); size_t nOutputApproximate = l1TracksHandle->size(); @@ -454,11 +551,17 @@ void L1TrackVertexAssociationProducer::produce(edm::StreamID, edm::Event& iEvent [track](const auto& ref) { return (*ref).getTrackWord() == track.getTrackWord(); }); bool passSelectionEmu = (itrEmu != l1SelectedTracksEmulationHandle->end()); // Associated tracks based on the bitwise accurate TTTrack_TrackWord - if (passLinkLimitEmu && passSelectionEmu && deltaZSelEmu(track, l1VerticesEmulationHandle->at(0))) { - vTTTrackAssociatedEmulationOutput->push_back(TTTrackRef(l1TracksHandle, i)); - } //end block for satisfying LinkLimitEmu and SelectionEmu criteria - } //end if (processEmulatedTracks_) - } //end loop over input converted tracks + if (useAssociationNetwork_) { + if (passLinkLimitEmu && passSelectionEmu && TTTrackNetworkSelector(track, l1VerticesEmulationHandle->at(0))) { + vTTTrackAssociatedEmulationOutput->push_back(TTTrackRef(l1TracksHandle, i)); + } + } else { + if (passLinkLimitEmu && passSelectionEmu && deltaZSelEmu(track, l1VerticesEmulationHandle->at(0))) { + vTTTrackAssociatedEmulationOutput->push_back(TTTrackRef(l1TracksHandle, i)); + } //end block for satisfying LinkLimitEmu and SelectionEmu criteria + } //end if use track association NN + } //end if (processEmulatedTracks_) + } //end loop over input converted tracks if (processSimulatedTracks_) { iEvent.put(std::move(vTTTrackAssociatedOutput), outputCollectionName_); @@ -507,6 +610,15 @@ void L1TrackVertexAssociationProducer::fillDescriptions(edm::ConfigurationDescri desc.add("processEmulatedTracks", true) ->setComment("return selected tracks after cutting on the bitwise emulated values"); desc.add("fwNTrackSetsTVA", 94)->setComment("firmware limit on processed tracks per GTT input link"); + desc.add("useAssociationNetwork", false)->setComment("Enable Association Network"); + desc.add("associationThreshold", 0)->setComment("Association Network threshold for PV tracks"); + desc.addOptional("associationGraph")->setComment("Location of Association Network model file"); + desc.add>("associationNetworkZ0binning", {}) + ->setComment("z0 binning used for setting the input feature digitisation"); + desc.add>("associationNetworkEtaBounds", {}) + ->setComment("Eta bounds used to set z0 resolution input feature"); + desc.add>("associationNetworkZ0ResBins", {})->setComment("z0 resolution input feature bins"); + desc.add("debug", 0)->setComment("Verbosity levels: 0, 1, 2, 3"); descriptions.addWithDefaultLabel(desc); } diff --git a/L1Trigger/L1TTrackMatch/python/l1tTrackVertexAssociationProducer_cfi.py b/L1Trigger/L1TTrackMatch/python/l1tTrackVertexAssociationProducer_cfi.py index 04f1cd320b3da..dfeb627994627 100644 --- a/L1Trigger/L1TTrackMatch/python/l1tTrackVertexAssociationProducer_cfi.py +++ b/L1Trigger/L1TTrackMatch/python/l1tTrackVertexAssociationProducer_cfi.py @@ -1,4 +1,5 @@ import FWCore.ParameterSet.Config as cms +from L1Trigger.VertexFinder.l1tVertexProducer_cfi import l1tVertexProducer l1tTrackVertexAssociationProducer = cms.EDProducer('L1TrackVertexAssociationProducer', l1TracksInputTag = cms.InputTag("l1tGTTInputProducer","Level1TTTracksConverted"), @@ -18,7 +19,18 @@ processSimulatedTracks = cms.bool(True), # return selected tracks after cutting on the floating point values processEmulatedTracks = cms.bool(True), # return selected tracks after cutting on the bitwise emulated values fwNTrackSetsTVA = cms.uint32(94), # firmware limit on number of GTT converted tracks considered for primary vertex association - debug = cms.int32(0) # Verbosity levels: 0, 1, 2, 3, 4 + debug = cms.int32(0), # Verbosity levels: 0, 1, 2, 3, 4 +) + +l1tTrackVertexNNAssociationProducer = l1tTrackVertexAssociationProducer.clone( + processSimulatedTracks = cms.bool(True), # return selected tracks after cutting on the floating point values + processEmulatedTracks = cms.bool(True), # return selected tracks after cutting on the bitwise emulated values + useAssociationNetwork = cms.bool(True), #Enable Association Network + associationThreshold = cms.double(0.1), #Association Network threshold for PV tracks + associationGraph = cms.FileInPath("L1Trigger/L1TTrackMatch/data/NNVtx_AssociationModelGraph.pb"), #Location of Association Network model file + associationNetworkZ0binning = l1tVertexProducer.VertexReconstruction.FH_HistogramParameters, #Z0 binning used for setting the input feature digitisation + associationNetworkEtaBounds = cms.vdouble(0.0, 0.01984126984126984, 0.03968253968253968, 0.05952380952380952, 0.07936507936507936, 0.0992063492063492, 0.11904761904761904, 0.1388888888888889, 0.15873015873015872, 0.17857142857142855, 0.1984126984126984, 0.21825396825396826, 0.23809523809523808, 0.2579365079365079, 0.2777777777777778, 0.2976190476190476, 0.31746031746031744, 0.33730158730158727, 0.3571428571428571, 0.376984126984127, 0.3968253968253968, 0.41666666666666663, 0.4365079365079365, 0.45634920634920634, 0.47619047619047616, 0.496031746031746, 0.5158730158730158, 0.5357142857142857, 0.5555555555555556, 0.5753968253968254, 0.5952380952380952, 0.615079365079365, 0.6349206349206349, 0.6547619047619048, 0.6746031746031745, 0.6944444444444444, 0.7142857142857142, 0.7341269841269841, 0.753968253968254, 0.7738095238095237, 0.7936507936507936, 0.8134920634920635, 0.8333333333333333, 0.8531746031746031, 0.873015873015873, 0.8928571428571428, 0.9126984126984127, 0.9325396825396824, 0.9523809523809523, 0.9722222222222222, 0.992063492063492, 1.0119047619047619, 1.0317460317460316, 1.0515873015873016, 1.0714285714285714, 1.0912698412698412, 1.1111111111111112, 1.130952380952381, 1.1507936507936507, 1.1706349206349205, 1.1904761904761905, 1.2103174603174602, 1.23015873015873, 1.25, 1.2698412698412698, 1.2896825396825395, 1.3095238095238095, 1.3293650793650793, 1.349206349206349, 1.369047619047619, 1.3888888888888888, 1.4087301587301586, 1.4285714285714284, 1.4484126984126984, 1.4682539682539681, 1.488095238095238, 1.507936507936508, 1.5277777777777777, 1.5476190476190474, 1.5674603174603174, 1.5873015873015872, 1.607142857142857, 1.626984126984127, 1.6468253968253967, 1.6666666666666665, 1.6865079365079365, 1.7063492063492063, 1.726190476190476, 1.746031746031746, 1.7658730158730158, 1.7857142857142856, 1.8055555555555554, 1.8253968253968254, 1.8452380952380951, 1.865079365079365, 1.8849206349206349, 1.9047619047619047, 1.9246031746031744, 1.9444444444444444, 1.9642857142857142, 1.984126984126984, 2.003968253968254, 2.0238095238095237, 2.0436507936507935, 2.0634920634920633, 2.083333333333333, 2.1031746031746033, 2.123015873015873, 2.142857142857143, 2.1626984126984126, 2.1825396825396823, 2.202380952380952, 2.2222222222222223, 2.242063492063492, 2.261904761904762, 2.2817460317460316, 2.3015873015873014, 2.321428571428571, 2.341269841269841, 2.361111111111111, 2.380952380952381, 2.4007936507936507, 2.4206349206349205, 2.4404761904761902, 2.46031746031746, 2.4801587301587302, 2.5), #Eta bounds used to set z0 resolution input feature + associationNetworkZ0ResBins = cms.vdouble(127.0, 126.0, 126.0, 126.0, 125.0, 124.0, 123.0, 122.0, 120.0, 119.0, 117.0, 115.0, 114.0, 112.0, 110.0, 107.0, 105.0, 103.0, 101.0, 98.0, 96.0, 94.0, 91.0, 89.0, 87.0, 85.0, 82.0, 80.0, 78.0, 76.0, 74.0, 72.0, 70.0, 68.0, 66.0, 64.0, 62.0, 61.0, 59.0, 57.0, 56.0, 54.0, 53.0, 51.0, 50.0, 48.0, 47.0, 46.0, 45.0, 43.0, 42.0, 41.0, 40.0, 39.0, 38.0, 37.0, 36.0, 35.0, 34.0, 33.0, 33.0, 32.0, 31.0, 30.0, 30.0, 29.0, 28.0, 28.0, 27.0, 26.0, 26.0, 25.0, 24.0, 24.0, 23.0, 23.0, 22.0, 22.0, 21.0, 21.0, 21.0, 20.0, 20.0, 19.0, 19.0, 18.0, 18.0, 18.0, 17.0, 17.0, 17.0, 16.0, 16.0, 16.0, 15.0, 15.0, 15.0, 15.0, 14.0, 14.0, 14.0, 14.0, 13.0, 13.0, 13.0, 13.0, 12.0, 12.0, 12.0, 12.0, 12.0, 11.0, 11.0, 11.0, 11.0, 11.0, 10.0, 10.0, 10.0, 10.0, 10.0, 10.0, 9.0, 9.0, 9.0, 9.0, 9.0, 0.0), #z0 resolution input feature bins ) l1tTrackVertexAssociationProducerExtended = l1tTrackVertexAssociationProducer.clone( diff --git a/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py b/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py index e2519cb5dff7e..270e5dbaf2e09 100644 --- a/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py +++ b/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py @@ -18,6 +18,8 @@ DISPLACED = '' + +runVtxNN = True ############################################################ # import standard configurations ############################################################ @@ -99,6 +101,34 @@ # Primary vertex ############################################################ process.pPV = cms.Path(process.l1tVertexFinder) + + +if runVtxNN: + process.l1tVertexFinderEmulator = process.l1tVertexProducer.clone() + process.l1tVertexFinderEmulator.VertexReconstruction.Algorithm = "NNEmulation" + + process.l1tTrackSelectionProducer.cutSet = cms.PSet(ptMin = cms.double(2.0), # pt must be greater than this value, [GeV] + absEtaMax = cms.double(2.4), # absolute value of eta must be less than this value + absZ0Max = cms.double(15.0), # z0 must be less than this value, [cm] + nStubsMin = cms.int32(4), # number of stubs must be greater than or equal to this value + nPSStubsMin = cms.int32(0), # the number of stubs in the PS Modules must be greater than or equal to this value + promptMVAMin = cms.double(-1.0), # MVA must be greater than this value + reducedBendChi2Max = cms.double(999.0), # bend chi2 must be less than this value + reducedChi2RZMax = cms.double(999.0), # chi2rz/dof must be less than this value + reducedChi2RPhiMax = cms.double(999.0), # chi2rphi/dof must be less than this value + reducedChi2MaxNstub4 = cms.double(999.9), # chi2/dof with nstub==4 must be less than this value + reducedChi2MaxNstub5 = cms.double(999.9), # chi2/dof with nstub>4 must be less than this value + reducedBendChi2MaxNstub4 = cms.double(999.9), # bend chi2 with nstub==4 must be less than this value + reducedBendChi2MaxNstub5 = cms.double(999.9), # bend chi2 with nstub>4 must be less than this value + ), + VertexAssociator = process.l1tTrackVertexNNAssociationProducer + AssociationName = "l1tTrackVertexNNAssociationProducer" +else: + process.l1tVertexFinderEmulator = process.l1tVertexProducer.clone() + process.l1tVertexFinderEmulator.VertexReconstruction.Algorithm = "FHEmulation" + VertexAssociator = process.l1tTrackVertexAssociationProducer + AssociationName = "l1tTrackVertexAssociationProducer" + process.pPVemu = cms.Path(process.l1tVertexFinderEmulator) # HYBRID: prompt tracking @@ -108,7 +138,7 @@ process.pL1TrackSelection = cms.Path(process.l1tTrackSelectionProducer * process.l1tTrackSelectionProducerForJets * process.l1tTrackSelectionProducerForEtMiss) - process.pL1TrackVertexAssociation = cms.Path(process.l1tTrackVertexAssociationProducer* + process.pL1TrackVertexAssociation = cms.Path(VertexAssociator * process.l1tTrackVertexAssociationProducerForJets* process.l1tTrackVertexAssociationProducerForEtMiss) process.pL1TrackJets = cms.Path(process.l1tTrackJets) @@ -148,7 +178,7 @@ process.pL1TrackSelection = cms.Path(process.l1tTrackSelectionProducer * process.l1tTrackSelectionProducerExtended * process.l1tTrackSelectionProducerForJets * process.l1tTrackSelectionProducerExtendedForJets * process.l1tTrackSelectionProducerForEtMiss * process.l1tTrackSelectionProducerExtendedForEtMiss) - process.pL1TrackVertexAssociation = cms.Path(process.l1tTrackVertexAssociationProducer * process.l1tTrackVertexAssociationProducerExtended * + process.pL1TrackVertexAssociation = cms.Path(VertexAssociator * process.l1tTrackVertexAssociationProducerExtended * process.l1tTrackVertexAssociationProducerForJets * process.l1tTrackVertexAssociationProducerExtendedForJets * process.l1tTrackVertexAssociationProducerForEtMiss * process.l1tTrackVertexAssociationProducerExtendedForEtMiss) process.pL1TrackJets = cms.Path(process.l1tTrackJets*process.l1tTrackJetsExtended) @@ -194,9 +224,9 @@ L1TrackExtendedGTTInputTag = cms.InputTag("l1tGTTInputProducerExtended","Level1TTTracksExtendedConverted"), # TTTracks, extended, GTT converted L1TrackSelectedInputTag = cms.InputTag("l1tTrackSelectionProducer", "Level1TTTracksSelected"), # TTTracks, prompt, selected L1TrackSelectedEmulationInputTag = cms.InputTag("l1tTrackSelectionProducer", "Level1TTTracksSelectedEmulation"), # TTTracks, prompt, emulation, selected - L1TrackSelectedAssociatedInputTag = cms.InputTag("l1tTrackVertexAssociationProducer", "Level1TTTracksSelectedAssociated"), # TTTracks, prompt, selected, associated - L1TrackSelectedAssociatedEmulationInputTag = cms.InputTag("l1tTrackVertexAssociationProducer", "Level1TTTracksSelectedAssociatedEmulation"), # TTTracks, prompt, emulation, selected, associated - + L1TrackSelectedAssociatedInputTag = cms.InputTag(AssociationName, "Level1TTTracksSelectedAssociated"), # TTTracks, prompt, selected, associated + L1TrackSelectedAssociatedEmulationInputTag = cms.InputTag(AssociationName, "Level1TTTracksSelectedAssociatedEmulation"), # TTTracks, prompt, emulation, selected, associated + L1TrackSelectedForJetsInputTag = cms.InputTag("l1tTrackSelectionProducerForJets", "Level1TTTracksSelected"), # TTTracks, prompt, selected L1TrackSelectedEmulationForJetsInputTag = cms.InputTag("l1tTrackSelectionProducerForJets", "Level1TTTracksSelectedEmulation"), # TTTracks, prompt, emulation, selected L1TrackSelectedAssociatedForJetsInputTag = cms.InputTag("l1tTrackVertexAssociationProducerForJets", "Level1TTTracksSelectedAssociated"), # TTTracks, prompt, selected, associated diff --git a/L1Trigger/TrackFindingTracklet/plugins/L1FPGATrackProducer.cc b/L1Trigger/TrackFindingTracklet/plugins/L1FPGATrackProducer.cc index 5f130d07dead1..b6381bd4221e9 100644 --- a/L1Trigger/TrackFindingTracklet/plugins/L1FPGATrackProducer.cc +++ b/L1Trigger/TrackFindingTracklet/plugins/L1FPGATrackProducer.cc @@ -739,9 +739,6 @@ void L1FPGATrackProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSe aTrack.setStubPtConsistency( StubPtConsistency::getConsistency(aTrack, theTrackerGeom, tTopo, settings_.bfield(), settings_.nHelixPar())); - // set TTTrack word - aTrack.setTrackWordBits(); - if (trackQuality_) { trackQualityModel_->setL1TrackQuality(aTrack); } @@ -751,6 +748,8 @@ void L1FPGATrackProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSe // trackQualityModel_->setBonusFeatures(hph.bonusFeatures()); // } + // set TTTrack word + aTrack.setTrackWordBits(); // test track word //aTrack.testTrackWordBits(); diff --git a/L1Trigger/VertexFinder/BuildFile.xml b/L1Trigger/VertexFinder/BuildFile.xml index 82b4a81c0ce38..93fedc05049fa 100644 --- a/L1Trigger/VertexFinder/BuildFile.xml +++ b/L1Trigger/VertexFinder/BuildFile.xml @@ -12,6 +12,7 @@ + diff --git a/L1Trigger/VertexFinder/interface/AlgoSettings.h b/L1Trigger/VertexFinder/interface/AlgoSettings.h index b3fb487c2765d..f66f92edb2cb5 100644 --- a/L1Trigger/VertexFinder/interface/AlgoSettings.h +++ b/L1Trigger/VertexFinder/interface/AlgoSettings.h @@ -6,6 +6,7 @@ #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/Utilities/interface/Exception.h" +#include "FWCore/ParameterSet/interface/FileInPath.h" namespace l1tVertexFinder { @@ -19,7 +20,8 @@ namespace l1tVertexFinder { PVR, adaptiveVertexReconstruction, HPV, - Kmeans + Kmeans, + NNEmulation }; enum class Precision { Simulation, Emulation }; @@ -54,6 +56,9 @@ namespace l1tVertexFinder { double vx_histogram_min() const { return vx_histogram_parameters_.at(0); } double vx_histogram_max() const { return vx_histogram_parameters_.at(1); } double vx_histogram_binwidth() const { return vx_histogram_parameters_.at(2); } + int vx_histogram_numbins() const { + return (vx_histogram_parameters_.at(1) - vx_histogram_parameters_.at(0)) / vx_histogram_parameters_.at(2); + } // fastHisto assumed vertex width float vx_width() const { return vx_width_; } // fastHisto track selection control @@ -73,6 +78,10 @@ namespace l1tVertexFinder { unsigned int vx_NStubMin() const { return vx_NStubMin_; } unsigned int vx_NStubPSMin() const { return vx_NStubPSMin_; } + // Functions for NN: + std::string vx_trkw_graph() const { return vx_trkw_graph_.fullPath(); } + std::string vx_pattrec_graph() const { return vx_pattrec_graph_.fullPath(); } + //=== Debug printout unsigned int debug() const { return debug_; } @@ -119,7 +128,8 @@ namespace l1tVertexFinder { float vx_dbscan_mintracks_; unsigned int vx_kmeans_iterations_; unsigned int vx_kmeans_nclusters_; - + edm::FileInPath vx_trkw_graph_; //For NNVtx (TrackWeight) + edm::FileInPath vx_pattrec_graph_; //For NNVtx (PatternRec) // Debug printout unsigned int debug_; }; diff --git a/L1Trigger/VertexFinder/interface/L1Track.h b/L1Trigger/VertexFinder/interface/L1Track.h index 00a669923fc0e..25b87283d5205 100644 --- a/L1Trigger/VertexFinder/interface/L1Track.h +++ b/L1Trigger/VertexFinder/interface/L1Track.h @@ -23,6 +23,9 @@ namespace l1tVertexFinder { float phi0() const { return track_->momentum().phi(); }; float pt() const { return track_->momentum().transverse(); }; float z0() const { return track_->POCA().z(); }; + float weight() const { return weight_; }; + void setWeight(float w) { weight_ = w; }; + double MVA1() const { return track_->trkMVA1(); }; //Track Quality BDT // FIXME: Double check nPar=4 is correct float chi2dof() const { return track_->chi2Red(); }; @@ -33,6 +36,7 @@ namespace l1tVertexFinder { private: edm::Ptr> track_; + float weight_; }; } // namespace l1tVertexFinder diff --git a/L1Trigger/VertexFinder/interface/VertexFinder.h b/L1Trigger/VertexFinder/interface/VertexFinder.h index 58b6de952c627..8a04c515a13dc 100644 --- a/L1Trigger/VertexFinder/interface/VertexFinder.h +++ b/L1Trigger/VertexFinder/interface/VertexFinder.h @@ -9,6 +9,7 @@ #include "FWCore/MessageLogger/interface/MessageLogger.h" #include "L1Trigger/VertexFinder/interface/AlgoSettings.h" #include "L1Trigger/VertexFinder/interface/RecoVertex.h" +#include "PhysicsTools/TensorFlow/interface/TensorFlow.h" #include #include @@ -100,6 +101,10 @@ namespace l1tVertexFinder { void fastHisto(const TrackerTopology* tTopo); /// Histogramming algorithm (emulation) void fastHistoEmulation(); + /// NNVtx algorithm + void NNVtxEmulation(tensorflow::Session* TrackWeightSesh = nullptr, + tensorflow::Session* PatternRecSesh = nullptr, + tensorflow::Session* AssociationSesh = nullptr); /// Sort vertices in pT void sortVerticesInPt(); diff --git a/L1Trigger/VertexFinder/interface/VertexProducer.h b/L1Trigger/VertexFinder/interface/VertexProducer.h index 5f4d93bcc74ee..599e404238b75 100644 --- a/L1Trigger/VertexFinder/interface/VertexProducer.h +++ b/L1Trigger/VertexFinder/interface/VertexProducer.h @@ -17,6 +17,7 @@ #include "L1Trigger/VertexFinder/interface/AlgoSettings.h" #include "L1Trigger/VertexFinder/interface/RecoVertex.h" #include "L1Trigger/VertexFinder/interface/VertexFinder.h" +#include "PhysicsTools/TensorFlow/interface/TensorFlow.h" #include #include @@ -46,6 +47,11 @@ class VertexProducer : public edm::global::EDProducer<> { const edm::ESGetToken tTopoToken; const std::string outputCollectionName_; + tensorflow::GraphDef* TrkWGraph_; + tensorflow::Session* TrkWSesh_; + tensorflow::GraphDef* PattRecGraph_; + tensorflow::Session* PattRecSesh_; + l1tVertexFinder::AlgoSettings settings_; }; diff --git a/L1Trigger/VertexFinder/plugins/VertexProducer.cc b/L1Trigger/VertexFinder/plugins/VertexProducer.cc index a300be0502b69..b9f679695c28b 100644 --- a/L1Trigger/VertexFinder/plugins/VertexProducer.cc +++ b/L1Trigger/VertexFinder/plugins/VertexProducer.cc @@ -56,14 +56,29 @@ VertexProducer::VertexProducer(const edm::ParameterSet& iConfig) case Algorithm::Kmeans: edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using a kmeans algorithm"; break; + case Algorithm::NNEmulation: + edm::LogInfo("VertexProducer") << "VertexProducer::Finding vertices using the Neural Network Emulation"; + break; } //--- Define EDM output to be written to file (if required) - if (settings_.vx_algo() == Algorithm::fastHistoEmulation) { + if (settings_.vx_algo() == Algorithm::fastHistoEmulation || settings_.vx_algo() == Algorithm::NNEmulation) { produces(outputCollectionName_ + "Emulation"); } else { produces(outputCollectionName_); } + + if (settings_.vx_algo() == Algorithm::NNEmulation) { + // load graphs, create a new session and add the graphDef + if (settings_.debug() > 1) { + edm::LogInfo("VertexProducer") << "loading TrkWeight graph from " << settings_.vx_trkw_graph() << std::endl; + edm::LogInfo("VertexProducer") << "loading PatternRec graph from " << settings_.vx_pattrec_graph() << std::endl; + } + TrkWGraph_ = tensorflow::loadGraphDef(settings_.vx_trkw_graph()); + TrkWSesh_ = tensorflow::createSession(TrkWGraph_); + PattRecGraph_ = tensorflow::loadGraphDef(settings_.vx_pattrec_graph()); + PattRecSesh_ = tensorflow::createSession(PattRecGraph_); + } } void VertexProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const { @@ -127,13 +142,16 @@ void VertexProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::Event case Algorithm::Kmeans: vf.Kmeans(); break; + case Algorithm::NNEmulation: + vf.NNVtxEmulation(TrkWSesh_, PattRecSesh_); + break; } vf.sortVerticesInPt(); vf.findPrimaryVertex(); // //=== Store output EDM track and hardware stub collections. - if (settings_.vx_algo() == Algorithm::fastHistoEmulation) { + if (settings_.vx_algo() == Algorithm::fastHistoEmulation || settings_.vx_algo() == Algorithm::NNEmulation) { std::unique_ptr product_emulation = std::make_unique(vf.verticesEmulation().begin(), vf.verticesEmulation().end()); iEvent.put(std::move(product_emulation), outputCollectionName_ + "Emulation"); diff --git a/L1Trigger/VertexFinder/python/l1tVertexProducer_cfi.py b/L1Trigger/VertexFinder/python/l1tVertexProducer_cfi.py index 682246189ed49..94e52373ccead 100644 --- a/L1Trigger/VertexFinder/python/l1tVertexProducer_cfi.py +++ b/L1Trigger/VertexFinder/python/l1tVertexProducer_cfi.py @@ -67,6 +67,10 @@ VxMinNStub = cms.uint32(4), # Minimum number of stubs in PS modules associated to a track VxMinNStubPS = cms.uint32(3), + # Track weight NN graph + TrackWeightGraph = cms.FileInPath("L1Trigger/VertexFinder/data/NNVtx_WeightModelGraph.pb"), + # Pattern recognition NN graph + PatternRecGraph = cms.FileInPath("L1Trigger/VertexFinder/data/NNVtx_PatternModelGraph.pb"), ), # Debug printout debug = cms.uint32(0) diff --git a/L1Trigger/VertexFinder/src/AlgoSettings.cc b/L1Trigger/VertexFinder/src/AlgoSettings.cc index f111a8bcbb63c..5022a345d2699 100644 --- a/L1Trigger/VertexFinder/src/AlgoSettings.cc +++ b/L1Trigger/VertexFinder/src/AlgoSettings.cc @@ -29,6 +29,8 @@ namespace l1tVertexFinder { vx_dbscan_mintracks_(vertex_.getParameter("DBSCANMinDensityTracks")), vx_kmeans_iterations_(vertex_.getParameter("KmeansIterations")), vx_kmeans_nclusters_(vertex_.getParameter("KmeansNumClusters")), + vx_trkw_graph_(vertex_.getParameter("TrackWeightGraph")), + vx_pattrec_graph_(vertex_.getParameter("PatternRecGraph")), // Debug printout debug_(iConfig.getParameter("debug")) { const std::string algoName(vertex_.getParameter("Algorithm")); @@ -64,7 +66,8 @@ namespace l1tVertexFinder { {"PVR", Algorithm::PVR}, {"adaptive", Algorithm::adaptiveVertexReconstruction}, {"HPV", Algorithm::HPV}, - {"K-means", Algorithm::Kmeans}}; + {"K-means", Algorithm::Kmeans}, + {"NNEmulation", Algorithm::NNEmulation}}; const std::map AlgoSettings::algoPrecisionMap = { {Algorithm::fastHisto, Precision::Simulation}, @@ -76,6 +79,7 @@ namespace l1tVertexFinder { {Algorithm::PVR, Precision::Simulation}, {Algorithm::adaptiveVertexReconstruction, Precision::Simulation}, {Algorithm::HPV, Precision::Simulation}, - {Algorithm::Kmeans, Precision::Simulation}}; + {Algorithm::Kmeans, Precision::Simulation}, + {Algorithm::NNEmulation, Precision::Emulation}}; } // end namespace l1tVertexFinder diff --git a/L1Trigger/VertexFinder/src/VertexFinder.cc b/L1Trigger/VertexFinder/src/VertexFinder.cc index d224bf00fd400..7778f4fe4fbfa 100644 --- a/L1Trigger/VertexFinder/src/VertexFinder.cc +++ b/L1Trigger/VertexFinder/src/VertexFinder.cc @@ -1099,4 +1099,126 @@ namespace l1tVertexFinder { pv_index_ = 0; } // end of fastHistoEmulation + void VertexFinder::NNVtxEmulation(tensorflow::Session* TrackWeightSesh, + tensorflow::Session* PatternRecSesh, + tensorflow::Session* AssociationSesh) { + // #### Weight Tracks: #### + // Loop over tracks -> weight the network -> set track weights + tensorflow::Tensor inputTrkWeight(tensorflow::DT_FLOAT, {1, 3}); //Single batch of 3 values + uint counter = 0; + + for (auto& track : fitTracks_) { + // Quantised Network: Use values from L1GTTInputProducer pT, MVA1, eta + auto& gttTrack = fitTracks_.at(counter); + + TTTrack_TrackWord::tanl_t etaEmulationBits = gttTrack.getTTTrackPtr()->getTanlWord(); + ap_fixed<16, 3> etaEmulation; + etaEmulation.V = (etaEmulationBits.range()); + + ap_uint<14> ptEmulationBits = gttTrack.getTTTrackPtr()->getTrackWord()( + TTTrack_TrackWord::TrackBitLocations::kRinvMSB - 1, TTTrack_TrackWord::TrackBitLocations::kRinvLSB); + ap_ufixed<14, 9> ptEmulation; + ptEmulation.V = (ptEmulationBits.range()); + + ap_ufixed<22, 9> ptEmulation_rescale; + ptEmulation_rescale = ptEmulation.to_double(); + + ap_ufixed<22, 9> etaEmulation_rescale; + etaEmulation_rescale = abs(etaEmulation.to_double()); + + ap_ufixed<22, 9> MVAEmulation_rescale; + MVAEmulation_rescale = gttTrack.getTTTrackPtr()->getMVAQualityBits(); + + inputTrkWeight.tensor()(0, 0) = ptEmulation_rescale.to_double(); + inputTrkWeight.tensor()(0, 1) = MVAEmulation_rescale.to_double(); + inputTrkWeight.tensor()(0, 2) = etaEmulation_rescale.to_double(); + + // CNN output: track weight + std::vector outputTrkWeight; + tensorflow::run(TrackWeightSesh, {{"weight:0", inputTrkWeight}}, {"Identity:0"}, &outputTrkWeight); + // Set track weight pack into tracks: + + ap_ufixed<16, 5> NNOutput; + NNOutput = (double)outputTrkWeight[0].tensor()(0, 0); + + //std::cout<<"NNOutput_weight_network = "<< NNOutput <vx_histogram_numbins(), 1}); //Single batch with 256 bins and 1 weight + std::vector outputPV; + RecoVertexCollection vertices(settings_->vx_histogram_numbins()); + std::map vertexMap; + std::map histogram; + std::map nnOutput; + + float binWidth = settings_->vx_histogram_binwidth(); + + int track_z = 0; + + for (int z = 0; z < settings_->vx_histogram_numbins(); z += 1) { + counter = 0; + double vxWeight = 0; + + for (const L1Track& track : fitTracks_) { + auto& gttTrack = fitTracks_.at(counter); + double temp_z0 = gttTrack.getTTTrackPtr()->z0(); + + track_z = std::floor((temp_z0 + settings_->vx_histogram_max()) / binWidth); + + if (track_z >= z && track_z < (z + 1)) { + vertices.at(z).insert(&track); + vxWeight += track.weight(); + } + ++counter; + } + // Get centre of bin before setting z0 + vertices.at(z).setZ0(((z + 0.5) * binWidth) - settings_->vx_histogram_max()); + + vertexMap[vxWeight] = z; + inputPV.tensor()(0, z, 0) = vxWeight; + //Fill histogram for 3 bin sliding window: + histogram[z] = vxWeight; + } + + // Run PV Network: + tensorflow::run(PatternRecSesh, {{"hist:0", inputPV}}, {"Identity:0"}, &outputPV); + // Threshold needed due to rounding differences in internal CNN layer emulation versus firmware + const float histogrammingThreshold_ = 0.0; + for (int i(0); i < settings_->vx_histogram_numbins(); ++i) { + if (outputPV[0].tensor()(0, i, 0) >= histogrammingThreshold_) { + nnOutput[i] = outputPV[0].tensor()(0, i, 0); + } + } + + //Find max then find all occurances of it in histogram and average their position -> python argmax layer + //Performance is not optimised for multiple peaks in histogram or spread peaks these are edge cases, need to revisit + int max_index = 0; + int num_maxes = 0; + float max_element = 0.0; + for (int i(0); i < settings_->vx_histogram_numbins(); ++i) { + if (nnOutput[i] > max_element) { + max_element = nnOutput[i]; + } + } + + for (int i(0); i < settings_->vx_histogram_numbins(); ++i) { + if (nnOutput[i] == max_element) { + num_maxes++; + max_index += i; + } + } + int PV_index = ceil((float)max_index / (float)num_maxes); + + edm::LogVerbatim("VertexFinder") << " NNVtxEmulation Chosen PV: prob: " << nnOutput[PV_index] + << " bin = " << PV_index << " z0 = " << vertices.at(PV_index).z0() << '\n'; + + verticesEmulation_.emplace_back(1, vertices.at(PV_index).z0(), 0, vertices.at(PV_index).pt(), 0, 0, 0); + + } // end of NNVtx Algorithm + } // namespace l1tVertexFinder