2

FFT振幅スペクトルがあり、それから周期的なノイズ源（正弦波スプリアスなど）を選択的に通過させ、ランダムなバックグラウンドノイズに関連する周波数ビンをゼロにするフィルターを作成したいと思います。このフィルターが時間領域にIFFTで戻ると、周波数領域での急激な遷移によってリンギングアーティファクトが発生することを理解しています...したがって、このようなフィルターで遷移をスムーズにして、このようなことを回避する方法があるかどうか疑問に思っています。リンギング。

たとえば、FFTに1Mの周波数ビンがあり、バックグラウンドノイズフロアから5つのスプリアスが突き出ている場合、5つのスプリアスのそれぞれに関連付けられたピークビンを除くすべてのビンをゼロにします。問題は、時間領域でのアーティファクトを防ぐために、隣接するスパービンをどのように処理するかです。たとえば、スパービンの両側のビンを50％の振幅に設定する必要がありますか？スパービンの両側にある2つのビンを使用する必要がありますか（最も近いものは50％、次に近いものは25％など）？どんな考えでも大歓迎です。ありがとう！

4

1 に答える 1

2

私は次の方法が好きです：

  • 理想的な振幅スペクトルを作成します（DCに対して対称にすることを忘れないでください）
  • 時間領域への逆変換
  • ブロックをブロックサイズの半分だけ回転させます
  • ハンウィンドウを適用する

私はそれが適度に滑らかな周波数領域の結果を生み出すことを発見しましたが、あなたが示唆しているほど鋭いものでそれを試したことはありません。カイザーベッセルウィンドウを使用すると、おそらくよりシャープなフィルターを作成できますが、パラメーターを適切に選択する必要があります。もっとシャープにすると、サイドローブを6dB程度減らすことができると思います。

これがMatlab/Octaveコードのサンプルです。結果をテストするために、を使用しfreqz(h, 1, length(h)*10);ました。

function [ht, htrot, htwin] = ArbBandPass(N, freqs)
%# N = desired filter length
%# freqs = array of frequencies, normalized by pi, to turn into passbands
%# returns raw, rotated, and rotated+windowed coeffs in time domain

if any(freqs >= 1) || any(freqs <= 0)
    error('0 < passband frequency < 1.0 required to fit within (DC,pi)')
end

hf = zeros(N,1); %# magnitude spectrum from DC to 2*pi is intialized to 0
%# In Matlabs FFT, idx 1 -> DC, idx 2 -> bin 1, idx N/2 -> Fs/2 - 1, idx N/2 + 1 -> Fs/2, idx N -> bin -1
idxs = round(freqs * N/2)+1; %# indeces of passband freqs between DC and pi
hf(idxs) = 1; %# set desired positive frequencies to 1
hf(N - (idxs-2)) = 1; %# make sure 2-sided spectrum is symmetric, guarantees real filter coeffs in time domain
ht = ifft(hf); %# this will have a small imaginary part due to numerical error
if any(abs(imag(ht)) > 2*eps(max(abs(real(ht)))))
    warning('Imaginary part of time domain signal surprisingly large - is the spectrum symmetric?')
end
ht = real(ht); %# discard tiny imag part from numerical error
htrot = [ht((N/2 + 1):end) ; ht(1:(N/2))]; %# circularly rotate time domain block by N/2 points
win = hann(N, 'periodic'); %# might want to use a window with a flatter mainlobe
htwin = htrot .* win;
htwin = htwin .* (N/sum(win)); %# normalize peak amplitude by compensating for width of window lineshape
于 2010-09-30T06:57:48.750 に答える