CHS and PUPPI track-vertex-association for Run3

CommonTools/ParticleFlow/plugins/PFPileUp.cc

@@ -76,6 +76,8 @@ class PFPileUp : public edm::stream::EDProducer<> {
   edm::EDGetTokenT<edm::ValueMap<int>> tokenVertexAssociationQuality_;
   bool fUseVertexAssociation;
   int vertexAssociationQuality_;
+  unsigned int fNumOfPUVtxsForCharged_;
+  double fDzCutForChargedFromPUVtxs_;
 };
 
 PFPileUp::PFPileUp(const edm::ParameterSet& iConfig) {
@@ -91,6 +93,8 @@ PFPileUp::PFPileUp(const edm::ParameterSet& iConfig) {
     tokenVertexAssociationQuality_ =
         consumes<edm::ValueMap<int>>(iConfig.getParameter<edm::InputTag>("vertexAssociation"));
   }
+  fNumOfPUVtxsForCharged_ = iConfig.getParameter<unsigned int>("NumOfPUVtxsForCharged");
+  fDzCutForChargedFromPUVtxs_ = iConfig.getParameter<double>("DzCutForChargedFromPUVtxs");
 
   enable_ = iConfig.getParameter<bool>("Enable");
 
@@ -172,7 +176,9 @@ void PFPileUp::produce(Event& iEvent, const EventSetup& iSetup) {
       for (auto& p : (*pfCandidatesRef)) {
         const reco::VertexRef& PVOrig = associatedPV[p];
         int quality = associationQuality[p];
-        if (PVOrig.isNonnull() && (PVOrig.key() > 0) && (quality >= vertexAssociationQuality_))
+        if ((PVOrig.isNonnull() && PVOrig.key() > 0 && PVOrig->ndof() > 4.0 && quality >= vertexAssociationQuality_) &&
+            (!(fNumOfPUVtxsForCharged_ > 0 && !vertices->empty() && PVOrig.key() <= fNumOfPUVtxsForCharged_ &&
+               std::abs(p->vertex().z() - vertices->at(0).z()) < fDzCutForChargedFromPUVtxs_)))
           pfCandidatesFromPU.push_back(p);
       }
       pOutput->insert(pOutput->end(), pfCandidatesFromPU.begin(), pfCandidatesFromPU.end());

CommonTools/ParticleFlow/python/pfCHS_cff.py

@@ -1,15 +1,10 @@
 import FWCore.ParameterSet.Config as cms
-from CommonTools.ParticleFlow.pfNoPileUpJME_cff import adapt, pfPileUpJME
-from CommonTools.RecoAlgos.sortedPackedPrimaryVertices_cfi import sortedPackedPrimaryVertices
+from CommonTools.ParticleFlow.pfNoPileUpJME_cff import primaryVertexAssociationJME
 
-packedPrimaryVertexAssociationJME = sortedPackedPrimaryVertices.clone(
-  produceSortedVertices = False,
-  producePileUpCollection  = False,
-  produceNoPileUpCollection = False
-)
-adapt(packedPrimaryVertexAssociationJME)
-
-from CommonTools.ParticleFlow.pfNoPileUpPacked_cfi import pfNoPileUpPacked as _pfNoPileUpPacked
-pfCHS = _pfNoPileUpPacked.clone(
-    vertexAssociationQuality=pfPileUpJME.vertexAssociationQuality
+pfCHS = cms.EDFilter("CandPtrSelector",
+              src = cms.InputTag("packedPFCandidates"),
+              cut = cms.string("fromPV(0)>0"+\
+	       " || (vertexRef().key<="+str(primaryVertexAssociationJME.assignment.NumOfPUVtxsForCharged.value())+" && "+\
+		    "abs(dz(0))<"+str(primaryVertexAssociationJME.assignment.DzCutForChargedFromPUVtxs.value())+")"
+                    )
 )

CommonTools/ParticleFlow/python/pfNoPileUpJME_cff.py

@@ -5,35 +5,36 @@
 from CommonTools.ParticleFlow.goodOfflinePrimaryVertices_cfi import goodOfflinePrimaryVertices
 from CommonTools.RecoAlgos.primaryVertexAssociation_cfi import primaryVertexAssociation
 
-def adapt(primaryVertexAssociationJME):
-  # options for quality PrimaryDz = 6 (used in PUPPI)
-  primaryVertexAssociationJME.assignment.maxDzSigForPrimaryAssignment = 1e10
-  primaryVertexAssociationJME.assignment.maxDzForPrimaryAssignment = 0.3
-  primaryVertexAssociationJME.assignment.maxDzErrorForPrimaryAssignment = 1e10
-  primaryVertexAssociationJME.assignment.NumOfPUVtxsForCharged = 2
-  primaryVertexAssociationJME.assignment.PtMaxCharged = 20.
-  primaryVertexAssociationJME.assignment.EtaMinUseDz = 2.4
-  primaryVertexAssociationJME.assignment.OnlyUseFirstDz = True
-  from Configuration.Eras.Modifier_phase2_common_cff import phase2_common
-  phase2_common.toModify(
-    primaryVertexAssociationJME.assignment,
-    maxDzForPrimaryAssignment=0.1,
-    EtaMinUseDz = 4.0
-    )
+# The following module implements vertex association for JME.
+# It is not run by default to save computing time (but can be run to investigate alternative vertex identification criteria).
+# Instead its parameters are used as input to faster implementations in pfPileUpJME and pfCHS and puppi.
 primaryVertexAssociationJME = primaryVertexAssociation.clone(vertices = 'goodOfflinePrimaryVertices')
-adapt(primaryVertexAssociationJME)
+primaryVertexAssociationJME.assignment.maxDzSigForPrimaryAssignment = 1e10
+primaryVertexAssociationJME.assignment.maxDzForPrimaryAssignment = 0.3
+primaryVertexAssociationJME.assignment.maxDzErrorForPrimaryAssignment = 1e10
+primaryVertexAssociationJME.assignment.NumOfPUVtxsForCharged = 2
+primaryVertexAssociationJME.assignment.DzCutForChargedFromPUVtxs = 0.2
+primaryVertexAssociationJME.assignment.PtMaxCharged = 20.
+primaryVertexAssociationJME.assignment.EtaMinUseDz = 2.4
+primaryVertexAssociationJME.assignment.OnlyUseFirstDz = True
+from Configuration.Eras.Modifier_phase2_common_cff import phase2_common
+phase2_common.toModify(
+  primaryVertexAssociationJME.assignment,
+  maxDzForPrimaryAssignment=0.1,
+  EtaMinUseDz = 4.0
+  )
 
 pfPileUpJME = _pfPileUp.clone(PFCandidates='particleFlowPtrs',
                               useVertexAssociation = True,
                               vertexAssociationQuality = 7,
-                              vertexAssociation = ('primaryVertexAssociationJME','original'),
+                              vertexAssociation = ('primaryVertexAssociation','original'),
+                              NumOfPUVtxsForCharged = primaryVertexAssociationJME.assignment.NumOfPUVtxsForCharged,
+                              DzCutForChargedFromPUVtxs = primaryVertexAssociationJME.assignment.DzCutForChargedFromPUVtxs,
                               )
 pfNoPileUpJME = _pfNoPileUp.clone(topCollection = 'pfPileUpJME',
                                   bottomCollection = 'particleFlowPtrs' )
 
 pfNoPileUpJMETask = cms.Task(
-    goodOfflinePrimaryVertices,
-    primaryVertexAssociationJME,
     pfPileUpJME,
     pfNoPileUpJME
     )

CommonTools/ParticleFlow/python/pfPileUp_cfi.py

@@ -12,4 +12,6 @@
     useVertexAssociation = cms.bool(False),
     vertexAssociationQuality = cms.int32(0),
     vertexAssociation = cms.InputTag(''),
+    NumOfPUVtxsForCharged = cms.uint32(0),
+    DzCutForChargedFromPUVtxs = cms.double(0.2),
     )

CommonTools/PileupAlgos/plugins/PuppiProducer.cc

@@ -72,6 +72,7 @@ class PuppiProducer : public edm::stream::EDProducer<> {
   bool fPuppiNoLep;
   bool fUseFromPVLooseTight;
   bool fUseDZ;
+  bool fUseDZforPileup;
   double fDZCut;
   double fEtaMinUseDZ;
   double fPtMaxCharged;
@@ -95,6 +96,7 @@ PuppiProducer::PuppiProducer(const edm::ParameterSet& iConfig) {
   fPuppiNoLep = iConfig.getParameter<bool>("puppiNoLep");
   fUseFromPVLooseTight = iConfig.getParameter<bool>("UseFromPVLooseTight");
   fUseDZ = iConfig.getParameter<bool>("UseDeltaZCut");
+  fUseDZforPileup = iConfig.getParameter<bool>("UseDeltaZCutForPileup");
   fDZCut = iConfig.getParameter<double>("DeltaZCut");
   fEtaMinUseDZ = iConfig.getParameter<double>("EtaMinUseDeltaZ");
   fPtMaxCharged = iConfig.getParameter<double>("PtMaxCharged");
@@ -281,8 +283,10 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
               pReco.id = 1;
             else if (std::abs(pReco.eta) > fEtaMaxCharged)
               pReco.id = 1;
-            else if ((fUseDZ) && (std::abs(pReco.eta) >= fEtaMinUseDZ))
-              pReco.id = (std::abs(pDZ) < fDZCut) ? 1 : 2;
+            else if ((fUseDZ) && (std::abs(pReco.eta) >= fEtaMinUseDZ) && (std::abs(pDZ) < fDZCut))
+              pReco.id = 1;
+            else if ((fUseDZforPileup) && (std::abs(pReco.eta) >= fEtaMinUseDZ) && (std::abs(pDZ) >= fDZCut))
+              pReco.id = 2;
             else if (fUseFromPVLooseTight && tmpFromPV == 1)
               pReco.id = 2;
             else if (fUseFromPVLooseTight && tmpFromPV == 2)
@@ -322,8 +326,10 @@ void PuppiProducer::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
               pReco.id = 1;
             else if (std::abs(pReco.eta) > fEtaMaxCharged)
               pReco.id = 1;
-            else if ((fUseDZ) && (std::abs(pReco.eta) >= fEtaMinUseDZ))
-              pReco.id = (std::abs(pDZ) < fDZCut) ? 1 : 2;
+            else if ((fUseDZ) && (std::abs(pReco.eta) >= fEtaMinUseDZ) && (std::abs(pDZ) < fDZCut))
+              pReco.id = 1;
+            else if ((fUseDZforPileup) && (std::abs(pReco.eta) >= fEtaMinUseDZ) && (std::abs(pDZ) >= fDZCut))
+              pReco.id = 2;
             else if (fUseFromPVLooseTight && lPack->fromPV() == (pat::PackedCandidate::PVLoose))
               pReco.id = 2;
             else if (fUseFromPVLooseTight && lPack->fromPV() == (pat::PackedCandidate::PVTight))
@@ -489,6 +495,7 @@ void PuppiProducer::fillDescriptions(edm::ConfigurationDescriptions& description
   desc.add<bool>("puppiNoLep", false);
   desc.add<bool>("UseFromPVLooseTight", false);
   desc.add<bool>("UseDeltaZCut", true);
+  desc.add<bool>("UseDeltaZCutForPileup", true);
   desc.add<double>("DeltaZCut", 0.3);
   desc.add<double>("EtaMinUseDeltaZ", 0.);
   desc.add<double>("PtMaxCharged", -1.);

CommonTools/PileupAlgos/python/Puppi_cff.py

@@ -25,15 +25,16 @@
                 )
 
 # from cfipython/, generated by PuppiProducer::fillDescriptions
+from CommonTools.ParticleFlow.pfNoPileUpJME_cff import primaryVertexAssociationJME
 import CommonTools.PileupAlgos.PuppiProducer_cfi as _mod
 puppi = _mod.PuppiProducer.clone(
-                       EtaMinUseDeltaZ = 2.4,
-                       NumOfPUVtxsForCharged = 2,
-                       PtMaxCharged   = 20.,
+                       UseDeltaZCutForPileup = False,
+                       DeltaZCut = primaryVertexAssociationJME.assignment.maxDzForPrimaryAssignment,
+                       EtaMinUseDeltaZ = primaryVertexAssociationJME.assignment.EtaMinUseDz,
+                       PtMaxCharged = primaryVertexAssociationJME.assignment.PtMaxCharged,
+                       NumOfPUVtxsForCharged = primaryVertexAssociationJME.assignment.NumOfPUVtxsForCharged,
+                       DeltaZCutForChargedFromPUVtxs = primaryVertexAssociationJME.assignment.DzCutForChargedFromPUVtxs,
                        PtMaxNeutralsStartSlope = 20.,
-                       useVertexAssociation = True,
-                       vertexAssociationQuality = 6,
-                       vertexAssociation = ('primaryVertexAssociationJME','original'),
                        clonePackedCands   = False, # should only be set to True for MiniAOD
                        algos          = { 
                         0: dict( 
@@ -99,8 +100,6 @@
 
 from Configuration.ProcessModifiers.pp_on_AA_cff import pp_on_AA
 pp_on_AA.toModify(puppi, algos = [])
-from Configuration.ProcessModifiers.run2_miniAOD_pp_on_AA_103X_cff import run2_miniAOD_pp_on_AA_103X
-run2_miniAOD_pp_on_AA_103X.toModify(puppi,useVertexAssociation = False) # because the association is only run on cleanedParticleFlow
 
 puppiNoLep = puppi.clone(
     puppiNoLep = True,

CommonTools/RecoAlgos/interface/PrimaryVertexAssignment.h

@@ -44,7 +44,7 @@ class PrimaryVertexAssignment {
         useTiming_(iConfig.getParameter<bool>("useTiming")),
         useVertexFit_(iConfig.getParameter<bool>("useVertexFit")),
         preferHighRanked_(iConfig.getParameter<bool>("preferHighRanked")),
-        fNumOfPUVtxsForCharged_(iConfig.getParameter<int>("NumOfPUVtxsForCharged")),
+        fNumOfPUVtxsForCharged_(iConfig.getParameter<unsigned int>("NumOfPUVtxsForCharged")),
         fDzCutForChargedFromPUVtxs_(iConfig.getParameter<double>("DzCutForChargedFromPUVtxs")),
         fPtMaxCharged_(iConfig.getParameter<double>("PtMaxCharged")),
         fEtaMinUseDz_(iConfig.getParameter<double>("EtaMinUseDz")),
@@ -154,7 +154,7 @@ class PrimaryVertexAssignment {
   bool useTiming_;
   bool useVertexFit_;
   bool preferHighRanked_;
-  int fNumOfPUVtxsForCharged_;
+  unsigned int fNumOfPUVtxsForCharged_;
   double fDzCutForChargedFromPUVtxs_;
   double fPtMaxCharged_;
   double fEtaMinUseDz_;

CommonTools/RecoAlgos/python/sortedPFPrimaryVertices_cfi.py

@@ -22,7 +22,7 @@
     preferHighRanked = cms.bool(False),
 
     #options used in PUPPI
-    NumOfPUVtxsForCharged = cms.int32(-1),
+    NumOfPUVtxsForCharged = cms.uint32(0),
     DzCutForChargedFromPUVtxs = cms.double(0.2),
     PtMaxCharged = cms.double(-1),
     EtaMinUseDz = cms.double(-1),

CommonTools/RecoAlgos/python/sortedPrimaryVertices_cfi.py

@@ -22,7 +22,7 @@
     preferHighRanked = cms.bool(False),
 
     #options used in PUPPI
-    NumOfPUVtxsForCharged = cms.int32(-1),
+    NumOfPUVtxsForCharged = cms.uint32(0),
     DzCutForChargedFromPUVtxs = cms.double(0.2),
     PtMaxCharged = cms.double(-1),
     EtaMinUseDz = cms.double(-1),

CommonTools/RecoAlgos/src/PrimaryVertexAssignment.cc

@@ -71,7 +71,7 @@ std::pair<int, PrimaryVertexAssignment::Quality> PrimaryVertexAssignment::charge
   if (!vertices.empty())
     firstVertexDz = std::abs(track->dz(vertices.at(0).position()));
   // recover cases where the primary vertex is split
-  if (useVertexFit_ && (iVertex > 0) && (iVertex <= fNumOfPUVtxsForCharged_) &&
+  if (useVertexFit_ && (iVertex > 0) && (iVertex <= int(fNumOfPUVtxsForCharged_)) &&
       firstVertexDz < fDzCutForChargedFromPUVtxs_)
     return {0, PrimaryVertexAssignment::PrimaryDz};
 
@@ -121,14 +121,12 @@ std::pair<int, PrimaryVertexAssignment::Quality> PrimaryVertexAssignment::charge
           (!useTime or dtmin / timeReso < maxDtSigForPrimaryAssignment_))
         iVertex = vtxIdMinSignif;
     } else {
-      // consider only distances to first vertex for association of pileup vertices (originally used in PUPPI)
-      if ((vtxIdMinSignif >= 0) && (std::abs(track->eta()) > fEtaMinUseDz_))
-        iVertex =
-            ((firstVertexDz < maxDzForPrimaryAssignment_ and firstVertexDz / dzE < maxDzSigForPrimaryAssignment_ and
-              track->dzError() < maxDzErrorForPrimaryAssignment_) and
-             (!useTime or std::abs(time - vertices.at(0).t()) / timeReso < maxDtSigForPrimaryAssignment_))
-                ? 0
-                : vtxIdMinSignif;
+      // consider only distances to first vertex for association (originally used in PUPPI)
+      if ((vtxIdMinSignif >= 0) && (std::abs(track->eta()) > fEtaMinUseDz_) and
+          ((firstVertexDz < maxDzForPrimaryAssignment_ and firstVertexDz / dzE < maxDzSigForPrimaryAssignment_ and
+            track->dzError() < maxDzErrorForPrimaryAssignment_) and
+           (!useTime or std::abs(time - vertices.at(0).t()) / timeReso < maxDtSigForPrimaryAssignment_)))
+        iVertex = 0;
     }
   }
 

PhysicsTools/PatAlgos/python/slimming/miniAOD_tools.py

@@ -199,7 +199,6 @@ def _add_slimmedMETsNoHF(process):
 
     #  ==================  CHSMET 
     process.load("CommonTools.ParticleFlow.pfCHS_cff")
-    task.add(process.packedPrimaryVertexAssociationJME)
     task.add(process.pfCHS)
 
     from RecoMET.METProducers.pfMet_cfi import pfMet

PhysicsTools/PatAlgos/python/slimming/puppiForMET_cff.py

@@ -12,26 +12,21 @@ def makePuppies( process ):
 
 def makePuppiesFromMiniAOD( process, createScheduledSequence=False ):
     task = getPatAlgosToolsTask(process)
-    from CommonTools.ParticleFlow.pfCHS_cff import packedPrimaryVertexAssociationJME
-    setattr(process, "packedPrimaryVertexAssociationJME", packedPrimaryVertexAssociationJME.clone())
-    task.add(process.packedPrimaryVertexAssociationJME)
     process.load('CommonTools.PileupAlgos.Puppi_cff')
     task.add(process.puppi)
     process.puppi.candName = 'packedPFCandidates'
     process.puppi.clonePackedCands = True
     process.puppi.vertexName = 'offlineSlimmedPrimaryVertices'
     process.puppi.useExistingWeights = True
-    process.puppi.vertexAssociation = 'packedPrimaryVertexAssociationJME:original'
     task.add(process.puppiNoLep)
     process.puppiNoLep.candName = 'packedPFCandidates'
     process.puppiNoLep.clonePackedCands = True
     process.puppiNoLep.vertexName = 'offlineSlimmedPrimaryVertices'
     process.puppiNoLep.useExistingWeights = True
-    process.puppiNoLep.vertexAssociation = 'packedPrimaryVertexAssociationJME:original'
 
     #making a sequence for people running the MET tool in scheduled mode
     if createScheduledSequence:
-        puppiMETTask = cms.Task(process.packedPrimaryVertexAssociationJME, process.puppi, process.puppiNoLep)
+        puppiMETTask = cms.Task(process.puppi, process.puppiNoLep)
         setattr(process, "puppiMETTask", puppiMETTask)
         puppiMETSequence = cms.Sequence(puppiMETTask)
         setattr(process, "puppiMETSequence", puppiMETSequence)