From 4b186f10f10c25ad4199e482083176948f2fd397 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 14:49:22 -0600
Subject: [PATCH 01/26] plumber2 scripts: Find csv file from anywhere.

---
 .../ctsm/site_and_regional/plumber2_shared.py | 25 +++++++++++++++++++
 .../plumber2_surf_wrapper.py                  | 10 ++++++--
 .../site_and_regional/plumber2_usermods.py    | 10 ++++++--
 3 files changed, 41 insertions(+), 4 deletions(-)
 create mode 100644 python/ctsm/site_and_regional/plumber2_shared.py

diff --git a/python/ctsm/site_and_regional/plumber2_shared.py b/python/ctsm/site_and_regional/plumber2_shared.py
new file mode 100644
index 0000000000..70ca2dd800
--- /dev/null
+++ b/python/ctsm/site_and_regional/plumber2_shared.py
@@ -0,0 +1,25 @@
+"""
+Things shared between plumber2 scripts
+"""
+
+import os
+import pandas as pd
+
+PLUMBER2_SITES_CSV = os.path.realpath(
+    os.path.join(
+        os.path.dirname(__file__),
+        os.pardir,
+        os.pardir,
+        os.pardir,
+        "tools",
+        "site_and_regional",
+        "PLUMBER2_sites.csv",
+    )
+)
+
+
+def read_plumber2_sites_csv():
+    """
+    Read PLUMBER2_sites.csv using pandas
+    """
+    return pd.read_csv(PLUMBER2_SITES_CSV, skiprows=4)
diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 022914d17e..b10068def8 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -23,10 +23,16 @@
 import argparse
 import logging
 import os
+import sys
 import subprocess
 import tqdm
 
-import pandas as pd
+# Get the ctsm tools
+_CTSM_PYTHON = os.path.abspath(os.path.join(os.path.dirname(__file__), "..", "..", "python"))
+sys.path.insert(1, _CTSM_PYTHON)
+
+# pylint:disable=wrong-import-position
+from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
 
 
 def get_parser():
@@ -90,7 +96,7 @@ def main():
     if args.verbose:
         logging.basicConfig(level=logging.DEBUG)
 
-    plumber2_sites = pd.read_csv("PLUMBER2_sites.csv", skiprows=4)
+    plumber2_sites = read_plumber2_sites_csv()
 
     for _, row in tqdm.tqdm(plumber2_sites.iterrows()):
         lat = row["Lat"]
diff --git a/python/ctsm/site_and_regional/plumber2_usermods.py b/python/ctsm/site_and_regional/plumber2_usermods.py
index 7b7f294a24..6fcd4a6224 100644
--- a/python/ctsm/site_and_regional/plumber2_usermods.py
+++ b/python/ctsm/site_and_regional/plumber2_usermods.py
@@ -11,9 +11,15 @@
 from __future__ import print_function
 
 import os
+import sys
 import tqdm
 
-import pandas as pd
+# Get the ctsm tools
+_CTSM_PYTHON = os.path.abspath(os.path.join(os.path.dirname(__file__), "..", "..", "python"))
+sys.path.insert(1, _CTSM_PYTHON)
+
+# pylint:disable=wrong-import-position
+from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
 
 
 # Big ugly function to create usermod_dirs for each site
@@ -155,7 +161,7 @@ def main():
     """
 
     # For now we can just run the 'main' program as a loop
-    plumber2_sites = pd.read_csv("PLUMBER2_sites.csv", skiprows=4)
+    plumber2_sites = read_plumber2_sites_csv()
 
     for _, row in tqdm.tqdm(plumber2_sites.iterrows()):
         lat = row["Lat"]

From 0bc1cce4e4145f1c5fdfdf205d586ea452ba8fb2 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 14:57:24 -0600
Subject: [PATCH 02/26] plumber2_surf_wrapper: Call subset_data directly.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 10 ++++------
 1 file changed, 4 insertions(+), 6 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index b10068def8..d2520a265f 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -24,7 +24,6 @@
 import logging
 import os
 import sys
-import subprocess
 import tqdm
 
 # Get the ctsm tools
@@ -33,6 +32,7 @@
 
 # pylint:disable=wrong-import-position
 from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
+from ctsm import subset_data
 
 
 def get_parser():
@@ -77,12 +77,10 @@ def execute(command):
     print("\n", " >>  ", *command, "\n")
 
     try:
-        subprocess.check_call(command, stdout=open(os.devnull, "w"), stderr=subprocess.STDOUT)
+        sys.argv = command
+        subset_data.main()
 
-    except subprocess.CalledProcessError as err:
-        # raise RuntimeError("command '{}' return with error
-        # (code {}): {}".format(e.cmd, e.returncode, e.output))
-        # print (e.ouput)
+    except Exception as err:  # pylint: disable=broad-exception-caught
         print(err)
 
 

From 55b97e6838cd96f405b1c19029208b8e7932450e Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 14:57:44 -0600
Subject: [PATCH 03/26] plumber2_surf_wrapper: Specify --lon-type 180.

---
 python/ctsm/site_and_regional/plumber2_shared.py       | 2 +-
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 4 ++++
 tools/site_and_regional/PLUMBER2_sites.csv             | 1 +
 3 files changed, 6 insertions(+), 1 deletion(-)

diff --git a/python/ctsm/site_and_regional/plumber2_shared.py b/python/ctsm/site_and_regional/plumber2_shared.py
index 70ca2dd800..d2232d0b2d 100644
--- a/python/ctsm/site_and_regional/plumber2_shared.py
+++ b/python/ctsm/site_and_regional/plumber2_shared.py
@@ -22,4 +22,4 @@ def read_plumber2_sites_csv():
     """
     Read PLUMBER2_sites.csv using pandas
     """
-    return pd.read_csv(PLUMBER2_SITES_CSV, skiprows=4)
+    return pd.read_csv(PLUMBER2_SITES_CSV, skiprows=5)
diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index d2520a265f..6b04b70a66 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -152,6 +152,8 @@ def main():
                 "--cap-saturation",
                 "--verbose",
                 "--overwrite",
+                "--lon-type",
+                "180",
             ]
         else:
             # use surface dataset with 78 pfts, and overwrite to 100% 1 dominant PFT
@@ -179,6 +181,8 @@ def main():
                 "--cap-saturation",
                 "--verbose",
                 "--overwrite",
+                "--lon-type",
+                "180",
             ]
         execute(subset_command)
 
diff --git a/tools/site_and_regional/PLUMBER2_sites.csv b/tools/site_and_regional/PLUMBER2_sites.csv
index f252fa1d61..8ace57ad7c 100644
--- a/tools/site_and_regional/PLUMBER2_sites.csv
+++ b/tools/site_and_regional/PLUMBER2_sites.csv
@@ -2,6 +2,7 @@
 #start_year and end_year will be used to define DATM_YR_ALIGH, DATM_YR_START and DATM_YR_END, and STOP_N in units of nyears.
 #RUN_STARTDATE and START_TOD are specified because we are starting at GMT corresponding to local midnight.
 #ATM_NCPL is specified so that the time step of the model matches the time interval specified by the atm forcing data.
+#longitudes must be in the range [-180,180]
 ,Site,Lat,Lon,pft1,pft1-%,pft1-cth,pft1-cbh,pft2,pft2-%,pft2-cth,pft2-cbh,start_year,end_year,RUN_STARTDATE,START_TOD,ATM_NCPL
 1,AR-SLu,-33.464802,-66.459808,5,50.00, 4.50, 0.13,7,50.00, 4.50, 2.59,2010,2010,2010-01-01,10800,48 
 2,AT-Neu,47.116669,11.317500,13,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2002,2012,2001-12-31,82800,48 

From 44d461d6b13751feb7e724ac0f8c7b177d362b41 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 15:03:05 -0600
Subject: [PATCH 04/26] plumber2_surf_wrapper: Respect --verbose.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 7 ++++---
 1 file changed, 4 insertions(+), 3 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 6b04b70a66..c494654334 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -51,7 +51,6 @@ def get_parser():
         help="Verbose mode will print more information. ",
         action="store_true",
         dest="verbose",
-        default=False,
     )
 
     parser.add_argument(
@@ -150,7 +149,6 @@ def main():
                 "--create-surface",
                 "--uniform-snowpack",
                 "--cap-saturation",
-                "--verbose",
                 "--overwrite",
                 "--lon-type",
                 "180",
@@ -179,11 +177,14 @@ def main():
                 "--create-surface",
                 "--uniform-snowpack",
                 "--cap-saturation",
-                "--verbose",
                 "--overwrite",
                 "--lon-type",
                 "180",
             ]
+
+        if args.verbose:
+            subset_command += ["--verbose"]
+
         execute(subset_command)
 
 

From 6d96107b1a768abdb435d87739515676d6c24962 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 15:04:34 -0600
Subject: [PATCH 05/26] plumber2_surf_wrapper: Stop on errors.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 8 ++------
 1 file changed, 2 insertions(+), 6 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index c494654334..9352eedb60 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -75,12 +75,8 @@ def execute(command):
     """
     print("\n", " >>  ", *command, "\n")
 
-    try:
-        sys.argv = command
-        subset_data.main()
-
-    except Exception as err:  # pylint: disable=broad-exception-caught
-        print(err)
+    sys.argv = command
+    subset_data.main()
 
 
 def main():

From 8781571ba1e33902edf424598155f7cbe0044e38 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 15:10:25 -0600
Subject: [PATCH 06/26] plumber2_surf_wrapper: Add optional
 --plumber2-sites-csv argument.

---
 python/ctsm/site_and_regional/plumber2_shared.py       |  4 ++--
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 10 ++++++++--
 2 files changed, 10 insertions(+), 4 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_shared.py b/python/ctsm/site_and_regional/plumber2_shared.py
