CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
List of all members | Public Member Functions | Private Member Functions | Private Attributes
TtSemiLepJetCombWMassDeltaTopMass Class Reference
Inheritance diagram for TtSemiLepJetCombWMassDeltaTopMass:
edm::global::EDProducer<> edm::global::EDProducerBase edm::ProducerBase edm::EDConsumerBase edm::ProductRegistryHelper

Public Member Functions

 TtSemiLepJetCombWMassDeltaTopMass (const edm::ParameterSet &)
 
- Public Member Functions inherited from edm::global::EDProducer<>
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
EDProduceroperator= (const EDProducer &)=delete
 
bool wantsGlobalLuminosityBlocks () const final
 
bool wantsGlobalRuns () const final
 
bool wantsInputProcessBlocks () const final
 
bool wantsProcessBlocks () const final
 
bool wantsStreamLuminosityBlocks () const final
 
bool wantsStreamRuns () const final
 
- Public Member Functions inherited from edm::global::EDProducerBase
 EDProducerBase ()
 
ModuleDescription const & moduleDescription () const
 
 ~EDProducerBase () override
 
- Public Member Functions inherited from edm::ProducerBase
void callWhenNewProductsRegistered (std::function< void(BranchDescription const &)> const &func)
 
std::vector
< edm::ProductResolverIndex >
const & 
indiciesForPutProducts (BranchType iBranchType) const
 
 ProducerBase ()
 
std::vector
< edm::ProductResolverIndex >
const & 
putTokenIndexToProductResolverIndex () const
 
void registerProducts (ProducerBase *, ProductRegistry *, ModuleDescription const &)
 
std::function< void(BranchDescription
const &)> 
registrationCallback () const
 used by the fwk to register list of products More...
 
void resolvePutIndicies (BranchType iBranchType, ModuleToResolverIndicies const &iIndicies, std::string const &moduleLabel)
 
 ~ProducerBase () noexcept(false) override
 
- Public Member Functions inherited from edm::EDConsumerBase
std::vector< ConsumesInfoconsumesInfo () const
 
void convertCurrentProcessAlias (std::string const &processName)
 Convert "@currentProcess" in InputTag process names to the actual current process name. More...
 
 EDConsumerBase ()
 
 EDConsumerBase (EDConsumerBase const &)=delete
 
 EDConsumerBase (EDConsumerBase &&)=default
 
ESProxyIndex const * esGetTokenIndices (edm::Transition iTrans) const
 
std::vector< ESProxyIndex > const & esGetTokenIndicesVector (edm::Transition iTrans) const
 
std::vector< ESRecordIndex >
const & 
esGetTokenRecordIndicesVector (edm::Transition iTrans) const
 
ProductResolverIndexAndSkipBit indexFrom (EDGetToken, BranchType, TypeID const &) const
 
void itemsMayGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
void itemsToGet (BranchType, std::vector< ProductResolverIndexAndSkipBit > &) const
 
std::vector
< ProductResolverIndexAndSkipBit >
const & 
itemsToGetFrom (BranchType iType) const
 
void labelsForToken (EDGetToken iToken, Labels &oLabels) const
 
void modulesWhoseProductsAreConsumed (std::array< std::vector< ModuleDescription const * > *, NumBranchTypes > &modulesAll, std::vector< ModuleProcessName > &modulesInPreviousProcesses, ProductRegistry const &preg, std::map< std::string, ModuleDescription const * > const &labelsToDesc, std::string const &processName) const
 
EDConsumerBase const & operator= (EDConsumerBase const &)=delete
 
EDConsumerBaseoperator= (EDConsumerBase &&)=default
 
bool registeredToConsume (ProductResolverIndex, bool, BranchType) const
 
bool registeredToConsumeMany (TypeID const &, BranchType) const
 
void selectInputProcessBlocks (ProductRegistry const &productRegistry, ProcessBlockHelperBase const &processBlockHelperBase)
 
