北京邮电大学《数字信号处理》课件-第8章 多采样率数字信号处理

第8章

多采样率数字信号处理

1北京邮电大学《数字信号处理》本章内容要点：8.1引言 8.2整数因子抽取 8.3整数因子内插 8.4按有理数因子I/D的采样率转换8.5整数倍抽取和内插在数字语音系统中的应用 8.6采样率转换滤波器的高效实现方法 8.6.1直接型FIR滤波器结构 8.6.2多相滤波器结构 8.7采样率转换器的MATLAB实现

2

8.1引言

前面所讨论的信号处理的各种方法都是把采样率Fs视为固定值，即在一个数字系统中只有一个采样频率。但在实际系统中，要求一个数字系统能工作在“多采样率”状态。例如： （1）在数字电视系统中，图像采集系统要将标准采集数字电视信号，转换成其它标准的数字信号。这就要求数字电视演播室系统在多采样率状态。 （2）在数字电话系统中，传输的信号既有语音信号，又有传真信号，甚至有视频信号。所以，该系统应具有多采样率功能。 （3）对一个非平稳随机信号作谱分析或编码时，对不同的信号段，可根据其频率成分的不同而采用不同的采样率，以达到既满足采样定理。 （4）如果以高采样率采集的数据存在冗余，这时就希望在该数字信号的基础上降低采样速率。

38.1引言一般概念：在满足采样定理的前提下，将以采样率F1采集的数字信号进行D/A转换变成模拟信号；按采样率F2进行A/D变换，从而实现从F1到F2的采样率转换。问题与缺点：较麻烦，且易使信号受到损伤实际实现方法（本节学习的内容）：在数字域直接改变采样率。根据采样率转换理论，对采样后的数字信号x(n)直接进行采样率转换，以得到新采样率下的采样数据。采样率转换的分类：

“抽取（Decimation）”：降低采样率以去掉多余数据的过程；“插值（Interpolation）”：提高采样率以增加数据的过程。本章主要内容：（1）讨论抽取和插值的一般概念；（2）整数倍抽取和内插在数字语音系统中的应用（3）讨论几种基本的高效实现方法。

本章的内容是语音及图像数据压缩新技术——子带编码的重要理论基础。

48.2整数因子抽取5设x(n1T1)是连续信号xa(t)的采样序列，采样率F1=1/T1(Hz)，T1称为采样间隔，单位为秒，即（8.2.1）希望将采样率降低到原来的1/D，D为大于1的整数，称为抽取因子。最简单的方法是对x(n1T1)每D点抽取1点，抽取的样点依次组成新序列y(n2T2)。y(n2T2)的采样间隔为T2，采样率为F2=1/T2(Hz)，T2与T1的关系为（8.2.2）为了后面叙述方便，我们将上述的抽取系统用图8.2.1(a)表示，x(n)和y(n)分别如图8.2.1（b）和（c）所示。n1和n2分别表示x(n)和y(n)序列的序号。（8.2.3）当n1=n2D时，8.2整数因子抽取为了后面叙述方便，我们将上述的抽取系统用图8.2.1(a)表示，x(n)和y(n)分别如图8.2.1（b）和（c）所示。n1和n2分别表示x(n)和y(n)序列的序号。

当n1=n2D时，

68.2整数因子抽取

上面直接抽取方法存在的问题：

直接每隔D–1个抽取一个样值，相当于降低了采样频率，可能会引起频谱混叠现象。下面讨论抽取过程中可能出现的频谱混叠及改进措施。78.2整数因子抽取如果x(n1T1)是连续信号xa(t)的采样信号，且8f为模拟频率变量，ω1为数字频率。

，亦称为采样频率。8.2整数因子抽取

为了对抽样前后的频谱进行比较，作图时均以模拟角频率Ω为自变量（横坐标），为此按（8.2.6）式将写成Ω的函数为98.2整数因子抽取10由图8.2.3可见，直接抽取确实产生频谱混叠，所以随意对进行抽取是不行的。只有在抽取后仍能满足采样定理时才能恢复出原来的信号，否则就必须另外采取措施。通常采取的措施是抗混叠滤波。所谓抗混叠滤波，就是在抽取之前先对信号进行低通滤波，把信号的频带限制在以下。对应的数字频率为这种抽取系统框图如图8.2.4所示。所以，在理想情况下，抗混叠低通滤波器的频率响应为8.2整数因子抽取图8.2.4中各点的信号在时域和频域中的示意图如图8.2.5所示。这种办法虽然把中的高频部分损失掉了，但由于抽取后避免了混叠，所以在中完好无损地保留了中的低频部分，可以从中恢复出的低频部分。118.2整数因子抽取12抽取前后的频域关系：

