ELEC6218 · Signal Processing

Lecture 2

The Transform Toolkit

This lecture connects time-domain equations with frequency-domain tools. You will use the Laplace transform to obtain transfer functions, poles and stability information. You will then see how the CTFT, CTFS, DTFT and DTFS describe different signal types, and how the DFT and FFT make frequency analysis practical on a computer.

Quiz

Podcast

InstructorProf. Abdulkadir Çelik (https://akadircelik.com) ModuleELEC6218 TopicsLaplace · poles & stability · CTFT · CTFS · DTFT · DTFS · DFT & FFT

▸ Laplace → Fourier ▸ Complex exponentials ▸ The Fourier ladder ▸ DFT & FFT

Learning outcomes

By the end of this lecture, you will be able to:

• Derive \(H(s)\) from a differential equation, invert it using partial fractions, and assess BIBO stability from its poles.
• Explain how the CTFT is obtained by evaluating the Laplace transform on the \(j\Omega\)-axis and why complex exponentials are eigenfunctions of LTI systems.
• Choose the correct transform—CTFT, CTFS, DTFT or DTFS—for a signal's time-domain properties.
• Use transform properties, calculate a Fourier series and a small DFT, and explain why the FFT is much faster than direct DFT evaluation.

Lecture 2

Contents

Sections 2.1–2.5 develop the Laplace toolkit: convergence, inversion, poles and stability. Sections 2.6–2.10 connect Laplace analysis to Fourier analysis and filtering. Sections 2.11–2.14 then introduce the four Fourier representations and finish with the DFT and FFT.

iNotation used in this lecture

\(H(s)\)

Laplace transfer function, \(H(s) = Y(s)/X(s)\)

\(s = \sigma + j\Omega\)

the complex Laplace variable

\(p_k\), \(z_k\)

poles and zeros of \(H(s)\)

\(X(\Omega)\)

CTFT of \(x(t)\); continuous spectrum, \(\Omega\) in rad/s

\(X_k\)

CTFS (or DTFS) coefficient of the \(k\)th harmonic

\(\Omega_0 = 2\pi/T_0\)

fundamental (radian) frequency of a periodic signal

\(X(e^{j\omega})\)

DTFT of \(x[n]\); continuous in \(\omega\), \(2\pi\)-periodic, \(\omega\) in rad/sample

\(X[k]\)

\(N\)-point DFT bin \(k\), \(k=0,\ldots,N-1\)

\(W_N = e^{-j2\pi/N}\)

the DFT twiddle factor; \([\mathbf{W}]_{k,n}=W_N^{kn}\)

Keep the frequency symbols separate. We use \(\Omega\) for continuous-time frequency in radians per second and \(\omega\) for discrete-time frequency in radians per sample. Therefore, \(X(\Omega)\) and \(H(j\Omega)\) describe continuous-time signals, while \(X(e^{j\omega})\) and \(H(e^{j\omega})\) describe sequences. Only one interval of length \(2\pi\), usually \((-\pi,\pi]\), is needed in discrete time. Lecture 3 connects the variables through \(\omega=\Omega T=\Omega/f_s\).

2.1  Differential LTI systems and the Laplace transform

Circuits, mechanical systems and control loops are often modelled by a linear constant-coefficient differential equation (LCCDE):

$$\sum_{k=0}^{N} a_k \frac{d^k}{dt^k}y(t) \;=\; \sum_{m=0}^{M} b_m \frac{d^m}{dt^m}x(t),$$

(2.1)

Here, \(N\) is the system order and the coefficients \(\{a_k\}\) and \(\{b_m\}\) are constant. The equation and the impulse response \(h(t)\) are two equivalent descriptions of the same LTI system. The total output contains two parts:

$$ y(t) = y_{\text{zero-input}}(t) + y_{\text{zero-state}}(t). $$

(2.2)

The zero-input response is caused by stored energy or initial conditions when \(x(t)=0\). The zero-state response is caused by the input when all initial conditions are zero. For an LTI system, the zero-state response is the convolution of \(x(t)\) with \(h(t)\).

2.1.1  The Laplace transform

The Laplace transform (LT) converts a causal time-domain function \(f(t)\) into a function of a complex variable \(s=\sigma+j\Omega\):

$$\mathcal{L}\{f(t)\} = F(s) = \int_{0}^{\infty} f(t)\,e^{-st}\,dt, \qquad \mathcal{L}^{-1}\{F(s)\} = f(t) = \frac{1}{2\pi j}\int_{\sigma-j\infty}^{\sigma+j\infty} F(s)\,e^{st}\,ds.$$

(2.3)

In the kernel \(e^{-st}=e^{-\sigma t}e^{-j\Omega t}\), the real part \(\sigma\) controls exponential decay or growth (convergence) while \(\Omega\) controls oscillation. In practice the inverse Bromwich integral is never evaluated directly: we use tables of Laplace pairs together with partial fraction expansion (Section 2.3).

Why leave the time domain at all?

Because in the time domain, cascading systems means convolving impulse responses — a messy integral every time. Transform domains trade that integral for ordinary multiplication of algebraic functions.

Case study. An RC low-pass filter obeys \(RC\,y'(t)+y(t)=x(t)\). Solved directly, a new differential equation must be re-solved for every input \(x(t)\). Laplace-transform both sides instead and the derivative becomes multiplication by \(s\): \((RCs+1)Y(s)=X(s)\), so \(H(s)=Y(s)/X(s)=1/(1+RCs)\) — one algebraic expression that describes the filter for every input at once. That is the whole motivation for this lecture: turn calculus into algebra, then read the system's behaviour straight off the resulting expression.

2.1.2  The Laplace transform of a convolution

The reason the Laplace transform is worth the extra abstraction is that it turns convolution into ordinary multiplication. For causal \(x,h\), with \(Y(s)=\mathcal{L}\{x(t)*h(t)\}\):

Step 1. Write out the double integral, noting the integration region is the triangle \(0\le\tau\le t<\infty\):

$$Y(s) = \int_{0}^{\infty} e^{-st}\left(\int_{0}^{t}x(\tau)h(t-\tau)\,d\tau\right)dt.$$

Step 2. Swap the order of integration (Fubini's theorem), so \(\tau\) runs from \(0\) to \(\infty\) and, for each \(\tau\), \(t\) runs from \(\tau\) to \(\infty\):

$$Y(s) = \int_{0}^{\infty}\int_{\tau}^{\infty} e^{-st}x(\tau)h(t-\tau)\,dt\,d\tau.$$

Step 3. Substitute \(\alpha=t-\tau\) (so \(t=\alpha+\tau\), \(dt=d\alpha\)), which decouples the two integrals:

$$Y(s) = \int_{0}^{\infty}\int_{0}^{\infty} e^{-s(\alpha+\tau)}x(\tau)h(\alpha)\,d\alpha\,d\tau = \left(\int_{0}^{\infty}e^{-s\tau}x(\tau)\,d\tau\right)\left(\int_{0}^{\infty}e^{-s\alpha}h(\alpha)\,d\alpha\right).$$

Step 4. Each factor is exactly a Laplace transform:

$$\boxed{\;\mathcal{L}\{x(t)*h(t)\} = X(s)H(s) = Y(s).\;}$$

(2.4)

Figure 2.1.a shows this result as a commutative diagram: it does not matter whether you convolve in the time domain then transform, or transform first and multiply, the two paths agree.

2.1.3  From LCCDE to transfer function

Apply the Laplace transform to the LCCDE, using the differentiation property \(\mathcal{L}\{d^rf/dt^r\} = s^rF(s) - s^{r-1}f(0^-) - \cdots - f^{(r-1)}(0^-)\). Assuming zero initial conditions (\(x^{(r)}(0^-)=y^{(r)}(0^-)=0\) for all \(r\), the standard assumption for system analysis) makes every initial-condition term vanish, converting the differential equation into a purely algebraic one:

$$\big(a_0+a_1s+\cdots+a_Ns^N\big)Y(s) = \big(b_0+b_1s+\cdots+b_Ms^M\big)X(s)$$

$$\boxed{H(s) = \frac{Y(s)}{X(s)} = \frac{b_0+b_1s+\cdots+b_Ms^M}{a_0+a_1s+\cdots+a_Ns^N} = K\,\frac{\prod_{k=1}^{M}(s-z_k)}{\prod_{\ell=1}^{N}(s-p_\ell)}.}$$

(2.5)

Factoring the numerator and denominator into roots exposes the zeros \(z_k\) (where \(H(s)=0\)) and poles \(p_\ell\) (where \(H(s)\to\infty\)); \(K\) is the overall gain. Zeros shape spectral notches; poles, as Section 2.4 shows, govern stability. A pole–zero plot (poles marked \(\times\), zeros marked \(\circ\), in the complex \(s\)-plane) is the single most useful sketch of a system, and is exactly the picture Lecture 5 draws over and over while designing filters.

MATLAB corner

MATLAB stores the numerator and denominator of Eq. (2.5) as two coefficient vectors. The tf() command turns them into a transfer-function model, and pole() and zero() find the roots for us. This is useful for checking hand calculations and for seeing stability directly on a pole–zero map.

Goal: build \(H(s)\) from the LCCDE's coefficients and read its poles/zeros straight off the object, instead of factoring by hand.

tf_from_lccde.m

MATLABCopy{{ l2ShowCodeTfLabel }}

% tf_from_lccde.m - transfer function from LCCDE, worked example with
% complex-conjugate poles  (Eq. 2.5, Examples 2.6 and 2.7)
num = [5 25]; den = [1 4 29];
H = tf(num,den);
disp('H(s) ='); disp(H)
p = pole(H); z = zero(H); disp('Poles:'); disp(p); disp('Zeros:'); disp(z)
pzmap(H); grid on; hold on
hp = plot(real(p), imag(p), 'rx', 'MarkerSize',14, 'LineWidth',2);    % poles
hz = plot(real(z), imag(z), 'go', 'MarkerSize',12, 'LineWidth',2);    % zeros
ylim([-6 6])
legend([hp hz], {'Poles','Zeros'}, 'Location','best')                 % legend for markers
% Poles at Re(s) = -2 < 0 (open LHP) => system is stable

2.2  The region of convergence: the half of the answer the algebra hides

A Laplace transform is more than an algebraic expression in \(s\). Its defining integral may converge for some values of \(s\) and diverge for others. The set of values for which it converges is the region of convergence (ROC). The pair “expression plus ROC” is needed to identify the time-domain signal uniquely.

The bilateral definition includes all time and therefore makes no assumption about what happened before \(t=0\):

$$F(s)=\int_{-\infty}^{\infty} f(t)\,e^{-st}\,dt, \qquad s=\sigma+j\Omega,$$

(2.6)

Because \(|e^{-st}|=e^{-\sigma t}\), convergence depends on \(\sigma=\operatorname{Re}\{s\}\), not on \(\Omega\). An ROC therefore appears as a vertical half-plane or strip in the \(s\)-plane. Its boundaries pass through poles, which are never included in the ROC.

2.2.1  Two signals, one algebraic expression

Transform the right-sided exponential \(f(t)=e^{-at}u(t)\). The integral runs from \(0\) to \(\infty\) and converges only when the integrand decays, i.e. when \(\sigma+a>0\):

$$e^{-at}u(t) \;\longleftrightarrow\; \frac{1}{s+a}, \qquad \text{ROC: } \sigma > -a.$$

(2.7)

Now transform the left-sided signal \(g(t)=-e^{-at}u(-t)\), which is zero for \(t>0\). The integral runs from \(-\infty\) to \(0\), and converges under the opposite condition:

$$-e^{-at}u(-t) \;\longleftrightarrow\; \frac{1}{s+a}, \qquad \text{ROC: } \sigma < -a.$$

(2.8)

Both signals produce the same expression, \(1/(s+a)\). The ROC tells us which signal is meant. An ROC to the right of the pole gives the right-sided, causal signal. An ROC to the left gives the left-sided, anti-causal signal. Lecture 4 uses the same idea for the \(z\)-transform.

Three ROC rules. First, an ROC never contains a pole. Second, a right-sided signal has an ROC to the right of its rightmost pole; a left-sided signal has an ROC to the left of its leftmost pole; and a two-sided signal has a strip between poles. Third, an absolutely integrable finite-duration signal has an ROC covering the whole \(s\)-plane.

EXAMPLE 2.1 · ROC decides the signal

A system has \(F(s)=\dfrac{1}{(s+1)(s-2)}\). Write down every signal \(f(t)\) consistent with this expression, and state which one is causal and which one is stable.

{{ l2Ex12Label }}

Step 1. Partial fractions give \(F(s)=\dfrac{-1/3}{s+1}+\dfrac{1/3}{s-2}\). The poles are at \(s=-1\) and \(s=+2\), so the plane is cut into three vertical strips, and each is a legitimate ROC.

Step 2. ROC \(\sigma>2\) (right of both poles): both terms are right-sided.

$$f(t) = \tfrac{1}{3}\big(e^{2t}-e^{-t}\big)u(t).$$

Step 3. ROC \(\sigma<-1\) (left of both poles): both terms are left-sided.

$$f(t) = \tfrac{1}{3}\big(e^{-t}-e^{2t}\big)u(-t).$$

Step 4. ROC \(-1<\sigma<2\) (the strip between them): the pole to the left contributes a right-sided term, the pole to the right a left-sided one.

$$f(t) = -\tfrac{1}{3}e^{-t}u(t) - \tfrac{1}{3}e^{2t}u(-t).$$

Step 5. The ROC \(\sigma>2\) is to the right of the rightmost pole, so its signal is causal. However, that ROC does not include the \(j\Omega\)-axis, so the signal is not stable. The middle strip \(-1<\sigma<2\) includes the \(j\Omega\)-axis, so it is stable, but it is two-sided and not causal. Because one pole lies in the right-half plane, no ROC can make this system both causal and stable.

2.2.2  Unilateral and bilateral: which one you are actually using

Many engineering problems begin at \(t=0\) and ask for the response afterwards. The unilateral Laplace transform is designed for this situation:

$$F(s)=\int_{0^-}^{\infty} f(t)\,e^{-st}\,dt.$$

(2.9)

The lower limit \(0^-\) includes an impulse at the origin, so \(\mathcal{L}\{\delta(t)\}=1\). A unilateral transform represents a right-sided signal, so its ROC is automatically to the right of the rightmost pole and is often omitted. Use the bilateral transform when behaviour before \(t=0\) matters or when the signal is two-sided.

The main advantage of the unilateral form is that it includes initial conditions directly. Its differentiation rules allow differential equations with a nonzero starting state to be solved algebraically:

$$\mathcal{L}\left\{\frac{df}{dt}\right\} = sF(s)-f(0^-), \qquad \mathcal{L}\left\{\frac{d^2f}{dt^2}\right\} = s^2F(s)-sf(0^-)-f'(0^-).$$

(2.10)

2.2.3  Initial and final value theorems

The initial and final values of a signal can sometimes be found directly from \(F(s)\), without calculating the full inverse transform:

$$f(0^+) = \lim_{s\to\infty} sF(s), \qquad \lim_{t\to\infty} f(t) = \lim_{s\to 0} sF(s).$$

(2.11)

For a strictly proper rational \(F(s)\), the initial value theorem is straightforward to use. The final value theorem requires an extra check: every pole of \(sF(s)\) must lie strictly in the left-half plane. Otherwise, the signal grows or continues to oscillate and has no final value. Always check the poles before using the second limit.

EXAMPLE 2.2 · Final value, applied and misapplied

For the step response of a stable first-order system, \(F(s)=\dfrac{5}{s(s+2)}\), find \(f(0^+)\) and \(f(\infty)\). Then try the same on \(G(s)=\dfrac{5}{s(s^2+9)}\) and say what goes wrong.

{{ l2Ex13Label }}

Step 1. Initial value: \(sF(s)=\dfrac{5}{s+2}\to 0\) as \(s\to\infty\), so \(f(0^+)=0\). The step response starts from rest, as it must.

Step 2. Final value: the poles of \(sF(s)=5/(s+2)\) are at \(s=-2\), strictly in the left-half plane, so the theorem applies:

$$\lim_{t\to\infty} f(t) = \lim_{s\to0}\frac{5}{s+2} = 2.5.$$

Inverting confirms it: \(f(t)=2.5\,(1-e^{-2t})u(t)\), which settles at \(2.5\).

Step 3. Now \(G(s)\). Mechanically, \(sG(s)=\dfrac{5}{s^2+9}\to 5/9\) as \(s\to0\) — an answer, and a plausible-looking one.

Step 4. But \(sG(s)\) has poles at \(s=\pm j3\), on the \(j\Omega\)-axis, not strictly left of it. Inverting gives \(g(t)=\tfrac{5}{9}\big(1-\cos 3t\big)u(t)\), which oscillates forever between \(0\) and \(10/9\) and has no limit. The theorem's precondition failed, so its output was meaningless. Always check the poles of \(sF(s)\) before quoting a final value.

2.3  Partial fractions: inverting the Laplace transform

Most transfer functions are ratios of polynomials in \(s\). Direct inversion may be difficult, so we first split the expression into simpler fractions. Each simple term then matches a standard Laplace-transform pair.

Distinct real poles. If \(F(s)=P(s)/Q(s)\) has \(N\) distinct poles \(p_1,\ldots,p_N\),

$$F(s) = \sum_{k=1}^{N}\frac{A_k}{s-p_k} \quad\Longrightarrow\quad f(t) = \sum_{k=1}^{N} A_k\,e^{p_kt}\,u(t).$$

(2.12)

Each pole \(p_k\) contributes one exponential term \(e^{p_kt}\): this is the key idea behind the whole method. Repeated poles use the table pair below, one term per power of the repeated factor:

$$ \mathcal{L}^{-1}\{1/(s+a)^n\} = \frac{t^{n-1}}{(n-1)!}e^{-at}u(t). $$

(2.13)

EXAMPLE 2.3 · Distinct real poles

Find \(f(t)=\mathcal{L}^{-1}\{F(s)\}\) for \(F(s) = \dfrac{4s+8}{s^2+5s+6}\).

{{ l2Ex2Label }}

Step 1. Factor: \(s^2+5s+6=(s+2)(s+3)\). Write \(F(s)=\dfrac{A}{s+2}+\dfrac{B}{s+3}\), so \(4s+8=A(s+3)+B(s+2)\).

Step 2. At \(s=-2\): \(0=A(1)\Rightarrow A=0\). At \(s=-3\): \(-4=B(-1)\Rightarrow B=4\).

$$F(s)=\frac{4}{s+3} \quad\Longrightarrow\quad \boxed{f(t)=4e^{-3t}u(t).}$$

EXAMPLE 2.4 · Repeated pole

Find \(f(t)\) for \(F(s)=\dfrac{3}{(s+1)^2}\).

{{ l2Ex3Label }}

Solution. Direct table lookup, \(\mathcal{L}^{-1}\{1/(s+a)^n\}=\frac{t^{n-1}}{(n-1)!}e^{-at}u(t)\) with \(n=2,a=1\):

$$\boxed{f(t) = 3t\,e^{-t}u(t).}$$

EXAMPLE 2.5 · Higher-order repeated pole

Find \(f(t)\) for \(F(s)=\dfrac{2s+5}{(s+2)^3}\).

{{ l2Ex4Label }}

Solution. Expand \(F(s)=\dfrac{A}{s+2}+\dfrac{B}{(s+2)^2}+\dfrac{C}{(s+2)^3}\); solving gives \(A=0,\ B=2,\ C=1\), so \(F(s)=\dfrac{2}{(s+2)^2}+\dfrac{1}{(s+2)^3}\). Inverting each term:

$$\boxed{f(t) = e^{-2t}\Big(2t+\tfrac{t^2}{2}\Big)u(t).}$$

2.3.1  Complex-conjugate poles

A complex-conjugate pole pair \(p_{1,2}=-\alpha\pm j\beta\) produces an oscillatory response. This occurs frequently in circuits, mechanical systems and control loops. A corresponding second-order term has the form \(H(s)=\dfrac{As+B}{(s+\alpha)^2+\beta^2}\). Two standard transform pairs are enough to invert it:

$$\mathcal{L}^{-1}\!\left\{\frac{s+\alpha}{(s+\alpha)^2+\beta^2}\right\} = e^{-\alpha t}\cos(\beta t), \qquad \mathcal{L}^{-1}\!\left\{\frac{\beta}{(s+\alpha)^2+\beta^2}\right\} = e^{-\alpha t}\sin(\beta t).$$

(2.14)

Step 1. Rewrite the numerator as a combination of \((s+\alpha)\) and \(\beta\): \(As+B = A(s+\alpha) + (B-A\alpha)\); writing \(C=(B-A\alpha)/\beta\) gives \(As+B = A(s+\alpha)+C\beta\), so

$$H(s) = A\,\frac{s+\alpha}{(s+\alpha)^2+\beta^2} + C\,\frac{\beta}{(s+\alpha)^2+\beta^2}.$$

Step 2. Apply the two standard pairs term by term:

$$\boxed{h(t) = e^{-\alpha t}\big[A\cos(\beta t)+C\sin(\beta t)\big]u(t) = D\,e^{-\alpha t}\sin(\beta t+\psi)\,u(t),}$$

(2.15)

where the amplitude–phase form uses \(D=\sqrt{A^2+C^2}\) and \(\tan\psi = A/C\) (both forms decay at the same rate \(e^{-\alpha t}\), only the bookkeeping differs).

EXAMPLE 2.6 · Complex-conjugate poles

Find \(f(t)\) for \(F(s) = \dfrac{5s+25}{s^2+4s+29}\).

{{ l2Ex5Label }}

Step 1. Complete the square: \(s^2+4s+29=(s+2)^2+25\), so \(\alpha=2,\ \beta=5\). Write the numerator as \(5s+25 = 5(s+2)+15\), so

$$F(s) = 5\,\frac{s+2}{(s+2)^2+25} + 15\,\frac{1}{(s+2)^2+25} = 5\,\frac{s+2}{(s+2)^2+25} + 3\,\frac{5}{(s+2)^2+25},$$

Step 2. matching the standard pairs with \(A=5\) and \(C=3\) (not \(15\): the second term's coefficient must be re-expressed in units of \(\beta=5\) before the pair applies). Then:

$$f(t) = e^{-2t}\big[5\cos(5t) + 3\sin(5t)\big]u(t).$$

Step 3. Combine the sine and cosine into one sinusoid. Expand \(R\sin(\beta t+\phi)\) and match coefficients with \(A\cos\beta t+C\sin\beta t\). This gives \(R\sin\phi=A=5\) and \(R\cos\phi=C=3\). Therefore:

$$R=\sqrt{A^2+C^2}=\sqrt{34}=5.83, \qquad \phi=\arctan\!\frac{A}{C}=\arctan\!\frac{5}{3}=1.03\ \text{rad},$$

$$\boxed{f(t) = e^{-2t}\big[5\cos(5t) + 3\sin(5t)\big]u(t) = 5.83\,e^{-2t}\sin(5t+1.03)\,u(t).}$$

The denominator is already in completed-square form, so the poles are \(s=-2\pm j5\). The imaginary part, \(5\), sets the oscillation frequency. The real part, \(-2\), produces the decaying envelope \(e^{-2t}\). Since both poles are in the open left-half plane, the causal system is stable.

MATLAB corner

Eq. (2.15) is the decaying oscillation \(h(t)\); its two exponential envelopes \(\pm\sqrt{A^2+C^2}\,e^{-\alpha t}\) pinch it at every peak, since \(\alpha\) alone sets the decay rate regardless of \(\beta\).

Goal: plot the actual decaying oscillation \(h(t)\) against its \(\pm e^{-\alpha t}\) envelope, to see the decay rate is set by \(\alpha\) alone.

% damped_sinusoid_demo.m - h(t) for complex-conjugate poles, Eq. (2.15)
% H(s) = (5s+25)/(s^2+4s+29): alpha=2, beta=5, A=5, C=3 (Examples 2.6 and 2.7)
alpha = 2; beta = 5; A = 5; C = 3;
t = 0:0.001:3;
h = exp(-alpha*t).*(A*cos(beta*t) + C*sin(beta*t));   % Eq. (2.15)
 
env = sqrt(A^2+C^2)*exp(-alpha*t);                    % +-e^(-alpha t) envelope
 
figure('Color','w'); hold on; grid on;
plot(t, h, 'LineWidth', 1.6);
plot(t, env, 'r--', t, -env, 'r--', 'LineWidth', 1.2);
xlabel('t (s)'); ylabel('h(t)');
legend('h(t) = e^{-2t}[5cos(5t)+3sin(5t)]', 'envelope \pm\surd(A^2+C^2)e^{-\alpha t}');
title('Damped oscillation decays inside its envelope, rate set by \alpha only');
 
% Sanity check against ilaplace, symbolically
syms s tt
Hs = (5*s+25)/(s^2+4*s+29);
ht = ilaplace(Hs, s, tt);
disp('Symbolic check via ilaplace:'); disp(simplify(ht))

2.4  Stability

All three reduce to one geometric test on the pole–zero plot: for a causal system, every pole must lie in the open left-half \(s\)-plane,

A pole on or to the right of the \(j\Omega\)-axis means \(h(t)\) does not decay, and BIBO stability fails.

2.4.1  Reading stability off a complex pole pair

For the important case \(p_{1,2}=-\alpha\pm j\beta\), the impulse response has the form \(h(t)=D e^{-\alpha t}\sin(\beta t+\phi)u(t)\). The two parts of each pole have different roles: \(\beta\) sets the oscillation frequency, while \(-\alpha\) determines whether the envelope decays, remains constant or grows.

So stability is decided by the sign of \(\alpha\) alone, and there are exactly three outcomes:

Marginal stability is a boundary case, not a comfortably stable one. Poles on the \(j\Omega\)-axis produce oscillations that do not decay. A bounded input at the same frequency can make the response grow through resonance, so BIBO stability fails. An ideal lossless \(LC\) circuit would behave this way. In a real \(RLC\) circuit, resistance moves the poles into the left-half plane and gives the resonance a finite height.

EXAMPLE 2.8 · An unstable system, from h(t) to the step response

A system has \(h(t)=\tfrac12 e^{-t}\big(e^{2t}-1\big)\). Find \(H(s)\), locate its poles, classify it, and find the unit-step response.

{{ l2Ex5bLabel }}

Step 1. Expanding first is the whole trick: \(h(t)=\tfrac12 e^{t}-\tfrac12 e^{-t}\), which is a sum of two exponentials, so Section 2.3 runs backwards and each term contributes one pole:

$$H(s)=\frac{1}{2(s-1)}-\frac{1}{2(s+1)}=\frac{1}{(s-1)(s+1)}.$$

Step 2. The poles are \(p=\pm1\). One of them, \(p=+1\), lies in the right half-plane, so Eq. (2.16) fails and the system is unstable — and the \(e^{+t}\) term in \(h(t)\) says so directly.

Step 3 · Step response. With \(x(t)=u(t)\), \(X(s)=1/s\), so \(Y(s)=H(s)/s\). Partial fractions on the three distinct poles \(\{0,+1,-1\}\) give

$$Y(s)=\frac{1}{2(s-1)}-\frac{1}{s}+\frac{1}{2(s+1)} \ \xrightarrow{\ \mathcal{L}^{-1}\ }\ \boxed{y(t)=-1+\tfrac12 e^{-t}+\tfrac12 e^{t}.}$$

The step input is bounded, but the output grows without limit because of the \(e^t\) term. The system therefore fails the BIBO stability test. The pole at \(+1\) predicted this result directly. No ROC for this causal system includes the \(j\Omega\)-axis, so a conventional frequency response \(H(j\Omega)\) does not exist.

How to read the result. The upper plots connect the impulse response with its poles. The pole at \(+1\) lies in the unstable half-plane and creates a growing exponential. The lower plot shows the consequence: a bounded step input produces an unbounded output. This confirms the pole test without needing to inspect every time-domain term.

MATLAB corner

MATLAB can apply the stability test in Eq. (2.16) directly: find the poles and inspect the signs of their real parts. The pzmap() command then makes the comparison visual. Poles in the left half-plane produce decaying modes; any pole in the right half-plane produces a growing mode.

Goal: verify the RHP/LHP stability test numerically on both a stable and an unstable system, and reproduce Figure 2.3 from real pole-zero maps.

% compare_systems_pz.m - increase font and emphasize poles/zeros
sys(1).num = [5 25]; sys(1).den = [1 4 29];
sys(2).num = [1];    sys(2).den = [1 0 -1];
figure('Color','w');
for k = 1:2
   p = roots(sys(k).den);
   stable = all(real(p) < 0);
   fprintf('System %d poles: %s -> stable: %d\n', k, mat2str(p,3), stable);
   ax = subplot(1,2,k);
   H = tf(sys(k).num, sys(k).den);
   pzmap(H); grid on; hold on
   plot(real(p), imag(p), 'rx', 'MarkerSize',20, 'LineWidth',2);   % emphasize poles
   z = zero(H);
   fontsize(gcf,16,"points")                % increase all figure font sizes
   plot(real(z), imag(z), 'go', 'MarkerSize',16, 'LineWidth',2);   % emphasize zeros
   ylim(ax,[-6 6])
   ax.FontSize = 18;
   if stable
      title(ax,'System ' + string(k) + ': stable');
   else
      title(ax,'System ' + string(k) + ': unstable');
   end
end

2.5  Reading a system off its pole-zero plot

A pole-zero plot collects several important properties in one picture. It shows natural modes, stability and the main features of the frequency response. For this reason, it is often the first plot used to inspect a system and the final plot checked after a design.

A rational transfer function is fixed, up to a gain, by where its poles and zeros sit:

At a chosen value of \(s\), each factor represents the distance from a pole or zero to that point. On the \(j\Omega\)-axis, the magnitude tends to increase near a pole and decrease near a zero. This geometric rule helps predict resonances and attenuation directly from the plot and becomes a central design idea in Lecture 5.

2.5.1  Invertibility, and why zeros have a side too

Inverting a system means undoing it: finding \(H_{\mathrm{inv}}(s)\) such that the cascade is the identity. Algebraically this is trivial —

The zeros of \(H\) become the poles of the inverse \(H_{\mathrm{inv}}\). Therefore, a right-half-plane zero in the original system becomes an unstable pole in its causal inverse. A stable causal system whose poles and zeros are all in the left-half plane is called minimum phase. Such a system can have a stable causal inverse.

The name comes from a second fact: among all systems sharing one magnitude response, the minimum-phase one has the least phase lag and the fastest energy build-up. Any other is that system cascaded with an all-pass factor,

which changes the phase while leaving the magnitude untouched. This is why a magnitude specification alone never pins down a filter, and why Lecture 5's design methods must say something about phase separately.

2.6  From Laplace to Fourier: the CTFT as a special case

Recall \(s=\sigma+j\Omega\), so the Laplace kernel splits into a damping factor and an oscillation: \(e^{-st}=e^{-\sigma t}e^{-j\Omega t}\). Setting \(\sigma=0\) strips away the damping and leaves a purely frequency-domain view:

This is the Continuous-Time Fourier Transform (CTFT). It exists when the ROC of \(X(s)\) includes the \(j\Omega\)-axis. For a causal stable LTI system, all poles are in the left-half plane, so this condition is satisfied. Its frequency response is therefore obtained by evaluating the transfer function at \(s=j\Omega\).

The transform runs backwards through the same axis, giving the inverse CTFT, and, for a real-valued signal, the conjugate symmetry already met in Lecture 1 (Eq. 1.7):

The analysis integral measures the amplitude and phase associated with each frequency. The synthesis integral adds those frequency components to reconstruct the signal. Because the transform is reversible, we can specify and design a system through its frequency response \(H(j\Omega)\), then return to the time domain when needed.

Removing the Laplace damping factor makes convergence more restrictive. In particular, the ordinary CTFT integral does not converge for a nonzero periodic signal:

This matches the distinction from Lecture 1. The ordinary CTFT is naturally suited to energy signals, while a nonzero periodic signal is a power signal. Periodic signals therefore need the Fourier series description introduced in Section 2.11. Discrete-time signals require the corresponding DTFT and DTFS forms developed later.

MATLAB corner

To obtain the frequency response, Eq. (2.19) evaluates the transfer function on the imaginary axis by setting \(s=j\Omega\). MATLAB's freqs() performs this calculation from the coefficient vectors. We compare its output with the closed-form expression from Example 2.10; the two should agree apart from numerical rounding.

Goal: evaluate \(H(s)=1/(s+a)\) at \(s=j\Omega\) via freqs() and confirm it matches Example 2.10's closed-form CTFT, including the \(-3\,\text{dB}\) point at \(\Omega=a\).

% ctft_from_laplace.m - Eq. (2.19): CTFT is H(s) evaluated at s = j*omega
a = 3;  num = 1;  den = [1 a];        % H(s) = 1/(s+a), Eq. (2.5) form
w = linspace(0, 15, 500);
 
Hjw = freqs(num, den, w);             % MATLAB sets s = j*w for you - Eq. (2.19)
Hclosed = 1 ./ (a + 1j*w);            % Example 2.10's closed-form answer
fprintf('max abs diff: %.2e\n', max(abs(Hjw - Hclosed)));
 
figure('Color','w');
subplot(2,1,1); plot(w,abs(Hjw),'LineWidth',1.8); hold on;
plot(w,abs(Hclosed),'go','MarkerSize',4); grid on;
legend('freqs(H(s)) at s=j\Omega', 'Example 2.10 closed form');
xlabel('\Omega (rad/s)'); ylabel('|X(\Omega)|');
 
subplot(2,1,2); plot(w, angle(Hjw)*180/pi, 'LineWidth', 1.8); grid on;
yline(-45,'r--'); xline(a,'r--');
xlabel('\Omega (rad/s)'); ylabel('phase (deg)');
title('phase = -45 deg exactly at \Omega = a');

2.6.1  Convergence: when does the CTFT actually exist?

The substitution \(s=j\Omega\) is valid only when the \(j\Omega\)-axis lies inside the ROC. If it does not, the Fourier integral diverges. Two common sufficient conditions help us decide when a CTFT exists.

The first is absolute integrability, which guarantees the integral converges for every \(\Omega\):

Together with the requirement that a finite interval contains only finitely many extrema and finite jumps, these form the Dirichlet conditions. At a continuous point, the inverse transform returns \(x(t)\). At a jump, it returns the midpoint of the two limiting values. Finite-energy signals may also converge in the mean-square sense, where the total error energy tends to zero even if individual points do not match exactly.

The sinc function decays as \(1/t\). It is square-integrable but not absolutely integrable, so its transform is interpreted through mean-square convergence. This issue appears again in the ideal brick-wall filter. Pure sinusoids and steps do not satisfy either ordinary condition; their Fourier transforms are represented using impulses in frequency. These impulse-valued transforms are useful, but they are generalized functions rather than ordinary curves.

2.7  Why complex exponentials? The basis an LTI system cannot break

Many sets of functions could be used to build a signal. Fourier analysis uses sinusoids, written compactly as complex exponentials, because they interact with LTI systems in an exceptionally simple way.

2.7.1  Sinusoids are the one shape an LTI system cannot change

Convolution adds scaled and delayed copies of the input. Delaying a sinusoid changes its phase but not its frequency. Adding several sinusoids at the same frequency produces another sinusoid at that frequency. Therefore, an LTI system can change a sinusoid's amplitude and phase, but not its frequency:

For example, a \(1\,\text{kHz}\) input to an LTI system produces a \(1\,\text{kHz}\) output. The output may be weaker, stronger or phase shifted, but it cannot become a different frequency or a square wave. A nonlinear system can create new frequencies; that behaviour is harmonic or intermodulation distortion.

2.7.2  Sinusoids of different frequencies are orthogonal

Closure alone would be useless if the components interfered with each other's bookkeeping. They do not: over a common period, distinct harmonics integrate against each other to exactly zero.

Different complex exponentials are orthogonal, so each Fourier coefficient can be measured independently. Multiply the signal by the conjugate of the required basis function, then integrate or sum. The formulas for \(X(\Omega)\), \(X_k\), \(X(e^{j\omega})\) and \(X[k]\) all follow this same pattern.

2.7.3  Euler's formula: one object instead of two

Sines and cosines work, but they are clumsy: a sinusoid carries an amplitude and a phase, and keeping them in separate columns makes every calculation twice as long. Euler's identity collapses the pair.

Three concrete payoffs follow, and each one shows up somewhere later in the module.

• Amplitude and phase become one number. The phasor \(Ae^{j\varphi}\) carries both, and complex multiplication applies a gain and a phase shift in one operation. That is exactly what a filter does to a tone, which is why a filter's whole effect at one frequency is the single complex number \(H(j\Omega)\).
• Calculus becomes arithmetic. \(\frac{d}{dt}e^{j\Omega t}=j\Omega\,e^{j\Omega t}\) and \(\int e^{j\Omega t}dt = \frac{1}{j\Omega}e^{j\Omega t}\): differentiation is multiplication by \(j\Omega\), integration is division by it. This is the Fourier-domain echo of the \(s^k\) that appeared when the LCCDE was Laplace-transformed in Section 2.1.3.
• Negative frequencies stop being mysterious. Eq. (2.25) says a real cosine is the sum of two counter-rotating phasors, one at \(+\Omega\) and one at \(-\Omega\). The negative half of every spectrum in this lecture is not a physical frequency; it is the other half of a real signal's conjugate pair, and it is why \(X(-\Omega)=X^{*}(\Omega)\).

2.7.4  The eigenfunction property: the reason the whole ladder exists

The three properties above converge on one statement. Push \(x[n]=e^{j\omega n}\) through the convolution sum of Eq. (1.14) and watch the exponential factor out:

The output is the input multiplied by a constant. In the language of linear algebra, complex exponentials are the eigenfunctions of every LTI/LSI system, and the eigenvalue is the frequency response, exactly as Lecture 1's Problem 1.7 proved for \(z^n\):

Now add linearity. If the input can be written as a sum of exponentials, the system acts on each one separately and never mixes them:

MATLAB corner

Eq. (2.26)–(2.27) are a strong claim: push a complex exponential through any LTI system and only a single complex number changes. That is easy to take on faith and easy to check.

Goal: filter a complex exponential through a real IIR system and verify the steady-state output is exactly the input multiplied by \(H(e^{j\omega_0})\).

% eigenfunction_demo_fixed.m - verify Eq. (2.26)-(2.27) numerically
b = [0.2 0.2]; a = [1 -0.6];
w0 = 0.4;                          % test frequency, rad/sample
n = (0:200)';
x = exp(1j*w0*n);
y = filter(b, a, x);
Hw0 = polyval(b, exp(-1j*w0)) ./ polyval(a, exp(-1j*w0));   % H(e^{j*w0}) evaluated directly
steady = 100:200;
predicted = Hw0 * x(steady+1);
fprintf('max steady-state error: %.2e\n', max(abs(y(steady+1) - predicted)));
figure('Color','w');
    fontsize(gcf,16,"points")                % increase all figure font sizes
subplot(2,1,1); plot(n, real(x), n, real(y), 'LineWidth', 1.4); grid on;
legend('Re{x[n]}', 'Re{y[n]}'); title('same frequency in and out - Eq. (2.26)');
subplot(2,1,2);
plot(n(steady+1), abs(y(steady+1)), 'LineWidth',1.4); hold on;
plot(n(steady+1), abs(predicted), 'ro', 'MarkerSize',4); grid on;
ylim([0.5 1]);                       % set y-axis between 0.5 and 1
legend('|y[n]|', '|H(e^{j\omega_0})| . |x[n]|');
title('eigenvalue confirmed');

2.8  Four kinds of signal, four transforms: the map

The Fourier formulas are not separate ideas to memorize. They all decompose a signal into complex exponentials. The correct formula depends on two questions: Is time continuous or discrete? Is the signal periodic or aperiodic?

Lecture 1 already supplied both axes, and neither has more than two values:

These two choices create four signal classes and four corresponding Fourier representations. The DFT is the finite, computational form associated with periodic discrete-time analysis. The FFT is not another transform; it is an efficient algorithm for calculating the DFT.

The fourth class is especially useful for computation. One period contains a finite number of time samples and the spectrum contains the same finite number of coefficients. Both lists fit in memory. This finite representation is the DFT.

2.9  The CTFT in use: pairs, and the properties that matter

The first case is an aperiodic continuous-time signal, described by the CTFT in Eq. (2.19). In practice, most transforms are found from a short table of standard pairs together with properties such as scaling, shifting and convolution.

2.9.1  Common CTFT pairs

Constant signals and pure sinusoids have infinite energy, so their CTFTs are represented by impulses in frequency. Also notice the duality between a rectangular pulse and a sinc: a rectangle in one domain corresponds to a sinc-shaped function in the other. This explains why rectangular windows and ideal rectangular filters produce closely related side-lobe behaviour.

2.9.2  Convolution becomes multiplication: this is filtering

This is the property the whole module is built on, and it is Eq. (2.4) again with \(s=j\Omega\):

A filter scales each frequency component by \(H(j\Omega)\). For a low-pass filter with \(f_c=150\,\text{Hz}\), the 50 Hz component of a two-tone input remains large, while the 300 Hz component is multiplied by a much smaller value. It is attenuated rather than mathematically erased. Filter design chooses \(H\) so that these frequency-dependent scale factors meet a specification.

2.9.3  Multiplication becomes convolution: windowing and leakage

Run the property the other way and it stops being a convenience and becomes a warning:

A measurement always has finite duration, which means the true signal has been multiplied by a window. Multiplication in time becomes convolution in frequency, so each ideal spectral line is spread by the window spectrum. This spreading is called spectral leakage; it is a mathematical consequence of finite observation, not a software error.

A rectangular window has a narrow main lobe, which helps separate nearby tones, but high side lobes cause strong leakage. Tapered windows such as Hann, Hamming and Blackman reduce the side lobes but widen the main lobe. Window choice is therefore a trade-off between frequency resolution and leakage.

MATLAB corner

The next script makes Eq. (2.29) and Eq. (2.30) visible in one sitting. Its first half low-pass-filters a two-tone signal purely by multiplying spectra; its second half shows why the window you truncate with decides how badly a clean tone smears.

Goal: show a filter removing a tone purely by spectral multiplication, then show how the truncation window you pick trades main-lobe width against leakage.

%% L2_7_fixed.m - Convolution/multiplication demos (non-ASCII spaces removed)
%% Part 1 - Convolution in time IS multiplication in frequency  (Eq. 2.29)
fs = 2000;                       % fine grid, so the FFT approximates the CTFT
t  = -0.1 : 1/fs : 0.1;
x  = sin(2*pi*50*t) + 0.5*sin(2*pi*300*t);
 
fc = 150;                        % low-pass cut-off in Hz
h  = 2*fc.*sinc(2*fc.*t) .* hann(numel(t))';   % windowed ideal LPF
 
N = 4096;  f = (-N/2:N/2-1)*fs/N;
X = fftshift(fft(x,N))/fs;
H = fftshift(fft(h,N))/fs;
Y = X .* H;                      % <-- the entire point of the section
 
figure('Color','w');
subplot(3,1,1); plot(f,abs(X),'LineWidth',1.2); xlim([-600 600]); grid on;
title('|X(j\Omega)|  -  50 Hz and 300 Hz'); xlabel('f (Hz)');
subplot(3,1,2); plot(f,abs(H),'LineWidth',1.2); xlim([-600 600]); grid on;
title('|H(j\Omega)|  -  low-pass, f_c = 150 Hz'); xlabel('f (Hz)');
subplot(3,1,3); plot(f,abs(Y),'LineWidth',1.2); xlim([-600 600]); grid on;
title('|Y| = |X||H|  -  the 300 Hz line has gone'); xlabel('f (Hz)');
 
%% Part 2 - Multiplication in time IS convolution in frequency  (Eq. 2.30)
% 0.1 cycles/sample over L = 128 samples is 12.8 cycles: NOT a whole
% number, so the assumed periodic extension has a jump -> leakage.
L = 128;  Nfft = 4096;  n = 0:L-1;
xs = sin(2*pi*0.1*n);
Wr = abs(fftshift(fft(xs .* rectwin(L)', Nfft)));
Wh = abs(fftshift(fft(xs .* hamming(L)', Nfft)));
w  = (-Nfft/2:Nfft/2-1)/Nfft*2*pi;
figure('Color','w'); hold on; grid on;
plot(w, 20*log10(Wr/max(Wr)), 'LineWidth', 1.1);
plot(w, 20*log10(Wh/max(Wh)), 'LineWidth', 1.1);
xlabel('\omega (rad/sample)'); ylabel('dB'); ylim([-100 5]);
legend('rectangular window', 'Hamming window');
title('Leakage: narrow main lobe vs low side lobes - pick one');

2.9.4  A delay is a linear phase ramp

A time delay leaves the magnitude unchanged and adds a phase term \(-\Omega t_0\), which is linear in frequency. Conversely, a filter with linear phase delays every frequency component by the same time, so the waveform shape is preserved. Nonlinear phase gives different delays to different components and can distort pulses even when the magnitude response is acceptable.

2.9.5  Parseval: the frequency domain is not a metaphor

One last property makes the spectrum quantitative rather than decorative. Taking inner products through the convolution theorem gives

Parseval's relation states that the signal energy calculated in time equals the energy calculated in frequency. The transform changes the representation, not the total energy. This allows filter performance to be measured by integrating \(|X(\Omega)|^2\) over passbands or stopbands.

2.10  Frequency-selective filters: the ideal, and why you cannot have it

Because convolution becomes multiplication in frequency, we can design a filter by choosing the shape of \(H(j\Omega)\). Most frequency-selective filters belong to one of four groups: lowpass, highpass, bandpass or bandstop.

The ideal lowpass filter is defined by nothing more than an indicator function on frequency, with cutoff \(\Omega_c\):

This ideal response is a useful design target, but it cannot be implemented exactly. Taking its inverse Fourier transform gives

The result is a sinc impulse response that exists for all time, including \(t<0\). A real-time filter cannot respond before the input arrives, so the ideal filter is non-causal. It is also infinitely long and decays slowly as \(1/t\). Cutting it to a finite length makes implementation possible, but introduces ripple in the frequency response. Lecture 7 shows how window functions control this trade-off.

A practical filter is therefore specified by four quantities: the passband edge, the stopband edge, the allowed passband ripple and the required stopband attenuation. Lectures 5 and 7 show how to meet these requirements with the lowest reasonable filter order. The brick-wall response remains a reference, not a realizable design.

2.11  Rung two: periodic signals and the CT Fourier series

The ordinary CTFT integral does not converge for a periodic signal. Instead, we describe the signal using a Fourier series. If \(x(t)\) repeats every \(T_0\), its frequency components occur only at integer multiples of the fundamental frequency \(\Omega_0=2\pi/T_0\). These frequencies \(k\Omega_0\) are called harmonics, so a sum is sufficient; we do not need an integral over every possible frequency.

The analysis formula isolates one harmonic at a time. Multiplying by \(e^{-jk\Omega_0t}\) and integrating over one period makes all other harmonics cancel because they are orthogonal. The complex coefficient \(X_k\) stores both the amplitude and phase of harmonic \(k\). For a real signal, \(X_{-k}=X_k^{*}\), so the magnitude spectrum is even and the phase spectrum is odd.

Parseval reappears too, now measuring power instead of energy, exactly as the energy/power split of Lecture 1 demands:

The average power of a periodic signal is the sum of the powers of its harmonics; a sinusoid of amplitude \(A\) has \(P=A^2/2\), split evenly between its two lines at \(\pm\Omega_0\).

2.11.1  The CTFS and the CTFT are not rivals

The CTFS and CTFT are complementary descriptions, not competing methods. Taking the CTFT of the Fourier-series expression term by term, and using \(e^{jk\Omega_0t}\leftrightarrow 2\pi\delta(\Omega-k\Omega_0)\), gives

In the generalized sense, the CTFT of a periodic signal is an impulse train. One impulse appears at each harmonic \(k\Omega_0\), weighted by the same coefficient \(X_k\) used in the Fourier series. The two descriptions therefore contain the same information. For calculations, it is usually simpler to list the coefficients \(X_k\) than to write the full impulse spectrum.

2.11.2  The periodic square wave, and a sinc that will come back

Take the standard test case: amplitude \(A\), pulse width \(\tau\), period \(T_0\). Only the pulse contributes to the integral, so

The middle step is just Euler's formula (2.25) read backwards, which is the third time in this lecture that a messy integral has been rescued by turning exponentials into a sine.

MATLAB corner

Example 2.13 found the Fourier-series coefficients by integration. Here we sample one period of the same square wave and calculate an FFT. Dividing the FFT by the number of samples gives coefficients that can be compared directly with the analytic values.

Goal: build Example 2.13's square wave, take its FFT, and check the bins against the analytic \(X_k\) of Eq. (2.38) and the power split of Eq. (2.36).

% square_wave_fft.m - Example 2.13's square wave, checked against the FFT
A = 1; tau = 1; T0 = 4; fs = 200;      % Example 2.13's numbers
t = (0:1/fs:T0-1/fs)';
x = A * (mod(t+tau/2, T0) < tau);     % one period, pulse centred at t=0
N = numel(x);
 
Xk_fft = fft(x)/N;                     % DFT/N approximates the CTFS, Eq. (2.47) vs (2.35)
k = (0:N-1)'; k(k>N/2) = k(k>N/2) - N;   % centre the harmonic index on 0
w0 = 2*pi/T0;
Xk_analytic = (A*tau/T0) * sinc(k*w0*tau/(2*pi));  % Eq. (2.38)
fprintf('max |FFT - analytic| for |k|<=10: %.2e\n', max(abs(Xk_fft(abs(k)<=10) - Xk_analytic(abs(k)<=10))));
 
P_total = A^2*tau/T0;                          % Eq. (2.36), direct power
P_partial = @(K) sum(abs(Xk_fft(abs(k)<=K)).^2);
fprintf('power in |k|<=7: %.1f%% of total\n', 100*P_partial(7)/P_total);
 
figure('Color','w');
subplot(2,1,1); plot(t, x, 'LineWidth', 1.6); grid on; ylim([-0.2 1.2]);
title('one period, A=1, \tau=1, T_0=4'); xlabel('t (s)');
 
subplot(2,1,2); stem(k(abs(k)<=10), abs(Xk_fft(abs(k)<=10)), 'filled'); hold on;
stem(k(abs(k)<=10), abs(Xk_analytic(abs(k)<=10)), 'ro'); grid on;
legend('|X_k| from fft(x)/N', '|X_k| from Eq. (2.38)'); xlabel('harmonic k');
title('FFT bins land exactly on the analytic sinc envelope');

This relationship will return in Section 2.14 with time and frequency exchanged. Making a signal periodic in time produces discrete frequency lines. In the opposite direction, sampling a spectrum at discrete frequencies makes the reconstructed time sequence periodic. This is the periodicity assumption built into the DFT.

2.11.3  Gibbs phenomenon: the overshoot that never leaves

Adding more harmonics gives a better square-wave approximation overall, because the total error energy decreases as \(N\) increases. Near each discontinuity, however, a visible overshoot remains.

Each partial sum overshoots the jump, and the overshoot does not shrink. Increasing \(N\) makes the ripple narrower and pushes it closer to the edge, but its height converges to a fixed fraction of the jump:

The limiting overshoot is approximately \(8.95\%\) of the jump height \(\Delta\) on each side. This is the Gibbs phenomenon. Increasing \(N\) makes the ripple narrower but does not remove its peak height. The ripple energy still tends to zero because it is concentrated into a smaller region.

% gibbs_demo.m - partial-sum reconstruction of a square wave, and the
% overshoot that will not go away  (Eq. 2.39)
T0 = 1; t = linspace(-1.5*T0, 1.5*T0, 4000);
Ns = [5 15 51 201];            % number of harmonics kept
% odd-harmonic Fourier series of a unit-amplitude square wave
figure; hold on
for N = Ns
    x = zeros(size(t));
    for k = 1:2:N
        x = x + (4/(pi*k)) * sin(2*pi*k*t/T0);
    end
    plot(t, x, 'LineWidth', 1.0)
end
yline(1, 'k--'); yline(1.0895, 'r:')   % ideal level and 8.95% overshoot
legend(compose('N = %d', Ns), 'Location', 'southeast')
xlabel('t / s'); ylabel('partial sum'); grid on
% peak overshoot, measured near the discontinuity at t = 0
for N = Ns
    x = zeros(size(t));
    for k = 1:2:N, x = x + (4/(pi*k))*sin(2*pi*k*t/T0); end
    fprintf('N = %3d  peak = %.4f\n', N, max(x));
end