diff --git a/.gitignore b/.gitignore
index 8877fff55..a374c129a 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,2 +1,3 @@
 build
 test_runs
+/src/MicroBooNE/scriptholder/
diff --git a/data/MicroBooNE/CC0pi_2025/2D/Extracted2DCrossSection_MicroBooNE_CC0PiSelection_BinningScheme1.root b/data/MicroBooNE/CC0pi_2025/2D/Extracted2DCrossSection_MicroBooNE_CC0PiSelection_BinningScheme1.root
new file mode 100644
index 000000000..a6b53a160
Binary files /dev/null and b/data/MicroBooNE/CC0pi_2025/2D/Extracted2DCrossSection_MicroBooNE_CC0PiSelection_BinningScheme1.root differ
diff --git a/data/MicroBooNE/CC0pi_2025/2D/bin_defsCC0pi_2025.txt b/data/MicroBooNE/CC0pi_2025/2D/bin_defsCC0pi_2025.txt
new file mode 100644
index 000000000..af92e66b0
--- /dev/null
+++ b/data/MicroBooNE/CC0pi_2025/2D/bin_defsCC0pi_2025.txt
@@ -0,0 +1,85 @@
+Muon2D
+stv_tree
+44
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.100 && mc_p3_mu.Mag() < 0.240 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < -0.550"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.100 && mc_p3_mu.Mag() < 0.240 && mc_p3_mu.CosTheta() >= -0.550 && mc_p3_mu.CosTheta() < 0.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.100 && mc_p3_mu.Mag() < 0.240 && mc_p3_mu.CosTheta() >= 0.000 && mc_p3_mu.CosTheta() < 0.450"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.100 && mc_p3_mu.Mag() < 0.240 && mc_p3_mu.CosTheta() >= 0.450 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < -0.550"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= -0.550 && mc_p3_mu.CosTheta() < -0.250"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= -0.250 && mc_p3_mu.CosTheta() < 0.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= 0.000 && mc_p3_mu.CosTheta() < 0.250"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= 0.250 && mc_p3_mu.CosTheta() < 0.450"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= 0.450 && mc_p3_mu.CosTheta() < 0.700"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.240 && mc_p3_mu.Mag() < 0.300 && mc_p3_mu.CosTheta() >= 0.700 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < -0.400"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= -0.400 && mc_p3_mu.CosTheta() < -0.100"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= -0.100 && mc_p3_mu.CosTheta() < 0.100"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= 0.100 && mc_p3_mu.CosTheta() < 0.350"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= 0.350 && mc_p3_mu.CosTheta() < 0.500"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= 0.500 && mc_p3_mu.CosTheta() < 0.700"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= 0.700 && mc_p3_mu.CosTheta() < 0.850"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.300 && mc_p3_mu.Mag() < 0.380 && mc_p3_mu.CosTheta() >= 0.850 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < 0.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= 0.000 && mc_p3_mu.CosTheta() < 0.500"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= 0.500 && mc_p3_mu.CosTheta() < 0.650"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= 0.650 && mc_p3_mu.CosTheta() < 0.800"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= 0.800 && mc_p3_mu.CosTheta() < 0.920"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.380 && mc_p3_mu.Mag() < 0.480 && mc_p3_mu.CosTheta() >= 0.920 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < 0.500"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= 0.500 && mc_p3_mu.CosTheta() < 0.650"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= 0.650 && mc_p3_mu.CosTheta() < 0.800"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= 0.800 && mc_p3_mu.CosTheta() < 0.875"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= 0.875 && mc_p3_mu.CosTheta() < 0.950"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.480 && mc_p3_mu.Mag() < 0.700 && mc_p3_mu.CosTheta() >= 0.950 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.700 && mc_p3_mu.Mag() < 0.850 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < 0.875"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.700 && mc_p3_mu.Mag() < 0.850 && mc_p3_mu.CosTheta() >= 0.875 && mc_p3_mu.CosTheta() < 0.950"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.700 && mc_p3_mu.Mag() < 0.850 && mc_p3_mu.CosTheta() >= 0.950 && mc_p3_mu.CosTheta() < 1.000"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.850 && mc_p3_mu.Mag() < 2.000 && mc_p3_mu.CosTheta() >= -1.000 && mc_p3_mu.CosTheta() < 0.900"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.850 && mc_p3_mu.Mag() < 2.000 && mc_p3_mu.CosTheta() >= 0.900 && mc_p3_mu.CosTheta() < 0.950"
+0 0 "mc_is_cc0pi_signal && mc_p3_mu.Mag() >= 0.850 && mc_p3_mu.Mag() < 2.000 && mc_p3_mu.CosTheta() >= 0.950 && mc_p3_mu.CosTheta() < 1.000"
+1 -1 "category == 0"
+1 -1 "category == 5"
+1 -1 "category == 7"
+1 -1 "category == 8"
+1 -1 "category == 9"
+1 -1 "category == 10"
+1 -1 "category == 11"
+37
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.100 && p3_mu.Mag() < 0.240 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < -0.550"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.100 && p3_mu.Mag() < 0.240 && p3_mu.CosTheta() >= -0.550 && p3_mu.CosTheta() < 0.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.100 && p3_mu.Mag() < 0.240 && p3_mu.CosTheta() >= 0.000 && p3_mu.CosTheta() < 0.450"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.100 && p3_mu.Mag() < 0.240 && p3_mu.CosTheta() >= 0.450 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < -0.550"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= -0.550 && p3_mu.CosTheta() < -0.250"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= -0.250 && p3_mu.CosTheta() < 0.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= 0.000 && p3_mu.CosTheta() < 0.250"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= 0.250 && p3_mu.CosTheta() < 0.450"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= 0.450 && p3_mu.CosTheta() < 0.700"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.240 && p3_mu.Mag() < 0.300 && p3_mu.CosTheta() >= 0.700 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < -0.400"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= -0.400 && p3_mu.CosTheta() < -0.100"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= -0.100 && p3_mu.CosTheta() < 0.100"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= 0.100 && p3_mu.CosTheta() < 0.350"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= 0.350 && p3_mu.CosTheta() < 0.500"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= 0.500 && p3_mu.CosTheta() < 0.700"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= 0.700 && p3_mu.CosTheta() < 0.850"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.300 && p3_mu.Mag() < 0.380 && p3_mu.CosTheta() >= 0.850 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < 0.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= 0.000 && p3_mu.CosTheta() < 0.500"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= 0.500 && p3_mu.CosTheta() < 0.650"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= 0.650 && p3_mu.CosTheta() < 0.800"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= 0.800 && p3_mu.CosTheta() < 0.920"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.380 && p3_mu.Mag() < 0.480 && p3_mu.CosTheta() >= 0.920 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < 0.500"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= 0.500 && p3_mu.CosTheta() < 0.650"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= 0.650 && p3_mu.CosTheta() < 0.800"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= 0.800 && p3_mu.CosTheta() < 0.875"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= 0.875 && p3_mu.CosTheta() < 0.950"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.480 && p3_mu.Mag() < 0.700 && p3_mu.CosTheta() >= 0.950 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.700 && p3_mu.Mag() < 0.850 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < 0.875"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.700 && p3_mu.Mag() < 0.850 && p3_mu.CosTheta() >= 0.875 && p3_mu.CosTheta() < 0.950"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.700 && p3_mu.Mag() < 0.850 && p3_mu.CosTheta() >= 0.950 && p3_mu.CosTheta() < 1.000"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.850 && p3_mu.Mag() < 2.000 && p3_mu.CosTheta() >= -1.000 && p3_mu.CosTheta() < 0.900"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.850 && p3_mu.Mag() < 2.000 && p3_mu.CosTheta() >= 0.900 && p3_mu.CosTheta() < 0.950"
+0 0 "sel_CC0pi && sel_muon_contained && p3_mu.Mag() >= 0.850 && p3_mu.Mag() < 2.000 && p3_mu.CosTheta() >= 0.950 && p3_mu.CosTheta() < 1.000"
diff --git a/src/FCN/SampleList.cxx b/src/FCN/SampleList.cxx
index 0eab78a45..3332912bb 100644
--- a/src/FCN/SampleList.cxx
+++ b/src/FCN/SampleList.cxx
@@ -186,6 +186,8 @@
 #include "MicroBooNE_CC1MuNp_XSec_1D_nu.h"
 #include "MicroBooNE_CC1ENp_XSec_1D_nu.h"
 #include "MicroBooNE_CCInc_XSec_2DPcos_nu.h"
