Note
Click here to download the full example code
import numpy as np import matplotlib.pyplot as plt def exponential_moving_average_filter_recursive(y, a, b=0, gamma=1): """ Moving average filter with sliding window borders a and b. Note: this is a naive (slow) implementation and for demonstration purpose only. :param y: signal :param a: left window boundary :param b: right window boundary (optional) :param gamma: exponential decay factor 0..1 :return: filtered signal """ y_avg = np.zeros_like(y) # create empty array K = len(y) if (b>0): # correct initialization of y_sum for non-causal windows pass # ... here whe should initialize y_avg[0] correctly ... mu = (gamma**(b)-gamma**(a-1))/(1-1/gamma) for k in range(K): y_sum = 0 if (k+a >= 0) and (k+a < K): y_sum = y_sum - (gamma**(a-1))*y[k+a] if (k+1+b >= 0) and (k+1+b < K): y_sum = y_sum + (gamma**(b))*y[k+1+b] if (k+1 >= 0) and (k+1 < K): y_avg[k+1] = y_avg[k]* gamma**(-1) + y_sum / mu return y_avg # Load Signal K = 200 # number of samples to load k = range(K) # generate rectangular pulase y = np.concatenate((np.zeros(int(K/4)), np.ones(int(K/2)), np.zeros(int(K/4)))) # --- moving average filter and peak detection L=20 # averaging window size gamma = 1.1 y_ma1 = exponential_moving_average_filter_recursive(y, a=-2*L, b=0, gamma=1.00001) # non-symmetric filter y_ma2 = exponential_moving_average_filter_recursive(y, a=-L, b=L, gamma=1.001) # symmetric filter y_ma3 = exponential_moving_average_filter_recursive(y, a=-200*L, b=0, gamma=gamma) # non-symmetric filter # Plots _, axs = plt.subplots(1, 3, sharex='all', figsize=(10, 2)) plt.subplots_adjust(hspace=0.3) ax = axs[0] ax.plot(k, y, lw=0.5, c='tab:gray', label='$y$') ax.plot(k, y_ma1, lw=1.0, c='b', label='${y}_{avg}$') ax.set(ylabel=r'$y_k$') ax.set(xlabel='k') ax.set_title('asymmetric, L='+str(2*L)) # ax.legend(loc=4) ax = axs[1] ax.plot(k, y, lw=0.5, c='tab:gray', label='$y$') ax.plot(k, y_ma2, lw=1.0, ls='-', c='b', label='${y}_{avg}$') ax.set(xlabel='k') ax.set_title('symmetric, L='+str(L)) # ax.legend(loc=4) ax = axs[2] ax.plot(k, y, lw=0.5, c='tab:gray', label='$y$') ax.plot(k, y_ma3, lw=1.0, ls='-', c='b', label='${y}_{avg}$') ax.set(xlabel='k') ax.set_title('exp. asymmetric, gamma='+str(gamma)) # ax.legend(loc=4) plt.show()
Total running time of the script: ( 0 minutes 0.183 seconds)
Gallery generated by Sphinx-Gallery