index d2232d0b2d..491b35e7d2 100644
--- a/python/ctsm/site_and_regional/plumber2_shared.py
+++ b/python/ctsm/site_and_regional/plumber2_shared.py
@@ -18,8 +18,8 @@
 )
 
 
-def read_plumber2_sites_csv():
+def read_plumber2_sites_csv(file=PLUMBER2_SITES_CSV):
     """
     Read PLUMBER2_sites.csv using pandas
     """
-    return pd.read_csv(PLUMBER2_SITES_CSV, skiprows=5)
+    return pd.read_csv(file, skiprows=5)
diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 9352eedb60..c0c287b356 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -31,7 +31,7 @@
 sys.path.insert(1, _CTSM_PYTHON)
 
 # pylint:disable=wrong-import-position
-from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
+from ctsm.site_and_regional.plumber2_shared import PLUMBER2_SITES_CSV, read_plumber2_sites_csv
 from ctsm import subset_data
 
 
@@ -61,6 +61,12 @@ def get_parser():
         default=True,
     )
 
+    parser.add_argument(
+        "--plumber2-sites-csv",
+        help=f"Comma-separated value (CSV) file with Plumber2 sites. Default: {PLUMBER2_SITES_CSV}",
+        default=PLUMBER2_SITES_CSV,
+    )
+
     return parser
 
 
@@ -89,7 +95,7 @@ def main():
     if args.verbose:
         logging.basicConfig(level=logging.DEBUG)
 
-    plumber2_sites = read_plumber2_sites_csv()
+    plumber2_sites = read_plumber2_sites_csv(args.plumber2_sites_csv)
 
     for _, row in tqdm.tqdm(plumber2_sites.iterrows()):
         lat = row["Lat"]

From 6943fe4055a9f021f8059f0f6ca7d2fcb6eb470d Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 16:04:17 -0600
Subject: [PATCH 07/26] plumber2_surf_wrapper: Fix handling of one-PFT sites.

Resolves ESCOMP/CTSM#3262.
---
 .../plumber2_surf_wrapper.py                  | 134 +++++++++---------
 .../site_and_regional/single_point_case.py    |   3 +-
 tools/site_and_regional/PLUMBER2_sites.csv    |  10 +-
 3 files changed, 77 insertions(+), 70 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index c0c287b356..c97f6772ab 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -85,6 +85,13 @@ def execute(command):
     subset_data.main()
 
 
