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

Fix handling of file attributes in seviri_l1b_nc reader #791

Merged
merged 2 commits into from
May 28, 2019
Merged

Fix handling of file attributes in seviri_l1b_nc reader #791

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
b144a0a
Fixes Issue 753 - handles nc file attrubutes correctly and reads/disp…
May 23, 2019
0188475
seviri level1b reader added comments for code clarification and link …
May 27, 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
31 changes: 14 additions & 17 deletions satpy/etc/readers/seviri_l1b_nc.yaml
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,12 +5,9 @@ reader:
reader: !!python/name:satpy.readers.yaml_reader.FileYAMLReader

file_types:
nc_seviri_l1b:
seviri_l1b_nc:
file_reader: !!python/name:satpy.readers.seviri_l1b_nc.NCSEVIRIFileHandler
file_patterns: ['W_XX-EUMETSAT-Darmstadt,VIS+IR+IMAGERY,{satid:4s}+SEVIRI_C_EUMG_{processing_time:%Y%m%d%H%M%S}.nc']
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Isn't the VIS+IR+IMAGERY pattern still needed if one de-selects HRV in the order (EUM data centre)?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

hi, when i ordered the data the Netcdf option didnt allow me to deselect any channels. But now I see that the Netcdf (GSICS) option does so i will update the yaml and test this.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

so there is an issue with the netcdf data ordered using netcdf (GSICS) as some attributes are missing, Also spacecraft id is missing from the filename. EUM has been informed

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh, too bad. Shall we wait for them to fix it or merge now and add support for GSICS nc-files later?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

after talking to @ColinDuff, I suggest we merge this now and deal with GSICS nc-files later. this way users can use the reader for "regular" nc-files.

nc_seviri_l1b_hrv:
file_reader: !!python/name:satpy.readers.seviri_l1b_nc.NCSEVIRIHRVFileHandler
file_patterns: ['W_XX-EUMETSAT-Darmstadt,HRV+IMAGERY,{satid:4s}+SEVIRI_C_EUMG_{processing_time:%Y%m%d%H%M%S}.nc']
file_patterns: ['W_XX-EUMETSAT-Darmstadt,VIS+IR+HRV+IMAGERY,{satid:4s}+SEVIRI_C_EUMG_{processing_time:%Y%m%d%H%M%S}.nc']

datasets:
HRV:
Expand All @@ -27,7 +24,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b_hrv
file_type: seviri_l1b_nc
nc_key: 'ch12'

IR_016:
Expand All @@ -44,7 +41,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch3'

IR_039:
Expand All @@ -61,7 +58,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch4'

IR_087:
Expand All @@ -78,7 +75,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch5'

IR_097:
Expand All @@ -95,7 +92,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch8'

IR_108:
Expand All @@ -112,7 +109,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch9'

IR_120:
Expand All @@ -129,7 +126,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch10'

IR_134:
Expand All @@ -146,7 +143,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch11'

VIS006:
Expand All @@ -163,7 +160,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch1'


Expand All @@ -181,7 +178,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch2'

WV_062:
Expand All @@ -198,7 +195,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch6'

WV_073:
Expand All @@ -215,7 +212,7 @@ datasets:
counts:
standard_name: counts
units: count
file_type: nc_seviri_l1b
file_type: seviri_l1b_nc
nc_key: 'ch7'


122 changes: 66 additions & 56 deletions satpy/readers/seviri_l1b_nc.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -20,16 +20,18 @@

# You should have received a copy of the GNU General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
"""SEVIRI netcdf format reader. """
"""SEVIRI netcdf format reader.

References:
MSG Level 1.5 Image Data Format Description
https://www.eumetsat.int/website/wcm/idc/idcplg?IdcService=GET_FILE&dDocName=PDF_TEN_05105_MSG_IMG_DATA&RevisionSelectionMethod=LatestReleased&Rendition=Web
"""

from satpy.readers.file_handlers import BaseFileHandler
from satpy.readers.seviri_base import (SEVIRICalibrationHandler,
CHANNEL_NAMES, CALIB, SATNUM)
import xarray as xr
# Needed for xarray netcdf workaround
from netCDF4 import Dataset as tmpDataset

from satpy import CHUNK_SIZE
from pyresample import geometry

import datetime
Expand Down Expand Up @@ -57,10 +59,7 @@ def _read_file(self):
self.nc = xr.open_dataset(self.filename,
decode_cf=True,
mask_and_scale=False,
chunks={'num_columns_vis_ir': CHUNK_SIZE,
'num_rows_vis_ir': CHUNK_SIZE})

self.nc = self.nc.rename({'num_columns_vis_ir': 'x', 'num_rows_vis_ir': 'y'})
chunks={})

# Obtain some area definition attributes
equatorial_radius = (self.nc.attrs['equatorial_radius'] * 1000.)
Expand All @@ -71,8 +70,11 @@ def _read_file(self):
'h': 35785831.00,
'ssp_longitude': ssp_lon}

self.mda['number_of_lines'] = int(self.nc.dims['y'])
self.mda['number_of_columns'] = int(self.nc.dims['x'])
self.mda['number_of_lines'] = int(self.nc.dims['num_rows_vis_ir'])
self.mda['number_of_columns'] = int(self.nc.dims['num_columns_vis_ir'])

self.mda['hrv_number_of_lines'] = int(self.nc.dims['num_rows_hrv'])
self.mda['hrv_number_of_columns'] = int(self.nc.dims['num_columns_hrv'])

self.deltaSt = self.reference + datetime.timedelta(
days=int(self.nc.attrs['true_repeat_cycle_start_day']),
Expand All @@ -82,18 +84,27 @@ def _read_file(self):
days=int(self.nc.attrs['planned_repeat_cycle_end_day']),
milliseconds=int(self.nc.attrs['planned_repeat_cycle_end_mi_sec']))

