From 7a9b6204ece7619cca657ec06f843699a63b2830 Mon Sep 17 00:00:00 2001
From: Steven Gardiner <gardiner@fnal.gov>
Date: Wed, 29 Oct 2025 21:28:30 -0500
Subject: [PATCH 1/2] Adjust the SpectralFunc class so that a missing spectral
 function can be optionally approximated by a spectral function for a
 different isotope of the same element. Enable this behavior in the AR25
 tunes.

---
 config/AR25_20i/ModelConfiguration.xml    | 10 +--
 config/AR25_22i/ModelConfiguration.xml    | 10 +--
 config/SpectralFunc.xml                   |  4 +
 src/Physics/NuclearState/SpectralFunc.cxx | 99 +++++++++++++++++++----
 src/Physics/NuclearState/SpectralFunc.h   |  5 ++
 5 files changed, 104 insertions(+), 24 deletions(-)

diff --git a/config/AR25_20i/ModelConfiguration.xml b/config/AR25_20i/ModelConfiguration.xml
index a91f7a616..63d8605d4 100644
--- a/config/AR25_20i/ModelConfiguration.xml
+++ b/config/AR25_20i/ModelConfiguration.xml
@@ -39,11 +39,11 @@ University of Liverpool
 
  <!--<param type="alg"  name="NuclearModel"> genie::FGMBodekRitchie/Default </param>-->
 
-  <!-- For 12C, 16O, and 40Ar, use the spectral function -->
-  <param type="alg" name="NuclearModel@Pdg=1000060120"> genie::SpectralFunc/Default </param>
-  <!--<param type="alg" name="NuclearModel@Pdg=1000260560"> genie::SpectralFunc/Default </param>-->
-  <param type="alg" name="NuclearModel@Pdg=1000080160"> genie::SpectralFunc/Default </param>
-  <param type="alg" name="NuclearModel@Pdg=1000180400"> genie::SpectralFunc/Default </param>
+  <!-- For C, O, and Ar targets, use the spectral function -->
+  <param type="alg" name="NuclearModel@Pdg=1000060120"> genie::SpectralFunc/ApproxElements </param>
+  <!--<param type="alg" name="NuclearModel@Pdg=1000260560"> genie::SpectralFunc/ApproxElements </param>-->
+  <param type="alg" name="NuclearModel@Pdg=1000080160"> genie::SpectralFunc/ApproxElements </param>
+  <param type="alg" name="NuclearModel@Pdg=1000180400"> genie::SpectralFunc/ApproxElements </param>
 
 <!--
   Example of specific model for specific nuclei
diff --git a/config/AR25_22i/ModelConfiguration.xml b/config/AR25_22i/ModelConfiguration.xml
index d83f2de98..3aa391911 100644
--- a/config/AR25_22i/ModelConfiguration.xml
+++ b/config/AR25_22i/ModelConfiguration.xml
@@ -39,11 +39,11 @@ University of Liverpool
 
  <!--<param type="alg"  name="NuclearModel"> genie::FGMBodekRitchie/Default </param>-->
 
-  <!-- For 12C, 16O, and 40Ar, use the spectral function -->
-  <param type="alg" name="NuclearModel@Pdg=1000060120"> genie::SpectralFunc/Default </param>
-  <!--<param type="alg" name="NuclearModel@Pdg=1000260560"> genie::SpectralFunc/Default </param>-->
-  <param type="alg" name="NuclearModel@Pdg=1000080160"> genie::SpectralFunc/Default </param>
-  <param type="alg" name="NuclearModel@Pdg=1000180400"> genie::SpectralFunc/Default </param>
+  <!-- For C, O, and Ar targets, use the spectral function -->
+  <param type="alg" name="NuclearModel@Pdg=1000060120"> genie::SpectralFunc/ApproxElements </param>
+  <!--<param type="alg" name="NuclearModel@Pdg=1000260560"> genie::SpectralFunc/ApproxElements </param>-->
+  <param type="alg" name="NuclearModel@Pdg=1000080160"> genie::SpectralFunc/ApproxElements </param>
+  <param type="alg" name="NuclearModel@Pdg=1000180400"> genie::SpectralFunc/ApproxElements </param>
 
 <!--
   Example of specific model for specific nuclei
