File indexing completed on 2021-02-14 14:20:32
0001 #include "HLTrigger/JetMET/interface/HLTJetCollectionsForLeptonPlusJets.h"
0002
0003 #include "FWCore/Framework/interface/MakerMacros.h"
0004
0005 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0006 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0007 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0008 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0009
0010 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0011
0012 #include "DataFormats/Common/interface/Handle.h"
0013
0014 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0015 #include "DataFormats/Math/interface/deltaR.h"
0016 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0017
0018 template <typename jetType>
0019 HLTJetCollectionsForLeptonPlusJets<jetType>::HLTJetCollectionsForLeptonPlusJets(const edm::ParameterSet& iConfig)
0020 : hltLeptonTag(iConfig.getParameter<edm::InputTag>("HltLeptonTag")),
0021 sourceJetTag(iConfig.getParameter<edm::InputTag>("SourceJetTag")),
0022 minDeltaR_(iConfig.getParameter<double>("minDeltaR")) {
0023 using namespace edm;
0024 using namespace std;
0025 typedef vector<RefVector<vector<jetType>, jetType, refhelper::FindUsingAdvance<vector<jetType>, jetType>>>
0026 JetCollectionVector;
0027 m_theLeptonToken = consumes<trigger::TriggerFilterObjectWithRefs>(hltLeptonTag);
0028 m_theJetToken = consumes<std::vector<jetType>>(sourceJetTag);
0029 produces<JetCollectionVector>();
0030 }
0031
0032 template <typename jetType>
0033 HLTJetCollectionsForLeptonPlusJets<jetType>::~HLTJetCollectionsForLeptonPlusJets() {
0034
0035
0036 }
0037
0038 template <typename jetType>
0039 void HLTJetCollectionsForLeptonPlusJets<jetType>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0040 edm::ParameterSetDescription desc;
0041 desc.add<edm::InputTag>("HltLeptonTag", edm::InputTag("triggerFilterObjectWithRefs"));
0042 desc.add<edm::InputTag>("SourceJetTag", edm::InputTag("caloJetCollection"));
0043 desc.add<double>("minDeltaR", 0.5);
0044 descriptions.add(defaultModuleLabel<HLTJetCollectionsForLeptonPlusJets<jetType>>(), desc);
0045 }
0046
0047
0048
0049
0050
0051
0052
0053 template <typename jetType>
0054 void HLTJetCollectionsForLeptonPlusJets<jetType>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0055 using namespace edm;
0056 using namespace std;
0057
0058 typedef vector<RefVector<vector<jetType>, jetType, refhelper::FindUsingAdvance<vector<jetType>, jetType>>>
0059 JetCollectionVector;
0060 typedef vector<jetType> JetCollection;
0061 typedef edm::RefVector<JetCollection> JetRefVector;
0062 typedef edm::Ref<JetCollection> JetRef;
0063
0064 Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0065 iEvent.getByToken(m_theLeptonToken, PrevFilterOutput);
0066
0067
0068
0069 vector<Ref<reco::RecoEcalCandidateCollection>> clusCands;
0070 PrevFilterOutput->getObjects(trigger::TriggerCluster, clusCands);
0071
0072 vector<Ref<reco::ElectronCollection>> eleCands;
0073 PrevFilterOutput->getObjects(trigger::TriggerElectron, eleCands);
0074
0075 trigger::VRphoton photonCands;
0076 PrevFilterOutput->getObjects(trigger::TriggerPhoton, photonCands);
0077
0078 vector<reco::RecoChargedCandidateRef> muonCands;
0079 PrevFilterOutput->getObjects(trigger::TriggerMuon, muonCands);
0080
0081 Handle<JetCollection> theJetCollectionHandle;
0082 iEvent.getByToken(m_theJetToken, theJetCollectionHandle);
0083
0084 const JetCollection& theJetCollection = *theJetCollectionHandle;
0085
0086 unique_ptr<JetCollectionVector> allSelections(new JetCollectionVector());
0087
0088 if (!clusCands.empty()) {
0089 for (auto& clusCand : clusCands) {
0090 JetRefVector refVector;
0091 for (unsigned int j = 0; j < theJetCollection.size(); j++) {
0092 if (deltaR(clusCand->superCluster()->position(), theJetCollection[j]) > minDeltaR_)
0093 refVector.push_back(JetRef(theJetCollectionHandle, j));
0094 }
0095 allSelections->push_back(refVector);
0096 }
0097 }
0098
0099 if (!eleCands.empty()) {
0100 for (auto& eleCand : eleCands) {
0101 JetRefVector refVector;
0102 for (unsigned int j = 0; j < theJetCollection.size(); j++) {
0103 if (deltaR(eleCand->superCluster()->position(), theJetCollection[j]) > minDeltaR_)
0104 refVector.push_back(JetRef(theJetCollectionHandle, j));
0105 }
0106 allSelections->push_back(refVector);
0107 }
0108 }
0109
0110 if (!photonCands.empty()) {
0111 for (auto& photonCand : photonCands) {
0112 JetRefVector refVector;
0113 for (unsigned int j = 0; j < theJetCollection.size(); j++) {
0114 if (deltaR(photonCand->superCluster()->position(), theJetCollection[j]) > minDeltaR_)
0115 refVector.push_back(JetRef(theJetCollectionHandle, j));
0116 }
0117 allSelections->push_back(refVector);
0118 }
0119 }
0120
0121 if (!muonCands.empty()) {
0122 for (auto& muonCand : muonCands) {
0123 JetRefVector refVector;
0124 for (unsigned int j = 0; j < theJetCollection.size(); j++) {
0125 if (deltaR(muonCand->p4(), theJetCollection[j]) > minDeltaR_)
0126 refVector.push_back(JetRef(theJetCollectionHandle, j));
0127 }
0128 allSelections->push_back(refVector);
0129 }
0130 }
0131
0132 iEvent.put(std::move(allSelections));
0133
0134 return;
0135 }