经过数学推导，可以得到如下频域关系式：

8.2整数因子抽取13（8.2.14）式就是与是的D个平移样本之和，相邻的平移样本在频率轴上相差2π/D，在模拟频率轴Ω上相差如图8.2.7和图8.2.8所示。的关系式，图8.2.7产生混叠失真图8.2.8未产生混叠失真8.3信号的整数倍内插1．整数倍内插的概念与内插方法

整数I倍内插是在已知的相邻两个原采样点之间插入I–1个新采样值点。由于这I–1个采样值并非已知的值，所以关键问题是如何求出这I–1个采样值。从理论上讲，可以对已知的采样序列进行D/A转换，得到原来的模拟信号xa(t)，然后再对xa(t)进行较高采样率的采样得到，这里（8.3.1）式中I为大于1的整数，称为内插因子。上述过程可用图8.3.1表示。但这样的插入方法是不经济的，且易对信号产生损伤。实际工作中采用下述内插方法。

148.3信号的整数倍内插实际工作中采用下述内插方法整数内插是先在已知采样序列的相邻两个样点之间等间隔插入I–1个0值点，然后进行低通滤波，即可求得I倍内插的结果。这种内插方案,如图8.3.2所示。图8.3.3中I=3

15图8.3.2零值内插方案原理框图图8.3.3内插过程时域波形

8.3信号的整数倍内插2．整数倍内插的频域解释

上述的零值内插方案中的频谱关系是怎样的？回答这一问题的过程就是解释为什么经低通滤波就能得出采样率升高I倍的。这样才能提出对图8.3.2中低通滤波器的技术要求。 为了回答上面的问题，设为模拟信号xa(t)的采样序列，并假定xa(t)及其傅里叶变换如图8.3.4所示。

16图8.3.4模拟信号及其频谱示意图

8.3信号的整数倍内插按照内插的概念：且满足T2=T1/I所以，及其傅里叶变换应该如图8.3.5所示。均为周期函数，若二者都用模拟频率Ω表示，则，，周期为；，周期为17图8.3.58.3信号的整数倍内插下面分析图8.3.2中的频谱，最后讨论为了得到满足插值要求的（如图8.3.5所示），对的技术要求。

18上式表明和的频谱相同，如图8.3.6所示。图8.3.6图8.3.6中，，。

8.3信号的整数倍内插镜像频谱及其滤波器：与图8.3.5中的相比较，图8.3.6中的多出了从的部分，通常将这部分频谱称为镜像频谱。由此可见，要想从得到如图8.3.5（b）所示的，就必须滤除这些镜像频谱。所以，要求滤波器的理想低通幅频特性如图8.3.7所示。实际工作中，所以允许有一定的过渡带，可用线性相位FIR滤波器实现。根据其功能，将称为镜像滤波器。

19图8.3.78.3信号的整数倍内插将理想镜像滤波器的阻带截止频率换算成数字频率为所以，理想情况下，镜像滤波器的频率响应特性为 （8.3.4）式中，C为定标系数。因此输出频谱为 （8.3.5）定标系数C的作用是，在时，确保输出序列。为了计算简单，取m=0来求解C的值。因为，所以,由此得出，定标系数C=I。208.3信号的整数倍内插3．内插器的输入、输出关系（1）时域输入、输出关系 由图8.3.2有21所以由（8.3.3）式知道，所以，因为内插器时域输入、输出关系。（2）频域输入、输出关系8.3信号的整数倍内插复频域输入与输出的关系：由图8.3.2可知：22

