Skip to content

Commit

Permalink
Move calibration algorithm to separate class
Browse files Browse the repository at this point in the history
  • Loading branch information
@sfinkens
sfinkens committed Dec 1, 2020
1 parent afb3d90 commit 6f51d30
Show file tree
Hide file tree
Showing 6 changed files with 280 additions and 354 deletions.
117 changes: 69 additions & 48 deletions satpy/readers/seviri_base.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -332,54 +332,15 @@ def images_used(self):
mpef_product_header = MpefProductHeader().get()


class SEVIRICalibrationHandler:
"""Calibration handler for SEVIRI HRIT-, native- and netCDF-formats."""
class SEVIRICalibrationAlgorithm:
"""SEVIRI calibration algorithms."""

def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time):
"""Initialize the calibration handler."""
def __init__(self, platform_id, scan_time):
"""Initialize the calibration algorithm."""
self.platform_id = platform_id
self.channel_name = channel_name
self.coefs = coefs
self.calib_mode = calib_mode.upper()
self.scan_time = scan_time

valid_modes = ('NOMINAL', 'GSICS')
if self.calib_mode not in valid_modes:
raise ValueError(
'Invalid calibration mode: {}. Choose one of {}'.format(
self.calib_mode, valid_modes)
)

def calibrate(self, data, calibration):
"""Calibrate the given data."""
if calibration == 'counts':
res = data
elif calibration in ['radiance', 'reflectance',
'brightness_temperature']:
# Convert to radiance
gain, offset = self._get_gain_offset()
data = data.where(data > 0)
res = self._convert_to_radiance(
data.astype(np.float32), gain, offset
)
else:
raise ValueError(
'Invalid calibration {} for channel {}'.format(
calibration, self.channel_name
)
)

if calibration == 'reflectance':
solar_irradiance = CALIB[self.platform_id][self.channel_name]["F"]
res = self._vis_calibrate(res, solar_irradiance)
elif calibration == 'brightness_temperature':
res = self._ir_calibrate(
res, self.channel_name, self.coefs['radiance_type']
)

return res

def _convert_to_radiance(self, data, gain, offset):
def convert_to_radiance(self, data, gain, offset):
"""Calibrate to radiance."""
return (data * gain + offset).clip(0.0, None)

Expand All @@ -392,7 +353,7 @@ def _erads2bt(self, data, channel_name):

return (self._tl15(data, wavenumber) - beta) / alpha

def _ir_calibrate(self, data, channel_name, cal_type):
def ir_calibrate(self, data, channel_name, cal_type):
"""Calibrate to brightness temperature."""
if cal_type == 1:
# spectral radiances
Expand All @@ -416,7 +377,7 @@ def _tl15(self, data, wavenumber):
return ((C2 * wavenumber) /
np.log((1.0 / data) * C1 * wavenumber ** 3 + 1.0))

