# 最低限必要なパッケージをimport using PyPlot using WAV # スペクトログラム可視化用の関数 function plotsp(X) figure(figsize=(16, 6), dpi=80, facecolor="w", edgecolor="k") imshow(X, origin="lower", aspect="auto") colorbar() end # 時間周波数解析のモジュール include("stft.jl") # 非負値行列因子分解の実装ですね include("nmf.jl") # 入力音声を読み込みます。ドレミの歌の冒頭です。 x, fs = wavread("test10k.wav") x = vec(x) println(fs) println(length(x)) # 時間周波数領域に変換します framelen = 512 hopsize = div(framelen, 2) X = stft(x, framelen, hopsize) # 振幅と位相に分けておく。位相は時間領域の信号を復元するのに使います。 Y, phase = abs(X), angle(X) @show size(Y) plotsp(20log(Y)) # ド、レ、ミ、ファの4つの音で構成されていることが今回はわかっている(普通はありえんのですが)ので、4つの基底に分解します K = 4 srand(98765) @time H, U = nmf_euc(Y, K) nothing # H: 学習された基底スペクトルを可視化します figure(figsize=(16, 10), dpi=80, facecolor="w", edgecolor="k") for k=1:K subplot(K, 1, k) plot(H[:,k], linewidth="2") end # 上手く行けば、ド、レ、ミ、ファに対応する基底スペクトルが得られます。 # U: 基底スペクトルに対応するアクティビティを可視化します figure(figsize=(16, 10), dpi=80, facecolor="w", edgecolor="k") for k=1:K subplot(K, 1, k) plot(vec(U[k,:]), linewidth="2") end # 基底ごとに時間領域の信号を構成し、可視化します figure(figsize=(16, 10), dpi=80, facecolor="w", edgecolor="k") for k=1:K # 時間領域の信号を復元 y = istft(H[:,k]*U[k,:] .* exp(im * phase), framelen, hopsize) # かしか subplot(K, 1, k) xlim(0, length(x)-1) ylim(-1.0, 1.0) plot(y, linewidth="4") # 結果をファイルに書き出したい場合は、以下のようにすればOK # wavwrite(float16(y), "適当なファイル名.wav", Fs=fs) end println("おわり")