python类convolve()的实例源码

kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 21 收藏 0 点赞 0 评论 0
def xcorr_offset(x1, x2):
    """
    Under MSR-LA License

    Based on MATLAB implementation from Spectrogram Inversion Toolbox

    References
    ----------
    D. Griffin and J. Lim. Signal estimation from modified
    short-time Fourier transform. IEEE Trans. Acoust. Speech
    Signal Process., 32(2):236-243, 1984.

    Malcolm Slaney, Daniel Naar and Richard F. Lyon. Auditory
    Model Inversion for Sound Separation. Proc. IEEE-ICASSP,
    Adelaide, 1994, II.77-80.

    Xinglei Zhu, G. Beauregard, L. Wyse. Real-Time Signal
    Estimation from Modified Short-Time Fourier Transform
    Magnitude Spectra. IEEE Transactions on Audio Speech and
    Language Processing, 08/2007.
    """
    x1 = x1 - x1.mean()
    x2 = x2 - x2.mean()
    frame_size = len(x2)
    half = frame_size // 2
    corrs = np.convolve(x1.astype('float32'), x2[::-1].astype('float32'))
    corrs[:half] = -1E30
    corrs[-half:] = -1E30
    offset = corrs.argmax() - len(x1)
    return offset
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 20 收藏 0 点赞 0 评论 0
def lsf_to_lpc(all_lsf):
    if len(all_lsf.shape) < 2:
        all_lsf = all_lsf[None]
    order = all_lsf.shape[1]
    all_lpc = np.zeros((len(all_lsf), order + 1))
    for i in range(len(all_lsf)):
        lsf = all_lsf[i]
        zeros = np.exp(1j * lsf)
        sum_zeros = zeros[::2]
        diff_zeros = zeros[1::2]
        sum_zeros = np.hstack((sum_zeros, np.conj(sum_zeros)))
        diff_zeros = np.hstack((diff_zeros, np.conj(diff_zeros)))
        sum_filt = np.poly(sum_zeros)
        diff_filt = np.poly(diff_zeros)

        if order % 2 != 0:
            deconv_diff = sg.convolve(diff_filt, [1, 0, -1])
            deconv_sum = sum_filt
        else:
            deconv_diff = sg.convolve(diff_filt, [1, -1])
            deconv_sum = sg.convolve(sum_filt, [1, 1])

        lpc = .5 * (deconv_sum + deconv_diff)
        # Last coefficient is 0 and not returned
        all_lpc[i] = lpc[:-1]
    return np.squeeze(all_lpc)
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 22 收藏 0 点赞 0 评论 0
def xcorr_offset(x1, x2):
    """
    Under MSR-LA License

    Based on MATLAB implementation from Spectrogram Inversion Toolbox

    References
    ----------
    D. Griffin and J. Lim. Signal estimation from modified
    short-time Fourier transform. IEEE Trans. Acoust. Speech
    Signal Process., 32(2):236-243, 1984.

    Malcolm Slaney, Daniel Naar and Richard F. Lyon. Auditory
    Model Inversion for Sound Separation. Proc. IEEE-ICASSP,
    Adelaide, 1994, II.77-80.

    Xinglei Zhu, G. Beauregard, L. Wyse. Real-Time Signal
    Estimation from Modified Short-Time Fourier Transform
    Magnitude Spectra. IEEE Transactions on Audio Speech and
    Language Processing, 08/2007.
    """
    x1 = x1 - x1.mean()
    x2 = x2 - x2.mean()
    frame_size = len(x2)
    half = frame_size // 2
    corrs = np.convolve(x1.astype('float32'), x2[::-1].astype('float32'))
    corrs[:half] = -1E30
    corrs[-half:] = -1E30
    offset = corrs.argmax() - len(x1)
    return offset
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 21 收藏 0 点赞 0 评论 0
def lsf_to_lpc(all_lsf):
    if len(all_lsf.shape) < 2:
        all_lsf = all_lsf[None]
    order = all_lsf.shape[1]
    all_lpc = np.zeros((len(all_lsf), order + 1))
    for i in range(len(all_lsf)):
        lsf = all_lsf[i]
        zeros = np.exp(1j * lsf)
        sum_zeros = zeros[::2]
        diff_zeros = zeros[1::2]
        sum_zeros = np.hstack((sum_zeros, np.conj(sum_zeros)))
        diff_zeros = np.hstack((diff_zeros, np.conj(diff_zeros)))
        sum_filt = np.poly(sum_zeros)
        diff_filt = np.poly(diff_zeros)

        if order % 2 != 0:
            deconv_diff = sg.convolve(diff_filt, [1, 0, -1])
            deconv_sum = sum_filt
        else:
            deconv_diff = sg.convolve(diff_filt, [1, -1])
            deconv_sum = sg.convolve(sum_filt, [1, 1])

        lpc = .5 * (deconv_sum + deconv_diff)
        # Last coefficient is 0 and not returned
        all_lpc[i] = lpc[:-1]
    return np.squeeze(all_lpc)
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 17 收藏 0 点赞 0 评论 0
def xcorr_offset(x1, x2):
    """
    Under MSR-LA License

    Based on MATLAB implementation from Spectrogram Inversion Toolbox

    References
    ----------
    D. Griffin and J. Lim. Signal estimation from modified
    short-time Fourier transform. IEEE Trans. Acoust. Speech
    Signal Process., 32(2):236-243, 1984.

    Malcolm Slaney, Daniel Naar and Richard F. Lyon. Auditory
    Model Inversion for Sound Separation. Proc. IEEE-ICASSP,
    Adelaide, 1994, II.77-80.

    Xinglei Zhu, G. Beauregard, L. Wyse. Real-Time Signal
    Estimation from Modified Short-Time Fourier Transform
    Magnitude Spectra. IEEE Transactions on Audio Speech and
    Language Processing, 08/2007.
    """
    x1 = x1 - x1.mean()
    x2 = x2 - x2.mean()
    frame_size = len(x2)
    half = frame_size // 2
    corrs = np.convolve(x1.astype('float32'), x2[::-1].astype('float32'))
    corrs[:half] = -1E30
    corrs[-half:] = -1E30
    offset = corrs.argmax() - len(x1)
    return offset
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 24 收藏 0 点赞 0 评论 0
def xcorr_offset(x1, x2):
    """
    Under MSR-LA License

    Based on MATLAB implementation from Spectrogram Inversion Toolbox

    References
    ----------
    D. Griffin and J. Lim. Signal estimation from modified
    short-time Fourier transform. IEEE Trans. Acoust. Speech
    Signal Process., 32(2):236-243, 1984.

    Malcolm Slaney, Daniel Naar and Richard F. Lyon. Auditory
    Model Inversion for Sound Separation. Proc. IEEE-ICASSP,
    Adelaide, 1994, II.77-80.

    Xinglei Zhu, G. Beauregard, L. Wyse. Real-Time Signal
    Estimation from Modified Short-Time Fourier Transform
    Magnitude Spectra. IEEE Transactions on Audio Speech and
    Language Processing, 08/2007.
    """
    x1 = x1 - x1.mean()
    x2 = x2 - x2.mean()
    frame_size = len(x2)
    half = frame_size // 2
    corrs = np.convolve(x1.astype('float32'), x2[::-1].astype('float32'))
    corrs[:half] = -1E30
    corrs[-half:] = -1E30
    offset = corrs.argmax() - len(x1)
    return offset