ProductResolverIndexAndSkipBit uncheckedIndexFrom (EDGetToken) const
 
void updateLookup (BranchType iBranchType, ProductResolverIndexHelper const &, bool iPrefetchMayGet)
 
void updateLookup (eventsetup::ESRecordsToProxyIndices const &)
 
virtual ~EDConsumerBase () noexcept(false)
 

Private Member Functions

bool isValid (const int &idx, const std::vector< pat::Jet > &jets) const
 
void produce (edm::StreamID, edm::Event &evt, const edm::EventSetup &setup) const override
 

Private Attributes

std::string bTagAlgorithm_
 
edm::EDGetTokenT< std::vector
< pat::Jet > > 
jetsToken_
 
edm::EDGetTokenT< edm::View
< reco::RecoCandidate > > 
lepsToken_
 
double maxBDiscLightJets_
 
int maxNJets_
 
edm::EDGetTokenT< std::vector
< pat::MET > > 
metsToken_
 
double minBDiscBJets_
 
int neutrinoSolutionType_
 
bool useBTagging_
 
double wMass_
 

Additional Inherited Members

- Public Types inherited from edm::global::EDProducerBase
typedef EDProducerBase ModuleType
 
- Public Types inherited from edm::ProducerBase
using ModuleToResolverIndicies = std::unordered_multimap< std::string, std::tuple< edm::TypeID const *, const char *, edm::ProductResolverIndex >>
 
typedef
ProductRegistryHelper::TypeLabelList 
TypeLabelList
 
- Public Types inherited from edm::EDConsumerBase
typedef ProductLabels Labels
 
- Static Public Member Functions inherited from edm::global::EDProducerBase
static const std::string & baseType ()
 
static void fillDescriptions (ConfigurationDescriptions &descriptions)
 
static void prevalidate (ConfigurationDescriptions &descriptions)
 
- Protected Member Functions inherited from edm::ProducerBase
ProducesCollector producesCollector ()
 
- Protected Member Functions inherited from edm::EDConsumerBase
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > consumes (edm::InputTag const &tag)
 
template<BranchType B = InEvent>
EDConsumerBaseAdaptor< Bconsumes (edm::InputTag tag) noexcept
 
