Fix numeric issues in PFCand scaling, add some debug output

Undo new line modification Remove nan check in pflinker Remove nan check from rescale
cms-sw · May 10, 2023 · f259152 · f259152
1 parent 9056c64
commit f259152
Show file tree

Hide file tree

Showing 3 changed files with 22 additions and 7 deletions.
diff --git a/DataFormats/ParticleFlowCandidate/src/PFCandidate.cc b/DataFormats/ParticleFlowCandidate/src/PFCandidate.cc
@@ -296,8 +296,10 @@ void PFCandidate::rescaleMomentum(double rescaleFactor) {
     throw cms::Exception(
         "NegativeScaling",
         "Scale factor " + std::to_string(rescaleFactor) + " is < 0. Cannot rescale momentum by this value");
-  float rescaleE = std::sqrt(p() * p() * (rescaleFactor * rescaleFactor - 1) / (energy() * energy()) + 1);
-  LorentzVector rescaledp4(rescaleFactor * px(), rescaleFactor * py(), rescaleFactor * pz(), rescaleE * energy());
+
+  float e = std::sqrt(p() * p() * rescaleFactor * rescaleFactor + mass() * mass());
+
+  LorentzVector rescaledp4(rescaleFactor * px(), rescaleFactor * py(), rescaleFactor * pz(), e);
   setP4(rescaledp4);
 }
 

diff --git a/RecoParticleFlow/PFProducer/plugins/PFLinker.cc b/RecoParticleFlow/PFProducer/plugins/PFLinker.cc
@@ -150,9 +150,6 @@ void PFLinker::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
       edm::Ptr<reco::PFCandidate> candPtr(pfCandidates, i);
       reco::PFCandidate cand(candPtr);
 
-      if (!(cand.energy() > 0))
-        continue;
-
       bool isphoton = cand.particleId() == reco::PFCandidate::gamma && cand.mva_nothing_gamma() > 0.;
       bool iselectron = cand.particleId() == reco::PFCandidate::e;
       // PFCandidates may have a valid MuonRef though they are not muons.

diff --git a/RecoParticleFlow/PFProducer/src/PFAlgo.cc b/RecoParticleFlow/PFProducer/src/PFAlgo.cc
@@ -2425,7 +2425,15 @@ void PFAlgo::createCandidatesHCAL(const reco::PFBlock& block,
       if (iTrack == corrTrack) {
         if (corrFact < 0.)
           corrFact = 0.;  // protect against negative scaling
-        (*pfCandidates_)[tmpi].rescaleMomentum(corrFact);
+        auto& candRescale = (*pfCandidates_)[tmpi];
+        LogTrace("PFAlgo|createCandidatesHCAL")
+            << "\tBefore rescaling: momentum " << candRescale.p() << " pT " << candRescale.pt() << " energy "
+            << candRescale.energy() << " mass " << candRescale.mass() << std::endl
+            << "\tTo rescale by " << corrFact << std::endl;
+        candRescale.rescaleMomentum(corrFact);
+        LogTrace("PFAlgo|createCandidatesHCAL")
+            << "\tRescaled candidate momentum " << candRescale.p() << " pT " << candRescale.pt() << " energy "
+            << candRescale.energy() << " mass " << candRescale.mass() << std::endl;
         trackMomentum *= corrFact;
       }
       chargedHadronsIndices.push_back(tmpi);
@@ -2508,7 +2516,15 @@ void PFAlgo::createCandidatesHCAL(const reco::PFBlock& block,
             double rescaleFactor = x(i) / hcalP[i];
             if (rescaleFactor < 0.)
               rescaleFactor = 0.;  // protect against negative scaling
-            (*pfCandidates_)[ich].rescaleMomentum(rescaleFactor);
+            auto& candRescale = (*pfCandidates_)[ich];
+            LogTrace("PFAlgo|createCandidatesHCAL")
+                << "\tBefore rescaling: momentum " << candRescale.p() << " pT " << candRescale.pt() << " energy "
+                << candRescale.energy() << " mass " << candRescale.mass() << std::endl
+                << "\tTo rescale by " << rescaleFactor << std::endl;
+            candRescale.rescaleMomentum(rescaleFactor);
+            LogTrace("PFAlgo|createCandidatesHCAL")
+                << "\tRescaled candidate momentum " << candRescale.p() << " pT " << candRescale.pt() << " energy "
+                << candRescale.energy() << " mass " << candRescale.mass() << std::endl;
 
             LogTrace("PFAlgo|createCandidatesHCAL")
                 << "\t\t\told p " << hcalP[i] << " new p " << x(i) << " rescale " << rescaleFactor;