+def is_valid_pft(pft_num):
+    """
+    Given a number, check whether it represents a valid PFT
+    """
+    return pft_num >= 1
+
+
 def main():
     """
     Read plumber2_sites from csv, iterate through sites, and add dominant PFT
@@ -101,89 +108,88 @@ def main():
         lat = row["Lat"]
         lon = row["Lon"]
         site = row["Site"]
+
+        clmsite = "1x1_PLUMBER2_" + site
+        print("Now processing site :", site)
+
+        # Set up part of subset_data command that is shared among all options
+        subset_command = [
+            "./subset_data",
+            "point",
+            "--lat",
+            str(lat),
+            "--lon",
+            str(lon),
+            "--site",
+            clmsite,
+            "--create-surface",
+            "--uniform-snowpack",
+            "--cap-saturation",
+            "--overwrite",
+            "--lon-type",
+            "180",
+        ]
+
+        # Read info for first PFT
         pft1 = row["pft1"]
+        if not is_valid_pft(pft1):
+            raise RuntimeError(f"pft1 must be a valid PFT; got {pft1}")
         pctpft1 = row["pft1-%"]
         cth1 = row["pft1-cth"]
         cbh1 = row["pft1-cbh"]
-        pft2 = row["pft2"]
-        pctpft2 = row["pft2-%"]
-        cth2 = row["pft2-cth"]
-        cbh2 = row["pft2-cbh"]
-        # overwrite missing values from .csv file
-        if pft1 == -999:
-            pft1 = 0
-            pctpft1 = 0
-            cth1 = 0
-            cbh1 = 0
-        if pft2 == -999:
-            pft2 = 0
-            pctpft2 = 0
-            cth2 = 0
-            cbh2 = 0
-        clmsite = "1x1_PLUMBER2_" + site
-        print("Now processing site :", site)
 
-        if args.pft_16:
-            # use surface dataset with 16 pfts, but overwrite to 100% 1 dominant PFT
-            # don't set crop flag
-            # set dominant pft
-            subset_command = [
-                "./subset_data",
-                "point",
-                "--lat",
-                str(lat),
-                "--lon",
-                str(lon),
-                "--site",
-                clmsite,
+        # Read info for second PFT, if a valid one is given in the .csv file
+        pft2 = row["pft2"]
+        if is_valid_pft(pft2):
+            pctpft2 = row["pft2-%"]
+            cth2 = row["pft2-cth"]
+            cbh2 = row["pft2-cbh"]
+
+        # Set dominant PFT(s)
+        if is_valid_pft(pft2):
+            subset_command += [
                 "--dompft",
                 str(pft1),
                 str(pft2),
                 "--pctpft",
                 str(pctpft1),
                 str(pctpft2),
-                "--cth",
-                str(cth1),
-                str(cth2),
-                "--cbh",
-                str(cbh1),
-                str(cbh2),
-                "--create-surface",
-                "--uniform-snowpack",
-                "--cap-saturation",
-                "--overwrite",
-                "--lon-type",
-                "180",
             ]
         else:
-            # use surface dataset with 78 pfts, and overwrite to 100% 1 dominant PFT
-            # NOTE: FATES will currently not run with a 78-PFT surface dataset
-            # set crop flag
-            # set dominant pft
-            subset_command = [
-                "./subset_data",
-                "point",
-                "--lat",
-                str(lat),
-                "--lon",
-                str(lon),
-                "--site",
-                clmsite,
-                "--crop",
+            subset_command += [
                 "--dompft",
                 str(pft1),
-                str(pft2),
                 "--pctpft",
                 str(pctpft1),
-                str(pctpft2),
-                "--create-surface",
-                "--uniform-snowpack",
-                "--cap-saturation",
-                "--overwrite",
-                "--lon-type",
-                "180",
             ]
 
+        if args.pft_16:
+            # use surface dataset with 16 pfts, but overwrite to 100% 1 dominant PFT
+            # don't set crop flag
+            # set canopy top and bottom heights
+            if is_valid_pft(pft2):
+                subset_command += [
+                    "--cth",
+                    str(cth1),
+                    str(cth2),
+                    "--cbh",
+                    str(cbh1),
+                    str(cbh2),
+                ]
+            else:
+                subset_command += [
+                    "--cth",
+                    str(cth1),
+                    "--cbh",
+                    str(cbh1),
+                ]
+        else:
+            # use surface dataset with 78 pfts, and overwrite to 100% 1 dominant PFT
+            # NOTE: FATES will currently not run with a 78-PFT surface dataset
+            # set crop flag
+            subset_command += ["--crop"]
+            # don't set canopy top and bottom heights
+
         if args.verbose:
             subset_command += ["--verbose"]
 
diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index ed8b4b5562..55813a423a 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -237,7 +237,8 @@ def check_dom_pft(self):
             if min_dom_pft < NAT_PFT <= max_dom_pft:
                 err_msg = (
                     "You are subsetting using mixed land units that have both "
-                    "natural pfts and crop cfts. Check your surface dataset. "
+                    "natural pfts and crop cfts. Check your surface dataset.\n"
+                    f"{min_dom_pft} < {NAT_PFT} <= {max_dom_pft}\n"
                 )
                 raise argparse.ArgumentTypeError(err_msg)
 
diff --git a/tools/site_and_regional/PLUMBER2_sites.csv b/tools/site_and_regional/PLUMBER2_sites.csv
index 8ace57ad7c..1097568051 100644
--- a/tools/site_and_regional/PLUMBER2_sites.csv
+++ b/tools/site_and_regional/PLUMBER2_sites.csv
@@ -74,7 +74,7 @@
 68,DK-Sor,55.485870,11.644640,7,100.00,25.00,14.37,-999,-999.00,-999.00,-999.00,1997,2014,1996-12-31,82800,48 
 69,DK-ZaH,74.473282,-20.550293,12,100.00, 0.47, 0.01,-999,-999.00,-999.00,-999.00,2000,2013,2000-01-01,0,48 
 70,ES-ES1,39.345970,-0.318817,1,100.00, 7.50, 3.75,-999,-999.00,-999.00,-999.00,1999,2006,1998-12-31,82800,48 
-71,ES-ES2,39.275558,-0.315277,-999,-999.00,-999.00,-999.00,16,100.00, 0.50, 0.01,2005,2006,2004-12-31,82800,48 
+71,ES-ES2,39.275558,-0.315277,16,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2005,2006,2004-12-31,82800,48 
 72,ES-LgS,37.097935,-2.965820,10,30.00, 0.20, 0.04,13,70.00, 0.50, 0.01,2007,2007,2006-12-31,82800,48 
 73,ES-LMa,39.941502,-5.773346,7,30.00, 8.00, 4.60,14,70.00, 0.50, 0.01,2004,2006,2003-12-31,82800,48 
 74,ES-VDA,42.152180, 1.448500,7,30.00, 0.50, 0.29,13,70.00, 0.50, 0.01,2004,2004,2003-12-31,82800,48 
@@ -95,7 +95,7 @@
 89,IE-Ca1,52.858791,-6.918152,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2004,2006,2004-01-01,0,48 
 90,IE-Dri,51.986691,-8.751801,13,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2003,2005,2003-01-01,0,48 
 91,IT-Amp,41.904099,13.605160,13,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2003,2006,2002-12-31,82800,48 
-92,IT-BCi,40.523800,14.957440,-999,-999.00,-999.00,-999.00,16,100.00, 0.50, 0.01,2005,2010,2004-12-31,82800,48 
+92,IT-BCi,40.523800,14.957440,16,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2005,2010,2004-12-31,82800,48 
 93,IT-CA1,42.380409,12.026560,7,100.00, 5.50, 3.16,-999,-999.00,-999.00,-999.00,2012,2013,2011-12-31,82800,48 
 94,IT-CA2,42.377220,12.026040,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2012,2013,2011-12-31,82800,48 
 95,IT-CA3,42.380001,12.022200,7,100.00, 3.50, 2.01,-999,-999.00,-999.00,-999.00,2012,2013,2011-12-31,82800,48 
@@ -152,8 +152,8 @@
 146,US-MMS,39.323200,-86.413086,7,100.00,27.00,15.52,-999,-999.00,-999.00,-999.00,1999,2014,1999-01-01,18000,24 
 147,US-MOz,38.744110,-92.200012,7,100.00,24.00,13.80,-999,-999.00,-999.00,-999.00,2005,2006,2005-01-01,21600,48 
 148,US-Myb,38.049801,-121.765106,13,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2011,2014,2011-01-01,28800,48 
-149,US-Ne1,41.165100,-96.476593,-999,-999.00,-999.00,-999.00,16,100.00, 0.50, 0.01,2002,2012,2002-01-01,21600,24 
-150,US-Ne2,41.164902,-96.470093,-999,-999.00,-999.00,-999.00,16,100.00, 0.50, 0.01,2002,2012,2002-01-01,21600,24 
+149,US-Ne1,41.165100,-96.476593,16,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2002,2012,2002-01-01,21600,24 
+150,US-Ne2,41.164902,-96.470093,16,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2002,2012,2002-01-01,21600,24 
 151,US-Ne3,41.179699,-96.439697,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2002,2012,2002-01-01,21600,24 
 152,US-NR1,40.032902,-105.546402,1,100.00,12.00, 6.00,-999,-999.00,-999.00,-999.00,1999,2014,1999-01-01,25200,48 
 153,US-PFa,45.945900,-90.272308,1, 8.18,30.00,15.00,7,91.82,30.00,17.25,1995,2014,1995-01-01,21600,24 
@@ -166,7 +166,7 @@
 160,US-Syv,46.242001,-89.347717,1, 4.91,27.00,13.50,7,95.09,27.00,15.53,2002,2008,2002-01-01,21600,48 
 161,US-Ton,38.431599,-120.966003,7,70.00, 7.10, 4.08,14,30.00, 0.50, 0.01,2001,2014,2001-01-01,28800,48 
 162,US-Tw4,38.103001,-121.641403,13,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2014,2014,2014-01-01,28800,48 
-163,US-Twt,38.108700,-121.653107,-999,-999.00,-999.00,-999.00,16,100.00, 0.50, 0.01,2010,2014,2010-01-01,28800,48 
+163,US-Twt,38.108700,-121.653107,16,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2010,2014,2010-01-01,28800,48 
 164,US-UMB,45.559799,-84.713806,7,100.00,20.00,11.50,-999,-999.00,-999.00,-999.00,2000,2014,2000-01-01,18000,24 
 165,US-Var,38.413300,-120.950729,14,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2001,2014,2001-01-01,28800,48 
 166,US-WCr,45.805901,-90.079895,7,100.00,24.00,13.80,-999,-999.00,-999.00,-999.00,1999,2006,1999-01-01,21600,48 

From d2c711825570bf2c8fcb546b6684badbab2b8f13 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Wed, 18 Jun 2025 16:43:07 -0600
Subject: [PATCH 08/26] Add Python system test for plumber2_surf_wrapper.

---
 .../test/test_sys_plumber2_surf_wrapper.py    | 70 +++++++++++++++++++
 1 file changed, 70 insertions(+)
 create mode 100755 python/ctsm/test/test_sys_plumber2_surf_wrapper.py

diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
new file mode 100755
index 0000000000..4bd2a6763a
--- /dev/null
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -0,0 +1,70 @@
+#!/usr/bin/env python3
+
+"""System tests for plumber2_surf_wrapper"""
+
+import glob
+import os
+import unittest
+import tempfile
+import shutil
+import sys
+
+from ctsm import unit_testing
+from ctsm.site_and_regional.plumber2_surf_wrapper import main
+from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
+from ctsm.path_utils import path_to_ctsm_root
+
+# Allow test names that pylint doesn't like; otherwise hard to make them
+# readable
+# pylint: disable=invalid-name
+
+
+class TestSysPlumber2SurfWrapper(unittest.TestCase):
+    """
+    System tests for plumber2_surf_wrapper
+    """
+
+    def setUp(self):
+        """
+        Make tempdir for use by these tests.
+        """
+        self._previous_dir = os.getcwd()
+        self._tempdir = tempfile.mkdtemp()
+        os.chdir(self._tempdir)  # cd to tempdir
+
+    def tearDown(self):
+        """
+        Remove temporary directory
+        """
+        os.chdir(self._previous_dir)
+        shutil.rmtree(self._tempdir, ignore_errors=True)
+
+    def test_plumber2_surf_wrapper(self):
+        """
+        Run the entire tool
+        """
+
+        tool_path = os.path.join(
+            path_to_ctsm_root(),
+            "tools",
+            "site_and_regional",
+            "plumber2_surf_wrapper",
+        )
+        sys.argv = [tool_path]
+        main()
+
+        # How many files do we expect?
+        plumber2_csv = read_plumber2_sites_csv()
+        n_files_expected = len(plumber2_csv)
+
+        # How many files did we get?
+        file_list = os.listdir("subset_data_single_point")
+        n_files = len(file_list)
+
+        # Check
+        self.assertEqual(n_files_expected, n_files)
+
+
+if __name__ == "__main__":
+    unit_testing.setup_for_tests()
+    unittest.main()

From 132d6f8c452255bb31c987c97dde7ef235f50ef2 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 15:28:58 -0600
Subject: [PATCH 09/26] Remove an unused import.

---
 python/ctsm/test/test_sys_plumber2_surf_wrapper.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index 4bd2a6763a..a35af29ad2 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -2,7 +2,6 @@
 
 """System tests for plumber2_surf_wrapper"""
 
-import glob
 import os
 import unittest
 import tempfile

From 194ce4397797fa15d937b09925a2c1007121155e Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 15:46:04 -0600
Subject: [PATCH 10/26] Test plumber2_surf_wrapper invalid-PFT error.

---
 .../test/test_sys_plumber2_surf_wrapper.py    | 34 +++++++++++++++----
 .../PLUMBER2_sites_invalid_pft.csv            |  8 +++++
 2 files changed, 35 insertions(+), 7 deletions(-)
 create mode 100644 python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv

diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index a35af29ad2..4cb34aa89a 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -31,6 +31,19 @@ def setUp(self):
         self._tempdir = tempfile.mkdtemp()
         os.chdir(self._tempdir)  # cd to tempdir
 
+        # Path to script
+        self.tool_path = os.path.join(
+            path_to_ctsm_root(),
+            "tools",
+            "site_and_regional",
+            "plumber2_surf_wrapper",
+        )
+
+        # Path to test inputs directory
+        self.test_inputs = os.path.join(
+            os.path.dirname(__file__), "testinputs", "plumber2_surf_wrapper"
+        )
+
     def tearDown(self):
         """
         Remove temporary directory
@@ -43,13 +56,7 @@ def test_plumber2_surf_wrapper(self):
         Run the entire tool
         """
 
-        tool_path = os.path.join(
-            path_to_ctsm_root(),
-            "tools",
-            "site_and_regional",
-            "plumber2_surf_wrapper",
-        )
-        sys.argv = [tool_path]
+        sys.argv = [self.tool_path]
         main()
 
         # How many files do we expect?
@@ -63,6 +70,19 @@ def test_plumber2_surf_wrapper(self):
         # Check
         self.assertEqual(n_files_expected, n_files)
 
+    def test_plumber2_surf_wrapper_invalid_pft(self):
+        """
+        plumber2_surf_wrapper should error if invalid PFT is given
+        """
+
+        sys.argv = [
+            self.tool_path,
+            "--plumber2-sites-csv",
+            os.path.join(self.test_inputs, "PLUMBER2_sites_invalid_pft.csv"),
+        ]
+        with self.assertRaisesRegex(RuntimeError, "must be a valid PFT"):
+            main()
+
 
 if __name__ == "__main__":
     unit_testing.setup_for_tests()
diff --git a/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv
new file mode 100644
index 0000000000..2d4b7dcb57
--- /dev/null
+++ b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv
@@ -0,0 +1,8 @@
+#pftX-cth and pftX-cbh are the site=specific canopy top and bottom heights
+#start_year and end_year will be used to define DATM_YR_ALIGH, DATM_YR_START and DATM_YR_END, and STOP_N in units of nyears.
+#RUN_STARTDATE and START_TOD are specified because we are starting at GMT corresponding to local midnight.
+#ATM_NCPL is specified so that the time step of the model matches the time interval specified by the atm forcing data.
+#longitudes must be in the range [-180,180]
+,Site,Lat,Lon,pft1,pft1-%,pft1-cth,pft1-cbh,pft2,pft2-%,pft2-cth,pft2-cbh,start_year,end_year,RUN_STARTDATE,START_TOD,ATM_NCPL
+26,Invalid-Pft,51.309166, 4.520560,0,19.22,21.00,10.50,7,80.78,21.00,12.08,2004,2014,2003-12-31,82800,48 
+27,BE-Lon,50.551590, 4.746130,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2005,2014,2004-12-31,82800,48 

From b9ed339e14e70a98e2efa80d2a2c0fa5257125a3 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 16:10:54 -0600
Subject: [PATCH 11/26] Replace plumber2_surf_wrapper args test with useful
 ones. Failing.

---
 .../plumber2_surf_wrapper.py                  |  8 +--
 .../test/test_unit_plumber2_surf_wrapper.py   | 56 +++++++++++++++++--
 2 files changed, 56 insertions(+), 8 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index c97f6772ab..46cdacf3ef 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -35,9 +35,9 @@
 from ctsm import subset_data
 
 
-def get_parser():
+def get_args():
     """
