CN1688104B - 数字信号处理方法及装置 - Google Patents

Abstract

本发明的信号处理装置包括：输入数据存储逻辑单元；FFT/FFT运算逻辑单元；FFT系数存储逻辑单元；乘加运算器；重叠缓冲器和参数与地址控制器。本发明的信号处理方法包括步骤：循环计数器的计数J置0，重叠缓冲器清零；选定FFT计算的平均长度N和运算长度FL；判断循环计数器的计数J和各滤波器冲击响应H_i值；将各个滤波器冲击响应H_i的FFT系数和输入信号的FFT系数点乘以及相加点乘；矢量相加M个点乘的结果；对矢量相加的结果进行IFFT运算并输出部分值；回到上述判断循环计数器的计数J和各滤波器冲击响应H_i值的步骤，然后重复进行以后的各个步骤。按照本发明，IFFT运算、时域加法和重叠缓冲器的操作减少了，因此，计算复杂度明显降低。

Description

数字信号处理方法及装置

技术领域

本发明涉及数字信号处理技术，特别是涉及利用滤波器的数字信号处理技术。

背景技术

在现有技术中，很多使用因果FIR(有限冲击响应)滤波器的数字信号处理都会使用如下形式的卷积求和运算，得到滤波结果Y：

$Y=\underset{i=0}{\overset{M}{\mathrm{\Σ}}}{X}_i*H_i$ (1)

Xi为第i路输入信号，H_i表示作用于第i路信号的FIR滤波器的冲击响应，表示为：

$H_{i,k}=\left\{\begin{array}{cc}0,& k<0\\ h_{i,k},& 0\&le;k<L_i\\ 0,& k\&GreaterEqual;L_i\end{array}\right\}---\left(2\right)$

即多个输入信号X₁，X₂，...，X_M分别经过多个不同的滤波器滤波之后相加并输出，式中H_i，k表示H_i在时间点k上的值。

图1是卷积求和运算的示意图。X₁，X₂，...，X_M是多个输入信号，每个信号X₁，X₂，...，X_M分别通过相应的滤波器滤波，输出已滤波的信号，各个滤波器的输出由加法器进行相加，得到一个滤波结果Y₁，完成一次卷积求和运算。然后，进行计算第二个滤波结果Y₂，依此类推，通过M次计算，可以得到M个滤波结果Y₁-Y_M。利用方程(1)计算信号的处理结果。

假设各滤波器的冲响应分别为H₁，H₂，...，H_M，阶数分别为L1，...Li，LM，平均长度为L，则针对这种形式的运算，一般的计算方法有两种：

1.根据方程(1)直接进行时域计算。这样计算单个滤波结果Y(n)所需要的计算复杂度为O(L×M)次乘加运算。

2.使用快速卷积算法Y_i＝X_i*H_i得到每个滤波器的滤波结果Y_i，然后再根据以下公式计算和，

$Y=\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{Y}_i.$

图2表示一种快速卷积算法的流程图。N点的输入信号X_i进行补零，形成P点的输入信号，P≥N+Li-1，P一般为2的整数次方。该输入信号进行FFT变换，输出P点的DET系数Xi，系数Xi和滤波器H_i的P点DFT系数进行乘法运算，输出P点IFFT结果，前L_i-1点的IFFT结果和LAP缓冲器(重叠缓冲器)中存储的值相加，把第N+1至N+L_i-1点的IFFT结果存储到LAP缓冲器，最后输出P点IFFT结果的前N个点IFFT结果。

假设平均滤波器长度L＜＜B。P≈B，则P点快速傅立叶变换(FFT)的运算量为O(Blog(B))。B为FFT计算的平均长度，P为快速卷积运算中平均的输入块长。假设滤波器H_i在整个滤波过程中的值不变，其离散傅立叶变换(DFT)系数可以保留在缓冲器中，那么整个快速卷积过程需要1个FFT和1个快速傅立叶反变换(IFFT)，DFT相乘需要2B次乘加运算，这样每次滤波器卷积运算的复杂度为O(Blog(B))。计算B个Y值输出需要M次滤波计算，复杂度为O(M(Nlog(N)))，计算每个Y值所需的复杂度为MO(log(B))。

