Lab 1 - Introduction to MATLAB, Signals and Digital Signal Processing
Theory:
MATLAB stands for Matrix Laboratory. It is a multipurpose software suite for numerical computing. MATLAB allows matrix manipulation, plotting of functions and data, implementation of algorithms and can interface with other programming languages.
MATLAB is an interactive, matrix-based system for scientific and engineering numeric computation and visualization. It allows for the computation of complex numerical problems and their solutions in a fraction of the time required in other programming languages.
Lab:
1. Create an vector.
2. Create a matrix.
3. Create a 5x1 vector of zeros. Create a 1x5 vector of random numbers.
4. Transpose a matrix.
5. Compute the inner product of two matrices. Compute the cross product of two matrices. Compute the inverse of a matrix.
6. Compute the element wise multiplication of a matrix and a scalar.
7. Concatenate two matrices.
8. Create a vector of complex numbers.
9. Multiply a row of a matrix with an element of that same matrix.
10. Generate a vector of values ranging from 0 to 500 with 100 elements.
11. Create a 2D plot of the sine function between 0 and 2Π.
12. Create a 3D plot of a surface by creating a grid along the X and Y axes and plotting the Z-coordinate according to the exponential function.
13. Write a script to plot a vector of random data. Draw a horizontal line at the mean. Save the script and run it from the command line.
14. Write another script that calculates the mean of five samples of data from a vector of random data. Calculate the overall mean. Use a for loop to perform the calculations. For each iteration of the loop print out the intermediate results. Use an if..else control block to display the results depending on whether the mean of the samples is less than, greater than or equal to the overall mean.
15. Create a function that calculates the sum of an arbitrary number of sinusoidal terms. Call this function from the MATLAB command line or in a MATLAB script (.m).
Lab 2 - Discrete-Time Signals
Theory:
Signals are broadly classified into two categories, analog and discrete signals. An analog signal is continuous and has values at any point. A discrete signal is non- continuous and only has values at certain points. A discrete signal can be created from an analog signal by sampling at intervals. Typically, these intervals are equally spaced, but this is not a requirement.
Prelab Assignment:
Lab:
Generate the following sequences using the basic MATLAB signal functions and the basic M ATLAB signal operations discussed in this chapter. Plot signal samples using the stem function.
1. x1(n) = 3δ(n+2) + 2δ(n) - δ(n-3) + 5δ(n-7), -5 ≤ n ≤ 15
2. x4(n) = e0.1n[u(n+20) - u(n-10)]
3. x7(n) = e-0.05n sin(0.1Πn + Π/3), 0 ≤ n ≤ 100.
Generate the following random sequences and obtain their histogram using the hist function with 100 bins. Use the bar function to plot each histogram.
1. x1(n) is a random sequence whose samples are independent and uniformly distributed over [0, 2] interval. Generate 100,000 samples.
2. x2(n) is a Gaussian random sequence whose samples are independent with mean 10 and variance 10. Generate 10,000 samples.
3. x4(n) = k=1Σ4 yk(n) where each random sequence yk(n) is independent of others with samples uniformly distributed over [-0.5, 0.5]. Comment on the shape of this histogram.
Generate the following periodic sequences and plot their samples (using the stem function) over the indicated number of periods.
1. x1‾(n) = {...., -2, -1, 0, 1, 2,.....}periodic. Plot 5 periods.
2. x3‾(n) = sin(0.1Πn)[u(n) - u(n-10)] Plot 4 periods.
Let x(n) = {2, 4, -3, 1, -5, 4, 7}
↑
Generate and plot the samples (use the stem function) of the following sequences.
1. x1(n) = 2x(n-3) + 3x (n+4) - x(n)
2. x4(n) = 2e0.5nx(n) + cos(0.1Πn)x(n+2), -10 ≤ n ≤ 10
The complex exponential sequence ejω0n or the sinusoidal sequence cos(ω0n) are periodic if the normalized frequency f0 Δ= ω0/2Π is a rational number; that is, f0 = K/N, where K and N are integers.
1. Prove the above the result.
2. Generate exp(j0.1Πn), -100 ≤ n ≤ 100. Plot its real and imaginary parts using the stem function. Is this sequence periodic? If it is, what is its fundamental period? From the examination of the plot, what interpretation can you give to the integers K and N above?
3. Generate and plot cos(0.1n), -20 ≤ n ≤ 20. Is this sequence periodic? What do you conclude from the plot? If necessary, examine the values of the sequence in MATLAB to arrive at your answer.
Lab 3 - Fourier Series and Fourier Transform
Theory:
By sampling the discrete in time Fourier transform in the frequency domain we can obtain the discrete Fourier series (DFS). The DFS can then be extended to finite- duration sequences which results in a new transform, the discrete Fourier transform (DFT). The numerical computation of the DFT for long sequences can be prohibitively time consuming. As a result, algorithms were developed to compute the DFT efficiently and is collectively known as the fast Fourier transform (FFT).
Prelab Assignment:
Lab:
1. Compute the DFS coefficients of the following periodic sequences using the DFS definition and then verify your answers using MATLAB.
a. x1~(n) = {4,1, -1, 1} N = 4
b. x2~(n) = {2, 0, 0, 0, -1, 0, 0, 0}, N = 8
c. x3~(n) = {1, 0, -1, -1, 0}, N = 5
d. x4~(n) = {0, 0, 2 j, 0, 2j, 0), N = 6
e. x5~(n) = (3, 2, 1), N = 3
2. Let X(ejω) be the DTFT of a finite lenght sequence
a. Let
y1(n) = IDFS10-point[X(ej0), X(ej2Π/10), X(ej4Π/10),......,X(ej18Π/10) ]
Determine y1(n) using the frequency sampling theorem. Verify your answer using MATLAB.
b. Let
y2(n) = IDFS200-point[X(ej0), X(ej2Π/200), X(ej4Π/200),......,X(ej398Π/200) ]
Determine y2(n) using the frequency sampling theorem. Verify your answer using MATLAB.
c. Comment on your results in parts (a) and (b).
3. Consider the following finite-length sequence:
x(n) = { sinc2{(n - 50)/2}, 0 ≤ n ≤ 100;
{ 0, else.
a. Determine the DFT X(k) of x(n). Plot (using the stem function) its magnitude and phase.
b. Plot the magnitude and phase of the DTFT X(ejω) of x(n) using MATLAB.
c. Verify that the above DFT is the sampled version of X(ejω). It might be helpful to combine the above two plots in one graph using the hold function.
d. Is it possible to reconstruct the DTFT X(ejω) from the DFT X(k)? If possible, give the necessary interpolation formula for reconstruction. If not possible, state why this reconstruction cannot be done.
4. Let a finite-length sequence be given by:
x(n) = { 2e-0.9|n|, -5 ≤ n ≤ 5;
{ 0, otherwise.
Plot the DTFT X(ejω) of the above sequence using DFT as a computation tool. Choose the length N of the DFT so that this plot appears to be a smooth graph.
5. Plot the DTFT magnitude and angle of each of the following sequences using the DFT as a computation tool. Make an educated guess about the length N so that your plots are meaningful.
a. x1 (n) = (0.6)|n|[u(n + 10) - u(n - 10)].
b. x2(n) = n(0.9)n, 0 ≤ n ≤ 20.
c. x3(n) = cos(0.5Πn) + j sin(0.5Πn), 0 ≤ n ≤ 50.
d. x(n) = {1, 2, 3, 4, 3, 2, 1}.
↑
e. x(n) = (-1,-2,-3,0,3,2,1).
↑
6. A 512-point DFT X(k) of a real-valued sequence x(n) has the DFT values:
X(0) = 20 + jα; X(5) = 20 + j30; X(k1) = -10 + j15; x(152) = 17 + j23;
X(k2) = 20 - j30; X(k3) = 17 - j23; X(480) = -10 - j15; x(256) = 30 + jβ;
And all other values are known to be zero.
a. Determine the real-valued coefficients α and β.
b. Determine the values of the integers k1, k2 and k3.
c. Determine the energy of the signal x(n).
d. Express the sequence x(n) in a closed form.
7. Compute the N-point circular convolution for the following sequences. Plot their samples.
a. x1(n) = sin(Πn/3)R6(n), x2(n) = cos(Πn/4)R8(n); N = 10
b. x1(n) = cos(2Πn/3)RN(n), x2(n) = sin(2Πn/4)RN(n); N = 32
c. x1(n) = (0.8)nRN(n), x2(n) = (-0.8)nRN(n); N = 10
d. x1(n) = nRN(n), x2(n) = (N - n)RN(n); N = 10
e. x1(n) = (0.8)nR20x2(n), x2(n) = u(n) -u(n-40); N = 10
8. Let x(n) be a uniformly distributed random number between [-1,1] for 0 ≤ n ≤ 106.
Let
h(n) = sin(0.4Πn), 0 ≤ n ≤ 100
a. Consider the conv function, determine the output sequence y(n) = x(n)*h(n)
b. Consider the overlap and save method of block convolution along with the FFT algorithm to implement high speed block convolution. Using this approach, determine y(n) with FFT sizes of 1024, 2048 and 4096.
c. Compare these approaches in terms of the convolution results and their execution times.
Lab 4 - Frequency and Spectral Analysis
Theory:
Many sinusoidal signals encode information. These sinusoidal signals have a period or frequency. This is true of naturally occurring signals, as well as those that have been created by humans. Speech is a result of vibration of the human vocal cords; stars and planets change their brightness as they rotate on their axes and revolve around each other; ship's propellers generate periodic displacement of the water, and so on. The time-domain representation of these waveforms is often complicated and difficult to understand. Analysis of these signals is often easier to perform in the frequency domain where the frequency, phase and amplitude of the component sinusoids is easily obtained.
Prelab Assignment:
Lab:
1. Let x1(n) = {1, 2, 2, 1} A new sequence x2(n) is formed using
a. Express X2(ejω) in terms of X1(ejω) without explicitly computing X1(ejω)
b. Verify your results using MATLAB by computing and plotting magnitudes of the respective DTFTs.
2. Using the definition of the DTFT, determine the sequences corresponding to the following DTFTs.
a. X(ejω) = 3 + 2cos(ω) + 4cos(2ω)
b. X(ejω) = [1- 6cos(3ω) + 8cos(5ω)]e-j3ω]
c. X(ejω) = 2 + j4sin(2ω) - 5cos(4ω)
d. X(ejω) = [ 1+ 2cos(ω) + 3cos(2ω)]cos(ω/2)e-j5ω/2
e. X(ejω) = j[3 + 2cos(ω) + 4cos(2ω)]sin(ω)e-j3ω
3. Using the definition of the inverse DTFT, determine the sequences corresponding to the following DTFTs.
4. For each of the linear, shift-invariant systems described by the impulse response, determine the frequency response function H(ejω). Plot the magnitude response |H(ejω)| and the phase response, H(ejω) over the interval [-Π,Π].
a. h(n) = (0.9)|n|
b. h(n) = sinc(0.2n)[u(n+20) - u(n -20)], where sinc 0 = 1
c. h(n) = sinc(0.2n)[u(n) - u(n -40)]
d. h(n) = [(0.5)n + (0.4)n]u(n)
e.h(n) = (0.5)|n| cos(0.1Πn) = 1/2.0.5|n|ej.0.1Πn + 1/2.0.5|n|e-j0.1Πn
5. Determine H(ejω)and plot its magnitude and phase for each of the following systems.
a. y(n) = 1/5 m=0∑4 x(n-m)
b. y(n) = x(n) - x(n-2) + 0.95y(n -1) - 0.9025y(n-2)
c. y(n) = x(n) - x(n-1) + x(n - 2) + 0.95y(n-1) - 0.9025y(n-2)
d. y(n) = x(n) - 1.7678x(n-1) + 1.5625x(n - 2) + 1.1314y(n-1) - 0.64y(n-2)
e. y(n) = x(n) - l=1∑5 (0.5)ly(n-l)
Lab 5 - FIR Filters
Theory:
In signal processing, a finite impulse response (FIR) filter is a filter whose impulse response is of finite length duration, because it settles to zero in finite time. FIR filters have linear phase, are stable and are typically easy to design and implement.
Prelab Assignment:
Lab:
1. Design a bandpass filter using the Hamming window design technique. The specifications are
Lower stopband edge: 0.3Π
Upper stopband edge: 0.6Π As = 50 dB
Lower passband edge: 0.4Π
Upper passband edge: 0.5Π Rp = 0.5 dB
Plot the impulse response and the magnitude response (in dB) of the designed filter. Do not use the fir1 function.
2. Design a highpass filter using one of the fixed window functions. The specifications are
Stopband edge: 0.4Π, As = 50 dB
Passband edge: 0.6Π, Rp = 0.004 dB
Plot the zoomed magnitude response (in dB) of the designed filter in the passband to verify the passband ripple Rp. Do not use the fir1 function.
3. Design a linear-phase bandpass filter using the Hann window design technique. The specifications are
Lower stopband edge: 0.2Π
Upper stopband edge: 0.75Π As = 40 dB
Lower passband edge: 0.35Π
Upper passband edge: 0.55Π Rp = 0.25 dB
Plot the impulse response and the magnitude response (in dB) of the designed filter. Do not use the fir1 function.
4. Design a bandstop filter using the Hamming window design technique. The specifications are
Lower stopband edge: 0.4Π
Upper stopband edge: 0.6Π As = 50 dB
Lower passband edge: 0.3Π
Upper passband edge: 0.7Π Rp = 0.2 dB
Plot the impulse response and the magnitude response (in dB) of the designed filter. Do not use the fir1 function.
5. Using the Kaiser window method, design a linear-phase FIR digital filter that meets the following specifications:
0.975 ≤ |H(ejω)| ≤ 1.025, 0 ≤ ω ≤ 0.25Π
0 ≤ |H(ejω)| ≤ 0.005, 0.35Π ≤ ω ≤ 0.65Π
0.975 ≤ |H(ejω)| ≤ 1.025, 0.75Π ≤ ω ≤ Π
Determine the minimum-length impulse response h(n) of such a filter. Provide a plot containing subplots of the amplitude response and the magnitude response in dB. Do not use the fir1 function.
6. We wish to use the Kaiser window method to design a linear-phase FIR digital filter that meets the following specifications:
0 ≤ |H(ejω)| ≤ 0.01, 0 ≤ ω ≤ 0.25Π
0.95 ≤ |H(ejω)| ≤ 1.05, 0.35Π ≤ ω ≤ 0.65Π
0 ≤ |H(ejω)| ≤ 0.01, 0.75Π ≤ ω ≤ Π
Determine the minimum length impulse response h(n) of such a filter. Provide a plot containing subpltos of the amplitude response and the magnitude response in dB. Do not use the fir1 function.
Lab 6 - IIR Filters
Theory:
IIR filters are digital filters with infinite impulse responses, and hence are called infinite impulse response filters. Unlike FIR filters, IIR filters have a feedback mechanism or recursive part. Because of their infinite impulse response, IIR filters can be matched to analog filters which also have an infinite impulse response. Typically, IIR filters can meet a given specification with a much lower filter order than corresponding FIR filters, thus making them more efficient and requiring less calculations.
Prelab Assignment:
Lab:
1. A digital resonator is to be designed with ω0 = Π/4 that has two zeros at z = 0.
a. Compute and plot the frequency response of this resonator for r = 0.8, 0.9 and 0.99.
b. For each case in part 1, determine the 3 dB bandwidth and the resonant frequency ωr from your magnitude plots.
c. Check if your results in part 2 are in agreement with the theoretical results.
2. A digital resonator is to be designed with ω0 = Π/4 that has two zeros at z = 1 and z = -1.
a. Compute and plot the frequency response of this resonator for r = 0.8, 0.9, and 0.99.
b. For each case in part 1 determine the 3 dB bandwidth and the resonant frequency ωr from your magnitude plots.
3. We want to design a digital resonator with the following requirements: a 3 dB bandwidth of 0.05 rad, a resonant frequency of 0.375 cycles/sam, and zeros at z = 1 and z = -1. Using trial and error approach, determine the difference equation of the resonator.
4. A notch filter is to be designed with a null at the frequency ω0 = Π/2.
a. Compute and plot the frequency response of this notch filter for r = 0.7, 0.9 and 0.99.
b. For each case in part 1, determine the 3 dB bandwidth from your magnitude plots.
c. By trial-and-error approach, determine the value of r if we want the 3 dB bandwidth to be 0.04 radians at the null frequency ω0 = Π/2.
5. Repeat 4 for a null of ω0 = Π/6.
6. A speech signal with bandwidth of 4 kHz is samples at 8 kHz. The signal is corrupted by sinusoids with frequencies 1kHz, 2kHz and 3kHz.
a. Design an IIR filter using notch filter components that eliminates these sinusoidal signals.
b. Choose the gain of the filter so that the maximum gain is equal to 1, and plot the log-magnitude response of your filter.
c. Load the handel sound file in MATLAB and add the preceding three sinusoidal signals to create a corrupted sound signal. Now filter the corrupted sound signal using your filter and comment on its performance.
Lab 7 - ADC and DAC
Theory:
Analog to digital and digital to analog conversion are extremely important in DSP applications. The representation of an analog signal in computer memory requires analog to digital conversion, requiring sampling and quantization. The playback or production of an analog signal from computer memory is known as digital to analog conversion. This requires reconstruction of a signal from sampled points and involves some form of interpolation. Depending on how the digital signal is stored, the analog signal can either be reconstructed completely (lossless), or some information must be created to fill in the gaps (lossy).
Prelab Assignment:
Lab:
1. An analog signal xa(t) = sin(1000Πt) is sampled using the following sampling intervals. In each case, plot the spectrum of the resulting discrete-time signal.
a. Ts = 0.1 msec
b. Ts = 1 msec
c. Ts = 0.01 sec
2. Consider an analog signal xa(t) = cos(20Πt), 0 ≤ t ≤ 1. It is sampled at Ts = 0.01, 0.05, and 0.1 sec intervals to obtain x(n).
a. For each Ts, plot x(n).
b. Reconstruct the analog signal ya(t) from the samples x(n) using the sinc interpolation (use ?t = 0.001) and determine the frequency in ya(t) from your plot. (Ignore the end effects.)
c. 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. (Again, ignore the end effects.)
d. Comment on your results.
3. Consider an analog signal xa(t) = cos(20Πt + θ), 0 ≤ t ≤ 1. It is sampled at Ts = 0.5sec intervals to obtain x(n). Let θ = 0, Π/6, Π/4, Π/3, Π/2. For each of these θ values, perform the following.
a. Plot xa(t) and superimpose x(n) on it using the plot(n,x,'o') function.
b. Reconstruct the analog signal ya(t) from the samples x(n) using the sinc interpolation (Use ?t = 0.001) and superimpose x(n) on it.
c. Reconstruct the analog signal ya(t) from the samples x(n) using the cubic spline interpolation and superimpose x(n) on it.
d. You should observe that the resultant reconstruction in each case has the correct frequency but a different amplitude. Explain this observation. Comment on the role of phase of xa(t) on the sampling and reconstruction of signals.
4. Consider a signal with spectrum
X(ω) = { non-zero, |ω| ≤ ω0;
{ 0, ω0 < |ω| ≤ Π
a. Show that the signal x(n) can be recovered from its samples x(mD) if the sampling frequency ωs = 2Π/D ≥ 2ω0.
b. Sketch the spectra of x(n) and x(mD) for D = 4.
c. Show that x(n) can be reconstructed from the bandlimited interpolation
k=-∞∑∞x(kD) sinc[fc(n - kD)]; fc = 1/D
For ω0 < 2Πfc < ωs - ω0
5. Using the function interp, study the operation of factor-of-4 interpolation on the following sequences. Use the stem function to plot the original and the interpolated sequences. Experiment, using the filter length parameter values equal to 3 and 5. Comment on any differences in performance of the interpolation.
a. x1(n) = sin(0.6Πn), 0 ≤ n ≤ 100.
b. x2(n) = sin(0.8Πn) + cos(0.5Πn), 0 ≤ n ≤ 100.
c. x3(n) = 1 + cos(Πn), 0 ≤ n ≤ 100.
d. x4(n) = 0.2n, 0 ≤ n ≤ 100.
e. x5(n) = {1, 1, 1, 1, 0, 0, 0, 0, 0, 0}PERIODIC, 0 ≤ n ≤ 100.
Lab 8 - Audio Processing
Theory:
Some interesting effects can be achieved by delaying audio samples. This introduces an echo effect to the audio samples. In listening spaces, the soundwaves reflect off surfaces and arrive at our ears with differing amounts of delay and gain. A single reflection or echo of a signal can be implemented by the following filter, which adds to the direct signal an attenuated and delayed copy of itself:
y(n) = x(n) + ax(n-D)
The following algorithm can accomplish the echo effect:
int delayMilliseconds = 500; // half a second
int delaySamples =
(int)((float)delayMilliseconds * 44.1f); // assumes 44100 Hz sample rate
float decay = 0.5f;
for (int i = 0; i < buffer.length - delaySamples; i++)
{
buffer[i + delaySamples] += (short)((float)buffer[i] * decay);
}
Lab:
Write an application to read a wav file into an array and introduce a delay. Reduce the gain of the delayed sample and feed it back into the input to create an echo effect.
Lab 9 - Audio Processing 2
Theory:
Flanging is an audio effect used in recording studios that can add a whoosing sound not unlike a jet. By allowing the delay time to vary, a flanging or chorus effect can be added.
y(n) = X(n) + ax(n-d(n))
Allowing the delay to vary periodically, a flanging effect is accomplished.
Lab:
Write an application to read a wav file into an array and introduce a periodic delay to create a flanging effect. The periodic delay can be accomplished by making the delay value modulated using the modulo operator %, or it can be a value retrieved from the sine function.
Lab 10 - Image Processing
Theory:
Digital images are prone to a variety of types of noise. Noise is the result of errors in the image acquisition process that result in pixel values that do not reflect the true intensities of the real scene. You can use linear filtering to remove certain types of noise. Certain filters, such as averaging or Gaussian filters, are appropriate for this purpose. For example, an averaging filter is useful for removing grain noise from a photograph. Because each pixel gets set to the average of the pixels in its neighborhood, local variations caused by grain are reduced.
Lab:
Write an application to read in a noisy image, filter the image and output the image to a new file for comparison with the original normal image.