function [sun]=main(n)
fplot('1/(x+2)',[-1,1],'r');
x=ones(n+2,1);
for j=0:n+1
x(j+1)=cos(pi*(n+1-j)/(n+1));
end
first=ones(n+2,1);
f=1./(x+2); %原函数 注意:: 这里出现了原函数。而我是没有函数只有数据,要作渐近线 也是为了最终辨识出系统模型,即它的函数模型
last=first;
for j=2:n+2
last(j)=(-1)*last(j-1);
end
A=ones(n+2,n+2);
A(:,1)=first;
A(:,n+2)=last;
for j=2:n+1
for t=2:j
A(:,j)=x.*A(:,j);
end
end
e=(1e-15)*first; %精度控制条件
sun=A\f;
while (1)
at='';
for i=1:n
for t=1:i
if (t==1)
at=strcat('(',num2str(sun(t+1)),')');
elseif (t>1)
xt='t';j=2;
while j<t
xt=strcat('t*',xt);j=j+1;
end
at=strcat(num2str(t),'*(',num2str(sun(t+1)),')*',xt,'+',at);
end
end
end
%以下得到逼近函数
ap1=sun(1:n+1,[1]);
for i=1:n+1
ap(i)=ap1(n+2-i);
end
yt=strcat('-1/(t+2)^2=',at);
[y]=solve(yt,'t');
y=numeric(y);
%以下得到一组新的交错点组
for i=1:n+1
if y(i) < 1 & y(i)>-1
for j=2:n+2
if y(i)<x(j)&y(i)>x(j-1)
if (1/(x(j-1)+2)-polyval(ap,x(j-1)))*(1/(y(i)+2)-polyval(ap,y(i)))> 0
x(j-1)=y(i);
elseif (1/(x(j-1)+2)-polyval(ap,x(j-1)))*(1/(y(i)+2)-polyval(ap,y(i)))< 0
x(j)=y(i);
end
end
end
end
end
A=ones(n+2,n+2);
A(:,1)=first;
A(:,n+2)=last;
for j=2:n+1
for t=2:j
A(:,j)=x.*A(:,j);
end
end
f=1./(x+2);
sun1=A\f;
if(abs(sun1-sun)<e)
break;
end
sun=sun1;
end
hold on;
funcion=poly2sym(ap);
ezplot(funcion,[-1,1]);
num=num2str(n);
legend('原函数曲线',strcat(num,'次逼近函数曲线'));
title('最佳逼近比较示意图');
xlabel('x的取值');
ylabel('f(x)的取值');
grid on;
end
happy说的是这个吧,我的问题根他的类似,但是有本质的区别,他是已经有原函数(在程序里面我已经标注出来了),再作图,再画渐近线,可我是没有原函数,只有原始作图的数据(各个点的数据),由数据得到的图再作渐近线,我的问题就是这样的!!
依然谢谢happy教授哦!! 另外说明的是,这个东西应该对大家做系统识别很重要,方法上完全符合系统识别最基础的理论,如果有做这方面的(控制系统分析及系统识别),欢迎大家一起探讨!!我的qq : 22007614
[ 本帖最后由 plsdd 于 2006-11-23 16:16 编辑 ] |