66 std::set<int> bHadronIdxs;
67 for (
size_t i=0,
n=genParticleHandle->size();
i<
n; ++
i ) {
70 if ( status == 1 )
continue;
78 std::vector<size_t> leptonIdxs;
79 for (
size_t i=0, n=finalStateHandle->size();
i<
n; ++
i ) {
82 if ( p.
status() != 1 )
continue;
91 leptonIdxs.push_back(
i);
93 case 12:
case 14:
case 16:
104 std::vector<fastjet::PseudoJet> fjLepInputs;
105 fjLepInputs.reserve(leptonIdxs.size());
106 for (
auto index : leptonIdxs ) {
110 fjLepInputs.push_back(fastjet::PseudoJet(p.
px(), p.
py(), p.
pz(), p.
energy()));
111 fjLepInputs.back().set_user_index(
index);
115 fastjet::ClusterSequence fjLepClusterSeq(fjLepInputs, *
fjLepDef_);
116 std::vector<fastjet::PseudoJet> fjLepJets = fastjet::sorted_by_pt(fjLepClusterSeq.inclusive_jets(
minLeptonPt_));
119 leptons->reserve(fjLepJets.size());
120 std::set<size_t> lepDauIdxs;
121 for (
auto& fjJet : fjLepJets ) {
125 const std::vector<fastjet::PseudoJet> fjConstituents = fastjet::sorted_by_pt(fjJet.constituents());
127 std::vector<reco::CandidatePtr> constituents;
129 for (
auto& fjConstituent : fjConstituents ) {
130 const size_t index = fjConstituent.user_index();
132 const int absPdgId =
abs(cand->
pdgId());
133 if ( absPdgId == 11
or absPdgId == 13 ) {
134 if ( lepCand.
isNonnull() and lepCand->
pt() > cand->
pt() )
continue;
137 constituents.push_back(cand);
139 if ( lepCand.
isNull() )
continue;
141 const LorentzVector jetP4(fjJet.px(), fjJet.py(), fjJet.pz(), fjJet.E());
145 lepJet.setPdgId(lepCand->
pdgId());
146 lepJet.setCharge(lepCand->
charge());
148 const double jetArea = fjJet.has_area() ? fjJet.area() : 0;
149 lepJet.setJetArea(jetArea);
154 for (
auto& fjConstituent : fjConstituents ) {
155 lepDauIdxs.insert(fjConstituent.user_index());
162 std::vector<fastjet::PseudoJet> fjJetInputs;
163 fjJetInputs.reserve(nStables);
164 for (
size_t i=0, n=finalStateHandle->size();
i<
n; ++
i ) {
166 if ( p.
status() != 1 )
continue;
170 if ( absId == 12
or absId == 14
or absId == 16 )
continue;
171 if ( lepDauIdxs.find(
i) != lepDauIdxs.end() )
continue;
173 fjJetInputs.push_back(fastjet::PseudoJet(p.
px(), p.
py(), p.
pz(), p.
energy()));
174 fjJetInputs.back().set_user_index(
i);
177 for (
auto index : bHadronIdxs ) {
181 const double scale = 1
e-20/p.
p();
183 fjJetInputs.back().set_user_index(index);
187 fastjet::ClusterSequence fjJetClusterSeq(fjJetInputs, *
fjJetDef_);
188 std::vector<fastjet::PseudoJet> fjJets = fastjet::sorted_by_pt(fjJetClusterSeq.inclusive_jets(
minJetPt_));
191 jets->reserve(fjJets.size());
192 std::vector<size_t> bjetIdxs, ljetIdxs;
193 for (
auto& fjJet : fjJets ) {
197 const std::vector<fastjet::PseudoJet> fjConstituents = fastjet::sorted_by_pt(fjJet.constituents());
199 std::vector<reco::CandidatePtr> constituents;
200 bool hasBHadron =
false;
201 for (
size_t j=0,
m=fjConstituents.size(); j<
m; ++j ) {
202 const int index = fjConstituents[j].user_index();
203 if ( bHadronIdxs.find(index) != bHadronIdxs.end() ) {
208 constituents.push_back(cand);
211 const LorentzVector jetP4(fjJet.px(), fjJet.py(), fjJet.pz(), fjJet.E());
215 const double jetArea = fjJet.has_area() ? fjJet.area() : 0;
216 genJet.setJetArea(jetArea);
219 bjetIdxs.push_back(
jets->size());
222 ljetIdxs.push_back(
jets->size());
225 jets->push_back(genJet);
231 if ( bjetIdxs.size() < 2 )
break;
238 if ( q1*q2 > 0 )
break;
248 int selNu1 = -1, selNu2 = -1;
249 for (
int i=0;
i<2; ++
i ) {
251 for (
int j=0; j<2; ++j ) {
252 if (
i == j )
continue;
254 const double newDm = dm1+dm2;
256 if ( newDm < dm ) { dm = newDm; selNu1 =
i; selNu2 = j; }
259 if ( dm >= 1e9 )
break;
263 const auto w1LVec = lepton1.p4()+nu1.p4();
264 const auto w2LVec = lepton2.p4()+nu2.p4();
268 int selB1 = -1, selB2 = -1;
269 for (
auto i : bjetIdxs ) {
271 for (
auto j : bjetIdxs ) {
272 if (
i == j )
continue;
274 const double newDm = dm1+dm2;
276 if ( newDm < dm ) { dm = newDm; selB1 =
i; selB2 = j; }
279 if ( dm >= 1e9 )
break;
283 const auto t1LVec = w1LVec +
bJet1.p4();
284 const auto t2LVec = w2LVec +
bJet2.p4();
287 const int lepQ1 = lepton1.charge();
288 const int lepQ2 = lepton2.charge();
316 else if (
leptons->size() >= 1 and
neutrinos->size() >= 1 and ljetIdxs.size() >= 2 ) {
318 const auto& lepton =
leptons->at(0);
322 bool hasVetoLepton =
false;
323 for (
auto vetoLeptonItr =
next(
leptons->begin()); vetoLeptonItr !=
leptons->end(); ++vetoLeptonItr ) {
326 hasVetoLepton =
true;
329 if ( hasVetoLepton )
break;
332 double metX = 0, metY = 0;
334 metX += neutrino.px();
335 metY += neutrino.py();
337 const double metPt = std::hypot(metX, metY);
340 const double mtW =
std::sqrt( 2*lepton.pt()*metPt*
cos(lepton.phi()-atan2(metX, metY)) );
349 const double lpz = lepton.pz(), le = lepton.energy(), lpt = lepton.pt();
350 const double cf = (
wMass_*
wMass_)/2 + lepton.px()*metX + lepton.py()*metY;
351 const double det = cf*cf*lpz*lpz - lpt*lpt*(le*le*metPt*metPt - cf*cf);
352 const double pz = (cf*lpz + (cf < 0 ? -
sqrt(det) : sqrt(det)))/lpt/lpt;
354 const auto w1LVec = lepton.p4()+nu1P4;
359 for (
auto b1Itr = bjetIdxs.begin(); b1Itr != bjetIdxs.end(); ++b1Itr ) {
360 const double dR =
deltaR(
jets->at(*b1Itr).p4(), w1LVec);
361 if ( dR < minDR ) { selB1 = *b1Itr; minDR =
dR; }
363 if ( selB1 == -1 )
break;
365 const auto t1LVec = w1LVec +
bJet1.p4();
368 const auto& wJet1 =
jets->at(ljetIdxs[0]);
369 const auto& wJet2 =
jets->at(ljetIdxs[1]);
370 const auto w2LVec = wJet1.p4() + wJet2.p4();
375 for (
auto i : bjetIdxs ) {
376 if (
int(
i) == selB1 )
continue;
378 if ( newDm < dm ) { dm = newDm; selB2 =
i; }
380 if ( dm >= 1e9 )
break;
383 const auto t2LVec = w2LVec +
bJet2.p4();
385 const int q = lepton.charge();
441 event.put(
std::move(neutrinos),
"neutrinos");
const edm::EDGetTokenT< edm::View< reco::Candidate > > genParticleToken_
bool isBHadron(const reco::Candidate *p) const
std::vector< GenParticle > GenParticleCollection
collection of GenParticles
reco::Particle::Point genVertex_
virtual double pz() const =0
z coordinate of momentum vector
const double maxLeptonEtaSemilepton_
common ppss p3p6s2 common epss epspn46 common const1 w2
const double maxLeptonEta_
const double minVetoLeptonPtSemilepton_
const double maxLeptonEtaDilepton_
std::shared_ptr< JetDef > fjJetDef_
std::vector< GenJet > GenJetCollection
collection of GenJet objects
bool isFromHadron(const reco::Candidate *p) const
virtual const Candidate * mother(size_type i=0) const =0
return pointer to mother
std::shared_ptr< JetDef > fjLepDef_
const double minMtWSemiLepton_
const double minMETSemiLepton_
const double minLeptonPtSemilepton_
virtual int status() const =0
status word
virtual double energy() const =0
energy
virtual double py() const =0
y coordinate of momentum vector
bJet1
do the razor with one or two b jets (medium CSV)
auto const T2 &decltype(t1.eta()) t2
virtual const LorentzVector & p4() const =0
four-momentum Lorentz vector
virtual int pdgId() const =0
PDG identifier.
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventID const &, edm::Timestamp const & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
virtual size_type numberOfMothers() const =0
number of mothers (zero or one in most of but not all the cases)
Cos< T >::type cos(const T &t)
bool isNull() const
Checks for null.
Abs< T >::type abs(const T &t)
Jets made from MC generator particles.
math::XYZTLorentzVector LorentzVector
virtual double p() const =0
magnitude of momentum vector
bool isNonnull() const
Checks for non-null.
double deltaR(double eta1, double eta2, double phi1, double phi2)
virtual double pt() const =0
transverse momentum
math::XYZTLorentzVector LorentzVector
Lorentz vector.
virtual int charge() const =0
electric charge
const double maxVetoLeptonEtaSemilepton_
const double minLeptonPt_
virtual const Point & vertex() const =0
vertex position
const double minDileptonMassDilepton_
virtual double px() const =0
x coordinate of momentum vector
const edm::EDGetTokenT< edm::View< reco::Candidate > > finalStateToken_
const double minLeptonPtDilepton_
void writeSpecific(reco::CaloJet &jet, reco::Particle::LorentzVector const &p4, reco::Particle::Point const &point, std::vector< reco::CandidatePtr > const &constituents, edm::EventSetup const &c)
Public Member Functions inherited from edm::stream::EDProducer<>
1.8.11
