# Can someone explain why signal.freqz() generates a different phase in the frequency fresponse than scipy.signal.dfreqresp()?

I have been using scipy to analyze filter performance for a single-pole IIR filter, and I noticed a disagreement between the phase in the outputs of the freqz() function compared to the output of the dfreqresp() function for the same system.

The transfer function for this filter is: $$H(z) = \frac{b}{1-az^{-1}}$$

According to the documentation that is the format freqz() expects: $$H(z) = \frac{b_0 + b_1z^{-1} + b_2z^{-2} + \dots + b_Nz^{-N}}{a_0 + a_1z^{-1} + a_2z^{-2} + \dots + a_Mz^{-M}}$$

but TransferFunction() expects the form with positive powers of z: $$H(z) = \frac{b_0z^{N} + b_1z^{N-1} + \dots + b_{N-1}z + b_N}{a_0z^{M} + a_1z^{M-1} + \dots + a_{M-1}z + a_M}$$

which I get by just multiplying by $$\frac{z}{z}$$: $$H(z) = \frac{b}{1-az^{-1}}\cdot\frac{z}{z} = \frac{bz}{z-a}$$

Luckily in this case the resulting coefficient arrays are identical ([b], [1, -a]), so that is not my problem.

Here is a short script to illustrate what I am talking about:

import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import freqz, TransferFunction, dfreqresp, bode

def plot_freq_response(title, w, h, fs):
"Plots the frequency response with magnitude in dB and phase in degrees"
# Convert frequency, magnitude, and phase
f = w * fs / (2*np.pi)
h_db = 20 * np.log10(np.abs(h))
h_phase = np.angle(h) * 180 / np.pi

# Plot info about the cutoff frequency
cutoff_idx = abs(h_db-(-3)).argmin()
cutoff_freq = f[cutoff_idx]
cutoff_mag = h_db[cutoff_idx]
cutoff_phase = h_phase[cutoff_idx]
cutoff_freq_str = f"Cutoff = {cutoff_freq:.2f} Hz"
cutoff_phase_str = f"Phase = {cutoff_phase:.2f} degrees"

# Magnitude
plt.figure(figsize=(10,7))
plt.subplot(2,1,1)
plt.title(title)
plt.xlabel("Frequency [Hz]")
plt.ylabel("Gain [dB]")
plt.semilogx(f, h_db)
plt.scatter(cutoff_freq, cutoff_mag, color='red', label=cutoff_freq_str)
plt.grid()
plt.legend()

# Phase
plt.subplot(2,1,2)
plt.ylabel("Phase [degrees]")
plt.xlabel("Frequency [Hz]")
plt.semilogx(f, h_phase)
plt.scatter(cutoff_freq, cutoff_phase, color='red', label=cutoff_phase_str)
plt.grid()
plt.legend()

plt.show()

# Filter parameters:
fs = 1000       # Sampling rate (Hz)
a = 0.75        # IIR filter decay value
b = 1 - a

# Use the same frequency range for both plots to make them comparable
w = np.linspace(0, np.pi, 1024)

# Generate the frequency response using freqz():
w, h_freqz = freqz([b], [1, -a], worN=w)

plot_freq_response(f"Low-Pass Frequency Response with decay={a:.2f} using freqz",w, h_freqz, fs)

# Now do the same with dfreqresp() which uses dlti() which uses TransferFunction():
s = TransferFunction([b], [1, -a], dt=1/fs)
w, h_dfreqresp = dfreqresp(s, w)

plot_freq_response(f"Low-Pass Frequency Response with decay={a:.2f} using dfreqresp", w, h_dfreqresp, fs)

# Plot Magnitude and phase errors
plt.figure(figsize=(10,4))
plt.xlabel("Frequency [Hz]")
plt.ylabel("Error")
f = w * fs / (2*np.pi)
mag_error = np.abs(np.abs(h_freqz) - np.abs(h_dfreqresp))
phase_error = np.abs(np.angle(h_freqz) - np.angle(h_dfreqresp)) * 180/np.pi
plt.plot(f, mag_error, label="Magnitude Error")
plt.plot(f, phase_error, label="Phase Error (degrees)")
plt.legend()
plt.show()


Here are the plots generated by the code above:   The phase of the freqz() plot clearly wraps back up to 0 degrees, while the dfreqresp() plot does not.

One interesting finding is that the magnitudes are the same (within error), but the difference between the phase in each plot is linear from 0 to 180 degrees as frequency goes from 0 to $$\frac{fs}{2}$$.

Can someone please explain why this is the case? I have run similar comparisons using both Octave and Matlab, and they behave the same way, so I think the error here is in my understanding of the differences between these two functions.

• Saluton DrEsperanto! Mi estas denaske Esperantisto. Bonvenon al DSP.SE! Feb 10, 2022 at 12:17
• @DanBoschen Ĉu vere? Dankon pro la bonveno! Mi neniam divenus ke mi hazarde renkontus denaskulon en forumo, kiu ne temas pri Esperanto mem ;) Mojose! Feb 10, 2022 at 19:27
• Mia patro: scua.library.umass.edu/boschen-allan-c Li parolis al mi nur en Esperanto gxis mi estis 13-jara, sed nun mi ne multe parolas gxin. Feb 10, 2022 at 19:39
• @Dan Boschen La mondo vere estas malgranda. Mi pensas ke vi renkontis mian koramikinon antaŭ kelkaj jaroj ĉe antaŭfunebra ceremonio en Beverly, MA. Ŝi rekonis vian verdan stelon. Kia eta mondo. Feb 13, 2022 at 3:34
• Ĉu vere! Jes ja, mi memoras tion kaj ŝin. La filino de Rikardo. Tio estas mirinda. Feb 13, 2022 at 6:32

Let's use $$H_1(z)$$ to denote the transfer function with negative powers of $$z$$, and $$H_2(z)$$ is the one with positive powers of $$z$$. If $$M=N$$ we have $$H_1(z)=H_2(z)$$. However, in general we have
$$H_2(z)=z^{N-M}H_1(z)\tag{1}$$
$$H_2(e^{j\omega})=e^{j\omega (N-M)}H_1(e^{j\omega})\tag{2}$$
I.e., their magnitudes are identical, and their phases differ by a linear phase term. In your example, $$N=0$$ and $$M=1$$, so $$H_2(z)$$ has an additional linear phase term $$e^{-j\omega}$$.