程序.m

给定一个音频范围内的模拟强噪声随机信号（强噪声+周期信号）

m=80000;
n=44100;
x=0:1/n:79999/n;
Signal=wavrecord(m,n,'double');
Signal=Signal';

figure(1)
plot(x,Signal);
axis([0,0.015,-2,2]);
title('源信号');

F1=fftshift(fft(Signal));%频谱函数数值
figure(2);
g = floor((-length(F1)/2+0.5):(length(F1)/2-0.5));
d=g*n/length(F1);
plot(d,abs(F1));
axis([0,10000,0,100]);
title('信号频谱');

GNoise=wgn(1,m,-20);
figure(3)
plot(x,GNoise);
axis([0,0.015,-2,2]);
title('高斯噪声');

F2=fftshift(fft(GNoise));%频谱函数数值
figure(4);
g = floor((-length(F2)/2+0.5):(length(F2)/2-0.5));
d=g*44100/length(F2);
plot(d,abs(F2));
axis([0,10000,0,100]);
title('噪声频谱');


V=Signal+GNoise;
figure(5)
plot(x,V);
axis([0,0.015,-2,2]);
title('加噪后的信号');

F3=fftshift(fft(V));%频谱函数数值
figure(6);
g = floor((-length(F3)/2+0.5):(length(F3)/2-0.5));
d=g*n/length(F3);
plot(d,abs(F3));
axis([0,10000,0,100]);
title('信号+噪声 频谱');

disp('自相关函数');
lv=xcorr(V,V,'unbiased');
figure(7);
plot(lv);axis([1000,1500,-0.3,0.3])
title('自相关函数示意图像');

%figure(20)
%lv3=xcorr(Signal,Signal,'unbiased')
%plot(lv3);axis([1000,2000,-0.3,0.3])
%title('测试信号');

figure(18);
lv1=10*diff(lv);
plot(lv1);
axis([1000,1500,-0.2,0.2]);
title('恢复信号');


F4=fftshift(fft(lv));%频谱函数数值
figure(8);
g = floor((-length(F4)/2+0.5):(length(F4)/2-0.5));
d=g*n/length(F4);
plot(d,abs(F4));
axis([0,10000,0,100]);
title('自相关函数恢复信号频谱');

%各函数功率谱
window=hamming(500);
noverlap=20;
range='half';
fs=n;

[Py,f]=pwelch(Signal,window,noverlap,m,fs,range);
[Py1,f]=pwelch(GNoise,window,noverlap,m,fs,range);
[Py2,f]=pwelch(V,window,noverlap,m,fs,range);
[Py3,f]=pwelch(lv,window,noverlap,m,fs,range);

DBPy=10*log10(Py);
DBPy1=10*log10(Py1);
DBPy2=10*log10(Py2);
DBPy3=10*log10(Py3);

figure(9)
plot(f,DBPy);
xlabel('F(Hz)');
ylabel('DB');
title('原信号功率谱');

figure(10)
plot(f,DBPy1);
xlabel('F(Hz)');
ylabel('DB');
title('噪声信号功率谱');

figure(11)
plot(f,DBPy2);
xlabel('F(Hz)');
ylabel('DB');
title('信号+噪声功率谱');

figure(12);
plot(f,DBPy3);
xlabel('F(Hz)');
ylabel('DB');
title('恢复信号功率谱');


%信号分析

%信号均值
disp('信号均值');
Ex1=mean(lv) 
%信号方差值
disp('信号方差值');
D=var(lv)
%信号标准差
disp('信号标准差');
S=std(lv)
%信号均方值
disp('信号均方值');
Ex2=Ex1^2+D



%恢复信号分析

%信号均值
disp('信号均值');
Ex1=mean(lv) 
%信号方差值
disp('信号方差值');
D=var(lv)
%信号标准差
disp('信号标准差');
S=std(lv)
%信号均方值
disp('信号均方值');
Ex2=Ex1^2+D