search

GeoscienceAustralia / anuga_core

AnuGA for the simulation of the shallow water equation
https://anuga.anu.edu.au
Other
192 stars 94 forks source link

GDAL generated input format #209

Open schoeller opened 5 years ago

schoeller commented 5 years ago

Dear all,

thanks for this OSS product.

I am trying to use GeoTIFF data for the terrain and convert it through GDAL to be usable by anuga.dem2pts directly. I took the cairns.dem file as target file format. gdalinfo delivers information as follows:

gdalinfo inputs\cairns.dem
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
Driver: netCDF/Network Common Data Format
Files: inputs\cairns.dem
Size is 1584, 1188
Coordinate System is `'
Metadata:
  NC_GLOBAL#cellsize=270.90250589988
  NC_GLOBAL#datum=WGS84
  NC_GLOBAL#description=NetCDF DEM format for compact and portable storage of spatial point data
  NC_GLOBAL#false_easting=500000
  NC_GLOBAL#false_northing=10000000
  NC_GLOBAL#institution=Geoscience Australia
  NC_GLOBAL#ncols=1584
  NC_GLOBAL#NODATA_value=-9999
  NC_GLOBAL#nrows=1188
  NC_GLOBAL#projection=UTM
  NC_GLOBAL#units=METERS
  NC_GLOBAL#xllcorner=285455.22987171
  NC_GLOBAL#yllcorner=7953325.2470893
  NC_GLOBAL#zone=55
Corner Coordinates:
Upper Left  (    0.0,    0.0)
Lower Left  (    0.0, 1188.0)
Upper Right ( 1584.0,    0.0)
Lower Right ( 1584.0, 1188.0)
Center      (  792.0,  594.0)
Band 1 Block=1584x1 Type=Float64, ColorInterp=Undefined
  NoData Value=9.969209968386869e+036
  Metadata:
    NETCDF_VARNAME=elevation

Based on above information the following parameters were applied to gdal_convert:

gdal_translate -of netCDF -co FORMAT=NC -ot Float64 N11E008.tif N11E008.dem

The information from gdalinfo for the conversion output now reads as follows:

gdalinfo inputs\N11E008.dem
Driver: netCDF/Network Common Data Format
Files: inputs\N11E008.dem
Size is 1201, 1201
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (7.999583333333334,12.000416666666666)
Pixel Size = (0.000833333333333,-0.000833333333333)
Metadata:
  Band1#grid_mapping=crs
  Band1#long_name=GDAL Band Number 1
  Band1#_FillValue=-32768
  crs#GeoTransform=7.999583333333334 0.000833333333333 0 12.00041666666667 0 -0.000833333333333
  crs#grid_mapping_name=latitude_longitude
  crs#inverse_flattening=298.257223563
  crs#longitude_of_prime_meridian=0
  crs#long_name=CRS definition
  crs#semi_major_axis=6378137
  crs#spatial_ref=GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]]
  lat#long_name=latitude
  lat#standard_name=latitude
  lat#units=degrees_north
  lon#long_name=longitude
  lon#standard_name=longitude
  lon#units=degrees_east
  NC_GLOBAL#Conventions=CF-1.5
  NC_GLOBAL#GDAL=GDAL 2.2.3, released 2017/11/20
  NC_GLOBAL#GDAL_AREA_OR_POINT=Area
  NC_GLOBAL#history=Tue Sep 17 15:16:39 2019: GDAL CreateCopy( N11E008.dem5, ... )
Corner Coordinates:
Upper Left  (   7.9995833,  12.0004167) (  7d59'58.50"E, 12d 0' 1.50"N)
Lower Left  (   7.9995833,  10.9995833) (  7d59'58.50"E, 10d59'58.50"N)
Upper Right (   9.0004167,  12.0004167) (  9d 0' 1.50"E, 12d 0' 1.50"N)
Lower Right (   9.0004167,  10.9995833) (  9d 0' 1.50"E, 10d59'58.50"N)
Center      (   8.5000000,  11.5000000) (  8d30' 0.00"E, 11d30' 0.00"N)
Band 1 Block=1201x1 Type=Float64, ColorInterp=Undefined
  NoData Value=-32768
  Metadata:
    grid_mapping=crs
    long_name=GDAL Band Number 1
    NETCDF_VARNAME=Band1
    _FillValue=-32768

When applying anuga.dem2pts(inputs_dir + "N11E008.dem", use_cache=False, verbose=True) from a Python script I receive the following error:

Traceback (most recent call last):
  File "runDamBreak.py", line 190, in <module>
    run_dambreak("2")
  File "runDamBreak.py", line 88, in run_dambreak
    anuga.dem2pts(inputs_dir + "N11E008.dem", use_cache=False, verbose=True)
  File "C:\ProgramData\Miniconda3\envs\python2\lib\site-packages\anuga\file_conversion\dem2pts.py", line 50, in dem2pts
    result = apply(_dem2pts, [name_in], kwargs)
  File "C:\ProgramData\Miniconda3\envs\python2\lib\site-packages\anuga\file_conversion\dem2pts.py", line 88, in _dem2pts
    ncols = int(infile.ncols)
  File "netCDF4\_netCDF4.pyx", line 2606, in netCDF4._netCDF4.Dataset.__getattr__
  File "netCDF4\_netCDF4.pyx", line 2552, in netCDF4._netCDF4.Dataset.getncattr
  File "netCDF4\_netCDF4.pyx", line 1188, in netCDF4._netCDF4._get_att
  File "netCDF4\_netCDF4.pyx", line 1638, in netCDF4._netCDF4._ensure_nc_success
AttributeError: NetCDF: Attribute not found

Any hints on how to proceed otherwise are highly welcome. Converting from tif->asc->dem I have not succeeded neither.

Best regards

Sebastian

schoeller commented 5 years ago

For the time being I tried as suggested at https://github.com/Hydrata/ChennaiFloodModel/blob/b13063198e991d5b93006be1539a1bfe24856773/runChennai.py#L85

stoiver commented 5 years ago

Hi Sebastian,

The dem format is not a standard format, it is just an anuga format which provides a binary (netcdf) copy of an asc type data file. So I'm not surprised that gdal had trouble with it. I would suggest using gdal on the asc, prj files instead.

But if you have a tiff file then your use of the utility function qs.composite_quantity_setting_function is the way to go (of course under the hood, this function uses gdal to read tiff file).

Cheers Steve

============================== Stephen Roberts Undergraduate Convenor Mathematical Sciences Institute Room 4.74 Hanna Neumann Building #145 The Australian National University Canberra, ACT 2600 AUSTRALIA Ph: +61 2 61254445 CRICOS: 00120C

From: Sebastian Schoeller notifications@github.com Sent: Tuesday, 17 September 2019 11:35 PM To: GeoscienceAustralia/anuga_core anuga_core@noreply.github.com Cc: Subscribed subscribed@noreply.github.com Subject: [GeoscienceAustralia/anuga_core] GDAL generated input format (#209)

Dear all,

thanks for this OSS product.

I am trying to use GeoTIFF data for the terrain and convert it through GDAL to be usable by anuga.dem2pts directly. I took the cairns.dem file as target file format. gdalinfo delivers information as follows:

gdalinfo inputs\cairns.dem Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute Driver: netCDF/Network Common Data Format Files: inputs\cairns.dem Size is 1584, 1188 Coordinate System is `' Metadata: NC_GLOBAL#cellsize=270.90250589988 NC_GLOBAL#datum=WGS84 NC_GLOBAL#description=NetCDF DEM format for compact and portable storage of spatial point data NC_GLOBAL#false_easting=500000 NC_GLOBAL#false_northing=10000000 NC_GLOBAL#institution=Geoscience Australia NC_GLOBAL#ncols=1584 NC_GLOBAL#NODATA_value=-9999 NC_GLOBAL#nrows=1188 NC_GLOBAL#projection=UTM NC_GLOBAL#units=METERS NC_GLOBAL#xllcorner=285455.22987171 NC_GLOBAL#yllcorner=7953325.2470893 NC_GLOBAL#zone=55 Corner Coordinates: Upper Left ( 0.0, 0.0) Lower Left ( 0.0, 1188.0) Upper Right ( 1584.0, 0.0) Lower Right ( 1584.0, 1188.0) Center ( 792.0, 594.0) Band 1 Block=1584x1 Type=Float64, ColorInterp=Undefined NoData Value=9.969209968386869e+036 Metadata: NETCDF_VARNAME=elevation

