４回目の記事では
MacBookAir(M1)のスピーカから１成分のサイン波形を再生して
それを録音してFFT解析を行います。
まずは、手持ちのMacBookのスピーカの周波数特性の解析を行います。

MacBook自体がモバイル製を重視しているので
そこまで周波数特性が良いとは思っていません。
ただし、部品への要求レベルが高いApple製のノートPCであり
それなりに良い周波数特性が期待できます。

どれぐらいのスペックなのでしょうか。　

以下に、コードを記載します。
再生と録音をスレッドを用いて並列処理で記載しています。

前処理

import threading
# import time
import wave
import pyaudio
import numpy as np
import matplotlib.pyplot as plt
# from datetime import datetime

nchannels = 1
f1 = 500

sampling_rate = 44100 # サンプリング周波数
sampwidth = 2 #16bit
total_secs = 5.0
out_wav_name = 'sin_out.wav'#作成する基準波形データ
in_wav_name = 'sin_in.wav' #録音された波形データ



chunk = 2**10       # チャンク（データ点数）
data_num = int(total_secs * sampling_rate) # 総データ点数
freq = np.linspace(0,sampling_rate*2, data_num)

解析用波形ファイルを作成する関数

# 波形ファイルを作成する
def gen_sound():
    t = np.linspace(0, total_secs, data_num)
    wav = 32767 * (np.sin(np.pi*f1*t))

    #窓関数の作成
    window_hn = np.hanning(data_num)
    wav_hnd = window_hn * wav
    plt.plot(t,wav_hnd)
    plt.show()

    f_hnd = np.fft.fft(wav_hnd)
    f_abs_hnd = np.abs(f_hnd)
    plt.plot(freq[:int(sampling_rate/100)+1], f_abs_hnd[:int(sampling_rate/100)+1])
    plt.show()

    data = wav_hnd.astype(np.int16)
    out_wave = wave.Wave_write(out_wav_name)
    out_wave.setnchannels(nchannels)
    out_wave.setsampwidth(sampwidth)
    out_wave.setframerate(sampling_rate)
    out_wave.writeframes(data)
    out_wave.close()

保存した音声ファイルを再生する関数

# 保存した音声ファイルを再生する
def playing():
    wf = wave.open(out_wav_name, "r")
    # ストリームを開く
    p = pyaudio.PyAudio()
    stream = p.open(format=p.get_format_from_width(wf.getsampwidth()),
                        channels=wf.getnchannels(),
                        rate=wf.getframerate(),
                        output=True)

    data = wf.readframes(chunk)
    while data != b'':
            stream.write(data)
            data = wf.readframes(chunk)
    stream.close()
    p.terminate()

保存した音声ファイルをFFT解析する関数

 #保存した音声ファイルをFFT解析する
def analize():
    in_wav = wave.Wave_read(in_wav_name)
    data = in_wav.readframes(in_wav.getnframes())
    data = np.frombuffer(data, dtype='int16')
    plt.plot(data)
    plt.show()

    freq = np.linspace(0,sampling_rate*2, data_num)
    f = np.fft.fft(data)
    f_abs = np.abs(f)
    plt.plot(freq[:int(sampling_rate/10)+1], f_abs[:int(sampling_rate/10)+1])
    plt.show()

録音する。

def recording():
    #rec_time = 5            # 録音時間[s]
    fmt = pyaudio.paInt16  # 音声のフォーマット
    ch = 1              # チャンネル1(モノラル)
    chunk = 2**11       # チャンク（データ点数）
    
    
    audio = pyaudio.PyAudio()
    index = 0 # 録音デバイスのインデックス番号（デフォルト1）

    stream = audio.open(format=fmt, channels=ch, rate=sampling_rate, input=True,
                        input_device_index = index,
                        frames_per_buffer=chunk)
    print("recording start...")
    
    # 録音処理
    frames = []
    for i in range(0, int(sampling_rate / chunk * total_secs)):
        data = stream.read(chunk)
        frames.append(data)

    print("recording  end...")

    # 録音終了処理
    stream.stop_stream()
    stream.close()
    audio.terminate()

    # 録音データをファイルに保存
    wav = wave.open(in_wav_name, 'wb')
    wav.setnchannels(ch)
    wav.setsampwidth(audio.get_sample_size(fmt))
    wav.setframerate(sampling_rate)
    wav.writeframes(b''.join(frames))
    wav.close()

メインの処理

# メインの処理
def main():
   """並列処理"""
   gen_sound()
   t1 = threading.Thread(target=playing)
   t2 = threading.Thread(target=recording)
   t1.start()
   t2.start()
   t1.join()
   t2.join()
   analize()


if __name__ == "__main__":
    main()

作成した解析用の波形（５００　[Hz]の１成分）

録音した波形

録音した波形のFFT結果

基準周波数の500[Hz]に対して、７００[Hz]あたりでノイズが発生しています。
これはPC本体が共振している可能性があります。
いわゆる、箱鳴と呼ばれる小型スピーカ特有の現象かもしれません。

今回は１成分のサイン波形でしたが、基準波形に複数の周波数を混ぜることで
音圧の周波数特性が得られます。

また、歪率の計算は

その波形に含まれる全高調波成分（ E2 ～ En ）の実効値の総和と基本波（ E1 ）の実効値との比として定義される。
https://www.nfcorp.co.jp/techinfo/dictionary/042/
とのことなので、次回、算出してみます。

続きます。