Python: высокочастотный FIR-фильтр с помощью Windowing
Я хотел бы создать базовый высокочастотный FIR-фильтр, используя Windowing в Python.
Мой код ниже и намеренно идиоматичен - я знаю, что вы можете (скорее всего) завершить это с помощью одной строки кода на Python, но я учусь. Я использовал базовую функцию sinc с прямоугольным окном: мой вывод работает для сигналов, которые являются аддитивными (f1 + f2), но не мультипликативными (f1 * f2), где f1 = 25 кГц и f2 = 1 МГц.
Мои вопросы: я неправильно понял что-то фундаментальное или мой код неверен? Таким образом, я хотел бы извлечь только сигнал верхних частот (f2 = 1 МГц) и отфильтровать все остальное. Я также включил снимки экрана того, что генерируется для (f1 + f2) и (f1 * f2):
import numpy as np
import matplotlib.pyplot as plt
# create an array of 1024 points sampled at 40MHz
# [each sample is 25ns apart]
Fs = 40e6
T = 1/Fs
t = np.arange(0,(1024*T),T)
# create an ip signal sampled at Fs, using two frequencies
F_low = 25e3 # 25kHz
F_high = 1e6 # 1MHz
ip = np.sin(2*np.pi*F_low*t) + np.sin(2*np.pi*F_high*t)
#ip = np.sin(2*np.pi*F_low*t) * np.sin(2*np.pi*F_high*t)
op = [0]*len(ip)
# Define -
# Fsample = 40MHz
# Fcutoff = 900kHz,
# this gives the normalised transition freq, Ft
Fc = 0.9e6
Ft = Fc/Fs
Length = 101
M = Length - 1
Weight = []
for n in range(0, Length):
if( n != (M/2) ):
Weight.append( -np.sin(2*np.pi*Ft*(n-(M/2))) / (np.pi*(n-(M/2))) )
else:
Weight.append( 1-2*Ft )
for n in range(len(Weight), len(ip)):
y = 0
for i in range(0, len(Weight)):
y += Weight[i]*ip[n-i]
op[n] = y
plt.subplot(311)
plt.plot(Weight,'ro', linewidth=3)
plt.xlabel( 'weight number' )
plt.ylabel( 'weight value' )
plt.grid()
plt.subplot(312)
plt.plot( ip,'r-', linewidth=2)
plt.xlabel( 'sample length' )
plt.ylabel( 'ip value' )
plt.grid()
plt.subplot(313)
plt.plot( op,'k-', linewidth=2)
plt.xlabel( 'sample length' )
plt.ylabel( 'op value' )
plt.grid()
plt.show()