频域实现维纳滤波，求助。。。

频域实现维纳滤波，求助。。。

% z为观测数据

% s为原始数据

% y为输出信号

% N为z，s的长度

clear;clc;

load chirp;

whos;% 显示变量y 、F s ,其中,向量y :表示创建的信号；标量Fs:表示以Hz为单位的频率。

s=y;

N=length(s);

z=awgn(s,10,'measured',1,'linear');% 加高斯白噪声，强度10dB 此为观测信号

window=kaiser(N);% 用beta=0.5时的凯泽窗

Gs_f=periodogram(s,window,'two sided',N);% 原始信号功率谱估计(周期图谱法)

Gz_f=periodogram(z,window,'two sided',N);% 观测信号功率谱估计(周期图谱法)

%频域法求维纳滤波传递函数

H_opt = Gs_f./Gz_f;

Y=H_opt.*fftn(z);

y=ifftn(Y);

Gy_f=periodogram(y,window,'two sided',N);% 输出信号功率谱估计(周期图谱法)

Z=fftn(z);

S=fftn(s);

figure,subplot(411);

plot(s);title('原始的语音信号'); axis([0 14000 -1 1]);

subplot(412);

plot(z-s);title('叠加的噪声信号'); axis([0 14000 -1 1]);

subplot(413);

plot(z);title('观测的输入信号'); axis([0 14000 -1 1]);

subplot(414);

plot(y);title('时域的输出信号'); axis([0 14000 -1 1]);

figure,subplot(311);

plot(10*log10(Gs_f));title('原始信号功率谱');

axis([0 N/2 -120 0]);

xlabel('f/Hz'),ylabel('Gs(f)/dB');

subplot(312);

plot(10*log10(Gz_f));title('观测信号功率谱');

axis([0 N/2 -120 0]);

xlabel('f/Hz'),ylabel('Gz(f)/dB');

subplot(313);

plot(10*log10(Gy_f));title('输出信号功率谱');

axis([0 N/2 -120 0]);

xlabel('f/Hz'),ylabel('Gy(f)/dB');

figure,subplot(411);

plot(abs(S));title('原始信号频谱'); axis([0 14000 0 150]);

subplot(412);

plot(abs(Z-S));title('噪声信号频谱'); axis([0 14000 0 150]);

subplot(413);

plot(abs(Z));title('输入信号频谱'); axis([0 14000 0 150]);

subplot(414);

plot(abs(Y));title('输出信号频谱'); axis([0 14000 0 150]);

wavplay(s);% 播放原始信号

pause;% 暂停，按任意键继续

wavplay(z);% 播放观测信号

pause;% 暂停，按任意键继续

wavplay(y);% 播放输出信号

谁能给个可以玩魂银的NDS中文模拟器或者下载地址。。。

493银魂+模拟器+金手指使用方法 邮箱