张春成
V2
2022/10/27阅读:36主题:默认主题
采样率越高越好
采样率越高越好
如果你的目的是寻找频率特征,那么对信号的采样率是越高越好。
这种好并不是无脑的好,而是有数字依据的好。
采样与奈奎斯特频率
信号处理的基本理论表明,采样会导致信号的数字功率谱发生变化,采样频率的一半称为奈奎斯特频率。只有奈氏频率之下的频率成分,才能在经过正、反傅里叶变换后“复原”回来。
但这是只是“可以”,而不是“必定”,因此这是一个必要条件。
事实上,并没有人保证,经过采样和频率变换的信号,一定能“复原”回来。
本文就通过计算模拟的方式,说明一个问题,那就是
即使在奈氏频率高于感兴趣频带最高频率的情况下,采样率过低也会导致信号失真。
因此,我认为,采样频率的追求越高越好,没有上限。
实验模拟
首先,生成一个随机信号,我们假设采样频率是 1000 Hz
.png)
原信号
之后计算它的频谱,不仅如此,还对信号进行降采样到 20、50、100 Hz,它们的频谱图如下。我们有理由认为蓝色线是真实曲线,而其他颜色分别是对它进行“降采样”妥协的结果。
.png)
多降采样的频谱特征
接下来,我们通过对应频率相减的方式量化降采样带来的绝对误差,可以看到,降采样越多,则频率误差越大。表现在频率特征上,就是会多出很多莫名其妙的峰值
.png)
绝对误差
最后,我们将减值除以能量值,得到的是“比例误差”,即真实信号失真的比例有多大。可以看到,误差相当大。
.png)
比例误差
比较代码
为了可追溯起见,本文的代码如下
import sklearn.preprocessing
import numpy as np
import pandas as pd
import plotly.express as px
import plotly.subplots as subplots
import plotly.graph_objects as go
from tqdm.auto import tqdm
from noise import pnoise1
# Generate data
def generate_data(repeat=50, num=1000):
data = []
seeds = [e for e in range(num)]
np.random.shuffle(seeds)
seeds = seeds[:repeat]
for seed in tqdm(seeds):
octaves = 0 * 4 + 1 * np.random.randint(5, 7)
persistence = 0 * 0.8 + 1 * np.random.uniform(0.7, 0.9)
s = np.array([
pnoise1(e, persistence=persistence, base=seed, octaves=octaves)
for e in np.linspace(0, 1, num)
])
data.append(s)
data = np.array(data)
data = sklearn.preprocessing.minmax_scale(data, axis=1)
for e in data:
e -= np.mean(e)
print(data.shape)
return data
data = generate_data(repeat=500)
px.imshow(data).show()
# Toolbox of comput frequency power different toolbox
def draw_diff(df):
a = df.query('type == "Sample at 1000.0Hz"')
b = df.query('type != "Sample at 1000.0Hz"')
c = pd.merge(b, a, on='freq')
c['diff'] = np.abs(c['engy_x'] - c['engy_y'])
c['diff_ratio'] = c['diff'] / c['engy_y']
c.index = range(len(c))
c.loc[0, 'type_x']='Sample at 1000.0Hz'
px.line(c, x='freq', y='diff', color='type_x', title='Diff', log_y=True).show()
px.line(c, x='freq', y='diff_ratio', color='type_x', range_y=[0, 1], title='Diff Ratio').show()
# Compute downsample
sample_rate = 1000
def mk_downsample(downsample, signal, offset=0, sample_rate=sample_rate):
if offset > downsample:
offset = downsample - 1
y = signal.copy()
y = y[offset::downsample]
fmax = sample_rate / downsample
n = len(y)
fft = np.fft.fftshift(np.fft.fft(y))
freqs = np.fft.fftshift(np.fft.fftfreq(n)) * fmax
df = pd.DataFrame()
df['freq'] = freqs
df['engy'] = np.abs(fft) / n
df['type'] = 'Sample at {}Hz'.format(fmax)
return df
# Try once
dfs = []
signal = np.mean(data, axis=0)
for down_sample in tqdm([1, 10, 20, 50], 'Simulation'):
dfs.append(mk_downsample(down_sample, signal))
df = pd.concat(dfs)
px.line(df, x='freq', y='engy', color='type', log_y=True).show()
draw_diff(df)
# Try twice
dfs = []
signal = np.mean(data, axis=0)
for down_sample in tqdm([1, 10, 20, 50], 'Simulation'):
dfs.append(mk_downsample(down_sample, signal, offset=5))
df = pd.concat(dfs)
px.line(df, x='freq', y='engy', color='type', log_y=True, ).show()
draw_diff(df)
# Origin signal
px.line(signal).show()
作者介绍
张春成
V2