在上述快速卷积算法中，计算最终输出Y值时，总共用了M次IFFT运算得到每个滤波器的输出值Y_i，然后把M个Y_i，进行时域相加。考虑到FFT运算为线性的，FFT运算和IFFT运算可以表示为：

$\mathrm{FFT}\left(\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{S}_i\right)=\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}\mathrm{FFT}\left(S_i\right)$

(3)

$\mathrm{IFFT}\left(\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{\mathrm{FS}}_i\right)=\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}\mathrm{IFFT}\left({\mathrm{FS}}_i\right)$

公式(3)中S_i为第i路时域信号，FS_i为第i路信号的FFT系数。

对于M路的信号输入，一般的算法每次进行一个块的计算需要M次的IFFT计算，另外还需要M次的重叠相加操作，计算复杂度较高。

发明内容

本发明的目的在于克服现有技术存在的问题，提供一种数字信号处理方法和装置，使得信号处理的运算复杂度明显降低。

本发明提供一种数字信号处理装置，该装置包括：

输入数据存储逻辑单元，用于存储输入的各路信号数据与滤波器数据，根据参数与地址控制端送来的补零个数对输入数据进行补零操作；

FFT/IFFT运算逻辑单元，根据输入数据存储逻辑单元提供的数据和参数与地址控制器提供的命令对输入数据进行FFT运算，然后把运算结果送到FFT系数存储逻辑单元；

FFT系数存储逻辑单元，在该参数与地址控制器的控制下把FFT系数存储在内存中；

乘加运算器，从FFT系数存储逻辑单元中取出各个输入数据的FFT系数以及各个滤波器的FFT系数，把相应的滤波器FFT系数和输入数据的FFT系数点乘，把各个点乘的结果矢量相加，最后的结果送回FFT系数存储逻辑单元；

所述FFT/IFFT运算逻辑单元，根据该参数与地址控制器发来的命令，对FFT系数存储逻辑单元保留的乘加运算结果进行IFFT运算，IFFT运算的结果由FFT运算存储逻辑进行保存；

重叠缓冲器，以每次IFFT运算结果的部分值更新其保存的值；

参数与地址控制器，控制输入数据存储逻辑单元、FFT/FFT运算逻辑单元和乘加运算器的操作。

所述乘加运算器把IFFT运算的部分结果和重叠缓冲器中保留的上一次运算的结果相加，并且保存回原来的位置。

本发明还提供一种数字信号处理方法，包括步骤：

循环计数器的计数J置0，重叠缓冲器清零；

选定FFT计算的平均长度N和运算长度FL；

判断循环计数器的计数J和各滤波器冲击响应H_i值，如果判断循环计数器的计数J＝0或者滤波器冲击响应H_i的值发生变化，把滤波器冲击响应H_i往后补零成为长度为FL的向量，并且对各个滤波器冲击响应H_i进行FFT运算，得到其FFT系数[h_i，0，h_i，1，...h_i，FL-1]并保留在滤波器冲击响应的缓冲器中；如果判断J≠0或者滤波器冲击响应H_i的值没有发生变化，则进入下一步骤；

对各路输入信号X_i按时间输入N个值，往后补零成为运算长度FL，并且对各路输入信号X_i进行FL点的FFT运算，得到其FFT系数[x_i，0，x_i，1，...x_i，FL-1]，i＝1，2，...M；

将各个滤波器冲击响应H_i的FFT系数[h_i，0，h_i，1，...h_i，FL-1]和各路输入信号X的FFT系数[x_i，0，x_i，1，...x_i，FL-1]点乘；

