信号処理

サウンドプログラミングで音色を変えるとき、いちばん基本になる道具がフィルターです。フィルターは、音を「加工する魔法の箱」ではありません。波形の中から、速く変化する部分を弱めたり、ゆっくり変化する部分を弱めたりする計算です。 この記事では、ローパスとハイパスに絞ります。まずは「ローパスのカットオフ周波数を、デジタルの1サンプルごとの計算へどう変換するのか」をはっきりさせます。 フィルターの音を曲で聴く フィルターは、単体で主役になるというより、音の明るさを開いたり閉じたりして曲の呼吸を作る道具です。特にサンプルを反復するダンスミュージックでは、同じループでも高域を削ると遠く、開くと一気に前へ出てきます。ここでは、その「明るさを動かす」感覚で聴きやすい曲を2つ挙げます。 Daft Punk / Discovery One More Time French house のフィルター感をつかむ入口として分かりやすい曲です。ループの明るさが前へ出たり奥へ引いたりする聴こえ方を、ローパスの「高域を削る」と結びつけて聴けます。 YouTubeで探す Amazonで探す Stardust / Single Music Sounds Better With You 短いサンプルを反復しながら、明るさと密度でグルーヴを育てる French house の代表例です。同じ素材でも、高域の残し方で近さや抜けが変わることを聴き取りやすいです。 YouTubeで探す Amazonで探す ジャケット画像はMusicBrainz / Cover Art Archive由来の画像を記事内に同梱して配信しています。 ローパスは波形をならす ローパスフィルターは、低い周波数を通し、高い周波数を落とすフィルターです。言い換えると、波形の急な変化をならします。矩形波やノコギリ波の角が丸くなり、音は明るさを失って、こもった方向へ変わります。 いちばん直感的なローパスは、移動平均です。いまの値だけでなく、少し前の値もまとめて平均します。たとえば4個の値を平均するなら、次のようになります。 $$ y[n] = \frac{x[n]+x[n-1]+x[n-2]+x[n-3]}{4} $$ $x[n]$ は入力、$y[n]$ は出力です。$n$ は「何番目の値か」を表す番号です。最後に4で割っているのは、4個の値を足して平均しているからです。3個を平均するなら3で割りますし、8個を平均するなら8で割ります。 この式は、波形の細かいギザギザをならします。細かいギザギザは高い周波数成分なので、高域が落ちます。まずはこれがローパスの出発点です。 実際のサウンドプログラミングでは、次のような形もよく使います。 $$ y[n] = ax[n] + (1-a)y[n-1] $$ これは、いまの入力 $x[n]$ と、ひとつ前の出力 $y[n-1]$ を混ぜる式です。$a$ が大きいほど入力に素早く追従し、$a$ が小さいほど前の出力を強く引きずります。前の値を引きずるほど、急な変化には追いつけません。だから高域が落ちます。 コードにすると、1サンプルごとの処理は次の形です。 let y = 0; function lowpass(x, a) { y = a * x + (1 - a) * y; return y; } ここではまだ、カットオフ周波数をどう計算するかは置いておきます。まず大事なのは、ローパスが「前後の値を混ぜて、急な変化をならす」処理だということです。 ...

