Adding in slimmed superclusters for EGM efficiency measurements to Mini

DataFormats/EgammaReco/interface/SlimmedSuperCluster.h

@@ -0,0 +1,51 @@
+#ifndef DataFormats_PatCandidates_SuperCluster_h
+#define DataFormats_PatCandidates_SuperCluster_h
+
+/***
+
+This is a slimmed data format to be able to store all EGamma superclusters 
+in MiniAOD as well as the necessary ID variables for EG workflows
+
+The goal is to catch the superclusters that dont make it to electrons/photons (of which there are many) to enable efficiency measurements of this stage. 
+
+So initially it was thought to make this a LeafCandidate but a 
+supercluster is really a point + energy not a momentum so while more awkward not to just treat it as a p4, its more accurate to treat it as a position + energy
+
+author: Sam Harper, Swagata Mukherjee
+
+***/
+
+#include "DataFormats/DetId/interface/DetId.h"
+#include "DataFormats/Math/interface/Point3D.h"
+#include <vector>
+
+namespace reco {
+  class SuperCluster;
+}
+
+namespace reco {
+
+  class SlimmedSuperCluster {
+  public:
+    SlimmedSuperCluster() : rawEnergy_(0.), preshowerEnergy_(0.), trkIso_(0.) {}
+    SlimmedSuperCluster(const reco::SuperCluster&, float trkIso = 0.);
+    float rawEnergy() const { return rawEnergy_; }
+    float preshowerEnergy() const { return preshowerEnergy_; }
+    DetId seedId() const { return clusterSeedIds_.empty() != 0 ? DetId(0) : clusterSeedIds_.front(); }
+    const std::vector<DetId>& clusterSeedIds() const { return clusterSeedIds_; }
+    math::XYZPoint position() const;
+    float trkIso() const { return trkIso_; }
+    void setTrkIso(float val) { trkIso_ = val; }
+
+  private:
+    float correctedEnergy_;
+    float rawEnergy_;
+    float preshowerEnergy_;
+    float rho_;
+    float eta_;
+    float phi_;
+    std::vector<DetId> clusterSeedIds_;  //the first one is always the seed
+    float trkIso_;
+  };
+}  // namespace reco
+#endif

DataFormats/EgammaReco/interface/SlimmedSuperClusterFwd.h

@@ -0,0 +1,25 @@
+#ifndef EgammaReco_SlimmedSuperClusterFwd_h
+#define EgammaReco_SlimmedSuperClusterFwd_h
+#include <vector>
+#include "DataFormats/Common/interface/Ref.h"
+#include "DataFormats/Common/interface/RefVector.h"
+#include "DataFormats/Common/interface/RefProd.h"
+
+namespace reco {
+  class SlimmedSuperCluster;
+
+  /// collection of SuperCluser objectr
+  typedef std::vector<SlimmedSuperCluster> SlimmedSuperClusterCollection;
+
+  /// reference to an object in a collection of SlimmedSuperCluster objects
+  typedef edm::Ref<SlimmedSuperClusterCollection> SlimmedSuperClusterRef;
+
+  /// reference to a collection of SlimmedSuperCluster objects
+  typedef edm::RefProd<SlimmedSuperClusterCollection> SlimmedSuperClusterRefProd;
+
+  /// vector of references to objects in the same colletion of SlimmedSuperCluster objects
+  typedef edm::RefVector<SlimmedSuperClusterCollection> SlimmedSuperClusterRefVector;
+
+}  // namespace reco
+
+#endif

DataFormats/EgammaReco/src/SlimmedSuperCluster.cc