由乘加运算器把M个点乘的结果[h_i，0x_i，0，h_i，1x_i，1，...h_i，FL-1x_i，FL-1]进行矢量相加；

对矢量相加的结果进行IFFT运算得到时域的FL点输出T_j；

该输出T_j中的第0点到第L_max-2点的值加上重叠缓冲器中第0点到第L_max-2点的值，L_max是滤波器的最大阶数；

将该输出T_j中的第0到第N-1点的值输出，并将该输出T_j中的第N点到第N+L_max-2点的值分别存储在重叠缓冲器中的0点到第L_max-2点的位置；

回到上述判断循环计数器的计数J和各滤波器冲击响应H_i值的步骤，然后重复进行以后的各个步骤，直到完成对所有的输入信号的处理。

本发明的算法和以往的快速卷积算法相比，IFFT运算减少了M-1次，减少了N-1次的时域加法，减少了M-1次的重叠缓冲器操作，只增加了(M-1)*FL次的复数加法，因此，计算复杂度明显降低。

附图说明

图1是卷积求和运算的示意图；

图2是快速卷积算法的流程图；

图3是根据本发明的信号处理方法的流程图；

图4是根据本发明的信号处理装置的结构方框图。

具体实施方式

考虑在频域进行M个DFT系数的相加得到Y对应的DFT系数，然后只进行一次IFFT运算得到时域Y信号。这样就要求所有单个滤波操作(DFT系数相乘)后的DTF系数长度必须相等，即所有快速卷积中的FFT运算长度必须相等。

每个滤波器的输入都可以表示为矩形加窗之后的信号之和，

$R\left(i\right)=\left\{\begin{array}{cc}1,& 0\&le;i<N\\ 0,& i<0\\ 0,& i\&GreaterEqual;N\end{array}\right\}$

$X_i=\underset{j=-\&infin;}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{X}_{i}R(i-\mathrm{jN})---\left(4\right)$

式中i，j是整数序号，N是矩形窗长。这样，由于卷积运算为线性运算，输出Y可以表示为：

$Y=\underset{j=-\&infin;}{\overset{\&infin;}{\mathrm{\&Sigma;}}}\underset{i=0}{\overset{M}{\mathrm{\&Sigma;}}}\left[X_{i}R(i-\mathrm{jN})\right]*H_i---\left(5\right)$

设定

$T_j=\underset{i=0}{\overset{M}{\mathrm{\&Sigma;}}}\left[X_{i}R(i-\mathrm{jN})\right]*H_i$

那么，方程式(5)可简化为

$Y=\underset{j=-\&infin;}{\overset{\&infin;}{\mathrm{\&Sigma;}}}{T}_j---\left(6\right)$

我们只需要计算各个时段的输出T_j，然后在时域上叠加各个T_j，就可以得到最终的输出Y。考虑T_j的计算，由于H_i为因果FIR滤波器的冲击响应，而各输入X_iR(i-jN)在时间jN之前全部为0，因此

T_jk＝0，k＜jN    (7)

假设滤波器H_i的最大长度为L_max，那么

T_jk＝0，k≥jN+N+L_max-1    (8)

因此T_j在区间jN＜k＜jN+N+L_max-1内的输出值可能不为0，区间的长度为N+L_max-1。各T_j之间重叠的部分为L_max-1。

考虑T_j的计算

$T_j=\underset{i=0}{\overset{M}{\mathrm{\&Sigma;}}}\left[X_{i}R(i-\mathrm{jN})\right]*H_i$

1)[X_iR(i-jN)]*Hⁱ的计算可以使用快速卷积算法，其快速卷积中FFT的长度必须大于等于N+L_i-1。

2)考虑到各个快速卷积运算中FFT运算长度必须一致，所选的长度必须大于等于N+L_max-1。而且为了FFT运算的方便，选择一个大于等于N+L_max-1的2的整数次方的一个数作为FFT的运算长度FL，即

FL≥N+L_max-1

FL＝2^k，k为整数。    (9)

