clear;close all;clc;
%x单调
x=[7;10.5;13;17.5;34;40.5;44.5;48;56;61;68.5;76.5;80.5;91;96;101;104;106.5;111.5;118;123.5;136.5;142;146;150;157;158];
y1=[44;45;47;50;50;38;30;30;34;36;34;41;45;46;43;37;33;28;32;65;55;54;52;50;66;66;68];
y2=[44;59;70;72;93;100;110;110;110;117;118;116;118;118;121;124;121;121;121;122;116;83;81;82;86;85;68];
x1=7:1e-1:158;
x=(10000000*x./4.5)/1000000;%km单位
x1=(10000000*x1./4.5)/1000000;%km单位
y1=(10000000*y1./4.5)/1000000;%km单位
y2=(10000000*y2./4.5)/1000000;%km单位
y11=spline(x,y1,x1);
y22=spline(x,y2,x1);
dx = x1(2: end) - x1(1: end-1);
dy1 = y11(2: end) - y11(1: end-1);
dy2 = y22(2: end) - y22(1: end-1);
l1 = sum(sqrt(dx.^2+dy1.^2));
l2 = sum(sqrt(dx.^2+dy2.^2));
figure('position', [100, 200, 800, 500]);
plot(x, y1, 'ro--');
hold on;
plot(x, y2, 'bo--');
hold on;
plot(x1,y11,'r');
hold on;
plot(x1,y22,'b');
legend({'原始数据y1', '原始数据y2', '插值数据y11', '插值数据y22',},...
'location', 'northwest');
grid on
title(sprintf('插值后曲线y11长度:%.3f\n插值后曲线y22长度:%.3f', l1, l2));