中插值与拟合

1、计算可视化1 插值与数据拟合1.1 一维数据的插值问题1.1.1一维插值问题的求解求解方法一维插值interp1()函数的调用格式为：说明函数根据x,y的值，计算函数在x1处的值。x,y是两个等长的已知向量，分别描述采样点和样本值，x1是一个向量或标量，描述欲插值的点，y1是一个与x1等长的插值结果。 插值方法linear，nearest，cubic，spline 注：X1的取值范围不能超出X的给定范围，否则，会给出“NaN”错误。例1-1：已知的数据点来自函数，根据生成的数据进行插值处理，得出较平滑的曲线直接生成数据。调用函数：Step1：x=0:.12:1;y=(x.2-3*x+5).*e

2、xp(-5*x).*sin(x);plot(x,y,x,y,'o')% 绘制原函数的图形%'linear'x1=0:.02:1; y0=(x1.2-3*x1+5).*exp(-5*x1).*sin(x1);y1=interp1(x,y,x1); %默认为'linear'，在样本点上斜率变化很大plot(x1, y1',':',x,y,'o',x1,y0)%'cubic'y2=interp1(x,y,x1, 'cubic'); %最占内存，输出结果与spline差不多plot(

3、x1,y2',':',x,y,'o',x1,y0)%'spline'y3=interp1(x,y,x1,'spline'); %最花时间，但输出结果也最平滑plot(x1,y3',':',x,y,'o',x1,y0)% 'nearest'y4=interp1(x,y,x1, 'nearest');%执行速度最快，输出结果为直角转折plot(x1,y4',':',x,y,'o',x1,y0)%'linear

4、', 'cubic', 'spline', 'nearest'plot(x1,y1',y2',y3',y4',':',x,y,'o',x1,y0);Step2：max(abs(y0(1:49)-y2(1:49),max(abs(y0-y3),max(abs(y0-y4)ans = 0.0177 0.0086 0.1598例1-2：编写一段程序，允许利用插值方法手工绘制一条光滑的曲线。function sketcher(vis) x=; y=; i=1; h=; axis(0,

5、1 0 1)while 1 x0,y0,but=ginput(1);if but=1, x=x,x0; y=y,y0; h(i)=line(x0,y0); set(h(i),'Marker','o'); i=i+1;else, break; endendif nargin=1, delete(h); end% 若需要，可以删除样本点标识 xx=x(1):(x(end)-x(1)/100: x(end); yy=interp1(x,y,xx,'spline'); line(xx,yy)1.1.2Lagrange插值算法及应用求解方法：已知xi，yi

6、点，可求出x向量上各点处的插值为：function y=lagrange(x0,y0,x)i1=1:length(x0);y=zeros(size(x);for i=ii ij=find(ii=i); y1=1; for j=1:length(ij),y1=y1.*(x-x0(ij(j); end y=y+y1*y0(i)/prod(x0(i)-x0(ij);end例1-3：对进行Lagrange插值。x0=-1+2*0:10/10; y0=1./(1+25*x0.2);x=-1:.01:1; y=lagrange(x0,y0,x); % Lagrange 插值ya=1./(1+25*x.2)

7、; plot(x,ya,x,y,':')y1=interp1(x0,y0,x,'cubic'); y2=interp1(x0,y0,x,'spline');plot(x,ya,x,y1,':',x,y2,'-')1.2已知样本点的定积分计算求解方法编写函数：function y=quadspln(x0,y0,a,b)f=inline('interp1(x0,y0,x,"spline")',.'x','x0','y0');y=quad

8、l(f,a,b,1e-8,x0,y0);%先用样条函数进行插值后再积分其中为样本点构成的横纵坐标向量，为积分区间例1-4：利用样条插值算法求解。求解：x0=0:pi/30:pi; y0=sin(x0);I=quadspln(x0,y0,0,pi)结果：I= 2.0000比较梯形法和插值法：梯形法：x0=0:pi/10:pi;y0=sin(x0);I1=trapz(x0,y0)结果：I1= 1.9835插值法：I=quadspln(x0,y0,0,pi)结果：I= 2.0000已知5个不均匀分布的样本点：x0=0,0.4,1 2,pi; y0=sin(x0); % 生成样本点plot(x0,y0

9、,x0,y0,'o') % 绘制出的样本点折线I=quadspln(x0,y0,0,pi)I=2.0191，大约有1%的相对误差，应该说是相当精确I1=trapz(x0,y0)I1=1.8416，用trapz() 函数将得出很大的误差(7.9%)样条插值的结果与理论之间的比较:x=0:0.01:pi,pi;y0a=sin(x); y=interp1(x0,y0,x,'spline');plot(x0,y0,x,y,':',x,y0a,x0,y0,'o')例1-5：已知其中的150个数据点，用quadspln()计算出该定积分的值。

10、求数值解：x=0:3*pi/2/200:3*pi/2; y=cos(15*x); I=quadspln(x,y,0,3*pi/2)结果：I= 0.0667绘制曲线：x0=0:3*pi/2/1000:3*pi/2; y0=cos(15*x0);y1=interp1(x,y,x0,'spline'); plot(x,y,x0,y1,':')1.3二维网格数据的插值问题求解方法二维插值interp2()函数的调用格式为：说明x0,y0,z0为已知数据，x1,y1为插值点构成的新的网格参数，z1矩阵为在所选插值网格点出的函数的近似值 插值方法linear，nearest

11、，cubic，spline 注：X1,Y1的取值范围不能超出X,Y的给定范围，否则，会给出“NaN”错误。例1-6：由可计算出一些较稀疏的网格数据，对整个函数曲面进行各种插值拟合，并比较拟合结果。绘制已知数据的网格图:x,y=meshgrid(-3:.6:3,-2:.4:2); z=(x.2-2*x).*exp(-x.2-y.2-x.*y);surf(x,y,z), axis(-3,3,-2,2,-0.7,1.5)默认插值算法进行插值:x1,y1=meshgrid(-3:.2:3, -2:.2:2);z1=interp2(x,y,z,x1,y1); surf(x1,y1,z1), axis(-3,3,-2,2,-0.7,1.5)立方插值z1=interp2(x,y,z,x1,y1,'cubic'); surf(x1,y1,z1), axis(-3,3,-2,2,-0.7,1.5)样条插值：z2=interp2(x,y,z,x1,y1,'spline');surf(x1,y1,z2), axis(-3,3,-2,2,