3)为了节省计算量，把计算T_j时的时域加法转换为频域的加法；

4)因为T_j与T_j+1之间有L_max-1的重叠，必须把T_j的FL点输出(其中第N+L_max-1到第FL-1点数值为0)中的第N点到第N+L_max-2点保存在缓冲器中，以便和T_j+1的值相加得到正确输出。为了避免T_j和T_j+2之间的重叠，可以设定FFT计算的平均长度

N＞L_max-2    (10)

下面结合图3说明信号处理的步骤。

假定有M个输入信号X₁至X_M。

1.置J＝0，重叠缓冲器清零，J是循环计数器的计数；

2.根据表示式(10)FFT计算的平均长度N＞L_max-2，根据运用中实时性的要求、缓冲器大小等情况确定FFT计算的平均长度B的大小；

3.根据表示式(9)，运算长度FL≥N+Lmax-1和运算长度FL必须为2的整数次方确定FFT的运算长度FL。

4.对各滤波器冲击响应H_i，如果判断J＝0或者滤波器冲击响应H_i的值发生变化，把滤波器冲击响应H_i往后补零成为长度为FL的向量，并且对各个滤波器冲击响应H_i进行FFT运算，得到其FFT系数[h_i，0，h_i，1，...h_i，FL-1]并保留在滤波器冲击响应的缓冲器中；如果判断J≠0或者滤波器冲击响应H_i的值没有发生变化，则进入下一步骤。

5.对各个输入信号X_i按时间输入N个值，往后补零成为长度为FL的向量，并且对X_i进行FL点的FFT运算，得到其FFT系数[x_i，0，x_i，1，...x_i，FL-1]。

6.对各个滤波器冲击响应H_i的FFT系数[h_i，0，h_i，1，...h_i，FL-1]和各个输入信号X_i的FFT系数[x_i，0，x_i，1，...x_i，FL-1]点乘，得到点乘结果[h_i，0x_i，0，h_i，1x_i，1，...h_i，FL-1x_i，FL-1]。此运算需要FL次复数乘法，但假如滤波器冲击响应H_i和输入信号X_i都为实数，根据DFT系数的对称性，则只需要进行FL/2+1次复数乘法运算。

7.把M个点乘结果[h_i，0x_i，0，h_i，1x_i，1，...h_i，FL-1x_i，FL-1]进行矢量相加，得到

$p=[\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{h}_{i,0}{x}_{i,0},...\underset{i=1}{\overset{M}{\mathrm{\&Sigma;}}}{h}_{i,\mathrm{FL}-1}{x}_{i,\mathrm{FL}-1}].$

此运算需要(M-1)*FL次复数加法，但是对于滤波器冲击响应H_i和输入信号X_i都为实数的情况，需要进行(M-1)*(FL/2+1)次复数加法。

8.对p进行IFFT运算得到时域的FL点输出T_j，T_j的第B+L_max-1到第FL-1点的值为零。

9.T_j的第0点到第L_max-2点的值加上重叠缓冲器中第0点到第L_max-2点的值。

10.T_j的第N点到第N+L_max-2点分别存储在重叠缓冲器中的0点到第L_max-2点。

11.输出T_j的第0到第N-1点。

12.循环计数器的计数J递增1，即J＝J+1，回到步骤4，重复以后的各个步骤，直到对M个输入信号的上述处理完成。

下面以具体的例子进行说明信号处理的过程。

假设共有3路输入信号X1，X2，X3，H1为34阶FIR滤波器，H2为99阶FIR滤波器，H3为59阶FIR滤波器，即L_max＝99，，计算滤波结果Y＝X1*H1+X2*H2+X3*H3，在计算过程种，H1，H2，H3的数值恒定。

