Wizard Notes

Pythonを使った音楽信号分析の技術録、作曲活動に関する雑記

LibROSA:音高分析のための定Q変換(CQT)と、逆定Q変換の性能評価

はじめに

f:id:Kurene:20210710233854p:plain

音信号を分析する時間周波数分析手法としてはSTFT(短時間フーリエ変換)が良く使われますが,特に音楽を対象にして音高、コード、メロディなどを分析する場合は定Q変換(CQT)という手法が便利です。

www.wizard-notes.com

www.wizard-notes.com

定Q変換を自分で実装するのは結構大変ですが,Pythonでは音楽信号分析モジュール LibROSA を使うことで、さくっと定Q変換を試すことができます

LibROSAでの定Q変換

librosa.cqt を呼び出すことで、数行で定Q変換を実装できます

import numpy as np
import librosa
import librosa.display
import soundfile as sf
import matplotlib.pyplot as plt


# 音源ファイル読込
sr = 44100
filepath = librosa.ex('trumpet')
y, sr = librosa.load(filepath, sr=sr)

# 定Q変換,振幅スペクトログラム
hop_length      = 128
n_octave        = 8
bins_per_octave = 12
C = librosa.cqt(y, sr=sr,
                hop_length=hop_length,
                n_bins=n_octave*bins_per_octave,
                bins_per_octave=bins_per_octave)  
M  = np.abs(C) 
librosa.display.specshow(librosa.amplitude_to_db(M, ref=np.max),
                         bins_per_octave=bins_per_octave,
                         sr=sr, x_axis='time', y_axis='cqt_note', cmap="jet")
plt.title('Constant-Q power spectrum')
plt.colorbar(format="%+2.0f dB")
plt.show()

librosa.cqt の重要な引数を列挙します.

  • hop_length
    • 信号に定Q変換をする際にどのくらいずらしながら分析するかを指定します
    • なるべく小さいサンプル数の方が時間変化を捉えやすくなりますが,計算量や出力されるCQTスペクトログラムのサイズが大きくなってしまいます
  • bins_per_octave
    • 1オクターブを何分割するかを指定します
    • 基本は12分割(12音平均律)でよいと思います
    • 目的によっては12分割で解像度が足りないことがあるため、その場合は24分割にするなど、分割数を増やします
  • n_bins
    • 全ての周波数ビンの数です。
    • 分割数12の場合,n_bins/12が何オクターブ分析するかに相当します
  • fmax
    • CQTの一番低い音高の周波数です。
    • librosaの場合、fmin=librosa.note_to_hz('C2')というようにlibrosa.note_to_hz()と音高名の文字列を使って指定するのがおすすめです
  • tuning

逆定Q変換

逆STFTのように、一応 逆定Q変換 もあります。

しかし,CQTは、

  • 周波数ビンによって分析する時間区間が異なる
  • 周波数ビンの中心周波数が高くなるにつれて,分析窓長を小さくして周波数分解能を下げる

という性質があるため、逆の変換時に品質劣化が起こることに注意が必要です *1

LibROSA では librosa.icqt として実装されています。

そこで,おそらく周波数解像度を上げれば逆変換の性能も上がると考え,bins_per_octaveを増やした時の性能変化を確認してみました.

sr = 44100
filepath = "miku.wav"
y, sr = librosa.load(filepath, sr=sr)

hop_length      = 128
n_octave        = 8
window          = "hann"
for k in range(8):
    bins_per_octave = 6 * (k+1)
    C = librosa.cqt(y, sr=sr,
                    hop_length=hop_length,
                    window=window,
                    n_bins=n_octave*bins_per_octave,
                    bins_per_octave=bins_per_octave) 
    y_hat = librosa.icqt(C=C, sr=sr, 
                         window=window,
                         hop_length=hop_length,
                         bins_per_octave=bins_per_octave)

    sf.write(f"cqt_icqt{bins_per_octave}.wav", y_hat, sr, format="WAV", subtype='PCM_16')
    error = 10 * (np.log10(np.sum(y**2)) - np.log10(np.sum((y[0:y_hat.shape[0]]-y_hat)**2)))
    print(f"{bins_per_octave},{error}")
    
sf.write("cqt_orig.wav", y,     sr, format="WAV", subtype='PCM_16')

SNR(品質)は、

f:id:Kurene:20210711002033p:plain

として計算しています。

結果は以下のようになりました。

f:id:Kurene:20210711002001p:plain

聴感的にも bins_per_octave=12は厳しく、 bins_per_octave=24 ~ 36は必要かなと思います。それでもオリジナルと比べると高い周波数は厳しいので実用では工夫が必要そうです。

以下, bins_per_octaveを変えた時の逆変換の音質比較デモ動画です。

youtu.be

*1:そもそもCQT変換行列の形状が(n_bins, n_fft) n_bins<<n_fftであり逆行列がない