3
我目前正在嘗試使用matlab將簡單逆濾波器與維納濾波器進行反捲積比較。我的起始信號是exp(-t^2)，這是用一個不等於-5的時間長度來卷積的。我正在引入振幅範圍爲-5到.5的噪音。將Naive逆濾波器與Wiener濾波器進行比較，以便在Matlab中解卷積
定義我的時域到頻域映射：
f = exp(-t^2) => F
s = rect => R
c = f*s => C
r = noise (see above) => R
with noise c becomes: c = f*s + n => C = FxS + N
對於第一種方法，我只是服用FT的c並通過f的FT劃分它，然後做逆FT。這相當於s = (approx.) ifft((FxS + N)/F)。
對於第二種方法，我使用維納濾波器W，並將其乘以C/R，然後進行逆FT。這相當於S = (approx.) ifft(CxW/R)。
維納過濾器是W = mag_squared(FxS)/(mag_squared(FxS) + mag_squared(N))。
我用'*'表示卷積，'x'表示乘法。
我想比較在時間間隔-3到3之間的矩形的兩個解卷積。 現在，我的解卷矩形的結果圖看起來與原始視圖沒有任何區別。
有人能指出我正確的方向，我做錯了什麼？我曾嘗試使用ifftshift和不同的順序，但似乎沒有任何工作。
感謝
我的MATLAB代碼如下：
%%using simple inverse filter
dt = 1/1000;
t = linspace(-3,3,1/dt); %time
s = zeros(1,length(t));
s(t>=-0.5 & t<=0.5) = 1; %rect
f = exp(-(t.^2)); %function
r = -.5 + rand(1,length(t)); %noise
S = fft(s);
F = fft(f);
R = fft(r);
C = F.*S + R;
S_temp = C./F;
s_recovered_1 = real(ifft(S_temp)); %correct?...works for signal without R (noise)
figure();
plot(t,s + r);
title('rect plus noise');
figure();
hold on;
plot(t,s,'r');
plot(t,f,'b');
legend('rect input','function');
title('inpute rect and exponential functions');
hold off;
figure();
plot(t,s_recovered_1,'black');
legend('recovered rect');
title('recovered rect using naive filter');
%% using wiener filter
N = length(s);
I_mag = abs(I).^2;
R_mag = abs(R).^2;
W = I_mag./(I_mag + R_mag);
S_temp = (C.*W)./F;
s_recovered_2 = abs(ifft(S_temp));
figure();
freq = -fs/2:fs/N:fs/2 - fs/N;
hold on;
plot(freq,10*log10(I_mag),'r');
plot(freq,10*log10(R_mag),'b');
grid on
legend('I_mag','R_mag');
title('Periodogram Using FFT')
xlabel('Frequency (Hz)')
ylabel('Power/Frequency (dB/Hz)')
figure();
plot(t,s_recovered_2);
legend('recovered rect');
title('recovered rect using wiener filter');
去除噪音和計算簡單逆濾波器揭示了原始矩形（假設我做's_recovered_1 =真實（ifft（S_temp））;'我已經改變了上述代碼來反映）。我期望簡單逆濾波器的輸出能夠給出很大的值，因爲我大部分被小的值所分割。這或多或少地匹配我實際得到的輸出。我認爲我現在的主要問題是計算維納濾波器。我已更新我的代碼以反映我現在的想法，但我對此非常不確定，它仍然不會產生任何類似原始矩形的東西。 – 2014-10-11 22:29:56
我也嘗試通過直接計算我認爲是I和R的雙側功率譜密度來計算維納濾波器。我更新了上面的代碼以反映這一點。我現在得到了一個像sinc的東西。所以這可能會更好，但它還沒有結束。 – 2014-10-11 23:42:49
我也嘗試過兩次關於維納濾波，這是愚蠢的，我知道，但它給出了正確的幅度矩（因爲第一個ifft是一個sinc），但寬度是錯誤的... 我已經更新了上面的代碼以顯示此行。 – 2014-10-11 23:48:50