From 5599279483abe0d257d5e04a18a9a7a069b00a2a Mon Sep 17 00:00:00 2001
From: ddobrigk <david.dobrigkeit.chinellato@cern.ch>
Date: Tue, 2 Sep 2025 15:09:12 +0200
Subject: [PATCH] [Common] Sanitize power-law mapping in corner cases

This PR addresses a corner case in the power-law mapping of data to MC used in the centrality framework that was leading to NaNs when specific sets of parameters were used.

Interestingly, even if the MC calibration mapping procedure has been power-law-based for many years, this corner case in which parameters led to NaNs in the intermediate evaluation has only appeared in light ion monte carlo. This adjustment will sanitize all instances of this effect also in the future, as the mapping function is adjusted in the evaluation code (no change needed in ccdb).
---
 Common/Tools/Multiplicity/MultModule.h | 33 ++++++++++++++++++--------
 1 file changed, 23 insertions(+), 10 deletions(-)

diff --git a/Common/Tools/Multiplicity/MultModule.h b/Common/Tools/Multiplicity/MultModule.h
index bed8817c124..f3122106de4 100644
--- a/Common/Tools/Multiplicity/MultModule.h
+++ b/Common/Tools/Multiplicity/MultModule.h
@@ -440,6 +440,7 @@ class MultModule
     internalOpts.mEnabledTables.resize(nTablesConst, 0);
 
     LOGF(info, "Configuring tables to generate");
+    LOGF(info, "Metadata information: isMC? %i", metadataInfo.isMC());
     const auto& workflows = context.services().template get<o2::framework::RunningWorkflowInfo const>();
 
     TString listOfRequestors[nTablesConst];
@@ -505,16 +506,6 @@ class MultModule
       listOfRequestors[kPVMults].Append(Form("%s ", "dependency check"));
     }
 
-    // capture the need for PYTHIA calibration in Pb-Pb runs
-    if (metadataInfo.isMC() && mRunNumber >= 544013 && mRunNumber <= 545367) {
-      internalOpts.generatorName.value = "PYTHIA";
-    }
-
-    // capture the need for PYTHIA calibration in light ion runs automatically
-    if (metadataInfo.isMC() && mRunNumber >= 564250 && mRunNumber <= 564472) {
-      internalOpts.generatorName.value = "PYTHIA";
-    }
-
     // list enabled tables
     for (int i = 0; i < nTablesConst; i++) {
       // printout to be improved in the future
@@ -1145,6 +1136,20 @@ class MultModule
   {
     if (bc.runNumber() != mRunNumberCentrality) {
       mRunNumberCentrality = bc.runNumber(); // mark that this run has been attempted already regardless of outcome
+      LOGF(info, "centrality loading procedure for timestamp=%llu, run number=%d", bc.timestamp(), bc.runNumber());
+
+      // capture the need for PYTHIA calibration in Pb-Pb runs
+      if (metadataInfo.isMC() && mRunNumber >= 544013 && mRunNumber <= 545367) {
+        LOGF(info, "This is MC for Pb-Pb. Setting generatorName automatically to PYTHIA");
+        internalOpts.generatorName.value = "PYTHIA";
+      }
+
+      // capture the need for PYTHIA calibration in light ion runs automatically
+      if (metadataInfo.isMC() && mRunNumber >= 564250 && mRunNumber <= 564472) {
+        LOGF(info, "This is MC for light ion runs. Setting generatorName automatically to PYTHIA");
+        internalOpts.generatorName.value = "PYTHIA";
+      }
+
       LOGF(info, "centrality loading procedure for timestamp=%llu, run number=%d", bc.timestamp(), bc.runNumber());
       TList* callst = nullptr;
       // Check if the ccdb path is a root file
@@ -1258,6 +1263,10 @@ class MultModule
       auto populateTable = [&](auto& table, struct CalibrationInfo& estimator, float multiplicity, bool isInelGt0) {
         const bool assignOutOfRange = internalOpts.embedINELgtZEROselection && !isInelGt0;
         auto scaleMC = [](float x, const float pars[6]) {
+          float core = ((pars[0] + pars[1] * std::pow(x, pars[2])) - pars[3]) / pars[4];
+          if (core < 0.0f) {
+            return 0.0f; // this should be marked as low multiplicity and not mapped, core^pars[5] would be NaN
+          }
           return std::pow(((pars[0] + pars[1] * std::pow(x, pars[2])) - pars[3]) / pars[4], 1.0f / pars[5]);
         };
 
@@ -1343,6 +1352,10 @@ class MultModule
       ConfigureCentralityRun2(ccdb, metadataInfo, firstbc);
 
       auto scaleMC = [](float x, const float pars[6]) {
+        float core = ((pars[0] + pars[1] * std::pow(x, pars[2])) - pars[3]) / pars[4];
+        if (core < 0.0f) {
+          return 0.0f; // this should be marked as low multiplicity and not mapped, core^pars[5] would be NaN
+        }
         return std::pow(((pars[0] + pars[1] * std::pow(x, pars[2])) - pars[3]) / pars[4], 1.0f / pars[5]);
       };