From 64198bc2b46592cfe943d236d66ed5d3784b1db3 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Tue, 27 Jul 2021 13:31:48 -0500 Subject: [PATCH 1/4] Fix awips_tiled writer producing invalid y coordinate variables --- satpy/tests/writer_tests/test_awips_tiled.py | 6 ++++++ satpy/writers/awips_tiled.py | 8 +++----- 2 files changed, 9 insertions(+), 5 deletions(-) diff --git a/satpy/tests/writer_tests/test_awips_tiled.py b/satpy/tests/writer_tests/test_awips_tiled.py index 83b8e04be9..d34f27ab3d 100644 --- a/satpy/tests/writer_tests/test_awips_tiled.py +++ b/satpy/tests/writer_tests/test_awips_tiled.py @@ -46,6 +46,7 @@ def check_required_common_attributes(ds): assert x_attrs.get('standard_name') == 'projection_x_coordinate' assert x_attrs.get('units') == 'meters' assert 'scale_factor' in x_attrs + assert x_attrs['scale_factor'] > 0 assert 'add_offset' in x_attrs assert 'y' in ds.coords @@ -55,6 +56,7 @@ def check_required_common_attributes(ds): assert y_attrs.get('standard_name') == 'projection_y_coordinate' assert y_attrs.get('units') == 'meters' assert 'scale_factor' in y_attrs + assert y_attrs['scale_factor'] < 0 assert 'add_offset' in y_attrs for attr_name in ('tile_row_offset', 'tile_column_offset', @@ -147,6 +149,10 @@ def test_basic_numbered_1_tile(self, use_save_dataset): scale_factor = output_ds['data'].encoding['scale_factor'] np.testing.assert_allclose(input_data_arr.values, output_ds['data'].data, atol=scale_factor / 2) + x_var = output_ds.coords['x'] + assert (np.diff(x_var.values) > 0).all() + y_var = output_ds.coords['y'] + assert (np.diff(y_var.values) < 0).all() @pytest.mark.parametrize( ("tile_count", "tile_size"), diff --git a/satpy/writers/awips_tiled.py b/satpy/writers/awips_tiled.py index 29bed7a321..476994c32e 100644 --- a/satpy/writers/awips_tiled.py +++ b/satpy/writers/awips_tiled.py @@ -303,7 +303,7 @@ def _get_tile_properties(self, tile_shape, tile_count): else: raise ValueError("Either 'tile_count' or 'tile_shape' must be provided") - # number of pixels per each tile + # number of pixels per each tile (rows, cols) self.tile_shape = tile_shape # number of tiles in each direction (rows, columns) self.tile_count = tile_count @@ -342,9 +342,7 @@ def _get_xy_arrays(self): new_extents, ) - x, y = imaginary_grid_def.get_proj_coords() - x = x[0].squeeze() # all rows should have the same coordinates - y = y[:, 0].squeeze() # all columns should have the same coordinates + x, y = imaginary_grid_def.get_proj_vectors() return x, y def _get_xy_scaling_parameters(self): @@ -352,7 +350,7 @@ def _get_xy_scaling_parameters(self): gd = self.area_definition bx = self.x.min() mx = gd.pixel_size_x - by = self.y.min() + by = self.y.max() my = -abs(gd.pixel_size_y) return mx, bx, my, by From 4612265c1099242429cb870e8f461520737a2fd0 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Tue, 27 Jul 2021 14:13:32 -0500 Subject: [PATCH 2/4] Refactor awips_tiled tests to be more specific --- satpy/tests/writer_tests/test_awips_tiled.py | 96 ++++++++++++-------- 1 file changed, 57 insertions(+), 39 deletions(-) diff --git a/satpy/tests/writer_tests/test_awips_tiled.py b/satpy/tests/writer_tests/test_awips_tiled.py index d34f27ab3d..71a2289e2b 100644 --- a/satpy/tests/writer_tests/test_awips_tiled.py +++ b/satpy/tests/writer_tests/test_awips_tiled.py @@ -37,8 +37,33 @@ def _check_production_location(ds): assert len(ds.attrs[prod_loc_name]) == 31 -def check_required_common_attributes(ds): +def check_required_properties(unmasked_ds, masked_ds): + """Check various aspects of coordinates and attributes for correctness.""" + _check_scaled_x_coordinate_variable(unmasked_ds, masked_ds) + _check_scaled_y_coordinate_variable(unmasked_ds, masked_ds) + _check_required_common_attributes(unmasked_ds) + + +def _check_required_common_attributes(ds): """Check common properties of the created AWIPS tiles for validity.""" + for attr_name in ('tile_row_offset', 'tile_column_offset', + 'product_tile_height', 'product_tile_width', + 'number_product_tiles', + 'product_rows', 'product_columns'): + assert attr_name in ds.attrs + _check_production_location(ds) + + for data_arr in ds.data_vars.values(): + if data_arr.ndim == 0: + # grid mapping variable + assert 'grid_mapping_name' in data_arr.attrs + continue + assert data_arr.encoding.get('zlib', False) + assert 'grid_mapping' in data_arr.attrs + assert data_arr.attrs['grid_mapping'] in ds + + +def _check_scaled_x_coordinate_variable(ds, masked_ds): assert 'x' in ds.coords x_coord = ds.coords['x'] np.testing.assert_equal(np.diff(x_coord), 1) @@ -49,6 +74,11 @@ def check_required_common_attributes(ds): assert x_attrs['scale_factor'] > 0 assert 'add_offset' in x_attrs + unscaled_x = masked_ds.coords['x'].values + assert (np.diff(unscaled_x) > 0).all() + + +def _check_scaled_y_coordinate_variable(ds, masked_ds): assert 'y' in ds.coords y_coord = ds.coords['y'] np.testing.assert_equal(np.diff(y_coord), 1) @@ -59,21 +89,8 @@ def check_required_common_attributes(ds): assert y_attrs['scale_factor'] < 0 assert 'add_offset' in y_attrs - for attr_name in ('tile_row_offset', 'tile_column_offset', - 'product_tile_height', 'product_tile_width', - 'number_product_tiles', - 'product_rows', 'product_columns'): - assert attr_name in ds.attrs - _check_production_location(ds) - - for data_arr in ds.data_vars.values(): - if data_arr.ndim == 0: - # grid mapping variable - assert 'grid_mapping_name' in data_arr.attrs - continue - assert data_arr.encoding.get('zlib', False) - assert 'grid_mapping' in data_arr.attrs - assert data_arr.attrs['grid_mapping'] in ds + unscaled_y = masked_ds.coords['y'].values + assert (np.diff(unscaled_y) < 0).all() class TestAWIPSTiledWriter: @@ -143,16 +160,12 @@ def test_basic_numbered_1_tile(self, use_save_dataset): assert len(all_files) == 1 assert os.path.basename(all_files[0]) == 'TESTS_AII_PLAT_SENSOR_test_ds_TEST_T001_20180101_1200.nc' for fn in all_files: - output_ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(output_ds) + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) output_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, output_ds) scale_factor = output_ds['data'].encoding['scale_factor'] np.testing.assert_allclose(input_data_arr.values, output_ds['data'].data, atol=scale_factor / 2) - x_var = output_ds.coords['x'] - assert (np.diff(x_var.values) > 0).all() - y_var = output_ds.coords['y'] - assert (np.diff(y_var.values) < 0).all() @pytest.mark.parametrize( ("tile_count", "tile_size"), @@ -190,12 +203,13 @@ def test_basic_numbered_tiles(self, tile_count, tile_size): expected_num_files = 0 if should_error else 9 assert len(all_files) == expected_num_files for fn in all_files: - ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(ds) - assert ds.attrs['my_global'] == 'TEST' - assert ds.attrs['sector_id'] == 'TEST' + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) + assert unmasked_ds.attrs['my_global'] == 'TEST' + assert unmasked_ds.attrs['sector_id'] == 'TEST' stime = input_data_arr.attrs['start_time'] - assert ds.attrs['start_date_time'] == stime.strftime('%Y-%m-%dT%H:%M:%S') + assert unmasked_ds.attrs['start_date_time'] == stime.strftime('%Y-%m-%dT%H:%M:%S') def test_basic_lettered_tiles(self): """Test creating a lettered grid.""" @@ -232,9 +246,10 @@ def test_basic_lettered_tiles(self): all_files = glob(os.path.join(self.base_dir, 'TESTS_AII*.nc')) assert len(all_files) == 16 for fn in all_files: - ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(ds) - assert ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) + assert masked_ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') def test_lettered_tiles_update_existing(self): """Test updating lettered tiles with additional data.""" @@ -379,9 +394,10 @@ def test_lettered_tiles_sector_ref(self): all_files = glob(os.path.join(self.base_dir, 'TESTS_AII*.nc')) assert len(all_files) == 16 for fn in all_files: - ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(ds) - assert ds.attrs['start_date_time'] == (now + timedelta(minutes=20)).strftime('%Y-%m-%dT%H:%M:%S') + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) + assert masked_ds.attrs['start_date_time'] == (now + timedelta(minutes=20)).strftime('%Y-%m-%dT%H:%M:%S') def test_lettered_tiles_no_fit(self): """Test creating a lettered grid with no data overlapping the grid.""" @@ -528,8 +544,9 @@ def test_basic_numbered_tiles_rgb(self): assert len(chan_files) == 9 all_files.extend(chan_files) for fn in all_files: - ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(ds) + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) @pytest.mark.parametrize( "sector", @@ -604,12 +621,13 @@ def test_multivar_numbered_tiles_glm(self, sector, extra_kwargs): all_files = glob(os.path.join(self.base_dir, fn_glob)) assert len(all_files) == 9 for fn in all_files: - ds = xr.open_dataset(fn, mask_and_scale=False) - check_required_common_attributes(ds) + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) if sector == 'C': - assert ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%S.%fZ') + assert masked_ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%S.%fZ') else: # 'F' - assert ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%SZ') + assert masked_ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%SZ') @staticmethod def _get_glm_glob_filename(extra_kwargs): From 7be8c2a57fe6db84b3e7e283dfba6146ac37800a Mon Sep 17 00:00:00 2001 From: David Hoese Date: Thu, 29 Jul 2021 11:31:10 -0500 Subject: [PATCH 3/4] Fix 'awips_tiled' writer not using data projection for lettered tile calculations --- satpy/tests/writer_tests/test_awips_tiled.py | 39 ++++++++++++++++++++ satpy/writers/awips_tiled.py | 28 ++++++++++---- 2 files changed, 59 insertions(+), 8 deletions(-) diff --git a/satpy/tests/writer_tests/test_awips_tiled.py b/satpy/tests/writer_tests/test_awips_tiled.py index 71a2289e2b..3c08019d61 100644 --- a/satpy/tests/writer_tests/test_awips_tiled.py +++ b/satpy/tests/writer_tests/test_awips_tiled.py @@ -251,6 +251,45 @@ def test_basic_lettered_tiles(self): check_required_properties(unmasked_ds, masked_ds) assert masked_ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') + def test_basic_lettered_tiles_diff_projection(self): + """Test creating a lettered grid from data with differing projection..""" + import xarray as xr + from satpy.writers.awips_tiled import AWIPSTiledWriter + from xarray import DataArray + from pyresample.geometry import AreaDefinition + w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) + area_def = AreaDefinition( + 'test', + 'test', + 'test', + '+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. +lat_0=45 +lat_1=45 +units=m +no_defs', + 1000, + 2000, + (-1000000., -1500000., 1000000., 1500000.), + ) + now = datetime(2018, 1, 1, 12, 0, 0) + ds = DataArray( + da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), + attrs=dict( + name='test_ds', + platform_name='PLAT', + sensor='SENSOR', + units='1', + area=area_def, + start_time=now, + end_time=now + timedelta(minutes=20)) + ) + # tile_count should be ignored since we specified lettered_grid + w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) + all_files = sorted(glob(os.path.join(self.base_dir, 'TESTS_AII*.nc'))) + assert len(all_files) == 24 + assert "TC02" in all_files[0] # the first tile should be TC02 + for fn in all_files: + unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) + masked_ds = xr.open_dataset(fn, mask_and_scale=True) + check_required_properties(unmasked_ds, masked_ds) + assert masked_ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') + def test_lettered_tiles_update_existing(self): """Test updating lettered tiles with additional data.""" import shutil diff --git a/satpy/writers/awips_tiled.py b/satpy/writers/awips_tiled.py index 476994c32e..e2cc824c73 100644 --- a/satpy/writers/awips_tiled.py +++ b/satpy/writers/awips_tiled.py @@ -227,7 +227,7 @@ from trollsift.parser import StringFormatter, Parser import numpy as np -from pyproj import Proj +from pyproj import Proj, CRS, Transformer import dask import dask.array as da import xarray as xr @@ -422,10 +422,20 @@ def __call__(self): class LetteredTileGenerator(NumberedTileGenerator): """Helper class to generate per-tile metadata for lettered tiles.""" - def __init__(self, area_definition, extents, + def __init__(self, area_definition, extents, sector_crs, cell_size=(2000000, 2000000), num_subtiles=None, use_sector_reference=False): - """Initialize tile information for later generation.""" + """Initialize tile information for later generation. + + Args: + area_definition (AreaDefinition): Area of the data being saved. + extents (tuple): Four element tuple of the configured lettered + area. + sector_crs (pyproj.CRS): CRS of the configured lettered sector + area. + cell_size (tuple): Two element tuple of resolution of each tile + in sector projection units (y, x). + """ # (row subtiles, col subtiles) self.num_subtiles = num_subtiles or (2, 2) self.cell_size = cell_size # (row tile height, col tile width) @@ -433,7 +443,8 @@ def __init__(self, area_definition, extents, self.ll_extents = extents[:2] # (x min, y min) self.ur_extents = extents[2:] # (x max, y max) self.use_sector_reference = use_sector_reference - super(LetteredTileGenerator, self).__init__(area_definition) + self._transformer = Transformer.from_crs(sector_crs, area_definition.crs) + super().__init__(area_definition) def _get_tile_properties(self, tile_shape, tile_count): """Calculate tile information for this particular sector/grid.""" @@ -445,8 +456,8 @@ def _get_tile_properties(self, tile_shape, tile_count): ad = self.area_definition x, y = ad.get_proj_vectors() - ll_xy = self.ll_extents - ur_xy = self.ur_extents + ll_xy = self._transformer.transform(*self.ll_extents) + ur_xy = self._transformer.transform(*self.ur_extents) cw = abs(ad.pixel_size_x) ch = abs(ad.pixel_size_y) st = self.num_subtiles @@ -1261,6 +1272,7 @@ def separate_init_kwargs(cls, kwargs): def _fill_sector_info(self): """Convert sector extents if needed.""" for sector_info in self.awips_sectors.values(): + sector_info['projection'] = CRS.from_user_input(sector_info['projection']) p = Proj(sector_info['projection']) if 'lower_left_xy' in sector_info: sector_info['lower_left_lonlat'] = p(*sector_info['lower_left_xy'], inverse=True) @@ -1297,6 +1309,7 @@ def _get_tile_generator(self, area_def, lettered_grid, sector_id, tile_gen = LetteredTileGenerator( area_def, sector_info['lower_left_xy'] + sector_info['upper_right_xy'], + sector_crs=sector_info['projection'], cell_size=sector_info['resolution'], num_subtiles=num_subtiles, use_sector_reference=use_sector_reference, @@ -1689,7 +1702,6 @@ def _create_debug_array(sector_info, num_subtiles, font_path='Verdana.ttf'): img.save("test.png") - from pyresample.utils import proj4_str_to_dict new_extents = ( ll_extent[0], ur_extent[1] - 1001. * meters_ppy, @@ -1700,7 +1712,7 @@ def _create_debug_array(sector_info, num_subtiles, font_path='Verdana.ttf'): 'debug_grid', 'debug_grid', 'debug_grid', - proj4_str_to_dict(sector_info['projection']), + sector_info['projection'], 1000, 1000, new_extents From 68162e192912ad61b855d0d5325c6f19edab54b5 Mon Sep 17 00:00:00 2001 From: David Hoese Date: Thu, 29 Jul 2021 12:10:28 -0500 Subject: [PATCH 4/4] Refactor awips_tiled writer tests to reuse test case creation --- satpy/tests/writer_tests/test_awips_tiled.py | 342 +++++-------------- 1 file changed, 81 insertions(+), 261 deletions(-) diff --git a/satpy/tests/writer_tests/test_awips_tiled.py b/satpy/tests/writer_tests/test_awips_tiled.py index 3c08019d61..ece67f660c 100644 --- a/satpy/tests/writer_tests/test_awips_tiled.py +++ b/satpy/tests/writer_tests/test_awips_tiled.py @@ -22,6 +22,7 @@ import numpy as np import dask.array as da +from pyproj import CRS import pytest @@ -100,6 +101,8 @@ def setup_method(self): """Create temporary directory to save files to.""" import tempfile self.base_dir = tempfile.mkdtemp() + self.start_time = datetime(2018, 1, 1, 12, 0, 0) + self.end_time = self.start_time + timedelta(minutes=20) def teardown_method(self): """Remove the temporary directory created for a test.""" @@ -114,32 +117,39 @@ def test_init(self): from satpy.writers.awips_tiled import AWIPSTiledWriter AWIPSTiledWriter(base_dir=self.base_dir) - def _get_test_lcc_data(self): - from xarray import DataArray + def _get_test_area(self, shape=(200, 100), crs=None, extents=None): from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict + if crs is None: + crs = CRS('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. +lat_0=25 +lat_1=25 +units=m +no_defs') + if extents is None: + extents = (-1000., -1500., 1000., 1500.) area_def = AreaDefinition( 'test', 'test', 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 100, - 200, - (-1000., -1500., 1000., 1500.), + crs, + shape[1], + shape[0], + extents, ) - now = datetime(2018, 1, 1, 12, 0, 0) - data = np.linspace(0., 1., 20000, dtype=np.float32).reshape((200, 100)) + return area_def + + def _get_test_data(self, shape=(200, 100), chunks=50): + data = np.linspace(0., 1., shape[0] * shape[1], dtype=np.float32).reshape(shape) + return da.from_array(data, chunks=chunks) + + def _get_test_lcc_data(self, dask_arr, area_def): + from xarray import DataArray ds = DataArray( - da.from_array(data, chunks=50), + dask_arr, attrs=dict( name='test_ds', platform_name='PLAT', sensor='SENSOR', units='1', area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) + start_time=self.start_time, + end_time=self.end_time) ) return ds @@ -149,7 +159,9 @@ def test_basic_numbered_1_tile(self, use_save_dataset): """Test creating a single numbered tile.""" import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - input_data_arr = self._get_test_lcc_data() + data = self._get_test_data() + area_def = self._get_test_area() + input_data_arr = self._get_test_lcc_data(data, area_def) w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) if use_save_dataset: w.save_dataset(input_data_arr, sector_id='TEST', source_name='TESTS') @@ -181,7 +193,9 @@ def test_basic_numbered_tiles(self, tile_count, tile_size): import dask from satpy.writers.awips_tiled import AWIPSTiledWriter from satpy.tests.utils import CustomScheduler - input_data_arr = self._get_test_lcc_data() + data = self._get_test_data() + area_def = self._get_test_area() + input_data_arr = self._get_test_lcc_data(data, area_def) w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) save_kwargs = dict( sector_id='TEST', @@ -215,32 +229,11 @@ def test_basic_lettered_tiles(self): """Test creating a lettered grid.""" import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = self._get_test_data(shape=(2000, 1000), chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) # tile_count should be ignored since we specified lettered_grid w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) all_files = glob(os.path.join(self.base_dir, 'TESTS_AII*.nc')) @@ -249,36 +242,18 @@ def test_basic_lettered_tiles(self): unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) masked_ds = xr.open_dataset(fn, mask_and_scale=True) check_required_properties(unmasked_ds, masked_ds) - assert masked_ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') + assert masked_ds.attrs['start_date_time'] == self.start_time.strftime('%Y-%m-%dT%H:%M:%S') def test_basic_lettered_tiles_diff_projection(self): """Test creating a lettered grid from data with differing projection..""" import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - '+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. +lat_0=45 +lat_1=45 +units=m +no_defs', - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + crs = CRS("+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. +lat_0=45 +lat_1=45 +units=m +no_defs") + data = self._get_test_data(shape=(2000, 1000), chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), crs=crs, + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) # tile_count should be ignored since we specified lettered_grid w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) all_files = sorted(glob(os.path.join(self.base_dir, 'TESTS_AII*.nc'))) @@ -288,44 +263,24 @@ def test_basic_lettered_tiles_diff_projection(self): unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) masked_ds = xr.open_dataset(fn, mask_and_scale=True) check_required_properties(unmasked_ds, masked_ds) - assert masked_ds.attrs['start_date_time'] == now.strftime('%Y-%m-%dT%H:%M:%S') + assert masked_ds.attrs['start_date_time'] == self.start_time.strftime('%Y-%m-%dT%H:%M:%S') def test_lettered_tiles_update_existing(self): """Test updating lettered tiles with additional data.""" import shutil import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict import dask first_base_dir = os.path.join(self.base_dir, 'first') w = AWIPSTiledWriter(base_dir=first_base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - data = np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)) + shape = (2000, 1000) + data = np.linspace(0., 1., shape[0] * shape[1], dtype=np.float32).reshape(shape) # pixels to be filled in later data[:, -200:] = np.nan - ds = DataArray( - da.from_array(data, chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = da.from_array(data, chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) # tile_count should be ignored since we specified lettered_grid w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) all_files = sorted(glob(os.path.join(first_base_dir, 'TESTS_AII*.nc'))) @@ -340,32 +295,15 @@ def test_lettered_tiles_update_existing(self): # Second writing/updating # Area is about 100 pixels to the right - area_def2 = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-800000., -1500000., 1200000., 1500000.), - ) + area_def2 = self._get_test_area(shape=(2000, 1000), + extents=(-800000., -1500000., 1200000., 1500000.)) data2 = np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)) # a gap at the beginning where old values remain data2[:, :200] = np.nan # a gap at the end where old values remain data2[:, -400:-300] = np.nan - ds2 = DataArray( - da.from_array(data2, chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def2, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data2 = da.from_array(data2, chunks=500) + ds2 = self._get_test_lcc_data(data2, area_def2) w = AWIPSTiledWriter(base_dir=second_base_dir, compress=True) # HACK: The _copy_to_existing function hangs when opening the output # file multiple times...sometimes. If we limit dask to one worker @@ -401,32 +339,11 @@ def test_lettered_tiles_sector_ref(self): """Test creating a lettered grid using the sector as reference.""" import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = self._get_test_data(shape=(2000, 1000), chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) w.save_datasets([ds], sector_id='LCC', source_name="TESTS", lettered_grid=True, use_sector_reference=True, use_end_time=True) @@ -436,37 +353,17 @@ def test_lettered_tiles_sector_ref(self): unmasked_ds = xr.open_dataset(fn, mask_and_scale=False) masked_ds = xr.open_dataset(fn, mask_and_scale=True) check_required_properties(unmasked_ds, masked_ds) - assert masked_ds.attrs['start_date_time'] == (now + timedelta(minutes=20)).strftime('%Y-%m-%dT%H:%M:%S') + expected_start = (self.start_time + timedelta(minutes=20)).strftime('%Y-%m-%dT%H:%M:%S') + assert masked_ds.attrs['start_date_time'] == expected_start def test_lettered_tiles_no_fit(self): """Test creating a lettered grid with no data overlapping the grid.""" from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (4000000., 5000000., 5000000., 6000000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = self._get_test_data(shape=(2000, 1000), chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), + extents=(4000000., 5000000., 5000000., 6000000.)) + ds = self._get_test_lcc_data(data, area_def) w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) # No files created all_files = glob(os.path.join(self.base_dir, 'TESTS_AII*.nc')) @@ -475,32 +372,11 @@ def test_lettered_tiles_no_fit(self): def test_lettered_tiles_no_valid_data(self): """Test creating a lettered grid with no valid data.""" from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.full((2000, 1000), np.nan, chunks=500, dtype=np.float32), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = da.full((2000, 1000), np.nan, chunks=500, dtype=np.float32), + area_def = self._get_test_area(shape=(2000, 1000), + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) w.save_datasets([ds], sector_id='LCC', source_name="TESTS", tile_count=(3, 3), lettered_grid=True) # No files created - all NaNs should result in no tiles being created all_files = glob(os.path.join(self.base_dir, 'TESTS_AII*.nc')) @@ -509,31 +385,11 @@ def test_lettered_tiles_no_valid_data(self): def test_lettered_tiles_bad_filename(self): """Test creating a lettered grid with a bad filename.""" from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True, filename="{Bad Key}.nc") - area_def = AreaDefinition( - 'test', - 'test', - 'test', - ('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 1000, - 2000, - (-1000000., -1500000., 1000000., 1500000.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 2000000, dtype=np.float32).reshape((2000, 1000)), chunks=500), - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = self._get_test_data(shape=(2000, 1000), chunks=500) + area_def = self._get_test_area(shape=(2000, 1000), + extents=(-1000000., -1500000., 1000000., 1500000.)) + ds = self._get_test_lcc_data(data, area_def) with pytest.raises(KeyError): w.save_datasets([ds], sector_id='LCC', @@ -545,33 +401,13 @@ def test_basic_numbered_tiles_rgb(self): """Test creating a multiple numbered tiles with RGB.""" from satpy.writers.awips_tiled import AWIPSTiledWriter import xarray as xr - from xarray import DataArray - from pyresample.geometry import AreaDefinition w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - ('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 100, - 200, - (-1000., -1500., 1000., 1500.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - ds = DataArray( - da.from_array(np.linspace(0., 1., 60000, dtype=np.float32).reshape((3, 200, 100)), chunks=50), - dims=('bands', 'y', 'x'), - coords={'bands': ['R', 'G', 'B']}, - attrs=dict( - name='test_ds', - platform_name='PLAT', - sensor='SENSOR', - units='1', - area=area_def, - start_time=now, - end_time=now + timedelta(minutes=20)) - ) + data = da.from_array(np.linspace(0., 1., 60000, dtype=np.float32).reshape((3, 200, 100)), chunks=50) + area_def = self._get_test_area() + ds = self._get_test_lcc_data(data, area_def) + ds = ds.rename(dict((old, new) for old, new in zip(ds.dims, ['bands', 'y', 'x']))) + ds.coords['bands'] = ['R', 'G', 'B'] + w.save_datasets([ds], sector_id='TEST', source_name="TESTS", tile_count=(3, 3)) chan_files = glob(os.path.join(self.base_dir, 'TESTS_AII*test_ds_R*.nc')) all_files = chan_files[:] @@ -604,34 +440,18 @@ def test_multivar_numbered_tiles_glm(self, sector, extra_kwargs): """Test creating a tiles with multiple variables.""" import xarray as xr from satpy.writers.awips_tiled import AWIPSTiledWriter - from xarray import DataArray - from pyresample.geometry import AreaDefinition - from pyresample.utils import proj4_str_to_dict import os os.environ['ORGANIZATION'] = '1' * 50 w = AWIPSTiledWriter(base_dir=self.base_dir, compress=True) - area_def = AreaDefinition( - 'test', - 'test', - 'test', - proj4_str_to_dict('+proj=lcc +datum=WGS84 +ellps=WGS84 +lon_0=-95. ' - '+lat_0=25 +lat_1=25 +units=m +no_defs'), - 100, - 200, - (-1000., -1500., 1000., 1500.), - ) - now = datetime(2018, 1, 1, 12, 0, 0) - end_time = now + timedelta(minutes=20) - ds1 = DataArray( - da.from_array(np.linspace(0., 1., 20000, dtype=np.float32).reshape((200, 100)), chunks=50), - attrs=dict( + data = self._get_test_data() + area_def = self._get_test_area() + ds1 = self._get_test_lcc_data(data, area_def) + ds1.attrs.update( + dict( name='total_energy', platform_name='GOES-17', sensor='SENSOR', units='1', - area=area_def, - start_time=now, - end_time=end_time, scan_mode='M3', scene_abbr=sector, platform_shortname="G17" @@ -664,9 +484,9 @@ def test_multivar_numbered_tiles_glm(self, sector, extra_kwargs): masked_ds = xr.open_dataset(fn, mask_and_scale=True) check_required_properties(unmasked_ds, masked_ds) if sector == 'C': - assert masked_ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%S.%fZ') + assert masked_ds.attrs['time_coverage_end'] == self.end_time.strftime('%Y-%m-%dT%H:%M:%S.%fZ') else: # 'F' - assert masked_ds.attrs['time_coverage_end'] == end_time.strftime('%Y-%m-%dT%H:%M:%SZ') + assert masked_ds.attrs['time_coverage_end'] == self.end_time.strftime('%Y-%m-%dT%H:%M:%SZ') @staticmethod def _get_glm_glob_filename(extra_kwargs):