Bagaimana cara menulis filter lowpass untuk sinyal sampel dalam Python?

16

Saya memiliki beberapa sinyal yang diambil sampelnya masing-masing 1 ns (1e-9 dtk) dan memiliki, katakanlah, 1e4 poin. Saya perlu memfilter frekuensi tinggi dari sinyal ini. Katakanlah saya perlu memfilter frekuensi lebih tinggi dari 10 MHz. Saya ingin agar frekuensi yang lebih rendah dari sinyal frekuensi cutoff akan dilewati tidak berubah. Ini berarti penguatan filter akan menjadi 1 untuk frekuensi yang lebih rendah dari frekuensi cutoff. Saya ingin dapat menentukan urutan filter. Maksud saya, filter urutan pertama memiliki 20 db / dekade kemiringan (power roll off) setelah frekuensi cutoff, filter urutan kedua memiliki 40 db / dec slope setelah frekuensi cutoff dan sebagainya. Kinerja kode yang tinggi itu penting.

Alex
sumber

Jawaban:

19

Respons frekuensi untuk filter yang dirancang menggunakan fungsi mentega adalah:

Respon Filter Butterworth

Tetapi tidak ada alasan untuk membatasi filter ke desain filter monotonik konstan. Jika Anda menginginkan atenuasi yang lebih tinggi di pita transisi stopband dan curam, ada opsi lain. Untuk informasi lebih lanjut tentang menentukan filter menggunakan iirdesing lihat ini . Seperti yang ditunjukkan oleh plot respons frekuensi untuk desain mentega frekuensi cutoff (-3dB point) jauh dari tujuan. Ini dapat dikurangi dengan pengambilan sampel sebelum disaring (fungsi desain akan mengalami kesulitan dengan filter yang sempit, 2% dari bandwidth). Mari kita lihat penyaringan laju sampel asli dengan cutoff yang ditentukan.

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

from scipy.signal import fir_filter_design as ffd
from scipy.signal import filter_design as ifd

# setup some of the required parameters
Fs = 1e9           # sample-rate defined in the question, down-sampled

# remez (fir) design arguements
Fpass = 10e6       # passband edge
Fstop = 11.1e6     # stopband edge, transition band 100kHz
Wp = Fpass/(Fs)    # pass normalized frequency
Ws = Fstop/(Fs)    # stop normalized frequency

# iirdesign agruements
Wip = (Fpass)/(Fs/2)
Wis = (Fstop+1e6)/(Fs/2)
Rp = 1             # passband ripple
As = 42            # stopband attenuation

# Create a FIR filter, the remez function takes a list of 
# "bands" and the amplitude for each band.
taps = 4096
br = ffd.remez(taps, [0, Wp, Ws, .5], [1,0], maxiter=10000) 

# The iirdesign takes passband, stopband, passband ripple, 
# and stop attenuation.
bc, ac = ifd.iirdesign(Wip, Wis, Rp, As, ftype='ellip')  
bb, ab = ifd.iirdesign(Wip, Wis, Rp, As, ftype='cheby2') 

Filter Nilai Sampel Asli

Seperti yang disebutkan, karena kami mencoba memfilter persentase kecil dari bandwidth, filter tidak akan memiliki cutoff yang tajam. Dalam hal ini, filter lowpass, kita dapat mengurangi bandwidth untuk mendapatkan filter yang terlihat lebih baik. Fungsi resample python / scipy.signal dapat digunakan untuk mengurangi bandwidth.

Perhatikan fungsi resample akan melakukan penyaringan untuk mencegah alias. Prefilter juga dapat dibuat (untuk mengurangi alias) dan dalam hal ini kita bisa dengan sederhana melakukan resample sebanyak 100 dan dilakukan , tetapi pertanyaannya adalah tentang membuat filter. Untuk contoh ini, kami akan menurunkan sampel sebanyak 25 dan membuat filter baru

R = 25;            # how much to down sample by
Fsr = Fs/25.       # down-sampled sample rate
xs = signal.resample(x, len(x)/25.)

Jika kami memperbarui parameter desain untuk filter FIR, respons yang baru adalah.

# Down sampled version, create new filter and plot spectrum
R = 25.             # how much to down sample by
Fsr = Fs/R          # down-sampled sample rate
Fstop = 11.1e6      # modified stopband
Wp = Fpass/(Fsr)    # pass normalized frequency
Ws = Fstop/(Fsr)    # stop normalized frequency
taps = 256
br = ffd.remez(taps, [0, Wp, Ws, .5], [1,0], maxiter=10000) 

Respons Filter Downsampled

Filter yang beroperasi pada data downsampled memiliki respons yang lebih baik. Manfaat lain menggunakan filter FIR adalah bahwa Anda akan memiliki respons fase linier.

