Inversion#

Classes#

Class

Description

NSBASMatrixFactory

Factory class to generate/format NSBAS matrix (d and G matrices)

NSBASInversion

a class used to operate NSBAS inversion.

PhaseDeformationConverter

A class to convert between phase and deformation (mm) for SAR interferometry.
Basis Functions#

Function

Description

censored_lstsq()

Solves least squares problem subject to missing data.

batch_lstsq()

calculates the least-squares solution for a batch of linear equations using the given G matrix and the data in d.

calculate_u()

Calculate correction matrix u by loop closure phase using least square.

Classes#

NSBASInversion#

class faninsar.NSBAS.inversion.NSBASInversion(matrix_factory: NSBASMatrixFactory, device: str | device | None = None, dtype: dtype = torch.float64, verbose=True)#

Bases: object

a class used to operate NSBAS inversion. The NSBAS inversion is usually expressed as: d = Gm, where d is the unwrapped interferograms matrix, G is the NSBAS matrix, and m is the model parameters, which is the combination of the deformation increment and the model parameters. see paper: TODO for more details.

Examples

__init__(matrix_factory: NSBASMatrixFactory, device: str | device | None = None, dtype: dtype = torch.float64, verbose=True)#

Initialize NSBASInversion

Parameters:

  • matrix_factory (NSBASMatrixFactory) – NSBASMatrixFactory instance

  • device (Optional[str | torch.device], optional) – device of torch.tensor used for computation. If None, use GPU if available, otherwise use CPU.

  • dtype (torch.dtype) – dtype of torch.tensor used for computation.

inverse() → Tuple[ndarray[tuple[int, ...], dtype[floating]], ndarray[tuple[int, ...], dtype[floating]], ndarray[tuple[int, ...], dtype[floating]], ndarray[tuple[int, ...], dtype[floating]]]#

Calculate increment displacement difference by NSBAS inversion.

Returns:

  • incs (np.ndarray (n_date - 1, n_pt)) – Incremental displacement

  • params (np.ndarray (n_param, n_pt)) – parameters of model in NSBAS inversion

  • residual_pair (np.ndarray (n_pair, n_pt)) – residual between interferograms and model result

  • residual_tsm (np.ndarray (n_date, n_pt)) – residual between time-series model and model result

NSBASMatrixFactory#

class faninsar.NSBAS.inversion.NSBASMatrixFactory(unw: ndarray[tuple[int, ...], dtype[floating]], pairs: Pairs | Sequence[str], model: TimeSeriesModels | None = None, gamma: float = 0.0001)#

Bases: object

Factory class to generate/format NSBAS matrix The NSBAS matrix is usually expressed as: d = Gm, where d is the unwrapped interferograms matrix, G is the NSBAS matrix, and m is the model parameters, which is the combination of the deformation increment and the model parameters. see paper: TODO for more details.

Note

After initialization, the d can still be updated by assigning a new unwrapped interferograms matrix to d. This is useful when the unwrapped interferograms is divided into multiple patches and the NSBAS matrix is calculated for each patch separately.

Examples

>>> import faninsar as fis
>>> import numpy as np

>>> names = ['20170111_20170204',
            '20170111_20170222',
            '20170111_20170318',
            '20170204_20170222',
            '20170204_20170318',
            '20170204_20170330',
            '20170222_20170318',
            '20170222_20170330',
            '20170222_20170411',
            '20170318_20170330']

>>> pairs = fis.Pairs.from_names(names)
>>> unw = np.random.randint(0, 255, (len(pairs),5))
>>> model = fis.AnnualSinusoidalModel(pairs.dates)
>>> nsbas_matrix = fis.NSBASMatrixFactory(unw, pairs, model)
>>> nsbas_matrix
NSBASMatrixFactory(
    pairs: Pairs(10)
    model: AnnualSinusoidalModel(dates: 6, unit: day)
    gamma: 0.0001
    G shape: (16, 9)
    d shape: (16, 5)
)

reset d by assigning a new unwrapped interferograms matrix with same pairs

>>> nsbas_matrix.d = np.random.randint(0, 255, (len(pairs), 10))
NSBASMatrixFactory(
    pairs: Pairs(10)
    model: AnnualSinusoidalModel(dates: 6, unit: day)
    gamma: 0.0001
    G shape: (16, 9)
    d shape: (16, 10)
)
slots = ['_pairs', '_model', '_gamma', '_G', '_d']#
__init__(unw: ndarray[tuple[int, ...], dtype[floating]], pairs: Pairs | Sequence[str], model: TimeSeriesModels | None = None, gamma: float = 0.0001)#

Initialize NSBASMatrixFactory

Parameters:

  • unw (NDArray (n_pairs, n_pixels)) – Unwrapped interferograms matrix

  • pairs (Pairs | Sequence[str]) – Pairs or Sequence of pair names

  • model (Optional[TimeSeriesModels], optional) – Time series model. If None, generate SBAS matrix rather than NSBAS matrix, by default None.

  • gamma (float, optional) – weight for the model component, by default 0.0001. This parameter will be ignored if model is None.

property pairs: Pairs#

Return pairs

property model: TimeSeriesModels | None#

Return model

property gamma: float#

Return gamma

property d: ndarray[tuple[int, ...], dtype[float32 | float64]]#

Return d matrix for NSBAS d = Gm

property G: ndarray[tuple[int, ...], dtype[float32]]#

Return G matrix for NSBAS d = Gm

