Models crash when "lat" array in forcing dataset is decreasing instead of increasing

[nikizadehgfdl]

This is a feature request to make models run when the coordinate "lat" in the forcing dataset is monotonically decreasing from +90 to -90.

A user's model run crashes as follows because the "lat" coordinate in the forcing file is monotonic but decreasing from +90 to -90,
unlike the forcing files we have been using till now which has "lat" increasing from -90 to +90.

 => Error: array must be monotonically increasing in "indp"
            when searching for nearest element to value=  -1.37252691311570
            array(i) < array(i-1) for i=           2
            array(i) for i=1..ia follows:
lat=89.82086, 89.58836, 89.35471, 89.12084, 88.88688, 88.65273, 88.41859,
    88.18436, 87.95024, 87.71603, 87.48183, 87.24764, 87.01344, 86.77924,
...
    -87.71603, -87.95024, -88.18436, -88.41859, -88.65273, -88.88688,
    -89.12084, -89.35471, -89.58836, -89.82086 ;
...
 lon = 0, 0.234375, 0.46875, 0.703125, 0.9375, 1.171875, 1.40625, 1.640625,
    1.875, 2.109375, 2.34375, 2.578125, 2.8125, 3.046875, 3.28125, 3.515625,
...
    356.9531, 357.1875, 357.4219, 357.6562, 357.8906, 358.125, 358.3594,
    358.5938, 358.8281, 359.0625, 359.2969, 359.5312, 359.7656 ;

forrtl: severe (174): SIGSEGV, segmentation fault occurred
Image              PC                Routine            Line        Source
MOM6               00000000019429BC  horiz_interp_bili        1385  horiz_interp_bilinear.F90
MOM6               00000000015EA759  horiz_interp_mod_         503  horiz_interp.F90
MOM6               00000000013C2D97  data_override_mod         990  data_override.F90
MOM6               00000000013C82E8  data_override_mod         636  data_override.F90
MOM6               00000000009AB80F  atm_land_ice_flux         892  atm_land_ice_flux_exchange.F90
MOM6               0000000000FE6BB2  MAIN__                    816  coupler_main.F90

Querying the new forcing file shows that "lat" array is decreasing:

ncdump /lustre/f2/dev/ncep/Hae-Cheol.Kim/MOM6_exercises/GLBb0.08_041/INPUT+09m11ob2/out_200001_GEFS2JRA.nc | grep lat

 lat = 89.82086, 89.58836, 89.35471, 89.12084, 88.88688, 88.65273, 88.41859, 
    88.18436, 87.95024, 87.71603, 87.48183, 87.24764, 87.01344, 86.77924, 
    86.54504, 86.31081, 86.07658, 85.84239, 85.60818, 85.37395, 85.13973, 
...
   -86.31081, -86.54504, -86.77924, -87.01344, -87.24764, -87.48183, 
    -87.71603, -87.95024, -88.18436, -88.41859, -88.65273, -88.88688, 
    -89.12084, -89.35471, -89.58836, -89.82086 ;

The forcing files in use at gfdl models have "lat" array increasing. E.g.,

 ncdump $PDATA/gfdl/gfdl_O/datasets/reanalysis/JRA55-do/v1.4.0/short_sample/tas_JRA_sample.nc | more

 lat = -89.5700895506066, -89.013176131022, -88.452973836713, 
...
    87.8920284453444, 88.452973836713, 89.013176131022, 89.5700895506066 ;

Describe the solution you'd like
Can we support forcing files with decreasing "lat"?

Describe alternatives you've considered
Suggest a way to the user to make the forcing file comply with the current restriction.

[nikizadehgfdl]

[email protected] says: The data is from NCEP/EMC. It is from Global Ensemble Forecast System (https://www.ncei.noaa.gov/products/weather-climate-models/global-ensemble-forecast). It is 1/4-th degree but the products were repackaged in CDEPS-compliant format by EMC folks.

@bensonr , @uramirez does it make sense to add support for this data. Do we have horiz_interp unit tests that can digest in this dataset? I do not have a test case for this but could probably come up with a small MOM6 test.

[mlee03]

@rem1776 your horz_interp unit tests could test monotonically decreasing latlon right?

[rem1776]

@mlee03 Yeah it could potentially, with mixed mode it added a test for the bilinear method that it's using here.

If we wanted a unit test to actually be reading in the data though, i think that would be in data_override.

[uramirez8707]

I think support for monotonically (increasing/decreasing) data can easily be added to horiz interp. The code that searches for the nearest data point to the point interpolating to requires the data array to be increasing but there is no reason for that restriction as long it is monotonic.

There’s some code in diag manager than does this already (to find the points of the subregion) that we can reuse.

@bensonr is this something that we want to support?

[bensonr]

@uramirez8707 - I think data_override should be able support both monotonically increasing and decreasing latitudinal representations. I hope the algorithm for weight generation isn't coded in such a way that only one representation is valid.

[StevePny]

The datasets from the ECMWF MARS server come as monotonically decreasing as well, so there is another use case. The ECMWF datasets also have a rotated longitude, which would be useful to support too. At the moment we are doing these as preprocessing steps with NCO:

# ----------------------------------------------------------
# 2. reorder longitudes
# ----------------------------------------------------------

docker run --rm --name=nco --volume=${root}/${datadir}/:/workdir \
  --workdir=/workdir -e YYYY -e MM -e DD -e HH ${nco_image} \
  ncks -O --msa -d longitude,0.,180. -d longitude,-180.,-0.125 \
  SFC_FCST_2D.nc IFS_SFC_FCST_${YYYY}${MM}${DD}_${HH}Z.nc

docker run --rm --name=nco --volume=${root}/${datadir}/:/workdir \
  --workdir=/workdir -e YYYY -e MM -e DD -e HH  ${nco_image} \
  ncap2 -O -s 'where(longitude<0) longitude=360+longitude' \
  IFS_SFC_FCST_${YYYY}${MM}${DD}_${HH}Z.nc IFS_SFC_FCST_${YYYY}${MM}${DD}_${HH}Z.nc

# ----------------------------------------------------------
# 3. flip the latitude axis
# ----------------------------------------------------------

docker run --rm --name=nco --volume=${root}/${datadir}/:/workdir \
  --workdir=/workdir -e YYYY -e MM -e DD -e HH ${nco_image} \
  ncpdq -O -a -latitude IFS_SFC_FCST_${YYYY}${MM}${DD}_${HH}Z.nc IFS_SFC_FCST_${YYYY}${MM}${DD}_${HH}Z.nc