interp_correlation1.m

信号提取是信号处理的一个重要技术

% analysis of discrete implementation of generalized corss correlator
% this algorithm is from IEEE Transaction on acoustic speech and signal
% processing, no 3, june, 1981.
% this program is writted by Tang Hong, 
% Feb 1, 2005.
% from 0.001---1

clear
LL=510;         % length of signal
L=500;
Fs=250e3;   % the sample frequency in Hz
M=200;         % numbers of interpolation
%D=11.45e-6        % set the TDOA in microsecond


%========filter by a low-pass filter=======
[b,a]=butter(10,0.3);
r=randn(1,LL);
s=filter(b,a,r);
clear r;

count=0;
for loop=1/M:1/M:1;
    count=count+1;
    D=(11+loop)/1e6;
    true(count)=D;
%========interpolation==========
interp_s=interp(s,M);
t1=(1/Fs)*(1/M);
delay=round(D/t1);

snr=0.1;
%======two receive signals========
x1=s+snr*randn(1,LL);
   tempt=interp_s(delay+1:length(interp_s));
x2=tempt(1:M:length(tempt));
x2=x2+snr*randn(1,length(x2));
%[m,i]=max(xcorr(x1(1:L),x2(1:L)));
%i-L

%=============

fx1=fft(x1(1:L));
fx2=fft(x2(1:L));
h=1./sqrt(abs(fx1).*abs(fx2));
gr=ifft(fx1.*conj(fx2).*h);
gr=abs(gr);
[m,i]=max((gr));
%i

es_d=(gr(i+1)-gr(i-1))/(2*(-gr(i+1)+2*gr(i)-gr(i-1)));

est(count)=(i+es_d-1)/Fs;

count
end   % loop
figure,hold on,grid on
plot(true)
plot(est,'r')