![]() |
![]() |
#include <TagProbeFitTreeProducer.cc>
Public Member Functions | |
TagProbeFitTreeProducer (const edm::ParameterSet &) | |
~TagProbeFitTreeProducer () | |
Private Member Functions | |
virtual void | analyze (const edm::Event &, const edm::EventSetup &) |
bool | checkMother (const reco::GenParticleRef &ref) const |
Return true if ref is not null and has an ancestor with pdgId inside 'motherPdgId_'. | |
virtual void | endJob () |
Private Attributes | |
edm::InputTag | allProbes_ |
InputTag to the collection of all probes. | |
bool | checkMotherInUnbiasEff_ |
Check mother pdgId in unbiased inefficiency measurement. | |
bool | isMC_ |
Is this sample MC? | |
bool | makeMCUnbiasTree_ |
Do we have to compute this. | |
std::auto_ptr < tnp::BaseTreeFiller > | mcFiller_ |
std::auto_ptr < tnp::BaseTreeFiller > | mcUnbiasFiller_ |
The object that actually computes variables and fills the tree for unbiased MC. | |
std::set< int32_t > | motherPdgId_ |
Possible pdgids for the mother. If empty, any truth-matched mu will be considered good. | |
std::auto_ptr < tnp::BaseTreeFiller > | oldTagFiller_ |
std::auto_ptr < tnp::BaseTreeFiller > | pairFiller_ |
edm::InputTag | probeMatches_ |
std::auto_ptr < tnp::BaseTreeFiller > | tagFiller_ |
edm::InputTag | tagMatches_ |
InputTag to an edm::Association<reco::GenParticle> from tags & probes to MC truth. | |
tnp::TagProbePairMaker | tagProbePairMaker_ |
The object that produces pairs of tags and probes, making any arbitration needed. | |
std::auto_ptr< tnp::TPTreeFiller > | treeFiller_ |
The object that actually computes variables and fills the tree for T&P. |
Description: <one line="" class="" summary>="">
Implementation: <Notes on="" implementation>="">
Definition at line 51 of file TagProbeFitTreeProducer.cc.
TagProbeFitTreeProducer::TagProbeFitTreeProducer | ( | const edm::ParameterSet & | iConfig | ) | [explicit] |
Definition at line 94 of file TagProbeFitTreeProducer.cc.
References edm::ParameterSet::addParameter(), allProbes_, edm::ParameterSet::empty(), edm::ParameterSet::existsAs(), edm::ParameterSet::getParameter(), isMC_, makeMCUnbiasTree_, mcFiller_, mcUnbiasFiller_, motherPdgId_, pairFiller_, probeMatches_, tagFiller_, tagMatches_, and treeFiller_.
: isMC_(iConfig.getParameter<bool>("isMC")), makeMCUnbiasTree_(isMC_ ? iConfig.getParameter<bool>("makeMCUnbiasTree") : false), checkMotherInUnbiasEff_(makeMCUnbiasTree_ ? iConfig.getParameter<bool>("checkMotherInUnbiasEff") : false), tagProbePairMaker_(iConfig), treeFiller_(new tnp::TPTreeFiller(iConfig)), oldTagFiller_((iConfig.existsAs<bool>("fillTagTree") && iConfig.getParameter<bool>("fillTagTree")) ? new tnp::BaseTreeFiller("tag_tree",iConfig) : 0) { if (isMC_) { // For mc efficiency we need the MC matches for tags & probes tagMatches_ = iConfig.getParameter<edm::InputTag>("tagMatches"); probeMatches_ = iConfig.getParameter<edm::InputTag>("probeMatches"); //.. and the pdgids of the possible mothers if (iConfig.existsAs<int32_t>("motherPdgId")) { motherPdgId_.insert(iConfig.getParameter<int32_t>("motherPdgId")); } else { std::vector<int32_t> motherIds = iConfig.getParameter<std::vector<int32_t> >("motherPdgId"); motherPdgId_.insert(motherIds.begin(), motherIds.end()); } // For unbiased efficiency we also need the collection of all probes if (makeMCUnbiasTree_) { allProbes_ = iConfig.getParameter<edm::InputTag>("allProbes"); mcUnbiasFiller_.reset(new tnp::BaseTreeFiller("mcUnbias_tree",iConfig)); } } edm::ParameterSet tagPSet; if (iConfig.existsAs<edm::ParameterSet>("tagVariables")) tagPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("tagVariables")); if (iConfig.existsAs<edm::ParameterSet>("tagFlags" )) tagPSet.addParameter<edm::ParameterSet>("flags", iConfig.getParameter<edm::ParameterSet>("tagFlags")); if (!tagPSet.empty()) { tagFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, tagPSet, "tag_")); } edm::ParameterSet mcPSet; if (iConfig.existsAs<edm::ParameterSet>("mcVariables")) mcPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("mcVariables")); if (iConfig.existsAs<edm::ParameterSet>("mcFlags" )) mcPSet.addParameter<edm::ParameterSet>("flags", iConfig.getParameter<edm::ParameterSet>("mcFlags")); if (!mcPSet.empty()) { mcFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, mcPSet, "mc_")); } edm::ParameterSet pairPSet; if (iConfig.existsAs<edm::ParameterSet>("pairVariables")) pairPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("pairVariables")); if (iConfig.existsAs<edm::ParameterSet>("pairFlags" )) pairPSet.addParameter<edm::ParameterSet>("flags", iConfig.getParameter<edm::ParameterSet>("pairFlags")); if (!pairPSet.empty()) { pairFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, pairPSet, "pair_")); } }
TagProbeFitTreeProducer::~TagProbeFitTreeProducer | ( | ) |
Definition at line 136 of file TagProbeFitTreeProducer.cc.
{ }
void TagProbeFitTreeProducer::analyze | ( | const edm::Event & | iEvent, |
const edm::EventSetup & | iSetup | ||
) | [private, virtual] |
Implements edm::EDAnalyzer.
Definition at line 147 of file TagProbeFitTreeProducer.cc.
References allProbes_, checkMother(), checkMotherInUnbiasEff_, edm::Event::getByLabel(), i, isMC_, edm::Ref< C, T, F >::isNonnull(), makeMCUnbiasTree_, mcFiller_, mcUnbiasFiller_, n, oldTagFiller_, pairFiller_, probeMatches_, tnp::TagProbePairMaker::run(), align_tpl::src, tagFiller_, tagMatches_, tagProbePairMaker_, and treeFiller_.
{ using namespace edm; using namespace std; Handle<reco::CandidateView> src, allProbes; Handle<Association<vector<reco::GenParticle> > > tagMatches, probeMatches; treeFiller_->init(iEvent); // read out info from the event if needed (external vars, list of passing probes, ...) if (oldTagFiller_.get()) oldTagFiller_->init(iEvent); if (tagFiller_.get()) tagFiller_->init(iEvent); if (pairFiller_.get()) pairFiller_->init(iEvent); if (mcFiller_.get()) mcFiller_->init(iEvent); // on mc we want to load also the MC match info if (isMC_) { iEvent.getByLabel(tagMatches_, tagMatches); iEvent.getByLabel(probeMatches_, probeMatches); } // get the list of (tag+probe) pairs, performing arbitration tnp::TagProbePairs pairs = tagProbePairMaker_.run(iEvent); // loop on them to fill the tree for (tnp::TagProbePairs::const_iterator it = pairs.begin(), ed = pairs.end(); it != ed; ++it) { // on mc, fill mc info (on non-mc, let it to 'true', the treeFiller will ignore it anyway bool mcTrue = false; if (isMC_) { reco::GenParticleRef mtag = (*tagMatches)[it->tag], mprobe = (*probeMatches)[it->probe]; mcTrue = checkMother(mtag) && checkMother(mprobe); if (mcTrue && mcFiller_.get()) mcFiller_->fill(reco::CandidateBaseRef(mprobe)); } // fill in the variables for this t+p pair if (tagFiller_.get()) tagFiller_->fill(it->tag); if (oldTagFiller_.get()) oldTagFiller_->fill(it->tag); if (pairFiller_.get()) pairFiller_->fill(it->pair); treeFiller_->fill(it->probe, it->mass, mcTrue); } if (isMC_ && makeMCUnbiasTree_) { // read full collection of probes iEvent.getByLabel(allProbes_, allProbes); // init the tree filler mcUnbiasFiller_->init(iEvent); // loop on probes for (size_t i = 0, n = allProbes->size(); i < n; ++i) { const reco::CandidateBaseRef & probe = allProbes->refAt(i); // check mc match, and possibly mother match reco::GenParticleRef probeMatch = (*probeMatches)[probe]; bool probeOk = checkMotherInUnbiasEff_ ? checkMother(probeMatch) : probeMatch.isNonnull(); // fill the tree only for good ones if (probeOk) mcUnbiasFiller_->fill(probe); } } }
bool TagProbeFitTreeProducer::checkMother | ( | const reco::GenParticleRef & | ref | ) | const [private] |
Return true if ref is not null and has an ancestor with pdgId inside 'motherPdgId_'.
Definition at line 202 of file TagProbeFitTreeProducer.cc.
References abs, edm::RefVector< C, T, F >::begin(), edm::RefVector< C, T, F >::end(), edm::Ref< C, T, F >::isNull(), m, and motherPdgId_.
Referenced by analyze().
{ if (ref.isNull()) return false; if (motherPdgId_.empty()) return true; if (motherPdgId_.find(abs(ref->pdgId())) != motherPdgId_.end()) return true; reco::GenParticle::mothers m = ref->motherRefVector(); for (reco::GenParticle::mothers::const_iterator it = m.begin(), e = m.end(); it != e; ++it) { if (checkMother(*it)) return true; } return false; }
void TagProbeFitTreeProducer::endJob | ( | void | ) | [private, virtual] |
Reimplemented from edm::EDAnalyzer.
Definition at line 216 of file TagProbeFitTreeProducer.cc.
References edm::getProcessParameterSet(), and treeFiller_.
{ // ask to write the current PSet info into the TTree header treeFiller_->writeProvenance(edm::getProcessParameterSet()); }
InputTag to the collection of all probes.
Definition at line 77 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
bool TagProbeFitTreeProducer::checkMotherInUnbiasEff_ [private] |
Check mother pdgId in unbiased inefficiency measurement.
Definition at line 75 of file TagProbeFitTreeProducer.cc.
Referenced by analyze().
bool TagProbeFitTreeProducer::isMC_ [private] |
Is this sample MC?
Definition at line 63 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
bool TagProbeFitTreeProducer::makeMCUnbiasTree_ [private] |
Do we have to compute this.
Definition at line 73 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
std::auto_ptr<tnp::BaseTreeFiller> TagProbeFitTreeProducer::mcFiller_ [private] |
Definition at line 88 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
std::auto_ptr<tnp::BaseTreeFiller> TagProbeFitTreeProducer::mcUnbiasFiller_ [private] |
The object that actually computes variables and fills the tree for unbiased MC.
Definition at line 84 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
std::set<int32_t> TagProbeFitTreeProducer::motherPdgId_ [private] |
Possible pdgids for the mother. If empty, any truth-matched mu will be considered good.
Definition at line 67 of file TagProbeFitTreeProducer.cc.
Referenced by checkMother(), and TagProbeFitTreeProducer().
std::auto_ptr<tnp::BaseTreeFiller> TagProbeFitTreeProducer::oldTagFiller_ [private] |
Definition at line 85 of file TagProbeFitTreeProducer.cc.
Referenced by analyze().
std::auto_ptr<tnp::BaseTreeFiller> TagProbeFitTreeProducer::pairFiller_ [private] |
Definition at line 87 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
Definition at line 65 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
std::auto_ptr<tnp::BaseTreeFiller> TagProbeFitTreeProducer::tagFiller_ [private] |
Definition at line 86 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
InputTag to an edm::Association<reco::GenParticle> from tags & probes to MC truth.
Definition at line 65 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), and TagProbeFitTreeProducer().
The object that produces pairs of tags and probes, making any arbitration needed.
Definition at line 80 of file TagProbeFitTreeProducer.cc.
Referenced by analyze().
std::auto_ptr<tnp::TPTreeFiller> TagProbeFitTreeProducer::treeFiller_ [private] |
The object that actually computes variables and fills the tree for T&P.
Definition at line 82 of file TagProbeFitTreeProducer.cc.
Referenced by analyze(), endJob(), and TagProbeFitTreeProducer().