周波数応答
離散時間で表されたインパルス応答から周波数応答を求める方法。
インパルス応答が
$$
h[n] =
\left\{
\begin{array}{ll}
\frac{1}{5}, & n = 0,1,2,3,4 \\
0, & \text{otherwise}
\end{array}
\right.
$$
であるようなシステムを考える。
伝達関数は、z変換を適用して
$$
\begin{eqnarray}
H(z) & = & \sum_{k=0}^4 \frac{1}{5} z^{-k} \
& = & \frac{1}{5} \frac{1-z^{-5}}{1-z^{-1}}
\end{eqnarray}
$$
となる。周波数応答は、$z=e^{j \omega}$ を代入して
$$ \begin{eqnarray} H(e^{j \omega}) & = & \frac{1}{5} \frac{1-e^{-j5\omega}}{1-e^{-j\omega}} \end{eqnarray} $$
絶対値をプロットすると以下のようになる。
ソースコード(Python)
import matplotlib.pyplot as plt
import numpy as np
omg = np.linspace(-np.pi, np.pi, 100)
alpha = 0.8
z = np.exp(1j * omg)
#H = 1 / (1 - alpha*(1/z))
H = (1/5) * (1 - z**(-5)) / (1- z**(-1))
print(plt.style.available)
plt.style.use('default') # use default style
fig,ax=plt.subplots()
ax.set_facecolor("white") # background color
ax.grid(color='lightgray') # grid color
#ax.spines["right"].set_color("blue") # right color = blue
#ax.spines["right"].set_color("none") # erase right line
#ax.spines["top"].set_color("none") # erase up line
#ax.spines["left"].set_color("m") # left color = magenta
#ax.spines["bottom"].set_color("c") # bottom color = cyan
plt.plot(omg, np.abs(H))
plt.xlim([-np.pi,np.pi])
plt.ylim([0,1])
plt.xticks([-np.pi,-0.5*np.pi,0,0.5*np.pi,np.pi],[r'$-\pi$',r'$-\frac{\pi}{2}$','0',r'$\frac{\pi}{2}$',r'$\pi$'])
plt.yticks([0.0, 0.2, 0.4, 0.6, 0.8, 1.0],['0.0', '0.2', '0.4', '0.6', '0.8', '1.0'])
plt.xlabel(r'$\omega$')
plt.ylabel(r'$|H(e^{j \omega})|$')
plt.title(r'Frequency Reponse of $H(z)$')
plt.savefig("fig1.png")
plt.savefig("fig1.eps")
plt.figure()
plt.plot(omg, 20*np.log10(np.abs(H)))
plt.grid(True)
#plt.yscale("log")
plt.savefig("fig2.png")
ソースコード(Mathematica)
h[n_] := (1/5)*(DiscreteDelta[n] + DiscreteDelta[n-1] + DiscreteDelta[n-2] + DiscreteDelta[n-3] + DiscreteDelta[n-4])
H[z_] := ZTransform[h[n], n, z]
H[w] := H[z] /. z-> E^(I w)
Plot[Abs[H[w]], {w, -Pi,Pi}, PlotRange->{0,1}]