Based on above information the following parameters were applied to gdal_convert:

gdal_translate -of netCDF -co FORMAT=NC -ot Float64 N11E008.tif N11E008.dem

The information from gdalinfo for the conversion output now reads as follows:

gdalinfo inputs\N11E008.dem Driver: netCDF/Network Common Data Format Files: inputs\N11E008.dem Size is 1201, 1201 Coordinate System is: GEOGCS["WGS 84", DATUM["WGS_1984", SPHEROID["WGS 84",6378137,298.257223563, AUTHORITY["EPSG","7030"]], AUTHORITY["EPSG","6326"]], PRIMEM["Greenwich",0], UNIT["degree",0.0174532925199433], AUTHORITY["EPSG","4326"]] Origin = (7.999583333333334,12.000416666666666) Pixel Size = (0.000833333333333,-0.000833333333333) Metadata: Band1#grid_mapping=crs Band1#long_name=GDAL Band Number 1 Band1#_FillValue=-32768 crs#GeoTransform=7.999583333333334 0.000833333333333 0 12.00041666666667 0 -0.000833333333333 crs#grid_mapping_name=latitude_longitude crs#inverse_flattening=298.257223563 crs#longitude_of_prime_meridian=0 crs#long_name=CRS definition crs#semi_major_axis=6378137 crs#spatial_ref=GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]] lat#long_name=latitude lat#standard_name=latitude lat#units=degrees_north lon#long_name=longitude lon#standard_name=longitude lon#units=degrees_east NC_GLOBAL#Conventions=CF-1.5 NC_GLOBAL#GDAL=GDAL 2.2.3, released 2017/11/20 NC_GLOBAL#GDAL_AREA_OR_POINT=Area NC_GLOBAL#history=Tue Sep 17 15:16:39 2019: GDAL CreateCopy( N11E008.dem5, ... ) Corner Coordinates: Upper Left ( 7.9995833, 12.0004167) ( 7d59'58.50"E, 12d 0' 1.50"N) Lower Left ( 7.9995833, 10.9995833) ( 7d59'58.50"E, 10d59'58.50"N) Upper Right ( 9.0004167, 12.0004167) ( 9d 0' 1.50"E, 12d 0' 1.50"N) Lower Right ( 9.0004167, 10.9995833) ( 9d 0' 1.50"E, 10d59'58.50"N) Center ( 8.5000000, 11.5000000) ( 8d30' 0.00"E, 11d30' 0.00"N) Band 1 Block=1201x1 Type=Float64, ColorInterp=Undefined NoData Value=-32768 Metadata: grid_mapping=crs long_name=GDAL Band Number 1 NETCDF_VARNAME=Band1 _FillValue=-32768

