2
$\begingroup$

I have a signal for which I need to calculate the magnitude and phase at 200 Hz frequency only. I would like to use Fourier transform for it. I am very new to signal processing. And this is my first time using a Fourier transform.

I found that I can use the scipy.fftpack.fft to calculate the FFT of the signal. Then use numpy.mag and numpyh.phase to calculate the magnitude and phases of the entire signal. But I would like to get the magnitude and phase value of the signal corresponding to 200 Hz frequency only. How can I do this using Python?

So far I have done.

from scipy.fftpack import fft
import numpy as np

fft_data = fft(signal)
magnitude = np.mag(fft_data)
phase = np.phase(fft_data)
$\endgroup$
4
  • 1
    $\begingroup$ Welcome to SE.SP! Do you know the sampling frequency of your data? $\endgroup$
    – Peter K.
    Dec 14 '20 at 19:12
  • $\begingroup$ the sampling frequency is 2 MHz $\endgroup$
    – thileepan
    Dec 15 '20 at 1:41
  • $\begingroup$ What do you have to separate it from? What is the rest of the signal that is not 200Hz? $\endgroup$
    – IanJ
    Dec 15 '20 at 3:08
  • $\begingroup$ I'm not interested in other frequencies. $\endgroup$
    – thileepan
    Dec 15 '20 at 3:25
1
$\begingroup$

You can find the index of the desired (or the closest one) frequency in the array of resulting frequency bins using np.fft.fftfreq function, then use np.abs and np.angle functions to get the magnitude and phase.

Here is an example using fft.fft function from numpy library for a synthetic signal.

import numpy as np
import matplotlib.pyplot as plt

# Number of sample points
N = 1000

# Sample spacing
T = 1.0 / 800.0     # f = 800 Hz

# Create a signal
x = np.linspace(0.0, N*T, N)
t0 = np.pi/6   # non-zero phase of the second sine
y = np.sin(50.0 * 2.0*np.pi*x) + 0.5*np.sin(200.0 * 2.0*np.pi*x + t0)
yf = np.fft.fft(y) # to normalize use norm='ortho' as an additional argument

# Where is a 200 Hz frequency in the results?
freq = np.fft.fftfreq(x.size, d=T)
index, = np.where(np.isclose(freq, 200, atol=1/(T*N)))

# Get magnitude and phase
magnitude = np.abs(yf[index[0]])
phase = np.angle(yf[index[0]])
print("Magnitude:", magnitude, ", phase:", phase)

# Plot a spectrum 
plt.plot(freq[0:N//2], 2/N*np.abs(yf[0:N//2]), label='amplitude spectrum')   # in a conventional form
plt.plot(freq[0:N//2], np.angle(yf[0:N//2]), label='phase spectrum')
plt.legend()
plt.grid()
plt.show()

And here is a useful manual with detailed explanations: reference.

$\endgroup$
3
  • $\begingroup$ thank you for the clear explanation. I have two doubts. 1. What should I do if I want to change the phase of the signals to 30 degrees? 2. Why are you not scaling the magnitude values with (2/N)? $\endgroup$
    – thileepan
    Dec 18 '20 at 11:11
  • 1
    $\begingroup$ @thileepan 1. The phase t0 would be an additional term in the argument of a sine: A*sin(wt+t0). t0 = np.pi/6 should shift the signal to 30 degrees. 2. The example shows the default fft results. You can normalize the magnitude by setting the "norm" parameter like this: yf = np.fft.fft(y, norm='ortho'). Btw, my bad, np.isclose does not work as intended. I will fix it in the answer. $\endgroup$
    – megasplash
    Dec 18 '20 at 14:03
  • $\begingroup$ @thileepan Regarding the normalization of the results, you may find this question to be helpful. $\endgroup$
    – megasplash
    Dec 18 '20 at 14:33

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.