L1TMuon/plugins/L1TMicroGMTInputProducerFromGen.cc

0001 // -*- C++ -*-
0002 //
0003 // Package:    L1TMicroGMTInputProducerFromGen
0004 // Class:      L1TMicroGMTInputProducerFromGen
0005 //
0006 /**\class L1TMicroGMTInputProducerFromGen L1TMicroGMTInputProducerFromGen.cc L1Trigger/L1TGlobalMuon/plugins/L1TMicroGMTInputProducerFromGen.cc
0007
0008  Description: takes generated muons and fills them in the expected collections for the MicroGMT
0009
0010  Implementation:
0011      [Notes on implementation]
0012 */
0013 //
0014 // Original Author:  Joschka Philip Lingemann,40 3-B01,+41227671598,
0015 //         Created:  Thu Oct  3 10:12:30 CEST 2013
0016 // $Id$
0017 //
0018 //
0019
0020 // system include files
0021 #include <memory>
0022 #include <fstream>
0023
0024 // user include files
0025 #include "FWCore/Framework/interface/Frameworkfwd.h"
0026 #include "FWCore/Framework/interface/stream/EDProducer.h"
0027
0028 #include "FWCore/Framework/interface/Event.h"
0029 #include "FWCore/Framework/interface/MakerMacros.h"
0030
0031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0032 #include "FWCore/Utilities/interface/Exception.h"
0033 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0034
0035 #include "DataFormats/L1TMuon/interface/RegionalMuonCandFwd.h"
0036 #include "DataFormats/L1TMuon/interface/RegionalMuonCand.h"
0037 #include "DataFormats/L1TMuon/interface/MuonCaloSumFwd.h"
0038 #include "DataFormats/L1TMuon/interface/MuonCaloSum.h"
0039
0040 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0041
0042 #include "TMath.h"
0043 #include "TRandom3.h"
0044
0045 //
0046 // class declaration
0047 //
0048 using namespace l1t;
0049
0050 class L1TMicroGMTInputProducerFromGen : public edm::stream::EDProducer<> {
0051 public:
0052   explicit L1TMicroGMTInputProducerFromGen(const edm::ParameterSet&);
0053
0054   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0055
0056 private:
0057   void produce(edm::Event&, const edm::EventSetup&) override;
0058
0059   static bool compareMuons(const RegionalMuonCand&, const RegionalMuonCand&);
0060
0061   // ----------member data ---------------------------
0062   edm::EDGetTokenT<reco::GenParticleCollection> genParticlesToken;
0063   int m_currEvt;
0064   const static int m_maxMuons = 108;
0065   TRandom3 m_rnd;
0066 };
0067
0068 //
0069 // constants, enums and typedefs
0070 //
0071
0072 //
0073 // static data member definitions
0074 //
0075
0076 //
0077 // constructors and destructor
0078 //
0079 L1TMicroGMTInputProducerFromGen::L1TMicroGMTInputProducerFromGen(const edm::ParameterSet& iConfig)
0080     : m_currEvt(0), m_rnd(0) {
0081   //register your inputs:
0082   genParticlesToken = consumes<reco::GenParticleCollection>(std::string("genParticles"));
0083   //register your products
0084   produces<RegionalMuonCandBxCollection>("BarrelTFMuons");
0085   produces<RegionalMuonCandBxCollection>("OverlapTFMuons");
0086   produces<RegionalMuonCandBxCollection>("ForwardTFMuons");
0087   produces<MuonCaloSumBxCollection>("TriggerTowerSums");
0088 }
0089
0090 //
0091 // member functions
0092 //
0093
0094 bool L1TMicroGMTInputProducerFromGen::compareMuons(const RegionalMuonCand& mu1, const RegionalMuonCand& mu2) {
0095   return mu1.processor() < mu2.processor();
0096 }
0097
0098 // ------------ method called to produce the data  ------------
0099 void L1TMicroGMTInputProducerFromGen::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0100   using namespace edm;
0101
0102   std::unique_ptr<RegionalMuonCandBxCollection> barrelMuons(new RegionalMuonCandBxCollection());
0103   std::unique_ptr<RegionalMuonCandBxCollection> overlapMuons(new RegionalMuonCandBxCollection());
0104   std::unique_ptr<RegionalMuonCandBxCollection> endcapMuons(new RegionalMuonCandBxCollection());
0105   std::unique_ptr<MuonCaloSumBxCollection> towerSums(new MuonCaloSumBxCollection());
0106
0107   std::vector<RegionalMuonCand> bmMuons;
0108   std::vector<RegionalMuonCand> omMuons;
0109   std::vector<RegionalMuonCand> emMuons;
0110
0111   std::vector<int> muIndices;
0112   edm::Handle<reco::GenParticleCollection> genParticles;
0113   // Make sure that you can get genParticles
0114   if (iEvent.getByToken(genParticlesToken, genParticles)) {
0115     int cntr = 0;
0116     for (auto it = genParticles->cbegin(); it != genParticles->cend(); ++it) {
0117       const reco::Candidate& mcParticle = *it;
0118       if (abs(mcParticle.pdgId()) == 13 && mcParticle.status() == 1)
0119         muIndices.push_back(cntr);
0120       cntr++;
0121     }
0122   } else {
0123     LogTrace("GlobalMuon") << " GenParticleCollection not found." << std::endl;
0124   }
0125
0126   RegionalMuonCand mu;
0127   MuonCaloSum tSum;
0128   // alternative scale (using full phi bit-width): 163.4521265553765f;
0129   const float phiToInt = 91.67324722093171f;
0130   // alternative scale: 100.0f;
0131   const float etaToInt = 90.9090909090f;
0132   const int maxPt = (1 << 9) - 1;
0133   int muCntr = 0;
0134
0135   double twoPi = TMath::Pi() * 2.;
0136
0137   for (auto it = muIndices.begin(); it != muIndices.end(); ++it) {
0138     // don't really care which muons are taken...
0139     // guess there ain't 108 generated anyways
0140     if (muCntr == m_maxMuons)
0141       break;
0142     int gen_idx = *it;
0143     const reco::Candidate& mcMuon = genParticles->at(gen_idx);
0144     double eta = mcMuon.eta();
0145     if (fabs(eta) > 2.45)
0146       continue;  // out of acceptance
0147     int hwPt = int(mcMuon.pt() * 2);
0148     hwPt = (hwPt < maxPt ? hwPt : maxPt);
0149     int hwEta = int(eta * etaToInt);
0150     double phi = mcMuon.phi();
0151     if (phi < 0)
0152       phi += twoPi;  // add 2*pi
0153     int hwPhi = (int(phi * phiToInt)) % 576;
0154     int hwQual = 8;
0155     int hwCharge = (mcMuon.charge() > 0) ? 0 : 1;
0156     int hwChargeValid = 1;
0157
0158     mu.setHwPt(hwPt);
0159
0160     tftype tf(tftype::bmtf);
0161     int globalWedgePhi = (hwPhi + 24) % 576;  // this sets CMS phi = 0 to -15 deg
0162     int localPhi = globalWedgePhi % 48;
0163     int processor = globalWedgePhi / 48 + 1;
0164     int globalSectorPhi = (hwPhi - 24);  // this sets CMS phi = 0 to +15 deg
0165     if (globalSectorPhi < 0) {
0166       globalSectorPhi += 576;
0167     }
0168
0169     if (fabs(eta) > 0.8) {
0170       if (fabs(eta) < 1.2) {
0171         tf = (eta > 0 ? tftype::omtf_pos : tftype::omtf_neg);
0172         processor = globalSectorPhi / 96 + 1;
0173         localPhi = globalSectorPhi % 96;
0174       } else {
0175         tf = (eta > 0 ? tftype::emtf_pos : tftype::emtf_neg);
0176         processor = globalSectorPhi / 96 + 1;
0177         localPhi = globalSectorPhi % 96;
0178       }
0179     }
0180     mu.setHwPhi(localPhi);
0181     mu.setTFIdentifiers(processor, tf);
0182
0183     mu.setHwEta(hwEta);
0184     mu.setHwSign(hwCharge);
0185     mu.setHwSignValid(hwChargeValid);
0186     mu.setHwQual(hwQual);
0187
0188     if (fabs(eta) < 0.8 && bmMuons.size() < 36) {
0189       bmMuons.push_back(mu);
0190       muCntr++;
0191     } else if (fabs(eta) < 1.2 && omMuons.size() < 36) {
0192       omMuons.push_back(mu);
0193       muCntr++;
0194     } else if (emMuons.size() < 36) {
0195       emMuons.push_back(mu);
0196       muCntr++;
0197     }
0198   }
0199
0200   std::sort(bmMuons.begin(), bmMuons.end(), L1TMicroGMTInputProducerFromGen::compareMuons);
0201   std::sort(omMuons.begin(), omMuons.end(), L1TMicroGMTInputProducerFromGen::compareMuons);
0202   std::sort(emMuons.begin(), emMuons.end(), L1TMicroGMTInputProducerFromGen::compareMuons);
0203
0204   for (const auto& mu : bmMuons) {
0205     barrelMuons->push_back(0, mu);
0206   }
0207
0208   for (const auto& mu : omMuons) {
0209     overlapMuons->push_back(0, mu);
0210   }
0211
0212   for (const auto& mu : emMuons) {
0213     endcapMuons->push_back(0, mu);
0214   }
0215
0216   for (int i = 0; i < 1008; ++i) {
0217     // from where could I take the tower energies?
0218     int energy = int(m_rnd.Gaus(12, 6));
0219     if (energy < 0)
0220       energy = 0;
0221     if (energy > 31)
0222       energy = 31;
0223     MuonCaloSum sum(energy, i / 28, i % 28, i);
0224     towerSums->push_back(0, sum);
0225   }
0226
0227   iEvent.put(std::move(barrelMuons), "BarrelTFMuons");
0228   iEvent.put(std::move(overlapMuons), "OverlapTFMuons");
0229   iEvent.put(std::move(endcapMuons), "ForwardTFMuons");
0230   iEvent.put(std::move(towerSums), "TriggerTowerSums");
0231   m_currEvt++;
0232 }
0233
0234 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0235 void L1TMicroGMTInputProducerFromGen::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0236   //The following says we do not know what parameters are allowed so do no validation
0237   // Please change this to state exactly what you do use, even if it is no parameters
0238   edm::ParameterSetDescription desc;
0239   desc.setUnknown();
0240   descriptions.addDefault(desc);
0241 }
0242
0243 //define this as a plug-in
0244 DEFINE_FWK_MODULE(L1TMicroGMTInputProducerFromGen);