splintfe.m

来自「可以方便的调用这个函数」· M 代码 · 共 23 行

%进行固定斜率端点条件的分段三阶样条插值
function yhat=splintFE(x,y,xhat,fp1,fpn)
x=x(:);y=y(:);xhat=xhat(:);
n=length(x);
dx=diff(x);
divdif=diff(y)./dx;
alpha=[0;dx(1:n-2);0];
beta=[1;2*(dx(1:n-2)+dx(2:n-1));1];
gamma=[0;dx(2:n-1);0];
A=tridiags(n,beta,alpha,gamma);
delta=[fp1;3*(divdif(2:n-1).*dx(1:n-2)+divdif(1:n-2).dx(2:n-1));fpn];
mmdflag=spparms('autommd',0);
b=A\delta;
spparms('autommd',mmdflag);
a=y(1:n-1);
c=(3*divdif-2*b(1:n-1)-b(2:n))./dx;
d=(b(1:n-1)-2*divdif+b(2:n))./dx.^2;
i=zeros(size(xhat));
for m=1:length(xhat)
    i(m)=binSearch(x,xhat(m));
end
xx=xhat-x(i);
yhat=a(i)+xx.*(b(i)+xx.*(c(i)+xx.*d(i)));