+#include "MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h"
+
 #endif
 
 #ifdef MINERvA_ENABLED
@@ -1146,7 +1148,9 @@ MeasurementBase *CreateSample(nuiskey samplekey) {
                !name.compare("MicroBooNE_CC1Mu1p_XSec_1DECal_nu") ||
                !name.compare("MicroBooNE_CC1Mu1p_XSec_1DEQE_nu")) {
       return (new MicroBooNE_CC1Mu1p_XSec_1D_nu(samplekey));
-    } else
+    }else if (!name.compare("MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D")){
+    return (new MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D(samplekey));
+    }else
 
 #endif
 #ifdef MINERvA_ENABLED
diff --git a/src/MicroBooNE/CMakeLists.txt b/src/MicroBooNE/CMakeLists.txt
index 71bc689f8..652b54906 100644
--- a/src/MicroBooNE/CMakeLists.txt
+++ b/src/MicroBooNE/CMakeLists.txt
@@ -22,6 +22,7 @@ set(MicroBooNE_Impl_Files
   MicroBooNE_CC1ENp_XSec_1D_nu.cxx
   MicroBooNE_CC1Mu2p_XSec_1D_nu.cxx
   MicroBooNE_CC1Mu1p_XSec_1D_nu.cxx
+  MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.cxx
   MicroBooNE_SignalDef.cxx
 )
 
