Python numpy.correlate() Examples

def _presample_varcov(self, params):
        """
        Returns the inverse of the presample variance-covariance.

        Notes
        -----
        See Hamilton p. 125
        """
        k = self.k_trend
        p = self.k_ar
        p1 = p+1

        # get inv(Vp) Hamilton 5.3.7
        params0 = np.r_[-1, params[k:]]

        Vpinv = np.zeros((p, p), dtype=params.dtype)
        for i in range(1, p1):
            Vpinv[i-1, i-1:] = np.correlate(params0, params0[:i],)[:-1]
            Vpinv[i-1, i-1:] -= np.correlate(params0[-i:], params0,)[:-1]

        Vpinv = Vpinv + Vpinv.T - np.diag(Vpinv.diagonal())
        return Vpinv 

def extractMseq(cover, stego, secret_length, m, tau=1):
	u"""Extract secret informations by spread spectrum using m-sequence.
	@param  cover         : cover data (2 dimensional np.ndarray)
	@param  stego         : stego data (2 dimension np.ndarray)
	@param  secret_length : length of secret information
	@param  m             : M-Sequence
	@param  tau           : embed shift interval
	@return secret        : extracted secret information
	"""

	cover = _image2vrctor(cover)
	stego = _image2vrctor(stego)

	m_length = len(m)

	data = stego - cover
	data = data[:m_length:tau]

	secret_data = correlate(m, data, cycle=CYCLE)
	center = ((m_length-1)*2+1)//2
	secret_data = secret_data[center:center+secret_length]
	secret_data = list(map(_checkData, secret_data))

	return secret_data 

def ar_generator(N=512, sigma=1.):
    # this generates a signal u(n) = a1*u(n-1) + a2*u(n-2) + ... + v(n)
    # where v(n) is a stationary stochastic process with zero mean
    # and variance = sigma
    # this sequence is shown to be estimated well by an order 8 AR system
    taps = np.array([2.7607, -3.8106, 2.6535, -0.9238])
    v = np.random.normal(size=N, scale=sigma**0.5)
    u = np.zeros(N)
    P = len(taps)
    for l in range(P):
        u[l] = v[l] + np.dot(u[:l][::-1], taps[:l])
    for l in range(P,N):
        u[l] = v[l] + np.dot(u[l-P:l][::-1], taps)
    return u, v, taps

#JP: small differences to using np.correlate, because assumes mean(s)=0
#    denominator is N, not N-k, biased estimator
#    misnomer: (biased) autocovariance not autocorrelation
#from nitime.utils 

