「インターホンの音を検知してLINE通知することで、どこにいてもイヤホンをしていても来客に気づきたい」と考えて作成したものをまとめています。
前回「(2)PyAudio録音時の警告・エラーに対処する」では、録音時に現れる警告やエラーへの対処方法をまとめました。録音自体は、「(1)インターホンの音を録音する」にまとめました。
今回は、録音したインターホンの音を解析して、「どういう音ならインターホンとみなすのか？」について検討します。
インターホン以外の生活音で誤作動しないようにしたいですよね。なお、最終的には運用しながら微調整していきます。

装置：Raspberry Pi 4
マイク：共立プロダクツ MI-305 [USBマイク]
プログラム言語：Python3 (PyAudioモジュールを使用。 PyAudio Documentation)

方針

「ピーンポーン　ピーンポーン」
という音が我が家のインターホンの音である。
音をFFT（高速フーリエ変換）して周波数分解したら、「ピーン」と「ポーン」の音の高さに相当する周波数にピークが出るはず。
これのピーク強度で判定したい。
また、通知までのラグを少なくしたいので、１つ目の「ピーン」だけで判定したい。
「ピーンポーン　ピーンポーン」が全体で４秒程度なので、音を１秒間拾ってFFTすることを想定する。

Pandasのインストール

参考：ラズパイにpandasをインストールする方法
データ解析と言えばpandasのdataframe。Raspberry Piの場合は、pipで入れると正常に動かないらしい。
（自分の場合、そもそもインストールが終わらなかった）
apt-getでインストール。後で使用する。

pi@raspberrypi:~ $ sudo apt-get install python3-pandas

音をプロットしてみる

参考：ラズパイにpandasをインストールする方法

plot_wave.py

import wave
import numpy as np
import matplotlib.pyplot as plt
 
file_name = "output.wav" # 録音ファイル
RATE = 44100 # 録音時に設定したRATE

wf = wave.open(file_name, "rb")
data = np.frombuffer(wf.readframes(wf.getnframes()), dtype='int16')
wf.close()

x = np.arange(data.shape[0]) / RATE
plt.plot(x, data)
plt.show() # 横軸:時間(sec)
x = [i for i in range(len(data))]
plt.plot(x, data)
plt.show() # 横軸:データ点インデックス

上記の２つの画像が表示される。１つ目は横軸が時間(秒)、２つ目は横軸がデータ点のインデックス番号になっている。
（なおRaspberry Pi上で実行した場合は、１つ目の画像ウインドウを閉じないと次が表示されない。）
「ピーンポーン　ピーンポーン」の波形が明確に見て取れる。

音をFFT（高速フーリエ変換）してみてる

参考：FFTでインターホンの音を検知する（Python）

plot_FFT_all.py

import wave
import numpy as np
import matplotlib.pyplot as plt
 
file_name = "output.wav" # 録音ファイル
RATE = 44100 # 録音時に設定したRATE

wf = wave.open(file_name, "rb")
data = np.frombuffer(wf.readframes(wf.getnframes()), dtype='int16')
wf.close()

fft_data = np.abs(np.fft.fft(data))    #FFTした信号の強度
freqList = np.fft.fftfreq(data.shape[0], d=1.0/RATE)    #周波数（グラフの横軸）の取得
plt.plot(freqList, fft_data)
plt.xlim(0, 5000)    #0～5000Hzまでとりあえず表示する
plt.show()

全体をFFTするとこんな感じ。横軸が周波数 (Hz)、縦軸が強度。
ちなみに人が聞き取れる音は20 Hzから20,000 Hz, 88鍵のピアノは27.5 Hz~4186 Hz, 男声は100 Hz, 女声は300 Hzらしい（参考）。
意外とピークが多い・・・。
続いて、１秒間だけ切り出してFFTしてみる。

plot_FFT_partial.py

import wave
import numpy as np
import matplotlib.pyplot as plt
 
file_name = "output.wav" # 録音ファイル
RATE = 44100 # 録音時に設定したRATE
CHUNK = 1024 * 8 # 録音時に設定したCHUNK
RECORD_SECONDS = 1 # 検出に使いたい秒数
pnts = int(RATE / CHUNK * RECORD_SECONDS) * CHUNK # dataが何点になるかを計算

start = 0 # ここをいろいろ変えてみる

wf = wave.open(file_name, "rb")
data = np.frombuffer(wf.readframes(wf.getnframes()), dtype='int16')
wf.close()

data = data[start:start+pnts]
fft_data = np.abs(np.fft.fft(data))    #FFTした信号の強度
freqList = np.fft.fftfreq(data.shape[0], d=1.0/RATE)    #周波数（グラフの横軸）の取得
plt.plot(freqList, fft_data)
plt.xlim(0, 5000)    #0～5000Hzまでとりあえず表示する
plt.show()

一部を切り取ってFFTしたもの（※）。start = 0が左、60000が右。startの値は、上で示した「横軸がデータ点インデックスのグラフ」を見ながら決めた。左が「ピーン」の音、右が「ポーン」の音らしい。左の音を検出したい。
※：１秒間の音を検出に使いたいので、１秒間を切り取るようにしている。

