LERP: fast n-dimensional linear interpolation & extrapolation

Syntax:
F=LERP(X,Y,XI)
F=LERP(X,Y,XI,EXTRAP)

LERP is akin to the Matlab function INTERP1Q and can be used with the same syntax as INTERP1Q with the same results at comparable speed. However, LERP has the following features and capabilities to set it apart from INTERP1Q and others:
•Vectorization allows LERP to be much faster that built-in MATLAB interpolation functions for some applications that would otherwise require many looped function calls.
•XI may be an n-dim array
•Y may be a column vector or a matrix (as in INTERP1Q)
•If Y is a matrix and XI is an array with N dimensions, F will have N+1 dimensions. The first N dimensions correspond to size(XI), and the N+1 dimension is size(Y,2).
•F=LERP(X,Y,XI,EXTRAP) will expand the function domain (defined in coordinates X) with linear extrapolation by a distance EXTRAP in either direction. EXTRAP may also be a 2-element vector with the first element defining the extrapolation distance below the beginning of X and the second element defining the extrapolation distance beyond the final value of X.

LERP has similar limitations to INTERP1Q

Example: interpolation for log and log10 functions with plot comparing extrapolation to reality
X = (linspace(.01, 2, 7)').^2;
Y = [log(X) log10(X)];
XI = magic(7)*.2 - .5;
F = lerp(X, Y, XI, [0 4])
x = linspace(4,8);
plot(X,Y,x,log(x),'r:',x,log10(x),'r:',XI,F(:,:,1),'bo',XI,F(:,:,2),'ko')