#01 Distortion編―MATLAB & Simulinkでオーディオエフェクト

#初めに
[前回](https://qiita.com/9W3R7Y/private/3a45ac71f8baf199f45d)の記事でMATLAB Simulinkでの音声処理の基礎を述べた．今回からは，実際に様々なエフェクトをMATLAB Simulinkで実装していく．今回はDistortionを実装する．

#Distortionとは
Distortionとは音を歪ませるエフェクトである．エレキギターの派手な音は，この歪によるものである．

Distortionでは，入力を歪ませた音が出力されるわけであるが，具体的には入力を非線形の関数に通したものが出力されている（線形の関数を用いるなら，それは単なるGainである）．この記事ではまず，Hard clipperとSoft clipperを実装する．これらのDistortionでは奇数倍音のみが生じるので，両者を実装した後に偶数倍音も生じさせるDistortionを作る．倍音については後述する．

#準備
以下のようなモデルを組む．入力のAudio Oscillatorでは40Hzのsin波を生成している．Audio Oscillatorの設定画面はダブルクリックで開ける．Audio Oscillatorが生成した音をdB Gain（これは前回作成したものではなく，標準で用意されたものである），MATLAB Function ブロックの順につなぎ，スピーカーから出力する．dB GainとMATLAB Function ブロックの出力はMatrix Concatenateで一つにまとめてTime ScopeとSpectrum Analyzerで可視化する．KnobはdB GainのGainと接続し，0 dBから24 dBの範囲で動かせるようにする．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/517d5bd7-175b-ca8a-939a-dbae681f79d4.png)

#Hard clipper
まずはハードクリッパーを実装する．ハードクリッパーは以下の式で表せる．

```math
f(x) = \left\{
\begin{array}{ll}
\,\,\,\,1 & (1 \leqq x) \\
\,\,\,\,x & (-1 \le x \le 1)\\
-1 & (x \leqq -1)
\end{array}
\right.
```
グラフの外形は以下のようになる．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/dbe07cb2-37b8-ac8f-8034-046df51d1743.png)

これを実装するには，単純にif文で場合分けすればよい．注意すべきことは，Simulinkではある程度の長さを持った信号が扱われるということである．例えば，Audio Oscillatorのデフォルトの設定では，512サンプルがまとめて出力されている．そのため，for文を回して全要素に対して計算を行う．MATLAB Funcitonを以下のように書き換える．

```matlab
function y = fcn(u)
    for i = 1:length(u)
        if (u(i) >= 1)
            u(i) = 1;
        elseif (-1 < u(i))
            u(i) = u(i);
        else
            u(i) = -1;
        end
    end
    y = u;
```
後は実行して，Gainをいじる．音と波形，スペクトルの変化に着目する．
<blockquote class="twitter-tweet"><p lang="ja" dir="ltr">こっちの方が色々分かりやすいな <a href="https://t.co/lMkQE3dWtq">pic.twitter.com/lMkQE3dWtq</a></p>&mdash; 9W3R7Y (@qwerty_16180339) <a href="https://twitter.com/qwerty_16180339/status/1339001591550988288?ref_src=twsrc%5Etfw">December 16, 2020</a></blockquote> <script async src="https://platform.twitter.com/widgets.js" charset="utf-8"></script>
入力のGainを上げると音が歪むのが分かる．また，波形も1以上と-1以下では真っ平にクリッピングされているのが分かる．スペクトルは，はじめは440Hz付近だけが盛り上がっていたが，Gainを上げるにしたがって440 Hz以外のピークが現れるのが分かる．

#Soft Clipper
つぎにSoft Clipperを実装する．Soft Clipperにはいくつか種類があるが，今回はCubic Nonlinear Distortionを実装する．このDistortionは以下の式で表せる．

```math
f(x) = \left\{
\begin{array}{ll}
\,\,\,\,\displaystyle \frac{2}{3} & (1 \leqq x) \\
\,\,\,\,x-\displaystyle \frac{x^3}{2} & (-1 \le x \le 1)\\
-\displaystyle \frac{2}{3} & (x \leqq -1)
\end{array}
\right.
```