kdllib.py 文件源码 项目:crikey 作者: kastnerkyle 项目源码 文件源码 阅读 27 收藏 0 点赞 0 评论 0
def lsf_to_lpc(all_lsf):
    if len(all_lsf.shape) < 2:
        all_lsf = all_lsf[None]
    order = all_lsf.shape[1]
    all_lpc = np.zeros((len(all_lsf), order + 1))
    for i in range(len(all_lsf)):
        lsf = all_lsf[i]
        zeros = np.exp(1j * lsf)
        sum_zeros = zeros[::2]
        diff_zeros = zeros[1::2]
        sum_zeros = np.hstack((sum_zeros, np.conj(sum_zeros)))
        diff_zeros = np.hstack((diff_zeros, np.conj(diff_zeros)))
        sum_filt = np.poly(sum_zeros)
        diff_filt = np.poly(diff_zeros)

        if order % 2 != 0:
            deconv_diff = sg.convolve(diff_filt, [1, 0, -1])
            deconv_sum = sum_filt
        else:
            deconv_diff = sg.convolve(diff_filt, [1, -1])
            deconv_sum = sg.convolve(sum_filt, [1, 1])

        lpc = .5 * (deconv_sum + deconv_diff)
        # Last coefficient is 0 and not returned
        all_lpc[i] = lpc[:-1]
    return np.squeeze(all_lpc)
