HDF error when trying to write Dataset read with rasterio to NetCDF

[loicdtx]

I'm getting an HDF error when trying to write a Dataset read from GeoTiff (rasterio backend) to NetCDF. See reproducible example below:

import urllib.request
import tempfile
import os

import xarray as xr

path = tempfile.gettempdir()
url = 'https://earthexplorer.usgs.gov/browse/gisready/landsat_8/LC08_L1TP_026047_20180110_20180119_01_T1.zip'
filename = os.path.join(path, url.split('/')[-1])
nc_name = os.path.join(path, 'landsat_rgb.nc')

# Download file if not exist (11 Mb)
if not os.path.isfile(filename):
    urllib.request.urlretrieve(url, filename)

# Read rgb file using rasterio backend
rgb_name = '/'.join(['/vsizip', filename,
                     os.path.basename(filename).split('.')[-2] + '.tif'])
ds = xr.open_rasterio(rgb_name)
ds = ds.to_dataset('band').rename({1:'blue', 2:'green', 3:'red'})
print(ds)

# <xarray.Dataset>
# Dimensions:  (x: 7611, y: 7761)
# Coordinates:
#   * y        (y) float64 2.193e+06 2.193e+06 2.193e+06 ... 1.961e+06 1.960e+06
#   * x        (x) float64 3.732e+05 3.732e+05 3.733e+05 ... 6.015e+05 6.015e+05
# Data variables:
#     blue     (y, x) uint8 ...
#     red      (y, x) uint8 ...
#     green    (y, x) uint8 ...
# Attributes:
#     transform:   (30.0, 0.0, 373185.0, 0.0, -30.0, 2193315.0)
#     crs:         +init=epsg:32614
#     res:         (30.0, 30.0)
#     is_tiled:    1
#     nodatavals:  (nan, nan, nan)


# Write to netcdf
ds.to_netcdf(nc_name)

Output of xr.show_versions()

python -c "import xarray as xr; xr.show_versions()"

INSTALLED VERSIONS

commit: None
python: 3.6.5.final.0
python-bits: 64
OS: Linux
OS-release: 4.15.0-36-generic
machine: x86_64
processor: x86_64
byteorder: little
LC_ALL: None
LANG: en_US.UTF-8
LOCALE: en_US.UTF-8

xarray: 0.10.9
pandas: 0.23.4
numpy: 1.15.3
scipy: None
netCDF4: 1.4.2
h5netcdf: None
h5py: None
Nio: None
zarr: 2.2.0
cftime: 1.0.2.1
PseudonetCDF: None
rasterio: 1.0.9
iris: None
bottleneck: None
cyordereddict: None
dask: 0.20.0
distributed: None
matplotlib: None
cartopy: None
seaborn: None
setuptools: 40.5.0
pip: 18.1
conda: None
pytest: None
IPython: 7.1.1
sphinx: None

[jhamman]

@loicdtx - thanks for the report. I just tried this using xarray master and didn't get an error. So, you could try that and see if that works. Since our last release, we've had a fairly significant refactor of the IO backends in xarray.

[loicdtx]

Hi @jhamman, just tried what you suggest. Apparently it's a scipy vs netCDF4 thing; it works with scipy but not with netCDF4 on both master HEAD and the latest stable release.

[jhamman]

Interesting. I'm not seeing any difference when using scipy or netcdf4 or h5netcdf.

(by the way, thank you for the reproducible example)

[loicdtx]

yes, that's interesting... are you using a different OS? Let me know if there's something else I can help with.

[shoyer]

Exactly what error message do you see?

[loicdtx]

@shoyer, below the full traceback; I installed everything with pip

mktmpenv -p python3
pip install xarray numpy netcdf4
pip install rasterio

Traceback (most recent call last):
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/api.py", line 724, in to_netcdf
    unlimited_dims=unlimited_dims, compute=compute)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/core/dataset.py", line 1179, in dump_to_store
    unlimited_dims=unlimited_dims)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/common.py", line 374, in store
    unlimited_dims=unlimited_dims)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/netCDF4_.py", line 406, in set_variables
    super(NetCDF4DataStore, self).set_variables(*args, **kwargs)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/common.py", line 413, in set_variables
    self.writer.add(source, target)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/common.py", line 272, in add
    target[...] = source
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/netCDF4_.py", line 48, in __setitem__
    data[key] = value
  File "netCDF4/_netCDF4.pyx", line 4648, in netCDF4._netCDF4.Variable.__setitem__
  File "netCDF4/_netCDF4.pyx", line 4913, in netCDF4._netCDF4.Variable._put
  File "netCDF4/_netCDF4.pyx", line 1754, in netCDF4._netCDF4._ensure_nc_success
RuntimeError: NetCDF: HDF error

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "rasterio_reprex.py", line 41, in <module>
    ds.to_netcdf(nc_name)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/core/dataset.py", line 1254, in to_netcdf
    compute=compute)
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/api.py", line 729, in to_netcdf
    store.close()
  File "/home/loic/.virtualenvs/tmp-3caa6b25124e2f31/lib/python3.6/site-packages/xarray/backends/netCDF4_.py", line 474, in close
    ds.close()
  File "netCDF4/_netCDF4.pyx", line 2276, in netCDF4._netCDF4.Dataset.close
  File "netCDF4/_netCDF4.pyx", line 2260, in netCDF4._netCDF4.Dataset._close
  File "netCDF4/_netCDF4.pyx", line 1754, in netCDF4._netCDF4._ensure_nc_success
RuntimeError: NetCDF: HDF error

[jhamman]

Is this the same error you get with xarray/master?

[fmaussion]

I can reproduce this on my machine - also linux. This is going to be hard to track down though.

[loicdtx]

@jhamman, yes, same error message when installing from master

I have similar problem, when importing rasterio in the same script (not even using it for anything). This fails with HDF error:

import xarray as xa
import numpy as np
#import netCDF4
import rasterio

ds = xa.Dataset()
ds['z'] = (('y', 'x'), np.zeros((100, 100), np.float32))
print(ds)
ds.to_netcdf('test.nc')
ds.close()

with xa.open_dataset('test.nc') as ds:
    print(ds)

If I import netCDF4 before rasterio it works fine (uncomment line 3). This is probably an issue with rasterio somehow.

I installed everything with pip:

$ pip install Cython
$ pip install netCDF4 xarray rasterio numpy

From pip freeze:

affine==2.2.1
attrs==18.2.0
cftime==1.0.3
Click==7.0
click-plugins==1.0.4
cligj==0.5.0
Cython==0.29.1
netCDF4==1.4.2
numpy==1.15.4
pandas==0.23.4
pyparsing==2.3.0
python-dateutil==2.7.5
pytz==2018.7
rasterio==1.0.11
six==1.11.0
snuggs==1.4.2
xarray==0.11.0

It seems that this is not a problem with xarray but only with rasterio and netCDF4. Also this fails:

import rasterio
import netCDF4

with netCDF4.Dataset('test.nc', mode='w') as ds:
    ds.createDimension('x')
    ds.createVariable('foo', float, dimensions=('x'))
    print(ds)

Commenting out import rasterio removes the HDF error. I’ll report this to rasterio.

[shoyer]

This looks like a binary incompatibility issue with wheels for rasterio and netCDF4 on PyPI.

Good alternatives (for now) would be building from source or using conda-forge.