1. 程式人生 > >【 MATLAB 】序列相關與序列卷積之間的關係

【 MATLAB 】序列相關與序列卷積之間的關係

關於序列卷積,之前寫了3篇博文:

這篇博文介紹的是MATLAB本身自帶的函式,但這個函式conv有個不如意的地方,就是求過卷積之後我們不知道各個卷積值的位置。

然後我們後面擴充套件了下這個函式,命名為conv_m,這個函式在這個博文的最後給出。

還有一篇博文:

這篇博文只是給出了一個例項,程式裡面隱藏了求卷積值位置的指令碼程式。

序列卷積和公式:

y(n) = \sum_{k=-\infty }^{+ \infty }x(k)h(n-k)

而序列的互相關公式為:

r_{x,y}(l)= \sum_{n=-\infty }^{+ \infty }x(n)y(n-l)

如果x等於y,那麼就得到自相關函式的公式:

r_{x,x}(l)= \sum_{n=-\infty }^{+ \infty }x(n)x(n-l)

比較卷積和公式和互相關函式的公式,我們可以發現二者之間的關係:

r_{x,y}(l)= x(n)*y(-n)

有了這個關係,我們就可以使用卷積的函式來求兩個序列的互相關了。

首先給出擴充套件後的卷積函式conv_m的指令碼:

function [y,ny] = conv_m(x,nx,h,nh)
% Modified convolution routine for signal processing
%___________________________________________________
% [y,ny] = conv_m(x,nx,h,nh)
% [y,ny] = convolution result
% [x,nx] = first signal
% [h,nh] = second signal
%
nyb = nx(1) + nh(1);
nye = nx(length(x)) + nh(length(h));
ny = nyb:nye;
y = conv(x,h);

兩個訊號相加的函式sigadd:

function [y,n] = sigadd(x1,n1,x2,n2)
% implements y(n) = x1(n) + x2(n)
% [y,n] = sigadd(x1,n1,x2,n2)
%____________________________________
% y = sum sequence over n, which includes n1 and n2
% x1 = first sequence over n1
% x2 = second sequence over n2( n2 can be different from n1)
%
n = min( min(n1), min(n2) ):max( max(n1), max(n2) ); %duration of y(n)
y1 = zeros(1,length(n)); y2 = y1; %initialization
y1( find( ( n >= min(n1) )&( n <= max(n1) ) == 1  )  ) = x1; %x1 with duration of y1
y2( find( ( n >= min(n2) )&( n <= max(n2) ) == 1  )  ) = x2; %x2 with duration of y2
y = y1 + y2;

訊號移位函式:

function [y,n] = sigshift(x,m,k)
%implements y(n) = x(n - k)
%_________________________
%[y,n] = sigshift(x,m,k)
%
n = m+k;
y = x;

預備工作做好了,下面給出一個例子:

設:

x(n)=[3, 11, 7, \underset{\Lambda }{0},-1,4,2]

是原序列,設 y(n) 是原 x(n) 受到噪聲汙染並移位了的序列

y(n) = x(n-2) + w(n)

這裡 w(n) 是均值為0,方差為1的高斯隨機序列。計算y(n)和x(n)之間的互相關。

題解:

clc
clear
close all

% noise sequence 1
nx = -3:3;
x = [3,11,7,0,-1,4,2];
% 
%  implements y(n) = x(n - k)
%  _________________________
%  [y,n] = sigshift(x,m,k)
[y1,n1]= sigshift(x,nx,2);

w = randn(1,length(y1));
nw = n1;

[y,ny] = sigadd(y1,n1,w,nw);

[x,nx]=sigfold(x,nx);

[rxy,nrxy]=conv_m(x,nx,y,ny);

subplot(2,1,1);
stem(nrxy,rxy);
title('noise sequence 1');


% noise sequence 2
nx = -3:3;
x = [3,11,7,0,-1,4,2];
% 
%  implements y(n) = x(n - k)
%  _________________________
%  [y,n] = sigshift(x,m,k)
[y1,n1]= sigshift(x,nx,2);

w = randn(1,length(y1));
nw = n1;

[y,ny] = sigadd(y1,n1,w,nw);

[x,nx]=sigfold(x,nx);

[rxy,nrxy]=conv_m(x,nx,y,ny);

subplot(2,1,2);
stem(nrxy,rxy);
title('noise sequence 2');

噪聲是隨機的,所以把互相關的計算執行了兩次,可見,兩幅圖的細節有一點點不同,但互相關的峰值都在l = 2上。