tsugaru’s blog

主に技術的なことをつぶやきます。プログラミング,ポケモンGO,音楽,アニメ等

楽器用チューナーアプリ作成の理論(第1回)

単にピッチを認識するだけだと思っていたが、いろいろな困難にぶち当たった。音程の認識は音高認識というようであるが、この検索ワードで調べても、音「声」認識のための分析であって、さらにリアルタイムではないので、似ている部分はあれど、精度もリアルタイム分析も難しいものがあった。様々な資料を読んで色々と調べてみたが、Philip McLeod のFast, Accurate Pitch Detection Tools for Music Analysisが最も要点をまとめ、網羅的であったので、これを読んで欲しい。自分なりに以下に概要をまとめる。

チューナーアプリとして求められるクオリティー

以下のように、測定精度、計算速度共に求められ、非常に困難が多い。

  • 人が区別でいる音程は区別できないといけない
    • 人が知覚できる違い(JND:just noticeable difference)は周波数や音の大きさ、長さ、突然の音の変更に依存するが、1000~4000HzではJNDはおよそ8セント程度
    • ちなみにYAMAHAのTDM-700は+-1セント
  • あらゆる楽器に対応しないといけない
    • コントラバスの最低音E1(41.2Hz)から、ピッコロのC8(4186Hz)まで
    • 声の場合はローパスフィルターを使って、音の間引きをして(ダウンサンプリング)効率化できたが、分析範囲の関係上使えない。(違うフィルターを用いる)
    • ナイキスト周波数を考えると、サンプリング周波数44100Hzは目的の周波数を測定可能である。
      • サンプリング周波数が大きいということはデータが多くなるので、効率の良いアルゴリズムが必須である。
    • ちなみにYAMAHAのTDM-700はC1(32.7Hz)からC8(4186Hz)
  • 音の入力から出力までの時間を最小限にしなければならない
    • その時刻での周波数を決定するための入力サンプル(窓の長さという)が少ないほどピッチ分析は不正確になるが、正確さを求めると長くなる
  • 音の入力間隔に対応できるだけの計算処理が求められる

使用する技術

当ブログ主の注

周波数波長に変換してピッチを分析する手法がある。 音声認識の世界ではたまに使われているようであるが、F0法である。周期的な周波数ピークのうち根音となるF0を取り出すというものである。これは音声認識の世界では有効であるが、二つ問題がある。FFTは一回で計算するデータ量に精度が依存する。例えばサンプリング周波数が44100Hzで窓の長さが1024とすると、FFTでの最も長い周期は1024/44100秒、つまり、周波数にして、43Hzである。これでは多少の音の工程がわかれど、とてもチューナーとして使えるものではない。さらにこれが使える分解能であったとしても、楽器の音程は周波数ドメインでのピークをみればわかるほど単純でもない。その理由は、別で記事を書く予定である。

本題

周波数ドメインにおいて、一般的な波を正確に分析するには3~5周期分のデータが必要だが、ここでは楽器の音に関連した波形を取り出したいので、その場合、自己相関法(autocorrelation)やSDF(Square Difference Function)を用いることで、2周期分程度の波形で分析ができる。細かいことは触れないが、どちらも時間ドメインでの手法で、音が一定であれば、どんな音色であっても波形に繰り返しがあるという想定で分析する。実用の際は、短い時間での分析であるので、音が一定であるという過程がある。

では自己相関とSDF、どちらを使うべきだろうか。長い時間一定の音においてはほぼ違いは出ないが、短時間においては大きな違いが出るようだ。以下でグラフをみてもらうために式を抑えておく。

窓の長さをWとする。自己相関方においては、波が周期\tauであるような度合い(大きいほどその周期であると言える)r'(\tau)は、r'(\tau)=\sum _ { j = 0 }^ { W-\tau -1 } x _ j *x _ {j+ \tau },ただし(0 \leq \tau \lt W)と表される。 SDFにおいては、波が周期\tauであるような度合い(小さいほどその周期であると言える)d'(\tau)は、d'(\tau)=\sum _ { j = 0 }^ { W- \tau -1 }( x _ j -x _ {j+ \tau })^ 2,ただし(0 \leq \tau \lt W)と表される。 いずれも0の時最大、最小(0)をとる。

ここで、同じ周期2.25で初期位相を変えた正弦波について、自己相関方とSDFで分析した結果が以下のようになる。

f:id:pika-bika:20200920160020p:plainf:id:pika-bika:20200920160042p:plainf:id:pika-bika:20200920160054p:plain
自己相関方とSDFの比較

上のように、自己相関方では初期位相によって、誤差が出ていることがわかる。そのため、以下ではSDFを用いることとし、どのようにして計算速度を上げるかに焦点を当てる。

SDF

式をみてわかるように、SDFにおいて時間計算量はO(W^ 2)となる。近年のプロセッサーが3GHz程度であることから、窓の大きさは1024か2048程度まで大きくしても実現可能であると言える。ここでd'は自己相関法で用いたr'を使って、d'(\tau)=\sum _ { j = 0 }^ { W- \tau -1 }( x _ j^ 2 +x _ {j+ \tau }^ 2)-2\sum _ { j = 0 }^ { W-\tau -1 } x _ j *x _ {j+ \tau }=m'(\tau)-2r'(\tau)とあらわあせるとすぐわかる。

