再议Uniform FFT modulated filterbank

最新推荐文章于 2024-07-04 09:58:48 发布

book_bbyuan

最新推荐文章于 2024-07-04 09:58:48 发布

阅读量688

点赞数 1

分类专栏：语音

本文链接：https://blog.csdn.net/book_bbyuan/article/details/83899119

版权

语音专栏收录该内容

12 篇文章 8 订阅

订阅专栏

前面写过一篇文章，但是感觉对于FFT 的filterbak还是没有说清楚，这边就再次把自己的一些理解说一下

1、这里面包含的几个概念

uniform，说明就是基础滤波器就一个，其他滤波器是它的一个变体，也就是调制（modulated）出来的，

polyphase分解，多相分解，主要目的是同一个多相组里面对应的“延迟“”是一样的，等会我后面会详细揭秘

2、FFT干嘛？只为了在analysis和 synthesis架设一个桥梁，实际上是analysis 就是简单滤波，通过滤波器

这每个channel出去的不是你理解的变成频域上的 X(0) X(1),...X(fft_len-1),它是每个子带的时域表达（直接滤波的结果），那为什么要求IDFT，是因为多相分解和调制的联合作用，“显得“而且“刚好” 可用一个IDFT来表达

3、fir滤波器的线性相位，一般设计的fir都是线性相位的，线性相位的 h(n),它的头尾是对称的，用科学的话就是偶对称或者奇对称，那为什么线性相位和对称有关系？具体可以看这个论文：

https://www.cnblogs.com/fall-li/p/4418646.html

4、我们举一个例子：

设计一个滤波器组，32个子带，它的原型滤波器截止到1/32 (对0~2*pi来说）

注意这个是全频带的，所谓的分32个子带是包括 fs/2~fs 的那一段的

滤波器调制公式

Hk = H0(z*exp(-j*2*pi*k/N)) k是第几个子带，从0开始，N是你一共有几个子带，（注意N不是说你滤波器tap数目）

或者说这么写 hk(i) = sigma i [exp(2*pi/M*k*i)*f(i)] , i=0~N-1 (因为H0（z）都是Z^-1的组合，z*exp(-j*2*pi*k/N代入后负数符号被干掉了，这里面的N又是滤波器长度

那么8个子带有了：它的顺序如下图：

实际上，Hk = H0(z*exp(j*2*pi*k/N)) 也是可以的，注意exp中没有负号了，它的子带滤波器的频带分配如下：

所以说上面那个图里面的分析部分的FFT和合成部分IFFT都是没有问题的，意思是你先分析用IDFT，合成用DFT或者先分析用DFT，后面再合成用IDFT都是可以的

原因上面说了

Hk = H0(z*exp(-j*2*pi*k/N)) 这是一种调制，对应IDFT

Hk = H0(z*exp(j*2*pi*k/N)) 这也是一种调制对应DFT

但是你要记住的是filterbank 的index是不一样的。他们对应的频带部分是不一样的。具体看上面的频带分配图

所以也不用奇怪纠结各种文献中的这个DFT和IDFT的先后次序问题了

5、另外说一个 FFT的 matrix变换方法：下面两个在matlab中是验证没有问题的，那么fft也就是可以用矩阵的方法来求，X是列向量

fft(x) = dftmtx(n)*x;
ifft(x) = conj(dftmtx(n))/n * x

6、举一个实际的例子：

input FIFO

filter已经分成 8组（或者8channel）多相滤波，每组 taps数目是3个，【不要认为这就是每个channel的滤波器系数，根本就不是，它只是假象而已】

他们对应相乘，每个ei都是多相part，另外如果所有的e求和，实际就是 x(n) 卷积 f(n),明眼人看出来就是滤波

% e7 = [ x16*f7 + x8*f15 + x0*f23
% e6=    x17*f6 + x9*f14 + x1*f22
% e5=    x18*f5 + x10*f13 + x2*f21
% e4=    x19*f4 + x11*f12 + x3*f20
% e3=    x20*f3 + x12*f11 + x4*f19
% e2 =   x21*f2 + x13*f10 + x5*f18
% e1 =   x22*f1 + x14*f9 + x6*f17
% e0 =   x23*f0 + x15*f8 + x7*f16 ]

容易求得 IDFT的变换矩阵是：

先算DFT的D（D是变换矩阵）

n = 8;
f = 2*pi/8;                 % Angular increment.
w = (0:f:2*pi-f/2).' * 1i; % Column.
x = 0:n-1;                  % Row.
D = exp(-w*x);              % Exponentiation of outer product.

IDFT的为 conj(dftmtx(n))/n = conj(D)/n

如果用 matrix 矩阵和多相矩阵相乘：也即（方阵）*（列阵列）

结果中间的每一行是：（每个滤波器或者每个通道的输出）

h0_part =sigma (1/8)*exp(2*pi/8*0*[0 1 2 3 4 5 6 7]).*[e7 e6 e5 e4 e3 e2 e1 e0];
h1_part =sigma (1/8)*exp(2*pi/8*1*[0 1 2 3 4 5 6 7]).*[e7 e6 e5 e4 e3 e2 e1 e0];

后面以此类推。。。

用h0_part不好说，因为都是0，在k=1的情况下，

h1_part = (1/8)* exp(2*pi/8*1*0) * (x23*f0 + x15*f8 + x7*f16)
   +(1/8)* exp(2*pi/8*1*1) * (x22*f1 + x14*f9 + x6*f17)
   +(1/8)* exp(2*pi/8*1*2) * (x21*f2 + x13*f10 + x5*f18)
   +(1/8)* exp(2*pi/8*1*3) * (x20*f3 + x12*f11 + x4*f19)
   +(1/8)* exp(2*pi/8*1*4) * (x19*f4 + x11*f12 + x3*f20)
   +(1/8)* exp(2*pi/8*1*5) * (x18*f5 + x10*f13 + x2*f21)
   +(1/8)* exp(2*pi/8*1*6) * (x17*f6 + x9*f14 + x1*f22)
   +(1/8)* exp(2*pi/8*1*7) * (x16*f7 + x8*f15 + x0*f23)

= (1/8)[ exp(2*pi/8*1*0) * (x23*f0) + exp(2*pi/8*1*(0+8)) * (x15*f8) +   exp(2*pi/8*1*(0+16))* x7*f16)
           exp(2*pi/8*1*1) * (x22*f1) + exp(2*pi/8*1*(1+8)) * (x14*f9) +   exp(2*pi/8*1*(1+16))* x6*f17)
           exp(2*pi/8*1*2) * (x21*f2) + exp(2*pi/8*1*(2+8)) * (x13*f10) + exp(2*pi/8*1*(2+16)*x5*f18)
+..

再把调制系数直接乘到滤波器系数上去（乘法结合律）.

= (1/8)[ (x23)*exp(2*pi/8*1*0) *(f0) + (x15)*exp(2*pi/8*1*(0+8)) *( f8) + (x7)*exp(2*pi/8*1*(0+16))*(f16)
(x22)*exp(2*pi/8*1*1) * (f1) + (x14)*exp(2*pi/8*1*(1+8)) *(f9) + (x6)*exp(2*pi/8*1*(1+16))*(f17)
(x21)*exp(2*pi/8*1*2) * f2) + (x13)*exp(2*pi/8*1*(2+8)) *(f10) + (x5)*exp(2*pi/8*1*(2+16)*(f18)
+..

是不是和一个调制滤波器相卷积？