[Git][debian-gis-team/pyspectral][master] 4 commits: New upstream version 0.9.0+ds
Antonio Valentino
gitlab at salsa.debian.org
Mon Sep 2 07:14:15 BST 2019
Antonio Valentino pushed to branch master at Debian GIS Project / pyspectral
Commits:
bdb76ea7 by Antonio Valentino at 2019-09-02T05:46:16Z
New upstream version 0.9.0+ds
- - - - -
a3beaafe by Antonio Valentino at 2019-09-02T05:46:17Z
Update upstream source from tag 'upstream/0.9.0+ds'
Update to upstream version '0.9.0+ds'
with Debian dir 7df66c97398a0f371925d08159b15eb1aef5c607
- - - - -
69437df5 by Antonio Valentino at 2019-09-02T05:46:53Z
New upstream release
- - - - -
a8df3a1b by Antonio Valentino at 2019-09-02T06:04:39Z
Set distribution to unstable
- - - - -
19 changed files:
- CHANGELOG.md
- README.md
- debian/changelog
- doc/37_reflectance.rst
- doc/rad_definitions.rst
- pyspectral/blackbody.py
- pyspectral/config.py
- pyspectral/etc/pyspectral.yaml
- pyspectral/rsr_reader.py
- pyspectral/tests/test_atm_correction_ir.py
- pyspectral/tests/test_blackbody.py
- pyspectral/tests/test_rad_tb_conversions.py
- pyspectral/tests/test_reflectance.py
- pyspectral/tests/test_utils.py
- pyspectral/utils.py
- pyspectral/version.py
- rsr_convert_scripts/README.rst
- + rsr_convert_scripts/agri_rsr.py
- + rsr_convert_scripts/virr_rsr.py
Changes:
=====================================
CHANGELOG.md
=====================================
@@ -1,3 +1,27 @@
+## Version <RELEASE_VERSION> (2019/08/30)
+
+### Issues Closed
+
+* [Issue 73](https://github.com/pytroll/pyspectral/issues/73) - Fix blackbody code to work with dask arrays ([PR 74](https://github.com/pytroll/pyspectral/pull/74))
+
+In this release 1 issue was closed.
+
+### Pull Requests Merged
+
+#### Bugs fixed
+
+* [PR 80](https://github.com/pytroll/pyspectral/pull/80) - Fix doc tests for python 2&3
+* [PR 79](https://github.com/pytroll/pyspectral/pull/79) - Fix rsr zenodo version
+* [PR 74](https://github.com/pytroll/pyspectral/pull/74) - Fix dask compatibility in blackbody functions ([73](https://github.com/pytroll/pyspectral/issues/73))
+
+#### Features added
+
+* [PR 78](https://github.com/pytroll/pyspectral/pull/78) - Add FY-3B VIRR and FY-3C VIRR RSRs
+* [PR 77](https://github.com/pytroll/pyspectral/pull/77) - Add FY-4A AGRI support
+
+In this release 5 pull requests were closed.
+
+
## Version <RELEASE_VERSION> (2019/06/07)
### Issues Closed
=====================================
README.md
=====================================
@@ -1,14 +1,11 @@
PySpectral
==========
-[](https://www.codacy.com/app/adybbroe/pyspectral?utm_source=github.com&utm_medium=referral&utm_content=pytroll/pyspectral&utm_campaign=badger)
[](https://travis-ci.org/pytroll/pyspectral)
[](https://ci.appveyor.com/project/pytroll/pyspectral)
[](https://coveralls.io/github/pytroll/pyspectral?branch=master)
-[](https://landscape.io/github/pytroll/pyspectral/master)
[](https://badge.fury.io/py/pyspectral)
[](https://codeclimate.com/github/pytroll/pyspectral)
-[](https://scrutinizer-ci.com/g/pytroll/pyspectral/?branch=master)
Given a passive sensor on a meteorological satellite PySpectral provides the
relative spectral response (rsr) function(s) and offer you some basic
=====================================
debian/changelog
=====================================
@@ -1,3 +1,9 @@
+pyspectral (0.9.0+ds-1) unstable; urgency=medium
+
+ * New upstream release.
+
+ -- Antonio Valentino <antonio.valentino at tiscali.it> Mon, 02 Sep 2019 06:04:23 +0000
+
pyspectral (0.8.9+ds-2) unstable; urgency=medium
* Bump Standards-Version to 4.4.0, no changes.
=====================================
doc/37_reflectance.rst
=====================================
@@ -46,7 +46,7 @@ expressed in :math:`W/m^2 sr^{-1} \mu m^{-1}`, or using SI units :math:`W/m^2 sr
>>> viirs = RadTbConverter('Suomi-NPP', 'viirs', 'M12')
>>> rad37 = viirs.tb2radiance(tb37)
>>> print([np.round(rad, 7) for rad in rad37['radiance']])
- [369717.4765726, 355110.5207853, 314684.2788726, 173143.5424898, 116408.0007877]
+ [369717.4972296, 355110.6414922, 314684.3507084, 173143.4836477, 116408.0022674]
>>> rad37['unit']
'W/m^2 sr^-1 m^-1'
@@ -59,7 +59,7 @@ In order to get the total radiance over the band one has to multiply with the eq
>>> viirs = RadTbConverter('Suomi-NPP', 'viirs', 'M12')
>>> rad37 = viirs.tb2radiance(tb37, normalized=False)
>>> print([np.round(rad, 8) for rad in rad37['radiance']])
- [0.07037968, 0.06759909, 0.05990352, 0.03295972, 0.02215951]
+ [0.07037968, 0.06759911, 0.05990353, 0.03295971, 0.02215951]
>>> rad37['unit']
'W/m^2 sr^-1'
@@ -218,17 +218,17 @@ We can try decompose equation :eq:`refl37` above using the example of VIIRS M12
>>> print(np.isnan(nomin))
[False False False False False]
>>> print([np.round(val, 8) for val in nomin])
- [0.05083677, 0.04805618, 0.0404157, 0.01279279, 0.00204485]
+ [0.05083677, 0.0480562, 0.04041571, 0.01279277, 0.00204485]
>>> denom = np.cos(np.deg2rad(sunz))/np.pi * sflux - rad11['radiance']
>>> print(np.isnan(denom))
[False False False False False]
>>> print([np.round(val, 8) for val in denom])
- [0.23646313, 0.23645682, 0.23650559, 0.23582015, 0.2358661]
+ [0.23646312, 0.23645681, 0.23650559, 0.23582014, 0.23586609]
>>> res = nomin/denom
>>> print(np.isnan(res))
[False False False False False]
>>> print([np.round(val, 8) for val in res])
- [0.21498817, 0.2032345, 0.17088689, 0.05424807, 0.00866955]
+ [0.21498817, 0.20323458, 0.17088693, 0.05424801, 0.00866952]
Derive the emissive part of the 3.7 micron band
@@ -255,5 +255,5 @@ Using the example of the VIIRS M12 band from above this gives the following spec
>>> ['{tb:6.3f}'.format(tb=np.round(t, 4)) for t in tb]
['266.996', '267.262', '267.991', '271.033', '271.927']
>>> rad = refl_m12.emissive_part_3x(tb=False)
- >>> ['{rad:6.3f}'.format(rad=np.round(r, 3)) for r in rad]
- ['80285.149', '81458.022', '84749.639', '99761.400', '104582.030']
+ >>> ['{rad:6.3f}'.format(rad=np.round(r, 3)) for r in rad.compute()]
+ ['80285.150', '81458.022', '84749.638', '99761.401', '104582.031']
=====================================
doc/rad_definitions.rst
=====================================
@@ -229,7 +229,7 @@ And using wavelength representation:
>>> wvl = 1./wavenumber
>>> rad = blackbody(wvl, [300., 301])
>>> print("{0:10.3f} {1:10.3f}".format(rad[0], rad[1]))
- 9573178.886 9714689.259
+ 9573177.494 9714687.157
Which are the spectral radiances in SI units around :math:`11 \mu m` at
temperatures 300 and 301 Kelvin. In units of :math:`mW/m^2\ m^{-1} sr^{-1}` this becomes:
=====================================
pyspectral/blackbody.py
=====================================
@@ -1,7 +1,7 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
-# Copyright (c) 2013-2018 Adam.Dybbroe
+# Copyright (c) 2013-2019 Adam.Dybbroe
# Author(s):
@@ -23,8 +23,13 @@
"""Planck radiation equation"""
import numpy as np
-
import logging
+
+try:
+ import dask.array as da
+except ImportError:
+ da = np
+
LOG = logging.getLogger(__name__)
H_PLANCK = 6.62606957 * 1e-34 # SI-unit = [J*s]
@@ -84,38 +89,39 @@ def blackbody_wn_rad2temp(wavenumber, radiance):
function. Wavenumber space"""
if np.isscalar(radiance):
- rad = np.array([radiance, ], dtype='float64')
- else:
- rad = np.array(radiance, dtype='float64')
+ radiance = np.array([radiance], dtype='float64')
+ elif isinstance(radiance, (list, tuple)):
+ radiance = np.array(radiance, dtype='float64')
if np.isscalar(wavenumber):
- wavnum = np.array([wavenumber, ], dtype='float64')
- else:
+ wavnum = np.array([wavenumber], dtype='float64')
+ elif isinstance(wavenumber, (list, tuple)):
wavnum = np.array(wavenumber, dtype='float64')
const1 = H_PLANCK * C_SPEED / K_BOLTZMANN
const2 = 2 * H_PLANCK * C_SPEED**2
- res = const1 * wavnum / np.log(np.divide(const2 * wavnum**3, rad) + 1.0)
+ res = const1 * wavnum / np.log(
+ np.divide(const2 * wavnum**3, radiance) + 1.0)
- shape = rad.shape
+ shape = radiance.shape
resshape = res.shape
if wavnum.shape[0] == 1:
- if rad.shape[0] == 1:
+ if radiance.shape[0] == 1:
return res[0]
else:
return res[::].reshape(shape)
else:
- if rad.shape[0] == 1:
+ if radiance.shape[0] == 1:
return res[0, :]
else:
if len(shape) == 1:
- return np.reshape(res, (shape[0], resshape[1]))
+ return res.reshape((shape[0], resshape[1]))
else:
- return np.reshape(res, (shape[0], shape[1], resshape[1]))
+ return res.reshape((shape[0], shape[1], resshape[1]))
-def planck(wave, temp, wavelength=True):
- """The Planck radiation or Blackbody radiation as a function of wavelength
+def planck(wave, temperature, wavelength=True):
+ """The Planck radiation or Blackbody radiation as a function of wavelength
or wavenumber. SI units.
_planck(wave, temperature, wavelength=True)
wave = Wavelength/wavenumber or a sequence of wavelengths/wavenumbers (m or m^-1)
@@ -136,12 +142,12 @@ def planck(wave, temp, wavelength=True):
units = ['wavelengths', 'wavenumbers']
if wavelength:
LOG.debug("Using {0} when calculating the Blackbody radiance".format(
- units[(wavelength == True) - 1]))
+ units[(wavelength is True) - 1]))
- if np.isscalar(temp):
- temperature = np.array([temp, ], dtype='float64')
- else:
- temperature = np.array(temp, dtype='float64')
+ if np.isscalar(temperature):
+ temperature = np.array([temperature, ], dtype='float64')
+ elif isinstance(temperature, (list, tuple)):
+ temperature = np.array(temperature, dtype='float64')
shape = temperature.shape
if np.isscalar(wave):
@@ -157,13 +163,19 @@ def planck(wave, temp, wavelength=True):
nom = 2 * H_PLANCK * (C_SPEED ** 2) * (wln ** 3)
arg1 = H_PLANCK * C_SPEED * wln / K_BOLTZMANN
- arg2 = np.where(np.greater(np.abs(temperature), EPSILON),
- np.array(1. / temperature), -9).reshape(-1, 1)
- arg2 = np.ma.masked_array(arg2, mask=arg2 == -9)
- LOG.debug("Max and min - arg1: %s %s", str(arg1.max()), str(arg1.min()))
- LOG.debug("Max and min - arg2: %s %s", str(arg2.max()), str(arg2.min()))
+ # use dask functions when needed
+ np_ = np if isinstance(temperature, np.ndarray) else da
+ arg2 = np_.where(np.greater(np.abs(temperature), EPSILON),
+ (1. / temperature), np.nan).reshape(-1, 1)
+ if isinstance(arg2, np.ndarray):
+ # don't compute min/max if we have dask arrays
+ LOG.debug("Max and min - arg1: %s %s",
+ str(np.nanmax(arg1)), str(np.nanmin(arg1)))
+ LOG.debug("Max and min - arg2: %s %s",
+ str(np.nanmax(arg2)), str(np.nanmin(arg2)))
+
try:
- exp_arg = np.multiply(arg1.astype('float32'), arg2.astype('float32'))
+ exp_arg = np.multiply(arg1.astype('float64'), arg2.astype('float64'))
except MemoryError:
LOG.warning(("Dimensions used in numpy.multiply probably reached "
"limit!\n"
@@ -171,9 +183,9 @@ def planck(wave, temp, wavelength=True):
"and try running again"))
raise
- LOG.debug("Max and min before exp: %s %s", str(exp_arg.max()),
- str(exp_arg.min()))
- if exp_arg.min() < 0:
+ if isinstance(exp_arg, np.ndarray) and exp_arg.min() < 0:
+ LOG.debug("Max and min before exp: %s %s",
+ str(exp_arg.max()), str(exp_arg.min()))
LOG.warning("Something is fishy: \n" +
"\tDenominator might be zero or negative in radiance derivation:")
dubious = np.where(exp_arg < 0)[0]
@@ -182,7 +194,6 @@ def planck(wave, temp, wavelength=True):
denom = np.exp(exp_arg) - 1
rad = nom / denom
- rad = np.where(rad.mask, np.nan, rad.data)
radshape = rad.shape
if wln.shape[0] == 1:
if temperature.shape[0] == 1:
@@ -194,9 +205,9 @@ def planck(wave, temp, wavelength=True):
return rad[0, :]
else:
if len(shape) == 1:
- return np.reshape(rad, (shape[0], radshape[1]))
+ return rad.reshape((shape[0], radshape[1]))
else:
- return np.reshape(rad, (shape[0], shape[1], radshape[1]))
+ return rad.reshape((shape[0], shape[1], radshape[1]))
def blackbody_wn(wavenumber, temp):
=====================================
pyspectral/config.py
=====================================
@@ -28,7 +28,12 @@ import os
from os.path import expanduser
from appdirs import AppDirs
import yaml
-from collections import Mapping
+try:
+ # python 3.3+
+ from collections.abc import Mapping
+except ImportError:
+ # deprecated (above can't be done in 2.7)
+ from collections import Mapping
import pkg_resources
try:
=====================================
pyspectral/etc/pyspectral.yaml
=====================================
@@ -209,6 +209,46 @@ download_from_internet: True
# ch15: GOES-R_ABI_FM2_SRF_CWG_ch15.txt
# ch16: GOES-R_ABI_FM2_SRF_CWG_ch16.txt
+# FY-4A-agri:
+# path: /path/to/original/fy4a/agri/data
+# ch1: FY4A_AGRI_SRF_CH01.txt
+# ch2: FY4A_AGRI_SRF_CH02.txt
+# ch3: FY4A_AGRI_SRF_CH03.txt
+# ch4: FY4A_AGRI_SRF_CH04.txt
+# ch5: FY4A_AGRI_SRF_CH05.txt
+# ch6: FY4A_AGRI_SRF_CH06.txt
+# ch7: FY4A_AGRI_SRF_CH07.txt
+# ch8: FY4A_AGRI_SRF_CH08.txt
+# ch9: FY4A_AGRI_SRF_CH09.txt
+# ch10: FY4A_AGRI_SRF_CH10.txt
+# ch11: FY4A_AGRI_SRF_CH11.txt
+# ch12: FY4A_AGRI_SRF_CH12.txt
+# ch13: FY4A_AGRI_SRF_CH13.txt
+# ch14: FY4A_AGRI_SRF_CH14.txt
+
+#FY-3B-virr:
+# path: /Users/davidh/repos/git/pyspectral/virr_srf/FY3B-VIRR
+# ch1: ch1.prn
+# ch2: ch2.prn
+# ch3: ch3.prn
+# ch4: ch4.prn
+# ch5: ch5.prn
+# ch6: ch6.prn
+# ch7: ch7.prn
+# ch8: ch8.prn
+# ch9: ch9.prn
+# ch10: ch10.prn
+
+#FY-3C-virr:
+# path: /Users/davidh/repos/git/pyspectral/virr_srf/FY3C_VIRR_SRF
+# ch1: FY3C_VIRR_CH01.txt
+# ch2: FY3C_VIRR_CH02.txt
+# ch6: FY3C_VIRR_CH06.txt
+# ch7: FY3C_VIRR_CH07.txt
+# ch8: FY3C_VIRR_CH08.txt
+# ch9: FY3C_VIRR_CH09.txt
+# ch10: FY3C_VIRR_CH10.txt
+
# FY-3D-mersi-2:
# path: /path/to/original/fy3d/mersi2/data
# ch1: FY3D_MERSI_SRF_CH01_Pub.txt
=====================================
pyspectral/rsr_reader.py
=====================================
@@ -27,16 +27,15 @@ import os
import numpy as np
from glob import glob
from os.path import expanduser
-
-import logging
-LOG = logging.getLogger(__name__)
-
from pyspectral.config import get_config
from pyspectral.utils import WAVE_NUMBER
from pyspectral.utils import WAVE_LENGTH
from pyspectral.utils import (INSTRUMENTS, download_rsr)
from pyspectral.utils import (RSR_DATA_VERSION_FILENAME, RSR_DATA_VERSION)
+import logging
+LOG = logging.getLogger(__name__)
+
class RSRDataBaseClass(object):
@@ -168,11 +167,18 @@ class RelativeSpectralResponse(RSRDataBaseClass):
no_detectors_message = False
with h5py.File(self.filename, 'r') as h5f:
- self.band_names = [b.decode('utf-8') for b in h5f.attrs['band_names'].tolist()]
- self.description = h5f.attrs['description'].decode('utf-8')
+ self.band_names = h5f.attrs['band_names'].tolist()
+ self.description = h5f.attrs['description']
+ if not isinstance(self.band_names[0], str):
+ # byte array in python 3
+ self.band_names = [x.decode('utf-8') for x in self.band_names]
+ self.description = self.description.decode('utf-8')
+
if not self.platform_name:
try:
- self.platform_name = h5f.attrs['platform_name'].decode('utf-8')
+ self.platform_name = h5f.attrs['platform_name']
+ if not isinstance(self.platform_name, str):
+ self.platform_name = self.platform_name.decode('utf-8')
except KeyError:
LOG.warning("No platform_name in HDF5 file")
try:
@@ -186,6 +192,8 @@ class RelativeSpectralResponse(RSRDataBaseClass):
if not self.instrument:
try:
self.instrument = h5f.attrs['sensor'].decode('utf-8')
+ if not isinstance(self.instrument, str):
+ self.instrument = self.instrument.decode('utf-8')
except KeyError:
LOG.warning("No sensor name specified in HDF5 file")
self.instrument = INSTRUMENTS.get(self.platform_name)
=====================================
pyspectral/tests/test_atm_correction_ir.py
=====================================
@@ -1,11 +1,11 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
-# Copyright (c) 2017 Adam.Dybbroe
+# Copyright (c) 2017 - 2019 Pytroll
# Author(s):
-# Adam.Dybbroe <a000680 at c20671.ad.smhi.se>
+# Adam.Dybbroe <adam.dybbroe at smhi.se>
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
@@ -19,27 +19,17 @@
# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
-
-"""Unit tests of the atmospherical correction in the ir spectral range
-"""
+"""Unit tests of the atmospherical correction in the ir spectral range."""
+import numpy as np
+from pyspectral.atm_correction_ir import AtmosphericalCorrection
import sys
if sys.version_info < (2, 7):
import unittest2 as unittest
else:
import unittest
-import numpy as np
-from pyspectral.atm_correction_ir import AtmosphericalCorrection
-#from mock import patch
-
-from pyspectral.tests.unittest_helpers import assertNumpyArraysEqual
-
-# Mock some modules, so we don't need them for tests.
-
-#sys.modules['pyresample'] = MagicMock()
-
SATZ = np.ma.array([[48.03, 48.03002, 48.03004, 48.03006, 48.03008, 48.0301,
48.03012, 48.03014, 48.03016, 48.03018],
[48.09, 48.09002, 48.09004, 48.09006, 48.09008, 48.0901,
@@ -125,28 +115,18 @@ RES = np.ma.array([[286.03159412, 286.03162417, 286.03165421, 286.03168426,
class TestAtmCorrection(unittest.TestCase):
-
- """Class for testing pyspectral.atm_correction_ir"""
-
- def setUp(self):
- """Setup the test"""
- pass
+ """Class for testing pyspectral.atm_correction_ir."""
def test_get_correction(self):
"""Test getting the atm correction"""
this = AtmosphericalCorrection('EOS-Terra', 'modis')
atm_corr = this.get_correction(SATZ, None, TBS)
- assertNumpyArraysEqual(TBS, atm_corr)
-
- def tearDown(self):
- """Clean up"""
- pass
+ np.testing.assert_almost_equal(TBS, atm_corr)
def suite():
- """The test suite for test_atm_correction_ir.
- """
+ """Create the test suite for test_atm_correction_ir."""
loader = unittest.TestLoader()
mysuite = unittest.TestSuite()
mysuite.addTest(loader.loadTestsFromTestCase(TestAtmCorrection))
=====================================
pyspectral/tests/test_blackbody.py
=====================================
@@ -1,11 +1,11 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
-# Copyright (c) 2013, 2014, 2015, 2016, 2017 Adam.Dybbroe
+# Copyright (c) 2013 - 2019 Pytroll
# Author(s):
-# Adam.Dybbroe <a000680 at c14526.ad.smhi.se>
+# Adam.Dybbroe <adam.dybbroe at smhi.se>
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
@@ -31,10 +31,8 @@ from pyspectral.tests.unittest_helpers import assertNumpyArraysEqual
import unittest
import numpy as np
-#RAD_11MICRON_300KELVIN = 9572498.1141643394
-RAD_11MICRON_300KELVIN = 9573177.8811719529
-#RAD_11MICRON_301KELVIN = 9713997.9623772576
-RAD_11MICRON_301KELVIN = 9714688.2959563732
+RAD_11MICRON_300KELVIN = 9573176.935507433
+RAD_11MICRON_301KELVIN = 9714686.576498277
# Radiances in wavenumber space (SI-units)
WN_RAD_11MICRON_300KELVIN = 0.00115835441353
@@ -43,18 +41,28 @@ WN_RAD_11MICRON_301KELVIN = 0.00117547716523
__unittest = True
-class TestBlackbody(unittest.TestCase):
+class CustomScheduler(object):
+ """Custom dask scheduler that raises an exception if dask is computed too many times."""
+
+ def __init__(self, max_computes=1):
+ """Set starting and maximum compute counts."""
+ self.max_computes = max_computes
+ self.total_computes = 0
+
+ def __call__(self, dsk, keys, **kwargs):
+ """Compute dask task and keep track of number of times we do so."""
+ import dask
+ self.total_computes += 1
+ if self.total_computes > self.max_computes:
+ raise RuntimeError("Too many dask computations were scheduled: {}".format(self.total_computes))
+ return dask.get(dsk, keys, **kwargs)
- """Unit testing the blackbody function"""
- def setUp(self):
- """Set up"""
- return
+class TestBlackbody(unittest.TestCase):
+ """Unit testing the blackbody function"""
def test_blackbody(self):
- """Calculate the blackbody radiation from wavelengths and
- temperatures
- """
+ """Calculate the blackbody radiation from wavelengths and temperatures."""
wavel = 11. * 1E-6
black = blackbody((wavel, ), [300., 301])
self.assertEqual(black.shape[0], 2)
@@ -71,14 +79,28 @@ class TestBlackbody(unittest.TestCase):
tb_therm = np.array([[300., 301], [299, 298], [279, 286]])
black = blackbody((10. * 1E-6, 11.e-6), tb_therm)
+ self.assertIsInstance(black, np.ndarray)
tb_therm = np.array([[300., 301], [0., 298], [279, 286]])
black = blackbody((10. * 1E-6, 11.e-6), tb_therm)
+ self.assertIsInstance(black, np.ndarray)
+
+ def test_blackbody_dask(self):
+ """Calculate the blackbody radiation from wavelengths and temperatures with dask arrays."""
+ import dask
+ import dask.array as da
+ tb_therm = da.from_array([[300., 301], [299, 298], [279, 286]], chunks=2)
+ with dask.config.set(scheduler=CustomScheduler(0)):
+ black = blackbody((10. * 1E-6, 11.e-6), tb_therm)
+ self.assertIsInstance(black, da.Array)
+
+ tb_therm = da.from_array([[300., 301], [0., 298], [279, 286]], chunks=2)
+ with dask.config.set(scheduler=CustomScheduler(0)):
+ black = blackbody((10. * 1E-6, 11.e-6), tb_therm)
+ self.assertIsInstance(black, da.Array)
def test_blackbody_wn(self):
- """Calculate the blackbody radiation from wavenumbers and
- temperatures
- """
+ """Calculate the blackbody radiation from wavenumbers and temperatures."""
wavenumber = 90909.1 # 11 micron band
black = blackbody_wn((wavenumber, ), [300., 301])
self.assertEqual(black.shape[0], 2)
@@ -106,9 +128,24 @@ class TestBlackbody(unittest.TestCase):
assertNumpyArraysEqual(t__, expected)
- def tearDown(self):
- """Clean up"""
- return
+ def test_blackbody_wn_dask(self):
+ """Test that blackbody rad2temp preserves dask arrays."""
+ import dask
+ import dask.array as da
+ wavenumber = 90909.1 # 11 micron band
+ radiances = da.from_array([0.001, 0.0009, 0.0012, 0.0018], chunks=2).reshape(2, 2)
+ with dask.config.set(scheduler=CustomScheduler(0)):
+ t__ = blackbody_wn_rad2temp(wavenumber, radiances)
+ self.assertIsInstance(t__, da.Array)
+ t__ = t__.compute()
+ expected = np.array([290.3276916, 283.76115441,
+ 302.4181330, 333.1414164]).reshape(2, 2)
+ self.assertAlmostEqual(t__[1, 1], expected[1, 1], 5)
+ self.assertAlmostEqual(t__[0, 0], expected[0, 0], 5)
+ self.assertAlmostEqual(t__[0, 1], expected[0, 1], 5)
+ self.assertAlmostEqual(t__[1, 0], expected[1, 0], 5)
+
+ assertNumpyArraysEqual(t__, expected)
def suite():
=====================================
pyspectral/tests/test_rad_tb_conversions.py
=====================================
@@ -1,7 +1,7 @@
#!/usr/bin/env python
# -*- coding: utf-8 -*-
-# Copyright (c) 2014-2018 Adam.Dybbroe
+# Copyright (c) 2014-2019 Adam.Dybbroe
# Author(s):
@@ -347,10 +347,10 @@ class TestRadTbConversions(unittest.TestCase):
self.assertTrue(np.allclose(TRUE_RADS * integral, res['radiance']))
res = self.modis.tb2radiance(237., lut=False)
- self.assertAlmostEqual(16570.592171157, res['radiance'])
+ self.assertAlmostEqual(16570.579551068, res['radiance'])
res = self.modis.tb2radiance(277., lut=False)
- self.assertAlmostEqual(167544.3823631, res['radiance'])
+ self.assertAlmostEqual(167544.39368663222, res['radiance'])
res = self.modis.tb2radiance(1.1, lut=False)
self.assertAlmostEqual(0.0, res['radiance'])
@@ -362,7 +362,7 @@ class TestRadTbConversions(unittest.TestCase):
self.assertAlmostEqual(5.3940515573e-06, res['radiance'])
res = self.modis.tb2radiance(200.1, lut=False)
- self.assertAlmostEqual(865.09776189, res['radiance'])
+ self.assertAlmostEqual(865.09759706, res['radiance'])
def tearDown(self):
"""Clean up"""
=====================================
pyspectral/tests/test_reflectance.py
=====================================
@@ -139,7 +139,10 @@ class TestReflectance(unittest.TestCase):
with patch('pyspectral.radiance_tb_conversion.RelativeSpectralResponse') as mymock:
instance = mymock.return_value
- instance.rsr = TEST_RSR
+ # VIIRS doesn't have a channel '20' like MODIS so the generic
+ # mapping this test will end up using will find 'ch20' for VIIRS
+ viirs_rsr = {'ch20': TEST_RSR['20'], '99': TEST_RSR['99']}
+ instance.rsr = viirs_rsr
instance.unit = '1e-6 m'
instance.si_scale = 1e-6
@@ -148,7 +151,7 @@ class TestReflectance(unittest.TestCase):
refl37 = Calculator('Suomi-NPP', 'viirs', 3.7)
self.assertEqual(refl37.bandwavelength, 3.7)
- self.assertEqual(refl37.bandname, '20')
+ self.assertEqual(refl37.bandname, 'ch20')
with patch('pyspectral.radiance_tb_conversion.RelativeSpectralResponse') as mymock:
instance = mymock.return_value
=====================================
pyspectral/tests/test_utils.py
=====================================
@@ -115,8 +115,7 @@ class TestUtils(unittest.TestCase):
self.rsr = RsrTestData()
def test_convert2wavenumber(self):
- """Testing the conversion of rsr from wavelength to wavenumber
- """
+ """Testing the conversion of rsr from wavelength to wavenumber."""
newrsr, info = utils.convert2wavenumber(TEST_RSR)
unit = info['unit']
self.assertEqual(unit, 'cm-1')
@@ -127,11 +126,7 @@ class TestUtils(unittest.TestCase):
self.assertTrue(np.allclose(wvn_res, wvn))
def test_get_bandname_from_wavelength(self):
- """Test that it is possible to get the right bandname provided the wavelength
- in micro meters
-
- """
-
+ """Test the right bandname is found provided the wavelength in micro meters."""
x = utils.get_bandname_from_wavelength('abi', 0.4, self.rsr.rsr)
self.assertEqual(x, 'ch1')
with self.assertRaises(AttributeError):
@@ -146,19 +141,16 @@ class TestUtils(unittest.TestCase):
x = utils.get_bandname_from_wavelength('abi', 1.0, self.rsr.rsr)
self.assertEqual(x, None)
+ # uses generic channel mapping where '20' -> 'ch20'
bandname = utils.get_bandname_from_wavelength('abi', 3.7, TEST_RSR)
- self.assertEqual(bandname, '20')
+ self.assertEqual(bandname, 'ch20')
bandname = utils.get_bandname_from_wavelength('abi', 3.0, TEST_RSR)
self.assertIsNone(bandname)
- def tearDown(self):
- """Clean up"""
- pass
-
def suite():
- """The suite for test_utils."""
+ """Create the suite for test_utils."""
loader = unittest.TestLoader()
mysuite = unittest.TestSuite()
mysuite.addTest(loader.loadTestsFromTestCase(TestUtils))
=====================================
pyspectral/utils.py
=====================================
@@ -78,9 +78,19 @@ BANDNAMES['generic'] = {'VIS006': 'VIS0.6',
'C15': 'ch15',
'C16': 'ch16',
}
+# handle arbitrary channel numbers
+for chan_num in range(1, 37):
+ BANDNAMES['generic'][str(chan_num)] = 'ch{:d}'.format(chan_num)
-BANDNAMES['avhrr-3'] = {'3b': 'ch3b',
- '3a': 'ch3a'}
+# MODIS RSR files were made before 'chX' became standard in pyspectral
+BANDNAMES['modis'] = {str(chan_num): str(chan_num) for chan_num in range(1, 37)}
+
+BANDNAMES['avhrr-3'] = {'1': 'ch1',
+ '2': 'ch2',
+ '3b': 'ch3b',
+ '3a': 'ch3a',
+ '4': 'ch4',
+ '5': 'ch5'}
BANDNAMES['ahi'] = {'B01': 'ch1',
'B02': 'ch2',
@@ -120,12 +130,16 @@ INSTRUMENTS = {'NOAA-19': 'avhrr/3',
'Suomi-NPP': 'viirs',
'NOAA-20': 'viirs',
'FY-3D': 'mersi-2',
- 'Feng-Yun 3D': 'mersi-2'
+ 'FY-3C': 'virr',
+ 'FY-3B': 'virr',
+ 'Feng-Yun 3D': 'mersi-2',
+ 'FY-4A': 'agri'
}
-HTTP_PYSPECTRAL_RSR = "https://zenodo.org/record/2653487/files/pyspectral_rsr_data.tgz"
+
+HTTP_PYSPECTRAL_RSR = "https://zenodo.org/record/3381130/files/pyspectral_rsr_data.tgz"
RSR_DATA_VERSION_FILENAME = "PYSPECTRAL_RSR_VERSION"
-RSR_DATA_VERSION = "v1.0.6"
+RSR_DATA_VERSION = "v1.0.9"
ATM_CORRECTION_LUT_VERSION = {}
ATM_CORRECTION_LUT_VERSION['antarctic_aerosol'] = {'version': 'v1.0.1',
=====================================
pyspectral/version.py
=====================================
@@ -46,9 +46,9 @@ def get_keywords():
# setup.py/versioneer.py will grep for the variable names, so they must
# each be defined on a line of their own. _version.py will just call
# get_keywords().
- git_refnames = " (tag: v0.8.9)"
- git_full = "b3a1f58f49e2ab7c81c331e379065f24eb21f1ba"
- git_date = "2019-06-07 09:27:42 +0200"
+ git_refnames = " (HEAD -> master, tag: v0.9.0)"
+ git_full = "2372396547ecace7fbcefc50d5dcbaa32e4a5177"
+ git_date = "2019-08-30 18:04:19 +0200"
keywords = {"refnames": git_refnames, "full": git_full, "date": git_date}
return keywords
=====================================
rsr_convert_scripts/README.rst
=====================================
@@ -176,3 +176,12 @@ the pyspectral.yaml file:
Adam Dybbroe
Sat Dec 1 17:39:48 2018
+.. code::
+
+ %> python virr_rsr.py
+
+Converting the FY-3B or FY-3C VIRR spectral responses to HDF5. Original files
+for FY-3B come as ``.prn`` text files for each channel (ex. ``ch1.prn``). For
+FY-3C they come as ``.txt`` text files for channels 1, 2, 6, 7, 8, 9, and 10
+only with names like ``FY3C_VIRR_CH01.txt``.
+
=====================================
rsr_convert_scripts/agri_rsr.py
=====================================
@@ -0,0 +1,108 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+
+# Copyright (c) 2018, 2019 Pytroll
+
+# Author(s):
+
+# Xin.Zhang <xinzhang1215 at gmail.com>
+# Adam.Dybbroe <adam.dybbroe at smhi.se>
+
+# This program 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.
+
+# This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
+
+"""Read the FY-4A AGRI relative spectral responses. Data from
+http://fy4.nsmc.org.cn/portal/cn/fycv/srf.html
+"""
+import os
+import numpy as np
+from pyspectral.utils import INSTRUMENTS
+from pyspectral.utils import convert2hdf5 as tohdf5
+from pyspectral.raw_reader import InstrumentRSR
+from pyspectral.utils import logging_on, get_logger
+
+FY4A_BAND_NAMES = ['ch1', 'ch2', 'ch3', 'ch4', 'ch5', 'ch6', 'ch7', 'ch8',
+ 'ch9', 'ch10', 'ch11', 'ch12', 'ch13', 'ch14']
+BANDNAME_SCALE2MICROMETERS = {'ch1': 0.001,
+ 'ch2': 0.001,
+ 'ch3': 0.001,
+ 'ch4': 1.0,
+ 'ch5': 1.0,
+ 'ch6': 1.0,
+ 'ch7': 1.0,
+ 'ch8': 1.0,
+ 'ch9': 1.0,
+ 'ch10': 1.0,
+ 'ch11': 1.0,
+ 'ch12': 1.0,
+ 'ch13': 1.0,
+ 'ch14': 1.0}
+
+
+class AGRIRSR(InstrumentRSR):
+ """Container for the FY-4 AGRI RSR data"""
+
+ def __init__(self, bandname, platform_name):
+ """Initialise the FY-4 AGRI relative spectral response data"""
+ super(AGRIRSR, self).__init__(bandname, platform_name, FY4A_BAND_NAMES)
+
+ self.instrument = INSTRUMENTS.get(platform_name, 'agri')
+
+ self._get_options_from_config()
+ self._get_bandfilenames()
+
+ LOG.debug("Filenames: %s", str(self.filenames))
+ if self.filenames[bandname] and os.path.exists(self.filenames[bandname]):
+ self.requested_band_filename = self.filenames[bandname]
+ scale = BANDNAME_SCALE2MICROMETERS.get(bandname)
+ if scale:
+ self._load(scale=scale)
+ else:
+ LOG.error(
+ "Failed determine the scale used to convert to wavelength in micrometers - channel = %s", bandname)
+ raise AttributeError('no scale for bandname %s', bandname)
+
+ else:
+ LOG.warning("Couldn't find an existing file for this band: %s",
+ str(self.bandname))
+
+ self.filename = self.requested_band_filename
+
+ def _load(self, scale=0.001):
+ """Load the AGRI RSR data for the band requested
+
+ Wavelength is given in nanometers.
+ """
+ data = np.genfromtxt(self.requested_band_filename,
+ unpack=True,
+ names=['wavelength',
+ 'response'],
+ skip_header=0)
+
+ wavelength = data['wavelength'] * scale
+ response = data['response']
+
+ self.rsr = {'wavelength': wavelength, 'response': response}
+
+
+def main():
+ """Main"""
+ for platform_name in ["FY-4A", ]:
+ tohdf5(AGRIRSR, platform_name, FY4A_BAND_NAMES)
+
+
+if __name__ == "__main__":
+ LOG = get_logger(__name__)
+ logging_on()
+
+ main()
=====================================
rsr_convert_scripts/virr_rsr.py
=====================================
@@ -0,0 +1,89 @@
+#!/usr/bin/env python
+# -*- coding: utf-8 -*-
+#
+# Copyright (c) 2019 Adam.Dybbroe
+#
+# Author(s):
+#
+# David Hoese <david.hoese at ssec.wisc.edu>
+#
+# This program 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.
+#
+# This program 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 this program. If not, see <http://www.gnu.org/licenses/>.
+"""Read the VIRR relative spectral responses.
+
+Data from http://gsics.nsmc.org.cn/portal/en/fycv/srf.html
+
+"""
+
+import os
+import numpy as np
+from pyspectral.utils import INSTRUMENTS
+from pyspectral.utils import convert2hdf5 as tohdf5
+from pyspectral.raw_reader import InstrumentRSR
+
+import logging
+LOG = logging.getLogger(__name__)
+
+VIRR_BAND_NAMES = {
+ 'FY-3B': ['ch{:d}'.format(x) for x in range(1, 11)],
+ 'FY-3C': ['ch1', 'ch2'] + ['ch{:d}'.format(x) for x in range(6, 11)],
+}
+
+
+class VirrRSR(InstrumentRSR):
+ """Container for the FY-3B/FY-3C VIRR RSR data."""
+
+ def __init__(self, bandname, platform_name):
+ """Verify that file exists and can be read."""
+ super(VirrRSR, self).__init__(bandname, platform_name, VIRR_BAND_NAMES[platform_name])
+
+ self.instrument = INSTRUMENTS.get(platform_name, 'virr')
+ self._get_options_from_config()
+ self._get_bandfilenames()
+
+ LOG.debug("Filenames: %s", str(self.filenames))
+ if self.filenames[bandname] and os.path.exists(self.filenames[bandname]):
+ self.requested_band_filename = self.filenames[bandname]
+ self._load()
+ else:
+ LOG.warning("Couldn't find an existing file for this band: %s",
+ str(self.bandname))
+
+ # To be compatible with VIIRS....
+ self.filename = self.requested_band_filename
+
+ def _load(self, scale=0.001):
+ """Load the VIRR RSR data for the band requested.
+
+ Wavelength is given in nanometers.
+ """
+ data = np.genfromtxt(self.requested_band_filename,
+ unpack=True,
+ names=['wavelength',
+ 'response'],
+ skip_header=0)
+
+ wavelength = data['wavelength'] * scale
+ response = data['response']
+
+ self.rsr = {'wavelength': wavelength, 'response': response}
+
+
+def main():
+ """Main"""
+ for platform_name, band_names in VIRR_BAND_NAMES.items():
+ tohdf5(VirrRSR, platform_name, band_names)
+
+
+if __name__ == "__main__":
+ main()
View it on GitLab: https://salsa.debian.org/debian-gis-team/pyspectral/compare/c9a3be279065fdb63807a9576f59a156996a9f51...a8df3a1bbb29fa6754bafa6b71920a3c9293f853
--
View it on GitLab: https://salsa.debian.org/debian-gis-team/pyspectral/compare/c9a3be279065fdb63807a9576f59a156996a9f51...a8df3a1bbb29fa6754bafa6b71920a3c9293f853
You're receiving this email because of your account on salsa.debian.org.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://alioth-lists.debian.net/pipermail/pkg-grass-devel/attachments/20190902/a751eaa2/attachment-0001.html>
More information about the Pkg-grass-devel
mailing list