非均匀B样条插值——怎样控制插值曲线的长度

hht123

如题，以下是我的程序
%------------------非均匀B样条拟合MATLAB程序-----------------
clear
clc
close all
k=3;                        %  k是k次B样条的意思，是个参数，得设置。如三次B样条，k=3。
T=1;
fs=1000;
t=1/fs:1/fs:T;
N=length(t);
Sig1=2*sin(60*pi*t+0.25*pi);
Sig2=2*sin(120*pi*t+0.25*pi);
xx=Sig1;
[num,T]=extrema_1(xx);

x=[T',xx(T)'];
[n,m]=size(x);
%% -----------弦长参数化--------------------------------------
u(k+n)=0;
for i=1:n-1
u(k+i+1)=u(k+i)+sqrt((x(i+1,1)-x(i,1))^2+(x(i+1,2)-x(i,2))^2);
end;
L=u(n+k);
for i=1:n
u(k+i)=u(k+i)/L;
end;
for i=1:3
u(k+i+n)=1;
end
%控制多边线
plot(x(:,1),x(:,2),'o');
hold on
%------------反求n+2个控制点--------------------
%首位重节点v1=v2
%首位与控制多边形相切
A=zeros(n+2);
A(1,1)=1;A(1,2)=-1;
A(2,2)=1;
A(n+2,n+1)=-1;A(n+2,n+2)=1;
A(n+1,n+1)=1;
for i=3:n
for j=0:2
   
A(i,i+j-1)=Base(i+j-1,k,u,u(i+2));
  
end
end
%e:方程右边.
e=0;
for i=1:m
   
e(n+2,i)=0;
end
for i=1:n
   
e(i+1,:)=x(i,:);
end
%求出控制点d
d=inv(A)*e;
plot(d(:,1),d(:,2),'g');

hold on
%------------插值并作出样条曲线-----------------
x=0;y=0;down=0;
for j=1:(n-1)
   
uu=u(j+3):1/N:u(j+4);
   
for kk=1:length(uu)
      
down=down+1;
x(down)=d(j,1)*Base(j,3,u,uu(kk))+d(j+1,1)*Base(j+1,3,u,uu(kk))+d(j+2,1)*Base(j+2,3,u,uu(kk))+d(j+3,1)*Base(j+3,3,u,uu(kk));
y(down)=d(j,2)*Base(j,3,u,uu(kk))+d(j+1,2)*Base(j+1,3,u,uu(kk))+d(j+2,2)*Base(j+2,3,u,uu(kk))+d(j+3,2)*Base(j+3,3,u,uu(kk));
   
end
end
axis('equal');
plot(xx);hold on
plot(x,y,'red');
xlabel('x');ylabel('y');
grid on

运行结果为y的长度为1019，信号原始长度为1000，请问怎样将插值曲线的长度控制为1000？


调用的函数为：
%% ************K次非均匀B样条基函数************
function result = Base(i,k,u,t)
%
%see also http://www.matlabsky.com
%Bbase()参见http://www.matlabsky.com/thread-665-1-1.html
%
%第i段k次B样条基,Deboor递推递归算法
%t为变量,u(i)<=t<u(i+1),k=0时result=1;
if(k==0)
    if(u(i)<=t && t<u(i+1))           %注意1=u(i)<=t<u(i+1)=1时的情况,这里要用t<=u(i+1);
        result=1;
        return;
    else
        result=0;
        return;
    end
else
    if(u(i+k)-u(i)==0)
        alpha=0;
    else
        alpha=(t-u(i))/(u(i+k)-u(i));
    end
    if(u(i+k+1)-u(i+1)==0)
         beta=0;
    else
        beta=(u(i+k+1)-t)/(u(i+k+1)-u(i+1));
    end
end
result=alpha*Base(i,k-1,u,t)+beta*Base(i+1,k-1,u,t);