Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 12:49:44

0001 /** \class L1TTkEleFilter
0002  *
0003  * See header file for documentation
0004  *
0005  *
0006  *  \author Martin Grunewald
0007  *  \author Simone Gennai
0008  *  \author Thiago Tomei
0009  *
0010  */
0011 
0012 #include "L1TTkEleFilter.h"
0013 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0014 
0015 #include "DataFormats/Common/interface/Handle.h"
0016 
0017 #include "DataFormats/Common/interface/Ref.h"
0018 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0019 #include "DataFormats/HLTReco/interface/TriggerTypeDefs.h"
0020 
0021 #include "FWCore/Framework/interface/EventSetup.h"
0022 #include "FWCore/Framework/interface/EventSetupRecord.h"
0023 #include "FWCore/Utilities/interface/Exception.h"
0024 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0025 
0026 //
0027 // constructors and destructor
0028 //
0029 
0030 L1TTkEleFilter::L1TTkEleFilter(const edm::ParameterSet& iConfig)
0031     : HLTFilter(iConfig),
0032       l1TkEleTag1_(iConfig.getParameter<edm::InputTag>("inputTag1")),
0033       l1TkEleTag2_(iConfig.getParameter<edm::InputTag>("inputTag2")),
0034       tkEleToken1_(consumes<TkEleCollection>(l1TkEleTag1_)),
0035       tkEleToken2_(consumes<TkEleCollection>(l1TkEleTag2_)) {
0036   min_Pt_ = iConfig.getParameter<double>("MinPt");
0037   min_N_ = iConfig.getParameter<int>("MinN");
0038   min_Eta_ = iConfig.getParameter<double>("MinEta");
0039   max_Eta_ = iConfig.getParameter<double>("MaxEta");
0040   scalings_ = iConfig.getParameter<edm::ParameterSet>("Scalings");
0041   barrelScalings_ = scalings_.getParameter<std::vector<double> >("barrel");
0042   endcapScalings_ = scalings_.getParameter<std::vector<double> >("endcap");
0043   etaBinsForIsolation_ = iConfig.getParameter<std::vector<double> >("EtaBinsForIsolation");
0044   trkIsolation_ = iConfig.getParameter<std::vector<double> >("TrkIsolation");
0045   quality1_ = iConfig.getParameter<int>("Quality1");
0046   quality2_ = iConfig.getParameter<int>("Quality2");
0047   qual1IsMask_ = iConfig.getParameter<bool>("Qual1IsMask");
0048   qual2IsMask_ = iConfig.getParameter<bool>("Qual2IsMask");
0049   applyQual1_ = iConfig.getParameter<bool>("ApplyQual1");
0050   applyQual2_ = iConfig.getParameter<bool>("ApplyQual2");
0051 
0052   if (etaBinsForIsolation_.size() != (trkIsolation_.size() + 1))
0053     throw cms::Exception("ConfigurationError")
0054         << "Vector of isolation values should have same size of vector of eta bins plus one.";
0055 }
0056 
0057 L1TTkEleFilter::~L1TTkEleFilter() = default;
0058 
0059 //
0060 // member functions
0061 //
0062 
0063 void L1TTkEleFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0064   edm::ParameterSetDescription desc;
0065 
0066   makeHLTFilterDescription(desc);
0067   desc.add<double>("MinPt", -1.0);
0068   desc.add<double>("MinEta", -5.0);
0069   desc.add<double>("MaxEta", 5.0);
0070   desc.add<int>("MinN", 1);
0071   desc.add<edm::InputTag>("inputTag1", edm::InputTag("L1TkElectrons1"));
0072   desc.add<edm::InputTag>("inputTag2", edm::InputTag("L1TkElectrons2"));
0073   desc.add<std::vector<double> >("EtaBinsForIsolation", {0.0, 1.479});
0074   desc.add<std::vector<double> >("TrkIsolation",
0075                                  {
0076                                      99999.9,
0077                                  });
0078   desc.add<int>("Quality1", 0);
0079   desc.add<int>("Quality2", 0);
0080   desc.add<bool>("Qual1IsMask", false);
0081   desc.add<bool>("Qual2IsMask", false);
0082   desc.add<bool>("ApplyQual1", false);
0083   desc.add<bool>("ApplyQual2", false);
0084 
0085   edm::ParameterSetDescription descScalings;
0086   descScalings.add<std::vector<double> >("barrel", {0.0, 1.0, 0.0});
0087   descScalings.add<std::vector<double> >("endcap", {0.0, 1.0, 0.0});
0088   desc.add<edm::ParameterSetDescription>("Scalings", descScalings);
0089 
0090   descriptions.add("L1TTkEleFilter", desc);
0091 }
0092 
0093 // ------------ method called to produce the data  ------------
0094 bool L1TTkEleFilter::hltFilter(edm::Event& iEvent,
0095                                const edm::EventSetup& iSetup,
0096                                trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0097   using namespace std;
0098   using namespace edm;
0099   using namespace reco;
0100   using namespace trigger;
0101 
0102   // All HLT filters must create and fill an HLT filter object,
0103   // recording any reconstructed physics objects satisfying (or not)
0104   // this HLT filter, and place it in the Event.
0105 
0106   // The filter object
0107   if (saveTags()) {
0108     filterproduct.addCollectionTag(l1TkEleTag1_);
0109     filterproduct.addCollectionTag(l1TkEleTag2_);
0110   }
0111 
0112   // Specific filter code
0113 
0114   // get hold of products from Event
0115 
0116   /// Barrel colleciton
0117   Handle<l1t::TkElectronCollection> tkEles1;
0118   iEvent.getByToken(tkEleToken1_, tkEles1);
0119 
0120   /// Endcap collection
0121   Handle<l1t::TkElectronCollection> tkEles2;
0122   iEvent.getByToken(tkEleToken2_, tkEles2);
0123 
0124   int ntrkEle(0);
0125   // Loop over first collection
0126   auto atrkEles(tkEles1->begin());
0127   auto otrkEles(tkEles1->end());
0128   TkEleCollection::const_iterator itkEle;
0129   for (itkEle = atrkEles; itkEle != otrkEles; itkEle++) {
0130     double offlinePt = this->TkEleOfflineEt(itkEle->pt(), itkEle->eta());
0131     bool passQuality(false);
0132     bool passIsolation(false);
0133 
0134     if (applyQual1_) {
0135       if (qual1IsMask_)
0136         passQuality = (itkEle->EGRef()->hwQual() & quality1_);
0137       else
0138         passQuality = (itkEle->EGRef()->hwQual() == quality1_);
0139     } else
0140       passQuality = true;
0141 
0142     // There has to be a better way to do this.
0143     for (unsigned int etabin = 1; etabin != etaBinsForIsolation_.size(); ++etabin) {
0144       if (std::abs(itkEle->eta()) < etaBinsForIsolation_.at(etabin) and
0145           std::abs(itkEle->eta()) > etaBinsForIsolation_.at(etabin - 1) and
0146           itkEle->trkIsol() < trkIsolation_.at(etabin - 1))
0147         passIsolation = true;
0148     }
0149 
0150     if (offlinePt >= min_Pt_ && itkEle->eta() <= max_Eta_ && itkEle->eta() >= min_Eta_ && passQuality &&
0151         passIsolation) {
0152       ntrkEle++;
0153       l1t::TkElectronRef ref1(l1t::TkElectronRef(tkEles1, distance(atrkEles, itkEle)));
0154       filterproduct.addObject(trigger::TriggerObjectType::TriggerL1TkEle, ref1);
0155     }
0156   }
0157 
0158   // Loop over second collection. Notice we don't reset ntrkEle
0159   atrkEles = tkEles2->begin();
0160   otrkEles = tkEles2->end();
0161   for (itkEle = atrkEles; itkEle != otrkEles; itkEle++) {
0162     double offlinePt = this->TkEleOfflineEt(itkEle->pt(), itkEle->eta());
0163     bool passQuality(false);
0164     bool passIsolation(false);
0165 
0166     if (applyQual2_) {
0167       if (qual2IsMask_)
0168         passQuality = (itkEle->EGRef()->hwQual() & quality2_);
0169       else
0170         passQuality = (itkEle->EGRef()->hwQual() == quality2_);
0171     } else
0172       passQuality = true;
0173 
0174     for (unsigned int etabin = 1; etabin != etaBinsForIsolation_.size(); ++etabin) {
0175       if (std::abs(itkEle->eta()) < etaBinsForIsolation_.at(etabin) and
0176           std::abs(itkEle->eta()) > etaBinsForIsolation_.at(etabin - 1) and
0177           itkEle->trkIsol() < trkIsolation_.at(etabin - 1))
0178         passIsolation = true;
0179     }
0180 
0181     if (offlinePt >= min_Pt_ && itkEle->eta() <= max_Eta_ && itkEle->eta() >= min_Eta_ && passQuality &&
0182         passIsolation) {
0183       ntrkEle++;
0184       l1t::TkElectronRef ref2(l1t::TkElectronRef(tkEles2, distance(atrkEles, itkEle)));
0185       filterproduct.addObject(trigger::TriggerObjectType::TriggerL1TkEle, ref2);
0186     }
0187   }
0188 
0189   // return with final filter decision
0190   const bool accept(ntrkEle >= min_N_);
0191   return accept;
0192 }
0193 
0194 double L1TTkEleFilter::TkEleOfflineEt(double Et, double Eta) const {
0195   if (std::abs(Eta) < 1.5)
0196     return (barrelScalings_.at(0) + Et * barrelScalings_.at(1) + Et * Et * barrelScalings_.at(2));
0197   else
0198     return (endcapScalings_.at(0) + Et * endcapScalings_.at(1) + Et * Et * endcapScalings_.at(2));
0199 }