Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

How to remove frequency from signal

I want to remove one frequency (one peak) from signal and plot my function without it. After fft I found frequency and amplitude and I am not sure what I need to do now. For example I want to remove my highest peak (marked with red dot on plot).

import numpy as np
import matplotlib.pyplot as plt

# create data
N = 4097
T = 100.0
t = np.linspace(-T/2,T/2,N)
f = np.sin(50.0 * 2.0*np.pi*t) + 0.5*np.sin(80.0 * 2.0*np.pi*t)

#plot function
plt.plot(t,f,'r')
plt.show()

# perform FT and multiply by dt
dt = t[1]-t[0]
ft = np.fft.fft(f) * dt      
freq = np.fft.fftfreq(N, dt)
freq = freq[:N/2+1]
amplitude = np.abs(ft[:N/2+1])
# plot results
plt.plot(freq, amplitude,'o-')
plt.legend(('numpy fft * dt'), loc='upper right')
plt.xlabel('f')
plt.ylabel('amplitude')
#plt.xlim([0, 1.4])


plt.plot(freq[np.argmax(amplitude)], max(amplitude), 'ro')
print "Amplitude: " + str(max(amplitude)) + "  Frequency: " + str(freq[np.argmax(amplitude)])

plt.show()
like image 894
Prefect Avatar asked Dec 14 '22 03:12

Prefect


1 Answers

One option is to transform the signal to the frequency domain then remove the selected frequency.

import numpy as np
import matplotlib.pyplot as plt
from scipy.fftpack import rfft, irfft, fftfreq, fft

# Number of samplepoints
N = 500
# sample spacing
T = 0.1

x = np.linspace(0.0, (N-1)*T, N)
# x = np.arange(0.0, N*T, T)  # alternate way to define x
y = 5*np.sin(x) + np.cos(2*np.pi*x) 

yf = fft(y)
xf = np.linspace(0.0, 1.0/(2.0*T), N//2)
#fft end

f_signal = rfft(y)
W = fftfreq(y.size, d=x[1]-x[0])

cut_f_signal = f_signal.copy()
cut_f_signal[(W>0.6)] = 0  # filter all frequencies above 0.6

cut_signal = irfft(cut_f_signal)

# plot results
f, axarr = plt.subplots(1, 3, figsize=(9, 3))
axarr[0].plot(x, y)
axarr[0].plot(x,5*np.sin(x),'g')

axarr[1].plot(xf, 2.0/N * np.abs(yf[:N//2]))
axarr[1].legend(('numpy fft * dt'), loc='upper right')
axarr[1].set_xlabel("f")
axarr[1].set_ylabel("amplitude")


axarr[2].plot(x,cut_signal)
axarr[2].plot(x,5*np.sin(x),'g')

plt.show()

enter image description here

enter image description here

like image 160
Prefect Avatar answered Dec 28 '22 07:12

Prefect