（8.3.12）式中所有变量都为，所以可去掉下标得到所以8.4按有理数因子I/D的采样率转换1、按有理数因子I/D采样率转换的一般原理由按整数因子I内插和整数因子D抽取的原理，显然，可以用图8.4.1所示方案实现有理数因子I/D采样率转换。图8.4.1按有理数因子I/D的采样率转换方法图8.4.2按有理数因子I/D采样率转换的实用原理方框图应当注意，先内插后抽取才能最大限度地保留输入序列的频谱成分。 用分别表示输入序列和输出序列的采样频率，则另外，图中镜像滤波器和抗混叠滤波器级联，而且工作在相同的采样频率，因此完全可以将它们合成为一个等效滤波器，所以，按有理数因子I/D采样率转换的实用原理方框图如图8.4.2所示。理想情况下，和均为理想低通滤波器，所以，的等效滤波器仍是理想低通滤波器，其等效带宽应当是和中最小的带宽。因此，的频率响应为 （8.4.1）238.4按有理数因子I/D的采样率转换2、输出序列的时域表达式

24零值内插器的输出序列为

线性滤波器输出序列为整数因子D抽取器最后输出序列为，其时域表达式为如果线性滤波器用FIR滤波器实现，则可以根据式（8.4.4）计算输出序列

除了前面采样率变换技术，在实际工作中还会遇到任意因子采样率转换（为任意有限数）。有兴趣的读者请参考文献[24]。

【24】宗孔德。多抽样率信号处理。北京：清华大学出版社,1996

。

8．5整数倍抽取和内插在数字语音系统中的应用

25本节对信号时域和频域的表示方法和描述符号:为模拟信号，表示对的采样序列，是对进行采样率转换（内插或抽取）后的序列。其中，数字频率与模拟频率的关系为的采样频率记为的采样频率记为

和在数字语音系统中的应用原理，本节全部以模拟角频率Ω为自变量（横坐标），并采用上面定义的符号，来绘制应用系统中各信号的频谱曲线。设并定义为了通过观察比较的频谱关系，理解采样率转换8.5.1数字语音系统中的信号采样过程及其

存在的问题在数字语音系统中，语音信号的采样过程如图8.5.1所示。26图8.5.1语音信号的一般采样过程示意图

存在的问题:1、对预滤波器的技术要求太高（要求过渡带宽度为0，即理想低通滤波器）。因而是难以设计与实现的。2、显然，在接收端D/A变换过程中同样会遇到此问题。3、如果简单地将采样率提高到Fs=16kHz，则预滤波器就容易实现（允许有4kHz的过渡带），但使采集信号的数据量加大一倍，传输带宽也加大一倍，极大地浪费带宽。8.5.2数字语音系统中改进的A/D转换方案27

采用整数因子抽取与整数因子内插来改进的A/D转换方案

，既解决上述问题，又不增加数据量和传输带宽。解决问题的思路：改进的A/D转换方案：工作原理：见下页图解1、用较高的采样率进行A/DC，如采样率2、在数字域抗混叠滤波，按因子D=2抽取，把采样率降至8kHz。这样，允许模拟预滤波器g(t)的过渡带为，

8.5.2数字语音系统中改进的A/D转换方案288.5.3接收端D/A转换器的改进方案

29改进的D/A转换方案：对进行理想变换，得到说明：这种D/A转换器难以实现，实际中常用零阶保持型D/A转换器代替，但其频响特性不理想，会引入幅频失真，这种失真可在数字域进行预处理补偿。

对模拟低通滤波器要求：阻带：

通带边缘频率：

过渡带：

由于模拟低通滤波器的过渡带较宽，所以设计与实现较容易！8.6采样率转换滤波器的高效实现方法

问题的提出：

高效实现的原则：在多采样率系统中，总是设法把乘法运算安排在低采样率一侧，以使每秒钟内的乘法次数最少。采样率转换基本方案存在的问题及解决途径：在前面介绍的二种采样率转换方案（见图）中，滤波器的卷积运算均在采样率较高的一侧。因此，必须对多采样率系统的网络结构进行研究，以便得到乘法次数最少的高效实现结构。

教学内容：本节只介绍FIR直接实现和多相结构，其它两种高效实现结构（多级实现和时变网络）请参阅其他文献。这是因为用FIR结构实现多采样率系统具有很大的优越性。FIR结构绝对稳定且很容易做成线性相位，特别是容易实现高效结构。所以在多采样率系统的实现中绝大多数采用FIR滤波器。308.6.1直接型FIR滤波器结构

