MESH_parameters now can be read as a NetCDF file

[kasra-keshavarz]

The MESH_parameters.nc file can now be prepared and fed into MESH >r1860 versions.

Reported by: @fuadyassin

[kasra-keshavarz]

An example shared by @zktesemma:

netcdf MESH_parameters {
dimensions:
    subbasin = 2122 ;
    nsol = 4 ;
    ngru = 21 ;
variables:
    float lon(subbasin) ;
        lon:units = "degrees_east" ;
        lon:_FillValue = -1.f ;
        lon:long_name = "longitude" ;
        lon:standard_name = "longitude" ;
        lon:axis = "X" ;
    float lat(subbasin) ;
        lat:units = "degrees_north" ;
        lat:_FillValue = -1.f ;
        lat:long_name = "latitude" ;
        lat:standard_name = "latitude" ;
        lat:axis = "Y" ;
    float time(subbasin) ;
        time:units = "days since 1980-10-01 00:00:00.0 -0:00" ;
        time:_FillValue = -1.f ;
        time:standard_name = "time" ;
        time:axis = "T" ;
    float clay(subbasin, nsol) ;
        clay:units = "%" ;
        clay:_FillValue = -1.f ;
        clay:long_name = "Clay Content of Soil Layer" ;
        clay:grid_mapping = "crs" ;
        clay:standard_name = "clay" ;
        clay:coordinates = "lat lon time" ;
    float sand(subbasin, nsol) ;
        sand:units = "%" ;
        sand:_FillValue = -1.f ;
        sand:long_name = "Sand Content of Soil Layer" ;
        sand:grid_mapping = "crs" ;
        sand:standard_name = "sand" ;
        sand:coordinates = "lat lon time" ;
    float orgm(subbasin, nsol) ;
        orgm:units = "%" ;
        orgm:_FillValue = -1.f ;
        orgm:long_name = "Organic Matter Content of Soil Layer" ;
        orgm:grid_mapping = "crs" ;
        orgm:standard_name = "orgm" ;
        orgm:coordinates = "lat lon time" ;
    float elevation(subbasin, ngru) ;
        elevation:units = "m" ;
        elevation:_FillValue = -1.f ;
        elevation:long_name = "Average Elevation of GRU" ;
        elevation:grid_mapping = "crs" ;
        elevation:standard_name = "elevation" ;
        elevation:coordinates = "lat lon time" ;
    float slope(subbasin, ngru) ;
        slope:units = "degree" ;
        slope:_FillValue = -1.f ;
        slope:long_name = "Average slope of GRU" ;
        slope:grid_mapping = "crs" ;
        slope:standard_name = "slope" ;
        slope:coordinates = "lat lon time" ;
    float aspect(subbasin, ngru) ;
        aspect:units = "degree" ;
        aspect:_FillValue = -1.f ;
        aspect:long_name = "Average Aspect of GRU" ;
        aspect:grid_mapping = "crs" ;
        aspect:standard_name = "aspect" ;
        aspect:coordinates = "lat lon time" ;
    float delta(subbasin, ngru) ;
        delta:units = "m" ;
        delta:_FillValue = -1.f ;
        delta:long_name = "Change in elevation between basin and of GRU" ;
        delta:grid_mapping = "crs" ;
        delta:standard_name = "delta" ;
        delta:coordinates = "lat lon time" ;
    float curvature(subbasin, ngru) ;
        curvature:units = "-" ;
        curvature:_FillValue = -1.f ;
        curvature:long_name = "Average curvature parameter of GRU" ;
        curvature:grid_mapping = "crs" ;
        curvature:standard_name = "curvature" ;
        curvature:coordinates = "lat lon time" ;
    float skyviewfactor(subbasin, ngru) ;
        skyviewfactor:units = "-" ;
        skyviewfactor:_FillValue = -1.f ;
        skyviewfactor:long_name = "Average skyviewfactor parameter of GRU" ;
        skyviewfactor:grid_mapping = "crs" ;
        skyviewfactor:standard_name = "skyviewfactor" ;
        skyviewfactor:coordinates = "lat lon time" ;
    float dd(subbasin) ;
        dd:units = "m m-2" ;
        dd:_FillValue = -1.f ;
        dd:long_name = "Subbasin average drainage density" ;
        dd:grid_mapping = "crs" ;
        dd:standard_name = "dd" ;
        dd:coordinates = "lat lon time" ;
    char GRUsClassName(ngru) ;
    int crs ;
        crs:grid_mapping_name = "latitude_longitude" ;
        crs:longitude_of_prime_meridian = 0 ;
        crs:semi_major_axis = 6378137 ;
        crs:inverse_flattening = 298.257223563 ;
// global attributes:
        :Conventions = "CF-1.0" ;
        :source = "MERIT Hydro for DEM and COE for land cover" ;
        :institution = "https://gwf.usask.ca/" ;
        :references = "xx et al. (xxxx) journal xx:xx-xx" ;
        :history = "Z.K. Tesemma, Thu Dec 14 20:49:09 2023" ;
        :featureType = "point"
}

[kasra-keshavarz]

The following is relevant: MESH-Model/MESH-Dev#55 (comment)

summarizing the email communication in response to the above:

for MESH_paramaters.nc, soil texture should be provided with an additional dimension called "nsol" corresponding to the number of levels/layers. This means 3D (lat, lon, level) for gridded setups, and 2D (subbasin, level) for vector setups.

The code expects though, that values are provided for all levels as per "MESH_input_soil_levels.txt" which means it behaves differently than using MESH_parameters.r2c or CLASS.ini. If the number of soil levels is large (as for permafrost applications) this is not very convenient. It also means experimenting with the number of soil levels comes requires recreating the texture parameters for each change.

In MESH_parameters.nc, %Sand, for example, is provided with parameters: sand, sand 1. sand 2, ... sand n. the first field is assigned to all layers then sand 1 overrides it for layer 1, sand 2 overrides it for layer 2; thus sand is used for all layers > n+1.

Similarly if texture is read from CLASS.ini, there is the NRSOILAYEREADFLAG to specify how many layers of data to be provided in CLASS.ini file. If there are more layers, they get assigned the value of the last layer as per NRSOILAYEREADFLAG.

In my opinion, behavior should be consistent weather soil texture data are supplied in CLASS.ini, MESH_parameters.r2c or MESH_parameters.nc. The three formats are not mutually exclusive. Special GRUs (with MID>100) keep their parameter values that were taken from CLASS,ini. One can supply INPUTPARAMSFORMFLAG r2c nc txt so it will read three additional parameter sets from those three files and will do all sorts of overrides.