Lab 1RT

Download as pdf or txt
Download as pdf or txt
You are on page 1of 8

EEE4114F - Lab Assignment 1

The aim of this assignment is to get you to be familiar with MATLAB,


and the functionality that it has for doing signal processing. MATLAB
also has excellent plotting facilities, which comes in useful when
preparing reports and theses.
Work through each of the following questions, and satisfy yourself that
you understand the content. In some cases part of a question may be
labelled as a Hand in task — you are expected to generate results for
these tasks, and compile them into a single document. You are welcome
to submit the report electronically, preferably as a PDF document.
Your hand-in solutions must be submitted by the due date.
MATLAB has a good help system. To get information on the usage of a
command, for example the plot command, type help plot at the
prompt. If you’re not sure what command you want, but know that it is
related to plotting, type lookfor plot at the prompt.
If you don’t want to use MATLAB you’re welcome to use any other
language. Python with numpy and matplotlib is a great option, and is
very close to MATLAB in syntax and functionality. Using this
combination through Jupyter notebook is convenient.

1
1. Basic MATLAB functions and plots
A basic cosine wave takes the form
x[n] = A cos(ω0 n + φ),
and is described completely by the parameters A (the amplitude),
ω0 (the frequency), and φ (the phase).
(a) Generate and plot each of the following sequences:
π 
x1 [n] = sin n , 0 ≤ n ≤ 25
17
π 
x2 [n] = sin n , −15 ≤ n ≤ 25
17
 π
x3 [n] = sin 3πn + , −10 ≤ n ≤ 10
 2

π
x4 [n] = cos √ n , 0 ≤ n ≤ 50.
23
In each case the horizontal axis should extend only over the
range indicated (see the axes command), and the axes should be
labelled accordingly (see xlabel and ylabel). Each sequence
should be displayed using stem.
It is important to be able to generate plots in MATLAB — they
are often required in reports, research papers, and theses. Other
MATLAB constructs that are useful are:
i. Variants of the plotting commands that allow multiple plots
on the same axes. For example, plot(n,s1,’k-’,n,s2,’g:’)
will plot the signal s1 as a solid black line, and s2 as a dotted
green line. Most other plot commands accept similar
arguments. The legend command is useful for labelling these
plots.
ii. The subplot command for placing multiple plot axes on the
same figure.
iii. The print command can be used for saving figures to file, for
later inclusion into reports. For example, the command

2
print(gcf, ’-depsc2’, ’outfile.eps’) will save an
Encapsulated PostScript representation of the current figure,
and print(gcf, ’-dpng’, ’outfile.png’) will save the
current figure in PNG format.
(b) Write a general MATLAB function that will generate a
finite-length sinusoid. The function will need five input
arguments: three for the parameters, and two more to indicate
the first and last n index of the signal. The function should
return a column vector that contains the values of the sinusoid.
Test the function by plotting the results for various choices of
the input parameters.
(c) Modify the function in part (b) to return two results: a vector of
indices over the range n, and the values of the signal.
(d) Hand in: Two plots, on the same figure (using subplot), of the
functions
 

x5 [n] = 3 cos n
10
 
2π π
x6 [n] = 2 cos n+
12 4
over the range −10 to 15.

2. A simple lowpass filter


The output y[n] of a linear filter is related to the input x[n] via the
relation
X∞
y[n] = h[k]x[n − k],
k=−∞

where h[n] is the filter impulse response. A very simple (and bad, as

3
you will discover later) lowpass filter has the impulse response

1 0≤n≤6
h[n] = 7
0 otherwise.

Because h[n] is non-zero only over a limited range, the limits of the
convolution sum can be narrowed in this case:
6
X
y[n] = h[k]x[n − k].
k=0

If we require the output at a given value of n, we can just calculate


the sum over seven product terms between elements of h and x.
This is easy (and instructive) to do using for loops in Matlab.
(a) Hand in: Suppose the input to the filter is the (one-sided)
sinusoidal sequence

x[n] = u[n] cos(ω0 n),

with ω0 = π/10. Calculate and plot the input and output


sequence on the same set of axes over the range −10, . . . , 70.
(There is a hold on command that will allow you to plot two
graphs on top of one another.) Note that after the initial
transient dies down, the output signal is also sinusoidal. What is
the magnitude and phase of this output relative to the input?
(b) Repeat the previous task for ω0 = 2π/5 and ω0 = 4π/5. Observe
that the filter attenuates high frequencies more than low
frequencies.
(c) If you didn’t already do so, explore the Matlab conv command
as a means of directly calculating the convolution between two
sequences. Note that the origin of the n-axis is nowhere
explicitly specified, so some interpretation is necessary.
3. Frequency response

