基于FPGA的低通滤波器,通过verilog实现并提供testbench测试文件

Posted 51matlab

tags:

篇首语:本文由小常识网(cha138.com)小编为大家整理,主要介绍了基于FPGA的低通滤波器,通过verilog实现并提供testbench测试文件相关的知识,希望对你有一定的参考价值。

1.算法仿真效果

matlab2022a仿真结果如下:

 

 

 

 

 

2.算法涉及理论知识概要

       FIR(Finite Impulse Response)滤波器:有限长单位冲激响应滤波器,又称为非递归型滤波器,是数字信号处理系统中最基本的元件,它可以在保证任意幅频特性的同时具有严格的线性相频特性,同时其单位抽样响应是有限长的,因而滤波器是稳定的系统。因此,FIR滤波器在通信、图像处理、模式识别等领域都有着广泛的应用。

        在进入FIR滤波器前,首先要将信号通过A/D器件进行模数转换,把模拟信号转化为数字信号;为了使信号处理能够不发生失真,信号的采样速度必须满足香农采样定理,一般取信号频率上限的4-5倍做为采样频率;一般可用速度较高的逐次逼进式A/D转换器,不论采用乘累加方法还是分布式算法设计FIR滤波器,滤波器输出的数据都是一串序列,要使它能直观地反应出来,还需经过数模转换,因此由FPGA构成的FIR滤波器的输出须外接D/A模块。FPGA有着规整的内部逻辑阵列和丰富的连线资源,特别适合于数字信号处理任务,相对于串行运算为主导的通用DSP芯片来说,其并行性和可扩展性更好,利用FPGA乘累加的快速算法,可以设计出高速的FIR数字滤波器。

(1) 系统的单位冲激响应h (n)在有限个n值处不为零

(2) 系统函数H(z)|z|>0处收敛,极点全部在z = 0处(因果系统)

(3) 结构上主要是非递归结构,没有输出到输入的反馈,但有些结构中(例如频率抽样结构)也包含有反馈的递归部分。

FIR滤波器的单位冲激响应h (n)为一个N点序列,0 n N 1,则滤波器的系统函数为

H(z)=h(n)*z^-k

就是说,它有(N1)阶极点在z = 0处,有(N1)个零点位于有限z平面的任何位置。

 

        对于FIR(有限长单位冲激响应)滤波器,其基本结构是一个分节的延时线,每一节的输出加权累加,得到滤波器的输出。其输出y就是输入x和系数h的内积:

 

 

 

 

3.MATLAB核心程序

 

module xfilter(clk, rst, 
               i_valid_new_pixel,
               i_new_pixel,
               i_valid_lpos, 
               i_valid_cpos, 
               i_valid_rpos, 
               i_rowM,
               o_valid_filt,
               o_filt_pixel,
               o_colN,
               o_rowM
               );
parameter XB = 10;
parameter PB = 8;
 
input clk;
input rst;
input i_valid_new_pixel;
input [(PB + 2) - 1:0] i_new_pixel;
input i_valid_lpos;
input i_valid_cpos;
input i_valid_rpos;
input i_rowM;
 
output o_valid_filt;
output [PB - 1:0] o_filt_pixel;
output o_colN;
output o_rowM;
 
reg [(PB + 4) - 1:0] r_pixel_sum = 0;
reg [PB - 1:0] rr_pixel_sum = 0;
 
reg r_valid_lpos = 0;
reg r_valid_cpos = 0;
reg r_valid_rpos = 0;
 
reg r_valid_filt_pixel = 0;
reg rr_valid_cpos = 0;
 
reg r_valid_pixel_in = 0;
reg rr_valid_pixel_in = 0;
........................................................................
 
always @(posedge clk)
begin
  //read in pixel
  r_valid_lpos <= i_valid_lpos;
  r_valid_cpos <= i_valid_cpos;
  r_valid_rpos <= i_valid_rpos;
  rr_valid_cpos <= r_valid_cpos;
 
  r_colN <= i_valid_rpos;
  rr_colN <= r_colN;
  rrr_colN <= rr_colN;
 
  r_rowM <= i_rowM;
  rr_rowM <= r_rowM;
  rrr_rowM <= rr_rowM;
 
  if (rst)
    begin
      r_valid_pixel_in <= 0;
      rr_valid_pixel_in <= 0;
      r_valid_filt_pixel <= 0;
    end
  else
    begin
      r_valid_pixel_in <= c_valid_pixel_in;
      rr_valid_pixel_in <= r_valid_pixel_in;
      r_valid_filt_pixel <= rr_valid_pixel_in;
    end
 
  r_pixel2 <= c_pixel2;
end
 
