Skip to content

Commit

Permalink
Merge pull request #2 from simonrp84/feature-ami-reader
Browse files Browse the repository at this point in the history 
Use internal calibration coefficients if available
  • Loading branch information
@djhoese
djhoese authored Oct 11, 2019
2 parents e6872b8 + 6088e5f commit 55798fb
Show file tree
Hide file tree
Showing 2 changed files with 78 additions and 54 deletions.
94 changes: 51 additions & 43 deletions satpy/readers/ami_l1b.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,36 +37,12 @@
'GK-2B': 'GEO-KOMPSAT-2B',
}

# Copied from 20190415_GK-2A_AMI_Conversion_Table_v3.0.xlsx
# Sheet: coeff.& equation_WN
# Visible channels
# channel_name -> (DN2Rad_Gain, DN2Rad_Offset, Rad. to Albedo)
# IR channels
# channel_name -> (DN2Rad_Gain, DN2Rad_Offset, c0, c1, c2)
CALIBRATION_COEFFS = {
"VI004": (0.363545805215835, -7.27090454101562, 0.001558245),
"VI005": (0.343625485897064, -6.87249755859375, 0.0016595767),
"VI006": (0.154856294393539, -6.19424438476562, 0.001924484),
"VI008": (0.0457241721451282, -3.65792846679687, 0.0032723873),
"NR013": (0.0346878096461296, -1.38751220703125, 0.0087081313),
"NR016": (0.0498007982969284, -0.996017456054687, 0.0129512876),
"SW038": (-0.00108296517282724, 17.699987411499, -0.447843939824124, 1.00065568090389, -0.0000000633824089912448),
"WV063": (-0.0108914673328399, 44.1777038574218, -1.76279494011147, 1.00414910562278, -0.000000983310914319385),
"WV069": (-0.00818779878318309, 66.7480773925781, -0.334311414359106, 1.00097359874468, -0.000000494603070252304),
"WV073": (-0.0096982717514038, 79.0608520507812, -0.0613124859696595, 1.00019008722941, -0.000000105863656750499),
"IR087": (-0.0144806550815701, 118.050903320312, -0.141418528203155, 1.00052232906885, -0.00000036287276076109),
"IR096": (-0.0178435463458299, 145.464874267578, -0.114017728158198, 1.00047380585402, -0.000000374931509928403),
"IR105": (-0.0198196955025196, 161.580139160156, -0.142866448475177, 1.00064069572049, -0.000000550443294960498),
"IR112": (-0.0216744858771562, 176.713439941406, -0.249111718496148, 1.00121166873756, -0.00000113167964011665),
"IR123": (-0.023379972204566, 190.649627685546, -0.458113885722738, 1.00245520975535, -0.00000253064314720476),
"IR133": (-0.0243037566542625, 198.224365234375, -0.0938521568527657, 1.00053982112966, -0.000000594913715312849),
}


class AMIL1bNetCDF(BaseFileHandler):
"""Base reader for AMI L1B NetCDF4 files."""

