open_mfdataset too many files

[rabernat]

I am very excited to try xray.

On my first attempt, I tried to use open_mfdataset on a set of ~8000 netcdf files. I hit a "RuntimeError: Too many open files". The ulimit on my system is 1024, so clearly that is the source of the error.

I am curious whether this is the desired behavior for open_mfdataset. Does xray have to keep all the files open? If so, I will work with my sysadmin to increase the ulimit.

It seems like the whole point of this function is to work with large collections of files, so this could be a significant limitation.

[rabernat]

Just a little follow up...I tried to work around the file limit by serializing the processing of the files and creating xray datasets with with fewer files in them. However, I still eventually hit this error, suggesting that the files are never being closed. For example

I would like to do

ds = xray.open_mfdataset(ddir + '*.nc' % yr, engine='scipy')
EKE = (ds.variables['u']**2 + ds.variables['v']**2).mean(dim='time').load()

This tries to open 8031 files and produces the error: [Errno 24] Too many open files

So then I try to create a new dataset for each year

EKE = []
for yr in xrange(1993,2015):
    print yr
    # this opens about 365 files
    ds = xray.open_mfdataset(ddir + '/dt_global_allsat_msla_uv_%04d*.nc' % yr, engine='scipy')
    EKE.append((ds.variables['u']**2 + ds.variables['v']**2).mean(dim='time').load())

This works okay for the first two years. However, by the third year, I still get the error: [Errno 24] Too many open files. This is when the ulimit of 1024 files is exceeded.

Using xray version 0.5.1 via conda module.

[shoyer]

Yes, this is a known issue, and I agree that it is annoying. We could work around this by opening up (and closing) netCDF files inside the __getitem__ call. If you're interested in possibly working on this, take a look at the netCDF4 backend for xray: https://github.com/xray/xray/blob/master/xray/backends/netCDF4_.py

[rabernat]

I am using the scipy backend because the netcdf4 backend doesn't work for me at all. It core dumps with the error

