Comparison PR for upgrade to 14_1

CommonTools/RecoAlgos/plugins/BuildFile.xml

@@ -16,4 +16,5 @@
   <use name="TrackingTools/Records"/>
   <use name="CommonTools/Utils"/>
   <use name="CommonTools/RecoAlgos"/>
+  <use name="JetMETCorrections/JetCorrector"/>
 </library>

CommonTools/RecoAlgos/plugins/HGCRecHitCollectionMerger.cc

@@ -0,0 +1,12 @@
+#include "FWCore/Framework/interface/MakerMacros.h"
+#include "CommonTools/UtilAlgos/interface/Merger.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHit.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
+#include "DataFormats/Common/interface/CloneTrait.h"
+#include "DataFormats/Common/interface/RefToBaseVector.h"
+#include "DataFormats/Common/interface/CopyPolicy.h"
+
+
+typedef Merger<HGCRecHitCollection, HGCRecHitCollection, edm::CopyPolicy<HGCRecHit>> HGCRecHitCollectionMerger;
+
+DEFINE_FWK_MODULE(HGCRecHitCollectionMerger);

CommonTools/RecoAlgos/plugins/RecHitToLayerClusterAssociationProducer.cc

@@ -0,0 +1,127 @@
+// system include files
+#include <memory>
+#include <string>
+
+// user include files
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+
+#include "FWCore/Framework/interface/ESHandle.h"
+
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
+#include "DataFormats/CaloRecHit/interface/CaloRecHit.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
+#include "DataFormats/CaloRecHit/interface/CaloCluster.h"
+
+#include "DataFormats/Common/interface/Association.h"
+#include "DataFormats/Common/interface/AssociationMap.h"
+#include "DataFormats/Common/interface/OneToManyWithQualityGeneric.h"
+
+#include "FWCore/Utilities/interface/transform.h"
+#include "FWCore/Utilities/interface/EDGetToken.h"
+#include <set>
+
+//
+// class decleration
+//
+typedef std::pair<size_t, float> IdxAndFraction;
+typedef edm::AssociationMap<edm::OneToManyWithQualityGeneric<
+    HGCRecHitCollection, std::vector<reco::CaloCluster>, float>> RecHitToLayerCluster;
+
+class RecHitToLayerClusterAssociationProducer : public edm::stream::EDProducer<> {
+public:
+  explicit RecHitToLayerClusterAssociationProducer(const edm::ParameterSet&);
+  ~RecHitToLayerClusterAssociationProducer() override;
+
+private:
+  void produce(edm::Event&, const edm::EventSetup&) override;
+
+  std::vector<edm::InputTag> caloRechitTags_;
+  std::vector<edm::EDGetTokenT<HGCRecHitCollection>> caloRechitCollectionTokens_;
+  edm::EDGetTokenT<std::vector<reco::CaloCluster>> layerClusterToken_;
+};
+
+RecHitToLayerClusterAssociationProducer::RecHitToLayerClusterAssociationProducer(const edm::ParameterSet& pset)
+    : caloRechitTags_(pset.getParameter<std::vector<edm::InputTag>>("caloRecHits")),
+      caloRechitCollectionTokens_(edm::vector_transform(
+          caloRechitTags_, [this](const edm::InputTag& tag) { return consumes<HGCRecHitCollection>(tag); })),
+      layerClusterToken_(consumes<std::vector<reco::CaloCluster>>(pset.getParameter<edm::InputTag>("layerClusters"))) {
+  for (auto& tag : caloRechitTags_) {
+    const std::string& label = !tag.instance().empty() ? tag.instance() : tag.label();
+    produces<edm::Association<std::vector<reco::CaloCluster>>>(label + "ToBestLayerCluster");
+    produces<RecHitToLayerCluster>(label + "ToLayerCluster");
+  }
+}
+
+RecHitToLayerClusterAssociationProducer::~RecHitToLayerClusterAssociationProducer() {}
+
+//
+// member functions
+//
+
+// ------------ method called to produce the data  ------------
+void RecHitToLayerClusterAssociationProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
+  edm::Handle<std::vector<reco::CaloCluster>> lcCollection;
+  iEvent.getByToken(layerClusterToken_, lcCollection);
+  std::unordered_map<size_t, std::vector<IdxAndFraction>> hitDetIdToIndices;
+
+  for (size_t j = 0; j < lcCollection->size(); j++) {
+      for (auto& idAndFrac : (*lcCollection)[j].hitsAndFractions()) {
+          hitDetIdToIndices[idAndFrac.first.rawId()].push_back({j, idAndFrac.second});
+      }
+  }
+
+  for (size_t i = 0; i < caloRechitCollectionTokens_.size(); i++) {
+    std::string label = caloRechitTags_.at(i).instance();
+	if (label.empty())
+		label = caloRechitTags_.at(i).label();
+    std::vector<int> rechitIndices;
+
+    edm::Handle<HGCRecHitCollection> caloRechitCollection;
+    iEvent.getByToken(caloRechitCollectionTokens_.at(i), caloRechitCollection);
+
+    auto assocMap = std::make_unique<RecHitToLayerCluster>(caloRechitCollection, lcCollection);
+
+    for (size_t h = 0; h < caloRechitCollection->size(); h++) {
+      HGCRecHitRef caloRh(caloRechitCollection, h);
+      size_t id = caloRh->detid().rawId();
+
+      // Need to sort before inserting into AssociationMap
+      auto match = hitDetIdToIndices.find(id);
+      if (match == std::end(hitDetIdToIndices)) {
+        rechitIndices.push_back(-1);
+        continue;
+      }
+      auto& lcIdxAndFrac = match->second;
+      // Sort by energy fraction
+      std::sort(std::begin(lcIdxAndFrac), std::end(lcIdxAndFrac), 
+            [](auto& a, auto& b) { return a.second > b.second; });
+
+      for (size_t m = 0; m < lcIdxAndFrac.size(); m++) {
+        float fraction = lcIdxAndFrac[m].second;
+        int lcIdx = lcIdxAndFrac[m].first;
+        // Best match is the layerCluster that carries the hit with the highest energy fraction
+        // (that is, the one responsible for the largest deposit in the detId)
+        if (m == 0)
+          rechitIndices.push_back(lcIdx);
+        edm::Ref<std::vector<reco::CaloCluster>> lc(lcCollection, lcIdx);
+        assocMap->insert(caloRh, std::make_pair(lc, fraction));
+      }
+    }
+
+    auto assoc = std::make_unique<edm::Association<std::vector<reco::CaloCluster>>>(lcCollection);
+    edm::Association<std::vector<reco::CaloCluster>>::Filler filler(*assoc);
+    filler.insert(caloRechitCollection, rechitIndices.begin(), rechitIndices.end());
+    filler.fill();
+    iEvent.put(std::move(assoc), label + "ToBestLayerCluster");
+    iEvent.put(std::move(assocMap), label + "ToLayerCluster");
+  }
+}
+
+// define this as a plug-in
+DEFINE_FWK_MODULE(RecHitToLayerClusterAssociationProducer);
+