-    Get parser object for this script.
+    Get arguments for this script.
     """
     parser = argparse.ArgumentParser(
         description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter
@@ -67,7 +67,7 @@ def get_parser():
         default=PLUMBER2_SITES_CSV,
     )
 
-    return parser
+    return parser.parse_args()
 
 
 def execute(command):
@@ -97,7 +97,7 @@ def main():
     Read plumber2_sites from csv, iterate through sites, and add dominant PFT
     """
 
-    args = get_parser().parse_args()
+    args = get_args()
 
     if args.verbose:
         logging.basicConfig(level=logging.DEBUG)
diff --git a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
index 66f5578caa..4e23fdd209 100755
--- a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
@@ -16,7 +16,7 @@
 
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
-from ctsm.site_and_regional.plumber2_surf_wrapper import get_parser
+from ctsm.site_and_regional.plumber2_surf_wrapper import get_args
 
 # pylint: disable=invalid-name
 
@@ -26,12 +26,60 @@ class TestPlumber2SurfWrapper(unittest.TestCase):
     Basic class for testing plumber2_surf_wrapper.py.
     """
 
-    def test_parser(self):
+    def setUp(self):
+        sys.argv = ["subset_data"]  # Could actually be anything
+
+    def test_parser_default_csv_exists(self):
+        """
+        Test that default PLUMBER2 sites CSV file exists
+        """
+
+        args = get_args()
+        self.assertTrue(os.path.exists(args.plumber2_sites_csv))
+
+    def test_parser_custom_csv(self):
+        """
+        Test that script accepts custom CSV file path
+        """
+
+        custom_path = "path/to/custom.csv"
+        sys.argv += ["--plumber2-sites-csv", custom_path]
+        args = get_args()
+        self.assertEqual(args.plumber2_sites_csv, custom_path)
+
+    def test_parser_verbose_false_default(self):
+        """
+        Test that script is not verbose by default
+        """
+
+        args = get_args()
+        self.assertFalse(args.verbose)
+
+    def test_parser_verbose_true(self):
+        """
+        Test that --verbose sets verbose to True
+        """
+
+        sys.argv += ["--verbose"]
+        args = get_args()
+        self.assertTrue(args.verbose)
+
+    def test_parser_16pft_false_default(self):
+        """
+        Test that script does not use 16pft mode by default
+        """
+
+        args = get_args()
+        self.assertFalse(args.pft_16)
+
+    def test_parser_16pft_true(self):
         """
-        Test that parser has same defaults as expected
+        Test that --16pft sets pft_16 to True
         """
 
-        self.assertEqual(get_parser().argument_default, None, "Parser not working as expected")
+        sys.argv += ["--16pft"]
+        args = get_args()
+        self.assertTrue(args.pft_16)
 
 
 if __name__ == "__main__":

From 3019c4eba0969dd7880de32e625a774cbf7fdfef Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 16:12:41 -0600
Subject: [PATCH 12/26] plumber2_surf_wrapper: Respect user not saying --16pft.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 1 -
 1 file changed, 1 deletion(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 46cdacf3ef..4a93859afe 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -58,7 +58,6 @@ def get_args():
         help="Create and/or modify 16-PFT surface datasets (e.g. for a FATES run) ",
         action="store_true",
         dest="pft_16",
-        default=True,
     )
 
     parser.add_argument(

From 4707843eb321f83c87806d1e4d427da1ef694e57 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 16:16:34 -0600
Subject: [PATCH 13/26] plumber2_surf_wrapper: Test full run with --16pft.

---
 .../test/test_sys_plumber2_surf_wrapper.py    | 23 ++++++++++++++++++-
 1 file changed, 22 insertions(+), 1 deletion(-)

diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index 4cb34aa89a..e405d487be 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -53,7 +53,8 @@ def tearDown(self):
 
     def test_plumber2_surf_wrapper(self):
         """
-        Run the entire tool
+        Run the entire tool with default settings.
+        CAN ONLY RUN ON SYSTEMS WITH INPUTDATA
         """
 
         sys.argv = [self.tool_path]
@@ -70,6 +71,26 @@ def test_plumber2_surf_wrapper(self):
         # Check
         self.assertEqual(n_files_expected, n_files)
 
+    def test_plumber2_surf_wrapper_16pft(self):
+        """
+        Run the entire tool with --16pft.
+        CAN ONLY RUN ON SYSTEMS WITH INPUTDATA
+        """
+
+        sys.argv = [self.tool_path, "--16pft"]
+        main()
+
+        # How many files do we expect?
+        plumber2_csv = read_plumber2_sites_csv()
+        n_files_expected = len(plumber2_csv)
+
+        # How many files did we get?
+        file_list = os.listdir("subset_data_single_point")
+        n_files = len(file_list)
+
+        # Check
+        self.assertEqual(n_files_expected, n_files)
+
     def test_plumber2_surf_wrapper_invalid_pft(self):
         """
         plumber2_surf_wrapper should error if invalid PFT is given

From 30f455806020d5e085bbff743ecc2d20f442f4e8 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 16:36:21 -0600
Subject: [PATCH 14/26] plumber2_surf_wrapper: Add --overwrite option.

---
 .../plumber2_surf_wrapper.py                  |  9 ++++-
 .../test/test_sys_plumber2_surf_wrapper.py    | 37 +++++++++++++++++++
 .../PLUMBER2_site_valid.csv                   |  7 ++++
 3 files changed, 52 insertions(+), 1 deletion(-)
 create mode 100644 python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_site_valid.csv

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 4a93859afe..4fd8436982 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -60,6 +60,12 @@ def get_args():
         dest="pft_16",
     )
 
+    parser.add_argument(
+        "--overwrite",
+        help="Overwrite any existing files",
+        action="store_true",
+    )
+
     parser.add_argument(
         "--plumber2-sites-csv",
         help=f"Comma-separated value (CSV) file with Plumber2 sites. Default: {PLUMBER2_SITES_CSV}",
@@ -124,7 +130,6 @@ def main():
             "--create-surface",
             "--uniform-snowpack",
             "--cap-saturation",
-            "--overwrite",
             "--lon-type",
             "180",
         ]
@@ -191,6 +196,8 @@ def main():
 
         if args.verbose:
             subset_command += ["--verbose"]
+        if args.overwrite:
+            subset_command += ["--overwrite"]
 
         execute(subset_command)
 
diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index e405d487be..de20fd11c3 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -104,6 +104,43 @@ def test_plumber2_surf_wrapper_invalid_pft(self):
         with self.assertRaisesRegex(RuntimeError, "must be a valid PFT"):
             main()
 
+    def test_plumber2_surf_wrapper_existing_no_overwrite_fails(self):
+        """
+        plumber2_surf_wrapper should fail if file exists but --overwrite isn't given
+        """
+
+        sys_argv_shared = [
+            self.tool_path,
+            "--plumber2-sites-csv",
+            os.path.join(self.test_inputs, "PLUMBER2_site_valid.csv"),
+        ]
+
+        # Run twice, expecting second to fail
+        sys.argv = sys_argv_shared
+        main()
+        sys.argv = sys_argv_shared
+        with self.assertRaisesRegex(SystemExit, "exists"):
+            main()
+
+    def test_plumber2_surf_wrapper_existing_overwrite_passes(self):
+        """
+        plumber2_surf_wrapper should pass if file exists and --overwrite is given
+        """
+
+        sys_argv_shared = [
+            self.tool_path,
+            "--plumber2-sites-csv",
+            os.path.join(self.test_inputs, "PLUMBER2_site_valid.csv"),
+        ]
+
+        # Run once to generate the files
+        sys.argv = sys_argv_shared
+        main()
+
+        # Run again with --overwrite, expecting pass
+        sys.argv = sys_argv_shared + ["--overwrite"]
+        main()
+
 
 if __name__ == "__main__":
     unit_testing.setup_for_tests()
diff --git a/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_site_valid.csv b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_site_valid.csv
new file mode 100644
index 0000000000..2c1580bc03
--- /dev/null
+++ b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_site_valid.csv
@@ -0,0 +1,7 @@
+#pftX-cth and pftX-cbh are the site=specific canopy top and bottom heights
+#start_year and end_year will be used to define DATM_YR_ALIGH, DATM_YR_START and DATM_YR_END, and STOP_N in units of nyears.
+#RUN_STARTDATE and START_TOD are specified because we are starting at GMT corresponding to local midnight.
+#ATM_NCPL is specified so that the time step of the model matches the time interval specified by the atm forcing data.
+#longitudes must be in the range [-180,180]
+,Site,Lat,Lon,pft1,pft1-%,pft1-cth,pft1-cbh,pft2,pft2-%,pft2-cth,pft2-cbh,start_year,end_year,RUN_STARTDATE,START_TOD,ATM_NCPL
+27,BE-Lon,50.551590, 4.746130,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2005,2014,2004-12-31,82800,48 

From 1c5c1e72944f0195790c1b3932869848c0bcec98 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 16:38:38 -0600
Subject: [PATCH 15/26] plumber2_surf_wrapper: Improve execute() comments.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 4fd8436982..2162a3830f 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -77,12 +77,12 @@ def get_args():
 
 def execute(command):
     """
-    Function for running a command on shell.
+    Runs subset_data with given arguments.
     Args:
-        command (str):
-            command that we want to run.
+        command (list):
+            list of args for command that we want to run.
     Raises:
-        Error with the return code from shell.
+        Whatever error subset_data gives, if any.
     """
     print("\n", " >>  ", *command, "\n")
 

