STFTとISTFTの秘密の関係

タイトルに深い意味はありません。
まあ、ここの記事で意味なんかあった試しがありませんけどね。

今回は、FFTW3ライブラリを用いたSTFTとISTFTについて、『変換後の逆変換で、元データを再現できるのか』という内容をお送りします。

ちなみに、今回実装は省きます。
需要もないですし、ソース貼り付けると重くなるので止めときます。
別に自分が書いたソースを見られるのが嫌なわけではありません。

そういうことなので、実装方面についてはあまり触れません。
面倒くさいので。
おや？ 何か見えましたか？
しかし、敢えて言いましょう。
世の中、気にした方が負けなのです。

……違うか。

STFTについて:

まあ、特にSTFTに関しては情報過多なくらいなので触れることもないでしょう。
ただ、僕の場合は元のデータの長さに合わせて、STFTする回数を動的に変える手法を使っていたので、その話でもしようと思います。

『元のデータ長に対して何回STFTが可能か』は、『いくつ紙テープを貼り付ければ期待する長さになるか』という問題に似ています。
STFTにおけるサンプリング数を紙テープの一枚あたりの長さに、サンプリングの間隔を紙テープの糊しろの長さに見立てれば、何となく分かってもらえると思います。
たとえば、サンプリング数を8、サンプリング間隔を4とすれば、それは8cmの長さの紙テープを糊しろ4cmでつなぎ合わせていくことと同じです。

さて、問題は元データの長さに合わせてSTFTの回数をどう変動させるかです。
元データ長に対してサンプリング数8、サンプリング間隔4のSTFTは何回実行できるか。
8cmの紙テープを4cmの糊しろで張り付けていき、cm以下で最も長くする場合に紙テープは何本必要か。
ここで、ならば、どうでしょう？

をそれぞれ回数と紙テープの数に見立てましょう。
紙テープは、よりです。
これは、最初の一枚を除いた結果なので、この場合、15枚でぴったり64cmになります。

STFTの回数も同じです。
すなわち、あるデータ長のデータに対して、サンプリング数、サンプリング間隔のSTFTを適用する場合の最大STFT回数は、

で求められます。
は、紙テープと同様に最初の一枚分です。

これで、元データ長に対して、STFTの回数を動的に求めることができます。

ISTFTについて:

ぶっちゃけ、今回の本題はここです。

まず、前提を思い出してください。
忘れた？ 困った人ですね。
これは、FFTW3ライブラリを使っている場合を想定しているということですよ。
まったく、僕が前置きで下らないことしか言わないと思ってましたね？
概ね、その認識で合ってます。
だからこそ、ここで再確認したわけで。

いいとして。
ここからは、長さ64のデータに対して、サンプリング数8、サンプリング間隔4でSTFTしたものをISTFTすると仮定して話を進めます。

FFTW3を使っている以上、注意しなければならないことがいくつかありますね。
まず、逆変換で得られる値をSTFTしたときのサンプル数で割らなければなりません。
これは、公式マニュアルの通りです。

さて、ISTFTの逆変換のイメージは大丈夫でしょうか？
まず、STFTなので得られた周波数特性も複数あると思います。
それをそのまま戻してつなぎ合わせると、どうなるか。

ここで、前節で話した紙テープの話を思い出してください。
ISTFTした結果をそのままつなぎ合わせることは、糊しろ0で貼り付けることと同じです。
要するに、間違いです。
紙テープの糊しろ、つまりサンプリング間隔は4なので、逆変換した結果も当然4ずつ重ねて繋げていかなければなりません。

ここで、別の問題が発生します。
重ねていくということは、重複するということです。
だからと言って、逆変換で得られた結果を勝手に削ってはいけません。
この問題は、該当する部分を重複した数で割ることで解決します。

ちょっとイメージし辛いですね。
まず、紙テープを思い浮かべてください。
長さ8cmの紙テープを8分割しましょう。
これを、Cの配列に見立ててください。
このとき、最初と最後の例外を除いて、配列の添字0から3は前の結果と、4から7の4つの値は後の結果とそれぞれ重複します。
なので、基本的には2で割ればよいのです。
繰りかえしますが、逆変換の結果は全て意味を持っているので、確固たる意味がある場合を除いて、勝手に削ってはいけません。

さて、ISTFTの面倒なところは、この『重複する数』がフレキシブルな点です。
つまり、サンプリング数とサンプリング間隔の組み合わせ次第では、単純に2で割っても期待する結果を得られないということです。
ここで、『重複する数』の法則は以下となります。
※いずれも最初と最後は例外的な扱いになります。

• サンプリング数 == サンプリング間隔*2
  • 重複は必ず2つ。
  • 中間は必ず重複する。
• サンプリング数 < サンプリング間隔*2
  • 重複は必ず2つ。
  • 中間に重複しない箇所がある。
• サンプリング数 > サンプリング間隔*2
  • 重複は複数。
  • 中間は必ず重なる。
  • 最初の方と最後の方で、重複数が段階的に変化する。

面倒くさいですよねー。
まあ、数式で求めることもできるんです。
できますが、とにかく面倒です。

この問題に対する僕のとった解決策は、カウンター配列を作ることでした。
要するに、元データ長と同じ長さの配列を作り、結果を格納する毎にカウントしていくんです。

ちょっと、想像の世界に飛び立ちましょうか。
言語はCでいきましょう。
そこ、嫌いとか言わない。

まず、ISTFTの結果格納配列resultを用意します。
次に、重複数をカウントする配列countも用意します。
これらの配列は、ともに元データの長さ分、ここでは64だけ確保します。
また、STFTの結果が15個のfftw_complex*ということはオーケーですね？
では、いきましょう。

まず、最初の逆変換で、8個の結果が返されます。これはSTFTのサンプリング数と同じです。
なので、resultの0から7番目にそれぞれ適切に結果を格納します。
合わせて、countの0から7番目もインクリメントします。

次の逆変換も同様に8個の結果を得られます。
ここで、サンプリング間隔が4であることを思い出してください。
したがって、resultの4から11番目に結果を格納することになります。
同様に、countも4から11番目をインクリメントします。

これを15回繰り返し、最後に

for(int i=0;i<64;i++)
  result[i] /= count[i];

みたいにすれば終わりです。

どうでしょう、イメージできました？
もちろん、最初に書いた通り、全ての結果はサンプリング数、つまりここでは8で割る必要があることも忘れてはいけません。
この方法はメモリを大量に消費しますが、一番楽な解決法だと思います。

最後に窓関数の扱いですね。
STFTの際に窓関数を用いていたのであれば、ISTFTでも同様の窓関数を使う必要があります。
STFTの際は乗算だったので、ISTFTでは結果に対して除算を適用すればよいのです。

僕の場合はWAVの音声ファイルに対して、これらの方法を使いました。
結果から言って、単純にSTFTしてISTFTしたものはオリジナルと同じでした。
diffをとろうが、cmpしようが、同じファイルという結果を得られました。
まあ、WAVに戻す処理を施す前はかなり細かいところで差異がありましたけどね。
e-30以下の領域で。

そんなわけで、『変換後の逆変換で、元データを再現できるのか』は、『ほぼできる』と結論づけていいと思います。

編集後記:

いやあ、久々に長文書いて疲れました。

まあ、今回も走り書きのメモみたいなものなので、もしかしら足りない部分があるかもしれません。
え？
『足りてる部分があるの？』……だと……？
……そうですね。
個人的には、これ以上の不足は必要ないですね。

こんな感じで。
今回は、終わりたいと思います。見つかっちゃったね！きらっ☆