File indexing completed on 2023-03-17 11:21:55
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014 #include <algorithm>
0015 #include <functional>
0016 #include <memory>
0017
0018 #include "FWCore/Framework/interface/stream/EDProducer.h"
0019 #include "FWCore/Framework/interface/EventSetup.h"
0020 #include "FWCore/Framework/interface/Event.h"
0021 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0022 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0023 #include "FWCore/Framework/interface/ConsumesCollector.h"
0024 #include <FWCore/ParameterSet/interface/ConfigurationDescriptions.h>
0025 #include <FWCore/ParameterSet/interface/ParameterSetDescription.h>
0026
0027 #include "RecoTauTag/RecoTau/interface/RecoTauPiZeroPlugins.h"
0028 #include "RecoTauTag/RecoTau/interface/RecoTauCleaningTools.h"
0029 #include "RecoTauTag/RecoTau/interface/RecoTauCommonUtilities.h"
0030 #include "RecoTauTag/RecoTau/interface/RecoTauQualityCuts.h"
0031
0032 #include "DataFormats/JetReco/interface/PFJetCollection.h"
0033 #include "DataFormats/TauReco/interface/JetPiZeroAssociation.h"
0034 #include "DataFormats/TauReco/interface/RecoTauPiZero.h"
0035 #include "DataFormats/Common/interface/Association.h"
0036
0037 #include "CommonTools/CandUtils/interface/AddFourMomenta.h"
0038 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0039
0040 class RecoTauPiZeroProducer : public edm::stream::EDProducer<> {
0041 public:
0042 typedef reco::tau::RecoTauPiZeroBuilderPlugin Builder;
0043 typedef reco::tau::RecoTauPiZeroQualityPlugin Ranker;
0044
0045 explicit RecoTauPiZeroProducer(const edm::ParameterSet& pset);
0046 ~RecoTauPiZeroProducer() override {}
0047 void produce(edm::Event& evt, const edm::EventSetup& es) override;
0048 void print(const std::vector<reco::RecoTauPiZero>& piZeros, std::ostream& out);
0049
0050 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0051
0052 private:
0053 typedef std::vector<std::unique_ptr<Ranker>> RankerList;
0054 typedef Builder::return_type PiZeroVector;
0055 typedef std::list<std::unique_ptr<reco::RecoTauPiZero>> PiZeroList;
0056
0057 typedef reco::tau::RecoTauLexicographicalRanking<RankerList, reco::RecoTauPiZero> PiZeroPredicate;
0058
0059 std::vector<std::unique_ptr<Builder>> builders_;
0060 RankerList rankers_;
0061 std::unique_ptr<PiZeroPredicate> predicate_;
0062 double piZeroMass_;
0063
0064
0065 std::unique_ptr<StringCutObjectSelector<reco::RecoTauPiZero>> outputSelector_;
0066
0067
0068 edm::EDGetTokenT<reco::JetView> cand_token;
0069
0070 double minJetPt_;
0071 double maxJetAbsEta_;
0072
0073 int verbosity_;
0074 };
0075
0076 RecoTauPiZeroProducer::RecoTauPiZeroProducer(const edm::ParameterSet& pset) {
0077 cand_token = consumes<reco::JetView>(pset.getParameter<edm::InputTag>("jetSrc"));
0078 minJetPt_ = pset.getParameter<double>("minJetPt");
0079 maxJetAbsEta_ = pset.getParameter<double>("maxJetAbsEta");
0080
0081 typedef std::vector<edm::ParameterSet> VPSet;
0082
0083 piZeroMass_ = pset.getParameter<double>("massHypothesis");
0084
0085
0086 const VPSet& builders = pset.getParameter<VPSet>("builders");
0087
0088 for (VPSet::const_iterator builderPSet = builders.begin(); builderPSet != builders.end(); ++builderPSet) {
0089
0090 const std::string& pluginType = builderPSet->getParameter<std::string>("plugin");
0091
0092 builders_.emplace_back(
0093 RecoTauPiZeroBuilderPluginFactory::get()->create(pluginType, *builderPSet, consumesCollector()));
0094 }
0095
0096
0097 const VPSet& rankers = pset.getParameter<VPSet>("ranking");
0098 for (VPSet::const_iterator rankerPSet = rankers.begin(); rankerPSet != rankers.end(); ++rankerPSet) {
0099 const std::string& pluginType = rankerPSet->getParameter<std::string>("plugin");
0100 rankers_.emplace_back(RecoTauPiZeroQualityPluginFactory::get()->create(pluginType, *rankerPSet));
0101 }
0102
0103
0104 predicate_ = std::make_unique<PiZeroPredicate>(rankers_);
0105
0106
0107 std::string selection = pset.getParameter<std::string>("outputSelection");
0108 if (!selection.empty()) {
0109 outputSelector_ = std::make_unique<StringCutObjectSelector<reco::RecoTauPiZero>>(selection);
0110 }
0111
0112 verbosity_ = pset.getParameter<int>("verbosity");
0113
0114 produces<reco::JetPiZeroAssociation>();
0115 }
0116
0117 void RecoTauPiZeroProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0118
0119 edm::Handle<reco::JetView> jetView;
0120 evt.getByToken(cand_token, jetView);
0121
0122
0123 for (auto& builder : builders_) {
0124 builder->setup(evt, es);
0125 }
0126
0127
0128 auto association = std::make_unique<reco::JetPiZeroAssociation>(reco::JetRefBaseProd(jetView));
0129
0130
0131 size_t nJets = jetView->size();
0132 for (size_t i = 0; i < nJets; ++i) {
0133 const reco::JetBaseRef jet(jetView->refAt(i));
0134
0135 if (jet->pt() - minJetPt_ < 1e-5)
0136 continue;
0137 if (std::abs(jet->eta()) - maxJetAbsEta_ > -1e-5)
0138 continue;
0139
0140 PiZeroList dirtyPiZeros;
0141
0142
0143 for (auto const& builder : builders_) {
0144 try {
0145 PiZeroVector result((*builder)(*jet));
0146 std::move(result.begin(), result.end(), std::back_inserter(dirtyPiZeros));
0147 } catch (cms::Exception& exception) {
0148 edm::LogError("BuilderPluginException")
0149 << "Exception caught in builder plugin " << builder->name() << ", rethrowing" << std::endl;
0150 throw exception;
0151 }
0152 }
0153
0154 dirtyPiZeros.sort([this](const auto& a, const auto& b) { return (*predicate_)(*a, *b); });
0155
0156
0157 std::vector<reco::RecoTauPiZero> cleanPiZeros;
0158 std::set<reco::CandidatePtr> photonsInCleanCollection;
0159 while (!dirtyPiZeros.empty()) {
0160
0161 std::unique_ptr<reco::RecoTauPiZero> toAdd(dirtyPiZeros.front().release());
0162 dirtyPiZeros.pop_front();
0163
0164 if (!(*outputSelector_)(*toAdd)) {
0165 continue;
0166 }
0167
0168 std::vector<reco::CandidatePtr> uniqueGammas;
0169 std::set_difference(toAdd->daughterPtrVector().begin(),
0170 toAdd->daughterPtrVector().end(),
0171 photonsInCleanCollection.begin(),
0172 photonsInCleanCollection.end(),
0173 std::back_inserter(uniqueGammas));
0174
0175
0176 if (uniqueGammas.empty()) {
0177 continue;
0178 } else if (uniqueGammas.size() == toAdd->daughterPtrVector().size()) {
0179
0180
0181 photonsInCleanCollection.insert(toAdd->daughterPtrVector().begin(), toAdd->daughterPtrVector().end());
0182 cleanPiZeros.push_back(*toAdd);
0183 } else {
0184
0185
0186 toAdd->clearDaughters();
0187
0188 for (auto const& gamma : uniqueGammas) {
0189 toAdd->addDaughter(gamma);
0190 }
0191
0192 AddFourMomenta p4Builder_;
0193 p4Builder_.set(*toAdd);
0194
0195 auto insertionPoint =
0196 std::lower_bound(dirtyPiZeros.begin(), dirtyPiZeros.end(), *toAdd, [this](const auto& a, const auto& b) {
0197 return (*predicate_)(*a, b);
0198 });
0199 dirtyPiZeros.insert(insertionPoint, std::move(toAdd));
0200 }
0201 }
0202
0203 if (piZeroMass_ >= 0) {
0204 for (auto& cleanPiZero : cleanPiZeros) {
0205 cleanPiZero.setMass(this->piZeroMass_);
0206 };
0207 }
0208
0209 if (verbosity_ >= 2) {
0210 print(cleanPiZeros, std::cout);
0211 }
0212 association->setValue(jet.key(), cleanPiZeros);
0213 }
0214 evt.put(std::move(association));
0215 }
0216
0217
0218 void RecoTauPiZeroProducer::print(const std::vector<reco::RecoTauPiZero>& piZeros, std::ostream& out) {
0219 const unsigned int width = 25;
0220 for (auto const& piZero : piZeros) {
0221 out << piZero;
0222 out << "* Rankers:" << std::endl;
0223 for (auto const& ranker : rankers_) {
0224 out << "* " << std::setiosflags(std::ios::left) << std::setw(width) << ranker->name() << " "
0225 << std::resetiosflags(std::ios::left) << std::setprecision(3) << (*ranker)(piZero);
0226 out << std::endl;
0227 }
0228 }
0229 }
0230
0231 void RecoTauPiZeroProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0232
0233 edm::ParameterSetDescription desc_ranking;
0234 desc_ranking.add<std::string>("selectionPassFunction", "Func");
0235 desc_ranking.add<double>("selectionFailValue", 1000);
0236 desc_ranking.add<std::string>("selection", "Sel");
0237 desc_ranking.add<std::string>("name", "name");
0238 desc_ranking.add<std::string>("plugin", "plugin");
0239 edm::ParameterSet pset_ranking;
0240 pset_ranking.addParameter<std::string>("selectionPassFunction", "");
0241 pset_ranking.addParameter<double>("selectionFailValue", 1000);
0242 pset_ranking.addParameter<std::string>("selection", "");
0243 pset_ranking.addParameter<std::string>("name", "");
0244 pset_ranking.addParameter<std::string>("plugin", "");
0245 std::vector<edm::ParameterSet> vpsd_ranking;
0246 vpsd_ranking.push_back(pset_ranking);
0247
0248 edm::ParameterSetDescription desc_qualityCuts;
0249 reco::tau::RecoTauQualityCuts::fillDescriptions(desc_qualityCuts);
0250
0251 edm::ParameterSet pset_builders;
0252 pset_builders.addParameter<std::string>("name", "");
0253 pset_builders.addParameter<std::string>("plugin", "");
0254 edm::ParameterSet qualityCuts;
0255 pset_builders.addParameter<edm::ParameterSet>("qualityCuts", qualityCuts);
0256 pset_builders.addParameter<int>("verbosity", 0);
0257
0258 {
0259
0260 edm::ParameterSetDescription desc;
0261 desc.add<double>("massHypothesis", 0.136);
0262 desc.addVPSet("ranking", desc_ranking, vpsd_ranking);
0263 desc.add<int>("verbosity", 0);
0264 desc.add<double>("maxJetAbsEta", 2.5);
0265 desc.add<std::string>("outputSelection", "pt > 0");
0266 desc.add<double>("minJetPt", 14.0);
0267 desc.add<edm::InputTag>("jetSrc", edm::InputTag("ak4PFJets"));
0268
0269 edm::ParameterSetDescription desc_builders;
0270 {
0271 edm::ParameterSetDescription psd0;
0272 psd0.add<std::string>("function", "TMath::Min(0.3, TMath::Max(0.05, [0]*TMath::Power(pT, -[1])))");
0273 psd0.add<double>("par1", 0.707716);
0274 psd0.add<double>("par0", 0.352476);
0275 desc_builders.addOptional<edm::ParameterSetDescription>("stripPhiAssociationDistanceFunc", psd0);
0276 }
0277 {
0278 edm::ParameterSetDescription psd0;
0279 psd0.add<std::string>("function", "TMath::Min(0.15, TMath::Max(0.05, [0]*TMath::Power(pT, -[1])))");
0280 psd0.add<double>("par1", 0.658701);
0281 psd0.add<double>("par0", 0.197077);
0282 desc_builders.addOptional<edm::ParameterSetDescription>("stripEtaAssociationDistanceFunc", psd0);
0283 }
0284 desc_builders.addOptional<double>("stripEtaAssociationDistance", 0.05);
0285 desc_builders.addOptional<double>("stripPhiAssociationDistance", 0.2);
0286
0287 desc_builders.add<edm::ParameterSetDescription>("qualityCuts", desc_qualityCuts);
0288
0289 desc_builders.add<std::string>("name");
0290 desc_builders.add<std::string>("plugin");
0291 desc_builders.add<int>("verbosity", 0);
0292
0293 desc_builders.addOptional<bool>("makeCombinatoricStrips");
0294 desc_builders.addOptional<int>("maxStripBuildIterations");
0295 desc_builders.addOptional<double>("minGammaEtStripAdd");
0296 desc_builders.addOptional<double>("minGammaEtStripSeed");
0297 desc_builders.addOptional<double>("minStripEt");
0298 desc_builders.addOptional<std::vector<int>>("stripCandidatesParticleIds");
0299 desc_builders.addOptional<bool>("updateStripAfterEachDaughter");
0300 desc_builders.addOptional<bool>("applyElecTrackQcuts");
0301
0302 std::vector<edm::ParameterSet> vpsd_builders;
0303 vpsd_builders.push_back(pset_builders);
0304 desc.addVPSet("builders", desc_builders, vpsd_builders);
0305
0306 descriptions.add("recoTauPiZeroProducer", desc);
0307 }
0308 }
0309
0310 #include "FWCore/Framework/interface/MakerMacros.h"
0311 DEFINE_FWK_MODULE(RecoTauPiZeroProducer);