def _autocorr(data):
    """
    Calculates the auto correlation of a given array.

    Parameters
    ----------
    data : array-like
        The array to calculate the autocorrelation of

    Returns
    -------
    u : ndarray
        The array of autocorrelations
    """
    u = np.correlate(data, data, mode='full')
    # Take upper half of correlation matrix
    return u[u.size // 2:] 

def test_rank0(self, dt):
        a = np.array(np.random.randn()).astype(dt)
        a += 1j * np.array(np.random.randn()).astype(dt)
        b = np.array(np.random.randn()).astype(dt)
        b += 1j * np.array(np.random.randn()).astype(dt)

        y_r = (correlate(a.real, b.real)
               + correlate(a.imag, b.imag)).astype(dt)
        y_r += 1j * (-correlate(a.real, b.imag) + correlate(a.imag, b.real))

        y = correlate(a, b, 'full')
        assert_array_almost_equal(y, y_r, decimal=self.decimal(dt) - 1)
        assert_equal(y.dtype, dt)

        assert_equal(correlate([1], [2j]), correlate(1, 2j))
        assert_equal(correlate([2j], [3j]), correlate(2j, 3j))
        assert_equal(correlate([3j], [4]), correlate(3j, 4)) 

def test_consistency_correlate_funcs(self):
        # Compare np.correlate, signal.correlate, signal.correlate2d
        a = np.arange(5)
        b = np.array([3.2, 1.4, 3])
        for mode in ['full', 'valid', 'same']:
            assert_almost_equal(np.correlate(a, b, mode=mode),
                                signal.correlate(a, b, mode=mode))
            assert_almost_equal(np.squeeze(signal.correlate2d([a], [b],
                                                              mode=mode)),
                                signal.correlate(a, b, mode=mode))

            # See gh-5897
            if mode == 'valid':
                assert_almost_equal(np.correlate(b, a, mode=mode),
                                    signal.correlate(b, a, mode=mode))
                assert_almost_equal(np.squeeze(signal.correlate2d([b], [a],
                                                                  mode=mode)),
                                    signal.correlate(b, a, mode=mode)) 

def plot(seq):
    """Plot the autocorrelation of the given sequence."""
    import matplotlib.pyplot as plt

    bipolar = np.where(seq, 1.0, -1.0)
    autocorr = np.correlate(bipolar, bipolar, 'same')

    plt.figure()
    plt.title("Length {} Gold code autocorrelation".format(len(seq)))
    xdata = np.arange(len(seq)) - len(seq) // 2
    plt.plot(xdata, autocorr, '.-')
    plt.show() 

def test_rank1_same(self, dt):
        a, b, y_r = self._setup_rank1(dt)
        y = correlate(a, b, 'same')
        assert_array_almost_equal(y, y_r[:-1])
        assert_equal(y.dtype, dt) 

def test_rank1_full(self, dt):
        a, b, y_r = self._setup_rank1(dt)
        y = correlate(a, b, 'full')
        assert_array_almost_equal(y, y_r)
        assert_equal(y.dtype, dt) 

def test_rank3_valid(self, dt):
        a, b, y_r = self._setup_rank3(dt)
        y = correlate(a, b, "valid")
        assert_array_almost_equal(y, y_r[1:2, 2:4, 3:5])
        assert_equal(y.dtype, dt)

        # See gh-5897
        y = correlate(b, a, "valid")
        assert_array_almost_equal(y, y_r[1:2, 2:4, 3:5][::-1, ::-1, ::-1])
        assert_equal(y.dtype, dt) 

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

def test_complex(self):
        x = np.array([1, 2, 3, 4+1j], dtype=complex)
        y = np.array([-1, -2j, 3+1j], dtype=complex)
        r_z = np.array([3-1j, 6, 8+1j, 11+5j, -5+8j, -4-1j], dtype=complex)
        r_z = r_z[::-1].conjugate()
        z = np.correlate(y, x, mode='full')
        assert_array_almost_equal(z, r_z) 

def _print_stats(seq):
    bipolar = np.where(seq, 1.0, -1.0)
    autocorr = np.correlate(bipolar, bipolar, 'same')

    peaks = np.sort(np.abs(autocorr))
    peak = peaks[-1]
    noise = np.sqrt(np.mean(peaks[:-1]**2))

    peak_to_peak2 = peak / peaks[-2]
    peak_to_noise = peak / noise

    print("Peak amplitude: {:.0f}".format(peak))
    print("Largest non-peak amplitude: {:.0f}".format(peaks[-2]))
    print("Peak-to-max: {:.2f}".format(peak_to_peak2))
    print("Peak-to-noise: {:.2f}".format(peak_to_noise)) 

def test_float(self):
        self._setup(float)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.x, self.y[:-1], 'full')
        assert_array_almost_equal(z, self.z1_4)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2)
        z = np.correlate(self.x[::-1], self.y, 'full')
        assert_array_almost_equal(z, self.z1r)
        z = np.correlate(self.y, self.x[::-1], 'full')
        assert_array_almost_equal(z, self.z2r)
        z = np.correlate(self.xs, self.y, 'full')
        assert_array_almost_equal(z, self.zs) 

def test_object(self):
        self._setup(Decimal)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2) 

def correlate(a, v, mode='valid'):
    """Returns the cross-correlation of two 1-dimensional sequences.

    Args:
        a (cupy.ndarray): first 1-dimensional input.
        v (cupy.ndarray): second 1-dimensional input.
        mode (str, optional): `valid`, `same`, `full`

    Returns:
        cupy.ndarray: Discrete cross-correlation of a and v.

    .. seealso:: :func:`numpy.correlate`

    """
    if a.size == 0 or v.size == 0:
        raise ValueError('Array arguments cannot be empty')
    if a.ndim != 1 or v.ndim != 1:
        raise ValueError('object too deep for desired array')
    # choose_conv_method does not choose from the values in
    # the input array, so no need to apply conj.
    method = cupyx.scipy.signal.choose_conv_method(a, v, mode)
    if method == 'direct':
        out = cupy.math.misc._dot_convolve(a, v.conj()[::-1], mode)
    elif method == 'fft':
        out = cupy.math.misc._fft_convolve(a, v.conj()[::-1], mode)
    else:
        raise ValueError('Unsupported method')
    return out 

def test_float(self):
        self._setup(float)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.x, self.y[:-1], 'full')
        assert_array_almost_equal(z, self.z1_4)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2)
        z = np.correlate(self.x[::-1], self.y, 'full')
        assert_array_almost_equal(z, self.z1r)
        z = np.correlate(self.y, self.x[::-1], 'full')
        assert_array_almost_equal(z, self.z2r)
        z = np.correlate(self.xs, self.y, 'full')
        assert_array_almost_equal(z, self.zs) 