はじめに 前回の記事では、RCローパスフィルタの伝達関数を計算し、ステップ応答を調べるところまでを行なった。今回はその続きとして、RCローパスフィルタの伝達関数を使って周波数応答（周波数振幅特性と周波数位相特性）を調べてみよう。 Pythonで伝達関数などを扱うにあたっては、こちらの書籍を大いに参考させてもらった。計算が複雑になりがちな制御工学の伝達関数なども、Pythonでらくらくシミュレーションできるようになるのでおすすめだ。 ステップ応答やインパルス応答などの過渡応答特性は時間関数f(t)の形であるが、周波数特性の場合は、$s$の代わりに$jω$となる。jは虚数（$j^2 = -1$）であり、ωは角周波数である$ω=2πf$。 実は伝達関数$G(s)$の代わりに$G(jω)$とすることで周波数振幅特性を調べることができる。また、$arg(G(jω))$とすることで周波数位相特性を調べることができる。 周波数振幅特性とは、さまざまな周波数の正弦波をシステムへ入力にした時、どの程度の振幅で出力されるか調べたもの である。たとえばシステムが単純なローパスフィルタなら、高周波になるほど出力は小さくるだろう。 また、 周波数位相特性とは、さまざまな周波数の正弦波をシステムへ入力にした時、どの程度の入力の位相とずれているかを調べたもの である。 これらのことが$G(jω)$や$arg(G(jω))$で調べられるわけだが、なぜ$s=jω$とするのかの理由はフーリエ変換やオイラーの公式を使って説明されるため非常に難しので本記事では省略する。しかし、 理屈を知らなくても伝達関数さえ分かっていれば簡単に周波数特性を調べられるので、本記事ではとにかく実践を大事に周波数応答の実験をしてみたいと思う。 具体的にはPythonでプログラミングして周波数特性をグラフ化していく。 周波数応答 さて、RCローパスフィルタの伝達関数は次の通りであった。 $$G(s)=\frac{1}{1+RCs}$$ よって、周波数振幅特性を調べるためには $$G(jω)=\frac{1}{1+jωRC}$$ とすれば良いことになる。 また、周波数位相特性を調べるためには$arg(G(jω))$であった。argとは複素数平面上での偏角のことである。 さて、$ω=2πf$であるので調べたい周波数$f$を変えて計算すれば良いわけだが、周波数は膨大にあるため手動で計算するのはほぼ不可能だ。そこでPythonを使って周波数応答をグラフ化していく。 次のプログラムは抵抗R=1kΩ、コンデンサの容量を1μF、入力のサイン波の電圧を1としたときの周波数応答をグラフ化させるプログラムである。 from control.matlab import * import matplotlib.pyplot as plt if __name__ == '__main__': R = 1 * pow(10, 3) # 1kΩ C = 1 * pow(10, -6) # 1μF K = 1 # ゲイン T = R * C # 時定数 G = tf([0, K], [T, 1]) # 伝達関数 W = logspace(1, 5) # 対数スケールの配列 gain, phase, w = bode(G, W, Hz=True) plt.show() このプログラムを実行すると、図のように周波数振幅特性（上）と周波数位相特性（下）が表示される。システムはローパスフィルタなので、周波数が高くなるにつれてゲインが下っていく。また、周波数が高くなるにつれて入力との位相ズレが出力にあらわれていく。 プログラムの説明 ここで、先ほどのプログラムをもう少し詳しく解説する。 ...

