File indexing completed on 2023-03-17 11:21:50
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include "FWCore/Framework/interface/stream/EDProducer.h"
0015 #include "FWCore/Framework/interface/EventSetup.h"
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0018
0019 #include "FWCore/Utilities/interface/Exception.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021
0022 #include <FWCore/ParameterSet/interface/ConfigurationDescriptions.h>
0023 #include <FWCore/ParameterSet/interface/ParameterSetDescription.h>
0024
0025 #include "RecoTauTag/RecoTau/interface/PFRecoTauChargedHadronPlugins.h"
0026 #include "RecoTauTag/RecoTau/interface/RecoTauCleaningTools.h"
0027 #include "RecoTauTag/RecoTau/interface/RecoTauCommonUtilities.h"
0028 #include "RecoTauTag/RecoTau/interface/RecoTauQualityCuts.h"
0029 #include "RecoTauTag/RecoTau/interface/pfRecoTauChargedHadronAuxFunctions.h"
0030
0031 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0032 #include "DataFormats/TauReco/interface/PFJetChargedHadronAssociation.h"
0033 #include "DataFormats/TauReco/interface/PFRecoTauChargedHadron.h"
0034 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0035 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0036 #include "DataFormats/TrackReco/interface/Track.h"
0037 #include "DataFormats/MuonReco/interface/Muon.h"
0038 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0039 #include "DataFormats/Common/interface/Association.h"
0040 #include "DataFormats/Math/interface/deltaR.h"
0041
0042 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0043
0044 #include <algorithm>
0045 #include <cmath>
0046 #include <functional>
0047 #include <list>
0048 #include <memory>
0049 #include <set>
0050 #include <string>
0051 #include <vector>
0052
0053 class PFRecoTauChargedHadronProducer : public edm::stream::EDProducer<> {
0054 public:
0055 typedef reco::tau::PFRecoTauChargedHadronBuilderPlugin Builder;
0056 typedef reco::tau::PFRecoTauChargedHadronQualityPlugin Ranker;
0057
0058 explicit PFRecoTauChargedHadronProducer(const edm::ParameterSet& cfg);
0059 ~PFRecoTauChargedHadronProducer() override {}
0060 void produce(edm::Event& evt, const edm::EventSetup& es) override;
0061 template <typename T>
0062 void print(const T& chargedHadron);
0063
0064 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0065
0066 private:
0067 typedef std::vector<std::unique_ptr<Ranker>> RankerList;
0068 typedef Builder::return_type ChargedHadronVector;
0069 typedef std::list<std::unique_ptr<reco::PFRecoTauChargedHadron>> ChargedHadronList;
0070
0071 typedef reco::tau::RecoTauLexicographicalRanking<RankerList, reco::PFRecoTauChargedHadron> ChargedHadronPredicate;
0072
0073 std::string moduleLabel_;
0074
0075
0076 edm::InputTag srcJets_;
0077 edm::EDGetTokenT<reco::JetView> Jets_token;
0078 double minJetPt_;
0079 double maxJetAbsEta_;
0080
0081
0082 std::vector<std::unique_ptr<Builder>> builders_;
0083 RankerList rankers_;
0084
0085 std::unique_ptr<ChargedHadronPredicate> predicate_;
0086
0087
0088 std::unique_ptr<StringCutObjectSelector<reco::PFRecoTauChargedHadron>> outputSelector_;
0089
0090
0091 int verbosity_;
0092 };
0093
0094 PFRecoTauChargedHadronProducer::PFRecoTauChargedHadronProducer(const edm::ParameterSet& cfg)
0095 : moduleLabel_(cfg.getParameter<std::string>("@module_label")) {
0096 srcJets_ = cfg.getParameter<edm::InputTag>("jetSrc");
0097 Jets_token = consumes<reco::JetView>(srcJets_);
0098 minJetPt_ = cfg.getParameter<double>("minJetPt");
0099 maxJetAbsEta_ = cfg.getParameter<double>("maxJetAbsEta");
0100 verbosity_ = cfg.getParameter<int>("verbosity");
0101
0102
0103 edm::VParameterSet psets_builders = cfg.getParameter<edm::VParameterSet>("builders");
0104 for (edm::VParameterSet::const_iterator pset = psets_builders.begin(); pset != psets_builders.end(); ++pset) {
0105 std::string pluginType = pset->getParameter<std::string>("plugin");
0106 edm::ParameterSet pset_modified = (*pset);
0107 pset_modified.addParameter<int>("verbosity", verbosity_);
0108 builders_.emplace_back(
0109 PFRecoTauChargedHadronBuilderPluginFactory::get()->create(pluginType, pset_modified, consumesCollector()));
0110 }
0111
0112
0113 edm::VParameterSet psets_rankers = cfg.getParameter<edm::VParameterSet>("ranking");
0114 for (edm::VParameterSet::const_iterator pset = psets_rankers.begin(); pset != psets_rankers.end(); ++pset) {
0115 std::string pluginType = pset->getParameter<std::string>("plugin");
0116 edm::ParameterSet pset_modified = (*pset);
0117 pset_modified.addParameter<int>("verbosity", verbosity_);
0118 rankers_.emplace_back(PFRecoTauChargedHadronQualityPluginFactory::get()->create(pluginType, pset_modified));
0119 }
0120
0121
0122 predicate_ = std::make_unique<ChargedHadronPredicate>(rankers_);
0123
0124
0125 std::string selection = cfg.getParameter<std::string>("outputSelection");
0126 if (!selection.empty()) {
0127 outputSelector_ = std::make_unique<StringCutObjectSelector<reco::PFRecoTauChargedHadron>>(selection);
0128 }
0129
0130 produces<reco::PFJetChargedHadronAssociation>();
0131 }
0132
0133 void PFRecoTauChargedHadronProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0134 if (verbosity_) {
0135 edm::LogPrint("PFRecoTauChHProducer") << "<PFRecoTauChargedHadronProducer::produce>:";
0136 edm::LogPrint("PFRecoTauChHProducer") << " moduleLabel = " << moduleLabel_;
0137 }
0138
0139
0140 for (auto& builder : builders_) {
0141 builder->setup(evt, es);
0142 }
0143
0144
0145 edm::Handle<reco::JetView> jets;
0146 evt.getByToken(Jets_token, jets);
0147
0148
0149 edm::RefToBaseVector<reco::Jet> pfJets;
0150 size_t nElements = jets->size();
0151 for (size_t i = 0; i < nElements; ++i) {
0152 pfJets.push_back(jets->refAt(i));
0153 }
0154
0155
0156 std::unique_ptr<reco::PFJetChargedHadronAssociation> pfJetChargedHadronAssociations;
0157
0158 if (!pfJets.empty()) {
0159 pfJetChargedHadronAssociations = std::make_unique<reco::PFJetChargedHadronAssociation>(reco::JetRefBaseProd(jets));
0160 } else {
0161 pfJetChargedHadronAssociations = std::make_unique<reco::PFJetChargedHadronAssociation>();
0162 }
0163
0164
0165 for (const auto& pfJet : pfJets) {
0166 if (pfJet->pt() - minJetPt_ < 1e-5)
0167 continue;
0168 if (std::abs(pfJet->eta()) - maxJetAbsEta_ > -1e-5)
0169 continue;
0170
0171
0172 ChargedHadronList uncleanedChargedHadrons;
0173
0174
0175 for (auto const& builder : builders_) {
0176 try {
0177 ChargedHadronVector result((*builder)(*pfJet));
0178 if (verbosity_) {
0179 edm::LogPrint("PFRecoTauChHProducer") << "result of builder = " << builder->name() << ":";
0180 for (auto const& had : result) {
0181 print(*had);
0182 }
0183 }
0184 std::move(result.begin(), result.end(), std::back_inserter(uncleanedChargedHadrons));
0185 } catch (cms::Exception& exception) {
0186 edm::LogError("BuilderPluginException")
0187 << "Exception caught in builder plugin " << builder->name() << ", rethrowing" << std::endl;
0188 throw exception;
0189 }
0190 }
0191
0192
0193 uncleanedChargedHadrons.sort([this](const auto& a, const auto& b) { return (*predicate_)(*a, *b); });
0194
0195
0196 std::vector<reco::PFRecoTauChargedHadron> cleanedChargedHadrons;
0197
0198
0199 typedef std::pair<double, double> etaPhiPair;
0200 std::list<etaPhiPair> tracksInCleanCollection;
0201 std::set<reco::CandidatePtr> neutralPFCandsInCleanCollection;
0202
0203 while (!uncleanedChargedHadrons.empty()) {
0204
0205 std::unique_ptr<reco::PFRecoTauChargedHadron> nextChargedHadron = std::move(uncleanedChargedHadrons.front());
0206 uncleanedChargedHadrons.pop_front();
0207 if (verbosity_) {
0208 edm::LogPrint("PFRecoTauChHProducer") << "processing nextChargedHadron:";
0209 edm::LogPrint("PFRecoTauChHProducer") << (*nextChargedHadron);
0210 }
0211
0212
0213 if (!(*outputSelector_)(*nextChargedHadron))
0214 continue;
0215
0216 const reco::Track* track = reco::tau::getTrackFromChargedHadron(*nextChargedHadron);
0217
0218
0219 bool isTrack_overlap = false;
0220 if (track) {
0221 double track_eta = track->eta();
0222 double track_phi = track->phi();
0223 for (std::list<etaPhiPair>::const_iterator trackInCleanCollection = tracksInCleanCollection.begin();
0224 trackInCleanCollection != tracksInCleanCollection.end();
0225 ++trackInCleanCollection) {
0226 double dR = deltaR(track_eta, track_phi, trackInCleanCollection->first, trackInCleanCollection->second);
0227 if (dR < 1.e-4)
0228 isTrack_overlap = true;
0229 }
0230 }
0231 if (verbosity_) {
0232 edm::LogPrint("PFRecoTauChHProducer") << "isTrack_overlap = " << isTrack_overlap;
0233 }
0234 if (isTrack_overlap)
0235 continue;
0236
0237
0238 bool isNeutralPFCand_overlap = false;
0239 if (nextChargedHadron->algoIs(reco::PFRecoTauChargedHadron::kPFNeutralHadron)) {
0240 for (std::set<reco::CandidatePtr>::const_iterator neutralPFCandInCleanCollection =
0241 neutralPFCandsInCleanCollection.begin();
0242 neutralPFCandInCleanCollection != neutralPFCandsInCleanCollection.end();
0243 ++neutralPFCandInCleanCollection) {
0244 if ((*neutralPFCandInCleanCollection) == nextChargedHadron->getChargedPFCandidate())
0245 isNeutralPFCand_overlap = true;
0246 }
0247 }
0248 if (verbosity_) {
0249 edm::LogPrint("PFRecoTauChHProducer") << "isNeutralPFCand_overlap = " << isNeutralPFCand_overlap;
0250 }
0251 if (isNeutralPFCand_overlap)
0252 continue;
0253
0254
0255 std::vector<reco::CandidatePtr> uniqueNeutralPFCands;
0256 std::set_difference(nextChargedHadron->getNeutralPFCandidates().begin(),
0257 nextChargedHadron->getNeutralPFCandidates().end(),
0258 neutralPFCandsInCleanCollection.begin(),
0259 neutralPFCandsInCleanCollection.end(),
0260 std::back_inserter(uniqueNeutralPFCands));
0261
0262 if (uniqueNeutralPFCands.size() ==
0263 nextChargedHadron->getNeutralPFCandidates()
0264 .size()) {
0265 if (track)
0266 tracksInCleanCollection.push_back(std::make_pair(track->eta(), track->phi()));
0267 neutralPFCandsInCleanCollection.insert(nextChargedHadron->getNeutralPFCandidates().begin(),
0268 nextChargedHadron->getNeutralPFCandidates().end());
0269 if (verbosity_) {
0270 edm::LogPrint("PFRecoTauChHProducer") << "--> adding nextChargedHadron to output collection.";
0271 }
0272 cleanedChargedHadrons.push_back(*nextChargedHadron);
0273 } else {
0274 nextChargedHadron->neutralPFCandidates_.clear();
0275 for (auto const& neutralPFCand : uniqueNeutralPFCands) {
0276 nextChargedHadron->neutralPFCandidates_.push_back(neutralPFCand);
0277 }
0278
0279 reco::tau::setChargedHadronP4(*nextChargedHadron);
0280
0281
0282 auto insertionPoint = std::lower_bound(uncleanedChargedHadrons.begin(),
0283 uncleanedChargedHadrons.end(),
0284 *nextChargedHadron,
0285 [this](const auto& a, const auto& b) { return (*predicate_)(*a, b); });
0286 if (verbosity_) {
0287 edm::LogPrint("PFRecoTauChHProducer") << "--> removing non-unique neutral PFCandidates and reinserting "
0288 "nextChargedHadron in uncleaned collection.";
0289 }
0290 uncleanedChargedHadrons.insert(insertionPoint, std::move(nextChargedHadron));
0291 }
0292 }
0293
0294 if (verbosity_) {
0295 for (auto const& had : cleanedChargedHadrons) {
0296 print(had);
0297 }
0298 }
0299
0300
0301 pfJetChargedHadronAssociations->setValue(pfJet.key(), cleanedChargedHadrons);
0302 }
0303
0304 evt.put(std::move(pfJetChargedHadronAssociations));
0305 }
0306
0307 template <typename T>
0308 void PFRecoTauChargedHadronProducer::print(const T& chargedHadron) {
0309 edm::LogPrint("PFRecoTauChHProducer") << chargedHadron;
0310 edm::LogPrint("PFRecoTauChHProducer") << "Rankers:";
0311 for (auto const& ranker : rankers_) {
0312 constexpr unsigned width = 25;
0313 edm::LogPrint("PFRecoTauChHProducer")
0314 << " " << std::setiosflags(std::ios::left) << std::setw(width) << ranker->name() << " "
0315 << std::resetiosflags(std::ios::left) << std::setprecision(3) << (*ranker)(chargedHadron) << std::endl;
0316 }
0317 }
0318
0319 void PFRecoTauChargedHadronProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0320
0321 edm::ParameterSetDescription desc;
0322 {
0323 edm::ParameterSetDescription desc_ranking;
0324 desc_ranking.add<std::string>("selectionPassFunction", "-pt");
0325 desc_ranking.add<double>("selectionFailValue", 1000.0);
0326 desc_ranking.add<std::string>("selection", "algoIs(\"kChargedPFCandidate\")");
0327 desc_ranking.add<std::string>("name", "ChargedPFCandidate");
0328 desc_ranking.add<std::string>("plugin", "PFRecoTauChargedHadronStringQuality");
0329
0330 edm::ParameterSet pset_ranking;
0331 pset_ranking.addParameter<std::string>("selectionPassFunction", "-pt");
0332 pset_ranking.addParameter<double>("selectionFailValue", 1000.0);
0333 pset_ranking.addParameter<std::string>("selection", "algoIs(\"kChargedPFCandidate\")");
0334 pset_ranking.addParameter<std::string>("name", "ChargedPFCandidate");
0335 pset_ranking.addParameter<std::string>("plugin", "PFRecoTauChargedHadronStringQuality");
0336 std::vector<edm::ParameterSet> vpsd_ranking;
0337 vpsd_ranking.push_back(pset_ranking);
0338
0339 desc.addVPSet("ranking", desc_ranking, vpsd_ranking);
0340 }
0341
0342 desc.add<int>("verbosity", 0);
0343 desc.add<double>("maxJetAbsEta", 2.5);
0344 desc.add<std::string>("outputSelection", "pt > 0.5");
0345 desc.add<double>("minJetPt", 14.0);
0346 desc.add<edm::InputTag>("jetSrc", edm::InputTag("ak4PFJets"));
0347
0348 {
0349 edm::ParameterSetDescription desc_builders;
0350 desc_builders.add<double>("minMergeChargedHadronPt");
0351 desc_builders.add<std::string>("name");
0352 desc_builders.add<std::string>("plugin");
0353 desc_builders.addOptional<double>("dRcone");
0354 desc_builders.addOptional<bool>("dRconeLimitedToJetArea");
0355 desc_builders.addOptional<double>("dRmergeNeutralHadron");
0356 desc_builders.addOptional<double>("dRmergePhoton");
0357 desc_builders.addOptional<edm::InputTag>("srcTracks");
0358
0359 edm::ParameterSetDescription desc_qualityCuts;
0360 reco::tau::RecoTauQualityCuts::fillDescriptions(desc_qualityCuts);
0361 desc_builders.add<edm::ParameterSetDescription>("qualityCuts", desc_qualityCuts);
0362
0363 desc_builders.add<double>("minMergeGammaEt");
0364 desc_builders.add<int>("verbosity", 0);
0365 desc_builders.add<double>("minMergeNeutralHadronEt");
0366
0367 desc_builders.addOptional<double>("dRmergePhotonWrtChargedHadron");
0368 desc_builders.addOptional<double>("dRmergePhotonWrtNeutralHadron");
0369 desc_builders.addOptional<int>("maxUnmatchedBlockElementsNeutralHadron");
0370 desc_builders.addOptional<double>("dRmergePhotonWrtElectron");
0371 desc_builders.addOptional<std::vector<int>>("chargedHadronCandidatesParticleIds");
0372 desc_builders.addOptional<int>("minBlockElementMatchesPhoton");
0373 desc_builders.addOptional<double>("dRmergeNeutralHadronWrtNeutralHadron");
0374 desc_builders.addOptional<int>("maxUnmatchedBlockElementsPhoton");
0375 desc_builders.addOptional<double>("dRmergeNeutralHadronWrtOther");
0376 desc_builders.addOptional<double>("dRmergeNeutralHadronWrtElectron");
0377 desc_builders.addOptional<int>("minBlockElementMatchesNeutralHadron");
0378 desc_builders.addOptional<double>("dRmergePhotonWrtOther");
0379 desc_builders.addOptional<double>("dRmergeNeutralHadronWrtChargedHadron");
0380
0381 edm::ParameterSet pset_builders;
0382 pset_builders.addParameter<std::string>("name", "");
0383 pset_builders.addParameter<std::string>("plugin", "");
0384 edm::ParameterSet qualityCuts;
0385 pset_builders.addParameter<edm::ParameterSet>("qualityCuts", qualityCuts);
0386 pset_builders.addParameter<int>("verbosity", 0);
0387 std::vector<edm::ParameterSet> vpsd_builders;
0388 vpsd_builders.push_back(pset_builders);
0389
0390 desc.addVPSet("builders", desc_builders, vpsd_builders);
0391 }
0392
0393 descriptions.add("pfRecoTauChargedHadronProducer", desc);
0394 }
0395
0396 #include "FWCore/Framework/interface/MakerMacros.h"
0397
0398 DEFINE_FWK_MODULE(PFRecoTauChargedHadronProducer);