Juliaが1.0になって、徐々に環境が整ってきました。ここ数日、LibSndFile.jlやPortAudio.jlを使ってゴニョゴニョやっていたことをまとめておきます。
入出力の同期とか考えずに最低限の音再生ができればいいという場合は、
using FileIO: load, save import LibSndFile using PortAudio snd = load("hoge.wav"); stream = PortAudioStream(0, 2); write(stream, snd);
とすれば、ファイルを読み込んで、デフォルトの2チャンネル再生デバイス（内蔵スピーカとか）から再生してくれます。
続きを読む上記のJulia 0.6でやったことと同じですが、Julia 1.0がリリースされたので変化があったか見てみました。
前回と同じくadd1()を関数定義します。これはJulia 0.6でも1.0でも同じ。
julia> function add1(x)
return x+1
end
add1 (generic function with 1 method)
Julia 0.6では@code_llvmと@code_nativeの出力結果は以下のようになっていました。
julia> @code_llvm add1(10)
define i64 @julia_add1_62318(i64) #0 !dbg !5 {
top:
%1 = add i64 %0, 1
ret i64 %1
}
julia> @code_native add1(10)
.section __TEXT,__text,regular,pure_instructions
Filename: REPL[1]
pushl %ebp
decl %eax
movl %esp, %ebp
Source line: 2
decl %eax
leal 1(%edi), %eax
popl %ebp
retl
nop
nop
nop
nop
nop
nop
Julia 1.0では以下のとおり。
julia> @code_llvm add1(10)
; Function add1
; Location: REPL[1]:2
define i64 @julia_add1_35956(i64) {
top:
; Function +; {
; Location: int.jl:53
%1 = add i64 %0, 1
;}
ret i64 %1
}
julia> @code_native add1(10)
.section __TEXT,__text,regular,pure_instructions
; Function add1 {
; Location: REPL[1]:2
; Function +; {
; Location: REPL[1]:2
decl %eax
leal 1(%edi), %eax
;}
retl
nopw %cs:(%eax,%eax)
;}
コメントが挟まるようになり、何をしているかが分かりやすくなりました。LLVMコードは同じままですが、ネイティブコードの見た目がだいぶ変わりました。最近のアセンブリはわからないのですが、なんとなく雰囲気からEAX（レジスタ名）、DEC（1だけ引き算）、LEA（当該アドレスを返す）、RET（プロシージャから帰る）、NOP（何もしない）あたりだと考えられます。末尾にLがついているのは64 bit（Long）、Wは任意長（Words）でしょうか？ NOPを続けるのは一連の命令の合計長さを64 bitの倍数にしたいから？
ちゃんと理解しようと思うとIntelのマニュアルとか読むといいのかもしれませんね。
Intel® 64 and IA-32 Architectures Software Developer Manuals | Intel® Software
ヨーロッパLISPコンファレンスで発表された上記のプレゼンの中に、Juliaコードをネイティブコードに変換して中身を確認するマクロが紹介されていました。
公式ドキュメントではEssentials · The Julia Languageに書かれていますが、面白そうだったので、ちょっとだけ試してみました。
まずは簡単な関数を定義します。
julia> function add1(x)
return x+1
end
add1 (generic function with 1 method)
@code_llvmや@code_nativeマクロを使ってあげると、LLVMバイトコードや、ネイティブコードのアセンブリが見られます。
julia> @code_llvm add1(10)
define i64 @julia_add1_62362(i64) #0 !dbg !5 {
top:
%1 = add i64 %0, 1
ret i64 %1
}
julia> @code_native add1(10)
.section __TEXT,__text,regular,pure_instructions
Filename: REPL[5]
pushl %ebp
decl %eax
movl %esp, %ebp
Source line: 2
decl %eax
leal 1(%edi), %eax
popl %ebp
retl
nop
nop
nop
nop
nop
nop
@code_native add1()のようにしてもエラーが出てしまうのは
多重ディスパッチの影響があるので、add1(10)と具体的な値をつけてあげることで関数が適用する型を定めてあげるのだと思います。たとえばadd1(10.0)と浮動小数点数を入れたり、add1(BigInt(10))と多倍長整数にすると、生成されるコードも異なるものになります。
julia> @code_llvm add1(10.0)
define double @julia_add1_62439(double) #0 !dbg !5 {
top:
%1 = fadd double %0, 1.000000e+00
ret double %1
}
julia> @code_native add1(10.0)
.section __TEXT,__text,regular,pure_instructions
Filename: REPL[5]
pushl %ebp
decl %eax
movl %esp, %ebp
decl %eax
movl $200081256, %eax ## imm = 0xBECFF68
addl %eax, (%eax)
addb %al, (%eax)
Source line: 2
addsd (%eax), %xmm0
popl %ebp
retl
nop
nop
nop
nop
nop
nop
nop
nop
nop
nop
nop
nop
ちなみにCommon Lispにも言語仕様に似たものがありますね。（以下はSBCLの場合）
* (defun add1 (x)
(+ x 1))
ADD1
* (disassemble #'add1)
; disassembly for ADD1
; Size: 39 bytes. Origin: #x1001AC6434
; 34: 498B4C2458 MOV RCX, [R12+88] ; no-arg-parsing entry point
; thread.binding-stack-pointer
; 39: 48894DF8 MOV [RBP-8], RCX
; 3D: BF02000000 MOV EDI, 2
; 42: 488BD3 MOV RDX, RBX
; 45: 41BB8004B021 MOV R11D, #x21B00480 ; GENERIC-+
; 4B: 41FFD3 CALL R11
; 4E: 488B5DF0 MOV RBX, [RBP-16]
; 52: 488BE5 MOV RSP, RBP
; 55: F8 CLC
; 56: 5D POP RBP
; 57: C3 RET
; 58: 0F0B10 BREAK 16 ; Invalid argument count trap
NIL
Lauridsenのモノラル→ステレオのアップミックス・エフェクタをサクッと作ってみました。簡単ですねー。参考文献は以下。（注：最近使われている手法とは違ってかなり原始的です！）
desc: maruware Lauridsen Decorrelator (mono-to-stereo) // The ring buffer code was taken from // http://www.auriculaonline.com/wp/?p=232 slider1:20<1,100,1>Delay [ms] slider2:0<0,100,1>Mix [%] @init buf = 0; // buffer exists at offset 0 bufposR = 0; delay = srate * 0.010; bufposW = delay; buflength = srate * 0.100; mix = 0.0; @slider delay = srate * slider1 / 1000.0; // # of samples bufposR = bufposW - delay; bufposR < 0 ? bufposR = bufPosR + buflength; mix = slider2 / 100.0; @sample x = spl0; spl0 = x + buf[bufposR] * mix; spl1 = x - buf[bufposR] * mix; buf[bufposW] = x; bufposR = bufposR + 1 ; bufposR > buflength ? bufposR = 0; bufposW = bufposW + 1 ; bufposW > buflength ? bufposW = 0;
REAPERのプラグインを簡単に書けるJSFXという仕組みがあります。それをつかってステレオ音源のオシロスコープ（リサジュー曲線）を描くプラグインを書いてみました。やっていることは
だけです。しかもディレイラインは Delay, Basic Ring Buffer | JSFX Plug-Ins に書いてあったリングバッファをほぼほぼ持ってきただけです。ただ、今回はリングバッファの特色を使っていないですし、もっと効率の良い書き方に工夫できそうですね。
desc: maruware Oscilloscope (stereo) @init bufL = 0; bufR = srate; bufposW = 0; buflength = srate * 0.025; @sample s0 = spl0; s1 = spl1; bufL[bufposW] = s0; bufR[bufposW] = s1; bufposW = bufposW + 1 ; bufposW > buflength ? bufposW = 0; @gfx 480 480 gfx_r = 0.4; gfx_g = 0.4; gfx_b = 0.4; gfx_x = 0; gfx_y = 0; gfx_lineto(gfx_w/2, gfx_h/2, 1); gfx_lineto(gfx_w, 0, 0); gfx_x = gfx_w/2; gfx_y = 0; gfx_lineto(gfx_w/2, gfx_h, 1); gfx_x = 0; gfx_y = gfx_h/2; gfx_lineto(gfx_w, gfx_h/2, 1); gfx_r = 0.7; gfx_g = 1.0; gfx_b = 0.4; gfx_x = gfx_w / 2; gfx_y = gfx_h / 2; n = 0; while (n < buflength) ( bx = -0.5 * bufL[n] + 0.5 * bufR[n]; by = -0.5 * bufL[n] - 0.5 * bufR[n]; px = (bx + 1.0) * gfx_w / 2; py = (by + 1.0) * gfx_h / 2; gfx_lineto(px, py, 1); n += 1; );
Julia言語でPinkTSP信号を作るプログラムを書いた。とは言っても、以前書いたMatlabのものをそのまま移植しただけ。そもそも参考にしたのは以下の文献。
最低でも2秒間の信号（4回の同期加算）を44100Hz/16bitで欲しいときには
julia> pinktsp_generate(2.0, 44100, 16, 4);
とすると、pinktsp_2972msec_44100Hz_16bits_4times.wavというファイルができる。1回あたりの信号長は2秒よりも長い2の累乗の長さになるようにしてあるので、結果として217÷44100=2.972秒となる（2の累乗の信号長にしておくとFFTが効率よく使えるので）。掃引方向を上向きにしたり（上記守谷文献参照）、スピーカ再生時にクリックノイズが出ないようにしたり、実務上の理由によりオリジナルの藤本文献とは違うことをいくつかしている。
この信号を使って録音したものをインパルス応答に戻すプログラムはまだ書いてない。また後日、時間があるときに。
コードは以下。
続きを読む| 給油日 | オドメーター (km) | 給油量 (L) | 単価 (円/L) | 燃費 (km/L) | 距離単価 (円/km) |
|---|---|---|---|---|---|
| 2017-03-18 | 8559.6 | 2.76 | 125.72 | 45.94 | 2.74 |
| 2017-04-29 | 8703.9 | 2.77 | 127.80 | 52.09 | 2.45 |
| 2017-04-29 | 8905.9 | 3.21 | 128.97 | 62.93 | 2.05 |
| 2017-04-30 | 9032.5 | 2.20 | 131.82 | 57.55 | 2.29 |
| 2017-05-01 | 9190.3 | 2.70 | 127.04 | 58.44 | 2.17 |
| 2017-05-01 | 9325.7 | 2.52 | 128.97 | 53.73 | 2.40 |
| 2017-05-01 | 9469.7 | 2.17 | 129.03 | 66.36 | 1.94 |
| 2017-05-02 | 9625.3 | 3.71 | 133.96 | 41.94 | 3.19 |
| 2017-05-02 | 9787.8 | 2.81 | 125.98 | 57.83 | 2.18 |
| 2017-09-15 | 10004.7 | 3.09 | 123.96 | 53.03 | 2.34 |
| 2017-10-04 | 10152.2 | 2.08 | 125.96 | 70.91 | 1.78 |
| 2017-10-18 | 10314.9 | 3.42 | 125.73 | 47.57 | 2.64 |
| 2018-05-05 | 10499.5 | 3.18 | 137.74 | 58.05 | 2.37 |
実に半年ぶりの給油。乗ってなかったのだなぁ。去年の同時期から1年間で1,000 kmも走ってない。
