1-3.m

累加弦长的参数三次样条曲线方法求插值曲线

x=[8.125 8.4 9 9.485 9.6 9.959 10.166 10.2];                                %%给所有变量赋初值
y=[0.0774 0.099 0.28 0.6 0.708 1.2 1.8 2.177];
n=length(x);                                                                %%求得插值点数
s=zeros(1,n);
a=0.01087;                                                                  %%边界斜率
b=100;                                                                  

for i=2:n
    s(i)=s(i-1)+sqrt((x(i)-x(i-1))^2+(y(i)-y(i-1))^2);                      %%求得各插值点对应的弦长值
end

x0=1/sqrt(1+a^2);
xn=1/sqrt(1+b^2);
y0=a/sqrt(1+a^2);
yn=b/sqrt(1+b^2);                                                           %%参数坐标下的边界条件
pp1=csape(s,x,[1,1],[x0,xn]);                                                   
pp2=csape(s,y,[1,1],[y0,yn]);                                               %%求得插值函数
x_values=fnval(pp1,0:0.0001:3.5);                                          
y_values=fnval(pp2,0:0.0001:3.5);                                       

plot(x_values,y_values,'LineWidth',2)
hold on
plot(x,y,'m:o')                                                             %%比较曲线
xlabel('x','FontSize',16)
ylabel('y','FontSize',16)
title('用累加弦长法求得的三次参数样条曲线','FontSize',12,'FontName','宋体')

dx_1=fnval(fnder(pp1),s)                                                    %%各插值点一阶导数
dy_1=fnval(fnder(pp2),s)
dx_2=fnval(fnder(pp1,2),s)                                                  %%各插值点二阶导数
dy_2=fnval(fnder(pp2,2),s)
for i=1:n
    k(i)=((dy_2(i))*(dx_1(i))+(dy_1(i))*(dx_2(i)))/...
         (sqrt(((dx_1(i))^2+(dy_1(i))^2)^3));                                %%各插值点曲率
end
k