ec-jrc / lisflood-usecases

Lisflood Test Catchments use cases
European Union Public License 1.2
16 stars 11 forks source link

Fixes in usecases settings #3

Closed domeniconappo closed 4 years ago

domeniconappo commented 5 years ago

Settings files still have paths of grid/hpc platforms, jrc network drives etc.

Need to use the new built-in variable like $(SettingsPath) to make paths relatives so that user don't need to edit anything.

domeniconappo commented 5 years ago

Sub-task in the context of #2

domeniconappo commented 5 years ago

Issue in forcings:

File "/DATA/virtualenvs/ttt/local/lib/python2.7/site-packages/lisflood/hydrological_modules/miscInitial.py", line 179, in initial raise Exception("If using projected coordinates (x, y), a variable with the 'proj4_params' attribute must be included in the precipitation file!") Exception: If using projected coordinates (x, y), a variable with the 'proj4_params' attribute must be included in the precipitation file!

domeniconappo commented 5 years ago

@valeriolorini could you produce again those files and add proj4_params attribute? Thanks!!

valeriolorini commented 5 years ago

@domeniconappo we should check with other cases too. The code has always worked fine with maps without spatial variable (which I am not sure it is a good practice as proj is specified in a metadata to main var). This could cause problem for EFAS GLOFAS and other users that produce input maps without proj variable (I saw ERA5 variables cause same problem). Why was this variable used? if it was a commit can't we use something else like metadata from variable?

domeniconappo commented 5 years ago

@valeriolorini @emiliano-gelati

Hi, proj4_params are used by Proj4 in this code by emiliano (not sure if it was there before..commit is end of April 2019):

https://github.com/ec-jrc/lisflood-code/blame/master/src/lisflood/hydrological_modules/miscInitial.py#L173-L184

I see that we use that convention in write netcdf functions (see lisvap, lisflood...I also use the same in pcr2nc) and it seems a standard to add a Proj4 string to a variable in netCDF. Code is looking for any variable in the file having that auxiliary attribute.

Maybe problem is that it looks for that attribute even if it's not really needed? @emiliano-gelati could you help?

valeriolorini commented 5 years ago

I couldn't find a specific recommendation from ucar about creating a variable about projection in the file. Maybe for latitude extraction we can use metadata?

emiliano-gelati commented 5 years ago

@valeriolorini @domeniconappo The code lines reading the "proj4" string were introduced in the snow/ice-melt bug-fix. The "proj4" string is necessary to compute the pixels' latitude for domains defined on projected coordinates (x,y). It is not used for geographic coordinates (lon,lat), so it should be not looked for when using ERA5: see https://github.com/ec-jrc/lisflood-code/blame/master/src/lisflood/hydrological_modules/miscInitial.py#L173-L184

The code was developed to use the metadata contained in the EFAS forcing NetCDF files. Has the format of these files changed?

valeriolorini commented 5 years ago

@emiliano-gelati @valeriolorini thanks, I thought the problem was relevant to Global lisflood run too. nevertheless maps which were cut for (x,y) use case didn't contain it neither (according to error above). shall we suggest users to recreate cookie-cut forcings maps including a variable for that? @ecCinziaMazzetti could you check the maps you sent for the use-case?

domeniconappo commented 5 years ago

Hi there,

we need to define a format to adopt (and to provide users with a validator to check their input files).

Forcings are not guaranteed to have a variable with that attribute so we could either:

  1. Enforce the adoption of a projection variable with all geo attributes within netcdf forcings (user might need rebuilding files to add that variable). Or, probably best thing, enforce adoption of CF-1.7 standard about projection variables (http://cfconventions.org/Data/cf-conventions/cf-conventions-1.7/build/cf-conventions.html#grid-mappings-and-projections)
  2. Use auxiliary maps as it was in the past (e.g. a lat.map/lat.nc file). User must compile his own latitude map.
  3. Let the user to define projection string (proj4/wkt format) inside settings xml
  4. Let the user to choose between a predefined set of well known and mostly used projections, simply by setting a variable in settings file (e.g. <textvar name="Projection" value="lambert_conformal_conic" />)

Real wkt/proj4 strings will be kept in lisflood code and user needs to know only projection name.

@emiliano-gelati @valeriolorini @ecCinziaMazzetti Can we have a briefing about this? We are stuck at this test and don't know how to proceed.

emiliano-gelati commented 5 years ago

Hi @domeniconappo, @ecCinziaMazzetti , @valeriolorini , @AdDeRoo

I propose going in the direction of proposal 1 by @domeniconappo, with a slight modification.

Current situation: Now the model looks for a NetCDF variable in the precipitation forcing file with the "'proj4_params" attribute if the simulated domain is defined on projected coordinates (x, y).

Proposal: We could enforce that a NetCDF file (still precipitation or any other) needs to have a global attribute (not attached to a specific variable, but to the whole file) holding the Proj4 string with the projection parameters. This is easy and computationally inexpensive to implement on any existing NetCDF file, e.g.:

import netCDF4
proj4_string = '+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs'
with netCDF4.Dataset('any_input.nc', 'a') as nc:
    nc.setncattr('proj4', proj4_string)

This would allow to keep using Proj4 strings and pyproj, which are extremely handy; and minimise code changes in https://github.com/ec-jrc/lisflood-code/blob/master/src/lisflood/hydrological_modules/miscInitial.py#L173-L184 This part would be simplified because it would not be necessary to scan all variables: instead, the model would look for the global attribute (properly named, "proj4" is an example), if the coordinates are x and y.

domeniconappo commented 5 years ago

@emiliano-gelati Solution is handy but I think it is worthy to adopt an established convention, like the CF-1.x. But in any case, we just need to decide and then clearly document the decision.

@gnrgomes In case we go for it, would it be hard for you to reproduce nc files to follow this convention? We would need an attribute at file level called proj4 holding a proj4 string. Just see the code snippet by Emiliano.

gnrgomes commented 5 years ago

@domeniconappo and @emiliano-gelati It isn't hard to change the files neither the script that produces them. I'm producing the netCDF files CF-1.6 compliant. I couldn't find if having proj4_string at the file attributes level is allowed or not.

I have this required info inside the file already:

variables: int lambert_azimuthal_equal_area; :grid_mapping_name = "lambert_azimuthal_equal_area"; :false_easting = 4321000.0; // double :false_northing = 3210000.0; // double :longitude_of_projection_origin = 10.0; // double :latitude_of_projection_origin = 52.0; // double :semi_major_axis = 6378137.0; // double :inverse_flattening = 298.257223563; // double :proj4_params = "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs";

emiliano-gelati commented 5 years ago

Hi @gnrgomes , @domeniconappo , @valeriolorini , @AdDeRoo If the current format, which we have used so far to force LISFLOOD on the Euro-Mediterranean domain, is CF-compliant, then I suggest we keep using it. Emiliano

domeniconappo commented 5 years ago

@emiliano-gelati @valeriolorini @gnrgomes @AdDeRoo

Yes, the only thing is that CF doesn't mention the proj4_params attribute. To be really compliant, we should compose proj4 string from other attributes (in case is not there) and it can be tricky because numbers and names of those attributes depend on the specific grid mapping (grid_mapping_name).

Any idea? For the moment we can keep the code as it is and keeping in mind that's not compliant and soon or later someone could complain :)

domeniconappo commented 5 years ago

Issue reported in https://github.com/ec-jrc/lisflood-usecases/issues/3#issuecomment-494450402 was fixed in lisflood (v. 2,8.22 https://github.com/ec-jrc/lisflood-code/issues/27).

There are more issues while trying to run this use case with current version of lisflood. All of these need to be addressed.

A new exception while reading ForestFraction maps/values

Traceback (most recent call last):
  File "/workarea/virtualenvs/lisflood27/bin/lisflood", line 7, in <module>
    exec(compile(f.read(), __file__, 'exec'))
  File "/workarea/Projects/lisflood-code/bin/lisflood", line 41, in <module>
    sys.exit(main())
  File "/workarea/Projects/lisflood-code/src/lisflood/main.py", line 224, in main
    Lisfloodexe(settings, optionxml)
  File "/workarea/Projects/lisflood-code/src/lisflood/main.py", line 78, in Lisfloodexe
    Lisflood = LisfloodModel()
  File "/workarea/Projects/lisflood-code/src/lisflood/Lisflood_initial.py", line 128, in __init__
    self.landusechange_module.initial()
  File "/workarea/Projects/lisflood-code/src/lisflood/hydrological_modules/landusechange.py", line 46, in initial
    self.var.ForestFraction = loadmap('ForestFraction',timestampflag='closest')
  File "/workarea/Projects/lisflood-code/src/lisflood/global_modules/add1.py", line 441, in loadmap
    mapC = compressArray(mapnp,pcr=False,name=filename)
  File "/workarea/Projects/lisflood-code/src/lisflood/global_modules/add1.py", line 253, in compressArray
    if np.max(np.isnan(mapC)):
  File "/workarea/virtualenvs/lisflood27/lib/python2.7/site-packages/numpy/core/fromnumeric.py", line 2505, in amax
    initial=initial)
  File "/workarea/virtualenvs/lisflood27/lib/python2.7/site-packages/numpy/core/fromnumeric.py", line 86, in _wrapreduction
    return ufunc.reduce(obj, axis, dtype, out, **passkwargs)
ValueError: zero-size array to reduction operation maximum which has no identity

@emiliano-gelati @valeriolorini @menta78 Any idea about this error?

domeniconappo commented 5 years ago

We have some issues on static dataset. Need to recover original data https://efascom.smhi.se/jira/browse/ECC-1354

domeniconappo commented 5 years ago

ETRS89 dataset has issues (to check: all forest fraction maps). We are going to re-cut files for that catchment from the entire EFAS domain as issue might be in how nc files were cut.

domeniconappo commented 4 years ago

both use cases run with no errors with latest OS lisflood model