%吳逢吉,無線定位TOA/TDOA參估計方法的研究,仿真問題(非高斯噪聲白信號+不相關白噪聲)(13-17)
clear;clc;
t=-1023:1023;
sx=-1+2.*rand(1,1024);%產生均值0,方差1/1.732,[-1,1]均勻分布的隨機序列.
sx1=1000.*sx';
sy1=1.5.*circshift(sx1,1004);%延時20單位sy1=6.*circshift(sx1,-1004);%延時20單位
sy=sy1';
x1=awgn(sx,0,'measured');%加入SNR為零的正態分布白噪聲隨機序列,'measured'表測量信號噪聲.x1=awgn(sx,0,'measured','db' )
y1=awgn(sy,0,'measured' );%第二種改法y1=awgn(sy,0,'measured','db' );
rxx=xcorr(x1,x1,'biased');%求信號x1自相關函數,度應是1024
ryy=xcorr(y1,y1,'biased');%求信號y1自相關函數
rxy=xcorr(x1,y1,'biased');%求信號x1,y1互相關函數[rxx,lags]=xcorr(x1,x1,'biased')
sxx=fft(rxx,2047);%求x1功率譜密度
syy=fft(ryy,2047);%求y1功率譜密度
sxy=fft(rxy,2047);%求互相關譜密度,
cxy=(sxy.* sxy)/(sxx.*syy);
Hs=cxy./((abs(sxy)).*(1-cxy));
sxyf=Hs.*sxy;
TDOA=ifft(sxyf,2047);%注2047必與t的長度一致,否則會報錯.
GCCA=(TDOA.*conj(TDOA)).^0.5;
t1=-30:0.01:30;
final=interp1(t,GCCA,t1,'spline');%插值
plot(t1,final)
axis([-40,40,-0.1,0.2])
grid on