shapely.errors.TopologicalError during InSAR water masking

[cirrusasf]

The bug

The insar workflow and INSAR_GAMMA jobs in HyP3 may fail with this stack trace when generating a water mask prior to phase unwrapping:

  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/__main__.py", line 38, in main
    process_entry_point.load()()
  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/__main__.py", line 185, in insar
    product_name = insar_sentinel_gamma(
  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/insar/ifm_sentinel.py", line 427, in insar_sentinel_gamma
    coords, ref_point_info = unwrapping_geocoding(reference, secondary, step="man", rlooks=rlooks, alooks=alooks,
  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/insar/unwrapping_geocoding.py", line 253, in unwrapping_geocoding
    get_water_mask(f"{ifgname}.adf.cc", width, lt, demw, demn, dempar)
  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/insar/unwrapping_geocoding.py", line 150, in get_water_mask
    create_water_mask(f'{temp_dir}/tmpgtiff_mask_geo.tif', 'water_mask.tif')
  File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/water_mask.py", line 71, in create_water_mask
    mask = gpd.clip(mask, envelope)
  File "/usr/local/lib/python3.10/dist-packages/geopandas/tools/clip.py", line 192, in clip
    clipped = _clip_gdf_with_mask(gdf, combined_mask)
  File "/usr/local/lib/python3.10/dist-packages/geopandas/tools/clip.py", line 69, in _clip_gdf_with_mask
    ] = gdf_sub.geometry.values[non_point_mask].intersection(mask)
  File "/usr/local/lib/python3.10/dist-packages/geopandas/array.py", line 671, in intersection
    self._binary_method("intersection", self, other), crs=self.crs
  File "/usr/local/lib/python3.10/dist-packages/geopandas/array.py", line 611, in _binary_method
    return getattr(vectorized, op)(left._data, right, **kwargs)
  File "/usr/local/lib/python3.10/dist-packages/geopandas/_vectorized.py", line 942, in intersection
    return _binary_geo("intersection", data, other)
  File "/usr/local/lib/python3.10/dist-packages/geopandas/_vectorized.py", line 306, in _binary_geo
    data[:] = [
  File "/usr/local/lib/python3.10/dist-packages/geopandas/_vectorized.py", line 307, in <listcomp>
    getattr(s, op)(right) if s is not None and right is not None else None
  File "/usr/lib/python3/dist-packages/shapely/geometry/base.py", line 689, in intersection
    return geom_factory(self.impl['intersection'](self, other))
  File "/usr/lib/python3/dist-packages/shapely/topology.py", line 72, in __call__
    self._check_topology(err, this, other)
  File "/usr/lib/python3/dist-packages/shapely/topology.py", line 37, in _check_topology
    raise TopologicalError(
shapely.errors.TopologicalError: The operation 'GEOSIntersection_r' could not be performed. Likely cause is invalidity of the geometry <shapely.geometry.multipolygon.MultiPolygon object at 0x7f8d9ef30100>

To date, we've seen this error for these InSAR pairs run with looks=10x2 and apply_water_mask=true:

S1A_IW_SLC__1SDV_20170728T051208_20170728T051234_017667_01D93F_8329,S1A_IW_SLC__1SDV_20180711T051213_20180711T051240_022742_0276FD_617A
S1A_IW_SLC__1SDV_20180711T051213_20180711T051240_022742_0276FD_617A,S1B_IW_SLC__1SDV_20190712T051136_20190712T051204_017096_0202A5_C8F2
S1A_IW_SLC__1SDV_20210707T051231_20210707T051258_038667_04901E_08FA,S1A_IW_SLC__1SDV_20220702T051237_20220702T051304_043917_053E27_9215
S1A_IW_SLC__1SDV_20220702T051237_20220702T051304_043917_053E27_9215,S1A_IW_SLC__1SDV_20230615T051241_20230615T051308_048992_05E43B_4288
S1B_IW_SLC__1SDV_20190712T051136_20190712T051204_017096_0202A5_C8F2,S1B_IW_SLC__1SDV_20200730T051144_20200730T051212_022696_02B139_7027
S1B_IW_SLC__1SDV_20200718T051143_20200718T051211_022521_02ABE5_9C44,S1A_IW_SLC__1SDV_20210719T051232_20210719T051259_038842_04955D_7838
S1B_IW_SLC__1SDV_20200718T051143_20200718T051211_022521_02ABE5_9C44,S1A_IW_SLC__1SDV_20210731T051233_20210731T051259_039017_049A8F_A8F7
S1B_IW_SLC__1SDV_20200718T051143_20200718T051211_022521_02ABE5_9C44,S1B_IW_SLC__1SDV_20210713T051149_20210713T051217_027771_03506A_A8C6
S1B_IW_SLC__1SDV_20200718T051143_20200718T051211_022521_02ABE5_9C44,S1B_IW_SLC__1SDV_20210725T051149_20210725T051217_027946_03559A_B4F3
S1B_IW_SLC__1SDV_20200730T051144_20200730T051212_022696_02B139_7027,S1A_IW_SLC__1SDV_20210707T051231_20210707T051258_038667_04901E_08FA

This issue was first observed on 2023-11-16 and was likely introduced with the water masking changes in hyp3-gamma v6.4.1.

Per @cirrusasf :
Version does not work for some pairs, for example, The error that hyp3-gamma throw outs is related to water_mask code. THe error: File "/usr/local/lib/python3.10/dist-packages/hyp3_gamma/water_mask.py", line 71, in create_water_mask
mask = gpd.clip(mask, envelope)
It turns out the that geopandas.clip() will fail if both mask and envelope are in projection coordinates and the mask is very large. If we change them into wgs84 coordinates, it works.