CMS 3D CMS Logo

TtSemiLepJetCombWMassMaxSumPt.cc
Go to the documentation of this file.
2 
4 
6  : jetsToken_(consumes<std::vector<pat::Jet> >(cfg.getParameter<edm::InputTag>("jets"))),
7  lepsToken_(consumes<edm::View<reco::RecoCandidate> >(cfg.getParameter<edm::InputTag>("leps"))),
8  maxNJets_(cfg.getParameter<int>("maxNJets")),
9  wMass_(cfg.getParameter<double>("wMass")),
10  useBTagging_(cfg.getParameter<bool>("useBTagging")),
11  bTagAlgorithm_(cfg.getParameter<std::string>("bTagAlgorithm")),
12  minBDiscBJets_(cfg.getParameter<double>("minBDiscBJets")),
13  maxBDiscLightJets_(cfg.getParameter<double>("maxBDiscLightJets")) {
14  if (maxNJets_ < 4 && maxNJets_ != -1)
15  throw cms::Exception("WrongConfig") << "Parameter maxNJets can not be set to " << maxNJets_ << ". \n"
16  << "It has to be larger than 4 or can be set to -1 to take all jets.";
17 
18  produces<std::vector<std::vector<int> > >();
19  produces<int>("NumberOfConsideredJets");
20 }
21 
23 
25  std::unique_ptr<std::vector<std::vector<int> > > pOut(new std::vector<std::vector<int> >);
26  std::unique_ptr<int> pJetsConsidered(new int);
27 
28  std::vector<int> match;
29  for (unsigned int i = 0; i < 4; ++i)
30  match.push_back(-1);
31 
32  // get jets
35 
36  // get leptons
39 
40  // skip events without lepton candidate or less than 4 jets or no MET
41  if (leps->empty() || jets->size() < 4) {
42  pOut->push_back(match);
43  evt.put(std::move(pOut));
44  *pJetsConsidered = jets->size();
45  evt.put(std::move(pJetsConsidered), "NumberOfConsideredJets");
46  return;
47  }
48 
49  unsigned maxNJets = maxNJets_;
50  if (maxNJets_ == -1 || (int)jets->size() < maxNJets_)
51  maxNJets = jets->size();
52  *pJetsConsidered = maxNJets;
53  evt.put(std::move(pJetsConsidered), "NumberOfConsideredJets");
54 
55  std::vector<bool> isBJet;
56  std::vector<bool> isLJet;
57  int cntBJets = 0;
58  if (useBTagging_) {
59  for (unsigned int idx = 0; idx < maxNJets; ++idx) {
60  isBJet.push_back(((*jets)[idx].bDiscriminator(bTagAlgorithm_) > minBDiscBJets_));
61  isLJet.push_back(((*jets)[idx].bDiscriminator(bTagAlgorithm_) < maxBDiscLightJets_));
62  if ((*jets)[idx].bDiscriminator(bTagAlgorithm_) > minBDiscBJets_)
63  cntBJets++;
64  }
65  }
66 
67  // -----------------------------------------------------
68  // associate those jets that get closest to the W mass
69  // with their invariant mass to the hadronic W boson
70  // -----------------------------------------------------
71  double wDist = -1.;
72  std::vector<int> closestToWMassIndices;
73  closestToWMassIndices.push_back(-1);
74  closestToWMassIndices.push_back(-1);
75  for (unsigned idx = 0; idx < maxNJets; ++idx) {
76  if (useBTagging_ && (!isLJet[idx] || (cntBJets <= 2 && isBJet[idx])))
77  continue;
78  for (unsigned jdx = (idx + 1); jdx < maxNJets; ++jdx) {
79  if (useBTagging_ &&
80  (!isLJet[jdx] || (cntBJets <= 2 && isBJet[jdx]) || (cntBJets == 3 && isBJet[idx] && isBJet[jdx])))
81  continue;
82  reco::Particle::LorentzVector sum = (*jets)[idx].p4() + (*jets)[jdx].p4();
83  if (wDist < 0. || wDist > fabs(sum.mass() - wMass_)) {
84  wDist = fabs(sum.mass() - wMass_);
85  closestToWMassIndices.clear();
86  closestToWMassIndices.push_back(idx);
87  closestToWMassIndices.push_back(jdx);
88  }
89  }
90  }
91 
92  // -----------------------------------------------------
93  // associate those jets with maximum pt of the vectorial
94  // sum to the hadronic decay chain
95  // -----------------------------------------------------
96  double maxPt = -1.;
97  int hadB = -1;
98  if (isValid(closestToWMassIndices[0], jets) && isValid(closestToWMassIndices[1], jets)) {
99  for (unsigned idx = 0; idx < maxNJets; ++idx) {
100  if (useBTagging_ && !isBJet[idx])
101  continue;
102  // make sure it's not used up already from the hadronic W
103  if ((int)idx != closestToWMassIndices[0] && (int)idx != closestToWMassIndices[1]) {
105  (*jets)[closestToWMassIndices[0]].p4() + (*jets)[closestToWMassIndices[1]].p4() + (*jets)[idx].p4();
106  if (maxPt < 0. || maxPt < sum.pt()) {
107  maxPt = sum.pt();
108  hadB = idx;
109  }
110  }
111  }
112  }
113 
114  // -----------------------------------------------------
115  // associate the remaining jet with maximum pt of the
116  // vectorial sum with the leading lepton with the
117  // leptonic b quark
118  // -----------------------------------------------------
119  maxPt = -1.;
120  int lepB = -1;
121  for (unsigned idx = 0; idx < maxNJets; ++idx) {
122  if (useBTagging_ && !isBJet[idx])
123  continue;
124  // make sure it's not used up already from the hadronic decay chain
125  if ((int)idx != closestToWMassIndices[0] && (int)idx != closestToWMassIndices[1] && (int)idx != hadB) {
126  reco::Particle::LorentzVector sum = (*jets)[idx].p4() + (*leps)[0].p4();
127  if (maxPt < 0. || maxPt < sum.pt()) {
128  maxPt = sum.pt();
129  lepB = idx;
130  }
131  }
132  }
133 
134  match[TtSemiLepEvtPartons::LightQ] = closestToWMassIndices[0];
135  match[TtSemiLepEvtPartons::LightQBar] = closestToWMassIndices[1];
138 
139  pOut->push_back(match);
140  evt.put(std::move(pOut));
141 }
electrons_cff.bool
bool
Definition: electrons_cff.py:372
mps_fire.i
i
Definition: mps_fire.py:355
TtSemiLepEvtPartons::LightQBar
Definition: TtSemiLepEvtPartons.h:25
TtSemiLepJetCombWMassMaxSumPt::TtSemiLepJetCombWMassMaxSumPt
TtSemiLepJetCombWMassMaxSumPt(const edm::ParameterSet &)
Definition: TtSemiLepJetCombWMassMaxSumPt.cc:5
TtSemiLepJetCombWMassMaxSumPt::isValid
bool isValid(const int &idx, const edm::Handle< std::vector< pat::Jet > > &jets)
Definition: TtSemiLepJetCombWMassMaxSumPt.h:20
sistrip::View
View
Definition: ConstantsForView.h:26
TtSemiLepJetCombWMassMaxSumPt::wMass_
double wMass_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:27
TtSemiLepEvtPartons::LepB
Definition: TtSemiLepEvtPartons.h:25
TtSemiLepJetCombWMassMaxSumPt::minBDiscBJets_
double minBDiscBJets_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:30
edm
HLT enums.
Definition: AlignableModifier.h:19
TtSemiLepJetCombWMassMaxSumPt.h
TtSemiLepJetCombWMassMaxSumPt::produce
void produce(edm::Event &evt, const edm::EventSetup &setup) override
Definition: TtSemiLepJetCombWMassMaxSumPt.cc:24
singleTopDQM_cfi.jets
jets
Definition: singleTopDQM_cfi.py:42
TtSemiLepJetCombWMassMaxSumPt::~TtSemiLepJetCombWMassMaxSumPt
~TtSemiLepJetCombWMassMaxSumPt() override
Definition: TtSemiLepJetCombWMassMaxSumPt.cc:22
reco
fixed size matrix
Definition: AlignmentAlgorithmBase.h:45
edm::Handle
Definition: AssociativeIterator.h:50
training_settings.idx
idx
Definition: training_settings.py:16
singleTopDQM_cfi.setup
setup
Definition: singleTopDQM_cfi.py:37
reco::Particle::LorentzVector
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:21
TtSemiLepJetCombWMassMaxSumPt::bTagAlgorithm_
std::string bTagAlgorithm_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:29
TtSemiLepHitFitProducer_Electrons_cfi.leps
leps
Definition: TtSemiLepHitFitProducer_Electrons_cfi.py:5
MuonErrorMatrixAnalyzer_cfi.maxPt
maxPt
Definition: MuonErrorMatrixAnalyzer_cfi.py:19
Jet
Definition: Jet.py:1
edm::Event::getByToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:528
TtSemiLepJetCombWMassMaxSumPt::useBTagging_
bool useBTagging_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:28
AlCaHLTBitMon_QueryRunRegistry.string
string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
TtSemiLepJetCombWMassMaxSumPt::jetsToken_
edm::EDGetTokenT< std::vector< pat::Jet > > jetsToken_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:22
HLT_2018_cff.InputTag
InputTag
Definition: HLT_2018_cff.py:79016
edm::ParameterSet
Definition: ParameterSet.h:36
TtSemiLepEvtPartons::HadB
Definition: TtSemiLepEvtPartons.h:25
match
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
TtSemiLepEvtPartons::LightQ
Definition: TtSemiLepEvtPartons.h:25
createfilelist.int
int
Definition: createfilelist.py:10
SUSYDQMAnalyzer_cfi.maxNJets
maxNJets
Definition: SUSYDQMAnalyzer_cfi.py:13
edm::Event::put
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:132
edm::EventSetup
Definition: EventSetup.h:57
pat
Definition: HeavyIon.h:7
looper.cfg
cfg
Definition: looper.py:297
TtSemiLepJetCombWMassMaxSumPt::maxNJets_
int maxNJets_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:26
eostools.move
def move(src, dest)
Definition: eostools.py:511
std
Definition: JetResolutionObject.h:76
TtSemiLepJetCombWMassMaxSumPt::maxBDiscLightJets_
double maxBDiscLightJets_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:31
Exception
Definition: hltDiff.cc:246
TtSemiLepJetCombWMassMaxSumPt::lepsToken_
edm::EDGetTokenT< edm::View< reco::RecoCandidate > > lepsToken_
Definition: TtSemiLepJetCombWMassMaxSumPt.h:25
ParameterSet.h
edm::Event
Definition: Event.h:73