CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_4_1_8_patch12/src/RecoTauTag/TauTagTools/plugins/RecoTauGenJetMatchSelector.cc

Go to the documentation of this file.
00001 /*
00002  * AssociationMatchRefSelector
00003  *
00004  * EDFilter template to select objects from a View<InputType> that have valid
00005  * matches in the given Association<MatchedType> matching, and produce a
00006  * RefToBaseVector pointing to the matched objects.
00007  *
00008  * Author: Evan K. Friis (UC Davis)
00009  *
00010  */
00011 
00012 #include "FWCore/Framework/interface/Frameworkfwd.h"
00013 #include "FWCore/Framework/interface/EDFilter.h"
00014 #include "FWCore/Framework/interface/Event.h"
00015 #include "FWCore/Utilities/interface/InputTag.h"
00016 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00017 #include "DataFormats/Common/interface/Handle.h"
00018 #include "DataFormats/Common/interface/Association.h"
00019 #include "DataFormats/Common/interface/RefToBaseVector.h"
00020 
00021 namespace reco { namespace tau {
00022 
00023 template<typename InputType, typename MatchedType>
00024 class AssociationMatchRefSelector : public edm::EDFilter {
00025   public:
00026     typedef typename edm::RefToBaseVector<InputType> OutputType;
00027     typedef typename edm::Association<MatchedType> AssocType;
00028     explicit AssociationMatchRefSelector(const edm::ParameterSet &pset) {
00029       src_ = pset.getParameter<edm::InputTag>("src");
00030       matching_ = pset.getParameter<edm::InputTag>("matching");
00031       filter_ = pset.getParameter<bool>("filter");
00032       produces<OutputType>();
00033     }
00034     bool filter(edm::Event& evt, const edm::EventSetup &es) {
00035       edm::Handle<edm::View<InputType> > input;
00036       evt.getByLabel(src_, input);
00037       edm::Handle<AssocType> match;
00038       evt.getByLabel(matching_, match);
00039       std::auto_ptr<OutputType> output(new OutputType);
00040       for (size_t i = 0; i < input->size(); ++i) {
00041         typename AssocType::reference_type matched = (*match)[input->refAt(i)];
00042         // Check if matched
00043         if (matched.isNonnull()) {
00044           output->push_back(input->refAt(i));
00045         }
00046       }
00047       bool notEmpty = output->size();
00048       evt.put(output);
00049       // Filter if no events passed
00050       return ( !filter_ || notEmpty );
00051     }
00052   private:
00053     edm::InputTag src_;
00054     edm::InputTag matching_;
00055     bool filter_;
00056 };
00057 }}  // end reco::tau
00058 
00059 #include "FWCore/Framework/interface/MakerMacros.h"
00060 #include "DataFormats/JetReco/interface/GenJetCollection.h"
00061 #include "DataFormats/JetReco/interface/PFJet.h"
00062 #include "DataFormats/TauReco/interface/PFTau.h"
00063 
00064 typedef reco::tau::AssociationMatchRefSelector<reco::Candidate,
00065           reco::GenJetCollection>  CandViewGenJetMatchRefSelector;
00066 
00067 
00068 DEFINE_FWK_MODULE(CandViewGenJetMatchRefSelector);