File indexing completed on 2021-02-14 12:50:03
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "DataFormats/PatCandidates/interface/Jet.h"
0020 #include "DataFormats/JetReco/interface/Jet.h"
0021 #include "DataFormats/JetReco/interface/PFJet.h"
0022 #include "DataFormats/JetReco/interface/CaloJet.h"
0023
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "CommonTools/Utils/interface/StringObjectFunction.h"
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "DataFormats/Common/interface/ValueMap.h"
0028 #include "DataFormats/Math/interface/deltaR.h"
0029
0030 template <class T, class C = T>
0031 class JetDeltaRValueMapProducer : public edm::stream::EDProducer<> {
0032 public:
0033 typedef edm::ValueMap<float> JetValueMap;
0034
0035 JetDeltaRValueMapProducer(edm::ParameterSet const& params)
0036 : srcToken_(consumes<typename edm::View<T> >(params.getParameter<edm::InputTag>("src"))),
0037 matchedToken_(consumes<typename edm::View<C> >(params.getParameter<edm::InputTag>("matched"))),
0038 distMax_(params.getParameter<double>("distMax")),
0039 value_(params.existsAs<std::string>("value") ? params.getParameter<std::string>("value") : ""),
0040 values_(params.existsAs<std::vector<std::string> >("values")
0041 ? params.getParameter<std::vector<std::string> >("values")
0042 : std::vector<std::string>()),
0043 valueLabels_(params.existsAs<std::vector<std::string> >("valueLabels")
0044 ? params.getParameter<std::vector<std::string> >("valueLabels")
0045 : std::vector<std::string>()),
0046 lazyParser_(params.existsAs<bool>("lazyParser") ? params.getParameter<bool>("lazyParser") : false),
0047 multiValue_(false) {
0048 if (!value_.empty()) {
0049 evaluationMap_.insert(std::make_pair(
0050 value_, std::unique_ptr<StringObjectFunction<C> >(new StringObjectFunction<C>(value_, lazyParser_))));
0051 produces<JetValueMap>();
0052 }
0053
0054 if (!valueLabels_.empty() || !values_.empty()) {
0055 if (valueLabels_.size() == values_.size()) {
0056 multiValue_ = true;
0057 for (size_t i = 0; i < valueLabels_.size(); ++i) {
0058 evaluationMap_.insert(std::make_pair(
0059 valueLabels_[i],
0060 std::unique_ptr<StringObjectFunction<C> >(new StringObjectFunction<C>(values_[i], lazyParser_))));
0061 produces<JetValueMap>(valueLabels_[i]);
0062 }
0063 } else
0064 edm::LogWarning("ValueLabelMismatch")
0065 << "The number of value labels does not match the number of values. Values will not be evaluated.";
0066 }
0067 }
0068
0069 ~JetDeltaRValueMapProducer() override {}
0070
0071 private:
0072 void produce(edm::Event& iEvent, const edm::EventSetup& iSetup) override {
0073 edm::Handle<typename edm::View<T> > h_jets1;
0074 iEvent.getByToken(srcToken_, h_jets1);
0075 edm::Handle<typename edm::View<C> > h_jets2;
0076 iEvent.getByToken(matchedToken_, h_jets2);
0077
0078 std::vector<float> values(h_jets1->size(), -99999);
0079 std::map<std::string, std::vector<float> > valuesMap;
0080 if (multiValue_) {
0081 for (size_t i = 0; i < valueLabels_.size(); ++i)
0082 valuesMap.insert(std::make_pair(valueLabels_[i], std::vector<float>(h_jets1->size(), -99999)));
0083 }
0084 std::vector<bool> jets1_locks(h_jets1->size(), false);
0085
0086 for (typename edm::View<C>::const_iterator ibegin = h_jets2->begin(), iend = h_jets2->end(), ijet = ibegin;
0087 ijet != iend;
0088 ++ijet) {
0089 float matched_dR2 = 1e9;
0090 int matched_index = -1;
0091
0092 for (typename edm::View<T>::const_iterator jbegin = h_jets1->begin(), jend = h_jets1->end(), jjet = jbegin;
0093 jjet != jend;
0094 ++jjet) {
0095 int index = jjet - jbegin;
0096
0097 if (jets1_locks.at(index))
0098 continue;
0099
0100 float temp_dR2 = reco::deltaR2(ijet->eta(), ijet->phi(), jjet->eta(), jjet->phi());
0101 if (temp_dR2 < matched_dR2) {
0102 matched_dR2 = temp_dR2;
0103 matched_index = index;
0104 }
0105 }
0106
0107 if (matched_index >= 0) {
0108 if (matched_dR2 > distMax_ * distMax_)
0109 edm::LogInfo("MatchedJetsFarApart") << "Matched jets separated by dR greater than distMax=" << distMax_;
0110 else {
0111 jets1_locks.at(matched_index) = true;
0112 if (!value_.empty())
0113 values.at(matched_index) = (*(evaluationMap_.at(value_)))(*ijet);
0114 if (multiValue_) {
0115 for (size_t i = 0; i < valueLabels_.size(); ++i)
0116 valuesMap.at(valueLabels_[i]).at(matched_index) = (*(evaluationMap_.at(valueLabels_[i])))(*ijet);
0117 }
0118 }
0119 }
0120 }
0121
0122 if (!value_.empty()) {
0123 std::unique_ptr<JetValueMap> jetValueMap(new JetValueMap());
0124
0125 JetValueMap::Filler filler(*jetValueMap);
0126 filler.insert(h_jets1, values.begin(), values.end());
0127 filler.fill();
0128
0129
0130 iEvent.put(std::move(jetValueMap));
0131 }
0132 if (multiValue_) {
0133 for (size_t i = 0; i < valueLabels_.size(); ++i) {
0134 std::unique_ptr<JetValueMap> jetValueMap(new JetValueMap());
0135
0136 JetValueMap::Filler filler(*jetValueMap);
0137 filler.insert(h_jets1, valuesMap.at(valueLabels_[i]).begin(), valuesMap.at(valueLabels_[i]).end());
0138 filler.fill();
0139
0140
0141 iEvent.put(std::move(jetValueMap), valueLabels_[i]);
0142 }
0143 }
0144 }
0145
0146 const edm::EDGetTokenT<typename edm::View<T> > srcToken_;
0147 const edm::EDGetTokenT<typename edm::View<C> > matchedToken_;
0148 const double distMax_;
0149 const std::string value_;
0150 const std::vector<std::string> values_;
0151 const std::vector<std::string> valueLabels_;
0152 const bool lazyParser_;
0153 bool multiValue_;
0154 std::map<std::string, std::unique_ptr<const StringObjectFunction<C> > > evaluationMap_;
0155 };
0156
0157 typedef JetDeltaRValueMapProducer<reco::Jet> RecoJetDeltaRValueMapProducer;
0158 typedef JetDeltaRValueMapProducer<reco::Jet, pat::Jet> RecoJetToPatJetDeltaRValueMapProducer;
0159 typedef JetDeltaRValueMapProducer<pat::Jet> PatJetDeltaRValueMapProducer;
0160
0161 DEFINE_FWK_MODULE(RecoJetDeltaRValueMapProducer);
0162 DEFINE_FWK_MODULE(RecoJetToPatJetDeltaRValueMapProducer);
0163 DEFINE_FWK_MODULE(PatJetDeltaRValueMapProducer);