CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/PhysicsTools/TagAndProbe/plugins/TagProbeFitTreeProducer.cc

Go to the documentation of this file.
00001 // -*- C++ -*-
00002 //
00003 // Package:    TagProbeFitTreeProducer
00004 // Class:      TagProbeFitTreeProducer
00005 // 
00013 //
00014 // Original Author:  Sep 15 09:45
00015 //         Created:  Mon Sep 15 09:49:08 CEST 2008
00016 // $Id: TagProbeFitTreeProducer.cc,v 1.7 2010/05/02 15:18:25 gpetrucc Exp $
00017 //
00018 //
00019 
00020 
00021 // system include files
00022 #include <memory>
00023 #include <ctype.h>
00024 
00025 // user include files
00026 #include "FWCore/Framework/interface/Frameworkfwd.h"
00027 #include "FWCore/Framework/interface/EDAnalyzer.h"
00028 
00029 #include "FWCore/Framework/interface/Event.h"
00030 #include "FWCore/Framework/interface/MakerMacros.h"
00031 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00032 #include "FWCore/Utilities/interface/InputTag.h"
00033 
00034 #include "DataFormats/Candidate/interface/CandidateFwd.h"
00035 #include "DataFormats/Candidate/interface/Candidate.h"
00036 
00037 #include "DataFormats/Common/interface/Association.h"
00038 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00039 
00040 #include "PhysicsTools/TagAndProbe/interface/TPTreeFiller.h"
00041 #include "PhysicsTools/TagAndProbe/interface/TagProbePairMaker.h"
00042 
00043 #include <set>
00044 
00045 #include "FWCore/ParameterSet/interface/Registry.h"
00046 
00047 //
00048 // class decleration
00049 //
00050 
00051 class TagProbeFitTreeProducer : public edm::EDAnalyzer {
00052    public:
00053       explicit TagProbeFitTreeProducer(const edm::ParameterSet&);
00054       ~TagProbeFitTreeProducer();
00055 
00056 
00057    private:
00058       virtual void analyze(const edm::Event&, const edm::EventSetup&);
00059       virtual void endJob() ;
00060 
00061       //---- MC truth information
00063       bool isMC_; 
00065       edm::InputTag tagMatches_, probeMatches_;
00067       std::set<int32_t> motherPdgId_;
00069       bool checkMother(const reco::GenParticleRef &ref) const ;
00070 
00071       //---- Unbiased MC truth information
00073       bool makeMCUnbiasTree_;
00075       bool checkMotherInUnbiasEff_;
00077       edm::InputTag allProbes_;
00078 
00080       tnp::TagProbePairMaker tagProbePairMaker_;
00082       std::auto_ptr<tnp::TPTreeFiller> treeFiller_; 
00084       std::auto_ptr<tnp::BaseTreeFiller> mcUnbiasFiller_;
00085       std::auto_ptr<tnp::BaseTreeFiller> oldTagFiller_;
00086       std::auto_ptr<tnp::BaseTreeFiller> tagFiller_;
00087       std::auto_ptr<tnp::BaseTreeFiller> pairFiller_;
00088       std::auto_ptr<tnp::BaseTreeFiller> mcFiller_;
00089 };
00090 
00091 //
00092 // constructors and destructor
00093 //
00094 TagProbeFitTreeProducer::TagProbeFitTreeProducer(const edm::ParameterSet& iConfig) :
00095     isMC_(iConfig.getParameter<bool>("isMC")),
00096     makeMCUnbiasTree_(isMC_ ? iConfig.getParameter<bool>("makeMCUnbiasTree") : false),
00097     checkMotherInUnbiasEff_(makeMCUnbiasTree_ ? iConfig.getParameter<bool>("checkMotherInUnbiasEff") : false),
00098     tagProbePairMaker_(iConfig),
00099     treeFiller_(new tnp::TPTreeFiller(iConfig)),
00100     oldTagFiller_((iConfig.existsAs<bool>("fillTagTree") && iConfig.getParameter<bool>("fillTagTree")) ? new tnp::BaseTreeFiller("tag_tree",iConfig) : 0)
00101 {
00102     if (isMC_) { 
00103         // For mc efficiency we need the MC matches for tags & probes
00104         tagMatches_ = iConfig.getParameter<edm::InputTag>("tagMatches");
00105         probeMatches_ = iConfig.getParameter<edm::InputTag>("probeMatches");
00106         //.. and the pdgids of the possible mothers
00107         if (iConfig.existsAs<int32_t>("motherPdgId")) {
00108             motherPdgId_.insert(iConfig.getParameter<int32_t>("motherPdgId"));
00109         } else {
00110             std::vector<int32_t> motherIds = iConfig.getParameter<std::vector<int32_t> >("motherPdgId");
00111             motherPdgId_.insert(motherIds.begin(), motherIds.end());
00112         }
00113 
00114         // For unbiased efficiency we also need the collection of all probes
00115         if (makeMCUnbiasTree_) {
00116             allProbes_ = iConfig.getParameter<edm::InputTag>("allProbes");
00117             mcUnbiasFiller_.reset(new tnp::BaseTreeFiller("mcUnbias_tree",iConfig));
00118         }
00119     }
00120 
00121     edm::ParameterSet tagPSet;
00122     if (iConfig.existsAs<edm::ParameterSet>("tagVariables")) tagPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("tagVariables"));
00123     if (iConfig.existsAs<edm::ParameterSet>("tagFlags"    )) tagPSet.addParameter<edm::ParameterSet>("flags",     iConfig.getParameter<edm::ParameterSet>("tagFlags"));
00124     if (!tagPSet.empty()) { tagFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, tagPSet, "tag_")); }
00125     edm::ParameterSet mcPSet;
00126     if (iConfig.existsAs<edm::ParameterSet>("mcVariables")) mcPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("mcVariables"));
00127     if (iConfig.existsAs<edm::ParameterSet>("mcFlags"    )) mcPSet.addParameter<edm::ParameterSet>("flags",     iConfig.getParameter<edm::ParameterSet>("mcFlags"));
00128     if (!mcPSet.empty()) { mcFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, mcPSet, "mc_")); }
00129     edm::ParameterSet pairPSet;
00130     if (iConfig.existsAs<edm::ParameterSet>("pairVariables")) pairPSet.addParameter<edm::ParameterSet>("variables", iConfig.getParameter<edm::ParameterSet>("pairVariables"));
00131     if (iConfig.existsAs<edm::ParameterSet>("pairFlags"    )) pairPSet.addParameter<edm::ParameterSet>("flags",     iConfig.getParameter<edm::ParameterSet>("pairFlags"));
00132     if (!pairPSet.empty()) { pairFiller_.reset(new tnp::BaseTreeFiller(*treeFiller_, pairPSet, "pair_")); }
00133 }
00134 
00135 
00136 TagProbeFitTreeProducer::~TagProbeFitTreeProducer()
00137 {
00138 }
00139 
00140 
00141 //
00142 // member functions
00143 //
00144 
00145 // ------------ method called to for each event  ------------
00146 void
00147 TagProbeFitTreeProducer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00148 {
00149     using namespace edm; using namespace std; 
00150     Handle<reco::CandidateView> src, allProbes;
00151     Handle<Association<vector<reco::GenParticle> > > tagMatches, probeMatches;
00152 
00153     treeFiller_->init(iEvent); // read out info from the event if needed (external vars, list of passing probes, ...)
00154     if (oldTagFiller_.get()) oldTagFiller_->init(iEvent);
00155     if (tagFiller_.get())    tagFiller_->init(iEvent);
00156     if (pairFiller_.get())   pairFiller_->init(iEvent);
00157     if (mcFiller_.get())     mcFiller_->init(iEvent);
00158 
00159     // on mc we want to load also the MC match info
00160     if (isMC_) {
00161         iEvent.getByLabel(tagMatches_, tagMatches);
00162         iEvent.getByLabel(probeMatches_, probeMatches);
00163     }
00164 
00165     // get the list of (tag+probe) pairs, performing arbitration 
00166     tnp::TagProbePairs pairs = tagProbePairMaker_.run(iEvent);
00167     // loop on them to fill the tree
00168     for (tnp::TagProbePairs::const_iterator it = pairs.begin(), ed = pairs.end(); it != ed; ++it) {
00169         // on mc, fill mc info (on non-mc, let it to 'true', the treeFiller will ignore it anyway
00170         bool mcTrue = false;
00171         if (isMC_) {
00172             reco::GenParticleRef mtag = (*tagMatches)[it->tag], mprobe = (*probeMatches)[it->probe];
00173             mcTrue = checkMother(mtag) && checkMother(mprobe);
00174             if (mcTrue && mcFiller_.get()) mcFiller_->fill(reco::CandidateBaseRef(mprobe));
00175         }
00176         // fill in the variables for this t+p pair
00177         if (tagFiller_.get())    tagFiller_->fill(it->tag);
00178         if (oldTagFiller_.get()) oldTagFiller_->fill(it->tag);
00179         if (pairFiller_.get())   pairFiller_->fill(it->pair);
00180         treeFiller_->fill(it->probe, it->mass, mcTrue);
00181     } 
00182 
00183     if (isMC_ && makeMCUnbiasTree_) {
00184         // read full collection of probes
00185         iEvent.getByLabel(allProbes_, allProbes);
00186         // init the tree filler
00187         mcUnbiasFiller_->init(iEvent);
00188         // loop on probes
00189         for (size_t i = 0, n = allProbes->size(); i < n; ++i) {
00190             const reco::CandidateBaseRef & probe = allProbes->refAt(i);
00191             // check mc match, and possibly mother match
00192             reco::GenParticleRef probeMatch = (*probeMatches)[probe];
00193             bool probeOk = checkMotherInUnbiasEff_ ? checkMother(probeMatch) : probeMatch.isNonnull();
00194             // fill the tree only for good ones
00195             if (probeOk) mcUnbiasFiller_->fill(probe);
00196         }
00197     }
00198    
00199 }
00200 
00201 bool
00202 TagProbeFitTreeProducer::checkMother(const reco::GenParticleRef &ref) const {
00203     if (ref.isNull()) return false;
00204     if (motherPdgId_.empty()) return true;
00205     if (motherPdgId_.find(abs(ref->pdgId())) != motherPdgId_.end()) return true;
00206     reco::GenParticle::mothers m = ref->motherRefVector();
00207     for (reco::GenParticle::mothers::const_iterator it = m.begin(), e = m.end(); it != e; ++it) {
00208         if (checkMother(*it)) return true;
00209     }
00210     return false;
00211 }
00212 
00213 
00214 // ------------ method called once each job just after ending the event loop  ------------
00215 void 
00216 TagProbeFitTreeProducer::endJob() {
00217     // ask to write the current PSet info into the TTree header
00218     treeFiller_->writeProvenance(edm::getProcessParameterSet());
00219 }
00220 
00221 //define this as a plug-in
00222 DEFINE_FWK_MODULE(TagProbeFitTreeProducer);