File indexing completed on 2021-03-05 03:57:50
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011 #include "SeedingRegionByL1.h"
0012
0013 #include <algorithm>
0014 #include <set>
0015 #include <vector>
0016
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018
0019 ticl::SeedingRegionByL1::SeedingRegionByL1(const edm::ParameterSet &conf, edm::ConsumesCollector &sumes)
0020 : SeedingRegionAlgoBase(conf, sumes),
0021 l1TkEmsToken_(sumes.consumes<std::vector<l1t::TkEm>>(conf.getParameter<edm::InputTag>("l1TkEmColl"))),
0022 algoVerbosity_(conf.getParameter<int>("algo_verbosity")),
0023 minPt_(conf.getParameter<double>("minPt")),
0024 minAbsEta_(conf.getParameter<double>("minAbsEta")),
0025 maxAbsEta_(conf.getParameter<double>("maxAbsEta")),
0026 endcapScalings_(conf.getParameter<std::vector<double>>("endcapScalings")),
0027 quality_(conf.getParameter<int>("quality")) {}
0028
0029 void ticl::SeedingRegionByL1::makeRegions(const edm::Event &ev,
0030 const edm::EventSetup &es,
0031 std::vector<TICLSeedingRegion> &result) {
0032 auto l1TrkEms = ev.getHandle(l1TkEmsToken_);
0033 edm::ProductID l1tkemsId = l1TrkEms.id();
0034
0035 for (size_t indx = 0; indx < (*l1TrkEms).size(); indx++) {
0036 const auto &l1TrkEm = (*l1TrkEms)[indx];
0037 double offlinePt = this->tkEmOfflineEt(l1TrkEm.pt());
0038 if ((offlinePt < minPt_) || (std::abs(l1TrkEm.eta()) < minAbsEta_) || (std::abs(l1TrkEm.eta()) > maxAbsEta_) ||
0039 (l1TrkEm.EGRef()->hwQual() != quality_)) {
0040 continue;
0041 }
0042
0043 int iSide = int(l1TrkEm.eta() > 0);
0044 result.emplace_back(GlobalPoint(l1TrkEm.p4().X(), l1TrkEm.p4().Y(), l1TrkEm.p4().Z()),
0045 GlobalVector(l1TrkEm.px(), l1TrkEm.py(), l1TrkEm.pz()),
0046 iSide,
0047 indx,
0048 l1tkemsId);
0049 }
0050
0051 std::sort(result.begin(), result.end(), [](const TICLSeedingRegion &a, const TICLSeedingRegion &b) {
0052 return a.directionAtOrigin.perp2() > b.directionAtOrigin.perp2();
0053 });
0054 }
0055
0056 double ticl::SeedingRegionByL1::tkEmOfflineEt(double et) const {
0057 return (endcapScalings_.at(0) + et * endcapScalings_.at(1) + et * et * endcapScalings_.at(2));
0058 }
0059
0060 void ticl::SeedingRegionByL1::fillPSetDescription(edm::ParameterSetDescription &desc) {
0061 desc.add<edm::InputTag>("l1TkEmColl", edm::InputTag("L1TkPhotonsHGC", "EG"));
0062 desc.add<double>("minPt", 10);
0063 desc.add<double>("minAbsEta", 1.479);
0064 desc.add<double>("maxAbsEta", 4.0);
0065 desc.add<std::vector<double>>("endcapScalings", {3.17445, 1.13219, 0.0});
0066 desc.add<int>("quality", 5);
0067 SeedingRegionAlgoBase::fillPSetDescription(desc);
0068 }