Signals and Systems ___

NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 01
OBJECT
Plotting of Basic Signals in MATLAB

THEORY
Signals

In the fields of communications, signal processing, and in electrical engineering more generally, a signal
is any timevarying or spatialvarying quantity. This variable (quantity) changes in time.

Speech or audio signal: Sound amplitude that varies in time

Temperature readings at different hours of a day
Stock price changes over days etc.

A discrete signal or discretetime signal is a time series, perhaps a signal that has been sampled
Signals can be classified by continuestime signal and discretetime signal:

A digital signal is a discretetime signal that takes on only a discrete set of values.
from continuous time signal.

MATLAB
MATrix LABoratory (MATLAB) is a powerful high-level programming language for scientific
computations. It supports a rich suite of mathematical, statistical and engineering functions and its
functionality is extended with interactive graphical capabilities for creating 2D as well as 3D plots. It
provides comprehensive toolboxes and various sets of algorithms.

Default Desktop Environment

a) Command Window
The main window in which commands are keyed in after the command prompt >>.
Results of most printing commands are displayed in this window.
b) Command History Window
This window records all of the executed commands as well as the date and time when these
commands were executed. This feature comes very handy when recalling previously executed
commands. Previously entered commands can also be re-invoked using up arrow key.
c) Current Directory Window
This window keeps track of the files in the current directory.

3
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

d) Workspace
This window is used to organize the loaded variables and displays the information such as size and
class of these variables.

Signals Representation
A signal in MATLAB is represented by a vector:
Examples:
x = [2, 3, -5, -3, 1]
n = 2:3:17 %( here step size is 3)
n = 2:17 %( Default Step size 1 is used)
Plotting in MATLAB
While plotting in MATLAB one must be careful that a vector is plotted against a vector and lengths

plot (for CT signals)

of vectors must match. Two functions are used for plotting:

stem (for DT signals )

Example:
x = 10sint

MATLAB Commands:
t = [-2:0.002:2]
x = 10 * sin (pi * t)
plot(t, x)
title(Example Sinusoid)
xlabel(time(sec))
ylabel(Amplitude)

Multiple Plots
For drawing multiple signals on the same graph, write first signals x and y axis vectors followed by
the next signal.
Syntax:
plot(X1,Y1,,Xn,Yn)
In order to differentiate them by colors, write line style specifier and color code.

4
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Example
plot(t, y, r-, t, x, g-);
legend(Sine curve, Cosine curve);

Generating Subplots
x=10*sin(-2*pi*t)
y=10*cos(-2*pi*t)
u=10*sin(-5*pi*t)
v=10*cos(-5*pi*t)
t = [-2:0.002:2]
subplot(2, 2, 1), plot(t, x);
xlabel(t),ylabel(x);
subplot(2, 2, 2), plot(t, y);
xlabel(t),ylabel(y);
subplot(2, 2, 4), plot(t, u);
xlabel(t),ylabel(u);
subplot(2, 2, 3), plot(t, v);
xlabel(t), ylabel(v);

5
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

DT Plots
Example:
Plot the DT sequences:
x = [2, 3, -1, 5, 4, 2, 3, 4, 6, 1]
x = [2, 3, -1, 5, 4, 2, 3, 4, 6, 1]
n = -6:3;
stem(n, x);

Zero & One Vectors

To generate zero or one vectors, use following statements:
zeros(1, 5)
Output: [0 0 0 0 0]
ones(1, 5)
Output: [1 1 1 1 1]

Some Common Types of Signals

For unit step:

x=(ones(size(n))).*(n>=0)
Note: size(n) returns dimension and number of elements in the array
For ramp:
x=n.*(n>=0)
For exponential:
x=(ones(size(n))).*(n>=0) %unit step
y=((0.5).^n).*x
stem(n,y)
For Rectangle:
t=-1:0.001:1;
y=rectpuls(t);
plot (t,y);
Triangle:
t=-1:0.001:1;
y=tripuls(t);
plot (t,y);
Sawtooth:
fs = 10000;
t = 0:1/fs:1.5;
x = sawtooth(2*pi*50*t);
plot(t,x), axis([0 0.2 -1 1]);
Square wave:
t=0:20;
y=square(t);
plot(t,y)
Sinc function:
t = -5:0.1:5;
y = sinc(t);
plot(t,y)

EXERCISES

1. Write MATLAB code to plot function x = sin(nx). Generate 8 subplots using for loop. Use step size of
0.05.

6
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

2. Write a sequence of MATLAB commands in the space below to plot the curves y = cosx and y
= x for 0 x 2 on the same figure. Then Zoom in to determine the point of intersection of the
two curves (and, hence, the root of x = cosx) to two significant figures. Your plot must be
properly labeled.

