Merge pull request #16 from 153hashimoto/khashimoto/chapter03

Add chapter03

Author: taishi-n

khashimoto/chapter03/q01.py

@@ -0,0 +1,15 @@
+# 線形畳み込み
+
+import numpy as np
+
+
+def LinearConv(x, h):
+    N = x.size
+    z = np.zeros(2 * N - 1)
+    for n in range(2 * N - 1):
+        for k in range(N):
+            if (n - k) < 0 or (n - k) > (N - 1):
+                z[n] += 0
+            else:
+                z[n] += x[k] * h[n - k]
+    return z

khashimoto/chapter03/q02.py

@@ -0,0 +1,11 @@
+# 巡回畳み込み
+
+import numpy as np
+
+
+def CircularConv(x, h):
+    N = x.size
+    z = np.zeros(N)
+    for k in range(N):
+        z += x[k] * h[(np.arange(N) - k) % N]
+    return z

khashimoto/chapter03/q03.py

@@ -0,0 +1,14 @@
+# 零詰め+巡回畳み込み
+
+import numpy as np
+
+
+def ZeroPaddingCircularConv(x, h):
+    N1 = x.size
+    x = np.hstack((x, np.zeros(N1 - 1)))
+    h = np.hstack((h, np.zeros(N1 - 1)))
+    N = 2 * N1 - 1
+    z = np.zeros(N)
+    for k in range(N):
+        z += x[k] * h[(np.arange(N) - k) % N]
+    return z

khashimoto/chapter03/q04.py

@@ -0,0 +1,15 @@
+# 各種畳み込みの関係
+
+import numpy as np
+from q01 import LinearConv
+from q02 import CircularConv
+from q03 import ZeroPaddingCircularConv
+
+x = np.array([4, 3, 2, 1])
+h = np.array([1, 0, -1, 0])
+print("線形畳み込み:")
+print(LinearConv(x, h))
+print("巡回畳み込み:")
+print(CircularConv(x, h))
+print("零詰めを行った巡回畳み込み:")
+print(ZeroPaddingCircularConv(x, h))

khashimoto/chapter03/q05.py

@@ -0,0 +1,12 @@
+# 差分方程式（再帰なし）
+
+import numpy as np
+import matplotlib.pyplot as plt
+
+x = np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0])
+y = np.zeros(10)
+n = np.arange(10)
+y[n] = 0.2 * x[n] + 0.2 * x[n - 1] + 0.2 * x[n - 2] + 0.2 * x[n - 3] + 0.2 * x[n - 4]
+
+plt.stem(y)
+plt.show()

khashimoto/chapter03/q06.py

@@ -0,0 +1,12 @@
+# 差分方程式（再帰あり）
+
+import numpy as np
+import matplotlib.pyplot as plt
+
+x = np.array([1, 0, 0, 0, 0, 0, 0, 0, 0, 0])
+y = np.zeros(10)
+for n in range(10):
+    y[n] = 0.3 * y[n - 1] + 0.4 * x[n]
+
+plt.stem(y)
+plt.show()

khashimoto/chapter03/q07.py

@@ -0,0 +1,17 @@
+# 差分方程式（一般系）
+
+import numpy as np
+
+
+def DifferenceEquation(a, b, x):
+    N = a.size - 1
+    M = b.size - 1
+    L = x.size
+
+    y = np.zeros(L)
+    ka = np.arange(1, N + 1)
+    kb = np.arange(M + 1)
+    for n in range(L):
+        y[n] = (-sum(a[ka] * y[n - ka]) + sum(b[kb] * x[n - kb])) / a[0]
+
+    return y

khashimoto/chapter03/q08.py

@@ -0,0 +1,15 @@
+# 周波数応答
+
+import numpy as np
+
+
+def FreqRes(a, b, fs, f):
+    N = a.size - 1
+    M = b.size - 1
+
+    ka = np.arange(1, N + 1)
+    kb = np.arange(M + 1)
+    H = sum(b[kb] * np.exp(-1j * 2 * np.pi * f / fs * kb)) / (
+        1 + sum(a[ka] * np.exp(-1j * 2 * np.pi * f / fs * ka))
+    )
+    return H

khashimoto/chapter03/q09.py

@@ -0,0 +1,22 @@
+# 周波数応答の確認（再帰なし）
+
+import numpy as np
+import matplotlib.pyplot as plt
+from q08 import FreqRes
+
+fs = 16000
+N = 10
+f = np.arange(0, fs, 1 / N)
+
+a = np.array([0])
+b = np.array([0.2, 0.2, 0.2, 0.2, 0.2])
+H = np.zeros(f.size, dtype=complex)
+for i in range(f.size):
+    H[i] = FreqRes(a, b, fs, f[i])
+
+plt.subplot(2, 1, 1)
+plt.plot(abs(H))
+
+plt.subplot(2, 1, 2)
+plt.plot(np.angle(H))
+plt.show()

khashimoto/chapter03/q10.py

@@ -0,0 +1,22 @@
+# 周波数応答の確認（再帰あり）
+
+import numpy as np
+import matplotlib.pyplot as plt
+from q08 import FreqRes
+
+fs = 16000
+N = 10
+f = np.arange(0, fs, 1 / N)
+
+a = np.array([0.3])
+b = np.array([0.4])
+H = np.zeros(f.size, dtype=complex)
+for i in range(f.size):
+    H[i] = FreqRes(a, b, fs, f[i])
+
+plt.subplot(2, 1, 1)
+plt.plot(abs(H))
+
+plt.subplot(2, 1, 2)
+plt.plot(np.angle(H))
+plt.show()