Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Use internal calibration coefficients if available #2

Merged
merged 8 commits into from
Oct 11, 2019
+79 −55
Merged

Use internal calibration coefficients if available #2

Changes from 2 commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
2b30720
Use internal calibration coefficients if available
simonrp84 Oct 8, 2019
e714ef2
Remove built-in slope/offset, use values in file. Add option to switc…
simonrp84 Oct 8, 2019
e6872b8
Fix AMI area extent test to match code
djhoese Oct 9, 2019
46cd7d5
Use internal calibration coefficients if available
simonrp84 Oct 8, 2019
35de058
Remove built-in slope/offset, use values in file. Add option to switc…
simonrp84 Oct 8, 2019
552b2e7
Add tests
simonrp84 Oct 11, 2019
796925e
Merge branch 'feature-ami-reader' of github.com:simonrp84/satpy into …
simonrp84 Oct 11, 2019
6088e5f
Move AMI IR calibration to separate method and add test for file calib
djhoese Oct 11, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
94 changes: 49 additions & 45 deletions satpy/readers/ami_l1b.py
View file
Open in desktop
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 @@ -110,9 +91,9 @@ def get_area_def(self, dsid):
bit_shift = 2**16
area_extent = (
h * np.deg2rad((0 - coff - 0.5) * bit_shift / cfac),
h * np.deg2rad((0 - loff - 0.5) * bit_shift / lfac),
-h * np.deg2rad((0 - loff - 0.5) * bit_shift / lfac),
h * np.deg2rad((cols - coff + 0.5) * bit_shift / cfac),
h * np.deg2rad((rows - loff + 0.5) * bit_shift / lfac),
-h * np.deg2rad((rows - loff + 0.5) * bit_shift / lfac),
)

proj_dict = {
Expand Down Expand Up @@ -182,32 +163,55 @@ 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':
# depends on the radiance calibration above
# Convert um to m^-1 (SI units for pyspectral)
wn = 1 / (dataset_id.wavelength[1] / 1e6)
# Convert cm^-1 (wavenumbers) and (mW/m^2)/(str/cm^-1) (radiance data)
# to SI units m^-1, mW*m^-3*str^-1.
bt_data = rad2temp(wn, data.data * 1e-5)
if isinstance(bt_data, np.ndarray):
# old versions of pyspectral produce numpy arrays
data.data = da.from_array(bt_data, chunks=data.data.chunks)
print(self.calib_mode)
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)
# Convert cm^-1 (wavenumbers) and (mW/m^2)/(str/cm^-1) (radiance data)
# to SI units m^-1, mW*m^-3*str^-1.
bt_data = rad2temp(wn, data.data * 1e-5)
if isinstance(bt_data, np.ndarray):
# old versions of pyspectral produce numpy arrays
data.data = da.from_array(bt_data, chunks=data.data.chunks)
else:
# new versions of pyspectral can do dask arrays
data.data = bt_data
else:
# new versions of pyspectral can do dask arrays
# 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
elif dataset_id.calibration not in ('counts', 'radiance'):
raise ValueError("Unknown calibration: '{}'".format(dataset_id.calibration))
Expand Down