diff --git a/config/SpectralFunc.xml b/config/SpectralFunc.xml
index 063a8f32c..e65543214 100644
--- a/config/SpectralFunc.xml
+++ b/config/SpectralFunc.xml
@@ -34,4 +34,8 @@ FermiMomentumTable      string   No         Table of Fermi momentum (kF) constan
 
   </param_set>
 
+  <param_set name="ApproxElements">
+    <param type="bool" name="AllowElementMatch"> true </param>
+  </param_set>
+
 </alg_conf>
diff --git a/src/Physics/NuclearState/SpectralFunc.cxx b/src/Physics/NuclearState/SpectralFunc.cxx
index 18ebfc222..b68086697 100644
--- a/src/Physics/NuclearState/SpectralFunc.cxx
+++ b/src/Physics/NuclearState/SpectralFunc.cxx
@@ -32,19 +32,9 @@
 #include "Framework/Messenger/Messenger.h"
 #include "Physics/NuclearState/SpectralFunc.h"
 #include "Framework/ParticleData/PDGCodes.h"
+#include "Framework/ParticleData/PDGUtils.h"
 #include "Framework/Numerical/RandomGen.h"
 
-namespace {
-
-  // TODO: Replace this with std::to_string when we switch to C++11
-  std::string replace_with_std_to_string(int an_integer) {
-    std::ostringstream oss;
-    oss << an_integer;
-    return oss.str();
-  }
-
-}
-
 using namespace genie;
 using namespace genie::constants;
 using namespace genie::controls;
@@ -198,6 +188,8 @@ void SpectralFunc::LoadConfig(void)
   }
 
   fDataPath = data_path;
+
+  this->GetParamDef( "AllowElementMatch", fAllowElementMatch, false );
 }
 //____________________________________________________________________________
 TH2D* SpectralFunc::SelectSpectralFunction(const Target& t) const
@@ -210,16 +202,95 @@ TH2D* SpectralFunc::SelectSpectralFunction(const Target& t) const
   if ( my_iter != fSpectralFunctionMap.end() ) return my_iter->second;
 
   // If not, attempt to build it
-  std::string target_pdg_string = replace_with_std_to_string( target_pdg );
-  std::string hitnuc_pdg_string = replace_with_std_to_string( hitnuc_pdg );
+  std::string target_pdg_string = std::to_string( target_pdg );
+  std::string hitnuc_pdg_string = std::to_string( hitnuc_pdg );
 
   RgKey sf_key( "SpectFuncTable@Pdg=" + target_pdg_string + "_" + hitnuc_pdg_string );
   std::string data_filename;
   this->GetParamDef( sf_key, data_filename, std::string() );
 
   if ( data_filename.empty() ) {
-    LOG("SpectralFunc", pERROR) << "** The spectral function for target "
+    LOG("SpectralFunc", pNOTICE) << "** The spectral function for target "
       << target_pdg << " and hit nuc PDG " << hitnuc_pdg << " isn't available";
+
+    // If element matching is enabled, then check for a spectral function with
+    // the same proton number Z as the requested nuclide.
+    if ( fAllowElementMatch ) {
+
+      // Get the proton number for the requested target nuclide
+      int tZ = genie::pdg::IonPdgCodeToZ( target_pdg );
+
+      // Not all configured spectral functions may have been built yet, so rely
+      // on the owned Registry to see if a suitable definition exists in the
+      // Algorithm configuration
+      const auto& conf_map = this->GetConfig().GetItemMap();
+
+      // Search the configuration for a matching proton number Z and matching
+      // hit nucleon PDG code. Store the nuclear PDG code for the result if
+      // we find a match.
+      int element_match_pdg;
+      auto iter = std::find_if( conf_map.cbegin(), conf_map.cend(),
+        [ tZ, hitnuc_pdg, &element_match_pdg ]( const auto& conf_pair ) -> bool {
+          const std::string prefix = "SpectFuncTable@Pdg=";
+          const std::size_t prefix_size = prefix.size();
+
+          // Check that the current config key begins with the exact prefix above
+          const auto& key = conf_pair.first;
+          // If not, then it cannot be a match
+          if ( key.compare(0, prefix_size, prefix) != 0 ) return false;
+
+          // Find an underscore following the prefix
+          size_t underscore_pos = key.find( '_', prefix_size );
+          // If there is no underscore, then the current key is not a match
+          if ( underscore_pos == std::string::npos ) return false;
+
+          // Extract the string between the prefix and the underscore. Assume
+          // it is a string representation of a nuclear PDG code
+          std::string tgt_pdg_str = key.substr( prefix_size,
+            underscore_pos - prefix_size );
+
+          // Also extract the remainder of the string following the underscore.
+          // Assume it is a string representation of a nucleon PDG code
+          std::string nucleon_pdg_str = key.substr( underscore_pos + 1u );
+
+          // Convert the string representations to integers
+          int tgt_pdg = std::stoi( tgt_pdg_str );
+          int nucleon_pdg = std::stoi( nucleon_pdg_str );
+
+          // Extract the proton number for the current nuclide
+          int target_Z = genie::pdg::IonPdgCodeToZ( tgt_pdg );
+
+          if ( hitnuc_pdg == nucleon_pdg && tZ == target_Z ) {
+            element_match_pdg = tgt_pdg;
+            return true;
+          }
+          return false;
+        }
+      );
+
+      // If a match was found, use it to approximate the spectral function for
+      // the current nuclide of interest
+      if ( iter != conf_map.end() ) {
+        LOG("SpectralFunc", pNOTICE) << "Approximating the spectral function"
+          << " for target " << target_pdg << " with the one for target "
+          << element_match_pdg;
+
+        // Use recursion to retrieve or build the spectral function for the
+        // nuclide with the matching proton number
+        genie::Target temp_tgt( element_match_pdg );
+        temp_tgt.SetHitNucPdg( hitnuc_pdg );
+        TH2D* sf_ptr = this->SelectSpectralFunction( temp_tgt );
+
+        // Add an entry for the current nuclide to the map so we don't have to do
+        // the Z-based search again
+        fSpectralFunctionMap[ std::make_pair( hitnuc_pdg, target_pdg ) ] = sf_ptr;
+        return sf_ptr;
+      }
+    }
+
+    // Fail with an error if the user requests a spectral function that
+    // is not available
+    LOG("SpectralFunc", pERROR) << "Failed to load spectral function data";
     std::exit( 1 );
   }
 
