IIR Filter: Lattice Form

IIR Filter: Lattice Form

Here, you will build a desktop application (UI, user interface) to simulate how is filtered using IIR filter with lattice form. A number of discrete signals are used as test signals are impulse, step, real exponential, sinusoidal, random, square, angled triangle, equilateral triangle, and trapezoidal signals. In the UI, this IIR filter is used to filter signals and audio file (wav) as well.

Follow these steps below:

1. Type this command on MATLAB prompt:

guide

The GUIDE Quick Start dialog box will show up, as shown in figure below:

2. Select the Blank GUI (Default) template, and then click OK:

3. Drag all components needed, as shown in figure below:

4. Right click on DRAW button and choose View Callbacks and then choose Callback. Define pbDraw_Callback() function below to read signals parameters and displays chosen signal in axes1 component:

% --- Executes on button press in pbDraw.
function pbDraw_Callback(hObject, eventdata, handles)
% hObject    handle to pbDraw (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

global x;

% Reads n0, n1, and n2
n0 = str2num(get(handles.n0,'String'));
n1 = str2num(get(handles.n1,'String'));
n2 = str2num(get(handles.n2,'String'));
exp = str2num(get(handles.expCoef,'String'));

% Reads parameters for sinusoidal
A = str2num(get(handles.A,'String'));
Freq = str2num(get(handles.Freq,'String'));
Fase = str2num(get(handles.Phase,'String'));
Fs = str2num(get(handles.Fs,'String'));

% Reads parameters for square and triangle signals
width = str2num(get(handles.width, 'String'));

% Choose option from popup menu popupmenusignal1
switch get(handles.popupmenusignal1,'Value')  
    case 1
       n = [n1:n2]; x = [(n-n0) == 0];
       
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); title('Discrete Impulse Signal');

    case 2
       n = [n1:n2]; x = [(n-n0) >= 0];
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); title('Discrete Step Signal');
       
    case 3
       n = [n1:n2]; x = [(n-n0) >= 0].*[(exp).^(n-n0)];
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); 
       title('Discrete Real Exponential Signal');
       
    case 4
       n = [n1:n2]; x = [(n-n0) >= 0].*[A*sin(2*pi*(Freq)*((n-n0)/Fs)+Fase)];
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); 
       title('Discrete Sinusoidal Signal');
       
    case 5
       n = [n1:n2]; x = [(n-n0) >= 0].*(rand(1,(n2-n1+1))-0.5);
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); 
       title('Discrete Random Signal');

    case 6
       n = [n1:n2]; x = [(n-n0) >= 0];
       ny = [n1:n2]; xy = [(ny-width-n0-1) >= 0];
       x = x - xy;
       
       % Display discrete signal
       axes(handles.axes1)
       stem(n,x,'linewidth',2,'color','b'); title('Discrete Square Signal');

    case 7
       n = [n1:n2]; x = n.*[n >= 0];
       ny = [n1:n2]; xy = ny.*[(ny-width-1) >= 0];
       x = (x - xy)/width;
       nb = [n1+n0:n2+n0];
       
       % Reads signal range
       set(handles.n1,'string',(n1+n0));
       set(handles.n2,'string',(n2+n0));

       % Display discrete signal
       axes(handles.axes1)
       stem(nb,x,'linewidth',2,'color','b'); 
       title('Discrete Angled Triangle Signal');
 
    case 8
       n = [n1:n2]; x = n.*[n >= 0];
       x2 = [zeros(1,width), x(1:end-width)];
       x1 = -x;
       x1 = [zeros(1,0.5*width), x1(1:end-0.5*width)];
       nb1 = n1+n0; nb2 = n2+n0;
       nb = [nb1:nb2];
       
       % Read signal range
       set(handles.n1,'string',(nb1));
       set(handles.n2,'string',(nb2));

       % Returns discrete equilateral triangle signal
       x = (x + 2*x1+x2)/(0.5*width);
       
       % Display discrete signal
       axes(handles.axes1)
       stem(nb,x,'linewidth',2,'color','b'); 
       title('Discrete Equilateral Triangle Signal');
       
    case 9
       n = [n1:n2]; x = n.*[n >= 0];
       x2 = [zeros(1,width), x(1:end-width)];
       x3 = [zeros(1,3*width), x(1:end-3*width)];
       x1 = x;
       x1 = [zeros(1,2*width), x1(1:end-2*width)];
       nb1 = n1+n0; nb2 = n2+n0;
       nb = [nb1:nb2];
       
       % Reads signal range
       set(handles.n1,'string',(nb1));
       set(handles.n2,'string',(nb2));

       % Returns discrete trapezoidal signal
       x = (x - x2 - x1 + x3)/width;
       
       % Display discrete signal
       axes(handles.axes1)
       stem(nb,x,'linewidth',2,'color','b'); 
       title('Discrete Trapezoidal Signal');
end  

5. You can run the UI. Choose one of signals and then click DRAW button:

6. Define latc2dir() function to convert lattice form to direct form:

function [b] = latc2dir(K)
% Convert from LATTICE to DIRECT
% ---------------------------------------------------
% [b] = latc2dir(K)
% b = Coefficent of DIRECT form
% K = Reflection coefficient of lattice form
%
M = length(K); J = 1; A = 1;
for m=2:1:M
A = [A,0]+conv([0,K(m)],J); J = fliplr(A);
end
b=A*K(1);

