Python numpy.convolve() Examples

def __box_filter_convolve(self, path, window_size):
        """
        An internal method that applies *normalized linear box filter* to path w.r.t averaging window
        
        Parameters:
        
        * path (numpy.ndarray): a cumulative sum of transformations
        * window_size (int): averaging window size
        """
        # pad path to size of averaging window
        path_padded = np.pad(path, (window_size, window_size), "median")
        # apply linear box filter to path
        path_smoothed = np.convolve(path_padded, self.__box_filter, mode="same")
        # crop the smoothed path to original path
        path_smoothed = path_smoothed[window_size:-window_size]
        # assert if cropping is completed
        assert path.shape == path_smoothed.shape
        # return smoothed path
        return path_smoothed 

def create_harmonic_mask(self, melody_signal):
        """
        Creates a harmonic mask from the melody signal. The mask is smoothed to reduce 
        the effects of discontinuities in the melody synthesizer.
        """
        stft = np.abs(melody_signal.stft())

        # Need to threshold the melody stft since the synthesized
        # F0 sequence overtones are at different weights.
        stft = stft ** self.compression
        stft /= np.maximum(np.max(stft, axis=1, keepdims=True), 1e-7)

        mask = np.empty(self.stft.shape)

        # Smoothing the mask row-wise using a low-pass filter to
        # get rid of discontuinities in the mask.
        kernel = np.full((1, self.smooth_length), 1 / self.smooth_length)
        for ch in range(self.audio_signal.num_channels):
            mask[..., ch] = convolve(stft[..., ch], kernel)
        return mask 

def testCausalConv(self):
        """Tests that the op is equivalent to a numpy implementation."""
        x1 = np.arange(1, 21, dtype=np.float32)
        x = np.append(x1, x1)
        x = np.reshape(x, [2, 20, 1])
        f = np.reshape(np.array([1, 1], dtype=np.float32), [2, 1, 1])
        out = causal_conv(x, f, 4)

        with self.test_session() as sess:
            result = sess.run(out)

        # Causal convolution using numpy
        ref = np.convolve(x1, [1, 0, 0, 0, 1], mode='valid')
        ref = np.append(ref, ref)
        ref = np.reshape(ref, [2, 16, 1])

        self.assertAllEqual(result, ref) 

def forecast3(self, step_ahead=1, start=None): #, end=None):
        '''another try for h-step ahead forecasting
        '''

        from .arima_process import arma2ma, ArmaProcess
        p,q = self.nar, self.nma
        k=0
        ar = self.params[k:k+p]
        ma = self.params[k+p:k+p+q]
        marep = arma2ma(ar,ma, start)[step_ahead+1:]  #truncated ma representation
        errors = self.error_estimate
        forecasts = np.convolve(errors, marep)
        return forecasts#[-(errors.shape[0] - start-5):] #get 5 overlapping for testing




    #copied from arima.ARIMA
    #TODO: is this needed as a method at all?
    #JP: not needed in this form, but can be replace with using the parameters 

def test_conv(mode, method):
    reload(convolution)
    # time vector for stimulus (long)
    stim_dur = 0.5  # seconds
    tsample = 0.001 / 1000
    t = np.arange(0, stim_dur, tsample)

    # stimulus (10 Hz anondic and cathodic pulse train)
    stim = np.zeros_like(t)
    stim[::1000] = 1
    stim[100::1000] = -1

    # kernel
    _, gg = gamma(1, 0.005, tsample)

    # make sure conv returns the same result as np.convolve for all modes:
    npconv = np.convolve(stim, gg, mode=mode)
    conv = convolution.conv(stim, gg, mode=mode, method=method)
    npt.assert_equal(conv.shape, npconv.shape)
    npt.assert_almost_equal(conv, npconv)

    with pytest.raises(ValueError):
        convolution.conv(gg, stim, mode="invalid")
    with pytest.raises(ValueError):
        convolution.conv(gg, stim, method="invalid") 