@@ -30,6 +31,7 @@ set(MicroBooNE_Hdr_Files
   MicroBooNE_CC1MuNp_XSec_1D_nu.h
   MicroBooNE_CC1Mu2p_XSec_1D_nu.h
   MicroBooNE_CC1Mu1p_XSec_1D_nu.h
+ MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h
   MicroBooNE_SignalDef.h
 )
 
@@ -41,4 +43,4 @@ install(TARGETS MicroBooNE
     LIBRARY DESTINATION lib/
     PUBLIC_HEADER DESTINATION include)
 
-add_library(NUIS::MicroBooNE ALIAS MicroBooNE)
\ No newline at end of file
+add_library(NUIS::MicroBooNE ALIAS MicroBooNE)
diff --git a/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.cxx b/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.cxx
new file mode 100644
index 000000000..3080d5317
--- /dev/null
+++ b/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.cxx
@@ -0,0 +1,414 @@
+// Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
+
+/*******************************************************************************
+*    This file is part of NUISANCE.
+*
+*    NUISANCE is free software: you can redistribute it and/or modify
+*    it under the terms of the GNU General Public License as published by
+*    the Free Software Foundation, either version 3 of the License, or
+*    (at your option) any later version.
+*
+*    NUISANCE is distributed in the hope that it will be useful,
+*    but WITHOUT ANY WARRANTY; without even the implied warranty of
+*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+*    GNU General Public License for more details.
+*
+*    You should have received a copy of the GNU General Public License
+*    along with NUISANCE.  If not, see <http://www.gnu.org/licenses/>.
+*******************************************************************************/
+
+#include <fstream>
+#include <set>
+
+#include "MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h"
+#include "MicroBooNE_SignalDef.h"
+#include "TH2D.h"
+#include "TMatrixD.h"
+
+namespace {
+  constexpr int MASS_NUMBER_40AR = 40;
+  constexpr int PROTON = 2212;
+  constexpr int MU_MINUS = 13;
+
+  constexpr double TARGET_MASS = 37.215526; // 40Ar, GeV
+  constexpr double NEUTRON_MASS = 0.93956541; // GeV
+  constexpr double PROTON_MASS = 0.93827208; // GeV
+  constexpr double MUON_MASS = 0.10565837; // GeV
+  constexpr double BINDING_ENERGY = 0.02478; // 40Ar, GeV
+
+/*
+  void compute_stvs( const TVector3& p3mu, const TVector3& p3p,
+    double& delta_pT, double& delta_phiT, double& delta_alphaT,
+    double& delta_pL, double& pn, double& delta_pTx, double& delta_pTy )
+  {
+    delta_pT = (p3mu + p3p).Perp();
+
+    delta_phiT = std::acos( (-p3mu.X()*p3p.X() - p3mu.Y()*p3p.Y())
+      / (p3mu.XYvector().Mod() * p3p.XYvector().Mod()) );
+
+    TVector2 delta_pT_vec = (p3mu + p3p).XYvector();
+    delta_alphaT = std::acos( (-p3mu.X()*delta_pT_vec.X()
+      - p3mu.Y()*delta_pT_vec.Y())
+      / (p3mu.XYvector().Mod() * delta_pT_vec.Mod()) );
+
+    double Emu = std::sqrt(std::pow(MUON_MASS, 2) + p3mu.Mag2());
+    double Ep = std::sqrt(std::pow(PROTON_MASS, 2) + p3p.Mag2());
+    double R = TARGET_MASS + p3mu.Z() + p3p.Z() - Emu - Ep;
+
+    // Estimated mass of the final remnant nucleus (CCQE assumption)
+    double mf = TARGET_MASS - NEUTRON_MASS + BINDING_ENERGY;
+    delta_pL = 0.5*R - (std::pow(mf, 2) + std::pow(delta_pT, 2)) / (2.*R);
+
+    pn = std::sqrt( std::pow(delta_pL, 2) + std::pow(delta_pT, 2) );
+
+    // Components of the 2D delta_pT vector (see arXiv:1910.08658)
+
+    // We assume that the neutrino travels along the +z direction (also done
+    // in the other expressions above)
+    TVector3 zUnit( 0., 0., 1. );
+
+    // Defines the x direction for the components of the delta_pT vector
+    TVector2 xTUnit = zUnit.Cross( p3mu ).XYvector().Unit();
+
+    delta_pTx = xTUnit.X()*delta_pT_vec.X() + xTUnit.Y()*delta_pT_vec.Y();
+
+    // Defines the y direction for the components of the delta_T vector
+    TVector2 yTUnit = ( -p3mu ).XYvector().Unit();
+
+    delta_pTy = yTUnit.X()*delta_pT_vec.X() + yTUnit.Y()*delta_pT_vec.Y();
+  }
+
+  double compute_stvs( FitEvent* ev, const std::string& var ) {
+    TVector3 p3mu = ev->GetHMFSParticle( MU_MINUS )->fP.Vect();
+    TVector3 p3p = ev->GetHMFSParticle( PROTON )->fP.Vect();
+
+    // Convert to GeV
+    p3mu *= 1e-3;
+    p3p *= 1e-3;
+
+    double delta_pT, delta_phiT, delta_alphaT, delta_pL;
+    double pn, delta_pTx, delta_pTy;
+
+    compute_stvs( p3mu, p3p, delta_pT, delta_phiT, delta_alphaT,
+      delta_pL, pn, delta_pTx, delta_pTy );
+
+    if ( var == "delta_pT" ) return delta_pT;
+    else if ( var == "delta_phiT" ) return delta_phiT;
+    else if ( var == "delta_alphaT" ) return delta_alphaT;
+    else if ( var == "delta_pL" ) return delta_pL;
+    else if ( var == "pn" ) return pn;
+    else if ( var == "delta_pTx" ) return delta_pTx;
+    else if ( var == "delta_pTy" ) return delta_pTy;
+    return 0.;
+  }
+*/
+
+
+}
+
+                                            
+  MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D( nuiskey samplekey )
+{
+  fSettings = LoadSampleSettings( samplekey );
+  std::string name = fSettings.GetS( "name" );
+
+  // The main histograms use the bin number on the x-axis
+  fSettings.SetXTitle( "bin number" );
+
+  if ( !name.compare("MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D") ) {
+    fSettings.SetYTitle( "#sigma (cm^{2}/^{40}Ar)" );
+  }
+  else {
+    assert( false );
+  }
+
+  // Sample overview ---------------------------------------------------
+  std::string descrip = name + " sample.\n" \
+                        "Target: Ar\n" \
+                        "Flux: BNB FHC numu\n" \
+                        "reference: https://arxiv.org/abs/2507.00921\n" \
+                        "Signal: CC0pi (yr 2025)\n";
+
+  fSettings.SetDescription( descrip );
+  fSettings.SetTitle( name );
+  fSettings.SetAllowedTypes( "FULL", "FIX/FULL" );
+  fSettings.SetEnuRange( 0.0, 10 ); /// using range we flux intergered 
+  fSettings.DefineAllowedTargets( "Ar" );
+  fSettings.DefineAllowedSpecies( "numu" );
+  FinaliseSampleSettings();
+
+  // Scale factor for the flux-averaged total cross section
+  // (10^{-38} cm^2 / Ar) in each bin
+  fScaleFactor = GetEventHistogram()->Integral( "width" )
+    * MASS_NUMBER_40AR / fNEvents / TotalIntegratedFlux();
+
+  // Get bin definitions
+  this->LoadBinDefinitions();
+
+  std::string data_file_name( FitPar::GetDataBase()
+    + "/MicroBooNE/CC0pi_2025/2D/Extracted2DCrossSection_MicroBooNE_CC0PiSelection_BinningScheme1.root" );
+
+  // Load the additional smearing matrix
+  TMatrixD* A_C = StatUtils::GetMatrixFromRootFile( data_file_name,
+    "AC_Matrix" );
+  fAddSmear.reset( A_C );
+
+  // Load the measured data points
+  fDataHist = PlotUtils::GetTH1DFromRootFile( data_file_name,
+    "h_ExtractedCrossSectionUnfoldedData" );
+
+  fDataHist->SetNameTitle( (fSettings.GetName() + "_data").c_str(),
+    fSettings.GetFullTitles().c_str() );
+
+  // Also retrieve the total covariance matrix and its inverse (pre-computed in
+  // the data release ROOT file for convenience)
+  //this->SetCovarFromRootFile( data_file_name, "cov_total" );
+  fFullCovar = StatUtils::GetCovarFromRootFile( data_file_name, "TotalCov_AfterUnfolding_Norm" );
+  //covar = StatUtils::GetCovarFromRootFile( data_file_name,
+  //  "inverse_cov_total" );
+
+  //fDecomp = StatUtils::GetDecomp( fFullCovar );
+  TDecompChol chol( *fFullCovar );
+  chol.Decompose();
+  fDecomp = new TMatrixDSym( fFullCovar->GetNrows(),
+    chol.GetU().GetMatrixArray(), "" );
+
+  // Push the diagonals of fFullCovar onto the data histogram
+  StatUtils::SetDataErrorFromCov( fDataHist, fFullCovar, 1.0, false );
+
+  // Setup fMCHist from data
+  fMCHist = dynamic_cast< TH1D* >( fDataHist->Clone() );
+  fMCHist->SetNameTitle( (fSettings.GetName() + "_MC").c_str(),
+   fSettings.GetFullTitles().c_str() );
+  fMCHist->Reset();
+
+  fMCStat = dynamic_cast< TH1D* >( fMCHist->Clone() );
+  fMCStat->Reset();
+
+  // Since we're using a 1D histogram with bin number along the x-axis, it
+  // doesn't make sense to subdivide the bins. Rather than doing that
+  // automatically, I just copy the original data histogram binning here. Thus,
+  // MCFine ends up using the same bins as regular MC.
+  fMCFine = dynamic_cast< TH1D* >( fMCHist->Clone() );
+  fMCFine->SetNameTitle( (fSettings.GetName() + "_MC_FINE").c_str(),
+    fSettings.GetFullTitles().c_str() );
+  fMCFine->Reset();
+
+  // Set up the MC modes histogram
+  fMCHist_Modes = new TrueModeStack( (fSettings.GetName() + "_MODES").c_str(),
+    "True Channels", fMCHist );
+  fMCHist_Modes->SetTitleX( fDataHist->GetXaxis()->GetTitle() );
+  fMCHist_Modes->SetTitleY( fDataHist->GetYaxis()->GetTitle() );
+  this->SetAutoProcessTH1( fMCHist_Modes, kCMD_Reset, kCMD_Norm, kCMD_Write );
+
+  //this->FinaliseMeasurement();
+}
+
+
+bool MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::isSignal( FitEvent* event ) {
+   return SignalDef::MicroBooNE::isCC1Mu0pi_2025(event, EnuMin, EnuMax);
+}
+
+void MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::FillEventVariables( FitEvent* event ) {
+
+  // Clear out the vector of passing bins, which may have already been filled
+  // for the previous event
+  fPassingBins.clear();
+
+  if ( event->NumFSParticle(MU_MINUS) == 0 ) return;
+
+  // Loop over each of the bin definitions. Keep track of the bins that
+  // pass all cuts (and should thus be filled in this event)
+  for ( size_t b = 0u; b < fBinDefinitions.size(); ++b ) {
+
+    // Start out by assuming that the event passes all cuts
+    bool passed_cuts = true;
+
+    // Apply all cuts
+    const auto& my_cuts = fBinDefinitions.at( b );
+    for ( const auto& cut : my_cuts ) {
+      passed_cuts &= cut.evaluate( event );
+    }
+
+    // If all cuts were passed, add this bin to the vector of indices for bins
+    // that should be filled
+    if ( passed_cuts ) fPassingBins.push_back( b );
+
+  } // loop over bins
+
+}
+/////////////////////////////////////////////////////////////////////
+//
+/////////////////////////////////////////////////////////////////////
+void MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::LoadBinDefinitions() {
+  std::string binning_file_name( FitPar::GetDataBase()
+    + "/MicroBooNE/CC0pi_2025/2D/bin_defsCC0pi_2025.txt" );
+
+  std::ifstream bin_file( binning_file_name );
+  std::string dummy_str;
+  int bin_type, dummy_int;
+  size_t num_true_bins;
+
+  // TODO: add I/O error handling hereg
+  // Skip the output TDirectoryFile name and the input TTree name
+  bin_file >> dummy_str >> dummy_str >> num_true_bins;
+
+  const std::string delimiter( "&&" );
+  for ( size_t tb = 0u; tb < num_true_bins; ++tb ) {
+
+    // Skip the true bin type and block index
+    bin_file >> bin_type >> dummy_int;
+
+    // Use two calls to std::getline using a double quote delimiter
+    // in order to get the contents of the next double-quoted string
+    std::string bin_def;
+    std::getline( bin_file, bin_def, '\"' );
+    std::getline( bin_file, bin_def, '\"' );
+
+    // Skip bins of type == 1 (background true bins)
+    if ( bin_type == 1 ) continue;
+
+    // Skip the text before the first "&&" (here we assume that it is a simple
+    // bool for the signal definition)
+    std::string all_cuts = bin_def.substr(
+      bin_def.find(delimiter) + delimiter.length() );
+
+    // Create an empty vector of MyCut objects to start defining the new bin
+    fBinDefinitions.emplace_back();
+
+    // Loop over the rest of the cuts and check each one
+    size_t pos = 0u;
+    do {
+      // Advance to the next individual cut, removing it from the temporary
+      // copy of the list of cuts
+      pos = all_cuts.find( delimiter );
+      std::string cut = all_cuts.substr( 0, pos );
+      all_cuts.erase( 0, pos + delimiter.length() );
+
+      // Parse the cut
+      std::string var_name;
+      std::string comp_op;
+      double cut_val;
+      std::stringstream cut_ss( cut );
+
+      cut_ss >> var_name >> comp_op >> cut_val;
+
+      // Interpret the variable for this cut and set up a function
+      // object that will calculate it from the input FitEvent
+      std::function< double(FitEvent*) > getter;
+
+      if ( var_name == "mc_p3_lead_p.Mag()" ) {
+        getter = [=]( FitEvent* ev ) -> double {
+          return ev->GetHMFSParticle( PROTON )->fP.Vect().Mag() / 1e3; // GeV
+        };
+      }
+      else if ( var_name == "mc_p3_lead_p.CosTheta()" ) {
+        getter = [=]( FitEvent* ev ) -> double {
+          return ev->GetHMFSParticle( PROTON )->fP.Vect().CosTheta();
+        };
+      }
+      else if ( var_name == "mc_p3_mu.Mag()" ) {
+        getter = [=]( FitEvent* ev ) -> double {
+          return ev->GetHMFSParticle( MU_MINUS )->fP.Vect().Mag() / 1e3; // GeV
+        };
+      }
+      else if ( var_name == "mc_p3_mu.CosTheta()" ) {
+        getter = [=]( FitEvent* ev ) -> double {
+          return ev->GetHMFSParticle( MU_MINUS )->fP.Vect().CosTheta();
+        };
+      }
+      else NUIS_ABORT( "Unrecognized cut variable " + var_name );
+     
+      // Test the cut based on the appropriate comparison operator
+      std::function< bool(double) > tester;
+      if ( comp_op == ">=" ) {
+        tester = [cut_val]( double x ) -> bool {
+          return x >= cut_val;
+        };
+      }
+      else if ( comp_op == "<" ) {
+        tester = [cut_val]( double x ) -> bool {
+          return x < cut_val;
+        };
+      }
+      else NUIS_ABORT( "Unrecognized comparison operator " + comp_op );
+
+      // Add the finished cut to the definition for the current bin
+      fBinDefinitions.back().emplace_back( getter, tester );
+
+    } while ( pos != std::string::npos );
+
+  } // loop over true bins
+
+}
+/////////////////////////////////////////////////////////////////////
+//
+/////////////////////////////////////////////////////////////////////
+void MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::FillHistograms() {
+
+  if ( !Signal ) return;
+
+  // The histograms have bin number as the x-axis variable, so just fill them
+  // using the bin index of each bin that passed all cuts
+  for ( const auto& bin : fPassingBins ) {
+    NUIS_LOG(DEB, "Fill MCHist: " << bin << ", " << Weight);
+
+    fMCHist->Fill( bin, Weight );
+    fMCStat->Fill( bin, 1.0 );
+    if ( fMCHist_Modes ) fMCHist_Modes->Fill( Mode, bin, Weight );
+
+    fMCFine->Fill( bin, Weight );
+    if ( fMCFine_Modes ) fMCFine_Modes->Fill( Mode, bin, Weight );
+  }
+
+}
+
+
+void MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::ConvertEventRates() {
+  // Do the standard conversion
+  Measurement1D::ConvertEventRates();
+
+  // TODO: Restore use of A_C once you make sure you aren't double-counting it
+  //// Build a column vector using the predicted cross sections
+  //int num_bins = fMCHist->GetNbinsX();
+  //TMatrixD pred( num_bins, 1 );
+  //for ( int b = 0; b < num_bins; ++b ) {
+  //  pred( b, 0 ) = fMCHist->GetBinContent( b + 1 );
+  //}
+
+  //// Apply the additional smearing matrix to create a new prediction
+  //TMatrixD new_pred( *fAddSmear, TMatrixD::kMult, pred );
+
+  //// Update the MC prediction histogram with the smeared version
+  //for ( int b = 0; b < num_bins; ++b ) {
+  //  double xsec = new_pred( b, 0 );
+  //  fMCHist->SetBinContent( b + 1, xsec );
+  //}
+}
+
+double MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D::GetLikelihood() {
+
+  if ( fNoData || !fDataHist ) return 0.;
+
+  // Apply Masking to MC if Required.
+  if ( fIsMask and fMaskHist ) {
+    NUIS_ERR(FTL, "Bin masks not yet supported by"
+      " the MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D sample" );
+    //PlotUtils::MaskBins(fMCHist, fMaskHist);
+  }
+
+  // Likelihood Calculation
+  double stat = 0.;
+  if ( fIsChi2 ) {
+    if ( fIsNS ) {
+      NUIS_ERR(FTL, "Norm-shape covariance not yet supported by"
+        " the MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D sample" );
+    }
+    stat = StatUtils::GetChi2FromCov( fDataHist, fMCHist, covar, NULL, 1.0,
+      1.0, fIsWriting ? fResidualHist : NULL );
+  }
+
+  fLikelihood = stat;
+
+  return stat;
+}
diff --git a/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h b/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h
new file mode 100644
index 000000000..7a8a22717
--- /dev/null
+++ b/src/MicroBooNE/MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D.h
@@ -0,0 +1,90 @@
+// Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret
+
+/*******************************************************************************
+*    This file is part of NUISANCE.
+*
+*    NUISANCE is free software: you can redistribute it and/or modify
+*    it under the terms of the GNU General Public License as published by
+*    the Free Software Foundation, either version 3 of the License, or
+*    (at your option) any later version.
+*
+*    NUISANCE is distributed in the hope that it will be useful,
+*    but WITHOUT ANY WARRANTY; without even the implied warranty of
+*    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+*    GNU General Public License for more details.
+*
+*    You should have received a copy of the GNU General Public License
+*    along with NUISANCE.  If not, see <http://www.gnu.org/licenses/>.
+*******************************************************************************/
+#ifndef MICROBOONE_BNB_NUMUCC0Pi_2025_XSEC_NU_2D_H_SEEN
+#define MICROBOONE_BNB_NUMUCC0Pi_2025_XSEC_NU_2D_H_SEEN
+
+
+
+#include <functional>
+#include <memory>
+#include <vector>
+
+#include "Measurement1D.h"
+
+class TH2D;
+
+struct MyCut {
+  MyCut( std::function< double(FitEvent*) > getter,
+    std::function< bool(double) > tester ) : getter_( getter ),
+    tester_( tester ) {}
+
+  std::function< double(FitEvent*) > getter_;
+  std::function< bool(double) > tester_;
+
+  inline bool evaluate( FitEvent* event ) const {
+    return tester_( getter_(event) );
+  }
+};
+
+class MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D : public Measurement1D {
+
+public:
+
+  /// Basic Constructor.
+  MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D( nuiskey samplekey );
+
+  /// Virtual Destructor
+  inline ~MicroBooNE_BNB_NumuCC0Pi_2025_XSec_nu_2D() {
+    fResidualHist = NULL;
+    fChi2LessBinHist = NULL;
+  };
+
+  /// Apply signal definition
+  virtual bool isSignal( FitEvent* nvect );
+
+  /// Fill kinematic distributions
+  void FillEventVariables( FitEvent* customEvent );
+
+  virtual void FillHistograms() override;
+
+  ///  Adjust the standard event rate conversion to include
+  ///  application of the additional smearing matrix A_C
+  void ConvertEventRates() override;
+
+  /// Work around some hard-coded assumptions in Measurement1D::GetLikelihood()
+  double GetLikelihood() override;
+
+private:
+
+  void LoadBinDefinitions();
+
+  // Each bin is defined as a series of cuts that are applied to a FitEvent to
+  // determine whether it belongs
+  std::vector< std::vector<MyCut> > fBinDefinitions;
+
+  // Temporary storage for the index of each bin that passed all cuts for any
+  // particular event
+  std::vector< size_t > fPassingBins;
+
+  // Additional smearing matrix used to transform the input MC predictions
+  std::unique_ptr< TMatrixD > fAddSmear;
+
+};
+
+#endif
diff --git a/src/MicroBooNE/MicroBooNE_SignalDef.cxx b/src/MicroBooNE/MicroBooNE_SignalDef.cxx
index b21f0c438..e1b63a92e 100644
--- a/src/MicroBooNE/MicroBooNE_SignalDef.cxx
+++ b/src/MicroBooNE/MicroBooNE_SignalDef.cxx
@@ -200,7 +200,40 @@ bool isCC1Mu1p(FitEvent* event, double EnuMin, double EnuMax) {
   return true;
 }
 
