<commit_before>"""
Functions for producing a time series in PyRate.
InSAR time series inversion without smoothing
Based on original Matlab code by <NAME>. Matlab 'tsinvnosm' function.

..codeauthor:: <NAME>, <NAME>
"""

from numpy import where, isnan, nan, diff, zeros, empty
from numpy import float32, cumsum, dot, delete, asarray
from numpy.linalg import matrix_rank, pinv
from scipy.linalg import qr
import numpy as np
import matplotlib.pyplot as plt
import itertools
import parmap

import pyrate.config as config
from algorithm import master_slave_ids, get_epochs
from osgeo import gdal
import pyrate.ifgconstants as ifc
from pyrate import config as cf


def time_series(ifgs, pthresh, params, vcmt, mst=None):
    """
    Returns time series data from the given ifgs.

    :param ifgs: Sequence of interferograms.
    :param pthresh: minimum number of coherent observations for a pixel.
    :param params: configuration parameters
    :param vcmt: Derived positive definite temporal variance covariance matrix
    :param mst: [optional] array of ifg indexes from the pixel by pixel MST.
    :param parallel: use multiprocessing or not.

    :return: Tuple with the elements:

        - incremental deformation time series (NB: not velocity),
        - cumulative deformation time series, and
        - velocity

        these outputs are multi-dimensional arrays of size(nrows, ncols, nepochs-1),
        where:

        - *nrows* is the number of rows in the ifgs,
        - *ncols* is the  number of columns in the ifgs, and
        - *nepochs* is the number of unique epochs (dates) covered by the ifgs.).
    """

    if len(ifgs) < 1:
        msg = 'Time series requires 2+ interferograms'
        raise TimeSeriesError(msg)

    # MULTIPROCESSING parameters
    parallel = params[cf.PARALLEL]
    processes = params[cf.PROCESSES]

    SMORDER = params[cf.TIME_SERIES_SM_ORDER]
    SMFACTOR = np.power(10, params[cf.TIME_SERIES_SM_FACTOR])
    head = ifgs[0]
    check_time_series_params(head, pthresh)
    epochlist = get_epochs(ifgs)

    nrows = ifgs[0].nrows
    ncols = ifgs[0].ncols
    nifgs = len(ifgs)
    span = diff(epochlist.spans)
    nepoch = len(epochlist.dates)  # epoch number
    nvelpar = nepoch - 1  # velocity parameters number
    nlap = nvelpar - SMORDER  # Laplacian observations number

    mast_slave_ids = master_slave_ids(epochlist.dates)
    imaster = [mast_slave_ids[ifg.master] for ifg in ifgs]
    islave = [mast_slave_ids[ifg.slave] for ifg in ifgs]
    imaster = min(imaster, islave)
    islave = max(imaster, islave)
    B0 = zeros((nifgs, nvelpar))
    for i in range(nifgs):
        B0[i, imaster[i]:islave[i]] = span[imaster[i]:islave[i]]

    # change the sign if slave is earlier than master
    isign = where(imaster > islave)
    B0[isign[0], :] = -B0[isign[0], :]

    #  Laplacian smoothing coefficient matrix
    BLap0 = np.zeros(shape=(nlap, nvelpar))

    for i in range(nlap):
      if SMORDER == 1:
        BLap0[i, i:i+2] = [-1, 1]
      else:
        BLap0[i, i:i+3] = [1, -2, 1]

    # smf: Laplacian smoothing factor
    BLap0 *= SMFACTOR

    tsvel_matrix = np.empty(shape=(nrows, ncols, nvelpar),
                            dtype=np.float32)

    ifg_data = np.zeros((nifgs, nrows, ncols), dtype=float32)
    for ifg_num in xrange(nifgs):
        # ifgs[ifg_num].convert_to_nans(0)
        ifg_data[ifg_num] = ifgs[ifg_num].phase_data

    if mst is None:
        mst = ~isnan(ifg_data)

    if parallel == 1:
        tsvel_matrix = parmap.map(time_series_by_rows, range(nrows), B0, BLap0,
                                  SMORDER, ifg_data, mst, ncols, nvelpar,
                                  pthresh, vcmt, processes=processes)
    elif parallel == 2:
        res = parmap.starmap(time_series_by_pixel,
                             itertools.product(range(nrows), range(ncols)),
                             B0, BLap0, SMORDER, ifg_data, mst, nvelpar,
                             pthresh, vcmt, processes=processes)
        res = np.array(res)
        tsvel_matrix = np.reshape(res, newshape=(nrows, ncols, res.shape[1]))
    else:
        for row in xrange(nrows):
            for col in xrange(ncols):
                tsvel_matrix[row, col] = time_series_by_pixel(
                    row, col, B0, BLap0, SMORDER, ifg_data, mst, nvelpar,
                    pthresh, vcmt)

    # SB: do the span multiplication as a numpy linalg operation, MUCH faster
    # not even this is necessary here, perform late for performance
    # tsincr = tsvel_matrix * span<commit_msg>[timeseries] following discussion with MG, removed  commented line "# ifgs[ifg_num].convert_to_nans(0)" as it was deemed not applicable here<commit_after>"""