1．整数倍抽取器的FIR直接实现

31整数（D）倍抽取器框图如图8.2.4所示。抗混叠低通滤波器用FIR结构时，抽取器的时域输入、输出关系为（设h(n1T1)长度为N）

滤波器用FIR直接型结构时，实现网络结构如图8.6.1（a）所示。最后将每隔D–1个取一个作为输出即中有（D–1）/D的样值都被舍弃了。的N次乘法和N–1次加法在一个所以这种结构是一种低效实现结构，而且要求计算每一个T1时间内完成。图8.6.1（a）图8.6.1（b）8.6.1直接型FIR滤波器结构

为了得到相应的高效FIR直接实现，对图8.6.1（a）进行等效变换。

显然，将图8.6.1（a）中的↓D移在N条乘法器支路中的乘法器之前，如图8.6.1（b）所示，所得与原结构输出相同，即图8.6.1（a）与图8.6.1（b）是等效的。 图8.6.1（b）中各条支路里的↓D同时在时开通，此刻开始计算N个支路的N次乘法和最后的N–1次加法，得到一个输出样值：由于在到来之前所有的↓D同时关闭，直到时，N个↓D才又同时开通，再开始计算下一个输出序列样值。

改进效果：将乘法运算移到低采样率一侧，使乘法运算速度要求降低到原来的1/D，即原来要在一个T1时间内完成的运算，现在只要在DT1时间之内完成就可以了。当然，也使计算量减少到原来的1/D。

故称之为高效结构。328.6.1直接型FIR滤波器结构

应当说明，图8.6.1（b）中将↓D放在h(0)，h(T1)，…，h[(N–1)T1]之前，减少了运算量，但这并不是把抗混叠滤波放到了抽取之后，而是与原来的滤波作用等效。对此作如下解释：

滤波和抽取的作用次序在FIR实现结构中体现在滤波器输入端及延迟链上所加的信号序列，如果所加信号是抽取以前的信号，则是先滤波后抽取，反之是先抽取后滤波。图8.6.1（b）中，所有↓D均安排在延迟链之后，即滤波器延迟链上各点的信号仍然是原序列x(n1T1)，x[(n1–1)T1]，…，x[(n1–N+1)T1]，而不是抽取过的信号。每当↓D开通时，进入左侧的信号是未抽取的原信号，即输出的与图8.6.1（a）中抽选的输出信号相同，而两次开通之间所阻挡的信号恰好就是图8.6.1（a）中将来要舍弃的部分，所以计算结果是正确的。但绝对不能将↓D提前到延迟链之前，那样才是真正的先抽取后滤波器，会产生严重的混叠现象。338.6.1直接型FIR滤波器结构线性相位的FIR滤波器：对线性相位的FIR滤波器，实现时可用FIR线性相位结构，这样又可以使乘法计算量减少一半。根据线性相位时域特性：可画出抽取器FIR结构的线性相位形式如图8.6.2所示。34图8.6.28.6.1直接型FIR滤波器结构2．整数倍内插器的FIR直接实现35整数倍内插系统框图如图8.3.2所示。镜像滤波器采用FIR结构，I倍内插器的FIR直接实现结构如图8.6.3所示。图8.3.2图8.6.3

图8.6.3中乘法是在高采样率一侧进行的，不是高效结构，应设法将乘法运算移到低采样率一侧以减少计算量。但不能直接将移到滤波器后面。得出相应的高效结构。为此，先将FIR滤波网络部分进行转置，得到图8.6.4所示的FIR转置型结构，再用其代替图8.6.3中的FIR滤波网络，得到图8.6.5所示的内插系统直接实现。

↑I8.6.1直接型FIR滤波器结构

先将FIR滤波网络部分进行转置，得到图8.6.4所示的FIR转置型结构，再用其代替图8.6.3中的FIR滤波网络，得到图8.6.5所示的内插系统直接实现。图8.6.5中，先零值内插后分支相乘与先分支相乘后零值内插等效。因此，可将图8.6.5中的分别移到FIR网络的各支路的乘法器之后，可得到图8.6.6所示的内插系统直接实现高效结构。由于延时链上所加的仍然是内插后的信号，所以等效变换后的高效结构仍是先内插后滤波。36图8.6.4