Christopher Felton
sumber
1
Terima kasih. Bagaimana Anda membuat grafik spektrum sinyal?
Alex
Terima kasih banyak atas jawaban yang sangat baik! Saya ingin tahu apakah Anda dapat menjelaskan cara menerapkan filter FIR berdasarkan koefisien yang dihitung menggunakan Remez? Saya mengalami kesulitan memahami apa yang filtfiltdiinginkan untuk aparameter.
ali_m
Setelah Anda memiliki koefisien dari desain filter, ( b untuk FIR b dan a untuk IIR), Anda dapat menggunakan beberapa fungsi yang berbeda untuk melakukan penyaringan: lfilter , convolve , filtfilt . Biasanya semua fungsi ini beroperasi serupa: y = filtfilt (b, a, x) Jika Anda memiliki filter FIR cukup atur a = 1 , x adalah sinyal input, b adalah koefisien FIR. Posting ini mungkin membantu juga.
Christopher Felton
5

Apakah ini berhasil?

from __future__ import division
from scipy.signal import butter, lfilter

fs = 1E9 # 1 ns -> 1 GHz
cutoff = 10E6 # 10 MHz
B, A = butter(1, cutoff / (fs / 2), btype='low') # 1st order Butterworth low-pass
filtered_signal = lfilter(B, A, signal, axis=0)

Namun, Anda benar, dokumentasinya tidak lengkap. Sepertinya butterini adalah pembungkus untuk iirfilter, yang lebih baik didokumentasikan :

N: int Urutan filter. Wn: array_like Urutan skalar atau panjang-2 yang memberikan frekuensi kritis.

Namun, sebagian besar barang ini diklon dari matlab, jadi Anda dapat melihat dokumentasi mereka juga:

frekuensi cutoff dinormalisasi Wn harus angka antara 0 dan 1, di mana 1 sesuai dengan frekuensi Nyquist, π radian per sampel.

Memperbarui:

Saya menambahkan dokumentasi untuk fungsi-fungsi ini. :) Github membuatnya mudah.

endolith
sumber
1

Tidak yakin apa aplikasi Anda, tetapi Anda mungkin ingin memeriksa Gnuradio: http://gnuradio.org/doc/doxygen/classgr__firdes.html

Blok pemrosesan sinyal ditulis dalam C ++ (walaupun grafik aliran Gnuradio menggunakan Python), tetapi Anda mengatakan kinerja tinggi itu penting.

bungkus bungkus
sumber
1

Saya mendapatkan hasil yang baik dengan filter FIR ini. Pemberitahuan itu menerapkan filter dua kali, maju "maju" dan "mundur", sehingga untuk mengimbangi offset sinyal ( filtfiltfungsi tidak bekerja, tidak tahu mengapa):

def firfilt(interval, freq, sampling_rate):
    nfreq = freq/(0.5*sampling_rate)
    taps =  sampling_rate + 1
    a = 1
    b = scipy.signal.firwin(taps, cutoff=nfreq)
    firstpass = scipy.signal.lfilter(b, a, interval)
    secondpass = scipy.signal.lfilter(b, a, firstpass[::-1])[::-1]
    return secondpass

Sumber yang bagus untuk memfilter desain dan penggunaan, dari mana saya mengambil kode ini, dan dari mana contoh band-pass dan hi-pass filter dapat diambil, adalah INI .

heltonbiker
sumber
Saya tidak percaya ada banyak manfaat memajukan dan membalikkan filter FIR. Filter IIR bisa mendapatkan keuntungan dari maju / mundur (filtfilt) karena Anda bisa mendapatkan fase linier dari filter fase non-linier dengan menyaring balik.
Christopher Felton
2
@ChristopherFelton Saya baru saja membalikkan untuk menyinkronkan sinyal elektromiografi RAW dengan versi itu sendiri. Saya tahu saya hanya bisa menggeser sinyalnya, tetapi menyaring dua kali akhirnya menjadi lebih sedikit masalah. Perlu diperhatikan bahwa pass kedua hampir tidak mengubah pass pertama yang sudah difilter ... Terima kasih telah mencatat!
heltonbiker
Ahh ya Untuk menghapus penundaan (keterlambatan grup), poin bagus.
Christopher Felton
1

Saya tidak punya hak komentar ...

@endolith: Saya menggunakan yang sama seperti Anda kecuali menggunakan scipy.signal.filtfilt (B, A, x) di mana x adalah vektor input untuk disaring - misalnya numpy.random.normal (size = (N)) . filtfilt membuat sinyal maju dan mundur. Demi kelengkapan (sebagian besar sama dengan @endolith):

import numpy as np
import scipy.signal as sps

input = np.random.normal(size=(N)) # Random signal as example
bz, az = sps.butter(FiltOrder, Bandwidth/(SamplingFreq/2)) # Gives you lowpass Butterworth as default
output = sps.filtfilt(bz, az, input) # Makes forward/reverse filtering (linear phase filter)

filtfilt seperti juga disarankan oleh @heltonbiker membutuhkan array koefisien yang saya percaya. Jika Anda perlu melakukan pemfilteran bandpass di baseband kompleks, konfigurasi yang lebih diperlukan diperlukan tetapi ini tampaknya tidak menjadi masalah di sini.

Lars1
sumber