update_z.py 文件源码 项目:alphacsc 作者: alphacsc 项目源码 文件源码 阅读 20 收藏 0 点赞 0 评论 0
def _fprime(ds, zi, Xi=None, sample_weights=None, reg=None, return_func=False):
    """np.dot(D.T, X[i] - np.dot(D, zi)) + reg

    Parameters
    ----------
    ds : array, shape (n_atoms, n_times_atom)
        The atoms
    zi : array, shape (n_atoms * n_times_valid)
        The activations
    Xi : array, shape (n_times, ) or None
        The data array for one trial
    sample_weights : array, shape (n_times, ) or None
        The sample weights for one trial
    reg : float or None
        The regularization constant
    return_func : boolean
        Returns also the objective function, used to speed up LBFGS solver

    Returns
    -------
    (func) : float
        The objective function
    grad : array, shape (n_atoms * n_times_valid)
        The gradient
    """
    n_atoms, n_times_atom = ds.shape
    zi_reshaped = zi.reshape((n_atoms, -1))
    Dzi = _choose_convolve(zi_reshaped, ds)
    if Xi is not None:
        Dzi -= Xi

    if sample_weights is not None:
        if return_func:
            # preserve Dzi, we don't want to apply the weights twice in func
            wDzi = sample_weights * Dzi
        else:
            Dzi *= sample_weights
            wDzi = Dzi
    else:
        wDzi = Dzi

    if return_func:
        func = 0.5 * np.dot(wDzi, Dzi.T)
        if reg is not None:
            func += reg * zi.sum()

    # Now do the dot product with the transpose of D (D.T) which is
    # the conv by the reversed filter (keeping valid mode)
    grad = np.concatenate(
        [signal.convolve(wDzi, d[::-1], 'valid') for d in ds])
    # grad = -np.dot(D.T, X[i] - np.dot(D, zi))
    if reg is not None:
        grad += reg

    if return_func:
        return func, grad
    else:
        return grad
getWeights.py 文件源码 项目:python-psignifit 作者: wichmann-lab 项目源码 文件源码 阅读 18 收藏 0 点赞 0 评论 0
def getWeights(X1D):

    d = len(X1D)

    ''' puts the X values in their respective dimensions to use bsxfun for
    evaluation'''
    Xreshape = []
    Xreshape.append(X1D[0].reshape(-1,1))
    if d >= 2:
        Xreshape.append(X1D[1].reshape([1,-1]))

    for idx in range (2,d):
        dims = [1]*idx
        dims.append(-1)
        Xreshape.append(reshape(X1D[idx], dims))


    # to find and handle singleton dimensions
    sizes = [X1D[ix].size for ix in range(0,d)]
    Xlength = array(sizes)

    ''' calculate weights '''
    #1) calculate mass/volume for each cuboid
    weight = 1
    for idx in range(0,d):
        if Xlength[idx] > 1:
            if idx > 1:
                weight = multiply(weight[...,newaxis],diff(Xreshape[idx], axis=idx))
            else:
                weight = multiply(weight,diff(Xreshape[idx], axis=idx))
        else:
            weight = weight[...,newaxis]
    #2) sum to get weight for each point
    if d > 1:
        dims = tile(2,[1,d])
        dims[0,where(Xlength == 1)] = 1
        d = sum(Xlength > 1)
        weight = (2.**(-d))*convn(weight, ones(dims.ravel()), mode='full')
    else:
        weight = (2.**(-1))*convolve(weight, array([[1],[1]]))

    return weight
ecg.py 文件源码 项目:CerebralCortex-2.0-legacy 作者: MD2Korg 项目源码 文件源码 阅读 20 收藏 0 点赞 0 评论 0
def compute_moving_window_int(sample: np.ndarray,
                              fs: float,
                              blackman_win_length: int,
                              filter_length: int = 257,
                              delta: float = .02) -> np.ndarray:
    """
    :param sample: ecg sample array
    :param fs: sampling frequency
    :param blackman_win_length: length of the blackman window on which to compute the moving window integration
    :param filter_length: length of the FIR bandpass filter on which filtering is done on ecg sample array
    :param delta: to compute the weights of each band in FIR filter
    :return: the Moving window integration of the sample array
    """
    # I believe these constants can be kept in a file

    # filter edges
    filter_edges = [0, 4.5 * 2 / fs, 5 * 2 / fs, 20 * 2 / fs, 20.5 * 2 / fs, 1]
    # gains at filter band edges
    gains = [0, 0, 1, 1, 0, 0]
    # weights
    weights = [500 / delta, 1 / delta, 500 / delta]
    # length of the FIR filter

    # FIR filter coefficients for bandpass filtering
    filter_coeff = signal.firls(filter_length, filter_edges, gains, weights)

    # bandpass filtered signal
    bandpass_signal = signal.convolve(sample, filter_coeff, 'same')
    bandpass_signal /= np.percentile(bandpass_signal, 90)

    # derivative array
    derivative_array = (np.array([-1.0, -2.0, 0, 2.0, 1.0])) * (1 / 8)
    # derivative signal (differentiation of the bandpass)
    derivative_signal = signal.convolve(bandpass_signal, derivative_array, 'same')
    derivative_signal /= np.percentile(derivative_signal, 90)

    # squared derivative signal
    derivative_squared_signal = derivative_signal ** 2
    derivative_squared_signal /= np.percentile(derivative_squared_signal, 90)

    # blackman window
    blackman_window = np.blackman(blackman_win_length)
    # moving window Integration of squared derivative signal
    mov_win_int_signal = signal.convolve(derivative_squared_signal, blackman_window, 'same')
    mov_win_int_signal /= np.percentile(mov_win_int_signal, 90)

    return mov_win_int_signal