グラフの外形は以下のようになる
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/54cd6b2b-0d9e-016b-42df-00f0eb6db1c5.png)
Hard clipper に比べると Softな見た目をしている．

これも先ほど同様に実装し，実行する．

```matlab
function y = fcn(u)
    for i = 1:length(u)
        if (u(i) >= 1)
            u(i) = 2/3;
        elseif (-1 < u(i))
            u(i) = u(i)-u(i)^3/3;
        else
            u(i) = -2/3;
        end
    end
    y = u;
```
<blockquote class="twitter-tweet"><p lang="en" dir="ltr">Cubic Nonlinear Distortion（Soft Clipper） <a href="https://t.co/sVtnV3ySgi">pic.twitter.com/sVtnV3ySgi</a></p>&mdash; 9W3R7Y (@qwerty_16180339) <a href="https://twitter.com/qwerty_16180339/status/1339006456486457344?ref_src=twsrc%5Etfw">December 16, 2020</a></blockquote> <script async src="https://platform.twitter.com/widgets.js" charset="utf-8"></script>

#倍音とは
倍音は，メインの周波数成分（基音）の整数倍の周波数の音のことである．sin波以外のすべての音は倍音成分を含んでいる．

#Distortionと倍音

さて，ここでHard clipperの倍音成分を確認してみる．Spectrum Analizerの設定で横軸を線形にし，拡大表示する．すると倍音成分は以下のグラフのようになっていることが分かる．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/447e593e-c9f2-047a-4466-393a09b43b9a.png)
一番左の成分は440 kHz程度であり，基音である．その次のピークは1.3 kHz程度となっている．ここで，$440 \times 3 = 1320$なので，この成分は3倍音の成分であることが分かる．同様に，その次の成分は5倍音，さらにその次の成分は7倍音の成分である．このように，Hard clipperでは奇数倍音しか生じないことが分かる．

Soft clipperの倍音成分は以下に示す．これも奇数倍音の成分しか生じていない．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/a6e2c33f-8e36-c9e6-d892-69906cbb638e.png)

調べると，Distortionの関数が奇関数である場合には，奇数倍音しか生じないようである．尚，奇関数とは$f(x)=-f(-x)$となるような関数のことである．偶数倍音を生じさせるには，奇関数の対称性を無くせばよいらしい．そのための簡単な方法としては，Distortionを掛ける前に適当な0でない定数を足すことが挙げられる（参考文献ではほかの手法も紹介されている）．

#偶数倍音が生じるDistiotionを作る
以下のようにブロックを組む．今までのモデルのOscillatorの後ろにAdd ブロックを置き，定数を加算する．その定数はKnobブロックで操作可能にする．GainのknobにはPre Gain,加算する定数のknobにはConstと名前を付けてある．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/db187db4-ed47-e533-16f8-e9b34260e263.png)
これを実行し，音と倍音構成の変化を確かめる．倍音構成は以下のようになる．
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/397506/08ca2a96-b882-f64e-653d-7c1d3eb7f325.png)
たしかに，今までは生じなかった偶数倍音が生じている．音も結構変わった．
<blockquote class="twitter-tweet"><p lang="ja" dir="ltr">Hard clipperの前に定数を加算して偶数倍音を生じさせる <a href="https://t.co/ZUyawvVhLA">pic.twitter.com/ZUyawvVhLA</a></p>&mdash; 9W3R7Y (@qwerty_16180339) <a href="https://twitter.com/qwerty_16180339/status/1339224686492717057?ref_src=twsrc%5Etfw">December 16, 2020</a></blockquote> <script async src="https://platform.twitter.com/widgets.js" charset="utf-8"></script>

#終わりに
今回はDistortionを作成した．[次回](https://qiita.com/9W3R7Y/private/e94c30d8bbddce5a3267)はDelayを作成する．

#GitHub
以下で完成したオーディオエフェクトを公開している．
https://github.com/qwerty16180339/matlab_audio_effects

#参考文献
* DSP Related.com NONLINEAR DISTORTION　https://www.dsprelated.com/freebooks/pasp/Nonlinear_Distortion.html