山区地貌问题

山区地貌问题：在某山区测得一些地点的高程如下表。（平面区域为：1200<=x<=4000,1200<=y<=3600)。试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。

解：x0=1200:400:4000;

y0=1200:400:3600;

z0=[1130,1250,1280,1230,1040,900,500,700;

1320,1450,1420,1400,1300,700,900,850;

1390,1500,1500,1400,900,1100,1060,950;

1500,1200,1100,1350,1450,1200,1150,1010;

1500,1200,1100,1550,1600,1550,1380,1070;

1500,1550,1600,1550,1600,1600,1600,1550;

1480,1500,1550,1510,1430,1300,1200,980];

meshc(x0,y0,z0)

xlabel('x');

ylabel('y');

zlabel('z')

title('原始图')

（1）分段线性插值

输入命令：

[xi2,yi2]=meshgrid(1200:10:4000,1200:10:3600); zi2=interp2(x0,y0,z0,xi2,yi2,'linear');

meshc(xi2,yi2,zi2)

xlabel('x');

ylabel('y');

zlabel('z')

title('插值后')

（2）三次样条插值

输入命令：

[xi,yi]=meshgrid(1200:10:4000,1200:10:3600);

zi=interp2(x0,y0,z0,xi,yi,'spline');

meshc(xi,yi,zi),

title('三次样条插值表面图')

xlabel('x');ylabel('y');zlabel('z')

（3）双三次插值

输入命令：

[xi,yi]=meshgrid(1200:10:4000,1200:10:3600); zi=interp2(x0,y0,z0,xi,yi,'cubic');

meshc(xi,yi,zi),

title('双三次插值后的表面图')

xlabel('x轴');

ylabel('y轴');

zlabel('z轴')

（4）画等高线图：

n=9;

pcolor(x0,y0,z0) shading interp

zmax=max(max(z0));zmin=min(min(z0)); caxis([zmin,zmax]) colorbar hold on

C=contour(x0,y0,z0,n,'k:'); clabel(C)

title('未拟合图形') hold off

1.49e+0031.38e+0031.27e+0

1.16e+940610

1.05

1.051.16e+500

600700800900

100011001200

13001400

1500

1600

pcolor(xi,yi,zi) shading interp

zmax1=max(max(zi));zmin1=min(min(zi)); caxis([zmin1,zmax1]) colorbar hold on

C=contour(xi,yi,zi,n,'k:'); clabel(C)

title('拟合图形')

hold off

1.53e+003

1.53e+003

1.53e+003 1.42e+003

1.42e+003

1.3e+003

1.3e+003

1.3e 1.18e+003

1.11.18e+0031.07e+003 1.0

1.07e+003954

838723

607

1.3e 1.18e+003 1.18e+0031.18e+003

6077238389541.07e+003 1.01.07e+003 1.11.3e+0031.3e+003 1.42e+003 1.42e+0031.53e+003 1.53e+0031.53e+003

相应的图形中都画出了地貌图和等高线，在以上的三种插值方法中，最为精确的要属双三次插值方法，可能在间隔取值较小时，它与三次样条插值所绘出的图形差异不大，不很明显。