Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix bug in OLCI reader that caused multiple error messages to print #945

Merged
merged 8 commits into from
Oct 21, 2019
92 changes: 52 additions & 40 deletions satpy/readers/olci_nc.py
Original file line number Diff line number Diff line change
Expand Up @@ -15,7 +15,31 @@
#
# You should have received a copy of the GNU General Public License along with
# satpy. If not, see <http://www.gnu.org/licenses/>.
"""Sentinel-3 OLCI reader."""
"""Sentinel-3 OLCI reader.

This reader supports an optional argument to choose the 'engine' for reading
simonrp84 marked this conversation as resolved.
Show resolved Hide resolved
OLCI netCDF4 files. By default, this reader uses the default xarray choice of
engine, as defined in the `xarray open_dataset documentation`_.

As an alternative, the user may wish to use the 'h5netcdf' engine, but that is
not default as it typically prints many non-fatal but confusing error messages
to the terminal.
To choose between engines the user can do as follows for the default::

scn = satpyScene(filenames=my_files, reader='olci_l1b')

or as follows for the h5netcdf engine::

scn = Scene(filenames=my_files,
reader='olci_l1b'), reader_kwargs={'engine': 'h5netcdf'})


References:
- `xarray open_dataset documentation`_
.. _xarray open_dataset: http://xarray.pydata.org/en/stable/generated/xarray.open_dataset.html

"""


import logging
from datetime import datetime
Expand Down Expand Up @@ -73,14 +97,15 @@ def __getitem__(self, item):
class NCOLCIBase(BaseFileHandler):
"""The OLCI reader base."""

def __init__(self, filename, filename_info, filetype_info):
def __init__(self, filename, filename_info, filetype_info,
engine=None):
"""Init the olci reader base."""
super(NCOLCIBase, self).__init__(filename, filename_info,
filetype_info)
self.nc = xr.open_dataset(self.filename,
decode_cf=True,
mask_and_scale=True,
engine='h5netcdf',
engine=engine,
chunks={'columns': CHUNK_SIZE,
simonrp84 marked this conversation as resolved.
Show resolved Hide resolved
'rows': CHUNK_SIZE})

Expand Down Expand Up @@ -109,6 +134,13 @@ def get_dataset(self, key, info):

return variable

def __del__(self):
"""Close the NetCDF file that may still be open."""
try:
self.nc.close()
except (IOError, OSError, AttributeError):
pass


class NCOLCICal(NCOLCIBase):
"""Dummy class for calibration."""
Expand All @@ -125,7 +157,8 @@ class NCOLCIGeo(NCOLCIBase):
class NCOLCIChannelBase(NCOLCIBase):
"""Base class for channel reading."""

def __init__(self, filename, filename_info, filetype_info):
def __init__(self, filename, filename_info, filetype_info,
engine=None):
"""Init the file handler."""
super(NCOLCIChannelBase, self).__init__(filename, filename_info,
filetype_info)
Expand All @@ -136,44 +169,13 @@ def __init__(self, filename, filename_info, filetype_info):
class NCOLCI1B(NCOLCIChannelBase):
"""File handler for OLCI l1b."""

def __init__(self, filename, filename_info, filetype_info, cal):
def __init__(self, filename, filename_info, filetype_info, cal,
engine=None):
"""Init the file handler."""
super(NCOLCI1B, self).__init__(filename, filename_info,
filetype_info)
self.cal = cal.nc

def _get_solar_flux_old(self, band):
"""Get the solar flux."""
# TODO: this could be replaced with vectorized indexing in the future.
from dask.base import tokenize
blocksize = CHUNK_SIZE

solar_flux = self.cal['solar_flux'].isel(bands=band).values
d_index = self.cal['detector_index'].fillna(0).astype(int)

shape = d_index.shape
vchunks = range(0, shape[0], blocksize)
hchunks = range(0, shape[1], blocksize)

token = tokenize(band, d_index, solar_flux)
name = 'solar_flux_' + token

def get_items(array, slices):
return solar_flux[d_index[slices].values]

dsk = {(name, i, j): (get_items,
d_index,
(slice(vcs, min(vcs + blocksize, shape[0])),
slice(hcs, min(hcs + blocksize, shape[1]))))
for i, vcs in enumerate(vchunks)
for j, hcs in enumerate(hchunks)
}

res = da.Array(dsk, name, shape=shape,
chunks=(blocksize, blocksize),
dtype=solar_flux.dtype)
return res

@staticmethod
def _get_items(idx, solar_flux):
"""Get items."""
Expand All @@ -184,7 +186,8 @@ def _get_solar_flux(self, band):
solar_flux = self.cal['solar_flux'].isel(bands=band).values
d_index = self.cal['detector_index'].fillna(0).astype(int)

return da.map_blocks(self._get_items, d_index.data, solar_flux=solar_flux, dtype=solar_flux.dtype)
return da.map_blocks(self._get_items, d_index.data,
solar_flux=solar_flux, dtype=solar_flux.dtype)

def get_dataset(self, key, info):
"""Load a dataset."""
Expand Down Expand Up @@ -247,7 +250,8 @@ class NCOLCIAngles(BaseFileHandler):
'solar_azimuth_angle': 'SAA',
'solar_zenith_angle': 'SZA'}

def __init__(self, filename, filename_info, filetype_info):
def __init__(self, filename, filename_info, filetype_info,
engine=None):
"""Init the file handler."""
super(NCOLCIAngles, self).__init__(filename, filename_info,
filetype_info)
Expand All @@ -258,6 +262,7 @@ def __init__(self, filename, filename_info, filetype_info):
self.cache = {}
self._start_time = filename_info['start_time']
self._end_time = filename_info['end_time']
self.engine = engine

def get_dataset(self, key, info):
"""Load a dataset."""
Expand All @@ -268,7 +273,7 @@ def get_dataset(self, key, info):
self.nc = xr.open_dataset(self.filename,
decode_cf=True,
mask_and_scale=True,
engine='h5netcdf',
engine=self.engine,
chunks={'tie_columns': CHUNK_SIZE,
'tie_rows': CHUNK_SIZE})

Expand Down Expand Up @@ -356,3 +361,10 @@ def start_time(self):
def end_time(self):
"""End the file handler."""
return self._end_time

def __del__(self):
"""Close the NetCDF file that may still be open."""
try:
self.nc.close()
except (IOError, OSError, AttributeError):
pass
57 changes: 48 additions & 9 deletions satpy/tests/reader_tests/test_olci_nc.py
Original file line number Diff line number Diff line change
Expand Up @@ -38,9 +38,19 @@ def test_instantiate(self, mocked_dataset):
from satpy.readers.olci_nc import (NCOLCIBase, NCOLCICal, NCOLCIGeo,
NCOLCIChannelBase, NCOLCI1B, NCOLCI2)
from satpy import DatasetID
import xarray as xr

cal_data = xr.Dataset(
{
'solar_flux': (('bands'), [0, 1, 2]),
'detector_index': (('bands'), [0, 1, 2]),
},
{'bands': [0, 1, 2], },
)

ds_id = DatasetID(name='foo')
filename_info = {'mission_id': 'S3A', 'dataset_name': 'foo', 'start_time': 0, 'end_time': 0}
ds_id = DatasetID(name='Oa01', calibration='reflectance')
ds_id2 = DatasetID(name='wsqf', calibration='reflectance')
filename_info = {'mission_id': 'S3A', 'dataset_name': 'Oa01', 'start_time': 0, 'end_time': 0}

test = NCOLCIBase('somedir/somefile.nc', filename_info, 'c')
test.get_dataset(ds_id, filename_info)
Expand All @@ -62,22 +72,19 @@ def test_instantiate(self, mocked_dataset):
mocked_dataset.assert_called()
mocked_dataset.reset_mock()

test = NCOLCI1B('somedir/somefile.nc', filename_info, 'c', mock.Mock())
cal = mock.Mock()
cal.nc = cal_data
test = NCOLCI1B('somedir/somefile.nc', filename_info, 'c', cal)
test.get_dataset(ds_id, filename_info)
mocked_dataset.assert_called()
mocked_dataset.reset_mock()

test = NCOLCI2('somedir/somefile.nc', filename_info, 'c')
test.get_dataset(ds_id, {'nc_key': 'the_key'})
test.get_dataset(ds_id2, {'nc_key': 'the_key'})
mocked_dataset.assert_called()
mocked_dataset.reset_mock()

# ds_id = DatasetID(name='solar_azimuth_angle')
# test = NCOLCIAngles('somedir/somefile.nc', filename_info, 'c')
# test.get_dataset(ds_id, filename_info)
# mocked_dataset.assert_called()
# mocked_dataset.reset_mock()

@mock.patch('xarray.open_dataset')
def test_get_dataset(self, mocked_dataset):
"""Test reading datasets."""
Expand All @@ -95,6 +102,38 @@ def test_get_dataset(self, mocked_dataset):
res = test.get_dataset(ds_id, {'nc_key': 'mask'})
self.assertEqual(res.dtype, np.dtype('bool'))

@mock.patch('xarray.open_dataset')
def test_olci_angles(self, mocked_dataset):
"""Test reading datasets."""
from satpy.readers.olci_nc import NCOLCIAngles
from satpy import DatasetID
import numpy as np
import xarray as xr
attr_dict = {
'ac_subsampling_factor': 1,
'al_subsampling_factor': 2,
}
mocked_dataset.return_value = xr.Dataset({'SAA': (['tie_rows', 'tie_columns'],
np.array([1 << x for x in range(30)]).reshape(5, 6)),
'SZA': (['tie_rows', 'tie_columns'],
np.array([1 << x for x in range(30)]).reshape(5, 6)),
'OAA': (['tie_rows', 'tie_columns'],
np.array([1 << x for x in range(30)]).reshape(5, 6)),
'OZA': (['tie_rows', 'tie_columns'],
np.array([1 << x for x in range(30)]).reshape(5, 6))},
coords={'rows': np.arange(5),
'columns': np.arange(6)},
attrs=attr_dict)
filename_info = {'mission_id': 'S3A', 'dataset_name': 'Oa01', 'start_time': 0, 'end_time': 0}

ds_id = DatasetID(name='solar_azimuth_angle')
ds_id2 = DatasetID(name='satellite_zenith_angle')
test = NCOLCIAngles('somedir/somefile.nc', filename_info, 'c')
test.get_dataset(ds_id, filename_info)
test.get_dataset(ds_id2, filename_info)
mocked_dataset.assert_called()
mocked_dataset.reset_mock()


class TestBitFlags(unittest.TestCase):
"""Test the bitflag reading."""
Expand Down