def __init__(self, filename, filename_info, filetype_info, allow_conditional_pixels=False):
def __init__(self, filename, filename_info, filetype_info,
calib_mode='PYSPECTRAL', allow_conditional_pixels=False):
"""Open the NetCDF file with xarray and prepare the Dataset for reading."""
super(AMIL1bNetCDF, self).__init__(filename, filename_info, filetype_info)
self.nc = xr.open_dataset(self.filename,
Expand All @@ -79,6 +55,11 @@ def __init__(self, filename, filename_info, filetype_info, allow_conditional_pix
self.platform_name = PLATFORM_NAMES.get(platform_shortname)
self.sensor = 'ami'
self.allow_conditional_pixels = allow_conditional_pixels
calib_mode_choices = ('FILE', 'PYSPECTRAL')
if calib_mode.upper() not in calib_mode_choices:
raise ValueError('Invalid calibration mode: {}. Choose one of {}'.format(
calib_mode, calib_mode_choices))
self.calib_mode = calib_mode.upper()

@property
def start_time(self):
Expand Down Expand Up @@ -182,21 +163,35 @@ def get_dataset(self, dataset_id, ds_info):
# only take "no error" pixels as valid
data = data.where(qf == 0)

channel_name = attrs.get('channel_name', dataset_id.name)
coeffs = CALIBRATION_COEFFS.get(channel_name)
if coeffs is None and dataset_id.calibration is not None:
raise ValueError("No coefficients configured for {}".format(dataset_id))
# Calibration values from file, fall back to built-in if unavailable
gain = self.nc.attrs['DN_to_Radiance_Gain']
offset = self.nc.attrs['DN_to_Radiance_Offset']

if dataset_id.calibration in ('radiance', 'reflectance', 'brightness_temperature'):
gain = coeffs[0]
offset = coeffs[1]
data = gain * data + offset
if dataset_id.calibration == 'reflectance':
# depends on the radiance calibration above
rad_to_alb = coeffs[2]
rad_to_alb = self.nc.attrs['Radiance_to_Albedo_c']
if ds_info.get('units') == '%':
rad_to_alb *= 100
data = data * rad_to_alb
elif dataset_id.calibration == 'brightness_temperature':
data = self._calibrate_ir(dataset_id, data)
elif dataset_id.calibration not in ('counts', 'radiance'):
raise ValueError("Unknown calibration: '{}'".format(dataset_id.calibration))

for attr_name in ('standard_name', 'units'):
attrs[attr_name] = ds_info[attr_name]
attrs.update(dataset_id.to_dict())
attrs['orbital_parameters'] = self.get_orbital_parameters()
attrs['platform_name'] = self.platform_name
attrs['sensor'] = self.sensor
data.attrs = attrs
return data

def _calibrate_ir(self, dataset_id, data):
"""Calibrate radiance data to BTs using either pyspectral or in-file coefficients."""
if self.calib_mode == 'PYSPECTRAL':
# depends on the radiance calibration above
# Convert um to m^-1 (SI units for pyspectral)
wn = 1 / (dataset_id.wavelength[1] / 1e6)
Expand All @@ -209,14 +204,27 @@ def get_dataset(self, dataset_id, ds_info):
else:
# new versions of pyspectral can do dask arrays
data.data = bt_data
elif dataset_id.calibration not in ('counts', 'radiance'):
raise ValueError("Unknown calibration: '{}'".format(dataset_id.calibration))

for attr_name in ('standard_name', 'units'):
attrs[attr_name] = ds_info[attr_name]
attrs.update(dataset_id.to_dict())
attrs['orbital_parameters'] = self.get_orbital_parameters()
attrs['platform_name'] = self.platform_name
attrs['sensor'] = self.sensor
data.attrs = attrs
else:
# IR coefficients from the file
# Channel specific
c0 = self.nc.attrs['Teff_to_Tbb_c0']
c1 = self.nc.attrs['Teff_to_Tbb_c1']
c2 = self.nc.attrs['Teff_to_Tbb_c2']

# These should be fixed, but load anyway
cval = self.nc.attrs['light_speed']
kval = self.nc.attrs['Boltzmann_constant_k']
hval = self.nc.attrs['Plank_constant_h']

# Compute wavenumber as cm-1
wn = (10000 / dataset_id.wavelength[1]) * 100

# Convert radiance to effective brightness temperature
e1 = (2 * hval * cval * cval) * np.power(wn, 3)
e2 = (data.data * 1e-5)
t_eff = ((hval * cval / kval) * wn) / np.log((e1 / e2) + 1)

# Now convert to actual brightness temperature
bt_data = c0 + c1 * t_eff + c2 * t_eff * t_eff
data.data = bt_data
return data
38 changes: 27 additions & 11 deletions satpy/tests/reader_tests/test_ami_l1b.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -115,6 +115,15 @@ def setUp(self, xr_, counts=None):
"number_of_lines": 22000,
"observation_mode": "FD",
"channel_spatial_resolution": "0.5",
"Radiance_to_Albedo_c": 1,
"DN_to_Radiance_Gain": -0.0144806550815701,
"DN_to_Radiance_Offset": 118.050903320312,
"Teff_to_Tbb_c0": -0.141418528203155,
"Teff_to_Tbb_c1": 1.00052232906885,
"Teff_to_Tbb_c2": -0.00000036287276076109,
"light_speed": 2.9979245800E+08,
"Boltzmann_constant_k": 1.3806488000E-23,
"Plank_constant_h": 6.6260695700E-34,
}
)

Expand Down Expand Up @@ -262,19 +271,26 @@ def setUp(self):
def test_ir_calibrate(self):
"""Test IR calibration."""
from satpy import DatasetID
res = self.reader.get_dataset(
DatasetID(name='IR087', wavelength=[8.415, 8.59, 8.765],
calibration='brightness_temperature'),
{
'file_key': 'image_pixel_values',
'wavelength': [8.415, 8.59, 8.765],
'standard_name': 'toa_brightness_temperature',
'units': 'K',
})

ds_id = DatasetID(name='IR087', wavelength=[8.415, 8.59, 8.765],
calibration='brightness_temperature')
ds_info = {
'file_key': 'image_pixel_values',
'wavelength': [8.415, 8.59, 8.765],
'standard_name': 'toa_brightness_temperature',
'units': 'K',
}
res = self.reader.get_dataset(ds_id, ds_info)
expected = np.array([[238.34385135, 238.31443527, 238.28500087, 238.25554813, 238.22607701],
[238.1965875, 238.16707956, 238.13755317, 238.10800829, 238.07844489]])
self.assertTrue(np.allclose(res.data.compute(), expected, equal_nan=True))
np.testing.assert_allclose(res.data.compute(), expected, equal_nan=True)
# make sure the attributes from the file are in the data array
self.assertEqual(res.attrs['standard_name'], 'toa_brightness_temperature')

# test builtin coefficients
self.reader.calib_mode = 'FILE'
res = self.reader.get_dataset(ds_id, ds_info)
# file coefficients are pretty close, give some wiggle room
np.testing.assert_allclose(res.data.compute(), expected, equal_nan=True, atol=0.04)
# make sure the attributes from the file are in the data array
self.assertEqual(res.attrs['standard_name'], 'toa_brightness_temperature')

Expand Down

0 comments on commit 55798fb

Please sign in to comment.