File indexing completed on 2021-02-14 14:20:32
0001 #include <string>
0002 #include <vector>
0003
0004 #include "TVector3.h"
0005
0006 #include "FWCore/Framework/interface/MakerMacros.h"
0007 #include "DataFormats/RecoCandidate/interface/RecoEcalCandidate.h"
0008 #include "DataFormats/EgammaCandidates/interface/Electron.h"
0009 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidate.h"
0010 #include "DataFormats/RecoCandidate/interface/RecoChargedCandidateFwd.h"
0011 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
0012 #include "DataFormats/Common/interface/Handle.h"
0013 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0014 #include "DataFormats/Math/interface/deltaR.h"
0015 #include "DataFormats/JetReco/interface/PFJet.h"
0016 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0017 #include "HLTrigger/JetMET/interface/HLTJetCollectionsForBoostedLeptonPlusJets.h"
0018 #include "HLTrigger/HLTcore/interface/defaultModuleLabel.h"
0019 #include "CommonTools/Utils/interface/PtComparator.h"
0020
0021 template <typename jetType>
0022 HLTJetCollectionsForBoostedLeptonPlusJets<jetType>::HLTJetCollectionsForBoostedLeptonPlusJets(
0023 const edm::ParameterSet& iConfig)
0024 : hltLeptonTag(iConfig.getParameter<edm::InputTag>("HltLeptonTag")),
0025 sourceJetTag(iConfig.getParameter<edm::InputTag>("SourceJetTag")),
0026 minDeltaR_(iConfig.getParameter<double>("minDeltaR")) {
0027 using namespace edm;
0028 using namespace std;
0029
0030 typedef vector<RefVector<vector<jetType>, jetType, refhelper::FindUsingAdvance<vector<jetType>, jetType>>>
0031 JetCollectionVector;
0032 typedef vector<jetType> JetCollection;
0033
0034 m_theLeptonToken = consumes<trigger::TriggerFilterObjectWithRefs>(hltLeptonTag);
0035 m_theJetToken = consumes<std::vector<jetType>>(sourceJetTag);
0036 produces<JetCollectionVector>();
0037 produces<JetCollection>();
0038 }
0039
0040 template <typename jetType>
0041 HLTJetCollectionsForBoostedLeptonPlusJets<jetType>::~HLTJetCollectionsForBoostedLeptonPlusJets() {
0042
0043
0044 }
0045
0046 template <typename jetType>
0047 void HLTJetCollectionsForBoostedLeptonPlusJets<jetType>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0048 edm::ParameterSetDescription desc;
0049 desc.add<edm::InputTag>("HltLeptonTag", edm::InputTag("triggerFilterObjectWithRefs"));
0050
0051 desc.add<edm::InputTag>("SourceJetTag", edm::InputTag("jetCollection"));
0052
0053 desc.add<double>("minDeltaR", 0.5);
0054 descriptions.add(defaultModuleLabel<HLTJetCollectionsForBoostedLeptonPlusJets<jetType>>(), desc);
0055 }
0056
0057
0058
0059
0060
0061
0062
0063 template <typename jetType>
0064 void HLTJetCollectionsForBoostedLeptonPlusJets<jetType>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0065 using namespace edm;
0066 using namespace std;
0067
0068 using namespace reco;
0069
0070 typedef vector<RefVector<vector<jetType>, jetType, refhelper::FindUsingAdvance<vector<jetType>, jetType>>>
0071 JetCollectionVector;
0072 typedef vector<jetType> JetCollection;
0073 typedef edm::Ref<JetCollection> JetRef;
0074 typedef edm::RefVector<JetCollection> JetRefVector;
0075
0076 Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0077 iEvent.getByToken(m_theLeptonToken, PrevFilterOutput);
0078
0079
0080
0081 vector<reco::RecoChargedCandidateRef> muonCands;
0082 PrevFilterOutput->getObjects(trigger::TriggerMuon, muonCands);
0083
0084 vector<Ref<reco::ElectronCollection>> eleCands;
0085 PrevFilterOutput->getObjects(trigger::TriggerElectron, eleCands);
0086
0087 trigger::VRphoton photonCands;
0088 PrevFilterOutput->getObjects(trigger::TriggerPhoton, photonCands);
0089
0090 vector<Ref<reco::RecoEcalCandidateCollection>> clusCands;
0091 PrevFilterOutput->getObjects(trigger::TriggerCluster, clusCands);
0092
0093 Handle<JetCollection> theJetCollectionHandle;
0094 iEvent.getByToken(m_theJetToken, theJetCollectionHandle);
0095
0096 typename JetCollection::const_iterator jet;
0097
0098 unique_ptr<JetCollection> allSelections(new JetCollection);
0099 unique_ptr<JetCollectionVector> product(new JetCollectionVector);
0100
0101 std::vector<size_t> usedCands;
0102
0103 if (!muonCands.empty()) {
0104 for (jet = theJetCollectionHandle->begin(); jet != theJetCollectionHandle->end(); jet++) {
0105
0106 jetType cleanedJet = *jet;
0107 for (size_t candNr = 0; candNr < muonCands.size(); candNr++) {
0108 if (std::find(usedCands.begin(), usedCands.end(), candNr) != usedCands.end())
0109 continue;
0110 if (deltaR((*muonCands[candNr]), cleanedJet) <= minDeltaR_) {
0111 std::vector<edm::Ptr<reco::PFCandidate>> pfConstituents = cleanedJet.getPFConstituents();
0112 for (std::vector<edm::Ptr<reco::PFCandidate>>::const_iterator i_candidate = pfConstituents.begin();
0113 i_candidate != pfConstituents.end();
0114 ++i_candidate) {
0115 if (deltaR((*muonCands[candNr]), (**i_candidate)) < 0.001) {
0116 cleanedJet.setP4(cleanedJet.p4() - muonCands[candNr]->p4());
0117 usedCands.push_back(candNr);
0118 break;
0119 }
0120 }
0121 }
0122 }
0123 allSelections->push_back(cleanedJet);
0124 }
0125 }
0126
0127 if (!eleCands.empty()) {
0128 for (jet = theJetCollectionHandle->begin(); jet != theJetCollectionHandle->end(); jet++) {
0129
0130 jetType cleanedJet = *jet;
0131 for (size_t candNr = 0; candNr < eleCands.size(); candNr++) {
0132 if (std::find(usedCands.begin(), usedCands.end(), candNr) != usedCands.end())
0133 continue;
0134 if (deltaR((*eleCands[candNr]), cleanedJet) <= minDeltaR_) {
0135 std::vector<edm::Ptr<reco::PFCandidate>> pfConstituents = cleanedJet.getPFConstituents();
0136 for (std::vector<edm::Ptr<reco::PFCandidate>>::const_iterator i_candidate = pfConstituents.begin();
0137 i_candidate != pfConstituents.end();
0138 ++i_candidate) {
0139 if (deltaR((*eleCands[candNr]), (**i_candidate)) < 0.001) {
0140 cleanedJet.setP4(cleanedJet.p4() - eleCands[candNr]->p4());
0141 usedCands.push_back(candNr);
0142 break;
0143 }
0144 }
0145 }
0146 }
0147 allSelections->push_back(cleanedJet);
0148 }
0149 }
0150
0151 if (!photonCands.empty()) {
0152 for (jet = theJetCollectionHandle->begin(); jet != theJetCollectionHandle->end(); jet++) {
0153
0154 jetType cleanedJet = *jet;
0155 for (size_t candNr = 0; candNr < photonCands.size(); candNr++) {
0156 if (std::find(usedCands.begin(), usedCands.end(), candNr) != usedCands.end())
0157 continue;
0158 if (deltaR((*photonCands[candNr]), cleanedJet) <= minDeltaR_) {
0159 std::vector<edm::Ptr<reco::PFCandidate>> pfConstituents = cleanedJet.getPFConstituents();
0160 for (std::vector<edm::Ptr<reco::PFCandidate>>::const_iterator i_candidate = pfConstituents.begin();
0161 i_candidate != pfConstituents.end();
0162 ++i_candidate) {
0163 if (deltaR((*photonCands[candNr]), (**i_candidate)) < 0.001) {
0164 cleanedJet.setP4(cleanedJet.p4() - photonCands[candNr]->p4());
0165 usedCands.push_back(candNr);
0166 break;
0167 }
0168 }
0169 }
0170 }
0171 allSelections->push_back(cleanedJet);
0172 }
0173 }
0174
0175 if (!clusCands.empty()) {
0176 for (jet = theJetCollectionHandle->begin(); jet != theJetCollectionHandle->end(); jet++) {
0177
0178 jetType cleanedJet = *jet;
0179 for (size_t candNr = 0; candNr < clusCands.size(); candNr++) {
0180 if (std::find(usedCands.begin(), usedCands.end(), candNr) != usedCands.end())
0181 continue;
0182 if (deltaR((*clusCands[candNr]), cleanedJet) <= minDeltaR_) {
0183 std::vector<edm::Ptr<reco::PFCandidate>> pfConstituents = cleanedJet.getPFConstituents();
0184 for (std::vector<edm::Ptr<reco::PFCandidate>>::const_iterator i_candidate = pfConstituents.begin();
0185 i_candidate != pfConstituents.end();
0186 ++i_candidate) {
0187 if (deltaR((*clusCands[candNr]), (**i_candidate)) < 0.001) {
0188 cleanedJet.setP4(cleanedJet.p4() - clusCands[candNr]->p4());
0189 usedCands.push_back(candNr);
0190 break;
0191 }
0192 }
0193 }
0194 }
0195 allSelections->push_back(cleanedJet);
0196 }
0197 }
0198
0199 NumericSafeGreaterByPt<jetType> compJets;
0200
0201 std::sort(allSelections->begin(), allSelections->end(), compJets);
0202 edm::OrphanHandle<JetCollection> cleanedJetHandle = iEvent.put(std::move(allSelections));
0203
0204 JetCollection const& jets = *cleanedJetHandle;
0205
0206 JetRefVector cleanedJetRefs;
0207
0208 for (unsigned iJet = 0; iJet < jets.size(); ++iJet) {
0209 cleanedJetRefs.push_back(JetRef(cleanedJetHandle, iJet));
0210 }
0211
0212 product->emplace_back(cleanedJetRefs);
0213 iEvent.put(std::move(product));
0214
0215 return;
0216 }