MATLAB自相关函数介绍

Matlab提供了计算互相关和自相关的函数xcorr函数
1.使用方法

c = xcorr(x,y)
c = xcorr(x)
c = xcorr(x,y,‘option’)
c = xcorr(x,‘option’)
c = xcorr(x,y,maxlags)
c = xcorr(x,maxlags)
c =
xcorr(x,y,maxlags,‘option’)
c =
xcorr(x,maxlags,‘option’)
[c,lags] = xcorr(…)
其中option为：
"biased"为有偏的互相关函数估计；
"unbiased"为无偏的互相关函数估计；
"coeff"为0延时的正规化序列的自相关计算；
"none"为原始的互相关计算
2.特别的：
c=
xcorr(x,‘option’)特指以上某个选项的自相关估计。
c = xcorr(x,y,maxlags)
返回一个延迟范围在[-maxlags,maxlags]的互相关函数序列，输出c的程度为2maxlags+1.
c = xcorr(x,maxlags)
返回一个延迟范围在[-maxlags,maxlags]的自相关函数序列，输出c的程度为2maxlags+1.
c = xcorr(x,y,maxlags,‘option’)
同时指定maxlags和option的互相关计算.
c = xcorr(x,maxlags,‘option’)
同时指定maxlags和option的自相关计算.
此部分转自
https://blog.csdn.net/lvsehaiyang1993/article/details/82913443?utm_medium=distribute.pc_relevant.none-task-blog-BlogCommendFromMachineLearnPai2-4.add_param_isCf&depth_1-utm_source=distribute.pc_relevant.none-task-blog-BlogCommendFromMachineLearnPai2-4.add_param_isCf

MATLAB功能代码实现

t=0:0.1:100;
noise=rand(1,1001);
input=sin(pit);
output=sin(pit)+oise;%添加大噪声
plot(output);

可以看到噪声信号已经淹没了输入信号

[x y]=xcorr(input,output,‘coeff’);%计算互相关性，存至数组
plot(y,x);%绘制互相关运算结果

这里得到了输出信号和输入信号的自相关函数，方便起见，假设此系统输入输出相等
[X Y]=xcorr(input,input,‘coeff’);%计算输入信号自相关函数
plot(Y,X);%绘制自相关运算结果

可以看到这个是标准的自相关函数求解结果
[a b]=xcorr(input,noise,‘coeff’);%计算输入与噪声互相关的结果
plot(b,a);%绘制互相关运算结果

可以看到经过一个互相关器运算，得到的增益是很小的，系统可以判断出两个信号相关性很差

试着继续加大噪声增益
output=sin(pit)+4noise;%添加大噪声
plot(output);
[x y]=xcorr(input,output,‘coeff’);%计算互相关性，存至数组
plot(y,x);%绘制互相关运算结果

可以看见随着噪声比例的增加，输出经过互相关器之后的运算结果有了明显的衰减

总结

通过相关运算可以很好的甄别有用信号，或者是延时信号，如果发现互相关结果出现明显衰减（本例均采用归一化手段了），就要注意对有用信号提供一定增益了，因为噪声以白噪声居多，频率含量丰富，有用信号也许就难以通过相关或者FFT手段提取了。

作者功力有限，字句、原理描述难免出现纰漏，若各位看官发现还请予以指正！

一、相关函数应用场景

求下面信号的 " 自相关函数 " :

$$x(n)=\sin (2\pi fn)+N(n)$$

其中 $N(n)$ 为 高斯白噪声 ;

高斯白噪声 符合 正态分布 特性 , 其 均值为 $0$ , 方差为 $1$ , 其功率谱密度是白的 , 在所有的频率上 , 其功率都相同 ;

在上一篇博客 【数字信号处理】相关函数应用 ( 相关函数应用场景 | 噪声中检测信号原理 ) 中 , 使用了公式推导的方法求相关函数 , 本篇博客使用 matlab 求相关函数 ;

下面开始使用 matlab 计算 $x(n)=\sin (2\pi fn)+N(n)$ 的相关函数 ;

1、生成高斯白噪声

生成 高斯白噪声 序列 , 均值 0 , 方差 1 ;

randn(1,200); 

上述代码 生成的 高斯白噪声 序列 是一个 有 200 个元素的行向量 ,

如果是 randn(2, 200) 代码 , 则生成的是一个 $2\times 200$ 的矩阵 ;

2、信噪比 SNR

SNR 是信噪比 , $P_s$ 是信号功率 , $P_N$ 是噪声功率 , 其关系是 :

$$SNR=10{\log }_{10}\frac{P_s}{P_N}$$

在 matlab 中 , 设置信噪比为 $7$ ;

% 设置 信噪比 
SNR = 7;

3、根据信噪比 SNR 求信号幅度

信号 $$A\sin \omega n$$ , 其功率是 $$P_s=\frac{A^2}{2}$$ , 噪声功率 是 $$P_N=1$$ , 将其带入到信噪比公式中 :

$$SNR=10{\log }_{10}\frac{P_s}{P_N}=10{\log }_{10}\frac{\frac{A^2}{2}}{1}=10{\log }_{10}\frac{A^2}{2}$$

根据信噪比 计算 幅度 $A$ :

$$\frac{SNR}{10}=lo{g}_{10}\frac{A^2}{2}$$

$$\frac{A^2}{2}=10^{\frac{SNR}{10}}$$

$$A^2=2\times 10^{\frac{SNR}{10}}$$

$$A=\sqrt{2\times 10^{\frac{SNR}{10}}}$$

写成 matlab 代码为 :

% 根据 信噪比 计算信号幅度
SignalAmplitude = sqrt( 2 * 10^(SNR/10) );

4、产生单载波信号及最终信号

根据 $$A\sin \omega n$$ 公式 , 产生 200 个信号 , 其中 $A$ 之前已经求出 ;

% 根据 Asinωn 产生 200 个信号 
x1 = A * sin( pi * 0.165 * (0:199) ); 

然后将 信号 与 高斯白噪声叠加 :

% 信号 + 高斯白噪声
x = x1 + xn; 

5、求自相关函数及功率

求自相关函数 , 使用 xcorr 函数 ;

% 求 x 的 自相关函数 , 长度为2N-1
y = xcorr(x, x); 

求信号功率 :

% 功率 : 自相关函数 幅度 求平均
% 自相关函数 就是 200 个功率之和
y = y / 200; 

6、matlab 完整代码

matlab 完整代码 :

% 清除之前的变量或内存
clear; 

% 生成 高斯白噪声 序列 , 均值 0 , 方差 1
% 生成的 高斯白噪声 序列 是一个 有 200 个元素的行向量
xn = randn(1,200); 

% 设置 信噪比 
SNR = 7;

% 根据 信噪比 计算信号幅度
A = sqrt( 2 * 10^(SNR/10) );

% 根据 Asinωn 产生 200 个信号 
x1 = A * sin( pi * 0.165 * (0:199) ); 

% 信号 + 高斯白噪声
x = x1 + xn; 

% 求 x 的 自相关函数 , 长度为2N-1
y = xcorr(x, x); 

% 功率 : 自相关函数 幅度 求平均
% 自相关函数 就是 200 个功率之和
y = y / 200; 

%建立幕布
figure;
%绘制 "输出序列" 图像 , 点用上三角表示
plot(y);

% 打开网格
grid on;

执行结果 :