File indexing completed on 2024-09-08 23:51:46
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 #include "DataFormats/Common/interface/Handle.h"
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/Framework/interface/Frameworkfwd.h"
0023 #include "FWCore/Framework/interface/global/EDFilter.h"
0024 #include "FWCore/Framework/interface/MakerMacros.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026 #include "FWCore/Utilities/interface/EDGetToken.h"
0027 #include "FWCore/Utilities/interface/InputTag.h"
0028 #include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
0029
0030 #include "TLorentzVector.h"
0031
0032 #include <algorithm>
0033 #include <cmath>
0034 #include <cstdlib>
0035
0036 class ZgammaMassFilter : public edm::global::EDFilter<> {
0037 public:
0038 explicit ZgammaMassFilter(const edm::ParameterSet&);
0039
0040 bool filter(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0041
0042 private:
0043 const edm::EDGetTokenT<edm::HepMCProduct> token_;
0044
0045 const double minPhotonPt;
0046 const double minLeptonPt;
0047
0048 const double minPhotonEta;
0049 const double minLeptonEta;
0050
0051 const double maxPhotonEta;
0052 const double maxLeptonEta;
0053
0054 const double minDileptonMass;
0055 const double minZgMass;
0056 };
0057
0058 namespace {
0059
0060 class orderByPt {
0061 public:
0062 bool operator()(TLorentzVector const& a, TLorentzVector const& b) {
0063 if (a.Pt() == b.Pt()) {
0064 return a.Pt() < b.Pt();
0065 } else {
0066 return a.Pt() > b.Pt();
0067 }
0068 }
0069 };
0070 }
0071
0072 using namespace edm;
0073 using namespace std;
0074
0075 ZgammaMassFilter::ZgammaMassFilter(const edm::ParameterSet& iConfig)
0076 : token_(consumes<edm::HepMCProduct>(iConfig.getParameter<InputTag>("HepMCProduct"))),
0077 minPhotonPt(iConfig.getParameter<double>("minPhotonPt")),
0078 minLeptonPt(iConfig.getParameter<double>("minLeptonPt")),
0079 minPhotonEta(iConfig.getParameter<double>("minPhotonEta")),
0080 minLeptonEta(iConfig.getParameter<double>("minLeptonEta")),
0081 maxPhotonEta(iConfig.getParameter<double>("maxPhotonEta")),
0082 maxLeptonEta(iConfig.getParameter<double>("maxLeptonEta")),
0083 minDileptonMass(iConfig.getParameter<double>("minDileptonMass")),
0084 minZgMass(iConfig.getParameter<double>("minZgMass")) {}
0085
0086 bool ZgammaMassFilter::filter(edm::StreamID, edm::Event& iEvent, const edm::EventSetup&) const {
0087 bool accepted = false;
0088 Handle<HepMCProduct> evt;
0089 iEvent.getByToken(token_, evt);
0090 const HepMC::GenEvent* myGenEvent = evt->GetEvent();
0091
0092 vector<TLorentzVector> Lepton;
0093 Lepton.clear();
0094 vector<TLorentzVector> Photon;
0095 Photon.clear();
0096 vector<float> Charge;
0097 Charge.clear();
0098
0099 for (HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin(); p != myGenEvent->particles_end();
0100 ++p) {
0101 if ((*p)->status() == 1 && (abs((*p)->pdg_id()) == 11 || abs((*p)->pdg_id()) == 13 || abs((*p)->pdg_id()) == 15)) {
0102 TLorentzVector LeptP((*p)->momentum().px(), (*p)->momentum().py(), (*p)->momentum().pz(), (*p)->momentum().e());
0103 if (LeptP.Pt() > minLeptonPt) {
0104 Lepton.push_back(LeptP);
0105 }
0106 }
0107
0108 if (abs((*p)->pdg_id()) == 22 && (*p)->status() == 1) {
0109 TLorentzVector PhotP((*p)->momentum().px(), (*p)->momentum().py(), (*p)->momentum().pz(), (*p)->momentum().e());
0110 if (PhotP.Pt() > minPhotonPt) {
0111 Photon.push_back(PhotP);
0112 }
0113 }
0114
0115 }
0116
0117
0118
0119
0120
0121
0122
0123
0124
0125
0126
0127 std::stable_sort(Photon.begin(), Photon.end(), orderByPt());
0128 std::stable_sort(Lepton.begin(), Lepton.end(), orderByPt());
0129
0130
0131
0132
0133
0134
0135
0136
0137
0138
0139
0140
0141 if (!Photon.empty() && Lepton.size() > 1 && Photon[0].Pt() > minPhotonPt && Lepton[0].Pt() > minLeptonPt &&
0142 Lepton[1].Pt() > minLeptonPt && Photon[0].Eta() > minPhotonEta && Lepton[0].Eta() > minLeptonEta &&
0143 Lepton[1].Eta() > minLeptonEta && Photon[0].Eta() < maxPhotonEta && Lepton[0].Eta() < maxLeptonEta &&
0144 Lepton[1].Eta() < maxLeptonEta && (Lepton[0] + Lepton[1]).M() > minDileptonMass &&
0145 (Lepton[0] + Lepton[1] + Photon[0]).M() >
0146 minZgMass) {
0147 accepted = true;
0148 }
0149
0150
0151
0152 return accepted;
0153 }
0154
0155 DEFINE_FWK_MODULE(ZgammaMassFilter);