PhaseDeformationConverter#

class faninsar.NSBAS.inversion.PhaseDeformationConverter(frequency: float = None, wavelength: float = None)#

Bases: object

A class to convert between phase and deformation (mm) for SAR interferometry.

__init__(frequency: float = None, wavelength: float = None) → None#

Initialize the converter. Either wavelength or frequency should be provided. If both are provided, wavelength will be recalculated by frequency.

Parameters:

  • frequency (float) – The frequency of the radar signal. Unit: GHz.

  • wavelength (float) – The wavelength of the radar signal. Unit: meter. this parameter will be ignored if frequency is provided.

phase2deformation(phase: ndarray[tuple[int, ...], dtype[floating]]) → ndarray[tuple[int, ...], dtype[floating]]#

Convert phase to deformation (mm)

deformation2phase(deformation: ndarray[tuple[int, ...], dtype[floating]]) → ndarray[tuple[int, ...], dtype[floating]]#

Convert deformation (mm) to phase (radian)

wrap_phase(phase: ndarray[tuple[int, ...], dtype[floating]]) → ndarray[tuple[int, ...], dtype[floating]]#

Wrap phase to [0, 2π]

Basis Functions#

batch_lstsq#

faninsar.NSBAS.inversion.batch_lstsq(G: ndarray | Tensor, d: ndarray | Tensor, dtype: dtype = torch.float64, device: str | device | None = None, verbose: bool = True, tqdm_args: dict = {}, return_numpy: bool = True) → ndarray[tuple[int, ...], dtype[floating]] | Tensor#

This function calculates the least-squares solution for a batch of linear equations using the given G matrix and the data in d.

Parameters:

  • G (np.ndarray | torch.Tensor) – model field matrix with shape of (n_im, n_param) or (n_pt, n_im, n_param). If G is 3D, the first dimension is the G matrix for every pixel.

  • d (np.ndarray | torch.Tensor) – data field matrix with shape of (n_im, n_pt).

  • dtype (torch.dtype, optional) – dtype of torch.tensor used for computation

  • device (Optional[str | torch.device], optional) – device of torch.tensor used for computation. If None, use GPU if available, otherwise use CPU.

  • verbose (bool, optional) – If True, show progress bar, by default True

  • tqdm_args (dict, optional) – Arguments to be passed to tqdm.tqdm Object for progress bar.

  • return_numpy (bool, optional) – If True, return a numpy array, otherwise return a torch tensor.

Returns:

X – (n_im x n_pt) matrix that minimizes norm(M*(GX - d)). If return_numpy is True, return a numpy array, otherwise return a torch tensor.

Return type:

np.ndarray | torch.Tensor

censored_lstsq#

faninsar.NSBAS.inversion.censored_lstsq(G: ndarray | Tensor, d: ndarray | Tensor, dtype: dtype = torch.float64, device: str | device | None = None, return_numpy: bool = True) → ndarray[tuple[int, ...], dtype[floating]] | Tensor#

Solves least squares problem subject to missing data. Reference: http://alexhwilliams.info/itsneuronalblog/2018/02/26/censored-lstsq/

Note

This function is used for solving the least squares problem with missing data. The missing data is represented by nan values in the data matrix d. If there are no nan values in d, you are recommended to use torch.linalg.lstsq() instead.

Parameters:

  • G (np.ndarray | torch.Tensor,) – model field matrix in shape of (n_im, n_param) or (n_pt, n_im, n_param). If G is 3D, the first dimension is the G matrix for each pixel.

  • d (np.ndarray | torch.Tensor, (n_im, n_pt) matrix) – data field matrix.

  • dtype (torch.dtype) – dtype of torch.tensor used for computation.

  • device (Optional[str | torch.device]) – device of torch.tensor used for computation. If None, use GPU if available, otherwise use CPU.

Returns:

X – (n_im x n_pt) matrix that minimizes norm(M*(GX - d)). If return_numpy is True, return a numpy array, otherwise return a torch tensor.

Return type:

np.ndarray | torch.Tensor

calculate_u#

faninsar.NSBAS.inversion.calculate_u(loops: Loops, unw_phases: ndarray, device: str | device | None = None, dtype: dtype = torch.float64) → ndarray[tuple[int, ...], dtype[floating]]#

Calculate correction matrix u by loop closure phase using least square. More details see paper:

Tip

The pairs in the loops may be fewer than the input pairs. To make sure the pairs in the loops are the same as the pairs in the unw_phases, you can use the Pairs.where() method to get the index/mask of the pairs in the loops from the input pairs.

Parameters:

  • loops (Loops) – Loops object. Loops are used to calculate the design matrix, which describes the relationship between the loops and pairs.

  • unw_phases (ndarray) – unwrapped interferograms phases with shape of (n_pair, n_pixel).

  • device (Optional[str | torch.device], optional) – device of torch.tensor used for computation. If None, use GPU if available, otherwise use CPU.

  • dtype (torch.dtype, optional) – dtype of torch.tensor used for computation.

Examples

get the loops from the pairs:

>>> loops = pairs.to_loops()
>>> idx = pairs.where(loops.pairs) # get the index of the pairs in the loops from the input pairs
>>> unw_used = unw[idx]

calculate u by loops and unwrapped interferometric phases:

>>> u = np.zeros_like(unw, dtype=np.float32)
>>> u[idx] = np.round(NSBAS.calculate_u(loops, unw_used))