Frequency response of discrete-time systems

The leaky integrator of Difference equations (6.1) gets four test sines, one after another. Watch the gain fall as the frequency $\Omega$ grows, and the output squares lag the input dots by less and less.

A test sine goes in, a scaled and shifted sine comes out

On the page Properties of LTI systems (5.4) you found that an arrow $e^{j\Omega n}$ goes through a linear time-invariant system and comes out as the same arrow times one number. That number is $H(e^{j\Omega})=\sum_kh[k]e^{-j\Omega k}$, the DTFT of the impulse response $h$ from The DTFT (12.2). It is the system’s frequency response. This page measures it, reads it off a difference equation, and then asks what its angle does to a signal.

A real signal is not a single arrow. From Properties of the DTFT (12.3), a cosine is two half arrows, at $+\Omega$ and at $-\Omega$, and a real $h$ gives $H(e^{-j\Omega})=H^*(e^{j\Omega})$. The two outputs are therefore a pair of conjugates, and they add to a cosine:

where $\lvert H\rvert$ and $\angle H$ are the size and angle of $H(e^{j\Omega})$ at that $\Omega$. This holds once the test sine has been playing for a long time, so that the start-up part of the output has died away. I call that the steady state. A negative $\angle H$ means the output crests come later than the input crests.

To get $H$ from a difference equation, put $x[n]=e^{j\Omega n}$ and $y[n]=He^{j\Omega n}$ into the leaky integrator, $y[n]=0.5y[n-1]+0.5x[n]$. The delayed output is $He^{j\Omega(n-1)}=e^{-j\Omega}He^{j\Omega n}$, so

In general, with $a_0=1$, the response is $\sum_kb_ke^{-j\Omega k}$ divided by $\sum_ka_ke^{-j\Omega k}$. That is what SciPy’s freqz(b, a) computes. Its size is $\lvert H\rvert=0.5/\sqrt{1.25-\cos\Omega}$, because $\lvert1-0.5e^{-j\Omega}\rvert^2=1+0.25-\cos\Omega$.

The four tests in the picture at the top lie on this $H(e^{j\Omega})$, and the picture ends by drawing it through them. The last test, at $\Omega=\pi$, is the sine $\cos(\pi n)=(-1)^n$, which alternates between $+1$ and $-1$. The phase is $0°$ at $\Omega=0$, dips to a minimum of $-30.0°$ at $\pi/3$, and returns to $0°$ at $\pi$.

The gain falls from 1 at $\Omega=0$ to $\tfrac13$ at $\Omega=\pi$, so the leaky integrator keeps slow waves and weakens fast ones. In decibels (1.3, How big is a signal) the four tests read $-3.37$, $-6.99$, $-8.94$ and $-9.54$ dB.

A moving average, probed

A moving average replaces each sample by the mean of the last few. The 4-point version is $y[n]=\tfrac14\bigl(x[n]+x[n-1]+x[n-2]+x[n-3]\bigr)$. If a sine fits a whole number of cycles into those four samples, the four values always cancel. So I expect the output to be exactly zero for some test sines, and nearly the input for slow ones. Watch the output go flat at 0 for two of the five tests.

Notice that the output is flat at 0 for $0.5\pi$ and $\pi$. When the clip has finished, drag across the gain panel to try any test frequency from 0 to $\pi$.

Now the formula. Each term is a delay, so $H=\tfrac14\bigl(1+e^{-j\Omega}+e^{-j2\Omega}+e^{-j3\Omega}\bigr)$. This is a short geometric sum, and it folds into

The fraction is the pulse spectrum of 12.2, and the factor $e^{-j1.5\Omega}$ is its centre shifted by 1.5 samples. The numerator is zero at $\Omega=0.5\pi$ and $\pi$, where the denominator is not. For $N$ points the zeros sit at $\Omega=2\pi k/N$.

That is why a 12-month rolling mean removes a yearly pattern completely. One cycle per year is $2\pi/12$ rad per month, which is the first zero. A 7-day mean puts its first zero at $2\pi/7$ rad per day, and removes the weekly cycle from daily counts.

One more system, the first difference $y[n]=x[n]-x[n-1]$ from Oversampling and noise shaping (11.3). There you fed it four test sines and plotted the dots. Its frequency response is

so $\lvert H\rvert=2\lvert\sin(\Omega/2)\rvert$. This is the curve that page promised. The four dots of 11.3 sit on it.

The squared curve is the power gain: $0.0246$, $1$, $2$, $3.9754$ and $4$ at $\pi/20$, $\pi/3$, $\pi/2$, $19\pi/20$ and $\pi$. The curve rises with $\Omega$, so the first difference is a high-pass system, the opposite of the moving average.

The ideal low-pass cannot be built

