diff --git a/satpy/readers/sar_c_safe.py b/satpy/readers/sar_c_safe.py index b6d84e8fb7..e8676c9b56 100644 --- a/satpy/readers/sar_c_safe.py +++ b/satpy/readers/sar_c_safe.py @@ -78,7 +78,7 @@ def _dictify(r): return int(r.text) except ValueError: try: - return float(r.text) + return np.float32(r.text) except ValueError: return r.text for x in r.findall("./*"): @@ -186,7 +186,7 @@ def get_dataset(self, key, info, chunks=None): def get_calibration_constant(self): """Load the calibration constant.""" - return float(self.root.find(".//absoluteCalibrationConstant").text) + return np.float32(self.root.find(".//absoluteCalibrationConstant").text) def _get_calibration_uncached(self, calibration, chunks=None): """Get the calibration array.""" @@ -341,7 +341,7 @@ def _get_array_pieces_for_current_line(self, current_line): current_blocks = self._find_blocks_covering_line(current_line) current_blocks.sort(key=(lambda x: x.coords["x"][0])) next_line = self._get_next_start_line(current_blocks, current_line) - current_y = np.arange(current_line, next_line) + current_y = np.arange(current_line, next_line, dtype=np.uint16) pieces = [arr.sel(y=current_y) for arr in current_blocks] return pieces @@ -389,7 +389,7 @@ def _get_padded_dask_pieces(self, pieces, chunks): @staticmethod def _fill_dask_pieces(dask_pieces, shape, chunks): if shape[1] > 0: - new_piece = da.full(shape, np.nan, chunks=chunks) + new_piece = da.full(shape, np.nan, chunks=chunks, dtype=np.float32) dask_pieces.append(new_piece) @@ -425,11 +425,10 @@ def expand(self, chunks): # corr = 1.5 data = self.lut * corr - x_coord = np.arange(self.first_pixel, self.last_pixel + 1) - y_coord = np.arange(self.first_line, self.last_line + 1) - - new_arr = (da.ones((len(y_coord), len(x_coord)), chunks=chunks) * - np.interp(y_coord, self.lines, data)[:, np.newaxis]) + x_coord = np.arange(self.first_pixel, self.last_pixel + 1, dtype=np.uint16) + y_coord = np.arange(self.first_line, self.last_line + 1, dtype=np.uint16) + new_arr = (da.ones((len(y_coord), len(x_coord)), dtype=np.float32, chunks=chunks) * + np.interp(y_coord, self.lines, data)[:, np.newaxis].astype(np.float32)) new_arr = xr.DataArray(new_arr, dims=["y", "x"], coords={"x": x_coord, @@ -438,29 +437,29 @@ def expand(self, chunks): @property def first_pixel(self): - return int(self.element.find("firstRangeSample").text) + return np.uint16(self.element.find("firstRangeSample").text) @property def last_pixel(self): - return int(self.element.find("lastRangeSample").text) + return np.uint16(self.element.find("lastRangeSample").text) @property def first_line(self): - return int(self.element.find("firstAzimuthLine").text) + return np.uint16(self.element.find("firstAzimuthLine").text) @property def last_line(self): - return int(self.element.find("lastAzimuthLine").text) + return np.uint16(self.element.find("lastAzimuthLine").text) @property def lines(self): lines = self.element.find("line").text.split() - return np.array(lines).astype(int) + return np.array(lines).astype(np.uint16) @property def lut(self): lut = self.element.find("noiseAzimuthLut").text.split() - return np.array(lut).astype(float) + return np.array(lut, dtype=np.float32) class XMLArray: @@ -487,7 +486,7 @@ def _read_xml_array(self): new_x = elt.find("pixel").text.split() y += [int(elt.find("line").text)] * len(new_x) x += [int(val) for val in new_x] - data += [float(val) + data += [np.float32(val) for val in elt.find(self.element_tag).text.split()] return np.asarray(data), (x, y) @@ -519,17 +518,17 @@ def interpolate_xarray_linear(xpoints, ypoints, values, shape, chunks=CHUNK_SIZE else: vchunks, hchunks = chunks, chunks - points = _ndim_coords_from_arrays(np.vstack((np.asarray(ypoints), - np.asarray(xpoints))).T) + points = _ndim_coords_from_arrays(np.vstack((np.asarray(ypoints, dtype=np.uint16), + np.asarray(xpoints, dtype=np.uint16))).T) interpolator = LinearNDInterpolator(points, values) - grid_x, grid_y = da.meshgrid(da.arange(shape[1], chunks=hchunks), - da.arange(shape[0], chunks=vchunks)) + grid_x, grid_y = da.meshgrid(da.arange(shape[1], chunks=hchunks, dtype=np.uint16), + da.arange(shape[0], chunks=vchunks, dtype=np.uint16)) # workaround for non-thread-safe first call of the interpolator: interpolator((0, 0)) - res = da.map_blocks(intp, grid_x, grid_y, interpolator=interpolator) + res = da.map_blocks(intp, grid_x, grid_y, interpolator=interpolator).astype(values.dtype) return DataArray(res, dims=("y", "x")) @@ -617,7 +616,7 @@ def _calibrate_and_denoise(self, data, key): def _get_digital_number(self, data): """Get the digital numbers (uncalibrated data).""" data = data.where(data > 0) - data = data.astype(np.float64) + data = data.astype(np.float32) dn = data * data return dn @@ -675,8 +674,8 @@ def get_gcps(self): for feature in gcps["features"]] gcp_array = np.array(gcp_list) - ypoints = np.unique(gcp_array[:, 0]) - xpoints = np.unique(gcp_array[:, 1]) + ypoints = np.unique(gcp_array[:, 0]).astype(np.uint16) + xpoints = np.unique(gcp_array[:, 1]).astype(np.uint16) gcp_lons = gcp_array[:, 2].reshape(ypoints.shape[0], xpoints.shape[0]) gcp_lats = gcp_array[:, 3].reshape(ypoints.shape[0], xpoints.shape[0]) @@ -686,6 +685,13 @@ def get_gcps(self): return (xpoints, ypoints), (gcp_lons, gcp_lats, gcp_alts), (rio_gcps, crs) + def get_bounding_box(self): + """Get the bounding box for the data coverage.""" + (xpoints, ypoints), (gcp_lons, gcp_lats, gcp_alts), (rio_gcps, crs) = self.get_gcps() + bblons = np.hstack((gcp_lons[0, :-1], gcp_lons[:-1, -1], gcp_lons[-1, :1:-1], gcp_lons[:1:-1, 0])) + bblats = np.hstack((gcp_lats[0, :-1], gcp_lats[:-1, -1], gcp_lats[-1, :1:-1], gcp_lats[:1:-1, 0])) + return bblons.tolist(), bblats.tolist() + @property def start_time(self): """Get the start time.""" @@ -733,7 +739,8 @@ def load(self, dataset_keys, **kwargs): gcps = get_gcps_from_array(val) from pyresample.future.geometry import SwathDefinition val.attrs["area"] = SwathDefinition(lonlats["longitude"], lonlats["latitude"], - attrs=dict(gcps=gcps)) + attrs=dict(gcps=gcps, + bounding_box=handler.get_bounding_box())) datasets[key] = val continue return datasets diff --git a/satpy/tests/reader_tests/test_sar_c_safe.py b/satpy/tests/reader_tests/test_sar_c_safe.py index 531305798b..0d88190a5d 100644 --- a/satpy/tests/reader_tests/test_sar_c_safe.py +++ b/satpy/tests/reader_tests/test_sar_c_safe.py @@ -215,22 +215,28 @@ def test_read_calibrated_natural(self, measurement_filehandler): calibration = Calibration.sigma_nought xarr = measurement_filehandler.get_dataset(DataQuery(name="measurement", polarization="vv", calibration=calibration, quantity="natural"), info=dict()) - expected = np.array([[np.nan, 0.02707529], [2.55858416, 3.27611055]]) + expected = np.array([[np.nan, 0.02707529], [2.55858416, 3.27611055]], dtype=np.float32) np.testing.assert_allclose(xarr.values[:2, :2], expected, rtol=2e-7) + assert xarr.dtype == np.float32 + assert xarr.compute().dtype == np.float32 def test_read_calibrated_dB(self, measurement_filehandler): """Test the calibration routines.""" calibration = Calibration.sigma_nought xarr = measurement_filehandler.get_dataset(DataQuery(name="measurement", polarization="vv", calibration=calibration, quantity="dB"), info=dict()) - expected = np.array([[np.nan, -15.674268], [4.079997, 5.153585]]) - np.testing.assert_allclose(xarr.values[:2, :2], expected) + expected = np.array([[np.nan, -15.674268], [4.079997, 5.153585]], dtype=np.float32) + np.testing.assert_allclose(xarr.values[:2, :2], expected, rtol=1e-6) + assert xarr.dtype == np.float32 + assert xarr.compute().dtype == np.float32 def test_read_lon_lats(self, measurement_filehandler): """Test reading lons and lats.""" query = DataQuery(name="longitude", polarization="vv") xarr = measurement_filehandler.get_dataset(query, info=dict()) np.testing.assert_allclose(xarr.values, expected_longitudes) + assert xarr.dtype == np.float64 + assert xarr.compute().dtype == np.float64 annotation_xml = b""" @@ -702,6 +708,8 @@ def test_get_noise_dataset(self, noise_filehandler): query = DataQuery(name="noise", polarization="vv") res = noise_filehandler.get_dataset(query, {}) np.testing.assert_allclose(res, self.expected_azimuth_noise * self.expected_range_noise) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_get_noise_dataset_has_right_chunk_size(self, noise_filehandler): """Test using get_dataset for the noise has right chunk size in result.""" @@ -724,12 +732,16 @@ def test_dn_calibration_array(self, calibration_filehandler): expected_dn = np.ones((10, 10)) * 1087 res = calibration_filehandler.get_calibration(Calibration.dn, chunks=5) np.testing.assert_allclose(res, expected_dn) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_beta_calibration_array(self, calibration_filehandler): """Test reading the beta calibration array.""" expected_beta = np.ones((10, 10)) * 1087 res = calibration_filehandler.get_calibration(Calibration.beta_nought, chunks=5) np.testing.assert_allclose(res, expected_beta) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_sigma_calibration_array(self, calibration_filehandler): """Test reading the sigma calibration array.""" @@ -737,18 +749,23 @@ def test_sigma_calibration_array(self, calibration_filehandler): 1277.968, 1277.968, 1277.968, 1277.968]]) * np.ones((10, 1)) res = calibration_filehandler.get_calibration(Calibration.sigma_nought, chunks=5) np.testing.assert_allclose(res, expected_sigma) - + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_gamma_calibration_array(self, calibration_filehandler): """Test reading the gamma calibration array.""" res = calibration_filehandler.get_calibration(Calibration.gamma, chunks=5) np.testing.assert_allclose(res, self.expected_gamma) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_get_calibration_dataset(self, calibration_filehandler): """Test using get_dataset for the calibration.""" query = DataQuery(name="gamma", polarization="vv") res = calibration_filehandler.get_dataset(query, {}) np.testing.assert_allclose(res, self.expected_gamma) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_get_calibration_dataset_has_right_chunk_size(self, calibration_filehandler): """Test using get_dataset for the calibration yields array with right chunksize.""" @@ -762,6 +779,7 @@ def test_get_calibration_constant(self, calibration_filehandler): query = DataQuery(name="calibration_constant", polarization="vv") res = calibration_filehandler.get_dataset(query, {}) assert res == 1 + assert type(res) is np.float32 def test_incidence_angle(annotation_filehandler): @@ -769,6 +787,8 @@ def test_incidence_angle(annotation_filehandler): query = DataQuery(name="incidence_angle", polarization="vv") res = annotation_filehandler.get_dataset(query, {}) np.testing.assert_allclose(res, 19.18318046) + assert res.dtype == np.float32 + assert res.compute().dtype == np.float32 def test_reading_from_reader(measurement_file, calibration_file, noise_file, annotation_file): @@ -787,7 +807,9 @@ def test_reading_from_reader(measurement_file, calibration_file, noise_file, ann array = dataset_dict["measurement"] np.testing.assert_allclose(array.attrs["area"].lons, expected_longitudes) expected_db = np.array([[np.nan, -15.674268], [4.079997, 5.153585]]) - np.testing.assert_allclose(array.values[:2, :2], expected_db) + np.testing.assert_allclose(array.values[:2, :2], expected_db, rtol=1e-6) + assert array.dtype == np.float32 + assert array.compute().dtype == np.float32 def test_filename_filtering_from_reader(measurement_file, calibration_file, noise_file, annotation_file, tmp_path): @@ -814,7 +836,7 @@ def test_filename_filtering_from_reader(measurement_file, calibration_file, nois pytest.fail(str(err)) -def test_swath_def_contains_gcps(measurement_file, calibration_file, noise_file, annotation_file): +def test_swath_def_contains_gcps_and_bounding_box(measurement_file, calibration_file, noise_file, annotation_file): """Test reading using the reader defined in the config.""" with open(Path(PACKAGE_CONFIG_PATH) / "readers" / "sar-c_safe.yaml") as fd: config = yaml.load(fd, Loader=yaml.UnsafeLoader) @@ -829,3 +851,4 @@ def test_swath_def_contains_gcps(measurement_file, calibration_file, noise_file, dataset_dict = reader.load([query]) array = dataset_dict["measurement"] assert array.attrs["area"].attrs["gcps"] is not None + assert array.attrs["area"].attrs["bounding_box"] is not None