我使用Matlab设计了一个非常简单的低通巴特沃斯滤波器。以下代码段演示了我所做的事情。

fs = 2.1e6;
flow = 44 * 1000;
fNorm =  flow / (fs / 2);
[b,a] = butter(10, fNorm, 'low');


在[b,a]中存储了滤波器系数。我想以整数形式获取[b,a],以便可以使用在线HDL代码生成器在Verilog中生成代码。

Matlab [b,a]的值似乎太小而无法与在线代码生成器一起使用(服务器端Perl脚本拒绝使用系数生成代码),我想知道是否有可能以可以用作适当输入的形式获取[b,a]。

我在Matlab中得到的a系数是:

1.0000
-9.1585
37.7780
-92.4225
148.5066
-163.7596
125.5009
-66.0030
22.7969
-4.6694
0.4307


我在Matlab中得到的b系数是:

1.0167e-012
1.0167e-011
4.5752e-011
1.2201e-010
2.1351e-010
2.5621e-010
2.1351e-010
1.2201e-010
4.5752e-011
1.0167e-011
1.0167e-012


我想使用在线生成器设计一个具有12位位宽和I或II滤波器形式的滤波器。我不知道上面链接中的“小数位”是什么意思。

使用以下命令运行代码生成器(http://www.spiral.net/hardware/filter.html)。上面列出的[b,a]系数,小数位设置为20,位宽为12,我收到以下运行错误:

Integer A constants: 1048576 -9603383 39613104 -96912015 155720456 -171714386 131597231 -69209161 23904282 -4896220 451621
Integer B constants: 0 0 0 0 0 0 0 0 0 0 0

Error: constants wider than 26 bits are not allowed, offending constant = -69209161, effective bitwidth = 7 mantissa + 20 fractional = 27 total.

An error has occurred - please revise the input parameters. 


如何更改我的设计,以便不会发生此错误?

更新:
使用Matlab生成6阶Butterworth滤波器,我得到以下系数:



1.0000
-5.4914
12.5848
-15.4051
10.6225
-3.9118
0.6010 
对于b:

0.0064e-005
0.0382e-005
0.0954e-005
0.1272e-005
0.0954e-005
0.0382e-005
0.0064e-005


运行在线代码生成器(http:// www。 spiral.net/hardware/filter.html),我现在收到以下错误(小数位数为8,位宽为20):

./iirGen.pl -A 256  '-1405' '3221' '-3943' '2719' '-1001' '153' -B  '0' '0' '0' '0' '0' '0' '0' -moduleName acm_filter -fractionalBits 8 -bitWidth 20 -inData inData  -inReg   -outReg  -outData outData -clk clk -reset reset -reset_edge negedge -filterForm 1  -debug  -outFile ../outputs/filter_1330617505.v 2>&1 
At least 1 non-zero-valued constant is required.  Please check the inputs and try again.


也许b-系数太小,或者也许代码生成器(http://www.spiral.net/hardware/filter.html)想要[b,a]采用另一种格式?

UPDATE:

也许我需要做的是sca根据小数位数对[b,a]系数进行取整,以获得整数形式的系数。

a .* 2^12
b .* 2^12


但是,我仍然认为b系数非常小。我在这里错了吗?

也许其他类型的过滤器(或过滤器设计方法)更合适?

更新:
正如Jason R和Christopher Felton在下面的评论中所建议的那样,SOS过滤器将更合适。我现在已经写了一些Matlab代码来获取SOS过滤器。

fs = 2.1e6;
flow = 44 * 1000;      
fNorm =  flow / (fs / 2);
[A,B,C,D] = butter(10, fNorm, 'low');
[sos,g] = ss2sos(A,B,C,D);


我得到的SOS矩阵是:

1.0000    3.4724    3.1253    1.0000   -1.7551    0.7705
1.0000    2.5057    1.9919    1.0000   -1.7751    0.7906
1.0000    1.6873    1.0267    1.0000   -1.8143    0.8301
1.0000    1.2550    0.5137    1.0000   -1.8712    0.8875
1.0000    1.0795    0.3046    1.0000   -1.9428    0.9598


是否仍可以使用Verilog代码生成工具(http://www.spiral.net/hardware/filter.html)来实现此SOS筛选器,还是我应该手工编写Verilog?

我想知道在这种情况下FIR滤波器是否会更好。

更多信息:递归IIR滤波器可以使用整数数学实现,具体方法如下:将系数表示为分数。 (有关更多详细信息,请参见Smith出色的DSP信号处理书:http://www.dspguide.com/ch19/5.htm)

以下Matlab程序使用Matlab将Butterworth滤波器系数转换为分数rat()函数。然后,如评论中所述,可以使用二阶部分来以数字方式实现过滤器(http://en.wikipedia.org/wiki/Digital_biquad_filter)。

% variables
% variables
fs = 2.1e6;                     % sampling frequency           
flow = 44 * 1000;               % lowpass filter


% pre-calculations
fNorm =  flow / (fs / 2);       % normalized freq for lowpass filter

% uncomment this to look at the coefficients in fvtool
% compute [b,a] coefficients
% [b,a] = butter(7, fNorm, 'low');
% fvtool(b,a)  

% compute SOS coefficients (7th order filter)
[z,p,k] = butter(7, fNorm, 'low');

% NOTE that we might have to scale things to make sure
% that everything works out well (see zp2sos help for 'up' and 'inf' options)
sos = zp2sos(z,p,k, 'up', 'inf'); 
[n,d] = rat(sos); 
sos_check = n ./ d;  % this should be the same as SOS matrix

% by here, n is the numerator and d is the denominator coefficients
% as an example, write the the coefficients into a C code header file
% for prototyping the implementation

 % write the numerator and denominator matices into a file
[rownum, colnum] = size(n);  % d should be the same
sections = rownum;           % the number of sections is the same as the number of rows
fid = fopen('IIR_coeff.h', 'w');

fprintf(fid, '#ifndef IIR_COEFF_H\n');
fprintf(fid, '#define IIR_COEFF_H\n\n\n');
for i = 1:rownum
   for j = 1:colnum

       if(j <= 3)  % b coefficients
            bn = ['b' num2str(j-1) num2str(i) 'n' ' = ' num2str(n(i,j))];
            bd = ['b' num2str(j-1) num2str(i) 'd' ' = ' num2str(d(i,j))];
            fprintf(fid, 'const int32_t %s;\n', bn);
            fprintf(fid, 'const int32_t %s;\n', bd);

       end
       if(j >= 5)  % a coefficients
            if(j == 5) 
                colstr = '1'; 
            end
            if(j == 6) 
                colstr = '2'; 
            end
            an = ['a' colstr num2str(i) 'n' ' = ' num2str(n(i,j))];
            ad = ['a' colstr num2str(i) 'd' ' = ' num2str(d(i,j))];
            fprintf(fid, 'const int32_t %s;\n', an);
            fprintf(fid, 'const int32_t %s;\n', ad);
       end
   end
end

% write the end of the file
fprintf(fid, '\n\n\n#endif');
fclose(fid);


评论

像这样的高阶IIR滤波器通常是使用二阶部分实现的。您可以通过级联多个第二阶(如果期望的阶为奇数,则使用一个第一阶)来获得所需的过滤器。与直接实现高阶滤波器相比,它通常是更健壮的实现。

如果您不执行@JasonR所建议的操作,则单词大小会非常大。当使用基本IIR结构实现时,此类过滤器可能会使单精度浮点运算失败,您需要SOS。

@JasonR:谢谢你的建议。我已通过上面的答案进行了更新。

@ChristopherFelton:谢谢您的帮助。

是的,您可以使用新的SOS矩阵从该站点创建3个过滤器。或者您可以在这里使用我的代码。它将与网站相同。我将很高兴将脚本更新为除SOS矩阵之外的内容。

#1 楼

如所讨论的,最好使用部分的总和,即将高阶滤波器分解为级联的二阶滤波器。更新后的问题具有SOS矩阵。使用此代码和此处的示例,Python对象可用于生成各个部分。

在matlab中

save SOS


在Python中

import shutil
import numpy
from scipy.io import loadmat
from siir import SIIR

matfile = loadmat('SOS.mat')  
SOS = matfile['SOS']
b = numpy.zeros((3,3))
a = numpy.zeros((3,3))
section = [None for ii in range(3)]
for ii in xrange(3):
    b[ii] = SOS[ii,0:3]
    a[ii] = SOS[ii,3:6]

    section[ii] = SIIR(b=b[ii], a=a[ii], W=(24,0))
    section[ii].Convert()  # Create the Verilog for the section
    shutil.copyfile('siir_hdl.v', 'iir_sos_section%d.v'%(ii))


有关定点的其他信息可在此处找到

评论


$ \ begingroup $
非常感谢您提供的所有有见地的链接以及Python代码;我希望您的答案(以及此处发布的其他答案)可以为其他许多人提供良好的参考。我希望我可以在这里将所有答复标记为已接受。
$ \ endgroup $
–尼古拉斯·基纳尔(Nicholas Kinar)
2012年3月2日14:38在

$ \ begingroup $
如果您有任何问题,请告诉我,如果它对您不起作用,我将更新/修复代码。我将对其进行修改(相对而言,很快),以直接接受SOS矩阵。
$ \ endgroup $
–克里斯托弗·费尔顿(Christopher Felton)
2012年3月2日15:52

$ \ begingroup $
我尝试从您的示例中实现自己的版本。在我的系统上,我必须使用“ from numpy import zeros”并将loatmat更改为loadmat()。 Matlab(mathworks.com/help/toolbox/signal/ref/ss2sos.html)给出的SOS矩阵是否与预期的格式相同?尝试访问SOS矩阵时收到以下错误:解释器到达“ b [ii] = SOS [0:3,ii]”行时,出现“ TypeError:unhashable type”
$ \ endgroup $
–尼古拉斯·基纳尔(Nicholas Kinar)
2012年3月2日在16:55

$ \ begingroup $
这取决于SOS.mat文件的格式。如果仅>>> print(matfile),它将显示已加载的.mat文件中的密钥。 scipy.io.loadmat始终作为字典(BOMK)加载。
$ \ endgroup $
–克里斯托弗·费尔顿(Christopher Felton)
2012年3月2日在18:05

$ \ begingroup $
是的,这是正确的,0的输出是1的输入,依此类推。需要对单词宽度稍加思考。默认值为两个使用24位小数(0整数,23小数,1个符号)。我相信您最初想使用较小的字宽。
$ \ endgroup $
–克里斯托弗·费尔顿(Christopher Felton)
2012年6月6日在2:33

#2 楼

“小数位”是总线上专用于表示数字小数部分的位数(例如,3.75中的.75)。

假设您拥有数字总线4位宽,1001代表什么数字?如果将其视为正整数(2 ^ 3 + 2 ^ 0 = 8 +1 = 9),则可能表示'9'。或者它可能表示二进制补码表示法为-7:(-2 ^ 3 + 2 ^ 0 = -8 + 1 = -7)。

带分数的数字,即“实数”呢?实数可以在硬件中表示为“定点”或“浮点”。看起来那些滤波器生成器使用定点。

回到我们的4位总线(1001)。让我们介绍一个二进制点,以便得到1.001。这意味着现在正在使用该点的RHS上的位构建整数,并使用LHS上的位构建分数。设置为1.001的数字总线表示的数字为1.125(1 * 2 ^ 0 + 0 * 2 ^ -1 + 0 * 2 ^ -2 + 1 * 2 ^ -3 = 1 + 0.125 = 1.125)。在这种情况下,在总线的4位中,我们使用3位来表示数字的小数部分。或者,我们有3个小数位。

所以,如果您有一个像上面一样的实数列表,则现在必须确定要表示多少个小数位。这是一个折衷方案:您使用的小数位数越多,则表示您想要的数字越接近,但是电路的体积就越大。而且,使用的小数位数越少,滤波器的实际频率响应与您最初设计的滤波器的频率响应就会越远!

更糟糕的是,您正在寻找建立无限脉冲响应(IIR)滤波器。如果您没有足够的小数和整数位,这些实际上可能变得不稳定!

评论


$ \ begingroup $
感谢您提供有见地的答案。我一直在尝试使用上面的小b系数运行代码生成器,但仍然收到一些错误。您能提出一些我可以做的正确运行发电机的建议吗?我将更新上面的答案以显示我已完成的工作。
$ \ endgroup $
–尼古拉斯·基纳尔(Nicholas Kinar)
2012年3月1日下午14:45

#3 楼

因此,马蒂(Marty)很好地解决了比特问题。在滤波器本身上,我认为您可能会从Matlab收到有关缩放系数不佳的警告或投诉?当我绘制滤波器时,从scipy到matlab都不是,但它可能非常相似。



通带下调了100 dB!因此,您可能需要确保使用较小的订单过滤器,无论如何这将有助于您的实施。当我到达6阶滤波器时,我就不再抱怨不良系数。也许尝试减少订单,看看它是否仍然符合您的要求。

评论


$ \ begingroup $
感谢您提出建议!我认为6阶滤波器同样适用。我认为使用matlab的fvtool,响应对于我的应用程序是很好的。我现在已经在上面更新了我的回复。但是,Verilog HDL代码生成器(spiral.net/hardware/filter.html)仍然出现问题。也许它想要[b,a]为另一种格式。另外,+ 1用于SciPy。
$ \ endgroup $
–尼古拉斯·基纳尔(Nicholas Kinar)
2012年3月1日于16:04