CommonTools/RecoAlgos/plugins/RecHitToPFCandAssociationProducer.cc

@@ -0,0 +1,133 @@
+// system include files
+#include <memory>
+#include <string>
+
+// user include files
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+
+#include "FWCore/Framework/interface/ESHandle.h"
+
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+
+#include "FWCore/MessageLogger/interface/MessageLogger.h"
+#include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
+#include "SimDataFormats/CaloAnalysis/interface/SimClusterFwd.h"
+#include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
+#include "DataFormats/CaloRecHit/interface/CaloRecHit.h"
+#include "SimDataFormats/CaloHit/interface/PCaloHit.h"
+#include "SimDataFormats/CaloHit/interface/PCaloHitContainer.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
+
+#include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
+#include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
+#include "DataFormats/ParticleFlowReco/interface/PFRecHit.h"
+#include "DataFormats/ParticleFlowReco/interface/PFRecHitFraction.h"
+#include "DataFormats/ParticleFlowReco/interface/PFBlock.h"
+#include "DataFormats/ParticleFlowReco/interface/PFBlockFwd.h"
+#include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
+#include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
+
+#include "DataFormats/Common/interface/Association.h"
+
+#include "FWCore/Utilities/interface/transform.h"
+#include "FWCore/Utilities/interface/EDGetToken.h"
+#include <set>
+
+//
+// class decleration
+//
+typedef std::pair<size_t, float> IdxAndFraction;
+
+class RecHitToPFCandAssociationProducer : public edm::stream::EDProducer<> {
+public:
+  explicit RecHitToPFCandAssociationProducer(const edm::ParameterSet&);
+  ~RecHitToPFCandAssociationProducer() override;
+
+private:
+  void produce(edm::Event&, const edm::EventSetup&) override;
+
+  std::vector<edm::InputTag> caloRechitTags_;
+  std::vector<edm::EDGetTokenT<edm::PCaloHitContainer>> caloSimhitCollectionTokens_;
+  std::vector<edm::EDGetTokenT<edm::View<CaloRecHit>>> caloRechitCollectionTokens_;
+  edm::EDGetTokenT<reco::PFCandidateCollection> pfCollectionToken_;
+};
+
+RecHitToPFCandAssociationProducer::RecHitToPFCandAssociationProducer(const edm::ParameterSet& pset)
+    : caloRechitTags_(pset.getParameter<std::vector<edm::InputTag>>("caloRecHits")),
+      caloRechitCollectionTokens_(edm::vector_transform(
+          caloRechitTags_, [this](const edm::InputTag& tag) { return consumes<edm::View<CaloRecHit>>(tag); })),
+      pfCollectionToken_(consumes<reco::PFCandidateCollection>(pset.getParameter<edm::InputTag>("pfCands"))) {
+  for (auto& tag : caloRechitTags_) {
+    const std::string& label = !tag.instance().empty() ? tag.instance() : tag.label();
+    produces<edm::Association<reco::PFCandidateCollection>>(label + "ToPFCand");
+  }
+}
+
+RecHitToPFCandAssociationProducer::~RecHitToPFCandAssociationProducer() {}
+
+//
+// member functions
+//
+
+// ------------ method called to produce the data  ------------
+void RecHitToPFCandAssociationProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
+  edm::Handle<reco::PFCandidateCollection> pfCollection;
+  iEvent.getByToken(pfCollectionToken_, pfCollection);
+  std::unordered_map<size_t, IdxAndFraction> hitDetIdToIndex;
+
+  for (size_t j = 0; j < pfCollection->size(); j++) {
+    const auto& pfCand = pfCollection->at(j);
+    const reco::PFCandidate::ElementsInBlocks& elements = pfCand.elementsInBlocks();
+    for (auto& element : elements) {
+      const reco::PFBlockRef blockRef = element.first;
+      if (!blockRef.isNonnull())
+          continue;
+      for (const auto& block : blockRef->elements()) {
+        // This seems to not work for PFTICL
+        //if (block.type() == reco::PFBlockElement::HGCAL) {
+        const reco::PFClusterRef cluster = block.clusterRef();
+        if (cluster.isNonnull()) {
+          const std::vector<reco::PFRecHitFraction>& rhf = cluster->recHitFractions();
+          for (const auto& hf : rhf) {
+            auto& hit = hf.recHitRef();
+            if (!hit)
+              throw cms::Exception("RecHitToPFCandAssociationProducer") << "Invalid RecHit ref";
+            size_t detId = hit->detId();
+            auto entry = hitDetIdToIndex.find(detId);
+            if (entry == hitDetIdToIndex.end() || entry->second.second < hf.fraction())
+              hitDetIdToIndex[detId] = {j, hf.fraction()};
+          }
+        }
+      }
+    }
+  }
+
+  for (size_t i = 0; i < caloRechitCollectionTokens_.size(); i++) {
+    std::string label = caloRechitTags_.at(i).instance();
+	if (label.empty())
+		label = caloRechitTags_.at(i).label();
+    std::vector<size_t> rechitIndices;
+
+    edm::Handle<edm::View<CaloRecHit>> caloRechitCollection;
+    iEvent.getByToken(caloRechitCollectionTokens_.at(i), caloRechitCollection);
+
+    for (size_t h = 0; h < caloRechitCollection->size(); h++) {
+      const CaloRecHit& caloRh = caloRechitCollection->at(h);
+      size_t id = caloRh.detid().rawId();
+      int match = hitDetIdToIndex.find(id) == hitDetIdToIndex.end() ? -1 : hitDetIdToIndex.at(id).first;
+      rechitIndices.push_back(match);
+    }
+
+    auto assoc = std::make_unique<edm::Association<reco::PFCandidateCollection>>(pfCollection);
+    edm::Association<reco::PFCandidateCollection>::Filler filler(*assoc);
+    filler.insert(caloRechitCollection, rechitIndices.begin(), rechitIndices.end());
+    filler.fill();
+    iEvent.put(std::move(assoc), label + "ToPFCand");
+  }
+}
+
+// define this as a plug-in
+DEFINE_FWK_MODULE(RecHitToPFCandAssociationProducer);