L1TCalorimeter/plugins/L1TPhysicalEtAdder.cc

0001 #include "L1Trigger/L1TCalorimeter/plugins/L1TPhysicalEtAdder.h"
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003
0004 #include "DataFormats/L1CaloTrigger/interface/L1CaloEmCand.h"
0005 #include "DataFormats/L1CaloTrigger/interface/L1CaloRegion.h"
0006 #include "DataFormats/L1CaloTrigger/interface/L1CaloCollections.h"
0007
0008 #include "DataFormats/L1Trigger/interface/BXVector.h"
0009 #include "DataFormats/L1Trigger/interface/EtSum.h"
0010
0011 #include "DataFormats/L1TCalorimeter/interface/CaloEmCand.h"
0012 #include "DataFormats/L1TCalorimeter/interface/CaloRegion.h"
0013
0014 #include "FWCore/Framework/interface/ESHandle.h"
0015 #include <vector>
0016
0017 namespace {
0018
0019   int getRegionEta(int gtEta, bool forward) {
0020     // backwards conversion is
0021     // unsigned rctEta = (iEta<11 ? 10-iEta : iEta-11);
0022     // return (((rctEta % 7) & 0x7) | (iEta<11 ? 0x8 : 0));
0023     int centralGtEta[] = {11, 12, 13, 14, 15, 16, 17, -100, 10, 9, 8, 7, 6, 5, 4};
0024     int forwardGtEta[] = {18, 19, 20, 21, -100, -100, -100, -100, 3, 2, 1, 0};
0025
0026     //printf("%i, %i\n",gtEta,forward);
0027
0028     int regionEta;
0029
0030     if (!forward) {
0031       regionEta = centralGtEta[gtEta];
0032     } else
0033       regionEta = forwardGtEta[gtEta];
0034
0035     if (regionEta == -100)
0036       edm::LogError("EtaIndexError") << "Bad eta index passed to L1TPhysicalEtAdder::getRegionEta, " << gtEta
0037                                      << std::endl;
0038
0039     return regionEta;
0040   }
0041
0042   // adapted these from the UCT2015 codebase.
0043   double getPhysicalEta(int gtEta, bool forward = false) {
0044     int etaIndex = getRegionEta(gtEta, forward);
0045
0046     const double rgnEtaValues[11] = {0.174,  // HB and inner HE bins are 0.348 wide
0047                                      0.522,
0048                                      0.870,
0049                                      1.218,
0050                                      1.566,
0051                                      1.956,  // Last two HE bins are 0.432 and 0.828 wide
0052                                      2.586,
0053                                      3.250,  // HF bins are 0.5 wide
0054                                      3.750,
0055                                      4.250,
0056                                      4.750};
0057     if (etaIndex < 11) {
0058       return -rgnEtaValues[-(etaIndex - 10)];  // 0-10 are negative eta values
0059     } else if (etaIndex < 22) {
0060       return rgnEtaValues[etaIndex - 11];  // 11-21 are positive eta values
0061     }
0062     return -9;
0063   }
0064
0065   double getPhysicalPhi(int phiIndex) {
0066     if (phiIndex < 10)
0067       return 2. * M_PI * phiIndex / 18.;
0068     if (phiIndex < 18)
0069       return -M_PI + 2. * M_PI * (phiIndex - 9) / 18.;
0070     return -9;
0071   }
0072
0073 }  // namespace
0074
0075 using namespace l1t;
0076
0077 L1TPhysicalEtAdder::L1TPhysicalEtAdder(const edm::ParameterSet& ps) {
0078   produces<EGammaBxCollection>();
0079   produces<TauBxCollection>("rlxTaus");
0080   produces<TauBxCollection>("isoTaus");
0081   produces<JetBxCollection>();
0082   produces<JetBxCollection>("preGtJets");
0083   produces<EtSumBxCollection>();
0084   produces<CaloSpareBxCollection>("HFRingSums");
0085   produces<CaloSpareBxCollection>("HFBitCounts");
0086
0087   EGammaToken_ = consumes<EGammaBxCollection>(ps.getParameter<edm::InputTag>("InputCollection"));
0088   RlxTauToken_ = consumes<TauBxCollection>(ps.getParameter<edm::InputTag>("InputRlxTauCollection"));
0089   IsoTauToken_ = consumes<TauBxCollection>(ps.getParameter<edm::InputTag>("InputIsoTauCollection"));
0090   JetToken_ = consumes<JetBxCollection>(ps.getParameter<edm::InputTag>("InputCollection"));
0091   preGtJetToken_ = consumes<JetBxCollection>(ps.getParameter<edm::InputTag>("InputPreGtJetCollection"));
0092   EtSumToken_ = consumes<EtSumBxCollection>(ps.getParameter<edm::InputTag>("InputCollection"));
0093   HfSumsToken_ = consumes<CaloSpareBxCollection>(ps.getParameter<edm::InputTag>("InputHFSumsCollection"));
0094   HfCountsToken_ = consumes<CaloSpareBxCollection>(ps.getParameter<edm::InputTag>("InputHFCountsCollection"));
0095   emScaleToken_ = esConsumes<L1CaloEtScale, L1EmEtScaleRcd>();
0096   jetScaleToken_ = esConsumes<L1CaloEtScale, L1JetEtScaleRcd>();
0097   htMissScaleToken_ = esConsumes<L1CaloEtScale, L1HtMissScaleRcd>();
0098 }
0099
0100 L1TPhysicalEtAdder::~L1TPhysicalEtAdder() {}
0101
0102 // ------------ method called to produce the data  ------------
0103 void L1TPhysicalEtAdder::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0104   // store new collections which include physical quantities
0105   std::unique_ptr<EGammaBxCollection> new_egammas(new EGammaBxCollection);
0106   std::unique_ptr<TauBxCollection> new_rlxtaus(new TauBxCollection);
0107   std::unique_ptr<TauBxCollection> new_isotaus(new TauBxCollection);
0108   std::unique_ptr<JetBxCollection> new_jets(new JetBxCollection);
0109   std::unique_ptr<JetBxCollection> new_preGtJets(new JetBxCollection);
0110   std::unique_ptr<EtSumBxCollection> new_etsums(new EtSumBxCollection);
0111   std::unique_ptr<CaloSpareBxCollection> new_hfsums(new CaloSpareBxCollection);
0112   std::unique_ptr<CaloSpareBxCollection> new_hfcounts(new CaloSpareBxCollection);
0113
0114   edm::Handle<EGammaBxCollection> old_egammas;
0115   edm::Handle<TauBxCollection> old_rlxtaus;
0116   edm::Handle<TauBxCollection> old_isotaus;
0117   edm::Handle<JetBxCollection> old_jets;
0118   edm::Handle<JetBxCollection> old_preGtJets;
0119   edm::Handle<EtSumBxCollection> old_etsums;
0120   edm::Handle<CaloSpareBxCollection> old_hfsums;
0121   edm::Handle<CaloSpareBxCollection> old_hfcounts;
0122
0123   iEvent.getByToken(EGammaToken_, old_egammas);
0124   iEvent.getByToken(RlxTauToken_, old_rlxtaus);
0125   iEvent.getByToken(IsoTauToken_, old_isotaus);
0126   iEvent.getByToken(JetToken_, old_jets);
0127   iEvent.getByToken(preGtJetToken_, old_preGtJets);
0128   iEvent.getByToken(EtSumToken_, old_etsums);
0129   iEvent.getByToken(HfSumsToken_, old_hfsums);
0130   iEvent.getByToken(HfCountsToken_, old_hfcounts);
0131
0132   //get the proper scales for conversion to physical et
0133   edm::ESHandle<L1CaloEtScale> emScale = iSetup.getHandle(emScaleToken_);
0134
0135   edm::ESHandle<L1CaloEtScale> jetScale = iSetup.getHandle(jetScaleToken_);
0136
0137   edm::ESHandle<L1CaloEtScale> htMissScale = iSetup.getHandle(htMissScaleToken_);
0138
0139   int firstBX = old_egammas->getFirstBX();
0140   int lastBX = old_egammas->getLastBX();
0141
0142   new_egammas->setBXRange(firstBX, lastBX);
0143   new_rlxtaus->setBXRange(firstBX, lastBX);
0144   new_isotaus->setBXRange(firstBX, lastBX);
0145   new_jets->setBXRange(firstBX, lastBX);
0146   new_preGtJets->setBXRange(firstBX, lastBX);
0147   new_etsums->setBXRange(firstBX, lastBX);
0148   new_hfsums->setBXRange(firstBX, lastBX);
0149   new_hfcounts->setBXRange(firstBX, lastBX);
0150
0151   for (int bx = firstBX; bx <= lastBX; ++bx) {
0152     for (EGammaBxCollection::const_iterator itEGamma = old_egammas->begin(bx); itEGamma != old_egammas->end(bx);
0153          ++itEGamma) {
0154       //const double pt = itEGamma->hwPt() * emScale->linearLsb();
0155       const double et = emScale->et(itEGamma->hwPt());
0156       const double eta = getPhysicalEta(itEGamma->hwEta());
0157       const double phi = getPhysicalPhi(itEGamma->hwPhi());
0158       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0159
0160       EGamma eg(*&p4, itEGamma->hwPt(), itEGamma->hwEta(), itEGamma->hwPhi(), itEGamma->hwQual(), itEGamma->hwIso());
0161       new_egammas->push_back(bx, *&eg);
0162     }
0163
0164     for (TauBxCollection::const_iterator itTau = old_rlxtaus->begin(bx); itTau != old_rlxtaus->end(bx); ++itTau) {
0165       // use the full-circle conversion to match l1extra, accounts for linearLsb and max value automatically
0166       //const uint16_t rankPt = jetScale->rank((uint16_t)itTau->hwPt());
0167       //const double et = jetScale->et( rankPt ) ;
0168
0169       // or use the emScale to get finer-grained et
0170       //const double et = itTau->hwPt() * emScale->linearLsb();
0171
0172       // we are now already in the rankPt
0173       const double et = jetScale->et(itTau->hwPt());
0174
0175       const double eta = getPhysicalEta(itTau->hwEta());
0176       const double phi = getPhysicalPhi(itTau->hwPhi());
0177       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0178
0179       Tau tau(*&p4, itTau->hwPt(), itTau->hwEta(), itTau->hwPhi(), itTau->hwQual(), itTau->hwIso());
0180       new_rlxtaus->push_back(bx, *&tau);
0181     }
0182
0183     for (TauBxCollection::const_iterator itTau = old_isotaus->begin(bx); itTau != old_isotaus->end(bx); ++itTau) {
0184       // use the full-circle conversion to match l1extra, accounts for linearLsb and max value automatically
0185       //const uint16_t rankPt = jetScale->rank((uint16_t)itTau->hwPt());
0186       //const double et = jetScale->et( rankPt ) ;
0187
0188       // or use the emScale to get finer-grained et
0189       //const double et = itTau->hwPt() * emScale->linearLsb();
0190
0191       // we are now already in the rankPt
0192       const double et = jetScale->et(itTau->hwPt());
0193
0194       const double eta = getPhysicalEta(itTau->hwEta());
0195       const double phi = getPhysicalPhi(itTau->hwPhi());
0196       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0197
0198       Tau tau(*&p4, itTau->hwPt(), itTau->hwEta(), itTau->hwPhi(), itTau->hwQual(), itTau->hwIso());
0199       new_isotaus->push_back(bx, *&tau);
0200     }
0201
0202     for (JetBxCollection::const_iterator itJet = old_jets->begin(bx); itJet != old_jets->end(bx); ++itJet) {
0203       // use the full-circle conversion to match l1extra, accounts for linearLsb and max value automatically
0204       //const uint16_t rankPt = jetScale->rank((uint16_t)itJet->hwPt());
0205       //const double et = jetScale->et( rankPt ) ;
0206
0207       // or use the emScale to get finer-grained et
0208       //const double et = itJet->hwPt() * emScale->linearLsb();
0209
0210       // we are now already in the rankPt
0211       const double et = jetScale->et(itJet->hwPt());
0212
0213       const bool forward = ((itJet->hwQual() & 0x2) != 0);
0214       const double eta = getPhysicalEta(itJet->hwEta(), forward);
0215       const double phi = getPhysicalPhi(itJet->hwPhi());
0216       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0217
0218       Jet jet(*&p4, itJet->hwPt(), itJet->hwEta(), itJet->hwPhi(), itJet->hwQual());
0219       new_jets->push_back(bx, *&jet);
0220     }
0221
0222     for (JetBxCollection::const_iterator itJet = old_preGtJets->begin(bx); itJet != old_preGtJets->end(bx); ++itJet) {
0223       // use the full-circle conversion to match l1extra, accounts for linearLsb and max value automatically
0224       //const uint16_t rankPt = jetScale->rank((uint16_t)itJet->hwPt());
0225       //const double et = jetScale->et( rankPt ) ;
0226
0227       // or use the emScale to get finer-grained et
0228       const double et = itJet->hwPt() * emScale->linearLsb();
0229
0230       // we are now already in the rankPt
0231       //const double et = jetScale->et( itJet->hwPt() );
0232
0233       const bool forward = ((itJet->hwQual() & 0x2) != 0);
0234       const double eta = getPhysicalEta(itJet->hwEta(), forward);
0235       const double phi = getPhysicalPhi(itJet->hwPhi());
0236       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0237
0238       Jet jet(*&p4, itJet->hwPt(), itJet->hwEta(), itJet->hwPhi(), itJet->hwQual());
0239       new_preGtJets->push_back(bx, *&jet);
0240     }
0241
0242     for (EtSumBxCollection::const_iterator itEtSum = old_etsums->begin(bx); itEtSum != old_etsums->end(bx); ++itEtSum) {
0243       double et = itEtSum->hwPt() * emScale->linearLsb();
0244       //hack while we figure out the right scales
0245       //double et = emScale->et( itEtSum->hwPt() );
0246       const EtSum::EtSumType sumType = itEtSum->getType();
0247
0248       const double eta = getPhysicalEta(itEtSum->hwEta());
0249       double phi = getPhysicalPhi(itEtSum->hwPhi());
0250       if (sumType == EtSum::EtSumType::kMissingHt) {
0251         et = htMissScale->et(itEtSum->hwPt());
0252         double regionPhiWidth = 2. * 3.1415927 / L1CaloRegionDetId::N_PHI;
0253         phi = phi + (regionPhiWidth / 2.);  // add the region half-width to match L1Extra MHT phi
0254       }
0255
0256       math::PtEtaPhiMLorentzVector p4(et, eta, phi, 0);
0257
0258       EtSum eg(*&p4, sumType, itEtSum->hwPt(), itEtSum->hwEta(), itEtSum->hwPhi(), itEtSum->hwQual());
0259       new_etsums->push_back(bx, *&eg);
0260     }
0261
0262     for (CaloSpareBxCollection::const_iterator itCaloSpare = old_hfsums->begin(bx); itCaloSpare != old_hfsums->end(bx);
0263          ++itCaloSpare) {
0264       //just pass through for now
0265       //a different scale is needed depending on the type
0266       new_hfsums->push_back(bx, *itCaloSpare);
0267     }
0268
0269     for (CaloSpareBxCollection::const_iterator itCaloSpare = old_hfcounts->begin(bx);
0270          itCaloSpare != old_hfcounts->end(bx);
0271          ++itCaloSpare) {
0272       //just pass through for now
0273       //a different scale is needed depending on the type
0274       new_hfcounts->push_back(bx, *itCaloSpare);
0275     }
0276   }
0277
0278   iEvent.put(std::move(new_egammas));
0279   iEvent.put(std::move(new_rlxtaus), "rlxTaus");
0280   iEvent.put(std::move(new_isotaus), "isoTaus");
0281   iEvent.put(std::move(new_jets));
0282   iEvent.put(std::move(new_preGtJets), "preGtJets");
0283   iEvent.put(std::move(new_etsums));
0284   iEvent.put(std::move(new_hfsums), "HFRingSums");
0285   iEvent.put(std::move(new_hfcounts), "HFBitCounts");
0286 }
0287
0288 // ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
0289 void L1TPhysicalEtAdder::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0290   //The following says we do not know what parameters are allowed so do no validation
0291   // Please change this to state exactly what you do use, even if it is no parameters
0292   edm::ParameterSetDescription desc;
0293   desc.setUnknown();
0294   descriptions.addDefault(desc);
0295 }
0296
0297 //define this as a plug-in
0298 DEFINE_FWK_MODULE(L1TPhysicalEtAdder);