File indexing completed on 2023-10-25 09:52:58
0001
0002
0003
0004
0005
0006
0007
0008 #include "HLTPMMassFilter.h"
0009
0010 #include <cstdlib>
0011 #include <cmath>
0012
0013
0014
0015
0016 HLTPMMassFilter::HLTPMMassFilter(const edm::ParameterSet& iConfig) : HLTFilter(iConfig), magFieldToken_(esConsumes()) {
0017 candTag_ = iConfig.getParameter<edm::InputTag>("candTag");
0018 beamSpot_ = iConfig.getParameter<edm::InputTag>("beamSpot");
0019 l1EGTag_ = iConfig.getParameter<edm::InputTag>("l1EGCand");
0020
0021 lowerMassCut_ = iConfig.getParameter<double>("lowerMassCut");
0022 upperMassCut_ = iConfig.getParameter<double>("upperMassCut");
0023 nZcandcut_ = iConfig.getParameter<int>("nZcandcut");
0024 reqOppCharge_ = iConfig.getUntrackedParameter<bool>("reqOppCharge", false);
0025 isElectron1_ = iConfig.getUntrackedParameter<bool>("isElectron1", true);
0026 isElectron2_ = iConfig.getUntrackedParameter<bool>("isElectron2", true);
0027
0028 candToken_ = consumes<trigger::TriggerFilterObjectWithRefs>(candTag_);
0029 beamSpotToken_ = consumes<reco::BeamSpot>(beamSpot_);
0030 }
0031
0032 HLTPMMassFilter::~HLTPMMassFilter() = default;
0033
0034 void HLTPMMassFilter::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0035 edm::ParameterSetDescription desc;
0036 makeHLTFilterDescription(desc);
0037 desc.add<edm::InputTag>("candTag", edm::InputTag("hltL1NonIsoDoublePhotonEt5UpsHcalIsolFilter"));
0038 desc.add<edm::InputTag>("beamSpot", edm::InputTag("hltOfflineBeamSpot"));
0039 desc.add<double>("lowerMassCut", 8.0);
0040 desc.add<double>("upperMassCut", 11.0);
0041 desc.add<int>("nZcandcut", 1);
0042 desc.addUntracked<bool>("reqOppCharge", true);
0043 desc.addUntracked<bool>("isElectron1", false);
0044 desc.addUntracked<bool>("isElectron2", false);
0045 desc.add<edm::InputTag>("l1EGCand", edm::InputTag("hltL1IsoRecoEcalCandidate"));
0046 descriptions.add("hltPMMassFilter", desc);
0047 }
0048
0049
0050 bool HLTPMMassFilter::hltFilter(edm::Event& iEvent,
0051 const edm::EventSetup& iSetup,
0052 trigger::TriggerFilterObjectWithRefs& filterproduct) const {
0053 using namespace std;
0054 using namespace edm;
0055 using namespace reco;
0056 using namespace trigger;
0057
0058
0059 if (saveTags()) {
0060 filterproduct.addCollectionTag(l1EGTag_);
0061 }
0062
0063 auto const& theMagField = iSetup.getData(magFieldToken_);
0064
0065 edm::Handle<trigger::TriggerFilterObjectWithRefs> PrevFilterOutput;
0066 iEvent.getByToken(candToken_, PrevFilterOutput);
0067
0068
0069 edm::Handle<reco::BeamSpot> recoBeamSpotHandle;
0070 iEvent.getByToken(beamSpotToken_, recoBeamSpotHandle);
0071
0072 const GlobalPoint vertexPos(
0073 recoBeamSpotHandle->position().x(), recoBeamSpotHandle->position().y(), recoBeamSpotHandle->position().z());
0074
0075 int n = 0;
0076
0077
0078
0079
0080
0081
0082 std::vector<TLorentzVector> pEleCh1;
0083 std::vector<TLorentzVector> pEleCh2;
0084 std::vector<double> charge;
0085
0086 if (isElectron1_ && isElectron2_) {
0087 Ref<ElectronCollection> refele;
0088
0089 vector<Ref<ElectronCollection> > electrons;
0090 PrevFilterOutput->getObjects(TriggerElectron, electrons);
0091
0092 for (auto& electron : electrons) {
0093 refele = electron;
0094
0095 TLorentzVector pThisEle(refele->px(), refele->py(), refele->pz(), refele->energy());
0096 pEleCh1.push_back(pThisEle);
0097 charge.push_back(refele->charge());
0098 }
0099
0100 std::vector<bool> save_cand(electrons.size(), true);
0101 for (unsigned int jj = 0; jj < electrons.size(); jj++) {
0102 for (unsigned int ii = jj + 1; ii < electrons.size(); ii++) {
0103 if (reqOppCharge_ && charge[jj] * charge[ii] > 0)
0104 continue;
0105 if (isGoodPair(pEleCh1[jj], pEleCh1[ii])) {
0106 n++;
0107 for (auto const idx : {jj, ii}) {
0108 if (save_cand[idx])
0109 filterproduct.addObject(TriggerElectron, electrons[idx]);
0110 save_cand[idx] = false;
0111 }
0112 }
0113 }
0114 }
0115
0116 } else {
0117 Ref<RecoEcalCandidateCollection> refsc;
0118
0119 vector<Ref<RecoEcalCandidateCollection> > scs;
0120 PrevFilterOutput->getObjects(TriggerCluster, scs);
0121 if (scs.empty())
0122 PrevFilterOutput->getObjects(TriggerPhoton, scs);
0123
0124 for (auto& i : scs) {
0125 refsc = i;
0126 const reco::SuperClusterRef sc = refsc->superCluster();
0127 TLorentzVector pscPos = approxMomAtVtx(theMagField, vertexPos, sc, 1);
0128 pEleCh1.push_back(pscPos);
0129
0130 TLorentzVector pscEle = approxMomAtVtx(theMagField, vertexPos, sc, -1);
0131 pEleCh2.push_back(pscEle);
0132 }
0133
0134 std::vector<bool> save_cand(scs.size(), true);
0135 for (unsigned int jj = 0; jj < scs.size(); jj++) {
0136 for (unsigned int ii = jj + 1; ii < scs.size(); ii++) {
0137 if (isGoodPair(pEleCh1[jj], pEleCh2[ii]) or isGoodPair(pEleCh2[jj], pEleCh1[ii])) {
0138 n++;
0139 for (auto const idx : {jj, ii}) {
0140 if (save_cand[idx])
0141 filterproduct.addObject(TriggerCluster, scs[idx]);
0142 save_cand[idx] = false;
0143 }
0144 }
0145 }
0146 }
0147 }
0148
0149
0150 bool accept(n >= nZcandcut_);
0151
0152 return accept;
0153 }
0154
0155 bool HLTPMMassFilter::isGoodPair(TLorentzVector const& v1, TLorentzVector const& v2) const {
0156 if (std::abs(v1.E() - v2.E()) < 0.00001)
0157 return false;
0158
0159 auto const mass = (v1 + v2).M();
0160 return (mass >= lowerMassCut_ and mass <= upperMassCut_);
0161 }
0162
0163 TLorentzVector HLTPMMassFilter::approxMomAtVtx(const MagneticField& magField,
0164 const GlobalPoint& xvert,
0165 const reco::SuperClusterRef sc,
0166 int charge) const {
0167 GlobalPoint xsc(sc->position().x(), sc->position().y(), sc->position().z());
0168 float energy = sc->energy();
0169 auto theFTS = trackingTools::ftsFromVertexToPoint(magField, xsc, xvert, energy, charge);
0170 float theApproxMomMod = theFTS.momentum().x() * theFTS.momentum().x() +
0171 theFTS.momentum().y() * theFTS.momentum().y() + theFTS.momentum().z() * theFTS.momentum().z();
0172 TLorentzVector theApproxMom(
0173 theFTS.momentum().x(), theFTS.momentum().y(), theFTS.momentum().z(), sqrt(theApproxMomMod + 2.61121E-7));
0174 return theApproxMom;
0175 }
0176
0177
0178 #include "FWCore/Framework/interface/MakerMacros.h"
0179 DEFINE_FWK_MODULE(HLTPMMassFilter);