@@ -0,0 +1,27 @@
+#include "DataFormats/EgammaReco/interface/SlimmedSuperCluster.h"
+#include "DataFormats/EgammaReco/interface/SuperCluster.h"
+#include "Math/GenVector/PositionVector3D.h"
+
+reco::SlimmedSuperCluster::SlimmedSuperCluster(const reco::SuperCluster& sc, float trkIso)
+    : correctedEnergy_(sc.correctedEnergy()),
+      rawEnergy_(sc.rawEnergy()),
+      preshowerEnergy_(sc.preshowerEnergy()),
+      rho_(sc.position().rho()),
+      eta_(sc.eta()),
+      phi_(sc.phi()),
+      trkIso_(trkIso) {
+  clusterSeedIds_.push_back(sc.seed()->seed());
+  for (const auto& clus : sc.clusters()) {
+    if (clus != sc.seed()) {
+      clusterSeedIds_.push_back(clus->seed());
+    }
+  }
+}
+
+math::XYZPoint reco::SlimmedSuperCluster::position() const {
+  ROOT::Math::PositionVector3D<ROOT::Math::CylindricalEta3D<float>> pos;
+  pos.SetRho(rho_);
+  pos.SetEta(eta_);
+  pos.SetPhi(phi_);
+  return math::XYZPoint(pos.x(), pos.y(), pos.z());
+}

DataFormats/EgammaReco/src/classes.h

@@ -7,6 +7,7 @@
 #include "DataFormats/EgammaReco/interface/BasicClusterFwd.h"
 #include "DataFormats/EgammaReco/interface/PreshowerClusterShape.h"
 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
+#include "DataFormats/EgammaReco/interface/SlimmedSuperCluster.h"
 #include "Rtypes.h"
 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
 #include "DataFormats/EgammaReco/interface/ClusterShape.h"

DataFormats/EgammaReco/src/classes_def.xml

@@ -159,5 +159,19 @@
    <class name="std::vector<reco::SuperClusterRef>" />
    <class name="edm::ValueMap<reco::SuperClusterRef>" />
    <class name="edm::Wrapper<edm::ValueMap<reco::SuperClusterRef> >" />
+
+   <class name="reco::SlimmedSuperCluster" ClassVersion="10">
+    <version ClassVersion="10" checksum="3805996953"/>
+   </class>
+   <class name="std::vector<reco::SlimmedSuperCluster>"/>
+   <class name="edm::Wrapper<std::vector<reco::SlimmedSuperCluster> >"/>
+   <class name="edm::Ref<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster,edm::refhelper::FindUsingAdvance<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster> >"/>
+   <class name="edm::RefProd<std::vector<reco::SlimmedSuperCluster> >"/>
+   <class name="edm::RefVector<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster,edm::refhelper::FindUsingAdvance<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster> >"/>
+   <class name="edm::Wrapper<edm::RefVector<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster,edm::refhelper::FindUsingAdvance<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster> > >"/>
+   <class name="std::vector<edm::Ref<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster,edm::refhelper::FindUsingAdvance<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster> > >"/>
+   <class name="edm::Wrapper<std::vector<edm::Ref<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster,edm::refhelper::FindUsingAdvance<std::vector<reco::SlimmedSuperCluster>,reco::SlimmedSuperCluster> > > >"/>
+
+
 </lcgdict>
 

RecoEgamma/EgammaIsolationAlgos/interface/EgammaHLTTrackIsolation.h

@@ -38,9 +38,18 @@
 #include "DataFormats/Math/interface/Point3D.h"
 
 #include "FWCore/MessageLogger/interface/MessageLogger.h"
