Merge pull request #35309 from swagata87/EGM_EtaExtendedEle

[EGM] Eta-Extended-Electrons (EEE)

HLTrigger/Configuration/python/customizeHLTforCMSSW.py

@@ -20,6 +20,18 @@
 #                     pset.minGoodStripCharge = cms.PSet(refToPSet_ = cms.string('HLTSiStripClusterChargeCutNone'))
 #     return process
 
+# Eta Extended Electrons 
+def customiseFor35309(process):
+    for pset in process._Process__psets.values():
+        if hasattr(pset,'ComponentType'):
+            if (pset.ComponentType == 'CkfBaseTrajectoryFilter'):
+                if not hasattr(pset, 'highEtaSwitch'):
+                    pset.highEtaSwitch = cms.double(5.0)
+                if not hasattr(pset, 'minHitsAtHighEta'):
+                    pset.minHitsAtHighEta = cms.int32(5)
+
+    return process
+
 def customiseHCALFor2018Input(process):
     """Customise the HLT to run on Run 2 data/MC using the old readout for the HCAL barel"""
 
@@ -155,6 +167,7 @@ def customizeHLTforCMSSW(process, menuType="GRun"):
     # add call to action function in proper order: newest last!
     # process = customiseFor12718(process)
 
+    process = customiseFor35309(process)
     process = customiseFor35315(process)
     process = customiseFor35269(process)
 

RecoParticleFlow/PFProducer/interface/PFEGammaFilters.h