6. Right click on FREQUENCY RESPONSE OF FILTER button and choose View Callbacks and then choose Callback. Define pbFreqResponse_Callback() function below to calculate FFT of designed IIR filter:

% --- Executes on button press in pbFrequencyResponse.
function pbFrequencyResponse_Callback(hObject, eventdata, handles)
% hObject    handle to pbFrequencyResponse (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

global x;

% Reads all filter coefficients
K2 = str2num(get(handles.K2,'String'));
K21 = str2num(get(handles.K21,'String'));
K1 = str2num(get(handles.K1,'String'));
K11 = str2num(get(handles.K11,'String'));
K0 = str2num(get(handles.K0,'String'));
K31 = str2num(get(handles.K31,'String'));
K3 = str2num(get(handles.K3,'String'));


% Convert from Lattice to Direct form
K = [K0 K1 K2 K3]

[b] = latc2dir(K);

% Calculates frequency response
[h,w] = freqz(b,1,'whole',2001);

axes(handles.axes2); 
plot(w/pi,20*log10(abs(h)),'linewidth',2,'color','r'); grid on;
set(gca,'color',[0,0,0]);
ax = gca;
ax.YLim = [-100 20];
ax.XTick = 0:.5:2;
xlabel('Normalized Frequency  (\times\pi rad/sample)')
ylabel('Magnitude (dB)'); title('Frequency Response of Filter')

7. You can run the UI. Choose one of signals and click on DRAW button. Then click on FREQUENCY RESPONSE OF FILTER button to calculate FFT of designed IIR filter:

8. Right click on FILTERING SIGNAL button and choose View Callbacks and then choose Callback. Define pbFilteringSignal_Callback() function below to read input signal and calculates filtered signal in axes3 component:

% --- Executes on button press in pbFilteringSignal.
function pbFilteringSignal_Callback(hObject, eventdata, handles)
% hObject    handle to pbFilteringSignal (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

global x;

% Reads all reflection coeffs
K2 = str2num(get(handles.K2,'String'));
K21 = str2num(get(handles.K21,'String'));
K1 = str2num(get(handles.K1,'String'));
K11 = str2num(get(handles.K11,'String'));
K0 = str2num(get(handles.K0,'String'));
K31 = str2num(get(handles.K31,'String'));
K3 = str2num(get(handles.K3,'String'));


% Conver to direct form
K = [K0 K1 K2 K3]

[b] = latc2dir(K);

% Calculates filtered signal
y = filter(b,1,x);

axes(handles.axes3); 
t = 0:length(x)-1;  

stem(t,y,'linewidth',1,'color','y'); 
title('Filtered Signal')
set(gca,'color',[0,0,0]);

9. Right click on LOAD AUDIO FILE button and choose View Callbacks and then choose Callback. Define pbLoadAudio_Callback() function below to read audio signal and display it in axes2 component:

% --- Executes on button press in pbLoadAudio.
function pbLoadAudio_Callback(hObject, eventdata, handles)
% hObject    handle to pbLoadAudio (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

global x_audio;

[filename,path] = uigetfile({'*.wav'},'load audio file');
[x,Fs] = audioread([path '/' filename]);
handles.x = x ./ max(abs(x));
handles.Fs = Fs;
axes(handles.axes2);
time = 0:1/Fs:(length(handles.x)-1)/Fs;
plot(time, handles.x);
axis([0 max(time) -1 1]); title('Audio Signal');

axes(handles.axes3);
specgram(handles.x, 1024, handles.Fs);
title('Signal Spectrum');
handles.fileLoad = 1;

handles.fileNoise = 0;
handles.fileFinal = 0;

set(handles.freqSamplingAudio, 'String', num2str(Fs));
set(handles.NSamples, 'String', num2str(length(handles.x)));

x_audio = handles.x;
guidata(hObject, handles);

10. You can run the UI. Choose one of audio file by clicking LOAD AUDIO FILE button. The audio signal will be displayed in axes2 component and signal spectrum will be displayed in axes3 component:

11. Right click on FILTERING AUDIO FILE button and choose View Callbacks and then choose Callback. Define pbFilteringAudio_Callback() function below to read audio signal and display it in axes2 component:

% --- Executes on button press in pbFilteringAudioFile.
function pbFilteringAudioFile_Callback(hObject, eventdata, handles)
% hObject    handle to pbFilteringAudioFile (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

global x_audio;

% Reads all reflection cofficients
K2 = str2num(get(handles.K2,'String'));
K21 = str2num(get(handles.K21,'String'));
K1 = str2num(get(handles.K1,'String'));
K11 = str2num(get(handles.K11,'String'));
K0 = str2num(get(handles.K0,'String'));
K31 = str2num(get(handles.K31,'String'));
K3 = str2num(get(handles.K3,'String'));


% Convert to direct form
K = [K0 K1 K2 K3]
[b] = latc2dir(K);

% Calculates filtered signal
y = filter(b,1,x_audio);

axes(handles.axes2); 
t = 0:length(y)-1; 
plot(t,y,'linewidth',1,'color','r');
title('Filtered Signal')
set(gca,'color',[0,0,0]);

axes(handles.axes3)
specgram(y, 1024, handles.freqSamplingAudio);
title('Signal Spectrum');
handles.fileLoad = 1;
set(gca,'color',[0,0,0]);
handles.x=y;
guidata(hObject, handles);

12. Now, you can run the UI and click FILTERING AUDIO FILE button and see the filtered audio and its spectrum: