tinterp1

function tinterp1

% Script illustrating an undesired behaviour in interp1
%   yi = interp1 (x, y, xi)
% One expects that if xi coincides with an x value, the output would be
% EXACTLY the element of y corresponding to the x value.

% This does not happen in the example below.

% Scenario
%   Piecewise linear mapping from x to y, some processing, then inverse
%   mapping. The processing leaves the end points unchanged.
%     % forward mapping with xi(end) == x(end)
%     yi = interp1(x, y, xi)
%     % modify yi (but leave end points unchanged)
%     yim = ...
%     % inverse mapping
%     xim = interp1(y, x, yim)
% This strategy fails since at the end point yi(end) ~= y(end). When
% I do the inverse mapping, yim(end) is outside the bounds and xim(end)
% returns NaN.

% In the code for interp1, the interpolation is carried out as
%  yi = yl + (yr-yl)*u    % where 0 <= u <= 1, and yl is the point to
%                         % the left and yr is the point to the right
% For u = 1, this becomes
%  yi = yl + (yr-yl)
% In the case below this does NOT equal yr as it should.

% There are other approaches to linear interpolation which do
% guarantee matching at the y values
%  yi = (1-u)*yl + u*yr

x = [0 0.0275 0.08];
y2 = 0.016603773584906 - 3.400058012914542e-016;
y = [0 y2 0.08];

xi = x;

% Forward mapping
yi = interp1(x, y, xi)
yi - y

% Inverse mapping, returns an NaN
xi = interp1(y, x, yi)
xi - x

% ----------------------------
% Try interp1q - same problem
yi = interp1q(x', y', xi')
yi - y'

return