1. 信号とシステムの基礎
1.1 信号の分類
信号(Signal)とは情報を伝達する物理量の変化です。主な分類を以下に示します。
**連続時間信号 vs 離散時間信号**
連続時間(CT)信号はすべての $t \in \mathbb{R}$ で定義されます。例: アナログ音声信号 $x(t) = A\cos(\omega_0 t + \phi)$
離散時間(DT)信号は整数インデックス $n \in \mathbb{Z}$ でのみ定義されます。例: デジタル音声サンプル列 $x[n]$
**エネルギー信号と電力信号**
信号のエネルギー $E$ と平均電力 $P$:
$$E = \int_{-\infty}^{\infty} |x(t)|^2 \, dt$$
$$P = \lim_{T \to \infty} \frac{1}{2T} \int_{-T}^{T} |x(t)|^2 \, dt$$
$E < \infty$ ならエネルギー信号、$0 < P < \infty$ なら電力信号です。
1.2 基本信号
**単位インパルス(Unit Impulse)**
$$\delta(t) = \begin{cases} \infty & t = 0 \\ 0 & t \neq 0 \end{cases}, \quad \int_{-\infty}^{\infty} \delta(t) \, dt = 1$$
サンプリング性質: $\int_{-\infty}^{\infty} x(t)\delta(t - t_0) \, dt = x(t_0)$
**単位ステップ(Unit Step)**
$$u(t) = \begin{cases} 1 & t \geq 0 \\ 0 & t < 0 \end{cases}$$
関係式 $u(t) = \int_{-\infty}^{t} \delta(\tau) \, d\tau$ が成立します。
1.3 線形時不変(LTI)システム
LTIシステムは2つの核心的な性質を満たします。
- **線形性(Linearity)**: $\mathcal{T}\{ax_1(t) + bx_2(t)\} = a\mathcal{T}\{x_1(t)\} + b\mathcal{T}\{x_2(t)\}$
- **時不変性(Time Invariance)**: $y(t) = \mathcal{T}\{x(t)\} \Rightarrow \mathcal{T}\{x(t-t_0)\} = y(t-t_0)$
1.4 畳み込み(Convolution)
LTIシステムの出力は入力とインパルス応答の畳み込みで表されます。
$$y(t) = x(t) * h(t) = \int_{-\infty}^{\infty} x(\tau) h(t - \tau) \, d\tau$$
離散時間の場合:
$$y[n] = x[n] * h[n] = \sum_{k=-\infty}^{\infty} x[k] h[n - k]$$
畳み込みは交換律・結合律・分配律を満たします。
2. フーリエ級数 (Fourier Series)
2.1 周期信号のフーリエ級数
周期 $T_0$ を持つ信号 $x(t)$ は基本周波数 $\omega_0 = 2\pi/T_0$ の複素指数関数の和で表せます。
$$x(t) = \sum_{k=-\infty}^{\infty} c_k e^{jk\omega_0 t}$$
**フーリエ係数**:
$$c_k = \frac{1}{T_0} \int_{T_0} x(t) e^{-jk\omega_0 t} \, dt$$
2.2 矩形波(Square Wave)の例
振幅 A、周期 T、デューティサイクル 50% の矩形波のフーリエ係数:
$$c_k = \begin{cases} A/2 & k = 0 \\ \frac{A}{k\pi} \sin\left(\frac{k\pi}{2}\right) & k \neq 0 \end{cases}$$
奇数次高調波のみが存在し、高次になるほど振幅は減少します。不連続点付近ではギブス現象(Gibbs Phenomenon)が現れます。
2.3 パーセバルの定理
$$\frac{1}{T_0} \int_{T_0} |x(t)|^2 \, dt = \sum_{k=-\infty}^{\infty} |c_k|^2$$
時間領域の平均電力と周波数領域の係数の二乗和が等しくなります。
3. フーリエ変換 (Fourier Transform)
3.1 連続時間フーリエ変換 (CTFT)
非周期信号のフーリエ変換ペア:
$$X(j\omega) = \int_{-\infty}^{\infty} x(t) e^{-j\omega t} \, dt \quad \text{(正変換)}$$
$$x(t) = \frac{1}{2\pi} \int_{-\infty}^{\infty} X(j\omega) e^{j\omega t} \, d\omega \quad \text{(逆変換)}$$
3.2 重要な変換ペア
| 信号 | フーリエ変換 |
| --------------------- | ------------------------------------------------------ |
| $\delta(t)$ | $1$ |
| $1$ | $2\pi\delta(\omega)$ |
| $e^{j\omega_0 t}$ | $2\pi\delta(\omega - \omega_0)$ |
| $\cos(\omega_0 t)$ | $\pi[\delta(\omega-\omega_0)+\delta(\omega+\omega_0)]$ |
| $u(t)$ | $\pi\delta(\omega) + 1/(j\omega)$ |
| $e^{-at}u(t), \, a>0$ | $1/(a+j\omega)$ |
| $\text{rect}(t/\tau)$ | $\tau \cdot \text{sinc}(\omega\tau/2)$ |
3.3 フーリエ変換の性質
- **線形性**: $ax_1(t) + bx_2(t) \leftrightarrow aX_1(j\omega) + bX_2(j\omega)$
- **時間シフト**: $x(t-t_0) \leftrightarrow e^{-j\omega t_0} X(j\omega)$
- **周波数シフト(変調)**: $e^{j\omega_0 t} x(t) \leftrightarrow X(j(\omega - \omega_0))$
- **畳み込み定理**: $x(t)*h(t) \leftrightarrow X(j\omega)H(j\omega)$
- **乗算定理**: $x(t)h(t) \leftrightarrow \frac{1}{2\pi}X(j\omega)*H(j\omega)$
- **微分**: $\frac{dx}{dt} \leftrightarrow j\omega X(j\omega)$
- **パーセバル**: $\int|x(t)|^2 dt = \frac{1}{2\pi}\int|X(j\omega)|^2 d\omega$
3.4 PythonによるFFT計算
from scipy.fft import fft, ifft, fftfreq, fftshift
矩形波信号の生成
fs = 1000 # サンプリング周波数 (Hz)
T = 1.0 # 信号長 (秒)
t = np.linspace(0, T, int(fs * T), endpoint=False)
from scipy.signal import square
x = square(2 * np.pi * 10 * t) # 10Hz 矩形波
FFT計算
N = len(x)
X = fft(x)
freqs = fftfreq(N, 1/fs)
両側スペクトル (正規化)
magnitude = np.abs(fftshift(X)) / N
freq_shifted = fftshift(freqs)
plt.figure(figsize=(12, 4))
plt.subplot(1, 2, 1)
plt.plot(t[:100], x[:100])
plt.title('Square Wave (time domain)')
plt.xlabel('Time (s)')
plt.subplot(1, 2, 2)
plt.plot(freq_shifted, magnitude)
plt.xlim(-100, 100)
plt.title('Fourier Transform (frequency domain)')
plt.xlabel('Frequency (Hz)')
plt.tight_layout()
plt.show()
4. 離散フーリエ変換 (DFT) と FFT
4.1 離散時間フーリエ変換 (DTFT)
離散時間信号 $x[n]$ の DTFT:
$$X(e^{j\omega}) = \sum_{n=-\infty}^{\infty} x[n] e^{-j\omega n}$$
逆変換:
$$x[n] = \frac{1}{2\pi} \int_{-\pi}^{\pi} X(e^{j\omega}) e^{j\omega n} \, d\omega$$
DTFTは連続な周波数 $\omega$ を持ち、$2\pi$ 周期です。
4.2 DFTの定義
N点DFTはDTFTをN個の周波数で等間隔サンプリングしたものです。
$$X[k] = \sum_{n=0}^{N-1} x[n] e^{-j2\pi kn/N}, \quad k = 0, 1, \ldots, N-1$$
逆DFT:
$$x[n] = \frac{1}{N} \sum_{k=0}^{N-1} X[k] e^{j2\pi kn/N}, \quad n = 0, 1, \ldots, N-1$$
4.3 FFTアルゴリズム
直接DFT計算の計算量は $O(N^2)$ ですが、FFT(クーリー・テューキーアルゴリズム)は $N = 2^m$ のとき分割統治法で $O(N \log N)$ を実現します。
旋回因子(Twiddle Factor): $W_N^k = e^{-j2\pi k/N}$
4.4 スペクトル解析の例
from scipy.fft import fft, fftfreq
複合信号の生成: 50Hz + 120Hz
fs = 1000 # サンプリング周波数 (Hz)
N = fs # サンプル数 (1秒)
t = np.linspace(0, 1, N, endpoint=False)
ノイズ付き複合信号
signal = (np.sin(2 * np.pi * 50 * t)
+ 0.5 * np.sin(2 * np.pi * 120 * t)
+ 0.2 * np.random.randn(N))
FFT計算
yf = fft(signal)
xf = fftfreq(N, 1/fs)
片側スペクトル (正の周波数のみ)
half = N // 2
amplitude = 2.0 / N * np.abs(yf[:half])
freq_pos = xf[:half]
plt.figure(figsize=(10, 4))
plt.plot(freq_pos, amplitude)
plt.xlabel('Frequency (Hz)')
plt.ylabel('Amplitude')
plt.title('Single-Sided Amplitude Spectrum')
plt.axvline(x=50, color='r', linestyle='--', label='50 Hz')
plt.axvline(x=120, color='g', linestyle='--', label='120 Hz')
plt.legend()
plt.grid(True)
plt.show()
print(f"周波数分解能: {fs/N} Hz")
print(f"ナイキスト周波数: {fs/2} Hz")
4.5 ナイキスト・シャノンのサンプリング定理
エイリアシング(折り返し歪み)なしに信号を復元するには、サンプリング周波数 $f_s$ が信号の最大周波数 $f_{max}$ の2倍以上必要です。
$$f_s \geq 2 f_{max} \quad \text{(ナイキスト条件)}$$
5. ラプラス変換 (Laplace Transform)
5.1 定義
両側ラプラス変換:
$$X(s) = \int_{-\infty}^{\infty} x(t) e^{-st} \, dt, \quad s = \sigma + j\omega$$
片側ラプラス変換:
$$X(s) = \int_{0^-}^{\infty} x(t) e^{-st} \, dt$$
5.2 重要な変換ペア
| 信号 $x(t)$ | 変換 $X(s)$ | 収束領域 |
| ---------------------- | --------------------------- | ------------------- |
| $\delta(t)$ | $1$ | 全s平面 |
| $u(t)$ | $1/s$ | $\text{Re}(s) > 0$ |
| $e^{-at}u(t)$ | $1/(s+a)$ | $\text{Re}(s) > -a$ |
| $te^{-at}u(t)$ | $1/(s+a)^2$ | $\text{Re}(s) > -a$ |
| $\sin(\omega_0 t)u(t)$ | $\omega_0/(s^2+\omega_0^2)$ | $\text{Re}(s) > 0$ |
| $\cos(\omega_0 t)u(t)$ | $s/(s^2+\omega_0^2)$ | $\text{Re}(s) > 0$ |
5.3 伝達関数 (Transfer Function)
LTIシステムの入出力関係をs領域で表現:
$$H(s) = \frac{Y(s)}{X(s)} = \frac{b_m s^m + b_{m-1} s^{m-1} + \cdots + b_0}{a_n s^n + a_{n-1} s^{n-1} + \cdots + a_0}$$
**極(Poles)**: $H(s) \to \infty$ となる $s$ の値(分母 = 0)
**零点(Zeros)**: $H(s) = 0$ となる $s$ の値(分子 = 0)
すべての極がs平面の左半面(LHP)にあればBIBO安定です。
5.4 2次標準系
$$H(s) = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2}$$
- $\omega_n$: 固有角周波数(Natural Frequency)
- $\zeta$: 減衰比(Damping Ratio)
- $\zeta < 1$: 不足減衰(Underdamped) — 振動応答
- $\zeta = 1$: 臨界減衰(Critically Damped) — 最速無振動
- $\zeta > 1$: 過減衰(Overdamped) — 緩慢無振動
6. Z変換
6.1 定義
離散時間信号 $x[n]$ のZ変換:
$$X(z) = \sum_{n=-\infty}^{\infty} x[n] z^{-n}, \quad z \in \mathbb{C}$$
$z = e^{j\omega}$ (単位円上)と置換するとDTFTに一致します。
6.2 Z変換とラプラス変換の関係
連続時間システムをサンプリング周期 $T_s$ で離散化すると:
$$z = e^{sT_s}$$
この写像によりs平面とz平面が対応します。
- s平面の虚軸 $\leftrightarrow$ z平面の単位円
- s平面の左半面 $\leftrightarrow$ z平面の単位円内部
6.3 安定性解析
離散時間システムが安定であるためには、すべての極が単位円の内部になければなりません。
$$|z_i| < 1 \quad \text{すべての極に対して}$$
6.4 重要なZ変換ペア
| 信号 $x[n]$ | Z変換 $X(z)$ |
| ---------------------- | ------------------------------------------ |
| $\delta[n]$ | $1$ |
| $u[n]$ | $z/(z-1)$ |
| $a^n u[n]$ | $z/(z-a)$ |
| $n \cdot a^n u[n]$ | $az/(z-a)^2$ |
| $\cos(\omega_0 n)u[n]$ | $z(z-\cos\omega_0)/(z^2-2z\cos\omega_0+1)$ |
7. デジタルフィルタ設計
7.1 FIRフィルタ vs IIRフィルタ
**FIR(有限インパルス応答)フィルタ**
$$y[n] = \sum_{k=0}^{M} b_k x[n-k]$$
- インパルス応答が有限長
- 常に安定(フィードバックなし)
- 線形位相を実現可能
- 急峻な遮断特性には高次数が必要
**IIR(無限インパルス応答)フィルタ**
$$y[n] = \sum_{k=0}^{M} b_k x[n-k] - \sum_{k=1}^{N} a_k y[n-k]$$
- フィードバックを含み、インパルス応答は無限
- 低次数で急峻な遮断特性
- 位相が非線形(位相歪みの可能性)
- 設計不良の場合は不安定になりうる
7.2 フィルタの種類
- **低域通過フィルタ(LPF)**: 低周波を通過、高周波を遮断
- **高域通過フィルタ(HPF)**: 高周波を通過、低周波を遮断
- **帯域通過フィルタ(BPF)**: 特定の周波数帯域のみ通過
- **帯域阻止フィルタ(BEF/Notch)**: 特定の周波数帯域を遮断
7.3 古典的フィルタの種類
**バタワース(Butterworth)**: 通過帯域で最大限平坦な応答
**チェビシェフ I 型**: 通過帯域にリップルを許容し、より急峻な遮断特性
**楕円(Elliptic)**: 通過帯域と阻止帯域の両方にリップルを持ち、最小次数で最大減衰
7.4 scipy.signalによるフィルタ設計
from scipy import signal
fs = 1000 # サンプリング周波数 (Hz)
バタワース低域通過フィルタ (4次、遮断周波数100Hz)
nyq = fs / 2
cutoff = 100 / nyq # 正規化周波数 (0〜1)
b_butter, a_butter = signal.butter(N=4, Wn=cutoff, btype='low')
チェビシェフ I 型 (4次、1dBリップル、遮断周波数100Hz)
b_cheby, a_cheby = signal.cheby1(N=4, rp=1, Wn=cutoff, btype='low')
周波数応答の計算
w, h_butter = signal.freqz(b_butter, a_butter, fs=fs)
w, h_cheby = signal.freqz(b_cheby, a_cheby, fs=fs)
plt.figure(figsize=(10, 5))
plt.subplot(1, 2, 1)
plt.plot(w, 20 * np.log10(np.abs(h_butter)), label='Butterworth')
plt.plot(w, 20 * np.log10(np.abs(h_cheby)), label='Chebyshev I')
plt.xlabel('Frequency (Hz)')
plt.ylabel('Magnitude (dB)')
plt.title('Filter Frequency Response')
plt.legend()
plt.grid(True)
plt.ylim(-80, 5)
ノイズ信号のフィルタリング
t = np.linspace(0, 1, fs, endpoint=False)
clean = np.sin(2 * np.pi * 50 * t)
noisy = clean + 0.5 * np.sin(2 * np.pi * 300 * t)
filtered = signal.filtfilt(b_butter, a_butter, noisy)
plt.subplot(1, 2, 2)
plt.plot(t[:200], noisy[:200], alpha=0.5, label='Noisy')
plt.plot(t[:200], filtered[:200], label='Filtered')
plt.xlabel('Time (s)')
plt.legend()
plt.title('Filtering Result')
plt.tight_layout()
plt.show()
7.5 FIRフィルタ設計 (窓関数法)
from scipy.signal import firwin, freqz
Hamming窓を用いたFIR LPF
タップ数 = 51、遮断周波数 = 100Hz
numtaps = 51
cutoff_hz = 100
b_fir = firwin(numtaps, cutoff_hz / (fs / 2), window='hamming')
w_fir, h_fir = freqz(b_fir, [1.0], fs=fs)
print(f"FIRフィルタ次数: {numtaps - 1}")
print(f"線形位相: True (対称係数)")
8. 制御システムの基礎
8.1 開ループ vs 閉ループ制御
**開ループ(Open-Loop)制御**: 出力が制御動作に影響を与えません。シンプルですが外乱やモデル誤差に弱いです。
**閉ループ(Closed-Loop)制御**: 出力を計測して基準入力と比較し、誤差を最小化します。フィードバックにより堅牢性が高まります。
基本的なフィードバックループの伝達関数:
$$\frac{Y(s)}{R(s)} = \frac{G(s)}{1 + G(s)H(s)}$$
ここで $G(s)$ はプラント、$H(s)$ はセンサの伝達関数です。
8.2 安定性: ラウス・フルビッツ基準
特性方程式 $1 + G(s)H(s) = 0$ のすべての根がs平面の左半面にあれば安定です。
ラウス配列(Routh Array)を用いると、係数だけから極の位置を判定できます。第1列の符号変化の回数が右半面の極の数と等しくなります。
8.3 ボード線図 (Bode Plot)
周波数応答 $H(j\omega)$ を対数スケールで表現します。
- **ゲイン線図**: $20\log_{10}|H(j\omega)|$ [dB] vs $\log_{10}\omega$
- **位相線図**: $\angle H(j\omega)$ [度] vs $\log_{10}\omega$
**安定余裕**:
- **位相余裕(Phase Margin, PM)**: ゲインが0 dBのときの位相が -180度からどれだけ余裕があるか
- **ゲイン余裕(Gain Margin, GM)**: 位相が -180度のときのゲインが 0 dB からどれだけ余裕があるか
PM > 45度、GM > 6 dB が一般的に良好な設計の目安です。
from scipy import signal
2次伝達関数: H(s) = 10 / (s^2 + 2s + 10)
num = [10]
den = [1, 2, 10]
sys = signal.TransferFunction(num, den)
w, mag, phase = signal.bode(sys)
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(8, 6))
ax1.semilogx(w, mag)
ax1.set_ylabel('Magnitude (dB)')
ax1.set_title('Bode Plot')
ax1.grid(True, which='both')
ax2.semilogx(w, phase)
ax2.set_ylabel('Phase (degrees)')
ax2.set_xlabel('Frequency (rad/s)')
ax2.grid(True, which='both')
plt.tight_layout()
plt.show()
9. PID制御器
9.1 PID制御動作
PID制御器は3つの制御動作の和です。
$$u(t) = K_p e(t) + K_i \int_0^t e(\tau) d\tau + K_d \frac{de(t)}{dt}$$
ここで $e(t) = r(t) - y(t)$ は誤差信号です。
- **比例(P)制御**: 現在の誤差に比例した制御出力。速い応答が得られますが定常偏差が残ります。
- **積分(I)制御**: 累積誤差を解消して定常偏差をゼロに。応答が遅くなります。
- **微分(D)制御**: 誤差の変化率に反応してオーバーシュートを抑制。ノイズに敏感。
伝達関数:
$$C(s) = K_p + \frac{K_i}{s} + K_d s = K_p\left(1 + \frac{1}{T_i s} + T_d s\right)$$
9.2 ジーグラー・ニコルス調整法
**方法1(開ループステップ応答)**:
プロセス反応曲線から無駄時間 $L$ と時定数 $T$ を求めます。
| 制御器 | $K_p$ | $T_i$ | $T_d$ |
| ------ | -------- | ----- | ------ |
| P | $T/L$ | - | - |
| PI | $0.9T/L$ | $3L$ | - |
| PID | $1.2T/L$ | $2L$ | $0.5L$ |
**方法2(閉ループ限界感度法)**:
$K_i = K_d = 0$ として $K_p$ を増加させ、持続振動が起きる限界ゲイン $K_u$ と振動周期 $T_u$ を求めます。
| 制御器 | $K_p$ | $T_i$ | $T_d$ |
| ------ | -------- | -------- | ---------- |
| PID | $0.6K_u$ | $0.5T_u$ | $0.125T_u$ |
9.3 PythonによるPIDシミュレーション
class PIDController:
def __init__(self, kp, ki, kd, dt):
self.kp = kp
self.ki = ki
self.kd = kd
self.dt = dt
self.integral = 0.0
self.prev_error = 0.0
def update(self, setpoint, measured):
error = setpoint - measured
self.integral += error * self.dt
derivative = (error - self.prev_error) / self.dt
output = self.kp * error + self.ki * self.integral + self.kd * derivative
self.prev_error = error
return output
1次プラント: dy/dt = -y + u (時定数1秒)
def plant_step(y, u, dt):
return y + dt * (-y + u)
dt = 0.01
t = np.arange(0, 10, dt)
setpoint = 1.0
pid = PIDController(kp=2.0, ki=1.0, kd=0.5, dt=dt)
y = 0.0
y_list = []
for _ in t:
u = pid.update(setpoint, y)
u = np.clip(u, -10, 10) # アクチュエータ飽和
y = plant_step(y, u, dt)
y_list.append(y)
plt.figure(figsize=(10, 4))
plt.plot(t, y_list, label='Plant Output')
plt.axhline(setpoint, color='r', linestyle='--', label='Setpoint')
plt.xlabel('Time (s)')
plt.ylabel('Output')
plt.title('PID Control Simulation')
plt.legend()
plt.grid(True)
plt.show()
10. 状態空間表現 (State Space)
10.1 状態方程式
連続時間LTIシステムの状態空間表現:
$$\dot{\mathbf{x}}(t) = \mathbf{A}\mathbf{x}(t) + \mathbf{B}\mathbf{u}(t)$$
$$\mathbf{y}(t) = \mathbf{C}\mathbf{x}(t) + \mathbf{D}\mathbf{u}(t)$$
- $\mathbf{x}$: 状態ベクトル ($n \times 1$)
- $\mathbf{u}$: 入力ベクトル ($m \times 1$)
- $\mathbf{y}$: 出力ベクトル ($p \times 1$)
- $\mathbf{A}$: システム行列 ($n \times n$)
- $\mathbf{B}$: 入力行列 ($n \times m$)
- $\mathbf{C}$: 出力行列 ($p \times n$)
- $\mathbf{D}$: 直達行列 ($p \times m$)
伝達関数との関係: $H(s) = \mathbf{C}(s\mathbf{I} - \mathbf{A})^{-1}\mathbf{B} + \mathbf{D}$
10.2 可制御性と可観測性
**可制御性(Controllability)**: 任意の初期状態から有限時間内に任意の目標状態へ移行できるか？
可制御性行列 $\mathcal{C} = [\mathbf{B} \; \mathbf{AB} \; \mathbf{A}^2\mathbf{B} \; \cdots \; \mathbf{A}^{n-1}\mathbf{B}]$ のランクが $n$ なら完全可制御。
**可観測性(Observability)**: 出力だけから初期状態を決定できるか？
可観測性行列 $\mathcal{O} = [\mathbf{C}^T \; (\mathbf{CA})^T \; \cdots \; (\mathbf{CA}^{n-1})^T]^T$ のランクが $n$ なら完全可観測。
10.3 線形二次レギュレータ (LQR)
二次コスト関数を最小化する最適状態フィードバック制御:
$$J = \int_0^{\infty} (\mathbf{x}^T \mathbf{Q} \mathbf{x} + \mathbf{u}^T \mathbf{R} \mathbf{u}) \, dt$$
最適制御入力: $\mathbf{u}^* = -\mathbf{K}\mathbf{x}$
$\mathbf{Q}$ は状態誤差の重み、$\mathbf{R}$ は制御入力エネルギーの重みです。
from scipy import linalg
二重積分器: x_ddot = u
A = np.array([[0, 1], [0, 0]])
B = np.array([[0], [1]])
LQR重み: 位置誤差を重視
Q = np.diag([10, 1])
R = np.array([[1]])
連続時間代数リカッチ方程式を解く
P = linalg.solve_continuous_are(A, B, Q, R)
K = np.linalg.inv(R) @ B.T @ P
print("LQRゲインK:", K)
print("閉ループ極:", np.linalg.eigvals(A - B @ K))
11. AI/MLと信号処理
11.1 1D CNNによる信号分類
1次元畳み込みニューラルネットワークはECG、地震波、音声などの時系列信号のパターン認識に優れています。
class Signal1DCNN(nn.Module):
def __init__(self, num_classes=5):
super().__init__()
self.conv_layers = nn.Sequential(
nn.Conv1d(1, 32, kernel_size=7, padding=3),
nn.BatchNorm1d(32),
nn.ReLU(),
nn.MaxPool1d(2),
nn.Conv1d(32, 64, kernel_size=5, padding=2),
nn.BatchNorm1d(64),
nn.ReLU(),
nn.MaxPool1d(2),
nn.Conv1d(64, 128, kernel_size=3, padding=1),
nn.ReLU(),
nn.AdaptiveAvgPool1d(8)
)
self.classifier = nn.Sequential(
nn.Flatten(),
nn.Linear(128 * 8, 256),
nn.ReLU(),
nn.Dropout(0.5),
nn.Linear(256, num_classes)
)
def forward(self, x):
x = self.conv_layers(x)
return self.classifier(x)
使用例
model = Signal1DCNN(num_classes=5)
入力形状: (batch_size, channels=1, sequence_length=1000)
x = torch.randn(8, 1, 1000)
output = model(x)
print("出力形状:", output.shape) # (8, 5)
11.2 STFTとメルスペクトログラム
音声AIシステムでは信号を時間-周波数表現に変換してからCNNで処理する手法が広く使われています。
STFT (Short-Time Fourier Transform)
D = librosa.stft(y)
S_db = librosa.amplitude_to_db(np.abs(D), ref=np.max)
メルスペクトログラム
mel_spec = librosa.feature.melspectrogram(y=y, sr=sr, n_mels=128)
mel_db = librosa.power_to_db(mel_spec, ref=np.max)
MFCC (メル周波数ケプストラム係数)
mfcc = librosa.feature.mfcc(y=y, sr=sr, n_mfcc=13)
print("STFT: ハン窓, hop_length=512, n_fft=2048")
print("メルフィルタバンク: 128バンド, 0Hz〜ナイキスト")
print("MFCC: 音声認識では13次元がよく使われる")
11.3 信号処理 + ディープラーニングの応用例
- **ECG不整脈分類**: 1D CNN + LSTMによる心拍リズム異常の検出
- **異常振動検知**: FFT特徴抽出後にオートエンコーダで異常検出
- **音声認識**: MFCC + Transformerアーキテクチャ(Whisper, wav2vec 2.0)
- **地震検知**: 地震計信号のP波・S波の自動分類
- **脳波(EEG)解析**: 周波数帯域別パワースペクトルで集中度・睡眠状態を分類
12. 確認問題
**解答**: 4000 Hz
**解説**: ナイキスト・シャノンのサンプリング定理により、エイリアシングなしに表現できる最大周波数はサンプリング周波数の半分であるナイキスト周波数(4000 Hz)です。電話音声コーデック(8 kHzサンプリング)が最大4 kHzの音声を扱える理由がこれです。
**解答**: 線形位相特性と無条件安定性
**解説**: FIRフィルタはフィードバックがないため常にBIBO安定です。また係数を対称に設計すれば線形位相(linear phase)が保証され、すべての周波数成分の群遅延が同一になります。これはデータ通信など位相歪みが許容されない場合に重要です。一方IIRフィルタは低次数で急峻な遮断特性が得られるという利点があります。
**解答**: 誤差の変化率に反応してオーバーシュートを抑制しますが、測定ノイズに非常に敏感です。
**解説**: D動作は誤差が急激に変化するときに強い制御出力を生成し、システムが目標値を行き過ぎる現象(オーバーシュート)を抑制します。しかしノイズを含む測定値を微分するとノイズが増幅されます。実用的なPIDでは微分項にローパスフィルタを追加するか、出力を直接微分するDerivative on Measurement方式が使われます。
**解答**: すべての極がz平面の単位円の内側に位置すること。
**解説**: 離散時間システムの極の位置と安定性の関係: 単位円内部(絶対値1未満)なら安定、単位円上なら限界安定、単位円外部なら不安定です。これは連続時間で極がs平面の左半面にある場合に安定することと対応しています。z = exp(sT)の写像により左半面が単位円内部に対応します。
**解答**: 不安定(Unstable)になります。
**解説**: 位相余裕はゲイン交差周波数(ゲインが0 dBになる周波数)における位相が -180度からどれだけ離れているかを示します。位相余裕が負であればゲインが0 dBのときにすでに位相が -180度を超えていることを意味し、閉ループシステムは不安定になります。一般に安定で良好な応答特性を得るためにPMは45度〜60度を目標に設計します。
参考文献
1. **Oppenheim, A. V. & Schafer, R. W.** — _Discrete-Time Signal Processing_ (3rd ed.), Prentice Hall, 2010
2. **Ogata, K.** — _Modern Control Engineering_ (5th ed.), Prentice Hall, 2010
3. **Proakis, J. G. & Manolakis, D. G.** — _Digital Signal Processing_, Prentice Hall, 2007
4. **SciPy Signal Processing Documentation** — [https://docs.scipy.org/doc/scipy/reference/signal.html](https://docs.scipy.org/doc/scipy/reference/signal.html)
5. **MIT OpenCourseWare 6.003** — Signals and Systems — [https://ocw.mit.edu/courses/6-003-signals-and-systems-fall-2011/](https://ocw.mit.edu/courses/6-003-signals-and-systems-fall-2011/)
6. **NumPy FFT ドキュメント** — [https://numpy.org/doc/stable/reference/routines.fft.html](https://numpy.org/doc/stable/reference/routines.fft.html)
현재 단락 (1/378)
信号(Signal)とは情報を伝達する物理量の変化です。主な分類を以下に示します。