-    //----------------------------------------//
+ //----------------------------------------//
+
+bool isCC1Mu0pi_2025(FitEvent* event, double EnuMin, double EnuMax) {
+  // Check CC inclusive
+  if (!SignalDef::isCCINC(event, 14, EnuMin, EnuMax)) return false;
+
+  // Veto events which don't have exactly 1 FS muon
+  if (event->NumFSMuon() != 1) return false;
+
+  // Muon momentum range
+  auto mu_mom = event->GetHMFSParticle(13)->p();
+  if ((mu_mom < 100) || (mu_mom > 2000)) return false;
+
+  // Reject events with neutral pions of any momenta
+  if (event->NumFSParticle(111) != 0) return false;
+
+  // Reject events with positively charged pions above 70 MeV/c
+  if (event->NumFSParticle(211) != 0) {
+    double ppiplus = event->GetHMFSParticle(211)->p();
+    if (ppiplus > 70) { return false; }
+  }
+
+  // Reject events with negatively charged pions above 70 MeV/c
+  if (event->NumFSParticle(-211) != 0) {
+    double ppiminus = event->GetHMFSParticle(-211)->p();
+    if (ppiminus > 70) { return false; }
+  }
+
+  return true;
+}
+
+
+
+
 
 
   }  // namespace MicroBooNE
diff --git a/src/MicroBooNE/MicroBooNE_SignalDef.h b/src/MicroBooNE/MicroBooNE_SignalDef.h
index 7513cb708..95e727eda 100644
--- a/src/MicroBooNE/MicroBooNE_SignalDef.h
+++ b/src/MicroBooNE/MicroBooNE_SignalDef.h
@@ -48,6 +48,9 @@ bool isCC1Mu2p(FitEvent* event, double EnuMin, double EnuMax);
 bool isCC1Mu1p(FitEvent* event, double EnuMin, double EnuMax);
 std::vector<FitParticle*> GetCC1Mu1pProtonsInPS(FitEvent* event);
 
+
+bool isCC1Mu0pi_2025(FitEvent* event, double EnuMin, double EnuMax);
+
   }  // namespace MicroBooNE
 }  // namespace SignalDef