fm.py 文件源码 项目:rtlsdr-rds-demod 作者: itdaniher 项目源码 文件源码 阅读 23 收藏 0 点赞 0 评论 0
def symbol_recovery_24(xdi, xdq):
    angles = numpy.where(xdi >= 0, numpy.arctan2(xdq, xdi), numpy.arctan2(-xdq, -xdi))
    theta = (signal.convolve(angles, smooth)) [-len(xdi):]

    xr = (xdi + 1j * xdq) * numpy.exp(-1j * theta)
    bi = (numpy.real(xr) >= 0) + 0
    # pll parameters
    period = 24
    halfPeriod = period / 2
    corr = period / 24.
    phase = 0

    res = []
    pin = 0

    stats = {0: 0, 1: 1}
    oddity = 0

    latestXrSquared = [0]*8
    lxsIndex = 0
    theta = [0]
    shift = 0

    # pll, system model, error calculation, estimate update
    for i in range(1, len(bi)):
        if bi[i-1] != bi[i]:
            if phase < halfPeriod-2:
                phase += corr
            elif phase > halfPeriod+2:
                phase -= corr
        if phase >= period:
            phase -= period
            latestXrSquared[lxsIndex] = (xdi[i] + 1j * xdq[i])**2
            lxsIndex += 1
            if lxsIndex >= len(latestXrSquared):
                lxsIndex = 0
            th = shift + cmath.phase(sum(latestXrSquared)) / 2
            if abs(th - theta[-1]) > 2:
                if th < theta[-1]:
                    shift += math.pi
                    th += math.pi
                else:
                    shift -= math.pi
                    th -= math.pi
            theta.append(th)
            oddity += 1
            if oddity == 2:
                oddity = 0
                yp = (xdi[i] + 1j * xdq[i])
                ypp = cmath.exp(-1j * th) * yp
                # bit decode
                nin = 1 * (ypp.real > 0)
                stats[nin] += 1
                res.append(pin ^ nin)
                pin = nin
        phase += 1
    return res
fm.py 文件源码 项目:rtlsdr-rds-demod 作者: itdaniher 项目源码 文件源码 阅读 25 收藏 0 点赞 0 评论 0
def demodulate_array(h, soft_lcd):
    # primary worker function, credit to all
    i = h[1:] * np.conj(h[:-1])
    j = np.angle(i)

    k = signal.convolve(j, basebandBP)

    # resample from 256kHz to 228kHz
    rdsBand = signal.resample(k, int(len(k)*228e3/256e3))

    # length modulo 4
    rdsBand = rdsBand[:(len(rdsBand)//4)*4]
    c57 = numpy.tile( [1., -1.], len(rdsBand)//4 )
    xi = rdsBand[::2] * c57
    xq = rdsBand[1::2] * (-c57)
    xfi = signal.convolve(xi, filtLP)
    xfq = signal.convolve(xq, filtLP)
    xsfi = signal.convolve(xfi, pulseFilt)
    xsfq = signal.convolve(xfq, pulseFilt)

    if len(xsfi) % 2 == 1:
        xsfi = xsfi[:-1]
        xsfq = xsfq[:-1]

    xdi = (xsfi[::2] + xsfi[1::2]) / 2
    xdq = xsfq[::2]

    res = symbol_recovery_24(xdi, xdq)
    hits = []
    for i in range(len(res)-26):
        h = rds_crc(res, i, 26)
        if h:
            hits.append( (i, h) )
    print(res,hits)
    packets = []
    print([decode_one(res, x[0]) for x in hits if x[1] == 'A'])
    for i in range(len(hits)-3):
        if hits[i][1] == "A":
            bogus = False
            for j,sp in enumerate("ABCD"):
                if 26*j != hits[i+j][0] - hits[i][0]:
                    bogus = True
                if hits[i+j][1] != sp:
                    bogus = True
            if not bogus:
                for j in range(4):
                    packets.append(decode_one(res, hits[i+j][0]))
    soft_lcd.update_state(packets)


问题


面经


文章

微信
公众号

扫码关注公众号