MyHDL и Scipy FIR

Еще один пример использования MyHDL в Python - реализация фильтра нижних частот. Пример кода:

from random import randint
import numpy as np
from myhdl import *
import matplotlib.pyplot as plt

class Cook_book():
    def PMC_sum_arg(self, i_arg1, i_arg2):
        o_len = 1 + max(len(i_arg1),len(i_arg2))
        o_len_signed = o_len-1
        return Signal(intbv(0,-2**o_len_signed, 2**o_len_signed))
    def PMC_mult_arg(self, i_arg1, i_arg2):
        o_len = len(i_arg1) + len(i_arg2)
        o_len_signed = o_len-1
        return Signal(intbv(0,-2**o_len_signed, 2**o_len_signed))
        
    def PMC_sum(self, i_arg1, i_arg2, o_out):
        w_out = self.PMC_sum_arg(i_arg1, i_arg2)
        @always_comb
        def exits():
            o_out.next = w_out
        @always_comb
        def comb():
            w_out.next = i_arg1 + i_arg2
        return exits, comb
####w_out
    def PMC_mult(self, i_arg1, i_arg2, o_out):
        w_out = self.PMC_mult_arg(i_arg1, i_arg2)
        @always_comb
        def exits():
            o_out.next = w_out
        @always_comb
        def comb():
            w_out.next = i_arg1*i_arg2
        return exits, comb
    def arg_ROM(self, i_adr, CONTENT):
        len_content = len(CONTENT)
        len_adr = 2**(len(i_adr))
        max_width = 0
        if len_content > len_adr:
            print "arg_ROM WARNING len_content > len_adr"
        else:
            for i in range(len(CONTENT)):
                if CONTENT[i].bit_length() > max_width:
                    max_width = CONTENT[i].bit_length()
        return Signal(intbv(0, -2**max_width, 2**max_width))
        
    def ROM(self, i_clk, i_adr, o_out, CONTENT):
        w_out = self.arg_ROM(i_adr, CONTENT)
        @always_comb
        def exits():
            o_out.next = w_out
        @always_comb
        def comb():
            w_out.next = CONTENT[int(i_adr)]
        return exits, comb
    def shift_arg(self, WIDTH, LEN):
        w_signed = WIDTH-1
        return [Signal(intbv(0,-2**w_signed, 2**w_signed)) for i in range(LEN)]

    def shift_reg(self, i_clk, i_rst, i_en, i_din, o_dout, WIDTH, LEN):
        w_signed = WIDTH-1
        #o_dout = self.shift_arg(WIDTH, LEN)
        r_shift = [Signal(intbv(0,-2**w_signed, 2**w_signed)) for i in range(LEN)]
        @always_comb
        def exits():
            for i in range(LEN):
                o_dout[i].next = r_shift[i]
        @always(i_clk.posedge, i_rst.negedge)
        def seq():
            if i_rst == 0:
                for i in range(LEN):
                    r_shift[i].next = 0
            elif i_en == 1:
                r_shift[0].next = i_din
                for i in range(1, LEN):
                    r_shift[i].next = r_shift[i-1]
        return exits, seq
    def serial_filter(self, i_clk, i_rst, i_newdata, i_din, o_dout, COEFFS, WIDTH, LEN):
        
        o_out_shift = self.shift_arg(16, len(COEFFS))
        uut_shift = self.shift_reg(i_clk, i_rst, i_newdata, i_din, o_out_shift, 16, LEN)
        
        MULT_MSB = 16+16
        MULT_WIDTH_SIGNED = 16+16-1
        w_mult = [Signal(intbv(0,-2**MULT_WIDTH_SIGNED, 2**MULT_WIDTH_SIGNED)) for i in range(LEN)]
        
        SUM_MSB = MULT_MSB + int(np.log2(LEN))
        SUM_WIDTH_SIGNED = MULT_WIDTH_SIGNED + int(np.log2(LEN))
        w_sum = Signal(intbv(0, -2**SUM_WIDTH_SIGNED, 2**SUM_WIDTH_SIGNED))
        
        @always_comb
        def comb_mult():
            for i in range(LEN):
                w_mult[i].next = o_out_shift[i] * COEFFS[i]
        @always_comb
        def comb_sum():
            w_sum.next = sum(w_mult)
        
        @always_comb
        def exits():
            o_dout.next = w_sum[SUM_MSB:SUM_MSB-16].signed()
        return uut_shift, comb_mult, comb_sum, exits
    
from scipy import signal

input_data = []
output_data = []
S_r = 16000.0
Nyq = S_r / 2
Cutoff = 2000.0
clk_rate = 16000.0*16

prescaller = int(clk_rate / S_r)


def tb():
    global input_data
    global output_data
    global S_r
    global Nyq
    global Cutoff
    global precaller
    
    
    clk = Signal(bool(0))
    i_en = Signal(bool(0))
    i_rst = Signal(bool(0))
    @instance
    def clk_gen():
        while True:
            yield delay(1)
            clk.next = not clk

    defines = Cook_book()

    N = 4
    N_signed = N - 1
    i_adr = Signal(modbv( 0, 0, 2**N))
    
    n = 2**len(i_adr)
    print n
    a = signal.firwin(n, cutoff = 0.1, window = "hamming") 
    fir_coeffs = signal.firwin(n, Cutoff/Nyq)
    
    M = 16-1
    CONTENT = [int((-1+2**M)*fir_coeffs[i]) for i in range(n)]
    plt.plot(CONTENT, 'r')
    print CONTENT
    i_din = Signal(intbv(0, -2**M, 2**M))
    o_dout = Signal(intbv(0, -2**M, 2**M))
    uut = defines.serial_filter(clk, i_rst, i_en, i_din, o_dout, CONTENT, 16, len(CONTENT))
    cnt = Signal(intbv(0,0,16))
    @instance
    def control():
        #global prescaller
        yield clk.posedge
        i_rst.next = 1
        yield clk.posedge
        yield clk.posedge
        yield clk.posedge
        i_en.next = 1
        while True:  
            
            yield clk.posedge
            cnt.next = (cnt + 1 ) %16
            if cnt.next == 0:
                i_en.next = 1
            else:
                i_en.next = 0
    
    
    @always(clk.posedge)
    def monitor():
        if i_en == 1:
            tmp = randint(1-2**15, -1+2**15)
            input_data.append(tmp)
            output_data.append(int(o_dout))
            i_din.next = tmp
            i_adr.next = i_adr + 1
    return clk_gen, monitor, uut, control

inst = traceSignals(tb)
sim = Simulation(inst)
sim.run(2000)

fig, ax = plt.subplots(4,1)
ax[0].plot(input_data)
X_in = np.fft.fft(input_data)
X_in = abs(X_in) / max(abs(X_in))
X_in_db = 20*np.log10(X_in)
ax[1].plot(X_in_db)
ax[2].plot(output_data)
print output_data
X_out = np.fft.fft(output_data)
X_out = abs(X_out) / max(abs(X_out))
X_out_db = 20*np.log10(X_out) 
ax[3].plot(X_out_db)

Полученные коэффициенты фильтра:

И картинка для привлечения внимания: входной сигнал, спектр входного сигнала, выходной сигнал, спектр выходного сигнала