Jet updates for Mini and (JME)Nano, add recipe for puppi JetMET reclustering

PhysicsTools/NanoAOD/plugins/JetPFConstituentVarProducer.cc

@@ -0,0 +1,264 @@
+#include <string>
+#include <memory>
+
+#include "FWCore/Framework/interface/Frameworkfwd.h"
+#include "FWCore/Framework/interface/stream/EDProducer.h"
+#include "FWCore/Framework/interface/Event.h"
+#include "FWCore/Framework/interface/MakerMacros.h"
+
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
+#include "FWCore/Utilities/interface/StreamID.h"
+
+#include "DataFormats/Common/interface/View.h"
+#include "DataFormats/Common/interface/ValueMap.h"
+#include "DataFormats/PatCandidates/interface/Jet.h"
+#include "DataFormats/PatCandidates/interface/PackedCandidate.h"
+#include "DataFormats/Candidate/interface/CandidateFwd.h"
+#include "DataFormats/NanoAOD/interface/FlatTable.h"
+
+template <typename T>
+class JetPFConstituentVarProducer : public edm::stream::EDProducer<> {
+public:
+  explicit JetPFConstituentVarProducer(const edm::ParameterSet&);
+  ~JetPFConstituentVarProducer() override;
+  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
+
+private:
+  void produce(edm::Event&, const edm::EventSetup&) override;
+  void PutValueMapInEvent(edm::Event&, const edm::Handle<edm::View<T>>&, const std::vector<float>&, std::string);
+
+  edm::EDGetTokenT<edm::View<T>> jet_token_;
+  edm::EDGetTokenT<edm::ValueMap<float>> puppi_value_map_token_;
+
+  const bool fallback_puppi_weight_;
+  bool use_puppi_value_map_;
+};
+
+//
+// constructors and destructor
+//
+template <typename T>
+JetPFConstituentVarProducer<T>::JetPFConstituentVarProducer(const edm::ParameterSet& iConfig)
+    : jet_token_(consumes<edm::View<T>>(iConfig.getParameter<edm::InputTag>("jets"))),
+      fallback_puppi_weight_(iConfig.getParameter<bool>("fallback_puppi_weight")),
+      use_puppi_value_map_(false) {
+  //
+  // Puppi Value Map
+  //
+  const auto& puppi_value_map_tag = iConfig.getParameter<edm::InputTag>("puppi_value_map");
+  if (!puppi_value_map_tag.label().empty()) {
+    puppi_value_map_token_ = consumes<edm::ValueMap<float>>(puppi_value_map_tag);
+    use_puppi_value_map_ = true;
+  }
+  produces<edm::ValueMap<float>>("leadConstNeHadEF");
+  produces<edm::ValueMap<float>>("leadConstChHadEF");
+  produces<edm::ValueMap<float>>("leadConstPhotonEF");
+  produces<edm::ValueMap<float>>("leadConstElectronEF");
+  produces<edm::ValueMap<float>>("leadConstMuonEF");
+  produces<edm::ValueMap<float>>("leadConstHFHADEF");
+  produces<edm::ValueMap<float>>("leadConstHFEMEF");
+  produces<edm::ValueMap<float>>("leadConstNeHadPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstChHadPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstPhotonPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstElectronPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstMuonPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstHFHADPuppiWeight");
+  produces<edm::ValueMap<float>>("leadConstHFEMPuppiWeight");
+}
+
+template <typename T>
+JetPFConstituentVarProducer<T>::~JetPFConstituentVarProducer() {}
+
+template <typename T>
+void JetPFConstituentVarProducer<T>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
+  // Input jets
+  auto jets = iEvent.getHandle(jet_token_);
+
+  // Get puppi value map
+  edm::Handle<edm::ValueMap<float>> puppi_value_map_;
+  if (use_puppi_value_map_) {
+    iEvent.getByToken(puppi_value_map_token_, puppi_value_map_);
+  }
+
+  std::vector<float> jet_leadConstNeHadEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstNeHadPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstChHadEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstChHadPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstPhotonEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstPhotonPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstElectronEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstElectronPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstMuonEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstMuonPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstHFHADEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstHFHADPuppiWeight(jets->size(), 0.f);
+
+  std::vector<float> jet_leadConstHFEMEF(jets->size(), 0.f);
+  std::vector<float> jet_leadConstHFEMPuppiWeight(jets->size(), 0.f);
+
+  // Loop over jet
+  for (std::size_t jet_idx = 0; jet_idx < jets->size(); jet_idx++) {
+    const auto& jet = (*jets)[jet_idx];
+
+    float jet_energy_raw = jet.energy();
+    if constexpr (std::is_same<T, pat::Jet>::value) {
+      jet_energy_raw = jet.correctedJet(0).energy();
+    }
+
+    reco::CandidatePtr leadConstNeHad;
+    reco::CandidatePtr leadConstChHad;
+    reco::CandidatePtr leadConstPhoton;
+    reco::CandidatePtr leadConstElectron;
+    reco::CandidatePtr leadConstMuon;
+    reco::CandidatePtr leadConstHFHAD;
+    reco::CandidatePtr leadConstHFEM;
+
+    float leadConstNeHadPuppiWeight = 1.f;
+    float leadConstChHadPuppiWeight = 1.f;
+    float leadConstPhotonPuppiWeight = 1.f;
+    float leadConstElectronPuppiWeight = 1.f;
+    float leadConstMuonPuppiWeight = 1.f;
+    float leadConstHFHADPuppiWeight = 1.f;
+    float leadConstHFEMPuppiWeight = 1.f;
+
+    //
+    // Loop over jet constituents
+    //
+    for (const reco::CandidatePtr& dau : jet.daughterPtrVector()) {
+      float puppiw = 1.f;
+
+      //
+      // Get Puppi weight from ValueMap, if provided.
+      //
+      if (use_puppi_value_map_) {
+        puppiw = (*puppi_value_map_)[dau];
+      } else if (!fallback_puppi_weight_) {
+        throw edm::Exception(edm::errors::InvalidReference, "PUPPI value map missing")
+            << "use fallback_puppi_weight option to use " << puppiw << " for cand as default";
+      }
+
+      //
+      // Find the highest energy constituents for each PF type
+      //
+      if (abs(dau->pdgId()) == 130) {
+        if (leadConstNeHad.isNull() ||
+            (puppiw * dau->energy() > leadConstNeHadPuppiWeight * leadConstNeHad->energy())) {
+          leadConstNeHad = dau;
+          leadConstNeHadPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 211) {
+        if (leadConstChHad.isNull() ||
+            (puppiw * dau->energy() > leadConstChHadPuppiWeight * leadConstChHad->energy())) {
+          leadConstChHad = dau;
+          leadConstChHadPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 22) {
+        if (leadConstPhoton.isNull() ||
+            (puppiw * dau->energy() > leadConstPhotonPuppiWeight * leadConstPhoton->energy())) {
+          leadConstPhoton = dau;
+          leadConstPhotonPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 11) {
+        if (leadConstElectron.isNull() ||
+            (puppiw * dau->energy() > leadConstElectronPuppiWeight * leadConstElectron->energy())) {
+          leadConstElectron = dau;
+          leadConstElectronPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 13) {
+        if (leadConstMuon.isNull() || (puppiw * dau->energy() > leadConstMuonPuppiWeight * leadConstMuon->energy())) {
+          leadConstMuon = dau;
+          leadConstMuonPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 1) {
+        if (leadConstHFHAD.isNull() ||
+            (puppiw * dau->energy() > leadConstHFHADPuppiWeight * leadConstHFHAD->energy())) {
+          leadConstHFHAD = dau;
+          leadConstHFHADPuppiWeight = puppiw;
+        }
+      } else if (abs(dau->pdgId()) == 2) {
+        if (leadConstHFEM.isNull() || (puppiw * dau->energy() > leadConstHFEMPuppiWeight * leadConstHFEM->energy())) {
+          leadConstHFEM = dau;
+          leadConstHFEMPuppiWeight = puppiw;
+        }
+      }
+    }  // End of Jet Constituents Loop
+
+    if (leadConstNeHad.isNonnull()) {
+      jet_leadConstNeHadEF[jet_idx] = (leadConstNeHad->energy() * leadConstNeHadPuppiWeight) / jet_energy_raw;
+      jet_leadConstNeHadPuppiWeight[jet_idx] = leadConstNeHadPuppiWeight;
+    }
+    if (leadConstChHad.isNonnull()) {
+      jet_leadConstChHadEF[jet_idx] = (leadConstChHad->energy() * leadConstChHadPuppiWeight) / jet_energy_raw;
+      jet_leadConstChHadPuppiWeight[jet_idx] = leadConstChHadPuppiWeight;
+    }
+    if (leadConstPhoton.isNonnull()) {
+      jet_leadConstPhotonEF[jet_idx] = (leadConstPhoton->energy() * leadConstPhotonPuppiWeight) / jet_energy_raw;
+      jet_leadConstPhotonPuppiWeight[jet_idx] = leadConstPhotonPuppiWeight;
+    }
+    if (leadConstElectron.isNonnull()) {
+      jet_leadConstElectronEF[jet_idx] = (leadConstElectron->energy() * leadConstElectronPuppiWeight) / jet_energy_raw;
+      jet_leadConstElectronPuppiWeight[jet_idx] = leadConstElectronPuppiWeight;
+    }
+    if (leadConstMuon.isNonnull()) {
+      jet_leadConstMuonEF[jet_idx] = (leadConstMuon->energy() * leadConstMuonPuppiWeight) / jet_energy_raw;
+      jet_leadConstMuonPuppiWeight[jet_idx] = leadConstMuonPuppiWeight;
+    }
+    if (leadConstHFHAD.isNonnull()) {
+      jet_leadConstHFHADEF[jet_idx] = (leadConstHFHAD->energy() * leadConstHFHADPuppiWeight) / jet_energy_raw;
+      jet_leadConstHFHADPuppiWeight[jet_idx] = leadConstHFHADPuppiWeight;
+    }
+    if (leadConstHFEM.isNonnull()) {
+      jet_leadConstHFEMEF[jet_idx] = (leadConstHFEM->energy() * leadConstHFEMPuppiWeight) / jet_energy_raw;
+      jet_leadConstHFEMPuppiWeight[jet_idx] = leadConstHFEMPuppiWeight;
+    }
+  }
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstNeHadEF, "leadConstNeHadEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstNeHadPuppiWeight, "leadConstNeHadPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstChHadEF, "leadConstChHadEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstChHadPuppiWeight, "leadConstChHadPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstPhotonEF, "leadConstPhotonEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstPhotonPuppiWeight, "leadConstPhotonPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstElectronEF, "leadConstElectronEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstElectronPuppiWeight, "leadConstElectronPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstMuonEF, "leadConstMuonEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstMuonPuppiWeight, "leadConstMuonPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstHFHADEF, "leadConstHFHADEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstHFHADPuppiWeight, "leadConstHFHADPuppiWeight");
+
+  PutValueMapInEvent(iEvent, jets, jet_leadConstHFEMEF, "leadConstHFEMEF");
+  PutValueMapInEvent(iEvent, jets, jet_leadConstHFEMPuppiWeight, "leadConstHFEMPuppiWeight");
+}
+template <typename T>
+void JetPFConstituentVarProducer<T>::PutValueMapInEvent(edm::Event& iEvent,
+                                                        const edm::Handle<edm::View<T>>& coll,
+                                                        const std::vector<float>& vec_var,
+                                                        std::string VMName) {
+  std::unique_ptr<edm::ValueMap<float>> VM(new edm::ValueMap<float>());
+  edm::ValueMap<float>::Filler fillerVM(*VM);
+  fillerVM.insert(coll, vec_var.begin(), vec_var.end());
+  fillerVM.fill();
+  iEvent.put(std::move(VM), VMName);
+}
+
+template <typename T>
+void JetPFConstituentVarProducer<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
+  edm::ParameterSetDescription desc;
+  desc.add<edm::InputTag>("jets", edm::InputTag("finalJetsPuppi"));
+  desc.add<edm::InputTag>("puppi_value_map", edm::InputTag("packedpuppi"));
+  desc.add<bool>("fallback_puppi_weight", false);
+  descriptions.addWithDefaultLabel(desc);
+}
+
+typedef JetPFConstituentVarProducer<pat::Jet> PatJetPFConstituentVarProducer;
+DEFINE_FWK_MODULE(PatJetPFConstituentVarProducer);

PhysicsTools/NanoAOD/python/autoNANO.py

@@ -37,7 +37,7 @@ def expandNanoMapping(seqList, mapping, key):
     'JME': {'sequence': '@PHYS',
             'customize': '@PHYS+PhysicsTools/NanoAOD/custom_jme_cff.PrepJMECustomNanoAOD'},
     'JMErePuppi': {'sequence': '@PHYS',
-                   'customize': '@PHYS+@JME+PhysicsTools/NanoAOD/custom_jme_cff.RecomputePuppiWeightsAndMET'},
+                   'customize': '@PHYS+@JME+PhysicsTools/NanoAOD/custom_jme_cff.RecomputePuppiWeightsMETAK8'},
     # L1 DPG (standalone with full calo TP info, L1T reemulation customization)
     'L1DPG' : {'sequence': 'DPGAnalysis/L1TNanoAOD/l1tNano_cff.l1tNanoSequence',
                'customize': ','.join(['PhysicsTools/NanoAOD/l1trig_cff.nanoL1TrigObjCustomizeFull',