Merge pull request #26616 from intrepid42/PLP_DressedLeptonUpdate

Fix dressed lepton constituents, add tau ghosts
cms-sw · May 4, 2019 · b990340 · b990340
2 parents 45b5000 + 61a9b1c
commit b990340
Show file tree

Hide file tree

Showing 2 changed files with 33 additions and 12 deletions.
diff --git a/GeneratorInterface/RivetInterface/interface/RivetAnalysis.h b/GeneratorInterface/RivetInterface/interface/RivetAnalysis.h
@@ -23,7 +23,7 @@ namespace Rivet {
   class RivetAnalysis : public Analysis {
 
     public:
-      std::vector<DressedLepton> leptons() const {return _leptons;}
+      ParticleVector leptons() const {return _leptons;}
       ParticleVector photons() const {return _photons;}
       ParticleVector neutrinos() const {return _neutrinos;}
       Jets jets() const {return _jets;}
@@ -40,8 +40,7 @@ namespace Rivet {
       double _jetConeSize, _jetMinPt, _jetMaxEta;
       double _fatJetConeSize, _fatJetMinPt, _fatJetMaxEta;
 
-      std::vector<DressedLepton> _leptons;
-      ParticleVector _photons, _neutrinos;
+      ParticleVector _leptons, _photons, _neutrinos;
       Jets _jets, _fatjets;
       Vector3 _met;
 
@@ -85,6 +84,7 @@ namespace Rivet {
         PromptFinalState prompt_leptons(charged_leptons);
         prompt_leptons.acceptMuonDecays(true);
         prompt_leptons.acceptTauDecays(true);
+        addProjection(prompt_leptons, "PromptLeptons");
 
         PromptFinalState prompt_photons(photons);
         prompt_photons.acceptMuonDecays(true);
@@ -151,7 +151,24 @@ namespace Rivet {
         Cut jet_cut    = (Cuts::abseta < _jetMaxEta)    and (Cuts::pT > _jetMinPt*GeV);
         Cut fatjet_cut = (Cuts::abseta < _fatJetMaxEta) and (Cuts::pT > _fatJetMinPt*GeV);
 
-        _leptons   = applyProjection<DressedLeptons>(event, "DressedLeptons").dressedLeptons();
+        _leptons   = applyProjection<DressedLeptons>(event, "DressedLeptons").particlesByPt();
+        // search tau ancestors
+        Particles promptleptons = applyProjection<PromptFinalState>(event, "PromptLeptons").particles();
+        for ( auto & lepton : _leptons ) {
+          const auto & cl = lepton.constituents().front(); // this has no ancestors anymore :(
+          for ( auto const & pl : promptleptons ) {
+            if (cl.momentum() == pl.momentum()) {
+              for ( auto & p : pl.ancestors()) {
+                if (p.abspid() == 15) {
+                  p.setMomentum(p.momentum()*10e-20);
+                  lepton.addConstituent(p, false);
+                }
+              }
+              break;
+            }
+          }
+        }
+
         _jets      = applyProjection<FastJets>(event, "Jets").jetsByPt(jet_cut);
         _fatjets   = applyProjection<FastJets>(event, "FatJets").jetsByPt(fatjet_cut);
         _photons   = applyProjection<FinalState>(event, "Photons").particlesByPt();

diff --git a/GeneratorInterface/RivetInterface/plugins/ParticleLevelProducer.cc b/GeneratorInterface/RivetInterface/plugins/ParticleLevelProducer.cc
@@ -123,28 +123,32 @@ void ParticleLevelProducer::produce(edm::Event& event, const edm::EventSetup& ev
 
   // Prompt leptons
   int iConstituent = -1;
+  int iTag = -1;
   for ( auto const & lepton : rivetAnalysis_->leptons() ) {
     reco::GenJet lepJet;
     lepJet.setP4(p4(lepton));
     lepJet.setVertex(genVertex_);
     lepJet.setPdgId(lepton.pdgId());
     lepJet.setCharge(lepton.charge());
 
-    const auto& cl = lepton.constituentLepton();
-    consts->push_back(reco::GenParticle(cl.charge(), p4(cl), genVertex_, cl.pdgId(), 1, true));
-    lepJet.addDaughter(edm::refToPtr(reco::GenParticleRef(constsRefHandle, ++iConstituent)));
-
-    for ( auto const & p : lepton.constituentPhotons()) {
-      consts->push_back(reco::GenParticle(p.charge(), p4(p), genVertex_, p.pdgId(), 1, true));
-      lepJet.addDaughter(edm::refToPtr(reco::GenParticleRef(constsRefHandle, ++iConstituent)));
+    for ( auto const & p : lepton.constituents()) {
+      // ghost taus (momentum scaled with 10e-20 in RivetAnalysis.h already)
+      if (p.abspid() == 15) {
+        tags->push_back(reco::GenParticle(p.charge(), p4(p), genVertex_, p.pdgId(), 2, true));
+        lepJet.addDaughter(edm::refToPtr(reco::GenParticleRef(tagsRefHandle, ++iTag)));
+      }
+      // electrons, muons, photons
+      else {
+        consts->push_back(reco::GenParticle(p.charge(), p4(p), genVertex_, p.pdgId(), 1, true));
+        lepJet.addDaughter(edm::refToPtr(reco::GenParticleRef(constsRefHandle, ++iConstituent)));
+      }
     }
 
     leptons->push_back(lepJet);
   }
   std::sort(leptons->begin(), leptons->end(), GreaterByPt<reco::GenJet>());
 
   // Jets with constituents and tag particles
-  int iTag = -1;
   for ( auto jet : rivetAnalysis_->jets() ) {
     addGenJet(jet, jets, consts, constsRefHandle, iConstituent, tags, tagsRefHandle, iTag);
   }