From 63d61f43157710e7376ce72c4217a3c51313c087 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Thu, 19 Jun 2025 17:00:42 -0600
Subject: [PATCH 16/26] plumber2_surf_wrapper: Switch --16pft to --78pft to
 preserve previous default.

---
 .../site_and_regional/plumber2_surf_wrapper.py     | 10 +++++-----
 python/ctsm/test/test_sys_plumber2_surf_wrapper.py |  6 +++---
 .../ctsm/test/test_unit_plumber2_surf_wrapper.py   | 14 +++++++-------
 3 files changed, 15 insertions(+), 15 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 2162a3830f..367fa97a77 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -54,10 +54,10 @@ def get_args():
     )
 
     parser.add_argument(
-        "--16pft",
-        help="Create and/or modify 16-PFT surface datasets (e.g. for a FATES run) ",
+        "--78pft",
+        help="Create and/or modify 78-PFT surface datasets (e.g. for a non-FATES run) ",
         action="store_true",
-        dest="pft_16",
+        dest="pft_78",
     )
 
     parser.add_argument(
@@ -167,8 +167,8 @@ def main():
                 str(pctpft1),
             ]
 
-        if args.pft_16:
-            # use surface dataset with 16 pfts, but overwrite to 100% 1 dominant PFT
+        if not args.pft_78:
+            # use surface dataset with 78 pfts, but overwrite to 100% 1 dominant PFT
             # don't set crop flag
             # set canopy top and bottom heights
             if is_valid_pft(pft2):
diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index de20fd11c3..a7dcf12821 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -71,13 +71,13 @@ def test_plumber2_surf_wrapper(self):
         # Check
         self.assertEqual(n_files_expected, n_files)
 
-    def test_plumber2_surf_wrapper_16pft(self):
+    def test_plumber2_surf_wrapper_78pft(self):
         """
-        Run the entire tool with --16pft.
+        Run the entire tool with --78pft.
         CAN ONLY RUN ON SYSTEMS WITH INPUTDATA
         """
 
-        sys.argv = [self.tool_path, "--16pft"]
+        sys.argv = [self.tool_path, "--78pft"]
         main()
 
         # How many files do we expect?
diff --git a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
index 4e23fdd209..e3e677a8a6 100755
--- a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
@@ -64,22 +64,22 @@ def test_parser_verbose_true(self):
         args = get_args()
         self.assertTrue(args.verbose)
 
-    def test_parser_16pft_false_default(self):
+    def test_parser_78pft_false_default(self):
         """
-        Test that script does not use 16pft mode by default
+        Test that script does not use 78pft mode by default
         """
 
         args = get_args()
-        self.assertFalse(args.pft_16)
+        self.assertFalse(args.pft_78)
 
-    def test_parser_16pft_true(self):
+    def test_parser_78pft_true(self):
         """
-        Test that --16pft sets pft_16 to True
+        Test that --78pft sets pft_78 to True
         """
 
-        sys.argv += ["--16pft"]
+        sys.argv += ["--78pft"]
         args = get_args()
-        self.assertTrue(args.pft_16)
+        self.assertTrue(args.pft_78)
 
 
 if __name__ == "__main__":

From 56a36f1adec07c4478d7e9b3de7c72838cb36661 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 16:38:08 -0600
Subject: [PATCH 17/26] plumber2 scripts: Remove addition of CTSM python dir to
 path.

---
 python/ctsm/site_and_regional/plumber2_surf_wrapper.py | 5 -----
 python/ctsm/site_and_regional/plumber2_usermods.py     | 5 -----
 2 files changed, 10 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 367fa97a77..86234aae9c 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -22,14 +22,9 @@
 
 import argparse
 import logging
-import os
 import sys
 import tqdm
 
-# Get the ctsm tools
-_CTSM_PYTHON = os.path.abspath(os.path.join(os.path.dirname(__file__), "..", "..", "python"))
-sys.path.insert(1, _CTSM_PYTHON)
-
 # pylint:disable=wrong-import-position
 from ctsm.site_and_regional.plumber2_shared import PLUMBER2_SITES_CSV, read_plumber2_sites_csv
 from ctsm import subset_data
diff --git a/python/ctsm/site_and_regional/plumber2_usermods.py b/python/ctsm/site_and_regional/plumber2_usermods.py
index 6fcd4a6224..7c8f37b1b5 100644
--- a/python/ctsm/site_and_regional/plumber2_usermods.py
+++ b/python/ctsm/site_and_regional/plumber2_usermods.py
@@ -11,13 +11,8 @@
 from __future__ import print_function
 
 import os
-import sys
 import tqdm
 
-# Get the ctsm tools
-_CTSM_PYTHON = os.path.abspath(os.path.join(os.path.dirname(__file__), "..", "..", "python"))
-sys.path.insert(1, _CTSM_PYTHON)
-
 # pylint:disable=wrong-import-position
 from ctsm.site_and_regional.plumber2_shared import read_plumber2_sites_csv
 

From cf86bfc0563d0cd089078b8be9cda464978afa24 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 16:43:18 -0600
Subject: [PATCH 18/26] Simplify PLUMBER2_SITES_CSV path.

---
 python/ctsm/site_and_regional/plumber2_shared.py | 16 ++++++----------
 1 file changed, 6 insertions(+), 10 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_shared.py b/python/ctsm/site_and_regional/plumber2_shared.py
index 491b35e7d2..d4ab9d00b3 100644
--- a/python/ctsm/site_and_regional/plumber2_shared.py
+++ b/python/ctsm/site_and_regional/plumber2_shared.py
@@ -4,17 +4,13 @@
 
 import os
 import pandas as pd
+from ctsm.path_utils import path_to_ctsm_root
 