def _zseries_mul(z1, z2):
    """Multiply two z-series.

    Multiply two z-series to produce a z-series.

    Parameters
    ----------
    z1, z2 : 1-D ndarray
        The arrays must be 1-D but this is not checked.

    Returns
    -------
    product : 1-D ndarray
        The product z-series.

    Notes
    -----
    This is simply convolution. If symmetric/anti-symmetric z-series are
    denoted by S/A then the following rules apply:

    S*S, A*A -> S
    S*A, A*S -> A

    """
    return np.convolve(z1, z2) 

def convolve_beam(self, current, wake):
        """
        convolve wake with beam current

        :param current: current[:, 0] - s in [m], current[:, 1] - current in [A]
        :param wake: wake function in form: wake(s, b, t, period)
        :return:
        """
        s_shift = current[0, 0]
        current[:, 0] -= s_shift
        s = current[:, 0]

        step = (s[-1] - s[0]) / (len(s) - 1)
        q = current[:, 1] / speed_of_light

        w = np.array(
            [wake(si, b=self.b, t=self.t, period=self.period) for si in s]) * 377 * speed_of_light / (
                    4 * np.pi)
        wake = np.convolve(q, w) * step
        s_new = np.cumsum(np.ones(len(wake))) * step
        wake_kick = np.vstack((s_new, wake))
        return wake_kick.T 

def compute_fitness(genome, net, episodes, min_reward, max_reward):
    m = int(round(np.log(0.01) / np.log(genome.discount)))
    discount_function = [genome.discount ** (m - i) for i in range(m + 1)]

    reward_error = []
    for score, data in episodes:
        # Compute normalized discounted reward.
        dr = np.convolve(data[:,-1], discount_function)[m:]
        dr = 2 * (dr - min_reward) / (max_reward - min_reward) - 1.0
        dr = np.clip(dr, -1.0, 1.0)

        for row, dr in zip(data, dr):
            observation = row[:8]
            action = int(row[8])
            output = net.activate(observation)
            reward_error.append(float((output[action] - dr) ** 2))

    return reward_error 

def fit(self, episode_data, resampling_factor=1):
        """
        Estimates `advised` actions probabilities distribution based on data received.

        Args:
            episode_data:           1D np.array of unscaled price values in OHL[CV] format
            resampling_factor:      factor by which to resample given data
                                    by taking min/max values inside every resampled bar

        Returns:
             Np.array of size [resampled_data_size, actions_space_size] of probabilities of advised actions, where
             resampled_data_size = int(len(episode_data) / resampling_factor) + 1/0

        """
        # Vector of advised actions:
        data = self.resample_data(episode_data, resampling_factor)
        signals = self.estimate_actions(data)
        signals = self.adjust_signals(signals)

        # One-hot actions encoding:
        actions_one_hot = np.zeros([signals.shape[0], len(self.action_space)])
        actions_one_hot[np.arange(signals.shape[0]), signals] = 1

        # Want a bit relaxed discrete distribution over actions instead of one hot (heuristic):
        actions_distr = np.zeros(actions_one_hot.shape)

        # For all actions except 'hold' (due to heuristic skewness):
        actions_distr[:, 0] = actions_one_hot[:, 0]

        # ...spread out actions probabilities by convolving with gaussian kernel :
        for channel in range(1, actions_one_hot.shape[-1]):
            actions_distr[:, channel] = np.convolve(actions_one_hot[:, channel], self.kernel, mode='same') + 0.1

        # Normalize:
        actions_distr /= actions_distr.sum(axis=-1)[..., None]

        return actions_distr 

def test_object(self):
        d = [1.] * 100
        k = [1.] * 3
        assert_array_almost_equal(np.convolve(d, k)[2:-2], np.full(98, 3)) 

def moving_average(values, window):
    """
    Smooth values by doing a moving average

    :param values: (numpy array)
    :param window: (int)
    :return: (numpy array)
    """
    weights = np.repeat(1.0, window) / window
    return np.convolve(values, weights, 'valid') 

