WAVファイルの音を解析する (1)
WAVファイルのデータをFFT解析するのに PowerShell でやる人が他にいるとも思えませんが、調べたことのまとめも兼ねて書いておきます。
コードはこちらに。
予めお断りしておきますが、FFTとか数値解析、音の波形などについてはまったくの素人です。FFTするとその時間帯の音が周波数でわかる、くらいのゆるい理解でやってます。ですのであちこちに試行錯誤が入ってますのでご了承ください。
ではコードを見ていきましょう。
PitchChecker.ps1
起動スクリプトになっている本体です。本体クラスが Run() で実行されると初期化処理の後、ExportTo が指定されていればデータの出力処理、それ以外の通常ケースはデータの解析処理をします。
Run() {
$this.Init()
if (!$script:ExportTo) {
$this.DecodeWave()
} else {
$this.ExportToCsv($script:ExportTo)
}
}
Init()
初期化関数 Init() では指定された WAV ファイルの読み込みと FFT 関連の初期化をします。
Init() {
if (!$script:File) {
throw "-File が指定されていません"
}
$this.LoadWaveFile($script:File)
$this.BufferSize = $this.BuffReader.Bytes.Length
$this.InitFFT()
}
LoadWaveFile()
WAVファイルの読み込み自体は別クラス WaveFile でやり、WAVファイルのパラメータを使いやすいようにクラスプロパティとして準備しておきます。
WAVファイルはいろいろな形式の音データが格納できるようになっていますがここでは 32bit Float と 16bit固定(PCM) だけを扱います。ネットで見つけたサンプルデータが 24bit だったのですが、24bit はプロ用機材とかで使われるとか。手元のスマホのPCM録音アプリは 16bit PCM、ハンディレコーダは 32bit Float で出力するので 24bit は後回しで。
ともあれ,ヘッダ情報に応じて 32bit float と 16bit データを取り扱うクラスのどちらかをインスタンス生成して BuffeReader として持っておきます。
SingnalScaleFactor は後で説明しますが、FFTの周波数解像度を上げるために0パディングすると信号強度が弱まってしまうので増幅しないといけないのですが、あれ? 設定だけして使っていないだと…。まぁいいか。将来使います。
LoadWaveFile([string]$filename) {
$wave = [WaveFile]::New()
$wave.Load($filename)
$this.SamplingRate = $wave.FmtChunk.SamplesPerSec
$this.BlockSize = $wave.FmtChunk.BlockAlign
$this.DataStartIndexInFile = $wave.DataChunk.DataStartIndex
$this.BytesPerSec = $this.SamplingRate * $this.BlockSize
switch ($wave.FmtChunk.BitsPerSample) {
32 {
if ($wave.FmtChunk.FormatTag -eq 3) {
$this.BuffReader = [Float32]::New($wave.FileBuffer, $wave.FmtChunk.BlockAlign)
$this.SignalScaleFactor = 10000
$this.PeakIgnoreThrethold = 0.1
}
break
}
16 {
if ($wave.FmtChunk.FormatTag -eq 1) {
$this.BuffReader = [PCM16]::New($wave.FileBuffer, $wave.FmtChunk.BlockAlign)
$this.SignalScaleFactor = 100
$this.PeakIgnoreThrethold = 0.0001 #??
}
}
}
if (!$this.BuffReader) {
throw "Wave fileformat not supported BitsPerSample=$($wave.FmtChunk.BitsPerSample), FormatTag=$($wave.FmtChunk.FormatTag)"
}
}
InitFFT()
FFT関連の初期化を行います。FFTは元データのサンプリングレートと1回で処理するデータの個数によって解析結果の解像度が決まるということなので、例えば 0.2 秒毎に解析すると 1 / 0.2 となって 5Hz 単位の結果にしかなりません。これでは話にならないのでそういう場合には0データを追加して見た目のデータ個数を増やすのだそうです。ただしこれをやると信号が「薄まる」ので予めデータを「増強」しておく必要がある、と。そのために SingaleScaleFactor を用意したのですが、使い忘れというか高速化のための変更をしているうちに漏れてしまったようです。まぁそれっぽく動くのでこの件は保留。
話を戻すと、もろもろの計算を簡素化するために、FFT処理するデータ数をサンプリングレートに合わせます。実データの数 WaveDataCount と実データを埋め込む場所 FFTDataStartIndex を計算しておきます。結果の周波数解像度 FreqResolution は 1 Hz になるはずですが、変えたくなった時のために変数で持っておきます。
また、FFT結果は 1Hz~データ数/2 までの周波数分が出力されますが、全部を調べるのは無駄なので最小MinFreq 、最大 MaxFreq を決めておきます。
最後にはまる原因になった Complex32[] の配列 FFTDataBlank を作っておきます。さすがにプロパティ名決め打ちで書き込みに行く関数はアレなので、指定のプロパティに書き込むようにしました。(見苦しくないとは言ってない)
InitFFT() {
$this.SampleCount = $this.SamplingRate
$this.WaveDataCount = [int]($script:Interval * $this.SamplingRate)
$this.FFTDataStartIndex = [int](($this.SampleCount - $this.WaveDataCount) / 2)
$this.FreqResolution = $this.SamplingRate / $this.SampleCount
$this.MinFreq = 170
$this.MaxFreq = 2000
SetComplexArray $this 'FFTDataBlank' $this.SampleCount
}
これで初期化関連は終わりです。
DecodeWave()
音データを解析する部分を見ていきましょう。
基本的な構造は、
- WavデータからComples32[]へデータを用意 (SetData())
- FFT実行 (外部関数 FFT_Forward)
- FFT結果から周波数を検出 (DetectPeak())
- 周波数を音に変換 (Pitch.Find())
- 音の良し悪しの評価 (ResultManager.Add())とその表示
をデータがあるだけ繰り返すとしています。
DecodeWave() {
$pitch = [Pitch]::New()
$resMan = [ResultManager]::New($script:PitchThrethold)
$lastIndex = $this.BufferSize - $script:Interval * $this.BytesPerSec
$preNote = ""
while (($pos = $this.BuffReader.GetCurrentPosition()) -le $lastIndex) {
$this.SetData()
FFT_Forward $this.FFTData
$this.DetectPeak()
if ($script:ShowPeak) {
log "DEBUG>> $pos $($this.Peaks |%{'[{0}hz,{1}]' -f $_.Freq,$_.Mag})"
}
$pk = ($this.Peaks |Sort Freq)[0]
if ($pk.Freq -le $this.MinFreq) { continue }
if ($pk.Mag -lt $this.PeakIgnoreThrethold) { continue }
$pitch.Find($pk.Freq)
$msg = $resMan.Add($pitch)
if ($script:ShowAll -or $preNote -ne $pitch.Name) {
log "$($this.PocToTimeSpan($pos)) $($pitch.Name) $($pk.Freq) ($('{0:+#;-#;0}' -f $pitch.FreqDiff)) $msg"
}
$preNote = $pitch.Name
}
$resMan.ShowAll()
}
SetData()
FFTするデータを準備する関数です。前述のように、解像度を上げるために0をダミーデータとして付け加えるのですが、FFTで使った配列領域は全域が更新されて返ってくるので毎回初期化が必要です。最初はバカ正直にループで0埋めしていましたが、メモリ管理のオーバヘッドを含んでも Array.Clone() の方が速いようなので初期化の際に作った配列(FFTDataBlank)を Clone() してデータ処理用の配列(FFTData)を毎回作成します。
実際のデータコピーはバイト列から数値をよろしく取ってきてくれる別クラス BUfferReader.CopyData() に任せます。この BufferReader は初期化時にWAVファイルのデータ形式に応じたクラスのインスタンスを準備してあるので、このループ内で都度フォーマットの違いを判断することは不要。毎回 If-Then で分岐するより速いはず。調べてないけど。
SetData() {
$this.FFTData = $this.FFTDataBlank.Clone()
$this.BuffReader.CopyData($this.FFTData, $this.FFTDataStartIndex, $this.WaveDataCount)
}
DetectPeak
さてFFT結果から周波数を検出するのですが、FFTの結果を見ると、
のように倍音? 3倍音? のところでも波を検出しています。っていうかむしろ倍音の 880Hz 近辺の方が音が大きい? どうなの? これが普通?
ともあれ 440Hz 近辺にもピークがあるのでこれを検出することに注力します。このデータで聞こえるのは下のラということと、上のラ(880Hz)を出そうとしていたとしてもこれだけ下のラ(440)が聞こえたらアカンでしょうということもあります。
そのあたりのグラフを拡大すると
となり、440だけでなく前後にもいろいろあることが判ります。
グラフのY軸は Magnitudeで複素数の絶対値、ここではその周波数で出ている音の大きさを表している、はず。実際の数値はこんな感じで、
Magnitude | 周波数 |
---|---|
612,744 | 437 |
1,058,301 | 438 |
1,404,385 | 439 |
1,573,249 | 440 |
1,529,941 | 441 |
1,289,699 | 442 |
912,674 | 443 |
このあたりをひとまとめとして判定したうえで一番 Magnitude の高いものが欲しいのです。
ということで、隣同士の Magnitude を比べて、ピークの左で上昇中($rising) から下りに転じたところ(ひとつ前 >= ココ) を頂点として検出します。ピークはいくつかあるので、ここでは大きさ順に上位 3 つを選んでいます。
DetectPeak() {
$nPeaks = 3
$this.Peaks = @(, [Peak]@{Freq = 9999; Mag = 0})
$cutline = 0.3
$preMag = -1
$rising = $false
$sIndex = $this.MinFreq / $this.FreqResolution - 1
$eIndex = $this.MaxFreq / $this.FreqResolution + 1
foreach ($i in $sIndex .. $eIndex) {
$freq = $i * $this.FreqResolution
$mag = $this.FFTData[$i].Magnitude
if ($mag -le $cutline) { continue }
if ($mag -le $preMag) {
if ($rising) {
$this.Peaks += ,[Peak]@{Freq = $freq - $this.FreqResolution; Mag = $preMag }
}
$rising = $false
}
else {
$rising = $true
}
$preMag = $mag
}
$this.Peaks = $this.Peaks | Sort Mag | Select -last $nPeaks
}
Pitch.Find(), ResultManager.Add()
さて、検出した周波数ピークから実際の音を求めていくのですが、音の大きさでみると倍音の方が大きかったりする例がある、というか最初に取った音がこうなっていて、WAVファイルを再生しても低いラで聞こえるので複数検出したピークの周波数が一番低いものを出ている音とすることにします。この例だと 440 Hzで、これを「ラ」に変換するのが別クラスの Pitch.Find()。
さらにそれの良し悪しを判断するのが、別クラス ResultManager.Add()。今回は「ラ」は 442 Hz に設定しているので 440 Hzは「2Hz 低いけどまぁ良いでしょう」とするかたち。
この周波数の誤差の良し悪し判定は適当です(開き直り)。音楽のプロの方は 1Hz ずれていると気持ち悪いそうなのですが、1Hzの誤差も許さん、とすると全滅しかねないので今のところは 1% の誤差はOKとする甘い評価です。
長くなってきたので別クラスは次回に。
(追記) 続きはこちら。
Discussion