[离散时间信号处理学习笔记] 13. 重采样

重采样常用于音频处理。在用麦克风对音频进行采集的时候，常见的采样率有8k（电话）、44.1k（CD）、48k（视频音轨）、96k/192k（Hi-Res），而某些系统会有默认固定的输出采样率（如Android的默认输出采样率为44.1k），此时就需要对输入音频数据进行重采样。

重采样的源样本序列为 $x[n]$

$x[n] = x_c(nT)$

重采样的目标序列为 $x'[n]$

$x'[n] = x_c(nT')$

如何通过 $x[n]$ 得到 $x'[n]$ 就是本文的讨论内容。

本文假设以采样周期为 $T$ 对 $x_c(t)$ 进行采样满足奈奎斯特采样定律。

减采样（downsampling）

减小采样率的过程被称为减采样，这一小节讨论的是按整数倍减小采样率。

按照我们一般的思维来说，按整数倍（倍数为 $M$ ）减少采样率应该是直接对源样本序列每隔 $M$ 个样本提取一个值

$x_d[n] = x[nM] = x_c(nMT)$

这种提取方法被称为采样率压缩器，简称压缩器（compressor）。可以看到所得的新序列是原始连续信号的一部分，并且新序列的采样周期为 $T_d = MT$ 。对于该新序列，我们可以分为两种情况进行讨论：

$T_d$ 符合奈奎斯特采样定理，即新序列能通过一个低通滤波器还原为原始的连续信号
$T_d$ 不符合奈奎斯特采样定理，即新序列发生混叠，无法还原为原始的连续信号

如下图假设信号在 $M=2$ 时恰好满足奈奎斯特采样定理，那么在 $M=3$ 时则会发生混叠

如果在采用了压缩率为 $M$ 的压缩器后，序列仍然符合奈奎斯特采样定理，我们可以直接进行使用 $x_d[n] = x[Mn]$ 来得到减采样序列。而发生混叠的情况则稍微复杂一点。观察混叠的频谱，可以发现只有低频部分保持了与原始信号频谱的一致性，而相当多的高频由于混叠而失去了原始频谱。

频谱丢失得越多说明信号的失真越大，因此为了减少失真，需要尽可能保留更多的原始信号频谱。我们可以先对元素信号进行低通滤波，然后再对滤波后的信号进行周期为 $MT$ 的采样即可得到失真更少的序列。按照这种思想，采样周期固定为 $MT$ ，如果一个被采样信号的采样周期为 $MT$ ，那么采样后不会混叠的条件就是该信号的截至频率为 $\frac{\pi}{MT}$ ，因此低通滤波的截至频率为 $\frac{\pi}{MT}$ 。

很明显，先低通滤波后采样的这种方法能最大限度地减低频谱的丢失，从而降低信号失真。

$\tilde{X}_d(j\Omega) = X_c(j\Omega)H_{M}(j\Omega)$

其中 $\tilde{X}_d(j\Omega)$ 就是对原始信号进行低通滤波后的信号 $\tilde{x}_d(t)$ 的傅里叶变换，低通滤波器的傅里叶变换为 $H_{M}(j\Omega)$ 。不过我们手中的并不是原始信号 $x_c(t)$ ，而是序列 $x[n]$ ，为了进入上述流程，我们需要先用 $x[n]$ 重构出 $x_c(t)$ ，才能进行低通滤波以及后续采集