def ngram_min_hash(string, ngram_range=(2, 4), seed=0, return_minmax=False):
    """ Compute the min/max hash of the ngrams of the string.

    Parameters
    ----------
    string : str
        String to encode.
    ngram_range : tuple (min_n, max_n), default=(2, 4)
        The lower and upper boundary of the range of n-values 
        for different n-grams to be extracted.
    seed : int, default=0
        Integer used to seed the hashing function.
    return_minmax : bool, default=False
        If True, returns both the minhash and maxhash of the string.
        Else, only returns the minhash. 
    Returns
    -------
    int or tuple 
        The min_hash or (min_hash, max_hash) of the n-grams of the string.

    """
    # Create a numerical 1D array from the string
    array = np.frombuffer(string.encode(), dtype='int8', count=len(string))

    max_hash = MININT32
    min_hash = MAXINT32
    for atom_len in range(ngram_range[0], ngram_range[1]):
        atom = gen_atom(atom_len, seed=seed)
        # np.correlate is faster than np.convolve
        # the convolution gives a hash for each ngram
        hashes = np.correlate(array, atom)
        min_hash = min(min_hash, hashes.min())
        if return_minmax:
            max_hash = max(max_hash, hashes.max())

    # We should check that longer windows do not have different
    # statistics from shorter ones
    if return_minmax:
        return min_hash, max_hash
    return min_hash 

def _generate(self, n: int) -> np.ndarray:
        res = np.array([0.0])
        for _ in range(n):
            r = np.random.normal(-res.sum() * self.pull, self.scale)
            res = np.insert(res, -1, r)

        return np.convolve(res, np.ones((self.mean_n,)) / self.mean_n)[
            (self.mean_n - 1) : -1
        ] 

def polymul(c1, c2):
    """
    Multiply one polynomial by another.

    Returns the product of two polynomials `c1` * `c2`.  The arguments are
    sequences of coefficients, from lowest order term to highest, e.g.,
    [1,2,3] represents the polynomial ``1 + 2*x + 3*x**2.``

    Parameters
    ----------
    c1, c2 : array_like
        1-D arrays of coefficients representing a polynomial, relative to the
        "standard" basis, and ordered from lowest order term to highest.

    Returns
    -------
    out : ndarray
        Of the coefficients of their product.

    See Also
    --------
    polyadd, polysub, polydiv, polypow

    Examples
    --------
    >>> from numpy.polynomial import polynomial as P
    >>> c1 = (1,2,3)
    >>> c2 = (3,2,1)
    >>> P.polymul(c1,c2)
    array([  3.,   8.,  14.,   8.,   3.])

    """
    # c1, c2 are trimmed copies
    [c1, c2] = pu.as_series([c1, c2])
    ret = np.convolve(c1, c2)
    return pu.trimseq(ret) 

def running_mean(signal, n, axis=0):

    signal = np.array(signal)
    if signal.ndim == 1:
        signal_sum = np.convolve(signal, np.ones((2*n+1,)), mode='same')
        signal_num = np.convolve(signal*0+1, np.ones((2*n+1,)), mode='same')
        return signal_sum/signal_num

    elif signal.ndim == 2:
        smoothed = np.empty(signal.shape)
        if axis == 0:
            for i, sig in enumerate(signal):
                sig_sum = np.convolve(sig, np.ones((2*n+1,)), mode='same')
                sig_num = np.convolve(sig*0+1, np.ones((2*n+1,)), mode='same')
                smoothed[i, :] = sig_sum / sig_num
        elif axis == 1:
            for i, sig in enumerate(signal.T):
                sig_sum = np.convolve(sig, np.ones((2*n+1,)), mode='same')
                sig_num = np.convolve(sig*0+1, np.ones((2*n+1,)), mode='same')
                smoothed[:, i] = sig_sum / sig_num
        else:
            print('wrong axis')
        return smoothed

    else:
        print('wrong dimensions')
        return None 

