%模拟风速的功率谱
function [svv]=svv(vv,N)
omegau=pi
omegak=0.002*pi
detaw=(omegau-omegak)/N;
%求相关系数
c=size(vv);
n=max(c(1),c(2))
x=-300:300
for m=1:c(1)
vvco(m,=xcorr(vv(m,1:600),300,'biased');
end
vvco1=zeros(1,c(2));
vvco1(1,=vvco(1,;
%求功率谱
for l=1:N
omegal(l)=omegak+(l-0.5)*detaw;
sv(l)=0;
for m=1:n
sv(l)=sv(l)+vvco1(1,m)*exp(-i*(m-1)*omegal(l));
end
svv(l)=abs(sv(l));
end
%绘曲线
subplot(211);
plot(x,vvco1);
xlabel('n');
ylabel('Rxx(t)');
title('相关函数');
grid on;
subplot(212);
plot(omegal,svv);
xlabel('w');
ylabel('Sxx(t)');
title('功率谱密度函数');
hold on
plotdaven(500);
grid on |