File indexing completed on 2023-10-25 10:01:05
0001 #include "RecoMET/METPUSubtraction/interface/NoPileUpMEtUtilities.h"
0002
0003 #include "DataFormats/Math/interface/deltaR.h"
0004
0005 #include <algorithm>
0006 #include <cmath>
0007
0008 NoPileUpMEtUtilities::NoPileUpMEtUtilities() {
0009 minPtDef_ = -1;
0010 maxPtDef_ = 1000000;
0011 }
0012
0013 NoPileUpMEtUtilities::~NoPileUpMEtUtilities() {}
0014
0015
0016
0017
0018
0019 void NoPileUpMEtUtilities::finalizeMEtData(CommonMETData& metData) {
0020 metData.met = sqrt(metData.mex * metData.mex + metData.mey * metData.mey);
0021 metData.mez = 0.;
0022 metData.phi = atan2(metData.mey, metData.mex);
0023 }
0024
0025
0026
0027
0028 reco::PUSubMETCandInfoCollection NoPileUpMEtUtilities::cleanJets(
0029 const reco::PUSubMETCandInfoCollection& jets,
0030 const std::vector<reco::Candidate::LorentzVector>& leptons,
0031 double dRoverlap,
0032 bool invert) {
0033 reco::PUSubMETCandInfoCollection retVal;
0034 for (reco::PUSubMETCandInfoCollection::const_iterator jet = jets.begin(); jet != jets.end(); ++jet) {
0035 bool isOverlap = false;
0036 for (std::vector<reco::Candidate::LorentzVector>::const_iterator lepton = leptons.begin(); lepton != leptons.end();
0037 ++lepton) {
0038 if (deltaR2(jet->p4(), *lepton) < dRoverlap * dRoverlap) {
0039 isOverlap = true;
0040 break;
0041 }
0042 }
0043 if ((!isOverlap && !invert) || (isOverlap && invert))
0044 retVal.push_back(*jet);
0045 }
0046 return retVal;
0047 }
0048
0049 CommonMETData NoPileUpMEtUtilities::computeCandidateSum(const reco::PUSubMETCandInfoCollection& cands,
0050 bool neutralFracOnly,
0051 double& sumAbsPx,
0052 double& sumAbsPy) {
0053 CommonMETData retVal;
0054 retVal.mex = 0.;
0055 retVal.mey = 0.;
0056 retVal.sumet = 0.;
0057 double retVal_sumAbsPx = 0.;
0058 double retVal_sumAbsPy = 0.;
0059 double pFrac = 1;
0060 for (reco::PUSubMETCandInfoCollection::const_iterator cand = cands.begin(); cand != cands.end(); ++cand) {
0061 pFrac = 1;
0062 if (neutralFracOnly)
0063 pFrac = (1 - cand->chargedEnFrac());
0064
0065 retVal.mex += cand->p4().px() * pFrac;
0066 retVal.mey += cand->p4().py() * pFrac;
0067 retVal.sumet += cand->p4().pt() * pFrac;
0068 retVal_sumAbsPx += std::abs(cand->p4().px());
0069 retVal_sumAbsPy += std::abs(cand->p4().py());
0070 }
0071 finalizeMEtData(retVal);
0072 sumAbsPx = retVal_sumAbsPx;
0073 sumAbsPy = retVal_sumAbsPy;
0074 return retVal;
0075 }
0076
0077
0078
0079
0080
0081 reco::PUSubMETCandInfoCollection NoPileUpMEtUtilities::cleanPFCandidates(
0082 const reco::PUSubMETCandInfoCollection& pfCandidates,
0083 const std::vector<reco::Candidate::LorentzVector>& leptons,
0084 double dRoverlap,
0085 bool invert) {
0086
0087
0088
0089 reco::PUSubMETCandInfoCollection retVal;
0090 for (reco::PUSubMETCandInfoCollection::const_iterator pfCandidate = pfCandidates.begin();
0091 pfCandidate != pfCandidates.end();
0092 ++pfCandidate) {
0093 bool isOverlap = false;
0094 for (std::vector<reco::Candidate::LorentzVector>::const_iterator lepton = leptons.begin(); lepton != leptons.end();
0095 ++lepton) {
0096 if (deltaR2(pfCandidate->p4(), *lepton) < dRoverlap * dRoverlap) {
0097 isOverlap = true;
0098 break;
0099 }
0100 }
0101 if ((!isOverlap && !invert) || (isOverlap && invert))
0102 retVal.push_back(*pfCandidate);
0103 }
0104 return retVal;
0105 }
0106
0107 reco::PUSubMETCandInfoCollection NoPileUpMEtUtilities::selectCandidates(const reco::PUSubMETCandInfoCollection& cands,
0108 double minPt,
0109 double maxPt,
0110 int type,
0111 bool isCharged,
0112 int isWithinJet) {
0113 reco::PUSubMETCandInfoCollection retVal;
0114 for (reco::PUSubMETCandInfoCollection::const_iterator cand = cands.begin(); cand != cands.end(); ++cand) {
0115 if (isCharged && cand->charge() == 0)
0116 continue;
0117 double jetPt = cand->p4().pt();
0118 if (jetPt < minPt || jetPt > maxPt)
0119 continue;
0120 if (type != reco::PUSubMETCandInfo::kUndefined && cand->type() != type)
0121 continue;
0122
0123
0124 if (isWithinJet != NoPileUpMEtUtilities::kAll && (cand->isWithinJet() != isWithinJet))
0125 continue;
0126
0127 retVal.push_back(*cand);
0128 }
0129 return retVal;
0130 }
0131
0132
0133
0134
0135
0136
0137
0138
0139
0140 void NoPileUpMEtUtilities::computeAllSums(const reco::PUSubMETCandInfoCollection& jets,
0141 const reco::PUSubMETCandInfoCollection& pfCandidates) {
0142 reco::PUSubMETCandInfoCollection pfcsCh = selectCandidates(
0143 pfCandidates, minPtDef_, maxPtDef_, reco::PUSubMETCandInfo::kUndefined, false, NoPileUpMEtUtilities::kAll);
0144
0145 reco::PUSubMETCandInfoCollection pfcsChHS = selectCandidates(
0146 pfCandidates, minPtDef_, maxPtDef_, reco::PUSubMETCandInfo::kChHS, true, NoPileUpMEtUtilities::kAll);
0147
0148 reco::PUSubMETCandInfoCollection pfcsChPU = selectCandidates(
0149 pfCandidates, minPtDef_, maxPtDef_, reco::PUSubMETCandInfo::kChPU, true, NoPileUpMEtUtilities::kAll);
0150
0151 reco::PUSubMETCandInfoCollection pfcsNUnclustered = selectCandidates(
0152 pfCandidates, minPtDef_, maxPtDef_, reco::PUSubMETCandInfo::kNeutral, false, NoPileUpMEtUtilities::kOutsideJet);
0153
0154 reco::PUSubMETCandInfoCollection jetsHS =
0155 selectCandidates(jets, 10.0, maxPtDef_, reco::PUSubMETCandInfo::kHS, false, NoPileUpMEtUtilities::kAll);
0156 reco::PUSubMETCandInfoCollection jetsPU =
0157 selectCandidates(jets, 10.0, maxPtDef_, reco::PUSubMETCandInfo::kPU, false, NoPileUpMEtUtilities::kAll);
0158
0159
0160
0161 chHSPfcSum_ = computeCandidateSum(pfcsChHS, false, chHSPfcSumAbsPx_, chHSPfcSumAbsPy_);
0162 chPUPfcSum_ = computeCandidateSum(pfcsChPU, false, chPUPfcSumAbsPx_, chPUPfcSumAbsPy_);
0163 nUncPfcSum_ = computeCandidateSum(pfcsNUnclustered, false, nUncPfcSumAbsPx_, nUncPfcSumAbsPy_);
0164
0165 nHSJetSum_ = computeCandidateSum(jetsHS, true, nHSJetSumAbsPx_, nHSJetSumAbsPy_);
0166 nPUJetSum_ = computeCandidateSum(jetsPU, true, nPUJetSumAbsPx_, nPUJetSumAbsPy_);
0167 }
0168
0169 CommonMETData NoPileUpMEtUtilities::computeRecoil(int metType, double& sumAbsPx, double& sumAbsPy) {
0170 CommonMETData retVal;
0171 double retSumAbsPx = 0.;
0172 double retSumAbsPy = 0.;
0173
0174
0175
0176
0177
0178
0179
0180
0181
0182
0183 if (metType == NoPileUpMEtUtilities::kChHSMET) {
0184 retVal.mex = chHSPfcSum_.mex;
0185 retVal.mey = chHSPfcSum_.mey;
0186 retVal.sumet = chHSPfcSum_.sumet;
0187 retSumAbsPx = chHSPfcSumAbsPx_;
0188 retSumAbsPy = chHSPfcSumAbsPy_;
0189 }
0190 if (metType == NoPileUpMEtUtilities::kChPUMET) {
0191 retVal.mex = chPUPfcSum_.mex;
0192 retVal.mey = chPUPfcSum_.mey;
0193 retVal.sumet = chPUPfcSum_.sumet;
0194 retSumAbsPx = chPUPfcSumAbsPx_;
0195 retSumAbsPy = chPUPfcSumAbsPy_;
0196 }
0197 if (metType == NoPileUpMEtUtilities::kNeutralUncMET) {
0198 retVal.mex = nUncPfcSum_.mex;
0199 retVal.mey = nUncPfcSum_.mey;
0200 retVal.sumet = nUncPfcSum_.sumet;
0201 retSumAbsPx = nUncPfcSumAbsPx_;
0202 retSumAbsPy = nUncPfcSumAbsPy_;
0203 }
0204 if (metType == NoPileUpMEtUtilities::kHadronicHSMET) {
0205 retVal.mex = chHSPfcSum_.mex + nHSJetSum_.mex;
0206 retVal.mey = chHSPfcSum_.mey + nHSJetSum_.mey;
0207 retVal.sumet = chHSPfcSum_.sumet + nHSJetSum_.sumet;
0208 retSumAbsPx = chHSPfcSumAbsPx_ + nHSJetSumAbsPx_;
0209 retSumAbsPy = chHSPfcSumAbsPy_ + nHSJetSumAbsPy_;
0210 }
0211 if (metType == NoPileUpMEtUtilities::kHadronicPUMET) {
0212 retVal.mex = chPUPfcSum_.mex + nHSJetSum_.mex;
0213 retVal.mey = chPUPfcSum_.mey + nHSJetSum_.mey;
0214 retVal.sumet = chPUPfcSum_.sumet + nHSJetSum_.sumet;
0215 retSumAbsPx = chPUPfcSumAbsPx_ + nHSJetSumAbsPx_;
0216 retSumAbsPy = chPUPfcSumAbsPy_ + nHSJetSumAbsPy_;
0217 }
0218
0219 sumAbsPx = retSumAbsPx;
0220 sumAbsPy = retSumAbsPy;
0221
0222 return retVal;
0223 }