def convolve(f, g):
    """
    FFT based convolution

    :param f: array
    :param g: array
    :return: array, (f * g)[n]
    """
    f_fft = fftpack.fftshift(fftpack.fftn(f))
    g_fft = fftpack.fftshift(fftpack.fftn(g))
    return fftpack.fftshift(fftpack.ifftn(fftpack.ifftshift(f_fft*g_fft))) 

def find_saturation(power, z, n_smooth=5):
    p = np.diff(np.log10(power))

    u = np.convolve(p, np.ones(n_smooth) / float(n_smooth), mode='same')
    um = np.max(u)
    
    ii = 0
    
    for i in range(len(u)):
        if u[i] < 0.0 * um and z[i] > 10: 
            ii = i
            break

    return z[ii+1], ii+1 

def polymul(c1, c2):
    """
    Multiply one polynomial by another.

    Returns the product of two polynomials `c1` * `c2`.  The arguments are
    sequences of coefficients, from lowest order term to highest, e.g.,
    [1,2,3] represents the polynomial ``1 + 2*x + 3*x**2.``

    Parameters
    ----------
    c1, c2 : array_like
        1-D arrays of coefficients representing a polynomial, relative to the
        "standard" basis, and ordered from lowest order term to highest.

    Returns
    -------
    out : ndarray
        Of the coefficients of their product.

    See Also
    --------
    polyadd, polysub, polydiv, polypow

    Examples
    --------
    >>> from numpy.polynomial import polynomial as P
    >>> c1 = (1,2,3)
    >>> c2 = (3,2,1)
    >>> P.polymul(c1,c2)
    array([  3.,   8.,  14.,   8.,   3.])

    """
    # c1, c2 are trimmed copies
    [c1, c2] = pu.as_series([c1, c2])
    ret = np.convolve(c1, c2)
    return pu.trimseq(ret) 

def polymul(c1, c2):
    """
    Multiply one polynomial by another.

    Returns the product of two polynomials `c1` * `c2`.  The arguments are
    sequences of coefficients, from lowest order term to highest, e.g.,
    [1,2,3] represents the polynomial ``1 + 2*x + 3*x**2.``

    Parameters
    ----------
    c1, c2 : array_like
        1-D arrays of coefficients representing a polynomial, relative to the
        "standard" basis, and ordered from lowest order term to highest.

    Returns
    -------
    out : ndarray
        Of the coefficients of their product.

    See Also
    --------
    polyadd, polysub, polymulx, polydiv, polypow

    Examples
    --------
    >>> from numpy.polynomial import polynomial as P
    >>> c1 = (1,2,3)
    >>> c2 = (3,2,1)
    >>> P.polymul(c1,c2)
    array([  3.,   8.,  14.,   8.,   3.])

    """
    # c1, c2 are trimmed copies
    [c1, c2] = pu.as_series([c1, c2])
    ret = np.convolve(c1, c2)
    return pu.trimseq(ret) 

def oral_cavity_filt(pole_amps, pole_freqs,fs):
    """
    This model for generating singing vowel sounds from sine tones comes
    from:

    https://simonl02.users.greyc.fr/files/Documents/Teaching/SignalProcessingLabs/lab3.pdf

    The above will be referred to throughout these docstrings as THE DOCUMENT.

    Original author: Fatemeh Pishdadian

    The arguments to the fucntions throughout this text follow the 
    signatures laid out in THE DOCUMENT.

    This is used in Melodia to generate the melody signal to produce a 
    mask.

    Solves "Q. Write a function to synthesize filter H(z)" in 
    THE DOCUMENT
    """
    num_pole_pair = len(pole_amps)
    poles = pole_amps * np.exp(1j * 2 * np.pi * pole_freqs / fs)
    poles_conj = np.conj(poles)

    denom_coeffs = 1
    for i in range(num_pole_pair):
        pole_temp = poles[i]
        pole_conj_temp = poles_conj[i]
        pole_pair_coeffs = np.convolve(np.array([1,-pole_temp]),np.array([1,-pole_conj_temp]))

        denom_coeffs = np.convolve(denom_coeffs, pole_pair_coeffs)
    return denom_coeffs 