3. Draw graphs of the functions for x =0:0.1:10 and label your graph properly.

i. y = sin(x)/x

ii. u = (1/(x-1)2)+x

iii. v = (x2+1)/(x2-4)

iv. z = ((10-x)1/3-1)/(4 - x2)1/2

7
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

4. Write MATLAB commands to plot following elementary DT signals. Properly label your graphs.
a) Unit Step
b) Unit Ramp
c) Real Exponential: x(n) = 2(0.25)n , 0 < n < 10

5. Generate multiple plots with the following data:

Suppose A=1, f=1Hz, t=0:0.01:1: y(t)=cos(2t) y(t)=cos(2t+/2) y(t)=cos(2t-/2) y(t)=sin(2t)
where A is the amplitude of signal. Use colors & line styles to distinguish the plots.

8
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 02
OBJECT
Implementing Elementary Operations on Signals

THEORY
Operation with signals means to add, subtract, multiply, divide, scale, exponentiation, shift, delay or
advance and to flip signals. MATLAB allows all these operations but we need to be careful in
computation because the vector representation of the signals should have the same time origins and the
same number of elements.

Basic Signal Operations:

Given the signals x1 and x2 perform the following operations:
y1 = x1 + x2; y2 = x1 - x2; y3 = x1 * x2; y4 = x1 / x2; y5 = 2 x1; y6 = (x1)3

Given is MATLAB code for first six basic operations on sinusoids:

x1 = 5*sin((pi/4)*[0:0.1:15]);
x2 = 3*cos((pi/7)*[0:0.1:15]);
% Plotting the signals
subplot(2,4,1), plot(x1)
title('x1 = 5 sin(pi/4)t ')
xlabel(' time (sec) ')
ylabel('x1 (volts) ')
subplot(2,4,2), plot(x2)
title('x2 = 3 cos(pi/7)t ')
xlabel(' time (sec) ')
ylabel('x2 (volts) ')
y1 = x1 + x2; % addition
y2 = x1 - x2; % subtraction
y3 = x1 .* x2; % multiplication
y4 = x1 ./ x2; % division
y5 = 2*x1; % scaling
y6 = x1 .^3; % exponentiation
% Plotting the signals
subplot(2,4,3), plot(y1)
title('y1 = x1 + x2 ')
xlabel(' time (sec) ')
ylabel('y1 (volts) ')
subplot(2,4,4), plot(y2)
title('y2 = x1 x2 ')
xlabel(' time (sec) ')
ylabel('y2 (volts) ')
subplot(2,4,5), plot(y3)
title('y3 = x1 * x2 ')
xlabel(' time (sec) ')
ylabel('y3 (volts)^2 ')
subplot(2,4,6), plot(y4)
title('y4 = x1 / x2 ')
xlabel(' time (sec) ')
ylabel('x1/x2 ')
subplot(2,4,7), plot(y5)
title('y5 = 2*x1 ')
xlabel(' time (sec) ')
ylabel('y5 (volts) ')
subplot(2,4,8), plot(y6)
title('y6 = x1 ^ 3 ')
xlabel(' time (sec) ')
ylabel('y6 (volts)^3 ')

9
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Plotting Piecewise Signals

Define a piecewise continuous function h(t)

t = [-1:0.01:2] ;
h=1*[t >= -1 & t <= 0];
h=h+[-1 * (t > 0 & t <= 2)];
plot(t,h,'linewidth',2)
grid
axis([-2 3 -1.2 1.2])

Exponentially Varying Sinusoids

Consider the following example:

Multiplication of two sinusoids

Consider the following signal

10
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

More on Signal Operations

Time Shifting
One of the most basic operations in a DSP system is to shift the time reference. The shift may be
considered either a delay or an advance.

MATLAB function for Signal Shifting:

function [y, n] = sigshift(x, n, k)
n = n + k;
y = x;

Example:

Shift x (t) by 0.5 units on the right

Time Scaling
If f(t) is compressed in time by a factor of a where (a>1) then the resulting signal is given by:
(t) = f(at)

Similarly if f(t) is expanded in time by a factor of a where (a>1) then the resulting signal is given by:
(t) = f(t/a)

If x = cos(t); for t = -0.5:0.01:0.5

Find: (t) = x(2*t) & (t) = x(t/2)

11
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Time Inversion (Folding)

Folding involves the reversal of the time axis.

(t) = f(-t)
(-t) = f(t)

To invert a signal we replace t by t.

Consider the following example:

Operations on DT Signals