はじめに この記事では、ステップ入力によるRCローパスフィルタの過渡応答を計算していく。 上の回路図のように抵抗Rと静電容量Cのコンデンサでローパスフィルタを組む。そしてこの回路に(Eu(t))のステップ電圧を入力した時の出力電圧の様子を調べていく。 ステップ入力によるRCローパスフィルタの過渡応答 ステップ電圧は、(t \geqq 0)で電圧が(E)ボルトとなるような信号のことである。ちょうど(t=0)で、電圧のスイッチをオンにしたと思ってもらえれば良い。 また、出力電圧はコンデンサの両端の電圧(V_c)である。 ここで電圧(V_c)は $$V_c(t)=\frac{q(t)}{C} \tag{1}$$ である。 また、out端子につながる機材はハイインピーダンスであると考えて、電流はoutへほとんど流れないものとする。するとコンデンサに蓄えられる電気量の変化量が回路へ流れる電流と考えることができるので $$\frac{dq(t)}{dt}=I(t) \tag{2}$$ である。よって、 $$q(t)=\int I(t)dt \tag{3}$$ となる。 式3を式1へ代入すると $$V_c(t)=\frac{1}{C}\int I(t)dt \tag{4}$$ が導き出される。 さて、ラプラス変換により電流(I(t))は次の通り導き出される。 $$I(t)=\frac{E}{R}e^{-\frac{1}{CR}t} \tag{5}$$ 式5を式4に代入すると $$V_c(t)=\frac{1}{C}\int \frac{E}{R}e^{-\frac{1}{CR}t} dt \tag{6}$$ である。 Tip 指数関数の積分は次の通り（ただしDは積分定数）。 $$\int e^{kx}dx = \frac{1}{k}e^{kx}+D \tag{7}$$ 式7より式6は、 $$V_c(t)=\frac{E}{CR}(\frac{1}{-\frac{1}{CR}}e^{-\frac{1}{CR}t} + D) \tag{8}$$ $$V_c(t)=-Ee^{-\frac{1}{CR}t}+\frac{E}{CR}D \tag{9}$$ となる。 (\frac{E}{CR}D)は定数なのでまとめてDと置き換えてしまう。すると式9は、 $$V_c(t)=-Ee^{-\frac{1}{CR}t}+D \tag{10}$$ となる。 ここでDについて考えてみよう。(t=0)の時、コンデンサの電圧は0Vであるから $$V_c(0)＝-Ee^{-\frac{1}{CR}0}+D=0 \tag{11}$$ である。よって、 $$D=E \tag{12}$$ である。 式12を式10へ代入すれば $$V_c(t)=-Ee^{-\frac{1}{CR}t}+E \tag{13}$$ となり、 $$V_c(t)=E(1 - e^{-\frac{1}{CR}t}) \tag{14}$$ となる。これがRCローパスフィルタにステップ入力したときの過渡応答である。 プログラミングでシミュレーション せっかくなので式14をPythonでプログラミングしてプロットしてみよう。前回と同様に、Eを5V、Rを1kΩ、Cを1μFとしてシミュレーションしてみた。 ...

はじめに この記事では、ラプラス変換の理屈は分からなくてもラプラス変換をとにかく使ってみようというテーマでお届けする。 電子回路に慣れ親しんでいる方なら、オペアンプの中身や仕組みが分からなくてもオペアンプを使うことはできるだろう。ラプラス変換も一緒で、 ラプラス変換そのものの理屈は分からなくても、ラプラス変換を扱うことは高校生レベルの数学で可能なので安心を！ ラプラス変換といえば制御工学で学ぶと思うが、教科書ではラプラス変換の仕組みから説明が入り、実際ラプラス変換のご利益に触れないまま工学分野を嫌いになってしまう方も多いのではないだろうか？実は私もその一人だ。 しかし、ラプラス変換を実際に使ってみるとこれが非常に面白い。いままで諦めていた微分方程式などもすらすら解けるようになるのだ。だからラプラス変換を知らないのは非常にもったいない。もう一度言おう。 ラプラス変換することは高校生レベルの数学を学んでいれば誰でもできる。 そこでこの記事では、ラプラス変換の難しい説明をできるだけ省き、とにかく手を動かしてラプラス変換を使ってその利便性に触れることを目的として進めていく。 ラプラス変換とは ラプラス変換をとにかく使ってみよう！ とはいえ、使う上で大体のイメージは必要である。 ラプラス変換とは、時間領域tの関数を複素数領域sへと置き換えて、微分・積分などの演算を簡単にするための手法である。 さっそく、複素数領域だとかなんのこっちゃ分からないかもしれないが、 大事なのは、微積分の面倒な計算を足し算に置き換えてくれる超便利なツールであるということ。 複素数が分からなくてもラプラス変換は使えるので安心してもらって良い。 実は、電気回路や制御工学などでは解に微分・積分が多く含まれてしまうため、ラプラス変換することで足し算・割り算レベルの簡単な数式に置き換えられるメリットがあるのだ。 しかも、ラプラス変換はすでに用意されている変換表と照らし合わせて変換すれば良いだけなので何も難しいことはない。 作業手順としては次のような流れになる。 時間関数 $f(t)$ を組み立てる ラプラス変換する s領域で計算する 逆ラプラス変換して再び時間領域tへ戻す ラプラス変換された関数は複素数領域の関数であるため、そのままでは現実世界に適応できない。そこで、再び時間領域へ戻すために逆ラプラス変換が必要である。しかしこの逆ラプラス変換というのも、変換表通りやれば高校生レベルでできる話なので安心して良い。 それでは、ラプラス変換対応表を見ていくことにしよう。 ラプラス変換対応表 下の表が、ラプラス変換対応表である。時間領域tの関数を複素数領域sへ置き換えることがラプラス変換である。sがなんのこっちゃ分からなくても、とにかく時間tをsの関数に置き換えるにはこの表に当てはめて変換すれば良い。 この表を使えば、複素数領域sから時間領域tへ置き換えること（逆ラプラス変換）も簡単にできる。 $$f(t)$$ $$F(s)$$ 説明 $$δ(t)$$ $$1$$ 単位インパルス関数 $$δ(t)=\begin{cases}\infty~(t=0)\0~~~(t\neq0)\end{cases}$$ $$u(t)$$ $$\frac{1}{s}$$ 単位ステップ関数 $$u(t)=\begin{cases} 0~(t \lt 0) \ 1~(t \geqq 0) \end{cases}$$ $$e^{-at}$$ $$\frac{1}{s+a}$$ 指数関数 $$sin(ωt)$$ $$\frac{ω}{s^2+ω^2}$$ 三角関数（正弦） $$cos(ωt)$$ $$\frac{s}{s^2+ω^2}$$ 三角関数（余弦） $$f’(x)$$ $$sF(s)-f(0)$$ 微分 $$f’’(x)$$ $$s^2F(s)-sf(0)-f’(0)$$ 二階微分 $$\int_0^tf(τ)dτ$$ $$\frac{1}{s}F(s)$$ 積分 この記事では、上記の表の単位ステップ関数 指数関数 積分だけを使うので、今これらの全部の変換を覚える必要はない。 ...

