File indexing completed on 2025-03-05 03:13:55
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include <memory>
0015
0016
0017 #include "DataFormats/Common/interface/Handle.h"
0018 #include "FWCore/Framework/interface/ESHandle.h"
0019 #include "FWCore/Framework/interface/Frameworkfwd.h"
0020 #include "FWCore/Framework/interface/MakerMacros.h"
0021
0022 #include "L1GctInternJetProducer.h"
0023
0024 L1GctInternJetProducer::L1GctInternJetProducer(const edm::ParameterSet& iConfig)
0025 : internalJetSource_(iConfig.getParameter<edm::InputTag>("internalJetSource")),
0026 caloGeomToken_(esConsumes<L1CaloGeometry, L1CaloGeometryRecord>()),
0027 jetScaleToken_(esConsumes<L1CaloEtScale, L1JetEtScaleRcd>()),
0028 centralBxOnly_(iConfig.getParameter<bool>("centralBxOnly")) {
0029 using namespace l1extra;
0030
0031 produces<L1JetParticleCollection>("Internal");
0032 }
0033
0034 void L1GctInternJetProducer::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0035
0036 using namespace edm;
0037 using namespace l1extra;
0038 using namespace std;
0039
0040 unique_ptr<L1JetParticleCollection> internJetColl(new L1JetParticleCollection);
0041
0042 ESHandle<L1CaloGeometry> caloGeomESH = iSetup.getHandle(caloGeomToken_);
0043 const L1CaloGeometry* caloGeom = &(*caloGeomESH);
0044
0045 ESHandle<L1CaloEtScale> jetScale = iSetup.getHandle(jetScaleToken_);
0046
0047 double etSumLSB = jetScale->linearLsb();
0048
0049
0050 Handle<L1GctInternJetDataCollection> hwIntJetCands;
0051 iEvent.getByLabel(internalJetSource_, hwIntJetCands);
0052
0053 if (!hwIntJetCands.isValid()) {
0054 std::cout << "These aren't the Jets you're looking for" << std::endl;
0055
0056 LogDebug("L1GctInternJetProducer") << "\nWarning: L1GctJetCandCollection with " << internalJetSource_
0057 << "\nrequested in configuration, but not found in the event." << std::endl;
0058
0059 } else {
0060
0061 L1GctInternJetDataCollection::const_iterator jetItr = hwIntJetCands->begin();
0062 L1GctInternJetDataCollection::const_iterator jetEnd = hwIntJetCands->end();
0063 int i;
0064 for (i = 0; jetItr != jetEnd; ++jetItr, ++i) {
0065
0066 if (!jetItr->empty() && (!centralBxOnly_ || jetItr->bx() == 0)) {
0067 double et = (jetItr->oflow() ? (double)0xfff : (double)jetItr->et()) * etSumLSB + 1.e-6;
0068
0069
0070
0071
0072 double eta = caloGeom->etaBinCenter(jetItr->regionId());
0073 double phi = caloGeom->emJetPhiBinCenter(jetItr->regionId());
0074
0075 internJetColl->push_back(
0076 L1JetParticle(math::PtEtaPhiMLorentzVector(et, eta, phi, 0.), Ref<L1GctJetCandCollection>(), jetItr->bx()));
0077 }
0078 }
0079 }
0080
0081 iEvent.put(std::move(internJetColl), "Internal");
0082 }
0083
0084
0085 DEFINE_FWK_MODULE(L1GctInternJetProducer);