请教频域信号转化到时域
下面的程序是关于陈果《铁路轨道不平顺随机过程的数值模拟》论文中将信号从频域转化到时域的代码,但是出来的结果和他在论文中给出的差距挺大的,主要是开始和结尾的值太大了,请高手帮忙看一下问题出在哪里,另外这个程序也是这个论坛别人发过的,当时讨论也没得出什么结论,我感觉也没什么问题,谢谢。% 频域离散点采样
%nextpow2() 取最接近的较大2次幂
clear;
clc;
Av=;
Aa=;
Ws=;
Wc=;
k=0.25;
rand('state',0);%%把正态随机数发生器置0
%yl=1.524;yu=304.8;
yl=0.5;yu=50;
pr1={'V='};
in1=inputdlg(pr1,'线路参数模拟信息',1,{'100'});
V=str2num(in1{1});
pr2={'选择路面等级i=','输入不平顺种类(1-高低不平顺,2-方向不平顺,3-水平及轨距不平顺j='};
in2=inputdlg(pr2,'线路参数模拟信息',1,{'6','1'});
i=str2num(in2{1});
j=str2num(in2{2});
v=V/3.6;
switch j
case 1
s=@(f)(k*Av(i)*Wc(i)^2./((f.*2*pi/v).^2.*((f*2*pi/v).^2+Wc(i)^2)))/v;%%s=k*Av(i)*Wc(i)^2/(w^2*(w^2+Wc(i)^2));
case 2
s=@(f)(k*Aa(i)*Wc(i)^2./((f*2*pi/v).^2.*((f*2*pi/v).^2+Wc(i)^2)))/v;
case 3
s=@(f)(4*k*Av(i)*Wc(i)^2./((f*2*pi/v).^2+Ws(i)^2)./((f*2*pi/v).^2+Wc(i)^2))/v;
end
fu=v/yl;fl=v/yu;
co=1;%判断点的个数和Nr/2的关系
while(co>0)
Ts=inputdlg('模拟时间Ts=','模拟的时间长度',1,{'10'});
Ts=str2num(Ts{1});
dt=0.0001;
N=Ts/dt;
Nr=2^(nextpow2(N));
Xm=zeros(1,Nr/2+1);
df=1/(Nr*dt);
n0=ceil(fl/df);
nf=ceil(fu/df);
if nf<=Nr/2
S=feval(s,(n0+1:nf)*df);
Xm(n0+1:nf)=Nr*sqrt(S*df);
break;
% N0+1~N0+Nf数值非零,其余为零。
end
end
n=nf-n0;
q=2*pi*rand(1,Nr/2+1);
X=Xm.*q;
temp=fliplr(X);
X=;
x=ifft(X,Nr);
t=(1:Nr)*dt;
%%-->频域法模拟结果
figure(1);
subplot(2,1,1);
plot(t,x);
grid on;
Fs=1/dt;
fn=(1:Nr/2)/Nr*Fs;
Xnn=fft(x,Nr);
Xn=Xnn(1:Nr/2);
Sn=(abs(Xn)*2/Nr).^2;
figure(1);
subplot(2,1,2);
loglog(fn,Sn,'* g',fn,s(fn),'r');
legend('模拟值','解析值');
xlim();
ylim(); 程序估计是没有什么问题,不过取值是不是有所不同?
现在国内文章都这样,很多取值都不给全,让人无法再现
页:
[1]