wire [PB + 4 - 1: 0] c_pixel_sum = (rr_valid_cpos) ? r_pixel_sum + r_pixel2 + 8: r_pixel_sum + 8;
reg r_valid_new_pixel = 0; //needed to check if new row\'s data was added
//add filter
always @(posedge clk)
begin
  r_valid_new_pixel <= i_valid_new_pixel;
  if (r_valid_rpos)
  begin
    if (r_valid_new_pixel) //new data added to queue
      r_pixel_sum <= c_pixel1, 1\'d0 + c_pixel2;
    else
      r_pixel_sum <= c_pixel0, 1\'d0 + c_pixel1;
  end
  else
  begin
    r_pixel_sum <= c_pixel0 + c_pixel1, 1\'d0;
  end
 
  rr_pixel_sum <= c_pixel_sum[PB + 4 -1: 4];
end
 
//3 pixel buffer
queue PIXBUF (.clk(clk), 
              .i_valid_pixel(i_valid_new_pixel),
              .i_pixel(i_new_pixel),
              .o_pixel0(c_pixel0),
              .o_pixel1(c_pixel1), 
              .o_pixel2(c_pixel2)
             );
endmodule
 
module queue(clk, i_valid_pixel, i_pixel, o_pixel0, o_pixel1, o_pixel2);
parameter PB = 8;
.....................................................
endmodule

  

基于冷冻电镜图像的低通滤波(Lowpass Filter)算法

目标:将冷冻电镜图像,通过低通滤波算法降噪,输入冷冻电镜图像,输出信息集中的降噪图像。

脚本位置:cryosparc/cryosparc_master/cryosparc_compute/jobs/imports/run.py

Particles:Lowpass Filtered Images:

Lowoass Filtered源码位置:

逻辑:

fig = plotutil.plot_images_simple(dispparts, rows=4, cols=8, radwn=6, figscale=1.5)
rc.log_plot(fig, 'Lowpass Filtered Images: ')

plotutil:cryosparc/cryosparc_master/cryosparc_compute/plotutil.py

解耦lowpass + fft逻辑:

img = images[imgi]
rD = img
N = rD.shape[0]
D = fourier.fft(rD)
Cval = -1.0
lowpass_mask = sigproc.lowpass_mask(N, 2, self.radwn, lowpassorder=self.lowpassorder)
dispim = fourier.ifft(lowpass_mask*D*Cval)
vmax = max(-dispim.min(), dispim.max())*vscale
vmin = -vmax
cmap='gray'
imshow(dispim , cmap, colorbar=False, vmin=vmin, vmax=vmax)

测试图像:

img_dir = "/nfs_baoding/chenlong/datasets/cryoEM/kongfang/cryo_dataset_4_imagenet_v1_1/train/1"
img_name = "018444013402718290241_stack_1348_cor2_DW_particles.mrc_7_1.png"
img_path = os.path.join(img_dir, img_name)
img_gray = cv2.imread(img_path, cv2.IMREAD_UNCHANGED)

dispim = img_gray

vscale = 1.5
vmax = max(-dispim.min(), dispim.max())*vscale
vmin = -vmax

cmap='gray'

# 原始图像
imshow(dispim, cmap, colorbar=False, vmin=vmin, vmax=vmax)

原始图像:

增加Conda环境,核心包:pyfftw

python -m ipykernel install --user --name cryosparc-master --display-name "cryosparc-master"

查看Conda环境,参考:Anaconda 查看、创建、管理和使用python环境

conda info --env

# conda environments:
#
base                  *  /home/chenlong/miniconda3
torch-def                /home/chenlong/miniconda3/envs/torch-def
                         /nfs_baoding/chenlong/workspace/cryosparc/cryosparc_master/deps/anaconda
                         /nfs_baoding/chenlong/workspace/cryosparc/cryosparc_worker/deps/anaconda

Conda重命名,参考:How can I rename a conda environment?

tmux attach -t download
conda create -n cryosparc-master --clone /nfs_baoding/chenlong/workspace/cryosparc/cryosparc_master/deps/anaconda

Source:      /nfs_baoding/chenlong/workspace/cryosparc/cryosparc_master/deps/anaconda
Destination: /home/chenlong/miniconda3/envs/cryosparc-master
The following packages cannot be cloned out of the root environment:
 - https://repo.anaconda.com/pkgs/main/linux-64::conda-4.8.3-py37_0
Packages: 34
Files: 1

测试脚本:

  • cryosparc/cryosparc_master/cryosparc_compute/20220802-lowpass_filter.ipynb
  • cryosparc/cryosparc_master/cryosparc_compute/20220802-lowpass_filter.py

matplotlib存储image,通过savefig,设置bbox_inchespad_inches,去除白边,但是,转换为numpy,无法去除左右的白边。

plt.savefig('match.png', bbox_inches='tight', transparent="True", pad_inches=0)