-PLUMBER2_SITES_CSV = os.path.realpath(
-    os.path.join(
-        os.path.dirname(__file__),
-        os.pardir,
-        os.pardir,
-        os.pardir,
-        "tools",
-        "site_and_regional",
-        "PLUMBER2_sites.csv",
-    )
+PLUMBER2_SITES_CSV = os.path.join(
+    path_to_ctsm_root(),
+    "tools",
+    "site_and_regional",
+    "PLUMBER2_sites.csv",
 )
 
 

From 602fe61fbfcd36b27fa1bdddd5bfb88adc23dea9 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:20:27 -0600
Subject: [PATCH 19/26] Move Python constants for PFT numbers to pft_utils.py.

---
 python/ctsm/pft_utils.py                             | 8 ++++++++
 python/ctsm/site_and_regional/single_point_case.py   | 5 +----
 python/ctsm/test/test_unit_singlept_data.py          | 3 ++-
 python/ctsm/test/test_unit_singlept_data_surfdata.py | 3 ++-
 python/ctsm/toolchain/gen_mksurfdata_namelist.py     | 3 ++-
 5 files changed, 15 insertions(+), 7 deletions(-)
 create mode 100644 python/ctsm/pft_utils.py

diff --git a/python/ctsm/pft_utils.py b/python/ctsm/pft_utils.py
new file mode 100644
index 0000000000..c7e36c9338
--- /dev/null
+++ b/python/ctsm/pft_utils.py
@@ -0,0 +1,8 @@
+"""
+Constants and functions relating to PFTs
+"""
+
+MIN_PFT = 0   # bare ground
+NAT_PFT = 15  # natural pfts
+NUM_PFT = 17  # for runs with generic crops
+MAX_PFT = 78  # for runs with explicit crops
diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index db33d875fc..d9e1a1513f 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -16,13 +16,10 @@
 from ctsm.site_and_regional.base_case import BaseCase, USRDAT_DIR, DatmFiles
 from ctsm.utils import add_tag_to_filename, ensure_iterable
 from ctsm.longitude import detect_lon_type
+from ctsm.pft_utils import NAT_PFT, NUM_PFT, MAX_PFT
 
 logger = logging.getLogger(__name__)
 
-NAT_PFT = 15  # natural pfts
-NUM_PFT = 17  # for runs with generic crops
-MAX_PFT = 78  # for runs with explicit crops
-
 
 class SinglePointCase(BaseCase):
     """
diff --git a/python/ctsm/test/test_unit_singlept_data.py b/python/ctsm/test/test_unit_singlept_data.py
index 644af82588..bc9bd1adb3 100755
--- a/python/ctsm/test/test_unit_singlept_data.py
+++ b/python/ctsm/test/test_unit_singlept_data.py
@@ -18,6 +18,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
+from ctsm.pft_utils import MAX_PFT
 
 # pylint: disable=invalid-name
 
@@ -223,7 +224,7 @@ def test_check_dom_pft_mixed_range(self):
             overwrite=self.overwrite,
         )
         single_point.dom_pft = [1, 5, 15]
-        single_point.num_pft = 78
+        single_point.num_pft = MAX_PFT
         with self.assertRaisesRegex(
             argparse.ArgumentTypeError, "You are subsetting using mixed land*"
         ):
diff --git a/python/ctsm/test/test_unit_singlept_data_surfdata.py b/python/ctsm/test/test_unit_singlept_data_surfdata.py
index 2106799a4b..71312c9db6 100755
--- a/python/ctsm/test/test_unit_singlept_data_surfdata.py
+++ b/python/ctsm/test/test_unit_singlept_data_surfdata.py
@@ -23,6 +23,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
+from ctsm.pft_utils import MAX_PFT
 
 # pylint: disable=invalid-name
 # pylint: disable=too-many-lines
@@ -667,7 +668,7 @@ class TestSinglePointCaseSurfaceCrop(unittest.TestCase):
     dom_pft = [17]
     evenly_split_cropland = False
     pct_pft = None
-    num_pft = 78
+    num_pft = MAX_PFT
     cth = 0.9
     cbh = 0.1
     include_nonveg = False
diff --git a/python/ctsm/toolchain/gen_mksurfdata_namelist.py b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
index 31fcbfe8ff..45e17bd504 100755
--- a/python/ctsm/toolchain/gen_mksurfdata_namelist.py
+++ b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
@@ -15,6 +15,7 @@
 
 from ctsm.path_utils import path_to_ctsm_root, path_to_cime
 from ctsm.ctsm_logging import setup_logging_pre_config, add_logging_args, process_logging_args
+from ctsm.pft_utils import MAX_PFT
 
 logger = logging.getLogger(__name__)
 
@@ -308,7 +309,7 @@ def main():
     if nocrop_flag:
         num_pft = "16"
     else:
-        num_pft = "78"
+        num_pft = str(MAX_PFT)
     logger.info("num_pft is %s", num_pft)
 
     # Write out if surface dataset will be created

From 8ca6e4556068f25155230bfc663aaaf64c21a8ff Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:28:14 -0600
Subject: [PATCH 20/26] Rename and replace some PFT number constants for
 clarity.

- Replace NAT_PFT=15 with MAX_NAT_PFT=14
- Replace NUM_PFT=17 with MAX_PFT_GENERICCROPS=16
- Rename MAX_PFT to MAX_PFT_MANAGEDCROPS
---
 python/ctsm/pft_utils.py                      |  6 ++--
 .../site_and_regional/single_point_case.py    | 28 +++++++++----------
 python/ctsm/test/test_unit_singlept_data.py   |  4 +--
 .../test/test_unit_singlept_data_surfdata.py  |  4 +--
 .../ctsm/toolchain/gen_mksurfdata_namelist.py |  4 +--
 5 files changed, 23 insertions(+), 23 deletions(-)

diff --git a/python/ctsm/pft_utils.py b/python/ctsm/pft_utils.py
index c7e36c9338..9588564c78 100644
--- a/python/ctsm/pft_utils.py
+++ b/python/ctsm/pft_utils.py
@@ -3,6 +3,6 @@
 """
 
 MIN_PFT = 0   # bare ground
-NAT_PFT = 15  # natural pfts
-NUM_PFT = 17  # for runs with generic crops
-MAX_PFT = 78  # for runs with explicit crops
+MAX_NAT_PFT = 14  # maximum natural pft
+MAX_PFT_GENERICCROPS = 16  # for runs with generic crops
+MAX_PFT_MANAGEDCROPS = 78  # for runs with explicit crops
diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index d9e1a1513f..c5777093ab 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -16,7 +16,7 @@
 from ctsm.site_and_regional.base_case import BaseCase, USRDAT_DIR, DatmFiles
 from ctsm.utils import add_tag_to_filename, ensure_iterable
 from ctsm.longitude import detect_lon_type
-from ctsm.pft_utils import NAT_PFT, NUM_PFT, MAX_PFT
+from ctsm.pft_utils import MAX_NAT_PFT, MAX_PFT_GENERICCROPS, MAX_PFT_MANAGEDCROPS
 
 logger = logging.getLogger(__name__)
 
@@ -187,20 +187,20 @@ def check_dom_pft(self):
           same range.
           e.g. If users specified multiple dom_pft, they should be
           either in :
-            - 0 - NAT_PFT-1 range
+            - 0 - MAX_NAT_PFT range
             or
-            - NAT_PFT - MAX_PFT range
+            - MAX_NAT_PFT+1 - MAX_PFT_MANAGEDCROPS range
             - give an error: mixed land units not possible
 
         -------------
         Raises:
             Error (ArgumentTypeError):
-                If any dom_pft is bigger than MAX_PFT.
+                If any dom_pft is bigger than MAX_PFT_MANAGEDCROPS.
             Error (ArgumentTypeError):
                 If any dom_pft is less than 1.
             Error (ArgumentTypeError):
                 If mixed land units are chosen.
-                dom_pft values are both in range of (0 - NAT_PFT-1) and (NAT_PFT - MAX_PFT).
+                dom_pft values are both in range of (0 - MAX_NAT_PFT) and (MAX_NAT_PFT+1 - MAX_PFT_MANAGEDCROPS).
 
 
         """
@@ -214,8 +214,8 @@ def check_dom_pft(self):
             min_dom_pft = min(self.dom_pft)
             max_dom_pft = max(self.dom_pft)
 
-            # -- check dom_pft values should be between 0-MAX_PFT
-            if min_dom_pft < 0 or max_dom_pft > MAX_PFT:
+            # -- check dom_pft values should be between 0-MAX_PFT_MANAGEDCROPS
+            if min_dom_pft < 0 or max_dom_pft > MAX_PFT_MANAGEDCROPS:
                 err_msg = "values for --dompft should be between 1 and 78."
                 raise argparse.ArgumentTypeError(err_msg)
 
@@ -225,17 +225,17 @@ def check_dom_pft(self):
                 raise argparse.ArgumentTypeError(err_msg)
 
             # -- check dom_pft vs MAX_pft
-            if self.num_pft - 1 < max_dom_pft < NUM_PFT:
+            if self.num_pft - 1 < max_dom_pft <= MAX_PFT_GENERICCROPS:
                 logger.info(
                     "WARNING, you trying to run with generic crops (16 PFT surface dataset)"
                 )
 
             # -- check if all dom_pft are in the same range:
-            if min_dom_pft < NAT_PFT <= max_dom_pft:
+            if min_dom_pft <= MAX_NAT_PFT < max_dom_pft:
                 err_msg = (
                     "You are subsetting using mixed land units that have both "
                     "natural pfts and crop cfts. Check your surface dataset.\n"
-                    f"{min_dom_pft} < {NAT_PFT} <= {max_dom_pft}\n"
+                    f"{min_dom_pft} <= {MAX_NAT_PFT} < {max_dom_pft}\n"
                 )
                 raise argparse.ArgumentTypeError(err_msg)
 
@@ -423,7 +423,7 @@ def modify_surfdata_atpoint(self, f_orig):
         if self.dom_pft is not None:
             max_dom_pft = max(self.dom_pft)
             # -- First initialize everything:
-            if max_dom_pft < NAT_PFT:
+            if max_dom_pft <= MAX_NAT_PFT :
                 f_mod["PCT_NAT_PFT"][:, :, :] = 0
             else:
                 f_mod["PCT_CFT"][:, :, :] = 0
@@ -442,10 +442,10 @@ def modify_surfdata_atpoint(self, f_orig):
                 if cth is not None:
                     f_mod["MONTHLY_HEIGHT_TOP"][:, :, :, dom_pft] = cth
                     f_mod["MONTHLY_HEIGHT_BOT"][:, :, :, dom_pft] = cbh
-                if dom_pft < NAT_PFT:
+                if dom_pft <= MAX_NAT_PFT:
                     f_mod["PCT_NAT_PFT"][:, :, dom_pft] = pct_pft
                 else:
-                    dom_pft = dom_pft - NAT_PFT
+                    dom_pft = dom_pft - (MAX_NAT_PFT + 1)
                     f_mod["PCT_CFT"][:, :, dom_pft] = pct_pft
 
         # -------------------------------
@@ -463,7 +463,7 @@ def modify_surfdata_atpoint(self, f_orig):
 
             if self.dom_pft is not None:
                 max_dom_pft = max(self.dom_pft)
-                if max_dom_pft < NAT_PFT:
+                if max_dom_pft <= MAX_NAT_PFT:
                     f_mod["PCT_NATVEG"][:, :] = 100
                     f_mod["PCT_CROP"][:, :] = 0
                 else:
diff --git a/python/ctsm/test/test_unit_singlept_data.py b/python/ctsm/test/test_unit_singlept_data.py
index bc9bd1adb3..dc6d655408 100755
--- a/python/ctsm/test/test_unit_singlept_data.py
+++ b/python/ctsm/test/test_unit_singlept_data.py
@@ -18,7 +18,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
-from ctsm.pft_utils import MAX_PFT
+from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
 
 # pylint: disable=invalid-name
 
@@ -224,7 +224,7 @@ def test_check_dom_pft_mixed_range(self):
             overwrite=self.overwrite,
         )
         single_point.dom_pft = [1, 5, 15]
-        single_point.num_pft = MAX_PFT
+        single_point.num_pft = MAX_PFT_MANAGEDCROPS
         with self.assertRaisesRegex(
             argparse.ArgumentTypeError, "You are subsetting using mixed land*"
         ):
diff --git a/python/ctsm/test/test_unit_singlept_data_surfdata.py b/python/ctsm/test/test_unit_singlept_data_surfdata.py
index 71312c9db6..fb6cc15720 100755
--- a/python/ctsm/test/test_unit_singlept_data_surfdata.py
+++ b/python/ctsm/test/test_unit_singlept_data_surfdata.py
@@ -23,7 +23,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
-from ctsm.pft_utils import MAX_PFT
+from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
 
 # pylint: disable=invalid-name
 # pylint: disable=too-many-lines
@@ -668,7 +668,7 @@ class TestSinglePointCaseSurfaceCrop(unittest.TestCase):
     dom_pft = [17]
     evenly_split_cropland = False
     pct_pft = None
-    num_pft = MAX_PFT
+    num_pft = MAX_PFT_MANAGEDCROPS
     cth = 0.9
     cbh = 0.1
     include_nonveg = False
diff --git a/python/ctsm/toolchain/gen_mksurfdata_namelist.py b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
index 45e17bd504..09bbb9c268 100755
--- a/python/ctsm/toolchain/gen_mksurfdata_namelist.py
+++ b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
@@ -15,7 +15,7 @@
 
 from ctsm.path_utils import path_to_ctsm_root, path_to_cime
 from ctsm.ctsm_logging import setup_logging_pre_config, add_logging_args, process_logging_args
-from ctsm.pft_utils import MAX_PFT
+from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
 
 logger = logging.getLogger(__name__)
 
@@ -309,7 +309,7 @@ def main():
     if nocrop_flag:
         num_pft = "16"
     else:
-        num_pft = str(MAX_PFT)
+        num_pft = str(MAX_PFT_MANAGEDCROPS)
     logger.info("num_pft is %s", num_pft)
 
     # Write out if surface dataset will be created

From f6a9afcd8b5624b14d61c93ff2b7688dcad1f2a4 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:40:44 -0600
Subject: [PATCH 21/26] Use MAX_PFT_ variables in more places.

---
 python/ctsm/site_and_regional/single_point_case.py   |  7 ++++---
 python/ctsm/subset_data.py                           |  5 +++--
 python/ctsm/test/test_unit_singlept_data.py          | 10 +++++-----
 python/ctsm/test/test_unit_singlept_data_surfdata.py |  4 ++--
 python/ctsm/toolchain/gen_mksurfdata_namelist.py     |  4 ++--
 5 files changed, 16 insertions(+), 14 deletions(-)

diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index c5777093ab..8c30b4c9a3 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -216,18 +216,19 @@ def check_dom_pft(self):
 
             # -- check dom_pft values should be between 0-MAX_PFT_MANAGEDCROPS
             if min_dom_pft < 0 or max_dom_pft > MAX_PFT_MANAGEDCROPS:
-                err_msg = "values for --dompft should be between 1 and 78."
+                err_msg = f"values for --dompft should be between 1 and {MAX_PFT_MANAGEDCROPS}."
                 raise argparse.ArgumentTypeError(err_msg)
 
             # -- check dom_pft vs num_pft
             if max_dom_pft > self.num_pft:
-                err_msg = "Please use --crop flag when --dompft is above 16."
+                err_msg = f"Please use --crop flag when --dompft is above {MAX_PFT_GENERICCROPS}."
                 raise argparse.ArgumentTypeError(err_msg)
 
             # -- check dom_pft vs MAX_pft
             if self.num_pft - 1 < max_dom_pft <= MAX_PFT_GENERICCROPS:
                 logger.info(
-                    "WARNING, you trying to run with generic crops (16 PFT surface dataset)"
+                    "WARNING, you are trying to run with generic crops (%s PFT surface dataset)",
+                    MAX_PFT_GENERICCROPS,
                 )
 
             # -- check if all dom_pft are in the same range:
diff --git a/python/ctsm/subset_data.py b/python/ctsm/subset_data.py
index 820391ff8b..de4e51db9b 100644
--- a/python/ctsm/subset_data.py
+++ b/python/ctsm/subset_data.py
@@ -70,6 +70,7 @@
 from ctsm.utils import abort
 from ctsm.config_utils import check_lon1_lt_lon2
 from ctsm.longitude import Longitude, detect_lon_type
+from ctsm.pft_utils import MAX_PFT_GENERICCROPS, MAX_PFT_MANAGEDCROPS
 
 # -- import ctsm logging flags
 from ctsm.ctsm_logging import (
@@ -597,9 +598,9 @@ def determine_num_pft(crop):
         num_pft (int) : number of pfts for surface dataset
     """
     if crop:
-        num_pft = "78"
+        num_pft = str(MAX_PFT_MANAGEDCROPS)
     else:
-        num_pft = "16"
+        num_pft = str(MAX_PFT_GENERICCROPS)
     logger.debug("crop_flag = %s => num_pft = %s", str(crop), num_pft)
     return num_pft
 
diff --git a/python/ctsm/test/test_unit_singlept_data.py b/python/ctsm/test/test_unit_singlept_data.py
index dc6d655408..bf29ced331 100755
--- a/python/ctsm/test/test_unit_singlept_data.py
+++ b/python/ctsm/test/test_unit_singlept_data.py
@@ -18,7 +18,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
-from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
+from ctsm.pft_utils import MAX_PFT_GENERICCROPS, MAX_PFT_MANAGEDCROPS
 
 # pylint: disable=invalid-name
 
@@ -39,7 +39,7 @@ class TestSinglePointCase(unittest.TestCase):
     dom_pft = [8]
     evenly_split_cropland = False
     pct_pft = None
-    num_pft = 16
+    num_pft = MAX_PFT_GENERICCROPS
     cth = [0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9]
     cbh = [0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1, 0.1]
     include_nonveg = False
@@ -132,7 +132,7 @@ def test_check_dom_pft_too_big(self):
             out_dir=self.out_dir,
             overwrite=self.overwrite,
         )
