Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-03-25 23:59:49

0001 /** \class HLTEgammaL1TMatchFilterRegional
0002  *
0003  *
0004  *  \author Monica Vazquez Acosta (CERN)
0005  *
0006  */
0007 #include "HLTEgammaL1TMatchFilterRegional.h"
0008 #include "DataFormats/Common/interface/Handle.h"
0009 #include "DataFormats/HLTReco/interface/TriggerFilterObjectWithRefs.h"
0010 #include "DataFormats/L1Trigger/interface/Jet.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/Framework/interface/EventSetup.h"
0015 
0016 #define TWOPI 6.283185308
0017 //
0018 // constructors and destructor
0019 //
0020 HLTEgammaL1TMatchFilterRegional::HLTEgammaL1TMatchFilterRegional(const edm::ParameterSet& iConfig)
0021     : HLTFilter(iConfig) {
0022   candIsolatedTag_ = iConfig.getParameter<edm::InputTag>("candIsolatedTag");
0023   l1EGTag_ = iConfig.getParameter<edm::InputTag>(
0024       "l1IsolatedTag");  //will be renamed l1EGTag for step 2 of the new L1 migration
0025   candNonIsolatedTag_ = iConfig.getParameter<edm::InputTag>("candNonIsolatedTag");
0026   L1SeedFilterTag_ = iConfig.getParameter<edm::InputTag>("L1SeedFilterTag");
0027   l1JetsTag_ = iConfig.getParameter<edm::InputTag>(
0028       "l1CenJetsTag");  //will be renamed l1JetsTag for step 2 of the new L1 migration
0029   l1TausTag_ =
0030       iConfig.getParameter<edm::InputTag>("l1TausTag");  //will be renamed l1JetsTag for step 2 of the new L1 migration
0031   ncandcut_ = iConfig.getParameter<int>("ncandcut");
0032   doIsolated_ = iConfig.getParameter<bool>("doIsolated");
0033   region_eta_size_ = iConfig.getParameter<double>("region_eta_size");
0034   region_eta_size_ecap_ = iConfig.getParameter<double>("region_eta_size_ecap");
0035   region_phi_size_ = iConfig.getParameter<double>("region_phi_size");
0036   barrel_end_ = iConfig.getParameter<double>("barrel_end");
0037   endcap_end_ = iConfig.getParameter<double>("endcap_end");
0038 
0039   candIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candIsolatedTag_);
0040   if (!doIsolated_)
0041     candNonIsolatedToken_ = consumes<reco::RecoEcalCandidateCollection>(candNonIsolatedTag_);
0042   L1SeedFilterToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(L1SeedFilterTag_);
0043 }
0044 
0045 HLTEgammaL1TMatchFilterRegional::~HLTEgammaL1TMatchFilterRegional() = default;
0046 
0047 void HLTEgammaL1TMatchFilterRegional::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0048   edm::ParameterSetDescription desc;
0049   makeHLTFilterDescription(desc);
0050   desc.add<edm::InputTag>("candIsolatedTag", edm::InputTag("hltRecoIsolatedEcalCandidate"));
0051   desc.add<edm::InputTag>("l1IsolatedTag",
0052                           edm::InputTag("hltCaloStage2Digis"));  //rename for step 2 of the L1 migration
0053   desc.add<edm::InputTag>("candNonIsolatedTag", edm::InputTag("hltRecoNonIsolatedEcalCandidate"));
0054   desc.add<edm::InputTag>("l1NonIsolatedTag",
0055                           edm::InputTag("l1extraParticles", "NonIsolated"));  //drop for step 2 of the L1 migration
0056   desc.add<edm::InputTag>("L1SeedFilterTag", edm::InputTag("theL1SeedFilter"));
0057   desc.add<edm::InputTag>("l1CenJetsTag", edm::InputTag("hltCaloStage2Digis"));      //rename for step 2 of L1 migration
0058   desc.add<edm::InputTag>("l1TausTag", edm::InputTag("hltCaloStage2Digis", "Tau"));  //rename for step 2 of L1 migration
0059   desc.add<int>("ncandcut", 1);
0060   desc.add<bool>("doIsolated", true);
0061   desc.add<double>("region_eta_size", 0.522);
0062   desc.add<double>("region_eta_size_ecap", 1.0);
0063   desc.add<double>("region_phi_size", 1.044);
0064   desc.add<double>("barrel_end", 1.4791);
0065   desc.add<double>("endcap_end", 2.65);
0066   descriptions.add("HLTEgammaL1TMatchFilterRegional", desc);
0067 }
0068 
0069 // ------------ method called to produce the data  ------------
0070 //configuration:
0071 //doIsolated=true, only isolated superclusters are allowed to match isolated L1 seeds
0072 //doIsolated=false, isolated superclusters are allowed to match either iso or non iso L1 seeds, non isolated superclusters are allowed only to match non-iso seeds. If no collection name is given for non-isolated superclusters, assumes the the isolated collection contains all (both iso + non iso) seeded superclusters.
0073 bool HLTEgammaL1TMatchFilterRegional::hltFilter(edm::Event& iEvent,
0074                                                 const edm::EventSetup& iSetup,
0075                                                 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0076   // std::cout <<"runnr "<<iEvent.id().run()<<" event "<<iEvent.id().event()<<std::endl;
0077   using namespace trigger;
0078   //using namespace l1extra;
0079 
0080   if (saveTags()) {
0081     filterproduct.addCollectionTag(l1EGTag_);
0082     filterproduct.addCollectionTag(l1JetsTag_);
0083     filterproduct.addCollectionTag(l1TausTag_);
0084     filterproduct.addCollectionTag(candIsolatedTag_);
0085     if (!doIsolated_ && !candNonIsolatedTag_.label().empty()) {
0086       filterproduct.addCollectionTag(candNonIsolatedTag_);
0087     }
0088   }
0089 
0090   edm::Ref<reco::RecoEcalCandidateCollection> ref;
0091 
0092   // look at all candidates,  check cuts and add to filter object
0093   int n(0);
0094 
0095   // Get the recoEcalCandidates
0096   edm::Handle<reco::RecoEcalCandidateCollection> recoIsolecalcands;
0097   iEvent.getByToken(candIsolatedToken_, recoIsolecalcands);
0098 
0099   edm::Handle<trigger::TriggerFilterObjectWithRefs> L1SeedOutput;
0100   iEvent.getByToken(L1SeedFilterToken_, L1SeedOutput);
0101 
0102   std::vector<l1t::EGammaRef> l1EGs;
0103   L1SeedOutput->getObjects(TriggerL1EG, l1EGs);
0104 
0105   std::vector<l1t::JetRef> l1Jets;
0106   L1SeedOutput->getObjects(TriggerL1Jet, l1Jets);
0107 
0108   std::vector<l1t::TauRef> l1Taus;
0109   L1SeedOutput->getObjects(TriggerL1Tau, l1Taus);
0110 
0111   int countCand = 0;
0112   for (auto recoecalcand = recoIsolecalcands->begin(); recoecalcand != recoIsolecalcands->end(); recoecalcand++) {
0113     countCand++;
0114     if (std::abs(recoecalcand->eta()) < endcap_end_) {
0115       //SC should be inside the ECAL fiducial volume
0116 
0117       //now EGamma is just one collection so automatically matches to Isolated and NonIsolated Seeds
0118       //it is assumed the HLTL1TSeed module fills it with the correct seeds
0119       bool matchedSCEG = matchedToL1Cand(l1EGs, recoecalcand->eta(), recoecalcand->phi());
0120       bool matchedSCJet = matchedToL1Cand(l1Jets, recoecalcand->eta(), recoecalcand->phi());
0121       bool matchedSCTau = matchedToL1Cand(l1Taus, recoecalcand->eta(), recoecalcand->phi());
0122 
0123       if (matchedSCEG || matchedSCJet || matchedSCTau) {
0124         n++;
0125         ref = edm::Ref<reco::RecoEcalCandidateCollection>(recoIsolecalcands,
0126                                                           distance(recoIsolecalcands->begin(), recoecalcand));
0127         filterproduct.addObject(TriggerCluster, ref);
0128       }  //end  matched check
0129 
0130     }  //end endcap fiduical check
0131 
0132   }  //end loop over all isolated RecoEcalCandidates
0133 
0134   //if doIsolated_ is false now run over the nonisolated superclusters and EG
0135   //however in the case we have a single collection of superclusters containing both iso L1 and non iso L1 seeded superclusters,
0136   //we do not have a non isolated collection of superclusters so we have to protect against that
0137   if (!doIsolated_ && !candNonIsolatedTag_.label().empty()) {
0138     edm::Handle<reco::RecoEcalCandidateCollection> recoNonIsolecalcands;
0139     iEvent.getByToken(candNonIsolatedToken_, recoNonIsolecalcands);
0140 
0141     for (auto recoecalcand = recoNonIsolecalcands->begin(); recoecalcand != recoNonIsolecalcands->end();
0142          recoecalcand++) {
0143       countCand++;
0144 
0145       if (std::abs(recoecalcand->eta()) < endcap_end_) {
0146         bool matchedSCEG = matchedToL1Cand(l1EGs, recoecalcand->eta(), recoecalcand->phi());
0147         bool matchedSCJet = matchedToL1Cand(l1Jets, recoecalcand->eta(), recoecalcand->phi());
0148         bool matchedSCTau = matchedToL1Cand(l1Taus, recoecalcand->eta(), recoecalcand->phi());
0149 
0150         if (matchedSCEG || matchedSCJet || matchedSCTau) {
0151           n++;
0152           ref = edm::Ref<reco::RecoEcalCandidateCollection>(recoNonIsolecalcands,
0153                                                             distance(recoNonIsolecalcands->begin(), recoecalcand));
0154           filterproduct.addObject(TriggerCluster, ref);
0155         }  //end  matched check
0156 
0157       }  //end endcap fiduical check
0158 
0159     }  //end loop over all isolated RecoEcalCandidates
0160   }    //end doIsolatedCheck
0161 
0162   // filter decision
0163   bool accept(n >= ncandcut_);
0164 
0165   return accept;
0166 }
0167 
0168 bool HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::EGammaRef>& l1Cands,
0169                                                       const float scEta,
0170                                                       const float scPhi) const {
0171   for (auto const& l1Cand : l1Cands) {
0172     //ORCA matching method
0173     double etaBinLow = 0.;
0174     double etaBinHigh = 0.;
0175     if (std::abs(scEta) < barrel_end_) {
0176       etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0177       etaBinHigh = etaBinLow + region_eta_size_;
0178     } else {
0179       etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0180       etaBinHigh = etaBinLow + region_eta_size_ecap_;
0181     }
0182 
0183     float deltaphi = std::abs(scPhi - l1Cand->phi());
0184     if (deltaphi > TWOPI)
0185       deltaphi -= TWOPI;
0186     if (deltaphi > TWOPI / 2.)
0187       deltaphi = TWOPI - deltaphi;
0188 
0189     if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0190       return true;
0191     }
0192   }
0193   return false;
0194 }
0195 
0196 bool
0197 //HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1extra::L1JetParticleRef >& l1Cands,const float scEta,const float scPhi) const
0198 HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::JetRef>& l1Cands,
0199                                                  const float scEta,
0200                                                  const float scPhi) const {
0201   for (auto const& l1Cand : l1Cands) {
0202     //ORCA matching method
0203     double etaBinLow = 0.;
0204     double etaBinHigh = 0.;
0205     if (std::abs(scEta) < barrel_end_) {
0206       etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0207       etaBinHigh = etaBinLow + region_eta_size_;
0208     } else {
0209       etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0210       etaBinHigh = etaBinLow + region_eta_size_ecap_;
0211     }
0212 
0213     float deltaphi = std::abs(scPhi - l1Cand->phi());
0214     if (deltaphi > TWOPI)
0215       deltaphi -= TWOPI;
0216     if (deltaphi > TWOPI / 2.)
0217       deltaphi = TWOPI - deltaphi;
0218 
0219     if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0220       return true;
0221     }
0222   }
0223   return false;
0224 }
0225 bool
0226 //HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1extra::L1JetParticleRef >& l1Cands,const float scEta,const float scPhi) const
0227 HLTEgammaL1TMatchFilterRegional::matchedToL1Cand(const std::vector<l1t::TauRef>& l1Cands,
0228                                                  const float scEta,
0229                                                  const float scPhi) const {
0230   for (auto const& l1Cand : l1Cands) {
0231     //ORCA matching method
0232     double etaBinLow = 0.;
0233     double etaBinHigh = 0.;
0234     if (std::abs(scEta) < barrel_end_) {
0235       etaBinLow = l1Cand->eta() - region_eta_size_ / 2.;
0236       etaBinHigh = etaBinLow + region_eta_size_;
0237     } else {
0238       etaBinLow = l1Cand->eta() - region_eta_size_ecap_ / 2.;
0239       etaBinHigh = etaBinLow + region_eta_size_ecap_;
0240     }
0241 
0242     float deltaphi = std::abs(scPhi - l1Cand->phi());
0243     if (deltaphi > TWOPI)
0244       deltaphi -= TWOPI;
0245     if (deltaphi > TWOPI / 2.)
0246       deltaphi = TWOPI - deltaphi;
0247 
0248     if (scEta < etaBinHigh && scEta > etaBinLow && deltaphi < region_phi_size_ / 2.) {
0249       return true;
0250     }
0251   }
0252   return false;
0253 }
0254 
0255 // declare this class as a framework plugin
0256 #include "FWCore/Framework/interface/MakerMacros.h"
0257 DEFINE_FWK_MODULE(HLTEgammaL1TMatchFilterRegional);