回答編集履歴
2
修正
test
CHANGED
|
@@ -1,6 +1,8 @@
|
|
|
1
|
-
|
|
1
|
+
データが無いので、てきとーに生成しています。
|
|
2
2
|
|
|
3
|
-
|
|
3
|
+
実空間・周波数空間の変換がミスっていました。
|
|
4
|
+
|
|
5
|
+
パラメータはダイレクトに入力して、付随して決まるものは全部計算によって出すと、ミスがすくなるかと。
|
|
4
6
|
|
|
5
7
|
|
|
6
8
|
|
1
修正
test
CHANGED
|
@@ -2,7 +2,7 @@
|
|
|
2
2
|
|
|
3
3
|
データが無いので、てきとーに生成していますが。
|
|
4
4
|
|
|
5
|
-
|
|
6
6
|
|
|
7
7
|
|
|
8
8
|
|
|
@@ -18,61 +18,69 @@
|
|
|
18
18
|
|
|
19
19
|
from scipy import signal
|
|
20
20
|
|
|
21
|
-
N = 10000
|
|
22
21
|
|
|
23
|
-
fs = 250 # サンプリング周波数
|
|
24
22
|
|
|
25
|
-
|
|
23
|
+
N = 1000
|
|
26
24
|
|
|
27
|
-
f
|
|
25
|
+
fs = 250
|
|
28
26
|
|
|
29
|
-
|
|
27
|
+
nyq = fs / 2
|
|
30
28
|
|
|
31
|
-
|
|
29
|
+
fc1 = 4.
|
|
32
30
|
|
|
33
|
-
|
|
31
|
+
fc2 = 30.
|
|
34
32
|
|
|
35
|
-
t =
|
|
33
|
+
numtaps = 255
|
|
36
34
|
|
|
35
|
+
dt = 1. / fs
|
|
36
|
+
|
|
37
|
+
t = dt*np.arange(N)
|
|
38
|
+
|
|
37
|
-
freq = np.
|
|
39
|
+
freq = np.fft.fftfreq(t.shape[-1])*fs
|
|
38
40
|
|
|
39
41
|
|
|
40
42
|
|
|
41
43
|
# フィルタ係数
|
|
42
44
|
|
|
43
|
-
b = signal.firwin(numtaps, [fc1, fc2], pass_zero=False)
|
|
45
|
+
b = signal.firwin(numtaps, [fc1, fc2], pass_zero=False, nyq=nyq)
|
|
44
46
|
|
|
45
47
|
#b = signal.firwin(numtaps, [fc1, fc2])
|
|
46
48
|
|
|
47
49
|
|
|
48
50
|
|
|
49
|
-
|
|
51
|
+
x = np.sum((np.random.random()*np.sin(i*t/(N*dt)) for i in range(N)), axis=1)
|
|
52
|
+
|
|
53
|
+
x /= np.max(np.abs(x))
|
|
50
54
|
|
|
51
55
|
# scipyによるFIRフィルタ
|
|
52
56
|
|
|
53
|
-
|
|
57
|
+
kx = signal.lfilter(b, 1, x)
|
|
54
58
|
|
|
55
59
|
|
|
56
60
|
|
|
57
|
-
F = np.fft.fft(
|
|
61
|
+
F = np.fft.fft(kx)
|
|
62
|
+
|
|
63
|
+
print(F.shape)
|
|
58
64
|
|
|
59
65
|
Amp = np.abs(F)
|
|
60
66
|
|
|
61
67
|
|
|
62
68
|
|
|
63
|
-
skip=
|
|
69
|
+
skip=1
|
|
64
70
|
|
|
65
71
|
plt.subplot(211)
|
|
66
72
|
|
|
67
|
-
plt.plot(
|
|
73
|
+
plt.plot(freq[::skip], F[::skip], marker='.', ls='none', label='t-F')
|
|
74
|
+
|
|
75
|
+
plt.xlim(0, 50)
|
|
68
76
|
|
|
69
77
|
plt.legend()
|
|
70
|
-
|
|
71
|
-
plt.xlim(0, 100)
|
|
72
78
|
|
|
73
79
|
plt.subplot(212)
|
|
74
80
|
|
|
75
81
|
plt.plot(freq[::skip], Amp[::skip], marker='.', ls='none', label='w-A')
|
|
82
|
+
|
|
83
|
+
plt.xlim(0, 50)
|
|
76
84
|
|
|
77
85
|
plt.legend()
|
|
78
86
|
|