diff --git a/doc/source/api/satpy.readers.rst b/doc/source/api/satpy.readers.rst
index ad79cd7dd4..06810f740b 100644
--- a/doc/source/api/satpy.readers.rst
+++ b/doc/source/api/satpy.readers.rst
@@ -364,6 +364,14 @@ satpy.readers.msi\_safe module
:undoc-members:
:show-inheritance:
+satpy.readers.mviri\_l1b\_fiduceo\_nc module
+--------------------------------------------
+
+.. automodule:: satpy.readers.mviri_l1b_fiduceo_nc
+ :members:
+ :undoc-members:
+ :show-inheritance:
+
satpy.readers.netcdf\_utils module
----------------------------------
diff --git a/doc/source/index.rst b/doc/source/index.rst
index 7f96b93106..3ec9dafdc1 100644
--- a/doc/source/index.rst
+++ b/doc/source/index.rst
@@ -93,6 +93,9 @@ the base Satpy installation.
* - MSG (Meteosat 8 to 11) L2 products in GRIB2 format
- `seviri_l2_grib`
- In development, CLM, OCA and FIR products supported
+ * - MFG (Meteosat 2 to 7) MVIRI data in netCDF format (FIDUCEO FCDR)
+ - `mviri_l1b_fiduceo_nc`
+ - Beta
* - Himawari 8 and 9 AHI data in HSD format
- `ahi_hsd`
- Nominal
diff --git a/satpy/etc/readers/mviri_l1b_fiduceo_nc.yaml b/satpy/etc/readers/mviri_l1b_fiduceo_nc.yaml
new file mode 100644
index 0000000000..bf17a427ee
--- /dev/null
+++ b/satpy/etc/readers/mviri_l1b_fiduceo_nc.yaml
@@ -0,0 +1,135 @@
+# References:
+# - MFG User Handbook
+# - FIDUCEO MVIRI FCDR Product User Guide
+
+reader:
+ name: mviri_l1b_fiduceo_nc
+ short_name: FIDUCEO MVIRI FCDR
+ long_name: >
+ Fundamental Climate Data Record of re-calibrated Level 1.5 Infrared, Water Vapour, and Visible radiances
+ from the Meteosat Visible Infra-Red Imager (MVIRI) instrument onboard the Meteosat First Generation satellites
+ description: >
+ Reader for FIDUCEO MVIRI FCDR data in netCDF format. For documentation see
+ http://doi.org/10.15770/EUM_SEC_CLM_0009 .
+ sensors: [mviri]
+ default_channels: [VIS, WV, IR]
+ reader: !!python/name:satpy.readers.yaml_reader.GEOFlippableFileYAMLReader
+
+file_types:
+ nc_easy:
+ file_reader: !!python/name:satpy.readers.mviri_l1b_fiduceo_nc.FiduceoMviriEasyFcdrFileHandler
+ file_patterns: [
+ 'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_EASY_{processor_version}_{format_version}.nc'
+ # Example: FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_EASY_v2.6_fv3.1.nc
+ ]
+ nc_full:
+ file_reader: !!python/name:satpy.readers.mviri_l1b_fiduceo_nc.FiduceoMviriFullFcdrFileHandler
+ file_patterns: [
+ 'FIDUCEO_FCDR_{level}_{sensor}_{platform}-{projection_longitude:f}_{start_time:%Y%m%d%H%M}_{end_time:%Y%m%d%H%M}_FULL_{processor_version}_{format_version}.nc'
+ # Example: FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_FULL_v2.6_fv3.1.nc
+ ]
+
+datasets:
+ VIS:
+ name: VIS
+ resolution: 2250
+ wavelength: [0.5, 0.7, 0.9]
+ calibration:
+ reflectance:
+ standard_name: toa_bidirectional_reflectance
+ units: "%"
+ radiance:
+ # Confirmed by EUM: No (1/wavenumber) here. Hence no standard name.
+ units: W m-2 sr-1
+ counts:
+ standard_name: counts
+ units: count
+ file_type: [nc_easy, nc_full]
+
+ WV:
+ name: WV
+ resolution: 4500
+ wavelength: [5.7, 6.4, 7.1]
+ calibration:
+ brightness_temperature:
+ standard_name: toa_brightness_temperature
+ units: K
+ radiance:
+ standard_name: toa_outgoing_radiance_per_unit_wavenumber
+ units: mW m-2 sr-1 (cm-1)-1
+ counts:
+ standard_name: counts
+ units: count
+ file_type: [nc_easy, nc_full]
+
+ IR:
+ name: IR
+ resolution: 4500
+ wavelength: [10.5, 11.5, 12.5]
+ calibration:
+ brightness_temperature:
+ standard_name: toa_brightness_temperature
+ units: K
+ radiance:
+ standard_name: toa_outgoing_radiance_per_unit_wavenumber
+ units: mW m-2 sr-1 (cm-1)-1
+ counts:
+ standard_name: counts
+ units: count
+ file_type: [nc_easy, nc_full]
+
+ quality_pixel_bitmask:
+ name: quality_pixel_bitmask
+ resolution: 2250
+ file_type: [nc_easy, nc_full]
+
+ data_quality_bitmask:
+ name: data_quality_bitmask
+ resolution: 2250
+ file_type: [nc_easy, nc_full]
+
+ u_independent_toa_bidirectional_reflectance:
+ name: u_independent_toa_bidirectional_reflectance
+ long_name: "independent uncertainty per pixel"
+ units: "%"
+ resolution: 2250
+ file_type: [nc_easy]
+
+ u_structured_toa_bidirectional_reflectance:
+ name: u_structured_toa_bidirectional_reflectance
+ long_name: "structured uncertainty per pixel"
+ units: "%"
+ resolution: 2250
+ file_type: [nc_easy]
+
+ solar_zenith_angle:
+ name: solar_zenith_angle
+ standard_name: solar_zenith_angle
+ long_name: "Solar zenith angle"
+ units: degree
+ resolution: [2250, 4500]
+ file_type: [nc_easy, nc_full]
+
+ solar_azimuth_angle:
+ name: solar_azimuth_angle
+ standard_name: solar_azimuth_angle
+ long_name: "Solar azimuth angle"
+ units: degree
+ resolution: [2250, 4500]
+ file_type: [nc_easy, nc_full]
+
+ satellite_zenith_angle:
+ name: satellite_zenith_angle
+ standard_name: sensor_zenith_angle
+ long_name: "Satellite zenith angle"
+ units: degree
+ resolution: [2250, 4500]
+ file_type: [nc_easy, nc_full]
+
+ satellite_azimuth_angle:
+ name: satellite_azimuth_angle
+ standard_name: sensor_azimuth_angle
+ long_name: "Satellite azimuth angle"
+ units: degree
+ resolution: [2250, 4500]
+ file_type: [nc_easy, nc_full]
diff --git a/satpy/readers/_geos_area.py b/satpy/readers/_geos_area.py
index b340f5fca0..51b39dd3c8 100644
--- a/satpy/readers/_geos_area.py
+++ b/satpy/readers/_geos_area.py
@@ -145,3 +145,24 @@ def get_area_definition(pdict, a_ext):
a_ext)
return a_def
+
+
+def sampling_to_lfac_cfac(sampling):
+ """Convert angular sampling to line/column scaling factor (aka LFAC/CFAC).
+
+ Reference: `MSG Ground Segment LRIT HRIT Mission Specific Implementation`_,
+ Appendix E.2.
+
+
+ .. _MSG Ground Segment LRIT HRIT Mission Specific Implementation: http://www.eumetsat.int/\
+website/wcm/idc/idcplg?IdcService=GET_FILE&dDocName=PDF_TEN_05057_SPE_MSG_LRIT_HRI\
+&RevisionSelectionMethod=LatestReleased&Rendition=Web
+
+ Args:
+ sampling: float
+ Angular sampling (rad)
+
+ Returns:
+ Line/column scaling factor (deg-1)
+ """
+ return 2.0 ** 16 / np.rad2deg(sampling)
diff --git a/satpy/readers/mviri_l1b_fiduceo_nc.py b/satpy/readers/mviri_l1b_fiduceo_nc.py
new file mode 100644
index 0000000000..c4fbcaf717
--- /dev/null
+++ b/satpy/readers/mviri_l1b_fiduceo_nc.py
@@ -0,0 +1,792 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+# Copyright (c) 2020 Satpy developers
+#
+# This file is part of satpy.
+#
+# satpy is free software: you can redistribute it and/or modify it under the
+# terms of the GNU General Public License as published by the Free Software
+# Foundation, either version 3 of the License, or (at your option) any later
+# version.
+#
+# satpy is distributed in the hope that it will be useful, but WITHOUT ANY
+# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
+# A PARTICULAR PURPOSE. See the GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License along with
+# satpy. If not, see .
+"""FIDUCEO MVIRI FCDR Reader.
+
+Introduction
+------------
+The FIDUCEO MVIRI FCDR is a Fundamental Climate Data Record (FCDR) of
+re-calibrated Level 1.5 Infrared, Water Vapour, and Visible radiances from
+the Meteosat Visible Infra-Red Imager (MVIRI) instrument onboard the
+Meteosat First Generation satellites. There are two variants of the dataset:
+The *full FCDR* and a simplified version called *easy FCDR*. Some datasets are
+only available in one of the two variants, see the corresponding YAML
+definition in ``satpy/etc/readers/``.
+
+Dataset Names
+-------------
+The FIDUCEO MVIRI readers use names ``VIS``, ``WV`` and ``IR`` for the visible,
+water vapor and infrared channels, respectively. These are different from
+the original netCDF variable names for the following reasons:
+
+- VIS channel is named differently in full FCDR (``counts_vis``) and easy FCDR
+ (``toa_bidirectional_reflectance_vis``)
+- netCDF variable names contain the calibration level (e.g. ``counts_...``),
+ which might be confusing for satpy users if a different calibration level
+ is chosen.
+
+Remaining datasets (such as quality flags and uncertainties) have the same
+name in the reader as in the netCDF file.
+
+
+Example
+-------
+This is how to read FIDUCEO MVIRI FCDR data in satpy:
+
+.. code-block:: python
+
+ from satpy import Scene
+
+ scn = Scene(filenames=['FIDUCEO_FCDR_L15_MVIRI_MET7-57.0...'],
+ reader='mviri_l1b_fiduceo_nc')
+ scn.load(['VIS', 'WV', 'IR'])
+
+Global netCDF attributes are available in the ``raw_metadata`` attribute of
+each loaded dataset.
+
+
+Image Orientation
+-----------------
+The images are stored in MVIRI scanning direction, that means South is up and
+East is right. This can be changed as follows:
+
+.. code-block:: python
+
+ scn.load(['VIS'], upper_right_corner='NE')
+
+
+Geolocation
+-----------
+In addition to the image data, FIDUCEO also provides so called *static FCDRs*
+containing latitude and longitude coordinates. In order to simplify their
+usage, the FIDUCEO MVIRI readers do not make use of these static files, but
+instead provide an area definition that can be used to compute longitude and
+latitude coordinates on demand.
+
+.. code-block:: python
+
+ area = scn['VIS'].attrs['area']
+ lons, lats = area.get_lonlats()
+
+Those were compared to the static FCDR and they agree very well, however there
+are small differences. The mean difference is < 1E3 degrees for all channels
+and projection longitudes.
+
+
+Huge VIS Reflectances
+---------------------
+You might encounter huge VIS reflectances (10^8 percent and greater) in
+situations where both radiance and solar zenith angle are small. The reader
+certainly needs some improvement in this regard. Maybe the corresponding
+uncertainties can be used to filter these cases before calculating reflectances.
+
+
+VIS Channel Quality Flags
+-------------------------
+Quality flags are available for the VIS channel only. A simple approach for
+masking bad quality pixels is to set the ``mask_bad_quality`` keyword argument
+to ``True``:
+
+.. code-block:: python
+
+ scn = Scene(filenames=['FIDUCEO_FCDR_L15_MVIRI_MET7-57.0...'],
+ reader='mviri_l1b_fiduceo_nc',
+ reader_kwargs={'mask_bad_quality': True})
+
+See :class:`FiduceoMviriBase` for an argument description. In some situations
+however the entire image can be flagged (look out for warnings). In that case
+check out the ``quality_pixel_bitmask`` and ``data_quality_bitmask`` datasets
+to find out why.
+
+
+Angles
+------
+The FIDUCEO MVIRI FCDR provides satellite and solar angles on a coarse tiepoint
+grid. By default these datasets will be interpolated to the higher VIS
+resolution. This can be changed as follows:
+
+.. code-block:: python
+
+ scn.load(['solar_zenith_angle'], resolution=4500)
+
+If you need the angles in both resolutions, use data queries:
+
+.. code-block:: python
+
+ from satpy import DataQuery
+
+ query_vis = DataQuery(
+ name='solar_zenith_angle',
+ resolution=2250
+ )
+ query_ir = DataQuery(
+ name='solar_zenith_angle',
+ resolution=4500
+ )
+ scn.load([query_vis, query_ir])
+
+ # Use the query objects to access the datasets as follows
+ sza_vis = scn[query_vis]
+
+
+References
+----------
+ - `[Handbook]`_ MFG User Handbook
+ - `[PUG]`_ FIDUCEO MVIRI FCDR Product User Guide
+
+.. _[Handbook]: http://www.eumetsat.int/\
+website/wcm/idc/idcplg?IdcService=GET_FILE&dDocName=PDF_TD06_MARF&\
+RevisionSelectionMethod=LatestReleased&Rendition=Web
+.. _[PUG]: http://doi.org/10.15770/EUM_SEC_CLM_0009
+"""
+
+import abc
+from functools import lru_cache
+import warnings
+
+import dask.array as da
+import numpy as np
+import xarray as xr
+
+from satpy import CHUNK_SIZE
+from satpy.readers._geos_area import (
+ sampling_to_lfac_cfac,
+ get_area_definition,
+ get_area_extent
+)
+from satpy.readers.file_handlers import BaseFileHandler
+
+EQUATOR_RADIUS = 6378140.0
+POLE_RADIUS = 6356755.0
+ALTITUDE = 42164000.0 - EQUATOR_RADIUS
+"""[Handbook] section 5.2.1."""
+
+MVIRI_FIELD_OF_VIEW = 18.0
+"""[Handbook] section 5.3.2.1."""
+
+CHANNELS = ['VIS', 'WV', 'IR']
+ANGLES = [
+ 'solar_zenith_angle',
+ 'solar_azimuth_angle',
+ 'satellite_zenith_angle',
+ 'satellite_azimuth_angle'
+]
+OTHER_REFLECTANCES = [
+ 'u_independent_toa_bidirectional_reflectance',
+ 'u_structured_toa_bidirectional_reflectance'
+]
+HIGH_RESOL = 2250
+
+
+class IRWVCalibrator:
+ """Calibrate IR & WV channels."""
+
+ def __init__(self, coefs):
+ """Initialize the calibrator.
+
+ Args:
+ coefs: Calibration coefficients.
+ """
+ self.coefs = coefs
+
+ def calibrate(self, counts, calibration):
+ """Calibrate IR/WV counts to the given calibration."""
+ if calibration == 'counts':
+ return counts
+ elif calibration in ('radiance', 'brightness_temperature'):
+ return self._calibrate_rad_bt(counts, calibration)
+ else:
+ raise KeyError(
+ 'Invalid IR/WV calibration: {}'.format(calibration.name)
+ )
+
+ def _calibrate_rad_bt(self, counts, calibration):
+ """Calibrate counts to radiance or brightness temperature."""
+ rad = self._counts_to_radiance(counts)
+ if calibration == 'radiance':
+ return rad
+ bt = self._radiance_to_brightness_temperature(rad)
+ return bt
+
+ def _counts_to_radiance(self, counts):
+ """Convert IR/WV counts to radiance.
+
+ Reference: [PUG], equations (4.1) and (4.2).
+ """
+ rad = self.coefs['a'] + self.coefs['b'] * counts
+ return rad.where(rad > 0, np.float32(np.nan))
+
+ def _radiance_to_brightness_temperature(self, rad):
+ """Convert IR/WV radiance to brightness temperature.
+
+ Reference: [PUG], equations (5.1) and (5.2).
+ """
+ bt = self.coefs['bt_b'] / (np.log(rad) - self.coefs['bt_a'])
+ return bt.where(bt > 0, np.float32(np.nan))
+
+
+class VISCalibrator:
+ """Calibrate VIS channel."""
+
+ def __init__(self, coefs, solar_zenith_angle=None):
+ """Initialize the calibrator.
+
+ Args:
+ coefs:
+ Calibration coefficients.
+ solar_zenith_angle (optional):
+ Solar zenith angle. Only required for calibration to
+ reflectance.
+ """
+ self.coefs = coefs
+ self.solar_zenith_angle = solar_zenith_angle
+
+ def calibrate(self, counts, calibration):
+ """Calibrate VIS counts."""
+ if calibration == 'counts':
+ return counts
+ elif calibration in ('radiance', 'reflectance'):
+ return self._calibrate_rad_refl(counts, calibration)
+ else:
+ raise KeyError(
+ 'Invalid VIS calibration: {}'.format(calibration.name)
+ )
+
+ def _calibrate_rad_refl(self, counts, calibration):
+ """Calibrate counts to radiance or reflectance."""
+ rad = self._counts_to_radiance(counts)
+ if calibration == 'radiance':
+ return rad
+ refl = self._radiance_to_reflectance(rad)
+ refl = self.update_refl_attrs(refl)
+ return refl
+
+ def _counts_to_radiance(self, counts):
+ """Convert VIS counts to radiance.
+
+ Reference: [PUG], equations (7) and (8).
+ """
+ years_since_launch = self.coefs['years_since_launch']
+ a_cf = (self.coefs['a0'] +
+ self.coefs['a1'] * years_since_launch +
+ self.coefs['a2'] * years_since_launch ** 2)
+ mean_count_space_vis = self.coefs['mean_count_space']
+ rad = (counts - mean_count_space_vis) * a_cf
+ return rad.where(rad > 0, np.float32(np.nan))
+
+ def _radiance_to_reflectance(self, rad):
+ """Convert VIS radiance to reflectance factor.
+
+ Note: Produces huge reflectances in situations where both radiance and
+ solar zenith angle are small. Maybe the corresponding uncertainties
+ can be used to filter these cases before calculating reflectances.
+
+ Reference: [PUG], equation (6).
+ """
+ sza = self.solar_zenith_angle.where(
+ da.fabs(self.solar_zenith_angle) < 90,
+ np.float32(np.nan)
+ ) # direct illumination only
+ cos_sza = np.cos(np.deg2rad(sza))
+ refl = (
+ (np.pi * self.coefs['distance_sun_earth'] ** 2) /
+ (self.coefs['solar_irradiance'] * cos_sza) *
+ rad
+ )
+ return self.refl_factor_to_percent(refl)
+
+ def update_refl_attrs(self, refl):
+ """Update attributes of reflectance datasets."""
+ refl.attrs['sun_earth_distance_correction_applied'] = True
+ refl.attrs['sun_earth_distance_correction_factor'] = self.coefs[
+ 'distance_sun_earth'].item()
+ return refl
+
+ @staticmethod
+ def refl_factor_to_percent(refl):
+ """Convert reflectance factor to percent."""
+ return refl * 100
+
+
+class Navigator:
+ """Navigate MVIRI images."""
+
+ def get_area_def(self, im_size, projection_longitude):
+ """Create MVIRI area definition."""
+ proj_params = self._get_proj_params(im_size, projection_longitude)
+ extent = get_area_extent(proj_params)
+ return get_area_definition(proj_params, extent)
+
+ def _get_proj_params(self, im_size, projection_longitude):
+ """Get projection parameters for the given settings."""
+ area_name = 'geos_mviri_{0}x{0}'.format(im_size)
+ lfac, cfac, loff, coff = self._get_factors_offsets(im_size)
+ return {
+ 'ssp_lon': projection_longitude,
+ 'a': EQUATOR_RADIUS,
+ 'b': POLE_RADIUS,
+ 'h': ALTITUDE,
+ 'units': 'm',
+ 'loff': loff - im_size,
+ 'coff': coff,
+ 'lfac': -lfac,
+ 'cfac': -cfac,
+ 'nlines': im_size,
+ 'ncols': im_size,
+ 'scandir': 'S2N', # Reference: [PUG] section 2.
+ 'p_id': area_name,
+ 'a_name': area_name,
+ 'a_desc': 'MVIRI Geostationary Projection'
+ }
+
+ def _get_factors_offsets(self, im_size):
+ """Determine line/column offsets and scaling factors."""
+ # For offsets see variables "asamp" and "aline" of subroutine
+ # "refgeo" in [Handbook] and in
+ # https://github.com/FIDUCEO/FCDR_MVIRI/blob/master/lib/nrCrunch/cruncher.f
+ loff = coff = im_size / 2 + 0.5
+ lfac = cfac = sampling_to_lfac_cfac(
+ np.deg2rad(MVIRI_FIELD_OF_VIEW) / im_size
+ )
+ return lfac, cfac, loff, coff
+
+
+class Interpolator:
+ """Interpolate datasets to another resolution."""
+
+ @staticmethod
+ def interp_tiepoints(ds, target_x, target_y):
+ """Interpolate dataset between tiepoints.
+
+ Uses linear interpolation.
+
+ FUTURE: [PUG] recommends cubic spline interpolation.
+
+ Args:
+ ds:
+ Dataset to be interpolated
+ target_x:
+ Target x coordinates
+ target_y:
+ Target y coordinates
+ """
+ # No tiepoint coordinates specified in the files. Use dimensions
+ # to calculate tiepoint sampling and assign tiepoint coordinates
+ # accordingly.
+ sampling = target_x.size // ds.coords['x'].size
+ ds = ds.assign_coords(x=target_x.values[::sampling],
+ y=target_y.values[::sampling])
+
+ return ds.interp(x=target_x.values, y=target_y.values)
+
+ @staticmethod
+ def interp_acq_time(time2d, target_y):
+ """Interpolate scanline acquisition time to the given coordinates.
+
+ The files provide timestamps per pixel for the low resolution
+ channels (IR/WV) only.
+
+ 1) Average values in each line to obtain one timestamp per line.
+ 2) For the VIS channel duplicate values in y-direction (as
+ advised by [PUG]).
+
+ Note that the timestamps do not increase monotonically with the
+ line number in some cases.
+
+ Returns:
+ Mean scanline acquisition timestamps
+ """
+ # Compute mean timestamp per scanline
+ time = time2d.mean(dim='x')
+
+ # If required, repeat timestamps in y-direction to obtain higher
+ # resolution
+ y = time.coords['y'].values
+ if y.size < target_y.size:
+ reps = target_y.size // y.size
+ y_rep = np.repeat(y, reps)
+ time_hires = time.reindex(y=y_rep)
+ time_hires = time_hires.assign_coords(y=target_y)
+ return time_hires
+ return time
+
+
+class VisQualityControl:
+ """Simple quality control for VIS channel."""
+
+ def __init__(self, mask):
+ """Initialize the quality control."""
+ self._mask = mask
+
+ def check(self):
+ """Check VIS channel quality and issue a warning if it's bad."""
+ use_with_caution = da.bitwise_and(self._mask, 2)
+ if use_with_caution.all():
+ warnings.warn(
+ 'All pixels of the VIS channel are flagged as "use with '
+ 'caution". Use datasets "quality_pixel_bitmask" and '
+ '"data_quality_bitmask" to find out why.'
+ )
+
+ def mask(self, ds):
+ """Mask VIS pixels with bad quality.
+
+ Pixels are considered bad quality if the "quality_pixel_bitmask" is
+ everything else than 0 (no flag set).
+ """
+ return ds.where(self._mask == 0, np.float32(np.nan))
+
+
+def is_high_resol(resolution):
+ """Identify high resolution channel."""
+ return resolution == HIGH_RESOL
+
+
+class DatasetWrapper:
+ """Helper class for accessing the dataset."""
+
+ def __init__(self, nc):
+ """Wrap the given dataset."""
+ self.nc = nc
+
+ @property
+ def attrs(self):
+ """Exposes dataset attributes."""
+ return self.nc.attrs
+
+ def __getitem__(self, item):
+ """Get a variable from the dataset."""
+ ds = self.nc[item]
+ if self._should_dims_be_renamed(ds):
+ ds = self._rename_dims(ds)
+ elif self._coordinates_not_assigned(ds):
+ ds = self._reassign_coords(ds)
+ self._cleanup_attrs(ds)
+ return ds
+
+ def _should_dims_be_renamed(self, ds):
+ """Determine whether dataset dimensions need to be renamed."""
+ return 'y_ir_wv' in ds.dims or 'y_tie' in ds.dims
+
+ def _rename_dims(self, ds):
+ """Rename dataset dimensions to match satpy's expectations."""
+ new_names = {
+ 'y_ir_wv': 'y',
+ 'x_ir_wv': 'x',
+ 'y_tie': 'y',
+ 'x_tie': 'x'
+ }
+ for old_name, new_name in new_names.items():
+ if old_name in ds.dims:
+ ds = ds.rename({old_name: new_name})
+ return ds
+
+ def _coordinates_not_assigned(self, ds):
+ return 'y' in ds.dims and 'y' not in ds.coords
+
+ def _reassign_coords(self, ds):
+ """Re-assign coordinates.
+
+ For some reason xarray doesn't assign coordinates to all high
+ resolution data variables.
+ """
+ return ds.assign_coords({'y': self.nc.coords['y'],
+ 'x': self.nc.coords['x']})
+
+ def _cleanup_attrs(self, ds):
+ """Cleanup dataset attributes."""
+ # Remove ancillary_variables attribute to avoid downstream
+ # satpy warnings.
+ ds.attrs.pop('ancillary_variables', None)
+
+ def get_time(self):
+ """Get time coordinate.
+
+ Variable is sometimes named "time" and sometimes "time_ir_wv".
+ """
+ try:
+ return self['time_ir_wv']
+ except KeyError:
+ return self['time']
+
+ def get_xy_coords(self, resolution):
+ """Get x and y coordinates for the given resolution."""
+ if is_high_resol(resolution):
+ return self.nc.coords['x'], self.nc.coords['y']
+ return self.nc.coords['x_ir_wv'], self.nc.coords['x_ir_wv']
+
+ def get_image_size(self, resolution):
+ """Get image size for the given resolution."""
+ if is_high_resol(resolution):
+ return self.nc.coords['y'].size
+ return self.nc.coords['y_ir_wv'].size
+
+
+class FiduceoMviriBase(BaseFileHandler):
+ """Baseclass for FIDUCEO MVIRI file handlers."""
+
+ nc_keys = {
+ 'WV': 'count_wv',
+ 'IR': 'count_ir'
+ }
+
+ def __init__(self, filename, filename_info, filetype_info,
+ mask_bad_quality=False):
+ """Initialize the file handler.
+
+ Args:
+ mask_bad_quality: Mask VIS pixels with bad quality, that means
+ any quality flag except "ok". If you need more control, use
+ the ``quality_pixel_bitmask`` and ``data_quality_bitmask``
+ datasets.
+ """
+ super(FiduceoMviriBase, self).__init__(
+ filename, filename_info, filetype_info)
+ self.mask_bad_quality = mask_bad_quality
+ nc_raw = xr.open_dataset(
+ filename,
+ chunks={'x': CHUNK_SIZE,
+ 'y': CHUNK_SIZE,
+ 'x_ir_wv': CHUNK_SIZE,
+ 'y_ir_wv': CHUNK_SIZE}
+ )
+ self.nc = DatasetWrapper(nc_raw)
+
+ # Projection longitude is not provided in the file, read it from the
+ # filename.
+ self.projection_longitude = float(filename_info['projection_longitude'])
+ self.calib_coefs = self._get_calib_coefs()
+
+ def get_dataset(self, dataset_id, dataset_info):
+ """Get the dataset."""
+ name = dataset_id['name']
+ resolution = dataset_id['resolution']
+ if name in ANGLES:
+ ds = self._get_angles(name, resolution)
+ elif name in CHANNELS:
+ ds = self._get_channel(name, resolution, dataset_id['calibration'])
+ else:
+ ds = self._get_other_dataset(name)
+ self._update_attrs(ds, dataset_info)
+ return ds
+
+ def get_area_def(self, dataset_id):
+ """Get area definition of the given dataset."""
+ im_size = self.nc.get_image_size(dataset_id['resolution'])
+ nav = Navigator()
+ return nav.get_area_def(
+ im_size=im_size,
+ projection_longitude=self.projection_longitude
+ )
+
+ def _get_channel(self, name, resolution, calibration):
+ """Get and calibrate channel data."""
+ ds = self.nc[self.nc_keys[name]]
+ ds = self._calibrate(
+ ds,
+ channel=name,
+ calibration=calibration
+ )
+ if name == 'VIS':
+ qc = VisQualityControl(self.nc['quality_pixel_bitmask'])
+ if self.mask_bad_quality:
+ ds = qc.mask(ds)
+ else:
+ qc.check()
+ ds['acq_time'] = ('y', self._get_acq_time(resolution))
+ return ds
+
+ @lru_cache(maxsize=8) # 4 angle datasets with two resolutions each
+ def _get_angles(self, name, resolution):
+ """Get angle dataset.
+
+ Files provide angles (solar/satellite zenith & azimuth) at a coarser
+ resolution. Interpolate them to the desired resolution.
+ """
+ angles = self.nc[name]
+ target_x, target_y = self.nc.get_xy_coords(resolution)
+ return Interpolator.interp_tiepoints(
+ angles,
+ target_x=target_x,
+ target_y=target_y
+ )
+
+ def _get_other_dataset(self, name):
+ """Get other datasets such as uncertainties."""
+ ds = self.nc[name]
+ if name in OTHER_REFLECTANCES:
+ ds = VISCalibrator.refl_factor_to_percent(ds)
+ return ds
+
+ def _update_attrs(self, ds, info):
+ """Update dataset attributes."""
+ ds.attrs.update(info)
+ ds.attrs.update({'platform': self.filename_info['platform'],
+ 'sensor': self.filename_info['sensor']})
+ ds.attrs['raw_metadata'] = self.nc.attrs
+ ds.attrs['orbital_parameters'] = self._get_orbital_parameters()
+
+ def _calibrate(self, ds, channel, calibration):
+ """Calibrate the given dataset."""
+ if channel == 'VIS':
+ return self._calibrate_vis(ds, channel, calibration)
+ calib = IRWVCalibrator(self.calib_coefs[channel])
+ return calib.calibrate(ds, calibration)
+
+ @abc.abstractmethod
+ def _calibrate_vis(self, ds, channel, calibration): # pragma: no cover
+ """Calibrate VIS channel. To be implemented by subclasses."""
+ raise NotImplementedError
+
+ def _get_calib_coefs(self):
+ """Get calibration coefficients for all channels.
+
+ Note: Only coefficients present in both file types.
+ """
+ coefs = {
+ 'VIS': {
+ 'distance_sun_earth': self.nc['distance_sun_earth'],
+ 'solar_irradiance': self.nc['solar_irradiance_vis']
+ },
+ 'IR': {
+ 'a': self.nc['a_ir'],
+ 'b': self.nc['b_ir'],
+ 'bt_a': self.nc['bt_a_ir'],
+ 'bt_b': self.nc['bt_b_ir']
+ },
+ 'WV': {
+ 'a': self.nc['a_wv'],
+ 'b': self.nc['b_wv'],
+ 'bt_a': self.nc['bt_a_wv'],
+ 'bt_b': self.nc['bt_b_wv']
+ },
+ }
+
+ # Convert coefficients to 32bit float to reduce memory footprint
+ # of calibrated data.
+ for ch in coefs:
+ for name in coefs[ch]:
+ coefs[ch][name] = np.float32(coefs[ch][name])
+
+ return coefs
+
+ @lru_cache(maxsize=3) # Three channels
+ def _get_acq_time(self, resolution):
+ """Get scanline acquisition time for the given resolution.
+
+ Note that the acquisition time does not increase monotonically
+ with the scanline number due to the scan pattern and rectification.
+ """
+ time2d = self.nc.get_time()
+ _, target_y = self.nc.get_xy_coords(resolution)
+ return Interpolator.interp_acq_time(time2d, target_y=target_y.values)
+
+ def _get_orbital_parameters(self):
+ """Get the orbital parameters."""
+ orbital_parameters = {
+ 'projection_longitude': self.projection_longitude,
+ 'projection_latitude': 0.0,
+ 'projection_altitude': ALTITUDE
+ }
+ ssp_lon, ssp_lat = self._get_ssp_lonlat()
+ if not np.isnan(ssp_lon) and not np.isnan(ssp_lat):
+ orbital_parameters.update({
+ 'satellite_actual_longitude': ssp_lon,
+ 'satellite_actual_latitude': ssp_lat,
+ # altitude not available
+ })
+ return orbital_parameters
+
+ def _get_ssp_lonlat(self):
+ """Get longitude and latitude at the subsatellite point.
+
+ Easy FCDR files provide satellite position at the beginning and
+ end of the scan. This method computes the mean of those two values.
+ In the full FCDR the information seems to be missing.
+
+ Returns:
+ Subsatellite longitude and latitude
+ """
+ ssp_lon = self._get_ssp('longitude')
+ ssp_lat = self._get_ssp('latitude')
+ return ssp_lon, ssp_lat
+
+ def _get_ssp(self, coord):
+ key_start = 'sub_satellite_{}_start'.format(coord)
+ key_end = 'sub_satellite_{}_end'.format(coord)
+ try:
+ sub_lonlat = np.nanmean(
+ [self.nc[key_start].values,
+ self.nc[key_end].values]
+ )
+ except KeyError:
+ # Variables seem to be missing in full FCDR
+ sub_lonlat = np.nan
+ return sub_lonlat
+
+
+class FiduceoMviriEasyFcdrFileHandler(FiduceoMviriBase):
+ """File handler for FIDUCEO MVIRI Easy FCDR."""
+
+ nc_keys = FiduceoMviriBase.nc_keys.copy()
+ nc_keys['VIS'] = 'toa_bidirectional_reflectance_vis'
+
+ def _calibrate_vis(self, ds, channel, calibration):
+ """Calibrate VIS channel.
+
+ Easy FCDR provides reflectance only, no counts or radiance.
+ """
+ if calibration == 'reflectance':
+ coefs = self.calib_coefs[channel]
+ cal = VISCalibrator(coefs)
+ refl = cal.refl_factor_to_percent(ds)
+ refl = cal.update_refl_attrs(refl)
+ return refl
+ elif calibration in ('counts', 'radiance'):
+ raise ValueError('Cannot calibrate to {}. Easy FCDR provides '
+ 'reflectance only.'.format(calibration.name))
+ else:
+ raise KeyError('Invalid calibration: {}'.format(calibration.name))
+
+
+class FiduceoMviriFullFcdrFileHandler(FiduceoMviriBase):
+ """File handler for FIDUCEO MVIRI Full FCDR."""
+
+ nc_keys = FiduceoMviriBase.nc_keys.copy()
+ nc_keys['VIS'] = 'count_vis'
+
+ def _get_calib_coefs(self):
+ """Add additional VIS coefficients only present in full FCDR."""
+ coefs = super()._get_calib_coefs()
+ coefs['VIS'].update({
+ 'years_since_launch': np.float32(self.nc['years_since_launch']),
+ 'a0': np.float32(self.nc['a0_vis']),
+ 'a1': np.float32(self.nc['a1_vis']),
+ 'a2': np.float32(self.nc['a2_vis']),
+ 'mean_count_space': np.float32(
+ self.nc['mean_count_space_vis']
+ )
+ })
+ return coefs
+
+ def _calibrate_vis(self, ds, channel, calibration):
+ """Calibrate VIS channel."""
+ sza = None
+ if calibration == 'reflectance':
+ sza = self._get_angles('solar_zenith_angle', HIGH_RESOL)
+ cal = VISCalibrator(self.calib_coefs[channel], sza)
+ return cal.calibrate(ds, calibration)
diff --git a/satpy/tests/reader_tests/test_geos_area.py b/satpy/tests/reader_tests/test_geos_area.py
index 53811e6500..87e0c2cf46 100644
--- a/satpy/tests/reader_tests/test_geos_area.py
+++ b/satpy/tests/reader_tests/test_geos_area.py
@@ -20,7 +20,8 @@
import unittest
from satpy.readers._geos_area import (get_xy_from_linecol,
get_area_extent,
- get_area_definition)
+ get_area_definition,
+ sampling_to_lfac_cfac)
import numpy as np
@@ -145,3 +146,9 @@ def test_get_area_definition(self):
self.assertEqual(a, 6378169)
self.assertEqual(b, 6356583.8)
self.assertEqual(a_def.proj_dict['h'], 35785831)
+
+ def test_sampling_to_lfac_cfac(self):
+ """Test conversion from angular sampling to line/column offset."""
+ lfac = 13642337 # SEVIRI LFAC
+ sampling = np.deg2rad(2 ** 16 / lfac)
+ np.testing.assert_allclose(sampling_to_lfac_cfac(sampling), lfac)
diff --git a/satpy/tests/reader_tests/test_mviri_l1b_fiduceo_nc.py b/satpy/tests/reader_tests/test_mviri_l1b_fiduceo_nc.py
new file mode 100644
index 0000000000..11bbc84c88
--- /dev/null
+++ b/satpy/tests/reader_tests/test_mviri_l1b_fiduceo_nc.py
@@ -0,0 +1,623 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+# Copyright (c) 2020 Satpy developers
+#
+# This file is part of satpy.
+#
+# satpy is free software: you can redistribute it and/or modify it under the
+# terms of the GNU General Public License as published by the Free Software
+# Foundation, either version 3 of the License, or (at your option) any later
+# version.
+#
+# satpy is distributed in the hope that it will be useful, but WITHOUT ANY
+# WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
+# A PARTICULAR PURPOSE. See the GNU General Public License for more details.
+#
+# You should have received a copy of the GNU General Public License along with
+# satpy. If not, see .
+"""Unit tests for the FIDUCEO MVIRI FCDR Reader."""
+
+import os
+from unittest import mock
+
+import dask.array as da
+import numpy as np
+import pytest
+import xarray as xr
+from pyresample.geometry import AreaDefinition
+from pyresample.utils import proj4_radius_parameters
+
+from satpy.readers.mviri_l1b_fiduceo_nc import (
+ ALTITUDE, EQUATOR_RADIUS, POLE_RADIUS, FiduceoMviriEasyFcdrFileHandler,
+ FiduceoMviriFullFcdrFileHandler, DatasetWrapper
+)
+from satpy.tests.utils import make_dataid
+
+attrs_exp = {
+ 'platform': 'MET7',
+ 'raw_metadata': {'foo': 'bar'},
+ 'sensor': 'MVIRI',
+ 'orbital_parameters': {
+ 'projection_longitude': 57.0,
+ 'projection_latitude': 0.0,
+ 'projection_altitude': 35785860.0,
+ 'satellite_actual_longitude': 57.1,
+ 'satellite_actual_latitude': 0.1,
+ }
+
+}
+attrs_refl_exp = attrs_exp.copy()
+attrs_refl_exp.update(
+ {'sun_earth_distance_correction_applied': True,
+ 'sun_earth_distance_correction_factor': 1.}
+)
+acq_time_vis_exp = [np.datetime64('1970-01-01 00:30'),
+ np.datetime64('1970-01-01 00:30'),
+ np.datetime64('1970-01-01 02:30'),
+ np.datetime64('1970-01-01 02:30')]
+vis_counts_exp = xr.DataArray(
+ np.array(
+ [[0., 17., 34., 51.],
+ [68., 85., 102., 119.],
+ [136., 153., np.nan, 187.],
+ [204., 221., 238., 255]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ 'acq_time': ('y', acq_time_vis_exp),
+ },
+ attrs=attrs_exp
+)
+vis_rad_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 18.56, 38.28, 58.],
+ [77.72, 97.44, 117.16, 136.88],
+ [156.6, 176.32, np.nan, 215.76],
+ [235.48, 255.2, 274.92, 294.64]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ 'acq_time': ('y', acq_time_vis_exp),
+ },
+ attrs=attrs_exp
+)
+vis_refl_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 23.440929, np.nan, np.nan],
+ [40.658744, 66.602233, 147.970867, np.nan],
+ [75.688217, 92.240733, np.nan, np.nan],
+ [np.nan, np.nan, np.nan, np.nan]],
+ dtype=np.float32
+ ),
+ # (0, 0) and (2, 2) are NaN because radiance is NaN
+ # (0, 2) is NaN because SZA >= 90 degrees
+ # Last row/col is NaN due to SZA interpolation
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ 'acq_time': ('y', acq_time_vis_exp),
+ },
+ attrs=attrs_refl_exp
+)
+u_vis_refl_exp = xr.DataArray(
+ np.array(
+ [[0.1, 0.2, 0.3, 0.4],
+ [0.5, 0.6, 0.7, 0.8],
+ [0.9, 1.0, 1.1, 1.2],
+ [1.3, 1.4, 1.5, 1.6]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ },
+ attrs=attrs_exp
+)
+acq_time_ir_wv_exp = [np.datetime64('1970-01-01 00:30'),
+ np.datetime64('1970-01-01 02:30')]
+wv_counts_exp = xr.DataArray(
+ np.array(
+ [[0, 85],
+ [170, 255]],
+ dtype=np.uint8
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+wv_rad_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 3.75],
+ [8, 12.25]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+wv_bt_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 230.461366],
+ [252.507448, 266.863289]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+ir_counts_exp = xr.DataArray(
+ np.array(
+ [[0, 85],
+ [170, 255]],
+ dtype=np.uint8
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+ir_rad_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 80],
+ [165, 250]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+ir_bt_exp = xr.DataArray(
+ np.array(
+ [[np.nan, 178.00013189],
+ [204.32955838, 223.28709913]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ 'acq_time': ('y', acq_time_ir_wv_exp),
+ },
+ attrs=attrs_exp
+)
+quality_pixel_bitmask_exp = xr.DataArray(
+ np.array(
+ [[0, 0, 0, 0],
+ [0, 0, 0, 0],
+ [0, 0, 1, 0],
+ [0, 0, 0, 0]],
+ dtype=np.uint8
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ },
+ attrs=attrs_exp
+)
+sza_vis_exp = xr.DataArray(
+ np.array(
+ [[45., 67.5, 90., np.nan],
+ [22.5, 45., 67.5, np.nan],
+ [0., 22.5, 45., np.nan],
+ [np.nan, np.nan, np.nan, np.nan]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ },
+ attrs=attrs_exp
+)
+sza_ir_wv_exp = xr.DataArray(
+ np.array(
+ [[45, 90],
+ [0, 45]],
+ dtype=np.float32
+ ),
+ dims=('y', 'x'),
+ coords={
+ 'y': [1, 2],
+ 'x': [1, 2],
+ },
+ attrs=attrs_exp
+)
+area_vis_exp = AreaDefinition(
+ area_id='geos_mviri_4x4',
+ proj_id='geos_mviri_4x4',
+ description='MVIRI Geostationary Projection',
+ projection={
+ 'proj': 'geos',
+ 'lon_0': 57.0,
+ 'h': ALTITUDE,
+ 'a': EQUATOR_RADIUS,
+ 'b': POLE_RADIUS
+ },
+ width=4,
+ height=4,
+ area_extent=[5621229.74392, 5621229.74392, -5621229.74392, -5621229.74392]
+)
+area_ir_wv_exp = area_vis_exp.copy(
+ area_id='geos_mviri_2x2',
+ proj_id='geos_mviri_2x2',
+ width=2,
+ height=2
+)
+
+
+@pytest.fixture(name='fake_dataset')
+def fixture_fake_dataset():
+ """Create fake dataset."""
+ count_ir = da.linspace(0, 255, 4, dtype=np.uint8).reshape(2, 2)
+ count_wv = da.linspace(0, 255, 4, dtype=np.uint8).reshape(2, 2)
+ count_vis = da.linspace(0, 255, 16, dtype=np.uint8).reshape(4, 4)
+ sza = da.from_array(
+ np.array(
+ [[45, 90],
+ [0, 45]],
+ dtype=np.float32
+ )
+ )
+ mask = da.from_array(
+ np.array(
+ [[0, 0, 0, 0],
+ [0, 0, 0, 0],
+ [0, 0, 1, 0], # 1 = "invalid"
+ [0, 0, 0, 0]],
+ dtype=np.uint8
+ )
+ )
+ time = np.arange(4).astype('datetime64[h]').reshape(2, 2)
+ ds = xr.Dataset(
+ data_vars={
+ 'count_vis': (('y', 'x'), count_vis),
+ 'count_wv': (('y_ir_wv', 'x_ir_wv'), count_wv),
+ 'count_ir': (('y_ir_wv', 'x_ir_wv'), count_ir),
+ 'toa_bidirectional_reflectance_vis': (
+ ('y', 'x'), vis_refl_exp / 100),
+ 'u_independent_toa_bidirectional_reflectance': (
+ ('y', 'x'), u_vis_refl_exp / 100),
+ 'quality_pixel_bitmask': (('y', 'x'), mask),
+ 'solar_zenith_angle': (('y_tie', 'x_tie'), sza),
+ 'time_ir_wv': (('y_ir_wv', 'x_ir_wv'), time),
+ 'a_ir': -5.0,
+ 'b_ir': 1.0,
+ 'bt_a_ir': 10.0,
+ 'bt_b_ir': -1000.0,
+ 'a_wv': -0.5,
+ 'b_wv': 0.05,
+ 'bt_a_wv': 10.0,
+ 'bt_b_wv': -2000.0,
+ 'years_since_launch': 20.0,
+ 'a0_vis': 1.0,
+ 'a1_vis': 0.01,
+ 'a2_vis': -0.0001,
+ 'mean_count_space_vis': 1.0,
+ 'distance_sun_earth': 1.0,
+ 'solar_irradiance_vis': 650.0,
+ 'sub_satellite_longitude_start': 57.1,
+ 'sub_satellite_longitude_end': np.nan,
+ 'sub_satellite_latitude_start': np.nan,
+ 'sub_satellite_latitude_end': 0.1,
+ },
+ coords={
+ 'y': [1, 2, 3, 4],
+ 'x': [1, 2, 3, 4],
+ 'y_ir_wv': [1, 2],
+ 'x_ir_wv': [1, 2],
+ 'y_tie': [1, 2],
+ 'x_tie': [1, 2]
+
+ },
+ attrs={'foo': 'bar'}
+ )
+ ds['count_ir'].attrs['ancillary_variables'] = 'a_ir b_ir'
+ ds['count_wv'].attrs['ancillary_variables'] = 'a_wv b_wv'
+ return ds
+
+
+@pytest.fixture(
+ name='file_handler',
+ params=[FiduceoMviriEasyFcdrFileHandler,
+ FiduceoMviriFullFcdrFileHandler]
+)
+def fixture_file_handler(fake_dataset, request):
+ """Create mocked file handler."""
+ marker = request.node.get_closest_marker("file_handler_data")
+ mask_bad_quality = True
+ if marker:
+ mask_bad_quality = marker.kwargs['mask_bad_quality']
+ fh_class = request.param
+ with mock.patch('satpy.readers.mviri_l1b_fiduceo_nc.xr.open_dataset') as open_dataset:
+ open_dataset.return_value = fake_dataset
+ return fh_class(
+ filename='filename',
+ filename_info={'platform': 'MET7',
+ 'sensor': 'MVIRI',
+ 'projection_longitude': '57.0'},
+ filetype_info={'foo': 'bar'},
+ mask_bad_quality=mask_bad_quality
+ )
+
+
+@pytest.fixture(name='reader')
+def fixture_reader():
+ """Return MVIRI FIDUCEO FCDR reader."""
+ from satpy.config import config_search_paths
+ from satpy.readers import load_reader
+
+ reader_configs = config_search_paths(
+ os.path.join("readers", "mviri_l1b_fiduceo_nc.yaml"))
+ reader = load_reader(reader_configs)
+ return reader
+
+
+class TestFiduceoMviriFileHandlers:
+ """Unit tests for FIDUCEO MVIRI file handlers."""
+
+ def test_init(self, file_handler):
+ """Test file handler initialization."""
+ assert file_handler.projection_longitude == 57.0
+ assert file_handler.mask_bad_quality is True
+
+ @pytest.mark.parametrize(
+ ('name', 'calibration', 'resolution', 'expected'),
+ [
+ ('VIS', 'counts', 2250, vis_counts_exp),
+ ('VIS', 'radiance', 2250, vis_rad_exp),
+ ('VIS', 'reflectance', 2250, vis_refl_exp),
+ ('WV', 'counts', 4500, wv_counts_exp),
+ ('WV', 'radiance', 4500, wv_rad_exp),
+ ('WV', 'brightness_temperature', 4500, wv_bt_exp),
+ ('IR', 'counts', 4500, ir_counts_exp),
+ ('IR', 'radiance', 4500, ir_rad_exp),
+ ('IR', 'brightness_temperature', 4500, ir_bt_exp),
+ ('quality_pixel_bitmask', None, 2250, quality_pixel_bitmask_exp),
+ ('solar_zenith_angle', None, 2250, sza_vis_exp),
+ ('solar_zenith_angle', None, 4500, sza_ir_wv_exp),
+ ('u_independent_toa_bidirectional_reflectance', None, 4500, u_vis_refl_exp)
+ ]
+ )
+ def test_get_dataset(self, file_handler, name, calibration, resolution,
+ expected):
+ """Test getting datasets."""
+ id_keys = {'name': name, 'resolution': resolution}
+ if calibration:
+ id_keys['calibration'] = calibration
+ dataset_id = make_dataid(**id_keys)
+ dataset_info = {'platform': 'MET7'}
+
+ is_easy = isinstance(file_handler, FiduceoMviriEasyFcdrFileHandler)
+ is_vis = name == 'VIS'
+ is_refl = calibration == 'reflectance'
+ if is_easy and is_vis and not is_refl:
+ # VIS counts/radiance not available in easy FCDR
+ with pytest.raises(ValueError):
+ file_handler.get_dataset(dataset_id, dataset_info)
+ else:
+ ds = file_handler.get_dataset(dataset_id, dataset_info)
+ xr.testing.assert_allclose(ds, expected)
+ assert ds.dtype == expected.dtype
+ assert ds.attrs == expected.attrs
+
+ def test_get_dataset_corrupt(self, file_handler):
+ """Test getting datasets with known corruptions."""
+ # Time may have different names and satellite position might be missing
+ file_handler.nc.nc = file_handler.nc.nc.rename(
+ {'time_ir_wv': 'time'}
+ )
+ file_handler.nc.nc = file_handler.nc.nc.drop_vars(
+ ['sub_satellite_longitude_start']
+ )
+
+ dataset_id = make_dataid(
+ name='VIS',
+ calibration='reflectance',
+ resolution=2250
+ )
+ ds = file_handler.get_dataset(dataset_id, {'platform': 'MET7'})
+ assert 'actual_satellite_longitude' not in ds.attrs['orbital_parameters']
+ assert 'actual_satellite_latitude' not in ds.attrs['orbital_parameters']
+ xr.testing.assert_allclose(ds, vis_refl_exp)
+
+ @mock.patch(
+ 'satpy.readers.mviri_l1b_fiduceo_nc.Interpolator.interp_acq_time'
+ )
+ def test_time_cache(self, interp_acq_time, file_handler):
+ """Test caching of acquisition times."""
+ dataset_id = make_dataid(
+ name='VIS',
+ resolution=2250,
+ calibration='reflectance'
+ )
+ info = {}
+ interp_acq_time.return_value = xr.DataArray([1, 2, 3, 4], dims='y')
+
+ # Cache init
+ file_handler.get_dataset(dataset_id, info)
+ interp_acq_time.assert_called()
+
+ # Cache hit
+ interp_acq_time.reset_mock()
+ file_handler.get_dataset(dataset_id, info)
+ interp_acq_time.assert_not_called()
+
+ # Cache miss
+ interp_acq_time.return_value = xr.DataArray([1, 2], dims='y')
+ another_id = make_dataid(
+ name='IR',
+ resolution=4500,
+ calibration='brightness_temperature'
+ )
+ interp_acq_time.reset_mock()
+ file_handler.get_dataset(another_id, info)
+ interp_acq_time.assert_called()
+
+ @mock.patch(
+ 'satpy.readers.mviri_l1b_fiduceo_nc.Interpolator.interp_tiepoints'
+ )
+ def test_angle_cache(self, interp_tiepoints, file_handler):
+ """Test caching of angle datasets."""
+ dataset_id = make_dataid(name='solar_zenith_angle',
+ resolution=2250)
+ info = {}
+
+ # Cache init
+ file_handler.get_dataset(dataset_id, info)
+ interp_tiepoints.assert_called()
+
+ # Cache hit
+ interp_tiepoints.reset_mock()
+ file_handler.get_dataset(dataset_id, info)
+ interp_tiepoints.assert_not_called()
+
+ # Cache miss
+ another_id = make_dataid(name='solar_zenith_angle',
+ resolution=4500)
+ interp_tiepoints.reset_mock()
+ file_handler.get_dataset(another_id, info)
+ interp_tiepoints.assert_called()
+
+ @pytest.mark.parametrize(
+ ('name', 'resolution', 'area_exp'),
+ [
+ ('VIS', 2250, area_vis_exp),
+ ('WV', 4500, area_ir_wv_exp),
+ ('IR', 4500, area_ir_wv_exp),
+ ('quality_pixel_bitmask', 2250, area_vis_exp),
+ ('solar_zenith_angle', 2250, area_vis_exp),
+ ('solar_zenith_angle', 4500, area_ir_wv_exp)
+ ]
+ )
+ def test_get_area_definition(self, file_handler, name, resolution,
+ area_exp):
+ """Test getting area definitions."""
+ dataset_id = make_dataid(name=name, resolution=resolution)
+ area = file_handler.get_area_def(dataset_id)
+ a, b = proj4_radius_parameters(area.proj_dict)
+ a_exp, b_exp = proj4_radius_parameters(area_exp.proj_dict)
+ assert a == a_exp
+ assert b == b_exp
+ assert area.width == area_exp.width
+ assert area.height == area_exp.height
+ for key in ['h', 'lon_0', 'proj', 'units']:
+ assert area.proj_dict[key] == area_exp.proj_dict[key]
+ np.testing.assert_allclose(area.area_extent, area_exp.area_extent)
+
+ def test_calib_exceptions(self, file_handler):
+ """Test calibration exceptions."""
+ with pytest.raises(KeyError):
+ file_handler.get_dataset(
+ make_dataid(name='solar_zenith_angle', calibration='counts'),
+ {}
+ )
+ with pytest.raises(KeyError):
+ file_handler.get_dataset(
+ make_dataid(
+ name='VIS',
+ resolution=2250,
+ calibration='brightness_temperature'),
+ {}
+ )
+ with pytest.raises(KeyError):
+ file_handler.get_dataset(
+ make_dataid(
+ name='IR',
+ resolution=4500,
+ calibration='reflectance'),
+ {}
+ )
+ if isinstance(file_handler, FiduceoMviriEasyFcdrFileHandler):
+ with pytest.raises(KeyError):
+ file_handler.get_dataset(
+ {'name': 'VIS', 'calibration': 'counts'},
+ {}
+ ) # not available in easy FCDR
+
+ @pytest.mark.file_handler_data(mask_bad_quality=False)
+ def test_bad_quality_warning(self, file_handler):
+ """Test warning about bad VIS quality."""
+ file_handler.nc.nc['quality_pixel_bitmask'] = 2
+ vis = make_dataid(name='VIS', resolution=2250,
+ calibration='reflectance')
+ with pytest.warns(UserWarning):
+ file_handler.get_dataset(vis, {})
+
+ def test_file_pattern(self, reader):
+ """Test file pattern matching."""
+ filenames = [
+ "FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_FULL_v2.6_fv3.1.nc",
+ "FIDUCEO_FCDR_L15_MVIRI_MET7-57.0_201701201000_201701201030_EASY_v2.6_fv3.1.nc",
+ "FIDUCEO_FCDR_L15_MVIRI_MET7-00.0_201701201000_201701201030_EASY_v2.6_fv3.1.nc",
+ "abcde",
+ ]
+
+ files = reader.select_files_from_pathnames(filenames)
+ # only 3 out of 4 above should match
+ assert len(files) == 3
+
+
+class TestDatasetWrapper:
+ """Unit tests for DatasetWrapper class."""
+
+ def test_reassign_coords(self):
+ """Test reassigning of coordinates.
+
+ For some reason xarray does not always assign (y, x) coordinates to
+ the high resolution datasets, although they have dimensions (y, x) and
+ coordinates y and x exist. A dataset with these properties seems
+ impossible to create (neither dropping, resetting or deleting
+ coordinates seems to work). Instead use mock as a workaround.
+ """
+ nc = mock.MagicMock(
+ coords={
+ 'y': [.1, .2],
+ 'x': [.3, .4]
+ },
+ dims=('y', 'x')
+ )
+ nc.__getitem__.return_value = xr.DataArray(
+ [[1, 2],
+ [3, 4]],
+ dims=('y', 'x')
+ )
+ foo_exp = xr.DataArray(
+ [[1, 2],
+ [3, 4]],
+ dims=('y', 'x'),
+ coords={
+ 'y': [.1, .2],
+ 'x': [.3, .4]
+ }
+ )
+ ds = DatasetWrapper(nc)
+ foo = ds['foo']
+ xr.testing.assert_equal(foo, foo_exp)