+#include "FWCore/ParameterSet/interface/ParameterSet.h"
 
 class EgammaHLTTrackIsolation {
 public:
+  EgammaHLTTrackIsolation(const edm::ParameterSet& iConfig)
+      : ptMin(iConfig.getParameter<double>("ptMin")),
+        conesize(iConfig.getParameter<double>("coneSize")),
+        zspan(iConfig.getParameter<double>("zspan")),
+        rspan(iConfig.getParameter<double>("rspan")),
+        vetoConesize(iConfig.getParameter<double>("vetoConeSize")),
+        stripBarrel(iConfig.getParameter<double>("stripBarrel")),
+        stripEndcap(iConfig.getParameter<double>("stripEndcap")) {}
   EgammaHLTTrackIsolation(double egTrkIso_PtMin,
                           double egTrkIso_ConeSize,
                           double egTrkIso_ZSpan,
@@ -54,45 +63,46 @@ class EgammaHLTTrackIsolation {
         rspan(egTrkIso_RSpan),
         vetoConesize(egTrkIso_VetoConeSize),
         stripBarrel(egTrkIso_stripBarrel),
-        stripEndcap(egTrkIso_stripEndcap) {
-    /*
-	std::cout << "EgammaHLTTrackIsolation instance:"
-	<< " ptMin=" << ptMin << " "
-	<< " conesize="<< conesize << " "
-	<< " zspan=" << zspan << " "
-	<< " rspan=" << rspan << " "
-	<< " vetoConesize="<< vetoConesize
-	<< std::endl;    
-      */
-  }
+        stripEndcap(egTrkIso_stripEndcap) {}
 
   /// Get number of tracks and Pt sum of tracks inside an isolation cone for electrons
-  std::pair<int, float> electronIsolation(const reco::Track* const tr, const reco::TrackCollection* isoTracks);
+  std::pair<int, float> electronIsolation(const reco::Track* const tr, const reco::TrackCollection* isoTracks) const;
   std::pair<int, float> electronIsolation(const reco::Track* const tr,
                                           const reco::ElectronCollection* allEle,
-                                          const reco::TrackCollection* isoTracks);
+                                          const reco::TrackCollection* isoTracks) const;
   std::pair<int, float> electronIsolation(const reco::Track* const tr,
                                           const reco::TrackCollection* isoTracks,
-                                          GlobalPoint vertex);
+                                          GlobalPoint vertex) const;
 
   /// Get number of tracks and Pt sum of tracks inside an isolation cone for photons
   /// set useVertex=true to use PhotonCandidate vertex from EgammaPhotonVtxFinder
   /// set useVertex=false to consider all tracks for isolation
   std::pair<int, float> photonIsolation(const reco::RecoCandidate* const recocand,
                                         const reco::TrackCollection* isoTracks,
-                                        bool useVertex);
+                                        bool useVertex) const;
   std::pair<int, float> photonIsolation(const reco::RecoCandidate* const recocand,
                                         const reco::TrackCollection* isoTracks,
-                                        GlobalPoint vertex);
+                                        GlobalPoint vertex) const;
+  std::pair<int, float> photonIsolation(const reco::RecoCandidate::Point& pos,
+                                        const reco::TrackCollection* isoTracks,
+                                        GlobalPoint vertex = GlobalPoint(0, 0, 0)) const;
+
   std::pair<int, float> photonIsolation(const reco::RecoCandidate* const recocand,
                                         const reco::ElectronCollection* allEle,
-                                        const reco::TrackCollection* isoTracks);
+                                        const reco::TrackCollection* isoTracks) const;
+
+  std::pair<int, float> photonIsolation(float phoEta,
+                                        float phiPhi,
+                                        const reco::ElectronCollection* allEle,
+                                        const reco::TrackCollection* isoTracks) const;
 
   /// Get number of tracks inside an isolation cone for electrons
