RecoTau/plugins/RecoTauProducer.cc

0001 /*
0002  * RecoTauProducer
0003  *
0004  * Interface between the various tau algorithms and the edm::Event.  The
0005  * RecoTauProducer takes as data input is a collection (view) of reco::PFJets,
0006  * and Jet-PiZero assoications that give the reco::RecoTauPiZeros for those
0007  * jets.  The actual building of taus is done by the list of builders - each of
0008  * which constructs a PFTau for each PFJet.  The output collection may have
0009  * multiple taus for each PFJet - these overlaps are to be resolved by the
0010  * RecoTauCleaner module.
0011  *
0012  * Additionally, there are "modifier" plugins, which can do things like add the
0013  * lead track significance, or electron rejection variables.
0014  *
0015  * Authors: Evan K. Friis (UC Davis),
0016  *          Christian Veelken (LLR)
0017  *
0018  */
0019
0020 #include <algorithm>
0021 #include <functional>
0022 #include <memory>
0023
0024 #include "FWCore/Framework/interface/stream/EDProducer.h"
0025 #include "FWCore/Framework/interface/EventSetup.h"
0026 #include "FWCore/Framework/interface/Event.h"
0027 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0028 #include <FWCore/ParameterSet/interface/ConfigurationDescriptions.h>
0029 #include <FWCore/ParameterSet/interface/ParameterSetDescription.h>
0030
0031 #include "RecoTauTag/RecoTau/interface/RecoTauBuilderPlugins.h"
0032 #include "RecoTauTag/RecoTau/interface/RecoTauCommonUtilities.h"
0033 #include "RecoTauTag/RecoTau/interface/RecoTauQualityCuts.h"
0034
0035 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0036 #include "DataFormats/TauReco/interface/PFJetChargedHadronAssociation.h"
0037 #include "DataFormats/TauReco/interface/PFRecoTauChargedHadron.h"
0038 #include "DataFormats/TauReco/interface/JetPiZeroAssociation.h"
0039 #include "DataFormats/TauReco/interface/PFTau.h"
0040 #include "DataFormats/JetReco/interface/JetCollection.h"
0041 #include "DataFormats/Common/interface/Association.h"
0042
0043 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0044
0045 class RecoTauProducer : public edm::stream::EDProducer<> {
0046 public:
0047   typedef reco::tau::RecoTauBuilderPlugin Builder;
0048   typedef reco::tau::RecoTauModifierPlugin Modifier;
0049   typedef std::vector<std::unique_ptr<Builder>> BuilderList;
0050   typedef std::vector<std::unique_ptr<Modifier>> ModifierList;
0051
0052   explicit RecoTauProducer(const edm::ParameterSet& pset);
0053   ~RecoTauProducer() override {}
0054   void produce(edm::Event& evt, const edm::EventSetup& es) override;
0055
0056   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0057
0058 private:
0059   edm::InputTag jetSrc_;
0060   edm::InputTag jetRegionSrc_;
0061   edm::InputTag chargedHadronSrc_;
0062   edm::InputTag piZeroSrc_;
0063
0064   double minJetPt_;
0065   double maxJetAbsEta_;
0066   //token definition
0067   edm::EDGetTokenT<reco::JetView> jet_token;
0068   edm::EDGetTokenT<edm::AssociationMap<edm::OneToOne<reco::JetView, reco::JetView>>> jetRegion_token;
0069   edm::EDGetTokenT<reco::PFJetChargedHadronAssociation> chargedHadron_token;
0070   edm::EDGetTokenT<reco::JetPiZeroAssociation> piZero_token;
0071
0072   BuilderList builders_;
0073   ModifierList modifiers_;
0074   // Optional selection on the output of the taus
0075   std::unique_ptr<StringCutObjectSelector<reco::PFTau>> outputSelector_;
0076   // Whether or not to add build a tau from a jet for which the builders
0077   // return no taus.  The tau will have no content, only the four vector of
0078   // the orginal jet.
0079   bool buildNullTaus_;
0080 };
0081
0082 RecoTauProducer::RecoTauProducer(const edm::ParameterSet& pset) {
0083   jetSrc_ = pset.getParameter<edm::InputTag>("jetSrc");
0084   jetRegionSrc_ = pset.getParameter<edm::InputTag>("jetRegionSrc");
0085   chargedHadronSrc_ = pset.getParameter<edm::InputTag>("chargedHadronSrc");
0086   piZeroSrc_ = pset.getParameter<edm::InputTag>("piZeroSrc");
0087
0088   minJetPt_ = pset.getParameter<double>("minJetPt");
0089   maxJetAbsEta_ = pset.getParameter<double>("maxJetAbsEta");
0090   //consumes definition
0091   jet_token = consumes<reco::JetView>(jetSrc_);
0092   jetRegion_token = consumes<edm::AssociationMap<edm::OneToOne<reco::JetView, reco::JetView>>>(jetRegionSrc_);
0093   chargedHadron_token = consumes<reco::PFJetChargedHadronAssociation>(chargedHadronSrc_);
0094   piZero_token = consumes<reco::JetPiZeroAssociation>(piZeroSrc_);
0095
0096   typedef std::vector<edm::ParameterSet> VPSet;
0097   // Get each of our tau builders
0098   const VPSet& builders = pset.getParameter<VPSet>("builders");
0099   for (VPSet::const_iterator builderPSet = builders.begin(); builderPSet != builders.end(); ++builderPSet) {
0100     // Get plugin name
0101     const std::string& pluginType = builderPSet->getParameter<std::string>("plugin");
0102     // Build the plugin
0103     builders_.emplace_back(RecoTauBuilderPluginFactory::get()->create(pluginType, *builderPSet, consumesCollector()));
0104   }
0105
0106   const VPSet& modfiers = pset.getParameter<VPSet>("modifiers");
0107   for (VPSet::const_iterator modfierPSet = modfiers.begin(); modfierPSet != modfiers.end(); ++modfierPSet) {
0108     // Get plugin name
0109     const std::string& pluginType = modfierPSet->getParameter<std::string>("plugin");
0110     // Build the plugin
0111     modifiers_.emplace_back(RecoTauModifierPluginFactory::get()->create(pluginType, *modfierPSet, consumesCollector()));
0112   }
0113
0114   // Check if we want to apply a final output selection
0115   std::string selection = pset.getParameter<std::string>("outputSelection");
0116   if (!selection.empty()) {
0117     outputSelector_ = std::make_unique<StringCutObjectSelector<reco::PFTau>>(selection);
0118   }
0119   buildNullTaus_ = pset.getParameter<bool>("buildNullTaus");
0120
0121   produces<reco::PFTauCollection>();
0122 }
0123
0124 void RecoTauProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0125   // Get the jet input collection via a view of Candidates
0126   edm::Handle<reco::JetView> jetView;
0127   evt.getByToken(jet_token, jetView);
0128
0129   // Get the jet region producer
0130   edm::Handle<edm::AssociationMap<edm::OneToOne<reco::JetView, reco::JetView>>> jetRegionHandle;
0131   evt.getByToken(jetRegion_token, jetRegionHandle);
0132
0133   // Get the charged hadron input collection
0134   edm::Handle<reco::PFJetChargedHadronAssociation> chargedHadronAssoc;
0135   evt.getByToken(chargedHadron_token, chargedHadronAssoc);
0136
0137   // Get the pizero input collection
0138   edm::Handle<reco::JetPiZeroAssociation> piZeroAssoc;
0139   evt.getByToken(piZero_token, piZeroAssoc);
0140
0141   // Update all our builders and modifiers with the event info
0142   for (auto& builder : builders_) {
0143     builder->setup(evt, es);
0144   }
0145   for (auto& modifier : modifiers_) {
0146     modifier->setup(evt, es);
0147   }
0148
0149   // Create output collection
0150   auto output = std::make_unique<reco::PFTauCollection>();
0151   output->reserve(jetView->size());
0152
0153   // Loop over the jets and build the taus for each jet
0154   for (size_t i_j = 0; i_j < jetView->size(); ++i_j) {
0155     const auto& jetRef = jetView->refAt(i_j);
0156     // Get the jet with extra constituents from an area around the jet
0157     if (jetRef->pt() - minJetPt_ < 1e-5)
0158       continue;
0159     if (std::abs(jetRef->eta()) - maxJetAbsEta_ > -1e-5)
0160       continue;
0161     reco::JetBaseRef jetRegionRef = (*jetRegionHandle)[jetRef];
0162     if (jetRegionRef.isNull()) {
0163       throw cms::Exception("BadJetRegionRef") << "No jet region can be found for the current jet: " << jetRef.id();
0164     }
0165     // Remove all the jet constituents from the jet extras
0166     std::vector<reco::CandidatePtr> jetCands = jetRef->daughterPtrVector();
0167     std::vector<reco::CandidatePtr> allRegionalCands = jetRegionRef->daughterPtrVector();
0168     // Sort both by ref key
0169     std::sort(jetCands.begin(), jetCands.end());
0170     std::sort(allRegionalCands.begin(), allRegionalCands.end());
0171     // Get the regional junk candidates not in the jet.
0172     std::vector<reco::CandidatePtr> uniqueRegionalCands;
0173
0174     // This can actually be less than zero, if the jet has really crazy soft
0175     // stuff really far away from the jet axis.
0176     if (allRegionalCands.size() > jetCands.size()) {
0177       uniqueRegionalCands.reserve(allRegionalCands.size() - jetCands.size());
0178     }
0179
0180     // Subtract the jet cands from the regional cands
0181     std::set_difference(allRegionalCands.begin(),
0182                         allRegionalCands.end(),
0183                         jetCands.begin(),
0184                         jetCands.end(),
0185                         std::back_inserter(uniqueRegionalCands));
0186
0187     // Get the charged hadrons associated with this jet
0188     const std::vector<reco::PFRecoTauChargedHadron>& chargedHadrons = (*chargedHadronAssoc)[jetRef];
0189
0190     // Get the pizeros associated with this jet
0191     const std::vector<reco::RecoTauPiZero>& piZeros = (*piZeroAssoc)[jetRef];
0192     // Loop over our builders and create the set of taus for this jet
0193     unsigned int nTausBuilt = 0;
0194     for (const auto& builder : builders_) {
0195       // Get a std::vector of std::unique_ptr to taus from the builder
0196       reco::tau::RecoTauBuilderPlugin::output_type taus(
0197           (*builder)(jetRef, chargedHadrons, piZeros, uniqueRegionalCands));
0198
0199       // Make sure all taus have their jetref set correctly
0200       std::for_each(taus.begin(), taus.end(), [&](auto& arg) { arg->setjetRef(reco::JetBaseRef(jetRef)); });
0201       // Copy without selection
0202       if (!outputSelector_.get()) {
0203         for (auto& tau : taus) {
0204           output->push_back(*tau);
0205         }
0206         nTausBuilt += taus.size();
0207       } else {
0208         // Copy only those that pass the selection.
0209         for (auto const& tau : taus) {
0210           if ((*outputSelector_)(*tau)) {
0211             nTausBuilt++;
0212             output->push_back(*tau);
0213           }
0214         }
0215       }
0216     }
0217     // If we didn't build *any* taus for this jet, build a null tau if desired.
0218     // The null PFTau has no content, but it's four vector is set to that of the
0219     // jet.
0220     if (!nTausBuilt && buildNullTaus_) {
0221       reco::PFTau nullTau(std::numeric_limits<int>::quiet_NaN(), jetRef->p4());
0222       nullTau.setjetRef(reco::JetBaseRef(jetRef));
0223       output->push_back(nullTau);
0224     }
0225   }
0226
0227   // Loop over the taus we have created and apply our modifiers to the taus
0228   for (auto& tau : *output) {
0229     for (const auto& modifier : modifiers_) {
0230       (*modifier)(tau);
0231     }
0232   }
0233
0234   for (auto& modifier : modifiers_) {
0235     modifier->endEvent();
0236   }
0237
0238   evt.put(std::move(output));
0239 }
0240
0241 void RecoTauProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0242   // combinatoricRecoTaus
0243   edm::ParameterSetDescription desc;
0244   desc.add<edm::InputTag>("piZeroSrc", edm::InputTag("ak4PFJetsRecoTauPiZeros"));
0245
0246   edm::ParameterSetDescription desc_qualityCuts;
0247   reco::tau::RecoTauQualityCuts::fillDescriptions(desc_qualityCuts);
0248
0249   {
0250     edm::ParameterSetDescription vpsd_modifiers;
0251     vpsd_modifiers.add<std::string>("name");
0252     vpsd_modifiers.add<std::string>("plugin");
0253     vpsd_modifiers.add<int>("verbosity", 0);
0254
0255     vpsd_modifiers.add<edm::ParameterSetDescription>("qualityCuts", desc_qualityCuts);
0256     vpsd_modifiers.addOptional<edm::InputTag>("ElectronPreIDProducer");
0257     vpsd_modifiers.addOptional<std::string>("DataType");
0258     vpsd_modifiers.addOptional<double>("maximumForElectrionPreIDOutput");
0259     vpsd_modifiers.addOptional<double>("ElecPreIDLeadTkMatch_maxDR");
0260     vpsd_modifiers.addOptional<double>("EcalStripSumE_minClusEnergy");
0261     vpsd_modifiers.addOptional<double>("EcalStripSumE_deltaPhiOverQ_minValue");
0262     vpsd_modifiers.addOptional<double>("EcalStripSumE_deltaPhiOverQ_maxValue");
0263     vpsd_modifiers.addOptional<double>("EcalStripSumE_deltaEta");
0264     vpsd_modifiers.addOptional<double>("dRaddNeutralHadron");
0265     vpsd_modifiers.addOptional<double>("minGammaEt");
0266     vpsd_modifiers.addOptional<double>("dRaddPhoton");
0267     vpsd_modifiers.addOptional<double>("minNeutralHadronEt");
0268     vpsd_modifiers.addOptional<edm::InputTag>("pfTauTagInfoSrc");
0269     vpsd_modifiers.addOptional<edm::InputTag>("trackSrc");
0270
0271     desc.addVPSet("modifiers", vpsd_modifiers);
0272   }
0273
0274   desc.add<edm::InputTag>("jetRegionSrc", edm::InputTag("recoTauAK4PFJets08Region"));
0275   desc.add<double>("maxJetAbsEta", 2.5);
0276   desc.add<std::string>("outputSelection", "leadPFChargedHadrCand().isNonnull()");
0277   desc.add<edm::InputTag>("chargedHadronSrc", edm::InputTag("ak4PFJetsRecoTauChargedHadrons"));
0278   desc.add<double>("minJetPt", 14.0);
0279   desc.add<edm::InputTag>("jetSrc", edm::InputTag("ak4PFJets"));
0280
0281   {
0282     edm::ParameterSetDescription desc_builders;
0283     desc_builders.add<std::string>("name");
0284     desc_builders.add<std::string>("plugin");
0285     desc_builders.add<int>("verbosity", 0);
0286
0287     desc_builders.add<edm::ParameterSetDescription>("qualityCuts", desc_qualityCuts);
0288     {
0289       edm::ParameterSetDescription desc_decayModes;
0290       desc_decayModes.add<unsigned int>("nPiZeros", 0);
0291       desc_decayModes.add<unsigned int>("maxPiZeros", 0);
0292       desc_decayModes.add<unsigned int>("nCharged", 1);
0293       desc_decayModes.add<unsigned int>("maxTracks", 6);
0294       desc_builders.addVPSetOptional("decayModes", desc_decayModes);
0295     }
0296     desc_builders.add<double>("minAbsPhotonSumPt_insideSignalCone", 2.5);
0297     desc_builders.add<double>("minRelPhotonSumPt_insideSignalCone", 0.1);
0298     desc_builders.add<edm::InputTag>("pfCandSrc", edm::InputTag("particleFlow"));
0299
0300     desc_builders.addOptional<std::string>("signalConeSize");
0301     desc_builders.addOptional<double>("isolationConeSize");
0302     desc_builders.addOptional<double>("minAbsPhotonSumPt_outsideSignalCone");
0303     desc_builders.addOptional<double>("minRelPhotonSumPt_outsideSignalCone");
0304     desc_builders.addOptional<std::string>("isoConeChargedHadrons");
0305     desc_builders.addOptional<std::string>("isoConeNeutralHadrons");
0306     desc_builders.addOptional<std::string>("isoConePiZeros");
0307     desc_builders.addOptional<double>("leadObjectPt");
0308     desc_builders.addOptional<std::string>("matchingCone");
0309     desc_builders.addOptional<int>("maxSignalConeChargedHadrons");
0310     desc_builders.addOptional<std::string>("signalConeChargedHadrons");
0311     desc_builders.addOptional<std::string>("signalConeNeutralHadrons");
0312     desc_builders.addOptional<std::string>("signalConePiZeros");
0313     desc_builders.addOptional<bool>("usePFLeptons");
0314
0315     std::vector<edm::ParameterSet> vpset_default;
0316     {
0317       edm::ParameterSet pset_default_builders;
0318       pset_default_builders.addParameter<std::string>("name", "");
0319       pset_default_builders.addParameter<std::string>("plugin", "");
0320       pset_default_builders.addParameter<int>("verbosity", 0);
0321       pset_default_builders.addParameter<double>("minAbsPhotonSumPt_insideSignalCone", 2.5);
0322       pset_default_builders.addParameter<double>("minRelPhotonSumPt_insideSignalCone", 0.1);
0323       pset_default_builders.addParameter<edm::InputTag>("pfCandSrc", edm::InputTag("particleFlow"));
0324       vpset_default.push_back(pset_default_builders);
0325     }
0326     desc.addVPSet("builders", desc_builders, vpset_default);
0327   }
0328
0329   desc.add<bool>("buildNullTaus", false);
0330   desc.add<int>("verbosity", 0);
0331   descriptions.add("combinatoricRecoTaus", desc);
0332 }
0333
0334 #include "FWCore/Framework/interface/MakerMacros.h"
0335 DEFINE_FWK_MODULE(RecoTauProducer);