【编程训练】山区地貌 - 插值

问题描述

山区地貌:在某山区测得一些地点的高程为一个矩阵。（平面区域\(1200<=x<=4000,1200<=y<=3600\)）

试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。

结果展示

原始数据

插值后

插值代码

注：以下采用的是两组线性插值来得到二维平面的插值数据。

对于相同的题型，仅需要修改下面x,y,z。 [MATLAB] main.m 主程序

clc;clear;
x = [1200	1600	2000	2400	2800	3200	3600	4000];
y = [3600	3200	2800	2400	2000	1600	1200]';
z = [1480	1500	1550	1510	1430	1300	1200	980
     1500	1550	1600	1550	1600	1600	1600	1550
     1500	1200	1100	1550	1600	1550	1380	1070
     1500	1200	1100	1350	1450	1200	1150	1010
     1390	1500	1500	1400	900	    1100	1060	950
     1320	1450	1420	1400	1300	700	    900	    850
     1130	1250	1280	1230	1040	900	    500	    700];

figure(2);
    [X,Y]=meshgrid(x,y);
    mesh(X,Y,z);
    title('原数据');

%% init
h = 1;
xco = [x(1):h:x(end)];
yco = [y(1):-h:y(end)]';
zco = [];zco2 = [];
%% spline
for i = 1:size(z,1)
    zco(i,:) = spline(x,z(i,:),xco);
end

for j = 1:size(zco,2)
    zco2(:,j)=spline(y,zco(:,j),yco);
end

figure(1);
    [X,Y]=meshgrid(xco,yco);
    Z=zco2;
    subplot(3,3,1); %三维图
        mesh(X,Y,Z);
    subtitle('插值后三维图');
    subplot(3,3,2); %三维等高线
        contour3(X,Y,Z,20);
    title('三次样条插值','Color','red');
    subtitle('插值后三维等高线');
    subplot(3,3,3); %平面等高线
        contour(X,Y,Z,20);
    subtitle('插值后平面等高线');

%% pchip
for i = 1:size(z,1)
    zco(i,:) = pchip(x,z(i,:),xco);
end

for j = 1:size(zco,2)
    zco2(:,j) = pchip(y,zco(:,j),yco);
end

figure(1);
    [X,Y]=meshgrid(xco,yco);
    Z=zco2;
    subplot(3,3,4); %三维图
        mesh(X,Y,Z);
    subtitle('插值后三维图');
    subplot(3,3,5); %三维等高线
        contour3(X,Y,Z,20);
    title('分段三次 Hermite 插值多项式','Color','red');
    subtitle('插值后三维等高线');
    subplot(3,3,6); %平面等高线
        contour(X,Y,Z,20);
    subtitle('插值后平面等高线');

%% makima
for i = 1:size(z,1)
    zco(i,:) = makima(x,z(i,:),xco);
end

for j = 1:size(zco,2)
    zco2(:,j) = makima(y,zco(:,j),yco);
end

figure(1);
    [X,Y]=meshgrid(xco,yco);
    Z=zco2;
    subplot(3,3,7); %三维图
        mesh(X,Y,Z);
    subtitle('插值后三维图');
    subplot(3,3,8); %三维等高线
        contour3(X,Y,Z,20);
    title('修正 Akima 分段三次 Hermite 插值','Color','red');
    subtitle('插值后三维等高线');
    subplot(3,3,9); %平面等高线
        contour(X,Y,Z,20);
    subtitle('插值后平面等高线');