通过以下方案转换numpy,去除白边:

  • 设置fig尺寸为输入图像尺寸:fig.set_size_inches(plt.figaspect(dispim))
  • 全部填充,没有pad:fig.tight_layout(pad=0, h_pad=None, w_pad=None)
fig = plt.gcf()
fig.set_size_inches(plt.figaspect(dispim))   # 设置fig尺寸为输入图像尺寸,例如dispim
fig.tight_layout(pad=0, h_pad=None, w_pad=None)   # 全部填充,没有pad
fig.canvas.draw()
img_plt = np.frombuffer(fig.canvas.tostring_rgb(), dtype=np.uint8)
img_plt = img_plt.reshape(fig.canvas.get_width_height()[::-1] + (3,))

存储matplotlib的图像,参考:

冷冻电镜图像的低通滤波,逻辑如下:

class LowpassFilter(object):
    """
    冷冻电镜图像的低通滤波
    """
    def __init__(self):
        self.radwn = 6  # 低通滤波mask参数
        self.lowpassorder = 1  # 低通滤波等级,mask参数

    @staticmethod
    def imshow(rD, cmap='turbo', colorbar=True, vmin=None, vmax=None):
        """
        展示图像
        """
        plt.imshow(rD, cmap=plt.get_cmap(cmap), vmin=vmin, vmax=vmax)
        if colorbar:
            plt.colorbar()
            # plt.colorbar(orientation='horizontal', shrink=0.6, fraction=0.05, pad=0.05);
        else:
            fig = plt.gcf()
            plt.axis('off')
            for ax in fig.axes:
                ax.get_xaxis().set_visible(False)
                ax.get_yaxis().set_visible(False)

    def lowpass_filter_for_cryoem(self, img_gray):
        """
        冷冻电镜图像的低通滤波算法
        """
        rD = img_gray
        N = rD.shape[0]

        # 计算低通滤波的mask
        lowpass_mask = sigproc.lowpass_mask(N, 2, self.radwn, lowpassorder=self.lowpassorder)
        D = fourier.fft(rD)  # 傅里叶变换
        Cval = -1.0
        dispim = fourier.ifft(lowpass_mask * D * Cval)  # 傅里叶变换

        vscale = 1.5
        vmax = max(-dispim.min(), dispim.max()) * vscale
        vmin = -vmax
        cmap = 'gray'

        # 原始图像
        LowpassFilter.imshow(dispim, cmap, colorbar=False, vmin=vmin, vmax=vmax)

        # 获取matplotlib显示的image,转换为numpy格式
        fig = plt.gcf()
        fig.set_size_inches(plt.figaspect(dispim))
        fig.tight_layout(pad=0, h_pad=None, w_pad=None)
        fig.canvas.draw()
        img_plt = np.frombuffer(fig.canvas.tostring_rgb(), dtype=np.uint8)
        img_plt = img_plt.reshape(fig.canvas.get_width_height()[::-1] + (3,))

        img_lowpass = cv2.cvtColor(img_plt, cv2.COLOR_RGB2GRAY)
        img_lowpass = cv2.resize(img_lowpass, img_gray.shape)

        return img_lowpass

低通滤波图像:

v1.0的lowpass filter算法,如下:

class LowpassFilter(object):
    """
    低通滤波类
    """
    def __init__(self):
        pass

    @staticmethod
    def gen_mask(image_size, radius=20):
        x, y = np.ogrid[0:image_size, 0:image_size]
        center_x = (image_size - 1) / 2
        center_y = (image_size - 1) / 2
        mask = (x - center_x) ** 2 + (y - center_y) ** 2 <= radius ** 2
        return np.array(mask, dtype=np.int32)

    @staticmethod
    def lowpass(img):
        """
        低通滤波
        """
        image_size = img.shape[0]
        img_dft = np.fft.fft2(img)
        dft_shift = np.fft.fftshift(img_dft)
        dft_mask = LowpassFilter.gen_mask(image_size, radius=20)
        dft_filter = dft_mask * dft_shift
        idft_shift = np.fft.ifftshift(dft_filter)
        img_back = np.abs(np.fft.ifft2(idft_shift))
        img_back = img_back.astype(np.uint8)

        return img_back

效果对比:

像素分布:

以上是关于基于FPGA的低通滤波器,通过verilog实现并提供testbench测试文件的主要内容,如果未能解决你的问题,请参考以下文章

FPGA教程案例46图像案例6——基于FPGA的图像高斯滤波verilog实现,通过MATLAB进行辅助验证

FPGA教程案例44图像案例4——基于FPGA的图像中值滤波verilog实现,通过MATLAB进行辅助验证

FPGA教程案例53语音案例2——基于FIR低通滤波器的语音信号降噪FPGA实现

基于冷冻电镜图像的低通滤波(Lowpass Filter)算法

基于FPGA的LMS自适应滤波器verilog实现,包括testbench

使用核心音频实现后处理低通滤波器