CODEWITHSUNDEEP
pythonpython-3.x
Jayson Herrera
Jayson Herrera

Reputation: 13

function to determine the frequency of a sinusoid

I am having trouble with my Digital Signal Processing homework. Using Python, I need to create a function that is able to determine the frequency of a sinusoid. I am given random frequencies form 0-4000 Hz with an Fs=8000. Can someone please help?

import numpy as np

def freqfinder(signal):
    """REPLACE"""
    x=np.fft.fft(signal)
    x=np.abs(x)
    x=np.max(x)
    return x


t=np.linspace(0,2*np.pi,8*8000)    
y=np.sin(2*t)

print(freqfinder(y))

z = np.fft.fft(y)
zz = np.abs(z)

plt.plot(zz)

I tried this as a test for the fft.

Upvotes: 1

Views: 2148

Answers (1)

SyntaxVoid
SyntaxVoid

Reputation: 2633

Your code is off to a good start. A few things to note:

  1. You should only look at the first half of your FFT -- For a REAL input, the output is symmetric around 0 and you only care about the frequencies greater than 0 (the first half of the fft output).

  2. You want the magnitude of each frequency - so you should then take the absolute value of the resulting fft.

  3. The max you are locating is NOT the frequency, but is related to the index of the frequency. It is the strength of the strongest frequency.

Here is a little script demonstrating these ideas:

import numpy as np
import matplotlib.pyplot as plt

fs = 8000
t = np.linspace(0, 2*np.pi, fs)
freqs = [  2, 152, 423, 2423, 3541] # Frequencies to test
amps  = [0.5, 0.5, 1.0, 0.8,   0.3] # Amplitude for each freq
y = np.zeros(len(t))
for freq, amp in zip(freqs, amps):
    y += amp*np.sin(freq*t)


fig, ax = plt.subplots(1, 2)
ax = ax.flatten()
ax[0].plot(t, y)
ax[0].set_title("Original signal")

y_fft = np.fft.fft(y)           # Original FFT
y_fft = y_fft[:round(len(t)/2)] # First half ( pos freqs )
y_fft = np.abs(y_fft)           # Absolute value of magnitudes
y_fft = y_fft/max(y_fft)        # Normalized so max = 1

freq_x_axis = np.linspace(0, fs/2, len(y_fft))
ax[1].plot(freq_x_axis, y_fft, "o-")
ax[1].set_title("Frequency magnitudes")
ax[1].set_xlabel("Frequency")
ax[1].set_ylabel("Magnitude")
plt.grid()
plt.tight_layout()
plt.show()

f_loc = np.argmax(y_fft) # Finds the index of the max
f_val = freq_x_axis[f_loc] # The strongest frequency value
print(f"The strongest frequency is f = {f_val}")

The output: enter image description here

The strongest frequency is f = 423.1057764441111

You can see on the right graph that there is a peak at each of the frequencies we specified in freqs, which is what is expected.

This kind of setup is fine if you only have one frequency you're looking for, but otherwise you may need to find and implement some peak finding algorithms to find all the indices of all the frequency peaks of y_fft and then correlate that with the frequencies in freq_x_axis

Upvotes: 1

Related Questions