home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 5.ddi / SPLINES.DI$ / SPAPI.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  1.5 KB  |  48 lines
  1. function sp=spapi(knots,x,y)
  2. % SPAPI    Spline interpolation.
  3. %
  4. %         spapi(knots,x,y)
  5. %
  6. %  returns in  sp  the spline  f  (if any) of order
  7. %               k := length(knots) - length(x)
  8. %  with knot sequence  knots  for which
  9. %                 y = f(x) .
  10. % This is taken in the osculatory sense in case some  x  are repeated, i.e.,
  11. % in the sense that  D^m(i) f(x(i))=y(i)  in case the  x  are in
  12. % nondecreasing order, with 
  13. %   m(i) := #{j<i: x(j) = x(i)}.
  14.  
  15. % C. de Boor / latest change: December 24, 1989
  16. % Copyright (c) 1990-92 by Carl de Boor and The MathWorks, Inc.
  17.  
  18. if (~isempty(find(diff(knots)<0))),
  19.    error('The knot sequence should be nondecreasing!')
  20. end
  21.  
  22. n=length(x);npk=length(knots);k=npk-n;
  23. if (k<1),
  24.    fprintf('number of data points, %.0f, should be less than ',n)
  25.    fprintf('the number, %.0f, of knots!\n',npk)
  26.    error('')
  27. end
  28.  
  29. % It is assumed that  y  has the same sense as  x  which is important in case
  30. % the entries of  y  are vectors rather than scalars
  31. [rx,cx]=size(x); [ry,cy]=size(y);
  32. if (rx>1),
  33.    if (ry~=rx), error(' y  should contain exactly as many entries as  x ')
  34.    else, y=y';end
  35. else,
  36.    if (cy~=cx), error(' y  should contain exactly as many entries as  x ')
  37.    end
  38. end
  39.  
  40. %  sort the given abscissae.
  41. [x,index]=sort(x);y=y(:,index);
  42.  
  43. %  Generate the collocation matrix and divide it into the possibly reordered
  44. %  sequence of given ordinates to generate the B-spline coefficients of the 
  45. %  interpolant, then put it all together into  sp .
  46.  
  47. sp=spmak(knots,slvblk(spcol(knots,k,x,1),y')');
  48.