Problems
Write programs (using Matlab or other software) to finish the exercises below.
Plot this signal and its frequency spectrum;
When the sampling period satisfies $T =1$，$T =π /2$，$T = 2$, respectively, please plot the sampling signal $f_p (n)$ and its frequency spectrum, respectively. Please give explanation of
these results;
Using lowpass filter with cutting frequency $ω_c =2.4$ to reconstruct signal $f_r(t)$ from
$f_p (n)$ . When the sampling period satisfies $T = 1$ ， $T = 2$ , respectively, please plot the reconstructed signal $f_r(t)$, and plot the absolute error between the reconstructed signal
$f_r (t)$ and the original signal $f (t)$. Please analyze these results.
Matlab functions potentially used:
plot;subplot;axis;exp;cos;sinc;ones;length;stem;abs
Original Signal & Frequency Spectrum
Experiment
- Sampling period $T = 0.1$
-
X label: $0\le t \le10$ - Highest frequency = 5 times the frequency of the reference wave
- Fourier transform:
Results
Implementation
% Initialization
T = 2 * pi;
f = 1 / T;
dt = 0.1;
x0 = -10;
xt = 10;
% Time series of discrete signals
n1 = x0 / dt : -pi / dt;
num_n1 = length(n1);
n2 = -pi / dt : pi / dt;
num_n2 = length(n2);
n3 = pi / dt : xt / dt;
num_n3 = length(n3);
% Original Signal
tn = [n1 * dt n2 * dt n3 * dt];
y1 = zeros(1,num_n1);
y2 = 0.5 * (1 + cos(2 * f * pi * n2 * dt));
y3 = zeros(1,num_n3);
y = [y1 y2 y3];
subplot(1,2,1);
plot(tn,y,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Original Signal']);
% Frequency Spectrum
N = length(tn);
fm = 5 * f;
Tm = 1 / fm;
wm = 2 * pi * fm;
k = 0 : N - 1;
w = k * wm / N;
F = y * exp(-1j*tn'*w)*dt;
subplot(1,2,2);
plot(w/(2*pi),abs(real(F)),'LineWidth',3);
axis([0 4*fm 1.2*min(abs(F)) 1.2*max(abs(F))]);
title(['Frequency Spectrum']);
figure;
Sampling Signals & Frequency Spectrums
Experiment
- Sampling period $T = 1, \pi/2, 2$
-
X label: $0\le t \le10$ - Highest frequency = 5 times the frequency of the reference wave
- Fourier transform:
Results
Insights
Suppose the highest frequency component for a given analog signal is $f_{max}$. According to the Nyquist Theorem, the sampling rate must be at least $2f_{max}$. The sampling in an analog-to-digital converter is actuated by a pulse generator (clock). If the sampling rate is less than $2f_{max}$, some of the highest frequency components in the analog input signal will not be correctly represented in the digitized output. When such a digital signal is converted back to analog form by a digital-to-analog converter, false frequency components appear that were not in the original analog signal. This undesirable condition is a form of distortion called aliasing.
Insight 1. In this experiment, $f_{max}=0.5$, $2\cdot f_{max}=1$. We can clearly see that when $T=2$, the highest frequency components in the analog input signal will not be correctly represented in the digitized output.
Insight 2. The spectrum of the sampled signal is the superposition of the shifted spectrum of the original signal. The spectrum of the sampled signal is the spectrum of the original analog signal along the frequency axis, and the sampling angular frequency (frequency of unit impulse signal) is repeated once at every interval, and the periodic function is formed by superposition. In other words, the spectrum of the ideal sampling signal is the spectrum of the original analog signal with period ($W_s$) and periodic extension.
Implementation
Ts = [1 pi/2 2];
fs = 1 / Ts(i);
% Time series of discrete signals
ns1 = x0 / Ts(i) : -pi / Ts(i);
num_ns1 = length(ns1);
ns2 = -pi / Ts(i) : pi / Ts(i);
num_ns2 = length(ns2);
ns3 = pi / Ts(i) : xt / Ts(i);
num_ns3 = length(ns3);
% Sampling Signal
tns = [ns1 * Ts(i) ns2 * Ts(i) ns3 * Ts(i)];
ts = 2 * f * pi * ns2 * Ts(i);
ys1 = zeros(1,num_ns1);
ys2 = 0.5 * (1 + cos(ts));
ys3 = zeros(1,num_ns3);
ys = [ys1 ys2 ys3];
subplot(5,2,i * 2 - 1);
plot(tns,ys,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Sampling Signal T = ',num2str(Ts(i))]);
% Sampling Frequency Spectrum
N = length(tns);
wm = 2 * pi * fs;
k = 0 : N - 1;
ws = k * wm / N;
Fs = fft(ys);
subplot(5,2,i * 2);
plot(ws/(2*pi),abs(Fs),'LineWidth',3);
axis([0 4*fm 1.2*min(abs(Fs)) 1.2*max(abs(Fs))]);
title(['Frequency Spectrum T = ',num2str(Ts(i))]);
Reconstructed Signals & Absolute Errors
Experiment
- Sampling period $T = 1, 2$
-
X label: $0\le t \le10$ - Highest frequency = 5 times the frequency of the reference wave
- Fourier transform:
Results
Insights
Insight 3. Since aliasing is not so accurate, the integrity of the new signal needs to be reconstructed, so there will be some errors between the original signal and the new signal.
Implementation
wc = 2.4;
if (i == 1) || (i == 3)
for j = 1 : N
if ws(j) > wc
Fs(j) = 0;
end
end
% Inverse Fourier transform
yr = ifft(Fs);
% Reconstructed Signal
subplot(5,2,6 + i);
plot(tns,yr,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Reconstructed Signal T = ',num2str(Ts(i))]);
% Absolute Error
err = yr - ys;
subplot(5,2,6 + i + 1);
stem(tns,err,'filled');
axis([x0 xt -1.5 1.5]);
title(['Absolute Error T = ',num2str(Ts(i))]);
end
Whole Implementation
% Initialization
T = 2 * pi;
f = 1 / T;
dt = 0.1;
x0 = -10;
xt = 10;
% Time series of discrete signals
n1 = x0 / dt : -pi / dt;
num_n1 = length(n1);
n2 = -pi / dt : pi / dt;
num_n2 = length(n2);
n3 = pi / dt : xt / dt;
num_n3 = length(n3);
% Original Signal
tn = [n1 * dt n2 * dt n3 * dt];
y1 = zeros(1,num_n1);
y2 = 0.5 * (1 + cos(2 * f * pi * n2 * dt));
y3 = zeros(1,num_n3);
y = [y1 y2 y3];
subplot(1,2,1);
plot(tn,y,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Original Signal']);
% Frequency Spectrum
N = length(tn);
fm = 5 * f;
Tm = 1 / fm;
wm = 2 * pi * fm;
k = 0 : N - 1;
w = k * wm / N;
F = y * exp(-1j*tn'*w)*dt;
subplot(1,2,2);
plot(w/(2*pi),abs(real(F)),'LineWidth',3);
axis([0 4*fm 1.2*min(abs(F)) 1.2*max(abs(F))]);
title(['Frequency Spectrum']);
figure;
for i = 1 : 3
Ts = [1 pi/2 2];
fs = 1 / Ts(i);
% Time series of discrete signals
ns1 = x0 / Ts(i) : -pi / Ts(i);
num_ns1 = length(ns1);
ns2 = -pi / Ts(i) : pi / Ts(i);
num_ns2 = length(ns2);
ns3 = pi / Ts(i) : xt / Ts(i);
num_ns3 = length(ns3);
% Sampling Signal
tns = [ns1 * Ts(i) ns2 * Ts(i) ns3 * Ts(i)];
ts = 2 * f * pi * ns2 * Ts(i);
ys1 = zeros(1,num_ns1);
ys2 = 0.5 * (1 + cos(ts));
ys3 = zeros(1,num_ns3);
ys = [ys1 ys2 ys3];
subplot(5,2,i * 2 - 1);
plot(tns,ys,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Sampling Signal T = ',num2str(Ts(i))]);
% Sampling Frequency Spectrum
N = length(tns);
wm = 2 * pi * fs;
k = 0 : N - 1;
ws = k * wm / N;
Fs = fft(ys);
subplot(5,2,i * 2);
plot(ws/(2*pi),abs(Fs),'LineWidth',3);
axis([0 4*fm 1.2*min(abs(Fs)) 1.2*max(abs(Fs))]);
title(['Frequency Spectrum T = ',num2str(Ts(i))]);
wc = 2.4;
if (i == 1) || (i == 3)
for j = 1 : N
if ws(j) > wc
Fs(j) = 0;
end
end
% Inverse Fourier transform
yr = ifft(Fs);
% Reconstructed Signal
subplot(5,2,6 + i);
plot(tns,yr,'LineWidth',3);
axis([x0 xt -0.5 1.5]);
title(['Reconstructed Signal T = ',num2str(Ts(i))]);
% Absolute Error
err = yr - ys;
subplot(5,2,6 + i + 1);
stem(tns,err,'filled');
axis([x0 xt -1.5 1.5]);
title(['Absolute Error T = ',num2str(Ts(i))]);
end
end
References
- Nyquist–Shannon sampling theorem
- Continuous Fourier Transform
- Inverse Fourier transform
- Matlab Tutorial