$\begin{align*} \tilde{x}_d(t) &= \mathcal{F}^{-1}\tilde{X}_d(j\Omega)\\ &= \mathcal{F}^{-1}\{X_c(j\Omega)H_{M}(j\Omega)\}\\ &= \mathcal{F}^{-1}\{X(e^{j\Omega T})H_r(j\Omega)H_{M}(j\Omega)\}\\ &= \mathcal{F}^{-1}\{X(e^{j\Omega T})TH_{M}(j\Omega)\} \qquad \left\{\begin{matrix} H_r(j\Omega)&=\left\{\begin{matrix}T, & |\Omega|<\pi/T\\ 0, & else\end{matrix}\right.\quad\\ H_M(j\Omega)&=\left\{\begin{matrix}1, &|\Omega|<\pi/MT \\ 0, & else\end{matrix}\right.\\ \end{matrix}\right .\\ &= \mathcal{F}^{-1}\{X_s(j\Omega)TH_M(j\Omega)\}\\ &= x_s(t)*[MT\cdot h_m(t)]/M\qquad fourier\ convolution\ theorem\\ &= \left\{\sum_{n=-\infty}^{\infty}x[n]\delta(t-nT)\right\}*\left\{\frac{sin(\pi t/MT)}{\pi t/MT} \right\}/M\\ &= \sum_{n=-\infty}^{\infty}x[n]\frac{sin[\pi(t-nT)/MT]}{\pi(t-nT)/MT}/M \end{align*}$

然后从 $\tilde{x}_d(t)$ 中以 $MT$ 为周期采集得到 $\tilde{x}_d[n]$

$\color{red}{\begin{align*} \tilde{x}_d[n] &= \tilde{x}_d(nMT)\\ &= \left .\sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(t-kT)/MT]}{\pi(t-kT)/MT}/M\right|_{t=nMT}\\ &= \sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(nMT-kT)/MT]}{\pi(nMT-kT)/MT}/M\\ &= \sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(nM-k)/M]}{\pi(nM-k)/M}/M \end{align*}}$

如果要从连续时间系统来理解的话，我们可以发现重建的连续信号 $\tilde{x}_d(t)$ 是由无数个 $x[k]$ 分别对相应位置的截至频率为 $\pi/MT$ 的sinc函数进行加权后叠加，然后对叠加得到的信号的幅度除以 $M$ 得到的。

如果要从离散时间系统来理解 $\tilde{x}_d[n]$ 的构造公式的话，可以发现上面公式中的sinc函数有如下规律（下二）

$\begin{matrix} H(e^{j\omega}) =\left\{\begin{matrix}M, &|\omega|<\pi/M \\0, &else \end{matrix}\right .& \Leftrightarrow & h[n] = \displaystyle{\frac{sin[\pi n/M]}{\pi n/M}}\\ H(e^{j\omega}) =\left\{\begin{matrix}1, &|\omega|<\pi/M \\0, &else \end{matrix}\right .& \Leftrightarrow & h[n]=\displaystyle{\frac{sin[\pi n/M]}{\pi n/M}/M} \end{matrix}$

可以发现该sinc函数是一个增益为1、截至频率为 $\pi/M$ 的低通滤波器，该滤波器后接因子为M的压缩器。那么，前面的公式可以理解为：用该滤波器对 $x[k]$ 进行滤波后再用因子为M的压缩器即可得到 $\tilde{x}_d[n]$ 。这一过程又被称为抽取（decimation）。

增采样（upsampling）

增加采样率的过程被称为增采样，这一小节讨论的是按整数倍增加采样率。

假设增加采样率的倍数为 $L$ ，那么增加采样率后的采样周期为 $T_i = LT$ ，有

$x_i[n] = x_c(nT_i) = x_c(nLT)$

不过实际上我们只有采样周期为 $T$ 的序列 $x[n]$ ，因此需要先进行 $x_c(t)$ 的重建，然后对 $x_c(t)$ 进行周期为 $T_i$ 的采样

$\color{red}{\begin{align*} x_i[n] &= x_c(nT/L)\\ &= \left .\sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(t-kT)/T]}{\pi(t-kT)/T}\right|_{t=nT/L}\\ &= \sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(nT/L-kT)/T]}{\pi(nT/L-kT)/T}\\ &= \sum_{k=-\infty}^{\infty}x[k]\frac{sin[\pi(n-kL)/L]}{\pi(n-kL)/L} \end{align*}}$

下面我们从频域来展开讨论

增采样即采样频率提高了，因此不会出现混叠的情况，即不会改变原始的连续时间信号。现假设用采样周期 $T$ 对 $x_c(t)$ 进行采样恰好满足奈奎斯特采样定理，即有

增采样后采样周期为 $T_i = T/L$ ，那么频谱将有如下变化

观察两张图右下角的 $X(e^{j\omega})$ 以及 $X_i(e^{j\omega})$ ，它们的时域表示分别就是我们的源序列 $x[n]$ 与目标序列 $x_i[n]$ 。从 $X(e^{j\omega})$ 变为 $X_i(e^{j\omega})$ 只需进行执行两个步骤：

