PCT_CLAY and PCT_SAND are truncated in the creation of regional surface datasets

[olyson]

Brief summary of bug

A couple of users have reported that PCT_CLAY and PCT_SAND fields are truncated in the creation of regional surface datasets. See below for what I mean by truncated.
The relevant DIscussCESM post is here:

https://bb.cgd.ucar.edu/cesm/threads/mksurfdata_esmf-tool-missing-percent-clay-and-sand.8597/

General bug information

CTSM version you are using: I'm using alpha-ctsm5.2.mksrf.25_ctsm5.1.dev171 to replicate the problem

Does this bug cause significantly incorrect results in the model's science? Yes

Configurations affected: Simulations that use regional surface datasets

Details of bug

The problem is illustrated by the PCT_SAND field for the top soil layer from a U.S. surface dataset I created following a user's approach:

Note the northern and southern parts of the domain where PCT_SAND is 43. The value 43 appears to be the value that is set when there is no soils data, according to the log file.
Maybe it has to do with the ocean/land mask? Or maybe we are missing a step in the process below? Suggestions for further debugging welcome.

Important details of your setup / configuration so we can reproduce the bug

These are the steps I took to create the surface dataset, following a user's approach:

cd /glade/work/oleson/ctsm5.2.mksurfdata_032124/tools/mksurfdata_esmf
./gen_mksurfdata_build.sh
Use scrip.nc from /glade/work/amans/data/MRB-CONUS_0.05 which is a 0.05x0.05 grid over the U.S.
/glade/u/apps/derecho/23.09/spack/opt/spack/esmf/8.6.0/cray-mpich/8.1.27/oneapi/2023.2.1/7haa/bin/ESMF_Scrip2Unstruct scrip.nc lnd_mesh.nc 0
ncap2 -s 'elementMask(:)=0' lnd_mesh.nc mask_mesh.nc
ncdump -v grid_dims scrip.nc (to get nx and ny)
conda activate ctsm_pylib
./gen_mksurfdata_namelist --start-year 2000 --end-year 2000 --model-mesh-nx 770 --model-mesh-ny 470 --model-mesh lnd_mesh.nc --res 0.05-conus-nomask
./gen_mksurfdata_jobscript_single --number-of-nodes 24 --tasks-per-node 12 --namelist-file surfdata_0.05-conus-nomask_hist_2000_78pfts_c240321.namelist
qsub ./mksurfdata_jobscript_single

[slevis-lmwg]

I remember discussing this with the first user who brought it up in the same thread. Erik and I decided at the time to recommend the ./subset_data workaround. When we decide that the workaround is insufficient, we should try to resolve the issue in mksurfdata_esmf.

[olyson]

I did try the workaround myself and was able to create a global surface dataset that had good fields of PCT_SAND/PCT_CLAY. The dataset at that desired resolution was nearly a TB in size, but I guess you can delete it after you've done the subsetting.

[wwieder]

In the SE meeting Erik thought this could be addressed in mksurfdata_esmf and wondered if it should abort for regional cases?

Sam R, wondered if this is caused by a slight mismatch between grids that xarray doesn't know how to handle, producing fill_values.

Sam L noted we're not using xarray, but calling esmf directly from fortran.

[olyson]

Another example of this problem at 0.01x0.01 resolution:

https://bb.cgd.ucar.edu/cesm/threads/creating-mapfiles-by-using-the-mkmapdata-sh.6546/#post-54853

[ekluzek]

@billsacks in mkmapdata we had to tell ESMF that the grids were regional -- I'm wondering if mksurfdata_esmf is just missing some type of flag to signal that the output grid is regional? Can you ask the ESMF team about this?

[billsacks]

As far as I can tell, there is no equivalent "regional" flag for the inline regridding used in mksurfdata_esmf. How do other fields look for a regional grid like this? Is it just the fields created by mksoiltexMod? If so, I think we should look into the details of what's going on there... note that it does something different from most mapping, in that it finds the appropriate mapunit and then assigns values based on that mapunit; this leverages ESMF dynamic masking in a way that I haven't come to understand.

One thing I'm curious about is whether you get the same problematic pattern with different processor counts. I'm not sure what I'd do with that information, but I'm wondering if there are somehow issues on certain processors of the decomposition for these fields.

[slevis-lmwg]

@olyson I am posting here because this issue may relate to the rest of the soil-texture issues that surfaced and got resolved last week. It's possible that this one is now resolved, too. We can discuss at our Wednesday meeting.

[billsacks]

Oh yes, I'd bet that the fix last Friday will resolve this one, too!

[slevis-lmwg]

By the way @billsacks and @mvertens:
BIG thanks for resolving #2744 so quickly!

[ekluzek]

Now, that ctsm5.3.0 is made it would be great if someone could test mksurfdat_esmf in it to see if this problem is still there. @olyson would you be up for that? Alternatively we could post about this on the Forum issue in question and ask them to try it out. Thoughts on doing that @olyson and @slevis-lmwg?

[olyson]

I have notes on this particular application, so I will try it out, thanks.

[olyson]

I generated a surface dataset for the U.S. as per the original application using mksurfdata_esmf in ctsm5.3.0. PCT_SAND, PCT_CLAY, and BULK look beautiful.

[olyson]

BULK

[ekluzek]

Awesome, thanks @olyson. Closing as fixed in ctsm5.3.0