txflt matlab,[转载]MATLAB以前版本里面的pchiptx.m函数的源代码

function

v = pchiptx(x,y,u)

%PCHIPTX Textbook piecewise cubic Hermite interpolation.

% v = pchiptx(x,y,u) finds the shape-preserving piecewise

% P(x), with P(x(j)) = y(j), and returns v(k) = P(u(k)).

% See PCHIP, SPLINETX.

% First

derivatives

h = diff(x);

delta = diff(y)./h;

d = pchipslopes(h,delta);

%

Piecewise polynomial coefficients

n = length(x);

c = (3*delta - 2*d(1:n-1) -

d(2:n))./h;

b = (d(1:n-1) - 2*delta +

d(2:n))./h.^2;

% Find

subinterval indices k so that x(k) <= u

< x(k+1)

k = ones(size(u));

for j = 2:n-1

k(x(j) <= u) = j;

end

%

Evaluate interpolant

s = u - x(k);

v = y(k) + s.*(d(k) + s.*(c(k) + s.*b(k)));

function d

= pchipslopes(h,delta)

% PCHIPSLOPES Slopes for shape-preserving Hermite cubic

% pchipslopes(h,delta) computes d(k) = P'(x(k)).

% Slopes at interior points

% delta = diff(y)./diff(x).

% d(k) = 0 if delta(k-1) and delta(k) have opposites signs or

either is zero.

% d(k) = weighted harmonic mean of delta(k-1) and delta(k) if they

have the same sign.

n = length(h)+1;

d = zeros(size(h));

k =

find(sign(delta(1:n-2)).*sign(delta(2:n-1))>0)+1;

w1 = 2*h(k)+h(k-1);

w2 = h(k)+2*h(k-1);

d(k) = (w1+w2)./(w1./delta(k-1) + w2./delta(k));

% Slopes at

endpoints

d(1) = pchipend(h(1),h(2),delta(1),delta(2));

d(n) = pchipend(h(n-1),h(n-2),delta(n-1),delta(n-2));

function d = pchipend(h1,h2,del1,del2)

%

Noncentered, shape-preserving, three-point formula.

d = ((2*h1+h2)*del1 - h1*del2)/(h1+h2);

if sign(d) ~= sign(del1)

d = 0;

elseif

(sign(del1)~=sign(del2))&(abs(d)>abs(3*del1))

d = 3*del1;

end