What would a perfect low-pass look like? It would keep every frequency below a cutoff $\Omega_c$ at gain 1 and remove all the others. In discrete time it is repeated every $2\pi$, so on $(-\pi,\pi]$ it is $H=1$ for $\lvert\Omega\rvert < \Omega_c$ and 0 elsewhere. Its impulse response is the inverse DTFT (12.2),

The stems are drawn for $\Omega_c=0.25\pi$. They shrink as $1/n$ but never stop, in either direction. So $h$ is not causal and not FIR (5.4): a system that responds before the input arrives, and with infinitely many taps, cannot be run. This is the discrete twin of the brick wall in Frequency response and Bode plots (8.4).

Delay $h$ and cut it off, and you get a filter you can run. The cut produces a ripple at the edge of the passband, the same overshoot as cutting a series short in Convergence and the Gibbs phenomenon (7.3). Designing with that trade is the subject of Window-method FIR design (19.1).

Straight phase, bent phase

Phase is a delay in disguise. In Frequency response and Bode plots (8.4, “Two delays”) you measured it in seconds. In samples, the phase delay is

the lateness of the crests inside a signal. The group delay is

the lateness of the envelope. Both use the unwrapped phase, the angle drawn without its jumps of $2\pi$ (12.2). The leaky integrator, for example, has $\tau_p=0.870$, $0.637$, $0.295$ and $0.108$ samples at $0.125\pi$, $0.25\pi$, $0.5\pi$ and $0.75\pi$.

A pure delay of $n_0$ samples has $H=e^{-j\Omega n_0}$. Its phase is the straight line $-\Omega n_0$, so $\tau_g=n_0$ at every $\Omega$ and the shape passes through unchanged. A filter whose phase is a straight line is linear phase.

Two short FIR filters show the difference. The triangle $(1,2,3,2,1)/9$ is the one from 12.3. Its response is

The bracket is squared, so it is never negative, and the phase is exactly $-2\Omega$. The ramp $(1,2,3)/6$ is not symmetric, and its phase bends. Watch the two group delays, minus the slopes of the phases: the triangle’s stays at 2.00 samples, and the ramp’s does not.

The ramp’s principal angle jumps from $-180°$ to $+180°$ at $0.608\pi$. Unwrapped, it runs on smoothly to $-255.4°$ at $0.75\pi$. When the clip has finished, drag the dotted line to read both delays at any $\Omega$; at $0.6\pi$ the ramp delays by 2.61 samples.

A bent phase delays each frequency by its own amount, so a signal made of several frequencies comes out reshaped. This is phase distortion. Picture runners who all keep the same pace: they cross the line in formation. Give each a different pace and the formation breaks up.

Linear phase comes from a symmetric $h$, and Linear-phase systems (17.3) builds on that. You can check it on the triangle. Its output for $\cos(0.25\pi n)$ equals $0.648\,x[n-2]$ exactly, from $n=4$ on. The gain $0.648$ is the bracket squared at $0.25\pi$, and the delay is the 2 samples of the straight line.

Worked example

Take the leaky integrator and read its response at $\Omega=0.5\pi$ by hand. Here $e^{-j0.5\pi}=-j$, so $H=0.5/(1+0.5j)$. Multiply top and bottom by $1-0.5j$: $H=0.5(1-0.5j)/1.25=0.4-0.2j$. Its size is $\sqrt{0.16+0.04}=0.4472$, and its angle is $-\arctan(0.5)=-26.57°$. This matches the second test in the picture at the top.

Now run the 4-point moving average from rest on $\cos(0.5\pi n)$. The samples are $1,0,-1,0,1,\dots$, so the outputs are $0.25$, $0.25$, then $0$ from $n=2$ on. The start-up part lasts only while the four-sample window is filling. After that the output is exactly the steady-state value, 0, which the gain of 0 predicted.

For the first difference at $\Omega=\pi/3$, $\lvert H\rvert=2\sin(\pi/6)=1$, so a sine at that frequency comes out the same size. At $\Omega=19\pi/20$ the gain is $2\sin(19\pi/40)=1.994$, close to the largest value, 2, at $\pi$.

The ramp $(1,2,3)/6$ at $\Omega=0.75\pi$ has gain $0.273$ and group delay $2.925$ samples. Its largest group delay, $3.061$ samples, is at $0.699\pi$. The triangle delays everything by exactly 2.

Where you will meet this

Every filter datasheet in discrete time is a plot of this page, a gain curve and a phase or delay curve against $\Omega$. The moving average is the simplest smoother, and Simple smoothing filters (18.2) compares it with others. The poles and zeros that shape such curves come in Transfer functions, poles & zeros (16.3). The ideal low-pass returns, cut down to a usable length, in Window-method FIR design (19.1).

Reference card