Signal Addition
This is sample-by-sample addition given by:
{x1(n)}+{x2(n)}={ x1(n) + x2(n)}
Implemented in MATLAB by + operator. However, the lengths of x1(n) and x2(n) must be same.

Signal Multiplication
This is sample-by-sample multiplication given by:
{x1(n)}.{x2(n)}={ x1(n) x2(n)}
It is implemented in MATLAB by * operator. However, the lengths of x1(n) and x2(n) must be same.

Delayed (Advanced) Impulse

Formula used is: {
MATLAB function
function [y, n] = impseq(n1, n2, k)
n = n1:n2;
y = [(n k) == 0];

Time Shifting, Inversion & scaling are also performed sample by sample on DT Signals.

EXERCISES

1. Write down MATLAB code to perform the signals operations discussed in lab.

12
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

2. The operations discussed in lab can also be performed on DT signals (sample by sample).
Write MATLAB scripts to generate and plot each of the following signals for the intervals
indicated.

a) [ ] [ ]

13
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

b) [ ] [ ]

14
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 03
OBJECTIVE

Sampling & reconstruction of CT sinusoids to understand aliasing phenomenon

THEORY

Sampling Principle

A CT sinusoid containing a maximum frequency of Fmax must be sampled at a sampling rate of

Fs > 2Fmax (Nyquist Rate) to avoid aliasing. If sampling rate is greater than Nyquist Rate, then
CT sinusoid can be uniquely recovered from its DT version.
Analog frequencies separated by integral multiple of a given sampling rate are alias of each other
Any CT sinusoid of frequency Fk when sampled at the sampling rate Fs will result in the same
DT sinusoid as does the CT sinusoid of frequency F0 sampled at Fs, where:
Fk = F0 + kFs where k=1, 2, 3
Assume two CT sinusoidal signals
x1(t) = cos(2t) F1 = 1 Hz
x2(t) = cos(6t) F2 = 3 Hz

These signals can be plotted using the MATLAB code shown:

t = -2:0.005:2;
x1 = cos(2*pi*t);
x2 = cos(6*pi*t);
subplot(3,2,1),
plot(t, x1);
axis([-2 2 -1 1]);
grid on;
xlabel('t'), ylabel('cos2\pit');
subplot(3,2,2),
plot(t, x2);
axis([-2 2 -1 1]);
grid on;
xlabel('t'), ylabel('cos6\pit');

Plotting DT Sinusoid

x1n = cos(2*pi*[-2:1/2:2]); %We will not plot it against n because this will depict
sampled signal incorrectly. Generate another vector from n as follows:
k = -2:length(n)-3;
subplot(3,2,3), stem(k, x1n);
axis([-2 length(n)-3 -1 1]);
grid on;
xlabel('n'),
ylabel('cos\pin');
X2n= cos(6*pi*n)
subplot(3,2,4),stem(k, x2n);
axis([-2 length(n)-3 -1 1]);
grid on;
xlabel('n'),
ylabel('cos3\pin');

15
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Reconstruction

D/A conversions are performed using interpolation. There are various approaches to interpolation
Zero-Order-Hold (ZOH)
Step interpolation
First-Order-Hold (FOH)
Samples are connected by straight lines
Cubic Spline Interpolation
It is invoked by spline (n, xn, t)

EXERCISES

1. Write a MATLAB script to carry out the following tasks:

a. Plot of two continuous time sinusoids x1(t) =cos2t and x2(t) = cos14t for 0 < t < 5 seconds.
Choose a suitable time step.
b. Sample them using a sampling rate of 3 samples/s.
c. Plot the resulting the resulting discrete time sinusoids.
d. Reconstruct the signals
e. You must divide the figure into six subplots so that reconstructed signal and DT version of
x1(t) lie in the same column and so are those of x2(t).

2. Does aliasing occur? Briefly explain why?

16
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

3. Consider an analog signal xa(t) = sin(20t), 0 < t < 1. It is sampled at sampling intervals of 0.01, 0.05,
and 0.1 seconds to obtain x(n).
a. Write MATLAB scripts to plot x(n) for each sampling interval
b. Reconstruct the analog signal ya(t) from the samples x(n) using the cubic spline interpolation and
determine the frequency in ya(t) from your plot. (Ignore the end effects)
c. Comment on your results.

17
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 04
OBJECTIVE
Understanding Fourier Series

THEORY
In mathematics, a Fourier series decomposes periodic functions or periodic signals into the sum of a
(possibly infinite) set of simple oscillating functions, namely sines and cosines (or complex exponentials).