对 $X(e^{j\omega})$ 的变量 $\omega$ 进行倍数为 $L$ 的扩展，得到 $X_e(e^{j\omega}) = X(e^{j\omega L})$
对所得的新频谱进行低通滤波，滤波器的增益为 $L$ ，截至频率为 $\frac{\pi}{L}$

对第一步有如下分析：

$\begin{align*} X(e^{j\omega}) &= \sum_{n=-\infty}^{\infty}x[n]e^{-j\omega n}\\ \Rightarrow X(e^{j\omega L}) &= \sum_{k=-\infty}^{\infty}x[k]e^{-j\omega Lk}\\ &= \sum_{n/L=-\infty}^{\infty}x[n/L]e^{-j\omega n} \qquad n=Lk\\ \Rightarrow x_e[n] &= \left\{\begin{matrix} x[n/L], &n = 0,\pm L, \pm 2L\cdot\cdot\cdot\\0, &else\end{matrix}\right. \end{align*}$

如下图所示，对序列 $x[n]$ 进行步长为 $L$ 的扩展，即可得到 $x_e[n]$ ，这种转换称为扩展器（expander）

第二步的低通滤波器的增益为 $L$ 、截至频率为 $\pi/L$ ，即其脉冲响应为（从前面的增采样公式同样也能得出该结论）

$h_i[n] = \frac{sin(\pi n/L)}{\pi n/L}$

因此增采样系统分解如下图，这一过程又被称为内插（interpolation）。

简单的内插滤波采样

在对序列进行内插时，某些情况下，我们并不用追求很准确的信号还原，此时用一些简单的滤波器即可达到不错的效果，如常见的线性内插、三次样条内插等。

$h_{lin}[n] =\left\{ \begin{matrix} 1-|n|/L , & |n|\leqslant L\\ 0, &else\end{matrix} \right. \qquad\qquad h_{cu}[n] = \left\{ \begin{matrix} (a+2)|n/L|^3 -(a+3)|n/L|^2+1, &0\leqslant |n|\leqslant L\\ a|n/L|^3-5a|n/L|^2+8a|n/L|-4a, &L\leqslant |n|\leqslant 2L\\ 0, &else \end{matrix} \right.$

线性内插只用到内插的位置两旁的两个样本，在图上表示的话，就是把两个样本用直线连接后，该直线在对应内插位置上的的值就是所求的内插值。
三次样条内插会用到内插位置两旁的四个样本，即左右各两个，曲线类似于对低通滤波器进行截取后的曲线（下面取 $L=5, a=-0.5$ ）。

观察上面两个脉冲响应，我们能发现这些简单内插器脉冲响应的一些规律：

零点的左右两边对称， $\tilde{h}[n] = \tilde{h}[-n]$
长度为 $2KL-1，K=1,2,\cdot\cdot\cdot$ ，也就是说左右两边距离零点超过KL时的值为0，即 $\tilde{h}[n] = 0,|n|\geqslant KL$ ·
零点处的值为1， $h[0] = 1$
KL处的值为0， $h[KL] = 0$

既然有这些简单但不够准确的内插器，那么就应该有相应的方法来判断这些内插器的内插效果。内插效果可以通过观察这些内插器的频谱来进行分析。

对比增采样所需要的低通滤波器跟简单的线性内插器频谱，可以发现在频谱的 $\omega<\frac{\pi}{L}$ 处本应做大小固定为 $L$ 的增益，但是线性内插在这部分有衰减，另外在 $\omega>\frac{\pi}{L}$ 处本应该是对原始序列频谱进行截断，但是由于线性内插器的频谱在 $\omega>\frac{\pi}{L}$ 仍然有较大的能量（幅度较高），因此可见线性内插可能不会获得很令人满意的内插效果。不过如果原本的采样频率就远大于原始信号的截至频率，则表明采样所得的序列的频谱较为集中，此时采用线性内插则会得到较好的效果。

对比线性内插以及三次样条内插，可以发现在 $\omega<\frac{\pi}{L}$ 处，三次样条的旁瓣较宽，而在 $\omega>\frac{\pi}{L}$ 处，三次样条的能量更小（幅度更低）,即三次样条更接近于低通滤波器，因此三次样条内插会比线性内插的效果更好。