From 785f1c5b68c14b986c022d18fce41e00a8c3bdbe Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Tue, 3 Oct 2023 23:04:50 +0200 Subject: [PATCH 01/18] add meirink calibration method --- satpy/readers/seviri_base.py | 69 +++++++++++++++++++++++++++++++- satpy/readers/seviri_l1b_hrit.py | 1 + 2 files changed, 68 insertions(+), 2 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 131fe39ad4..c1d288f4df 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -166,7 +166,7 @@ """ import warnings -from datetime import timedelta +from datetime import datetime, timedelta import dask.array as da import numpy as np @@ -353,6 +353,62 @@ 'ALPHA': 0.9981, 'BETA': 0.5635}} +# Calibration coefficients from Meirink, J.F., R.A. Roebeling and P. Stammes, 2013: +# Inter-calibration of polar imager solar channels using SEVIRI, Atm. Meas. Tech., 6, +# 2495-2508, doi:10.5194/amt-6-2495-2013 + +# To obtain the slope for the calibration, one should use the routine get_seviri_meirink_slope + +# Epoch for the MEIRINK re-calibration +DATE_2000 = datetime(2000, 1, 1) + +MEIRINK_COEFS = {} + +# Meteosat-8 + +MEIRINK_COEFS[321] = {'VIS006': (24.346, 0.3739), + 'VIS008': (30.989, 0.3111), + 'IR_016': (22.869, 0.0065) + } + +# Meteosat-9 + +MEIRINK_COEFS[322] = {'VIS006': (21.026, 0.3739), + 'VIS008': (26.875, 0.3111), + 'IR_016': (21.394, 0.0065) + } + +# Meteosat-10 + +MEIRINK_COEFS[323] = {'VIS006': (19.829, 0.5856), + 'VIS008': (25.284, 0.6787), + 'IR_016': (23.066, -0.0286) + } + +# Meteosat-11 + +MEIRINK_COEFS[324] = {'VIS006': (20.515, 0.3600), + 'VIS008': (25.803, 0.4844), + 'IR_016': (22.354, -0.0187) + } + + +def get_meirink_slope(meirink_coefs, acquisition_time): + """Compute the slope for the visible channel calibration according to Meirink 2013. + + S = A + B * 1.e-3* Day + + S is here in µW m-2 sr-1 (cm-1)-1 + + EUMETSAT calibration is given in mW m-2 sr-1 (cm-1)-1, so an extra factor of 1/1000 must + be applied. + """ + A = meirink_coefs[0] + B = meirink_coefs[1] + delta_t = (acquisition_time - DATE_2000).total_seconds() + S = A + B * delta_t / (3600*24) / 1000. + return S/1000 + def get_cds_time(days, msecs): """Compute timestamp given the days since epoch and milliseconds of the day. @@ -559,6 +615,11 @@ def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time): self._platform_id = platform_id self._channel_name = channel_name self._coefs = coefs + if channel_name in ['VIS006', 'VIS008', 'IR_016']: + self._coefs['coefs']['MEIRINK'] = MEIRINK_COEFS[platform_id][channel_name] + else: + self._coefs['coefs']['MEIRINK'] = None + self._calib_mode = calib_mode.upper() self._scan_time = scan_time self._algo = SEVIRICalibrationAlgorithm( @@ -566,7 +627,7 @@ def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time): scan_time=self._scan_time ) - valid_modes = ('NOMINAL', 'GSICS') + valid_modes = ('NOMINAL', 'GSICS', 'MEIRINK') if self._calib_mode not in valid_modes: raise ValueError( 'Invalid calibration mode: {}. Choose one of {}'.format( @@ -622,6 +683,10 @@ def get_gain_offset(self): internal_gain = gsics_gain internal_offset = gsics_offset + if self._calib_mode == 'MEIRINK': + if coefs['MEIRINK'] is not None: + internal_gain = get_meirink_slope(coefs['MEIRINK'], self._scan_time) + # Override with external coefficients, if any. gain = coefs['EXTERNAL'].get('gain', internal_gain) offset = coefs['EXTERNAL'].get('offset', internal_offset) diff --git a/satpy/readers/seviri_l1b_hrit.py b/satpy/readers/seviri_l1b_hrit.py index 2b153edfcc..4480afdbfb 100644 --- a/satpy/readers/seviri_l1b_hrit.py +++ b/satpy/readers/seviri_l1b_hrit.py @@ -247,6 +247,7 @@ mask_bad_quality, pad_data_horizontally, round_nom_time, + MEIRINK_CALIB, ) from satpy.readers.seviri_l1b_native_hdr import hrit_epilogue, hrit_prologue, impf_configuration from satpy.utils import get_legacy_chunk_size From 984a7331ecccabb5fd976d8f5899e7eb91ea5e7a Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Thu, 5 Oct 2023 15:00:56 +0200 Subject: [PATCH 02/18] add test for the slope of the Meirink coefficients --- satpy/tests/reader_tests/test_seviri_base.py | 32 ++++++++++++++++++++ 1 file changed, 32 insertions(+) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index 157ed88bbf..f7c8525ab7 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -37,6 +37,9 @@ pad_data_horizontally, pad_data_vertically, round_nom_time, + SEVIRICalibrationHandler, + MEIRINK_COEFS, + DATE_2000, ) from satpy.utils import get_legacy_chunk_size @@ -358,3 +361,32 @@ def test_get_orbit_polynomial_exceptions(self, orbit_polynomials, time): finder = OrbitPolynomialFinder(orbit_polynomials) with pytest.raises(NoValidOrbitParams): finder.get_orbit_polynomial(time=time) + + +class TestMeirinkSlope: + """Unit tests for the slope of Meirink calibration.""" + + @pytest.mark.parametrize('platform_id', [321, 322, 323, 324]) + @pytest.mark.parametrize('channel_name', ['VIS006', 'VIS008', 'IR_016']) + def test_get_meirink_slope_epoch(self, platform_id, channel_name): + """Test the value of the slope of the Meirink calibration on 2000-01-01.""" + coefs = {'coefs': {}} + coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} + coefs['coefs']['EXTERNAL'] = {} + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', DATE_2000) + assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS[platform_id][channel_name][0]/1000. + + @pytest.mark.parametrize('platform_id', [321, 322, 323, 324]) + @pytest.mark.parametrize('channel_name', ['VIS006', 'VIS008', 'IR_016']) + def test_get_meirink_slope_2020(self, platform_id, channel_name): + """Test the value of the slope of the Meirink calibration on 2020-01-01.""" + DATE_2020 = datetime(2020, 1, 1) + coefs = {'coefs': {}} + coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} + coefs['coefs']['EXTERNAL'] = {} + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', DATE_2020) + A, B = MEIRINK_COEFS[platform_id][channel_name] + delta_t = (DATE_2020 - DATE_2000).total_seconds() + S = A + B * delta_t / (3600*24) / 1000. + S = S/1000 + assert calibration_handler.get_gain_offset()[0] == S From a9521ef3f9af36a7bd7aa8686443649c8f2e12a5 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Thu, 5 Oct 2023 16:45:28 +0200 Subject: [PATCH 03/18] remove un-necessary import --- satpy/readers/seviri_l1b_hrit.py | 1 - 1 file changed, 1 deletion(-) diff --git a/satpy/readers/seviri_l1b_hrit.py b/satpy/readers/seviri_l1b_hrit.py index 4480afdbfb..2b153edfcc 100644 --- a/satpy/readers/seviri_l1b_hrit.py +++ b/satpy/readers/seviri_l1b_hrit.py @@ -247,7 +247,6 @@ mask_bad_quality, pad_data_horizontally, round_nom_time, - MEIRINK_CALIB, ) from satpy.readers.seviri_l1b_native_hdr import hrit_epilogue, hrit_prologue, impf_configuration from satpy.utils import get_legacy_chunk_size From 688eb795d67490655aacd1f1cad76b317c158b6d Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Fri, 6 Oct 2023 11:18:56 +0200 Subject: [PATCH 04/18] fix MEIRINK calibration coefficients (typos from copying the coefficients over) --- satpy/readers/seviri_base.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index c1d288f4df..c54a7586dc 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -373,9 +373,9 @@ # Meteosat-9 -MEIRINK_COEFS[322] = {'VIS006': (21.026, 0.3739), - 'VIS008': (26.875, 0.3111), - 'IR_016': (21.394, 0.0065) +MEIRINK_COEFS[322] = {'VIS006': (21.026, 0.2556), + 'VIS008': (26.875, 0.1835), + 'IR_016': (21.394, 0.0498) } # Meteosat-10 From 0abbb9a10376271c228839f53e59d70037aec413 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Fri, 6 Oct 2023 14:04:45 +0200 Subject: [PATCH 05/18] used fixed values to test meirink calibration --- satpy/tests/reader_tests/test_seviri_base.py | 26 +++++++++++--------- 1 file changed, 15 insertions(+), 11 deletions(-) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index f7c8525ab7..b924929d84 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -376,17 +376,21 @@ def test_get_meirink_slope_epoch(self, platform_id, channel_name): calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', DATE_2000) assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS[platform_id][channel_name][0]/1000. - @pytest.mark.parametrize('platform_id', [321, 322, 323, 324]) - @pytest.mark.parametrize('channel_name', ['VIS006', 'VIS008', 'IR_016']) - def test_get_meirink_slope_2020(self, platform_id, channel_name): - """Test the value of the slope of the Meirink calibration on 2020-01-01.""" - DATE_2020 = datetime(2020, 1, 1) + @pytest.mark.parametrize('platform_id,time,expected', ( + (321, datetime(2005, 1, 18, 0, 0), [0.0250354716, 0.0315626684, 0.022880986]), + (321, datetime(2010, 12, 31, 0, 0), [0.0258479563, 0.0322386887, 0.022895110500000003]), + (322, datetime(2010, 1, 18, 0, 0), [0.021964051999999998, 0.027548445, 0.021576766]), + (322, datetime(2015, 6, 1, 0, 0), [0.022465028, 0.027908105, 0.021674373999999996]), + (323, datetime(2005, 1, 18, 0, 0), [0.0209088464, 0.0265355228, 0.0230132616]), + (323, datetime(2010, 12, 31, 0, 0), [0.022181355200000002, 0.0280103379, 0.0229511138]), + (324, datetime(2010, 1, 18, 0, 0), [0.0218362, 0.027580748, 0.022285370999999998]), + (324, datetime(2015, 6, 1, 0, 0), [0.0225418, 0.028530172, 0.022248718999999997]), + )) + def test_get_meirink_slope_2020(self, platform_id, time, expected): + """Test the value of the slope of the Meirink calibration.""" coefs = {'coefs': {}} coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} coefs['coefs']['EXTERNAL'] = {} - calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', DATE_2020) - A, B = MEIRINK_COEFS[platform_id][channel_name] - delta_t = (DATE_2020 - DATE_2000).total_seconds() - S = A + B * delta_t / (3600*24) / 1000. - S = S/1000 - assert calibration_handler.get_gain_offset()[0] == S + for i, channel_name in enumerate(['VIS006', 'VIS008', 'IR_016']): + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', time) + assert abs(calibration_handler.get_gain_offset()[0] - expected[i]) < 1e-6 From 909fb05bac0235ef80b58b11646495c12f44540e Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Fri, 6 Oct 2023 14:24:41 +0200 Subject: [PATCH 06/18] isort --- satpy/tests/reader_tests/test_seviri_base.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index b924929d84..f5ba38cea5 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -26,9 +26,12 @@ import xarray as xr from satpy.readers.seviri_base import ( + DATE_2000, + MEIRINK_COEFS, NoValidOrbitParams, OrbitPolynomial, OrbitPolynomialFinder, + SEVIRICalibrationHandler, chebyshev, dec10216, get_cds_time, @@ -37,9 +40,6 @@ pad_data_horizontally, pad_data_vertically, round_nom_time, - SEVIRICalibrationHandler, - MEIRINK_COEFS, - DATE_2000, ) from satpy.utils import get_legacy_chunk_size From 5c2234223f9e9fb44739830ffb24ce1ac9cf520f Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Fri, 6 Oct 2023 15:22:46 +0200 Subject: [PATCH 07/18] fix "multiple spaces after ','" for codefactor --- satpy/tests/reader_tests/test_seviri_base.py | 16 ++++++++-------- 1 file changed, 8 insertions(+), 8 deletions(-) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index f5ba38cea5..252da43e75 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -377,14 +377,14 @@ def test_get_meirink_slope_epoch(self, platform_id, channel_name): assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS[platform_id][channel_name][0]/1000. @pytest.mark.parametrize('platform_id,time,expected', ( - (321, datetime(2005, 1, 18, 0, 0), [0.0250354716, 0.0315626684, 0.022880986]), - (321, datetime(2010, 12, 31, 0, 0), [0.0258479563, 0.0322386887, 0.022895110500000003]), - (322, datetime(2010, 1, 18, 0, 0), [0.021964051999999998, 0.027548445, 0.021576766]), - (322, datetime(2015, 6, 1, 0, 0), [0.022465028, 0.027908105, 0.021674373999999996]), - (323, datetime(2005, 1, 18, 0, 0), [0.0209088464, 0.0265355228, 0.0230132616]), - (323, datetime(2010, 12, 31, 0, 0), [0.022181355200000002, 0.0280103379, 0.0229511138]), - (324, datetime(2010, 1, 18, 0, 0), [0.0218362, 0.027580748, 0.022285370999999998]), - (324, datetime(2015, 6, 1, 0, 0), [0.0225418, 0.028530172, 0.022248718999999997]), + (321, datetime(2005, 1, 18, 0, 0), [0.0250354716, 0.0315626684, 0.022880986]), + (321, datetime(2010, 12, 31, 0, 0), [0.0258479563, 0.0322386887, 0.022895110500000003]), + (322, datetime(2010, 1, 18, 0, 0), [0.021964051999999998, 0.027548445, 0.021576766]), + (322, datetime(2015, 6, 1, 0, 0), [0.022465028, 0.027908105, 0.021674373999999996]), + (323, datetime(2005, 1, 18, 0, 0), [0.0209088464, 0.0265355228, 0.0230132616]), + (323, datetime(2010, 12, 31, 0, 0), [0.022181355200000002, 0.0280103379, 0.0229511138]), + (324, datetime(2010, 1, 18, 0, 0), [0.0218362, 0.027580748, 0.022285370999999998]), + (324, datetime(2015, 6, 1, 0, 0), [0.0225418, 0.028530172, 0.022248718999999997]), )) def test_get_meirink_slope_2020(self, platform_id, time, expected): """Test the value of the slope of the Meirink calibration.""" From 09efba8ef643ae2a0af9ba61778a9bf4c6912a5d Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Mon, 9 Oct 2023 14:18:58 +0200 Subject: [PATCH 08/18] add comment with the units of the Meirink coefficients --- satpy/readers/seviri_base.py | 5 ++++- 1 file changed, 4 insertions(+), 1 deletion(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index c54a7586dc..3f548aa8ea 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -356,7 +356,10 @@ # Calibration coefficients from Meirink, J.F., R.A. Roebeling and P. Stammes, 2013: # Inter-calibration of polar imager solar channels using SEVIRI, Atm. Meas. Tech., 6, # 2495-2508, doi:10.5194/amt-6-2495-2013 - +# +# The coefficients are stored in pairs of A, B (see function `get_meirink_slope`) where the +# units of A are µW m-2 sr-1 (cm-1)-1 and those of B are µW m-2 sr-1 (cm-1)-1 (86400 s)-1 +# # To obtain the slope for the calibration, one should use the routine get_seviri_meirink_slope # Epoch for the MEIRINK re-calibration From 8193298786f788c81b3fddc8b161679fb982c4d1 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Mon, 9 Oct 2023 16:11:57 +0200 Subject: [PATCH 09/18] Use MeirinkCalibrationHandler to manage the coefficients Co-authored-by: Stephan Finkensieper --- satpy/readers/seviri_base.py | 61 +++++++++++++++++++++--------------- 1 file changed, 36 insertions(+), 25 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 3f548aa8ea..394540cb8c 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -366,34 +366,35 @@ DATE_2000 = datetime(2000, 1, 1) MEIRINK_COEFS = {} +MEIRINK_COEFS['2013'] = {} # Meteosat-8 -MEIRINK_COEFS[321] = {'VIS006': (24.346, 0.3739), - 'VIS008': (30.989, 0.3111), - 'IR_016': (22.869, 0.0065) - } +MEIRINK_COEFS['2013'][321] = {'VIS006': (24.346, 0.3739), + 'VIS008': (30.989, 0.3111), + 'IR_016': (22.869, 0.0065) + } # Meteosat-9 -MEIRINK_COEFS[322] = {'VIS006': (21.026, 0.2556), - 'VIS008': (26.875, 0.1835), - 'IR_016': (21.394, 0.0498) - } +MEIRINK_COEFS['2013'][322] = {'VIS006': (21.026, 0.2556), + 'VIS008': (26.875, 0.1835), + 'IR_016': (21.394, 0.0498) + } # Meteosat-10 -MEIRINK_COEFS[323] = {'VIS006': (19.829, 0.5856), - 'VIS008': (25.284, 0.6787), - 'IR_016': (23.066, -0.0286) - } +MEIRINK_COEFS['2013'][323] = {'VIS006': (19.829, 0.5856), + 'VIS008': (25.284, 0.6787), + 'IR_016': (23.066, -0.0286) + } # Meteosat-11 -MEIRINK_COEFS[324] = {'VIS006': (20.515, 0.3600), - 'VIS008': (25.803, 0.4844), - 'IR_016': (22.354, -0.0187) - } +MEIRINK_COEFS['2013'][324] = {'VIS006': (20.515, 0.3600), + 'VIS008': (25.803, 0.4844), + 'IR_016': (22.354, -0.0187) + } def get_meirink_slope(meirink_coefs, acquisition_time): @@ -413,6 +414,21 @@ def get_meirink_slope(meirink_coefs, acquisition_time): return S/1000 +class MeirinkCalibrationHandler: + """Re-calibration of the SEVIRI visible channels slop (see Meirink 2013).""" + + def __init__(self, coefs=MEIRINK_COEFS, calib_mode=None): + """Initialize the calibration handler.""" + if calib_mode is None: + raise ValueError("Missing calib_mode") + self.coefs = coefs[calib_mode.split('-')[1]] + + def get_slope(self, platform, channel, time): + """Return the slope using the provided calibration coefficients.""" + coefs = self.coefs[platform][channel] + return get_meirink_slope(coefs, time) + + def get_cds_time(days, msecs): """Compute timestamp given the days since epoch and milliseconds of the day. @@ -618,11 +634,6 @@ def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time): self._platform_id = platform_id self._channel_name = channel_name self._coefs = coefs - if channel_name in ['VIS006', 'VIS008', 'IR_016']: - self._coefs['coefs']['MEIRINK'] = MEIRINK_COEFS[platform_id][channel_name] - else: - self._coefs['coefs']['MEIRINK'] = None - self._calib_mode = calib_mode.upper() self._scan_time = scan_time self._algo = SEVIRICalibrationAlgorithm( @@ -630,7 +641,7 @@ def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time): scan_time=self._scan_time ) - valid_modes = ('NOMINAL', 'GSICS', 'MEIRINK') + valid_modes = ('NOMINAL', 'GSICS', 'MEIRINK-2013') if self._calib_mode not in valid_modes: raise ValueError( 'Invalid calibration mode: {}. Choose one of {}'.format( @@ -686,9 +697,9 @@ def get_gain_offset(self): internal_gain = gsics_gain internal_offset = gsics_offset - if self._calib_mode == 'MEIRINK': - if coefs['MEIRINK'] is not None: - internal_gain = get_meirink_slope(coefs['MEIRINK'], self._scan_time) + if "MEIRINK" in self._calib_mode and self._channel_name in ['VIS006', 'VIS008', 'IR_016']: + meirink = MeirinkCalibrationHandler(calib_mode=self._calib_mode) + internal_gain = meirink.get_slope(self._platform_id, self._channel_name, self._scan_time) # Override with external coefficients, if any. gain = coefs['EXTERNAL'].get('gain', internal_gain) From ad1b690d01609b2b50692dda3cebca14de2adf12 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Mon, 9 Oct 2023 23:22:17 +0200 Subject: [PATCH 10/18] Move criterion for application of meirink calibration in method Also fix typo --- satpy/readers/seviri_base.py | 10 ++++++++-- 1 file changed, 8 insertions(+), 2 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 394540cb8c..f41b3cec40 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -414,8 +414,14 @@ def get_meirink_slope(meirink_coefs, acquisition_time): return S/1000 +def should_apply_meirink(calib_mode, channel_name): + """Decide whether to use the Meirink calibration coefficients.""" + + return "MEIRINK" in calib_mode and channel_name in ['VIS006', 'VIS008', 'IR_016'] + + class MeirinkCalibrationHandler: - """Re-calibration of the SEVIRI visible channels slop (see Meirink 2013).""" + """Re-calibration of the SEVIRI visible channels slope (see Meirink 2013).""" def __init__(self, coefs=MEIRINK_COEFS, calib_mode=None): """Initialize the calibration handler.""" @@ -697,7 +703,7 @@ def get_gain_offset(self): internal_gain = gsics_gain internal_offset = gsics_offset - if "MEIRINK" in self._calib_mode and self._channel_name in ['VIS006', 'VIS008', 'IR_016']: + if should_apply_meirink(self._calib_mode, self._channel_name): meirink = MeirinkCalibrationHandler(calib_mode=self._calib_mode) internal_gain = meirink.get_slope(self._platform_id, self._channel_name, self._scan_time) From 65503a00ae7cfa7492f811fd6eac7a5de2012676 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Tue, 10 Oct 2023 10:12:56 +0200 Subject: [PATCH 11/18] flake8 --- satpy/readers/seviri_base.py | 1 - 1 file changed, 1 deletion(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index f41b3cec40..bd4ce962ad 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -416,7 +416,6 @@ def get_meirink_slope(meirink_coefs, acquisition_time): def should_apply_meirink(calib_mode, channel_name): """Decide whether to use the Meirink calibration coefficients.""" - return "MEIRINK" in calib_mode and channel_name in ['VIS006', 'VIS008', 'IR_016'] From dbefc0fbe67270ba04a42aae6f71426c95a010b6 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Tue, 10 Oct 2023 16:00:08 +0200 Subject: [PATCH 12/18] rename DATE_2000 to MEIRINK_EPOCH --- satpy/readers/seviri_base.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index bd4ce962ad..88994d3f18 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -363,7 +363,7 @@ # To obtain the slope for the calibration, one should use the routine get_seviri_meirink_slope # Epoch for the MEIRINK re-calibration -DATE_2000 = datetime(2000, 1, 1) +MEIRINK_EPOCH = datetime(2000, 1, 1) MEIRINK_COEFS = {} MEIRINK_COEFS['2013'] = {} @@ -409,7 +409,7 @@ def get_meirink_slope(meirink_coefs, acquisition_time): """ A = meirink_coefs[0] B = meirink_coefs[1] - delta_t = (acquisition_time - DATE_2000).total_seconds() + delta_t = (acquisition_time - MEIRINK_EPOCH).total_seconds() S = A + B * delta_t / (3600*24) / 1000. return S/1000 From 60bb30f76fb13fa95b957e6168a1ec1e1c938403 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Tue, 10 Oct 2023 16:01:26 +0200 Subject: [PATCH 13/18] remove coef argument in MeirinkCalibrationHandler --- satpy/readers/seviri_base.py | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 88994d3f18..6471c42639 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -422,11 +422,9 @@ def should_apply_meirink(calib_mode, channel_name): class MeirinkCalibrationHandler: """Re-calibration of the SEVIRI visible channels slope (see Meirink 2013).""" - def __init__(self, coefs=MEIRINK_COEFS, calib_mode=None): + def __init__(self, calib_mode): """Initialize the calibration handler.""" - if calib_mode is None: - raise ValueError("Missing calib_mode") - self.coefs = coefs[calib_mode.split('-')[1]] + self.coefs = MEIRINK_COEFS[calib_mode.split('-')[1]] def get_slope(self, platform, channel, time): """Return the slope using the provided calibration coefficients.""" From 69c5751dacd76501569ddef56cc37bf7c52faa33 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Tue, 10 Oct 2023 22:13:01 +0200 Subject: [PATCH 14/18] add documentation about Meirink calibration --- satpy/readers/seviri_base.py | 34 +++++++++++++++++++++++++++------- 1 file changed, 27 insertions(+), 7 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 6471c42639..73a1454cb8 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -61,14 +61,22 @@ reader_kwargs={'calib_mode': 'GSICS'}) scene.load(['VIS006', 'IR_108']) -Furthermore, it is possible to specify external calibration coefficients -for the conversion from counts to radiances. External coefficients take -precedence over internal coefficients, but you can also mix internal and -external coefficients: If external calibration coefficients are specified -for only a subset of channels, the remaining channels will be calibrated -using the chosen file-internal coefficients (nominal or GSICS). +In addition, two other calibration methods are available: + +1. It is possible to specify external calibration coefficients for the + conversion from counts to radiances. External coefficients take + precedence over internal coefficients and over the Meirink + coefficients, but you can also mix internal and external coefficients: + If external calibration coefficients are specified for only a subset + of channels, the remaining channels will be calibrated using the + chosen file-internal coefficients (nominal or GSICS). Calibration + coefficients must be specified in [mW m-2 sr-1 (cm-1)-1]. + +2. The calibration mode ``meirink-2013`` uses coefficients based on an + intercalibration with Aqua-MODIS for the visible channels, as found in + `Inter-calibration of polar imager solar channels using SEVIRI`_ + (2013) by J. F. Meirink, R. A. Roebeling, and P. Stammes. -Calibration coefficients must be specified in [mW m-2 sr-1 (cm-1)-1]. In the following example we use external calibration coefficients for the ``VIS006`` & ``IR_108`` channels, and nominal coefficients for the @@ -93,6 +101,15 @@ 'ext_calib_coefs': coefs}) scene.load(['VIS006', 'VIS008', 'IR_108', 'IR_120']) +In the next example we use the mode ``meirink-2013`` calibration +coefficients for all visible channels and nominal coefficients for the +rest:: + + scene = satpy.Scene(filenames, + reader='seviri_l1b_...', + reader_kwargs={'calib_mode': 'meirink-2013'}) + scene.load(['VIS006', 'VIS008', 'IR_016']) + Calibration to reflectance ^^^^^^^^^^^^^^^^^^^^^^^^^^ @@ -163,6 +180,9 @@ .. _Radiometric Calibration of MSG SEVIRI Level 1.5 Image Data in Equivalent Spectral Blackbody Radiance: https://www-cdn.eumetsat.int/files/2020-04/pdf_ten_msg_seviri_rad_calib.pdf +.. _Inter-calibration of polar imager solar channels using SEVIRI: + http://dx.doi.org/10.5194/amt-6-2495-2013 + """ import warnings From 714130cba99751d3b604e9ff8252662b9334dddb Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Wed, 11 Oct 2023 09:59:30 +0200 Subject: [PATCH 15/18] Update test for changes in Meirink calibration --- satpy/tests/reader_tests/test_seviri_base.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index 252da43e75..3b357512e4 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -26,7 +26,7 @@ import xarray as xr from satpy.readers.seviri_base import ( - DATE_2000, + MEIRINK_EPOCH, MEIRINK_COEFS, NoValidOrbitParams, OrbitPolynomial, @@ -373,8 +373,8 @@ def test_get_meirink_slope_epoch(self, platform_id, channel_name): coefs = {'coefs': {}} coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} coefs['coefs']['EXTERNAL'] = {} - calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', DATE_2000) - assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS[platform_id][channel_name][0]/1000. + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2013', MEIRINK_EPOCH) + assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS['2013'][platform_id][channel_name][0]/1000. @pytest.mark.parametrize('platform_id,time,expected', ( (321, datetime(2005, 1, 18, 0, 0), [0.0250354716, 0.0315626684, 0.022880986]), @@ -392,5 +392,5 @@ def test_get_meirink_slope_2020(self, platform_id, time, expected): coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} coefs['coefs']['EXTERNAL'] = {} for i, channel_name in enumerate(['VIS006', 'VIS008', 'IR_016']): - calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK', time) + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2013', time) assert abs(calibration_handler.get_gain_offset()[0] - expected[i]) < 1e-6 From 51a3331206c158cdf464eab219da121da34cd869 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Wed, 11 Oct 2023 11:56:39 +0200 Subject: [PATCH 16/18] isort --- satpy/tests/reader_tests/test_seviri_base.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index 3b357512e4..910623bcd5 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -26,8 +26,8 @@ import xarray as xr from satpy.readers.seviri_base import ( - MEIRINK_EPOCH, MEIRINK_COEFS, + MEIRINK_EPOCH, NoValidOrbitParams, OrbitPolynomial, OrbitPolynomialFinder, From 3b465a0dcbac71f3e007a8e825e1b0c1aa5b6910 Mon Sep 17 00:00:00 2001 From: Pierre de Buyl Date: Wed, 11 Oct 2023 14:06:12 +0200 Subject: [PATCH 17/18] Change version of Meirink coefficients to 2023. The coefficients on the webpage https://msgcpp.knmi.nl/solar-channel-calibration.html were updated in place. The current set of coefficients were obtained in 2023, the code now reflects this correctly. --- satpy/readers/seviri_base.py | 21 +++++++++++--------- satpy/tests/reader_tests/test_seviri_base.py | 6 +++--- 2 files changed, 15 insertions(+), 12 deletions(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 73a1454cb8..83977bf40c 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -72,7 +72,7 @@ chosen file-internal coefficients (nominal or GSICS). Calibration coefficients must be specified in [mW m-2 sr-1 (cm-1)-1]. -2. The calibration mode ``meirink-2013`` uses coefficients based on an +2. The calibration mode ``meirink-2023`` uses coefficients based on an intercalibration with Aqua-MODIS for the visible channels, as found in `Inter-calibration of polar imager solar channels using SEVIRI`_ (2013) by J. F. Meirink, R. A. Roebeling, and P. Stammes. @@ -101,13 +101,13 @@ 'ext_calib_coefs': coefs}) scene.load(['VIS006', 'VIS008', 'IR_108', 'IR_120']) -In the next example we use the mode ``meirink-2013`` calibration +In the next example we use the mode ``meirink-2023`` calibration coefficients for all visible channels and nominal coefficients for the rest:: scene = satpy.Scene(filenames, reader='seviri_l1b_...', - reader_kwargs={'calib_mode': 'meirink-2013'}) + reader_kwargs={'calib_mode': 'meirink-2023'}) scene.load(['VIS006', 'VIS008', 'IR_016']) @@ -377,6 +377,9 @@ # Inter-calibration of polar imager solar channels using SEVIRI, Atm. Meas. Tech., 6, # 2495-2508, doi:10.5194/amt-6-2495-2013 # +# The coeffients in the 2023 entry have been obtained from the webpage +# https://msgcpp.knmi.nl/solar-channel-calibration.html on 2023-10-11. +# # The coefficients are stored in pairs of A, B (see function `get_meirink_slope`) where the # units of A are µW m-2 sr-1 (cm-1)-1 and those of B are µW m-2 sr-1 (cm-1)-1 (86400 s)-1 # @@ -386,32 +389,32 @@ MEIRINK_EPOCH = datetime(2000, 1, 1) MEIRINK_COEFS = {} -MEIRINK_COEFS['2013'] = {} +MEIRINK_COEFS['2023'] = {} # Meteosat-8 -MEIRINK_COEFS['2013'][321] = {'VIS006': (24.346, 0.3739), +MEIRINK_COEFS['2023'][321] = {'VIS006': (24.346, 0.3739), 'VIS008': (30.989, 0.3111), 'IR_016': (22.869, 0.0065) } # Meteosat-9 -MEIRINK_COEFS['2013'][322] = {'VIS006': (21.026, 0.2556), +MEIRINK_COEFS['2023'][322] = {'VIS006': (21.026, 0.2556), 'VIS008': (26.875, 0.1835), 'IR_016': (21.394, 0.0498) } # Meteosat-10 -MEIRINK_COEFS['2013'][323] = {'VIS006': (19.829, 0.5856), +MEIRINK_COEFS['2023'][323] = {'VIS006': (19.829, 0.5856), 'VIS008': (25.284, 0.6787), 'IR_016': (23.066, -0.0286) } # Meteosat-11 -MEIRINK_COEFS['2013'][324] = {'VIS006': (20.515, 0.3600), +MEIRINK_COEFS['2023'][324] = {'VIS006': (20.515, 0.3600), 'VIS008': (25.803, 0.4844), 'IR_016': (22.354, -0.0187) } @@ -664,7 +667,7 @@ def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time): scan_time=self._scan_time ) - valid_modes = ('NOMINAL', 'GSICS', 'MEIRINK-2013') + valid_modes = ('NOMINAL', 'GSICS', 'MEIRINK-2023') if self._calib_mode not in valid_modes: raise ValueError( 'Invalid calibration mode: {}. Choose one of {}'.format( diff --git a/satpy/tests/reader_tests/test_seviri_base.py b/satpy/tests/reader_tests/test_seviri_base.py index 910623bcd5..32918ea45b 100644 --- a/satpy/tests/reader_tests/test_seviri_base.py +++ b/satpy/tests/reader_tests/test_seviri_base.py @@ -373,8 +373,8 @@ def test_get_meirink_slope_epoch(self, platform_id, channel_name): coefs = {'coefs': {}} coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} coefs['coefs']['EXTERNAL'] = {} - calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2013', MEIRINK_EPOCH) - assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS['2013'][platform_id][channel_name][0]/1000. + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2023', MEIRINK_EPOCH) + assert calibration_handler.get_gain_offset()[0] == MEIRINK_COEFS['2023'][platform_id][channel_name][0]/1000. @pytest.mark.parametrize('platform_id,time,expected', ( (321, datetime(2005, 1, 18, 0, 0), [0.0250354716, 0.0315626684, 0.022880986]), @@ -392,5 +392,5 @@ def test_get_meirink_slope_2020(self, platform_id, time, expected): coefs['coefs']['NOMINAL'] = {'gain': -1, 'offset': -1} coefs['coefs']['EXTERNAL'] = {} for i, channel_name in enumerate(['VIS006', 'VIS008', 'IR_016']): - calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2013', time) + calibration_handler = SEVIRICalibrationHandler(platform_id, channel_name, coefs, 'MEIRINK-2023', time) assert abs(calibration_handler.get_gain_offset()[0] - expected[i]) < 1e-6 From 056a3a9da5c8de9d811655441c234f735482f1a8 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Wed, 11 Oct 2023 10:47:00 -0500 Subject: [PATCH 18/18] Fix type annotations in seviri_base.py --- satpy/readers/seviri_base.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/satpy/readers/seviri_base.py b/satpy/readers/seviri_base.py index 83977bf40c..5e6d69ea68 100644 --- a/satpy/readers/seviri_base.py +++ b/satpy/readers/seviri_base.py @@ -184,6 +184,7 @@ http://dx.doi.org/10.5194/amt-6-2495-2013 """ +from __future__ import annotations import warnings from datetime import datetime, timedelta @@ -388,7 +389,7 @@ # Epoch for the MEIRINK re-calibration MEIRINK_EPOCH = datetime(2000, 1, 1) -MEIRINK_COEFS = {} +MEIRINK_COEFS: dict[str, dict[int, dict[str, tuple[float, float]]]] = {} MEIRINK_COEFS['2023'] = {} # Meteosat-8