EDGetToken consumes (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken consumes (TypeToGet const &id, edm::InputTag const &tag)
 
ConsumesCollector consumesCollector ()
 Use a ConsumesCollector to gather consumes information from helper functions. More...
 
template<typename ProductType , BranchType B = InEvent>
void consumesMany ()
 
void consumesMany (const TypeToGet &id)
 
template<BranchType B>
void consumesMany (const TypeToGet &id)
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes ()
 
template<typename ESProduct , typename ESRecord , Transition Tr = Transition::Event>
auto esConsumes (ESInputTag const &tag)
 
template<Transition Tr = Transition::Event>
constexpr auto esConsumes ()
 
template<Transition Tr = Transition::Event>
auto esConsumes (ESInputTag tag)
 
template<Transition Tr = Transition::Event>
ESGetTokenGeneric esConsumes (eventsetup::EventSetupRecordKey const &iRecord, eventsetup::DataKey const &iKey)
 Used with EventSetupRecord::doGet. More...
 
template<typename ProductType , BranchType B = InEvent>
EDGetTokenT< ProductType > mayConsume (edm::InputTag const &tag)
 
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
template<BranchType B>
EDGetToken mayConsume (const TypeToGet &id, edm::InputTag const &tag)
 
void resetItemsToGetFrom (BranchType iType)
 

Detailed Description

Definition at line 13 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Constructor & Destructor Documentation

TtSemiLepJetCombWMassDeltaTopMass::TtSemiLepJetCombWMassDeltaTopMass ( const edm::ParameterSet cfg)
explicit

Definition at line 36 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

References Exception, and maxNJets_.

37  : jetsToken_(consumes<std::vector<pat::Jet>>(cfg.getParameter<edm::InputTag>("jets"))),
39  metsToken_(consumes<std::vector<pat::MET>>(cfg.getParameter<edm::InputTag>("mets"))),
40  maxNJets_(cfg.getParameter<int>("maxNJets")),
41  wMass_(cfg.getParameter<double>("wMass")),
42  useBTagging_(cfg.getParameter<bool>("useBTagging")),
43  bTagAlgorithm_(cfg.getParameter<std::string>("bTagAlgorithm")),
44  minBDiscBJets_(cfg.getParameter<double>("minBDiscBJets")),
45  maxBDiscLightJets_(cfg.getParameter<double>("maxBDiscLightJets")),
46  neutrinoSolutionType_(cfg.getParameter<int>("neutrinoSolutionType")) {
47  if (maxNJets_ < 4 && maxNJets_ != -1)
48  throw cms::Exception("WrongConfig") << "Parameter maxNJets can not be set to " << maxNJets_ << ". \n"
49  << "It has to be larger than 4 or can be set to -1 to take all jets.";
50 
51  produces<std::vector<std::vector<int>>>();
52  produces<int>("NumberOfConsideredJets");
53 }
edm::EDGetTokenT< std::vector< pat::Jet > > jetsToken_
EDGetTokenT< ProductType > consumes(edm::InputTag const &tag)
edm::EDGetTokenT< edm::View< reco::RecoCandidate > > lepsToken_
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
edm::EDGetTokenT< std::vector< pat::MET > > metsToken_

Member Function Documentation

bool TtSemiLepJetCombWMassDeltaTopMass::isValid ( const int &  idx,
const std::vector< pat::Jet > &  jets 
) const
inlineprivate

Definition at line 20 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by ntupleDataFormat._Object::_checkIsValid(), produce(), and core.AutoHandle.AutoHandle::ReallyLoad().

20  {
21  return (0 <= idx && idx < (int)jets.size());
22  };
vector< PseudoJet > jets
void TtSemiLepJetCombWMassDeltaTopMass::produce ( edm::StreamID  ,
edm::Event evt,
const edm::EventSetup setup 
) const
overrideprivatevirtual

Implements edm::global::EDProducerBase.

Definition at line 55 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

References bTagAlgorithm_, MEzCalculator::Calculate(), Exception, edm::Event::get(), TtSemiLepEvtPartons::HadB, mps_fire::i, isValid(), fwrapper::jets, jetsToken_, TtSemiLepEvtPartons::LepB, lepsToken_, TtSemiLepEvtPartons::LightQ, TtSemiLepEvtPartons::LightQBar, match(), maxBDiscLightJets_, maxNJets_, metsToken_, minBDiscBJets_, eostools::move(), neutrinoSolutionType_, edm::Event::put(), MEzCalculator::SetLepton(), MEzCalculator::SetMET(), mathSSE::sqrt(), useBTagging_, and wMass_.

55  {
56  auto pOut = std::make_unique<std::vector<std::vector<int>>>();
57  auto pJetsConsidered = std::make_unique<int>(0);
58 
59  std::vector<int> match;
60  for (unsigned int i = 0; i < 4; ++i)
61  match.push_back(-1);
62 
63  // get jets
64  const auto& jets = evt.get(jetsToken_);
65 
66  // get leptons
67  const auto& leps = evt.get(lepsToken_);
68 
69  // get MET
70  const auto& mets = evt.get(metsToken_);
71 
72  // skip events without lepton candidate or less than 4 jets or no MET
73  if (leps.empty() || jets.size() < 4 || mets.empty()) {
74  pOut->push_back(match);
75  evt.put(std::move(pOut));
76  *pJetsConsidered = jets.size();
77  evt.put(std::move(pJetsConsidered), "NumberOfConsideredJets");
78  return;
79  }
80 
81  unsigned maxNJets = maxNJets_;
82  if (maxNJets_ == -1 || (int)jets.size() < maxNJets_)
83  maxNJets = jets.size();
84  *pJetsConsidered = maxNJets;
85  evt.put(std::move(pJetsConsidered), "NumberOfConsideredJets");
86 
87  std::vector<bool> isBJet;
88  std::vector<bool> isLJet;
89  int cntBJets = 0;
90  if (useBTagging_) {
91  for (unsigned int idx = 0; idx < maxNJets; ++idx) {
92  isBJet.push_back((jets[idx].bDiscriminator(bTagAlgorithm_) > minBDiscBJets_));
93  isLJet.push_back((jets[idx].bDiscriminator(bTagAlgorithm_) < maxBDiscLightJets_));
94  if (jets[idx].bDiscriminator(bTagAlgorithm_) > minBDiscBJets_)
95  cntBJets++;
96  }
97  }
98 
99  // -----------------------------------------------------
100  // associate those jets that get closest to the W mass
101  // with their invariant mass to the hadronic W boson
102  // -----------------------------------------------------
103  double wDist = -1.;
104  std::vector<int> closestToWMassIndices;
105  closestToWMassIndices.push_back(-1);
106  closestToWMassIndices.push_back(-1);
107  for (unsigned idx = 0; idx < maxNJets; ++idx) {
108  if (useBTagging_ && (!isLJet[idx] || (cntBJets <= 2 && isBJet[idx])))
109  continue;
110  for (unsigned jdx = (idx + 1); jdx < maxNJets; ++jdx) {
111  if (useBTagging_ &&
112  (!isLJet[jdx] || (cntBJets <= 2 && isBJet[jdx]) || (cntBJets == 3 && isBJet[idx] && isBJet[jdx])))
113  continue;
114  reco::Particle::LorentzVector sum = jets[idx].p4() + jets[jdx].p4();
115  if (wDist < 0. || wDist > fabs(sum.mass() - wMass_)) {
116  wDist = fabs(sum.mass() - wMass_);
117  closestToWMassIndices.clear();
118  closestToWMassIndices.push_back(idx);
119  closestToWMassIndices.push_back(jdx);
120  }
121  }
122  }
123 
124  // -----------------------------------------------------
125  // build a leptonic W boson from the lepton and the MET
126  // -----------------------------------------------------
127  double neutrino_pz = 0.;
128  if (neutrinoSolutionType_ != -1) {
129  MEzCalculator mez;
130  mez.SetMET(*mets.begin());
131  if (dynamic_cast<const reco::Muon*>(&(leps.front())))
132  mez.SetLepton(leps[0], true);
133  else if (dynamic_cast<const reco::GsfElectron*>(&(leps.front())))
134  mez.SetLepton(leps[0], false);
135  else
136  throw cms::Exception("UnimplementedFeature")
137  << "Type of lepton given together with MET for solving neutrino kinematics is neither muon nor electron.\n";
138  neutrino_pz = mez.Calculate(neutrinoSolutionType_);
139  }
140  const math::XYZTLorentzVector neutrino(mets.begin()->px(),
141  mets.begin()->py(),
142  neutrino_pz,
143  sqrt(mets.begin()->px() * mets.begin()->px() +
144  mets.begin()->py() * mets.begin()->py() + neutrino_pz * neutrino_pz));
145  const reco::Particle::LorentzVector lepW = neutrino + leps.front().p4();
146 
147  // -----------------------------------------------------
148  // associate those jets to the hadronic and the leptonic
149  // b quark that minimize the difference between
150  // hadronic and leptonic top mass
151  // -----------------------------------------------------
152  double deltaTop = -1.;
153  int hadB = -1;
154  int lepB = -1;
155  if (isValid(closestToWMassIndices[0], jets) && isValid(closestToWMassIndices[1], jets)) {
156  const reco::Particle::LorentzVector hadW =
157  jets[closestToWMassIndices[0]].p4() + jets[closestToWMassIndices[1]].p4();
158  // find hadronic b candidate
159  for (unsigned idx = 0; idx < maxNJets; ++idx) {
160  if (useBTagging_ && !isBJet[idx])
161  continue;
162  // make sure it's not used up already from the hadronic W
163  if ((int)idx != closestToWMassIndices[0] && (int)idx != closestToWMassIndices[1]) {
164  reco::Particle::LorentzVector hadTop = hadW + jets[idx].p4();
165  // find leptonic b candidate
166  for (unsigned jdx = 0; jdx < maxNJets; ++jdx) {
167  if (useBTagging_ && !isBJet[jdx])
168  continue;
169  // make sure it's not used up already from the hadronic branch
170  if ((int)jdx != closestToWMassIndices[0] && (int)jdx != closestToWMassIndices[1] && jdx != idx) {
171  reco::Particle::LorentzVector lepTop = lepW + jets[jdx].p4();
172  if (deltaTop < 0. || deltaTop > fabs(hadTop.mass() - lepTop.mass())) {
173  deltaTop = fabs(hadTop.mass() - lepTop.mass());
174  hadB = idx;
175  lepB = jdx;
176  }
177  }
178  }
179  }
180  }
181  }
182 
183  match[TtSemiLepEvtPartons::LightQ] = closestToWMassIndices[0];
184  match[TtSemiLepEvtPartons::LightQBar] = closestToWMassIndices[1];
185  match[TtSemiLepEvtPartons::HadB] = hadB;
186  match[TtSemiLepEvtPartons::LepB] = lepB;
187 
188  pOut->push_back(match);
189  evt.put(std::move(pOut));
190 }
edm::EDGetTokenT< std::vector< pat::Jet > > jetsToken_
void SetLepton(const pat::Particle &lepton, bool isMuon=true)
Set lepton.
Definition: MEzCalculator.h:31
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:133
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
void SetMET(const pat::MET &MET)
Set MET.
Definition: MEzCalculator.h:25
double Calculate(int type=1)
member functions
T sqrt(T t)
Definition: SSEVec.h:19
vector< PseudoJet > jets
def move
Definition: eostools.py:511
bool get(ProductID const &oid, Handle< PROD > &result) const
Definition: Event.h:346
edm::EDGetTokenT< edm::View< reco::RecoCandidate > > lepsToken_
bool isValid(const int &idx, const std::vector< pat::Jet > &jets) const
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:10
edm::EDGetTokenT< std::vector< pat::MET > > metsToken_
math::XYZTLorentzVector LorentzVector
Lorentz vector.
Definition: Particle.h:21

Member Data Documentation

std::string TtSemiLepJetCombWMassDeltaTopMass::bTagAlgorithm_
private

Definition at line 30 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

edm::EDGetTokenT<std::vector<pat::Jet> > TtSemiLepJetCombWMassDeltaTopMass::jetsToken_
private

Definition at line 22 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

edm::EDGetTokenT<edm::View<reco::RecoCandidate> > TtSemiLepJetCombWMassDeltaTopMass::lepsToken_
private

Definition at line 25 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

double TtSemiLepJetCombWMassDeltaTopMass::maxBDiscLightJets_
private

Definition at line 32 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

int TtSemiLepJetCombWMassDeltaTopMass::maxNJets_
private
edm::EDGetTokenT<std::vector<pat::MET> > TtSemiLepJetCombWMassDeltaTopMass::metsToken_
private

Definition at line 26 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

double TtSemiLepJetCombWMassDeltaTopMass::minBDiscBJets_
private

Definition at line 31 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

int TtSemiLepJetCombWMassDeltaTopMass::neutrinoSolutionType_
private

Definition at line 33 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

bool TtSemiLepJetCombWMassDeltaTopMass::useBTagging_
private

Definition at line 29 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().

double TtSemiLepJetCombWMassDeltaTopMass::wMass_
private

Definition at line 28 of file TtSemiLepJetCombWMassDeltaTopMass.cc.

Referenced by produce().