diff --git a/src/Physics/NuclearState/SpectralFunc.h b/src/Physics/NuclearState/SpectralFunc.h
index d7c57da4a..21006c0c0 100644
--- a/src/Physics/NuclearState/SpectralFunc.h
+++ b/src/Physics/NuclearState/SpectralFunc.h
@@ -81,6 +81,11 @@ class SpectralFunc : public NuclearModelI {
   /// The number of removal energy bins to use when making spectral function
   /// histograms
   int fNumYBins;
+
+  /// If this flag is true, then a cached spectral function with the same
+  /// proton number Z will be used when a perfect match is not available
+  /// for the requested nuclide.
+  bool fAllowElementMatch = false;
 };
 
 }      // genie namespace

From d577187ed833c5af5092eee8d142f7244c441b33 Mon Sep 17 00:00:00 2001
From: kjplows <kom.plows@gmail.com>
Date: Thu, 30 Oct 2025 14:51:03 -0500
Subject: [PATCH 2/2] Update destructor to protect against double free

---
 src/Physics/NuclearState/SpectralFunc.cxx | 8 +-------
 1 file changed, 1 insertion(+), 7 deletions(-)

diff --git a/src/Physics/NuclearState/SpectralFunc.cxx b/src/Physics/NuclearState/SpectralFunc.cxx
index b68086697..10d4395c4 100644
--- a/src/Physics/NuclearState/SpectralFunc.cxx
+++ b/src/Physics/NuclearState/SpectralFunc.cxx
@@ -53,13 +53,7 @@ SpectralFunc::SpectralFunc(string config) :
 //____________________________________________________________________________
 SpectralFunc::~SpectralFunc()
 {
-  // Delete the TH2D objects from the spectral function map
-  std::map<std::pair<int,int>, TH2D*>::iterator begin = fSpectralFunctionMap.begin();
-  std::map<std::pair<int,int>, TH2D*>::iterator end = fSpectralFunctionMap.end();
-  for (std::map<std::pair<int,int>, TH2D*>::iterator iter = begin; iter != end; ++iter) {
-    TH2D* hist = iter->second;
-    if ( hist ) delete hist;
-  }
+  fSpectralFunctionMap.clear();
 }
 //____________________________________________________________________________
 bool SpectralFunc::GenerateNucleon(const Target& target) const