-  int electronTrackCount(const reco::Track* const tr, const reco::TrackCollection* isoTracks) {
+  int electronTrackCount(const reco::Track* const tr, const reco::TrackCollection* isoTracks) const {
     return electronIsolation(tr, isoTracks).first;
   }
-  int electronTrackCount(const reco::Track* const tr, const reco::TrackCollection* isoTracks, GlobalPoint vertex) {
+  int electronTrackCount(const reco::Track* const tr,
+                         const reco::TrackCollection* isoTracks,
+                         GlobalPoint vertex) const {
     return electronIsolation(tr, isoTracks, vertex).first;
   }
 
@@ -101,72 +111,81 @@ class EgammaHLTTrackIsolation {
   /// set useVertex=false to consider all tracks for isolation
   int photonTrackCount(const reco::RecoCandidate* const recocand,
                        const reco::TrackCollection* isoTracks,
-                       bool useVertex) {
+                       bool useVertex) const {
     return photonIsolation(recocand, isoTracks, useVertex).first;
   }
   int photonTrackCount(const reco::RecoCandidate* const recocand,
                        const reco::TrackCollection* isoTracks,
-                       GlobalPoint vertex) {
+                       GlobalPoint vertex) const {
     return photonIsolation(recocand, isoTracks, vertex).first;
   }
   int photonTrackCount(const reco::RecoCandidate* const recocand,
                        const reco::ElectronCollection* allEle,
-                       const reco::TrackCollection* isoTracks) {
+                       const reco::TrackCollection* isoTracks) const {
     return photonIsolation(recocand, allEle, isoTracks).first;
   }
 
   /// Get Pt sum of tracks inside an isolation cone for electrons
-  float electronPtSum(const reco::Track* const tr, const reco::TrackCollection* isoTracks) {
+  float electronPtSum(const reco::Track* const tr, const reco::TrackCollection* isoTracks) const {
     return electronIsolation(tr, isoTracks).second;
   }
-  float electronPtSum(const reco::Track* const tr, const reco::TrackCollection* isoTracks, GlobalPoint vertex) {
+  float electronPtSum(const reco::Track* const tr, const reco::TrackCollection* isoTracks, GlobalPoint vertex) const {
     return electronIsolation(tr, isoTracks, vertex).second;
   }
   float electronPtSum(const reco::Track* const tr,
                       const reco::ElectronCollection* allEle,
-                      const reco::TrackCollection* isoTracks) {
+                      const reco::TrackCollection* isoTracks) const {
     return electronIsolation(tr, allEle, isoTracks).second;
   }
 
   /// Get Pt sum of tracks inside an isolation cone for photons
   /// set useVertex=true to use Photon vertex from EgammaPhotonVtxFinder
   /// set useVertex=false to consider all tracks for isolation
-  float photonPtSum(const reco::RecoCandidate* const recocand, const reco::TrackCollection* isoTracks, bool useVertex) {
+  float photonPtSum(const reco::RecoCandidate* const recocand,
+                    const reco::TrackCollection* isoTracks,
+                    bool useVertex) const {
     return photonIsolation(recocand, isoTracks, useVertex).second;
   }
   float photonPtSum(const reco::RecoCandidate* const recocand,
                     const reco::TrackCollection* isoTracks,
-                    GlobalPoint vertex) {
+                    GlobalPoint vertex) const {
     return photonIsolation(recocand, isoTracks, vertex).second;
   }
   float photonPtSum(const reco::RecoCandidate* const recocand,
                     const reco::ElectronCollection* allEle,
-                    const reco::TrackCollection* isoTracks) {
+                    const reco::TrackCollection* isoTracks) const {
     return photonIsolation(recocand, allEle, isoTracks).second;
   }
 
   /// Get pt cut for itracks.
-  double getPtMin() { return ptMin; }
+  double getPtMin() const { return ptMin; }
   /// Get isolation cone size.
-  double getConeSize() { return conesize; }
+  double getConeSize() const { return conesize; }
   /// Get maximum ivertex z-coordinate spread.
-  double getZspan() { return zspan; }
+  double getZspan() const { return zspan; }
   /// Get maximum transverse distance of ivertex from beam line.
-  double getRspan() { return rspan; }
+  double getRspan() const { return rspan; }
   /// Get veto cone size
-  double getvetoConesize() { return vetoConesize; }
+  double getvetoConesize() const { return vetoConesize; }
 
 private:
   // Call track reconstruction
   std::pair<int, float> findIsoTracks(GlobalVector mom,
                                       GlobalPoint vtx,
                                       const reco::TrackCollection* isoTracks,
                                       bool isElectron,
-                                      bool useVertex = true);
+                                      bool useVertex = true) const;
   std::pair<int, float> findIsoTracksWithoutEle(GlobalVector mom,
                                                 GlobalPoint vtx,
                                                 const reco::ElectronCollection* allEle,
-                                                const reco::TrackCollection* isoTracks);
+                                                const reco::TrackCollection* isoTracks) const {
+    return findIsoTracksWithoutEle(mom.eta(), mom.phi(), vtx, allEle, isoTracks);
+  }
+  std::pair<int, float> findIsoTracksWithoutEle(float centreEta,
+                                                float centrePhi,
+                                                GlobalPoint vtx,
+                                                const reco::ElectronCollection* allEle,
+                                                const reco::TrackCollection* isoTracks) const;
 
   // Parameters of isolation cone geometry.
   double ptMin;