%%
% by Ljiljana Milic
% Supplemental material for Chapter VI:
%% 6. Sampling Rate Conversion by a Fractional Factor
%%
% Table of Contents
% 6. Sampling Rate Conversion by a Fractional Factor
%     6.1 Sampling Rate Conversion by a Rational Factor
%         Example 6.1
%     6.2 Efficient Implementation of Polynomial Interpolation Filters Using Farrow Structure
%         Example 6.2
%     6.3 Fractional Delay Filters
%         Example 6.3
%         Example 6.4
%% 6.1 Sampling Rate Conversion by a Rational Factor
%% Example 6.1
%
% In this example, the MATLAB function resample is used for sampling-rate conversion by a rational factor.
%
% Sampling-rate conversion by L/M = 3/2

close all, clear all
Fx = 20;   Tx=1/Fx; % Original sampling frequency in Hz
tx = 0:Tx:1; % Time vector tx
x = 0.9*sin(2*pi*tx); % Input sequence
y = resample(x,3,2); % Re-sampling
ty = (0:(length(y)-1))*2*Tx/3; % New time vector ty
figure 
subplot(2,1,1), stem(tx,x,'*k')
hold on
stem(ty,y,'-.k')
legend('original','resampled')
xlabel('Time'), ylabel('Amplitude')
axis([0,1,-1,1])
disp('END OF EXAMPLE 6.1)')

%
%% 6.2 Efficient Implementation of Polynomial Interpolation Filters Using Farrow Structure
%
% Remainder
% Using the polynomial-based approach and the Farrow  structure, the sampling rate conversion by an arbitrary factor can be implemented efficiently. The above mentioned structure is indicated in the Figure.
%
%
% Farrow structure for polynomial-based interpolation filters:
% (a) The overall structure. (b) FIR filter details.
% In the following example, we demonstrate the application of the Farrow structure for converting the sampling rate of the input signal by the factor 16/15.
%
%% Example 6.2
%
% This example demonstrates sampling rate alteration by the fractional factor R = 16/15.
%

clear all, close all 

% Generating the input signal x[n]

n = -2:32;
f1 = 0.1; f2 = 0.16;
x = n+2*sin(pi*f1*n)+5*cos(pi*f2*n)+1; 
n=0:length(x)-5;
figure (1)
subplot(2,1,1),stem(n,x(3:length(x)-2),'k');
xlabel('Time index [n]'), ylabel('x[n]')
title('Original signal')

R = 16/15; %  Re-sampling factor
s = 1:15;
mu = [0,1-(1-1/R)*s];
mu = repmat(mu,1,3); % Vector of fractional intervals

% Farrow filter, coefficients:

C0 = [1/6,-1/2,1/2,-1/6]; % Vector of Lagrange coefficients for C_0(k)
C1 = [0,1/2,-1,1/2]; % Vector of Lagrange coefficients for C_1(k)   
C2 = ([-1/6,1,-1/2,-1/3]); % Vector of Lagrange coefficients for C_2(k)
C3 = [0,0,1];  % Vector of Lagrange coefficients for C_3(k)

% Farrow filter, initial states:

xs0 = [0,0,0]; xs1 = [0,0,0]; xs2 = [0,0,0]; xs3 = [0,0]; 

% Farrow filtering:

l=1;
for n=1:1:31
    xnl = x(n+3);
[v0,xs0] = filter(C0,1,xnl,xs0);
[v1,xs1] = filter(C1,1,xnl,xs1);
[v2,xs2] = filter(C2,1,xnl,xs2);
[v3,xs3] = filter(C3,1,xnl,xs3);
if mu(l)==0
     y(l) = ((v0*mu(l)+v1)*mu(l)+v2)*mu(l)+x(n+2);
    y(l+1) = ((v0*mu(l+1)+v1)*mu(l+1)+v2)*mu(l+1)+x(n+2);
    l = l + 2;
else
    y(l) = ((v0*mu(l)+v1)*mu(l)+v2)*mu(l)+x(n+2);
    l = l + 1;
end
end

% Displaying the results

figure 
subplot(2,1,2)
stem(0:length(y) -1,y,'sk'), xlabel ('Time index [l]'),ylabel('y[l]')
title('Re-sampled signal')
axis([0,33,0,40])
disp('END OF EXAMPLE 6.2')

%
%% 6.3 Fractional Delay Filters
%% Example 6.3
% We generate impulse responses of a family of fractional delay filters using the MATLAB object dfilt.farrowlinearfd.

clear all, close all

% Designing family of fractional delay filters

for k = 0:9
    mu = k/10;                         
    h(k+1)= dfilt.farrowlinearfd(mu);
end

% Compute and plot phase delay characteristics

[phi,f] = phasedelay(h,512,2); % Computing the phase delay
figure 
plot(f,phi/pi,'k')
xlabel('\omega/\pi'), ylabel('Phase delay (samples)'), grid

% Compute and plot magnitude responses

[H,f] = freqz(h,512,2);  % Computing the frequency response
figure 
plot(f,abs(H),'k')
xlabel('\omega/\pi'), ylabel('Magnitude'), axis([0,1,0,1.1]),grid

% Aplying the fractional delay filter

t=0:9;
x=0.5*sin(pi*0.2*t)+(sin(pi*0.1*t));  % Original sequence
y1=filter(h(3),x); % 0.2-sample delay
y2=filter(h(6),x); % 0.5-sample delay
figure 
stem(t,x,'k'); hold on; stem(t-0.2,y1,'kv--','filled')
hold on; stem(t-0.5,y2,'ks-.','filled')
xlabel('n'),ylabel('Amplitude')
legend('\mu=0','\mu=0.2','\mu=0.5'), grid
disp('END OF EXAMPLE 6.3')

%
%% Example 6.4
%
% This example illustrates processing of a sinusoidal sequence when the fractional delay factor is changed from \mu =0\ldotp 4 to \mu =0\ldotp 2. Program generates the 1st order Lagrange fractional-delay filter, creates 10 samples of the sinusoidal sequence, and performs filtering with the sudden change of the delay factor from \mu =0\ldotp 4 to \mu =0\ldotp 2 at the time-index n = 5. The Matlab object dfilt.farrowlinearfd has been used.

clear all, close all
t = 0:9;
x = sin(2*pi*0.1*t);
for i=1:10
    if i <= 5, mu(i)=0.4;
    else
        mu(i)=0.2;
    end
end
h = dfilt.farrowlinearfd(mu);                           
h.PersistentMemory = true;
for i=1:10
    h.Fracdelay = mu(i);
    ty(i) = t(i)-mu(i);
    y(i) = filter(h,x(i));
end
figure 
stem(t,x); hold on; stem(ty,y,'rs')
xlabel('n'),ylabel('Amplitude')
legend('Original Signal','Delayed Signal')
grid
disp('END OF EXAMPLE 6.4')
disp('                     END OF CHAPRER VI')

%