検出に使うデータ点を決める

録音したデータが「ピーン」の音かどうかを判定するには、上記のピークの周波数を覚えておき、録音データがその周波数に十分な強度を持つかどうかを判定すれば良い。
今回は周波数ではなく「dataの中の何番目の点か」を覚えておいて、そこでの強度を足し合わせることにした。

print_indices.py

import wave
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
 
file_name = "output.wav" # 録音ファイル
RATE = 44100 # 録音時に設定したRATE
CHUNK = 1024 * 8 # 録音時に設定したCHUNK
RECORD_SECONDS = 1 # 検出に使いたい秒数
pnts = int(RATE / CHUNK * RECORD_SECONDS) * CHUNK # dataが何点になるかを計算

start = 0 # ここをいろいろ変えてみる

wf = wave.open(file_name, "rb")
data = np.frombuffer(wf.readframes(wf.getnframes()), dtype='int16')
wf.close()

data = data[start:start+pnts]
fft_data = np.abs(np.fft.fft(data))
freqList = np.fft.fftfreq(data.shape[0], d=1.0/RATE)

df = pd.DataFrame(dict(freq = freqList, amp = fft_data))
df = df[df['freq']>500] # 500 Hz以下は無視する。
df = df[df['amp']>0.5e7] # 0.5e7以上の強度を持つ点を覚える。
print(list(df.index))
# [610, 611, 612, 613, 615, 616, 1831, 1832, 1833, 1834, 1835, 1836, 3056, 3057, 3058, 3059, 4277, 4278, 4280, 4281, 4282, 4283, 4285]
print(list(df['freq']))
# [656.7626953125, 657.83935546875, 658.916015625, 659.99267578125, 662.14599609375, 663.22265625, 1971.36474609375, 1972.44140625, 1973.51806640625, 1974.5947265625, 1975.67138671875, 1976.748046875, 3290.2734375, 3291.35009765625, 3292.4267578125, 3293.50341796875, 4604.87548828125, 4605.9521484375, 4608.10546875, 4609.18212890625, 4610.2587890625, 4611.33544921875, 4613.48876953125]
print(len(df))
# 23

最後のところでpandasを使ってデータをフィルタリングしている。今回は全部で23点残った。
点数を調整したければ0.5e7になっているところを調整する（0.5e7 = 0.5*10**7）
ちなみに周波数としては658 Hzあたりと、その3, 5, 7倍のところになっている。
（本当にこう鳴っているのか、あるいは、FFTのときだけ奇数倍の周波数が現れる理由があるのか・・・？）

ここでprintされたindexを使って、

print_amp.py

import wave
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
 
file_name = "output.wav" # 録音ファイル
RATE = 44100 # 録音時に設定したRATE
CHUNK = 1024 * 8 # 録音時に設定したCHUNK
RECORD_SECONDS = 1 # 検出に使いたい秒数
pnts = int(RATE / CHUNK * RECORD_SECONDS) * CHUNK # dataが何点になるかを計算

# ここに先ほどの結果を入れる
freq_indices = [610, 611, 612, 613, 615, 616, 1831, 1832, 1833, 1834, 1835, 1836, 3056, 3057, 3058, 3059, 4277, 4278, 4280, 4281, 4282, 4283, 4285]

start = 0 # ここをいろいろ変えてみる

wf = wave.open(file_name, "rb")
data = np.frombuffer(wf.readframes(wf.getnframes()), dtype='int16')
wf.close()

data = data[start:start+pnts]
fft_data = np.abs(np.fft.fft(data))    #FFTした信号の強度

amp = 0
for i in freq_indices:
    amp += fft_data[i]

print('{:.2e}'.format(amp))
# 3.24e+08

を実行すれば、ピーク位置での強度を足し合わせた値が出てくる。
これが閾値（運用しながら調整）よりも大きいか小さいかで判定できそう。

なお、おそらく、「FFTデータの何番目の点が何Hzか」は元データの点数やRATEに依存する。
録音条件や、録音時間を変更した場合は再度freq_indicesを調べる必要があるかも。

まとめ

今回は、インターホンの録音データを調べながら、どのように検知するかについて検討しました。
「FFTしたとき、特定の点（周波数）の合計値が閾値以上かどうか」で検知できそうだということがわかりました。
次回は、この方式を使って実際に検知してみます。

2022/03/10: 検知基準に「他の周波数での強度との比」を加えると精度が向上しました。第4回に記載。

その他の記事：
Raspberry Piでインターホンの音を検知してLINEに通知する (1)インターホンの音を録音する
 Raspberry Piでインターホンの音を検知してLINEに通知する (2)PyAudio録音時の警告・エラーに対処する
 Raspberry Piでインターホンの音を検知してLINEに通知する (4)検知してLINEに通知する

参考

ラズパイにpandasをインストールする方法
 FFTでインターホンの音を検知する（Python）
AGC Glass Plaza > ８．音−１（「音」の基本から考えましょう）