@@ -32,6 +32,7 @@ class PFEGammaFilters {
 
 private:
   bool passGsfElePreSelWithOnlyConeHadem(const reco::GsfElectron &) const;
+  bool thisEleIsNotAllowedInPF(const reco::GsfElectron &, bool) const;
 
   // Photon selections
   const float ph_Et_;
@@ -52,6 +53,7 @@ class PFEGammaFilters {
   const int ele_missinghits_;
   const float ele_ecalDrivenHademPreselCut_;
   const float ele_maxElePtForOnlyMVAPresel_;
+  const bool allowEEEinPF_;
   float ele_maxNtracks_;
   float ele_maxHcalE_;
   float ele_maxTrackPOverEele_;

RecoParticleFlow/PFProducer/python/particleFlow_cff.py

@@ -10,6 +10,10 @@
 
 particleFlowTmp = particleFlow.clone()
 
+## temporary for 12_1; EtaExtendedEles do not enter PF because ID/regression of EEEs are not ready yet
+## In 12_2, we expect to have EEE's ID/regression, then this switch can flip to True
+particleFlowTmp.PFEGammaFiltersParameters.allowEEEinPF = cms.bool(False)
+
 from Configuration.Eras.Modifier_pf_badHcalMitigationOff_cff import pf_badHcalMitigationOff
 pf_badHcalMitigationOff.toModify(particleFlowTmp.PFEGammaFiltersParameters,
                                  electron_protectionsForBadHcal = dict(enableProtections = False),

RecoParticleFlow/PFProducer/src/PFEGammaFilters.cc

@@ -36,7 +36,8 @@ PFEGammaFilters::PFEGammaFilters(const edm::ParameterSet& cfg)
       ele_noniso_mva_(cfg.getParameter<double>("electron_noniso_mvaCut")),
       ele_missinghits_(cfg.getParameter<unsigned int>("electron_missinghits")),
       ele_ecalDrivenHademPreselCut_(cfg.getParameter<double>("electron_ecalDrivenHademPreselCut")),
-      ele_maxElePtForOnlyMVAPresel_(cfg.getParameter<double>("electron_maxElePtForOnlyMVAPresel")) {
+      ele_maxElePtForOnlyMVAPresel_(cfg.getParameter<double>("electron_maxElePtForOnlyMVAPresel")),
+      allowEEEinPF_(cfg.getParameter<bool>("allowEEEinPF")) {
   auto const& eleProtectionsForBadHcal = cfg.getParameter<edm::ParameterSet>("electron_protectionsForBadHcal");
   auto const& eleProtectionsForJetMET = cfg.getParameter<edm::ParameterSet>("electron_protectionsForJetMET");
   auto const& phoProtectionsForBadHcal = cfg.getParameter<edm::ParameterSet>("photon_protectionsForBadHcal");
@@ -164,6 +165,14 @@ bool PFEGammaFilters::passElectronSelection(const reco::GsfElectron& electron,
     }
   }
 
+  //TEMPORARY hack for 12_1.
+  //Do not allow new EtaExtendedEle to enter PF, until ID, regression of EtaExtendedEle are in place.
+  //In 12_2, we expect to have EtaExtendedEle's ID/regression, then this switch can flip to True
+  //this is to be taken care of by EGM POG
+  //https://github.com/cms-sw/cmssw/issues/35374
+  if (thisEleIsNotAllowedInPF(electron, allowEEEinPF_))
+    passEleSelection = false;
+
   return passEleSelection && passGsfElePreSelWithOnlyConeHadem(electron);
 }
 
@@ -323,6 +332,14 @@ bool PFEGammaFilters::isElectronSafeForJetMET(const reco::GsfElectron& electron,
     isSafeForJetMET = false;
   }
 
+  //TEMPORARY hack for 12_1.
+  //Do not allow new EtaExtendedEle to be SafeForJetMET, until ID, regression of EtaExtendedEle are in place.
+  //In 12_2, we expect to have EtaExtendedEle's ID/regression, then this switch can flip to True
+  //this is to be taken care of by EGM POG
+  //https://github.com/cms-sw/cmssw/issues/35374
+  if (thisEleIsNotAllowedInPF(electron, allowEEEinPF_))
+    isSafeForJetMET = false;
+
   return isSafeForJetMET;
 }
 bool PFEGammaFilters::isPhotonSafeForJetMET(const reco::Photon& photon, const reco::PFCandidate& pfcand) const {
@@ -400,6 +417,18 @@ bool PFEGammaFilters::passGsfElePreSelWithOnlyConeHadem(const reco::GsfElectron&
     return passCutBased || passMVA;
 }
 
+bool PFEGammaFilters::thisEleIsNotAllowedInPF(const reco::GsfElectron& electron, bool allowEtaExtEleinPF) const {
+  bool returnVal = false;
+  if (!allowEtaExtEleinPF) {
+    const auto nHitGsf = electron.gsfTrack()->numberOfValidHits();
+    const auto absEleEta = std::abs(electron.eta());
+    if ((absEleEta > 2.5) && (nHitGsf < 5)) {
+      returnVal = true;
+    }
+  }
+  return returnVal;
+}
+
 void PFEGammaFilters::fillPSetDescription(edm::ParameterSetDescription& iDesc) {
   // Electron selection cuts
   iDesc.add<double>("electron_iso_pt", 10.0);
@@ -411,6 +440,7 @@ void PFEGammaFilters::fillPSetDescription(edm::ParameterSetDescription& iDesc) {
   iDesc.add<unsigned int>("electron_missinghits", 1);
   iDesc.add<double>("electron_ecalDrivenHademPreselCut", 0.15);
   iDesc.add<double>("electron_maxElePtForOnlyMVAPresel", 50.0);
+  iDesc.add<bool>("allowEEEinPF", false);
   {
     edm::ParameterSetDescription psd;
     psd.add<double>("maxNtracks", 3.0)->setComment("Max tracks pointing at Ele cluster");

TrackingTools/GsfTracking/python/CkfElectronCandidateMaker_cff.py

@@ -21,6 +21,8 @@
     maxConsecLostHits = 1,
     nSigmaMinPt = 5.0,
     minimumNumberOfHits = 5,
+    highEtaSwitch = 2.5,
+    minHitsAtHighEta = 3,
     maxCCCLostHits = 9999,
     minGoodStripCharge = dict(refToPSet_ = 'SiStripClusterChargeCutNone')
 )

TrackingTools/GsfTracking/python/GsfElectronFittingSmoother_cfi.py

@@ -4,5 +4,7 @@
 GsfElectronFittingSmoother = TrackingTools.TrackFitters.KFFittingSmoother_cfi.KFFittingSmoother.clone(
     ComponentName = 'GsfElectronFittingSmoother',
     Fitter        = 'GsfTrajectoryFitter',
-    Smoother      = 'GsfTrajectorySmoother'
+    Smoother      = 'GsfTrajectorySmoother',
+    MinNumberOfHitsHighEta = 3,
+    HighEtaSwitch = 2.5
 )

TrackingTools/TrackFitters/plugins/KFFittingSmootherESProducer.cc

@@ -28,6 +28,8 @@ namespace {
           theNoOutliersBeginEnd(conf.getParameter<bool>("NoOutliersBeginEnd")),
           theMinDof(conf.getParameter<int>("MinDof")),
           theMinNumberOfHits(conf.getParameter<int>("MinNumberOfHits")),
+          theMinNumberOfHitsHighEta(conf.getParameter<int>("MinNumberOfHitsHighEta")),
+          theHighEtaSwitch(conf.getParameter<double>("HighEtaSwitch")),
           rejectTracksFlag(conf.getParameter<bool>("RejectTracks")),
           breakTrajWith2ConsecutiveMissing(conf.getParameter<bool>("BreakTrajWith2ConsecutiveMissing")),
           noInvalidHitsBeginEnd(conf.getParameter<bool>("NoInvalidHitsBeginEnd")) {}
@@ -40,6 +42,8 @@ namespace {
     int theMinDof;
 
     int theMinNumberOfHits;
+    int theMinNumberOfHitsHighEta;
+    double theHighEtaSwitch;
     bool rejectTracksFlag;
     bool breakTrajWith2ConsecutiveMissing;
     bool noInvalidHitsBeginEnd;
@@ -62,12 +66,34 @@ namespace {
       desc.add<int>("MinDof", 2);
       desc.add<bool>("NoOutliersBeginEnd", false);
       desc.add<int>("MinNumberOfHits", 5);
+      desc.add<int>("MinNumberOfHitsHighEta", 5);
+      desc.add<double>("HighEtaSwitch", 5.0);
       desc.add<bool>("RejectTracks", true);
       desc.add<bool>("BreakTrajWith2ConsecutiveMissing", true);
       desc.add<bool>("NoInvalidHitsBeginEnd", true);
       desc.add<double>("LogPixelProbabilityCut", 0);
     }
 
+    int getNhitCutValue(const Trajectory& t,
+                        double theHighEtaSwitch,
+                        int theMinNumberOfHitsHighEta,
+                        int theMinNumberOfHits) const {
+      double sinhTrajEta2 = std::numeric_limits<double>::max();
+      if (!t.empty() && t.isValid()) {
+        /* in principle we can access eta() and check it w.r.t theHighEtaSwitch.
+	   but eta() is expensive, so we are making use of the following relation
+	   sinh(eta) = pz/pt (will square on both side to get rid of sign)
+	 */
+        double pt = t.lastMeasurement().updatedState().freeTrajectoryState()->momentum().perp();
+        double pz = t.lastMeasurement().updatedState().freeTrajectoryState()->momentum().z();
+        sinhTrajEta2 = (pz * pz) / (pt * pt);
+      }
+      double myEtaSwitch = sinh(theHighEtaSwitch);
+      const auto thisHitCut =
+          sinhTrajEta2 > (myEtaSwitch * myEtaSwitch) ? theMinNumberOfHitsHighEta : theMinNumberOfHits;
+      return thisHitCut;
+    }
+
     Trajectory fitOne(const Trajectory& t, fitType type) const override;
     Trajectory fitOne(const TrajectorySeed& aSeed,
                       const RecHitContainer& hits,
@@ -93,14 +119,16 @@ namespace {
         : KFFittingSmootherParam(other), theFitter(aFitter.clone()), theSmoother(aSmoother.clone()) {}
 
     Trajectory smoothingStep(Trajectory&& fitted) const {
+      const auto thisHitCut = getNhitCutValue(fitted, theHighEtaSwitch, theMinNumberOfHitsHighEta, theMinNumberOfHits);
+
       if (theEstimateCut > 0) {
         // remove "outlier" at the end of Traj
         while (
-            !fitted.empty() && fitted.foundHits() >= theMinNumberOfHits &&
+            !fitted.empty() && fitted.foundHits() >= thisHitCut &&
             (!fitted.lastMeasurement().recHitR().isValid() || (fitted.lastMeasurement().recHitR().det() != nullptr &&
                                                                fitted.lastMeasurement().estimate() > theEstimateCut)))
           fitted.pop();
-        if (fitted.foundHits() < theMinNumberOfHits)
+        if (fitted.foundHits() < thisHitCut)
           return Trajectory();
       }
       return theSmoother->trajectory(fitted);
@@ -199,11 +227,14 @@ namespace {
 #endif
 
       bool hasNaN = false;
-      if (!smoothed.isValid() || (hasNaN = !checkForNans(smoothed)) || (smoothed.foundHits() < theMinNumberOfHits)) {
+      const auto thisHitCut =
+          getNhitCutValue(smoothed, theHighEtaSwitch, theMinNumberOfHitsHighEta, theMinNumberOfHits);
+
+      if (!smoothed.isValid() || (hasNaN = !checkForNans(smoothed)) || (smoothed.foundHits() < thisHitCut)) {
         if (hasNaN)
           edm::LogWarning("TrackNaN") << "Track has NaN or the cov is not pos-definite";
-        if (smoothed.foundHits() < theMinNumberOfHits)
-          LogTrace("TrackFitters") << "smoothed.foundHits()<theMinNumberOfHits";
+        if (smoothed.foundHits() < thisHitCut)
+          LogTrace("TrackFitters") << "smoothed.foundHits()<thisHitCut";
         DPRINT("TrackFitters") << "smoothed invalid => trajectory rejected with nhits/chi2 " << smoothed.foundHits()
                                << '/' << smoothed.chiSquared() << "\n";
         if (rejectTracksFlag) {

TrackingTools/TrajectoryFiltering/interface/MinHitsTrajectoryFilter.h

@@ -12,11 +12,19 @@
 
 class MinHitsTrajectoryFilter final : public TrajectoryFilter {
 public:
-  explicit MinHitsTrajectoryFilter(int minHits = 5, int seedPairPenalty = 0)
-      : theMinHits(minHits), theSeedPairPenalty(seedPairPenalty) {}
+  explicit MinHitsTrajectoryFilter(int minHits = 5,
+                                   double highEtaSwitch = 5.0,
+                                   int minHitsAtHighEta = 5,
+                                   int seedPairPenalty = 0)
+      : theMinHits(minHits),
+        theHighEtaSwitch(highEtaSwitch),
+        theMinHitsAtHighEta(minHitsAtHighEta),
+        theSeedPairPenalty(seedPairPenalty) {}
 
   MinHitsTrajectoryFilter(const edm::ParameterSet& pset, edm::ConsumesCollector& iC)
       : theMinHits(pset.getParameter<int>("minimumNumberOfHits")),
+        theHighEtaSwitch(pset.getParameter<double>("highEtaSwitch")),
+        theMinHitsAtHighEta(pset.getParameter<int>("minHitsAtHighEta")),
         theSeedPairPenalty(pset.getParameter<int>("seedPairPenalty")) {}
 
   bool qualityFilter(const Trajectory& traj) const override { return QF<Trajectory>(traj); }
@@ -30,6 +38,8 @@ class MinHitsTrajectoryFilter final : public TrajectoryFilter {
   inline edm::ParameterSetDescription getFilledConfigurationDescription() {
     edm::ParameterSetDescription desc;
     desc.add<int>("minimumNumberOfHits", 5);
+    desc.add<double>("highEtaSwitch", 5.0);
+    desc.add<int>("minHitsAtHighEta", 5);
     desc.add<int>("seedPairPenalty", 0);
     return desc;
   }
@@ -38,10 +48,24 @@ class MinHitsTrajectoryFilter final : public TrajectoryFilter {
   template <class T>
   bool QF(const T& traj) const {
     int seedPenalty = (2 == traj.seedNHits()) ? theSeedPairPenalty : 0;  // increase by one if seed-doublet...
-    return (traj.foundHits() >= theMinHits + seedPenalty);
+    bool passed = false;
+    double pt = traj.lastMeasurement().updatedState().freeTrajectoryState()->momentum().perp();
+    double pz = traj.lastMeasurement().updatedState().freeTrajectoryState()->momentum().z();
+    double sinhTrajEta2 = (pz * pz) / (pt * pt);
+    double myEtaSwitch = sinh(theHighEtaSwitch);
+    if (sinhTrajEta2 < (myEtaSwitch * myEtaSwitch)) {
+      if (traj.foundHits() >= theMinHits + seedPenalty)
+        passed = true;
+    } else {  //absTrajEta>theHighEtaSwitch, so apply relaxed cuts
+      if (traj.foundHits() >= theMinHitsAtHighEta + seedPenalty)
+        passed = true;
+    }
+    return passed;
   }
 
   int theMinHits;
+  double theHighEtaSwitch;
+  int theMinHitsAtHighEta;
   int theSeedPairPenalty;
 };
 

TrackingTools/TrajectoryFiltering/python/TrajectoryFilter_cff.py

@@ -8,6 +8,8 @@
 #--- Cuts applied to completed trajectory
 # At least this many hits (counting matched hits as 1)
     minimumNumberOfHits = cms.int32(5),
+    highEtaSwitch = cms.double(5.0),
+    minHitsAtHighEta = cms.int32(5),
 # add this if seed is a Pair  (opposed to a triplet)
     seedPairPenalty = cms.int32(0),
 # What is this ?
@@ -68,7 +70,9 @@
 )
 MinHitsTrajectoryFilter_block = cms.PSet(
     ComponentType = cms.string('MinHitsTrajectoryFilter'),
-    minimumNumberOfHits = cms.int32(5)
+    minimumNumberOfHits = cms.int32(5),
+    highEtaSwitch = cms.double(5.0),
+    minHitsAtHighEta = cms.int32(5)
 )
 MinPtTrajectoryFilter_block = cms.PSet(
     ComponentType = cms.string('MinPtTrajectoryFilter'),