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/',mode='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/',mode='r')[:]
n = zarr.open('../../data/n.zarr/',mode='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 numpy as np
from holoviews import opts
from bokeh.models import WheelZoomTool

hv.extension('bokeh')

ⓘ

unw_plot = mr.pc_plot(unw,ds_n, ds_e)
ph_plot = mr.pc_plot(np.angle(ds_ph_), ds_n, ds_e)

unw_plot = unw_plot.redim(
    x=hv.Dimension('lon', label='Longitude'), y=hv.Dimension('lat',label='Latitude'),
    z=hv.Dimension('unw',label='Unwrapped phase')
)
unw_plot.opts(
    opts.Image(cmap='colorwheel',width=600, height=400, colorbar=True,
               default_tools=['pan',WheelZoomTool(zoom_on_axis=False),'save','reset','hover'],
               active_tools=['wheel_zoom']
              ),
    opts.Points(color='unw', cmap='colorwheel',width=600, height=400, colorbar=True,
                default_tools=['pan',WheelZoomTool(zoom_on_axis=False),'save','reset','hover'],
                active_tools=['wheel_zoom']
               )
)

ph_plot = ph_plot.redim(
    x=hv.Dimension('lon', label='Longitude'), y=hv.Dimension('lat',label='Latitude'),
    z=hv.Dimension('ph',label='Wrapped phase')
)
ph_plot.opts(
    opts.Image(cmap='colorwheel',width=600, height=400, colorbar=True,
               default_tools=['pan',WheelZoomTool(zoom_on_axis=False),'save','reset','hover'],
               active_tools=['wheel_zoom']
              ),
    opts.Points(color='ph', cmap='colorwheel',width=600, height=400, colorbar=True,
                default_tools=['pan',WheelZoomTool(zoom_on_axis=False),'save','reset','hover'],
                active_tools=['wheel_zoom']
               )
)