Sampling Frequency offset impairment in MATLAB

SCO = 50e-6; % value of SCO
out = RX_SCO(out,SCO);
figure;plot(real(out),imag(out),'b+');title('constellation with Sampling Frequency Offset impairment');

function out = RX_SCO(in,SCO);
[M N] = size(in);
if(SCO ~= 0)
indx = (1:M)*(1+SCO);
if SCO > 0
in = [in; zeros(fix(indx(end))-M,N)]; % zero-padding: extra
% samples are needed
for kk=1:N
in2(:,kk) = lrange_interp(in(:,kk),indx,'cubic');

end
else
in = [zeros(1,N); in];
indx = indx+1;
for kk=1:N
in2(:,kk) = lrange_interp(in(:,kk),indx,'cubic');
end
end
out = in2;
else
out = in;
end

function yi = lrange_interp(y,xi,interptype)
% Check arguments
if(nargin ~= 3)
error('fast_interp expects 3 arguments');
end
if( (xi > length(y)) | (xi < 1) )
error('xi must be between 1 and the length of y');
end
if(~ischar(interptype))
error('Interpolation type must be a string');
end
% For linear interpolation, just use matlab's linear
% interpolation function

if(strcmp(interptype,'linear') == 1) yi = interp1(1:length(y),y,xi,'*linear');
% For cubic interpolation, calculate piecewise lagrange
% polynomial at specified points
elseif(strcmp(interptype,'cubic') == 1)
yi = cubic_lrange(y,xi);
% Otherwise print an error to the screen
else
error('interptype must be either linear or cubic');
end

function [coeffs] = lrange_coeff(order,mu)
% Calculates the coefficients of an interpolator filter with
% a given interpolation order and a time shift 'mu'.
% mu is centered around a symmetrical interval [-1,1]
mu = (mu-0.5)*2/order;
% coefficients are calculated following Lagrange interpolation
if (order == 1) % linear
coeffs = [0.5+0.5.*mu 0.5-0.5.*mu];
elseif (order == 2) % quadratic
coeffs = [ 0.5.*mu + 0.5.*mu.^2 ...
1 - mu.^2 ...
-0.5.*mu + 0.5.*mu.^2 ];
elseif (order == 3) % cubic
mu2 = mu.*mu; % mu squared
mu3 = mu2.*mu;% mu cubed
coeffs = [ -0.0625 - 0.0625.*mu + 0.5625.*mu2 + 0.5625.*mu3 ...
0.5625 + 1.6875.*mu - 0.5625.*mu2 - 1.6875.*mu3 ...
0.5625 - 1.6875.*mu - 0.5625.*mu2 + 1.6875.*mu3 ...
-0.0625 + 0.0625.*mu + 0.5625.*mu2 - 0.5625.*mu3 ];
end;

function yi = cubic_lrange(y,xi);
y = reshape(y,length(y),1); % make sure y is a column vector
y(end+1:end+2) = 0;
y = [0;y];
xi = xi + 1;
% Get fractional part of indices
mu = xi - floor(xi);
% Get integer part of indices
xi = floor(xi);
% Get values from y used for cubic interpolation
xi = reshape([y(xi+2) y(xi+1) y(xi) y(xi-1)],length(mu),4);
% Perform interpolation using piecewise Lagrange polynomials
yi = sum((xi.*lrange_coeff(3,mu.')).');