Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:01:11

0001 /** \class LeptonFSRProducer
0002  * Search for FSR photons for muons and electrons.
0003  *
0004  * Photon candidates are searched among the "packedPFCandidates" collection with the specified cuts, and are required to be isolated 
0005  * (relIso, with a cone of 0.3) and not to be in the footprint of all electrons in the "electrons" collection.
0006  * Each photon is matched by DeltaR to the closest among all muons and electrons and stored if passing dR/ET^2<deltaROverEt2Max.
0007  * In addition ValueMaps are stored, with links to one photon per muon/electron. For this purpose, if more than a photon
0008  * is matched to a lepton, the lowest-DR/ET^2 is chosen.
0009  *
0010  */
0011 
0012 #include <memory>
0013 #include "FWCore/Framework/interface/Frameworkfwd.h"
0014 #include "FWCore/Framework/interface/global/EDProducer.h"
0015 
0016 #include "FWCore/Framework/interface/Event.h"
0017 #include "FWCore/Framework/interface/MakerMacros.h"
0018 
0019 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0020 #include "FWCore/Utilities/interface/StreamID.h"
0021 
0022 #include "DataFormats/Candidate/interface/Candidate.h"
0023 #include "DataFormats/PatCandidates/interface/PackedCandidate.h"
0024 #include "DataFormats/PatCandidates/interface/GenericParticle.h"
0025 #include "DataFormats/Math/interface/LorentzVector.h"
0026 
0027 #include "DataFormats/PatCandidates/interface/Muon.h"
0028 #include "DataFormats/PatCandidates/interface/Electron.h"
0029 #include "DataFormats/Common/interface/ValueMap.h"
0030 
0031 class LeptonFSRProducer : public edm::global::EDProducer<> {
0032 public:
0033   explicit LeptonFSRProducer(const edm::ParameterSet& iConfig)
0034       : pfcands_{consumes<pat::PackedCandidateCollection>(iConfig.getParameter<edm::InputTag>("packedPFCandidates"))},
0035         electronsForVeto_{consumes<pat::ElectronCollection>(iConfig.getParameter<edm::InputTag>("slimmedElectrons"))},
0036         muons_{consumes<edm::View<reco::Muon>>(iConfig.getParameter<edm::InputTag>("muons"))},
0037         electrons_{consumes<edm::View<reco::GsfElectron>>(iConfig.getParameter<edm::InputTag>("electrons"))},
0038         ptCutMu(iConfig.getParameter<double>("muonPtMin")),
0039         etaCutMu(iConfig.getParameter<double>("muonEtaMax")),
0040         ptCutE(iConfig.getParameter<double>("elePtMin")),
0041         etaCutE(iConfig.getParameter<double>("eleEtaMax")),
0042         photonPtCut(iConfig.getParameter<double>("photonPtMin")),
0043         drEtCut(iConfig.getParameter<double>("deltaROverEt2Max")),
0044         isoCut(iConfig.getParameter<double>("isolation")),
0045         drSafe(0.0001) {
0046     produces<std::vector<pat::GenericParticle>>();
0047     produces<edm::ValueMap<int>>("muFsrIndex");
0048     produces<edm::ValueMap<int>>("eleFsrIndex");
0049   }
0050   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0051     edm::ParameterSetDescription desc;
0052     desc.add<edm::InputTag>("packedPFCandidates", edm::InputTag("packedPFCandidates"))
0053         ->setComment("packed pf candidates where to look for photons");
0054     desc.add<edm::InputTag>("slimmedElectrons", edm::InputTag("slimmedElectrons"))
0055         ->setComment(
0056             "electrons to check for footprint, the electron collection must have proper linking with the "
0057             "packedCandidate collection");
0058     desc.add<edm::InputTag>("muons", edm::InputTag("slimmedMuons"))
0059         ->setComment("collection of muons to match with FSR ");
0060     desc.add<edm::InputTag>("electrons", edm::InputTag("slimmedElectrons"))
0061         ->setComment("collection of electrons to match with FSR ");
0062     desc.add<double>("muonPtMin", 3.)->setComment("minimum pt of the muon to look for a near photon");
0063     desc.add<double>("muonEtaMax", 2.4)->setComment("max eta of the muon to look for a near photon");
0064     desc.add<double>("elePtMin", 5.)->setComment("minimum pt of the electron to look for a near photon");
0065     desc.add<double>("eleEtaMax", 2.5)->setComment("max eta of the electron to look for a near photon");
0066     desc.add<double>("photonPtMin", 2.0)->setComment("minimum photon Pt");
0067     desc.add<double>("deltaROverEt2Max", 0.05)->setComment("max ratio of deltaR(lep,photon) over et2 of the photon");
0068     desc.add<double>("isolation", 2.0)->setComment("photon relative isolation cut");
0069 
0070     descriptions.addWithDefaultLabel(desc);
0071   }
0072   ~LeptonFSRProducer() override = default;
0073 
0074 private:
0075   void produce(edm::StreamID, edm::Event&, const edm::EventSetup&) const override;
0076 
0077   double computeRelativeIsolation(const pat::PackedCandidate& photon,
0078                                   const pat::PackedCandidateCollection& pfcands,
0079                                   const double& isoConeMax2,
0080                                   const double& isoConeMin2) const;
0081 
0082   bool electronFootprintVeto(pat::PackedCandidateRef& pfcandRef,
0083                              edm::Handle<pat::ElectronCollection> electronsForVeto) const;
0084 
0085   // ----------member data ---------------------------
0086   const edm::EDGetTokenT<pat::PackedCandidateCollection> pfcands_;
0087   const edm::EDGetTokenT<pat::ElectronCollection> electronsForVeto_;
0088   const edm::EDGetTokenT<edm::View<reco::Muon>> muons_;
0089   const edm::EDGetTokenT<edm::View<reco::GsfElectron>> electrons_;
0090   const double ptCutMu;
0091   const double etaCutMu;
0092   const double ptCutE;
0093   const double etaCutE;
0094   const double photonPtCut;
0095   const double drEtCut;
0096   const double isoCut;
0097   const double drSafe;
0098 };
0099 
0100 void LeptonFSRProducer::produce(edm::StreamID streamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0101   using namespace std;
0102 
0103   edm::Handle<pat::PackedCandidateCollection> pfcands;
0104   iEvent.getByToken(pfcands_, pfcands);
0105   edm::Handle<edm::View<reco::Muon>> muons;
0106   iEvent.getByToken(muons_, muons);
0107   edm::Handle<edm::View<reco::GsfElectron>> electrons;
0108   iEvent.getByToken(electrons_, electrons);
0109   edm::Handle<pat::ElectronCollection> electronsForVeto;
0110   iEvent.getByToken(electronsForVeto_, electronsForVeto);
0111 
0112   // The output collection of FSR photons
0113   auto fsrPhotons = std::make_unique<std::vector<pat::GenericParticle>>();
0114 
0115   std::vector<int> muPhotonIdxs(muons->size(), -1);
0116   std::vector<double> muPhotonDRET2(muons->size(), 1e9);
0117 
0118   std::vector<int> elePhotonIdxs(electrons->size(), -1);
0119   std::vector<double> elePhotonDRET2(electrons->size(), 1e9);
0120 
0121   //----------------------
0122   // Loop on photon candidates
0123   //----------------------
0124 
0125   for (auto pc = pfcands->begin(); pc != pfcands->end(); pc++) {
0126     // consider only photons, with pT and eta cuts
0127     if (abs(pc->pdgId()) != 22 || pc->pt() < photonPtCut || abs(pc->eta()) > 2.5)
0128       continue;
0129 
0130     //------------------------------------------------------
0131     // Get the closest lepton
0132     //------------------------------------------------------
0133     double dRMin(0.5);
0134     int closestMu = -1;
0135     int closestEle = -1;
0136     double photon_relIso03 = 1e9;  // computed only if necessary
0137     bool skipPhoton = false;
0138 
0139     for (auto muon = muons->begin(); muon != muons->end(); ++muon) {
0140       if (muon->pt() < ptCutMu || std::abs(muon->eta()) > etaCutMu)
0141         continue;
0142 
0143       int muonIdx = muon - muons->begin();
0144       double dR = deltaR(muon->eta(), muon->phi(), pc->eta(), pc->phi());
0145       if (dR < dRMin && dR > drSafe && dR < drEtCut * pc->pt() * pc->pt()) {
0146         // Check if photon is isolated
0147         photon_relIso03 = computeRelativeIsolation(*pc, *pfcands, 0.3 * 0.3, drSafe * drSafe);
0148         if (photon_relIso03 > isoCut) {
0149           skipPhoton = true;
0150           break;  // break loop on muons -> photon will be skipped
0151         }
0152         // Check that photon is not in footprint of an electron
0153         pat::PackedCandidateRef pfcandRef = pat::PackedCandidateRef(pfcands, pc - pfcands->begin());
0154         skipPhoton = electronFootprintVeto(pfcandRef, electronsForVeto);
0155         if (skipPhoton)
0156           break;  // break loop on muons -> photon will be skipped
0157 
0158         // Candidate matching
0159         dRMin = dR;
0160         closestMu = muonIdx;
0161       }
0162     }  // end of loop on muons
0163 
0164     if (skipPhoton)
0165       continue;  // photon does not pass iso or ele footprint veto; do not look for electrons
0166 
0167     for (auto ele = electrons->begin(); ele != electrons->end(); ++ele) {
0168       if (ele->pt() < ptCutE || std::abs(ele->eta()) > etaCutE)
0169         continue;
0170 
0171       int eleIdx = ele - electrons->begin();
0172       double dR = deltaR(ele->eta(), ele->phi(), pc->eta(), pc->phi());
0173       if (dR < dRMin && dR > drSafe && dR < drEtCut * pc->pt() * pc->pt()) {
0174         // Check if photon is isolated (no need to recompute iso if already done for muons above)
0175         if (photon_relIso03 > 1e8) {
0176           photon_relIso03 = computeRelativeIsolation(*pc, *pfcands, 0.3 * 0.3, drSafe * drSafe);
0177         }
0178         if (photon_relIso03 > isoCut) {
0179           break;  // break loop on electrons -> photon will be skipped
0180         }
0181         // Check that photon is not in footprint of an electron
0182         pat::PackedCandidateRef pfcandRef = pat::PackedCandidateRef(pfcands, pc - pfcands->begin());
0183         if (electronFootprintVeto(pfcandRef, electronsForVeto)) {
0184           break;  // break loop on electrons -> photon will be skipped
0185         }
0186 
0187         // Candidate matching
0188         dRMin = dR;
0189         closestEle = eleIdx;
0190         closestMu = -1;  // reset match to muons
0191       }
0192     }  // end loop on electrons
0193 
0194     if (closestMu >= 0 || closestEle >= 0) {
0195       // Add FSR photon to the output collection
0196       double dRET2 = dRMin / pc->pt() / pc->pt();
0197       int iPhoton = fsrPhotons->size();
0198       fsrPhotons->push_back(pat::GenericParticle(*pc));
0199       fsrPhotons->back().addUserFloat("relIso03", photon_relIso03);
0200       fsrPhotons->back().addUserFloat("dROverEt2", dRET2);
0201 
0202       if (closestMu >= 0) {
0203         fsrPhotons->back().addUserCand("associatedMuon", reco::CandidatePtr(muons, closestMu));
0204         // Store the backlink to the photon: choose the lowest-dRET2 photon for each mu...
0205         if (dRET2 < muPhotonDRET2[closestMu]) {
0206           muPhotonDRET2[closestMu] = dRET2;
0207           muPhotonIdxs[closestMu] = iPhoton;
0208         }
0209       } else if (closestEle >= 0) {
0210         // ...and same for eles
0211         fsrPhotons->back().addUserCand("associatedElectron", reco::CandidatePtr(electrons, closestEle));
0212         if (dRET2 < elePhotonDRET2[closestEle]) {
0213           elePhotonDRET2[closestEle] = dRET2;
0214           elePhotonIdxs[closestEle] = iPhoton;
0215         }
0216       }
0217     }
0218   }  // end of loop over pfCands
0219 
0220   edm::OrphanHandle<std::vector<pat::GenericParticle>> oh = iEvent.put(std::move(fsrPhotons));
0221 
0222   {
0223     std::unique_ptr<edm::ValueMap<int>> bareIdx(new edm::ValueMap<int>());
0224     edm::ValueMap<int>::Filler fillerBareIdx(*bareIdx);
0225     fillerBareIdx.insert(muons, muPhotonIdxs.begin(), muPhotonIdxs.end());
0226     fillerBareIdx.fill();
0227     iEvent.put(std::move(bareIdx), "muFsrIndex");
0228   }
0229 
0230   {
0231     std::unique_ptr<edm::ValueMap<int>> bareIdx(new edm::ValueMap<int>());
0232     edm::ValueMap<int>::Filler fillerBareIdx(*bareIdx);
0233     fillerBareIdx.insert(electrons, elePhotonIdxs.begin(), elePhotonIdxs.end());
0234     fillerBareIdx.fill();
0235     iEvent.put(std::move(bareIdx), "eleFsrIndex");
0236   }
0237 }
0238 
0239 double LeptonFSRProducer::computeRelativeIsolation(const pat::PackedCandidate& photon,
0240                                                    const pat::PackedCandidateCollection& pfcands,
0241                                                    const double& isoConeMax2,
0242                                                    const double& isoConeMin2) const {
0243   double ptsum = 0;
0244 
0245   for (const auto& pfcand : pfcands) {
0246     // Isolation cone
0247     double dRIsoCone2 = deltaR2(photon.eta(), photon.phi(), pfcand.eta(), pfcand.phi());
0248     if (dRIsoCone2 > isoConeMax2 || dRIsoCone2 < isoConeMin2)
0249       continue;
0250 
0251     // Charged hadrons
0252     if (pfcand.charge() != 0 && abs(pfcand.pdgId()) == 211 && pfcand.pt() > 0.2 && dRIsoCone2 > drSafe * drSafe) {
0253       ptsum += pfcand.pt();
0254       // Neutral hadrons + photons
0255     } else if (pfcand.charge() == 0 && (abs(pfcand.pdgId()) == 22 || abs(pfcand.pdgId()) == 130) && pfcand.pt() > 0.5 &&
0256                dRIsoCone2 > 0.01 * 0.01) {
0257       ptsum += pfcand.pt();
0258     }
0259   }
0260 
0261   return ptsum / photon.pt();
0262 }
0263 
0264 bool LeptonFSRProducer::electronFootprintVeto(pat::PackedCandidateRef& pfcandRef,
0265                                               edm::Handle<pat::ElectronCollection> electronsForVeto) const {
0266   bool skipPhoton = false;
0267   for (auto electrons_iter = electronsForVeto->begin(); electrons_iter != electronsForVeto->end(); ++electrons_iter) {
0268     for (auto const& cand : electrons_iter->associatedPackedPFCandidates()) {
0269       if (!cand.isAvailable())
0270         continue;
0271       if (cand.id() != pfcandRef.id())
0272         throw cms::Exception("Configuration")
0273             << "The electron associatedPackedPFCandidates item does not have "
0274             << "the same ID of packed candidate collection used for cleaning the electron footprint: " << cand.id()
0275             << " (" << pfcandRef.id() << ")\n";
0276       if (cand.key() == pfcandRef.key()) {
0277         skipPhoton = true;
0278         break;
0279       }
0280     }
0281   }
0282   return skipPhoton;
0283 }
0284 
0285 //define this as a plug-in
0286 DEFINE_FWK_MODULE(LeptonFSRProducer);