Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-12-12 03:12:30

0001 #include "FWCore/Framework/interface/stream/EDProducer.h"
0002 #include "FWCore/Framework/interface/Event.h"
0003 #include "FWCore/Framework/interface/EventSetup.h"
0004 #include "FWCore/Utilities/interface/InputTag.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "DataFormats/PatCandidates/interface/Tau.h"
0007 #include "DataFormats/PatCandidates/interface/Jet.h"
0008 #include "FWCore/Utilities/interface/transform.h"
0009 #include "DataFormats/JetReco/interface/GenJetCollection.h"
0010 #include "RecoTauTag/RecoTau/interface/RecoTauCommonUtilities.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "CommonTools/Utils/interface/PtComparator.h"
0014 
0015 class PATTauHybridProducer : public edm::stream::EDProducer<> {
0016 public:
0017   explicit PATTauHybridProducer(const edm::ParameterSet&);
0018   ~PATTauHybridProducer() override {}
0019 
0020   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0021   void produce(edm::Event&, const edm::EventSetup&) override;
0022 
0023 private:
0024   void fillTauFromJet(reco::PFTau& pfTau, const reco::JetBaseRef& jet);
0025   //--- configuration parameters
0026   edm::EDGetTokenT<pat::TauCollection> tausToken_;
0027   edm::EDGetTokenT<pat::JetCollection> jetsToken_;
0028   bool addGenJetMatch_;
0029   edm::EDGetTokenT<edm::Association<reco::GenJetCollection>> genJetMatchToken_;
0030   const float dR2Max_, jetPtMin_, jetEtaMax_;
0031   const std::string UTagLabel_;
0032   const std::string tagPrefix_;
0033   std::vector<std::string> utagTauScoreNames_;
0034   std::vector<std::string> utagJetScoreNames_;
0035   std::vector<std::string> utagLepScoreNames_;
0036   std::string UtagPtCorrName_;
0037   const float tauScoreMin_, vsJetMin_, chargeAssignmentProbMin_;
0038   const bool checkTauScoreIsBest_;
0039   const bool usePFLeptonsAsChargedHadrons_;
0040 
0041   GreaterByPt<pat::Tau> pTTauComparator_;
0042   const std::map<std::string, int> tagToDM_;
0043   enum class tauId_utag_idx : size_t { dm = 0, vsjet, vse, vsmu, ptcorr, qconf, pdm0, pdm1, pdm2, pdm10, pdm11, last };
0044   enum class tauId_min_idx : size_t { hpsnew = 0, last };
0045 };
0046 PATTauHybridProducer::PATTauHybridProducer(const edm::ParameterSet& cfg)
0047     : tausToken_(consumes<pat::TauCollection>(cfg.getParameter<edm::InputTag>("src"))),
0048       jetsToken_(consumes<pat::JetCollection>(cfg.getParameter<edm::InputTag>("jetSource"))),
0049       addGenJetMatch_(cfg.getParameter<bool>("addGenJetMatch")),
0050       dR2Max_(std::pow(cfg.getParameter<double>("dRMax"), 2)),
0051       jetPtMin_(cfg.getParameter<double>("jetPtMin")),
0052       jetEtaMax_(cfg.getParameter<double>("jetEtaMax")),
0053       UTagLabel_(cfg.getParameter<std::string>("UTagLabel")),
0054       tagPrefix_(cfg.getParameter<std::string>("tagPrefix")),
0055       UtagPtCorrName_(cfg.getParameter<std::string>("UtagPtCorrName")),
0056       tauScoreMin_(cfg.getParameter<double>("tauScoreMin")),
0057       vsJetMin_(cfg.getParameter<double>("vsJetMin")),
0058       chargeAssignmentProbMin_(cfg.getParameter<double>("chargeAssignmentProbMin")),
0059       checkTauScoreIsBest_(cfg.getParameter<bool>("checkTauScoreIsBest")),
0060       usePFLeptonsAsChargedHadrons_(cfg.getParameter<bool>("usePFLeptonsAsChargedHadrons")),
0061       tagToDM_({{"1h0p", 0}, {"1h1or2p", 1}, {"1h1p", 1}, {"1h2p", 2}, {"3h0p", 10}, {"3h1p", 11}}) {
0062   // Read the different Unified Tagger score names
0063   std::vector<std::string> UTagScoreNames = cfg.getParameter<std::vector<std::string>>("UTagScoreNames");
0064   for (const auto& scoreName : UTagScoreNames) {
0065     // Check that discriminator matches tagger specified
0066     if (scoreName.find(UTagLabel_) == std::string::npos)
0067       continue;
0068     size_t labelLength = scoreName.find(':') == std::string::npos ? 0 : scoreName.find(':') + 1;
0069     std::string name = scoreName.substr(labelLength);
0070     if (name.find("prob") == std::string::npos)
0071       continue;
0072     if (name.find("probtau") != std::string::npos)
0073       utagTauScoreNames_.push_back(name);
0074     else if (name == "probele" || name == "probmu")
0075       utagLepScoreNames_.push_back(name);
0076     else if (name.find("data") == std::string::npos && name.find("mc") == std::string::npos)
0077       utagJetScoreNames_.push_back(name);
0078     if (UtagPtCorrName_.find(':') != std::string::npos)
0079       UtagPtCorrName_ = UtagPtCorrName_.substr(UtagPtCorrName_.find(':') + 1);
0080   }
0081   // GenJet matching
0082   if (addGenJetMatch_) {
0083     genJetMatchToken_ =
0084         consumes<edm::Association<reco::GenJetCollection>>(cfg.getParameter<edm::InputTag>("genJetMatch"));
0085   }
0086 
0087   produces<std::vector<pat::Tau>>();
0088   //FIXME: produce a separate collection for PNet-recovered taus?
0089 }
0090 
0091 void PATTauHybridProducer::produce(edm::Event& evt, const edm::EventSetup& es) {
0092   // Get the vector of taus
0093   edm::Handle<pat::TauCollection> inputTaus;
0094   evt.getByToken(tausToken_, inputTaus);
0095 
0096   auto outputTaus = std::make_unique<std::vector<pat::Tau>>();
0097   outputTaus->reserve(inputTaus->size());
0098 
0099   // Get the vector of jets
0100   edm::Handle<pat::JetCollection> jets;
0101   evt.getByToken(jetsToken_, jets);
0102 
0103   // Switch off gen-matching for real data
0104   if (evt.isRealData()) {
0105     addGenJetMatch_ = false;
0106   }
0107   edm::Handle<edm::Association<reco::GenJetCollection>> genJetMatch;
0108   if (addGenJetMatch_)
0109     evt.getByToken(genJetMatchToken_, genJetMatch);
0110 
0111   // Minimal HPS-like tauID list
0112   std::vector<pat::Tau::IdPair> tauIds_minimal((size_t)tauId_min_idx::last);
0113   tauIds_minimal[(size_t)tauId_min_idx::hpsnew] = std::make_pair("decayModeFindingNewDMs", -1);
0114 
0115   // Unified Tagger tauID list
0116   std::vector<pat::Tau::IdPair> tauIds_utag((size_t)tauId_utag_idx::last);
0117   tauIds_utag[(size_t)tauId_utag_idx::dm] = std::make_pair(tagPrefix_ + "DecayMode", reco::PFTau::kNull);
0118   tauIds_utag[(size_t)tauId_utag_idx::vsjet] = std::make_pair(tagPrefix_ + "VSjetraw", -1);
0119   tauIds_utag[(size_t)tauId_utag_idx::vse] = std::make_pair(tagPrefix_ + "VSeraw", -1);
0120   tauIds_utag[(size_t)tauId_utag_idx::vsmu] = std::make_pair(tagPrefix_ + "VSmuraw", -1);
0121   tauIds_utag[(size_t)tauId_utag_idx::ptcorr] = std::make_pair(tagPrefix_ + "PtCorr", 1);
0122   tauIds_utag[(size_t)tauId_utag_idx::qconf] = std::make_pair(tagPrefix_ + "QConf", 0);
0123   tauIds_utag[(size_t)tauId_utag_idx::pdm0] = std::make_pair(tagPrefix_ + "Prob1h0pi0", -1);
0124   tauIds_utag[(size_t)tauId_utag_idx::pdm1] = std::make_pair(tagPrefix_ + "Prob1h1pi0", -1);
0125   tauIds_utag[(size_t)tauId_utag_idx::pdm2] = std::make_pair(tagPrefix_ + "Prob1h2pi0", -1);
0126   tauIds_utag[(size_t)tauId_utag_idx::pdm10] = std::make_pair(tagPrefix_ + "Prob3h0pi0", -1);
0127   tauIds_utag[(size_t)tauId_utag_idx::pdm11] = std::make_pair(tagPrefix_ + "Prob3h1pi0", -1);
0128 
0129   std::set<unsigned int> matched_taus;
0130   size_t jet_idx = 0;
0131   for (const auto& jet : *jets) {
0132     jet_idx++;
0133     if (jet.pt() < jetPtMin_)
0134       continue;
0135     if (std::abs(jet.eta()) > jetEtaMax_)
0136       continue;
0137     size_t tau_idx = 0;
0138     bool matched = false;
0139 
0140     // Analyse Tagger scores
0141     std::pair<std::string, float> bestUtagTauScore("probtauundef", -1);
0142     float sumOfUtagTauScores = 0;
0143     std::vector<float> tauPerDMScores(5);
0144     float plusChargeProb = 0;
0145     for (const auto& scoreName : utagTauScoreNames_) {
0146       float score = jet.bDiscriminator(UTagLabel_ + ":" + scoreName);
0147       sumOfUtagTauScores += score;
0148       if (scoreName.find("taup") != std::string::npos)
0149         plusChargeProb += score;
0150       if (score > bestUtagTauScore.second) {
0151         bestUtagTauScore.second = score;
0152         bestUtagTauScore.first = scoreName;
0153       }
0154       if (scoreName.find("1h0p") != std::string::npos)
0155         tauPerDMScores[0] += score;
0156       else if (scoreName.find("1h1") !=
0157                std::string::
0158                    npos)  //Note: final "p" in "1p" ommited to enble matching also with "1h1or2p" from early trainings
0159         tauPerDMScores[1] += score;
0160       else if (scoreName.find("1h2p") != std::string::npos)
0161         tauPerDMScores[2] += score;
0162       else if (scoreName.find("3h0p") != std::string::npos)
0163         tauPerDMScores[3] += score;
0164       else if (scoreName.find("3h1p") != std::string::npos)
0165         tauPerDMScores[4] += score;
0166     }
0167     if (sumOfUtagTauScores > 0)
0168       plusChargeProb /= sumOfUtagTauScores;
0169 
0170     float sumOfUtagEleScores = 0, sumOfUtagMuScores = 0;
0171     bool isTauScoreBest = (sumOfUtagTauScores > 0);
0172     for (const auto& scoreName : utagLepScoreNames_) {
0173       float score = jet.bDiscriminator(UTagLabel_ + ":" + scoreName);
0174       if (scoreName.find("ele") != std::string::npos)
0175         sumOfUtagEleScores += score;
0176       else if (scoreName.find("mu") != std::string::npos)
0177         sumOfUtagMuScores += score;
0178       if (score > bestUtagTauScore.second)
0179         isTauScoreBest = false;
0180     }
0181     if (checkTauScoreIsBest_ && isTauScoreBest) {  //if needed iterate over jet scores
0182       for (const auto& scoreName : utagJetScoreNames_)
0183         if (jet.bDiscriminator(UTagLabel_ + ":" + scoreName) > bestUtagTauScore.second)
0184           isTauScoreBest = false;
0185     }
0186 
0187     // Unified Tagger discriminants vs jets, electrons and muons
0188     tauIds_utag[(size_t)tauId_utag_idx::vsjet].second =
0189         sumOfUtagTauScores /
0190         (1. - sumOfUtagEleScores -
0191          sumOfUtagMuScores);  //vsJet: tau scores by sum of tau and jet scores or equally by  1 - sum of lepton scores
0192     tauIds_utag[(size_t)tauId_utag_idx::vse].second =
0193         sumOfUtagTauScores /
0194         (sumOfUtagTauScores + sumOfUtagEleScores);  //vsEle: tau scores by sum of tau and ele scores
0195     tauIds_utag[(size_t)tauId_utag_idx::vsmu].second =
0196         sumOfUtagTauScores / (sumOfUtagTauScores + sumOfUtagMuScores);  //vsMu: tau scores by sum of tau and mu scores
0197 
0198     // Decay mode and charge of the highest tau score
0199     int bestCharge = 0;
0200     size_t pos =
0201         bestUtagTauScore.first.find("tau") + 3;  //this is well defined by constraction as name is "probtauXXXX"
0202     const char q = (pos < bestUtagTauScore.first.size()) ? bestUtagTauScore.first[pos] : 'u';
0203     if (q == 'm') {  //minus
0204       pos++;
0205       bestCharge = -1;
0206     } else if (q == 'p') {  //plus
0207       pos++;
0208       bestCharge = 1;
0209     }
0210     auto UtagDM = tagToDM_.find(bestUtagTauScore.first.substr(pos));
0211     if (UtagDM != tagToDM_.end())
0212       tauIds_utag[(size_t)tauId_utag_idx::dm].second = UtagDM->second;
0213 
0214     // Unified tagger Pt correction
0215     float ptcorr = jet.bDiscriminator(UTagLabel_ + ":" + UtagPtCorrName_);
0216     if (ptcorr > -1000.)  // -1000. is default for not found discriminantor
0217       tauIds_utag[(size_t)tauId_utag_idx::ptcorr].second = ptcorr;
0218 
0219     // Unified Tagger charge confidence
0220     tauIds_utag[(size_t)tauId_utag_idx::qconf].second = (plusChargeProb - 0.5);
0221 
0222     // Unified Tagger per decay mode normalised score
0223     tauIds_utag[(size_t)tauId_utag_idx::pdm0].second = tauPerDMScores[0] / sumOfUtagTauScores;
0224     tauIds_utag[(size_t)tauId_utag_idx::pdm1].second = tauPerDMScores[1] / sumOfUtagTauScores;
0225     tauIds_utag[(size_t)tauId_utag_idx::pdm2].second = tauPerDMScores[2] / sumOfUtagTauScores;
0226     tauIds_utag[(size_t)tauId_utag_idx::pdm10].second = tauPerDMScores[3] / sumOfUtagTauScores;
0227     tauIds_utag[(size_t)tauId_utag_idx::pdm11].second = tauPerDMScores[4] / sumOfUtagTauScores;
0228 
0229     // Search for matching tau
0230     for (const auto& inputTau : *inputTaus) {
0231       tau_idx++;
0232       if (matched_taus.count(tau_idx - 1) > 0)
0233         continue;
0234       float dR2 = deltaR2(jet, inputTau);
0235       // select 1st found match rather than best match (both should be equivalent for reasonable dRMax)
0236       if (dR2 < dR2Max_) {
0237         matched_taus.insert(tau_idx - 1);
0238         pat::Tau outputTau(inputTau);
0239         const size_t nTauIds = inputTau.tauIDs().size();
0240         std::vector<pat::Tau::IdPair> tauIds(nTauIds + tauIds_utag.size());
0241         for (size_t i = 0; i < nTauIds; ++i)
0242           tauIds[i] = inputTau.tauIDs()[i];
0243         for (size_t i = 0; i < tauIds_utag.size(); ++i) {
0244           if ((tauIds_utag[i].first.find("PtCorr") != std::string::npos) &&
0245               (inputTau.tauID("decayModeFindingNewDMs") == -1)) {
0246             // if jet is matched to a recovered tau (i.e. non-HPS) then the Pt Correction
0247             // should be adjusted so that it can still be applied as PtCorr * TauPt
0248             // (as the original PtCorr will be w.r.t the jet pT, but recovered tau
0249             // pT is not necessarily set by same jet algorithm if adding both CHS
0250             // and PUPPI based taggers)
0251             tauIds[nTauIds + i].first = tauIds_utag[i].first;
0252             tauIds[nTauIds + i].second = tauIds_utag[i].second * jet.correctedP4("Uncorrected").pt() / inputTau.pt();
0253           } else {
0254             tauIds[nTauIds + i] = tauIds_utag[i];
0255           }
0256         }
0257         outputTau.setTauIDs(tauIds);
0258         matched = true;
0259         outputTaus->push_back(outputTau);
0260 
0261         break;
0262       }
0263     }  // end of tau loop
0264     if (matched)
0265       continue;
0266 
0267     // Accept only jets passing minimal tau-like selection, i.e. with one of the tau score being globally the best and above some threshold, and with good quality of charge assignment
0268     if ((checkTauScoreIsBest_ && !isTauScoreBest) || bestUtagTauScore.second < tauScoreMin_ ||
0269         tauIds_utag[(size_t)tauId_utag_idx::vsjet].second < vsJetMin_ ||
0270         std::abs(0.5 - plusChargeProb) < chargeAssignmentProbMin_)
0271       continue;
0272 
0273     // Build taus from non-matched jets
0274     // "Null" pftau with raw (uncorrected) jet kinematics
0275     reco::PFTau pfTauFromJet(bestCharge, jet.correctedP4("Uncorrected"));
0276     // Set PDGid
0277     pfTauFromJet.setPdgId(bestCharge < 0 ? 15 : -15);
0278     // and decay mode predicted by unified Tagger
0279     pfTauFromJet.setDecayMode(
0280         static_cast<const reco::PFTau::hadronicDecayMode>(int(tauIds_utag[(size_t)tauId_utag_idx::dm].second)));
0281     // Fill tau content using only jet consistunets within cone around leading
0282     // charged particle
0283     // FIXME: more sophisticated finding of tau constituents will be considered later
0284     pfTauFromJet.setSignalConeSize(
0285         std::clamp(3.6 / jet.correctedP4("Uncorrected").pt(), 0.08, 0.12));  // shrinking cone in function of jet-Pt
0286     const edm::Ref<pat::JetCollection> jetRef(jets, jet_idx - 1);
0287     fillTauFromJet(pfTauFromJet, reco::JetBaseRef(jetRef));
0288 
0289     // PATTau
0290     pat::Tau outputTauFromJet(pfTauFromJet);
0291     // Add tauIDs
0292     std::vector<pat::Tau::IdPair> newtauIds(tauIds_minimal.size() + tauIds_utag.size());
0293     for (size_t i = 0; i < tauIds_minimal.size(); ++i)
0294       newtauIds[i] = tauIds_minimal[i];
0295     for (size_t i = 0; i < tauIds_utag.size(); ++i)
0296       newtauIds[tauIds_minimal.size() + i] = tauIds_utag[i];
0297     outputTauFromJet.setTauIDs(newtauIds);
0298     // Add genTauJet match
0299     if (addGenJetMatch_) {
0300       reco::GenJetRef genJetTau = (*genJetMatch)[jetRef];
0301       if (genJetTau.isNonnull() && genJetTau.isAvailable()) {
0302         outputTauFromJet.setGenJet(genJetTau);
0303       }
0304     }
0305     outputTaus->push_back(outputTauFromJet);
0306 
0307   }  // end of jet loop
0308 
0309   // Taus non-matched to jets (usually at pt-threshold or/and eta boundaries)
0310   if (matched_taus.size() < inputTaus->size()) {
0311     for (size_t iTau = 0; iTau < inputTaus->size(); ++iTau) {
0312       if (matched_taus.count(iTau) > 0)
0313         continue;
0314       const pat::Tau& inputTau = inputTaus->at(iTau);
0315       pat::Tau outputTau(inputTau);
0316       const size_t nTauIds = inputTau.tauIDs().size();
0317       std::vector<pat::Tau::IdPair> tauIds(nTauIds + tauIds_utag.size());
0318       for (size_t i = 0; i < nTauIds; ++i)
0319         tauIds[i] = inputTau.tauIDs()[i];
0320       for (size_t i = 0; i < tauIds_utag.size(); ++i) {
0321         tauIds[nTauIds + i] = tauIds_utag[i];
0322         tauIds[nTauIds + i].second =
0323             (i != (size_t)tauId_utag_idx::ptcorr ? (i != (size_t)tauId_utag_idx::qconf ? -1 : 0) : 1);
0324       }
0325       outputTau.setTauIDs(tauIds);
0326       outputTaus->push_back(outputTau);
0327     }
0328   }  //non-matched taus
0329 
0330   // sort taus in pT
0331   std::sort(outputTaus->begin(), outputTaus->end(), pTTauComparator_);
0332 
0333   evt.put(std::move(outputTaus));
0334 }
0335 
0336 void PATTauHybridProducer::fillTauFromJet(reco::PFTau& pfTau, const reco::JetBaseRef& jet) {
0337   // Use tools as in PFTau builders to select tau decay products and isolation candidates
0338   typedef std::vector<reco::CandidatePtr> CandPtrs;
0339 
0340   // Get the charged hadron candidates
0341   CandPtrs pfChs, pfChsSig;
0342   // Check if we want to include electrons and muons in "charged hadron"
0343   // collection. This is the preferred behavior, as the PF lepton selections
0344   // are very loose.
0345   if (!usePFLeptonsAsChargedHadrons_) {
0346     pfChs = reco::tau::pfCandidatesByPdgId(*jet, 211);
0347   } else {
0348     pfChs = reco::tau::pfChargedCands(*jet);
0349   }
0350   // take 1st charged candidate with charge as of tau (collection is pt-sorted)
0351   if (pfTau.charge() == 0 || pfChs.size() == 1) {
0352     pfTau.setleadChargedHadrCand(pfChs[0]);
0353     pfTau.setleadCand(pfChs[0]);
0354     pfChsSig.push_back(pfChs[0]);
0355     pfChs.erase(pfChs.begin());
0356   } else {
0357     for (CandPtrs::iterator it = pfChs.begin(); it != pfChs.end();) {
0358       if ((*it)->charge() == pfTau.charge()) {
0359         pfTau.setleadChargedHadrCand(*it);
0360         pfTau.setleadCand(*it);
0361         pfChsSig.push_back(*it);
0362         it = pfChs.erase(it);
0363         break;
0364       } else {
0365         ++it;
0366       }
0367     }
0368     // In case of lack of candidate with charge same as of tau use leading charged candidate
0369     if (pfTau.leadChargedHadrCand().isNull() && !pfChs.empty()) {
0370       pfTau.setleadChargedHadrCand(pfChs[0]);
0371       pfTau.setleadCand(pfChs[0]);
0372       pfChsSig.push_back(pfChs[0]);
0373       pfChs.erase(pfChs.begin());
0374     }
0375   }
0376   // if more than one charged decay product is expected add all inside signal
0377   // cone around the leading track
0378   const double dR2Max = std::pow(pfTau.signalConeSize(), 2);
0379   if (pfTau.decayMode() >= reco::PFTau::kThreeProng0PiZero && pfTau.leadChargedHadrCand().isNonnull()) {
0380     for (CandPtrs::iterator it = pfChs.begin(); it != pfChs.end();) {
0381       if (deltaR2((*it)->p4(), pfTau.leadChargedHadrCand()->p4()) < dR2Max) {
0382         pfChsSig.push_back(*it);
0383         it = pfChs.erase(it);
0384       } else {
0385         ++it;
0386       }
0387     }
0388   }
0389   // Clean isolation candidates from low-pt and leptonic ones
0390   pfChs.erase(std::remove_if(pfChs.begin(),
0391                              pfChs.end(),
0392                              [](auto const& cand) { return cand->pt() < 0.5 || std::abs(cand->pdgId()) != 211; }),
0393               pfChs.end());
0394   // Set charged candidates
0395   pfTau.setsignalChargedHadrCands(pfChsSig);
0396   pfTau.setisolationChargedHadrCands(pfChs);
0397 
0398   // Get the gamma candidates (pi0 decay products)
0399   CandPtrs pfGammas, pfGammasSig;
0400   pfGammas = reco::tau::pfCandidatesByPdgId(*jet, 22);
0401   // In case of lack of leading charged candidate substiute it with leading gamma candidate
0402   if (pfTau.leadChargedHadrCand().isNull() && !pfGammas.empty()) {
0403     pfTau.setleadChargedHadrCand(pfGammas[0]);
0404     pfTau.setleadCand(pfGammas[0]);
0405     pfGammasSig.push_back(pfGammas[0]);
0406     pfGammas.erase(pfGammas.begin());
0407   }
0408   // Clean gamma candidates from low-pt ones
0409   pfGammas.erase(std::remove_if(pfGammas.begin(), pfGammas.end(), [](auto const& cand) { return cand->pt() < 0.5; }),
0410                  pfGammas.end());
0411   // if decay mode with pi0s is expected look for signal gamma candidates
0412   // within eta-phi strips around leading track
0413   if (pfTau.decayMode() % 5 != 0 && pfTau.leadChargedHadrCand().isNonnull()) {
0414     for (CandPtrs::iterator it = pfGammas.begin(); it != pfGammas.end();) {
0415       if (std::abs((*it)->eta() - pfTau.leadChargedHadrCand()->eta()) <
0416               std::clamp(0.2 * std::pow((*it)->pt(), -0.66), 0.05, 0.15) &&
0417           deltaPhi((*it)->phi(), pfTau.leadChargedHadrCand()->phi()) <
0418               std::clamp(0.35 * std::pow((*it)->pt(), -0.71), 0.05, 0.3)) {
0419         pfGammasSig.push_back(*it);
0420         it = pfGammas.erase(it);
0421       } else {
0422         ++it;
0423       }
0424     }
0425   }
0426   // Set gamma candidates
0427   pfTau.setsignalGammaCands(pfGammasSig);
0428   pfTau.setisolationGammaCands(pfGammas);
0429 }
0430 
0431 void PATTauHybridProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0432   // patTauHybridProducer
0433   edm::ParameterSetDescription desc;
0434 
0435   desc.add<edm::InputTag>("src", edm::InputTag("slimmedTaus"));
0436   desc.add<edm::InputTag>("jetSource", edm::InputTag("slimmedJetsUpdated"));
0437   desc.add<double>("dRMax", 0.4);
0438   desc.add<double>("jetPtMin", 20.0);
0439   desc.add<double>("jetEtaMax", 2.5);
0440   desc.add<std::string>("UTagLabel", "pfParticleNetAK4JetTags");
0441   desc.add<std::string>("tagPrefix", "byUTag")->setComment("Prefix to be set for PUPPI or CHS tagger");
0442   desc.add<std::vector<std::string>>("UTagScoreNames", {})
0443       ->setComment("Output nodes for Unified Tagger (different for PNet vs UParT)");
0444   desc.add<std::string>("UtagPtCorrName", "ptcorr");
0445   desc.add<double>("tauScoreMin", -1)->setComment("Minimal value of the best tau score to built recovery tau");
0446   desc.add<double>("vsJetMin", -1)->setComment("Minimal value of UTag tau-vs-jet discriminant to built recovery tau");
0447   desc.add<bool>("checkTauScoreIsBest", false)
0448       ->setComment("If true, recovery tau is built only if one of tau scores is the highest");
0449   desc.add<double>("chargeAssignmentProbMin", 0.2)
0450       ->setComment("Minimal value of charge assignment probability to built recovery tau (0,0.5)");
0451   desc.add<bool>("addGenJetMatch", true)->setComment("add MC genTauJet matching");
0452   desc.add<edm::InputTag>("genJetMatch", edm::InputTag("tauGenJetMatch"));
0453   desc.add<bool>("usePFLeptonsAsChargedHadrons", true)
0454       ->setComment("If true, all charged particles are used as charged hadron candidates");
0455 
0456   descriptions.addWithDefaultLabel(desc);
0457 }
0458 
0459 #include "FWCore/Framework/interface/MakerMacros.h"
0460 
0461 DEFINE_FWK_MODULE(PATTauHybridProducer);