Phase Unwrapping

phase unwrapping

import zarr
import moraine as mr

---

source

gamma_mcf_pt

 gamma_mcf_pt (pc_x:numpy.ndarray, pc_y:numpy.ndarray, ph:numpy.ndarray,
               ph_weight:numpy.ndarray=None, ref_point:int=1)

A simple wrapper for mcf_pt in GAMMA software.

	Type	Default	Details
pc_x	ndarray		x coordinate, shape of (N,)
pc_y	ndarray		y coordinate, shape of (N,)
ph	ndarray		wrapped phase, shape of (N,) or (N,M)
ph_weight	ndarray	None	point weight, shape of (N,) or (N,M), optional
ref_point	int	1	reference point, the first point by default
Returns	ndarray		unwrapped phase, shape of (N,) or (N,M)

usage:

rslc = zarr.open('../../data/rslc.zarr/','r')[:]

# SHP selection
az_half_win = 5; r_half_win = 5
az_win = 2*az_half_win+1; r_win = 2*r_half_win+1

rmli = np.abs(rslc)**2
p = mr.ks_test(rmli,az_half_win=az_half_win,r_half_win=r_half_win)
is_shp = p < 0.05

# Select DS candidate
shp_num = np.count_nonzero(is_shp,axis=(-2,-1))
is_ds_can = shp_num >= 50

ds_can_is_shp = is_shp[is_ds_can]
ds_can_gix = np.stack(np.where(is_ds_can),axis=-1)
ds_can_ph, ds_can_emi_quality, ds_can_t_coh = mr.emperical_co_emi_temp_coh_pc(rslc,ds_can_gix,ds_can_is_shp,batch_size=1000)

_is_ds_can_refined = (ds_can_emi_quality>=1.0) & (ds_can_emi_quality <1.2) & (ds_can_t_coh > 0.7) & (ds_can_t_coh <= 1.0)

ds_gix = ds_can_gix[_is_ds_can_refined]
ds_ph = ds_can_ph[_is_ds_can_refined]

ds_hix = mr.pc_hix(ds_gix,shape=rslc.shape[:2])
hix_sort_key = mr.pc_sort(ds_hix)

ds_hix = ds_hix[hix_sort_key]
ds_gix = ds_gix[hix_sort_key]
ds_ph = ds_ph[hix_sort_key]

e = zarr.open('../../data/e.zarr/','r')[:]
n = zarr.open('../../data/n.zarr/','r')[:]

ds_e = mr.ras2pc(e,ds_gix)
ds_n = mr.ras2pc(n,ds_gix)

ds_ph = ds_ph*ds_ph[:,0:1].conj() # set the first image as secondary

Unwrap one single interferogram:

ds_ph_ = ds_ph[:,14]*ds_ph[:,9].conj()

unw = gamma_mcf_pt(ds_e, ds_n, ds_ph_)

def wrap(unw):
    return np.mod(unw + np.pi, 2 * np.pi) - np.pi

np.testing.assert_array_almost_equal(wrap(unw),np.angle(ds_ph_),decimal=4)

Unwrap multiple interferograms:

unw2 = gamma_mcf_pt(ds_e, ds_n, ds_ph[:,2:4])
np.testing.assert_array_almost_equal(wrap(unw2),np.angle(ds_ph[:,2:4]),decimal=4)

import holoviews as hv
import datashader as ds
import spatialpandas
import zarr
import numpy as np
from holoviews import opts
from bokeh.models import WheelZoomTool

hv.extension('bokeh')

coordinates = spatialpandas.geometry.PointArray((ds_e, ds_n))
data = spatialpandas.GeoDataFrame({'geometry':coordinates,'phase':unw})
unw_plot = mr.points(data,kdims=['e','n'],pdim='phase',google_earth=True)

coordinates = spatialpandas.geometry.PointArray((ds_e, ds_n))
data = spatialpandas.GeoDataFrame({'geometry':coordinates,'phase':np.angle(ds_ph_)})
ph_plot = mr.points(data,kdims=['e','n'],pdim='phase',prange=(-np.pi,np.pi),google_earth=True)

unw_plot.opts(opts.Image(cmap='colorwheel',width=600, height=400, colorbar=True, xlabel='Longitude', ylabel='Latitude'),
          opts.Points(marker='o',size=10,tools=['hover']))

ph_plot.opts(opts.Image(cmap='colorwheel',width=600, height=400, colorbar=True, xlabel='Longitude', ylabel='Latitude'),
          opts.Points(marker='o',size=10,tools=['hover']))