-        single_point.dom_pft = [16, 36, 79]
+        single_point.dom_pft = [MAX_PFT_GENERICCROPS, 36, 79]
         with self.assertRaisesRegex(argparse.ArgumentTypeError, "values for --dompft should*"):
             single_point.check_dom_pft()
 
@@ -162,7 +162,7 @@ def test_check_dom_pft_too_small(self):
             out_dir=self.out_dir,
             overwrite=self.overwrite,
         )
-        single_point.dom_pft = [16, 36, -1]
+        single_point.dom_pft = [MAX_PFT_GENERICCROPS, 36, -1]
         with self.assertRaisesRegex(argparse.ArgumentTypeError, "values for --dompft should*"):
             single_point.check_dom_pft()
 
@@ -193,7 +193,7 @@ def test_check_dom_pft_numpft(self):
             overwrite=self.overwrite,
         )
         single_point.dom_pft = [15, 53]
-        single_point.num_pft = 16
+        single_point.num_pft = MAX_PFT_GENERICCROPS
         with self.assertRaisesRegex(argparse.ArgumentTypeError, "Please use --crop*"):
             single_point.check_dom_pft()
 
diff --git a/python/ctsm/test/test_unit_singlept_data_surfdata.py b/python/ctsm/test/test_unit_singlept_data_surfdata.py
index fb6cc15720..d163c29e4f 100755
--- a/python/ctsm/test/test_unit_singlept_data_surfdata.py
+++ b/python/ctsm/test/test_unit_singlept_data_surfdata.py
@@ -23,7 +23,7 @@
 # pylint: disable=wrong-import-position
 from ctsm import unit_testing
 from ctsm.site_and_regional.single_point_case import SinglePointCase
-from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
+from ctsm.pft_utils import MAX_PFT_GENERICCROPS, MAX_PFT_MANAGEDCROPS
 
 # pylint: disable=invalid-name
 # pylint: disable=too-many-lines
@@ -47,7 +47,7 @@ class TestSinglePointCaseSurfaceNoCrop(unittest.TestCase):
     dom_pft = [8]
     evenly_split_cropland = False
     pct_pft = None
-    num_pft = 16
+    num_pft = MAX_PFT_GENERICCROPS
     cth = 0.9
     cbh = 0.1
     include_nonveg = False
diff --git a/python/ctsm/toolchain/gen_mksurfdata_namelist.py b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
index 09bbb9c268..3a405bf5fa 100755
--- a/python/ctsm/toolchain/gen_mksurfdata_namelist.py
+++ b/python/ctsm/toolchain/gen_mksurfdata_namelist.py
@@ -15,7 +15,7 @@
 
 from ctsm.path_utils import path_to_ctsm_root, path_to_cime
 from ctsm.ctsm_logging import setup_logging_pre_config, add_logging_args, process_logging_args
-from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS
+from ctsm.pft_utils import MAX_PFT_GENERICCROPS, MAX_PFT_MANAGEDCROPS
 
 logger = logging.getLogger(__name__)
 
@@ -307,7 +307,7 @@ def main():
 
     # Determine num_pft
     if nocrop_flag:
-        num_pft = "16"
+        num_pft = str(MAX_PFT_GENERICCROPS)
     else:
         num_pft = str(MAX_PFT_MANAGEDCROPS)
     logger.info("num_pft is %s", num_pft)

From ddca30ff7254b628937a4ced60f1d099e701bc1d Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:44:37 -0600
Subject: [PATCH 22/26] plumber2_surf_wrapper: Avoid mentioning 78.

---
 .../ctsm/site_and_regional/plumber2_surf_wrapper.py  | 12 ++++++++----
 python/ctsm/test/test_sys_plumber2_surf_wrapper.py   |  4 ++--
 python/ctsm/test/test_unit_plumber2_surf_wrapper.py  |  8 ++++----
 3 files changed, 14 insertions(+), 10 deletions(-)

diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 86234aae9c..117ddb4a29 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -28,6 +28,7 @@
 # pylint:disable=wrong-import-position
 from ctsm.site_and_regional.plumber2_shared import PLUMBER2_SITES_CSV, read_plumber2_sites_csv
 from ctsm import subset_data
+from ctsm import pft_utils
 
 
 def get_args():
@@ -49,10 +50,13 @@ def get_args():
     )
 
     parser.add_argument(
-        "--78pft",
-        help="Create and/or modify 78-PFT surface datasets (e.g. for a non-FATES run) ",
+        "--crop",
+        help=(
+            f"Create and/or modify {pft_utils.MAX_PFT_MANAGEDCROPS}-PFT ",
+            "surface datasets (e.g. for a non-FATES run)",
+        ),
         action="store_true",
-        dest="pft_78",
+        dest="use_managed_crops",
     )
 
     parser.add_argument(
@@ -162,7 +166,7 @@ def main():
                 str(pctpft1),
             ]
 
-        if not args.pft_78:
+        if not args.use_managed_crops:
             # use surface dataset with 78 pfts, but overwrite to 100% 1 dominant PFT
             # don't set crop flag
             # set canopy top and bottom heights
