Julia 1.0でオーディオ入出力(とインパルス応答測定)

Juliaが1.0になって、徐々に環境が整ってきました。ここ数日、LibSndFile.jlPortAudio.jlを使ってゴニョゴニョやっていたことをまとめておきます。

TLDR。入出力の同期とか考えずに最低限の音再生ができればいいという場合は、

using FileIO: load, save
import LibSndFile
using PortAudio
snd = load("hoge.wav");
stream = PortAudioStream(0, 2);
write(stream, snd);

とすれば、ファイルを読み込んで、デフォルトの2チャンネル再生デバイス(内蔵スピーカとか)から再生してくれます。

続きを読む

Juliaコードのディスアセンブル(Julia 1.0)

marui.hatenablog.com

上記の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

Juliaコードのディスアセンブル(Julia 0.6)

www.youtube.com

ヨーロッパ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

REAPERのJSFXでモノラル→ステレオのアップミックスやってみた

Lauridsenのモノラル→ステレオのアップミックス・エフェクタをサクッと作ってみました。簡単ですねー。参考文献は以下。(注:最近使われている手法とは違ってかなり原始的です!)

www.aes.org

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でステレオ・オシロスコープ作ってみた

REAPERのプラグインを簡単に書けるJSFXという仕組みがあります。それをつかってステレオ音源のオシロスコープ(リサジュー曲線)を描くプラグインを書いてみました。やっていることは

  1. ディレイライン2ch分のオーディオ信号を溜める
  2. それを使ってグラフ描画

だけです。しかもディレイラインDelay, Basic Ring Buffer | JSFX Plug-Ins に書いてあったリングバッファをほぼほぼ持ってきただけです。ただ、今回はリングバッファの特色を使っていないですし、もっと効率の良い書き方に工夫できそうですね。

f:id:amarui:20180626124640g:plain

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;
);