4
To characterise linear time-invariant systems, we are interested in
their response to the complex exponential sequences ejωn for
different values of ω. Since Matlab deals with complex numbers as a
matter of course, we can generate portions of such sequences easily
using the exp command and a given value of ω.
A problem that arises is that such sequences take on non-zero
values for all n, both positive and negative. For computational
purposes some truncation is necessary, so instead we consider input
sequences of the form

x[n] = u[n]ejωn .

The sudden onset of the input at time n = 0 will cause transients,


but for a finite impulse response (FIR) filter with p taps you can
convince yourself that these will die within at most p samples.
(This should be clear from the expression for the convolution sum
in the previous question.)
(a) Use the exp function to generate the sequence x[n] above for
ω = π/10, over the range −10 to 70. Now calculate the output
sequence y[n] for the filter h[n] used in the previous question.
Convince yourself that the result (once the transients have died
out) is also a complex exponential sequence with frequency ω.
(One way to do this is to plot the real and imaginary parts of
y[n] — you may think of a more compelling argument.)
(b) Hand in: Based in the last result, we conclude that (once the
transients have died out) the output sequence can be written in
the form
y[n] = Aejφ ejωn ,


for some positive value of A and some value −π < φ ≤ π. Use


your results to estimate the corresponding values of A and φ for
the specific value of ω under investigation. (One way to do this

5
is to consider a value of n large enough so that the transients
have settled — now the ratio y[n]/x[n] should give you the
required complex value.)
(c) Repeat for ω = 2π/5 and ω = 4π/5. In total you should now
have three values that describe the frequency response of the
filter to three different frequencies. Convince yourself again that
the filter exhibits lowpass behaviour.
4. More frequency response
Formally, the response of the system to the input x[n] = ejωn is
given by

X ∞
X
y[n] = h[k]x[n − k] = h[k]ejω(n−k)
k=−∞ k=−∞
∞ ∞
!
X X
= h[k]ejωn e−jωk = h[k]e−jωk ejωn
k=−∞ k=−∞

= H(ejω )ejωn .

For each frequency ω the term in brackets is just a complex


number, designated H(ejω ), which specifies how the system changes
the magnitude and phase of the input complex exponential at the
specified frequency.
(The notation H(ejω ) above may seem strange to you, but makes
sense when you know about the z-transform — for now just think of
it as a function that takes a value ω and returns a complex scalar.)
(a) Hand in: Use the formula derived above, namely

X

H(e ) = h[k]e−jωk ,
k=−∞

to directly calculate the frequency response at ω = π/10,


ω = 2π/5, and ω = 4π/5. Note that h[n] is only nonzero over a

6
limited range, so the limits of the sum in this expression become
finite. Compare the results to the values obtained for A and φ
for all cases in the previous question.
(b) Hand in: If you didn’t already do so, write a Matlab function
that takes in a filter impulse response h[n] and a frequency ω,
and returns the value of H(ejω ) under the assumptions
described above. The function declaration may look as follows:
function H = freqresp(h, omega)
where h is a vector of values h = (h[1], h[2], . . . , h[p])T , ω is a
scalar frequency, and H is a scalar return value. Use this
function to find the frequency response H for a range of
frequencies ω = 0, 0.01, . . . , 4π. Plot the magnitude abs(H) as a
function of ω. On a separate set of axes, plot the phase
angle(H) as a function of ω. Interpret the results you obtain,
and note the periodicity of the frequency response. What effect
will the system have on an input complex exponential with
frequency ω = π/2?
(c) Do you see from the previous plots why the system is a bad
lowpass filter? Note the large sidelobes. Now plot the magnitude
and phase response (also over the range 0 to 4π) of the filter
with modified impulse response

1
 (n + 1) 0≤n≤3

 16
1
ht [n] = 16 (7 − n) 4≤n≤6


0 otherwise.

What do you notice?


5. Noise filtering (Optional and simple)
We don’t really deal in detail with random signals in this course,
but they form the basis of advanced DSP.
(a) Create a zero-mean white Gaussian noise (WGN) input sequence

7
with 100 elements, using the command
n = randn(100,1);
Plot the signal on a set of axes.
(b) Now filter the signal with one of the lowpass filters used before.
Plot the result on another set of axes, and compare your two
plots. The input noise signal has equal power at all frequencies,
while the output signal has a spectrum determined by the filter
frequency response.
(c) Clearly, if we have a low-frequency signal corrupted by noise, a
lowpass filter can increase our signal-to-noise ratio. Generate a
noisy signal
sn = sin((1:100)’/10) + randn(100,1);
and plot it both before and after filtering. Notice the reduction
in noise?

You might also like