python: posixio.c:366: px_rel: Assertion `pxp->bf_offset <= offset && offset < pxp->bf_offset + (off_t) pxp->bf_extent' failed.
Aborted (core dumped)

Are you suggesting I work on the scipy backend?

[shoyer]

Sure, you could do this on the scipy backend -- the logic will be essentially the same on both backends.

I believe your issue with netCDF4 backend is the same as this one: #444. This will be fixed in the next release.

[rabernat]

Ok, I will have a look at this. I would be happy to contribute to this awesome project.

By the way, by monitoring /proc, I was able to see that the scipy backend actually opens each file TWICE, exacerbating the problem.

[rabernat]

I came up with a solution for this, but it is so slow that it is useless.

[shoyer]

Hmm. How big are each of your netCDF files?

[rabernat]

8 MB. This is daily satellite data, with one file per time point. (Most satellite data is distributed this way.)

There are many other workarounds to this problem. You can try to increase your ulimits. Or you can join these small netcdf files together into a big one. I had daily data files, and I used NCO to concatentate them into monthly files. That basically solved my problem. But of course that involves going out of xray.

[cpaulik]

I've run into the same problem and have been looking at the netCDF backend. A solution does not seem to be so easy as to open and close the file in the __getitem__ method since this closes the file also for any other access e.g. attributes like shape or dtype.

Short of decorating all the functions of the netCDF4 package I can not think of a workable solution to this. But maybe I'm overlooking something fundamental.

[shoyer]

I think we can actually read all the variable metadata (shape and dtype) in when we open the file -- we already do that for reading in attributes. Something like this prototype, which would also be useful for reading compressed netCDF4 files with multiprocessing: dask/dask#457 (comment)

[cpaulik]

I've pushed a few commits trying this out to https://github.com/cpaulik/xray/tree/closing_netcdf_backend . I can open a WIP PR if this would be easier to discuss there.

There are however a few tests that keep failing and I can not figure out why.

e.g.: test_backends.py::NetCDF4ViaDaskDataTest::test_compression_encoding:

If I set a breakpoint at line 941 of dataset.py and just continue the test fails.

If I however evaluate self.variables.items() or even self.variables at the breakpoint I get the correct output and the test passes when continued. I can not really see the difference between me evaluating this in ipdb and the code that is on the line.

The error I get when running the test without interference is:

test_backends.py::NetCDF4ViaDaskDataTest::test_compression_encoding FAILED

====================================================== FAILURES =======================================================
__________________________________ NetCDF4ViaDaskDataTest.test_compression_encoding ___________________________________

self = <xray.test.test_backends.NetCDF4ViaDaskDataTest testMethod=test_compression_encoding>

    def test_compression_encoding(self):
        data = create_test_data()
        data['var2'].encoding.update({'zlib': True,
                                      'chunksizes': (5, 5),
                                      'fletcher32': True})
>       with self.roundtrip(data) as actual:

test_backends.py:502: 
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
/usr/lib/python2.7/contextlib.py:17: in __enter__
    return self.gen.next()
test_backends.py:596: in roundtrip
    yield ds.chunk()
../core/dataset.py:942: in chunk
    for k, v in self.variables.items()])
../core/dataset.py:935: in maybe_chunk
    token2 = tokenize(name, token if token else var._data)
/home/cpa/.virtualenvs/xray/local/lib/python2.7/site-packages/dask/base.py:152: in tokenize
    return md5(str(tuple(map(normalize_token, args))).encode()).hexdigest()
../core/indexing.py:301: in __repr__
    (type(self).__name__, self.array, self.key))
../core/utils.py:377: in __repr__
    return '%s(array=%r)' % (type(self).__name__, self.array)
../core/indexing.py:301: in __repr__
    (type(self).__name__, self.array, self.key))
../core/utils.py:377: in __repr__
    return '%s(array=%r)' % (type(self).__name__, self.array)
netCDF4/_netCDF4.pyx:2931: in netCDF4._netCDF4.Variable.__repr__ (netCDF4/_netCDF4.c:25068)
    ???
netCDF4/_netCDF4.pyx:2938: in netCDF4._netCDF4.Variable.__unicode__ (netCDF4/_netCDF4.c:25243)
    ???
netCDF4/_netCDF4.pyx:3059: in netCDF4._netCDF4.Variable.dimensions.__get__ (netCDF4/_netCDF4.c:27486)
    ???
_ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _

>   ???
E   RuntimeError: NetCDF: Not a valid ID

netCDF4/_netCDF4.pyx:2994: RuntimeError
============================================== 1 failed in 0.50 seconds ===============================================

[shoyer]

@cpaulik I wonder if the issue is this section in your __getitem__ method:

        data = getitem(self.array, key)
        try:
            self.store.ensure_open()
            data = getitem(self.array, key)
        except RuntimeError as e:
            raise e
            pass
        if self.ndim == 0:
            # work around for netCDF4-python's broken handling of 0-d
            # arrays (slicing them always returns a 1-dimensional array):
            # https://github.com/Unidata/netcdf4-python/pull/220
            data = np.asscalar(data)
        self.store.close()
        return data

I would put self.store.close() in a finally clause following the getitem clause.

Actually, you probably want to put this in a context manager that automatically closes the file, something like:

with self.store.opened():
    data = getitem(self.array, key)

[cpaulik]

I've only put the try - except there to conditionally set the breakpoint. How does it make a difference if the self.store.close is called? It it is not called then the dataset remains opened which should not cause the weird behaviour reported above?

Nevertheless I have updated my branch to use a contextmanager because it is a better solution but I still have this strange behaviour of only printing the variable altering the test outcome.

[shoyer]

OK, so the problem is that self.array on NetCDF4ArrayWrapper is retaining a reference to netCDF4.Variable object on the closed dataset. It's not enough to merely ensure that a netCDF4 dataset is opened -- you also need to ensure that no references to variables on the old dataset are still around. So get_variables/open_store_variable may need a refactor to deal with this.

[cpaulik]

OK, I'll try. Thanks.

But I originally tested if netCDF4 can work with a closed/reopened variable like this:

In [1]: import netCDF4

In [2]: a = netCDF4.Dataset("temp.nc", mode="w")

In [3]: a.createDimension("lon")
Out[3]: <class 'netCDF4._netCDF4.Dimension'> (unlimited): name = 'lon', size = 0

In [4]: a.createVariable("lon", "f8", dimensions=("lon"))
Out[4]:
    <class 'netCDF4._netCDF4.Variable'>
    float64 lon(lon)
    unlimited dimensions: lon
    current shape = (0,)
    filling on, default _FillValue of 9.969209968386869e+36 used

In [5]: v = a.variables['lon']

In [6]: v
Out[6]:
    <class 'netCDF4._netCDF4.Variable'>
    float64 lon(lon)
    unlimited dimensions: lon
    current shape = (0,)
    filling on, default _FillValue of 9.969209968386869e+36 used

In [7]: a.close()

In [8]: v
Out[8]: ---------------------------------------------------------------------------
    RuntimeError                              Traceback (most recent call last)
    /home/cp/.pyenv/versions/miniconda3-3.16.0/envs/xray-3.5.0/lib/python3.5/site-packages/IPython/core/formatters.py in __call__(self, obj)
    695                 type_pprinters=self.type_printers,
    696                 deferred_pprinters=self.deferred_printers)
    --> 697             printer.pretty(obj)
    698             printer.flush()
    699             return stream.getvalue()
    /home/cp/.pyenv/versions/miniconda3-3.16.0/envs/xray-3.5.0/lib/python3.5/site-packages/IPython/lib/pretty.py in pretty(self, obj)
    381                             if callable(meth):
    382                                 return meth(obj, self, cycle)
    --> 383             return _default_pprint(obj, self, cycle)
    384         finally:
    385             self.end_group()
    /home/cp/.pyenv/versions/miniconda3-3.16.0/envs/xray-3.5.0/lib/python3.5/site-packages/IPython/lib/pretty.py in _default_pprint(obj, p, cycle)
    501     if _safe_getattr(klass, '__repr__', None) not in _baseclass_reprs:
    502         # A user-provided repr. Find newlines and replace them with p.break_()
    --> 503         _repr_pprint(obj, p, cycle)
    504         return
    505     p.begin_group(1, '<')
    /home/cp/.pyenv/versions/miniconda3-3.16.0/envs/xray-3.5.0/lib/python3.5/site-packages/IPython/lib/pretty.py in _repr_pprint(obj, p, cycle)
    683     """A pprint that just redirects to the normal repr function."""
    684     # Find newlines and replace them with p.break_()
    --> 685     output = repr(obj)
    686     for idx,output_line in enumerate(output.splitlines()):
    687         if idx:
    netCDF4/_netCDF4.pyx in netCDF4._netCDF4.Variable.__repr__ (netCDF4/_netCDF4.c:25045)()
    netCDF4/_netCDF4.pyx in netCDF4._netCDF4.Variable.__unicode__ (netCDF4/_netCDF4.c:25243)()
    netCDF4/_netCDF4.pyx in netCDF4._netCDF4.Variable.dimensions.__get__ (netCDF4/_netCDF4.c:27486)()
    netCDF4/_netCDF4.pyx in netCDF4._netCDF4.Variable._getdims (netCDF4/_netCDF4.c:26297)()
    RuntimeError: NetCDF: Not a valid ID
In [9]: a = netCDF4.Dataset("temp.nc")
In [10]: v
Out[10]:
    class 'netCDF4._netCDF4.Variable'>
    lon(lon)
    dimensions: lon
    shape = (0,)
    on, default _FillValue of 9.969209968386869e+36 used