CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_3/src/RecoTauTag/TauTagTools/plugins/RecoTauDecayModeTruthMatchPlugin.cc

Go to the documentation of this file.
00001 /*
00002  * RecoTauDecayModeTruthMatchPlugin
00003  *
00004  * Author: Evan K. Friis, UC Davis
00005  *
00006  * Implements a RecoTauCleaner plugin that returns the difference
00007  * between the reconstructed decay mode and true decay mode index.
00008  *
00009  * By requiring the return value to be zero one can select reco taus
00010  * that have the decay mode correctly reconstructed.
00011  *
00012  * $Id. $
00013  */
00014 
00015 #include "RecoTauTag/RecoTau/interface/RecoTauBuilderPlugins.h"
00016 
00017 #include "DataFormats/TauReco/interface/PFTau.h"
00018 #include "DataFormats/TauReco/interface/PFTauFwd.h"
00019 #include "DataFormats/JetReco/interface/GenJet.h"
00020 #include "DataFormats/JetReco/interface/GenJetCollection.h"
00021 #include "DataFormats/Common/interface/Association.h"
00022 
00023 #include "RecoTauTag/RecoTau/interface/PFTauDecayModeTools.h"
00024 
00025 #include <boost/foreach.hpp>
00026 
00027 
00028 namespace tautools {
00029 
00030 class RecoTauDecayModeTruthMatchPlugin : public reco::tau::RecoTauCleanerPlugin
00031 {
00032   public:
00033     explicit RecoTauDecayModeTruthMatchPlugin(const edm::ParameterSet& pset);
00034     virtual ~RecoTauDecayModeTruthMatchPlugin() {}
00035     double operator()(const reco::PFTauRef&) const;
00036     void beginEvent();
00037 
00038   private:
00039     edm::InputTag matchingSrc_;
00040     typedef edm::Association<reco::GenJetCollection> GenJetAssociation;
00041     edm::Handle<GenJetAssociation> genTauMatch_;
00042 };
00043 
00044 // ctor
00045 RecoTauDecayModeTruthMatchPlugin::RecoTauDecayModeTruthMatchPlugin(
00046     const edm::ParameterSet& pset): RecoTauCleanerPlugin(pset),
00047   matchingSrc_(pset.getParameter<edm::InputTag>("matching")) {}
00048 
00049 // Called by base class at the beginning of each event
00050 void RecoTauDecayModeTruthMatchPlugin::beginEvent() {
00051   // Load the matching information
00052   evt()->getByLabel(matchingSrc_, genTauMatch_);
00053 }
00054 
00055 // Determine a number giving the quality of the input tau.  Lower numbers are
00056 // better - zero indicates that the reco decay mode matches the truth.
00057 double RecoTauDecayModeTruthMatchPlugin::operator()(const reco::PFTauRef& tau)
00058   const {
00059   GenJetAssociation::reference_type truth = (*genTauMatch_)[tau];
00060   // Check if the matching exists, if not return +infinity
00061   if (truth.isNull())
00062     return std::numeric_limits<double>::infinity();
00063   // Get the difference in decay mode.  The closer to zero, the more the decay
00064   // mode is matched.
00065   return std::abs(
00066       reco::tau::getDecayMode(truth.get()) - tau->decayMode());
00067 }
00068 
00069 } // end tautools namespace
00070 
00071 #include "FWCore/Framework/interface/MakerMacros.h"
00072 DEFINE_EDM_PLUGIN(RecoTauCleanerPluginFactory, tautools::RecoTauDecayModeTruthMatchPlugin, "RecoTauDecayModeTruthMatchPlugin");