图8.6.5↑I图8.6.68.6.1直接型FIR滤波器结构

当满足线性相位条件时，可用线性相位结构实，将乘法次数再减少一半。取N=9，画出内插器的线性相位FIR直接高效实现如图8.6.7所示。37图8.6.7按有理数因子I/D的采样率转换系统的高效FIR滤波器结构

为了叙述方便，先由图8.4.2画出按有理数因子I/D采样率转换系统的直接型FIR结构，如图8.6.8所示。38图8.6.8指导思想

基于内插系统的高效FIR滤波器结构与抽取系统的高效FIR滤波器结构进行设计。

（1）当I>D时，Fy

>Fx，将图8.6.8中的直接型FIR结构与前面的用图8.6.6所示的整数因子I内插器的高效FIR滤波器结构代替即可。↑I（2）当I<D时，Fy

<Fx，将图8.6.8中的直接型FIR结构与后面的用图8.6.1（b）所示的整数因子D抽取器的高效FIR滤波器结构代替即可。

↓D

8.6.2多相滤波器结构

1、内插系统的多相结构：可以证明，图8.6.6所示的按整数因子I内插系统的高效FIR滤波器结构可以用一组较短的多相滤波器组实现。如果FIR滤波器总长度为M=NI，则多相滤波器组由I个长度为N=M/I的短滤波器构成，且I个短滤波器轮流分时工作。

证明：

