gradient_search leaves part of the destination empty for a partial geostationary disk

[lorenzo-cingano]

Describe the bug

When the source is a partial geostationary disk (for example MSG SEVIRI rapid scan data, which covers only the northern part of the disk, or a region of interest), resampling to an area in a different projection with gradient_search leaves a large part of the destination empty even where the source clearly has data. The populated region is bounded by a smooth arc.

The cause is in AreaSlicer.get_polygon_to_contain (pyresample/slicer.py). For a geostationary source it builds the bounding polygon of the disk and reprojects it into the destination projection to compute the source slice. When the disk is clipped to a partial area extent, the clip introduces straight chord edges (for instance the southern boundary of a rapid scan strip), each represented by only its two end points. Those chords are not straight in the destination projection, but with only two vertices the reprojected polygon keeps them straight and cuts the corner. The resulting source slice is too small, so resample_blocks fills the uncovered destination pixels with the fill value. A full disk source is unaffected, because its whole boundary is the densely sampled disk arc with no straight chord.

To Reproduce

from pyresample import AreaDefinition
from pyresample.slicer import create_slicer

# Partial geostationary disk: the northern 1392 lines of an MSG rapid scan image
source = AreaDefinition(
    "rss", "rss", None,
    {"a": 6378169.0, "b": 6356583.8, "h": 35785831.0, "lon_0": 9.5,
     "proj": "geos", "units": "m"},
    3712, 1392,
    (5568748.2758, 5568748.2758, -5568748.2758, 1392187.0689))

# Polar stereographic destination over Europe, fully inside the disk
euro4 = AreaDefinition(
    "euro4", "euro4", None,
    {"proj": "stere", "ellps": "bessel", "lat_0": 90.0, "lon_0": 14.0,
     "lat_ts": 60.0, "units": "m"},
    1024, 1024,
    (-2717181.7304994687, -5571048.14031214,
     1378818.2695005313, -1475048.1403121399))

x_slice, y_slice = create_slicer(source, euro4).get_slices()
print(y_slice)

Expected behavior

The whole euro4 area falls inside the rapid scan strip, so the source slice should cover every source row the destination needs (down to about row 694), and a gradient_search resample should fill the whole destination.

Actual results

slice(1151, 1335, None)

Only the northern rows are kept (the slice should start near row 694). A gradient_search resample of this source to euro4 then fills only about 44 percent of the destination, with the populated area bounded by an arc.

Environment Info:

• pyresample: reproduced with 1.34.2 and 1.35.0 (the relevant code is unchanged on main)
• shapely: 2.x
• pyproj: 3.x

Additional context

Reported downstream in satpy as pytroll/satpy#3403. A fix that densifies the bounding polygon with segmentize before reprojecting it makes the example above fill the whole destination. Happy to open a pull request.