THEOREM:
Let x(t) be a bounded periodic signal with period T. Then x(t) can be expanded as a weighted sum of
sinusoids with angular frequencies that are integer multiples of 0 =2/T:
X (t) = a0 + a1 cos(0t) + a2 cos(20t) + a3 cos(30t) + . . .+ b1 sin(0t) + b2 sin(20t) + b3 sin(30t) +
...
This is the trigonometric Fourier series expansion of x(t).

Consider the given signal. Find the Trigonometric Fourier series coefficients and plot the magnitude
and phase spectra for the periodic signals shown below in MATLAB:

Fourier series coefficients for the periodic signal shown above are:

( ) ( )

The code for plotting the magnitude and phase spectra for the given co-efficients is given below
using for loop:

n=1:7;
a0=0.504;
b0=0.504*(8*0/(1+16*0^2)); % b0=0;
Cn=a0;

18
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

theta0=atan(-b0/a0);
thetan=theta0;
den=(1+16*n.^2);
N=length(den);
for i=1:N
an(i)=0.504*2/den(i);
bn(i)=0.504*8*n(i)/den(i);
cn=sqrt(an(i)^2+bn(i)^2);
Cn=[Cn cn];
theta=atan(-bn(i)/an(i));
thetan=[thetan theta];
end
n=0:7;
subplot(211),plot(n,,'o'),grid, xlabel('n'),ylabel(C_n),title(Fourier Series)
subplot(212),plot(n,thetan,'o'),grid,xlabel('n'),ylabel('\theta_n (rad)')

EXERCISE

1. If f (t) is defined below:

Find the trigonometric Fourier series coefficients for the periodic signal shown below and use them
to plot the magnitude and phase spectra.

19
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

20
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 05
OBJECT
To study Convolution and its Applications

THEORY
What is Convolution?

Convolution is a mathematical method for combining two signals to form a third signal. It is one of
the most important techniques in Digital Signal Processing. Convolution is important because it
relates the three signals of interest: the input signal, the output signal, and the impulse response.

Expression for discrete time convolution;

Where; y: output signal

x: input signal
h: impulse response of system

Impulse Response of a System:

Impulse response of a system can be regarded as the transfer function of a system. If we know a
system's impulse response, then we can calculate what the output will be for any possible input
signal. This means we know everything about the system. There is nothing more that can be learned
about a linear system's characteristics.

The impulse response goes by a different name in some applications. If the system being considered
is a filter, the impulse response is called the filter kernel, the convolution kernel, or simply,
the kernel. In image processing, the impulse response is called the point spread function. While
these terms are used in slightly different ways, they all mean the same thing, the signal produced by a
system when the input is a delta function.

How to evaluate the impulse response of a system?

First, the input signal can be decomposed into a set of impulses, each of which can be viewed as a
scaled and shifted delta function. Second, the output resulting from each impulse is a scaled and
shifted version of the impulse response. Third, the overall output signal can be found by adding these
scaled and shifted impulse responses.

21
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

An input signal, x[n], enters a linear system with an impulse response, h[n], resulting in an output
signal, y[n]. In equation form: x[n] * h[n] = y[n]. Expressed in words, the input signal convolved with
the impulse response is equal to the output signal. Just as addition is represented by the plus, +, and
multiplication by the cross, , convolution is represented by the star, *.

22
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

LAB PERFORMANCE
Convolution in MATLAB

The MATLAB command 'conv' computes the convolution of two vectors.

A single vector can be used to represent a discrete time signal with `t` distance between each sample.
Two vectors represent two distinct discrete time signals.

Let two digital signals of different lengths are given as v1 = 1 0 1 1 0 1 1 0 0 and v2= 0 1 0 apply
convolution and write down the resultant signal.

MATLAB Code:

v1 = [1 0 1 1 0 1 1 0 0];
v2 = [0 1 0];
u = conv(v1, v2);

In the above program the value of the resultant signal should be passed to the variable mentioned on
the left hand side of equality sign where we have used the conv function.

Task: Create subplots for v1, v2 and u signals.

EXERCISES

1. Write MATLAB code to perform convolution between a square pulse [0 0 0 0 4 4 4 4 0 0 0 0] and a

triangular wave [0 0 0 1 2 3 4 3 2 1 0 0]. Using stem instruction, generate subplots for the two given
vectors and the resultant vector (Attach subplots to this lab session)

2. Write a sequence of MATLAB statements to compute convolution of signal x(t) and impulse
response h(t) defined for times t = [-10,+10];

23
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

3. Perform and write down MATLAB code of any program that employs convolution operations
by the advice of lab instructor.

4. Write MATLAB commands to perform convolution between signal a & b for time interval [-20,
+20]. Also do attach a printout of resultant signal.