首先，根据(10)式和(9)式选定FFT计算的平均长度N＝128，运算长度FL＝256，N和FL的数值由控制信号输入端进入参数与地址控制器确定，重叠缓冲器的大小不小于L_max-1＝98，初始值为0。数据信号分别由输入端输入h1的滤波器系数，输入数据存储逻辑根据控制器送来的补零数分别把系数补零成256位，FFT/IFFT运算逻辑单元从数据存储逻辑取数后进行FFT运算，由FFT系数存储逻辑单元存储其FFT系数；分别对h2，h3进行类似的运算，h1，h2，h3的FFT系数分别记为FH1，FH2，FH3。输入128位信号X1的值，补零成256位之后，由FFT/IFFT运算逻辑单元进行FFT运算之后，由FFT系数存储逻辑单元存放其FFT系数值。对信号X2，X3进行相应的操作，信号X1，X2，X3的FFT系数分别记为FX1，FX2，FX3。由乘加运算器计算FFT系数的乘加结果FS＝FX1*FH2+FX2*FH2+FX3*FH3，由FFT系数存储逻辑单元存储运算结果FS(256位复系数)。FFT/IFFT运算逻辑对FS进行IFFT运算，由FFT系数存储逻辑对结果S进行保存。通过乘加运算器把结果S的第0个值到第97个值与重叠缓冲器中的第0个值到第97个值相加，结果保存在结果S的相应位置，即S[0]＝S[0]+重叠缓冲器[0]，...，S[97]＝S[97]+重叠缓冲器[97]。S[0]，...，S[127]位计算结果进行输出。把S[128]，...，S[128+(98-1)]保存到重叠缓冲器中的重叠缓冲器[0]到重叠缓冲器[97]的位置。然后回到对各个滤波器冲击响应H_i和各路输入信号X_i补零以及FFT运算的步骤，并重复进行以后的各个步骤，直到信号X1，X2，X3的输入结束为止。

图4是根据本发明的信号处理装置的结构方框图。

输入数据存储逻辑单元存储输入的各路信号数据与滤波器数据，根据参数与地址控制端送来的补零个数对输入数据进行补零操作。FFT/IFFT运算逻辑单元，根据输入数据存储逻辑单元提供的数据和参数与地址控制器提供的命令对输入数据进行FFT运算，然后把运算结果送到FFT系数存储逻辑单元。FFT存储逻辑单元在参数与地址控制器的控制下把FFT系数存储到内存的某个位置。当各个输入数据的FFT运算及存储完成时，由乘加运算器从FFT系数存储逻辑单元中取出各个输入数据的FFT系数以及各个滤波器的FFT系数，把相应的滤波器FFT系数和输入数据的FFT系数点乘，把各个点乘的结果矢量相加，最后的结果送回FFT系数存储逻辑单元。此后，FFT/IFFT运算逻辑单元根据参数与地址控制器的命令对FFT系数存储逻辑单元保留的乘加运算结果进行IFFT运算，IFFT运算的结果由FFT运算存储逻辑单元进行保存。最后乘加运算器把IFFT运算的部分结果和重叠缓冲器中保留的上一次运算的结果相加，并且保存在回原来的位置。重叠缓冲器的值更新为IFFT运算结果的部分值。输出的运算结果为IFFT运算结果的一部分值。这里没有对硬件执行此种计算的各个步骤进行详细说明，因为硬件的步骤和前面已经详述的算法步骤是一一对应的。

本发明以具体的实施例结合附图进行了说明，但是具体的实施例不是用来限定本发明，本领域的技术人员可以在不脱离本发明的精神和范围的情况下进行修改、变更或替换。

Claims (5)

1.一种信号处理方法，包括以下步骤：

循环计数器的计数J置0，重叠缓冲器清零；

选定FFT计算的平均长度N和运算长度FL；

判断循环计数器的计数J和各滤波器冲击响应H_i值，如果判断循环计数器的计数J＝0或者滤波器冲击响应H_i的值发生变化，把滤波器冲击响应H_i往后补零成为长度为FL的向量，并且对各个滤波器冲击响应H_i进行FFT运算，得到其FFT系数[h_i，0，h_i，1，...h_i，FL-1]并保留在滤波器冲击响应的缓冲器中；如果判断J≠0或者滤波器冲击响应H_i的值没有发生变化，则进入下一步骤；

