clc
clear all
t=-1023:1023;
fs=80;
tt=(1:1024)/fs;
lent=length(tt);
T=lent/fs;
B=80;
k=B/T;
sx=sin(pi*k*tt.*tt);
% plot(real(sx))
sx1=1.*sx';
sy1=0.5.*circshift(sx1,1004);%延时20单位
sy=sy1';
for snr2=-20:2:0
n=0;
for i=1:1000
x1=awgn(sx,-5,'measured');
y1=awgn(sy,snr2,'measured');
rxx=xcorr(x1,x1,'biased');%求信号x1自相关函数
ryy=xcorr(y1,y1,'biased');%求信号y1自相关函数
rxy=xcorr(x1,y1,'biased');%求信号x1,y1互相关函数
sxx=fft(rxx,2047);%求x1功率谱密度
syy=fft(ryy,2047);%求y1功率谱密度
sxy=fft(rxy,2047);%求互相关谱密度,
Hs=1;
sxyf=Hs.*sxy;
TDOA=ifft(sxyf,2047);
CC=(TDOA.*conj(TDOA)).^0.5;
shiyan(i)=abs(find(max(CC)==CC)-1024);
if shiyan(i)==20
n=n+1;
end
end
a((snr2+20)/2+1)=n/1000
end
评论0