Juliaで楽器を作ろう（コナミ効果）！

# (#JuliaLang) [Julia Advent Calendar 2014](http://qiita.com/advent-calendar/2014/julialang)
# 12月25日担当 kimrin(twitter:kimrin)

楽しかった[Julia Advent Calendar](http://qiita.com/advent-calendar/2014/julialang)も今回が最後になりました(^^)
見ていてくださったオーディエンスの皆様、ありがとうございました。
[Julia Advent Calendar](http://qiita.com/advent-calendar/2014/julialang)の執筆者一同に代わって、お礼申し上げます。

## Juliaで楽器を作る
### Juliaでソフトシンセを改造しましょう！！

まずソフトシンセとは何かについて説明します。
ソフトシンセとは、シンセサイザー（あのビョーンとかヒュイーンとか鳴るあれです）
をソフトウェアで実現したものです。僕のMacだと、Logic Pro X
という所謂音楽編集ソフト（DAWといいます）から、外部の楽器として
Audio Unitという規格のソフトシンセを使うことができます。

で、Audio Unitを作ってしまおう、という猛者のためにAppleはサンプルの
SDKを用意しているのですが、これを変更してJuliaをembeddedしてしまうでね、
というわけです。

タイトルほど大げさではないのですが（ぉぃ）Juliaで正弦波を出す
Audio Unit モジュール（要するにソフトシンセ）を作ってみました。
今回は波形を作っている箇所をJuliaで書いて、あとはC++のコード内に
Juliaをembeddedして、Juliaのコードを呼べる仕組みを作ります。

楽しみですね、楽しみですね！

### Audio Unitとは

詳しいことを説明し始めると本が一冊書けるので（ぉぃ）
要点を説明します。

Audio Unitがやる仕事は、

* UIにボリュームなどのパラメーターを見せる
* DAW側から要求されたサンプル数の音データをレンダリング（実際のデータにして転送）する

です。で、UIのパラメーターを変更したら、レンダリングのデータに変更を反映する必要がありますです。

具体的にはC++のRenderというメソッドが細かい周期で呼ばれます。で、引数にサンプル数が指定されています。このメソッドRenderはその指定されたサンプル数の音データをleftとrightという二つのfloatの配列に書き込み、転送します。

SinSynthというサンプルのAudio Unit音源（正弦波音源）があり、Appleから提供されています。
これをハックして、JuliaにRenderの機能を持って行って、スクリプトの変更だけで音生成のロジックを変えてしまおう、というわけです。

## コナミ効果とは
往年のコナミのゲームのBGMでは、「コナミ効果」と呼ばれる
特殊な効果を使って基本的な波形の楽音に厚みを付けていたと
言われています。

これは具体的には微妙に周波数の異なる波形を二つ生成して合成し、
「うなり」を伴う音を生成するというものです。

物理的現象として、二つの正弦波などを微妙に違う周波数で一緒に鳴らすと、
周期 $f$で強弱が付くことがわかっています。これを「うなり」
などと言います。

二つの周波数をそれぞれ $f_1$、 $f_2$とすると、うなりの周波数 $f$は、
$|f_1 - f_2|$となることが知られています。つまり今、$f_1 = 440[Hz]$
（これは中央のラの音で、A（アー）などといいます）に対して $f_2 = 441[Hz]$
の音を鳴らすと、440.5[Hz]の音が周期 $f = 1 [Hz]$ でうねります。つまり１秒ごとに
弱くなったり、強くなったりします。

今回はこの往年の効果をAudio Unitで実現してみましょう！！

## とっかかりから

とりあえず`JULIALANG_HACK`という`#define`を切って、コードに変更を入れていきます。

```cpp:SinSynth.h
#define JULIALANG_HACK (1)
#ifdef JULIALANG_HACK
#include <julia.h>
#endif

```
で、このとき`julia.h`をインクルードします。そしてライブラリ参照に`libjulia.dylib`を追加し、`/usr/local/include/julia`にインクルードのサーチパスを通します。

これでまず、SinSynthに`libjulia.dylib`が結びつけられ、`libjulia`内の関数を呼べるようになりましたー！

---

次に音のレンダリングを実現する構造体である、TestNoteを変更します。
ついでにいくつかあるパラメータの説明をします。

コナミ効果のためにphase2という変数をここで追加します。

基本的に正弦波は$\sin(phase)$という感じで表現して、毎サンプルごとにこのphaseに対して鳴らしている音の周波数から計算した小さな値、freqを足し算していきます。

またampとmaxampがそれぞれ振幅（音の大きさ）と最大振幅です。
まず鍵盤が打鍵されると急激に音が大きくなり、これがmaxampに達すると
こんどは少しずつ小さくなるか、音を維持するかのどちらかになります。

で、鍵盤から手を離すと、今度はmaxampから急激に音が小さくなります。
このときそれぞれ打鍵時のプラス成分、手を離したときのマイナス成分を
それぞれ定数として、
up_slope, dn_slopeという変数に入れておきます。
これは毎サンプルごとに打鍵時はup_slope分だけampが大きくなる、
という感じです。

```cpp:SinSynth.h
struct TestNote : public SynthNote
{
	virtual					~TestNote() {}

	virtual bool			Attack(const MusicDeviceNoteParams &inParams)
								{ 
#if DEBUG_PRINT
									printf("TestNote::Attack %p %d\n", this, GetState());
#endif
									double sampleRate = SampleRate();
									phase = 0.;
#ifdef JULIALANG_HACK
                                    phase2 = 0.;
#endif
                                    amp = 0.;
									maxamp = 0.4 * pow(inParams.mVelocity/127., 3.); 
									up_slope = maxamp / (0.1 * sampleRate);
									dn_slope = -maxamp / (0.9 * sampleRate);
									fast_dn_slope = -maxamp / (0.005 * sampleRate);
									return true;
								}
	virtual void			Kill(UInt32 inFrame); // voice is being stolen.
	virtual void			Release(UInt32 inFrame);
	virtual void			FastRelease(UInt32 inFrame);
	virtual Float32			Amplitude() { return amp; } // used for finding quietest note for voice stealing.
	virtual OSStatus		Render(UInt64 inAbsoluteSampleFrame, UInt32 inNumFrames, AudioBufferList** inBufferList, UInt32 inOutBusCount);

    double phase, amp, maxamp;
#ifdef JULIALANG_HACK
    double phase2;
#endif
	double up_slope, dn_slope, fast_dn_slope;
};
```

これでヘッダーファイルへの変更は終わりです。

## SinSynth.cppへの変更点

さて、では実装を行っているSinSynth.cppに、Julia言語組み込みを行って行きましょう。

このうち、UIから見えるパラメータに0から2のKONAMI Effectというパラメータを追加したのですが、これについては省略します。具体的にはKonamiHzという変数に0から2のうなりの振動数が入っていると思ってください。KonamiHzが1のときは、周波数を1Hzずらして鳴らそうというわけです。

```cpp:SinSynth.cpp
OSStatus SinSynth::Initialize()
{	
#if DEBUG_PRINT
	printf("->SinSynth::Initialize\n");
#endif
	AUMonotimbralInstrumentBase::Initialize();
	
	SetNotes(kNumNotes, kMaxActiveNotes, mTestNotes, sizeof(TestNote));
#ifdef JULIALANG_HACK
    char jlpath[] = "/usr/local/lib/";
    jl_init(jlpath);
    JL_SET_STACK_BASE;
    char str[] = "include(\"/Users/kimrin/Synth.jl\")";
    jl_eval_string(str);
#endif
#if DEBUG_PRINT
	printf("<-SinSynth::Initialize\n");
#endif
	
	return noErr;
}
```

`SinSynth::Initialize()`で、Juliaコンテキストの初期化を行います。
まず`jl_init()`でコンテキストを初期化します。
このときイメージファイルを見つけられるようにJULIA_HOMEのパスを
指定してあげます。brewでJuliaをインストールした場合は、
```bash:bash
$ ls /usr/local/lib/../lib/julia/sys.ji
```
のように、`../lib/julia/sys.jl` でイメージファイルが見つけられるパス、
すなわち`/usr/local/lib`を指定します。

`JL_SET_STACK_BASE`はおまじないだと思ってください。

そしておもむろに僕のホームにある、`Synth.jl`をインクルードします。
この`Synth.jl`の中に必要な関数を定義しておきます。詳しくは後述します。

これでJuliaのファイルが無事evalされました。
あとはC++の世界から必要に応じてJuliaの関数を呼べばいいだけです。

---
## Synth.jl

ここで先に先ほど`include`した`Synth.jl`を見ていきましょー

```jl:Synth.jl
twopi = 2.0 * 3.14159265358979

function getKonamiFreq(frequency,sampleRate,konamiConstant::Float32)
    float64((frequency+konamiConstant)*(twopi/sampleRate))
end

function attacked(left,right,amp,maxamp,phase,phase2,slope,globalVol,num,freq,freq2,ch)
    for frame = 1:num
        if amp < maxamp
            amp = amp + slope
        end
        #out = (sin(phase)^5) * amp * globalVol
        out = (((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        phase = phase + freq
        phase2 = phase2 + freq2
        if phase > twopi
            phase = phase - twopi
        end
        if phase2 > twopi
            phase2 = phase2 - twopi
        end
        left[frame] = left[frame] + out
        if ch == 2
            right[frame] = right[frame] + out
        end
    end
    Float64[amp,phase,phase2]
end

function released(left,right,amp,maxamp,phase,phase2,slope,globalVol,num,freq,freq2,ch)
    endFrame = 0xFFFFFFFF
    for frame = 1:num
        if amp > 0.0
            amp = amp + slope
        elseif endFrame == 0xFFFFFFFF
            endFrame = frame
        end
        #out = (sin(phase)^5) * amp * globalVol
        out = (((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        phase = phase + freq
        phase2 = phase2 + freq2
        left[frame] = left[frame] + out
        if ch == 2
            right[frame] = right[frame] + out
        end
    end
    Float64[amp,phase,float64(endFrame),phase2]
end

```
まず一番目の関数、`getKonamiFreq`について説明します。
先ほど`phase`は毎サンプルごと、`freq`だけ加算されるという
話をしました。で、うなりを出すためにもう一つ正弦波を作って、
それを元の正弦波と加算しようというわけです。で、その２番目の
`freq`にあたる`freq2`(C++の世界では）を計算するのが、この関数です。

```jl:Synth.jl
twopi = 2.0 * 3.14159265358979

function getKonamiFreq(frequency,sampleRate,konamiConstant::Float32)
    float64((frequency+konamiConstant)*(twopi/sampleRate))
end

```
`frequency`に打鍵したノートの周波数が入ってきます。で、これに$[0,2]$の$konamiConstant$を加算してあげます。で、$\frac{2\pi}{sampleRate}$で
一回のサンプルの前進で`sin()`関数の引数に加算する量を計算します。

次に`attacked`関数の説明をします。ちょっと難しい感じになっていますが、
要点はここです。

```jl:Synth.jl
        #out = (sin(phase)^5) * amp * globalVol
        out = (((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        phase = phase + freq
        phase2 = phase2 + freq2
```

上のコメントアウトしてある`out`が従来の正弦波です。
これに対して二つの変数`phase`と`phase2`による正弦波を平均しているところが、
コナミ効果です。5乗しているのは、音の大きさを整えるためと思ってください。
それに対して現在の振幅`amp*globalVol`を掛けています。
で、それが終わったら、`phase`と`phase2`にそれぞれ`freq`と`freq2`を
足してあげます。このときの`freq2`が、先ほどの関数`getKonamiFreq`で
求めた値になります（ここではC++ワールドから渡してもらっています）。

そして`left`と`right`の配列にデータを`num`個分積み、帰ります。
`attacked`は打鍵されて音が大きくなるところと、ampmaxに達して
音が持続するところを担当しています。

`released`も同様ですが、`slope`に負の値が入ってくるところが違います。
あとサンプルが0個の場合の処理（つまりイベントの終端）があるのが異なります。

## SinSynth.cppの変更箇所

まず始めに、サンプルを積んでいる、C++のコードを掲載します。

```cpp:SinSynth.cpp
OSStatus		TestNote::Render(UInt64 inAbsoluteSampleFrame, UInt32 inNumFrames, AudioBufferList** inBufferList, UInt32 inOutBusCount)
{
	float *left, *right;

	//（中略）
	switch (GetState())
	{
		case kNoteState_Attacked :
		case kNoteState_Sostenutoed :
		case kNoteState_ReleasedButSostenutoed :
		case kNoteState_ReleasedButSustained :
			{
                
				for (UInt32 frame=0; frame<inNumFrames; ++frame)
				{
					if (amp < maxamp) amp += up_slope;
					float out = pow5(sin(phase)) * amp * globalVol;
					phase += freq;
					if (phase > twopi) phase -= twopi;
					left[frame] += out;
					if (right) right[frame] += out;
				}
            }
			break;

```
先ほどの`Synth.jl`の`attacked`を、この`for`ブロックの代わりに呼びたいと思います。

```cpp:SinSynth.cpp
#ifdef JULIALANG_HACK
    double freq2;
    jl_function_t *func0 = jl_get_function(jl_main_module, "getKonamiFreq");
    jl_value_t *argument = jl_box_float64(Frequency());
    jl_value_t *argument2 = jl_box_float64(sampleRate);
    jl_value_t *argument3 = jl_box_float32(konamiHz);
    jl_value_t *ret = jl_call3(func0,argument,argument2,argument3);
    freq2 = jl_unbox_float64(ret);
    
    jl_value_t* array_type = jl_apply_array_type(jl_float32_type, 1);
    jl_array_t *left_jl;
    jl_array_t *right_jl;
    
    left_jl = jl_ptr_to_array_1d(array_type, left, inNumFrames, 0);
    if (right) {
        right_jl = jl_ptr_to_array_1d(array_type, right, inNumFrames, 0);
    }
    jl_function_t *func1 = jl_get_function(jl_main_module, "attacked");
    jl_function_t *func2 = jl_get_function(jl_main_module, "released");
    
#endif
    
	
	switch (GetState())
	{
		case kNoteState_Attacked :
		case kNoteState_Sostenutoed :
		case kNoteState_ReleasedButSostenutoed :
		case kNoteState_ReleasedButSustained :
			{
                
#ifdef JULIALANG_HACK
                jl_value_t **args;
                JL_GC_PUSHARGS(args, 12); // args can now hold 2 `jl_value_t*` objects
                args[0] = (jl_value_t *)left_jl;
                args[1] = (jl_value_t *)right_jl;
                args[2] = jl_box_float64(amp);
                args[3] = jl_box_float64(maxamp);
                args[4] = jl_box_float64(phase);
                args[5] = jl_box_float64(phase2);
                args[6] = jl_box_float64(up_slope);
                args[7] = jl_box_float32(globalVol);
                args[8] = jl_box_int32(inNumFrames);
                args[9] = jl_box_float64(freq);
                args[10] = jl_box_float64(freq2);
                UInt32 ch = (right)?2:1;
                args[11] = jl_box_int32(ch);
                // Do something with args (e.g. call jl_... functions)
                jl_array_t *y_out = (jl_array_t*)jl_call(func1, args, 12);
                double *xData = (double*)jl_array_data(y_out);
                amp = xData[0];
                phase = xData[1];
                phase2 = xData[2];
                
                JL_GC_POP();
                    
#else
				for (UInt32 frame=0; frame<inNumFrames; ++frame)
				{
					if (amp < maxamp) amp += up_slope;
					float out = pow5(sin(phase)) * amp * globalVol;
					phase += freq;
					if (phase > twopi) phase -= twopi;
					left[frame] += out;
					if (right) right[frame] += out;
				}
#endif
            }
			break;
```

長々していますが、これで`Julia`の関数を呼べます。具体的には`jl_call`がそれを担当しています。`Julia`に渡すときはboxして渡し、`Julia`から`C++`に戻すときはunboxします。

`released`についても掲載します。

```cpp:SinSynth.cpp
		case kNoteState_Released :
			{
                UInt32 endFrame = 0xFFFFFFFF;
                
#ifdef JULIALANG_HACK
                jl_value_t **args;
                JL_GC_PUSHARGS(args, 12); // args can now hold 2 `jl_value_t*` objects
                args[0] = (jl_value_t *)left_jl;
                args[1] = (jl_value_t *)right_jl;
                args[2] = jl_box_float64(amp);
                args[3] = jl_box_float64(maxamp);
                args[4] = jl_box_float64(phase);
                args[5] = jl_box_float64(phase2);
                args[6] = jl_box_float64(dn_slope);
                args[7] = jl_box_float32(globalVol);
                args[8] = jl_box_int32(inNumFrames);
                args[9] = jl_box_float64(freq);
                args[10] = jl_box_float64(freq2);
                UInt32 ch = (right)?2:1;
                args[11] = jl_box_int32(ch);
                // Do something with args (e.g. call jl_... functions)
                jl_array_t *y_out = (jl_array_t*)jl_call(func2, args, 12);
                double *xData = (double*)jl_array_data(y_out);
                amp = xData[0];
                phase = xData[1];
                endFrame = static_cast<UInt32>(xData[2]);
                phase2 = xData[3];
                
                JL_GC_POP();
                
#else
				for (UInt32 frame=0; frame<inNumFrames; ++frame)
				{
					if (amp > 0.0) amp += dn_slope;
					else if (endFrame == 0xFFFFFFFF) endFrame = frame;
					float out = pow5(sin(phase)) * amp * globalVol;
					phase += freq;
					left[frame] += out;
					if (right) right[frame] += out;
				}
#endif
                if (endFrame != 0xFFFFFFFF) {
#if DEBUG_PRINT
					printf("TestNote::NoteEnded  %p %d %g %g\n", this, GetState(), phase, amp);
#endif
					NoteEnded(endFrame);
                }
			}
			break;
			
```

実際にはもう一つslopeの違う`fastRelease`というのもあります。

---
## インストールと実行♡

Xcodeでコンパイルした結果を所定の場所に置くのがAudio Unitの流儀です。

```bash:bash
kimuraken-no-MacBook-Prox:~ kimrin$ sudo rm -rf /Library/Audio/Plug-Ins/Components/SinSynth.component
Password:
kimuraken-no-MacBook-Prox:~ kimrin$ sudo cp -R Library/Developer/Xcode/DerivedData/SinSynth-cjkaayhphtgvfjdlzbbxdzumybjk/Build/Products/Development/SinSynth.component /Library/Audio/Plug-Ins/Components/SinSynth.component
kimuraken-no-MacBook-Prox:~ kimrin$ sudo chmod -R 775 /Library/Audio/Plug-Ins/Components/SinSynth.component
kimuraken-no-MacBook-Prox:~ kimrin$ 
```

`/Library/Audio/Plug-Ins/Components/`に、`SinSynth.component`一式を置きます。

同様に`/Users/kimrin/Synth.jl`がちゃんとあることを確認します。本来は配布物に含めるべきですが。。。

さて、Logic Pro X を開きます！！


![Logic1.png](https://qiita-image-store.s3.amazonaws.com/0/61379/c7e577d1-e9d4-8e2c-63ce-b05bb3c25651.png)

Audio UnitをInst1(ソフト音源）に割り当ててみましょう。

![Logic2.png](https://qiita-image-store.s3.amazonaws.com/0/61379/bd447aa8-6922-d20c-997d-9fc05c566396.png)

おおーー！！　KONAMI EffectがちゃんとUIに表示されてる！！


ではバウンスしてmp3にした結果ですー！！

---

## 結果


# 結構大きな音ですので、注意してください！！


### f=0（同じ周波数）
[One Sine Wave](https://soundcloud.com/kimrinkimrin/sinsynthjulia-languageone-sin-wave)

### f=1（違う周波数）

[Two Sine Waves](https://soundcloud.com/kimrinkimrin/sinsynthjulia-languagetwo-sin-wave)

### f=0.02 （かなり近い周波数の二つの正弦波）

[Two Sine Waves(f=0.02Hz)](https://soundcloud.com/kimrinkimrin/julialongfnear-two-sine-waves-that-differ-002hz)

### ノコギリ波（f=0.01Hz）

[Two SawTooth Waves(f=0.01Hz)](https://soundcloud.com/kimrinkimrin/sawtooth001)

（サンクラがエンベデッドできなかったんじゃーw）

## おわりに

最終的なノコギリ波の`Synth.jl`を引用掲載します。

```jl:Synth.jl(最終盤)
twopi = 2.0 * 3.14159265358979

function getKonamiFreq(frequency,sampleRate,konamiConstant::Float32)
    float64((frequency+konamiConstant)*(twopi/sampleRate))
end

function attacked(left,right,amp,maxamp,phase,phase2,slope,globalVol,num,freq,freq2,ch)
    for frame = 1:num
        if amp < maxamp
            amp = amp + slope
        end
        #out = (sin(phase)^5) * amp * globalVol
        saw0 = 2.0 * (phase/twopi - floor(phase/twopi + 0.5))
        saw1 = 2.0 * (phase2/twopi - floor(phase2/twopi + 0.5))
        saw2 = (((saw0+saw1)*0.5)^5) * amp * globalVol
        houkei = (sign((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        out = (((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        phase = phase + freq
        phase2 = phase2 + freq2
        if phase > twopi
            phase = phase - twopi
        end
        if phase2 > twopi
            phase2 = phase2 - twopi
        end
        left[frame] = left[frame] + saw2 #out
        if ch == 2
            right[frame] = right[frame] + saw2 #out
        end
    end
    Float64[amp,phase,phase2]
end

function released(left,right,amp,maxamp,phase,phase2,slope,globalVol,num,freq,freq2,ch)
    endFrame = 0xFFFFFFFF
    for frame = 1:num
        if amp > 0.0
            amp = amp + slope
        elseif endFrame == 0xFFFFFFFF
            endFrame = frame
        end
        #out = (sin(phase)^5) * amp * globalVol
        saw0 = 2.0 * (phase/twopi - floor(phase/twopi + 0.5))
        saw1 = 2.0 * (phase2/twopi - floor(phase2/twopi + 0.5))
        saw2 = (((saw0+saw1)*0.5)^5) * amp * globalVol
        houkei = (sign((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        out = (((sin(phase)+sin(phase2))*0.5)^5) * amp * globalVol
        phase = phase + freq
        phase2 = phase2 + freq2
        if phase > twopi
            phase = phase - twopi
        end
        if phase2 > twopi
            phase2 = phase2 - twopi
        end
        left[frame] = left[frame] + saw2 #out
        if ch == 2
            right[frame] = right[frame] + saw2 #out
        end
    end
    Float64[amp,phase,float64(endFrame),phase2]
end
```

いかがだったでしょうか。ノコギリ波でなんとかコナミ効果っぽい感じでした。
市販のソフトシンセには遠く及びませんね(^^;)

そういうわけで、Juliaのembeddedの参考になればと思いますー。
ありがとうございましたー

また来年お会いしましょうー　ではー(^o^)/

（文責 kimrin）