Matlab插值拟合（专题三）

 

一维插值

（1）分段线性插值

计算简单，处理速度快，适用于光滑性要求不高的插值问题

缺点是对于海量数据及非线性问题，处理的误差较大，所得到的曲线并不是光滑的。

废话不多说，直接上代码

下面展示

format bank
xx=1:1:17;
yy=rand(1,17)*10;
xxi=1:0.3:17;
y=interp1(xx,yy,xxi,'linear');
plot(xxi,y,'r',xx,yy,'g')
legend('线性插值','原始数值')

需要注意的是：咱们这次数据点较少，且原数据不是曲线，所以线性插值拟合的值基本重合。

为了更直观的展示，我们把插值的线变成圈圈来表示

plot(xxi,y,'ro',xx,yy,'g')
legend('线性插值','原始数值')

接下来我们演示函数曲线的插值效果。

x=-4:4;y=1./(x.^2+1);%原始数据产生的结点
xx=-4:0.1:3;yy=1./(xx.^2+1);%从节点中产生数据插值点
y1=interp1(x,y,xx,'linear');%线性插值
plot(xx,yy,xx,y1,'r')%作图对比
legend('原始数值','线性插值')

可以看出插值的结果并不是很理想，我们可以用三次样条插值对比下效果。

（2）三次样条插值

 x=-4:4;y=1./(x.^2+1);%原始数据产生的结点
xx=-4:0.1:3;yy=1./(xx.^2+1);%从节点中产生数据插值点
y1=interp1(x,y,xx,'spline');%三次样条插值
plot(xx,yy,xx,y1,'r')%作图对比
legend('原始数值','三次样条插值')

可以看出三次样条插值对于曲线来说效果更好，几乎重合。

（3）三次多项式插值

精确度更高更好，插值的曲线较平滑，可用作对于插值精确度要求较高的计算中

同时计算需要的内存也多，时间也较长。

x=-4:4;y=1./(x.^2+1);%原始数据产生的结点
xx=-4:0.1:3;yy=1./(xx.^2+1);%从节点中产生数据插值点
y1=interp1(x,y,xx,'cubic');%分段三次多项式插值
plot(xx,yy,xx,y1,'r')%作图对比
legend('原始数值','三次多项式插值')

可以看出，重合度比三次样条更高。

（4）最邻近插值

通常用于图像放缩中，在数值分析中较为少用

x=-4:4;y=1./(x.^2+1);%原始数据产生的结点
xx=-4:0.1:3;yy=1./(xx.^2+1);
y1=interp1(x,y,xx,'nearest');
plot(xx,yy,xx,y1,'r')%作图对比
legend('原始数值','最邻近区域插值')

二维插值

interp2

zi=interp2(x,y,z,xi,yi,method)

x,y,z为原始节点

xi,yi表示产生的数据插值点。zi表示插值后得到的近似值。

对于[xi,yi]与zi间的对应关系，则[x,y]与z的关系相同

常用的几种method如下：

‘linear’：双线性插值，（默认值）

‘cubic’：双三次插值

‘nearest’:最邻近插值

下面来个例子

下表是气象学家测量得到的气象资料，它们分别表示在南半球地区按不同维度、不同月份的平均气旋数字。根据这些数据，绘制出气旋分布曲面图形.

具体代码如下 