a) Gaussian curve
b) Unit Ramp signal

24
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

25
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 06
OBJECT
Performing statistical measurements on signals

THEORY

This lab session shows how to perform statistical measurements on an input data stream using
DSP System Toolbox functionality available at the MATLAB command line. You will compute
the signal statistics minimum, maximum, mean and the spectrum and plot them.
Introduction
This lab session computes signal statistics using DSP System Toolbox System objects. These
objects handle their states automatically reducing the amount of hand code needed to update
states reducing the possible chance of coding errors.
These System objects pre-compute many values used in the processing. This is very useful when
you are processing signals of same properties in a loop. For example, in computing an FFT, the
values of sine and cosine can be computed and stored once you know the properties of the input
and these values can be reused for subsequent calls. Also the objects check only whether the
input properties are of same type as previous inputs in each call.
Initialization
Here you initialize some of the variables used in the code and instantiate the System objects used
in your processing. These objects also pre-compute any necessary variables or tables resulting in
efficient processing calls later inside a loop.
frameSize = 1024; % Size of one chunk of signal to be processed
in one loop

% Here you create a System object to read from a specified audio

file and
% set its output data type.
hfileIn = dsp.AudioFileReader(which('speech_dft.mp3'), ...
'SamplesPerFrame', frameSize, ...
'OutputDataType', 'double');

fileInfo = info(hfileIn);
Fs = fileInfo.SampleRate;

Create an FFT System object to compute the FFT of the input.

hfft = dsp.FFT;
Create System objects to calculate mean, minimum and maximum and set them to running mode.
In running mode, you compute the statistics of the input for its entire length in the past rather
than the statistics for just the current input.
hmean = dsp.Mean('RunningMean', true);

26
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

hmin = dsp.Minimum('RunningMinimum', true);

hmax = dsp.Maximum('RunningMaximum', true);

Create audio output System object. Note that audio output to speakers from small chunks of data
in a loop is possible using the AudioPlayer System object. Using sound or audioplayer in
MATLAB, either overlaps or introduces gaps in audio playback.

haudioOut = dsp.AudioPlayer('SampleRate', Fs, ...

'QueueDuration', 1);

% Initialize figures for plotting

s = hfigsstats(frameSize, Fs);

Stream Processing Loop

Here you call your processing loop which will calculate the mean, min, max, FFT, and filter the
data using the System objects.
Note that inside the loop you are reusing the same FFT System object twice. Since the input data
properties do not change, this enables reuse of objects here. This reduces memory usage. The
loop is stopped when you reach the end of input file, which is detected by the AudioFileReader
object.
while ~isDone(hfileIn)
% Audio input from file
sig = step(hfileIn);

% Compute FFT of the input audio data

fftoutput = step(hfft, sig);
fftoutput = fftoutput(1:512); % Store for plotting

% The hmean System object keeps track of the information

about past
% samples and gives you the mean value reached until now.
The same is
% true for hmin and hmax System objects.
meanval = step(hmean, sig);
minimum = step(hmin, sig);
maximum = step(hmax, sig);

% Play output audio

step(haudioOut, sig);

% Plot the data you have processed

s = plotstatsdata(s, minimum, maximum, sig.', meanval,
fftoutput);
end
pause(haudioOut.QueueDuration); % Wait for audio to finish

27
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Release
Here you call the release method on the System objects to close any open files and devices.
release(hfileIn);
release(haudioOut);

Conclusion
You have seen visually that the code involves just calling successive System objects with
appropriate input arguments and does not involve maintaining any more variables like indices or
counters to compute the statistics. This helps in quicker and error free coding. Pre-computation
of constant variables inside the objects generally lead to faster processing time.

EXERCISE:

1. For the same data signal and frame size used in lab session, determine the following statistical
measurements: standard deviation, variance, histogarm and peak to peak value.
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________

28
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________

29
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 07
OBJECTIVE

Understanding Z-transform for the analysis of LTI Systems

THEORY
Z Transform

Z transform is a fundamental tool for the analysis of LTI systems. Especially, it plays very significant
role in digital filter design. Hence, it becomes vitally important to understand its use and draw
conclusions about systems using z transform analysis.

MATLAB Functions for Z Transform

Poles & Zeroes

The poles and zeroes of a system can be found by applying MATLAB function roots to the
denominator and numerator polynomials of its system function H(z). The arguments to the function
roots are the coefficients of the respective polynomial in the ascending powers of z-1.

Example
To find the roots of the polynomial 1 + 4z-1 6z-2.