When applying anuga.dem2pts(inputs_dir + "N11E008.dem", use_cache=False, verbose=True) from a Python script I receive the following error:

Traceback (most recent call last): File "runDamBreak.py", line 190, in run_dambreak("2") File "runDamBreak.py", line 88, in run_dambreak anuga.dem2pts(inputs_dir + "N11E008.dem", use_cache=False, verbose=True) File "C:\ProgramData\Miniconda3\envs\python2\lib\site-packages\anuga\file_conversion\dem2pts.py", line 50, in dem2pts result = apply(_dem2pts, [name_in], kwargs) File "C:\ProgramData\Miniconda3\envs\python2\lib\site-packages\anuga\file_conversion\dem2pts.py", line 88, in _dem2pts ncols = int(infile.ncols) File "netCDF4_netCDF4.pyx", line 2606, in netCDF4._netCDF4.Dataset.getattr File "netCDF4_netCDF4.pyx", line 2552, in netCDF4._netCDF4.Dataset.getncattr File "netCDF4_netCDF4.pyx", line 1188, in netCDF4._netCDF4._get_att File "netCDF4_netCDF4.pyx", line 1638, in netCDF4._netCDF4._ensure_nc_success AttributeError: NetCDF: Attribute not found

Any hints on how to proceed otherwise are highly welcome. Converting from tif->asc->dem I have not succeeded neither.

Best regards

Sebastian

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHubhttps://github.com/GeoscienceAustralia/anuga_core/issues/209?email_source=notifications&email_token=ABKOKJGARJ5NNCTNXRJ7CDTQKDMKBA5CNFSM4IXQBD72YY3PNVWWK3TUL52HS4DFUVEXG43VMWVGG33NNVSW45C7NFSM4HL3N2AQ, or mute the threadhttps://github.com/notifications/unsubscribe-auth/ABKOKJFOGOXFZUHQ4UXUU3DQKDMKBANCNFSM4IXQBD7Q.