复件 zhengshu.m

来自「上传一个自己编的计算系统功率谱的程序」· M 代码 · 共 58 行

function fangcheng
[t,y]=ode45(@weifen,[0:0.0025:50],[0.5 0.20 5]);


%figure,plot(t,y(:,1));
figure,plot(y(:,1),y(:,3));
%figure,plot(y(:,1),y(:,2));
%figure,plot(y(:,2),y(:,3));
%figure,plot(y(:,3),y(:,1));
%figure,plot3(y(:,1),y(:,2),y(:,3));
 
%Fs = 10;   t = 0:1/Fs:500;
% h = spectrum.welch;    
  %d=log(abs(y(:,1)))   ;     % Instantiate a welch object. 
 %psd(h,y(:,1),'Fs',Fs);
 %  psd(h,y(:,1)); 
   
   
  % b = firgr(27,[0 .4 .6 1],[1 1 0 0]);
%h = dfilt.dffir(b); % Create the filter object.
%h = dfilt.dffir(y(:,1));
%h.arithmetic = 'fixed'; % Quantize the filter to fixed-point.
%hpsd = noisepsd(h);
%plot(hpsd)



%x=y(:,2);
%N=length(x);
%fx=fft(x);
%df=10/N;
%n=0:N/2;
%f=n*df;
%figure,subplot(211);
%plot(t,x);grid;


%A=log(abs(fx(n+1)));
%A=fx(n+1);
%figure,plot(f,A); grid; %(abs(fx(n+1))

%b = firgr(27,[0 .4 .6 1],[1 1 0 0]);
b = firgr(27,y(:,2));
h = dfilt.dffir(b); % Create the filter object.
h.arithmetic = 'fixed'; % Quantize the filter to fixed-point.
hpsd = noisepsd(h);
plot(hpsd)


function dy=weifen(t,y)
dy=zeros(3,1);
a=10;b=40;c=2.5;k=1;h=4;

dy(1)=a*(y(2)-y(1));
dy(2)=b*y(1)-k*y(1)*y(3);
dy(3)=-c*y(3)+h*y(1)^2;