HcalTrigPrimProducers/plugins/HcalTTPDigiProducer.cc

0001 #include "CalibFormats/HcalObjects/interface/HcalTPGRecord.h"
0002 #include "CalibFormats/HcalObjects/interface/HcalTPGCoder.h"
0003 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0004 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0005 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
0006 #include "DataFormats/HcalDigi/interface/HFDataFrame.h"
0007 #include "DataFormats/Common/interface/Handle.h"
0008 #include "FWCore/Framework/interface/stream/EDProducer.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "FWCore/Framework/interface/EventSetup.h"
0012 #include "FWCore/Framework/interface/ESHandle.h"
0013 #include "FWCore/Utilities/interface/ESGetToken.h"
0014 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0015
0016 #include <cstdio>
0017
0018 class HcalTTPDigiProducer : public edm::stream::EDProducer<> {
0019 public:
0020   explicit HcalTTPDigiProducer(const edm::ParameterSet& ps);
0021   ~HcalTTPDigiProducer() override = default;
0022
0023   void produce(edm::Event& e, const edm::EventSetup& c) override;
0024
0025 private:
0026   bool isMasked(HcalDetId id);
0027   bool decision(int nP, int nM, int bit);
0028
0029   edm::EDGetTokenT<HFDigiCollection> tok_hf_;
0030   edm::ESGetToken<HcalTPGCoder, HcalTPGRecord> tok_tpgCoder_;
0031   std::vector<unsigned int> maskedChannels_;
0032   std::string bit_[4];
0033   int calc_[4];
0034   int nHits_[4], nHFp_[4], nHFm_[4];
0035   char pReq_[4], mReq_[4], pmLogic_[4];
0036   int id_, samples_, presamples_;
0037   int fwAlgo_;
0038   int iEtaMin_, iEtaMax_;
0039   unsigned int threshold_;
0040
0041   int SoI_;
0042
0043   static const int inputs_[];
0044 };
0045
0046 // DO NOT MODIFY: Mapping between iphi (array index) and TTP input (value) for HF
0047 const int HcalTTPDigiProducer::inputs_[] = {30, 66, 4,  44, 4,  44, 0,  68, 0,  68, 16, 48, 16, 48, 6,  46, 6,  46,
0048                                             2,  70, 2,  70, 18, 50, 18, 50, 12, 40, 12, 40, 8,  52, 8,  52, 20, 36,
0049                                             20, 36, 14, 42, 14, 42, 10, 54, 10, 54, 22, 38, 22, 38, 24, 56, 24, 56,
0050                                             32, 60, 32, 60, 28, 64, 28, 64, 26, 58, 26, 58, 34, 62, 34, 62, 30, 66};
0051
0052 HcalTTPDigiProducer::HcalTTPDigiProducer(const edm::ParameterSet& ps) {
0053   tok_hf_ = consumes<HFDigiCollection>(ps.getParameter<edm::InputTag>("HFDigiCollection"));
0054   tok_tpgCoder_ = esConsumes<HcalTPGCoder, HcalTPGRecord>();
0055   maskedChannels_ = ps.getParameter<std::vector<unsigned int> >("maskedChannels");
0056   bit_[0] = ps.getParameter<std::string>("defTT8");
0057   bit_[1] = ps.getParameter<std::string>("defTT9");
0058   bit_[2] = ps.getParameter<std::string>("defTT10");
0059   bit_[3] = ps.getParameter<std::string>("defTTLocal");
0060
0061   for (int i = 0; i < 4; i++) {
0062     nHits_[i] = -1;
0063     nHFp_[i] = -1;
0064     nHFm_[i] = -1;
0065     pReq_[i] = ' ';
0066     mReq_[i] = ' ';
0067     pmLogic_[i] = ' ';
0068     calc_[i] = sscanf(bit_[i].c_str(),
0069                       "hits>=%d:hfp%c=%d%chfm%c=%d",
0070                       &(nHits_[i]),
0071                       &(pReq_[i]),
0072                       &(nHFp_[i]),
0073                       &(pmLogic_[i]),
0074                       &(mReq_[i]),
0075                       &(nHFm_[i]));
0076     if (calc_[i] == 1) {
0077       if (nHits_[i] < 0)
0078         throw cms::Exception("HcalTTPDigiProducer") << "Unable to read logic for technical trigger";
0079     } else if (calc_[i] == 6) {
0080       if (nHits_[i] < 0 || nHFp_[i] < 0 || nHFm_[i] < 0)
0081         throw cms::Exception("HcalTTPDigiProducer") << "Unable to read logic for technical trigger";
0082       if ((pReq_[i] != '>' && pReq_[i] != '<') || (mReq_[i] != '>' && mReq_[i] != '<') ||
0083           (pmLogic_[i] != ':' && pmLogic_[i] != '|'))
0084         throw cms::Exception("HcalTTPDigiProducer") << "Technical Trigger logic must obey the following format:\n"
0085                                                        "\"hits>=[A1]:hfp[B1]=[A2][C]hfm[B2]=[A3]\",\n"
0086                                                        "or \"hits>=[A1]\",\n"
0087                                                        "with A# >= 0, B# = (</>) and C = (:/|)";
0088     } else {
0089       throw cms::Exception("HcalTTPDigiProducer") << "Unable to read logic for technical trigger";
0090     }
0091   }
0092
0093   id_ = ps.getUntrackedParameter<int>("id", -1);
0094   samples_ = ps.getParameter<int>("samples");
0095   presamples_ = ps.getParameter<int>("presamples");
0096   iEtaMin_ = ps.getParameter<int>("iEtaMin");
0097   iEtaMax_ = ps.getParameter<int>("iEtaMax");
0098   threshold_ = ps.getParameter<unsigned int>("threshold");
0099   fwAlgo_ = ps.getParameter<int>("fwAlgorithm");
0100
0101   SoI_ = ps.getParameter<int>("HFSoI");
0102
0103   if (samples_ > 8) {
0104     samples_ = 8;
0105     edm::LogWarning("HcalTTPDigiProducer") << "Samples forced to maximum value of 8";
0106   }
0107   if (presamples_ - SoI_ > 0) {  // Too many presamples
0108     presamples_ = SoI_;
0109     edm::LogWarning("HcalTTPDigiProducer") << "Presamples reset to HF SoI value";
0110   }
0111
0112   produces<HcalTTPDigiCollection>();
0113 }
0114
0115 bool HcalTTPDigiProducer::isMasked(HcalDetId id) {
0116   for (unsigned int i = 0; i < maskedChannels_.size(); i++)
0117     if (id.rawId() == maskedChannels_.at(i))
0118       return true;
0119   return false;
0120 }
0121
0122 bool HcalTTPDigiProducer::decision(int nP, int nM, int bit) {
0123   bool pOK = false;
0124   bool mOK = false;
0125   if ((nP + nM) < nHits_[bit])
0126     return false;
0127   if (calc_[bit] == 1)
0128     return ((nP + nM) >= nHits_[bit]);
0129
0130   if (pReq_[bit] == '>')
0131     pOK = (nP >= nHFp_[bit]);
0132   else if (pReq_[bit] == '<')
0133     pOK = (nP <= nHFp_[bit]);
0134
0135   if (mReq_[bit] == '>')
0136     mOK = (nM >= nHFm_[bit]);
0137   else if (mReq_[bit] == '<')
0138     mOK = (nM <= nHFm_[bit]);
0139
0140   if (pmLogic_[bit] == ':')
0141     return (pOK && mOK);
0142   else if (pmLogic_[bit] == '|')
0143     return (pOK || mOK);
0144
0145   // Should not ever get here...need to create a warning message
0146   edm::LogWarning("HcalTTPDigiProducer") << "Trigger logic exhausted.  Returning false";
0147   return false;
0148 }
0149
0150 void HcalTTPDigiProducer::produce(edm::Event& e, const edm::EventSetup& eventSetup) {
0151   // Step A: Get Inputs
0152   edm::Handle<HFDigiCollection> hfDigiCollection;
0153   e.getByToken(tok_hf_, hfDigiCollection);
0154   edm::ESHandle<HcalTPGCoder> inputCoder = eventSetup.getHandle(tok_tpgCoder_);
0155
0156   // Step B: Create empty output
0157   std::unique_ptr<HcalTTPDigiCollection> ttpResult(new HcalTTPDigiCollection());
0158
0159   // Step C: Compute TTP inputs
0160   uint16_t trigInputs[40];
0161   int nP[8];
0162   int nM[8];
0163   for (int i = 0; i < 8; i++) {
0164     nP[i] = 0;
0165     nM[i] = 0;
0166     for (int j = 0; j < 5; j++)
0167       trigInputs[j * 8 + i] = 0;
0168   }
0169   for (HFDigiCollection::const_iterator theDigi = hfDigiCollection->begin(); theDigi != hfDigiCollection->end();
0170        theDigi++) {
0171     HcalDetId id = HcalDetId(theDigi->id());
0172     if (isMasked(id))
0173       continue;
0174     if (id.ietaAbs() < iEtaMin_ || id.ietaAbs() > iEtaMax_)
0175       continue;
0176
0177     IntegerCaloSamples samples(id, theDigi->size());
0178     inputCoder->adc2Linear(*theDigi, samples);
0179
0180     for (int relSample = -presamples_; relSample < (samples_ - presamples_); relSample++) {
0181       if (samples[SoI_ + relSample] >= threshold_) {
0182         int linSample = presamples_ + relSample;
0183         int offset = (-1 + id.zside()) / 2;
0184         int shift = inputs_[id.iphi() + offset];
0185         int group = 0;
0186         while (shift >= 16) {
0187           shift -= 16;
0188           group++;
0189         }
0190         if (!(trigInputs[(linSample * 8) + group] & (1 << shift)))
0191           (id.ieta() > 0) ? (nP[linSample]++) : (nM[linSample]++);
0192         trigInputs[(linSample * 8) + group] |= (1 << shift);
0193       }
0194     }
0195   }
0196
0197   // Step D: Compute trigger decision and fill TTP digi
0198   uint8_t trigOutput[8];
0199   uint32_t algoDepBits[8];
0200   HcalTTPDigi ttpDigi(id_, samples_, presamples_, 0, fwAlgo_, 0);
0201   for (int linSample = 0; linSample < 8; linSample++) {
0202     trigOutput[linSample] = 0;
0203     algoDepBits[linSample] = 0;
0204     if (linSample < samples_) {
0205       for (int j = 0; j < 4; j++)
0206         trigOutput[linSample] |= (decision(nP[linSample], nM[linSample], j) << j);
0207       int nT = nP[linSample] + nM[linSample];
0208
0209       // Algorithm Dependent bits for FW flavor = 1
0210       // NOTE: this disagrees with the fw var. names that implies (LSB) T,M,P (MSB)
0211       if (fwAlgo_ == 1)
0212         algoDepBits[linSample] = (nT & 0x7F) | ((nP[linSample] & 0x3F) << 7) | ((nM[linSample] & 0x3F) << 13);
0213       ttpDigi.setSample(
0214           (linSample - presamples_), &trigInputs[linSample * 8], algoDepBits[linSample], trigOutput[linSample]);
0215     }
0216   }
0217   ttpResult->push_back(ttpDigi);
0218
0219   // Step E: Put outputs into event
0220   e.put(std::move(ttpResult));
0221 }
0222
0223 #include "FWCore/PluginManager/interface/ModuleDef.h"
0224 #include "FWCore/Framework/interface/MakerMacros.h"
0225
0226 DEFINE_FWK_MODULE(HcalTTPDigiProducer);