File indexing completed on 2024-09-07 04:37:32
0001 #include "RecoJets/JetAlgorithms/interface/CATopJetHelper.h"
0002
0003 struct GreaterByPtCandPtr {
0004 bool operator()(const edm::Ptr<reco::Candidate>& t1, const edm::Ptr<reco::Candidate>& t2) const {
0005 return t1->pt() > t2->pt();
0006 }
0007 };
0008
0009 reco::CATopJetProperties CATopJetHelper::operator()(reco::Jet const& ihardJet) const {
0010 reco::CATopJetProperties properties;
0011
0012 reco::Jet::Constituents subjets = ihardJet.getJetConstituents();
0013 properties.nSubJets = subjets.size();
0014 properties.topMass = ihardJet.mass();
0015 properties.wMass = 99999.;
0016 properties.minMass = 999999.;
0017
0018
0019 if (properties.nSubJets >= 3) {
0020
0021 sort(subjets.begin(), subjets.end(), GreaterByPtCandPtr());
0022
0023
0024 for (int isub = 0; isub < 2; ++isub) {
0025
0026 reco::Jet::Constituent icandJet = subjets[isub];
0027
0028
0029
0030 for (int jsub = isub + 1; jsub < 3; ++jsub) {
0031
0032 reco::Jet::Constituent jcandJet = subjets[jsub];
0033
0034 reco::Candidate::LorentzVector wCand = icandJet->p4() + jcandJet->p4();
0035
0036
0037 double imw = wCand.mass();
0038
0039
0040 if (fabs(imw - WMass_) < fabs(properties.wMass - WMass_)) {
0041 properties.wMass = imw;
0042 }
0043
0044 if (fabs(imw) < properties.minMass) {
0045 properties.minMass = imw;
0046 }
0047 }
0048 }
0049 }
0050
0051 if (properties.minMass == 999999) {
0052 properties.minMass = -1;
0053 }
0054
0055 return properties;
0056 }