diff --git a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
index a7dcf12821..12ca561150 100755
--- a/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_sys_plumber2_surf_wrapper.py
@@ -73,11 +73,11 @@ def test_plumber2_surf_wrapper(self):
 
     def test_plumber2_surf_wrapper_78pft(self):
         """
-        Run the entire tool with --78pft.
+        Run the entire tool with --crop.
         CAN ONLY RUN ON SYSTEMS WITH INPUTDATA
         """
 
-        sys.argv = [self.tool_path, "--78pft"]
+        sys.argv = [self.tool_path, "--crop"]
         main()
 
         # How many files do we expect?
diff --git a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
index e3e677a8a6..4b84752edb 100755
--- a/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
+++ b/python/ctsm/test/test_unit_plumber2_surf_wrapper.py
@@ -70,16 +70,16 @@ def test_parser_78pft_false_default(self):
         """
 
         args = get_args()
-        self.assertFalse(args.pft_78)
+        self.assertFalse(args.use_managed_crops)
 
     def test_parser_78pft_true(self):
         """
-        Test that --78pft sets pft_78 to True
+        Test that --crop sets use_managed_crops to True
         """
 
-        sys.argv += ["--78pft"]
+        sys.argv += ["--crop"]
         args = get_args()
-        self.assertTrue(args.pft_78)
+        self.assertTrue(args.use_managed_crops)
 
 
 if __name__ == "__main__":

From b7f7af386b02a5d10ac14eb27e6e2fa32843891f Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:54:29 -0600
Subject: [PATCH 23/26] plumber2_surf_wrapper: Better is_valid_pft() check.

---
 python/ctsm/pft_utils.py                      | 15 +++++++++++-
 .../plumber2_surf_wrapper.py                  | 23 ++++++-------------
 .../PLUMBER2_sites_invalid_pft.csv            |  2 +-
 3 files changed, 22 insertions(+), 18 deletions(-)

diff --git a/python/ctsm/pft_utils.py b/python/ctsm/pft_utils.py
index 9588564c78..40ab8b9f23 100644
--- a/python/ctsm/pft_utils.py
+++ b/python/ctsm/pft_utils.py
@@ -2,7 +2,20 @@
 Constants and functions relating to PFTs
 """
 
-MIN_PFT = 0   # bare ground
+MIN_PFT = 0  # bare ground
+MIN_NAT_PFT = 1  # minimum natural pft (not including bare ground)
 MAX_NAT_PFT = 14  # maximum natural pft
 MAX_PFT_GENERICCROPS = 16  # for runs with generic crops
 MAX_PFT_MANAGEDCROPS = 78  # for runs with explicit crops
+
+
+def is_valid_pft(pft_num, managed_crops):
+    """
+    Given a number, check whether it represents a valid PFT (bare ground OK)
+    """
+    if managed_crops:
+        max_allowed_pft = MAX_PFT_MANAGEDCROPS
+    else:
+        max_allowed_pft = MAX_PFT_GENERICCROPS
+
+    return MIN_PFT <= pft_num <= max_allowed_pft
diff --git a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
index 117ddb4a29..cedc6b25e0 100755
--- a/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
+++ b/python/ctsm/site_and_regional/plumber2_surf_wrapper.py
@@ -28,7 +28,7 @@
 # pylint:disable=wrong-import-position
 from ctsm.site_and_regional.plumber2_shared import PLUMBER2_SITES_CSV, read_plumber2_sites_csv
 from ctsm import subset_data
-from ctsm import pft_utils
+from ctsm.pft_utils import MAX_PFT_MANAGEDCROPS, is_valid_pft
 
 
 def get_args():
@@ -51,10 +51,8 @@ def get_args():
 
     parser.add_argument(
         "--crop",
-        help=(
-            f"Create and/or modify {pft_utils.MAX_PFT_MANAGEDCROPS}-PFT ",
-            "surface datasets (e.g. for a non-FATES run)",
-        ),
+        help=f"Create and/or modify {MAX_PFT_MANAGEDCROPS}-PFT "
+        "surface datasets (e.g. for a non-FATES run)",
         action="store_true",
         dest="use_managed_crops",
     )
@@ -89,13 +87,6 @@ def execute(command):
     subset_data.main()
 
 
-def is_valid_pft(pft_num):
-    """
-    Given a number, check whether it represents a valid PFT
-    """
-    return pft_num >= 1
-
-
 def main():
     """
     Read plumber2_sites from csv, iterate through sites, and add dominant PFT
@@ -135,7 +126,7 @@ def main():
 
         # Read info for first PFT
         pft1 = row["pft1"]
-        if not is_valid_pft(pft1):
+        if not is_valid_pft(pft1, args.use_managed_crops):
             raise RuntimeError(f"pft1 must be a valid PFT; got {pft1}")
         pctpft1 = row["pft1-%"]
         cth1 = row["pft1-cth"]
@@ -143,13 +134,13 @@ def main():
 
         # Read info for second PFT, if a valid one is given in the .csv file
         pft2 = row["pft2"]
-        if is_valid_pft(pft2):
+        if is_valid_pft(pft2, args.use_managed_crops):
             pctpft2 = row["pft2-%"]
             cth2 = row["pft2-cth"]
             cbh2 = row["pft2-cbh"]
 
         # Set dominant PFT(s)
-        if is_valid_pft(pft2):
+        if is_valid_pft(pft2, args.use_managed_crops):
             subset_command += [
                 "--dompft",
                 str(pft1),
@@ -170,7 +161,7 @@ def main():
             # use surface dataset with 78 pfts, but overwrite to 100% 1 dominant PFT
             # don't set crop flag
             # set canopy top and bottom heights
-            if is_valid_pft(pft2):
+            if is_valid_pft(pft2, args.use_managed_crops):
                 subset_command += [
                     "--cth",
                     str(cth1),
diff --git a/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv
index 2d4b7dcb57..e8f0eb8fbb 100644
--- a/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv
+++ b/python/ctsm/test/testinputs/plumber2_surf_wrapper/PLUMBER2_sites_invalid_pft.csv
@@ -4,5 +4,5 @@
 #ATM_NCPL is specified so that the time step of the model matches the time interval specified by the atm forcing data.
 #longitudes must be in the range [-180,180]
 ,Site,Lat,Lon,pft1,pft1-%,pft1-cth,pft1-cbh,pft2,pft2-%,pft2-cth,pft2-cbh,start_year,end_year,RUN_STARTDATE,START_TOD,ATM_NCPL
-26,Invalid-Pft,51.309166, 4.520560,0,19.22,21.00,10.50,7,80.78,21.00,12.08,2004,2014,2003-12-31,82800,48 
+26,Invalid-Pft,51.309166, 4.520560,-1,19.22,21.00,10.50,7,80.78,21.00,12.08,2004,2014,2003-12-31,82800,48 
 27,BE-Lon,50.551590, 4.746130,15,100.00, 0.50, 0.01,-999,-999.00,-999.00,-999.00,2005,2014,2004-12-31,82800,48 

From 75db098206b064b8b7b2a0604d3f0bf8fdb950cc Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:55:30 -0600
Subject: [PATCH 24/26] Reformat with black.

---
 python/ctsm/site_and_regional/single_point_case.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index 8c30b4c9a3..9f709be648 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -424,7 +424,7 @@ def modify_surfdata_atpoint(self, f_orig):
         if self.dom_pft is not None:
             max_dom_pft = max(self.dom_pft)
             # -- First initialize everything:
-            if max_dom_pft <= MAX_NAT_PFT :
+            if max_dom_pft <= MAX_NAT_PFT:
                 f_mod["PCT_NAT_PFT"][:, :, :] = 0
             else:
                 f_mod["PCT_CFT"][:, :, :] = 0

From 0eb376d61cc0ba3cc91af5c75a9cdbbc4314f0d1 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:55:52 -0600
Subject: [PATCH 25/26] Add previous commit to .git-blame-ignore-revs.

---
 .git-blame-ignore-revs | 1 +
 1 file changed, 1 insertion(+)

diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs
index 6cffe9dd35..518de3672c 100644
--- a/.git-blame-ignore-revs
+++ b/.git-blame-ignore-revs
@@ -68,3 +68,4 @@ cdf40d265cc82775607a1bf25f5f527bacc97405
 3dd489af7ebe06566e2c6a1c7ade18550f1eb4ba
 742cfa606039ab89602fde5fef46458516f56fd4
 4ad46f46de7dde753b4653c15f05326f55116b73
+75db098206b064b8b7b2a0604d3f0bf8fdb950cc

From c3258057f33952d1929a57e7baab803b460325a3 Mon Sep 17 00:00:00 2001
From: Sam Rabin <sam.rabin@gmail.com>
Date: Fri, 20 Jun 2025 17:57:08 -0600
Subject: [PATCH 26/26] Resolve pylint complaint.

---
 python/ctsm/site_and_regional/single_point_case.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/python/ctsm/site_and_regional/single_point_case.py b/python/ctsm/site_and_regional/single_point_case.py
index 9f709be648..d71d014f36 100644
--- a/python/ctsm/site_and_regional/single_point_case.py
+++ b/python/ctsm/site_and_regional/single_point_case.py
@@ -200,7 +200,8 @@ def check_dom_pft(self):
                 If any dom_pft is less than 1.
             Error (ArgumentTypeError):
                 If mixed land units are chosen.
-                dom_pft values are both in range of (0 - MAX_NAT_PFT) and (MAX_NAT_PFT+1 - MAX_PFT_MANAGEDCROPS).
+                dom_pft values are both in range of
+                (0 - MAX_NAT_PFT) and (MAX_NAT_PFT+1 - MAX_PFT_MANAGEDCROPS).
 
 
         """