File indexing completed on 2023-07-18 08:35:33
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #include <memory>
0011
0012
0013 #include "FWCore/Framework/interface/Frameworkfwd.h"
0014 #include "FWCore/Framework/interface/stream/EDProducer.h"
0015 #include "FWCore/Framework/interface/Event.h"
0016 #include "FWCore/Framework/interface/MakerMacros.h"
0017 #include "FWCore/Framework/interface/EventSetup.h"
0018 #include "DataFormats/Common/interface/Handle.h"
0019 #include "FWCore/Utilities/interface/InputTag.h"
0020 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0021 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0022 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0023 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0024 #include "DataFormats/Math/interface/LorentzVector.h"
0025
0026
0027 #include "DataFormats/L1TrackTrigger/interface/TTTypes.h"
0028 #include "DataFormats/L1TCorrelator/interface/TkJet.h"
0029 #include "DataFormats/L1TCorrelator/interface/TkJetFwd.h"
0030 #include "DataFormats/L1Trigger/interface/Vertex.h"
0031
0032
0033 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0034 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0035
0036
0037 #include "SimTracker/TrackTriggerAssociation/interface/TTTrackAssociationMap.h"
0038 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0039
0040 #include <fastjet/JetDefinition.hh>
0041
0042 #include <string>
0043 #include "TMath.h"
0044 #include "TH1.h"
0045
0046 using namespace l1t;
0047 using namespace edm;
0048 using namespace std;
0049
0050
0051
0052
0053
0054
0055
0056 class L1FastTrackingJetProducer : public edm::stream::EDProducer<> {
0057 public:
0058 typedef TTTrack<Ref_Phase2TrackerDigi_> L1TTTrackType;
0059 typedef std::vector<L1TTTrackType> L1TTTrackCollectionType;
0060
0061 explicit L1FastTrackingJetProducer(const edm::ParameterSet&);
0062 ~L1FastTrackingJetProducer() override;
0063
0064 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0065
0066 private:
0067 virtual void beginJob();
0068 void produce(edm::Event&, const edm::EventSetup&) override;
0069 virtual void endJob();
0070
0071
0072 const float trkZMax_;
0073 const float trkChi2dofMax_;
0074 const double trkBendChi2Max_;
0075 const float trkPtMin_;
0076 const float trkEtaMax_;
0077 const int trkNStubMin_;
0078 const int trkNPSStubMin_;
0079 const double deltaZ0Cut_;
0080 const double coneSize_;
0081 const bool doTightChi2_;
0082 const float trkPtTightChi2_;
0083 const float trkChi2dofTightChi2_;
0084 const bool displaced_;
0085 bool selectTrkMatchGenTight_;
0086 bool selectTrkMatchGenLoose_;
0087 bool selectTrkMatchGenOrPU_;
0088
0089 const edm::EDGetTokenT<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>> trackToken_;
0090 edm::EDGetTokenT<std::vector<l1t::Vertex>> pvToken_;
0091 const edm::EDGetTokenT<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>> genToken_;
0092 edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0093 edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tGeomToken_;
0094 };
0095
0096
0097 L1FastTrackingJetProducer::L1FastTrackingJetProducer(const edm::ParameterSet& iConfig)
0098 : trkZMax_((float)iConfig.getParameter<double>("trk_zMax")),
0099 trkChi2dofMax_((float)iConfig.getParameter<double>("trk_chi2dofMax")),
0100 trkBendChi2Max_(iConfig.getParameter<double>("trk_bendChi2Max")),
0101 trkPtMin_((float)iConfig.getParameter<double>("trk_ptMin")),
0102 trkEtaMax_((float)iConfig.getParameter<double>("trk_etaMax")),
0103 trkNStubMin_((int)iConfig.getParameter<int>("trk_nStubMin")),
0104 trkNPSStubMin_((int)iConfig.getParameter<int>("trk_nPSStubMin")),
0105 deltaZ0Cut_((float)iConfig.getParameter<double>("deltaZ0Cut")),
0106 coneSize_((float)iConfig.getParameter<double>("coneSize")),
0107 doTightChi2_(iConfig.getParameter<bool>("doTightChi2")),
0108 trkPtTightChi2_((float)iConfig.getParameter<double>("trk_ptTightChi2")),
0109 trkChi2dofTightChi2_((float)iConfig.getParameter<double>("trk_chi2dofTightChi2")),
0110 displaced_(iConfig.getParameter<bool>("displaced")),
0111 selectTrkMatchGenTight_(iConfig.getParameter<bool>("selectTrkMatchGenTight")),
0112 selectTrkMatchGenLoose_(iConfig.getParameter<bool>("selectTrkMatchGenLoose")),
0113 selectTrkMatchGenOrPU_(iConfig.getParameter<bool>("selectTrkMatchGenOrPU")),
0114 trackToken_(consumes<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>>(
0115 iConfig.getParameter<edm::InputTag>("L1TrackInputTag"))),
0116 pvToken_(consumes<std::vector<l1t::Vertex>>(iConfig.getParameter<edm::InputTag>("L1PrimaryVertexTag"))),
0117 genToken_(
0118 consumes<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>>(iConfig.getParameter<edm::InputTag>("GenInfo"))),
0119 tTopoToken_(esConsumes<TrackerTopology, TrackerTopologyRcd>(edm::ESInputTag("", ""))),
0120 tGeomToken_(esConsumes<TrackerGeometry, TrackerDigiGeometryRecord>(edm::ESInputTag("", ""))) {
0121 if (displaced_)
0122 produces<TkJetCollection>("L1FastTrackingJetsExtended");
0123 else
0124 produces<TkJetCollection>("L1FastTrackingJets");
0125 }
0126
0127
0128 L1FastTrackingJetProducer::~L1FastTrackingJetProducer() {}
0129
0130
0131 void L1FastTrackingJetProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0132 std::unique_ptr<TkJetCollection> L1FastTrackingJets(new TkJetCollection);
0133
0134
0135 edm::Handle<std::vector<TTTrack<Ref_Phase2TrackerDigi_>>> TTTrackHandle;
0136 iEvent.getByToken(trackToken_, TTTrackHandle);
0137 std::vector<TTTrack<Ref_Phase2TrackerDigi_>>::const_iterator iterL1Track;
0138
0139
0140 edm::Handle<TTTrackAssociationMap<Ref_Phase2TrackerDigi_>> MCTrkAssociation;
0141 if (!iEvent.isRealData())
0142 iEvent.getByToken(genToken_, MCTrkAssociation);
0143
0144
0145 const TrackerTopology& tTopo = iSetup.getData(tTopoToken_);
0146
0147
0148 edm::Handle<std::vector<l1t::Vertex>> L1VertexHandle;
0149 iEvent.getByToken(pvToken_, L1VertexHandle);
0150
0151 fastjet::JetDefinition jet_def(fastjet::antikt_algorithm, coneSize_);
0152 std::vector<fastjet::PseudoJet> JetInputs;
0153
0154 float recoVtx = L1VertexHandle->begin()->z0();
0155 unsigned int this_l1track = 0;
0156 for (iterL1Track = TTTrackHandle->begin(); iterL1Track != TTTrackHandle->end(); iterL1Track++) {
0157 this_l1track++;
0158 float trk_pt = iterL1Track->momentum().perp();
0159 float trk_z0 = iterL1Track->z0();
0160 float trk_chi2dof = iterL1Track->chi2Red();
0161 float trk_bendchi2 = iterL1Track->stubPtConsistency();
0162 std::vector<edm::Ref<edmNew::DetSetVector<TTStub<Ref_Phase2TrackerDigi_>>, TTStub<Ref_Phase2TrackerDigi_>>>
0163 theStubs = iterL1Track->getStubRefs();
0164 int trk_nstub = (int)theStubs.size();
0165
0166 if (std::abs(trk_z0) > trkZMax_)
0167 continue;
0168 if (std::abs(iterL1Track->momentum().eta()) > trkEtaMax_)
0169 continue;
0170 if (trk_pt < trkPtMin_)
0171 continue;
0172 if (trk_nstub < trkNStubMin_)
0173 continue;
0174 if (trk_chi2dof > trkChi2dofMax_)
0175 continue;
0176 if (trk_bendchi2 > trkBendChi2Max_)
0177 continue;
0178 if (doTightChi2_ && (trk_pt > trkPtTightChi2_ && trk_chi2dof > trkChi2dofTightChi2_))
0179 continue;
0180
0181 int trk_nPS = 0;
0182 for (int istub = 0; istub < trk_nstub; istub++) {
0183 DetId detId(theStubs.at(istub)->getDetId());
0184 bool tmp_isPS = false;
0185 if (detId.det() == DetId::Detector::Tracker) {
0186 if (detId.subdetId() == StripSubdetector::TOB && tTopo.tobLayer(detId) <= 3)
0187 tmp_isPS = true;
0188 else if (detId.subdetId() == StripSubdetector::TID && tTopo.tidRing(detId) <= 9)
0189 tmp_isPS = true;
0190 }
0191 if (tmp_isPS)
0192 trk_nPS++;
0193 }
0194 if (trk_nPS < trkNPSStubMin_)
0195 continue;
0196 if (std::abs(recoVtx - trk_z0) > deltaZ0Cut_)
0197 continue;
0198 if (!iEvent.isRealData()) {
0199 edm::Ptr<TTTrack<Ref_Phase2TrackerDigi_>> trk_ptr(TTTrackHandle, this_l1track);
0200
0201 if (!(MCTrkAssociation->isGenuine(trk_ptr)) && selectTrkMatchGenTight_)
0202 continue;
0203 if (!(MCTrkAssociation->isLooselyGenuine(trk_ptr) || MCTrkAssociation->isGenuine(trk_ptr)) &&
0204 selectTrkMatchGenLoose_)
0205 continue;
0206 if (!(MCTrkAssociation->isLooselyGenuine(trk_ptr) || MCTrkAssociation->isGenuine(trk_ptr) ||
0207 MCTrkAssociation->isCombinatoric(trk_ptr)) &&
0208 selectTrkMatchGenOrPU_)
0209 continue;
0210 }
0211
0212 fastjet::PseudoJet psuedoJet(iterL1Track->momentum().x(),
0213 iterL1Track->momentum().y(),
0214 iterL1Track->momentum().z(),
0215 iterL1Track->momentum().mag());
0216 JetInputs.push_back(psuedoJet);
0217 JetInputs.back().set_user_index(this_l1track - 1);
0218 }
0219
0220 fastjet::ClusterSequence cs(JetInputs, jet_def);
0221 std::vector<fastjet::PseudoJet> JetOutputs =
0222 fastjet::sorted_by_pt(cs.inclusive_jets(0));
0223
0224 for (unsigned int ijet = 0; ijet < JetOutputs.size(); ++ijet) {
0225 math::XYZTLorentzVector jetP4(
0226 JetOutputs[ijet].px(), JetOutputs[ijet].py(), JetOutputs[ijet].pz(), JetOutputs[ijet].modp());
0227 float sumpt = 0;
0228 float avgZ = 0;
0229 std::vector<edm::Ptr<L1TTTrackType>> L1TrackPtrs;
0230 std::vector<fastjet::PseudoJet> fjConstituents = fastjet::sorted_by_pt(cs.constituents(JetOutputs[ijet]));
0231
0232 for (unsigned int i = 0; i < fjConstituents.size(); ++i) {
0233 auto index = fjConstituents[i].user_index();
0234 edm::Ptr<L1TTTrackType> trkPtr(TTTrackHandle, index);
0235 L1TrackPtrs.push_back(trkPtr);
0236 sumpt = sumpt + trkPtr->momentum().perp();
0237 avgZ = avgZ + trkPtr->momentum().perp() * trkPtr->z0();
0238 }
0239 avgZ = avgZ / sumpt;
0240 edm::Ref<JetBxCollection> jetRef;
0241 TkJet trkJet(jetP4, jetRef, L1TrackPtrs, avgZ);
0242 L1FastTrackingJets->push_back(trkJet);
0243 }
0244
0245 if (displaced_)
0246 iEvent.put(std::move(L1FastTrackingJets), "L1FastTrackingJetsExtended");
0247 else
0248 iEvent.put(std::move(L1FastTrackingJets), "L1FastTrackingJets");
0249 }
0250
0251 void L1FastTrackingJetProducer::beginJob() {}
0252
0253 void L1FastTrackingJetProducer::endJob() {}
0254
0255 void L1FastTrackingJetProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0256 {
0257
0258 edm::ParameterSetDescription desc;
0259 desc.add<edm::InputTag>("L1TrackInputTag", edm::InputTag("l1tTTTracksFromTrackletEmulation", "Level1TTTracks"));
0260 desc.add<std::string>("L1PrimaryVertexTag", "L1Vertices");
0261 desc.add<edm::InputTag>("GenInfo", edm::InputTag("TTTrackAssociatorFromPixelDigis", "Level1TTTracks"));
0262 desc.add<double>("trk_zMax", 15.0);
0263 desc.add<double>("trk_chi2dofMax", 10.0);
0264 desc.add<double>("trk_bendChi2Max", 2.2);
0265 desc.add<double>("trk_ptMin", 2.0);
0266 desc.add<double>("trk_etaMax", 2.5);
0267 desc.add<int>("trk_nStubMin", 4);
0268 desc.add<int>("trk_nPSStubMin", -1);
0269 desc.add<double>("deltaZ0Cut", 0.5);
0270 desc.add<bool>("doTightChi2", true);
0271 desc.add<double>("trk_ptTightChi2", 20.0);
0272 desc.add<double>("trk_chi2dofTightChi2", 5.0);
0273 desc.add<double>("coneSize", 0.4);
0274 desc.add<bool>("displaced", false);
0275 desc.add<bool>("selectTrkMatchGenTight", true);
0276 desc.add<bool>("selectTrkMatchGenLoose", false);
0277 desc.add<bool>("selectTrkMatchGenOrPU", false);
0278 descriptions.add("L1FastTrackingJets", desc);
0279
0280
0281 }
0282
0283
0284
0285
0286
0287
0288
0289
0290
0291
0292
0293
0294
0295
0296
0297
0298
0299
0300
0301
0302 }
0303
0304
0305 DEFINE_FWK_MODULE(L1FastTrackingJetProducer);