对各路输入信号X_i按时间输入N个值，往后补零成为运算长度FL，并且对各路输入信号X_i进行FL点的FFT运算，得到其FFT系数[x_i，0，x_i，1，...x_i，FL-1]，i＝1，2，...M；

将各个滤波器冲击响应H_i的FFT系数[h_i，0，h_i，1，...h_i，FL-1]和各路输入信号X_i的FFT系数[x_i，0，x_i，1，...x_i，FL-1]点乘；

由乘加运算器把M个点乘的结果[h_i，0x_i，0，h_i，1x_i，1，...h_i，FL-1x_i，FL-1]进行矢量相加；

对矢量相加的结果进行IFFT运算得到时域的FL点输出T_j；其中j为整数序号；

该输出T_j中的第0点到第L_max-2点的值加上重叠缓冲器中第0点到第L_max-2点的值，L_max是滤波器的最大阶数；

将该输出T_j中的第0到第N-1点的值输出，并将该输出T_j中的第N点到第N+L_max-2点的值分别存储在重叠缓冲器中的0点到第L_max-2点的位置；

回到上述判断循环计数器的计数J和各滤波器冲击响应H_i值的步骤，然后重复进行以后的各个步骤，直到完成对所有的输入信号的处理；

其中所述选定FFT计算的平均长度是根据下式确定的：

N＞L_max-2

和运算长度根据以下公式确定：

FL≥N+L_max-1，FL＝2^k，k为整数，L_max是滤波器的最大阶数。

2.根据权利要求1的信号处理方法，其中所述各路输入信号X_i的系数与滤波器冲击响应H_i的FFT系数的点乘步骤进行FL次；其中所述矢量相加运算步骤需要进行(M-1)*FL次。

3.根据权利要求2的信号处理方法，其中在滤波器冲击响应H_i的FFT系数和所述各路输入信号X_i的FFT系数都是实数的情况下，需要进行FL/2+1次复数乘法运算和(M-1)*(FL/2+1)次复数加法运算。

4.一种信号处理装置，用于权利要求1所述的信号处理方法，包括：

输入数据存储逻辑单元，用于存储输入的各路信号数据与滤波器数据，根据参数与地址控制端送来的补零个数对输入数据进行补零操作；

FFT/IFFT运算逻辑单元，根据输入数据存储逻辑单元提供的数据和参数与地址控制器提供的命令对输入数据进行FFT运算，然后把运算结果送到FFT系数存储逻辑单元；

FFT系数存储逻辑单元，在该参数与地址控制器的控制下把FFT系数存储在内存中；

乘加运算器，从FFT系数存储逻辑单元中取出各个输入数据的FFT系数以及各个滤波器的FFT系数，把相应的滤波器冲击响应的FFT系数和输入数据的FFT系数点乘，把各个点乘的结果矢量相加，最后的结果送回FFT系数存储逻辑单元；

所述FFT/IFFT运算逻辑单元，根据该参数与地址控制器发来的命令，对FFT系数存储逻辑单元保留的乘加运算结果进行IFFT运算，IFFT运算的结果由FFT系数存储逻辑单元进行保存；

重叠缓冲器，保存T_j的FL点输出中的第N点到第N+L_max-2点，以每次IFFT运算结果的部分值更新其保存的值，所述部分值为所述运算结果中的第0个值到第L_max-2个值；

参数与地址控制器，控制输入数据存储逻辑单元、FFT/IFFT运算逻辑单元和乘加运算器的操作。

5.根据权利要求4的信号处理装置，其特征在于所述乘加运算器把IFFT运算结果的所述部分值和重叠缓冲器中保留的上一次运算的结果相加，并且保存回所述部分值在所述FFT系数存储逻辑单元中的位置。

Images