|
|
验证方法是：
读取 piano_binaural_10hz.wav
截取中间一段，例如 10 秒附近的 2 秒
分别分析左声道和右声道
对信号加 Hann 窗
做 FFT
绘制 150–250 Hz 局部频谱
在图中明确标出 200 Hz 和 210 Hz
同时打印 200 Hz、210 Hz 附近的幅度
最后比较左右声道：
左声道应该重点看到 200 Hz
右声道应该重点看到 210 Hz
import numpy as np
import soundfile as sf
import matplotlib.pyplot as plt
# ============================================================
# 1. 参数设置
# ============================================================
# 我们要验证的 WAV 文件
audio_file = "piano_binaural_10hz.wav"
# ------------------------------------------------------------
# FFT 分析开始时间
# ------------------------------------------------------------
#
# 不建议从 0 秒开始。
#
# 因为我们之前给双耳信号设置了 3 秒淡入，
# 所以这里从 10 秒开始分析比较合适。
#
start_time = 10.0
# ------------------------------------------------------------
# FFT 分析时长
# ------------------------------------------------------------
#
# 这里使用 2 秒。
#
# 分析时间越长：
#
# FFT 频率分辨率越高
#
# 例如：
#
# Fs = 44100
# T = 2 秒
#
# N = 88200
#
# Δf = Fs / N
# = 44100 / 88200
# = 0.5 Hz
#
# 这样就可以比较清楚地区分：
#
# 200 Hz
# 210 Hz
#
analysis_duration = 2.0
# ------------------------------------------------------------
# 重点观察的频率范围
# ------------------------------------------------------------
#
# 我们不看整个频谱。
#
# 只观察：
#
# 150 Hz ~ 250 Hz
#
# 这样 200 Hz / 210 Hz 会非常明显。
#
min_freq = 150
max_freq = 250
# ============================================================
# 2. 读取 WAV
# ============================================================
audio, sample_rate = sf.read(
audio_file
)
print("========================================")
print("音频基本信息")
print("========================================")
print(
"采样率:",
sample_rate,
"Hz"
)
print(
"音频数据形状:",
audio.shape
)
# ============================================================
# 3. 检查是否为立体声
# ============================================================
if audio.ndim != 2:
raise ValueError(
"这个 WAV 文件不是立体声文件。"
)
if audio.shape[1] != 2:
raise ValueError(
"这个 WAV 文件必须包含左右两个声道。"
)
print(
"声道数量:",
audio.shape[1]
)
# ============================================================
# 4. 计算音频总时长
# ============================================================
total_samples = len(audio)
duration = (
total_samples
/ sample_rate
)
print(
f"音频时长: {duration:.2f} 秒"
)
# ============================================================
# 5. 检查分析时间
# ============================================================
if start_time >= duration:
raise ValueError(
"start_time 超过了音频总长度"
)
# ============================================================
# 6. 截取分析片段
# ============================================================
start_sample = int(
start_time
* sample_rate
)
end_sample = int(
(start_time + analysis_duration)
* sample_rate
)
end_sample = min(
end_sample,
total_samples
)
segment = audio[
start_sample:end_sample
]
actual_duration = (
len(segment)
/ sample_rate
)
print()
print("========================================")
print("FFT 分析参数")
print("========================================")
print(
"开始时间:",
start_time,
"秒"
)
print(
"分析时间:",
actual_duration,
"秒"
)
print(
"分析采样点:",
len(segment)
)
# ============================================================
# 7. 分离左右声道
# ============================================================
# 左声道
left_signal = segment[:, 0]
# 右声道
right_signal = segment[:, 1]
# ============================================================
# 8. 对左右声道加 Hann 窗
# ============================================================
# 为什么要加 Hann 窗？
#
# 我们截取的是一段有限长度的音频。
#
# 如果直接截断：
#
# 原始信号 ──────────────┐
# │
# └── 突然变成 0
#
# FFT 会产生比较严重的频谱泄漏。
#
# Hann 窗可以降低这个问题。
#
window = np.hanning(
len(left_signal)
)
left_windowed = (
left_signal
* window
)
right_windowed = (
right_signal
* window
)
# ============================================================
# 9. FFT
# ============================================================
left_fft = np.fft.rfft(
left_windowed
)
right_fft = np.fft.rfft(
right_windowed
)
# ============================================================
# 10. 计算频率轴
# ============================================================
freqs = np.fft.rfftfreq(
len(left_windowed),
d=1 / sample_rate
)
# ============================================================
# 11. 计算 FFT 幅度
# ============================================================
left_magnitude = np.abs(
left_fft
)
right_magnitude = np.abs(
right_fft
)
# ============================================================
# 12. 转换为 dB
# ============================================================
# 使用各自声道的最大值进行归一化。
#
# 这样主要是为了方便比较局部频谱。
#
left_magnitude = (
left_magnitude
/ (left_magnitude.max() + 1e-12)
)
right_magnitude = (
right_magnitude
/ (right_magnitude.max() + 1e-12)
)
left_db = (
20
* np.log10(
left_magnitude + 1e-12
)
)
right_db = (
20
* np.log10(
right_magnitude + 1e-12
)
)
# ============================================================
# 13. 只提取 150 ~ 250 Hz
# ============================================================
freq_mask = (
(freqs >= min_freq)
&
(freqs <= max_freq)
)
display_freqs = freqs[
freq_mask
]
display_left_db = left_db[
freq_mask
]
display_right_db = right_db[
freq_mask
]
# ============================================================
# 14. 找到最接近 200 Hz 的 FFT 点
# ============================================================
target_200_index = np.argmin(
np.abs(freqs - 200)
)
target_210_index = np.argmin(
np.abs(freqs - 210)
)
# 实际 FFT 频率
actual_200_freq = freqs[
target_200_index
]
actual_210_freq = freqs[
target_210_index
]
# ============================================================
# 15. 打印 200 Hz / 210 Hz 的结果
# ============================================================
print()
print("========================================")
print("200 Hz / 210 Hz 检测结果")
print("========================================")
print()
print(
f"目标 200 Hz"
)
print(
f"实际 FFT 频率: "
f"{actual_200_freq:.2f} Hz"
)
print(
f"左声道幅度: "
f"{left_db[target_200_index]:.2f} dB"
)
print(
f"右声道幅度: "
f"{right_db[target_200_index]:.2f} dB"
)
print()
print(
f"目标 210 Hz"
)
print(
f"实际 FFT 频率: "
f"{actual_210_freq:.2f} Hz"
)
print(
f"左声道幅度: "
f"{left_db[target_210_index]:.2f} dB"
)
print(
f"右声道幅度: "
f"{right_db[target_210_index]:.2f} dB"
)
# ============================================================
# 16. 绘制左右声道频谱
# ============================================================
plt.figure(
figsize=(14, 7)
)
# ------------------------------------------------------------
# 左声道
# ------------------------------------------------------------
plt.plot(
display_freqs,
display_left_db,
label="Left channel"
)
# ------------------------------------------------------------
# 右声道
# ------------------------------------------------------------
plt.plot(
display_freqs,
display_right_db,
label="Right channel"
)
# ------------------------------------------------------------
# 标记 200 Hz
# ------------------------------------------------------------
plt.axvline(
200,
linestyle="--",
label="200 Hz"
)
# ------------------------------------------------------------
# 标记 210 Hz
# ------------------------------------------------------------
plt.axvline(
210,
linestyle="--",
label="210 Hz"
)
# ------------------------------------------------------------
# 图形设置
# ------------------------------------------------------------
plt.xlabel(
"Frequency (Hz)"
)
plt.ylabel(
"Magnitude (dB)"
)
plt.title(
"Piano + 10 Hz Binaural Beat: 150–250 Hz Spectrum"
)
plt.xlim(
min_freq,
max_freq
)
plt.ylim(
-80,
5
)
plt.grid(
True,
alpha=0.3
)
plt.legend()
plt.tight_layout()
plt.show()
# ============================================================
# 17. 单独放大左声道
# ============================================================
plt.figure(
figsize=(14, 6)
)
plt.plot(
display_freqs,
display_left_db
)
# 200 Hz 参考线
plt.axvline(
200,
linestyle="--",
label="Target 200 Hz"
)
# 210 Hz 参考线
plt.axvline(
210,
linestyle="--",
label="Target 210 Hz"
)
plt.xlabel(
"Frequency (Hz)"
)
plt.ylabel(
"Magnitude (dB)"
)
plt.title(
"Left Channel: 200 Hz Verification"
)
plt.xlim(
150,
250
)
plt.ylim(
-80,
5
)
plt.grid(
True,
alpha=0.3
)
plt.legend()
plt.tight_layout()
plt.show()
# ============================================================
# 18. 单独放大右声道
# ============================================================
plt.figure(
figsize=(14, 6)
)
plt.plot(
display_freqs,
display_right_db
)
# 200 Hz
plt.axvline(
200,
linestyle="--",
label="200 Hz"
)
# 210 Hz
plt.axvline(
210,
linestyle="--",
label="Target 210 Hz"
)
plt.xlabel(
"Frequency (Hz)"
)
plt.ylabel(
"Magnitude (dB)"
)
plt.title(
"Right Channel: 210 Hz Verification"
)
plt.xlim(
150,
250
)
plt.ylim(
-80,
5
)
plt.grid(
True,
alpha=0.3
)
plt.legend()
plt.tight_layout()
plt.show()
print()
print("========================================")
print("FFT 分析完成")
print("========================================")
分析结果如下：
在左声道的200Hz处出现波峰，在右声道的210Hz处出现波峰
|
|
