To get autocorrelation function of periodic signal with xcorr function

17 次查看(过去 30 天)
periodic_signal.JPG
I couldn't get autocorrelation function of my periodic signal with xcorr() function as I have expected. My signal is periodic and its autocorrelation function should be periodic with same amplitudes, I got transient signal when I used xcorr() function. Where is the mistake?
autocorrelation_funct.JPG
clc,clear,close all
T=2;A=2; %period and amplitude of signal
t=linspace(-4*T,4*T,1000);
x=mod(t,T)
x=A*(-x+A/2)%Periodic signal which it's autcorr. will be determined
plot(t,x,'linewidth',2);
ylim([-4 4]);
grid on
xlabel('t');ylabel('x(t)');title('Real Function, T=2, A=2');
%return
a=xcorr(x,((length(t))/2));
dt=median(diff(t)); % diff() calculates distance between adjacent elements in vector, median() avaraging
a=a*dt; %xcorr function calculate corelation as a discrete time signal, we converted to continous time signal
t(1001)=t(1000)+dt;
plot(t,a,'r');hold on

回答(1 个)

Luis J Gilarranz
Luis J Gilarranz 2019-5-26
No idea about xcor function, but you can write a lower level code that does the job:
Autocorrelation is just the pearson correlation coefficient between a time series and the same time series but with a lag. In this code I do the lag-1 autocorrelation inside a moving window of 100 datapoints.
clc,clear,close all
T=2;A=2; %period and amplitude of signal
t=linspace(-4*T,4*T,1000);
x=mod(t,T)
x=A*(-x+A/2)%Periodic signal which it's autcorr. will be determined
%I transpose x so that it works with the code I had
x=x'
HalfMW=50;%half of the moving window
cont=0;
AutCorVsTime=zeros(cont,1);
for d=HalfMW+1:1000-HalfMW
cont=cont+1;
Chunk = x(d-HalfMW:d+HalfMW,1);
SizeChunk = size(Chunk,1);
AutCorVsTime(d,1)=corr(Chunk(2:SizeChunk,1),Chunk(1:SizeChunk-1,1),'rows','complete'); %Pearson correlation coefficient ignoring NaN
end
plot(AutCorVsTime)
And as you can see the resulting figure is indeed periodic.
untitled.png
If you change the size of the moving window the plot obviously changes.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by