ここで、自己相関法で求めるr'は定義通りであれば計算量はO(W^ 2)であるが、パワースペクトル密度関数の逆フーリエ変換と一致することで知られている(Wiener-Khintchineの定理)。逆フーリエ変換に高速離散フーリエ逆変換(IFTT:inverse fast Fourier transform)を使うことで計算量をO(W \rm {log} W)に抑えられる。

さらに、m'について、m'(\tau)=m'(\tau -1) - x^ 2 _ {\tau -1} - x^ 2 _ {W-\tau}とあらわせるので、m'(0)の計算をすれば、帰納的にすぐにもとまることがわかる。

結果として、全体で時間計算量はO(W \rm {log} W)となる。FFTによって、丸め誤差が発生することもあるが、32bitの不動小数点を使えば十分である。

SDFの改良

SDFによって求めるためには、0に近いもの周期を取り出せばよかった。しかし、それはどのくらい近ければ良いのだろうか。最小値を求めると、もしかしたら、取り出すべき周期の整数倍の波長を取り出すかもしれない。そのように、基準が難しい。それに対し、自己相関法では初めて0を通ってから次の0までの最大値を求めるようにすれば、基準が明確である(原文:a zero centre-reference)。そこで、SDFの自己相関法に勝る正確性を生かしつつ、取り出すべき周期の基準を明確にできるような工夫をしたい。

唐突ではあるが、ここでn'(\tau)=\dfrac{2\sum _ { j = 0 }^ { W-\tau -1 } x _ j x _ {j+ \tau }}{\sum _ { j = 0 }^ { W- \tau -1 }( x _ j^ 2 +x _ {j+ \tau }^ 2)}=\dfrac{2r'(\tau)}{m'(\tau)}なるn'を定義する。これは計算するとわかるが自己相関方の値について-1から1になるように正規化したもので、このような正規化の仕方を特にSNAC(specially-normalised autocorrelation)と呼ぶ。この式はこれまでに示した式を用いてn'(\tau)=\dfrac{m'(\tau)- d'(\tau)}{m'(\tau)}=1-\dfrac{d'(\tau)}{m'(\tau)}と表せる。したがって、n'はSDFを利用しつつ、自己相関法での周期特性表記のメリットを生かした関数であると言える。さらに正規化しているため、ピークについての相対的な基準も決めやすい(あまりにゼロで囲まれた最初のピークが小さい時、飛ばしたりできる)。この経緯から、SNACのことを、NSDF(Normalized SDF)と呼ぶこともある。この論文では発音がしやすいでしょ、ってことでSNACを採用しているのでそれに従う。

SNACと共に窓関数を利用する

計算量

上での分析の主役d'の代わりに\hat{d'}(\tau)=\sum _ { j = 0 }^ { W- \tau -1 }w _ j w _ {j+\tau}( x _ j -x _ {j+ \tau })^ 2を分析の主役とする。wは窓関数である。この時、\hat{d'}(\tau)=\sum _ { j = 0 }^ { W- \tau -1 }w _ j w _ {j+\tau}( x _ j^ 2 +x _ {j+ \tau }^ 2)-2\sum _ { j = 0 }^ { W-\tau -1 }w _ j w _ {j+\tau} x _ j x _ {j+ \tau }=\hat{m'}(\tau)-2\hat { r'}(\tau)、\hat{m'}(\tau)=\sum _ { j = 0 }^ { W-\tau -1 } (x _ j ^ 2 w _ j * w _ {j+\tau}) +\sum _ { j = 0 }^ { W-\tau -1 } (x _ {j+\tau} ^ 2 w _ {j+\tau} * w _ {j}) 、\hat{n'}(\tau)=\dfrac{2\hat{r'}(\tau)}{\hat{m'}(\tau)}=\dfrac{2\sum _ { j = 0 }^ { W-\tau -1 } w _ j w _ {j+\tau} x _ j x _ {j+ \tau }}{\sum _ { j = 0 }^ { W- \tau -1 }( w _ j w _ {j+\tau} x _ j^ 2 +x _ {j+ \tau }^ 2)}というように対応させる。この時、先ほどのように、\hat{d'}の代わりに\hat{n'}を用いて分析する。ここで、 \hat{n'}の分子\hat{m'}は上の式で分解したものをみてもらうと、jとj+\tauで分離して表記しているのがわかる。それぞれを関数として、相互相関関数とみてFFTを使うことにより計算時間をO(W \rm {log} W)に抑えられる。分母についてはこれまでのようにFFTを用いる。このようにして求めることをWSNAC(windowed SNAC)と名付けている。

窓関数

ハン窓を使うとある。 窓関数については筆者が勉強不足なので略す。

放物線補完

\tauはこれまで整数であるので、精度を上げるために放物線補完をする。ピークの点と両サイドの点を通る放物線の頂点を計算することで求め、精度をより上げる。

精度の比較

これらは原文のChapter 5にグラフ付きで詳しく述べてある。 SNAC,Unbiased Autocorrelation,Autocorrelation,Windowed Unviassed Autocorrelation,WSNACで比較をしている。 まず、サイン波を周波数を変えて測定する。次に、複数周波数のサイン波を混ぜて測定する。いずれにおいてもSNACが最も精度が良くWSNACが次によかった。SNACはずれても2セントほどであったが、自己相関法は10セント程度なら平気でずれてしまうほど、誤差に違いが出ていた。次に、グリッサンドやビブラートの影響、音量の影響を調べるため、定常波でなく、時間変動する波の実験をする。すると、SNACよりWSNACの方が良い精度を得られていることが多くなった。

ピッチと基本周波数

これについては後日記事を書く予定。 原文においては6章に相当。

続きを別で書く予定。