From dd4abdb4d40e46982de5180d50a82d12f9b2fd84 Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Mon, 9 Jan 2023 13:09:03 +0000 Subject: [PATCH 01/16] Update L1T sum branches in downstream modules. --- .../python/HLT_75e33/modules/l1tPFPuppiHT400offMaxEta2p4_cfi.py | 2 +- .../python/HLT_75e33/modules/l1tPFPuppiHT450off_cfi.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT400offMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT400offMaxEta2p4_cfi.py index 648da465f075f..d387a3640748b 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT400offMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT400offMaxEta2p4_cfi.py @@ -6,5 +6,5 @@ theScalings = cms.vdouble(50.0182, 1.0961, 0) ), TypeOfSum = cms.string('HT'), - inputTag = cms.InputTag("l1tPhase1JetSumsProducer","Sums") + inputTag = cms.InputTag("l1tPhase1JetSumsProducer9x9trimmed","Sums") ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT450off_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT450off_cfi.py index 090b3d17a160b..7d8cad4520804 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT450off_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT450off_cfi.py @@ -6,5 +6,5 @@ theScalings = cms.vdouble(50.0182, 1.0961, 0) ), TypeOfSum = cms.string('HT'), - inputTag = cms.InputTag("l1tPhase1JetSumsProducer","Sums") + inputTag = cms.InputTag("l1tPhase1JetSumsProducer9x9trimmed","Sums") ) From 13e52047b3d7649475e056d514df1a5be2d5521c Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Wed, 18 Jan 2023 21:02:10 +0000 Subject: [PATCH 02/16] Update L1T jets consumed by HLT modules from 9x9 to 9x9trimmed. --- .../HLT_75e33/modules/l1t1PFPuppiJet70offMaxEta2p4_cfi.py | 2 +- .../HLT_75e33/modules/l1t2PFPuppiJet55offMaxEta2p4_cfi.py | 2 +- .../HLT_75e33/modules/l1t4PFPuppiJet25OnlineMaxEta2p4_cfi.py | 2 +- .../HLT_75e33/modules/l1t4PFPuppiJet40offMaxEta2p4_cfi.py | 2 +- .../modules/l1tDoublePFPuppiJet112offMaxEta2p4_cfi.py | 2 +- .../modules/l1tDoublePFPuppiJets112offMaxDeta1p6_cfi.py | 4 ++-- .../python/HLT_75e33/modules/l1tPFPuppiHT_cfi.py | 2 +- .../python/HLT_75e33/modules/l1tSinglePFPuppiJet230off_cfi.py | 2 +- 8 files changed, 9 insertions(+), 9 deletions(-) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t1PFPuppiJet70offMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t1PFPuppiJet70offMaxEta2p4_cfi.py index 3eeda4877538e..af21009718979 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t1PFPuppiJet70offMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t1PFPuppiJet70offMaxEta2p4_cfi.py @@ -10,5 +10,5 @@ endcap = cms.vdouble(42.4039, 1.33052, 0), overlap = cms.vdouble(24.8375, 1.4152, 0) ), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates") + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates") ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t2PFPuppiJet55offMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t2PFPuppiJet55offMaxEta2p4_cfi.py index 9d17f392e98a7..4ce4aad79115b 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t2PFPuppiJet55offMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t2PFPuppiJet55offMaxEta2p4_cfi.py @@ -10,5 +10,5 @@ endcap = cms.vdouble(42.4039, 1.33052, 0), overlap = cms.vdouble(24.8375, 1.4152, 0) ), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates") + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates") ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet25OnlineMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet25OnlineMaxEta2p4_cfi.py index 77cf6d8138274..2c4f2b65f3f3f 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet25OnlineMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet25OnlineMaxEta2p4_cfi.py @@ -5,5 +5,5 @@ MinEta = cms.double(-2.4), MinN = cms.int32(4), MinPt = cms.double(25.0), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates") + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates") ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet40offMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet40offMaxEta2p4_cfi.py index 9e0c55b29708a..0ba10ab1e6e92 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet40offMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1t4PFPuppiJet40offMaxEta2p4_cfi.py @@ -10,5 +10,5 @@ endcap = cms.vdouble(42.4039, 1.33052, 0), overlap = cms.vdouble(24.8375, 1.4152, 0) ), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates") + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates") ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJet112offMaxEta2p4_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJet112offMaxEta2p4_cfi.py index 8bc28b6c44b02..22bef3d0647f0 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJet112offMaxEta2p4_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJet112offMaxEta2p4_cfi.py @@ -10,6 +10,6 @@ endcap = cms.vdouble(42.4039, 1.33052, 0), overlap = cms.vdouble(24.8375, 1.4152, 0) ), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates"), + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates"), saveTags = cms.bool(True) ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJets112offMaxDeta1p6_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJets112offMaxDeta1p6_cfi.py index f8c6ea593444d..79b52b0ae52a4 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJets112offMaxDeta1p6_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tDoublePFPuppiJets112offMaxDeta1p6_cfi.py @@ -14,8 +14,8 @@ MinPt = cms.double(0.0), inputTag1 = cms.InputTag("l1tDoublePFPuppiJet112offMaxEta2p4"), inputTag2 = cms.InputTag("l1tDoublePFPuppiJet112offMaxEta2p4"), - originTag1 = cms.VInputTag(cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates")), - originTag2 = cms.VInputTag(cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates")), + originTag1 = cms.VInputTag(cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates")), + originTag2 = cms.VInputTag(cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates")), saveTags = cms.bool(True), triggerType1 = cms.int32(-116), triggerType2 = cms.int32(-116) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT_cfi.py index 0f7a28e9ea4d9..7c37e7576eab5 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tPFPuppiHT_cfi.py @@ -1,7 +1,7 @@ import FWCore.ParameterSet.Config as cms l1tPFPuppiHT = cms.EDProducer("HLTHtMhtProducer", - jetsLabel = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates"), + jetsLabel = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates"), maxEtaJetHt = cms.double(2.4), minPtJetHt = cms.double(30.0) ) diff --git a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tSinglePFPuppiJet230off_cfi.py b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tSinglePFPuppiJet230off_cfi.py index 27463576e938a..25694aea4859c 100644 --- a/HLTrigger/Configuration/python/HLT_75e33/modules/l1tSinglePFPuppiJet230off_cfi.py +++ b/HLTrigger/Configuration/python/HLT_75e33/modules/l1tSinglePFPuppiJet230off_cfi.py @@ -10,5 +10,5 @@ endcap = cms.vdouble(42.4039, 1.33052, 0), overlap = cms.vdouble(24.8375, 1.4152, 0) ), - inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9","Phase1L1TJetFromPfCandidates") + inputTag = cms.InputTag("l1tPhase1JetCalibrator9x9trimmed","Phase1L1TJetFromPfCandidates") ) From bb2ce17f41e6e2dad2be4bc577a6a6961ecc52ee Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Mon, 6 Feb 2023 21:13:45 +0100 Subject: [PATCH 03/16] Add CTL1 producers using extended tracks to task. --- L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py b/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py index ee6fc05d8726a..733a56803047f 100644 --- a/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py +++ b/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py @@ -537,7 +537,9 @@ L1TLayer1Task = cms.Task( l1tLayer1Barrel, + l1tLayer1BarrelExtended, l1tLayer1HGCal, + l1tLayer1HGCalExtended, l1tLayer1HGCalNoTK, l1tLayer1HF, l1tLayer1, From eb6a0b9f92b4a81b33bfd57c1628d6bc8db7a1df Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Mon, 6 Feb 2023 21:14:15 +0100 Subject: [PATCH 04/16] Fix name of extended tracks input. --- .../Phase2L1ParticleFlow/python/l1tPFTracksFromL1Tracks_cfi.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/L1Trigger/Phase2L1ParticleFlow/python/l1tPFTracksFromL1Tracks_cfi.py b/L1Trigger/Phase2L1ParticleFlow/python/l1tPFTracksFromL1Tracks_cfi.py index 35f44e8d32e2e..7ebf7d8236d1c 100644 --- a/L1Trigger/Phase2L1ParticleFlow/python/l1tPFTracksFromL1Tracks_cfi.py +++ b/L1Trigger/Phase2L1ParticleFlow/python/l1tPFTracksFromL1Tracks_cfi.py @@ -26,7 +26,7 @@ ) l1tPFTracksFromL1TracksExtended = l1tPFTracksFromL1Tracks.clone( - L1TrackTag = cms.InputTag("TTTracksFromExtendedTrackletEmulation", "Level1TTTracks"), + L1TrackTag = cms.InputTag("l1tTTTracksFromExtendedTrackletEmulation", "Level1TTTracks"), nParam = 5, qualityBits = cms.vstring( "momentum.perp > 2 && getStubRefs.size >= 4 && chi2Red < 15 && POCA.x < 1.0 && POCA.x > -1.0 && POCA.y < 1.0 && POCA.y > -1.0", From 1270c84cad03ea8582b4cb746674001e79e9bae6 Mon Sep 17 00:00:00 2001 From: Gianluca Date: Fri, 10 Mar 2023 12:15:31 +0100 Subject: [PATCH 05/16] Add Phase2 Correlator Layer2 e/g objects to the event content --- L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py | 1 + 1 file changed, 1 insertion(+) diff --git a/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py b/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py index 559144d3dfaf0..c412f98d4f517 100644 --- a/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py +++ b/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py @@ -211,6 +211,7 @@ def _appendPhase2Digis(obj): 'keep *_l1tLayer1HF_*_*', 'keep *_l1tLayer1_*_*', 'keep *_l1tLayer1EG_*_*', + 'keep *_l1tLayer2EG_*_*', 'keep *_l1tMETPFProducer_*_*', 'keep *_l1tNNTauProducer_*_*', 'keep *_l1tNNTauProducerPuppi_*_*', From 4acc8e50bb4701434de3736b96e9bd5f3aa3168a Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Thu, 30 Mar 2023 15:36:37 +0100 Subject: [PATCH 06/16] Add SC jets from extended tracks to event content. --- L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py | 5 +++++ 1 file changed, 5 insertions(+) diff --git a/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py b/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py index c412f98d4f517..6399f341f7984 100644 --- a/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py +++ b/L1Trigger/Configuration/python/L1Trigger_EventContent_cff.py @@ -199,6 +199,8 @@ def _appendPhase2Digis(obj): 'keep *_l1tPFTracksFromL1TracksHGCal_*_*', 'keep *_l1tSCPFL1PuppiCorrectedEmulator_*_*', 'keep *_l1tSCPFL1PuppiCorrectedEmulatorMHT_*_*', + 'keep *_l1tSCPFL1PuppiExtendedCorrectedEmulator_*_*', + 'keep *_l1tSCPFL1PuppiExtendedCorrectedEmulatorMHT_*_*', 'keep *_l1tPhase1JetProducer9x9_*_*', 'keep *_l1tPhase1JetCalibrator9x9_*_*', 'keep *_l1tPhase1JetSumsProducer9x9_*_*', @@ -210,6 +212,9 @@ def _appendPhase2Digis(obj): 'keep *_l1tLayer1HGCalNoTK_*_*', 'keep *_l1tLayer1HF_*_*', 'keep *_l1tLayer1_*_*', + 'keep *_l1tLayer1BarrelExtended_*_*', + 'keep *_l1tLayer1HGCalExtended_*_*', + 'keep *_l1tLayer1Extended_*_*', 'keep *_l1tLayer1EG_*_*', 'keep *_l1tLayer2EG_*_*', 'keep *_l1tMETPFProducer_*_*', From a9b5a0d299493dc7fe1634289aee73e4f27c45af Mon Sep 17 00:00:00 2001 From: Zhenbin Wu Date: Mon, 20 Mar 2023 10:40:10 -0500 Subject: [PATCH 07/16] Add ap_data type for GT interface Update GMT emulator for new ap format Also clean up the config file Update Z0/D0 according to GT Using 0.05 for Z0 and 0.03 for d0. update code format Bug fix in tracker muon matching Thanks to Santi for spotting this. This line was commented out by someone, likely by accident. Update the reference type to SA muons --- .../L1TMuonPhase2/interface/Constants.h | 46 +++++++++++++++---- DataFormats/L1TMuonPhase2/interface/SAMuon.h | 10 ++++ .../L1TMuonPhase2/interface/TrackerMuon.h | 18 ++++++-- DataFormats/L1TMuonPhase2/src/classes_def.xml | 4 +- .../plugins/Phase2L1TGMTSAMuonProducer.cc | 14 +++--- .../plugins/TrackMuonMatchAlgorithm.h | 6 +-- L1Trigger/Phase2L1GMT/python/gmt_cfi.py | 2 +- L1Trigger/Phase2L1GMT/test/runGMT.py | 20 -------- 8 files changed, 76 insertions(+), 44 deletions(-) diff --git a/DataFormats/L1TMuonPhase2/interface/Constants.h b/DataFormats/L1TMuonPhase2/interface/Constants.h index d3f80a3421948..ea8c2f4c0c5cb 100644 --- a/DataFormats/L1TMuonPhase2/interface/Constants.h +++ b/DataFormats/L1TMuonPhase2/interface/Constants.h @@ -54,7 +54,7 @@ namespace Phase2L1GMT { // Bitwidth for standalone muons to CL1 and GT const int BITSSAZ0 = 5; const int BITSSAD0 = 7; - const int BITSSAQUALITY = 4; + const int BITSSAQUAL = 4; // Bitwidth for dataformat to GT const int BITSGTPT = 16; @@ -62,14 +62,24 @@ namespace Phase2L1GMT { const int BITSGTETA = 14; const int BITSGTZ0 = 10; const int BITSGTD0 = 10; - const int BITSGTQUALITY = 8; + const int BITSGTQUAL = 8; const int BITSGTISO = 4; + const int BITSGTBETA = 4; + + // Bitwidth for Tau->3mu object + const int BITSTMPT = 8; + const int BITSTMPHI = 8; + const int BITSTMETA = 8; + const int BITSTMMASS2 = 8; + const int BITSTMTYPE = 6; + const int BITSTMIDX = 4; + const int BITSTMQUAL = 4; const float maxCurv_ = 0.00855; // 2 GeV pT Rinv is in cm const float maxPhi_ = 1.026; // relative to the center of the sector const float maxTanl_ = 8.0; - const float maxZ0_ = 30.; - const float maxD0_ = 15.4; + const float maxZ0_ = 25.6; + const float maxD0_ = 15.36; // Updated barrelLimit according to Karol, https://indico.cern.ch/event/1113802/#1-phase2-gmt-performance-and-i const int barrelLimit0_ = 1.4 / 0.00076699039 / 8; const int barrelLimit1_ = 1.1 / 0.00076699039 / 8; @@ -81,12 +91,32 @@ namespace Phase2L1GMT { const float LSBpt = 0.03125; const float LSBphi = 2. * M_PI / pow(2, BITSPHI); const float LSBeta = 2. * M_PI / pow(2, BITSETA); - const float LSBGTz0 = 2. * maxZ0_ / pow(2, BITSZ0); - const float LSBGTd0 = 2. * maxD0_ / pow(2, BITSD0); - const float LSBSAz0 = 1.875; - const float LSBSAd0 = 3.85; + const float LSBGTz0 = 0.05; // 0.5mm, in sync with GTT and Correlator + const float LSBGTd0 = 0.03; // from GT interface doc + const float LSBSAz0 = 1.6; // 0.05 * 32 cm, with range +- 25.6 + const float LSBSAd0 = 3.84; // 0.03 * 128 cm, with range +- 245.76 typedef ap_uint<64> wordtype; + typedef ap_uint<1> valid_gt_t; //valid + typedef ap_uint<1> q_gt_t; //charge + typedef ap_uint pt_gt_t; //pt of tracker muon + typedef ap_int phi_gt_t; //phi of tracker muon + typedef ap_int eta_gt_t; //eta of tracker muon + typedef ap_int z0_gt_t; //z0 of tracker muon + typedef ap_int d0_gt_t; //d0 of tracker muon + typedef ap_uint iso_gt_t; //isolation of tracker muon + typedef ap_uint beta_gt_t; //beta of tracker muon + typedef ap_uint qual_gt_t; //quality of tracker muon + + //Standalone muon datatype + typedef ap_uint<1> valid_sa_t; //valid + typedef ap_uint pt_sa_t; //pt of standalone muon + typedef ap_int phi_sa_t; //phi of standalone muon + typedef ap_int eta_sa_t; //eta of standalone muon + typedef ap_int z0_sa_t; //z0 of standalone muon + typedef ap_int d0_sa_t; //d0 of standalone muon + typedef ap_uint<1> q_sa_t; //charge of standalone muon + typedef ap_uint qual_sa_t; //quality of standalone muon inline uint64_t twos_complement(long long int v, uint bits) { uint64_t mask = (1 << bits) - 1; diff --git a/DataFormats/L1TMuonPhase2/interface/SAMuon.h b/DataFormats/L1TMuonPhase2/interface/SAMuon.h index d46d6cab3610c..12ea802dae04f 100644 --- a/DataFormats/L1TMuonPhase2/interface/SAMuon.h +++ b/DataFormats/L1TMuonPhase2/interface/SAMuon.h @@ -31,6 +31,16 @@ namespace l1t { const uint hwBeta() const { return hwBeta_; } void setBeta(uint beta) { hwBeta_ = beta; } + // For GT, returning ap_ type + const Phase2L1GMT::valid_sa_t apValid() const { return Phase2L1GMT::valid_sa_t(hwPt() > 0); }; + const Phase2L1GMT::pt_sa_t apPt() const { return Phase2L1GMT::pt_sa_t(hwPt()); }; + const Phase2L1GMT::phi_sa_t apPhi() const { return Phase2L1GMT::phi_sa_t(hwPhi()); }; + const Phase2L1GMT::eta_sa_t apEta() const { return Phase2L1GMT::eta_sa_t(hwEta()); }; + const Phase2L1GMT::z0_sa_t apZ0() const { return Phase2L1GMT::z0_sa_t(hwZ0()); }; + const Phase2L1GMT::d0_sa_t apD0() const { return Phase2L1GMT::d0_sa_t(hwD0()); }; + const Phase2L1GMT::q_sa_t apCharge() const { return Phase2L1GMT::q_sa_t(hwCharge()); }; + const Phase2L1GMT::qual_sa_t apQual() const { return Phase2L1GMT::qual_sa_t(hwQual()); }; + // For HLT const double phZ0() const { return Phase2L1GMT::LSBSAz0 * hwZ0(); } const double phD0() const { return Phase2L1GMT::LSBSAd0 * hwD0(); } diff --git a/DataFormats/L1TMuonPhase2/interface/TrackerMuon.h b/DataFormats/L1TMuonPhase2/interface/TrackerMuon.h index 78bb288685d37..1dd8cc90bf344 100644 --- a/DataFormats/L1TMuonPhase2/interface/TrackerMuon.h +++ b/DataFormats/L1TMuonPhase2/interface/TrackerMuon.h @@ -32,7 +32,7 @@ namespace l1t { ~TrackerMuon() override; const edm::Ptr& trkPtr() const { return trkPtr_; } - const edm::Ref& muonRef() const { return muRef_; } + const std::vector& muonRef() const { return muRef_; } const bool hwCharge() const { return hwCharge_; } const int hwZ0() const { return hwZ0_; } @@ -41,10 +41,22 @@ namespace l1t { const int hwIsoSumAp() const { return hwIsoSumAp_; } const uint hwBeta() const { return hwBeta_; } void setBeta(uint beta) { hwBeta_ = beta; } - void setMuonRef(const edm::Ref& p) { muRef_ = p; } + void setMuonRef(const std::vector& p) { muRef_ = p; } void setHwIsoSum(int isoSum) { hwIsoSum_ = isoSum; } void setHwIsoSumAp(int isoSum) { hwIsoSumAp_ = isoSum; } + // For GT, returning ap_ type + const Phase2L1GMT::valid_gt_t apValid() const { return Phase2L1GMT::valid_gt_t(hwPt() > 0); }; + const Phase2L1GMT::pt_gt_t apPt() const { return Phase2L1GMT::pt_gt_t(hwPt()); }; + const Phase2L1GMT::phi_gt_t apPhi() const { return Phase2L1GMT::phi_gt_t(hwPhi()); }; + const Phase2L1GMT::eta_gt_t apEta() const { return Phase2L1GMT::eta_gt_t(hwEta()); }; + const Phase2L1GMT::z0_gt_t apZ0() const { return Phase2L1GMT::z0_gt_t(hwZ0()); }; + const Phase2L1GMT::d0_gt_t apD0() const { return Phase2L1GMT::d0_gt_t(hwD0()); }; + const Phase2L1GMT::q_gt_t apCharge() const { return Phase2L1GMT::q_gt_t(hwCharge()); }; + const Phase2L1GMT::qual_gt_t apQual() const { return Phase2L1GMT::qual_gt_t(hwQual()); }; + const Phase2L1GMT::iso_gt_t apIso() const { return Phase2L1GMT::iso_gt_t(hwIso()); }; + const Phase2L1GMT::beta_gt_t apBeta() const { return Phase2L1GMT::beta_gt_t(hwBeta()); }; + // For HLT const double phZ0() const { return Phase2L1GMT::LSBGTz0 * hwZ0(); } const double phD0() const { return Phase2L1GMT::LSBGTd0 * hwD0(); } @@ -76,7 +88,7 @@ namespace l1t { //Store the eneryg sum for isolation with ap_type int hwIsoSumAp_; - edm::Ref muRef_; + std::vector muRef_; MuonStubRefVector stubs_; }; } // namespace l1t diff --git a/DataFormats/L1TMuonPhase2/src/classes_def.xml b/DataFormats/L1TMuonPhase2/src/classes_def.xml index b20feacd4c091..d847c7d613ec8 100644 --- a/DataFormats/L1TMuonPhase2/src/classes_def.xml +++ b/DataFormats/L1TMuonPhase2/src/classes_def.xml @@ -8,8 +8,8 @@ - - + + diff --git a/L1Trigger/Phase2L1GMT/plugins/Phase2L1TGMTSAMuonProducer.cc b/L1Trigger/Phase2L1GMT/plugins/Phase2L1TGMTSAMuonProducer.cc index 5f64b43b10ccc..c6284873729c5 100644 --- a/L1Trigger/Phase2L1GMT/plugins/Phase2L1TGMTSAMuonProducer.cc +++ b/L1Trigger/Phase2L1GMT/plugins/Phase2L1TGMTSAMuonProducer.cc @@ -127,17 +127,17 @@ void Phase2L1TGMTSAMuonProducer::produce(edm::Event& iEvent, const edm::EventSet // Description: // =========================================================================== SAMuon Phase2L1TGMTSAMuonProducer::Convertl1tMuon(const l1t::Muon& mu, const int bx_) { - ap_uint qual = mu.hwQual(); + qual_sa_t qual = mu.hwQual(); int charge = mu.charge() > 0 ? 0 : 1; - ap_uint pt = round(mu.pt() / LSBpt); - ap_int phi = round(mu.phi() / LSBphi); - ap_int eta = round(mu.eta() / LSBeta); + pt_sa_t pt = round(mu.pt() / LSBpt); + phi_sa_t phi = round(mu.phi() / LSBphi); + eta_sa_t eta = round(mu.eta() / LSBeta); // FIXME: Below are not well defined in phase1 GMT // Using the version from Correlator for now - ap_int z0 = 0; // No tracks info in Phase 1 + z0_sa_t z0 = 0; // No tracks info in Phase 1 // Use 2 bits with LSB = 30cm for BMTF and 25cm for EMTF currently, but subjet to change - ap_int d0 = mu.hwDXY(); + d0_sa_t d0 = mu.hwDXY(); int bstart = 0; wordtype word(0); @@ -148,7 +148,7 @@ SAMuon Phase2L1TGMTSAMuonProducer::Convertl1tMuon(const l1t::Muon& mu, const int bstart = wordconcat(word, bstart, z0, BITSSAZ0); bstart = wordconcat(word, bstart, d0, BITSSAD0); bstart = wordconcat(word, bstart, charge, 1); - bstart = wordconcat(word, bstart, qual, BITSSAQUALITY); + bstart = wordconcat(word, bstart, qual, BITSSAQUAL); SAMuon samuon(mu, charge, pt.to_uint(), eta.to_int(), phi.to_int(), z0.to_int(), d0.to_int(), qual.to_uint()); samuon.setWord(word); diff --git a/L1Trigger/Phase2L1GMT/plugins/TrackMuonMatchAlgorithm.h b/L1Trigger/Phase2L1GMT/plugins/TrackMuonMatchAlgorithm.h index 3df98ee6ed63e..30f3ad6d8a3f5 100644 --- a/L1Trigger/Phase2L1GMT/plugins/TrackMuonMatchAlgorithm.h +++ b/L1Trigger/Phase2L1GMT/plugins/TrackMuonMatchAlgorithm.h @@ -105,7 +105,7 @@ namespace Phase2L1GMT { if (out.size() == maximum) break; l1t::TrackerMuon muon(mu.trkPtr(), mu.charge(), mu.pt(), mu.eta(), mu.phi(), mu.z0(), mu.d0(), mu.quality()); - //muon.setMuonRef(mu.muonRef()); + muon.setMuonRef(mu.muonRef()); for (const auto& stub : mu.stubs()) muon.addStub(stub); out.push_back(muon); @@ -131,9 +131,9 @@ namespace Phase2L1GMT { bstart = 0; bstart = wordconcat(word2, bstart, mu.hwCharge(), 1); - bstart = wordconcat(word2, bstart, mu.hwQual(), BITSGTQUALITY); + bstart = wordconcat(word2, bstart, mu.hwQual(), BITSGTQUAL); bstart = wordconcat(word2, bstart, mu.hwIso(), BITSGTISO); - bstart = wordconcat(word2, bstart, mu.hwBeta(), BITSMUONBETA); + bstart = wordconcat(word2, bstart, mu.hwBeta(), BITSGTBETA); std::array wordout = {{word1, word2}}; mu.setWord(wordout); diff --git a/L1Trigger/Phase2L1GMT/python/gmt_cfi.py b/L1Trigger/Phase2L1GMT/python/gmt_cfi.py index 6e8bcde812de1..e5d24cf7bd955 100644 --- a/L1Trigger/Phase2L1GMT/python/gmt_cfi.py +++ b/L1Trigger/Phase2L1GMT/python/gmt_cfi.py @@ -70,7 +70,7 @@ RelIsoThresholdM = cms.double(0.05), RelIsoThresholdT = cms.double(0.01), verbose = cms.int32(0), - IsodumpForHLS = cms.int32(1), + IsodumpForHLS = cms.int32(0), ), tauto3mu = cms.PSet() diff --git a/L1Trigger/Phase2L1GMT/test/runGMT.py b/L1Trigger/Phase2L1GMT/test/runGMT.py index 4a870ad535b43..26c0aec35ae64 100644 --- a/L1Trigger/Phase2L1GMT/test/runGMT.py +++ b/L1Trigger/Phase2L1GMT/test/runGMT.py @@ -115,30 +115,10 @@ process.dtTriggerPhase2PrimitiveDigis.scenario = 0 process.load("L1Trigger.Phase2L1GMT.gmt_cff") -process.l1tGMTMuons.trackMatching.verbose=1 -process.l1tGMTMuons.verbose=0 -process.l1tGMTMuons.trackConverter.verbose=0 - - - -#process.schedule = cms.Schedule(process.L1TrackTrigger_step,process.pL1TMuonTPS,process.endjob_step,process.e) # Adding MuonTPS - # Path and EndPath definitions process.raw2digi_step = cms.Path(process.RawToDigi) process.L1TrackTrigger_step = cms.Path(process.L1TrackTrigger) -#process.pL1TkPhotonsCrystal = cms.Path(process.L1TkPhotonsCrystal) -#process.pL1TkIsoElectronsCrystal = cms.Path(process.L1TkIsoElectronsCrystal) -#process.pL1TkElectronsLooseCrystal = cms.Path(process.L1TkElectronsLooseCrystal) -#process.pL1TkElectronsLooseHGC = cms.Path(process.L1TkElectronsLooseHGC) -#process.pL1TkElectronsHGC = cms.Path(process.L1TkElectronsHGC) -process.pL1TkMuon = cms.Path(process.L1TkMuons+process.L1TkMuonsTP) -#process.pL1TkElectronsEllipticMatchHGC = cms.Path(process.L1TkElectronsEllipticMatchHGC) -#process.pL1TkElectronsCrystal = cms.Path(process.L1TkElectronsCrystal) -#process.pL1TkPhotonsHGC = cms.Path(process.L1TkPhotonsHGC) -#process.pL1TkIsoElectronsHGC = cms.Path(process.L1TkIsoElectronsHGC) -#process.pL1TkElectronsEllipticMatchCrystal = cms.Path(process.L1TkElectronsEllipticMatchCrystal) -#process.L1simulation_step = cms.Path(process.SimL1Emulator) process.testpath=cms.Path(process.CalibratedDigis*process.dtTriggerPhase2PrimitiveDigis*process.phase2GMT) process.endjob_step = cms.EndPath(process.endOfProcess) process.FEVTDEBUGHLToutput_step = cms.EndPath(process.FEVTDEBUGHLToutput) From 50aab0099bb41d2be7a8086ab0118e0db3139d2e Mon Sep 17 00:00:00 2001 From: Andrew David Loeliger Date: Thu, 6 Apr 2023 05:47:40 +0200 Subject: [PATCH 08/16] Add back removed data version --- DataFormats/L1TMuonPhase2/src/classes_def.xml | 1 + 1 file changed, 1 insertion(+) diff --git a/DataFormats/L1TMuonPhase2/src/classes_def.xml b/DataFormats/L1TMuonPhase2/src/classes_def.xml index d847c7d613ec8..b3cbf00fcf430 100644 --- a/DataFormats/L1TMuonPhase2/src/classes_def.xml +++ b/DataFormats/L1TMuonPhase2/src/classes_def.xml @@ -10,6 +10,7 @@ + From d6e23fdbed4902eb8a681a3bdefec6e4cafd1c3f Mon Sep 17 00:00:00 2001 From: Sioni Summers <14807534+thesps@users.noreply.github.com> Date: Fri, 31 Mar 2023 13:41:44 +0200 Subject: [PATCH 09/16] Provide method to access GT and CT HW jet objects from l1t::PFJet, update consumers to pick the appropriate HW object where needed Address review comments Fix HT emulator input collection and output encoding --- DataFormats/L1TParticleFlow/interface/PFJet.h | 19 ++++++++++++++++--- DataFormats/L1TParticleFlow/interface/sums.h | 2 +- .../L1TParticleFlow/src/classes_def.xml | 3 ++- .../plugins/L1MHtPFProducer.cc | 13 +++++++------ .../plugins/L1SeedConePFJetProducer.cc | 3 ++- .../python/l1pfJetMet_cff.py | 2 +- 6 files changed, 29 insertions(+), 13 deletions(-) diff --git a/DataFormats/L1TParticleFlow/interface/PFJet.h b/DataFormats/L1TParticleFlow/interface/PFJet.h index e1cfc282659c0..5dd73fdc0bdec 100644 --- a/DataFormats/L1TParticleFlow/interface/PFJet.h +++ b/DataFormats/L1TParticleFlow/interface/PFJet.h @@ -5,6 +5,8 @@ #include "DataFormats/L1Trigger/interface/L1Candidate.h" #include "DataFormats/L1TParticleFlow/interface/PFCandidate.h" #include "DataFormats/Common/interface/Ptr.h" +#include "DataFormats/L1TParticleFlow/interface/jets.h" +#include "DataFormats/L1TParticleFlow/interface/gt_datatypes.h" namespace l1t { @@ -38,13 +40,24 @@ namespace l1t { edm::Ptr daughterPtr(size_type i) const { return constituents_[i]; } // Get and set the encodedJet_ bits. The Jet is encoded in 128 bits as a 2-element array of uint64_t - std::array encodedJet() { return encodedJet_; } - void setEncodedJet(std::array jet) { encodedJet_ = jet; } + // We store encodings both for Correlator internal usage and for Global Trigger + enum class HWEncoding { CT, GT }; + typedef std::array PackedJet; + const PackedJet& encodedJet(const HWEncoding encoding = HWEncoding::GT) const { + return encodedJet_[static_cast(encoding)]; + } + void setEncodedJet(const HWEncoding encoding, const PackedJet jet) { + encodedJet_[static_cast(encoding)] = jet; + } + + // Accessors to HW objects with ap_* types from encoded words + l1gt::Jet getHWJetGT() const { return l1gt::Jet::unpack(encodedJet(HWEncoding::GT)); } + l1ct::Jet getHWJetCT() const { return l1ct::Jet::unpack(encodedJet(HWEncoding::CT)); } private: float rawPt_; Constituents constituents_; - std::array encodedJet_ = {{0, 0}}; + std::array encodedJet_ = {{{{0, 0}}, {{0, 0}}}}; }; typedef std::vector PFJetCollection; diff --git a/DataFormats/L1TParticleFlow/interface/sums.h b/DataFormats/L1TParticleFlow/interface/sums.h index 195566f61afbe..b91aeec3b93d6 100644 --- a/DataFormats/L1TParticleFlow/interface/sums.h +++ b/DataFormats/L1TParticleFlow/interface/sums.h @@ -52,7 +52,7 @@ namespace l1ct { sum.valid = (hwPt != 0) || (hwSumPt != 0); sum.vector_pt = CTtoGT_pt(hwPt); sum.vector_phi = CTtoGT_phi(hwPhi); - sum.scalar_pt = CTtoGT_phi(hwSumPt); + sum.scalar_pt = CTtoGT_pt(hwSumPt); return sum; } }; diff --git a/DataFormats/L1TParticleFlow/src/classes_def.xml b/DataFormats/L1TParticleFlow/src/classes_def.xml index 5f90ab44bb63d..eaf2e1bff560b 100644 --- a/DataFormats/L1TParticleFlow/src/classes_def.xml +++ b/DataFormats/L1TParticleFlow/src/classes_def.xml @@ -33,7 +33,8 @@ - + + diff --git a/L1Trigger/Phase2L1ParticleFlow/plugins/L1MHtPFProducer.cc b/L1Trigger/Phase2L1ParticleFlow/plugins/L1MHtPFProducer.cc index 741b7a84e550c..1a970da6e8ad0 100644 --- a/L1Trigger/Phase2L1ParticleFlow/plugins/L1MHtPFProducer.cc +++ b/L1Trigger/Phase2L1ParticleFlow/plugins/L1MHtPFProducer.cc @@ -65,7 +65,7 @@ void L1MhtPfProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::Even std::vector L1MhtPfProducer::convertEDMToHW(std::vector edmJets) const { std::vector hwJets; std::for_each(edmJets.begin(), edmJets.end(), [&](l1t::PFJet jet) { - l1ct::Jet hwJet = l1ct::Jet::unpack(jet.encodedJet()); + l1ct::Jet hwJet = jet.getHWJetCT(); hwJets.push_back(hwJet); }); return hwJets; @@ -75,17 +75,18 @@ std::vector L1MhtPfProducer::convertHWToEDM(l1ct::Sum hwSums) const std::vector edmSums; reco::Candidate::PolarLorentzVector htVector; - htVector.SetPt(hwSums.hwSumPt.to_double()); + l1gt::Sum gtSum = hwSums.toGT(); + htVector.SetPt(l1gt::Scales::floatPt(gtSum.scalar_pt)); htVector.SetPhi(0); htVector.SetEta(0); reco::Candidate::PolarLorentzVector mhtVector; - mhtVector.SetPt(hwSums.hwPt.to_double()); - mhtVector.SetPhi(l1ct::Scales::floatPhi(hwSums.hwPhi)); + mhtVector.SetPt(l1gt::Scales::floatPt(gtSum.vector_pt)); + mhtVector.SetPhi(l1gt::Scales::floatPhi(gtSum.vector_phi)); mhtVector.SetEta(0); - l1t::EtSum ht(htVector, l1t::EtSum::EtSumType::kTotalHt); - l1t::EtSum mht(mhtVector, l1t::EtSum::EtSumType::kMissingHt); + l1t::EtSum ht(htVector, l1t::EtSum::EtSumType::kTotalHt, gtSum.scalar_pt.bits_to_uint64()); + l1t::EtSum mht(mhtVector, l1t::EtSum::EtSumType::kMissingHt, gtSum.vector_pt.bits_to_uint64(), 0, gtSum.vector_phi); edmSums.push_back(ht); edmSums.push_back(mht); diff --git a/L1Trigger/Phase2L1ParticleFlow/plugins/L1SeedConePFJetProducer.cc b/L1Trigger/Phase2L1ParticleFlow/plugins/L1SeedConePFJetProducer.cc index 4f64649d48103..5cf255385f433 100644 --- a/L1Trigger/Phase2L1ParticleFlow/plugins/L1SeedConePFJetProducer.cc +++ b/L1Trigger/Phase2L1ParticleFlow/plugins/L1SeedConePFJetProducer.cc @@ -201,7 +201,8 @@ std::vector L1SeedConePFJetProducer::convertHWToEDM( gtJet.v3.pt.V, gtJet.v3.eta.V, gtJet.v3.phi.V); - edmJet.setEncodedJet(jet.toGT().pack()); + edmJet.setEncodedJet(l1t::PFJet::HWEncoding::CT, jet.pack()); + edmJet.setEncodedJet(l1t::PFJet::HWEncoding::GT, jet.toGT().pack()); // get back the references to the constituents std::vector> constituents; std::for_each(jet.constituents.begin(), jet.constituents.end(), [&](auto constituent) { diff --git a/L1Trigger/Phase2L1ParticleFlow/python/l1pfJetMet_cff.py b/L1Trigger/Phase2L1ParticleFlow/python/l1pfJetMet_cff.py index a795b5f590952..60d58d2ac23dd 100644 --- a/L1Trigger/Phase2L1ParticleFlow/python/l1pfJetMet_cff.py +++ b/L1Trigger/Phase2L1ParticleFlow/python/l1pfJetMet_cff.py @@ -25,7 +25,7 @@ phase2_hgcalV11.toModify(_correctedJets, correctorFile = "L1Trigger/Phase2L1ParticleFlow/data/jecs/jecs_20220308.root") from L1Trigger.Phase2L1ParticleFlow.l1tMHTPFProducer_cfi import l1tMHTPFProducer -l1tSCPFL1PuppiCorrectedEmulatorMHT = l1tMHTPFProducer.clone() +l1tSCPFL1PuppiCorrectedEmulatorMHT = l1tMHTPFProducer.clone(jets = 'l1tSCPFL1PuppiCorrectedEmulator') L1TPFJetsTask = cms.Task( l1tLayer2Deregionizer, l1tSCPFL1PF, l1tSCPFL1Puppi, l1tSCPFL1PuppiEmulator, l1tSCPFL1PuppiCorrectedEmulator, l1tSCPFL1PuppiCorrectedEmulatorMHT From 7b66861abbddbe5291e71cf58f9ede266aaebf1d Mon Sep 17 00:00:00 2001 From: Gianluca Date: Thu, 13 Apr 2023 17:25:00 +0200 Subject: [PATCH 10/16] fix isolation type for GT objects --- DataFormats/L1TParticleFlow/interface/gt_datatypes.h | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/DataFormats/L1TParticleFlow/interface/gt_datatypes.h b/DataFormats/L1TParticleFlow/interface/gt_datatypes.h index 12cec406492f7..178e332dbedf4 100644 --- a/DataFormats/L1TParticleFlow/interface/gt_datatypes.h +++ b/DataFormats/L1TParticleFlow/interface/gt_datatypes.h @@ -30,7 +30,7 @@ namespace l1gt { typedef ap_uint<1> valid_t; // E/gamma fields - typedef ap_fixed<11, 9> iso_t; + typedef ap_ufixed<11, 9> iso_t; typedef ap_uint<4> egquality_t; // tau fields From 49c8205cd99327f38d9703246b3e34c7cd43cbc1 Mon Sep 17 00:00:00 2001 From: Georgios Karathanasis Date: Wed, 3 May 2023 14:30:44 +0200 Subject: [PATCH 11/16] L1track jet emulator update (backport of PR40563) --- DataFormats/L1Trigger/interface/TkJetWord.h | 77 +- DataFormats/L1Trigger/src/TkJetWord.cc | 35 +- DataFormats/L1Trigger/src/classes_def.xml | 3 +- .../src/classes_def_others.xml | 1 + .../Configuration/python/SimL1Emulator_cff.py | 4 +- .../DemonstratorTools/src/codecs_tkjets.cc | 20 +- .../L1TTrackMatch/interface/L1Clustering.h | 252 ------- .../interface/L1TrackJetEmulationProducer.h | 79 -- .../plugins/L1TrackJetClustering.h | 365 ++++++++++ .../plugins/L1TrackJetEmulationProducer.cc | 681 ------------------ .../plugins/L1TrackJetEmulatorProducer.cc | 456 ++++++++++++ .../plugins/L1TrackJetProducer.cc | 412 ++++++----- .../python/l1tTrackJetsEmulation_cfi.py | 65 +- .../L1TTrackMatch/python/l1tTrackJets_cfi.py | 18 +- .../test/L1TrackObjectNtupleMaker_cfg.py | 33 +- .../test/L1TrackNtupleMaker_cfg.py | 83 ++- 16 files changed, 1240 insertions(+), 1344 deletions(-) delete mode 100644 L1Trigger/L1TTrackMatch/interface/L1Clustering.h delete mode 100644 L1Trigger/L1TTrackMatch/interface/L1TrackJetEmulationProducer.h create mode 100644 L1Trigger/L1TTrackMatch/plugins/L1TrackJetClustering.h delete mode 100644 L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulationProducer.cc create mode 100644 L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulatorProducer.cc diff --git a/DataFormats/L1Trigger/interface/TkJetWord.h b/DataFormats/L1Trigger/interface/TkJetWord.h index a4f12dd57b627..8f45482c4e355 100644 --- a/DataFormats/L1Trigger/interface/TkJetWord.h +++ b/DataFormats/L1Trigger/interface/TkJetWord.h @@ -1,8 +1,10 @@ + #ifndef FIRMWARE_TkJetWord_h #define FIRMWARE_TkJetWord_h -// Class to store the 96-bit TkJet word for L1 Track Trigger. +// Class to store the 128-bit TkJet word for L1 Track Trigger. // Author: Emily MacDonald, updated by Benjamin Radburn-Smith (September 2022) +// 2nd update: George Karathanasis (Oct 2022) #include #include @@ -10,17 +12,15 @@ #include #include #include +#include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h" namespace l1t { class TkJetWord { public: // ----------constants, enums and typedefs --------- - int INTPHI_PI = 720; - int INTPHI_TWOPI = 2 * INTPHI_PI; - float INTPT_LSB = 1 >> 5; - float ETAPHI_LSB = M_PI / (1 << 12); - float Z0_LSB = 0.05; + static constexpr double MAX_Z0 = 30.; + static constexpr double MAX_ETA = 8.; enum TkJetBitWidths { kPtSize = 16, @@ -30,33 +30,38 @@ namespace l1t { kZ0Size = 10, kNtSize = 5, kXtSize = 4, - kUnassignedSize = 8, - kTkJetWordSize = kPtSize + kGlbEtaSize + kGlbPhiSize + kZ0Size + kNtSize + kXtSize + kUnassignedSize, + kDispFlagSize = 1, + kUnassignedSize = 65, + kTkJetWordSize = + kPtSize + kGlbEtaSize + kGlbPhiSize + kZ0Size + kNtSize + kXtSize + kDispFlagSize + kUnassignedSize, }; enum TkJetBitLocations { kPtLSB = 0, kPtMSB = kPtLSB + TkJetBitWidths::kPtSize - 1, - kGlbEtaLSB = kPtMSB + 1, - kGlbEtaMSB = kGlbEtaLSB + TkJetBitWidths::kGlbEtaSize - 1, - kGlbPhiLSB = kGlbEtaMSB + 1, + kGlbPhiLSB = kPtMSB + 1, kGlbPhiMSB = kGlbPhiLSB + TkJetBitWidths::kGlbPhiSize - 1, - kZ0LSB = kGlbPhiMSB + 1, + kGlbEtaLSB = kGlbPhiMSB + 1, + kGlbEtaMSB = kGlbEtaLSB + TkJetBitWidths::kGlbEtaSize - 1, + kZ0LSB = kGlbEtaMSB + 1, kZ0MSB = kZ0LSB + TkJetBitWidths::kZ0Size - 1, kNtLSB = kZ0MSB + 1, kNtMSB = kNtLSB + TkJetBitWidths::kNtSize - 1, kXtLSB = kNtMSB + 1, kXtMSB = kXtLSB + TkJetBitWidths::kXtSize - 1, - kUnassignedLSB = kXtMSB + 1, + kDispFlagLSB = kXtMSB + 1, + kDispFlagMSB = kDispFlagLSB + TkJetBitWidths::kDispFlagSize - 1, + kUnassignedLSB = kDispFlagMSB + 1, kUnassignedMSB = kUnassignedLSB + TkJetBitWidths::kUnassignedSize - 1, }; typedef ap_ufixed pt_t; typedef ap_int glbeta_t; typedef ap_int glbphi_t; - typedef ap_int z0_t; // 40cm / 0.1 - typedef ap_uint nt_t; //number of tracks - typedef ap_uint nx_t; //number of tracks with xbit = 1 + typedef ap_int z0_t; // 40cm / 0.1 + typedef ap_uint nt_t; //number of tracks + typedef ap_uint nx_t; //number of tracks with xbit = 1 + typedef ap_uint dispflag_t; typedef ap_uint tkjetunassigned_t; // Unassigned bits typedef std::bitset tkjetword_bs_t; typedef ap_uint tkjetword_t; @@ -64,7 +69,14 @@ namespace l1t { public: // ----------Constructors -------------------------- TkJetWord() {} - TkJetWord(pt_t pt, glbeta_t eta, glbphi_t phi, z0_t z0, nt_t nt, nx_t nx, tkjetunassigned_t unassigned); + TkJetWord(pt_t pt, + glbeta_t eta, + glbphi_t phi, + z0_t z0, + nt_t nt, + nx_t nx, + dispflag_t dispflag, + tkjetunassigned_t unassigned); ~TkJetWord() {} @@ -109,6 +121,12 @@ namespace l1t { ret.V = tkJetWord()(TkJetBitLocations::kXtMSB, TkJetBitLocations::kXtLSB); return ret; } + dispflag_t dispFlagWord() const { + dispflag_t ret; + ret.V = tkJetWord()(TkJetBitLocations::kDispFlagMSB, TkJetBitLocations::kDispFlagLSB); + return ret; + } + tkjetunassigned_t unassignedWord() const { return tkJetWord()(TkJetBitLocations::kUnassignedMSB, TkJetBitLocations::kUnassignedLSB); } @@ -122,20 +140,37 @@ namespace l1t { unsigned int z0Bits() const { return z0Word().to_uint(); } unsigned int ntBits() const { return ntWord().to_uint(); } unsigned int xtBits() const { return xtWord().to_uint(); } + unsigned int dispFlagBits() const { return dispFlagWord().to_uint(); } unsigned int unassignedBits() const { return unassignedWord().to_uint(); } // These functions return the unpacked and converted values // These functions return real numbers converted from the digitized quantities by unpacking the 64-bit vertex word float pt() const { return ptWord().to_float(); } - float glbeta() const { return glbEtaWord().to_float() * ETAPHI_LSB; } - float glbphi() const { return glbPhiWord().to_float() * ETAPHI_LSB; } - float z0() const { return z0Word().to_float() * Z0_LSB; } + float glbeta() const { + return unpackSignedValue( + glbEtaWord(), TkJetBitWidths::kGlbEtaSize, (MAX_ETA) / (1 << TkJetBitWidths::kGlbEtaSize)); + } + float glbphi() const { + return unpackSignedValue( + glbPhiWord(), TkJetBitWidths::kGlbPhiSize, (2. * std::abs(M_PI)) / (1 << TkJetBitWidths::kGlbPhiSize)); + } + float z0() const { + return unpackSignedValue(z0Word(), TkJetBitWidths::kZ0Size, MAX_Z0 / (1 << TkJetBitWidths::kZ0Size)); + } int nt() const { return (ap_ufixed(ntWord())).to_int(); } int xt() const { return (ap_ufixed(xtWord())).to_int(); } + int dispflag() const { return (ap_ufixed(dispFlagWord())).to_int(); } unsigned int unassigned() const { return unassignedWord().to_uint(); } // ----------member functions (setters) ------------ - void setTkJetWord(pt_t pt, glbeta_t eta, glbphi_t phi, z0_t z0, nt_t nt, nx_t nx, tkjetunassigned_t unassigned); + void setTkJetWord(pt_t pt, + glbeta_t eta, + glbphi_t phi, + z0_t z0, + nt_t nt, + nx_t nx, + dispflag_t dispflag, + tkjetunassigned_t unassigned); private: // ----------private member functions -------------- diff --git a/DataFormats/L1Trigger/src/TkJetWord.cc b/DataFormats/L1Trigger/src/TkJetWord.cc index 2464c4a2db941..87897de7b07b0 100644 --- a/DataFormats/L1Trigger/src/TkJetWord.cc +++ b/DataFormats/L1Trigger/src/TkJetWord.cc @@ -4,26 +4,39 @@ #include "DataFormats/L1Trigger/interface/TkJetWord.h" namespace l1t { - TkJetWord::TkJetWord(pt_t pt, glbeta_t eta, glbphi_t phi, z0_t z0, nt_t nt, nx_t nx, tkjetunassigned_t unassigned) { - setTkJetWord(pt, eta, phi, z0, nt, nx, unassigned); + TkJetWord::TkJetWord(pt_t pt, + glbeta_t eta, + glbphi_t phi, + z0_t z0, + nt_t nt, + nx_t nx, + dispflag_t dispflag, + tkjetunassigned_t unassigned) { + setTkJetWord(pt, eta, phi, z0, nt, nx, dispflag, unassigned); } - void TkJetWord::setTkJetWord( - pt_t pt, glbeta_t eta, glbphi_t phi, z0_t z0, nt_t nt, nx_t nx, tkjetunassigned_t unassigned) { + void TkJetWord::setTkJetWord(pt_t pt, + glbeta_t eta, + glbphi_t phi, + z0_t z0, + nt_t nt, + nx_t nx, + dispflag_t dispflag, + tkjetunassigned_t unassigned) { // pack the TkJet word unsigned int offset = 0; for (unsigned int b = offset; b < (offset + TkJetBitWidths::kPtSize); b++) { tkJetWord_.set(b, pt[b - offset]); } offset += TkJetBitWidths::kPtSize; - for (unsigned int b = offset; b < (offset + TkJetBitWidths::kGlbEtaSize); b++) { - tkJetWord_.set(b, eta[b - offset]); - } - offset += TkJetBitWidths::kGlbEtaSize; for (unsigned int b = offset; b < (offset + TkJetBitWidths::kGlbPhiSize); b++) { tkJetWord_.set(b, phi[b - offset]); } offset += TkJetBitWidths::kGlbPhiSize; + for (unsigned int b = offset; b < (offset + TkJetBitWidths::kGlbEtaSize); b++) { + tkJetWord_.set(b, eta[b - offset]); + } + offset += TkJetBitWidths::kGlbEtaSize; for (unsigned int b = offset; b < (offset + TkJetBitWidths::kZ0Size); b++) { tkJetWord_.set(b, z0[b - offset]); } @@ -36,9 +49,13 @@ namespace l1t { tkJetWord_.set(b, nx[b - offset]); } offset += TkJetBitWidths::kXtSize; + for (unsigned int b = offset; b < (offset + TkJetBitWidths::kDispFlagSize); b++) { + tkJetWord_.set(b, nx[b - offset]); + } + offset += TkJetBitWidths::kDispFlagSize; for (unsigned int b = offset; b < (offset + TkJetBitWidths::kUnassignedSize); b++) { tkJetWord_.set(b, unassigned[b - offset]); } } -} //namespace l1t \ No newline at end of file +} //namespace l1t diff --git a/DataFormats/L1Trigger/src/classes_def.xml b/DataFormats/L1Trigger/src/classes_def.xml index 6b4f9444aaf04..b57276c8ebd03 100644 --- a/DataFormats/L1Trigger/src/classes_def.xml +++ b/DataFormats/L1Trigger/src/classes_def.xml @@ -304,7 +304,8 @@ - + + diff --git a/DataFormats/StdDictionaries/src/classes_def_others.xml b/DataFormats/StdDictionaries/src/classes_def_others.xml index e2c769fea5a1d..962d89e945c9b 100644 --- a/DataFormats/StdDictionaries/src/classes_def_others.xml +++ b/DataFormats/StdDictionaries/src/classes_def_others.xml @@ -12,6 +12,7 @@ + diff --git a/L1Trigger/Configuration/python/SimL1Emulator_cff.py b/L1Trigger/Configuration/python/SimL1Emulator_cff.py index 17600bf2e63a9..5dc7ba09569f6 100644 --- a/L1Trigger/Configuration/python/SimL1Emulator_cff.py +++ b/L1Trigger/Configuration/python/SimL1Emulator_cff.py @@ -151,9 +151,9 @@ from L1Trigger.L1TTrackMatch.l1tTrackerEtMiss_cfi import * from L1Trigger.L1TTrackMatch.l1tTrackerHTMiss_cfi import * # make the input tags consistent with the choice L1VertexFinder above -l1tTrackJets.L1PVertexCollection = ("l1tVertexFinder", "l1vertices") +l1tTrackJets.L1PVertexInputTag = ("l1tVertexFinderEmulator","l1verticesEmulation") l1tTrackFastJets.L1PrimaryVertexTag = ("l1tVertexFinder", "l1vertices") -l1tTrackJetsExtended.L1PVertexCollection = ("l1tVertexFinder", "l1vertices") +l1tTrackJetsExtended.L1PVertexInputTag = ("l1tVertexFinderEmulator","l1verticesEmulation") #l1tTrackerEtMiss.L1VertexInputTag = ("l1tVertexFinder", "l1vertices") #l1tTrackerEtMissExtended.L1VertexInputTag = ("l1tVertexFinder", "l1vertices") diff --git a/L1Trigger/DemonstratorTools/src/codecs_tkjets.cc b/L1Trigger/DemonstratorTools/src/codecs_tkjets.cc index 75de4a5c52e0b..069f23f5107da 100644 --- a/L1Trigger/DemonstratorTools/src/codecs_tkjets.cc +++ b/L1Trigger/DemonstratorTools/src/codecs_tkjets.cc @@ -34,13 +34,19 @@ namespace l1t::demo::codecs { //if (not x.test(TkJetWord::kValidLSB)) // break; - TkJetWord j(TkJetWord::pt_t(frames[f](TkJetWord::kPtMSB, TkJetWord::kPtLSB)), - TkJetWord::glbeta_t(frames[f](TkJetWord::kGlbEtaMSB, TkJetWord::kGlbEtaLSB)), - TkJetWord::glbphi_t(frames[f](TkJetWord::kGlbPhiMSB, TkJetWord::kGlbPhiLSB)), - TkJetWord::z0_t(frames[f](TkJetWord::kZ0MSB, TkJetWord::kZ0LSB)), - TkJetWord::nt_t(frames[f](TkJetWord::kNtMSB, TkJetWord::kNtLSB)), - TkJetWord::nx_t(frames[f](TkJetWord::kXtMSB, TkJetWord::kXtLSB)), - TkJetWord::tkjetunassigned_t(frames[f](TkJetWord::kUnassignedMSB, TkJetWord::kUnassignedLSB))); + TkJetWord j( + TkJetWord::pt_t(frames[f](TkJetWord::TkJetBitLocations::kPtMSB, TkJetWord::TkJetBitLocations::kPtLSB)), + TkJetWord::glbphi_t( + frames[f](TkJetWord::TkJetBitLocations::kGlbPhiMSB, TkJetWord::TkJetBitLocations::kGlbPhiLSB)), + TkJetWord::glbeta_t( + frames[f](TkJetWord::TkJetBitLocations::kGlbEtaMSB, TkJetWord::TkJetBitLocations::kGlbEtaLSB)), + TkJetWord::z0_t(frames[f](TkJetWord::TkJetBitLocations::kZ0MSB, TkJetWord::TkJetBitLocations::kZ0LSB)), + TkJetWord::nt_t(frames[f](TkJetWord::TkJetBitLocations::kNtMSB, TkJetWord::TkJetBitLocations::kNtLSB)), + TkJetWord::nx_t(frames[f](TkJetWord::TkJetBitLocations::kXtMSB, TkJetWord::TkJetBitLocations::kXtLSB)), + TkJetWord::dispflag_t( + frames[f](TkJetWord::TkJetBitLocations::kDispFlagMSB, TkJetWord::TkJetBitLocations::kDispFlagLSB)), + TkJetWord::tkjetunassigned_t( + frames[f](TkJetWord::TkJetBitLocations::kUnassignedMSB, TkJetWord::TkJetBitLocations::kUnassignedLSB))); tkJets.push_back(j); } diff --git a/L1Trigger/L1TTrackMatch/interface/L1Clustering.h b/L1Trigger/L1TTrackMatch/interface/L1Clustering.h deleted file mode 100644 index f97a11e8431ae..0000000000000 --- a/L1Trigger/L1TTrackMatch/interface/L1Clustering.h +++ /dev/null @@ -1,252 +0,0 @@ -#ifndef L1Trigger_L1TTrackMatch_L1Clustering_HH -#define L1Trigger_L1TTrackMatch_L1Clustering_HH -#include -#include -#include -#include - -//Each individual box in the eta and phi dimension. -// Also used to store final cluster data for each zbin. -struct EtaPhiBin { - float pTtot = 0; - int numtracks = 0; - int numttrks = 0; - int numtdtrks = 0; - int numttdtrks = 0; - bool used = false; - float phi; //average phi value (halfway b/t min and max) - float eta; //average eta value - std::vector trackidx; -}; - -//store important information for plots -struct MaxZBin { - int znum; //Numbered from 0 to nzbins (16, 32, or 64) in order - int nclust; //number of clusters in this bin - float zbincenter; - std::vector clusters; //list of all the clusters in this bin - float ht; //sum of all cluster pTs--only the zbin with the maximum ht is stored -}; - -inline std::vector L1_clustering(EtaPhiBin *phislice, int etaBins_, float etaStep_) { - std::vector clusters; - // Find eta-phibin with maxpT, make center of cluster, add neighbors if not already used - int nclust = 0; - - // get tracks in eta bins in increasing eta order - for (int etabin = 0; etabin < etaBins_; ++etabin) { - float my_pt = 0, previousbin_pt = 0; //, nextbin_pt=0, next2bin_pt=0; - float nextbin_pt = 0, nextbin2_pt = 0; - - // skip (already) used tracks - if (phislice[etabin].used) - continue; - my_pt = phislice[etabin].pTtot; - if (my_pt == 0) - continue; - //get previous bin pT - if (etabin > 0 && !phislice[etabin - 1].used) - previousbin_pt = phislice[etabin - 1].pTtot; - - // get next bins pt - if (etabin < etaBins_ - 1 && !phislice[etabin + 1].used) { - nextbin_pt = phislice[etabin + 1].pTtot; - if (etabin < etaBins_ - 2 && !phislice[etabin + 2].used) { - nextbin2_pt = phislice[etabin + 2].pTtot; - } - } - // check if pT of current cluster is higher than neighbors - if (my_pt < previousbin_pt || my_pt <= nextbin_pt) { - // if unused pT in the left neighbor, spit it out as a cluster - if (previousbin_pt > 0) { - clusters.push_back(phislice[etabin - 1]); - phislice[etabin - 1].used = true; - nclust++; - } - continue; //if it is not the local max pT skip - } - // here reach only unused local max clusters - clusters.push_back(phislice[etabin]); - phislice[etabin].used = true; //if current bin a cluster - if (previousbin_pt > 0) { - clusters[nclust].pTtot += previousbin_pt; - clusters[nclust].numtracks += phislice[etabin - 1].numtracks; - clusters[nclust].numtdtrks += phislice[etabin - 1].numtdtrks; - for (unsigned int itrk = 0; itrk < phislice[etabin - 1].trackidx.size(); itrk++) - clusters[nclust].trackidx.push_back(phislice[etabin - 1].trackidx[itrk]); - } - - if (my_pt >= nextbin2_pt && nextbin_pt > 0) { - clusters[nclust].pTtot += nextbin_pt; - clusters[nclust].numtracks += phislice[etabin + 1].numtracks; - clusters[nclust].numtdtrks += phislice[etabin + 1].numtdtrks; - for (unsigned int itrk = 0; itrk < phislice[etabin + 1].trackidx.size(); itrk++) - clusters[nclust].trackidx.push_back(phislice[etabin + 1].trackidx[itrk]); - phislice[etabin + 1].used = true; - } - - nclust++; - - } // for each etabin - - // Merge close-by clusters - for (int m = 0; m < nclust - 1; ++m) { - if (std::abs(clusters[m + 1].eta - clusters[m].eta) < 1.5 * etaStep_) { - if (clusters[m + 1].pTtot > clusters[m].pTtot) { - clusters[m].eta = clusters[m + 1].eta; - } - clusters[m].pTtot += clusters[m + 1].pTtot; - clusters[m].numtracks += clusters[m + 1].numtracks; // total ntrk - clusters[m].numtdtrks += clusters[m + 1].numtdtrks; // total ndisp - for (unsigned int itrk = 0; itrk < clusters[m + 1].trackidx.size(); itrk++) - clusters[m].trackidx.push_back(clusters[m + 1].trackidx[itrk]); - - // if remove the merged cluster - all the others must be closer to 0 - for (int m1 = m + 1; m1 < nclust - 1; ++m1) { - clusters[m1] = clusters[m1 + 1]; - //clusters.erase(clusters.begin()+m1); - } - // clusters[m1] = clusters[m1 + 1]; - clusters.erase(clusters.begin() + nclust); - nclust--; - m = -1; - } // end if clusters neighbor in eta - } // end for (m) loop - - return clusters; -} - -inline void Fill_L2Cluster(EtaPhiBin &bin, float pt, int ntrk, int ndtrk, std::vector trkidx) { - bin.pTtot += pt; - bin.numtracks += ntrk; - bin.numtdtrks += ndtrk; - for (unsigned int itrk = 0; itrk < trkidx.size(); itrk++) - bin.trackidx.push_back(trkidx[itrk]); -} - -inline float DPhi(float phi1, float phi2) { - float x = phi1 - phi2; - if (x >= M_PI) - x -= 2 * M_PI; - if (x < -1 * M_PI) - x += 2 * M_PI; - return x; -} - -inline std::vector L2_clustering(std::vector> &L1clusters, - int phiBins_, - float phiStep_, - float etaStep_) { - std::vector clusters; - for (int phibin = 0; phibin < phiBins_; ++phibin) { //Find eta-phibin with highest pT - if (L1clusters[phibin].empty()) - continue; - - // sort L1 clusters max -> min - sort(L1clusters[phibin].begin(), L1clusters[phibin].end(), [](struct EtaPhiBin &a, struct EtaPhiBin &b) { - return a.pTtot > b.pTtot; - }); - for (unsigned int imax = 0; imax < L1clusters[phibin].size(); ++imax) { - if (L1clusters[phibin][imax].used) - continue; - float pt_current = L1clusters[phibin][imax].pTtot; //current cluster (pt0) - float pt_next = 0; // next phi bin (pt1) - float pt_next2 = 0; // next to next phi bin2 (pt2) - int trk1 = 0; - int trk2 = 0; - int tdtrk1 = 0; - int tdtrk2 = 0; - std::vector trkidx1; - std::vector trkidx2; - clusters.push_back(L1clusters[phibin][imax]); - - L1clusters[phibin][imax].used = true; - - // if we are in the last phi bin, dont check phi+1 phi+2 - if (phibin == phiBins_ - 1) - continue; - std::vector used_already; //keep phi+1 clusters that have been used - for (unsigned int icluster = 0; icluster < L1clusters[phibin + 1].size(); ++icluster) { - if (L1clusters[phibin + 1][icluster].used) - continue; - if (std::abs(L1clusters[phibin + 1][icluster].eta - L1clusters[phibin][imax].eta) > 1.5 * etaStep_) - continue; - pt_next += L1clusters[phibin + 1][icluster].pTtot; - trk1 += L1clusters[phibin + 1][icluster].numtracks; - tdtrk1 += L1clusters[phibin + 1][icluster].numtdtrks; - for (unsigned int itrk = 0; itrk < L1clusters[phibin + 1][icluster].trackidx.size(); itrk++) - trkidx1.push_back(L1clusters[phibin + 1][icluster].trackidx[itrk]); - used_already.push_back(icluster); - } - - if (pt_next < pt_current) { // if pt1 used_already2; //keep used clusters in phi+2 - for (unsigned int icluster = 0; icluster < L1clusters[phibin + 2].size(); ++icluster) { - if (L1clusters[phibin + 2][icluster].used) - continue; - if (std::abs(L1clusters[phibin + 2][icluster].eta - L1clusters[phibin][imax].eta) > 1.5 * etaStep_) - continue; - pt_next2 += L1clusters[phibin + 2][icluster].pTtot; - trk2 += L1clusters[phibin + 2][icluster].numtracks; - tdtrk2 += L1clusters[phibin + 2][icluster].numtdtrks; - for (unsigned int itrk = 0; itrk < L1clusters[phibin + 2][icluster].trackidx.size(); itrk++) - trkidx2.push_back(L1clusters[phibin + 2][icluster].trackidx[itrk]); - used_already2.push_back(icluster); - } - if (pt_next2 < pt_next) { - std::vector trkidx_both; - trkidx_both.reserve(trkidx1.size() + trkidx2.size()); - trkidx_both.insert(trkidx_both.end(), trkidx1.begin(), trkidx1.end()); - trkidx_both.insert(trkidx_both.end(), trkidx2.begin(), trkidx2.end()); - Fill_L2Cluster(clusters[clusters.size() - 1], pt_next + pt_next2, trk1 + trk2, tdtrk1 + tdtrk2, trkidx_both); - clusters[clusters.size() - 1].phi = L1clusters[phibin + 1][used_already[0]].phi; - for (unsigned int iused : used_already) - L1clusters[phibin + 1][iused].used = true; - for (unsigned int iused : used_already2) - L1clusters[phibin + 2][iused].used = true; - } - } // for each L1 cluster - } // for each phibin - - int nclust = clusters.size(); - - // merge close-by clusters - for (int m = 0; m < nclust - 1; ++m) { - for (int n = m + 1; n < nclust; ++n) { - if (clusters[n].eta != clusters[m].eta) - continue; - if (std::abs(DPhi(clusters[n].phi, clusters[m].phi)) > 1.5 * phiStep_) - continue; - - if (clusters[n].pTtot > clusters[m].pTtot) - clusters[m].phi = clusters[n].phi; - - clusters[m].pTtot += clusters[n].pTtot; - clusters[m].numtracks += clusters[n].numtracks; - clusters[m].numtdtrks += clusters[n].numtdtrks; - for (unsigned int itrk = 0; itrk < clusters[n].trackidx.size(); itrk++) - clusters[m].trackidx.push_back(clusters[n].trackidx[itrk]); - for (int m1 = n; m1 < nclust - 1; ++m1) - clusters[m1] = clusters[m1 + 1]; - clusters.erase(clusters.begin() + nclust); - - nclust--; - m = -1; - } // end of n-loop - } // end of m-loop - return clusters; -} -#endif diff --git a/L1Trigger/L1TTrackMatch/interface/L1TrackJetEmulationProducer.h b/L1Trigger/L1TTrackMatch/interface/L1TrackJetEmulationProducer.h deleted file mode 100644 index dd46370d13af4..0000000000000 --- a/L1Trigger/L1TTrackMatch/interface/L1TrackJetEmulationProducer.h +++ /dev/null @@ -1,79 +0,0 @@ -#ifndef L1Trigger_L1TTrackMatch_L1TrackJetEmulationProducer_HH -#define L1Trigger_L1TTrackMatch_L1TrackJetEmulationProducer_HH -#include -#include -#include -#include -#include -#include -#include "DataFormats/L1Trigger/interface/TkJetWord.h" - -//For precision studies -const int PT_EXTRABITS = 0; -const int ETA_EXTRABITS = 0; -const int PHI_EXTRABITS = 0; -const int Z0_EXTRABITS = 0; - -typedef ap_ufixed<16 + PT_EXTRABITS, 11, AP_TRN, AP_SAT> pt_intern; -typedef ap_int<14 + ETA_EXTRABITS> glbeta_intern; -typedef ap_int<14 + PHI_EXTRABITS> glbphi_intern; -typedef ap_int<10 + Z0_EXTRABITS> z0_intern; // 40cm / 0.1 - -namespace convert { - const int INTPHI_PI = 720; - const int INTPHI_TWOPI = 2 * INTPHI_PI; - static const float INTPT_LSB_POW = pow(2.0, -5 - PT_EXTRABITS); - static const float INTPT_LSB = INTPT_LSB_POW; - static const float ETA_LSB_POW = pow(2.0, -1 * ETA_EXTRABITS); - static const float ETA_LSB = M_PI / pow(2.0, 12) * ETA_LSB_POW; - static const float PHI_LSB_POW = pow(2.0, -1 * PHI_EXTRABITS); - static const float PHI_LSB = M_PI / pow(2.0, 12) * PHI_LSB_POW; - static const float Z0_LSB_POW = pow(2.0, -1 * Z0_EXTRABITS); - static const float Z0_LSB = 0.05 * Z0_LSB_POW; - inline float floatPt(pt_intern pt) { return pt.to_float(); } - inline int intPt(pt_intern pt) { return (ap_ufixed<18 + PT_EXTRABITS, 13 + PT_EXTRABITS>(pt)).to_int(); } - inline float floatEta(glbeta_intern eta) { return eta.to_float() * ETA_LSB; } - inline float floatPhi(glbphi_intern phi) { return phi.to_float() * PHI_LSB; } - inline float floatZ0(z0_intern z0) { return z0.to_float() * Z0_LSB; } - - inline pt_intern makePt(int pt) { return ap_ufixed<18 + PT_EXTRABITS, 13 + PT_EXTRABITS>(pt); } - inline pt_intern makePtFromFloat(float pt) { return pt_intern(INTPT_LSB_POW * round(pt / INTPT_LSB_POW)); } - inline z0_intern makeZ0(float z0) { return z0_intern(round(z0 / Z0_LSB)); } - - inline ap_uint ptToInt(pt_intern pt) { - // note: this can be synthethized, e.g. when pT is used as intex in a LUT - ap_uint ret = 0; - ret(pt_intern::width - 1, 0) = pt(pt_intern::width - 1, 0); - return ret; - } - - inline glbeta_intern makeGlbEta(float eta) { return round(eta / ETA_LSB); } - inline glbeta_intern makeGlbEtaRoundEven(float eta) { - glbeta_intern ghweta = round(eta / ETA_LSB); - return (ghweta % 2) ? glbeta_intern(ghweta + 1) : ghweta; - } - - inline glbphi_intern makeGlbPhi(float phi) { return round(phi / PHI_LSB); } - -}; // namespace convert - -//Each individual box in the eta and phi dimension. -// Also used to store final cluster data for each zbin. -struct TrackJetEmulationEtaPhiBin { - pt_intern pTtot; - l1t::TkJetWord::nt_t ntracks; - l1t::TkJetWord::nx_t nxtracks; - bool used; - glbphi_intern phi; //average phi value (halfway b/t min and max) - glbeta_intern eta; //average eta value -}; - -//store important information for plots -struct TrackJetEmulationMaxZBin { - int znum; //Numbered from 0 to nzbins (16, 32, or 64) in order - int nclust; //number of clusters in this bin - z0_intern zbincenter; - TrackJetEmulationEtaPhiBin *clusters; //list of all the clusters in this bin - pt_intern ht; //sum of all cluster pTs--only the zbin with the maximum ht is stored -}; -#endif diff --git a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetClustering.h b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetClustering.h new file mode 100644 index 0000000000000..c4fdf3dfd1c02 --- /dev/null +++ b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetClustering.h @@ -0,0 +1,365 @@ +#ifndef L1Trigger_L1TTrackMatch_L1TrackJetClustering_HH +#define L1Trigger_L1TTrackMatch_L1TrackJetClustering_HH +#include +#include +#include +#include +#include +#include +#include "DataFormats/L1Trigger/interface/TkJetWord.h" +#include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h" +#include "DataFormats/L1TrackTrigger/interface/TTTypes.h" +#include "FWCore/MessageLogger/interface/MessageLogger.h" + +namespace l1ttrackjet { + //For precision studies + const unsigned int PT_INTPART_BITS{9}; + const unsigned int ETA_INTPART_BITS{3}; + const unsigned int kExtraGlobalPhiBit{4}; + + typedef ap_ufixed pt_intern; + typedef ap_fixed glbeta_intern; + typedef ap_int glbphi_intern; + typedef ap_int z0_intern; // 40cm / 0.1 + typedef ap_uint d0_intern; + + inline const unsigned int DoubleToBit(double value, unsigned int maxBits, double step) { + unsigned int digitized_value = std::floor(std::abs(value) / step); + unsigned int digitized_maximum = (1 << (maxBits - 1)) - 1; // The remove 1 bit from nBits to account for the sign + if (digitized_value > digitized_maximum) + digitized_value = digitized_maximum; + if (value < 0) + digitized_value = (1 << maxBits) - digitized_value; // two's complement encoding + return digitized_value; + } + inline const double BitToDouble(unsigned int bits, unsigned int maxBits, double step) { + int isign = 1; + unsigned int digitized_maximum = (1 << maxBits) - 1; + if (bits & (1 << (maxBits - 1))) { // check the sign + isign = -1; + bits = (1 << (maxBits + 1)) - bits; + } + return (double(bits & digitized_maximum) + 0.5) * step * isign; + } + + // eta/phi clusters - simulation + struct EtaPhiBin { + float pTtot; + int ntracks; + int nxtracks; + bool used; + float phi; //average phi value (halfway b/t min and max) + float eta; //average eta value + std::vector trackidx; + }; + // z bin struct - simulation (used if z bin are many) + struct MaxZBin { + int znum; //Numbered from 0 to nzbins (16, 32, or 64) in order + int nclust; //number of clusters in this bin + float zbincenter; + std::vector clusters; //list of all the clusters in this bin + float ht; //sum of all cluster pTs--only the zbin with the maximum ht is stored + }; + + // eta/phi clusters - emulation + struct TrackJetEmulationEtaPhiBin { + pt_intern pTtot; + l1t::TkJetWord::nt_t ntracks; + l1t::TkJetWord::nx_t nxtracks; + bool used; + glbphi_intern phi; //average phi value (halfway b/t min and max) + glbeta_intern eta; //average eta value + std::vector trackidx; + }; + + // z bin struct - emulation (used if z bin are many) + struct TrackJetEmulationMaxZBin { + int znum; //Numbered from 0 to nzbins (16, 32, or 64) in order + int nclust; //number of clusters in this bin + z0_intern zbincenter; + std::vector clusters; //list of all the clusters in this bin + pt_intern ht; //sum of all cluster pTs--only the zbin with the maximum ht is stored + }; + + // track quality cuts + inline bool TrackQualitySelection(int trk_nstub, + double trk_chi2, + double trk_bendchi2, + double nStubs4PromptBend_, + double nStubs5PromptBend_, + double nStubs4PromptChi2_, + double nStubs5PromptChi2_, + double nStubs4DisplacedBend_, + double nStubs5DisplacedBend_, + double nStubs4DisplacedChi2_, + double nStubs5DisplacedChi2_, + bool displaced_) { + bool PassQuality = false; + if (!displaced_) { + if (trk_nstub == 4 && trk_bendchi2 < nStubs4PromptBend_ && + trk_chi2 < nStubs4PromptChi2_) // 4 stubs are the lowest track quality and have different cuts + PassQuality = true; + if (trk_nstub > 4 && trk_bendchi2 < nStubs5PromptBend_ && + trk_chi2 < nStubs5PromptChi2_) // above 4 stubs diffent selection imposed (genrally looser) + PassQuality = true; + } else { + if (trk_nstub == 4 && trk_bendchi2 < nStubs4DisplacedBend_ && + trk_chi2 < nStubs4DisplacedChi2_) // 4 stubs are the lowest track quality and have different cuts + PassQuality = true; + if (trk_nstub > 4 && trk_bendchi2 < nStubs5DisplacedBend_ && + trk_chi2 < nStubs5DisplacedChi2_) // above 4 stubs diffent selection imposed (genrally looser) + PassQuality = true; + } + return PassQuality; + } + + // L1 clustering (in eta) + template + inline std::vector L1_clustering(T *phislice, int etaBins_, Eta etaStep_) { + std::vector clusters; + // Find eta bin with maxpT, make center of cluster, add neighbors if not already used + int nclust = 0; + + // get tracks in eta bins in increasing eta order + for (int etabin = 0; etabin < etaBins_; ++etabin) { + Pt my_pt = 0; + Pt previousbin_pt = 0; + Pt nextbin_pt = 0; + Pt nextbin2_pt = 0; + + // skip (already) used tracks + if (phislice[etabin].used) + continue; + + my_pt = phislice[etabin].pTtot; + if (my_pt == 0) + continue; + + //get previous bin pT + if (etabin > 0 && !phislice[etabin - 1].used) + previousbin_pt = phislice[etabin - 1].pTtot; + + // get next bins pt + if (etabin < etaBins_ - 1 && !phislice[etabin + 1].used) { + nextbin_pt = phislice[etabin + 1].pTtot; + if (etabin < etaBins_ - 2 && !phislice[etabin + 2].used) { + nextbin2_pt = phislice[etabin + 2].pTtot; + } + } + // check if pT of current cluster is higher than neighbors + if (my_pt < previousbin_pt || my_pt <= nextbin_pt) { + // if unused pT in the left neighbor, spit it out as a cluster + if (previousbin_pt > 0) { + clusters.push_back(phislice[etabin - 1]); + phislice[etabin - 1].used = true; + nclust++; + } + continue; //if it is not the local max pT skip + } + // here reach only unused local max clusters + clusters.push_back(phislice[etabin]); + phislice[etabin].used = true; //if current bin a cluster + if (previousbin_pt > 0) { + clusters[nclust].pTtot += previousbin_pt; + clusters[nclust].ntracks += phislice[etabin - 1].ntracks; + clusters[nclust].nxtracks += phislice[etabin - 1].nxtracks; + for (unsigned int itrk = 0; itrk < phislice[etabin - 1].trackidx.size(); itrk++) + clusters[nclust].trackidx.push_back(phislice[etabin - 1].trackidx[itrk]); + } + + if (my_pt >= nextbin2_pt && nextbin_pt > 0) { + clusters[nclust].pTtot += nextbin_pt; + clusters[nclust].ntracks += phislice[etabin + 1].ntracks; + clusters[nclust].nxtracks += phislice[etabin + 1].nxtracks; + for (unsigned int itrk = 0; itrk < phislice[etabin + 1].trackidx.size(); itrk++) + clusters[nclust].trackidx.push_back(phislice[etabin + 1].trackidx[itrk]); + phislice[etabin + 1].used = true; + } + + nclust++; + + } // for each etabin + + // Merge close-by clusters + for (int m = 0; m < nclust - 1; ++m) { + if (((clusters[m + 1].eta - clusters[m].eta) < (3 * etaStep_) / 2) && + (-(3 * etaStep_) / 2 < (clusters[m + 1].eta - clusters[m].eta))) { + if (clusters[m + 1].pTtot > clusters[m].pTtot) { + clusters[m].eta = clusters[m + 1].eta; + } + clusters[m].pTtot += clusters[m + 1].pTtot; + clusters[m].ntracks += clusters[m + 1].ntracks; // total ntrk + clusters[m].nxtracks += clusters[m + 1].nxtracks; // total ndisp + for (unsigned int itrk = 0; itrk < clusters[m + 1].trackidx.size(); itrk++) + clusters[m].trackidx.push_back(clusters[m + 1].trackidx[itrk]); + + // if remove the merged cluster - all the others must be closer to 0 + for (int m1 = m + 1; m1 < nclust - 1; ++m1) + clusters[m1] = clusters[m1 + 1]; + + clusters.erase(clusters.begin() + nclust); + nclust--; + } // end if for cluster merging + } // end for (m) loop + + return clusters; + } + + // Fill L2 clusters (helper function) + template + inline void Fill_L2Cluster(T &bin, Pt pt, int ntrk, int ndtrk, std::vector trkidx) { + bin.pTtot += pt; + bin.ntracks += ntrk; + bin.nxtracks += ndtrk; + for (unsigned int itrk = 0; itrk < trkidx.size(); itrk++) + bin.trackidx.push_back(trkidx[itrk]); + } + + inline glbphi_intern DPhi(glbphi_intern phi1, glbphi_intern phi2) { + glbphi_intern x = phi1 - phi2; + if (x >= DoubleToBit( + M_PI, TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, TTTrack_TrackWord::stepPhi0)) + x -= DoubleToBit( + 2 * M_PI, TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, TTTrack_TrackWord::stepPhi0); + if (x < DoubleToBit(-1 * M_PI, + TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, + TTTrack_TrackWord::stepPhi0)) + x += DoubleToBit( + 2 * M_PI, TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, TTTrack_TrackWord::stepPhi0); + return x; + } + + inline float DPhi(float phi1, float phi2) { + float x = phi1 - phi2; + if (x >= M_PI) + x -= 2 * M_PI; + if (x < -1 * M_PI) + x += 2 * M_PI; + return x; + } + + // L2 clustering (in phi) + template + inline std::vector L2_clustering(std::vector > &L1clusters, + int phiBins_, + Phi phiStep_, + Eta etaStep_) { + std::vector clusters; + for (int phibin = 0; phibin < phiBins_; ++phibin) { //Find eta-phibin with highest pT + if (L1clusters[phibin].empty()) + continue; + + // sort L1 clusters max -> min + sort(L1clusters[phibin].begin(), L1clusters[phibin].end(), [](T &a, T &b) { return a.pTtot > b.pTtot; }); + + for (unsigned int imax = 0; imax < L1clusters[phibin].size(); ++imax) { + if (L1clusters[phibin][imax].used) + continue; + Pt pt_current = L1clusters[phibin][imax].pTtot; //current cluster (pt0) + Pt pt_next = 0; // next phi bin (pt1) + Pt pt_next2 = 0; // next to next phi bin2 (pt2) + int trk1 = 0; + int trk2 = 0; + int tdtrk1 = 0; + int tdtrk2 = 0; + std::vector trkidx1; + std::vector trkidx2; + clusters.push_back(L1clusters[phibin][imax]); + + L1clusters[phibin][imax].used = true; + + // if we are in the last phi bin, dont check phi+1 phi+2 + if (phibin == phiBins_ - 1) + continue; + + std::vector used_already; //keep phi+1 clusters that have been used + for (unsigned int icluster = 0; icluster < L1clusters[phibin + 1].size(); ++icluster) { + if (L1clusters[phibin + 1][icluster].used) + continue; + + if (((L1clusters[phibin + 1][icluster].eta - L1clusters[phibin][imax].eta) > (3 * etaStep_) / 2) || + ((L1clusters[phibin + 1][icluster].eta - L1clusters[phibin][imax].eta) < -(3 * etaStep_) / 2)) + continue; + + pt_next += L1clusters[phibin + 1][icluster].pTtot; + trk1 += L1clusters[phibin + 1][icluster].ntracks; + tdtrk1 += L1clusters[phibin + 1][icluster].nxtracks; + + for (unsigned int itrk = 0; itrk < L1clusters[phibin + 1][icluster].trackidx.size(); itrk++) + trkidx1.push_back(L1clusters[phibin + 1][icluster].trackidx[itrk]); + used_already.push_back(icluster); + } + + if (pt_next < pt_current) { // if pt1(clusters[clusters.size() - 1], pt_next, trk1, tdtrk1, trkidx1); + for (unsigned int iused : used_already) + L1clusters[phibin + 1][iused].used = true; + continue; + } + // if phi = next to last bin there is no "next to next" + if (phibin == phiBins_ - 2) { + Fill_L2Cluster(clusters[clusters.size() - 1], pt_next, trk1, tdtrk1, trkidx1); + clusters[clusters.size() - 1].phi = L1clusters[phibin + 1][used_already[0]].phi; + for (unsigned int iused : used_already) + L1clusters[phibin + 1][iused].used = true; + continue; + } + std::vector used_already2; //keep used clusters in phi+2 + for (unsigned int icluster = 0; icluster < L1clusters[phibin + 2].size(); ++icluster) { + if (L1clusters[phibin + 2][icluster].used) + continue; + if (((L1clusters[phibin + 2][icluster].eta - L1clusters[phibin][imax].eta) > (3 * etaStep_) / 2) || + ((L1clusters[phibin + 2][icluster].eta - L1clusters[phibin][imax].eta) < -(3 * etaStep_) / 2)) + continue; + pt_next2 += L1clusters[phibin + 2][icluster].pTtot; + trk2 += L1clusters[phibin + 2][icluster].ntracks; + tdtrk2 += L1clusters[phibin + 2][icluster].nxtracks; + + for (unsigned int itrk = 0; itrk < L1clusters[phibin + 2][icluster].trackidx.size(); itrk++) + trkidx2.push_back(L1clusters[phibin + 2][icluster].trackidx[itrk]); + used_already2.push_back(icluster); + } + if (pt_next2 < pt_next) { + std::vector trkidx_both; + trkidx_both.reserve(trkidx1.size() + trkidx2.size()); + trkidx_both.insert(trkidx_both.end(), trkidx1.begin(), trkidx1.end()); + trkidx_both.insert(trkidx_both.end(), trkidx2.begin(), trkidx2.end()); + Fill_L2Cluster( + clusters[clusters.size() - 1], pt_next + pt_next2, trk1 + trk2, tdtrk1 + tdtrk2, trkidx_both); + clusters[clusters.size() - 1].phi = L1clusters[phibin + 1][used_already[0]].phi; + for (unsigned int iused : used_already) + L1clusters[phibin + 1][iused].used = true; + for (unsigned int iused : used_already2) + L1clusters[phibin + 2][iused].used = true; + } + } // for each L1 cluster + } // for each phibin + + int nclust = clusters.size(); + + // merge close-by clusters + for (int m = 0; m < nclust - 1; ++m) { + for (int n = m + 1; n < nclust; ++n) { + if (clusters[n].eta != clusters[m].eta) + continue; + if ((DPhi(clusters[n].phi, clusters[m].phi) > (3 * phiStep_) / 2) || + (DPhi(clusters[n].phi, clusters[m].phi) < -(3 * phiStep_) / 2)) + continue; + if (clusters[n].pTtot > clusters[m].pTtot) + clusters[m].phi = clusters[n].phi; + clusters[m].pTtot += clusters[n].pTtot; + clusters[m].ntracks += clusters[n].ntracks; + clusters[m].nxtracks += clusters[n].nxtracks; + for (unsigned int itrk = 0; itrk < clusters[n].trackidx.size(); itrk++) + clusters[m].trackidx.push_back(clusters[n].trackidx[itrk]); + for (int m1 = n; m1 < nclust - 1; ++m1) + clusters[m1] = clusters[m1 + 1]; + clusters.erase(clusters.begin() + nclust); + + nclust--; + } // end of n-loop + } // end of m-loop + return clusters; + } +} // namespace l1ttrackjet +#endif diff --git a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulationProducer.cc b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulationProducer.cc deleted file mode 100644 index 701be4b196efa..0000000000000 --- a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulationProducer.cc +++ /dev/null @@ -1,681 +0,0 @@ -/* Software to emulate the hardware 2-layer jet-finding algorithm (fixed-point). *Layers 1 and 2* - * - * 2021 - * - * Created based on Rishi Patel's L1TrackJetProducer.cc file. - * Authors: Samuel Edwin Leigh, Tyler Wu - * Rutgers, the State University of New Jersey - * Revolutionary for 250 years - */ - -//Holds data from tracks, converted from their integer versions. - -// system include files - -#include "DataFormats/Common/interface/Ref.h" -#include "DataFormats/L1TCorrelator/interface/TkJet.h" -#include "DataFormats/L1TCorrelator/interface/TkJetFwd.h" -#include "DataFormats/L1TrackTrigger/interface/TTTypes.h" -#include "DataFormats/L1TrackTrigger/interface/TTTrack.h" -#include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h" -#include "DataFormats/L1Trigger/interface/TkJetWord.h" -#include "DataFormats/L1Trigger/interface/VertexWord.h" -// user include files -#include "FWCore/Framework/interface/Frameworkfwd.h" -#include "FWCore/Framework/interface/stream/EDProducer.h" -#include "FWCore/Framework/interface/Event.h" -#include "FWCore/Framework/interface/EventSetup.h" -#include "FWCore/Framework/interface/MakerMacros.h" -#include "FWCore/ParameterSet/interface/ParameterSet.h" -#include "FWCore/Utilities/interface/StreamID.h" -#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" -#include "Geometry/CommonTopologies/interface/PixelGeomDetUnit.h" -#include "Geometry/CommonTopologies/interface/PixelGeomDetType.h" -#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" - -#include "DataFormats/L1Trigger/interface/Vertex.h" - -#include -#include -#include -#include -#include -#include -#include -#include "TH1D.h" -#include "TH2D.h" -#include -#include -#include "L1Trigger/L1TTrackMatch/interface/L1TrackJetEmulationProducer.h" - -using namespace std; -using namespace edm; -class L1TrackJetEmulationProducer : public stream::EDProducer<> { -public: - explicit L1TrackJetEmulationProducer(const ParameterSet &); - ~L1TrackJetEmulationProducer() override; - typedef TTTrack L1TTTrackType; - typedef vector L1TTTrackCollectionType; - - static void fillDescriptions(ConfigurationDescriptions &descriptions); - bool trackQualityCuts(int trk_nstub, float trk_chi2, float trk_bendchi2); - void L2_cluster(vector> L1TrkPtrs_, TrackJetEmulationMaxZBin &mzb); - virtual TrackJetEmulationEtaPhiBin *L1_cluster(TrackJetEmulationEtaPhiBin *phislice); - -private: - void beginStream(StreamID) override; - void produce(Event &, const EventSetup &) override; - void endStream() override; - - // ----------member data --------------------------- - - vector> L1TrkPtrs_; - vector zBinCount_; - vector tdtrk_; - const float MaxDzTrackPV; - const float trkZMax_; - const float trkPtMax_; - const float trkPtMin_; - const float trkEtaMax_; - const float trkChi2dofMax_; - const float trkBendChi2Max_; - const int trkNPSStubMin_; - const double minTrkJetpT_; - int etaBins_; - int phiBins_; - int zBins_; - float d0CutNStubs4_; - float d0CutNStubs5_; - int lowpTJetMinTrackMultiplicity_; - float lowpTJetMinpT_; - int highpTJetMinTrackMultiplicity_; - float highpTJetMinpT_; - bool displaced_; - float nStubs4DisplacedChi2_; - float nStubs5DisplacedChi2_; - float nStubs4DisplacedBend_; - float nStubs5DisplacedBend_; - int nDisplacedTracks_; - - float PVz; - z0_intern zStep_; - glbeta_intern etaStep_; - glbphi_intern phiStep_; - - const EDGetTokenT>> trackToken_; - const EDGetTokenT PVtxToken_; - ESGetToken tTopoToken_; -}; - -L1TrackJetEmulationProducer::L1TrackJetEmulationProducer(const ParameterSet &iConfig) - : MaxDzTrackPV((float)iConfig.getParameter("MaxDzTrackPV")), - trkZMax_((float)iConfig.getParameter("trk_zMax")), - trkPtMax_((float)iConfig.getParameter("trk_ptMax")), - trkPtMin_((float)iConfig.getParameter("trk_ptMin")), - trkEtaMax_((float)iConfig.getParameter("trk_etaMax")), - trkChi2dofMax_((float)iConfig.getParameter("trk_chi2dofMax")), - trkBendChi2Max_((float)iConfig.getParameter("trk_bendChi2Max")), - trkNPSStubMin_((int)iConfig.getParameter("trk_nPSStubMin")), - minTrkJetpT_(iConfig.getParameter("minTrkJetpT")), - etaBins_((int)iConfig.getParameter("etaBins")), - phiBins_((int)iConfig.getParameter("phiBins")), - zBins_((int)iConfig.getParameter("zBins")), - d0CutNStubs4_((float)iConfig.getParameter("d0_cutNStubs4")), - d0CutNStubs5_((float)iConfig.getParameter("d0_cutNStubs5")), - lowpTJetMinTrackMultiplicity_((int)iConfig.getParameter("lowpTJetMinTrackMultiplicity")), - lowpTJetMinpT_((float)iConfig.getParameter("lowpTJetMinpT")), - highpTJetMinTrackMultiplicity_((int)iConfig.getParameter("highpTJetMinTrackMultiplicity")), - highpTJetMinpT_((float)iConfig.getParameter("highpTJetMinpT")), - displaced_(iConfig.getParameter("displaced")), - nStubs4DisplacedChi2_((float)iConfig.getParameter("nStubs4DisplacedChi2")), - nStubs5DisplacedChi2_((float)iConfig.getParameter("nStubs5DisplacedChi2")), - nStubs4DisplacedBend_((float)iConfig.getParameter("nStubs4Displacedbend")), - nStubs5DisplacedBend_((float)iConfig.getParameter("nStubs5Displacedbend")), - nDisplacedTracks_((int)iConfig.getParameter("nDisplacedTracks")), - trackToken_(consumes>>(iConfig.getParameter("L1TrackInputTag"))), - PVtxToken_(consumes(iConfig.getParameter("VertexInputTag"))), - tTopoToken_(esConsumes(edm::ESInputTag("", ""))) { - zStep_ = convert::makeZ0(2.0 * trkZMax_ / (zBins_ + 1)); - etaStep_ = convert::makeGlbEta(2.0 * trkEtaMax_ / etaBins_); //etaStep is the width of an etabin - phiStep_ = convert::makeGlbPhi(2.0 * M_PI / phiBins_); ////phiStep is the width of a phibin - - if (displaced_) - produces("L1TrackJetsExtended"); - else - produces("L1TrackJets"); -} - -L1TrackJetEmulationProducer::~L1TrackJetEmulationProducer() {} - -void L1TrackJetEmulationProducer::produce(Event &iEvent, const EventSetup &iSetup) { - unique_ptr L1L1TrackJetProducer(new l1t::TkJetWordCollection); - - // For TTStubs - const TrackerTopology &tTopo = iSetup.getData(tTopoToken_); - - edm::Handle>> TTTrackHandle; - iEvent.getByToken(trackToken_, TTTrackHandle); - vector>::const_iterator iterL1Track; - - edm::Handle PVtx; - iEvent.getByToken(PVtxToken_, PVtx); - float PVz = (PVtx->at(0)).z0(); - - L1TrkPtrs_.clear(); - zBinCount_.clear(); - tdtrk_.clear(); - unsigned int this_l1track = 0; - for (iterL1Track = TTTrackHandle->begin(); iterL1Track != TTTrackHandle->end(); iterL1Track++) { - edm::Ptr trkPtr(TTTrackHandle, this_l1track); - this_l1track++; - float trk_pt = trkPtr->momentum().perp(); - int trk_nstubs = (int)trkPtr->getStubRefs().size(); - float trk_chi2dof = trkPtr->chi2Red(); - float trk_d0 = trkPtr->d0(); - float trk_bendchi2 = trkPtr->stubPtConsistency(); - float trk_z0 = trkPtr->z0(); - - int trk_nPS = 0; - for (int istub = 0; istub < trk_nstubs; istub++) { // loop over the stubs - DetId detId(trkPtr->getStubRefs().at(istub)->getDetId()); - if (detId.det() == DetId::Detector::Tracker) { - if ((detId.subdetId() == StripSubdetector::TOB && tTopo.tobLayer(detId) <= 3) || - (detId.subdetId() == StripSubdetector::TID && tTopo.tidRing(detId) <= 9)) - trk_nPS++; - } - } - - if (trk_nPS < trkNPSStubMin_) - continue; - if (!trackQualityCuts(trk_nstubs, trk_chi2dof, trk_bendchi2)) - continue; - if (std::abs(trk_z0 - PVz) > MaxDzTrackPV) - continue; - if (std::abs(trk_z0) > trkZMax_) - continue; - if (std::abs(trkPtr->momentum().eta()) > trkEtaMax_) - continue; - if (trk_pt < trkPtMin_) - continue; - L1TrkPtrs_.push_back(trkPtr); - zBinCount_.push_back(0); - - if ((std::abs(trk_d0) > d0CutNStubs5_ && trk_nstubs >= 5 && d0CutNStubs5_ >= 0) || - (trk_nstubs == 4 && std::abs(trk_d0) > d0CutNStubs4_ && d0CutNStubs4_ >= 0)) - tdtrk_.push_back(1); //displaced track - else - tdtrk_.push_back(0); // not displaced track - } - - if (!L1TrkPtrs_.empty()) { - TrackJetEmulationMaxZBin mzb; - - L2_cluster(L1TrkPtrs_, mzb); - if (mzb.clusters != nullptr) { - for (int j = 0; j < mzb.nclust; ++j) { - //FILL Two Layer Jets for Jet Collection - if (mzb.clusters[j].pTtot < convert::makePtFromFloat(trkPtMin_)) - continue; //protects against reading bad memory - if (mzb.clusters[j].ntracks < 1) - continue; - if (mzb.clusters[j].ntracks > 5000) - continue; - l1t::TkJetWord::glbeta_t jetEta = mzb.clusters[j].eta * convert::ETA_LSB_POW; - l1t::TkJetWord::glbphi_t jetPhi = mzb.clusters[j].phi * convert::PHI_LSB_POW; - l1t::TkJetWord::z0_t jetZ0 = mzb.zbincenter * convert::Z0_LSB_POW; - l1t::TkJetWord::pt_t jetPt = mzb.clusters[j].pTtot; - l1t::TkJetWord::nt_t totalntracks_ = mzb.clusters[j].ntracks; - l1t::TkJetWord::nx_t totalxtracks_ = mzb.clusters[j].nxtracks; - l1t::TkJetWord::tkjetunassigned_t unassigned_ = 0; - - l1t::TkJetWord trkJet(jetPt, jetEta, jetPhi, jetZ0, totalntracks_, totalxtracks_, unassigned_); - //trkJet.setDispCounters(DispCounters); - L1L1TrackJetProducer->push_back(trkJet); - } - } else if (mzb.clusters == nullptr) { - edm::LogWarning("L1TrackJetEmulationProducer") << "mzb.clusters Not Assigned!\n"; - } - - if (displaced_) - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended"); - else - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets"); - delete[] mzb.clusters; - } else if (L1TrkPtrs_.empty()) { - edm::LogWarning("L1TrackJetEmulationProducer") << "L1TrkPtrs Not Assigned!\n"; - if (displaced_) - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended"); - else - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets"); - } -} - -void L1TrackJetEmulationProducer::L2_cluster(vector> L1TrkPtrs_, TrackJetEmulationMaxZBin &mzb) { - enum TrackBitWidths { - kEtaSize = 16, // Width of z-position (40cm / 0.1) - kEtaMagSize = 3, // Width of z-position magnitude (signed) - kPtSize = 14, // Width of pt - kPtMagSize = 9, // Width of pt magnitude (unsigned) - kPhiSize = 12, - kPhiMagSize = 1, - }; - - const int nz = zBins_ + 1; - TrackJetEmulationMaxZBin all_zBins[nz]; - static TrackJetEmulationMaxZBin mzbtemp; - for (int z = 0; z < nz; ++z) - all_zBins[z] = mzbtemp; - - z0_intern zmin = convert::makeZ0(-1.0 * trkZMax_); - z0_intern zmax = zmin + 2 * zStep_; - - TrackJetEmulationEtaPhiBin epbins[phiBins_][etaBins_]; // create grid of phiBins - glbphi_intern phi = convert::makeGlbPhi(-1.0 * M_PI); - glbeta_intern eta; - glbeta_intern etamin, etamax, phimin, phimax; - for (int i = 0; i < phiBins_; ++i) { - eta = convert::makeGlbEta(-1.0 * trkEtaMax_); - for (int j = 0; j < etaBins_; ++j) { - phimin = phi; - phimax = phi + phiStep_; - etamin = eta; - eta = eta + etaStep_; - etamax = eta; - epbins[i][j].phi = (phimin + phimax) / 2; - epbins[i][j].eta = (etamin + etamax) / 2; - epbins[i][j].pTtot = 0; - epbins[i][j].ntracks = 0; - epbins[i][j].nxtracks = 0; - } // for each etabin - phi = phi + phiStep_; - } // for each phibin (finished creating epbins) - - mzb = all_zBins[0]; - mzb.ht = 0; - int ntracks = L1TrkPtrs_.size(); - // uninitalized arrays - TrackJetEmulationEtaPhiBin *L1clusters[phiBins_]; - TrackJetEmulationEtaPhiBin L2cluster[ntracks]; - - for (int zbin = 0; zbin < zBins_; ++zbin) { - for (int i = 0; i < phiBins_; ++i) { //First initialize pT, numtracks, used to 0 (or false) - for (int j = 0; j < etaBins_; ++j) { - epbins[i][j].pTtot = 0; - epbins[i][j].used = false; - epbins[i][j].ntracks = 0; - epbins[i][j].nxtracks = 0; - } //for each etabin - L1clusters[i] = epbins[i]; - } //for each phibin - - for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) { - ap_ufixed inputTrkPt = 0; - inputTrkPt.V = L1TrkPtrs_[k]->getTrackWord()(TTTrack_TrackWord::TrackBitLocations::kRinvMSB - 1, - TTTrack_TrackWord::TrackBitLocations::kRinvLSB); - pt_intern trkpt = inputTrkPt; - - ap_fixed trketainput = 0; - trketainput.V = L1TrkPtrs_[k]->getTrackWord()(TTTrack_TrackWord::TrackBitLocations::kTanlMSB, - TTTrack_TrackWord::TrackBitLocations::kTanlLSB); - ap_ufixed<64 + ETA_EXTRABITS, 32 + ETA_EXTRABITS> eta_conv = - 1.0 / convert::ETA_LSB; //conversion factor from input eta format to internal format - glbeta_intern trketa = eta_conv * trketainput; - - int Sector = L1TrkPtrs_[k]->phiSector(); - glbphi_intern trkphiSector = 0; - if (Sector < 5) { - trkphiSector = Sector * convert::makeGlbPhi(2.0 * M_PI / 9.0); - } else { - trkphiSector = - convert::makeGlbPhi(-1.0 * M_PI + M_PI / 9.0) + (Sector - 5) * convert::makeGlbPhi(2.0 * M_PI / 9.0); - } - - ap_fixed trkphiinput = 0; - trkphiinput.V = L1TrkPtrs_[k]->getTrackWord()(TTTrack_TrackWord::TrackBitLocations::kPhiMSB, - TTTrack_TrackWord::TrackBitLocations::kPhiLSB); - ap_ufixed<64 + PHI_EXTRABITS, 32 + PHI_EXTRABITS> phi_conv = - TTTrack_TrackWord::stepPhi0 / convert::PHI_LSB * - (1 << (TrackBitWidths::kPhiSize - TrackBitWidths::kPhiMagSize)); - //phi_conv is a conversion factor from input phi format to the internal format - glbeta_intern localphi = phi_conv * trkphiinput; - glbeta_intern trkphi = localphi + trkphiSector; - - ap_int inputTrkZ0 = L1TrkPtrs_[k]->getTrackWord()( - TTTrack_TrackWord::TrackBitLocations::kZ0MSB, TTTrack_TrackWord::TrackBitLocations::kZ0LSB); - ap_ufixed<32 + Z0_EXTRABITS, 1 + Z0_EXTRABITS> z0_conv = - TTTrack_TrackWord::stepZ0 / convert::Z0_LSB; //conversion factor from input z format to internal format - z0_intern trkZ = z0_conv * inputTrkZ0; - - for (int i = 0; i < phiBins_; ++i) { - for (int j = 0; j < etaBins_; ++j) { - L2cluster[k] = epbins[i][j]; - if ((zmin <= trkZ && zmax >= trkZ) && - ((epbins[i][j].eta - etaStep_ / 2 < trketa && epbins[i][j].eta + etaStep_ / 2 >= trketa) && - epbins[i][j].phi - phiStep_ / 2 < trkphi && epbins[i][j].phi + phiStep_ / 2 >= trkphi && - (zBinCount_[k] != 2))) { - zBinCount_.at(k) = zBinCount_.at(k) + 1; - - if (trkpt < convert::makePtFromFloat(trkPtMax_)) - epbins[i][j].pTtot += trkpt; - else - epbins[i][j].pTtot += convert::makePtFromFloat(trkPtMax_); - ++epbins[i][j].ntracks; - //x-bit is currently not used in firmware, so we leave nxtracks = 0 for now - } // if right bin - } // for each phibin: j loop - } // for each phibin: i loop - } // end loop over tracks - - for (int phislice = 0; phislice < phiBins_; ++phislice) { - L1clusters[phislice] = L1_cluster(epbins[phislice]); - for (int ind = 0; L1clusters[phislice][ind].pTtot != 0; ++ind) { - L1clusters[phislice][ind].used = false; - } - } - - //Create clusters array to hold output cluster data for Layer2; can't have more clusters than tracks. - //Find eta-phibin with maxpT, make center of cluster, add neighbors if not already used. - pt_intern hipT = 0; - int nclust = 0; - int phibin = 0; - int imax = -1; - int index1; //index of clusters array for each phislice - pt_intern E1 = 0; - pt_intern E0 = 0; - pt_intern E2 = 0; - l1t::TkJetWord::nt_t ntrk1, ntrk2; - l1t::TkJetWord::nx_t nxtrk1, nxtrk2; - int used1, used2, used3, used4; - - for (phibin = 0; phibin < phiBins_; ++phibin) { //Find eta-phibin with highest pT - while (true) { - hipT = 0; - for (index1 = 0; L1clusters[phibin][index1].pTtot > 0; ++index1) { - if (!L1clusters[phibin][index1].used && L1clusters[phibin][index1].pTtot >= hipT) { - hipT = L1clusters[phibin][index1].pTtot; - imax = index1; - } - } // for each index within the phibin - - if (hipT == 0) - break; // If highest pT is 0, all bins are used - E0 = hipT; // E0 is pT of first phibin of the cluster - E1 = 0; - E2 = 0; - ntrk1 = 0; - ntrk2 = 0; - nxtrk1 = 0; - nxtrk2 = 0; - L2cluster[nclust] = L1clusters[phibin][imax]; - - L1clusters[phibin][imax].used = true; - // Add pT of upper neighbor - // E1 is pT of the middle phibin (should be highest pT) - if (phibin != phiBins_ - 1) { - used1 = -1; - used2 = -1; - for (index1 = 0; L1clusters[phibin + 1][index1].pTtot != 0; ++index1) { - if (L1clusters[phibin + 1][index1].used) - continue; - if (L1clusters[phibin + 1][index1].eta - L1clusters[phibin][imax].eta <= 3 * etaStep_ / 2 && - L1clusters[phibin][imax].eta - L1clusters[phibin + 1][index1].eta <= 3 * etaStep_ / 2) { - E1 += L1clusters[phibin + 1][index1].pTtot; - ntrk1 += L1clusters[phibin + 1][index1].ntracks; - nxtrk1 += L1clusters[phibin + 1][index1].nxtracks; - if (used1 < 0) - used1 = index1; - else - used2 = index1; - } // if cluster is within one phibin - } // for each cluster in above phibin - - if (E1 < E0) { // if E1 isn't higher, E0 and E1 are their own cluster - L2cluster[nclust].pTtot += E1; - L2cluster[nclust].ntracks += ntrk1; - L2cluster[nclust].nxtracks += nxtrk1; - if (used1 >= 0) - L1clusters[phibin + 1][used1].used = true; - if (used2 >= 0) - L1clusters[phibin + 1][used2].used = true; - nclust++; - continue; - } - - if (phibin != phiBins_ - 2) { // E2 will be the pT of the third phibin (should be lower than E1) - used3 = -1; - used4 = -1; - for (index1 = 0; L1clusters[phibin + 2][index1].pTtot != 0; ++index1) { - if (L1clusters[phibin + 2][index1].used) - continue; - if (L1clusters[phibin + 2][index1].eta - L1clusters[phibin][imax].eta <= 3 * etaStep_ / 2 && - L1clusters[phibin][imax].eta - L1clusters[phibin + 2][index1].eta <= 3 * etaStep_ / 2) { - E2 += L1clusters[phibin + 2][index1].pTtot; - ntrk2 += L1clusters[phibin + 2][index1].ntracks; - nxtrk2 += L1clusters[phibin + 2][index1].nxtracks; - if (used3 < 0) - used3 = index1; - else - used4 = index1; - } - } - // if indeed E2 < E1, add E1 and E2 to E0, they're all a cluster together - // otherwise, E0 is its own cluster - if (E2 < E1) { - L2cluster[nclust].pTtot += E1 + E2; - L2cluster[nclust].ntracks += ntrk1 + ntrk2; - L2cluster[nclust].nxtracks += nxtrk1 + nxtrk2; - L2cluster[nclust].phi = L1clusters[phibin + 1][used1].phi; - if (used1 >= 0) - L1clusters[phibin + 1][used1].used = true; - if (used2 >= 0) - L1clusters[phibin + 1][used2].used = true; - if (used3 >= 0) - L1clusters[phibin + 2][used3].used = true; - if (used4 >= 0) - L1clusters[phibin + 2][used4].used = true; - } - nclust++; - continue; - } // end Not phiBins-2 - else { - L2cluster[nclust].pTtot += E1; - L2cluster[nclust].ntracks += ntrk1; - L2cluster[nclust].nxtracks += nxtrk1; - L2cluster[nclust].phi = L1clusters[phibin + 1][used1].phi; - if (used1 >= 0) - L1clusters[phibin + 1][used1].used = true; - if (used2 >= 0) - L1clusters[phibin + 1][used2].used = true; - nclust++; - continue; - } - } //End not last phibin(23) - else { //if it is phibin 23 - L1clusters[phibin][imax].used = true; - nclust++; - } - } // while hipT not 0 - } // for each phibin - - for (phibin = 0; phibin < phiBins_; ++phibin) - delete[] L1clusters[phibin]; - - // Now merge clusters, if necessary - for (int m = 0; m < nclust - 1; ++m) { - for (int n = m + 1; n < nclust; ++n) { - if (L2cluster[n].eta == L2cluster[m].eta && - ((L2cluster[n].phi - L2cluster[m].phi < 3 * phiStep_ / 2 && - L2cluster[m].phi - L2cluster[n].phi < 3 * phiStep_ / 2) || - (L2cluster[n].phi - L2cluster[m].phi > 2 * convert::makeGlbPhi(M_PI) - phiStep_ || - L2cluster[m].phi - L2cluster[n].phi > 2 * convert::makeGlbPhi(M_PI) - phiStep_))) { - if (L2cluster[n].pTtot > L2cluster[m].pTtot) { - L2cluster[m].phi = L2cluster[n].phi; - } - L2cluster[m].pTtot += L2cluster[n].pTtot; - L2cluster[m].ntracks += L2cluster[n].ntracks; - L2cluster[m].nxtracks += L2cluster[n].nxtracks; - for (int m1 = n; m1 < nclust - 1; ++m1) { - L2cluster[m1] = L2cluster[m1 + 1]; - } - nclust--; - m = -1; - break; //????? - } // end if clusters neighbor in eta - } - } // end for (m) loop - // sum up all pTs in this zbin to find ht - - pt_intern ht = 0; - for (int k = 0; k < nclust; ++k) { - if (L2cluster[k].pTtot > convert::makePtFromFloat(lowpTJetMinpT_) && - L2cluster[k].ntracks < lowpTJetMinTrackMultiplicity_) - continue; - if (L2cluster[k].pTtot > convert::makePtFromFloat(highpTJetMinpT_) && - L2cluster[k].ntracks < highpTJetMinTrackMultiplicity_) - continue; - if (L2cluster[k].pTtot > convert::makePtFromFloat(minTrkJetpT_)) - ht += L2cluster[k].pTtot; - } - // if ht is larger than previous max, this is the new vertex zbin - all_zBins[zbin].znum = zbin; - all_zBins[zbin].clusters = new TrackJetEmulationEtaPhiBin[nclust]; - all_zBins[zbin].nclust = nclust; - all_zBins[zbin].zbincenter = (zmin + zmax) / 2; - for (int k = 0; k < nclust; ++k) { - all_zBins[zbin].clusters[k].phi = L2cluster[k].phi; - all_zBins[zbin].clusters[k].eta = L2cluster[k].eta; - all_zBins[zbin].clusters[k].pTtot = L2cluster[k].pTtot; - all_zBins[zbin].clusters[k].ntracks = L2cluster[k].ntracks; - all_zBins[zbin].clusters[k].nxtracks = L2cluster[k].nxtracks; - } - all_zBins[zbin].ht = ht; - if (ht >= mzb.ht) { - mzb = all_zBins[zbin]; - mzb.zbincenter = (zmin + zmax) / 2; - } - // Prepare for next zbin! - zmin = zmin + zStep_; - zmax = zmax + zStep_; - } // for each zbin - - for (int zbin = 0; zbin < zBins_; ++zbin) { - if (zbin == mzb.znum) { - continue; - } - delete[] all_zBins[zbin].clusters; - } -} - -TrackJetEmulationEtaPhiBin *L1TrackJetEmulationProducer::L1_cluster(TrackJetEmulationEtaPhiBin *phislice) { - TrackJetEmulationEtaPhiBin *clusters = new TrackJetEmulationEtaPhiBin[etaBins_ / 2]; - for (int etabin = 0; etabin < etaBins_ / 2; ++etabin) { - clusters[etabin].pTtot = 0; - clusters[etabin].ntracks = 0; - clusters[etabin].nxtracks = 0; - clusters[etabin].phi = 0; - clusters[etabin].eta = 0; - clusters[etabin].used = false; - } - - if (clusters == nullptr) - edm::LogWarning("L1TrackJetEmulationProducer") << "Clusters memory not assigned!\n"; - - // Find eta-phibin with maxpT, make center of cluster, add neighbors if not already used - pt_intern my_pt, left_pt, right_pt, right2pt; - int nclust = 0; - right2pt = 0; - for (int etabin = 0; etabin < etaBins_; ++etabin) { - // assign values for my pT and neighbors' pT - if (phislice[etabin].used) - continue; - my_pt = phislice[etabin].pTtot; - if (etabin > 0 && !phislice[etabin - 1].used) { - left_pt = phislice[etabin - 1].pTtot; - } else - left_pt = 0; - if (etabin < etaBins_ - 1 && !phislice[etabin + 1].used) { - right_pt = phislice[etabin + 1].pTtot; - if (etabin < etaBins_ - 2 && !phislice[etabin + 2].used) { - right2pt = phislice[etabin + 2].pTtot; - } else - right2pt = 0; - } else - right_pt = 0; - - // if I'm not a cluster, move on - if (my_pt < left_pt || my_pt <= right_pt) { - // if unused pT in the left neighbor, spit it out as a cluster - if (left_pt > 0) { - clusters[nclust] = phislice[etabin - 1]; - phislice[etabin - 1].used = true; - nclust++; - } - continue; - } - - // I guess I'm a cluster-- should I use my right neighbor? - // Note: left neighbor will definitely be used because if it - // didn't belong to me it would have been used already - clusters[nclust] = phislice[etabin]; - phislice[etabin].used = true; - if (left_pt > 0) { - clusters[nclust].pTtot += left_pt; - clusters[nclust].ntracks += phislice[etabin - 1].ntracks; - clusters[nclust].nxtracks += phislice[etabin - 1].nxtracks; - } - if (my_pt >= right2pt && right_pt > 0) { - clusters[nclust].pTtot += right_pt; - clusters[nclust].ntracks += phislice[etabin + 1].ntracks; - clusters[nclust].nxtracks += phislice[etabin + 1].nxtracks; - phislice[etabin + 1].used = true; - } - nclust++; - } // for each etabin - - // Now merge clusters, if necessary - for (int m = 0; m < nclust - 1; ++m) { - if (clusters[m + 1].eta - clusters[m].eta < 3 * etaStep_ / 2 && - clusters[m].eta - clusters[m + 1].eta < 3 * etaStep_ / 2) { - if (clusters[m + 1].pTtot > clusters[m].pTtot) { - clusters[m].eta = clusters[m + 1].eta; - } - clusters[m].pTtot += clusters[m + 1].pTtot; - clusters[m].ntracks += clusters[m + 1].ntracks; // Previous version didn't add tracks when merging - clusters[m].nxtracks += clusters[m + 1].nxtracks; - for (int m1 = m + 1; m1 < nclust - 1; ++m1) - clusters[m1] = clusters[m1 + 1]; - nclust--; - m = -1; - } // end if clusters neighbor in eta - } // end for (m) loop - - for (int i = nclust; i < etaBins_ / 2; ++i) // zero out remaining unused clusters - clusters[i].pTtot = 0; - return clusters; -} - -void L1TrackJetEmulationProducer::beginStream(StreamID) {} - -void L1TrackJetEmulationProducer::endStream() {} - -bool L1TrackJetEmulationProducer::trackQualityCuts(int trk_nstub, float trk_chi2, float trk_bendchi2) { - bool PassQuality = false; - if (trk_bendchi2 < trkBendChi2Max_ && trk_chi2 < trkChi2dofMax_ && trk_nstub >= 4 && !displaced_) - PassQuality = true; - if (displaced_ && trk_bendchi2 < nStubs4DisplacedBend_ && trk_chi2 < nStubs4DisplacedChi2_ && trk_nstub == 4) - PassQuality = true; - if (displaced_ && trk_bendchi2 < nStubs5DisplacedBend_ && trk_chi2 < nStubs5DisplacedChi2_ && trk_nstub > 4) - PassQuality = true; - return PassQuality; -} - -void L1TrackJetEmulationProducer::fillDescriptions(ConfigurationDescriptions &descriptions) { - //The following says we do not know what parameters are allowed so do no validation - // Please change this to state exactly what you do use, even if it is no parameters - ParameterSetDescription desc; - desc.setUnknown(); - descriptions.addDefault(desc); -} - -//define this as a plug-in -DEFINE_FWK_MODULE(L1TrackJetEmulationProducer); diff --git a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulatorProducer.cc b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulatorProducer.cc new file mode 100644 index 0000000000000..96ce847d7efa5 --- /dev/null +++ b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetEmulatorProducer.cc @@ -0,0 +1,456 @@ +// Original Author: Samuel Edwin Leigh, Tyler Wu, +// Rutgers, the State University of New Jersey +// +// Rewritting/Improvements: George Karathanasis, +// georgios.karathanasis@cern.ch, CU Boulder +// +// Created: Wed, 01 Aug 2018 14:01:41 GMT +// Latest update: Nov 2022 (by GK) +// +// Track jets are clustered in a two-layer process, first by clustering in phi, +// then by clustering in eta. The code proceeds as following: putting all tracks// in a grid of eta vs phi space, and then cluster them. Finally we merge the cl +// usters when needed. The code is improved to use the same module between emula +// tion and simulation was also improved, with bug fixes and being faster. +// Introduction to object (p10-13): +// https://indico.cern.ch/event/791517/contributions/3341650/attachments/1818736/2973771/TrackBasedAlgos_L1TMadrid_MacDonald.pdf +// New and improved version: https://indico.cern.ch/event/1203796/contributions/5073056/attachments/2519806/4333006/trackjet_emu.pdf + +// L1T include files +#include "DataFormats/L1TCorrelator/interface/TkJet.h" +#include "DataFormats/L1TCorrelator/interface/TkJetFwd.h" +#include "DataFormats/L1TrackTrigger/interface/TTTypes.h" +#include "DataFormats/L1TrackTrigger/interface/TTTrack.h" +#include "DataFormats/L1TrackTrigger/interface/TTTrack_TrackWord.h" +#include "DataFormats/L1Trigger/interface/TkJetWord.h" +#include "DataFormats/L1Trigger/interface/VertexWord.h" + +// system include files +#include "FWCore/Framework/interface/Frameworkfwd.h" +#include "FWCore/Framework/interface/stream/EDProducer.h" +#include "FWCore/Framework/interface/Event.h" +#include "FWCore/Framework/interface/EventSetup.h" +#include "FWCore/Framework/interface/MakerMacros.h" +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#include "FWCore/Utilities/interface/StreamID.h" +#include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" +#include "DataFormats/Common/interface/Ref.h" + +//own headers +#include "L1TrackJetClustering.h" + +//general +#include + +using namespace std; +using namespace edm; +using namespace l1t; +using namespace l1ttrackjet; + +class L1TrackJetEmulatorProducer : public stream::EDProducer<> { +public: + explicit L1TrackJetEmulatorProducer(const ParameterSet &); + ~L1TrackJetEmulatorProducer() override = default; + typedef TTTrack L1TTTrackType; + typedef vector L1TTTrackCollectionType; + static void fillDescriptions(ConfigurationDescriptions &descriptions); + +private: + void produce(Event &, const EventSetup &) override; + + // ----------member data --------------------------- + + std::vector> L1TrkPtrs_; + vector tdtrk_; + const float trkZMax_; + const float trkPtMax_; + const float trkPtMin_; + const float trkEtaMax_; + const float nStubs4PromptChi2_; + const float nStubs5PromptChi2_; + const float nStubs4PromptBend_; + const float nStubs5PromptBend_; + const int trkNPSStubMin_; + const int lowpTJetMinTrackMultiplicity_; + const float lowpTJetThreshold_; + const int highpTJetMinTrackMultiplicity_; + const float highpTJetThreshold_; + const int zBins_; + const int etaBins_; + const int phiBins_; + const double minTrkJetpT_; + const bool displaced_; + const float d0CutNStubs4_; + const float d0CutNStubs5_; + const float nStubs4DisplacedChi2_; + const float nStubs5DisplacedChi2_; + const float nStubs4DisplacedBend_; + const float nStubs5DisplacedBend_; + const int nDisplacedTracks_; + const float dzPVTrk_; + + float PVz; + float zStep_; + glbeta_intern etaStep_; + glbphi_intern phiStep_; + + TTTrack_TrackWord trackword; + + edm::ESGetToken tTopoToken_; + const EDGetTokenT>> trackToken_; + const EDGetTokenT PVtxToken_; +}; + +//constructor +L1TrackJetEmulatorProducer::L1TrackJetEmulatorProducer(const ParameterSet &iConfig) + : trkZMax_(iConfig.getParameter("trk_zMax")), + trkPtMax_(iConfig.getParameter("trk_ptMax")), + trkPtMin_(iConfig.getParameter("trk_ptMin")), + trkEtaMax_(iConfig.getParameter("trk_etaMax")), + nStubs4PromptChi2_(iConfig.getParameter("nStubs4PromptChi2")), + nStubs5PromptChi2_(iConfig.getParameter("nStubs5PromptChi2")), + nStubs4PromptBend_(iConfig.getParameter("nStubs4PromptBend")), + nStubs5PromptBend_(iConfig.getParameter("nStubs5PromptBend")), + trkNPSStubMin_(iConfig.getParameter("trk_nPSStubMin")), + lowpTJetMinTrackMultiplicity_(iConfig.getParameter("lowpTJetMinTrackMultiplicity")), + lowpTJetThreshold_(iConfig.getParameter("lowpTJetThreshold")), + highpTJetMinTrackMultiplicity_(iConfig.getParameter("highpTJetMinTrackMultiplicity")), + highpTJetThreshold_(iConfig.getParameter("highpTJetThreshold")), + zBins_(iConfig.getParameter("zBins")), + etaBins_(iConfig.getParameter("etaBins")), + phiBins_(iConfig.getParameter("phiBins")), + minTrkJetpT_(iConfig.getParameter("minTrkJetpT")), + displaced_(iConfig.getParameter("displaced")), + d0CutNStubs4_(iConfig.getParameter("d0_cutNStubs4")), + d0CutNStubs5_(iConfig.getParameter("d0_cutNStubs5")), + nStubs4DisplacedChi2_(iConfig.getParameter("nStubs4DisplacedChi2")), + nStubs5DisplacedChi2_(iConfig.getParameter("nStubs5DisplacedChi2")), + nStubs4DisplacedBend_(iConfig.getParameter("nStubs4DisplacedBend")), + nStubs5DisplacedBend_(iConfig.getParameter("nStubs5DisplacedBend")), + nDisplacedTracks_(iConfig.getParameter("nDisplacedTracks")), + dzPVTrk_(iConfig.getParameter("MaxDzTrackPV")), + tTopoToken_(esConsumes(edm::ESInputTag("", ""))), + trackToken_(consumes>>(iConfig.getParameter("L1TrackInputTag"))), + PVtxToken_(consumes(iConfig.getParameter("L1PVertexInputTag"))) { + zStep_ = 2.0 * trkZMax_ / (zBins_ + 1); // added +1 in denom + etaStep_ = glbeta_intern(2.0 * trkEtaMax_ / etaBins_); //etaStep is the width of an etabin + phiStep_ = DoubleToBit(2.0 * (M_PI) / phiBins_, + TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, + TTTrack_TrackWord::stepPhi0); ///phiStep is the width of a phibin + + if (displaced_) + produces("L1TrackJetsExtended"); + else + produces("L1TrackJets"); +} + +void L1TrackJetEmulatorProducer::produce(Event &iEvent, const EventSetup &iSetup) { + unique_ptr L1TrackJetContainer(new l1t::TkJetWordCollection); + // Read inputs + const TrackerTopology &tTopo = iSetup.getData(tTopoToken_); + + edm::Handle>> TTTrackHandle; + iEvent.getByToken(trackToken_, TTTrackHandle); + + edm::Handle PVtx; + iEvent.getByToken(PVtxToken_, PVtx); + float PVz = (PVtx->at(0)).z0(); + + L1TrkPtrs_.clear(); + tdtrk_.clear(); + // track selection + for (unsigned int this_l1track = 0; this_l1track < TTTrackHandle->size(); this_l1track++) { + edm::Ptr trkPtr(TTTrackHandle, this_l1track); + float trk_pt = trkPtr->momentum().perp(); + int trk_nstubs = (int)trkPtr->getStubRefs().size(); + float trk_chi2dof = trkPtr->chi2Red(); + float trk_bendchi2 = trkPtr->stubPtConsistency(); + int trk_nPS = 0; + for (int istub = 0; istub < trk_nstubs; istub++) { + DetId detId(trkPtr->getStubRefs().at(istub)->getDetId()); + if (detId.det() == DetId::Detector::Tracker) { + if ((detId.subdetId() == StripSubdetector::TOB && tTopo.tobLayer(detId) <= 3) || + (detId.subdetId() == StripSubdetector::TID && tTopo.tidRing(detId) <= 9)) + trk_nPS++; + } + } + // selection tracks - supposed to happen on seperate module (kept for legacy/debug reasons) + if (trk_nPS < trkNPSStubMin_) + continue; + if (!TrackQualitySelection(trk_nstubs, + trk_chi2dof, + trk_bendchi2, + nStubs4PromptBend_, + nStubs5PromptBend_, + nStubs4PromptChi2_, + nStubs5PromptChi2_, + nStubs4DisplacedBend_, + nStubs5DisplacedBend_, + nStubs4DisplacedChi2_, + nStubs5DisplacedChi2_, + displaced_)) + continue; + if (std::abs(PVz - trkPtr->z0()) > dzPVTrk_ && dzPVTrk_ > 0) + continue; + if (std::abs(trkPtr->z0()) > trkZMax_) + continue; + if (std::abs(trkPtr->momentum().eta()) > trkEtaMax_) + continue; + if (trk_pt < trkPtMin_) + continue; + L1TrkPtrs_.push_back(trkPtr); + } + + // if no tracks pass selection return empty containers + if (L1TrkPtrs_.empty()) { + if (displaced_) + iEvent.put(std::move(L1TrackJetContainer), "L1TrackJetsExtended"); + else + iEvent.put(std::move(L1TrackJetContainer), "L1TrackJets"); + return; + } + + TrackJetEmulationMaxZBin mzb; + mzb.znum = 0; + mzb.nclust = 0; + mzb.ht = 0; + + TrackJetEmulationEtaPhiBin epbins_default[phiBins_][etaBins_]; // create grid of phiBins + glbphi_intern phi = DoubleToBit( + -1.0 * M_PI, TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, TTTrack_TrackWord::stepPhi0); + for (int i = 0; i < phiBins_; ++i) { + glbeta_intern eta = -1 * trkEtaMax_; + for (int j = 0; j < etaBins_; ++j) { + epbins_default[i][j].phi = (phi + (phi + phiStep_)) / 2; // phi coord of bin center + epbins_default[i][j].eta = (eta + (eta + etaStep_)) / 2; // eta coord of bin center + eta += etaStep_; + } // for each etabin + phi += phiStep_; + } // for each phibin (finished creating epbins) + + // bins for z if multibin option is selected + std::vector zmins, zmaxs; + for (int zbin = 0; zbin < zBins_; zbin++) { + zmins.push_back(DoubleToBit( + -1.0 * trkZMax_ + zStep_ * zbin, TTTrack_TrackWord::TrackBitWidths::kZ0Size, TTTrack_TrackWord::stepZ0)); + zmaxs.push_back(DoubleToBit(-1.0 * trkZMax_ + zStep_ * zbin + 2 * zStep_, + TTTrack_TrackWord::TrackBitWidths::kZ0Size, + TTTrack_TrackWord::stepZ0)); + } + + // create vectors that hold clusters + std::vector> L1clusters; + L1clusters.reserve(phiBins_); + std::vector L2clusters; + + // default (empty) grid + for (int i = 0; i < phiBins_; ++i) { + for (int j = 0; j < etaBins_; ++j) { + epbins_default[i][j].pTtot = 0; + epbins_default[i][j].used = false; + epbins_default[i][j].ntracks = 0; + epbins_default[i][j].nxtracks = 0; + epbins_default[i][j].trackidx.clear(); + } + } + + // logic: loop over z bins find tracks in this bin and arrange them in a 2D eta-phi matrix + for (unsigned int zbin = 0; zbin < zmins.size(); ++zbin) { + // initialize matrices for every z bin + z0_intern zmin = zmins[zbin]; + z0_intern zmax = zmaxs[zbin]; + TrackJetEmulationEtaPhiBin epbins[phiBins_][etaBins_]; + std::copy(&epbins_default[0][0], &epbins_default[0][0] + phiBins_ * etaBins_, &epbins[0][0]); + + //clear containers + L1clusters.clear(); + L2clusters.clear(); + + // fill grid + for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) { + //// conversions + //-z0 + z0_intern trkZ = L1TrkPtrs_[k]->getZ0Word(); + + if (zmax < trkZ) + continue; + if (zmin > trkZ) + continue; + if (zbin == 0 && zmin == trkZ) + continue; + + //-pt + ap_uint ptEmulationBits = L1TrkPtrs_[k]->getRinvWord(); + pt_intern trkpt; + trkpt.V = ptEmulationBits.range(); + + //-eta + TTTrack_TrackWord::tanl_t etaEmulationBits = L1TrkPtrs_[k]->getTanlWord(); + glbeta_intern trketa; + trketa.V = etaEmulationBits.range(); + + //-phi + int Sector = L1TrkPtrs_[k]->phiSector(); + double sector_phi_value = 0; + if (Sector < 5) { + sector_phi_value = 2.0 * M_PI * Sector / 9.0; + } else { + sector_phi_value = (-1.0 * M_PI + M_PI / 9.0 + (Sector - 5) * 2.0 * M_PI / 9.0); + } + + glbphi_intern trkphiSector = DoubleToBit(sector_phi_value, + TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, + TTTrack_TrackWord::stepPhi0); + glbphi_intern local_phi = 0; + local_phi.V = L1TrkPtrs_[k]->getPhiWord(); + glbphi_intern local_phi2 = + DoubleToBit(BitToDouble(local_phi, TTTrack_TrackWord::TrackBitWidths::kPhiSize, TTTrack_TrackWord::stepPhi0), + TTTrack_TrackWord::TrackBitWidths::kPhiSize + kExtraGlobalPhiBit, + TTTrack_TrackWord::stepPhi0); + glbphi_intern trkphi = local_phi2 + trkphiSector; + + //-d0 + d0_intern abs_trkD0 = L1TrkPtrs_[k]->getD0Word(); + + //-nstub + int trk_nstubs = (int)L1TrkPtrs_[k]->getStubRefs().size(); + + // now fill the 2D grid with tracks + for (int i = 0; i < phiBins_; ++i) { + for (int j = 0; j < etaBins_; ++j) { + glbeta_intern eta_min = epbins[i][j].eta - etaStep_ / 2; //eta min + glbeta_intern eta_max = epbins[i][j].eta + etaStep_ / 2; //eta max + glbphi_intern phi_min = epbins[i][j].phi - phiStep_ / 2; //phi min + glbphi_intern phi_max = epbins[i][j].phi + phiStep_ / 2; //phi max + if ((trketa < eta_min) && j != 0) + continue; + if ((trketa > eta_max) && j != etaBins_ - 1) + continue; + if ((trkphi < phi_min) && i != 0) + continue; + if ((trkphi > phi_max) && i != phiBins_ - 1) + continue; + + if (trkpt < pt_intern(trkPtMax_)) + epbins[i][j].pTtot += trkpt; + else + epbins[i][j].pTtot += pt_intern(trkPtMax_); + if ((abs_trkD0 > + DoubleToBit(d0CutNStubs5_, TTTrack_TrackWord::TrackBitWidths::kD0Size, TTTrack_TrackWord::stepD0) && + trk_nstubs >= 5 && d0CutNStubs5_ >= 0) || + (abs_trkD0 > + DoubleToBit(d0CutNStubs4_, TTTrack_TrackWord::TrackBitWidths::kD0Size, TTTrack_TrackWord::stepD0) && + trk_nstubs == 4 && d0CutNStubs4_ >= 0)) + epbins[i][j].nxtracks += 1; + + epbins[i][j].trackidx.push_back(k); + ++epbins[i][j].ntracks; + } // for each etabin + } // for each phibin + } //end loop over tracks + + // first layer clustering - in eta using grid + for (int phibin = 0; phibin < phiBins_; ++phibin) { + L1clusters.push_back(L1_clustering( + epbins[phibin], etaBins_, etaStep_)); + } + + //second layer clustering in phi for all eta clusters + L2clusters = L2_clustering( + L1clusters, phiBins_, phiStep_, etaStep_); + + // sum all cluster pTs in this zbin to find max + pt_intern sum_pt = 0; + for (unsigned int k = 0; k < L2clusters.size(); ++k) { + if (L2clusters[k].pTtot > pt_intern(highpTJetThreshold_) && L2clusters[k].ntracks < lowpTJetMinTrackMultiplicity_) + continue; + if (L2clusters[k].pTtot > pt_intern(highpTJetThreshold_) && + L2clusters[k].ntracks < highpTJetMinTrackMultiplicity_) + continue; + + if (L2clusters[k].pTtot > pt_intern(minTrkJetpT_)) + sum_pt += L2clusters[k].pTtot; + } + if (sum_pt < mzb.ht) + continue; + // if ht is larger than previous max, this is the new vertex zbin + mzb.ht = sum_pt; + mzb.znum = zbin; + mzb.clusters = L2clusters; + mzb.nclust = L2clusters.size(); + mzb.zbincenter = (zmin + zmax) / 2.0; + } //zbin loop + + vector> L1TrackAssocJet; + for (unsigned int j = 0; j < mzb.clusters.size(); ++j) { + if (mzb.clusters[j].pTtot < pt_intern(trkPtMin_)) + continue; + + l1t::TkJetWord::glbeta_t jetEta = DoubleToBit(double(mzb.clusters[j].eta), + TkJetWord::TkJetBitWidths::kGlbEtaSize, + TkJetWord::MAX_ETA / (1 << TkJetWord::TkJetBitWidths::kGlbEtaSize)); + l1t::TkJetWord::glbphi_t jetPhi = DoubleToBit( + BitToDouble(mzb.clusters[j].phi, TTTrack_TrackWord::TrackBitWidths::kPhiSize + 4, TTTrack_TrackWord::stepPhi0), + TkJetWord::TkJetBitWidths::kGlbPhiSize, + (2. * std::abs(M_PI)) / (1 << TkJetWord::TkJetBitWidths::kGlbPhiSize)); + l1t::TkJetWord::z0_t jetZ0 = 0; + l1t::TkJetWord::pt_t jetPt = mzb.clusters[j].pTtot; + l1t::TkJetWord::nt_t total_ntracks = mzb.clusters[j].ntracks; + l1t::TkJetWord::nx_t total_disptracks = mzb.clusters[j].nxtracks; + l1t::TkJetWord::dispflag_t dispflag = 0; + l1t::TkJetWord::tkjetunassigned_t unassigned = 0; + + if (total_disptracks > nDisplacedTracks_ || total_disptracks == nDisplacedTracks_) + dispflag = 1; + L1TrackAssocJet.clear(); + for (unsigned int itrk = 0; itrk < mzb.clusters[j].trackidx.size(); itrk++) + L1TrackAssocJet.push_back(L1TrkPtrs_[mzb.clusters[j].trackidx[itrk]]); + + l1t::TkJetWord trkJet(jetPt, jetEta, jetPhi, jetZ0, total_ntracks, total_disptracks, dispflag, unassigned); + + L1TrackJetContainer->push_back(trkJet); + } + + std::sort(L1TrackJetContainer->begin(), L1TrackJetContainer->end(), [](auto &a, auto &b) { return a.pt() > b.pt(); }); + if (displaced_) + iEvent.put(std::move(L1TrackJetContainer), "L1TrackJetsExtended"); + else + iEvent.put(std::move(L1TrackJetContainer), "L1TrackJets"); +} + +void L1TrackJetEmulatorProducer::fillDescriptions(ConfigurationDescriptions &descriptions) { + //The following says we do not know what parameters are allowed so do no validation + // Please change this to state exactly what you do use, even if it is no parameters + ParameterSetDescription desc; + desc.add("L1TrackInputTag", edm::InputTag("l1tTTTracksFromTrackletEmulation", "Level1TTTracks")); + desc.add("L1PVertexInputTag", edm::InputTag("l1tVertexFinderEmulator", "l1verticesEmulation")); + desc.add("MaxDzTrackPV", 1.0); + desc.add("trk_zMax", 15.0); + desc.add("trk_ptMax", 200.0); + desc.add("trk_ptMin", 3.0); + desc.add("trk_etaMax", 2.4); + desc.add("nStubs4PromptChi2", 5.0); + desc.add("nStubs4PromptBend", 1.7); + desc.add("nStubs5PromptChi2", 2.75); + desc.add("nStubs5PromptBend", 3.5); + desc.add("trk_nPSStubMin", -1); + desc.add("minTrkJetpT", -1.0); + desc.add("etaBins", 24); + desc.add("phiBins", 27); + desc.add("zBins", 1); + desc.add("d0_cutNStubs4", -1); + desc.add("d0_cutNStubs5", -1); + desc.add("lowpTJetMinTrackMultiplicity", 2); + desc.add("lowpTJetThreshold", 50.0); + desc.add("highpTJetMinTrackMultiplicity", 3); + desc.add("highpTJetThreshold", 100.0); + desc.add("displaced", false); + desc.add("nStubs4DisplacedChi2", 5.0); + desc.add("nStubs4DisplacedBend", 1.7); + desc.add("nStubs5DisplacedChi2", 2.75); + desc.add("nStubs5DisplacedBend", 3.5); + desc.add("nDisplacedTracks", 2); + descriptions.add("l1tTrackJetsEmulator", desc); +} + +//define this as a plug-in +DEFINE_FWK_MODULE(L1TrackJetEmulatorProducer); diff --git a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetProducer.cc b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetProducer.cc index 32c40ba9d37a4..9a93ddb686ea6 100644 --- a/L1Trigger/L1TTrackMatch/plugins/L1TrackJetProducer.cc +++ b/L1Trigger/L1TTrackMatch/plugins/L1TrackJetProducer.cc @@ -1,20 +1,28 @@ -// Original Author: Rishi Patel -// Modifications: George Karathanasis, georgios.karathanasis@cern.ch, CU Boulder +// Original simulation Author: Rishi Patel +// +// Rewritting/improvements: George Karathanasis, +// georgios.karathanasis@cern.ch, CU Boulder +// // Created: Wed, 01 Aug 2018 14:01:41 GMT -// Latest update: Nov 2021 (by GK) +// Latest update: Nov 2022 (by GK) // // Track jets are clustered in a two-layer process, first by clustering in phi, -// then by clustering in eta +// then by clustering in eta. The code proceeds as following: putting all tracks +// in a grid of eta vs phi space, and then cluster them. Finally we merge the cl +// usters when needed. The code is improved to use the same module between emula +// tion and simulation was also improved, with bug fixes and being faster. // Introduction to object (p10-13): // https://indico.cern.ch/event/791517/contributions/3341650/attachments/1818736/2973771/TrackBasedAlgos_L1TMadrid_MacDonald.pdf +// New and improved version: https://indico.cern.ch/event/1203796/contributions/5073056/attachments/2519806/4333006/trackjet_emu.pdf // system include files - #include "DataFormats/Common/interface/Ref.h" #include "DataFormats/L1TCorrelator/interface/TkJet.h" #include "DataFormats/L1TCorrelator/interface/TkJetFwd.h" #include "DataFormats/L1TrackTrigger/interface/TTTypes.h" #include "DataFormats/L1TrackTrigger/interface/TTTrack.h" +#include "DataFormats/L1Trigger/interface/VertexWord.h" + // user include files #include "FWCore/Framework/interface/Frameworkfwd.h" #include "FWCore/Framework/interface/stream/EDProducer.h" @@ -24,42 +32,25 @@ #include "FWCore/ParameterSet/interface/ParameterSet.h" #include "FWCore/Utilities/interface/StreamID.h" #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h" -#include "Geometry/CommonTopologies/interface/PixelGeomDetUnit.h" -#include "Geometry/CommonTopologies/interface/PixelGeomDetType.h" -#include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h" - -#include "DataFormats/L1Trigger/interface/Vertex.h" -#include "L1Trigger/L1TTrackMatch/interface/L1Clustering.h" - -#include "TH1D.h" -#include "TH2D.h" -#include -#include -#include -#include -#include -#include -#include +//own headers +#include "L1TrackJetClustering.h" using namespace std; using namespace edm; using namespace l1t; +using namespace l1ttrackjet; class L1TrackJetProducer : public stream::EDProducer<> { public: explicit L1TrackJetProducer(const ParameterSet &); - ~L1TrackJetProducer() override; + ~L1TrackJetProducer() override = default; typedef TTTrack L1TTTrackType; typedef vector L1TTTrackCollectionType; - static void fillDescriptions(ConfigurationDescriptions &descriptions); - bool trackQualityCuts(int trk_nstub, float trk_chi2, float trk_bendchi2); private: - void beginStream(StreamID) override; void produce(Event &, const EventSetup &) override; - void endStream() override; // ----------member data --------------------------- @@ -97,39 +88,39 @@ class L1TrackJetProducer : public stream::EDProducer<> { edm::ESGetToken tTopoToken_; const EDGetTokenT>> trackToken_; - const edm::EDGetTokenT> PVtxToken_; + const EDGetTokenT PVtxToken_; }; L1TrackJetProducer::L1TrackJetProducer(const ParameterSet &iConfig) - : trkZMax_((float)iConfig.getParameter("trk_zMax")), - trkPtMax_((float)iConfig.getParameter("trk_ptMax")), - trkPtMin_((float)iConfig.getParameter("trk_ptMin")), - trkEtaMax_((float)iConfig.getParameter("trk_etaMax")), - nStubs4PromptChi2_((float)iConfig.getParameter("nStubs4PromptChi2")), - nStubs5PromptChi2_((float)iConfig.getParameter("nStubs5PromptChi2")), - nStubs4PromptBend_((float)iConfig.getParameter("nStubs4PromptBend")), - nStubs5PromptBend_((float)iConfig.getParameter("nStubs5PromptBend")), - trkNPSStubMin_((int)iConfig.getParameter("trk_nPSStubMin")), - lowpTJetMinTrackMultiplicity_((int)iConfig.getParameter("lowpTJetMinTrackMultiplicity")), - lowpTJetThreshold_((float)iConfig.getParameter("lowpTJetThreshold")), - highpTJetMinTrackMultiplicity_((int)iConfig.getParameter("highpTJetMinTrackMultiplicity")), - highpTJetThreshold_((float)iConfig.getParameter("highpTJetThreshold")), - zBins_((int)iConfig.getParameter("zBins")), - etaBins_((int)iConfig.getParameter("etaBins")), - phiBins_((int)iConfig.getParameter("phiBins")), + : trkZMax_(iConfig.getParameter("trk_zMax")), + trkPtMax_(iConfig.getParameter("trk_ptMax")), + trkPtMin_(iConfig.getParameter("trk_ptMin")), + trkEtaMax_(iConfig.getParameter("trk_etaMax")), + nStubs4PromptChi2_(iConfig.getParameter("nStubs4PromptChi2")), + nStubs5PromptChi2_(iConfig.getParameter("nStubs5PromptChi2")), + nStubs4PromptBend_(iConfig.getParameter("nStubs4PromptBend")), + nStubs5PromptBend_(iConfig.getParameter("nStubs5PromptBend")), + trkNPSStubMin_(iConfig.getParameter("trk_nPSStubMin")), + lowpTJetMinTrackMultiplicity_(iConfig.getParameter("lowpTJetMinTrackMultiplicity")), + lowpTJetThreshold_(iConfig.getParameter("lowpTJetThreshold")), + highpTJetMinTrackMultiplicity_(iConfig.getParameter("highpTJetMinTrackMultiplicity")), + highpTJetThreshold_(iConfig.getParameter("highpTJetThreshold")), + zBins_(iConfig.getParameter("zBins")), + etaBins_(iConfig.getParameter("etaBins")), + phiBins_(iConfig.getParameter("phiBins")), minTrkJetpT_(iConfig.getParameter("minTrkJetpT")), displaced_(iConfig.getParameter("displaced")), - d0CutNStubs4_((float)iConfig.getParameter("d0_cutNStubs4")), - d0CutNStubs5_((float)iConfig.getParameter("d0_cutNStubs5")), - nStubs4DisplacedChi2_((float)iConfig.getParameter("nStubs4DisplacedChi2")), - nStubs5DisplacedChi2_((float)iConfig.getParameter("nStubs5DisplacedChi2")), - nStubs4DisplacedBend_((float)iConfig.getParameter("nStubs4DisplacedBend")), - nStubs5DisplacedBend_((float)iConfig.getParameter("nStubs5DisplacedBend")), - nDisplacedTracks_((int)iConfig.getParameter("nDisplacedTracks")), - dzPVTrk_((float)iConfig.getParameter("MaxDzTrackPV")), + d0CutNStubs4_(iConfig.getParameter("d0_cutNStubs4")), + d0CutNStubs5_(iConfig.getParameter("d0_cutNStubs5")), + nStubs4DisplacedChi2_(iConfig.getParameter("nStubs4DisplacedChi2")), + nStubs5DisplacedChi2_(iConfig.getParameter("nStubs5DisplacedChi2")), + nStubs4DisplacedBend_(iConfig.getParameter("nStubs4DisplacedBend")), + nStubs5DisplacedBend_(iConfig.getParameter("nStubs5DisplacedBend")), + nDisplacedTracks_(iConfig.getParameter("nDisplacedTracks")), + dzPVTrk_(iConfig.getParameter("MaxDzTrackPV")), tTopoToken_(esConsumes(edm::ESInputTag("", ""))), trackToken_(consumes>>(iConfig.getParameter("L1TrackInputTag"))), - PVtxToken_(consumes>(iConfig.getParameter("L1PVertexCollection"))) { + PVtxToken_(consumes(iConfig.getParameter("L1PVertexInputTag"))) { zStep_ = 2.0 * trkZMax_ / (zBins_ + 1); // added +1 in denom etaStep_ = 2.0 * trkEtaMax_ / etaBins_; //etaStep is the width of an etabin phiStep_ = 2 * M_PI / phiBins_; ////phiStep is the width of a phibin @@ -140,8 +131,6 @@ L1TrackJetProducer::L1TrackJetProducer(const ParameterSet &iConfig) produces("L1TrackJets"); } -L1TrackJetProducer::~L1TrackJetProducer() {} - void L1TrackJetProducer::produce(Event &iEvent, const EventSetup &iSetup) { unique_ptr L1L1TrackJetProducer(new TkJetCollection); @@ -151,13 +140,14 @@ void L1TrackJetProducer::produce(Event &iEvent, const EventSetup &iSetup) { edm::Handle>> TTTrackHandle; iEvent.getByToken(trackToken_, TTTrackHandle); - edm::Handle> PVtx; + edm::Handle PVtx; iEvent.getByToken(PVtxToken_, PVtx); float PVz = (PVtx->at(0)).z0(); L1TrkPtrs_.clear(); tdtrk_.clear(); + // track selection for (unsigned int this_l1track = 0; this_l1track < TTTrackHandle->size(); this_l1track++) { edm::Ptr trkPtr(TTTrackHandle, this_l1track); float trk_pt = trkPtr->momentum().perp(); @@ -165,7 +155,6 @@ void L1TrackJetProducer::produce(Event &iEvent, const EventSetup &iSetup) { float trk_chi2dof = trkPtr->chi2Red(); float trk_d0 = trkPtr->d0(); float trk_bendchi2 = trkPtr->stubPtConsistency(); - float trk_z0 = trkPtr->z0(); int trk_nPS = 0; for (int istub = 0; istub < trk_nstubs; istub++) { // loop over the stubs @@ -179,11 +168,22 @@ void L1TrackJetProducer::produce(Event &iEvent, const EventSetup &iSetup) { // select tracks if (trk_nPS < trkNPSStubMin_) continue; - if (!trackQualityCuts(trk_nstubs, trk_chi2dof, trk_bendchi2)) + if (!TrackQualitySelection(trk_nstubs, + trk_chi2dof, + trk_bendchi2, + nStubs4PromptBend_, + nStubs5PromptBend_, + nStubs4PromptChi2_, + nStubs5PromptChi2_, + nStubs4DisplacedBend_, + nStubs5DisplacedBend_, + nStubs4DisplacedChi2_, + nStubs5DisplacedChi2_, + displaced_)) continue; - if (std::abs(PVz - trk_z0) > dzPVTrk_ && dzPVTrk_ > 0) + if (std::abs(PVz - trkPtr->z0()) > dzPVTrk_ && dzPVTrk_ > 0) continue; - if (std::abs(trk_z0) > trkZMax_) + if (std::abs(trkPtr->z0()) > trkZMax_) continue; if (std::abs(trkPtr->momentum().eta()) > trkEtaMax_) continue; @@ -198,167 +198,165 @@ void L1TrackJetProducer::produce(Event &iEvent, const EventSetup &iSetup) { tdtrk_.push_back(0); // not displaced track } - if (!L1TrkPtrs_.empty()) { - MaxZBin mzb; - mzb.znum = -1; - mzb.nclust = -1; - mzb.ht = -1; - EtaPhiBin epbins_default[phiBins_][etaBins_]; // create grid of phiBins - - float phi = -1.0 * M_PI; - for (int i = 0; i < phiBins_; ++i) { - float eta = -1.0 * trkEtaMax_; - for (int j = 0; j < etaBins_; ++j) { - epbins_default[i][j].phi = (phi + (phi + phiStep_)) / 2.0; // phimin=phi; phimax= phimin+step - epbins_default[i][j].eta = (eta + (eta + etaStep_)) / 2.0; // phimin=phi; phimax phimin+step - eta += etaStep_; - } // for each etabin - phi += phiStep_; - } // for each phibin (finished creating epbins) - - std::vector zmins, zmaxs; - for (int zbin = 0; zbin < zBins_; zbin++) { - zmins.push_back(-1.0 * trkZMax_ + zStep_ * zbin); - zmaxs.push_back(-1.0 * trkZMax_ + zStep_ * zbin + zStep_ * 2); + // if no tracks pass selection return empty containers + if (L1TrkPtrs_.empty()) { + if (displaced_) + iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended"); + else + iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets"); + return; + } + + MaxZBin mzb; + mzb.znum = -1; + mzb.nclust = -1; + mzb.ht = -1; + + // create 2D grid of eta phi bins + EtaPhiBin epbins_default[phiBins_][etaBins_]; + float phi = -1.0 * M_PI; + for (int i = 0; i < phiBins_; ++i) { + float eta = -1.0 * trkEtaMax_; + for (int j = 0; j < etaBins_; ++j) { + epbins_default[i][j].phi = (phi + (phi + phiStep_)) / 2.0; + epbins_default[i][j].eta = (eta + (eta + etaStep_)) / 2.0; + eta += etaStep_; + } // for each etabin + phi += phiStep_; + } // for each phibin (finished creating bins) + + // create z-bins (might be useful for displaced if we run w/o dz cut) + std::vector zmins, zmaxs; + for (int zbin = 0; zbin < zBins_; zbin++) { + zmins.push_back(-1.0 * trkZMax_ + zStep_ * zbin); + zmaxs.push_back(-1.0 * trkZMax_ + zStep_ * zbin + zStep_ * 2); + } + + // create vectors of clusters + std::vector> L1clusters; + L1clusters.reserve(phiBins_); + std::vector L2clusters; + + // default (empty) grid + for (int i = 0; i < phiBins_; ++i) { + for (int j = 0; j < etaBins_; ++j) { + epbins_default[i][j].pTtot = 0; + epbins_default[i][j].used = false; + epbins_default[i][j].ntracks = 0; + epbins_default[i][j].nxtracks = 0; + epbins_default[i][j].trackidx.clear(); } + } - // create vectors that hold data - std::vector> L1clusters; - L1clusters.reserve(phiBins_); - std::vector L2clusters; - - for (unsigned int zbin = 0; zbin < zmins.size(); ++zbin) { - // initialize matrices - float zmin = zmins[zbin]; - float zmax = zmaxs[zbin]; - EtaPhiBin epbins[phiBins_][etaBins_]; - std::copy(&epbins_default[0][0], &epbins_default[0][0] + phiBins_ * etaBins_, &epbins[0][0]); - - //clear containers - L1clusters.clear(); - L2clusters.clear(); - - // fill grid - for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) { - float trkZ = L1TrkPtrs_[k]->z0(); - if (zmax < trkZ) - continue; - if (zmin > trkZ) - continue; - if (zbin == 0 && zmin == trkZ) - continue; - float trkpt = L1TrkPtrs_[k]->momentum().perp(); - float trketa = L1TrkPtrs_[k]->momentum().eta(); - float trkphi = L1TrkPtrs_[k]->momentum().phi(); - for (int i = 0; i < phiBins_; ++i) { - for (int j = 0; j < etaBins_; ++j) { - float eta_min = epbins[i][j].eta - etaStep_ / 2.0; //eta min - float eta_max = epbins[i][j].eta + etaStep_ / 2.0; //eta max - float phi_min = epbins[i][j].phi - phiStep_ / 2.0; //phi min - float phi_max = epbins[i][j].phi + phiStep_ / 2.0; //phi max - - if ((trketa < eta_min) || (trketa > eta_max) || (trkphi < phi_min) || (trkphi > phi_max)) - continue; - - if (trkpt < trkPtMax_) - epbins[i][j].pTtot += trkpt; - else - epbins[i][j].pTtot += trkPtMax_; - epbins[i][j].numtdtrks += tdtrk_[k]; - epbins[i][j].trackidx.push_back(k); - ++epbins[i][j].numtracks; - } // for each phibin - } // for each phibin - } //end loop over tracks - - // first layer clustering - in eta for all phi bins - for (int phibin = 0; phibin < phiBins_; ++phibin) { - L1clusters.push_back(L1_clustering(epbins[phibin], etaBins_, etaStep_)); - } + for (unsigned int zbin = 0; zbin < zmins.size(); ++zbin) { + // initialize grid + float zmin = zmins[zbin]; + float zmax = zmaxs[zbin]; + EtaPhiBin epbins[phiBins_][etaBins_]; + std::copy(&epbins_default[0][0], &epbins_default[0][0] + phiBins_ * etaBins_, &epbins[0][0]); + + //clear cluster containers + L1clusters.clear(); + L2clusters.clear(); + + // fill grid with tracks + for (unsigned int k = 0; k < L1TrkPtrs_.size(); ++k) { + float trkZ = L1TrkPtrs_[k]->z0(); + if (zmax < trkZ) + continue; + if (zmin > trkZ) + continue; + if (zbin == 0 && zmin == trkZ) + continue; + float trkpt = L1TrkPtrs_[k]->momentum().perp(); + float trketa = L1TrkPtrs_[k]->momentum().eta(); + float trkphi = L1TrkPtrs_[k]->momentum().phi(); + for (int i = 0; i < phiBins_; ++i) { + for (int j = 0; j < etaBins_; ++j) { + float eta_min = epbins[i][j].eta - etaStep_ / 2.0; //eta min + float eta_max = epbins[i][j].eta + etaStep_ / 2.0; //eta max + float phi_min = epbins[i][j].phi - phiStep_ / 2.0; //phi min + float phi_max = epbins[i][j].phi + phiStep_ / 2.0; //phi max + + if ((trketa < eta_min) || (trketa > eta_max) || (trkphi < phi_min) || (trkphi > phi_max)) + continue; + + if (trkpt < trkPtMax_) + epbins[i][j].pTtot += trkpt; + else + epbins[i][j].pTtot += trkPtMax_; + epbins[i][j].nxtracks += tdtrk_[k]; + epbins[i][j].trackidx.push_back(k); + ++epbins[i][j].ntracks; + } // for each etabin + } // for each phibin + } //end loop over tracks + + // cluster tracks in eta (first layer) using grid + for (int phibin = 0; phibin < phiBins_; ++phibin) { + L1clusters.push_back(L1_clustering(epbins[phibin], etaBins_, etaStep_)); + } - //second layer clustering in phi for all eta clusters - L2clusters = L2_clustering(L1clusters, phiBins_, phiStep_, etaStep_); - - // sum all cluster pTs in this zbin to find max - float sum_pt = 0; - for (unsigned int k = 0; k < L2clusters.size(); ++k) { - if (L2clusters[k].pTtot > lowpTJetThreshold_ && L2clusters[k].numtracks < lowpTJetMinTrackMultiplicity_) - continue; - if (L2clusters[k].pTtot > highpTJetThreshold_ && L2clusters[k].numtracks < highpTJetMinTrackMultiplicity_) - continue; - if (L2clusters[k].pTtot > minTrkJetpT_) - sum_pt += L2clusters[k].pTtot; - } + // second layer clustering in phi for using eta clusters + L2clusters = L2_clustering(L1clusters, phiBins_, phiStep_, etaStep_); - if (sum_pt < mzb.ht) + // sum all cluster pTs in this zbin to find max + float sum_pt = 0; + for (unsigned int k = 0; k < L2clusters.size(); ++k) { + if (L2clusters[k].pTtot > lowpTJetThreshold_ && L2clusters[k].ntracks < lowpTJetMinTrackMultiplicity_) + continue; + if (L2clusters[k].pTtot > highpTJetThreshold_ && L2clusters[k].ntracks < highpTJetMinTrackMultiplicity_) continue; - // if ht is larger than previous max, this is the new vertex zbin - mzb.ht = sum_pt; - mzb.znum = zbin; - mzb.clusters = L2clusters; - mzb.nclust = L2clusters.size(); - mzb.zbincenter = (zmin + zmax) / 2.0; - } //zbin loop - - vector> L1TrackAssocJet; - for (unsigned int j = 0; j < mzb.clusters.size(); ++j) { - float jetEta = mzb.clusters[j].eta; - float jetPhi = mzb.clusters[j].phi; - float jetPt = mzb.clusters[j].pTtot; - float jetPx = jetPt * cos(jetPhi); - float jetPy = jetPt * sin(jetPhi); - float jetPz = jetPt * sinh(jetEta); - float jetP = jetPt * cosh(jetEta); - int totalDisptrk = mzb.clusters[j].numtdtrks; - bool isDispJet = (totalDisptrk > nDisplacedTracks_ || totalDisptrk == nDisplacedTracks_); - - math::XYZTLorentzVector jetP4(jetPx, jetPy, jetPz, jetP); - L1TrackAssocJet.clear(); - for (unsigned int itrk = 0; itrk < mzb.clusters[j].trackidx.size(); itrk++) - L1TrackAssocJet.push_back(L1TrkPtrs_[mzb.clusters[j].trackidx[itrk]]); - - TkJet trkJet(jetP4, L1TrackAssocJet, mzb.zbincenter, mzb.clusters[j].numtracks, 0, totalDisptrk, 0, isDispJet); - - if (!L1TrackAssocJet.empty()) - L1L1TrackJetProducer->push_back(trkJet); + if (L2clusters[k].pTtot > minTrkJetpT_) + sum_pt += L2clusters[k].pTtot; } - if (displaced_) - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended"); - else - iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets"); - } -} + if (sum_pt < mzb.ht) + continue; -void L1TrackJetProducer::beginStream(StreamID) {} - -void L1TrackJetProducer::endStream() {} - -bool L1TrackJetProducer::trackQualityCuts(int trk_nstub, float trk_chi2, float trk_bendchi2) { - bool PassQuality = false; - if (!displaced_) { - if (trk_nstub == 4 && trk_bendchi2 < nStubs4PromptBend_ && - trk_chi2 < nStubs4PromptChi2_) // 4 stubs are the lowest track quality and have different cuts - PassQuality = true; - if (trk_nstub > 4 && trk_bendchi2 < nStubs5PromptBend_ && - trk_chi2 < nStubs5PromptChi2_) // above 4 stubs diffent selection imposed (genrally looser) - PassQuality = true; - } else { - if (trk_nstub == 4 && trk_bendchi2 < nStubs4DisplacedBend_ && - trk_chi2 < nStubs4DisplacedChi2_) // 4 stubs are the lowest track quality and have different cuts - PassQuality = true; - if (trk_nstub > 4 && trk_bendchi2 < nStubs5DisplacedBend_ && - trk_chi2 < nStubs5DisplacedChi2_) // above 4 stubs diffent selection imposed (genrally looser) - PassQuality = true; + // if ht is larger than previous max, this is the new vertex zbin + mzb.ht = sum_pt; + mzb.znum = zbin; + mzb.clusters = L2clusters; + mzb.nclust = L2clusters.size(); + mzb.zbincenter = (zmin + zmax) / 2.0; + } //zbin loop + + // output + vector> L1TrackAssocJet; + for (unsigned int j = 0; j < mzb.clusters.size(); ++j) { + float jetEta = mzb.clusters[j].eta; + float jetPhi = mzb.clusters[j].phi; + float jetPt = mzb.clusters[j].pTtot; + float jetPx = jetPt * cos(jetPhi); + float jetPy = jetPt * sin(jetPhi); + float jetPz = jetPt * sinh(jetEta); + float jetP = jetPt * cosh(jetEta); + int totalDisptrk = mzb.clusters[j].nxtracks; + bool isDispJet = (totalDisptrk > nDisplacedTracks_ || totalDisptrk == nDisplacedTracks_); + + math::XYZTLorentzVector jetP4(jetPx, jetPy, jetPz, jetP); + L1TrackAssocJet.clear(); + for (unsigned int itrk = 0; itrk < mzb.clusters[j].trackidx.size(); itrk++) + L1TrackAssocJet.push_back(L1TrkPtrs_[mzb.clusters[j].trackidx[itrk]]); + + TkJet trkJet(jetP4, L1TrackAssocJet, mzb.zbincenter, mzb.clusters[j].ntracks, 0, totalDisptrk, 0, isDispJet); + + L1L1TrackJetProducer->push_back(trkJet); } - return PassQuality; + + std::sort( + L1L1TrackJetProducer->begin(), L1L1TrackJetProducer->end(), [](auto &a, auto &b) { return a.pt() > b.pt(); }); + if (displaced_) + iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJetsExtended"); + else + iEvent.put(std::move(L1L1TrackJetProducer), "L1TrackJets"); } void L1TrackJetProducer::fillDescriptions(ConfigurationDescriptions &descriptions) { - // l1tTrackJets - edm::ParameterSetDescription desc; + ParameterSetDescription desc; desc.add("L1TrackInputTag", edm::InputTag("l1tTTTracksFromTrackletEmulation", "Level1TTTracks")); - desc.add("L1PVertexCollection", edm::InputTag("l1tVertexProducer", "l1vertices")); + desc.add("L1PVertexInputTag", edm::InputTag("l1tVertexFinderEmulator", "l1verticesEmulation")); desc.add("MaxDzTrackPV", 1.0); desc.add("trk_zMax", 15.0); desc.add("trk_ptMax", 200.0); diff --git a/L1Trigger/L1TTrackMatch/python/l1tTrackJetsEmulation_cfi.py b/L1Trigger/L1TTrackMatch/python/l1tTrackJetsEmulation_cfi.py index 97b8651af4bdd..9276870a56a70 100644 --- a/L1Trigger/L1TTrackMatch/python/l1tTrackJetsEmulation_cfi.py +++ b/L1Trigger/L1TTrackMatch/python/l1tTrackJetsEmulation_cfi.py @@ -1,59 +1,46 @@ import FWCore.ParameterSet.Config as cms -l1tTrackJetsEmulation = cms.EDProducer('L1TrackJetEmulationProducer', +l1tTrackJetsEmulation = cms.EDProducer('L1TrackJetEmulatorProducer', L1TrackInputTag= cms.InputTag("l1tGTTInputProducer", "Level1TTTracksConverted"), - VertexInputTag=cms.InputTag("l1tVertexFinderEmulator", "l1verticesEmulation"), - MaxDzTrackPV = cms.double(0.5), + L1PVertexInputTag= cms.InputTag("l1tVertexFinderEmulator","l1verticesEmulation"), + MaxDzTrackPV = cms.double(1.0), trk_zMax = cms.double (15.) , # maximum track z trk_ptMax = cms.double(200.), # maximumum track pT before saturation [GeV] trk_ptMin = cms.double(2.0), # minimum track pt [GeV] trk_etaMax = cms.double(2.4), # maximum track eta - trk_chi2dofMax=cms.double(10.), # maximum track chi2/dof - trk_bendChi2Max=cms.double(2.2), # maximum track bendchi2 + nStubs4PromptChi2=cms.double(10.0), #Prompt track quality flags for loose/tight + nStubs4PromptBend=cms.double(2.2), + nStubs5PromptChi2=cms.double(10.0), + nStubs5PromptBend=cms.double(2.2), trk_nPSStubMin=cms.int32(-1), # minimum PS stubs, -1 means no cut - minTrkJetpT=cms.double(5.), # minimum track pt to be considered for track jet + minTrkJetpT=cms.double(-1.), # minimum track pt to be considered for track jet etaBins=cms.int32(24), phiBins=cms.int32(27), zBins=cms.int32(1), - d0_cutNStubs4=cms.double(0.15), - d0_cutNStubs5=cms.double(0.5), + d0_cutNStubs4=cms.double(-1), + d0_cutNStubs5=cms.double(-1), lowpTJetMinTrackMultiplicity=cms.int32(2), - lowpTJetMinpT=cms.double(50.), + lowpTJetThreshold=cms.double(50.), highpTJetMinTrackMultiplicity=cms.int32(3), - highpTJetMinpT=cms.double(100.), + highpTJetThreshold=cms.double(100.), displaced=cms.bool(False), #Flag for displaced tracks nStubs4DisplacedChi2=cms.double(5.0), #Displaced track quality flags for loose/tight - nStubs4Displacedbend=cms.double(1.7), + nStubs4DisplacedBend=cms.double(1.7), nStubs5DisplacedChi2=cms.double(2.75), - nStubs5Displacedbend=cms.double(3.5), + nStubs5DisplacedBend=cms.double(3.5), nDisplacedTracks=cms.int32(2) #Number of displaced tracks required per jet ) -l1tTrackJetsExtendedEmulation = cms.EDProducer('L1TrackJetEmulationProducer', - L1TrackInputTag= cms.InputTag("l1tGTTInputProducerExtended", "Level1TTTracksExtendedConverted"), - VertexInputTag=cms.InputTag("l1tVertexFinderEmulator", "l1verticesEmulation"), - MaxDzTrackPV = cms.double(4.0), - trk_zMax = cms.double (15.) , # maximum track z - trk_ptMax = cms.double(200.), # maximumum track pT before saturation [GeV] - trk_ptMin = cms.double(3.0), # minimum track pt [GeV] - trk_etaMax = cms.double(2.4), # maximum track eta - trk_chi2dofMax=cms.double(40.), # maximum track chi2/dof - trk_bendChi2Max=cms.double(40.), # maximum track bendchi2 - trk_nPSStubMin=cms.int32(-1), # minimum # PS stubs, -1 means no cut - minTrkJetpT=cms.double(5.), # minimum track pt to be considered for track jet - etaBins=cms.int32(24), - phiBins=cms.int32(27), - zBins=cms.int32(10), - d0_cutNStubs4=cms.double(-1), # -1 excludes nstub=4 from disp tag - d0_cutNStubs5=cms.double(0.22), - lowpTJetMinTrackMultiplicity=cms.int32(2), - lowpTJetMinpT=cms.double(50.), - highpTJetMinTrackMultiplicity=cms.int32(3), - highpTJetMinpT=cms.double(100.), - displaced=cms.bool(True), #Flag for displaced tracks - nStubs4DisplacedChi2=cms.double(3.3), #Disp tracks selection [trkcut excluded diff --git a/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py b/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py index f8cffe057f9dd..f4f0ea8e8d126 100644 --- a/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py +++ b/L1Trigger/L1TTrackMatch/test/L1TrackObjectNtupleMaker_cfg.py @@ -25,8 +25,8 @@ process.load('Configuration.StandardSequences.Services_cff') process.load('Configuration.EventContent.EventContent_cff') process.load('Configuration.StandardSequences.MagneticField_cff') -process.load('Configuration.Geometry.GeometryExtended2026D77Reco_cff') -process.load('Configuration.Geometry.GeometryExtended2026D77_cff') +process.load('Configuration.Geometry.GeometryExtended2026D95Reco_cff') +process.load('Configuration.Geometry.GeometryExtended2026D95_cff') process.load('Configuration.StandardSequences.EndOfProcess_cff') process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff') @@ -43,15 +43,7 @@ process.maxEvents = cms.untracked.PSet(input = cms.untracked.int32(10)) readFiles = cms.untracked.vstring( - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/c6df2819-ed05-4b98-8f92-81b7d1b1092e.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/3f476d95-1ef7-4be6-977b-6bcd1a7c5678.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/68d651da-4cb7-4bf4-b002-66aecc57a2bc.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/db0e0ce2-4c5a-4988-9dbd-52066e40b9d2.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/257a9712-0a96-47b7-897e-f5d980605e46.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/bee31399-8559-4243-b539-cae1ea897def.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/24629540-2377-4168-9ae5-518ddd4c43a9.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/e31ba8f0-332a-4a1a-8bc0-91a12a5fe3db.root', - '/store/relval/CMSSW_12_3_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_123X_mcRun4_realistic_v3_2026D77PU200-v1/2580000/17902198-4db6-4fcc-9e8c-787991b4db32.root', + '/store/relval/CMSSW_13_0_0/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/130X_mcRun4_realistic_v2_2026D95noPU-v1/00000/16f6615d-f98c-475f-ad33-0e89934b6c7f.root' ) secFiles = cms.untracked.vstring() @@ -67,7 +59,7 @@ useJobReport = cms.untracked.bool(True) ) -process.TFileService = cms.Service("TFileService", fileName = cms.string('GTTObjects_ttbar200PU.root'), closeFileFast = cms.untracked.bool(True)) +process.TFileService = cms.Service("TFileService", fileName = cms.string('GTTObjects_ttbar200PU_v2p2.root'), closeFileFast = cms.untracked.bool(True)) ############################################################ @@ -86,9 +78,8 @@ # DTC emulation -process.load('L1Trigger.TrackerDTC.ProducerES_cff') process.load('L1Trigger.TrackerDTC.ProducerED_cff') -process.dtc = cms.Path(process.TrackerDTCProducer)#*process.TrackerDTCAnalyzer) +process.dtc = cms.Path(process.TrackerDTCProducer) process.load("L1Trigger.TrackFindingTracklet.L1HybridEmulationTracks_cff") process.load("L1Trigger.L1TTrackMatch.l1tTrackSelectionProducer_cfi") @@ -116,8 +107,8 @@ process.l1tTrackFastJets.L1PrimaryVertexTag = cms.InputTag("l1tVertexFinder", "l1vertices") process.l1tTrackFastJetsExtended.L1PrimaryVertexTag = cms.InputTag("l1tVertexFinder", "l1vertices") -process.l1tTrackJets.L1PVertexCollection = cms.InputTag("l1tVertexFinder", "l1vertices") -process.l1tTrackJetsExtended.L1PVertexCollection = cms.InputTag("l1tVertexFinder", "l1vertices") +process.l1tTrackJets.L1PVertexInputTag = cms.InputTag("l1tVertexFinderEmulator","l1verticesEmulation") +process.l1tTrackJetsExtended.L1PVertexInputTag = cms.InputTag("l1tVertexFinderEmulator","l1verticesEmulation") process.l1tTrackerEtMiss.L1VertexInputTag = cms.InputTag("l1tVertexFinder", "l1vertices") process.l1tTrackerHTMiss.L1VertexInputTag = cms.InputTag("l1tVertexFinder", "l1vertices") process.l1tTrackerEtMissExtended.L1VertexInputTag = cms.InputTag("l1tVertexFinder", "l1vertices") @@ -150,7 +141,6 @@ process.pL1GTTInput = cms.Path(process.l1tGTTInputProducerExtended) process.pL1TrackJetsEmu = cms.Path(process.l1tTrackJetsExtendedEmulation) process.pTkMET = cms.Path(process.l1tTrackerEtMissExtended) - #process.pTkMETEmu = cms.Path(process.L1TrackerEmuEtMissExtended) #Doesn't exist process.pTkMHT = cms.Path(process.l1tTrackerHTMissExtended) process.pTkMHTEmulator = cms.Path(process.l1tTrackerEmuHTMissExtended) DISPLACED = 'Displaced'# @@ -235,9 +225,10 @@ process.ntuple = cms.Path(process.L1TrackNtuple) process.out = cms.OutputModule( "PoolOutputModule", - fastCloning = cms.untracked.bool( False ), + outputCommands = process.RAWSIMEventContent.outputCommands, fileName = cms.untracked.string("test.root" ) ) +process.out.outputCommands.append('keep *TTTrack*_*_*_*') process.pOut = cms.EndPath(process.out) @@ -247,5 +238,7 @@ # use this if cluster/stub associators not available # process.schedule = cms.Schedule(process.TTClusterStubTruth,process.TTTracksEmuWithTruth,process.ntuple) -process.schedule = cms.Schedule(process.TTClusterStub, process.TTClusterStubTruth, process.dtc, process.TTTracksEmuWithTruth, process.pL1GTTInput, process.pPV, process.pPVemu, process.pL1TrackSelection, process.pL1TrackJets, process.pL1TrackJetsEmu, -process.pL1TrackFastJets, process.pTkMET, process.pTkMETEmu, process.pTkMHT, process.pTkMHTEmulator, process.ntuple) +process.schedule = cms.Schedule(process.TTClusterStub, process.TTClusterStubTruth, process.dtc, process.TTTracksEmuWithTruth, process.pL1GTTInput, process.pPV, process.pPVemu, process.pL1TrackSelection, process.pL1TrackJets, process.pL1TrackJetsEmu,process.pL1TrackFastJets, process.pTkMET, process.pTkMETEmu, process.pTkMHT, process.pTkMHTEmulator, process.ntuple) + + + diff --git a/L1Trigger/TrackFindingTracklet/test/L1TrackNtupleMaker_cfg.py b/L1Trigger/TrackFindingTracklet/test/L1TrackNtupleMaker_cfg.py index 17a89466a4c2f..eeb0ea60cd8ac 100644 --- a/L1Trigger/TrackFindingTracklet/test/L1TrackNtupleMaker_cfg.py +++ b/L1Trigger/TrackFindingTracklet/test/L1TrackNtupleMaker_cfg.py @@ -11,9 +11,14 @@ # edit options here ############################################################ -GEOMETRY = "D49" -# Set L1 tracking algorithm: -# 'HYBRID' (baseline, 4par fit) or 'HYBRID_DISPLACED' (extended, 5par fit). +# D76 used for old CMSSW_11_3 MC datasets. D88 used for CMSSW_12_6 datasets. +#GEOMETRY = "D76" +GEOMETRY = "D88" + +# Set L1 tracking algorithm: +# 'HYBRID' (baseline, 4par fit) or 'HYBRID_DISPLACED' (extended, 5par fit). +# 'HYBRID_NEWKF' (baseline, 4par fit, with bit-accurate KF emulation), +# 'HYBRID_REDUCED' to use the "Summer Chain" configuration with reduced inputs. # (Or legacy algos 'TMTT' or 'TRACKLET'). L1TRKALGO = 'HYBRID' @@ -31,14 +36,12 @@ process.MessageLogger.L1track = dict(limit = -1) process.MessageLogger.Tracklet = dict(limit = -1) -if GEOMETRY == "D49": - print "using geometry " + GEOMETRY + " (tilted)" - process.load('Configuration.Geometry.GeometryExtended2026D49Reco_cff') - process.load('Configuration.Geometry.GeometryExtended2026D49_cff') -elif GEOMETRY == "D76": - print "using geometry " + GEOMETRY + " (tilted)" - process.load('Configuration.Geometry.GeometryExtended2026D76Reco_cff') - process.load('Configuration.Geometry.GeometryExtended2026D76_cff') +if GEOMETRY == "D76" or GEOMETRY == "D88": + # Use D88 for both, as both correspond to identical CMS Tracker design, and D76 + # unavailable in CMSSW_12_6_0. + print("using geometry " + GEOMETRY + " (tilted)") + process.load('Configuration.Geometry.GeometryExtended2026D88Reco_cff') + process.load('Configuration.Geometry.GeometryExtended2026D88_cff') else: print "this is not a valid geometry!!!" @@ -61,13 +64,13 @@ #from MCsamples.Scripts.getCMSdata_cfi import * #from MCsamples.Scripts.getCMSlocaldata_cfi import * -if GEOMETRY == "D49": +if GEOMETRY == "D76": # Read data from card files (defines getCMSdataFromCards()): #from MCsamples.RelVal_1120.PU200_TTbar_14TeV_cfi import * #inputMC = getCMSdataFromCards() # Or read .root files from directory on local computer: - #dirName = "$myDir/whatever/" + #dirName = "$scratchmc/MCsamples1130_D76/RelVal/TTbar/PU200/" #inputMC=getCMSlocaldata(dirName) # Or read specified dataset (accesses CMS DB, so use this method only occasionally): @@ -79,11 +82,30 @@ elif GEOMETRY == "D76": inputMC = ["/store/relval/CMSSW_11_3_0_pre6/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_113X_mcRun4_realistic_v6_2026D76PU200-v1/00000/00026541-6200-4eed-b6f8-d3a1fd720e9c.root"] + +elif GEOMETRY == "D88": + + # Read specified .root file: + inputMC = ["/store/relval/CMSSW_12_6_0_pre4/RelValTTbar_14TeV/GEN-SIM-DIGI-RAW/PU_125X_mcRun4_realistic_v2_2026D88PU200-v1/2590000/00b3d04b-4c7b-4506-8d82-9538fb21ee19.root"] + else: print "this is not a valid geometry!!!" process.source = cms.Source("PoolSource", fileNames = cms.untracked.vstring(*inputMC)) +if GEOMETRY == "D76": + # If reading old MC dataset, drop incompatible EDProducts. + process.source.dropDescendantsOfDroppedBranches = cms.untracked.bool(False) + process.source.inputCommands = cms.untracked.vstring() + process.source.inputCommands.append('keep *_*_*Level1TTTracks*_*') + process.source.inputCommands.append('keep *_*_*StubAccepted*_*') + process.source.inputCommands.append('keep *_*_*ClusterAccepted*_*') + process.source.inputCommands.append('keep *_*_*MergedTrackTruth*_*') + process.source.inputCommands.append('keep *_genParticles_*_*') + +# Use skipEvents to select particular single events for test vectors +#process.source.skipEvents = cms.untracked.uint32(11) + process.TFileService = cms.Service("TFileService", fileName = cms.string('TTbar_PU200_'+GEOMETRY+'.root'), closeFileFast = cms.untracked.bool(True)) process.Timing = cms.Service("Timing", summaryOnly = cms.untracked.bool(True)) @@ -93,6 +115,8 @@ ############################################################ process.load('L1Trigger.TrackTrigger.TrackTrigger_cff') +process.load('L1Trigger.TrackerTFP.ProducerES_cff') +process.load('L1Trigger.TrackerTFP.ProducerLayerEncoding_cff') # remake stubs? #from L1Trigger.TrackTrigger.TTStubAlgorithmRegister_cfi import * @@ -142,11 +166,32 @@ L1TRK_NAME = "l1tTTTracksFromExtendedTrackletEmulation" L1TRK_LABEL = "Level1TTTracks" L1TRUTH_NAME = "TTTrackAssociatorFromPixelDigisExtended" - -# LEGACY ALGORITHM (EXPERTS ONLY): TRACKLET + +# HYBRID_NEWKF: prompt tracking or reduced +elif (L1TRKALGO == 'HYBRID_NEWKF' or L1TRKALGO == 'HYBRID_REDUCED'): + process.load( 'L1Trigger.TrackFindingTracklet.Producer_cff' ) + NHELIXPAR = 4 + L1TRK_NAME = process.TrackFindingTrackletProducer_params.LabelTT.value() + L1TRK_LABEL = process.TrackFindingTrackletProducer_params.BranchAcceptedTracks.value() + L1TRUTH_NAME = "TTTrackAssociatorFromPixelDigis" + process.TTTrackAssociatorFromPixelDigis.TTTracks = cms.VInputTag( cms.InputTag(L1TRK_NAME, L1TRK_LABEL) ) + process.HybridNewKF = cms.Sequence(process.L1THybridTracks + process.TrackFindingTrackletProducerTBout + process.TrackFindingTrackletProducerKFin + process.TrackFindingTrackletProducerKF + process.TrackFindingTrackletProducerTT + process.TrackFindingTrackletProducerAS + process.TrackFindingTrackletProducerKFout) + process.TTTracksEmulation = cms.Path(process.HybridNewKF) + #process.TTTracksEmulationWithTruth = cms.Path(process.HybridNewKF + process.TrackTriggerAssociatorTracks) + # Optionally include code producing performance plots & end-of-job summary. + process.load( 'SimTracker.TrackTriggerAssociation.StubAssociator_cff' ) + process.load( 'L1Trigger.TrackFindingTracklet.Analyzer_cff' ) + process.TTTracksEmulationWithTruth = cms.Path(process.HybridNewKF + process.TrackTriggerAssociatorTracks + process.StubAssociator + process.TrackFindingTrackletAnalyzerTracklet + process.TrackFindingTrackletAnalyzerTBout + process.TrackFindingTrackletAnalyzerKFin + process.TrackFindingTrackletAnalyzerKF + process.TrackFindingTrackletAnalyzerKFout) + from L1Trigger.TrackFindingTracklet.Customize_cff import * + if (L1TRKALGO == 'HYBRID_NEWKF'): + fwConfig( process ) + if (L1TRKALGO == 'HYBRID_REDUCED'): + reducedConfig( process ) + +# LEGACY ALGORITHM (EXPERTS ONLY): TRACKLET elif (L1TRKALGO == 'TRACKLET'): - print "\n WARNING: This is not the baseline algorithm! Prefer HYBRID or HYBRID_DISPLACED!" - print "\n To run the Tracklet-only algorithm, ensure you have commented out 'CXXFLAGS=-DUSEHYBRID' in BuildFile.xml & recompiled! \n" + print("\n WARNING: This is not the baseline algorithm! Prefer HYBRID or HYBRID_DISPLACED!") + print("\n To run the Tracklet-only algorithm, ensure you have commented out 'CXXFLAGS=-DUSEHYBRID' in BuildFile.xml & recompiled! \n") process.TTTracksEmulation = cms.Path(process.L1THybridTracks) process.TTTracksEmulationWithTruth = cms.Path(process.L1THybridTracksWithAssociators) NHELIXPAR = 4 @@ -248,4 +293,6 @@ process.pd = cms.EndPath(process.writeDataset) process.schedule.append(process.pd) - + + + From 9097cfdfe458f5c7b27103d8605724dcf85e0a16 Mon Sep 17 00:00:00 2001 From: Emyr Clement Date: Mon, 15 May 2023 16:46:26 +0100 Subject: [PATCH 12/16] Ensure cluster ecal_eta/phi are set when not running CombinedCaloLinker. --- L1Trigger/Phase2L1ParticleFlow/src/CaloClusterer.cc | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/L1Trigger/Phase2L1ParticleFlow/src/CaloClusterer.cc b/L1Trigger/Phase2L1ParticleFlow/src/CaloClusterer.cc index 269bcc7de53c0..c42d8614a0128 100644 --- a/L1Trigger/Phase2L1ParticleFlow/src/CaloClusterer.cc +++ b/L1Trigger/Phase2L1ParticleFlow/src/CaloClusterer.cc @@ -589,6 +589,8 @@ void l1tpf_calo::SimpleCaloLinker::run() { float wecal = cluster.ecal_et / cluster.et, whcal = 1.0 - wecal; cluster.eta = ecal.eta * wecal + hcal.eta * whcal; cluster.phi = ecal.phi * wecal + hcal.phi * whcal; + cluster.ecal_eta = cluster.eta; + cluster.ecal_phi = cluster.phi; // wrap around phi cluster.phi = reco::reduceRange(cluster.phi); cluster.constituents.emplace_back(-i - 1, 1); @@ -617,6 +619,8 @@ void l1tpf_calo::SimpleCaloLinker::run() { cluster.et = myet + etot; cluster.eta = hcal.eta + avg_eta / cluster.et; cluster.phi = hcal.phi + avg_phi / cluster.et; + cluster.ecal_eta = cluster.eta; + cluster.ecal_phi = cluster.phi; // wrap around phi cluster.phi = reco::reduceRange(cluster.phi); } @@ -676,6 +680,8 @@ void l1tpf_calo::FlatCaloLinker::run() { dst.et = src.et; dst.eta = src.eta; dst.phi = src.phi; + dst.ecal_eta = src.eta; + dst.ecal_phi = src.phi; dst.ecal_et = 0; dst.hcal_et = 0; for (const auto &pair : src.constituents) { From 032ebd46537667cb20543ae2d8c896d376cc8212 Mon Sep 17 00:00:00 2001 From: SimranjitSingh Date: Tue, 20 Sep 2022 15:01:13 +0200 Subject: [PATCH 13/16] pftkegsorter_barrel emulator added --- .../egamma/pftkegsorter_barrel_ref.h | 202 ++++++++++++++++++ 1 file changed, 202 insertions(+) create mode 100644 L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h diff --git a/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h new file mode 100644 index 0000000000000..a1aa8fa4ee90d --- /dev/null +++ b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h @@ -0,0 +1,202 @@ +#ifndef REF_PFTKEGSORTER_BARREL_REF_H +#define REF_PFTKEGSORTER_BARREL_REF_H + +#include +#include + +#include "DataFormats/L1TParticleFlow/interface/layer1_emulator.h" +#include "L1Trigger/Phase2L1ParticleFlow/interface/common/bitonic_hybrid_sort_ref.h" + +#ifdef CMSSW_GIT_HASH +#include "FWCore/ParameterSet/interface/ParameterSet.h" +#endif + +namespace edm { + class ParameterSet; +} + +namespace l1ct { + class PFTkEGSorterBarrelEmulator { + public: + PFTkEGSorterBarrelEmulator(const unsigned int nObjToSort = 10, const unsigned int nObjSorted = 16) + : nObjToSort_(nObjToSort), nObjSorted_(nObjSorted), debug_(false) {} + + PFTkEGSorterBarrelEmulator(const edm::ParameterSet& iConfig); + + virtual ~PFTkEGSorterBarrelEmulator() {} + + void setDebug(bool debug = true) { debug_ = debug; }; + + void toFirmware_pho(const OutputRegion& outregions, EGIsoObj (&photons_in)[nObjSorted_]) const { + for (unsigned int i = 0; i < nObjToSort_; i++) { + EGIsoObj pho; + if (i < outregions.egphoton.size()) { + pho = outregions.egphoton[i]; + } else + pho.clear(); + + photons_in[i] = pho; + } + } + + void toFirmware_ele(const OutputRegion& outregions, EGIsoEleObj (&eles_in)[nObjSorted_]) const { + for (unsigned int i = 0; i < nObjToSort_; i++) { + EGIsoEleObj ele; + if (i < outregions.egelectron.size()) { + ele = outregions.egelectron[i]; + } else + ele.clear(); + + eles_in[i] = ele; + } + } + + template + void run(const std::vector& pfregions, + const std::vector& outregions, + const std::vector& region_index, + std::vector& eg_sorted_inBoard) { + // we copy to be able to resize them + std::vector> objs_in; + objs_in.reserve(nObjToSort_); + for (unsigned int i : region_index) { + std::vector objs; + extractEGObjEmu(pfregions[i].region, outregions[i], objs); + if (debug_) + dbgCout() << "objs size " << objs.size() << "\n"; + resize_input(objs); + objs_in.push_back(objs); + if (debug_) + dbgCout() << "objs (re)size and total objs size " << objs.size() << " " << objs_in.size() << "\n"; + } + + merge(objs_in, eg_sorted_inBoard); + + if (debug_) { + dbgCout() << "objs.size() size " << eg_sorted_inBoard.size() << "\n"; + for (const auto& out : eg_sorted_inBoard) + dbgCout() << "kinematics of sorted objects " << out.hwPt << " " << out.hwEta << " " << out.hwPhi << "\n"; + } + } + + private: + void extractEGObjEmu(const PFRegionEmu& region, + const l1ct::OutputRegion& outregion, + std::vector& eg) { + extractEGObjEmu(region, outregion.egphoton, eg); + } + void extractEGObjEmu(const PFRegionEmu& region, + const l1ct::OutputRegion& outregion, + std::vector& eg) { + extractEGObjEmu(region, outregion.egelectron, eg); + } + template + void extractEGObjEmu(const PFRegionEmu& region, + const std::vector& regional_objects, + std::vector& global_objects) { + for (const auto& reg_obj : regional_objects) { + global_objects.emplace_back(reg_obj); + global_objects.back().hwEta = region.hwGlbEta(reg_obj.hwEta); //<=========== uncomment this + global_objects.back().hwPhi = region.hwGlbPhi(reg_obj.hwPhi); //<=========== uncomment this + } + } + + template + void resize_input(std::vector& in) const { + if (in.size() > nObjToSort_) { + in.resize(nObjToSort_); + } else if (in.size() < nObjToSort_) { + for (unsigned int i = 0, diff = (nObjToSort_ - in.size()); i < diff; ++i) { + in.push_back(T()); + in.back().clear(); + } + } + } + + template + void merge_regions(const std::vector& in_region1, + const std::vector& in_region2, + std::vector& out, + unsigned int nOut) const { + // we crate a bitonic list + out = in_region1; + if (debug_) + for (const auto& tmp : out) + dbgCout() << "out " << tmp.hwPt << " " << tmp.hwEta << " " << tmp.hwPhi << "\n"; + std::reverse(out.begin(), out.end()); + if (debug_) + for (const auto& tmp : out) + dbgCout() << "out reverse " << tmp.hwPt << " " << tmp.hwEta << " " << tmp.hwPhi << "\n"; + std::copy(in_region2.begin(), in_region2.end(), std::back_inserter(out)); + if (debug_) + for (const auto& tmp : out) + dbgCout() << "out inserted " << tmp.hwPt << " " << tmp.hwEta << " " << tmp.hwPhi << "\n"; + + hybridBitonicMergeRef(&out[0], out.size(), 0, false); + + if (out.size() > nOut) { + out.resize(nOut); + if (debug_) + for (const auto& tmp : out) + dbgCout() << "final " << tmp.hwPt << " " << tmp.hwEta << " " << tmp.hwPhi << "\n"; + } + } + + template + void merge(const std::vector>& in_objs, std::vector& out) const { + if (in_objs.size() == 1) { //size is 1, fine! + std::copy(in_objs[0].begin(), in_objs[0].end(), std::back_inserter(out)); + if (out.size() > nObjSorted_) + out.resize(nObjSorted_); //size 16 + } else if (in_objs.size() == 2) { //size is 2, fine! + merge_regions(in_objs[0], in_objs[1], out, nObjSorted_); //10, 10, 16 + } else { + std::vector pair_merge_01; //size is >2, merge 0 and 1 regions always + merge_regions(in_objs[0], in_objs[1], pair_merge_01, nObjSorted_); //10, 10, 16 + + std::vector> to_merge; + if (in_objs.size() == 3) + to_merge.push_back(pair_merge_01); //push 01 only if size is 3 //and then in_objs[id] will be pushed into it + + std::vector pair_merge_tmp = pair_merge_01; // + for (unsigned int id = 2, idn = 3; id < in_objs.size(); id += 2, idn = id + 1) { + if (idn >= in_objs.size()) { //if size is odd number starting from 3 + to_merge.push_back(in_objs[id]); //size 10 + } else { + std::vector pair_merge; + merge_regions(in_objs[id], + in_objs[idn], + pair_merge, + nObjSorted_); //10, 10, 16 // merge two regions: 23, 45, 67, and so on + + //pair_merge_tmp.resize(nObjToSort_); + //pair_merge.resize(nObjToSort_); + merge_regions( + pair_merge_tmp, + pair_merge, + pair_merge_tmp, + nObjSorted_); //16, 16, 16 //merge 23 with 01 for first time // then merge 45, and then 67, and then 89, and so on + to_merge.push_back( + pair_merge_tmp); //push back 0123, 012345, 01234567, and so on (remember 01 is the 0th element) + } + } + if (in_objs.size() % 2 == 1) { + //to_merge[to_merge.size()-2].resize(nObjToSort_); + merge_regions( + to_merge[to_merge.size() - 2], + to_merge[to_merge.size() - 1], + out, + nObjSorted_); //16, 16, 16 //e.g. is size is 3, we merge 01 and 2, if size is 7, we merge 012345 and 6 + } else + out = + pair_merge_tmp; //if size is even number, e.g. 6, out is 012345 (or we can say to_merge[to_merge.size()-1]) + } + } + + unsigned int nObjToSort_; + unsigned int nObjSorted_; + bool debug_; + }; +} // namespace l1ct + +#endif From 13744ec7cdb96ac161bff3eb1ab76f9ae8e0932a Mon Sep 17 00:00:00 2001 From: Giovanni Date: Mon, 20 Mar 2023 14:50:25 +0100 Subject: [PATCH 14/16] Add track linking info in egamma debug printouts --- .../Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp b/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp index 2fc881879455d..89ee83ee6b481 100644 --- a/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp +++ b/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp @@ -135,6 +135,16 @@ void PFTkEGAlgoEmulator::link_emCalo2tk(const PFRegionEmu &r, float dEtaMax = cfg.dEtaValues[eta_index]; float dPhiMax = cfg.dPhiValues[eta_index]; + if (debug_ > 2 && calo.hwPt > 0) { + dbgCout() << "[REF] tried to link calo " << ic << " (pt " << calo.intPt() << ", eta " << calo.intEta() + << ", phi " << calo.intPhi() << ") " + << " to tk " << itk << " (pt " << tk.intPt() << ", eta " << tk.intEta() << ", phi " << tk.intPhi() + << "): " + << " eta_index " << eta_index << ", " + << " dEta " << d_eta << " (max " << dEtaMax << "), dPhi " << d_phi << " (max " << dPhiMax << ") " + << " ellipse = " + << (((d_phi / dPhiMax) * (d_phi / dPhiMax)) + ((d_eta / dEtaMax) * (d_eta / dEtaMax))) << "\n"; + } if ((((d_phi / dPhiMax) * (d_phi / dPhiMax)) + ((d_eta / dEtaMax) * (d_eta / dEtaMax))) < 1.) { // NOTE: for now we implement only best pt match. This is NOT what is done in the L1TkElectronTrackProducer if (fabs(tk.floatPt() - calo.floatPt()) < dPtMin) { From 3aa9d3c49370381ead2386dc4ba5a7da72a9295e Mon Sep 17 00:00:00 2001 From: Giovanni Date: Mon, 20 Mar 2023 11:26:40 +0100 Subject: [PATCH 15/16] Improve implementation, integrate with CMSSW --- .../egamma/pftkegsorter_barrel_ref.h | 120 ++++++++---------- .../interface/egamma/pftkegsorter_ref.h | 25 ++-- .../regionizer/tdr_regionizer_elements_ref.h | 4 - .../plugins/L1TCorrelatorLayer1Producer.cc | 16 ++- .../python/l1ctLayer1_cff.py | 4 + .../test/make_l1ctLayer1_dumpFiles_cfg.py | 25 +++- .../test/make_l1ct_binaryFiles_cfg.py | 24 +++- 7 files changed, 131 insertions(+), 87 deletions(-) diff --git a/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h index a1aa8fa4ee90d..2895848dd11fe 100644 --- a/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h +++ b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h @@ -1,33 +1,28 @@ -#ifndef REF_PFTKEGSORTER_BARREL_REF_H -#define REF_PFTKEGSORTER_BARREL_REF_H +#ifndef L1Trigger_Phase2L1ParticleFlow_egamma_pftkegsorter_barrel_ref_h +#define L1Trigger_Phase2L1ParticleFlow_egamma_pftkegsorter_barrel_ref_h #include #include +#include "pftkegsorter_ref.h" #include "DataFormats/L1TParticleFlow/interface/layer1_emulator.h" #include "L1Trigger/Phase2L1ParticleFlow/interface/common/bitonic_hybrid_sort_ref.h" -#ifdef CMSSW_GIT_HASH -#include "FWCore/ParameterSet/interface/ParameterSet.h" -#endif - -namespace edm { - class ParameterSet; -} - namespace l1ct { - class PFTkEGSorterBarrelEmulator { + class PFTkEGSorterBarrelEmulator : public PFTkEGSorterEmulator { public: PFTkEGSorterBarrelEmulator(const unsigned int nObjToSort = 10, const unsigned int nObjSorted = 16) - : nObjToSort_(nObjToSort), nObjSorted_(nObjSorted), debug_(false) {} + : PFTkEGSorterEmulator(nObjToSort, nObjSorted) {} - PFTkEGSorterBarrelEmulator(const edm::ParameterSet& iConfig); - - virtual ~PFTkEGSorterBarrelEmulator() {} +#ifdef CMSSW_GIT_HASH + PFTkEGSorterBarrelEmulator(const edm::ParameterSet& iConfig) + : PFTkEGSorterEmulator(iConfig.getParameter("nObjToSort"), + iConfig.getParameter("nObjSorted")) {} +#endif - void setDebug(bool debug = true) { debug_ = debug; }; + ~PFTkEGSorterBarrelEmulator() override {} - void toFirmware_pho(const OutputRegion& outregions, EGIsoObj (&photons_in)[nObjSorted_]) const { + void toFirmware_pho(const OutputRegion& outregions, EGIsoObj photons_in[/*nObjSorted_*/]) const { for (unsigned int i = 0; i < nObjToSort_; i++) { EGIsoObj pho; if (i < outregions.egphoton.size()) { @@ -39,7 +34,7 @@ namespace l1ct { } } - void toFirmware_ele(const OutputRegion& outregions, EGIsoEleObj (&eles_in)[nObjSorted_]) const { + void toFirmware_ele(const OutputRegion& outregions, EGIsoEleObj eles_in[/*nObjSorted_*/]) const { for (unsigned int i = 0; i < nObjToSort_; i++) { EGIsoEleObj ele; if (i < outregions.egelectron.size()) { @@ -51,6 +46,19 @@ namespace l1ct { } } + void runPho(const std::vector& pfregions, + const std::vector& outregions, + const std::vector& region_index, + std::vector& eg_sorted_inBoard) override { + run(pfregions, outregions, region_index, eg_sorted_inBoard); + } + void runEle(const std::vector& pfregions, + const std::vector& outregions, + const std::vector& region_index, + std::vector& eg_sorted_inBoard) override { + run(pfregions, outregions, region_index, eg_sorted_inBoard); + } + template void run(const std::vector& pfregions, const std::vector& outregions, @@ -96,8 +104,8 @@ namespace l1ct { std::vector& global_objects) { for (const auto& reg_obj : regional_objects) { global_objects.emplace_back(reg_obj); - global_objects.back().hwEta = region.hwGlbEta(reg_obj.hwEta); //<=========== uncomment this - global_objects.back().hwPhi = region.hwGlbPhi(reg_obj.hwPhi); //<=========== uncomment this + global_objects.back().hwEta = region.hwGlbEta(reg_obj.hwEta); + global_objects.back().hwPhi = region.hwGlbPhi(reg_obj.hwPhi); } } @@ -144,58 +152,32 @@ namespace l1ct { template void merge(const std::vector>& in_objs, std::vector& out) const { - if (in_objs.size() == 1) { //size is 1, fine! - std::copy(in_objs[0].begin(), in_objs[0].end(), std::back_inserter(out)); - if (out.size() > nObjSorted_) - out.resize(nObjSorted_); //size 16 - } else if (in_objs.size() == 2) { //size is 2, fine! - merge_regions(in_objs[0], in_objs[1], out, nObjSorted_); //10, 10, 16 - } else { - std::vector pair_merge_01; //size is >2, merge 0 and 1 regions always - merge_regions(in_objs[0], in_objs[1], pair_merge_01, nObjSorted_); //10, 10, 16 - - std::vector> to_merge; - if (in_objs.size() == 3) - to_merge.push_back(pair_merge_01); //push 01 only if size is 3 //and then in_objs[id] will be pushed into it - - std::vector pair_merge_tmp = pair_merge_01; // - for (unsigned int id = 2, idn = 3; id < in_objs.size(); id += 2, idn = id + 1) { - if (idn >= in_objs.size()) { //if size is odd number starting from 3 - to_merge.push_back(in_objs[id]); //size 10 - } else { - std::vector pair_merge; - merge_regions(in_objs[id], - in_objs[idn], - pair_merge, - nObjSorted_); //10, 10, 16 // merge two regions: 23, 45, 67, and so on - - //pair_merge_tmp.resize(nObjToSort_); - //pair_merge.resize(nObjToSort_); - merge_regions( - pair_merge_tmp, - pair_merge, - pair_merge_tmp, - nObjSorted_); //16, 16, 16 //merge 23 with 01 for first time // then merge 45, and then 67, and then 89, and so on - to_merge.push_back( - pair_merge_tmp); //push back 0123, 012345, 01234567, and so on (remember 01 is the 0th element) + unsigned int nregions = in_objs.size(); + std::vector pair_merge(nObjSorted_); + if (nregions == 18) { // merge pairs, and accumulate pairs + for (unsigned int i = 0; i < nregions; i += 2) { + merge_regions(in_objs[i + 0], in_objs[i + 1], pair_merge, nObjSorted_); + if (i == 0) + out = pair_merge; + else + merge_regions(out, pair_merge, out, nObjSorted_); + } + } else if (nregions == 9) { // simple accumulation + for (unsigned int i = 0; i < nregions; ++i) { + for (unsigned int j = 0, nj = in_objs[i].size(); j < nObjSorted_; ++j) { + if (j < nj) + pair_merge[j] = in_objs[i][j]; + else + pair_merge[j].clear(); } + if (i == 0) + out = pair_merge; + else + merge_regions(out, pair_merge, out, nObjSorted_); } - if (in_objs.size() % 2 == 1) { - //to_merge[to_merge.size()-2].resize(nObjToSort_); - merge_regions( - to_merge[to_merge.size() - 2], - to_merge[to_merge.size() - 1], - out, - nObjSorted_); //16, 16, 16 //e.g. is size is 3, we merge 01 and 2, if size is 7, we merge 012345 and 6 - } else - out = - pair_merge_tmp; //if size is even number, e.g. 6, out is 012345 (or we can say to_merge[to_merge.size()-1]) - } + } else + throw std::runtime_error("This sorter requires 18 or 9 regions"); } - - unsigned int nObjToSort_; - unsigned int nObjSorted_; - bool debug_; }; } // namespace l1ct diff --git a/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h index 4a417fd48d04a..fffb556552b0d 100644 --- a/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h +++ b/L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h @@ -1,5 +1,5 @@ -#ifndef FIRMWARE_PFTKEGSORTER_REF_H -#define FIRMWARE_PFTKEGSORTER_REF_H +#ifndef L1Trigger_Phase2L1ParticleFlow_egamma_pftkegsorter_ref_h +#define L1Trigger_Phase2L1ParticleFlow_egamma_pftkegsorter_ref_h #include #include @@ -11,10 +11,6 @@ #include "FWCore/ParameterSet/interface/ParameterSet.h" #endif -namespace edm { - class ParameterSet; -} - namespace l1ct { class PFTkEGSorterEmulator { public: @@ -28,10 +24,23 @@ namespace l1ct { #endif - ~PFTkEGSorterEmulator(){}; + virtual ~PFTkEGSorterEmulator() {} void setDebug(bool debug = true) { debug_ = debug; }; + virtual void runPho(const std::vector& pfregions, + const std::vector& outregions, + const std::vector& region_index, + std::vector& eg_sorted_inBoard) { + run(pfregions, outregions, region_index, eg_sorted_inBoard); + } + virtual void runEle(const std::vector& pfregions, + const std::vector& outregions, + const std::vector& region_index, + std::vector& eg_sorted_inBoard) { + run(pfregions, outregions, region_index, eg_sorted_inBoard); + } + template void run(const std::vector& pfregions, const std::vector& outregions, @@ -64,7 +73,7 @@ namespace l1ct { } } - private: + protected: unsigned int nObjToSort_, nObjSorted_; bool debug_; diff --git a/L1Trigger/Phase2L1ParticleFlow/interface/regionizer/tdr_regionizer_elements_ref.h b/L1Trigger/Phase2L1ParticleFlow/interface/regionizer/tdr_regionizer_elements_ref.h index ae76062083fea..ec194e127ef9c 100644 --- a/L1Trigger/Phase2L1ParticleFlow/interface/regionizer/tdr_regionizer_elements_ref.h +++ b/L1Trigger/Phase2L1ParticleFlow/interface/regionizer/tdr_regionizer_elements_ref.h @@ -8,11 +8,7 @@ #include #include -#ifdef CMSSW_GIT_HASH #include "L1Trigger/Phase2L1ParticleFlow/interface/dbgPrintf.h" -#else -#include "../../utils/dbgPrintf.h" -#endif namespace l1ct { namespace tdr_regionizer { diff --git a/L1Trigger/Phase2L1ParticleFlow/plugins/L1TCorrelatorLayer1Producer.cc b/L1Trigger/Phase2L1ParticleFlow/plugins/L1TCorrelatorLayer1Producer.cc index f5599d121d0b9..c27b0187ef408 100644 --- a/L1Trigger/Phase2L1ParticleFlow/plugins/L1TCorrelatorLayer1Producer.cc +++ b/L1Trigger/Phase2L1ParticleFlow/plugins/L1TCorrelatorLayer1Producer.cc @@ -31,6 +31,7 @@ #include "L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegalgo_ref.h" #include "L1Trigger/Phase2L1ParticleFlow/interface/pf/pfalgo_common_ref.h" #include "L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_ref.h" +#include "L1Trigger/Phase2L1ParticleFlow/interface/egamma/pftkegsorter_barrel_ref.h" #include "L1Trigger/Phase2L1ParticleFlow/interface/L1TCorrelatorLayer1PatternFileWriter.h" #include "DataFormats/L1TCorrelator/interface/TkElectron.h" @@ -236,8 +237,15 @@ L1TCorrelatorLayer1Producer::L1TCorrelatorLayer1Producer(const edm::ParameterSet l1tkegalgo_ = std::make_unique( l1ct::PFTkEGAlgoEmuConfig(iConfig.getParameter("tkEgAlgoParameters"))); - l1tkegsorter_ = - std::make_unique(iConfig.getParameter("tkEgSorterParameters")); + const std::string &egsortalgo = iConfig.getParameter("tkEgSorterAlgo"); + if (egsortalgo == "Barrel") { + l1tkegsorter_ = std::make_unique( + iConfig.getParameter("tkEgSorterParameters")); + } else if (egsortalgo == "Endcap") { + l1tkegsorter_ = + std::make_unique(iConfig.getParameter("tkEgSorterParameters")); + } else + throw cms::Exception("Configuration", "Unsupported tkEgSorterAlgo"); if (l1tkegalgo_->writeEgSta()) produces>("L1Eg"); @@ -439,8 +447,8 @@ void L1TCorrelatorLayer1Producer::produce(edm::Event &iEvent, const edm::EventSe // l1tkegsorter_->setDebug(true); for (auto &board : event_.board_out) { - l1tkegsorter_->run(event_.pfinputs, event_.out, board.region_index, board.egphoton); - l1tkegsorter_->run(event_.pfinputs, event_.out, board.region_index, board.egelectron); + l1tkegsorter_->runPho(event_.pfinputs, event_.out, board.region_index, board.egphoton); + l1tkegsorter_->runEle(event_.pfinputs, event_.out, board.region_index, board.egelectron); } // save PF into the event diff --git a/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py b/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py index ee6fc05d8726a..aa5ebd4147e39 100644 --- a/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py +++ b/L1Trigger/Phase2L1ParticleFlow/python/l1ctLayer1_cff.py @@ -106,6 +106,7 @@ nEMCALO_EGIN = 10, nEM_EGOUT = 10, ), + tkEgSorterAlgo = cms.string("Barrel"), tkEgSorterParameters=tkEgSorterParameters.clone( nObjToSort = 10 ), @@ -262,6 +263,7 @@ doEndcapHwQual=True, writeBeforeBremRecovery=False, writeEGSta=True), + tkEgSorterAlgo = cms.string("Endcap"), tkEgSorterParameters=tkEgSorterParameters.clone( nObjToSort = 5 ), @@ -349,6 +351,7 @@ doEndcapHwQual=True, writeBeforeBremRecovery=False, writeEGSta=True), + tkEgSorterAlgo = cms.string("Endcap"), tkEgSorterParameters=tkEgSorterParameters.clone( nObjToSort=5 ), @@ -432,6 +435,7 @@ nEM_EGOUT = 5, # to be defined doBremRecovery=True, writeEGSta=True), + tkEgSorterAlgo = cms.string("Endcap"), tkEgSorterParameters=tkEgSorterParameters.clone(), caloSectors = cms.VPSet( cms.PSet( diff --git a/L1Trigger/Phase2L1ParticleFlow/test/make_l1ctLayer1_dumpFiles_cfg.py b/L1Trigger/Phase2L1ParticleFlow/test/make_l1ctLayer1_dumpFiles_cfg.py index c8eddfb3b029c..72c77480e8d2f 100644 --- a/L1Trigger/Phase2L1ParticleFlow/test/make_l1ctLayer1_dumpFiles_cfg.py +++ b/L1Trigger/Phase2L1ParticleFlow/test/make_l1ctLayer1_dumpFiles_cfg.py @@ -38,6 +38,28 @@ from L1Trigger.Phase2L1GMT.gmt_cfi import l1tStandaloneMuons process.l1tSAMuonsGmt = l1tStandaloneMuons.clone() +process.l1tLayer1BarrelSerenity = process.l1tLayer1Barrel.clone() +process.l1tLayer1BarrelSerenity.regionizerAlgo = "MultififoBarrel" +process.l1tLayer1BarrelSerenity.regionizerAlgoParameters = cms.PSet( + barrelSetup = cms.string("Full54"), + useAlsoVtxCoords = cms.bool(True), + nClocks = cms.uint32(54), + nHCalLinks = cms.uint32(2), + nECalLinks = cms.uint32(1), + nTrack = cms.uint32(22), + nCalo = cms.uint32(15), + nEmCalo = cms.uint32(12), + nMu = cms.uint32(2)) +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nTrack = 22 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nSelCalo = 15 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nCalo = 15 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nAllNeutral = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nTrack = 22 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nIn = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nOut = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.finalSortAlgo = "FoldedHybrid" + + process.l1tLayer1Barrel9 = process.l1tLayer1Barrel.clone() process.l1tLayer1Barrel9.puAlgo.nFinalSort = 32 process.l1tLayer1Barrel9.regions[0].etaBoundaries = [ -1.5, -0.5, 0.5, 1.5 ] @@ -55,6 +77,7 @@ process.l1tGTTInputProducer + process.l1tVertexFinderEmulator + process.l1tLayer1Barrel + + process.l1tLayer1BarrelSerenity + process.l1tLayer1Barrel9 + process.l1tLayer1HGCal + process.l1tLayer1HGCalNoTK + @@ -63,7 +86,7 @@ process.runPF.associate(process.L1TLayer1TaskInputsTask) -for det in "Barrel", "Barrel9", "HGCal", "HGCalNoTK", "HF": +for det in "Barrel", "BarrelSerenity", "Barrel9", "HGCal", "HGCalNoTK", "HF": l1pf = getattr(process, 'l1tLayer1'+det) l1pf.dumpFileName = cms.untracked.string("TTbar_PU200_"+det+".dump") diff --git a/L1Trigger/Phase2L1ParticleFlow/test/make_l1ct_binaryFiles_cfg.py b/L1Trigger/Phase2L1ParticleFlow/test/make_l1ct_binaryFiles_cfg.py index 62125c730587a..1ad6414ee65bb 100644 --- a/L1Trigger/Phase2L1ParticleFlow/test/make_l1ct_binaryFiles_cfg.py +++ b/L1Trigger/Phase2L1ParticleFlow/test/make_l1ct_binaryFiles_cfg.py @@ -82,6 +82,27 @@ regions=cms.vuint32(*[6+9*ie+i for ie in range(3) for i in range(3)])), ) +process.l1tLayer1BarrelSerenity = process.l1tLayer1Barrel.clone() +process.l1tLayer1BarrelSerenity.regionizerAlgo = "MultififoBarrel" +process.l1tLayer1BarrelSerenity.regionizerAlgoParameters = cms.PSet( + barrelSetup = cms.string("Full54"), + useAlsoVtxCoords = cms.bool(True), + nClocks = cms.uint32(54), + nHCalLinks = cms.uint32(2), + nECalLinks = cms.uint32(1), + nTrack = cms.uint32(22), + nCalo = cms.uint32(15), + nEmCalo = cms.uint32(12), + nMu = cms.uint32(2)) +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nTrack = 22 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nSelCalo = 15 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nCalo = 15 +process.l1tLayer1BarrelSerenity.pfAlgoParameters.nAllNeutral = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nTrack = 22 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nIn = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.nOut = 27 +process.l1tLayer1BarrelSerenity.puAlgoParameters.finalSortAlgo = "FoldedHybrid" + from L1Trigger.Phase2L1ParticleFlow.l1ctLayer1_patternWriters_cff import * if not args.patternFilesOFF: process.l1tLayer1Barrel.patternWriters = cms.untracked.VPSet(*barrelWriterConfigs) @@ -95,6 +116,7 @@ process.l1tGTTInputProducer + process.l1tVertexFinderEmulator + process.l1tLayer1Barrel + + process.l1tLayer1BarrelSerenity + process.l1tLayer1Barrel9 + process.l1tLayer1HGCal + process.l1tLayer1HGCalNoTK + @@ -124,7 +146,7 @@ process.l1tLayer2SeedConeJetWriter.maxLinesPerFile = eventsPerFile_*54 if not args.dumpFilesOFF: - for det in "Barrel", "Barrel9", "HGCal", "HGCalNoTK", "HF": + for det in "Barrel", "BarrelSerenity", "Barrel9", "HGCal", "HGCalNoTK", "HF": l1pf = getattr(process, 'l1tLayer1'+det) l1pf.dumpFileName = cms.untracked.string("TTbar_PU200_"+det+".dump") From 952496ecb304dc03ede9c45e334fdf23b8831b1a Mon Sep 17 00:00:00 2001 From: Gianluca Date: Tue, 16 May 2023 18:18:24 +0200 Subject: [PATCH 16/16] implement review comments --- DataFormats/L1TCorrelator/src/classes_def.xml | 5 ++--- L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp | 4 ++-- 2 files changed, 4 insertions(+), 5 deletions(-) diff --git a/DataFormats/L1TCorrelator/src/classes_def.xml b/DataFormats/L1TCorrelator/src/classes_def.xml index 723849343a5a9..6b5acaf05cbc1 100644 --- a/DataFormats/L1TCorrelator/src/classes_def.xml +++ b/DataFormats/L1TCorrelator/src/classes_def.xml @@ -44,9 +44,8 @@ - - - + + diff --git a/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp b/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp index c149f57006096..984dca69f36a5 100644 --- a/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp +++ b/L1Trigger/Phase2L1ParticleFlow/src/egamma/pftkegalgo_ref.cpp @@ -185,9 +185,9 @@ void PFTkEGAlgoEmulator::link_emCalo2tk_composite(const PFRegionEmu &r, float d_phi = deltaPhi(tk.floatPhi(), calo.floatPhi()); float d_eta = tk.floatEta() - calo.floatEta(); // We only use it squared - float dR = std::sqrt((d_phi * d_phi) + (d_eta * d_eta)); + float dR2 = (d_phi * d_phi) + (d_eta * d_eta); - if (dR < 0.2) { + if (dR2 < 0.04) { // Only store indices, dR and dpT for now. The other quantities are computed only for the best nCandPerCluster. CompositeCandidate cand; cand.cluster_idx = ic;