def test_rank1_valid(self, dt):
        a, b, y_r = self._setup_rank1(dt)
        y = correlate(a, b, 'valid')
        assert_array_almost_equal(y, y_r[1:4])
        assert_equal(y.dtype, dt)

        # See gh-5897
        y = correlate(b, a, 'valid')
        assert_array_almost_equal(y, y_r[1:4][::-1])
        assert_equal(y.dtype, dt) 

def test_method(self, dt):
        if dt == Decimal:
            method = choose_conv_method([Decimal(4)], [Decimal(3)])
            assert_equal(method, 'direct')
        else:
            a, b, y_r = self._setup_rank3(dt)
            y_fft = correlate(a, b, method='fft')
            y_direct = correlate(a, b, method='direct')

            assert_array_almost_equal(y_r, y_fft, decimal=self.equal_tolerance(y_fft.dtype))
            assert_array_almost_equal(y_r, y_direct, decimal=self.equal_tolerance(y_fft.dtype))
            assert_equal(y_fft.dtype, dt)
            assert_equal(y_direct.dtype, dt) 

def test_invalid_params(self):
        a = [3, 4, 5]
        b = [1, 2, 3]
        assert_raises(ValueError, correlate, a, b, mode='spam')
        assert_raises(ValueError, correlate, a, b, mode='eggs', method='fft')
        assert_raises(ValueError, correlate, a, b, mode='ham', method='direct')
        assert_raises(ValueError, correlate, a, b, mode='full', method='bacon')
        assert_raises(ValueError, correlate, a, b, mode='same', method='bacon') 

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

def test_object(self):
        self._setup(Decimal)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2) 

def test_float(self):
        self._setup(float)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.x, self.y[:-1], 'full')
        assert_array_almost_equal(z, self.z1_4)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2)
        z = np.correlate(self.x[::-1], self.y, 'full')
        assert_array_almost_equal(z, self.z1r)
        z = np.correlate(self.y, self.x[::-1], 'full')
        assert_array_almost_equal(z, self.z2r)
        z = np.correlate(self.xs, self.y, 'full')
        assert_array_almost_equal(z, self.zs) 

def rolling_mean(X, window_size):
    w = 1.0 / window_size * np.ones((window_size))
    return np.correlate(X, w, 'valid') 

def autocorrelation(x):
    n = x.shape[0]
    x = x - x.mean()
    r = np.correlate(x, x, mode='full')[-n:]
    return r/(x.var()*np.arange(n, 0, -1)) 

def correlate(data1, data2, cycle=NON_CYCLE):
	u"""Calculate correlate function.
	@param  data1     : data1
	@param  data2     : data2
	@param  cycle     : CYCLE or NON_CYCLE (Default is NON_CYCLE)
	@return correlate : correlate list
	"""
	if cycle == NON_CYCLE:
		correlate = np.correlate(data1, data2, 'full')

	elif cycle == CYCLE:
		N = len(data1)
		correlate = list()

		for i in np.arange(-N+1, N):
			ans = 0
			data1_shifted = np.roll(data1,i)
			for j in range(0,N):
				ans += data1_shifted[j]*data2[j]
			correlate.append(ans)

		correlate = np.array(correlate)

	else:
		print('Invalid argument.')
		print('Please review parameter in correlate function.')
		sys.exit()

	return correlate 

def test_complex(self):
        x = np.array([1, 2, 3, 4+1j], dtype=complex)
        y = np.array([-1, -2j, 3+1j], dtype=complex)
        r_z = np.array([3-1j, 6, 8+1j, 11+5j, -5+8j, -4-1j], dtype=complex)
        r_z = r_z[::-1].conjugate()
        z = np.correlate(y, x, mode='full')
        assert_array_almost_equal(z, r_z) 

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

def test_object(self):
        self._setup(Decimal)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2) 

def test_float(self):
        self._setup(float)
        z = np.correlate(self.x, self.y, 'full')
        assert_array_almost_equal(z, self.z1)
        z = np.correlate(self.x, self.y[:-1], 'full')
        assert_array_almost_equal(z, self.z1_4)
        z = np.correlate(self.y, self.x, 'full')
        assert_array_almost_equal(z, self.z2)
        z = np.correlate(self.x[::-1], self.y, 'full')
        assert_array_almost_equal(z, self.z1r)
        z = np.correlate(self.y, self.x[::-1], 'full')
        assert_array_almost_equal(z, self.z2r)
        z = np.correlate(self.xs, self.y, 'full')
        assert_array_almost_equal(z, self.zs)