time   calls  line
                 92 function d = pchipslopes(x,y,del)
                 93 %PCHIPSLOPES  Derivative values for shape-preserving Piecewise Cubic Hermite
                 94 % Interpolation.
                 95 % d = pchipslopes(x,y,del) computes the first derivatives, d(k) = P'(x(k)).
                 96 
                 97 %  Special case n=2, use linear interpolation.
                 98 
            18   99    n = length(x); 
            18  100    if n==2   
                101       d = repmat(del(1),size(y));
                102       return
                103    end
                104 
                105 %  Slopes at interior points.
                106 %  d(k) = weighted average of del(k-1) and del(k) when they have the same sign.
                107 %  d(k) = 0 when del(k-1) and del(k) have opposites signs or either is zero.
                108 
            18  109    d = zeros(size(y)); 
                110   
  0.01      18  111    k = find(sign(del(1:n-2)).*sign(del(2:n-1)) > 0); 
                112   
            18  113    h = diff(x); 
            18  114    hs = h(k)+h(k+1); 
  0.04      18  115    w1 = (h(k)+hs)./(3*hs); 
            18  116    w2 = (hs+h(k+1))./(3*hs); 
            18  117    dmax = max(abs(del(k)), abs(del(k+1))); 
  0.01      18  118    dmin = min(abs(del(k)), abs(del(k+1))); 
  0.03      18  119    d(k+1) = dmin./conj(w1.*(del(k)./dmax) + w2.*(del(k+1)./dmax)); 
                120 
                121 %  Slopes at end points.
                122 %  Set d(1) and d(n) via non-centered, shape-preserving three-point formulae.
                123 
            18  124    d(1) = ((2*h(1)+h(2))*del(1) - h(1)*del(2))/(h(1)+h(2)); 
            18  125    if sign(d(1)) ~= sign(del(1)) 
                126       d(1) = 0;
            18  127    elseif (sign(del(1)) ~= sign(del(2))) && (abs(d(1)) > abs(3*del(1))) 
             1  128       d(1) = 3*del(1); 
             1  129    end 
            18  130    d(n) = ((2*h(n-1)+h(n-2))*del(n-1) - h(n-1)*del(n-2))/(h(n-1)+h(n-2)); 
            18  131    if sign(d(n)) ~= sign(del(n-1)) 
             2  132       d(n) = 0; 
            16  133    elseif (sign(del(n-1)) ~= sign(del(n-2))) && (abs(d(n)) > abs(3*del(n-1))) 
             1  134       d(n) = 3*del(n-1); 
             1  135    end