我们从Python开源项目中,提取了以下42个代码示例,用于说明如何使用scipy.signal.convolve()。
def _signal_template(self, signal_period, frequency_sampling, frequency_cut, mean, std, filter_std=25, filter_window_size=500, seed=None, n: int=1000): randgen = np.random.RandomState(seed=seed) noisy = randgen.normal(mean, std, size=self._size) high_band_filter = firwin(filter_window_size, frequency_cut, nyq=frequency_sampling, window=('gaussian', filter_std)) filtered_noisy = convolve(noisy, high_band_filter, mode='same') return filtered_noisy
def boxcar(y, window_size=3): """ Smooth the input vector using the mean of the neighboring values, where neighborhood size is defined by the window. Parameters ========== y : array The vector to be smoothed. window_size : int An odd integer describing the window size. Returns ======= : array The smoothed array. """ filt = np.ones(window_size) / window_size return Series(np.convolve(y, filt, mode='same'), index=y.index)
def gaussian(y, window_size=3, sigma=2): """ Apply a gaussian filter to smooth the input vector Parameters ========== y : array The input array window_size : int An odd integer describing the size of the filter. sigma : float The numver of standard deviation """ filt = signal.gaussian(window_size, sigma) return Series(signal.convolve(y, filt, mode='same'), index=y.index)
def test_convolve_generalization(): ag_convolve = autograd.scipy.signal.convolve A_35 = R(3, 5) A_34 = R(3, 4) A_342 = R(3, 4, 2) A_2543 = R(2, 5, 4, 3) A_24232 = R(2, 4, 2, 3, 2) for mode in ['valid', 'full']: assert npo.allclose(ag_convolve(A_35, A_34, axes=([1], [0]), mode=mode)[1, 2], sp_convolve(A_35[1,:], A_34[:, 2], mode)) assert npo.allclose(ag_convolve(A_35, A_34, axes=([],[]), dot_axes=([0], [0]), mode=mode), npo.tensordot(A_35, A_34, axes=([0], [0]))) assert npo.allclose(ag_convolve(A_35, A_342, axes=([1],[2]), dot_axes=([0], [0]), mode=mode)[2], sum([sp_convolve(A_35[i, :], A_342[i, 2, :], mode) for i in range(3)])) assert npo.allclose(ag_convolve(A_2543, A_24232, axes=([1, 2],[2, 4]), dot_axes=([0, 3], [0, 3]), mode=mode)[2], sum([sum([sp_convolve(A_2543[i, :, :, j], A_24232[i, 2, :, j, :], mode) for i in range(2)]) for j in range(3)]))
def main(): """Initialize kernel, apply it to an image (via crosscorrelation, convolution).""" img = data.camera() kernel = np.array([ [-1, -2, -1], [0, 0, 0], [1, 2, 1] ]) cc_response = crosscorrelate(img, kernel) cc_gt = signal.correlate(img, kernel, mode="same") conv_response = convolve(img, kernel) conv_gt = signal.convolve(img, kernel, mode="same") util.plot_images_grayscale( [img, cc_response, cc_gt, conv_response, conv_gt], ["Image", "Cross-Correlation", "Cross-Correlation (Ground Truth)", "Convolution", "Convolution (Ground Truth)"] )
def conv4d(x, weights, bias, output): # print 'called' assert len(x.shape) == 4 and len(output.shape) == 4 batch_size, input_channel = x.shape[:2] output_batch_size, output_channel = output.shape[:2] num_filters, filter_channel = weights.shape[:2] assert batch_size == output_batch_size, '%d vs %d' % (batch_size, output_batch_size) assert output_channel == num_filters assert filter_channel == input_channel # func = convolve if true_conv else correlate for img_idx in range(batch_size): for c in range(output_channel): output[img_idx][c] = (correlate(x[img_idx], weights[c], mode='valid') + bias[c].reshape((1, 1, 1))) # if img_idx == 0 and c == 0: # print output[img_idx][c] # print bias[c].reshape((1, 1, 1))
def dietrich_baseline(bands, intensities, half_window=16, num_erosions=10): ''' Fast and precise automatic baseline correction of ... NMR spectra, 1991. http://www.sciencedirect.com/science/article/pii/002223649190402F http://www.inmr.net/articles/AutomaticBaseline.html ''' # Step 1: moving-window smoothing w = half_window * 2 + 1 window = np.ones(w) / float(w) Y = intensities.copy() if Y.ndim == 2: window = window[None] Y[..., half_window:-half_window] = convolve(Y, window, mode='valid') # Step 2: Derivative. dY = np.diff(Y) ** 2 # Step 3: Iterative thresholding. is_baseline = np.ones(Y.shape, dtype=bool) is_baseline[..., 1:] = iterative_threshold(dY) # Step 3: Binary erosion, to get rid of peak-tops. mask = np.zeros_like(is_baseline) mask[..., half_window:-half_window] = True s = np.ones(3, dtype=bool) if Y.ndim == 2: s = s[None] is_baseline = binary_erosion(is_baseline, structure=s, iterations=num_erosions, mask=mask) # Step 4: Reconstruct baseline via interpolation. if Y.ndim == 2: return np.row_stack([np.interp(bands, bands[m], y[m]) for y, m in zip(intensities, is_baseline)]) return np.interp(bands, bands[is_baseline], intensities[is_baseline])
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)
def polyphase_core(x, m, f): # x = input data # m = decimation rate # f = filter # Hack job - append zeros to match decimation rate if x.shape[0] % m != 0: x = np.append(x, np.zeros((m - x.shape[0] % m,))) if f.shape[0] % m != 0: f = np.append(f, np.zeros((m - f.shape[0] % m,))) polyphase = p = np.zeros((m, (x.shape[0] + f.shape[0]) / m), dtype=x.dtype) p[0, :-1] = np.convolve(x[::m], f[::m]) # Invert the x values when applying filters for i in range(1, m): p[i, 1:] = np.convolve(x[m - i::m], f[i::m]) return p
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
def _dense_convolve(Zi, ds): """Convolve Zi[k] and ds[k] for each atom k, and return the sum.""" return sum([signal.convolve(zik, dk) for zik, dk in zip(Zi, ds)], 0)
def convolve(img, kernel): """Apply a kernel/filter via convolution to an image. Args: img The image kernel The kernel/filter to apply Returns: New image """ return crosscorrelate(img, np.flipud(np.fliplr(kernel)))
def bandpassfilter(x,fs): """ :param x: a list of samples :param fs: sampling frequency :return: filtered list """ x = signal.detrend(x) b = signal.firls(129,[0,0.6*2/fs,0.7*2/fs,3*2/fs,3.5*2/fs,1],[0,0,1,1,0,0],[100*0.02,0.02,0.02]) return signal.convolve(x,b,'valid')
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
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
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
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
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)