用于数据平滑的卷积算法

convolution algorithm for data smoothing

所以我决定编写自己的卷积来平滑我的数据,其作用与 np.convolve 相同。唯一的问题是我得到的振幅比我预期的要高一点。我不知道我应该把我的数据分成哪一部分才能得到正确的结果。

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

def dotProduct(a,b):
    dot = 0
    if (len(a) == len(b)):
        for i in range(len(a)):
            dot += a[i] * b[i]
    elif (len(a) > len(b)):
        for i in range(len(b)):
            dot += a[i] * b[i]
    elif (len(b) > len(a)):
        for i in range(len(a)):
            dot += a[i] * b[i]

    return dot


def convolution(signal,magnitude=3):
    kernel = np.linspace(0,1,4) * magnitude
    convolution = np.zeros(len(signal))
    sameSizeKernel = np.zeros(len(signal))
    for i in range(len(convolution)):
        if ((len(kernel)+i) <= len(convolution)):
            sameSizeKernel[i:len(kernel)+i] = kernel
            convolution[i+int(len(kernel)/2)] = dotProduct(signal,sameSizeKernel)/4
        else:
            kernel = kernel[:len(kernel)-1]
            sameSizeKernel[i:len(kernel)+i] = kernel
            convolution[i+int(len(kernel)/2)] = dotProduct(signal,sameSizeKernel)/4
    return convolution


t = np.linspace(0, 1, 500)
triangle = signal.sawtooth(2 * np.pi * 5 * t, 0.5)
conv_ = convolution(triangle,10)

plt.figure(figsize=(15,8))
plt.plot(triangle, label='signal')
plt.plot(conv_, label='convolution')
plt.legend(loc=1)

所以我自己找到了答案:) 问题是,如果您想看到正确的结果,内核数组的总和应该不超过 1。因此,为了规范化,我们应该将每个点除以内核索引的总和。

convolution[i+int(len(kernel)/2)] = dotProduct(signal,sameSizeKernel)/np.sum(kernel)

完整代码在这里:

import numpy as np
import matplotlib.pyplot as plt
from scipy import signal

def dotProduct(a,b):
    dot = 0
    if (len(a) == len(b)):
        for i in range(len(a)):
            dot += a[i] * b[i]
    elif (len(a) > len(b)):
        for i in range(len(b)):
            dot += a[i] * b[i]
    elif (len(b) > len(a)):
        for i in range(len(a)):
            dot += a[i] * b[i]

    return dot

dotProduct(a,b)

def convolution(signal,magnitude=3):
    kernel = np.linspace(0,magnitude,10)
    convolution = np.zeros(len(signal))
    sameSizeKernel = np.zeros(len(signal))
    for i in range(len(convolution)):
        if ((len(kernel)+i) <= len(convolution)):
            sameSizeKernel[i:len(kernel)+i] = kernel
            convolution[i+int(len(kernel)/2)] = dotProduct(signal,sameSizeKernel)/np.sum(kernel)
        else:
            kernel = kernel[:len(kernel)-1]
            sameSizeKernel[i:len(kernel)+i] = kernel
            convolution[i+int(len(kernel)/2)] = dotProduct(signal,sameSizeKernel)/np.sum(kernel)
    return convolution


t = np.linspace(0, 1, 50)
triangle = np.zeros(len(t))
triangle[10:20] = 10
conv_ = convolution(triangle,10)

plt.figure(figsize=(15,8))
plt.plot(triangle, label='signal')
plt.plot(conv_, label='convolution')
plt.legend(loc=1)