こんにちは。すごい改善のAitoです。
今日は少し変わったネタで書いていきます。それも「ExcelでWAVファイルから周波数分析を行う」です!一体何を言っているんだ?という方も、これやったことあるという方も様々な反応が出ると思います。理系で学んでいると専攻によってはどこかで音声処理を学ぶ機会があります。ここでは最も身近な分析ソフトであるExcelを使って解析をしてみましょう。ロジックの勉強にもなって一石二鳥です。少し重ための内容になりますが、是非お付き合いください!
では始めていきましょう!
まずは、事前知識から。
さて、音とは何でしょうか。音は空気の振動で、疎密波と呼ばれる波です。例えばこんな感じのものです。
実際に私たちが聞く音は、連続した波形になります。いわゆるアナログデータというやつですね。一方、コンピュータで扱う値はごく短い時間で区切って音の強弱を数値で表現した離散値を取ります。こちらはいわゆるデジタルデータと呼ばれるものですね。上記の図はデジタルデータをAudacityというソフトで可視化したものになります。
では、この音を分析していくにはどうしましょう?巷ではよく「データ分析!」と叫ばれることが増えましたが、ではどうやって分析するの?となると途端に「わからない!」となる方も多いのではないでしょうか。一度やり方を見てしまえば、自分でもできるようになると思いますので、今回はこんな風にやるんだという一例として見てみてください。話がそれました、分析手法についてですが、今回は「(離散)フーリエ変換」を用います。フーリエ変換とは簡単に言うとある音(の波形)の中にどのような周波数が含まれているかを調べる手法です。入力信号(いわゆるこちらを表示したものが波形になります)の種類で分析手法が変化します。アナログで周期的な信号の場合はフーリエ級数展開を用い、アナログで周期的な信号の場合はフーリエ変換を用い、デジタル信号の場合は離散フーリエ変換を用います。今回はデジタル値を扱いますので、離散フーリエ変換ですね。
離散フーリエ変換について少し説明しておきます。波を表す関数は
と表すことができます。関数は三角関数の和の形で表すことができると言えます(オイラーの公式より、三角関数に変換可能)。係数は各周波数における周波数成分の大きさを表しているので、2番目の式に関してkの関数F(k)とみなすと、
となります。デジタル値は離散値をとるので、積分を総和に置き換えてあげると、
と表せます。これで各に対応する周波数成分の割合を特定することができます。
今回プログラムで実装するにあたって、この離散フーリエ変換を限定的な場面で高速に実行できる高速フーリエ変換を用いますが、原理は同じだと考えてもらって大丈夫です。この実装はこちらの記事を参考にして使わせていただきます。
【FFT】
http://blog.livedoor.jp/sce_info3-craft/archives/7774799.html
次はコンピュータにおいて、音の取り扱いについて解説します。マイクに入力された音はごく短い時間で音量が正規化されて振幅が保存されます。このデータを保存しておくのがWAVファイルになります。WAVファイルは音の波形が何も加工されずに直接書かれている形のファイルです。最近はあまり見なくなってしまった気がしないでもないですが、一番音質が良い形式の一つとなります。(いわゆる無劣化保存ができます)
Excelでこのファイルを読み込むにはVBAを使わないといけません。こちらも先人がありがたいことに読み込みのための記事を書いてくださっていたので、少し手を加えて使用します。
【WAVファイル読み込み】
http://otktake.blogspot.com/2015/01/excelwave.html
さて、基礎的なところは解説したので、次は実際に読み込んで解析するところに入っていきましょう。
まずはExcelを起動し、Alt + F11を押してVBEを起動します。起動したら上の挿入メニューから標準モジュールとクラスモジュールを追加し、次の図のように名前を変更します。
ここまで出来たら次はWAVファイルの読み込みから作っていきます。今回解析する曲は「カエルの歌」にでもしておきます(笑)
ここからダウンロードしてください。この音源を使用します。
先ほどの記事を参考に、書きます。VarModuleに次のコードを書きます。
'共通変数
'**************************************************************WAVEファイル***************************************************************************************
Public Type RIFF_CHUNK
ChunkID As String * 4 '"RIFF"(0x52494646)で固定
ChunkSize As Long 'チャンクサイズ。ファイル全体のサイズからファイル全体サイズからRIFFとWAVEのバイト数(8Byte)を引いた数。この情報をもとにWavファイルのファイルサイズを算出できる。
Format As String * 4 '"WAVE"(0x57415645)で固定
End Type
Public Type fmt_chunk
SubchunkID As String * 4 '"fmt"(0x666D7420)で固定
SubchunkSize As Long 'fmtチャンクのバイト数。リニアPCMならば16(0x10000000)。その他は16+拡張パラメータ
AudioFormat As Integer '音声フォーマット。非圧縮のリニアPCMフォーマットは1(0x0100)。
NumChannels As Integer 'チャンネル数。モノラルは1、ステレオは2
SamplingRate As Long 'サンプリング周波数(Hz)。44.1kHzの場合なら(0x44AC0000)
ByteRate As Long '一秒あたりのバイト数の平均。サンプリング周波数*ブロックサイズで求める。44.1kHz、16bit、モノラルならば44100x2x1=88,200(0x88580100)
BlockAlign As Integer 'ブロックサイズ。チャンネル数 * 1サンプルあたりのビット数 / 8で求める。モノラル16bitなら1*16bit = 16bit = 2byte(0x0200)、ステレオ16bitなら4(0x0400)
BitsPerSample As Integer 'ビット/サンプル。1サンプルに必要なビット数。8ビットの場合は8(0x0800)、16ビットの場合は16(0x1000)など。
End Type
Public Type DATA_CHUNK
SubchunkID As String * 4 '"data"(0x64617461)で固定
SubchunkSize As Long '波形データのバイト数。
End Type
'******************************************************************************************************************************************************************
そして、WAVクラスモジュールに次のコードを書き込みます。
'**************************************************再定義**************************************************************************
Private riffChunk As RIFF_CHUNK
Private fmtChunk As fmt_chunk
Private dataChunk As DATA_CHUNK
Private data() As Integer
'**********************************************************************************************************************************
'シートに音階と周波数の対応表がある
Private Pitch_Hz(87) As Double
Private Pitch_Code(87) As String
Const PI = 3.14159265358979
Private Sampling As Long
Private Peek_Hz_Before As Long
'WAVファイルを読み込む(44.1kHz/16bit 44.1kHz/8bit 両方に対応)
Public Function Read_WaveFile(ByVal WaveFilePath As String) As Integer()
Open WaveFilePath For Binary As #1 'WAVファイルをバイナリ形式で開く
Get #1, , riffChunk 'RIFF識別子の読み込み
Get #1, , fmtChunk 'fmt識別子の読み込み
Get #1, , dataChunk '波形データに関する情報の読み込み
ReDim data(dataChunk.SubchunkSize / 2) As Integer 'dataの配列数をファイル形式に合わせて再定義
Get #1, , data '波形データの読み込み
Close #1
Sampling = fmtChunk.SamplingRate
Read_WaveFile = data '波形データを返す
End Function
これで読み込みができるようになりました。試しに読み込んで表示するプログラムを書いてみましょうか。
MainModuleに次のように書きます。
Sub ReadTest()
On Error Resume Next
Dim CWAV As New WAV 'WAVクラスモジュールを使う
Const FilePath As String = "D:\Song.wav" 'ファイルパス
Dim data() As Integer: data = CWAV.Read_WaveFile(FilePath) 'データの読み込み
Dim i As Long '表示用
For i = 0 To 4410 '0.1秒分表示(サンプリング周波数が44100Hzなので)
Cells(i, 1).Value = data(i) 'とりあえずA列に書き出す
Next
End Sub