SinglePhotonJetPlusHOFilter

Line Code
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163
// Package:    SinglePhotonJetPlusHOFilter
// Class:      SinglePhotonJetPlusHOFilter
//
/*
 Description: [one line class summary]
Skimming of SinglePhoton data set for the study of HO absolute weight calculation
* Skimming Efficiency : ~ 2 %
 Implementation:
     [Notes on implementation]
     For Secondary Datasets (SD)
*/
//
// Original Author:  Gobinda Majumder & Suman Chatterjee
//         Created:  Fri July 29 14:52:17 IST 2016
// $Id$
//
//

// system include files
#include <memory>

// class declaration
//

#include "FWCore/Framework/interface/global/EDFilter.h"

#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"

#include "DataFormats/JetReco/interface/PFJet.h"
#include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
#include "DataFormats/ParticleFlowReco/interface/PFCluster.h"

#include "DataFormats/EgammaCandidates/interface/Photon.h"
#include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
#include "DataFormats/Math/interface/deltaR.h"

class SinglePhotonJetPlusHOFilter : public edm::global::EDFilter<> {
public:
  explicit SinglePhotonJetPlusHOFilter(const edm::ParameterSet&);

  static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);

private:
  bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;

  // ----------member data ---------------------------
  const double jtptthr_;
  const double jtetath_;
  const double hothres_;
  const double pho_Ptcut_;

  const edm::EDGetTokenT<reco::PFJetCollection> tok_PFJets_;
  const edm::EDGetTokenT<reco::PFClusterCollection> tok_hoht_;
  const edm::EDGetTokenT<edm::View<reco::Photon> > tok_photons_;
};

SinglePhotonJetPlusHOFilter::SinglePhotonJetPlusHOFilter(const edm::ParameterSet& iConfig)
    : jtptthr_{iConfig.getParameter<double>("Ptcut")},
      jtetath_{iConfig.getParameter<double>("Etacut")},
      hothres_{iConfig.getParameter<double>("HOcut")},
      pho_Ptcut_{iConfig.getParameter<double>("Pho_Ptcut")},
      tok_PFJets_{consumes<reco::PFJetCollection>(iConfig.getParameter<edm::InputTag>("PFJets"))},
      tok_hoht_{consumes<reco::PFClusterCollection>(iConfig.getParameter<edm::InputTag>("particleFlowClusterHO"))},
      tok_photons_{consumes<edm::View<reco::Photon> >(iConfig.getParameter<edm::InputTag>("Photons"))} {}

// ------------ method called on each new Event  ------------
bool SinglePhotonJetPlusHOFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
  using namespace std;
  using namespace edm;
  using namespace reco;

  bool passed = false;
  vector<pair<double, double> > jetdirection;
  vector<double> jetspt;
  if (auto PFJets = iEvent.getHandle(tok_PFJets_)) {
    for (const auto& jet : *PFJets) {
      if ((jet.pt() < jtptthr_) || (abs(jet.eta()) > jtetath_))
        continue;

      jetdirection.emplace_back(jet.eta(), jet.phi());
      jetspt.push_back(jet.pt());
      passed = true;
    }
  }
  if (!passed)
    return passed;

  bool pho_passed = false;
  vector<pair<double, double> > phodirection;
  vector<double> phopT;
  if (auto photons = iEvent.getHandle(tok_photons_)) {
    for (const auto& gamma1 : *photons) {
      if (gamma1.pt() < pho_Ptcut_)
        continue;
      phodirection.emplace_back(gamma1.eta(), gamma1.phi());
      phopT.push_back(gamma1.pt());

      pho_passed = true;
    }
  }

  if (!pho_passed)
    return false;

  bool isJetDir = false;
  if (auto hhoht = iEvent.getHandle(tok_hoht_)) {
    const auto& hoht = *hhoht;
    if (!hoht.empty()) {
      for (const auto& jet : jetdirection) {
        bool matched = false;
        for (const auto& ph : phodirection) {
          if (abs(deltaPhi(ph.second, jet.second)) > 2.0) {
            matched = true;
            break;
          }
        }
        if (matched) {
          for (const auto& ij : hoht) {
            double hoenr = ij.energy();
            if (hoenr < hothres_)
              continue;

            const math::XYZPoint& cluster_pos = ij.position();

            double hoeta = cluster_pos.eta();
            double hophi = cluster_pos.phi();

            double delta = deltaR2(jet.first, jet.second, hoeta, hophi);
            if (delta < 0.5) {
              isJetDir = true;
              break;
            }
          }
        }
        if (isJetDir) {
          break;
        }
      }
    }
  }

  return isJetDir;
}

// ------------ method fills 'descriptions' with the allowed parameters for the module  ------------
void SinglePhotonJetPlusHOFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
  //The following says we do not know what parameters are allowed so do no validation
  // Please change this to state exactly what you do use, even if it is no parameters
  edm::ParameterSetDescription desc;
  desc.add<double>("Ptcut", 90.0);
  desc.add<double>("Etacut", 1.5);
  desc.add<double>("HOcut", 5);
  desc.add<double>("Pho_Ptcut", 120);
  desc.add<edm::InputTag>("PFJets");
  desc.add<edm::InputTag>("particleFlowClusterHO");
  desc.add<edm::InputTag>("Photons");
  descriptions.addWithDefaultLabel(desc);
}

//define this as a plug-in
DEFINE_FWK_MODULE(SinglePhotonJetPlusHOFilter);