From e1e1aa6275eae88e8bad470c2ee50092817cfbc0 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Thu, 5 May 2022 18:10:24 +0000 Subject: [PATCH 01/11] Refactor calibration tests --- level1c4pps/tests/test_seviri2pps.py | 128 ++++++++++++++++----------- 1 file changed, 74 insertions(+), 54 deletions(-) diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index 0acfee8..d18f710 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -24,6 +24,7 @@ import datetime as dt import numpy as np +import pytest import unittest try: from unittest import mock @@ -490,7 +491,7 @@ def test_set_nominal_scan_time(self): self.assertEqual(arr.attrs['end_time'], end_time) -class TestCalibration(unittest.TestCase): +class TestCalibration: """Test SEVIRI calibration.""" def test_get_calibration_for_date(self): @@ -502,64 +503,83 @@ def test_get_calibration_for_date(self): 'VIS008': {'gain': 0.029757990399999996, 'offset': -1.5176575103999999}, 'IR_016': {'gain': 0.0228774688, 'offset': -1.1667509087999999}} - for channel in REF.keys(): - self.assertEqual(coefs[channel]['gain'], REF[channel]['gain']) - self.assertEqual(coefs[channel]['offset'], REF[channel]['offset']) - - def test_get_calibration_for_time(self): + self._assert_coefs_close(coefs, REF) + + def _assert_coefs_close(self, coefs, expected): + for channel in expected.keys(): + np.testing.assert_allclose(coefs[channel]['gain'], + expected[channel]['gain']) + np.testing.assert_allclose(coefs[channel]['offset'], + expected[channel]['offset']) + + @pytest.mark.parametrize( + "platform,timestamp,expected", + ( + [ + ( + "MSG1", + dt.datetime(2005, 1, 18, 0, 0), + { + 'VIS006': {'gain': 0.0250354716, + 'offset': -1.2768090516000001}, + 'VIS008': {'gain': 0.0315626684, + 'offset': -1.6096960884}, + 'IR_016': {'gain': 0.022880986, + 'offset': -1.166930286} + } + ), + ( + "MSG2", + dt.datetime(2010, 1, 18, 0, 0), + { + 'VIS006': {'gain': 0.021964051999999998, + 'offset': -1.120166652}, + 'VIS008': {'gain': 0.027548445, + 'offset': -1.404970695}, + 'IR_016': {'gain': 0.021576766, + 'offset': -1.100415066} + }, + ), + ( + 'MSG3', + dt.datetime(2018, 1, 18, 0, 0), + { + 'VIS006': {'gain': 0.023689275200000002, + 'offset': -1.2081530352}, + 'VIS008': {'gain': 0.029757990399999996, + 'offset': -1.5176575103999999}, + 'IR_016': {'gain': 0.0228774688, + 'offset': -1.1667509087999999} + } + ), + ( + 'MSG4', + dt.datetime(2019, 1, 18, 0, 0), + { + 'VIS006': {'gain': 0.0230415289, + 'offset': -1.1751179739}, + 'VIS008': {'gain': 0.0291916818, + 'offset': -1.4887757718}, + 'IR_016': {'gain': 0.022223894, + 'offset': -1.1334185940000001} + } + ) + ] + ) + ) + def test_get_calibration_for_time(self, platform, timestamp, expected): """Test MODIS-intercalibrated gain and offset for specific time.""" - REF = { - ('MSG1', dt.datetime(2005, 1, 18, 0, 0)): { - 'VIS006': {'gain': 0.0250354716, - 'offset': -1.2768090516000001}, - 'VIS008': {'gain': 0.0315626684, - 'offset': -1.6096960884}, - 'IR_016': {'gain': 0.022880986, - 'offset': -1.166930286}}, - ('MSG2', dt.datetime(2010, 1, 18, 0, 0)): { - 'VIS006': {'gain': 0.021964051999999998, - 'offset': -1.120166652}, - 'VIS008': {'gain': 0.027548445, - 'offset': -1.404970695}, - 'IR_016': {'gain': 0.021576766, - 'offset': -1.100415066}}, - ('MSG3', dt.datetime(2018, 1, 18, 0, 0)): { - 'VIS006': {'gain': 0.023689275200000002, - 'offset': -1.2081530352}, - 'VIS008': {'gain': 0.029757990399999996, - 'offset': -1.5176575103999999}, - 'IR_016': {'gain': 0.0228774688, - 'offset': -1.1667509087999999}}, - ('MSG4', dt.datetime(2019, 1, 18, 0, 0)): { - 'VIS006': {'gain': 0.0230415289, - 'offset': -1.1751179739}, - 'VIS008': {'gain': 0.0291916818, - 'offset': -1.4887757718}, - 'IR_016': {'gain': 0.022223894, - 'offset': -1.1334185940000001}} - } - for (platform, time), ref in REF.items(): - coefs = calib.get_calibration_for_time(platform=platform, - time=time) - for channel in ref.keys(): - np.testing.assert_allclose(coefs[channel]['gain'], - ref[channel]['gain']) - np.testing.assert_allclose(coefs[channel]['offset'], - ref[channel]['offset']) - - def test_get_calibration(self): - """Test MODIS-intercalibrated for date and time.""" + coefs = calib.get_calibration_for_time(platform=platform, + time=timestamp) + self._assert_coefs_close(coefs, expected) + + def test_calibration_is_smooth(self): + """Test that calibration is smooth in time.""" coefs1 = calib.get_calibration_for_time( platform='MSG3', time=dt.datetime(2018, 1, 18, 23, 59)) coefs2 = calib.get_calibration_for_date( platform='MSG3', date=dt.date(2018, 1, 19)) - for channel in coefs1.keys(): - self.assertAlmostEqual(coefs1[channel]['gain'], - coefs2[channel]['gain'], - delta=10e-8) - self.assertAlmostEqual(coefs1[channel]['offset'], - coefs2[channel]['offset'], - delta=10e-8) + self._assert_coefs_close(coefs1, coefs2) class TestSEVIRIFilenameParser(unittest.TestCase): From 80e98391e18d862c0a72229c3adfa825b3fff620 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Thu, 5 May 2022 18:31:26 +0000 Subject: [PATCH 02/11] Merge date and time calibration methods --- level1c4pps/calibration_coefs.py | 36 +++------------------------- level1c4pps/seviri2pps_lib.py | 4 ++-- level1c4pps/tests/test_seviri2pps.py | 29 ++++++++++++---------- 3 files changed, 22 insertions(+), 47 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 587262a..23859c4 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -49,7 +49,6 @@ ) ) -REF_DATE = datetime.date(2000, 1, 1) REF_TIME = datetime.datetime(2000, 1, 1, 0, 0) @@ -60,10 +59,11 @@ def calib_meirink(platform, channel, time): :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] """ + if isinstance(time, datetime.date): + time = datetime.datetime.combine(time, datetime.time(0)) if time < REF_TIME: raise ValueError('Given time ({0}) is < reference time ({1})'.format( time, REF_TIME)) - a = COEFS_MEIRINK[platform][channel]['a'] b = COEFS_MEIRINK[platform][channel]['b'] delta_days = (time - REF_TIME).total_seconds() / 3600.0 / 24.0 @@ -73,26 +73,7 @@ def calib_meirink(platform, channel, time): return gain, offset -def calib_meirink_date(platform, channel, date): - """Get MODIS-intercalibrated gain and offset for SEVIRI VIS channels. - - Reference: http://msgcpp.knmi.nl/mediawiki/index.php/MSG-SEVIRI_solar_channel_calibration - - :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] - """ - if date < REF_DATE: - raise ValueError('Given date ({0}) is < reference date ({1})'.format( - date, REF_DATE)) - - a = COEFS_MEIRINK[platform][channel]['a'] - b = COEFS_MEIRINK[platform][channel]['b'] - gain = (b + a*(date - REF_DATE).days) / 1000.0 # micro Watts -> milli Watts - offset = -51.0 * gain # Space count is 51 - - return gain, offset - - -def get_calibration_for_time(platform, time): +def get_calibration(platform, time): """Get MODIS-intercalibrated gain and offset for specific time.""" coefs = {} for channel in ('VIS006', 'VIS008', 'IR_016'): @@ -103,17 +84,6 @@ def get_calibration_for_time(platform, time): return coefs -def get_calibration_for_date(platform, date): - """Get MODIS-intercalibrated gain and offset for specific date.""" - coefs = {} - for channel in ('VIS006', 'VIS008', 'IR_016'): - gain, offset = calib_meirink_date(platform=platform, channel=channel, - date=date) - coefs[channel] = {'gain': gain, 'offset': offset} - - return coefs - - if __name__ == '__main__': time = datetime.datetime(2018, 1, 18, 12, 0) platform = 'MSG3' diff --git a/level1c4pps/seviri2pps_lib.py b/level1c4pps/seviri2pps_lib.py index 6a94f17..c8e66aa 100644 --- a/level1c4pps/seviri2pps_lib.py +++ b/level1c4pps/seviri2pps_lib.py @@ -43,7 +43,7 @@ from pyorbital.astronomy import get_alt_az, sun_zenith_angle from pyorbital.orbital import get_observer_look -from level1c4pps.calibration_coefs import get_calibration_for_time, CALIB_MODE +from level1c4pps.calibration_coefs import get_calibration, CALIB_MODE from level1c4pps import make_azidiff_angle, get_encoding, compose_filename, update_angle_attributes @@ -107,7 +107,7 @@ def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate): parser = SEVIRIFilenameParser() file_format, info = parser.parse(os.path.basename(filenames[0])) - calib_coefs = get_calibration_for_time( + calib_coefs = get_calibration( platform=info['platform_shortname'], time=info['start_time'] ) diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index d18f710..a431ab5 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -494,10 +494,14 @@ def test_set_nominal_scan_time(self): class TestCalibration: """Test SEVIRI calibration.""" - def test_get_calibration_for_date(self): + @pytest.mark.parametrize( + "time", + (dt.date(2018, 1, 18), dt.datetime(2018, 1, 18)) + ) + def test_get_calibration_for_date(self, time): """Test MODIS-intercalibrated gain and offset for specific date.""" - coefs = calib.get_calibration_for_date( - platform='MSG3', date=dt.date(2018, 1, 18)) + coefs = calib.get_calibration( + platform='MSG3', time=time) REF = { 'VIS006': {'gain': 0.023689275200000002, 'offset': -1.2081530352}, 'VIS008': {'gain': 0.029757990399999996, @@ -505,12 +509,14 @@ def test_get_calibration_for_date(self): 'IR_016': {'gain': 0.0228774688, 'offset': -1.1667509087999999}} self._assert_coefs_close(coefs, REF) - def _assert_coefs_close(self, coefs, expected): + def _assert_coefs_close(self, coefs, expected, **tol): for channel in expected.keys(): np.testing.assert_allclose(coefs[channel]['gain'], - expected[channel]['gain']) + expected[channel]['gain'], + **tol) np.testing.assert_allclose(coefs[channel]['offset'], - expected[channel]['offset']) + expected[channel]['offset'], + **tol) @pytest.mark.parametrize( "platform,timestamp,expected", @@ -569,17 +575,16 @@ def _assert_coefs_close(self, coefs, expected): ) def test_get_calibration_for_time(self, platform, timestamp, expected): """Test MODIS-intercalibrated gain and offset for specific time.""" - coefs = calib.get_calibration_for_time(platform=platform, - time=timestamp) + coefs = calib.get_calibration(platform=platform, time=timestamp) self._assert_coefs_close(coefs, expected) def test_calibration_is_smooth(self): """Test that calibration is smooth in time.""" - coefs1 = calib.get_calibration_for_time( + coefs1 = calib.get_calibration( platform='MSG3', time=dt.datetime(2018, 1, 18, 23, 59)) - coefs2 = calib.get_calibration_for_date( - platform='MSG3', date=dt.date(2018, 1, 19)) - self._assert_coefs_close(coefs1, coefs2) + coefs2 = calib.get_calibration( + platform='MSG3', time=dt.datetime(2018, 1, 19)) + self._assert_coefs_close(coefs1, coefs2, atol=1e-4) class TestSEVIRIFilenameParser(unittest.TestCase): From dae200c76c02f8736dfecd42a5f58cdf61d39cc9 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 09:24:18 +0000 Subject: [PATCH 03/11] Refactor calibration --- level1c4pps/calibration_coefs.py | 36 +++++++++++++++++++++++++------- 1 file changed, 29 insertions(+), 7 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 23859c4..97c6430 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -50,6 +50,7 @@ ) REF_TIME = datetime.datetime(2000, 1, 1, 0, 0) +SPACE_COUNT = -51.0 def calib_meirink(platform, channel, time): @@ -59,20 +60,41 @@ def calib_meirink(platform, channel, time): :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] """ - if isinstance(time, datetime.date): - time = datetime.datetime.combine(time, datetime.time(0)) + time = _convert_to_datetime(time) + _check_time(time) + a = COEFS_MEIRINK[platform][channel]['a'] + b = COEFS_MEIRINK[platform][channel]['b'] + days_since_ref_time = _get_days_since_ref_time(time) + return _calc_gain_offset(a, b, days_since_ref_time) + + +def _check_time(time): if time < REF_TIME: raise ValueError('Given time ({0}) is < reference time ({1})'.format( time, REF_TIME)) - a = COEFS_MEIRINK[platform][channel]['a'] - b = COEFS_MEIRINK[platform][channel]['b'] - delta_days = (time - REF_TIME).total_seconds() / 3600.0 / 24.0 - gain = (b + a * delta_days) / 1000.0 # micro Watts -> milli Watts - offset = -51.0 * gain # Space count is 51 + +def _convert_to_datetime(date_or_time): + if isinstance(date_or_time, datetime.date): + return datetime.datetime.combine(date_or_time, datetime.time(0)) + return date_or_time + + +def _get_days_since_ref_time(time): + return (time - REF_TIME).total_seconds() / 3600.0 / 24.0 + + +def _calc_gain_offset(a, b, days_since_ref_time): + gain = (b + a * days_since_ref_time) + gain = _microwatts_to_milliwatts(gain) + offset = SPACE_COUNT * gain return gain, offset +def _microwatts_to_milliwatts(microwatts): + return microwatts / 1000.0 + + def get_calibration(platform, time): """Get MODIS-intercalibrated gain and offset for specific time.""" coefs = {} From 4141ba21d121eea14f3618f0ef449f36706f813d Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 09:25:02 +0000 Subject: [PATCH 04/11] Fix broken link to reference --- level1c4pps/calibration_coefs.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 97c6430..9b580e8 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -56,7 +56,7 @@ def calib_meirink(platform, channel, time): """Get MODIS-intercalibrated gain and offset for SEVIRI VIS channels. - Reference: http://msgcpp.knmi.nl/mediawiki/index.php/MSG-SEVIRI_solar_channel_calibration + Reference: https://msgcpp.knmi.nl/solar-channel-calibration.html :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] """ From 9fdc642b9e4f3844bd12977b6bce39a92e61b82b Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 12:06:11 +0000 Subject: [PATCH 05/11] Add clipping for SEVIRI calibration coefficients --- level1c4pps/calibration_coefs.py | 147 +++++++++++++++++---------- level1c4pps/seviri2pps_lib.py | 12 ++- level1c4pps/tests/test_seviri2pps.py | 18 ++++ 3 files changed, 118 insertions(+), 59 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 9b580e8..f251ceb 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -24,54 +24,75 @@ """Module with calibration coefficients for SEVIRI.""" import datetime - -CALIB_MODE = 'Nominal' -COEFS_MEIRINK = dict( - MSG1=dict( - VIS006=dict(b=24.346, a=0.3739E-3), - VIS008=dict(b=30.989, a=0.3111E-3), - IR_016=dict(b=22.869, a=0.0065E-3) - ), - MSG2=dict( - VIS006=dict(b=21.026, a=0.2556E-3), - VIS008=dict(b=26.875, a=0.1835E-3), - IR_016=dict(b=21.394, a=0.0498E-3) - ), - MSG3=dict( - VIS006=dict(b=19.829, a=0.5856E-3), - VIS008=dict(b=25.284, a=0.6787E-3), - IR_016=dict(b=23.066, a=-0.0286E-3) - ), - MSG4=dict( - VIS006=dict(b=21.040, a=0.2877E-3), - VIS008=dict(b=24.966, a=0.6074E-3), - IR_016=dict(b=21.236, a=0.1420E-3) +from enum import Enum + + +class CalibrationData(Enum): + COEFS = dict( + MSG1=dict( + VIS006=dict(b=24.346, a=0.3739E-3), + VIS008=dict(b=30.989, a=0.3111E-3), + IR_016=dict(b=22.869, a=0.0065E-3) + ), + MSG2=dict( + VIS006=dict(b=21.026, a=0.2556E-3), + VIS008=dict(b=26.875, a=0.1835E-3), + IR_016=dict(b=21.394, a=0.0498E-3) + ), + MSG3=dict( + VIS006=dict(b=19.829, a=0.5856E-3), + VIS008=dict(b=25.284, a=0.6787E-3), + IR_016=dict(b=23.066, a=-0.0286E-3) + ), + MSG4=dict( + VIS006=dict(b=21.040, a=0.2877E-3), + VIS008=dict(b=24.966, a=0.6074E-3), + IR_016=dict(b=21.236, a=0.1420E-3) + ) ) -) + SPACE_COUNT = -51.0 + REF_TIME = datetime.datetime(2000, 1, 1, 0, 0) + TIME_COVERAGE = { + "start": datetime.datetime(2004, 1, 1, 0, 0), + "end": datetime.datetime(2021, 1, 1, 0, 0) + } + SATPY_CALIB_MODE = 'Nominal' + + +def get_calibration(platform, time, clip_time=True): + """Get MODIS-intercalibrated gain and offset for specific time. + + Args: + platform: Platform name. + time: Date or time of observations to be calibrated. + clip_time: If True, do not extrapolate calibration + coefficients beyond the time coverage of the calibration dataset. + Instead, clip at the boundaries, that means return the boundary + coefficients for timestamps outside the coverage. + """ + coefs = {} + for channel in ('VIS006', 'VIS008', 'IR_016'): + coefs[channel] = _get_single_channel_calibration( + platform=platform, + channel=channel, + time=time, + clip_time=clip_time + ) + return coefs -REF_TIME = datetime.datetime(2000, 1, 1, 0, 0) -SPACE_COUNT = -51.0 +def _get_single_channel_calibration(platform, channel, time, clip_time): + time = _prepare_time(time, clip_time) + gain, offset = calib_meirink(platform, channel, time) + return {'gain': gain, 'offset': offset} -def calib_meirink(platform, channel, time): - """Get MODIS-intercalibrated gain and offset for SEVIRI VIS channels. - - Reference: https://msgcpp.knmi.nl/solar-channel-calibration.html - :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] - """ +def _prepare_time(time, clip_time): time = _convert_to_datetime(time) _check_time(time) - a = COEFS_MEIRINK[platform][channel]['a'] - b = COEFS_MEIRINK[platform][channel]['b'] - days_since_ref_time = _get_days_since_ref_time(time) - return _calc_gain_offset(a, b, days_since_ref_time) - - -def _check_time(time): - if time < REF_TIME: - raise ValueError('Given time ({0}) is < reference time ({1})'.format( - time, REF_TIME)) + if clip_time: + time = _clip_time_at_coverage_bounds(time) + return time def _convert_to_datetime(date_or_time): @@ -80,14 +101,43 @@ def _convert_to_datetime(date_or_time): return date_or_time +def _check_time(time): + ref_time = CalibrationData.REF_TIME.value + if time < ref_time: + raise ValueError('Given time ({0}) is < reference time ({1})'.format( + time, ref_time)) + + +def _clip_time_at_coverage_bounds(time): + time_cov = CalibrationData.TIME_COVERAGE.value + time = max(time, time_cov["start"]) + time = min(time, time_cov["end"]) + return time + + +def calib_meirink(platform, channel, time): + """Get MODIS-intercalibrated gain and offset for SEVIRI VIS channels. + + Reference: https://msgcpp.knmi.nl/solar-channel-calibration.html + + :returns: gain, offset [mW m-2 sr-1 (cm-1)-1] + """ + coefs = CalibrationData.COEFS.value + a = coefs[platform][channel]['a'] + b = coefs[platform][channel]['b'] + days_since_ref_time = _get_days_since_ref_time(time) + return _calc_gain_offset(a, b, days_since_ref_time) + + def _get_days_since_ref_time(time): - return (time - REF_TIME).total_seconds() / 3600.0 / 24.0 + ref_time = CalibrationData.REF_TIME.value + return (time - ref_time).total_seconds() / 3600.0 / 24.0 def _calc_gain_offset(a, b, days_since_ref_time): gain = (b + a * days_since_ref_time) gain = _microwatts_to_milliwatts(gain) - offset = SPACE_COUNT * gain + offset = CalibrationData.SPACE_COUNT.value * gain return gain, offset @@ -95,17 +145,6 @@ def _microwatts_to_milliwatts(microwatts): return microwatts / 1000.0 -def get_calibration(platform, time): - """Get MODIS-intercalibrated gain and offset for specific time.""" - coefs = {} - for channel in ('VIS006', 'VIS008', 'IR_016'): - gain, offset = calib_meirink(platform=platform, channel=channel, - time=time) - coefs[channel] = {'gain': gain, 'offset': offset} - - return coefs - - if __name__ == '__main__': time = datetime.datetime(2018, 1, 18, 12, 0) platform = 'MSG3' diff --git a/level1c4pps/seviri2pps_lib.py b/level1c4pps/seviri2pps_lib.py index c8e66aa..8fc5f14 100644 --- a/level1c4pps/seviri2pps_lib.py +++ b/level1c4pps/seviri2pps_lib.py @@ -43,7 +43,7 @@ from pyorbital.astronomy import get_alt_az, sun_zenith_angle from pyorbital.orbital import get_observer_look -from level1c4pps.calibration_coefs import get_calibration, CALIB_MODE +from level1c4pps.calibration_coefs import get_calibration, CalibrationData from level1c4pps import make_azidiff_angle, get_encoding, compose_filename, update_angle_attributes @@ -123,10 +123,12 @@ def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate): def _create_scene(file_format, filenames, calib_coefs): - return Scene(reader=file_format, - filenames=filenames, - reader_kwargs={'calib_mode': CALIB_MODE, - 'ext_calib_coefs': calib_coefs}) + return Scene( + reader=file_format, + filenames=filenames, + reader_kwargs={'calib_mode': CalibrationData.SATPY_CALIB_MODE.value, + 'ext_calib_coefs': calib_coefs} + ) def _check_is_seviri_data(scene): diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index a431ab5..3e2b482 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -586,6 +586,24 @@ def test_calibration_is_smooth(self): platform='MSG3', time=dt.datetime(2018, 1, 19)) self._assert_coefs_close(coefs1, coefs2, atol=1e-4) + @pytest.mark.parametrize( + "coverage_boundary,outside_coverage", + [ + (dt.datetime(2004, 1, 1), dt.datetime(2003, 1, 1)), + (dt.datetime(2021, 1, 1), dt.datetime(2022, 1, 1)), + ] + ) + def test_clip_at_time_coverage_bounds(self, coverage_boundary, outside_coverage): + coefs1 = calib.get_calibration( + platform='MSG1', + time=coverage_boundary + ) + coefs2 = calib.get_calibration( + platform='MSG1', + time=outside_coverage + ) + self._assert_coefs_close(coefs1, coefs2) + class TestSEVIRIFilenameParser(unittest.TestCase): def test_parse_native(self): From a26bfe0808963abf2e96ca41b3f04d2785fcda58 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 12:42:02 +0000 Subject: [PATCH 06/11] Expose option in high-level methods --- level1c4pps/calibration_coefs.py | 24 ++++++++++++------------ level1c4pps/seviri2pps_lib.py | 15 +++++++++++---- 2 files changed, 23 insertions(+), 16 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index f251ceb..74ac880 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -59,16 +59,16 @@ class CalibrationData(Enum): SATPY_CALIB_MODE = 'Nominal' -def get_calibration(platform, time, clip_time=True): +def get_calibration(platform, time, clip=True): """Get MODIS-intercalibrated gain and offset for specific time. Args: platform: Platform name. time: Date or time of observations to be calibrated. - clip_time: If True, do not extrapolate calibration - coefficients beyond the time coverage of the calibration dataset. - Instead, clip at the boundaries, that means return the boundary - coefficients for timestamps outside the coverage. + clip: If True, do not extrapolate calibration coefficients beyond the + time coverage of the calibration dataset. Instead, clip at the + boundaries, that means return the boundary coefficients for + timestamps outside the coverage. """ coefs = {} for channel in ('VIS006', 'VIS008', 'IR_016'): @@ -76,22 +76,22 @@ def get_calibration(platform, time, clip_time=True): platform=platform, channel=channel, time=time, - clip_time=clip_time + clip=clip ) return coefs -def _get_single_channel_calibration(platform, channel, time, clip_time): - time = _prepare_time(time, clip_time) +def _get_single_channel_calibration(platform, channel, time, clip): + time = _prepare_time(time, clip) gain, offset = calib_meirink(platform, channel, time) return {'gain': gain, 'offset': offset} -def _prepare_time(time, clip_time): +def _prepare_time(time, clip): time = _convert_to_datetime(time) _check_time(time) - if clip_time: - time = _clip_time_at_coverage_bounds(time) + if clip: + time = _clip_at_coverage_bounds(time) return time @@ -108,7 +108,7 @@ def _check_time(time): time, ref_time)) -def _clip_time_at_coverage_bounds(time): +def _clip_at_coverage_bounds(time): time_cov = CalibrationData.TIME_COVERAGE.value time = max(time, time_cov["start"]) time = min(time, time_cov["end"]) diff --git a/level1c4pps/seviri2pps_lib.py b/level1c4pps/seviri2pps_lib.py index 8fc5f14..6c8707b 100644 --- a/level1c4pps/seviri2pps_lib.py +++ b/level1c4pps/seviri2pps_lib.py @@ -89,7 +89,8 @@ class UnexpectedSatpyVersion(Exception): # MSG4-SEVI-MSG15-1234-NA-20190409121243.927000000Z -def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate): +def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate, + clip_calib): """Load and calibrate data. Uses inter-calibration coefficients from Meirink et al. @@ -99,6 +100,9 @@ def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate): apply_sun_earth_distance_correction: If True, apply sun-earth-distance- correction to visible channels. rotate: Rotate image so that pixel (0, 0) is N-W. + clip_calib: If True, do not extrapolate calibration coefficients beyond + the time coverage of the calibration dataset. Instead, clip at the + boundaries. Returns: Satpy scene holding calibrated channels @@ -109,7 +113,8 @@ def load_and_calibrate(filenames, apply_sun_earth_distance_correction, rotate): calib_coefs = get_calibration( platform=info['platform_shortname'], - time=info['start_time'] + time=info['start_time'], + clip=clip_calib ) scn_ = _create_scene(file_format, filenames, calib_coefs) _check_is_seviri_data(scn_) @@ -558,7 +563,8 @@ def _postproc_hrit(self, parsed): def process_one_scan(tslot_files, out_path, rotate=True, engine='h5netcdf', use_nominal_time_in_filename=False, - apply_sun_earth_distance_correction=True): + apply_sun_earth_distance_correction=True, + clip_calib=False): """Make level 1c files in PPS-format.""" for fname in tslot_files: if not os.path.isfile(fname): @@ -568,7 +574,8 @@ def process_one_scan(tslot_files, out_path, rotate=True, engine='h5netcdf', scn_ = load_and_calibrate( tslot_files, apply_sun_earth_distance_correction=apply_sun_earth_distance_correction, - rotate=rotate + rotate=rotate, + clip_calib=clip_calib ) # Find lat/lon data From 9d17d6a515f7c172cf5746dc90558061afe0bee3 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 12:49:45 +0000 Subject: [PATCH 07/11] Fix tests --- level1c4pps/tests/test_seviri2pps.py | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index 3e2b482..3455d6b 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -76,7 +76,8 @@ def test_load_and_calibrate(self, mocked_scene): res = seviri2pps.load_and_calibrate( filenames, apply_sun_earth_distance_correction=False, - rotate=False + rotate=False, + clip_calib=False ) # Compare results and expectations @@ -105,7 +106,8 @@ def test_load_and_calibrate_with_rotation(self, mocked_scene): res = seviri2pps.load_and_calibrate( filenames, apply_sun_earth_distance_correction=False, - rotate=True + rotate=True, + clip_calib=False ) scene.load.assert_called_with(mock.ANY, upper_right_corner='NE') self.assertTrue(res.attrs['image_rotated']) From 45eee9a2fea6abfcb4a3eb281fdb0633abf87e36 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 12:53:20 +0000 Subject: [PATCH 08/11] Improve naming of test methods --- level1c4pps/tests/test_seviri2pps.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index 3455d6b..c088ece 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -500,7 +500,7 @@ class TestCalibration: "time", (dt.date(2018, 1, 18), dt.datetime(2018, 1, 18)) ) - def test_get_calibration_for_date(self, time): + def test_get_calibration_can_handle_both_date_and_time(self, time): """Test MODIS-intercalibrated gain and offset for specific date.""" coefs = calib.get_calibration( platform='MSG3', time=time) @@ -575,7 +575,7 @@ def _assert_coefs_close(self, coefs, expected, **tol): ] ) ) - def test_get_calibration_for_time(self, platform, timestamp, expected): + def test_get_calibration(self, platform, timestamp, expected): """Test MODIS-intercalibrated gain and offset for specific time.""" coefs = calib.get_calibration(platform=platform, time=timestamp) self._assert_coefs_close(coefs, expected) From e5c49fc1c3df8997240f7503e8ddbee5eef548d8 Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Mon, 9 May 2022 13:15:56 +0000 Subject: [PATCH 09/11] Fix test misses --- level1c4pps/calibration_coefs.py | 11 ++++++++--- level1c4pps/tests/test_seviri2pps.py | 7 +++++++ 2 files changed, 15 insertions(+), 3 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 74ac880..aa6f8d2 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -89,19 +89,24 @@ def _get_single_channel_calibration(platform, channel, time, clip): def _prepare_time(time, clip): time = _convert_to_datetime(time) - _check_time(time) + _check_is_valid_time(time) if clip: time = _clip_at_coverage_bounds(time) return time def _convert_to_datetime(date_or_time): - if isinstance(date_or_time, datetime.date): + if _is_date(date_or_time): return datetime.datetime.combine(date_or_time, datetime.time(0)) return date_or_time -def _check_time(time): +def _is_date(date_or_time): + # datetime is a subclass of date, therefore we cannot use isinstance here + return type(date_or_time) == datetime.date + + +def _check_is_valid_time(time): ref_time = CalibrationData.REF_TIME.value if time < ref_time: raise ValueError('Given time ({0}) is < reference time ({1})'.format( diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index c088ece..b1f7482 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -606,6 +606,13 @@ def test_clip_at_time_coverage_bounds(self, coverage_boundary, outside_coverage) ) self._assert_coefs_close(coefs1, coefs2) + def test_fails_with_invalid_time(self): + with pytest.raises(ValueError): + calib.get_calibration( + platform='MSG1', + time=dt.datetime(1999, 1, 1) + ) + class TestSEVIRIFilenameParser(unittest.TestCase): def test_parse_native(self): From 9c1ea90c41de78f5e6baa361a2c2fb03417fb7ce Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Tue, 17 May 2022 13:57:49 +0000 Subject: [PATCH 10/11] Disable clipping by default in lower level methods. --- level1c4pps/calibration_coefs.py | 2 +- level1c4pps/tests/test_seviri2pps.py | 16 +++++++++++++++- 2 files changed, 16 insertions(+), 2 deletions(-) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index aa6f8d2..7101603 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -59,7 +59,7 @@ class CalibrationData(Enum): SATPY_CALIB_MODE = 'Nominal' -def get_calibration(platform, time, clip=True): +def get_calibration(platform, time, clip=False): """Get MODIS-intercalibrated gain and offset for specific time. Args: diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index b1f7482..471eb4f 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -571,6 +571,18 @@ def _assert_coefs_close(self, coefs, expected, **tol): 'IR_016': {'gain': 0.022223894, 'offset': -1.1334185940000001} } + ), + ( + 'MSG4', + dt.datetime(2024, 1, 18, 0, 0), + { + 'VIS006': {'gain': 0.0235668691, + 'offset': -1.2019103241}, + 'VIS008': {'gain': 0.0303007942, + 'offset': -1.5453405042000001}, + 'IR_016': {'gain': 0.022483186, + 'offset': -1.146642486} + } # extrapolation beyond time coverage of calib. dataset ) ] ) @@ -596,13 +608,15 @@ def test_calibration_is_smooth(self): ] ) def test_clip_at_time_coverage_bounds(self, coverage_boundary, outside_coverage): + """Test clipping beyond time coverage of the calibration dataset.""" coefs1 = calib.get_calibration( platform='MSG1', time=coverage_boundary ) coefs2 = calib.get_calibration( platform='MSG1', - time=outside_coverage + time=outside_coverage, + clip=True ) self._assert_coefs_close(coefs1, coefs2) From 61a886ee4ced3d2697b88af16dbee4496d8ff5fb Mon Sep 17 00:00:00 2001 From: Stephan Finkensieper Date: Tue, 17 May 2022 13:58:20 +0000 Subject: [PATCH 11/11] Add docstrings. --- level1c4pps/calibration_coefs.py | 8 ++++++++ level1c4pps/tests/test_seviri2pps.py | 1 + 2 files changed, 9 insertions(+) diff --git a/level1c4pps/calibration_coefs.py b/level1c4pps/calibration_coefs.py index 7101603..d7dce62 100644 --- a/level1c4pps/calibration_coefs.py +++ b/level1c4pps/calibration_coefs.py @@ -82,6 +82,14 @@ def get_calibration(platform, time, clip=False): def _get_single_channel_calibration(platform, channel, time, clip): + """Get calibration coefficients for a single channel. + + Args: + platform: Platform name. + channel: Channel name. + time: Observation date or time. + clip: Clip at time coverage boundaries of the calibration dataset. + """ time = _prepare_time(time, clip) gain, offset = calib_meirink(platform, channel, time) return {'gain': gain, 'offset': offset} diff --git a/level1c4pps/tests/test_seviri2pps.py b/level1c4pps/tests/test_seviri2pps.py index 471eb4f..75e7c82 100644 --- a/level1c4pps/tests/test_seviri2pps.py +++ b/level1c4pps/tests/test_seviri2pps.py @@ -621,6 +621,7 @@ def test_clip_at_time_coverage_bounds(self, coverage_boundary, outside_coverage) self._assert_coefs_close(coefs1, coefs2) def test_fails_with_invalid_time(self): + """Test that calibration fails with timestamps < reference time.""" with pytest.raises(ValueError): calib.get_calibration( platform='MSG1',