def _vis_calibrate(self, data, solar_irradiance):
def vis_calibrate(self, data, solar_irradiance):
"""Calibrate to reflectance.
This uses the method described in Conversion from radiances to
Expand All @@ -425,8 +386,68 @@ def _vis_calibrate(self, data, solar_irradiance):
reflectance = np.pi * data * 100.0 / solar_irradiance
return apply_earthsun_distance_correction(reflectance, self.scan_time)

def _get_gain_offset(self):
"""Get gain & offset for calibration from counts to radiance."""

class SEVIRICalibrationHandler:
"""Calibration handler for SEVIRI HRIT-, native- and netCDF-formats.
Handles selection of calibration coefficients and calls the appropriate
calibration algorithm.
"""

def __init__(self, platform_id, channel_name, coefs, calib_mode, scan_time):
"""Initialize the calibration handler."""
self.platform_id = platform_id
self.channel_name = channel_name
self.coefs = coefs
self.calib_mode = calib_mode.upper()
self.scan_time = scan_time
self.algo = SEVIRICalibrationAlgorithm(
platform_id=self.platform_id,
scan_time=self.scan_time
)

valid_modes = ('NOMINAL', 'GSICS')
if self.calib_mode not in valid_modes:
raise ValueError(
'Invalid calibration mode: {}. Choose one of {}'.format(
self.calib_mode, valid_modes)
)

def calibrate(self, data, calibration):
"""Calibrate the given data."""
if calibration == 'counts':
res = data
elif calibration in ['radiance', 'reflectance',
'brightness_temperature']:
# Convert to radiance
gain, offset = self.get_gain_offset()
data = data.where(data > 0)
res = self.algo.convert_to_radiance(
data.astype(np.float32), gain, offset
)
else:
raise ValueError(
'Invalid calibration {} for channel {}'.format(
calibration, self.channel_name
)
)

if calibration == 'reflectance':
solar_irradiance = CALIB[self.platform_id][self.channel_name]["F"]
res = self.algo.vis_calibrate(res, solar_irradiance)
elif calibration == 'brightness_temperature':
res = self.algo.ir_calibrate(
res, self.channel_name, self.coefs['radiance_type']
)

return res

def get_gain_offset(self):
"""Get gain & offset for calibration from counts to radiance.
Choices for internal coefficients are nominal or GSICS. External
coefficients take precedence.
"""
coefs = self.coefs['coefs']

# Select internal coefficients for the given calibration mode
Expand Down
197 changes: 1 addition & 196 deletions satpy/tests/reader_tests/test_seviri_base.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -17,16 +17,11 @@
# satpy. If not, see <http://www.gnu.org/licenses/>.
"""Test the MSG common (native and hrit format) functionionalities."""

from datetime import datetime
import unittest

import numpy as np
import xarray as xr
import pytest

from satpy.readers.seviri_base import (
dec10216, chebyshev, get_cds_time, SEVIRICalibrationHandler
)
from satpy.readers.seviri_base import dec10216, chebyshev, get_cds_time


def chebyshev4(c, x, domain):
Expand Down Expand Up @@ -71,193 +66,3 @@ def test_get_cds_time(self):
msecs = 12*3600*1000
expected = np.datetime64('2016-03-03 12:00:00.000')
np.testing.assert_equal(get_cds_time(days=days, msecs=msecs), expected)


class TestCalibrationBase:
"""Base class for calibration tests."""

platform_id = 324
gains_nominal = np.arange(1, 13)
offsets_nominal = np.arange(-1, -13, -1)
gains_gsics = np.arange(0.1, 1.3, 0.1)
offsets_gsics = np.arange(-0.1, -1.3, -0.1)
radiance_types = 2 * np.ones(12)
scan_time = datetime(2020, 1, 1)
external_coefs = {
'VIS006': {'gain': 10, 'offset': -10},
'IR_108': {'gain': 20, 'offset': -20}
}
spectral_channel_ids = {'VIS006': 1, 'IR_108': 9}
expected = {
'VIS006': {
'counts': {
'NOMINAL': xr.DataArray(
[[0, 10],
[100, 255]],
dims=('y', 'x')
)
},
'radiance': {
'NOMINAL': xr.DataArray(
[[np.nan, 9],
[99, 254]],
dims=('y', 'x')
),
'GSICS': xr.DataArray(
[[np.nan, 9],
[99, 254]],
dims=('y', 'x')
),
'EXTERNAL': xr.DataArray(
[[np.nan, 90],
[990, 2540]],
dims=('y', 'x')
)
},
'reflectance': {
'NOMINAL': xr.DataArray(
[[np.nan, 40.47923],
[445.27155, 1142.414]],
dims=('y', 'x')
),
'EXTERNAL': xr.DataArray(
[[np.nan, 404.7923],
[4452.7153, 11424.14]],
dims=('y', 'x')
)
}
},
'IR_108': {
'counts': {
'NOMINAL': xr.DataArray(
[[0, 10],
[100, 255]],
dims=('y', 'x')
)
},
'radiance': {
'NOMINAL': xr.DataArray(
[[np.nan, 81],
[891, 2286]],
dims=('y', 'x')
),
'GSICS': xr.DataArray(
[[np.nan, 8.19],
[89.19, 228.69]],
dims=('y', 'x')
),
'EXTERNAL': xr.DataArray(
[[np.nan, 180],
[1980, 5080]],
dims=('y', 'x')
)
},
'brightness_temperature': {
'NOMINAL': xr.DataArray(
[[np.nan, 279.82318],
[543.2585, 812.77167]],
dims=('y', 'x')
),
'GSICS': xr.DataArray(
[[np.nan, 189.20985],
[285.53293, 356.06668]],
dims=('y', 'x')
),
'EXTERNAL': xr.DataArray(
[[np.nan, 335.14236],
[758.6249, 1262.7567]],
dims=('y', 'x')
),
}
}
}

@pytest.fixture(name='counts')
def counts(self):
"""Provide fake image counts."""
return xr.DataArray(
[[0, 10],
[100, 255]],
dims=('y', 'x')
)

def _get_expected(
self, channel, calibration, calib_mode, use_ext_coefs
):
if use_ext_coefs:
return self.expected[channel][calibration]['EXTERNAL']
return self.expected[channel][calibration][calib_mode]


class TestSeviriCalibrationHandler:
"""Unit tests for calibration handler."""

def test_init(self):
"""Test initialization of the calibration handler."""
with pytest.raises(ValueError):
SEVIRICalibrationHandler(
platform_id=None,
channel_name=None,
coefs=None,
calib_mode='invalid',
scan_time=None
)

@pytest.fixture(name='counts')
def counts(self):
"""Provide fake counts."""
return xr.DataArray(
[[1, 2],
[3, 4]],
dims=('y', 'x')
)

@pytest.fixture(name='calib')
def calib(self):
"""Provide a calibration handler."""
return SEVIRICalibrationHandler(
platform_id=324,
channel_name='IR_108',
coefs={
'coefs': {
'NOMINAL': {
'gain': 10,
'offset': -1
},
'GSICS': {
'gain': 20,
'offset': -2
},
'EXTERNAL': {}
},
'radiance_type': 1
},
calib_mode='NOMINAL',
scan_time=None
)

def test_calibrate_exceptions(self, counts, calib):
"""Test exception raised by the calibration handler."""
with pytest.raises(ValueError):
# Invalid calibration
calib.calibrate(counts, 'invalid')
with pytest.raises(NotImplementedError):
# Invalid radiance type
calib.coefs['radiance_type'] = 999
calib.calibrate(counts, 'brightness_temperature')

@pytest.mark.parametrize(
('calib_mode', 'ext_coefs', 'expected'),
[
('NOMINAL', {}, (10, -1)),
('GSICS', {}, (20, -40)),
('GSICS', {'gain': 30, 'offset': -3}, (30, -3)),
('NOMINAL', {'gain': 30, 'offset': -3}, (30, -3))
]
)
def test_get_gain_offset(self, calib, calib_mode, ext_coefs, expected):
"""Test selection of gain and offset."""
calib.calib_mode = calib_mode
calib.coefs['coefs']['EXTERNAL'] = ext_coefs
coefs = calib._get_gain_offset()
assert coefs == expected

0 comments on commit 6f51d30

Please sign in to comment.