Code:
roots([1 4 -6])
Output:
-5.1623 1.1623
The function tf2zp(num, den) can also be used to compute poles and zeros of a system, where num
and den have their usual meanings.
[z, p, k] = tf2zp(num, den)
Zeros and poles are respectively returned in the column vectors z and p, whereas, the gain constant is
returned in the variable k. However, the length of vectors num and den must be made equal by zero
padding.

Example
Determine the poles and zeros of the following systems function using tf2zp function.

MATLAB code:
num=[1 0.2 -0.15]
den= [1 -0.6 -0.19 0.144 -.018]
[z,p,k]=tf2zp(num,den)
Output:
Zeros: 0.3 -0.5 0 0
Poles: 0.3 -0.5 0.2 0.6
The function [num, den] = zp2tf(z, p, k) does the reverse of the function tf2zp.

30
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Zplane
The pole-zero plot can be obtained using the function zplane. The arguments to this function can be
specified in two ways:

1. Passing coefficients of numerator and denominator polynomials of system function: they are
specified as two row vectors
Example:

MATLAB code:
b = 0.094 * [1, 4, 6, 4, 1];
a = [1, 0, 0.4860, 0, 0.0177];
zplane(b, a);

Output:

2. Passing zeros and poles

The function zplane can also be passed zeros and poles of the system. However, the difference must
be appreciated between passing (zeros, poles) and (num, den) as arguments. That is, zeros and poles
are entered as column vectors while num and den are entered as row vectors.

Factored Form of H(z)

From the pole-zero description, the factored form of the transfer function can be obtained using the
function sos = zp2sos(z, p, k), where sos stands for second-order sections. The function
computes the coefficients of each second-order section (factor) given as an L x 6 matrix sos, where,

The lth row of sos contains the coefficients of the numerator and the denominator of the lth second-
order factor of the z-transform G(z):

31
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

EXERCISES
1. Write MATLAB command to draw pole-zero plot of the following system.

2. Use tf2zp(num, den)to determine poles, zeros and gain constant of the system function in Q.1

3.Write MATLAB command(s) to determine numerator and denominator polynomials of a rational

system function using the following information:
zeros
-3.0407
0.9211 + j0.8873
0.9211 - j0.8873
-0.4008 + j0.9190
-0.4008 - j0.9190
poles
-0.3987 + j0.9142
-0.3987 - j0.9142
0.5631 + j0.5424
0.5631 - j0.5424
-0.3289
gain constant
0.2

32
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

4. Use MATLAB to determine the poles, zeros and factored form of the following rational z transform:

Report your findings.

33
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 08
OBJECT
Designing low pass FIR Filters

THEORY

FIR filters are widely used due to the powerful design algorithms that exist for them, their
inherent stability when implemented in non-recursive form, the ease with which one can attain
linear phase, and the ample hardware support that exists for them among other reasons. This
example showcases functionality in the DSP System Toolbox for the design of low pass FIR
filters with a variety of characteristics.

A Simple Low Pass Filter Design

An ideal low pass filter requires an infinite impulse response. Truncating (or windowing) the
impulse response results in the so-called window method of FIR filter design. We will use filter
design objects (fdesign) throughout this lab session. Consider a simple design of a low pass filter
with a cutoff frequency of 0.4*pi radians per sample:

Fc = 0.4;
N = 100; % FIR filter order
Hf = fdesign.lowpass('N,Fc',N,Fc);

We can design this low pass filter using the window method. For example, we can use a
Hamming window or a Dolph-Chebyshev window:

Hd1 = design(Hf,'window','window'
@hamming,'SystemObject',true);
Hd2 =
design(Hf,'window','window',{@chebwin,50},'SystemObject',true);
hfvt = fvtool(Hd1,Hd2,'Color','White');
legend(hfvt,'Hamming window design','Dolph-Chebyshev window
design')

34
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

The choice of filter order was arbitrary. Since ideally the order should be infinite, in general, a
larger order results in a better approximation to ideal at the expense of a more costly
implementation. For instance, with a Dolph-Chebyshev window, we can decrease the transition
region by increasing the filter order:

Hf.FilterOrder = 200;
Hd3 =
design(Hf,'window','window',{@chebwin,50},'SystemObject',true);
hfvt = fvtool(Hd2,Hd3,'Color','White');
legend(hfvt,'Dolph-Chebyshev window design. Order = 100',...
'Dolph-Chebyshev window design. Order = 200')

35
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Minimum-Order Low Pass Filter Design

In order to determine a suitable filter order, it is necessary to specify the amount of passband
ripple and stopband attenuation that will be tolerated. It is also necessary to specify the width of
the transition region around the ideal cutoff frequency. The latter is done by setting the passband
edge frequency and the stopband edge frequency. The difference between the two determines the
transition width.

Fp = 0.38;
Fst = 0.42;
% Fc = (Fp+Fst)/2; Transition Width = Fst - Fp
Ap = 0.06;
Ast = 60;
setspecs(Hf,'Fp,Fst,Ap,Ast',Fp,Fst,Ap,Ast);

We can still use the window method, along with a Kaiser window, to design the low pass filter.
Hd4 = design(Hf,'kaiserwin','SystemObject',true);
measure(Hd4)

36
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

ans =

Sample Rate : N/A (normalized frequency)

Passband Edge : 0.38

3-dB Point : 0.39539

6-dB Point : 0.4

Stopband Edge : 0.42

Passband Ripple : 0.016058 dB

Stopband Atten. : 60.092 dB

Transition Width : 0.04

One thing to note is that the transition width as specified is centered around the cutoff frequency
of 0.4 pi. This will become the point at which the gain of the low pass filter is half the passband
gain (or the point at which the filter reaches 6 dB of attenuation).

EXERCISE:

1. Follow the steps explained in lab session, design a minimum order ideal low pass FIR filter
with cutoff frequency 0.2*pi radians per sample. Use Kaiser window to design the filter.
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________

37
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 09
OBJECT
Designing low pass IIR Filters

THEORY

In this lab session, we will learn how to design classic IIR filters. The initial focus is on the
situation for which the critical design parameter is the cutoff frequency at which the filter's power
decays to half (-3 dB) the nominal passband value.
This lab illustrates how easy it is to replace a Butterworth design with either a Chebyshev or an
elliptic filter of the same order and obtain a steeper rolloff at the expense of some ripple in the
passband and/or stopband of the filter. After this, minimum-order designs are explored.

Lowpass Filters
Let's design an 8th order filter with a normalized cutoff frequency of 0.4pi. First, we design a
Butterworth filter which is maximally flat (no ripple in the passband or in the stopband):

N = 8; F3dB = .4;
d = fdesign.lowpass('N,F3dB',N,F3dB);

Hbutter = design(d,'butter','SystemObject',true);

A Chebyshev Type I design allows for the control of ripples in the passband. There are still no
ripples in the stopband. Larger ripples enable a steeper rolloff. Here, we specify peak-to-peak
ripples of 0.5dB:

Ap = .5;
setspecs(d,'N,F3dB,Ap',N,F3dB,Ap);
Hcheby1 = design(d,'cheby1','SystemObject',true);
hfvt = fvtool(Hbutter,Hcheby1,'Color','white');
axis([0 .44 -5 .1])
legend(hfvt,'Butterworth','Chebyshev Type I');

38
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

A Chebyshev Type II design allows for the control of the stopband attenuation. There are no
ripples in the passband. A smaller stopband attenuation enables a steeper rolloff. Here we specify
a stopband attenuation of 80 dB:

Ast = 80;
setspecs(d,'N,F3dB,Ast',N,F3dB,Ast);
Hcheby2 = design(d,'cheby2','SystemObject',true);
hfvt = fvtool(Hbutter,Hcheby2,'Color','white');
axis([0 1 -90 2])
legend(hfvt,'Butterworth','Chebyshev Type II');

39
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Finally, an elliptic filter can provide the steeper rolloff compared to previous designs by allowing
ripples both in the stopband and the passband. To illustrate that, we reuse the same passband and
stopband characteristic as above:

setspecs(d,'N,F3dB,Ap,Ast',N,F3dB,Ap,Ast);
Hellip = design(d,'ellip','SystemObject',true);
hfvt = fvtool(Hbutter,Hcheby1,Hcheby2,Hellip,'Color','white');
axis([0 1 -90 2])
legend(hfvt, ...
'Butterworth','Chebyshev Type I','Chebyshev Type
II','Elliptic');

40
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

By zooming in the passband, we verify that all filters have the same -3dB frequency point and
that only Butterworth and Chebyshev Type II designs have a perfectly flat passband:
axis([0 .44 -5 .1])

Minimum Order Designs

In cases where the 3dB cutoff frequency is not of primary interest but instead both the passband
and stopband are fully specified in terms of frequencies and the amount of tolerable ripples, we
can use a minimum order design technique:

41
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Fp = .1; Fst = .3; Ap = 1; Ast = 60;

setspecs(d,'Fp,Fst,Ap,Ast',Fp,Fst,Ap,Ast);
Hbutter = design(d,'butter','SystemObject',true);
Hcheby1 = design(d,'cheby1','SystemObject',true);
Hcheby2 = design(d,'cheby2','SystemObject',true);
Hellip = design(d,'ellip','SystemObject',true);
hfvt = fvtool(Hbutter,Hcheby1,Hcheby2,Hellip, 'DesignMask',
'on',...
'Color','white');
axis([0 1 -70 2])
legend(hfvt, ...
'Butterworth','Chebyshev Type I','Chebyshev Type
II','Elliptic');

A 7th order filter is necessary to meet the specification with a Butterworth design whereas a 5th
order is sufficient with either Chebyshev techniques. The order of the filter can even be reduced
to 4 with an elliptic design:
order(Hbutter)

order(Hcheby1)

order(Hcheby2)

42
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

order(Hellip)

ans =

ans =

ans =

ans =

EXERCISE:
1. Follow the steps explained in lab session, compare Butterworth, Chebyshev Type I & II and
Ellliptic low pass IIR filter design with normalized cutoff frequency 0.2*pi radians per
sample.
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________
_________________________________________________________________________

43
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Lab Session 10
OBJECT
Exploring Amplitude Modulation Detection

THEORY

Envelope Detection

This lab session shows how to implement envelope detection using squaring and lowpass
filtering method.

Introduction
The signal's envelope is equivalent to its outline and an envelope detector connects all the peaks
in this signal. Envelope detection has numerous applications in the fields of signal processing and
communications, one of which is amplitude modulation (AM) detection. The following block
diagram shows the implementation of the envelope detection using squaring and lowpass
filtering method.

Squaring and Lowpass Filtering Method

This envelope detection method involves squaring the input signal and sending this signal
through a lowpass filter. Squaring the signal demodulates the input by using the input as its own
carrier wave. This means that half the energy of the signal is pushed up to higher frequencies and
half is shifted down toward DC. You then downsample this signal to reduce the sampling
frequency. You can do downsampling if the signal does not have any high frequencies which
could cause aliasing. Otherwise an FIR decimation should be used which applies a low pass filter
before downsampling the signal. After this, pass the signal through a minimum-phase, lowpass
filter to eliminate the high frequency energy. Finally you are left with only the envelope of the
signal.
To maintain the correct scale, you must perform two additional operations. First, you must
amplify the signal by a factor of two. Since you are keeping only the lower half of the signal
energy, this gain matches the final energy to its original energy. Second, you must take the
square root of the signal to reverse the scaling distortion that resulted from squaring the signal.
This envelope detection method is easy to implement and can be done with a low-order filter,
which minimizes the lag of the output.

44
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

Initialization
Initialize required variables such as those for the frame size and file name, see lab session 6.
Creating and initializing your System objects before they are used in a processing loop is critical
for getting optimal performance.
Fs = 22050;
numSamples = 10000;
DownsampleFactor = 15;
frameSize = 10*DownsampleFactor;
Create a sine wave System object and set its properties to generate two sine waves. One sine
wave will act as the message signal and the other sine wave will be the carrier signal to produce
Amplitude Modulation.
hsin = dsp.SineWave( [0.4 1], [10 200], ...
'SamplesPerFrame', frameSize, ...
'SampleRate', Fs);

Create a lowpass FIR filter for filtering the squared signal to detect its envelope.
hlowpass1 = dsp.FIRFilter(...
'Numerator', firpm(20, [0 0.03 0.1 1], [1 1 0 0]));

Create and configure time scope System object to plot the input signal and its envelope.
hts1 = dsp.TimeScope(...
'NumInputPorts', 2, ...
'Name', 'Envelope detection using Amplitude Modulation', ...
'SampleRate', [Fs Fs/DownsampleFactor], ...
'TimeDisplayOffset', [(N/2+frameSize)/Fs 0], ...
'TimeSpanSource', 'Property', ...
'TimeSpan', 0.45, ...
'YLimits', [-2.5 2.5]);
pos = hts1.Position;

Stream Processing Loop

Create the processing loop to perform envelope detection on the input signal. This loop uses the
System objects you instantiated.
for i=1:numSamples/frameSize
sig = step(hsin);
sig = (1+sig(:,1)).*sig(:, 2); % Amplitude modulation

% Envelope detector by squaring the signal and lowpass

filtering
sigsq = 2*sig.*sig;
sigenv1 = sqrt(step(hlowpass1, downsample(sigsq,
DownsampleFactor)));

45
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

% Plot the signals and envelopes

step(hts1, sig, sigenv1);
end

EXERCISE:

1. Perform the squaring and low pass filtering method of envelop detection for Fs=11025 and
numsamples=5000. Use appropriate downsample factor.
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________

46
Signals and Systems ___
NED University of Engineering & Technology Department of Computer & Information Systems Engineering

______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________
______________________________________________________________________________

47