Merge pull request #5 from pytroll/master

Bring my fork up to date

CHANGELOG.md

@@ -1,3 +1,27 @@
+## Version <RELEASE_VERSION> (2019/04/29)
+
+### Issues Closed
+
+* [Issue 70](https://github.com/pytroll/pyspectral/issues/70) - Update yaml usage to work with pyyaml 5.1+ ([PR 71](https://github.com/pytroll/pyspectral/pull/71))
+* [Issue 66](https://github.com/pytroll/pyspectral/issues/66) - Throws a warning about non-existing directory - storing/reading cached radiance-tb look-up-tables ([PR 67](https://github.com/pytroll/pyspectral/pull/67))
+* [Issue 61](https://github.com/pytroll/pyspectral/issues/61) - can this program be used for user-defined sensor or rsr？ ([PR 62](https://github.com/pytroll/pyspectral/pull/62))
+* [Issue 58](https://github.com/pytroll/pyspectral/issues/58) - Use dask instead of numpy masked arrays ([PR 59](https://github.com/pytroll/pyspectral/pull/59))
+
+In this release 4 issues were closed.
+
+### Pull Requests Merged
+
+#### Features added
+
+* [PR 71](https://github.com/pytroll/pyspectral/pull/71) - Fix yaml 5.1+ support with unsafe loading ([70](https://github.com/pytroll/pyspectral/issues/70), [70](https://github.com/pytroll/pyspectral/issues/70))
+* [PR 69](https://github.com/pytroll/pyspectral/pull/69) - Feature rayleigh catch exception
+* [PR 68](https://github.com/pytroll/pyspectral/pull/68) - Feaure metimage multiple detectors
+* [PR 65](https://github.com/pytroll/pyspectral/pull/65) - Feature add metimage
+
+
+In this release 4 pull requests were closed.
+
+
 ## Version <RELEASE_VERSION> (2019/04/09)
 
 ### Issues Closed
@@ -12,8 +36,10 @@ In this release 3 issues were closed.
 
 #### Bugs fixed
 
+* [PR 67](https://github.com/pytroll/pyspectral/pull/67) - Bugfix tb2rad lut caching ([66](https://github.com/pytroll/pyspectral/issues/66))
 * [PR 64](https://github.com/pytroll/pyspectral/pull/64) - Fix interp function in rayleigh correction to be serializable
 
+
 #### Features added
 
 * [PR 59](https://github.com/pytroll/pyspectral/pull/59) - Daskify NIR reflectance calculations ([58](https://github.com/pytroll/pyspectral/issues/58))

doc/usage.rst

@@ -18,7 +18,7 @@ Now, you can work with the data as you wish, make some simple plot for instance:
   >>> [str(b) for b in olci.band_names]
   ['Oa01', 'Oa02', 'Oa03', 'Oa04', 'Oa05', 'Oa06', 'Oa07', 'Oa08', 'Oa09', 'Oa10', 'Oa11', 'Oa12', 'Oa13', 'Oa14', 'Oa15', 'Oa16', 'Oa17', 'Oa18', 'Oa19', 'Oa20', 'Oa21']
   >>> print("Central wavelength = {wvl:7.6f}".format(wvl=olci.rsr['Oa01']['det-1']['central_wavelength']))
-  Central wavelength = 0.400123
+  Central wavelength = 0.400303
   >>> import matplotlib.pyplot as plt
   >>> dummy = plt.figure(figsize=(10, 5))
   >>> import numpy as np

pyspectral/config.py

@@ -31,6 +31,11 @@
 from collections import Mapping
 import pkg_resources
 
+try:
+    from yaml import UnsafeLoader
+except ImportError:
+    from yaml import Loader as UnsafeLoader
+
 
 LOG = logging.getLogger(__name__)
 
@@ -72,7 +77,7 @@ def get_config():
 
     config = {}
     with open(configfile, 'r') as fp_:
-        config = recursive_dict_update(config, yaml.load(fp_))
+        config = recursive_dict_update(config, yaml.load(fp_, Loader=UnsafeLoader))
 
     app_dirs = AppDirs('pyspectral', 'pytroll')
     user_datadir = app_dirs.user_data_dir

pyspectral/utils.py

@@ -123,9 +123,9 @@
                'Feng-Yun 3D': 'mersi-2'
                }
 
-HTTP_PYSPECTRAL_RSR = "https://zenodo.org/record/2617441/files/pyspectral_rsr_data.tgz"
+HTTP_PYSPECTRAL_RSR = "https://zenodo.org/record/2653487/files/pyspectral_rsr_data.tgz"
 RSR_DATA_VERSION_FILENAME = "PYSPECTRAL_RSR_VERSION"
-RSR_DATA_VERSION = "v1.0.5"
+RSR_DATA_VERSION = "v1.0.6"
 
 ATM_CORRECTION_LUT_VERSION = {}
 ATM_CORRECTION_LUT_VERSION['antarctic_aerosol'] = {'version': 'v1.0.1',

rsr_convert_scripts/metimage_rsr.py

@@ -28,7 +28,7 @@
 
 import os
 import numpy as np
-from pyspectral.utils import convert2hdf5 as tohdf5
+from pyspectral.utils import get_central_wave
 from pyspectral.raw_reader import InstrumentRSR
 import logging
 
@@ -80,13 +80,84 @@ def _load(self, scale=1.0):
         wavelength = 1. / data['wavenumber'] * 10000.
         response = data['response']
 
-        self.rsr = {'wavelength': wavelength, 'response': response}
+        # The real MetImage has 24 detectors. However, for now we store the
+        # single rsr as 'detector-1', indicating that there will be multiple
+        # detectors in the future:
+        detectors = {}
+        detectors['det-1'] = {'wavelength': wavelength, 'response': response}
+        self.rsr = detectors
+
+
+def generate_metimage_file(platform_name):
+    """Retrieve original RSR data and convert to internal hdf5 format.
+    """
+    import h5py
+
+    bandnames = METIMAGE_BAND_NAMES
+    instr = MetImageRSR(bandnames[0], platform_name)
+    instr_name = instr.instrument.replace('/', '')
+    filename = os.path.join(instr.output_dir,
+                            "rsr_{0}_{1}.h5".format(instr_name,
+                                                    platform_name))
+
+    with h5py.File(filename, "w") as h5f:
+        h5f.attrs['description'] = ('Relative Spectral Responses for ' +
+                                    instr.instrument.upper())
+        h5f.attrs['platform_name'] = platform_name
+        h5f.attrs['band_names'] = bandnames
+
+        for chname in bandnames:
+            metimage = MetImageRSR(chname, platform_name)
+            grp = h5f.create_group(chname)
+            grp.attrs['number_of_detectors'] = len(metimage.rsr.keys())
+            # Loop over each detector to check if the sampling wavelengths are
+            # identical:
+            det_names = list(metimage.rsr.keys())
+            wvl = metimage.rsr[det_names[0]]['wavelength']
+            wvl, idx = np.unique(wvl, return_index=True)
+            wvl_is_constant = True
+            for det in det_names[1:]:
+                det_wvl = np.unique(metimage.rsr[det]['wavelength'])
+                if not np.alltrue(wvl == det_wvl):
+                    LOG.warning(
+                        "Wavelngth arrays are not the same among detectors!")
+                    wvl_is_constant = False
+
+            if wvl_is_constant:
+                arr = wvl
+                dset = grp.create_dataset('wavelength', arr.shape, dtype='f')
+                dset.attrs['unit'] = 'm'
+                dset.attrs['scale'] = 1e-06
+                dset[...] = arr
+
+            # Loop over each detector:
+            for det in metimage.rsr:
+                det_grp = grp.create_group(det)
+                wvl = metimage.rsr[det]['wavelength'][
+                    ~np.isnan(metimage.rsr[det]['wavelength'])]
+                rsp = metimage.rsr[det]['response'][
+                    ~np.isnan(metimage.rsr[det]['wavelength'])]
+                wvl, idx = np.unique(wvl, return_index=True)
+                rsp = np.take(rsp, idx)
+                LOG.debug("wvl.shape: %s", str(wvl.shape))
+                det_grp.attrs[
+                    'central_wavelength'] = get_central_wave(wvl, rsp)
+                if not wvl_is_constant:
+                    arr = wvl
+                    dset = det_grp.create_dataset(
+                        'wavelength', arr.shape, dtype='f')
+                    dset.attrs['unit'] = 'm'
+                    dset.attrs['scale'] = 1e-06
+                    dset[...] = arr
+
+                dset = det_grp.create_dataset('response', rsp.shape, dtype='f')
+                dset[...] = rsp
 
 
 def main():
     """Main"""
     for platform_name in ["Metop-SG-A1", ]:
-        tohdf5(MetImageRSR, platform_name, METIMAGE_BAND_NAMES)
+        generate_metimage_file(platform_name)
 
 
 if __name__ == "__main__":

rsr_convert_scripts/olci_rsr.py

@@ -1,7 +1,7 @@
 #!/usr/bin/env python
 # -*- coding: utf-8 -*-
 
-# Copyright (c) 2016, 2017, 2018 Adam.Dybbroe
+# Copyright (c) 2016 - 2019 Adam.Dybbroe
 
 # Author(s):
 
@@ -23,7 +23,10 @@
 """
 Reading the Sentinel-3 OLCI relative spectral responses
 
-https://sentinel.esa.int/documents/247904/322304/OLCI+SRF+%28NetCDF%29/15cfd7a6-b7bc-4051-87f8-c35d765ae43a
+https://sentinel.esa.int/documents/247904/322304/OLCI+SRF+%28NetCDF%29/
+
+https://sentinel.esa.int/web/sentinel/technical-guides/sentinel-3-olci/olci-instrument/spectral-response-function-data
+
 """
 
 from netCDF4 import Dataset
@@ -34,7 +37,8 @@
 
 LOG = logging.getLogger(__name__)
 
-RSRFILE = '/home/a000680/data/SpectralResponses/olci/OLCISRFNetCDF.nc4'
+RSRFILE = {'Sentinel-3A': '/home/a000680/data/SpectralResponses/olci/S3A_OL_SRF_20160713_mean_rsr.nc4',
+           'Sentinel-3B': '/home/a000680/data/SpectralResponses/olciB/S3B_OL_SRF_0_20180109_mean_rsr.nc4'}
 
 
 OLCI_BAND_NAMES = ['Oa01', 'Oa02', 'Oa03', 'Oa04',
@@ -73,18 +77,24 @@ def _load(self, scale=0.001):
         ncf = Dataset(self.path, 'r')
 
         bandnum = OLCI_BAND_NAMES.index(self.bandname)
-        cam = 0
-        view = 0
+        # cam = 0
+        # view = 0
+        # resp = ncf.variables[
+        #     'spectral_response_function'][bandnum, cam, view, :]
+        # wvl = ncf.variables[
+        #     'spectral_response_function_wavelength'][bandnum, cam, view, :] * scale
         resp = ncf.variables[
-            'spectral_response_function'][bandnum, cam, view, :]
+            'mean_spectral_response_function'][bandnum, :]
         wvl = ncf.variables[
-            'spectral_response_function_wavelength'][bandnum, cam, view, :] * scale
+            'mean_spectral_response_function_wavelength'][bandnum, :] * scale
+
         self.rsr = {'wavelength': wvl, 'response': resp}
 
 
 def main():
     """Main"""
-    for platform_name in ['Sentinel-3A', ]:
+    # for platform_name in ['Sentinel-3A', 'Sentinel-3B']:
+    for platform_name in ['Sentinel-3B', ]:
         tohdf5(OlciRSR, platform_name, OLCI_BAND_NAMES)
 
 
@@ -119,5 +129,6 @@ def testplot():
     fig.savefig('olci_band1.png')
     # pylab.show()
 
+
 if __name__ == "__main__":
     main()