File indexing completed on 2023-03-17 10:45:37
0001 #ifndef UtilAlgos_NewMatcher_h
0002 #define UtilAlgos_NewMatcher_h
0003
0004
0005
0006
0007
0008 #include "FWCore/Framework/interface/global/EDProducer.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "CommonTools/UtilAlgos/interface/DeltaR.h"
0012 #include "CommonTools/UtilAlgos/interface/MasterCollectionHelper.h"
0013 #include "FWCore/Framework/interface/Event.h"
0014 #include "DataFormats/Common/interface/Handle.h"
0015 #include "DataFormats/Common/interface/Association.h"
0016 #include "DataFormats/Common/interface/getRef.h"
0017 #include "DataFormats/Common/interface/View.h"
0018
0019 namespace reco {
0020 namespace modulesNew {
0021
0022 template <typename C1, typename C2, typename S, typename D = DeltaR<typename C1::value_type, typename C2::value_type> >
0023 class Matcher : public edm::global::EDProducer<> {
0024 public:
0025 Matcher(const edm::ParameterSet& cfg);
0026 ~Matcher() override;
0027
0028 private:
0029 typedef typename C1::value_type T1;
0030 typedef typename C2::value_type T2;
0031 typedef edm::Association<C2> MatchMap;
0032 void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0033 edm::EDGetTokenT<C1> srcToken_;
0034 edm::EDGetTokenT<C2> matchedToken_;
0035 double distMin_;
0036 double matchDistance(const T1& c1, const T2& c2) const { return distance_(c1, c2); }
0037 bool select(const T1& c1, const T2& c2) const { return select_(c1, c2); }
0038 S select_;
0039 D distance_;
0040 };
0041
0042 namespace helper {
0043 typedef std::pair<size_t, double> MatchPair;
0044
0045 struct SortBySecond {
0046 bool operator()(const MatchPair& p1, const MatchPair& p2) const { return p1.second < p2.second; }
0047 };
0048 }
0049
0050 template <typename C1, typename C2, typename S, typename D>
0051 Matcher<C1, C2, S, D>::Matcher(const edm::ParameterSet& cfg)
0052 : srcToken_(consumes<C1>(cfg.template getParameter<edm::InputTag>("src"))),
0053 matchedToken_(consumes<C2>(cfg.template getParameter<edm::InputTag>("matched"))),
0054 distMin_(cfg.template getParameter<double>("distMin")),
0055 select_(reco::modules::make<S>(cfg)),
0056 distance_(reco::modules::make<D>(cfg)) {
0057 produces<MatchMap>();
0058 }
0059
0060 template <typename C1, typename C2, typename S, typename D>
0061 Matcher<C1, C2, S, D>::~Matcher() {}
0062
0063 template <typename C1, typename C2, typename S, typename D>
0064 void Matcher<C1, C2, S, D>::produce(edm::StreamID, edm::Event& evt, const edm::EventSetup&) const {
0065 using namespace edm;
0066 using namespace std;
0067 Handle<C2> matched;
0068 evt.getByToken(matchedToken_, matched);
0069 Handle<C1> cands;
0070 evt.getByToken(srcToken_, cands);
0071 unique_ptr<MatchMap> matchMap(new MatchMap(matched));
0072 size_t size = cands->size();
0073 if (size != 0) {
0074 typename MatchMap::Filler filler(*matchMap);
0075 ::helper::MasterCollection<C1> master(cands, evt);
0076 std::vector<int> indices(master.size(), -1);
0077 for (size_t c = 0; c != size; ++c) {
0078 const T1& cand = (*cands)[c];
0079 vector<helper::MatchPair> v;
0080 for (size_t m = 0; m != matched->size(); ++m) {
0081 const T2& match = (*matched)[m];
0082 if (select(cand, match)) {
0083 double dist = matchDistance(cand, match);
0084 if (dist < distMin_)
0085 v.push_back(make_pair(m, dist));
0086 }
0087 }
0088 if (!v.empty()) {
0089 size_t idx = master.index(c);
0090 assert(idx < indices.size());
0091 indices[idx] = min_element(v.begin(), v.end(), helper::SortBySecond())->first;
0092 }
0093 }
0094 filler.insert(master.get(), indices.begin(), indices.end());
0095 filler.fill();
0096 }
0097 evt.put(std::move(matchMap));
0098 }
0099
0100 }
0101 }
0102
0103 #endif