はじめに デジタル信号におけるRCハイパスフィルタは次の式となる。 $$y_i=ay_{i-1}+a(x_i-x_{i-1}) \tag{1}$$ なぜこの式がRCハイパスフィルタになるのか、ここで詳しく解説してみたいと思う。 基本的にローパスフィルタの記事で説明した方法と同じやり方でハイパスフィルタは導き出すことができる。 アナログ回路のRCハイパスフィルタ 電圧を分かりやすくするために最初の回路図を少し変形させた。 回路図より、 $$E_R(t)=I(t)R \tag{1}$$ である。 また、 $$q(t)=CE_C=C(E(t)-E_R(t)) \tag{2}$$ $$I(t)=\frac{dq(t)}{dt} \tag{3}$$ なので、式1は $$E_R(t)=\frac{dq(t)}{dt}R=C(\frac{dE(t)}{dt}-\frac{dE_R(t)}{dt})R\tag{4}$$ となる。よって出力電圧(E_R(t))は次式で導き出される。 $$E_R(t)=RC(\frac{dE(t)}{dt}-\frac{dE_R(t)}{dt})\tag{5}$$ 電圧をデジタルデータで考える ここで、電圧をサンプリングデータに置き換える。(E(t))を(x_i)とし、(E_R(t))を(y_i)とすると、式5は次のようになる。ただし、(ΔT)はサンプリング周期である。 $$y_i=RC(\frac{x_i-x_{i-1}}{ΔT}-\frac{y_i-y_{i-1}}{ΔT}) \tag{6}$$ さて、式6を(y_i)について解くと、 $$y_i=\frac{RC}{RC+ΔT}y_{i-1}+\frac{RC}{RC+ΔT}(x_i-x_{i-1})\tag{7}$$ となる。 ここで(a=\frac{RC}{RC+ΔT})に置き換えると、 $$y_i=ay_{i-1}+a(x_i-x_{i-1}) \tag{8}$$ となり、RCハイパスフィルタのデジタルフィルタ計算式が導き出された。 カットオフ周波数 最後にカットオフ周波数と(a)の関係を導き出してみよう。 係数(a)は、 $$a=\frac{RC}{RC+ΔT}\tag{9}$$ であった。 また、RC回路におけるカットオフ周波数(f_c)は次式で求められる。 $$f_c=\frac{1}{2πRC}\tag{10}$$ これを展開すると、 $$RC=\frac{1}{2πf_c}\tag{11}$$ となる。式11を式9へ代入して係数(a)について展開すれば、 $$a=\frac{1}{2πΔTf_c+1} \tag{12}$$ の関係となる。