def test_convolve_empty(self):
        # Convolve should raise an error for empty input array.
        assert_raises(ValueError, np.convolve, [], [1])
        assert_raises(ValueError, np.convolve, [1], []) 

def test_convolve_empty(self, level=rlevel):
        # Convolve should raise an error for empty input array.
        self.assertRaises(ValueError, np.convolve, [], [1])
        self.assertRaises(ValueError, np.convolve, [1], []) 

def test_no_overwrite(self):
        d = np.ones(100)
        k = np.ones(3)
        np.convolve(d, k)
        assert_array_equal(d, np.ones(100))
        assert_array_equal(k, np.ones(3)) 

def test_object(self):
        d = [1.] * 100
        k = [1.] * 3
        assert_array_almost_equal(np.convolve(d, k)[2:-2], np.full(98, 3)) 

def exp_MA(array, window):
    weights = numpy.exp(numpy.linspace(-1,0,window))
    weights /= sum(weights)
    return numpy.convolve(array,weights,'valid') 

def MA(array, window):
    weights = numpy.repeat(1, window)/window
    return numpy.convolve(array,weights,'valid') 

def polymul(c1, c2):
    """
    Multiply one polynomial by another.

    Returns the product of two polynomials `c1` * `c2`.  The arguments are
    sequences of coefficients, from lowest order term to highest, e.g.,
    [1,2,3] represents the polynomial ``1 + 2*x + 3*x**2.``

    Parameters
    ----------
    c1, c2 : array_like
        1-D arrays of coefficients representing a polynomial, relative to the
        "standard" basis, and ordered from lowest order term to highest.

    Returns
    -------
    out : ndarray
        Of the coefficients of their product.

    See Also
    --------
    polyadd, polysub, polydiv, polypow

    Examples
    --------
    >>> from numpy.polynomial import polynomial as P
    >>> c1 = (1,2,3)
    >>> c2 = (3,2,1)
    >>> P.polymul(c1,c2)
    array([  3.,   8.,  14.,   8.,   3.])

    """
    # c1, c2 are trimmed copies
    [c1, c2] = pu.as_series([c1, c2])
    ret = np.convolve(c1, c2)
    return pu.trimseq(ret) 

def smooth_curve(x):
    """用于使损失函数的图形变圆滑
    参考：http://glowingpython.blogspot.jp/2012/02/convolution-with-numpy.html
    """
    window_len = 11
    s = np.r_[x[window_len-1:0:-1], x, x[-1:-window_len:-1]]
    w = np.kaiser(window_len, 2)
    y = np.convolve(w/w.sum(), s, mode='valid')
    return y[5:len(y)-5] 

def test_convmtx(par):
    """Compare convmtx with np.convolve
    """
    x = np.random.normal(0, 1, par['nt']) + \
        par['imag'] * np.random.normal(0, 1, par['nt'])

    nh = 7
    h = np.hanning(7)
    H = convmtx(h, par['nt'])
    H = H[:, nh//2:-nh//2+1]

    y = np.convolve(x, h, mode='same')
    y1 = np.dot(H, x)
    assert_array_almost_equal(y, y1, decimal=4) 

def smooth_reward_curve(x, y):
    halfwidth = int(np.ceil(len(x) / 60))  # Halfwidth of our smoothing convolution
    k = halfwidth
    xsmoo = x
    ysmoo = np.convolve(y, np.ones(2 * k + 1), mode='same') / np.convolve(np.ones_like(y), np.ones(2 * k + 1),
        mode='same')
    return xsmoo, ysmoo