Lowpass : What is wrong with average filtering?

Question

The question is : What is wrong with averaging as low pass filter ?

The details : I want to lowpass filter a signal to downsample it. The constraints are : I have no RAM available and I work in streaming, therefore I cannot use a weighting window to filter. Consequently, I have to do either an average on the last X samples, or a classic filter like :

signal_t+1 = signal_t * (1-factor) + new_measurement * factor

with factor a value adapted to the sampling. What is more, I know there are signals around a frequency just above the downsampling frequency.

Is it better to use the average or to use the 1-order filter described above ?

More details and tests : I work with real samples, not complex. Sampling frequency = 200 Hz, frequency after downsampling = 4Hz.

Let's compare a 1-order filter, with a cutoff frequency of 1Hz to limitate aliasing, with two averaging filters. The first averaging filter is the average of 50 samples to downsample from 200Hz to 4Hz. The second averaging filter is an average on 66 samples to get as much rejection as with the 1-order filter.

I think that the filter "average on 66 samples" is the best : The 50-samples average filter has a rejection higher than the 1-order filter, and I really need to limitate aliasing. But the 66-samples filter has a cutoff frequency larger than the 1-order filter and a rejection at least equal after 2Hz. What is more, the 66-samples filter cuts the frequencies around 3Hz and I know there will be signals on 3Hz that will make my signal noisy after the 4Hz downsampling - I work with real samples, not complex.

Then, I would chose the 66-samples averaging filter.

Is it something wrong in the logic above ? Would you have suggestions ?

Here is the code :

import numpy as np
import matplotlib.pyplot as plt

nfft = 8192
fsampl = 200
dtsampl = 1 / 200
fcutoff = 1
dtcutoff = 1 / (2 * np.pi * fcutoff)
faverage1 = 4
faverage2 = 3

print("generate signal")
raw_signal = np.zeros(nfft)
raw_signal[int(nfft/2)]=1

print("filter signal")
filtered_signal = 0.*np.copy(raw_signal)
for i in range(1, np.size(raw_signal)):
    filtered_signal[i] = filtered_signal[i-1]+ (raw_signal[i-1]-filtered_signal[i-1]) * dtsampl / dtcutoff

averaged_signal_1 = np.copy( raw_signal )
num_average = fsampl / faverage1
for i in range(1, np.size(raw_signal)):
    averaged_signal_1[i] = np.average( raw_signal[int( max( 0, i - num_average ) ):i] )

averaged_signal_2 = np.copy( raw_signal )
num_average = fsampl / faverage2
for i in range(1, np.size(raw_signal)):
    averaged_signal_2[i] = np.average( raw_signal[int( max( 0, i - num_average ) ):i] )

print("get FFT")
raw_signal_fft = np.fft.rfft(raw_signal)
filtered_signal_fft = np.fft.rfft(filtered_signal)
averaged_signal1_fft = np.fft.rfft( averaged_signal_1)
averaged_signal2_fft = np.fft.rfft( averaged_signal_2)

raw_signal_fftlog = 20*np.log10(np.abs(raw_signal_fft))
filtered_signal_fftlog = 20*np.log10(np.abs(filtered_signal_fft))
averaged_signal1_fftlog = 20*np.log10(np.abs(averaged_signal1_fft))
averaged_signal2_fftlog = 20*np.log10(np.abs(averaged_signal2_fft))

print("plot")
freq_axe= np.linspace( 0, fsampl / 2, int( nfft / 2 ) + 1 )
plt.plot(freq_axe, raw_signal_fftlog, '-')
plt.plot(freq_axe, filtered_signal_fftlog, 'r-')
plt.plot(freq_axe, averaged_signal1_fftlog, '-')
plt.plot(freq_axe, averaged_signal2_fftlog, '-')
plt.plot([0.02, 100], [-3, -3], 'k:')
plt.legend(['raw signal (dirac)', '1-order 1Hz filter', 'average on 50', 'average on 66'],  loc="upper right")
plt.ylabel('amplitude (dB)')
plt.xlabel('frequencies (Hz)')
plt.xlim([0,8])
plt.ylim([-30,2])

plt.show()

Answer 1

The question was : What is best between mooving average and 1st order filter when no RAM is available ?

The solution would be : None of them, except for a specific need : rejection of a particular frequency. Otherwise, it is better to use filter of higher orders than 1st order.