观察图8.6.3给出的整数因子I内插系统的直接型FIR滤波器结构。为了下面描述简单，定义一下符号：39用x(n)表示x(n1T1)用v(m)表示v(n2T2）用y(m)表示y(n2T2)8.6.2多相滤波器结构40所以，m=jI时刻：

输出序列为因为m=jI+1时刻，式中右移1位，N个x(n)的非零值与h(n)的对应关系也右移1位，所以，8.6.2多相滤波器结构

依此类推，当m=jI+I=(j+1)I时刻，

N个x(n)的值与h(n)的对应关系，只是x(n)又移进1位，所以，41综上所述，当时，有把上式中的看作长度N=M/I的子滤波器的单位脉冲响应，表示：

并用这样，从m=0开始，整数因子I内插系统的输出序列计算如下：

8.6.2多相滤波器结构式中，显然，当从0开始增大时，k从0开始以I为周期循环取值；j表示循环周期数。所以，实现上式y(m)的多相滤波器结构如下图42所示。图9.5.6特点：（1）I个子滤波器均运行于低采样率Fx下（2）系数少，计算量小。所以多相滤波器结构是一种高效结构。8.6.2多相滤波器结构工作过程：

输入端的x(n)每移入一个样值，I个子滤波器分别计算出y(m)的I个样值，选择电子开关以高采样率Fy=IFx，依次逆时针循环选取I个子滤波器的输出，形成输出序列。实现了整数因子I内插功能。43“多相滤波器”的解释：对低通滤波器h(n)按整数因子I抽取得到子滤波器h(n)是截止频率为的理想低通滤波器，所以的截止频率必然是，即I个子滤波器都是全通滤波器，幅度特性相同，它们的唯一区别是相位特性不同，故称为“多相滤波器”结构。

形成多相特性的机理：由h(n)的I个不同的起始点抽取得到I个子滤波器。8.6.2多相滤波器结构抽取系统的多相结构：44

思想：整数因子I内插器的实现结构与整数因子D=I抽取器的实现结构互为转置关系，将图9.5.6给出的整数因子I内插系统的多相滤波器结构进行转置，则得到图9.5.7所示的整数因子D抽取系统的多相滤波器结构。图9.5.7多相滤波器的单位脉冲响应：

式中，N为抗混叠FIR滤波器的总长度

M=DN，N=M/D

的长度。一般选择8.6.2多相滤波器结构抽取系统多相结构正确性的验证：45以N=D=2，M=DN=4为例，验证图9.5.7所示的抽取系统多相结构的正确性。首先根据图9.5.2（a）计算出抽取器的正确输出y(m):假设x(n)为因果信号，则8.6.2多相滤波器结构根据图9.5.7计算多相实现结构的输出y(m):

开始k=0,n=0，只有x(0)进入p0(n)，p1(n)中无信号，所以总输出y(0)=p0(0)x(0)=h(0)x(0)。

逆时针旋转开始下一周期:

k=D-1=1时，电子开关转到p1(n),x(1)进入p1(n)，p1(n)的输出为p1(0)x(1)=h(1)x(1)；

k=0时，电子开关又转到p0(n)，此时，x(2)进入p0(n)第一节，上一周期中进入p0(n)的x(0)移位到p0(n)的第二节，所以p0(n)的输出为46p0(0)x(2)+

p0(1)x(0)=h(0)x(2)+h(2)x(0)总的输出y(1)为p0(n)与p1(n)输出之和，即8.6.2多相滤波器结构同样道理，可求出下一旋转周期得到的输出所求y(0),y(1)和y(2)与式（9.5.10）相同，所以，图9.5.7所给结构是正确的。47例9.5.1设计一个按因子I=5的内插器，要求镜像滤波器通带最大衰减为0.1dB，阻带最小衰减为30dB，过渡带宽度不大于π/20。设计FIR滤波器系数h(n)，并求出多相滤波器实现结构中的5个多相滤波器系数。解:由式（9.3.5）知道FIR滤波器h(n)的阻带截止频率为π/5，根据题意可知滤波器其他指标参数:通带截止频率为π/5-π/20=3π/20，通带最大衰减为0.1dB，阻带最小衰减为30dB。调用remezord函数求得h(n)长度M=47，为了满足5的整数倍，取M=50。调用remez函数求得h(n)如下8.6.2多相滤波器结构h(0)=6.684246e-002=h(49) h(13)=-1.800562e-003=h(36)h(1)=-3.073256e-002=h(48) h(14)=-7.220485e-002=h(35)h(2)=-4.303671e-002=h(47) h(15)=-1.370181e-001=h(34)h(3)=-5.803096e-002=h(46) h(16)=-1.740193e-001=h(33)h(4)=-6.759203e-002=h(45) h(17)=-1.631924e-001=h(32)h(5)=-6.493009e-002=h(44) h(18)=-9.215300e-002=h(31)h(6)=-4.657608e-002=h(43) h(19)=4.004513e-002=h(30)h(7)=-1.386252e-002=h(42) h(20)=2.202029e-001=h(29)h(8)=2.674276e-002=h(41) h(21)=4.239994e-001=h(28)h(9)=6.463158e-002=h(40) h(22)=6.191918e-001=h(27)h(10)=8.776083e-002=h(39) h(23)=7.725483e-001=h(26)h(11)=8.607506e-002=h(38) h(24)=8.568808e-001=h(25)h(12)=5.500303e-002=h(37)488.6.2多相滤波器结构根据式（9.5.7）确定多相滤波器实现结构中的5个多相滤波器系数如下498.7采样率转换器的MATLAB实现MATLAB信号处理工具箱提供的采样率转换函数:upfirdn,interp,decimate,resample功能:Y=upfirdn(X,H,I,D)

先对输入信号向量X进行I倍零值内插，再用H提供的FIRDF对内插结果滤波，其中H为FIR数字滤波器的单位脉冲向量，FIR数字滤波器采用高效的多相实现结构。最后按因子D抽取得到输出信号向量Y。Y=interp(X,I)采用低通滤波插值法实现对序列向量X的I倍插值，其中的插值滤波器让原序列无失真通过，并在X的两个相邻样值之间按照最小均方误差准则插入I-1个序列值。得到的输出信号向量Y的长度为X长度的I倍。508.7采样率转换器的MATLAB实现Y=decimate(X,D,N)

先对序列X抗混叠滤波，再按整数因子D对序列X抽取。输出序列Y的长度是X长度的1/D。抗混叠滤波用N阶切比雪夫Ⅰ型低通滤波器，阻带截止频率为0.8Fs/(2D)，如果省略N，则默认用8阶切比雪夫Ⅰ型低通滤波；Y=decimate(X,D,N,'FIR')

用长度为N的FIR滤波器，FIR滤波器是抽取函数decimate自动调用fir1(N,1/D)设计的。省略N，则默认用30点FIR数字滤波器。其中，1/D为归一化-6dB截止频率，（存在问题？会有