[Feature]: OpenStreetMap Roads Mask

[ClassiCcOwl]

In S1A_Stack_CPGF_T173_roads notebook there is OpenStreetMap Roads , when in try to achive same result with PyGMTSAR 2 i get stuck in part where

sbas.unwrap_parallel(pairs, mask=osmmask)

can we have a updated version for that notebook with PYGMTSAR

[AlexeyPechnikov]

You’re trying to use outdated functions. I previously shared an example for sparse data processing, see the explanation at https://www.patreon.com/posts/pygmtsar-106488442 However, things are much simpler now.

[ClassiCcOwl]

I have phases and corrs gecoded using sbas.ra2ll() which is is provided in v2, and then i was able to get osmmask = sbas.ll2ra(osmmask_ll) ,i only need help on how to apply them together, in older version you have used masked paramater in unwrap_parallel function, but i couldn't find the equal function or parameter in new version.

[AlexeyPechnikov]

You can apply the mask in radar coordinates to your phase and correlation grids. After that, process them as usual.

[ClassiCcOwl]

You can apply the mask in radar coordinates to your phase and correlation grids. After that, process them as usual.

How can i apply that mask?

this is my corr and phase code.
data = sbas.open_data()
intensity = sbas.multilooking(np.square(np.abs(data)), wavelength=400, coarsen=(1,4))
phase = sbas.multilooking(sbas.phasediff(pairs), wavelength=400, coarsen=(1,4))
corr = sbas.correlation(phase, intensity)
intf_filt = sbas.interferogram(sbas.goldstein(phase, corr, 32))
decimator = sbas.decimator()
tqdm_dask(result := dask.persist(decimator(corr), decimator(intf_filt)), desc='Compute Phase and Correlation')
corr60m, intf60m = result

tqdm_dask(phasefilts_ll := sbas.ra2ll(intf60m), desc='Geocoding')
tqdm_dask(corrs_ll := sbas.ra2ll(corr60m), desc='Geocoding')

and this is my mask code

def dload_osm_roads(lon_min, lat_min, lon_max, lat_max):
# OpenStreetMap graph calculations
import osmnx
# graph calculations
import networkx as nx

#north, south, east, west
G = osmnx.graph_from_bbox(lat_max, lat_min, lon_min, lon_max, network_type='drive')
osmroads = osmnx.graph_to_gdfs(G, nodes=False, edges=True)
return osmroads

workaround to run osmnx code in a separate process and fix the related Dask library error:

TypeError: _config_dns.._getaddrinfo() got an unexpected keyword argument 'family'

pathos isolate it well when joblib and multiprocessing do not

import pathos.multiprocessing as mp
with mp.Pool(processes=1) as pool:
osmroads = pool.apply(dload_osm_roads, args=(lon_min, lat_min, lon_max, lat_max))

add 90m buffer

osmroads_buffer = osmroads.to_crs("EPSG:32611").buffer(90).to_crs("EPSG:4326")

osmroads_buffer = osmroads.buffer(90)

osmmask_ll = sbas.as_geo(xr.ones_like(phasefilts_ll[0]).compute().rename({'lon': 'x', 'lat': 'y'}))
osmmask_ll = osmmask_ll.rio.clip(osmroads_buffer).rename({'x': 'lon', 'y': 'lat'})

convert to radar coordinates

osmmask = sbas.ll2ra(osmmask_ll)

[AlexeyPechnikov]

Your osmmask should be the raster mask in radar coordinates. You might want to plot it to check. Simply multiply it by your interferogram or correlation raster, or use a condition like intf60m.where(np.isfinite(osmmask)) or intf60m.where(osmmask > 0) depending on the mask datatype (for float or boolean/binary masks).

[ClassiCcOwl]

when i use

sbas.plot_interferograms(intf60m.where(osmmask > 0), cols=3, size=3, caption='Phase, [rad]')
plt.savefig('Phase masked, [rad].jpg')

it never stops calculating

can you update the notebook with new version please?

[AlexeyPechnikov]

If your osmmask is a lazy object, materialize it with osmmask.compute() or by using persist() with a progress bar, or by using PyGMTSAR’s sync_* functions. Additionally, you can use a materialized intf60m if you still encounter performance issues.

[ClassiCcOwl]

Thanks i will try them

[ClassiCcOwl]

still wasn't able to fix it

[AlexeyPechnikov]

You don’t need to focus on ‘INFO’ messages, as these are typically logged for the benefit of Dask developers and are meant for informational purposes only. If seeing these messages is bothersome, you can opt to run the Dask cluster from the command line outside of the notebook and connect to it from within your notebook. This approach helps prevent the display of such messages.

[ClassiCcOwl]

When i use phase.where(osmask) i get no result in plots because it says there is no numeric to plot, apparently the osmmask that i have and phases do not fit or work with each other like they were used to. What alternative method do we have for sbas.unwrap_parallel(pairs, mask=osmmask) in new version of pygmtsar, sir?

[AlexeyPechnikov]

As I mentioned earlier, you need to verify whether your mask is in float or boolean format and apply it accordingly. Please refer to the examples provided above for guidance.