Skip to content

Fix bugs in Ionosphere Phase Uncertainty & connected components mask #220

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
oberonia78 wants to merge 8 commits into
base: develop
Choose a base branch
from
Open

Fix bugs in Ionosphere Phase Uncertainty & connected components mask #220

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
8 commits
Select commit Hold shift + click to select a range
e1f1404
fix uncertainty values
Feb 12, 2026
eea9ed0
simplify unwrap error
Feb 12, 2026
5ff452c
fix bug in reading connected components
Feb 12, 2026
adc9c63
update unitest for unwrapping error
Feb 12, 2026
c821d03
revert the changes
Mar 5, 2026
8a1c775
add log instead of print
Mar 5, 2026
e95c1cd
revert changes for unittest
Mar 5, 2026
a9b2c3a
cleaning unittest
Mar 5, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
73 changes: 29 additions & 44 deletions python/packages/isce3/atmosphere/main_band_estimation.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -421,6 +421,7 @@ def estimate_iono_std(
self,
main_coh=None,
side_coh=None,
diff_ms_coh=None,
low_band_coh=None,
high_band_coh=None,
diff_low_high_coh=None,
Expand All @@ -437,6 +438,8 @@ def estimate_iono_std(
coherence of main-band interferogram
side_coh : numpy.ndarray
coherence of side-band interferogram
diff_ms_coh : numpy.ndarray
coherence of difference (main-side) interferogram
low_band_coh : numpy.ndarray
coherence of low subband interferogram
high_band_coh : numpy.ndarray
Expand All @@ -458,9 +461,22 @@ def estimate_iono_std(
sig_nondisp : numpy.ndarray
phase standard deviation of non-dispersive
"""
error_channel = journal.error(
'MainBandIonosphereEstimation.estimate_iono_std')
if side_coh is None and diff_ms_coh is None:
err_str = "estimate_iono_std requires at least one of "\
"`diff_ms_coh` or `side_coh` to be provided."
error_channel.log(err_str)
raise ValueError(err_str)
if main_coh is None:
err_str = "estimate_iono_std requires `main_coh` to be provided."
error_channel.log(err_str)
raise ValueError(err_str)
side_used_coh = diff_ms_coh if diff_ms_coh is not None else side_coh

# resample coherences array of frequency A to
# frequency B grid
if (side_coh is not None) and (resample_flag):
if resample_flag:
if slant_main is None:
slant_main = self.slant_main
if slant_side is None:
Expand All @@ -471,18 +487,18 @@ def estimate_iono_std(
slant_side,
main_coh)

sig_phi_main = np.divide(
np.sqrt(1 - main_coh**2),
main_coh * np.sqrt(2 * number_looks),
out=np.zeros_like(main_coh),
where=main_coh != 0)

# estimate sigma from main- and side- band coherences
if (main_coh is not None) & (side_coh is not None):
sig_phi_main = np.divide(
np.sqrt(1 - main_coh**2),
main_coh / np.sqrt(2 * number_looks),
out=np.zeros_like(main_coh),
where=main_coh != 0)
sig_phi_side = np.divide(
np.sqrt(1 - side_coh**2),
side_coh / np.sqrt(2 * number_looks),
out=np.zeros_like(side_coh),
where=side_coh != 0)
sig_phi_side = np.divide(
np.sqrt(1 - side_used_coh**2),
side_used_coh * np.sqrt(2 * number_looks),
out=np.zeros_like(side_used_coh),
where=side_used_coh != 0)

sig_phi_iono, sig_nondisp = \
self.estimate_sigma(
Expand All @@ -491,37 +507,6 @@ def estimate_iono_std(

return sig_phi_iono, sig_nondisp

def estimate_sigma_main_diff(
self,
sig_phi0,
sig_phi1):
"""Estimate sigma from coherence for main_diff_ms method

Parameters
----------
sig_phi0 : numpy.ndarray
phase standard deviation of main-band interferogram
sig_phi1 : numpy.ndarray
phase standard deviation of side-band interferogram

Returns
-------
sig_iono : numpy.ndarray
2D array of phase standard deviation of dispersive
sig_nondisp : numpy.ndarray
2D array of phase standard deviation of non-dispersive
"""

a = self.f1 / (self.f1 + self.f0)
b = (self.f0 * self.f1) / (self.f0**2 - self.f1**2)
sig_phi01 = np.sqrt(sig_phi0**2 + sig_phi1**2)
sig_iono = np.sqrt(a**2 * sig_phi0**2 + b**2 * sig_phi01**2)

c = self.f0**2 / (self.f0**2 - self.f1**2)
d = self.f0 * self.f1 / (self.f0**2 - self.f1**2)
sig_nondisp = np.sqrt(c**2 * sig_phi0**2 + d**2 * sig_phi01**2)

return sig_iono, sig_nondisp

def estimate_sigma_main_side(
self,
Expand Down Expand Up @@ -669,7 +654,7 @@ def __init__(self,
high_center_freq, method)

self.estimate_iono = estimate_iono_main_diff
self.estimate_sigma = self.estimate_sigma_main_diff
self.estimate_sigma = self.estimate_sigma_main_side
self.compute_unwrap_err = compute_unwrapp_error_main_diff_ms_band


Expand Down
19 changes: 9 additions & 10 deletions python/packages/isce3/atmosphere/split_band_estimation.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -461,6 +461,7 @@ def estimate_iono_std(
self,
main_coh=None,
side_coh=None,
diff_ms_coh=None,
low_band_coh=None,
high_band_coh=None,
diff_low_high_coh=None,
Expand All @@ -477,6 +478,8 @@ def estimate_iono_std(
coherence of main-band interferogram
side_coh : numpy.ndarray
coherence of side-band interferogram
diff_ms_coh : numpy.ndarray
coherence of difference (main-side) interferogram
low_band_coh : numpy.ndarray
coherence of low subband interferogram
high_band_coh : numpy.ndarray
Expand Down Expand Up @@ -817,19 +820,15 @@ def compute_unwrapp_error_main_diff_low_high(
"""

freq_diff = freq_high - freq_low
freq_sum = freq_high + freq_low
freq_multi = freq_high * freq_low

x_coeff = freq_multi / f0 / freq_sum
z_coeff = - freq_multi / 2 / f0 / freq_diff
diff_phase = diff_low_high_runw - \
(freq_diff / f0 * nondisp_array + f0 *
(1/freq_high - 1/freq_low) * disp_array)

diff_unw_coeff = np.round((2 * nondisp_array / z_coeff -
disp_array / z_coeff * x_coeff -
diff_low_high_runw) / 2 / np.pi)
com_unw_coeff = np.round((
(nondisp_array + disp_array) / (2 * x_coeff + 1) -
main_runw / 2) / np.pi)
comm_phase = main_runw - (nondisp_array + disp_array)

com_unw_coeff = np.rint(comm_phase / (2*np.pi)).astype(np.int32)
diff_unw_coeff = np.rint(diff_phase / (2*np.pi)).astype(np.int32)
return com_unw_coeff, diff_unw_coeff


Expand Down
111 changes: 77 additions & 34 deletions python/packages/nisar/workflows/ionosphere.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1262,6 +1262,7 @@ def run(cfg: dict, runw_hdf5: str):
main_conn_image = None
side_conn_image = None
diff_ms_conn_image = None
diff_subband_conn_image = None

# ionosphere method using sub-bands uses subswath mask in freq A
# ionosphere method using side-band uses subswath masks in freq A and B
Expand Down Expand Up @@ -1416,7 +1417,7 @@ def run(cfg: dict, runw_hdf5: str):
np.s_[
row_start:row_start + block_rows_data,
:])
src_high_h5[rcom_path_freq_a].read_direct(
src_diff_h5[rcom_path_freq_a].read_direct(
diff_subband_conn_image,
np.s_[
row_start:row_start + block_rows_data,
Expand Down Expand Up @@ -1446,9 +1447,10 @@ def run(cfg: dict, runw_hdf5: str):
side_conn_image = np.empty(
[block_rows_data, cols_side],
dtype=float)
diff_ms_conn_image = np.empty(
[block_rows_data, cols_side],
dtype=float)
if iono_method == 'main_diff_ms_band':
diff_ms_conn_image = np.empty(
[block_rows_data, cols_side],
dtype=float)

if "subswath_mask" in mask_type:
subswath_mask_main_image = np.empty(
Expand Down Expand Up @@ -1525,7 +1527,7 @@ def run(cfg: dict, runw_hdf5: str):
:]
)
if "connected_components" in mask_type:
src_side_h5[rcom_path_freq_b].read_direct(
src_diff_h5[rcom_path_freq_b].read_direct(
diff_ms_conn_image,
np.s_[row_start:row_start + block_rows_data, :])

Expand All @@ -1537,13 +1539,23 @@ def run(cfg: dict, runw_hdf5: str):
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)
side_image = bridge_unwrapped_phase(
side_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)

if iono_method == 'main_side_band':
side_image = bridge_unwrapped_phase(
side_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)
elif iono_method == 'main_diff_ms_band':
diff_ms_image = bridge_unwrapped_phase(
diff_ms_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)

if unwrap_correction_bool:
main_image = unwrapping_correction_with_filter(
Expand All @@ -1554,14 +1566,25 @@ def run(cfg: dict, runw_hdf5: str):
sig_kernel_y=kernel_sigma_azimuth,
iterations=filter_iterations,
filter_method='convolution')
side_image = unwrapping_correction_with_filter(
side_image,
kernel_width=kernel_range_size,
kernel_length=kernel_azimuth_size,
sig_kernel_x=kernel_sigma_range,
sig_kernel_y=kernel_sigma_azimuth,
iterations=filter_iterations,
filter_method='convolution')

if iono_method == 'main_side_band':
side_image = unwrapping_correction_with_filter(
side_image,
kernel_width=kernel_range_size,
kernel_length=kernel_azimuth_size,
sig_kernel_x=kernel_sigma_range,
sig_kernel_y=kernel_sigma_azimuth,
iterations=filter_iterations,
filter_method='convolution')
elif iono_method == 'main_diff_ms_band':
diff_ms_image = unwrapping_correction_with_filter(
diff_ms_image,
kernel_width=kernel_range_size,
kernel_length=kernel_azimuth_size,
sig_kernel_x=kernel_sigma_range,
sig_kernel_y=kernel_sigma_azimuth,
iterations=filter_iterations,
filter_method='convolution')

# Estimate dispersive and non-dispersive phase
dispersive, non_dispersive = iono_phase_obj.compute_disp_nondisp(
Expand Down Expand Up @@ -1611,6 +1634,7 @@ def run(cfg: dict, runw_hdf5: str):
iono_std, nondisp_std = iono_phase_obj.estimate_iono_std(
main_coh=main_coh_image,
side_coh=side_coh_image,
diff_ms_coh=diff_ms_coh_image,
low_band_coh=sub_low_coh_image,
high_band_coh=sub_high_coh_image,
diff_low_high_coh=diff_coh_image,
Expand Down Expand Up @@ -1687,6 +1711,7 @@ def run(cfg: dict, runw_hdf5: str):
diff_ms_array=diff_ms_conn_image,
low_band_array=sub_low_conn_image,
high_band_array=sub_high_conn_image,
diff_low_high_band_array=diff_subband_conn_image,
slant_main=main_slant,
slant_side=side_slant)
mask_array = mask_array & mask_image
Expand All @@ -1701,14 +1726,14 @@ def run(cfg: dict, runw_hdf5: str):
mask_array = mask_array & mask_image

if "subswath_mask" in mask_type:
mask_array = iono_phase_obj.get_subswath_mask_array(
mask_subswath = iono_phase_obj.get_subswath_mask_array(
main_array=subswath_mask_main_image,
side_array=subswath_mask_side_image,
low_band_array=subswath_mask_image,
high_band_array=subswath_mask_image,
slant_main=main_slant,
slant_side=side_slant)
mask_array = mask_array & mask_image
mask_array &= mask_subswath

if "water" in mask_type:
# Extract preprocessing dictionary and open arrays
Expand Down Expand Up @@ -1844,7 +1869,7 @@ def run(cfg: dict, runw_hdf5: str):
zip(block_params, block_params_side)):

block_rows_data = block_parm.read_length
row_start = block * blocksize
row_start = block_ind * blocksize

if iono_method in iono_method_sideband:
block_parm_iono = block_parm_side
Expand Down Expand Up @@ -1977,25 +2002,43 @@ def run(cfg: dict, runw_hdf5: str):
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)
side_image = bridge_unwrapped_phase(
side_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)

if iono_method == 'main_side_band':
side_image = bridge_unwrapped_phase(
side_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)

elif iono_method == 'main_diff_ms_band':
diff_ms_image = bridge_unwrapped_phase(
diff_ms_image,
radius=bridge_radius,
min_num_pixel=bridge_minimum_samples,
erosion_size=bridge_erosion_size,
ramp_type=bridge_deramp_type,
deramp_max_num_sample=bridge_ramp_maximum_pixel)

if block_ind > 0:
main_image, _ = compute_phase_jump(
previous_main_with_pad,
main_image,
half_pad_length)
side_image, _ = compute_phase_jump(
previous_side_with_pad,
side_image,
half_pad_length)
if iono_method == 'main_side_band':
side_image, _ = compute_phase_jump(
previous_side_with_pad,
side_image,
half_pad_length)
elif iono_method == 'main_diff_ms_band':
diff_ms_image, _ = compute_phase_jump(
previous_diff_ms_with_pad,
diff_ms_image,
half_pad_length)
previous_main_with_pad = main_image
previous_side_with_pad = side_image
previous_diff_ms_with_pad = diff_ms_image

# Estimating phase unwrapping errors
com_unw_err, diff_unw_err = \
Expand Down
Loading
Loading