Projektowanie filtra Butterwortha w Matlabie i uzyskiwanie współczynników filtra [ab] jako liczb całkowitych dla internetowego generatora kodów HDL Verilog

Za pomocą Matlaba zaprojektowałem bardzo prosty dolnoprzepustowy filtr Butterwortha. Poniższy fragment kodu pokazuje, co zrobiłem.

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

W [b, a] są przechowywane współczynniki filtra. Chciałbym uzyskać [b, a] jako liczby całkowite, dzięki czemu będę mógł używać internetowego generatora kodów HDL do generowania kodu w Verilog.

Wartości Matlab [b, a] wydają się zbyt małe, aby można je było stosować z generatorem kodu online (skrypt Perla po stronie serwera odmawia generowania kodu ze współczynnikami) i zastanawiam się, czy byłoby możliwe uzyskanie [b, a] w formie, która może być użyta jako właściwy wkład.

Współczynniki, które otrzymuję w Matlabie to:

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

Współczynniki b, które otrzymuję w Matlabie to:

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

Korzystając z generatora online, chciałbym zaprojektować filtr o 12-bitowej przepustowości i formie filtra I lub II. Nie wiem, co należy rozumieć przez „bity ułamkowe” w powyższym linku.

Po uruchomieniu generatora kodu (http://www.spiral.net/hardware/filter.html) ze współczynnikami [b, a] wymienionymi powyżej, z bitami ułamkowymi ustawionymi na 20 i przepustowością 12, pojawia się następujący błąd uruchamiania :

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. 

Jak mogę zmienić swój projekt, aby ten błąd nie wystąpił?

AKTUALIZACJA: Używając Matlaba do wygenerowania filtra Butterwortha 6. rzędu, otrzymuję następujące współczynniki:

Dla:

1.0000
-5.4914
12.5848
-15.4051
10.6225
-3.9118
0.6010 

dla b:

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

Po uruchomieniu generatora kodów online (http://www.spiral.net/hardware/filter.html) otrzymuję następujący błąd (z bitami ułamkowymi jako 8 i przepustowością 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.

Być może współczynniki b są zbyt małe, a może generator kodu (http://www.spiral.net/hardware/filter.html) chce [b, a] w innym formacie?

AKTUALIZACJA:

Być może to, co muszę zrobić, to skalować współczynniki [b, a] według liczby bitów ułamkowych, aby uzyskać współczynniki jako liczby całkowite.

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

Jednak nadal uważam, że współczynniki b są wyjątkowo małe. Co robię tutaj źle?

Być może inny rodzaj filtra (lub metoda projektowania filtrów) byłaby bardziej odpowiednia? Czy ktoś może coś zasugerować?

AKTUALIZACJA: Jak sugerują Jason R i Christopher Felton w komentarzach poniżej, bardziej odpowiedni byłby filtr SOS. Napisałem teraz kod Matlaba, aby uzyskać filtr 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);

Matryca SOS, którą otrzymuję, to:

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

Czy nadal można użyć narzędzia do generowania kodu Verilog (http://www.spiral.net/hardware/filter.html), aby zaimplementować ten filtr SOS, czy powinienem po prostu napisać Verilog ręcznie? Czy dostępne jest dobre referencje?

Zastanawiam się, czy lepiej byłoby użyć filtra FIR w takiej sytuacji.

WIĘCEJ: Rekurencyjne filtry IIR można zaimplementować za pomocą matematyki liczb całkowitych, wyrażając współczynniki jako ułamki. (Więcej informacji można znaleźć w doskonałej książce przetwarzania sygnałów DSP Smitha: http://www.dspguide.com/ch19/5.htm )

Poniższy program Matlab konwertuje współczynniki filtra Butterwortha na części ułamkowe za pomocą funkcji rat () Matlaba. Następnie, jak wspomniano w komentarzach, sekcje drugiego rzędu można wykorzystać do numerycznej implementacji filtra (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);

Jak omówiono, najlepiej jest użyć sumy sekcji, czyli rozbić filtr wyższego rzędu na kaskadowe filtry drugiego rzędu. Zaktualizowane pytanie ma matrycę SOS. Korzystając z tego kodu i przykładu tutaj, obiekt Python może zostać użyty do wygenerowania poszczególnych sekcji.

W Matlabie

save SOS

W Pythonie

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))

Dodatkowe informacje na temat punktu stałego można znaleźć tutaj

„Bity ułamkowe” to liczba bitów w magistrali, które zostały przeznaczone do reprezentowania ułamkowej części liczby (np. 0,75 w 3,75).

Powiedzmy, że masz cyfrową magistralę o szerokości 4 bitów, jaką liczbę 1001reprezentuje? Może to oznaczać „9”, jeśli potraktujesz go jako dodatnią liczbę całkowitą (2 ^ 3 + 2 ^ 0 = 8 + 1 = 9). Lub może oznaczać -7 w notacji dopełniacza dwóch: (-2 ^ 3 + 2 ^ 0 = -8 + 1 = -7).

Co z liczbami z niektórymi ułamkami, tj. Liczbami „rzeczywistymi”? Liczby rzeczywiste mogą być reprezentowane sprzętowo jako „stały punkt” lub „zmiennoprzecinkowy”. Wygląda na to, że te generatory filtrów używają punktu stałego.

Powrót do naszej 4-bitowej magistrali ( 1001). Wprowadźmy punkt binarny, abyśmy otrzymali 1.001. Oznacza to, że używali teraz bitów na RHS punktu do budowania liczb całkowitych i bitów na LHS do budowania ułamka. Liczba reprezentowana przez cyfrową magistralę 1.001to 1,125 ( 1* 2 ^ 0 + 0* 2 ^ -1 + 0* 2 ^ -2 + 1* 2 ^ -3 = 1 + 0,125 = 1,125). W tym przypadku spośród 4 bitów w magistrali używamy 3 z nich do reprezentowania ułamkowej części liczby. Lub mamy 3 bity ułamkowe.

Więc jeśli masz listę liczb rzeczywistych, tak jak wcześniej, musisz teraz zdecydować, ile bitów ułamkowych chcesz je reprezentować. I tutaj jest kompromis: im więcej bitów ułamkowych użyjesz, tym bliżej możesz przedstawić żądaną liczbę, ale im większy będzie twój obwód. Co więcej, im mniej bitów ułamkowych używasz, tym bardziej rzeczywista odpowiedź częstotliwościowa filtra będzie odbiegać od tej, którą zaprojektowałeś na początku!

Co gorsza, chcesz zbudować filtr Infinite Impulse Response (IIR). Mogą one faktycznie stać się niestabilne, jeśli nie masz wystarczającej liczby bitów ułamkowych i całkowitych!

Więc Marty dobrze zadbał o to pytanie. Jeśli chodzi o sam filtr, myślę, że prawdopodobnie otrzymujesz ostrzeżenie lub skargę od Matlaba na temat źle skalowanych współczynników? Kiedy kreślę filtr, od scipy nie matlab, ale prawdopodobnie jest bardzo podobny.

To jest 100 dB w dół w paśmie pasma! Możesz więc chcieć upewnić się, że potrzebujesz mniejszego filtra zamówień, który i tak pomoże ci we wdrożeniu. Kiedy dochodzę do filtra 6. rzędu, przestaję narzekać na złe współczynniki. Może spróbuj zmniejszyć zamówienie i sprawdź, czy nadal spełnia twoje wymagania.