白噪声和高斯白噪声

1.白噪声

什么是白噪声

对于时间序列{wt} 若满足下面三个条件,该序列为一个离散的白噪声(white noise):

  1. 每个时间点均值为0:E(wt)=0
  2. 每个时间点方差为σ2Var(wt)=σ2
  3. 对于任意k≥1,自相关系ρk=0Cov(wt,wt+1)=0

使用matlab产生白噪声

fs=10000;                       %采样频率
Ns=10000;                       %采样点数
t=0:1/fs:Ns/fs;
s3=rand(1, Ns+1);%生成一个白噪声数组

1‑1白噪声信号

nfft=length(s3);
cxn=xcorr(s3,'unbiased');                  % 计算序列的自相关函数
CXk=fft(cxn,nfft);
Pxx=abs(CXk);
index=0:round(nfft/2-1);
k=index*fs/nfft;
%plot_Pxx=10*log10(Pxx(index+1)); 
plot_Pxx=(Pxx(index+1)); 

 先进行自相关求解,然后进行fft,就能得到功率谱。

 1‑2白噪声信号功率谱

histogram(s3,25) %直方图

1‑3白噪声信号直方图

 2.高斯白噪声

什么是高斯噪声

高斯白噪声(white Gaussian noise; WGN):均匀分布于给定频带上的高斯噪声

  1. 如果一个噪声,它的幅度服从高斯分布,而它的功率谱密度又是均匀分布的,则称它为高斯白噪声。
  2. 高斯白噪声中的高斯是指:概率分布是正态函数,而白噪声是指:它的二阶矩不相关,一阶矩为常数,是指先后信号在时间上的相关性。这是考察一个信号的两个不同方面的问题。
  3. 热噪声和散粒噪声是高斯白噪声。

 使用matlab产生白噪声

s3=randn(1, Ns+1);%生成一个高斯白噪声数组

方法与白噪声同理,产生高斯白噪声使用randn函数。

图 2‑1高斯白噪声信号

 图 2‑2高斯白噪声信号功率谱

 高斯白噪声直方图

图 2‑3高斯白噪声信号直方图

 matlab滤波器使用

 

3‑1正弦信号和噪声信号

图 1‑4 红色信号为频率为5Hz的正弦波信号幅值为10,绿色正弦波信号幅值为1,频率为1KHz的信号,黑色信号为白噪声信号。

3‑2信号进行相加得到的结果

图 1‑5 为上面三个信号相加得到的结果,1KHz正弦波和白噪声为噪声信号,下面使用FIR滤波器将噪声滤掉。

3‑3使用FIR滤波器设计

图 1‑6 使用了matlab中的滤波器设计模块进行设计,设计的采样频率为10K,截至频率为10HZ,阶数为100,得到了上图的结果,最后将它生成为代码,进行程序的编写。

 

 3‑4滤波前傅里叶变换

图 3-4可以看出,改信号中包含了10Hz、1KHz以及白噪声的幅值谱。

图 3‑5滤波后信号

图 3-5 可以看出通过滤波器后得到了一个较好的正弦波信号。

 

 

3‑6滤波后傅里叶变换

图 3-6 可以明显的看出幅值谱已经没有了频率为1K的成分,已经为滤波器滤掉。

clear;
close all;
clc;
f1=10;f2=300;     %原始信号频率
fs=10000;                       %采样频率
Ns=10000;                       %采样点数
t=0:1/fs:Ns/fs;
s1=10*sin(2*pi*f1*t);
s2=sin(2*pi*f2*t);
s3=randn(1, Ns+1);%生成一个白噪声数组
nfft=length(s3);
cxn=xcorr(s3,'unbiased');                  % 计算序列的自相关函数
CXk=fft(cxn,nfft);
Pxx=abs(CXk);
index=0:round(nfft/2-1);
k=index*fs/nfft;
%plot_Pxx=10*log10(Pxx(index+1)); 
plot_Pxx=(Pxx(index+1)); 
figure(1);
subplot(211)
plot(t,s1,'r',t,s2,'g',t,s3,'k');
xlim([0 5000/fs])
title('正弦信号');xlabel('t/s');
y=s1+s2+s3;
set(gcf,'color','w');
subplot(212)
plot(t,y);
title('原始信号');xlabel('t/s');
set(gcf,'color','w');
N=4096;
f=(0:N-1)*fs/N;
yy=abs(fft(y,N)*2/N);
figure(2);
subplot(311)
plot(f(1:N/2),yy(1:N/2));
title('原始信号的FFT');xlabel('f/Hz');
xlim([0 1500]);
h=myfilter;
ylp=filter(h,y);
subplot(312)
plot(t,ylp);
xlim([0.2 0.8]);
title('低通滤波之后的信号')
xlabel('t/s');
yyy=abs(fft(ylp,N)*2/N);
subplot(313)
plot(f(1:N/2),yyy(1:N/2));
xlabel('f/Hz');
title('低通滤波之后的信号频谱');
xlim([0 1500]);
figure(3);
subplot(311)
plot(t,s3);
xlim([0 500/fs])
title('噪声信号');xlabel('t/s');
subplot(312);plot(k,plot_Pxx);title('间接法原始噪声信号的功率谱');     %绘制分贝形式的功率谱
subplot(313);histogram(s3,25);title('噪声直方图');     %直方图
function Hd = myfilter
%MYFILTER 返回离散时间滤波器对象。

% MATLAB Code
% Generated by MATLAB(R) 9.8 and DSP System Toolbox 9.10.
% Generated on: 17-Sep-2022 20:01:49

% FIR Window Lowpass filter designed using the FIR1 function.

% All frequency values are in Hz.
Fs = 10000;  % Sampling Frequency

N             = 100;      % Order
Fc            = 10;       % Cutoff Frequency
flag          = 'scale';  % Sampling Flag
SidelobeAtten = 100;      % Window Parameter

% Create the window vector for the design algorithm.
win = chebwin(N+1, SidelobeAtten);

% Calculate the coefficients using the FIR1 function.
b  = fir1(N, Fc/(Fs/2), 'low', win, flag);
Hd = dfilt.dffir(b);

% [EOF]

  • 8
    点赞
  • 98
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值