|
|
|
#include <TopQuarkAnalysis/TopPairBSM/plugins/BoostedTopProducer.cc>
Public Member Functions | |
| BoostedTopProducer (const edm::ParameterSet &) | |
| ~BoostedTopProducer () | |
Private Member Functions | |
| virtual void | beginJob (const edm::EventSetup &) |
| virtual void | endJob () |
| virtual void | produce (edm::Event &, const edm::EventSetup &) |
| double | Psi (TLorentzVector p1, TLorentzVector p2, double mass) |
Private Attributes | |
| double | caloIsoCut_ |
| edm::Handle< std::vector < pat::Electron > > | electronHandle_ |
| edm::InputTag | eleLabel_ |
| edm::Handle< std::vector < pat::Jet > > | jetHandle_ |
| edm::InputTag | jetLabel_ |
| edm::Handle< std::vector < pat::MET > > | METHandle_ |
| edm::InputTag | metLabel_ |
| double | mTop_ |
| edm::InputTag | muoLabel_ |
| edm::Handle< std::vector < pat::Muon > > | muonHandle_ |
| edm::InputTag | solLabel_ |
Description: Class to examine boosted ttbar pairs in multiple mass regions.
This will produce a ttbar solution, which will take one of two forms: a) lv jj jj Full reconstruction.
ttbar-> (hadt -> (hadW -> hadp + hadq) + hadb) + (lept -> (lepW -> lepton + neutrino) + lepb)
b) lv jj (j) Partial reconstruction, associate at least 1 jet to the lepton hemisphere, and at least one jet in the opposite hemisphere.
ttbar-> (hadt -> (hadJet1 [+ hadJet2] ) ) + (lept -> (lepW -> lepton + neutrino) + lepJet1 )
There will also be two subcategories of (b) that will correspond to physics cases:
b1) Lepton is isolated: Moderate ttbar mass. b2) Lepton is nonisolated: High ttbar mass.
Implementation: To implement this, we use the NamedCompositeCandidate structures from the Candidate model. This provides flexibility in the definition of the output objects and allows automatic plotting in the Starter Kit. We use the PAT objects to construct the ttbar solutions in the different ranges as follows: a) Full reconstruction: We use TtSemiEventSolutions made upstream of this module. b) Partial reconstruction: Association of variables using the "Psi" variable, which is a more rapidity-invariant version of deltaR.
Definition at line 84 of file BoostedTopProducer.h.
| BoostedTopProducer::BoostedTopProducer | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 23 of file BoostedTopProducer.cc.
: eleLabel_ (iConfig.getParameter<edm::InputTag> ("electronLabel")), muoLabel_ (iConfig.getParameter<edm::InputTag> ("muonLabel")), jetLabel_ (iConfig.getParameter<edm::InputTag> ("jetLabel")), metLabel_ (iConfig.getParameter<edm::InputTag> ("metLabel")), solLabel_ (iConfig.getParameter<edm::InputTag> ("solLabel")), caloIsoCut_ (iConfig.getParameter<double> ("caloIsoCut") ), mTop_ (iConfig.getParameter<double> ("mTop") ) { //register products produces<std::vector<reco::CompositeCandidate> > (); }
| BoostedTopProducer::~BoostedTopProducer | ( | ) |
Definition at line 37 of file BoostedTopProducer.cc.
{
}
| void BoostedTopProducer::beginJob | ( | const edm::EventSetup & | ) | [private, virtual] |
Definition at line 375 of file BoostedTopProducer.cc.
{
}
| void BoostedTopProducer::endJob | ( | void | ) | [private, virtual] |
| void BoostedTopProducer::produce | ( | edm::Event & | iEvent, |
| const edm::EventSetup & | iSetup | ||
| ) | [private, virtual] |
Implements edm::EDProducer.
Definition at line 48 of file BoostedTopProducer.cc.
References reco::CompositeCandidate::addDaughter(), METzCalculator::Calculate(), caloIsoCut_, gather_cfg::cout, debug, metsig::electron, HI_PhotonSkim_cff::electrons, eleLabel_, reco::LeafCandidate::energy(), edm::Event::getByLabel(), METzCalculator::IsComplex(), reco::isMuon(), edm::HandleBase::isValid(), jetLabel_, fwrapper::jets, TtEvent::kMVADisc, metLabel_, mTop_, muoLabel_, metsig::muon, patZpeak::muons, reco::CompositeCandidate::numberOfDaughters(), Pi, Psi(), edm::Event::put(), reco::LeafCandidate::px(), reco::LeafCandidate::py(), reco::LeafCandidate::pz(), alignCSCRings::s, METzCalculator::SetMET(), METzCalculator::SetMuon(), reco::LeafCandidate::setPz(), solLabel_, and TablePrint::write.
{
using namespace edm;
bool debug = false;
// -----------------------------------------------------
// get the bare PAT objects
// -----------------------------------------------------
edm::Handle<std::vector<pat::Muon> > muonHandle;
iEvent.getByLabel(muoLabel_,muonHandle);
std::vector<pat::Muon> const & muons = *muonHandle;
edm::Handle<std::vector<pat::Jet> > jetHandle;
iEvent.getByLabel(jetLabel_,jetHandle);
std::vector<pat::Jet> const & jets = *jetHandle;
edm::Handle<std::vector<pat::Electron> > electronHandle;
iEvent.getByLabel(eleLabel_,electronHandle);
std::vector<pat::Electron> const & electrons = *electronHandle;
edm::Handle<std::vector<pat::MET> > metHandle;
iEvent.getByLabel(metLabel_,metHandle);
std::vector<pat::MET> const & mets = *metHandle;
// -----------------------------------------------------
// Event Preselection:
// <= 1 isolated electron or muon
// >= 1 electron or muon
// >= 2 jets
// >= 1 missing et
//
// To explain:
// We want to look at leptons within "top jets" in some
// cases. This means the isolation will kill those events.
// However, if there IS an isolated lepton, we want only
// one of them.
//
// So to select the prompt W lepton, the logic is:
// 1. If there is an isolated lepton, accept it as the W lepton.
// 2. Else, take the highest Pt lepton (possibly non-isolated)
//
// -----------------------------------------------------
bool preselection = true;
// This will hold the prompt W lepton candidate, and a
// maximum pt decision variable
double maxWLeptonPt = -1;
reco::Candidate const * Wlepton = 0;
// ----------------------
// Find isolated muons, and highest pt lepton
// ----------------------
std::vector<pat::Muon>::const_iterator isolatedMuon = muons.end();
std::vector<pat::Muon>::const_iterator muon = muons.end();
bool nIsolatedMuons = 0;
std::vector<pat::Muon>::const_iterator muonIt = muons.begin(),
muonEnd = muons.end();
for (; muonIt != muonEnd; ++muonIt ) {
// Find highest pt lepton
double pt = muonIt->pt();
if ( pt > maxWLeptonPt ) {
maxWLeptonPt = pt;
muon = muonIt;
}
// Find any isolated muons
double caloIso = muonIt->caloIso();
if ( caloIso >= 0 && caloIso < caloIsoCut_ ) {
nIsolatedMuons++;
isolatedMuon = muonIt;
}
}
// ----------------------
// Find isolated electrons, and highest pt lepton
// ----------------------
std::vector<pat::Electron>::const_iterator isolatedElectron = electrons.end();
std::vector<pat::Electron>::const_iterator electron = electrons.end();
bool nIsolatedElectrons = 0;
std::vector<pat::Electron>::const_iterator electronIt = electrons.begin(),
electronEnd = electrons.end();
for (; electronIt != electronEnd; ++electronIt ) {
// Find highest pt lepton
double pt = electronIt->pt();
if ( pt > maxWLeptonPt ) {
maxWLeptonPt = pt;
electron = electronIt;
}
// Find any isolated electrons
double caloIso = electronIt->caloIso();
if ( caloIso >= 0 && caloIso < caloIsoCut_ ) {
nIsolatedElectrons++;
isolatedElectron = electronIt;
}
}
// ----------------------
// Now decide on the "prompt" lepton from the W:
// Choose isolated leptons over all, and if no isolated,
// then take highest pt lepton.
// ----------------------
bool isMuon = true;
if ( isolatedMuon != muonEnd ) { muon = isolatedMuon; isMuon = true; }
else if ( isolatedElectron != electronEnd ) { electron = isolatedElectron; isMuon = false; }
else {
// Set to the highest pt lepton
if ( muon != muonEnd && electron == electronEnd ) isMuon = true;
else if ( muon == muonEnd && electron != electronEnd ) isMuon = false;
else if ( muon != muonEnd && electron != electronEnd ) {
isMuon = muon->pt() > electron->pt();
}
}
// ----------------------
// Veto events that have more than one isolated lepton
// ----------------------
int nIsolatedLeptons = nIsolatedMuons + nIsolatedElectrons;
if ( nIsolatedLeptons > 1 ) {
preselection = false;
}
// ----------------------
// Veto events that have no prompt lepton candidates
// ----------------------
if ( muon == muonEnd && electron == electronEnd ) {
preselection = false;
}
// ----------------------
// Veto events with < 2 jets or no missing et
// ----------------------
if ( jets.size() < 2 ||
mets.size() == 0 ) {
preselection = false;
}
bool write = false;
// -----------------------------------------------------
//
// CompositeCandidates to store the event solution.
// This will take one of two forms:
// a) lv jj jj Full reconstruction.
//
// ttbar->
// (hadt -> (hadW -> hadp + hadq) + hadb) +
// (lept -> (lepW -> lepton + neutrino) + lepb)
//
// b) lv jj (j) Partial reconstruction, associate
// at least 1 jet to the lepton
// hemisphere, and at least one jet in
// the opposite hemisphere.
//
// ttbar->
// (hadt -> (hadJet1 [+ hadJet2] ) ) +
// (lept -> (lepW -> lepton + neutrino) + lepJet1 )
//
// There will also be two subcategories of (b) that
// will correspond to physics cases:
//
// b1) Lepton is isolated: Moderate ttbar mass.
// b2) Lepton is nonisolated: High ttbar mass.
//
// -----------------------------------------------------
reco::CompositeCandidate ttbar("ttbar");
AddFourMomenta addFourMomenta;
// Main decisions after preselection
if ( preselection ) {
if ( debug ) cout << "Preselection is satisfied" << endl;
if ( debug ) cout << "Jets.size() = " << jets.size() << endl;
// This will be modified for the z solution, so make a copy
pat::MET neutrino( mets[0] );
// 1. First examine the low mass case with 4 jets and widely separated
// products. We take out the TtSemiLeptonicEvent from the TQAF and
// form the ttbar invariant mass.
if ( jets.size() >= 4 ) {
if ( debug ) cout << "Getting ttbar semileptonic solution" << endl;
// get the ttbar semileptonic event solution if there are more than 3 jets
edm::Handle< TtSemiLeptonicEvent > eSol;
iEvent.getByLabel(solLabel_, eSol);
// Have solution, continue
if ( eSol.isValid() ) {
if ( debug ) cout << "Got a nonzero size solution vector" << endl;
// Just set the ttbar solution to the best ttbar solution from
// TtSemiEvtSolutionMaker
ttbar = eSol->eventHypo(TtSemiLeptonicEvent::kMVADisc);
write = true;
}
// No ttbar solution with 4 jets, something is weird, print a warning
else {
edm::LogWarning("DataNotFound") << "BoostedTopProducer: Cannot find TtSemiEvtSolution\n";
}
}
// 2. With 2 or 3 jets, we decide based on the separation between
// the lepton and the closest jet in that hemisphere whether to
// consider it "moderate" or "high" mass.
else if ( jets.size() == 2 || jets.size() == 3 ) {
// ------------------------------------------------------------------
// First create a leptonic W candidate
// ------------------------------------------------------------------
reco::CompositeCandidate lepW("lepW");
if ( isMuon ) {
if ( debug ) cout << "Adding muon as daughter" << endl;
lepW.addDaughter( *muon, "muon" );
} else {
if ( debug ) cout << "Adding electron as daughter" << endl;
lepW.addDaughter( *electron, "electron" );
}
if ( debug ) cout << "Adding neutrino as daughter" << endl;
lepW.addDaughter ( neutrino, "neutrino");
addFourMomenta.set( lepW );
bool nuzHasComplex = false;
METzCalculator zcalculator;
zcalculator.SetMET( neutrino );
if ( isMuon )
zcalculator.SetMuon( *muon );
else
zcalculator.SetMuon( *electron ); // This name is misleading, should be setLepton
double neutrinoPz = zcalculator.Calculate(1);// closest to the lepton Pz
if (zcalculator.IsComplex()) nuzHasComplex = true;
// Set the neutrino pz
neutrino.setPz( neutrinoPz );
if ( debug ) cout << "Set neutrino pz to " << neutrinoPz << endl;
// ------------------------------------------------------------------
// Next ensure that there is a jet within the hemisphere of the
// leptonic W, and one in the opposite hemisphere
// ------------------------------------------------------------------
reco::CompositeCandidate hadt("hadt");
reco::CompositeCandidate lept("lept");
if ( debug ) cout << "Adding lepW as daughter" << endl;
lept.addDaughter( lepW, "lepW" );
std::string hadName("hadJet");
std::string lepName("lepJet");
// Get the W momentum
TLorentzVector p4_W (lepW.px(), lepW.py(), lepW.pz(), lepW.energy() );
// Loop over the jets
std::vector<pat::Jet>::const_iterator jetit = jets.begin(),
jetend = jets.end();
unsigned long ii = 1; // Count by 1 for naming histograms
for ( ; jetit != jetend; ++jetit, ++ii ) {
// Get this jet's momentum
TLorentzVector p4_jet (jetit->px(), jetit->py(), jetit->pz(), jetit->energy() );
// Calculate psi (like DeltaR, only more invariant under Rapidity)
double psi = Psi( p4_W, p4_jet, mTop_ );
// Get jets that are in the leptonic hemisphere
if ( psi < TMath::Pi() ) {
// Add this jet to the leptonic top
std::stringstream s;
s << lepName << ii;
if ( debug ) cout << "Adding daughter " << s.str() << endl;
lept.addDaughter( *jetit, s.str() );
}
// Get jets that are in the hadronic hemisphere
if ( psi > TMath::Pi() ) {
// Add this jet to the hadronic top. We don't
// make any W hypotheses in this case, since
// we cannot determine which of the three
// jets are merged.
std::stringstream s;
s << hadName << ii;
if ( debug ) cout << "Adding daughter " << s.str() << endl;
hadt.addDaughter( *jetit, s.str() );
}
} // end loop over jets
addFourMomenta.set( lept );
addFourMomenta.set( hadt );
bool lepWHasJet = lept.numberOfDaughters() >= 2; // W and >= 1 jet
bool hadWHasJet = hadt.numberOfDaughters() >= 1; // >= 1 jet
if ( lepWHasJet && hadWHasJet ) {
if ( debug ) cout << "Adding daughters lept and hadt" << endl;
ttbar.addDaughter( lept, "lept");
ttbar.addDaughter( hadt, "hadt");
addFourMomenta.set( ttbar );
write = true;
} // end of hadronic jet and leptonic jet
} // end if there are 2 or 3 jets
} // end if preselection is satisfied
// Write the solution to the event record
std::vector<reco::CompositeCandidate> ttbarList;
if ( write ) {
if ( debug ) cout << "Writing out" << endl;
ttbarList.push_back( ttbar );
}
std::auto_ptr<std::vector<reco::CompositeCandidate> > pTtbar ( new std::vector<reco::CompositeCandidate>(ttbarList) );
iEvent.put( pTtbar );
}
| double BoostedTopProducer::Psi | ( | TLorentzVector | p1, |
| TLorentzVector | p2, | ||
| double | mass | ||
| ) | [private] |
Definition at line 385 of file BoostedTopProducer.cc.
Referenced by produce().
{
TLorentzVector ptot = p1 + p2;
Double_t theta1 = TMath::ACos( (p1.Vect().Dot(ptot.Vect()))/(p1.P()*ptot.P()) );
Double_t theta2 = TMath::ACos( (p2.Vect().Dot(ptot.Vect()))/(p2.P()*ptot.P()) );
//Double_t sign = 1.;
//if ( (theta1+theta2) > (TMath::Pi()/2) ) sign = -1.;
double th1th2 = theta1 + theta2;
double psi = (p1.P()+p2.P())*TMath::Abs(TMath::Sin(th1th2))/(2.* mass );
if ( th1th2 > (TMath::Pi()/2) )
psi = (p1.P()+p2.P())*( 1. + TMath::Abs(TMath::Cos(th1th2)))/(2.* mass );
return psi;
}
double BoostedTopProducer::caloIsoCut_ [private] |
Definition at line 110 of file BoostedTopProducer.h.
Referenced by produce().
edm::Handle<std::vector<pat::Electron> > BoostedTopProducer::electronHandle_ [private] |
Definition at line 105 of file BoostedTopProducer.h.
edm::InputTag BoostedTopProducer::eleLabel_ [private] |
Definition at line 97 of file BoostedTopProducer.h.
Referenced by produce().
edm::Handle<std::vector<pat::Jet> > BoostedTopProducer::jetHandle_ [private] |
Definition at line 106 of file BoostedTopProducer.h.
edm::InputTag BoostedTopProducer::jetLabel_ [private] |
Definition at line 99 of file BoostedTopProducer.h.
Referenced by produce().
edm::Handle<std::vector<pat::MET> > BoostedTopProducer::METHandle_ [private] |
Definition at line 107 of file BoostedTopProducer.h.
edm::InputTag BoostedTopProducer::metLabel_ [private] |
Definition at line 100 of file BoostedTopProducer.h.
Referenced by produce().
double BoostedTopProducer::mTop_ [private] |
Definition at line 111 of file BoostedTopProducer.h.
Referenced by produce().
edm::InputTag BoostedTopProducer::muoLabel_ [private] |
Definition at line 98 of file BoostedTopProducer.h.
Referenced by produce().
edm::Handle<std::vector<pat::Muon> > BoostedTopProducer::muonHandle_ [private] |
Definition at line 104 of file BoostedTopProducer.h.
edm::InputTag BoostedTopProducer::solLabel_ [private] |
Definition at line 101 of file BoostedTopProducer.h.
Referenced by produce().
1.7.3