• Mooving average :
  • Drawbacks : delay and poor rejection.
  • Assets : huge rejection at several particular frequencies
• 1-order filter
  • Drawbacks : rejection which is poor whatever frequency
  • Assets : no delay
• Higher order filters. Higher order filters include : cascaded 1st order filters (here under an exemple with 5 cascaded 1st order filters), butterworth filters, etc.
  • Assets : better rejection than 1st order filter and no need of much RAM. Small delay.

View of the result :

Thank you @user883521 and @Florent Ecochard for your help !

Here is the Pyhton code for a 5th order filter (for other orders, cf links here under. for other filters, see signal.lfilter)

import numpy as np
import matplotlib.pyplot as plt

nfft = 8192
fsampl = 200
dtsampl = 1 / 200

fcutoff_1st_order = 1
dtcutoff_1storder = 1 / (2 * np.pi * fcutoff_1st_order)
fcutoff_5thorder =2.635
faverage2 = 3

print("generate signal")
raw_signal = np.zeros(nfft)
raw_signal[int(nfft/2)]=1

print("filter signal")
filtered_signal = 0.*np.copy(raw_signal)
for i in range(1, np.size(raw_signal)):
    filtered_signal[i] = filtered_signal[i-1]+ (raw_signal[i-1]-filtered_signal[i-1]) * dtsampl / dtcutoff_1storder
filtered_signal5= np.copy(raw_signal)
x = np.exp( -2 * np.pi * fcutoff_5thorder / fsampl )
for i in range(5, np.size(raw_signal)):
    filtered_signal5[i] = np.power( 1 - x, 5 ) * raw_signal[i - 1] \
                          + 5 * np.power( x, 1 ) * filtered_signal5[i - 1] + \
                          -10 * np.power( x, 2 ) * filtered_signal5[i - 2] + \
                          +10 * np.power( x, 3 ) * filtered_signal5[i - 3] + \
                          -5 * np.power( x, 4 ) * filtered_signal5[i - 4] + \
                          +np.power( x, 5 ) * filtered_signal5[i-5]
averaged_signal_66 = np.copy( raw_signal )
num_average = fsampl / faverage2
for i in range(1, np.size(raw_signal)):
    averaged_signal_66[i] = np.average( raw_signal[int( max( 0, i - num_average ) ):i] )

print("get FFT")
raw_signal_fft = np.fft.rfft(raw_signal)
filtered_signal_fft = np.fft.rfft(filtered_signal)
filtered_signal5_fft = np.fft.rfft(filtered_signal5)
averaged_signal66_fft = np.fft.rfft( averaged_signal_66 )
raw_signal_fftlog = 20*np.log10(np.abs(raw_signal_fft))
filtered_signal_fftlog = 20*np.log10(np.abs(filtered_signal_fft))
filtered_signal5_fftlog = 20*np.log10(np.abs(filtered_signal5_fft))
averaged_signal66_fftlog = 20 * np.log10( np.abs( averaged_signal66_fft ) )

print("plot")
freq_axe= np.linspace( 0, fsampl / 2, int( nfft / 2 ) + 1 )
plt.plot(freq_axe, raw_signal_fftlog, '-')
plt.plot(freq_axe, filtered_signal_fftlog, '-')
plt.plot(freq_axe, filtered_signal5_fftlog, '-')
plt.plot( freq_axe, averaged_signal66_fftlog, '-' )
plt.plot([0.02, 100], [-3, -3], 'b--')
plt.plot([2, 2],[-50, 0],  'k:')
plt.plot([3, 3],[-50, 0],  'k--')
plt.legend(['raw signal (dirac)', '1-order 1Hz filter', '5-order 1Hz filter', 'mooving average 66 samples',
            '-3dB', 'aliasing freq', 'unwanted signals frequency'],  loc="upper right")
plt.ylabel('amplitude (dB)')
plt.xlabel('frequency (Hz)')
plt.xlim([0,8])
plt.ylim([-30,2])

plt.show()

useful links : http://www.analog.com/media/en/technical-documentation/dsp-book/dsp_book_Ch14.pdf http://www.analog.com/media/en/technical-documentation/dsp-book/dsp_book_Ch15.pdf http://www.analog.com/media/en/technical-documentation/dsp-book/dsp_book_Ch19.pdf