x=1:12;y=5:10:85;
z =[
2.4 18.7 20.8 22.1 37.3 48.2 25.6 5.3 0.3
1.6 21.4 18.5 20.1 28.8 36.6 24.2 5.3 0
2.4 16.2 18.2 20.5 27.8 35.5 25.5 5.4 0
3.2 9.2 16.6 25.1 37.2 40 24.6 4.9 0.3
1.0 2.8 12.9 29.2 40.3 37.6 21.1 4.9 0
0.5 1.7 10.1 32.6 41.7 35.4 22.2 7.1 0
0.4 1.4 8.3 33.0 46.2 35 20.2 5.3 0.1
0.2 2.4 11.2 31.0 39.9 34.7 21.2 7.3 0.2
0.5 5.8 12.5 28.6 35.9 35.7 22.6 7 0.3
0.8 9.2 21.1 32.0 40.3 39.5 28.5 8.6 0
2.4 10.3 23.9 28.1 38.2 40 25.3 6.3 0.1
3.6 16 25.5 25.6 43.4 41.9 24.3 6.6 0.3];
[x1,y1]=meshgrid(x,y);
mesh(x1,y1,z')%z转置是为了保持维度一致性。
xlabel('月份');ylabel('维度');zlabel('平均气旋数字')

 

接着上面的代码进行插值处理

xi=1:1:12;yi=5:1:85;
[xj,yj]=meshgrid(xi,yi);
zi=interp2(x,y,z',xj,yj);%第六个参数未赋值，默认用linear方法。
mesh(xi,yi,zi)

成功插值后的图形 

griddata

语法和interp2一样。咱们直接上题目。

在某海域测得一些点(x，y)处的水深，在矩形区域(75，200)×(-50，150)内画出海地曲面图形。数据下表所示。试着画出海底地貌图和等高线图。

x	129	140	103.5	88	185.5	195	105	157.5	107.5	77	81	162	162	117.5
y	7.5	141.5	23	147	22.5	137.5	85.5	-6.5	-81	3	56.5	-66.5	84	-33.5
z	4	8	6	8	6	8	8	9	9	8	8	9	4	9

x=[129 140 103.5 88 185.5 195 105 157.5 107.5 77 81 162 162 117.5];
y=[7.5 141.5 23 147 22.5 137.5 85.5 -6.5 -81 3 56.5 -66.5 84 -33.5]; 
z=-[4 8 6 8 6 8 8 9 9 8 8 9 4 9];
cx=75:5:200;%插值x
cy=-70:8:150;%插值y
[CX,CY]=meshgrid(cx,cy);
cz=griddata(x,y,z,CX,CY,'cubic');%为了使图像更加精确和光滑，我们采用cubic方法。
mesh(CX,CY,cz)

接着上面的代码，画出等高线图。 

contour3(CX,CY,cz,50)

 

 再例如在海底某个矩形内测得的海底深度z,如下，绘制海底深度图及等值线图。

z=[100,99,100,99,100,99,99,99,100;
100,99,99,99,100,99,100,99,99;
99,99,98,98,100,99,100,100,100;
100,98,97,97,99,100,100,100,99;
101,100,98,98,100,102,103,100,100;
102,103, 101, 100,102,106,104,101,100;
99,102,100,100,102,106,104,101,100;
97,99,100,100,102,105,103,101,100;
100,102,103,101,102,103,102,100,99;
100,102,103,102,101,101,100,99,99;
100,100,101,101,100,100,100,99,99;
100,100,100,100,100,99,99,99,99;
100 100,100,99,99,100,99,100,99];
x=0:8;
y=0:12;
mesh(x,y,z)
hold on
[xx,yy]=meshgrid(x,y);
plot3(xx,yy,z+0.1,'ro')

得到的原始图形（未插值。）

为了更直观深度，我们采用邻近插值法。

 

z=[100,99,100,99,100,99,99,99,100;
100,99,99,99,100,99,100,99,99;
99,99,98,98,100,99,100,100,100;
100,98,97,97,99,100,100,100,99;
101,100,98,98,100,102,103,100,100;
102,103, 101, 100,102,106,104,101,100;
99,102,100,100,102,106,104,101,100;
97,99,100,100,102,105,103,101,100;
100,102,103,101,102,103,102,100,99;
100,102,103,102,101,101,100,99,99;
100,100,101,101,100,100,100,99,99;
100,100,100,100,100,99,99,99,99;
100 100,100,99,99,100,99,100,99];
x=0:8;
y=0:12;
cx=0:0.4:8;
cy=0:0.4:12;
[xx,yy]=meshgrid(cx,cy);
[x1,y1]=meshgrid(x,y);
cz=interp2(x,y,z,xx,yy,'nearest');
mesh(xx,yy,cz)
hold on
plot3(x1,y1,(z+0.1)','ro')
hold off

 下面绘制海底深度等值线图。

z=[100,99,100,99,100,99,99,99,100;
100,99,99,99,100,99,100,99,99;
99,99,98,98,100,99,100,100,100;
100,98,97,97,99,100,100,100,99;
101,100,98,98,100,102,103,100,100;
102,103, 101, 100,102,106,104,101,100;
99,102,100,100,102,106,104,101,100;
97,99,100,100,102,105,103,101,100;
100,102,103,101,102,103,102,100,99;
100,102,103,102,101,101,100,99,99;
100,100,101,101,100,100,100,99,99;
100,100,100,100,100,99,99,99,99;
100 100,100,99,99,100,99,100,99];
x=0:8;
y=0:12;
cx=0:0.4:8;
cy=0:0.4:12;
[xx,yy]=meshgrid(cx,cy);
cz=interp2(x,y,z,xx,yy,'bicubic');
pcolor(xx,yy,cz)%做颜色图
shading interp%通过在每个线条或面中对颜色图索引或真彩色值进行插值来改变该线条或面中的颜色。
hold on
contour(xx,yy,cz,15,'k')
colormap(cool)%着色系
colorbar('vert')%s竖直颜色条
hold on