# Netcdf xrray dimensions workaround - these 4 dimensions not found in the dataset
# but ncdump confirms they are there
tmp = tmpDataset(self.filename)
self.north = int(tmp.dimensions['north_most_line'].size)
self.east = int(tmp.dimensions['east_most_pixel'].size)
self.west = int(tmp.dimensions['west_most_pixel'].size)
self.south = int(tmp.dimensions['south_most_line'].size)
self.north = int(self.nc.attrs['north_most_line'])
self.east = int(self.nc.attrs['east_most_pixel'])
self.west = int(self.nc.attrs['west_most_pixel'])
self.south = int(self.nc.attrs['south_most_line'])

def get_dataset(self, dataset_id, dataset_info):

channel = dataset_id.name
i = list(CHANNEL_NAMES.values()).index(channel)

if (channel == 'HRV'):
self.nc = self.nc.rename({'num_columns_hrv': 'x', 'num_rows_hrv': 'y'})
else:
# the first channel of a composite will rename the dimension variable
# but the later channels will raise a value error as its already been renamed
# we can just ignore these exceptions
try:
self.nc = self.nc.rename({'num_columns_vis_ir': 'x', 'num_rows_vis_ir': 'y'})
except ValueError:
pass

dataset = self.nc[dataset_info['nc_key']]

dataset.attrs.update(dataset_info)
Expand Down Expand Up @@ -147,6 +158,7 @@ def get_area_def(self, dataset_id):
else:
nlines = self.mda['number_of_lines']
ncols = self.mda['number_of_columns']

area = geometry.AreaDefinition(
'some_area_name',
"On-the-fly area",
Expand All @@ -159,49 +171,47 @@ def get_area_def(self, dataset_id):
return area

def get_area_extent(self, dsid):
if dsid.name != 'HRV':

# following calculations assume grid origin is south-east corner
# section 7.2.4 of MSG Level 1.5 Image Data Format Description
origins = {0: 'NW', 1: 'SW', 2: 'SE', 3: 'NE'}
grid_origin = self.nc.attrs['vis_ir_grid_origin']
grid_origin = int(grid_origin, 16)
if grid_origin != 2:
raise NotImplementedError(
'Grid origin not supported number: {}, {} corner'
.format(grid_origin, origins[grid_origin])
)

center_point = 3712/2

column_step = self.nc.attrs['vis_ir_column_dir_grid_step'] * 1000.0

line_step = self.nc.attrs['vis_ir_line_dir_grid_step'] * 1000.0
# section 3.1.4.2 of MSG Level 1.5 Image Data Format Description
earth_model = int(self.nc.attrs['type_of_earth_model'], 16)
if earth_model == 2:
ns_offset = 0 # north +ve
we_offset = 0 # west +ve
elif earth_model == 1:
ns_offset = -0.5 # north +ve
we_offset = 0.5 # west +ve
else:
raise NotImplementedError(
'unrecognised earth model: {}'.format(earth_model)
)
# section 3.1.5 of MSG Level 1.5 Image Data Format Description
ll_c = (center_point - self.west - 0.5 + we_offset) * column_step
ll_l = (self.south - center_point - 0.5 + ns_offset) * line_step
ur_c = (center_point - self.east + 0.5 + we_offset) * column_step
ur_l = (self.north - center_point + 0.5 + ns_offset) * line_step
area_extent = (ll_c, ll_l, ur_c, ur_l)

# following calculations assume grid origin is south-east corner
# section 7.2.4 of MSG Level 1.5 Image Data Format Description
origins = {0: 'NW', 1: 'SW', 2: 'SE', 3: 'NE'}
grid_origin = self.nc.attrs['vis_ir_grid_origin']
grid_origin = int(grid_origin, 16)
if grid_origin != 2:
raise NotImplementedError(
'Grid origin not supported number: {}, {} corner'
.format(grid_origin, origins[grid_origin])
)

center_point = 3712/2

column_step = self.nc.attrs['vis_ir_column_dir_grid_step'] * 1000.0

line_step = self.nc.attrs['vis_ir_line_dir_grid_step'] * 1000.0

# check for Earth model as this affects the north-south and
# west-east offsets
# section 3.1.4.2 of MSG Level 1.5 Image Data Format Description
earth_model = int(self.nc.attrs['type_of_earth_model'], 16)
if earth_model == 2:
ns_offset = 0 # north +ve
we_offset = 0 # west +ve
elif earth_model == 1:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe add a short comment that the data is shifted to N-W and to correct for it, we have to adjust the corners in the opposite direction, i.e. S-E

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

the comment can be copied from the seviri_l1b_native-reader.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Working on creating tests for this reader now

I will add the comment as suggested by sjoro

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

comment added

ns_offset = -0.5 # north +ve
we_offset = 0.5 # west +ve
else:

raise NotImplementedError('HRV not supported!')
raise NotImplementedError(
'unrecognised earth model: {}'.format(earth_model)
)
# section 3.1.5 of MSG Level 1.5 Image Data Format Description
ll_c = (center_point - self.west - 0.5 + we_offset) * column_step
ll_l = (self.south - center_point - 0.5 + ns_offset) * line_step
ur_c = (center_point - self.east + 0.5 + we_offset) * column_step
ur_l = (self.north - center_point + 0.5 + ns_offset) * line_step
area_extent = (ll_c, ll_l, ur_c, ur_l)

return area_extent


class NCSEVIRIHRVFileHandler():
class NCSEVIRIHRVFileHandler(BaseFileHandler, SEVIRICalibrationHandler):
pass