home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 5.ddi / SPLINES.DI$ / NEWKNT.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  2.7 KB  |  66 lines
  1. function [newbrk,distfn]=newknt(pp,newl)
  2. % NEWKNT    Improve knot sequence distribution.
  3. %
  4. %        [newbrk,distfn]=newknt(pp,newl)
  5. %
  6. % returns the breakpoint sequence  newbrk  which cuts the basic interval
  7. % for  pp  into  newl  pieces in such a way that a certain p.linear monotone
  8. % function (whose pp-form is in  distfn ) related to the high derivative of  
  9. %  pp  is equidistributed.
  10. %
  11. %   The intent is to choose a break sequence suitable to the fine 
  12. % approximation of a function  f  whose rough approximation in  pp  is assumed
  13. % to contain enough information about  f  to make this feasible.
  14.  
  15. % C. de Boor / latest change: May 24, 1989
  16. % C. de Boor / latest change: Nov.14, 1990  (add missing final shift)
  17. % Copyright (c) 1990-92 by Carl de Boor and The MathWorks, Inc.
  18.  
  19. [breaks,coefs,l,k,d] = ppbrk(pp);
  20. if d>1, fprintf('newknt only functions for functions.\n'),end
  21. %  it would be feasible to have it function for curves as well by choosing
  22. %   distfn  to have as its derivative the maximum over all component
  23. %  functions.
  24.  
  25. if l==1,   % if there is only one piece, return uniform break sequence
  26.    newbrk=breaks(1) + [0:newl]*(breaks(l+1)-breaks(1))/newl;
  27.    return
  28. end
  29.  
  30. % The distribution function  distfn  is constructed as the integral of the
  31. % k-th root of the absolute value of the k-th derivative of  pp .
  32. %  Since  pp  is of degree < k , this requires some approximation. The ap-
  33. % proximation to the absolute value of the k-th derivative is found as the 
  34. % piecewise constant on the same break sequence whose value on  
  35. % breaks(i)..breaks(i+1)  is the derivative at  
  36. %   (break(i-1)+3break(i)+3break(i+1)+break(i+2))/8  
  37. % of the parabola which agrees with the variation of  the (k-1)st derivative 
  38. % of  pp  at the three points breaks(i-1/2), breaks(i+1/2), breaks(i+3/2) .
  39.  
  40. abs(diff(coefs(:,1))'./(breaks(3:l+1)-breaks(1:l-1)));
  41. (ans([1 1:l-1])+ans([1:l-1 l-1])).^(1/k);
  42. distfn = fnint(ppmak(breaks,ans(:)'));
  43.  
  44. % The total variation of  distfn  is its value at  break(l+1) , i.e,
  45. var = fnval(distfn,breaks(l+1));
  46. if var==0 ,
  47.    newbrk=breaks(1) + [0:newl]*(breaks(l+1)-breaks(1))/newl;
  48.    return
  49. end
  50.  
  51. % The break sequence is to be chosen so that  
  52. %        distfn(newbrk) = 
  53. steps = [0:newl]*(var/newl);
  54. newbrk=zeros(1,newl+1);newbrk([1 newl+1])=breaks([1 l+1]);
  55.  
  56. % For this, do the inverse interpolation:
  57. [junk,coefs]=ppbrk(distfn);
  58. flats=find(coefs(:,1)==0);
  59. if (~isempty(flats)),coefs(flats,1)=coefs(flats,1)-1;end
  60. breaks=breaks-breaks(1);
  61. invbrk=[coefs(:,2)' var];
  62. coefs=[ones(l,1)./coefs(:,1),breaks(1:l)'];
  63. if (~isempty(flats)),coefs(flats,1)=(breaks(flats)+breaks(flats+1))'/2;end
  64.  
  65. newbrk(2:newl)=newbrk(1) + fnval(ppmak(invbrk,coefs,1),steps(2:newl));
  66.