Functions for producing a time series in PyRate.
InSAR time series inversion without smoothing
Based on original Matlab code by H<NAME>. Matlab 'tsinvnosm' function.

..codeauthor:: <NAME>, <NAME>
"""

from numpy import where, isnan, nan, diff, zeros, empty
from numpy import float32, cumsum, dot, delete, asarray
from numpy.linalg import matrix_rank, pinv
from scipy.linalg import qr
import numpy as np
import matplotlib.pyplot as plt
import itertools
import parmap

import pyrate.config as config
from algorithm import master_slave_ids, get_epochs
from osgeo import gdal
import pyrate.ifgconstants as ifc
from pyrate import config as cf


def time_series(ifgs, pthresh, params, vcmt, mst=None):
    """
    Returns time series data from the given ifgs.

    :param ifgs: Sequence of interferograms.
    :param pthresh: minimum number of coherent observations for a pixel.
    :param params: configuration parameters
    :param vcmt: Derived positive definite temporal variance covariance matrix
    :param mst: [optional] array of ifg indexes from the pixel by pixel MST.
    :param parallel: use multiprocessing or not.

    :return: Tuple with the elements:

        - incremental deformation time series (NB: not velocity),
        - cumulative deformation time series, and
        - velocity

        these outputs are multi-dimensional arrays of size(nrows, ncols, nepochs-1),
        where:

        - *nrows* is the number of rows in the ifgs,
        - *ncols* is the  number of columns in the ifgs, and
        - *nepochs* is the number of unique epochs (dates) covered by the ifgs.).
    """

    if len(ifgs) < 1:
        msg = 'Time series requires 2+ interferograms'
        raise TimeSeriesError(msg)

    # MULTIPROCESSING parameters
    parallel = params[cf.PARALLEL]
    processes = params[cf.PROCESSES]

    SMORDER = params[cf.TIME_SERIES_SM_ORDER]
    SMFACTOR = np.power(10, params[cf.TIME_SERIES_SM_FACTOR])
    head = ifgs[0]
    check_time_series_params(head, pthresh)
    epochlist = get_epochs(ifgs)

    nrows = ifgs[0].nrows
    ncols = ifgs[0].ncols
    nifgs = len(ifgs)
    span = diff(epochlist.spans)
    nepoch = len(epochlist.dates)  # epoch number
    nvelpar = nepoch - 1  # velocity parameters number
    nlap = nvelpar - SMORDER  # Laplacian observations number

    mast_slave_ids = master_slave_ids(epochlist.dates)
    imaster = [mast_slave_ids[ifg.master] for ifg in ifgs]
    islave = [mast_slave_ids[ifg.slave] for ifg in ifgs]
    imaster = min(imaster, islave)
    islave = max(imaster, islave)
    B0 = zeros((nifgs, nvelpar))
    for i in range(nifgs):
        B0[i, imaster[i]:islave[i]] = span[imaster[i]:islave[i]]

    # change the sign if slave is earlier than master
    isign = where(imaster > islave)
    B0[isign[0], :] = -B0[isign[0], :]

    #  Laplacian smoothing coefficient matrix
    BLap0 = np.zeros(shape=(nlap, nvelpar))

    for i in range(nlap):
      if SMORDER == 1:
        BLap0[i, i:i+2] = [-1, 1]
      else:
        BLap0[i, i:i+3] = [1, -2, 1]

    # smf: Laplacian smoothing factor
    BLap0 *= SMFACTOR

    tsvel_matrix = np.empty(shape=(nrows, ncols, nvelpar),
                            dtype=np.float32)

    ifg_data = np.zeros((nifgs, nrows, ncols), dtype=float32)
    for ifg_num in xrange(nifgs):
        ifg_data[ifg_num] = ifgs[ifg_num].phase_data

    if mst is None:
        mst = ~isnan(ifg_data)

    if parallel == 1:
        tsvel_matrix = parmap.map(time_series_by_rows, range(nrows), B0, BLap0,
                                  SMORDER, ifg_data, mst, ncols, nvelpar,
                                  pthresh, vcmt, processes=processes)
    elif parallel == 2:
        res = parmap.starmap(time_series_by_pixel,
                             itertools.product(range(nrows), range(ncols)),
                             B0, BLap0, SMORDER, ifg_data, mst, nvelpar,
                             pthresh, vcmt, processes=processes)
        res = np.array(res)
        tsvel_matrix = np.reshape(res, newshape=(nrows, ncols, res.shape[1]))
    else:
        for row in xrange(nrows):
            for col in xrange(ncols):
                tsvel_matrix[row, col] = time_series_by_pixel(
                    row, col, B0, BLap0, SMORDER, ifg_data, mst, nvelpar,
                    pthresh, vcmt)

    # SB: do the span multiplication as a numpy linalg operation, MUCH faster
    # not even this is necessary here, perform late for performance
    # tsincr = tsvel_matrix * span
{'fraction_non_alphanumeric': 0.09215835969026066, 'fraction_numerical': 0.009379430690369724, 'mean_word_length': 2.721590909090909, 'pattern_counts': {'":': 0, '<': 11, '<?xml version=': 0, '>': 11, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '2601212', 'n_tokens_mistral': 3142, 'n_tokens_neox': 3028, 'n_words': 901}