ISTA.py 文件源码

def ISTA(
    fmatA,
    arrB,
    numLambda=0.1,
    numMaxSteps=100
):
    '''
    Wrapper around the ISTA algrithm to allow processing of arrays of signals
        fmatA         - input system matrix
        arrB          - input data vector (measurements)
        numLambda     - balancing parameter in optimization problem
                        between data fidelity and sparsity
        numMaxSteps   - maximum number of steps to run
        numL          - step size during the conjugate gradient step
    '''

    if len(arrB.shape) > 2:
        raise ValueError("Only n x m arrays are supported for ISTA")

    # calculate the largest singular value to get the right step size
    numL = 1.0 / (fmatA.largestSV ** 2)

    arrX = np.zeros(
        (fmatA.numM, arrB.shape[1]),
        dtype=np.promote_types(np.float32, arrB.dtype)
    )

    # start iterating
    for numStep in range(numMaxSteps):
        # do the gradient step and threshold

        arrStep = arrX -  numL * fmatA.backward(fmatA.forward(arrX) - arrB)
        arrX = _softThreshold(arrStep, numL * numLambda * 0.5)

    # return the unthresholded values for all non-zero support elements
    return np.where(arrX != 0, arrStep, arrX)


################################################################################
###  Maintenance and Documentation
################################################################################

################################################## inspection interface