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.

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 spline interpolation for faster processing #2927

Merged
merged 6 commits into from
Oct 17, 2024
Merged
Show file tree
Hide file tree
Changes from 5 commits
Commits
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
13 changes: 9 additions & 4 deletions satpy/readers/sar_c_safe.py
Original file line number Diff line number Diff line change
Expand Up @@ -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
Expand Down Expand Up @@ -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
mraspaud marked this conversation as resolved.
Show resolved Hide resolved

interpolator = MultipleSplineInterpolator((ypoints, xpoints), x, y, z, gcp_alts, kx=kx, ky=ky)
mraspaud marked this conversation as resolved.
Show resolved Hide resolved
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"])
Expand Down
4 changes: 2 additions & 2 deletions satpy/readers/sgli_l1b.py
Original file line number Diff line number Diff line change
Expand Up @@ -175,15 +175,15 @@ 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"])

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

Expand Down
131 changes: 28 additions & 103 deletions satpy/tests/reader_tests/test_sar_c_safe.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,10 +26,12 @@
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")

Expand Down Expand Up @@ -174,102 +176,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):
Expand Down Expand Up @@ -304,8 +230,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"""<?xml version="1.0" encoding="UTF-8"?>
Expand Down Expand Up @@ -860,7 +785,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)

Expand Down
Loading