From 879893eb2cda176d48a1ec1f728a6f3690686348 Mon Sep 17 00:00:00 2001 From: Martin Raspaud Date: Fri, 11 Oct 2024 12:44:06 +0200 Subject: [PATCH 1/5] Use spline interpolation for faster processing This requires https://github.com/pytroll/python-geotiepoints/pull/85 to be merged and released. --- satpy/readers/sar_c_safe.py | 13 ++- satpy/readers/sgli_l1b.py | 4 +- satpy/tests/reader_tests/test_sar_c_safe.py | 121 ++++---------------- 3 files changed, 33 insertions(+), 105 deletions(-) diff --git a/satpy/readers/sar_c_safe.py b/satpy/readers/sar_c_safe.py index 986440759a..0d42491cbd 100644 --- a/satpy/readers/sar_c_safe.py +++ b/satpy/readers/sar_c_safe.py @@ -51,7 +51,7 @@ import xarray as xr from dask import array as da from geotiepoints.geointerpolator import lonlat2xyz, xyz2lonlat -from geotiepoints.interpolator import MultipleGridInterpolator +from geotiepoints.interpolator import MultipleSplineInterpolator from xarray import DataArray from satpy.dataset.data_dict import DatasetDict @@ -636,10 +636,15 @@ def _get_lonlatalts_uncached(self): fine_points = [np.arange(size) for size in shape] x, y, z = lonlat2xyz(gcp_lons, gcp_lats) - interpolator = MultipleGridInterpolator((ypoints, xpoints), x, y, z, gcp_alts) - hx, hy, hz, altitudes = interpolator.interpolate(fine_points, method="cubic", chunks=self.chunks) - longitudes, latitudes = xyz2lonlat(hx, hy, hz) + kx = 2 + ky = 2 + + interpolator = MultipleSplineInterpolator((ypoints, xpoints), x, y, z, gcp_alts, kx=kx, ky=ky) + hx, hy, hz, altitudes = interpolator.interpolate(fine_points, chunks=self.chunks) + + + longitudes, latitudes = xyz2lonlat(hx, hy, hz) altitudes = xr.DataArray(altitudes, dims=["y", "x"]) longitudes = xr.DataArray(longitudes, dims=["y", "x"]) latitudes = xr.DataArray(latitudes, dims=["y", "x"]) diff --git a/satpy/readers/sgli_l1b.py b/satpy/readers/sgli_l1b.py index 079f93d2f3..f22f77b03a 100644 --- a/satpy/readers/sgli_l1b.py +++ b/satpy/readers/sgli_l1b.py @@ -175,7 +175,7 @@ def get_lon_lats(self, key): def interpolate_spherical(self, azimuthal_angle, polar_angle, resampling_interval): """Interpolate spherical coordinates.""" - from geotiepoints.geointerpolator import GeoGridInterpolator + from geotiepoints.geointerpolator import GeoSplineInterpolator full_shape = (self.h5file["Image_data"].attrs["Number_of_lines"], self.h5file["Image_data"].attrs["Number_of_pixels"]) @@ -183,7 +183,7 @@ def interpolate_spherical(self, azimuthal_angle, polar_angle, resampling_interva tie_lines = np.arange(0, polar_angle.shape[0] * resampling_interval, resampling_interval) tie_cols = np.arange(0, polar_angle.shape[1] * resampling_interval, resampling_interval) - interpolator = GeoGridInterpolator((tie_lines, tie_cols), azimuthal_angle, polar_angle, method="slinear") + interpolator = GeoSplineInterpolator((tie_lines, tie_cols), azimuthal_angle, polar_angle, kx=2, ky=2) new_azi, new_pol = interpolator.interpolate_to_shape(full_shape, chunks="auto") return new_azi, new_pol diff --git a/satpy/tests/reader_tests/test_sar_c_safe.py b/satpy/tests/reader_tests/test_sar_c_safe.py index 9e24c00c4e..26ff603b10 100644 --- a/satpy/tests/reader_tests/test_sar_c_safe.py +++ b/satpy/tests/reader_tests/test_sar_c_safe.py @@ -174,102 +174,26 @@ def measurement_filehandler(measurement_file, noise_filehandler, calibration_fil -expected_longitudes = np.array([[3.79492915e-16, 5.91666667e-01, 9.09722222e-01, - 1.00000000e+00, 9.08333333e-01, 6.80555556e-01, - 3.62500000e-01, 8.32667268e-17, -3.61111111e-01, - -6.75000000e-01, -8.95833333e-01, -9.77777778e-01, - -8.75000000e-01, -5.41666667e-01, 6.80555556e-02, - 1.00000000e+00], - [1.19166667e+00, 1.32437500e+00, 1.36941964e+00, - 1.34166667e+00, 1.25598214e+00, 1.12723214e+00, - 9.70282738e-01, 8.00000000e-01, 6.31250000e-01, - 4.78898810e-01, 3.57812500e-01, 2.82857143e-01, - 2.68898810e-01, 3.30803571e-01, 4.83437500e-01, - 7.41666667e-01], - [1.82638889e+00, 1.77596726e+00, 1.72667765e+00, - 1.67757937e+00, 1.62773172e+00, 1.57619402e+00, - 1.52202558e+00, 1.46428571e+00, 1.40203373e+00, - 1.33432894e+00, 1.26023065e+00, 1.17879819e+00, - 1.08909084e+00, 9.90167942e-01, 8.81088790e-01, - 7.60912698e-01], - [2.00000000e+00, 1.99166667e+00, 1.99305556e+00, - 2.00000000e+00, 2.00833333e+00, 2.01388889e+00, - 2.01250000e+00, 2.00000000e+00, 1.97222222e+00, - 1.92500000e+00, 1.85416667e+00, 1.75555556e+00, - 1.62500000e+00, 1.45833333e+00, 1.25138889e+00, - 1.00000000e+00], - [1.80833333e+00, 2.01669643e+00, 2.18011267e+00, - 2.30119048e+00, 2.38253827e+00, 2.42676446e+00, - 2.43647747e+00, 2.41428571e+00, 2.36279762e+00, - 2.28462160e+00, 2.18236607e+00, 2.05863946e+00, - 1.91605017e+00, 1.75720663e+00, 1.58471726e+00, - 1.40119048e+00], - [1.34722222e+00, 1.89627976e+00, 2.29940830e+00, - 2.57341270e+00, 2.73509779e+00, 2.80126842e+00, - 2.78872945e+00, 2.71428571e+00, 2.59474206e+00, - 2.44690334e+00, 2.28757440e+00, 2.13356009e+00, - 2.00166525e+00, 1.90869473e+00, 1.87145337e+00, - 1.90674603e+00], - [7.12500000e-01, 1.67563988e+00, 2.36250177e+00, - 2.80892857e+00, 3.05076318e+00, 3.12384850e+00, - 3.06402742e+00, 2.90714286e+00, 2.68903770e+00, - 2.44555485e+00, 2.21253720e+00, 2.02582766e+00, - 1.92126913e+00, 1.93470451e+00, 2.10197669e+00, - 2.45892857e+00], - [5.55111512e-16, 1.40000000e+00, 2.38095238e+00, - 3.00000000e+00, 3.31428571e+00, 3.38095238e+00, - 3.25714286e+00, 3.00000000e+00, 2.66666667e+00, - 2.31428571e+00, 2.00000000e+00, 1.78095238e+00, - 1.71428571e+00, 1.85714286e+00, 2.26666667e+00, - 3.00000000e+00], - [-6.94444444e-01, 1.11458333e+00, 2.36631944e+00, - 3.13888889e+00, 3.51041667e+00, 3.55902778e+00, - 3.36284722e+00, 3.00000000e+00, 2.54861111e+00, - 2.08680556e+00, 1.69270833e+00, 1.44444444e+00, - 1.42013889e+00, 1.69791667e+00, 2.35590278e+00, - 3.47222222e+00], - [-1.27500000e+00, 8.64613095e-01, 2.33016227e+00, - 3.21785714e+00, 3.62390731e+00, 3.64452239e+00, - 3.37591199e+00, 2.91428571e+00, 2.35585317e+00, - 1.79682398e+00, 1.33340774e+00, 1.06181406e+00, - 1.07825255e+00, 1.47893282e+00, 2.36006448e+00, - 3.81785714e+00], - [-1.64583333e+00, 6.95312500e-01, 2.28404018e+00, - 3.22916667e+00, 3.63950893e+00, 3.62388393e+00, - 3.29110863e+00, 2.75000000e+00, 2.10937500e+00, - 1.47805060e+00, 9.64843750e-01, 6.78571429e-01, - 7.28050595e-01, 1.22209821e+00, 2.26953125e+00, - 3.97916667e+00], - [-1.71111111e+00, 6.51904762e-01, 2.23951247e+00, - 3.16507937e+00, 3.54197279e+00, 3.48356009e+00, - 3.10320862e+00, 2.51428571e+00, 1.83015873e+00, - 1.16419501e+00, 6.29761905e-01, 3.40226757e-01, - 4.08956916e-01, 9.49319728e-01, 2.07468254e+00, - 3.89841270e+00], - [-1.37500000e+00, 7.79613095e-01, 2.20813846e+00, - 3.01785714e+00, 3.31605017e+00, 3.20999858e+00, - 2.80698342e+00, 2.21428571e+00, 1.53918651e+00, - 8.88966837e-01, 3.70907738e-01, 9.22902494e-02, - 1.60395408e-01, 6.82504252e-01, 1.76589782e+00, - 3.51785714e+00], - [-5.41666667e-01, 1.12366071e+00, 2.20147747e+00, - 2.77976190e+00, 2.94649235e+00, 2.78964711e+00, - 2.39720451e+00, 1.85714286e+00, 1.25744048e+00, - 6.86075680e-01, 2.31026786e-01, -1.97278912e-02, - 2.17899660e-02, 4.43558673e-01, 1.33355655e+00, - 2.77976190e+00], - [8.84722222e-01, 1.72927083e+00, 2.23108879e+00, - 2.44305556e+00, 2.41805060e+00, 2.20895337e+00, - 1.86864335e+00, 1.45000000e+00, 1.00590278e+00, - 5.89231151e-01, 2.52864583e-01, 4.96825397e-02, - 3.25644841e-02, 2.54389881e-01, 7.68038194e-01, - 1.62638889e+00], - [3.00000000e+00, 2.64166667e+00, 2.30853175e+00, - 2.00000000e+00, 1.71547619e+00, 1.45436508e+00, - 1.21607143e+00, 1.00000000e+00, 8.05555556e-01, - 6.32142857e-01, 4.79166667e-01, 3.46031746e-01, - 2.32142857e-01, 1.36904762e-01, 5.97222222e-02, - 0.00000000e+00]]) +expected_longitudes = np.array([[-0., 0.54230055, 0.87563228, 1., 0.91541479, + 0.62184442, 0.26733714, -0., -0.18015287, -0.27312165], + [1.0883956 , 1.25662247, 1.34380634, 1.34995884, 1.2750712 , + 1.11911385, 0.9390845 , 0.79202785, 0.67796547, 0.59691204], + [1.75505196, 1.74123364, 1.71731849, 1.68330292, 1.63918145, + 1.58494674, 1.52376394, 1.45880655, 1.39007883, 1.31758574], + [2., 1.99615628, 1.99615609, 2., 2.00768917, + 2.0192253 , 2.02115051, 2. , 1.95576762, 1.88845002], + [1.82332931, 2.02143515, 2.18032829, 2.30002491, 2.38053511, + 2.4218612 , 2.43113105, 2.41546985, 2.37487052, 2.3093278 ], + [1.22479001, 1.81701462, 2.26984318, 2.58335874, 2.75765719, + 2.79279164, 2.75366973, 2.70519769, 2.64737395, 2.58019762], + [0.51375081, 1.53781389, 2.3082042 , 2.82500549, 3.0885147 , + 3.09893859, 2.98922885, 2.89232293, 2.8082302 , 2.7369586 ], + [0., 1.33889733, 2.33891557, 3., 3.32266837, + 3.30731797, 3.1383157 , 3., 2.8923933 , 2.81551297], + [-0.31638932, 1.22031759, 2.36197571, 3.10836734, 3.46019271, + 3.41800603, 3.20098223, 3.02826595, 2.89989242, 2.81588745], + [-0.43541441, 1.18211505, 2.37738272, 3.1501186 , 3.50112948, + 3.43104055, 3.17724665, 2.97712796, 2.83072911, 2.73808164]]) class Calibration(Enum): @@ -304,8 +228,7 @@ 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()) - expected = expected_longitudes - np.testing.assert_allclose(xarr.values, expected[:10, :10], atol=1e-3) + np.testing.assert_allclose(xarr.values, expected_longitudes) annotation_xml = b""" @@ -860,7 +783,7 @@ def test_reading_from_reader(measurement_file, calibration_file, noise_file, ann query = DataID(reader._id_keys, **query.to_dict()) dataset_dict = reader.load([query]) array = dataset_dict["measurement"] - np.testing.assert_allclose(array.attrs["area"].lons, expected_longitudes[:10, :10], atol=1e-3) + 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) From 98762e4375368757f2d4d496fe3038439d49e942 Mon Sep 17 00:00:00 2001 From: Martin Raspaud Date: Fri, 11 Oct 2024 12:47:29 +0200 Subject: [PATCH 2/5] Fix style --- satpy/tests/reader_tests/test_sar_c_safe.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/satpy/tests/reader_tests/test_sar_c_safe.py b/satpy/tests/reader_tests/test_sar_c_safe.py index 26ff603b10..3ec3aa9577 100644 --- a/satpy/tests/reader_tests/test_sar_c_safe.py +++ b/satpy/tests/reader_tests/test_sar_c_safe.py @@ -174,8 +174,8 @@ def measurement_filehandler(measurement_file, noise_filehandler, calibration_fil -expected_longitudes = np.array([[-0., 0.54230055, 0.87563228, 1., 0.91541479, - 0.62184442, 0.26733714, -0., -0.18015287, -0.27312165], +expected_longitudes = np.array([[-0., 0.54230055, 0.87563228, 1., 0.91541479, + 0.62184442, 0.26733714, -0., -0.18015287, -0.27312165], [1.0883956 , 1.25662247, 1.34380634, 1.34995884, 1.2750712 , 1.11911385, 0.9390845 , 0.79202785, 0.67796547, 0.59691204], [1.75505196, 1.74123364, 1.71731849, 1.68330292, 1.63918145, From 63de5e07965b7ae3d421722baa52e07d83b9cbad Mon Sep 17 00:00:00 2001 From: Martin Raspaud Date: Mon, 14 Oct 2024 11:20:56 +0200 Subject: [PATCH 3/5] Run sar tests only with geotiepoints >= 1.7.5 --- satpy/tests/reader_tests/test_sar_c_safe.py | 1 + 1 file changed, 1 insertion(+) diff --git a/satpy/tests/reader_tests/test_sar_c_safe.py b/satpy/tests/reader_tests/test_sar_c_safe.py index 3ec3aa9577..f7191a951a 100644 --- a/satpy/tests/reader_tests/test_sar_c_safe.py +++ b/satpy/tests/reader_tests/test_sar_c_safe.py @@ -32,6 +32,7 @@ from satpy.readers.sar_c_safe import Calibrator, Denoiser, SAFEXMLAnnotation rasterio = pytest.importorskip("rasterio") +geotiepoints = pytest.importorskip("geotiepoints", "1.7.5") dirname_suffix = "20190201T024655_20190201T024720_025730_02DC2A_AE07" From ee7edefbe57be93c56db5ed1d32a4a9e170b8359 Mon Sep 17 00:00:00 2001 From: Martin Raspaud Date: Mon, 14 Oct 2024 11:45:30 +0200 Subject: [PATCH 4/5] Do importskip before importing satpy --- satpy/tests/reader_tests/test_sar_c_safe.py | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/satpy/tests/reader_tests/test_sar_c_safe.py b/satpy/tests/reader_tests/test_sar_c_safe.py index f7191a951a..531305798b 100644 --- a/satpy/tests/reader_tests/test_sar_c_safe.py +++ b/satpy/tests/reader_tests/test_sar_c_safe.py @@ -26,13 +26,14 @@ import pytest import yaml -from satpy._config import PACKAGE_CONFIG_PATH -from satpy.dataset import DataQuery -from satpy.dataset.dataid import DataID -from satpy.readers.sar_c_safe import Calibrator, Denoiser, SAFEXMLAnnotation +geotiepoints = pytest.importorskip("geotiepoints", "1.7.5") + +from satpy._config import PACKAGE_CONFIG_PATH # noqa: E402 +from satpy.dataset import DataQuery # noqa: E402 +from satpy.dataset.dataid import DataID # noqa: E402 +from satpy.readers.sar_c_safe import Calibrator, Denoiser, SAFEXMLAnnotation # noqa: E402 rasterio = pytest.importorskip("rasterio") -geotiepoints = pytest.importorskip("geotiepoints", "1.7.5") dirname_suffix = "20190201T024655_20190201T024720_025730_02DC2A_AE07" From 12ebe3feab2d9f7b92bca8eabe059776291500cf Mon Sep 17 00:00:00 2001 From: Martin Raspaud Date: Tue, 15 Oct 2024 22:02:13 +0200 Subject: [PATCH 5/5] Update satpy/readers/sar_c_safe.py --- satpy/readers/sar_c_safe.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/satpy/readers/sar_c_safe.py b/satpy/readers/sar_c_safe.py index 0d42491cbd..b6d84e8fb7 100644 --- a/satpy/readers/sar_c_safe.py +++ b/satpy/readers/sar_c_safe.py @@ -637,10 +637,8 @@ def _get_lonlatalts_uncached(self): fine_points = [np.arange(size) for size in shape] x, y, z = lonlat2xyz(gcp_lons, gcp_lats) - kx = 2 - ky = 2 - interpolator = MultipleSplineInterpolator((ypoints, xpoints), x, y, z, gcp_alts, kx=kx, ky=ky) + interpolator = MultipleSplineInterpolator((ypoints, xpoints), x, y, z, gcp_alts, kx=2, ky=2) hx, hy, hz, altitudes = interpolator.interpolate(fine_points, chunks=self.chunks)