pythonscipyfilteringsignal-processinglowpass-filter

scipy.signal.firwin lowpass filter acts like highpass filter


When designing low/high/bandpass filters I encountered the problem that they do not work as expected (see code and output). I want to isolate one of the three frequency peaks (the one with the lowest frequency with lowpass filter, etc.). However, the wrong peak is isolated... Does anyone know what my mistake is?

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

    def create_single_freq_state(N, w):
        """
        Creates a state that only contains the frequencies in w. Outputs array with N entries.
        """
        res = np.exp(2.0*np.pi*w[:,np.newaxis]*np.arange(N)/N*1j)
        return np.sum(res, axis=0)

    def band_pass_filter(f1, f2, numtaps, window="hamming"):

        fil = signal.firwin(numtaps, [f1, f2], pass_zero=False, fs=2*numtaps, window=window, scale=False)
        return fil

    def high_pass_filter(f, numtaps, window="hamming"):

        fil = signal.firwin(numtaps, f, pass_zero=False, fs=2*numtaps, window=window, scale=False)
        return fil

    def low_pass_filter(f, numtaps, window="hamming"):

        fil = signal.firwin(numtaps, f, pass_zero=True, fs=2*numtaps, window=window, scale=False)
        return fil

    def plot_freq_response(fil):

        fs = 2*len(fil)
        w, h = signal.freqz(fil, 1, worN=2000)
        plt.plot((fs * 0.5 / np.pi) * w, abs(h), label="Hamming window")
        plt.show()

    def plot_func(f):

        x = np.arange(1, len(f)+1)  
        plt.plot(x, f, color ="red")  
        plt.show()

    #create state with several frequencies (those contained in the variable ws)
    N = 60
    ws = np.array([15, 30, 50])
    f = create_single_freq_state(N, ws)


    numtaps = 60-1
    freq1 = 19.9
    freq2 = 40.0

    #choose filter
    fil = low_pass_filter(freq1, numtaps)
    #fil = band_pass_filter(freq1, freq2, numtaps)
    #fil = high_pass_filter(freq2, numtaps)

    #plot signal in fourierspace
    plot_func(np.fft.fft(f))
    #plot frequency response
    plot_freq_response(fil)
    #plot result of the application of the filter
    plot_func(np.absolute(np.fft.fft(np.convolve(fil, f, mode="same"))))

Output:The output of my code. The wrong frequency is isolated.


Solution

  • Two of your frequencies are so high that they appear as DC and a low frequency (this is called aliasing, and is a consequence of the Nyquist sampling theorem). The fact that you're working with a complex input signal, and that your signal is short relative to your filter, make this hard to see.

    I changed ws to be ws = np.array([7.5, 15, 25]), and changed your first call to plot_func to be plot_func(np.abs(np.fft.fft(f))) (because one of the peaks is entirely in the complex part of the transform). Now the input and output graphs look as follows: input with three peaks output with one peak