はじめに デジタル信号におけるRCローパスフィルタは次の式となる。 $$y_i=ax_i+(1-a)y_{i-1} \tag{1}$$ なぜこの式がRCローパスフィルタになるのか、ここで詳しく解説してみたいと思う。 アナログ回路のRCローパスフィルタ 話を分かりやすくするために、最初の回路図を少し変形させてみよう。 信号源を(E(t))とし、時間(t)によって信号の電圧は変化する。 抵抗(R)に流れる電流は(I(t))なので、抵抗の両端にかかる電圧は(RI(t))である。また、コンデンサの電圧は(Ec(t))とすると、(Ec(t)=\frac{q(t)}{C})が成り立つ。 よって、キルヒホッフの法則より、 $$E(t)=RI(t)+\frac{q(t)}{C} \tag{2}$$ となる。 また、ある時間における電流(I(t))はコンデンサの電気量(q)の変化量であるから $$I(t)=\frac{dq(t)}{dt} \tag{3}$$ で表せる。 ここでポイントなのは、図のOUT端子（抵抗とコンデンサの間）へ接続される負荷はハイインピーダンスであるため、出力端子へほとんど電流が流れていかないものとしている。 式3を式2へ代入すると、 $$E(t)=R\frac{dq(t)}{dt}+\frac{q(t)}{C}\tag{4}$$ となる。 コンデンサの電圧は(CE_c(t)=q(t))なので $$E_c(t)=\frac{q(t)}{C}\tag{5}$$ である。式5を式4へ代入すると、 $$RC\frac{dE_c(t)}{dt}+E_c(t)=E(t)\tag{6}$$ である。 電圧をデジタルデータで考える ここで、電圧をデジタルのサンプリングデータとして考えると(x_i=E(t))、(y_i=E_c(t))で置き換えることができる。 よって、式6は次のようになる。ただし(ΔT)はサンプリング周期である。 $$RC\frac{y_i-y_{i-1}}{ΔT}+y_i=x_i\tag{7}$$ これを(y_i)について解くと、 $$y_i=\frac{ΔT}{RC+ΔT}x_i+\frac{RC}{RC+ΔT}y_{i-1}\tag{8}$$ となる。 また、 $$\frac{RC}{RC+ΔT}=\frac{ΔT-ΔT+RC}{RC+ΔT}=\frac{RC+ΔT}{RC+ΔT}-\frac{ΔT}{RC+ΔT}$$ なので、ここで(a=\frac{ΔT}{RC+ΔT})とおけば、式8は $$y_i=ax_i+(1-a)y_{i-1}\tag{9}$$ となり、デジタル信号におけるRCローパスフィルタが導き出された。 カットオフ周波数 最後にカットオフ周波数と(a)の関係を導き出してみよう。 係数(a)は、 $$a=\frac{ΔT}{RC+ΔT}\tag{10}$$ であった。 また、RC回路におけるカットオフ周波数(f_c)は次式で求められる。 $$f_c=\frac{1}{2πRC}\tag{11}$$ これを展開すると、 $$RC=\frac{1}{2πf_c}\tag{12}$$ となる。式12を式10へ代入して計算すれば、 $$a=\frac{2πΔTf_c}{2πΔTf_c+1}\tag{13}$$ の関係が導き出される。