简介
为了入门Verilog图像处理,熟悉modelsim的仿真是必不可少的,因此选择将直方图均衡作为一个上手项目,使用vivado对zynq进行开发,配合modelsim进行仿真,在此过程中需要借助matlab进行数据处理和结果展示。
原理
数据预处理
预处理借助matlab实现,因为modelsim可以读取txt文件,可以利用matlab将图像转为txt文件,再利用modelsim读取txt文件进行仿真。
首先,如果使用的图像不是灰度图,就需要将图像转换为灰度图
%将RGB888格式的图像转为灰度图像
clear;clc;
image=imread("lena.png");
%imshow(image),title('原始图片');
gray=rgb2gray(image);
imwrite(gray,"./gray.jpg");
imshow(gray)
然后将灰度图转为txt文件
%将灰度格式的图像转为8bit的hex格式
clear;clc;
file_path="./image_gray.txt";
file=fopen(file_path,"w");
image=imread("gray.jpg");
imshow(image),title('图片');
[image_width,image_height,channel]=size(image);
for x=1:image_height
for y=1 : image_width
fprintf(file,'%d\n',image(x,y));
end
end
fclose(file);
转换后可以打开image_gray.txt文件查看数据,可以发现其实际上就是将图像中的灰度值每一行存一个的方式存入txt文件中。
代码
顶层模块
顶层模块中使用了三个RAM ip核,用于存储数据:u_transform256x8存储映射关系、Ram32x256存储灰度直方图、u_image_mem存储原始图像,均衡化之后的图像并没有存储在FPGA里面,而是完成映射关系的计算之后直接输出了。这里有两个自定义模块:histogram负责计算直方图,并将直方图存储到Ram32x256中、reflect实现映射关系的计算并将映射关系存储到u_transform256x8中。顶层模块的主体是一个三个状态的状态机:
- 状态000:将外部输入的灰度值存入u_image_mem中,对图像进行存储,同时计算灰度直方图。其内部也是一个状态机:
整个top模块的忙信号wen管脚如果为1,则根据histogram的wea管脚进行拉低,因为该管脚意味着histogram模块已经完成计算,整个模块进入空闲状态。- 状态0:等待外部输入的灰度值,如果外部输入灰度值并且pin_clken为高电平,则将灰度值作为u_image_mem写入数值,并且自增u_image_mem的写入地址,并且灰度值作为histogram的输入和Ram32x256的读出地址,拉高top模块的wen管脚进入忙状态,进入状态001;
- 状态1:检查histogram模块的忙信号,如果处于空闲状态就拉高其写入信号pin_clken_hist,使能u_image_mem的写入,将灰度值写入对应的位置上,返回状态0。
Ram32x256的写入受控于histogram模块,histogram根据输入的灰度值读取Ram32x256中的数值,然后加1再次写入,实现直方图的计算。
- 状态001:在这个状态下,top模块的wen一直拉高,因为接来都在进行映射计算,因此不能接收下一帧图像。在这里面也有也给状态机:
- 状态0:此时Ram32x256输出地址为addr的内容,u_transform256x8的地址跳转为addr,将addr作为reflect模块输入,拉高reflect的输入时钟,进入状态1;
- 状态1:对addr进行自增,用于上一个周期拉高了reflect的输入时钟,所以在这个周期里reflec模块会将addr作为灰度和Ram32x256的输出作为参数进行计算。
u_transform256x8的写入受控于reflect模块,因为Ram32x256存储的是灰度直方图,所以其addr就是灰度值,RAM中的内容是灰度值对应的像素数量,因此需要将两者输入到reflect模块中计算灰度分布函数,计算映射关系,最后将原灰度需要映射的灰度值写入到u_transform256x8中。
- 状态100:在这个状态下,需要完成将原图像的灰度根据映射关系输出对应的像素,因此这里也是有一个状态机:
- 读取u_image_mem中存储的原图像灰度值;
- 将灰度值作为地址输入到u_transform256x8;
- 输出u_transform256x8的输出
u_transform256x8的地址就是原灰度值,而存储的内容则是原灰度值映射过去的灰度值。
module top_histogram(
input clk,
input rst_n,
input [7:0] gray,
input pin_clken,
output reg wen,//高电平表示此模块正忙,无法接受数据
output wire [7:0] hist_gray,
output reg pout_clken
);
parameter image_wdith = 10'd512;
parameter image_height = 10'd512;
parameter total_pixel = 20'h40000;
wire wen_hist;//灰度直方图统计模块的忙标志
reg wen_hist_r;
reg [7:0] addr;
wire [7:0] addr_w;
wire [31:0] dout;
wire [31:0] din;
wire wea;//灰度直方图的写入使能
reg [17:0] addr_image;//image RAM的地址
wire [7:0] dout_image;//输出的图像灰度
(* DONT_TOUCH = "1" *)
reg [7:0] din_image;//输入的图像灰度
reg wea_image;//image RAM的写入使能
reg [7:0] addr_trans;//转换表的地址
wire [7:0] addr_trans_w;//转换表的地址
wire [7:0] dout_trans;//转换表的输出数据
wire [7:0] din_trans;//转换表的输入数据
wire wea_trans;//转换表的写入使能
reg [31:0] pixel_index;//记录直方图已经统计的像素数量
reg [2:0] state;//直方图均衡的状态机,0:统计灰度值 1:采用灰度数量分布代替概率分布进行映射转换
reg [1:0] state_reflct;//建立映射关系时使用的状态机
reg pin_clken_reflect;//映射关系使用的像素时钟
reg pin_clken_hist;
wire pout_clken_reflect;
wire bussy_reflect;//映射模块忙标志
reg [7:0] gray_in_reflect;//映射模块的输入灰度值
reg [2:0] state_hiostogram;//直方图t统计的状态机
reg [2:0] state_grayout;//完成灰度变换的状态机
always @(posedge clk or negedge rst_n) begin
if(!rst_n) begin
wea_image<=1'b0;
pixel_index<=32'd0;
state<=3'b001;
state_reflct<=2'd0;
pin_clken_reflect<=1'b0;
state_grayout<=3'd0;
state_hiostogram<=3'd0;
pin_clken_hist<=1'b0;
wen<=1'b0;
pout_clken<=1'b0;
addr_image<=18'd0;
din_image<=8'd0;
end
else begin
case(state)
//统计灰度值和将灰度图像存入RAM
3'b001:begin
wen<=(wen==1'b1)?(~wea):wen;//如果wen为1了,就要根据wea进行拉低,因为wea拉高了,证明histogram已经处理完成了。
case(state_hiostogram)
3'd0:begin
pin_clken_hist<=1'b0;
wea_image<=1'b0;//禁止写入图像存储RAM
if(pixel_index<total_pixel) begin
if(pin_clken) begin
wen<=1'b1;//此时模块的忙于灰度直方图统计模块一致
pixel_index<=pixel_index+1'b1;
din_image<=gray;//准备将灰度值写入存储图像的RAM
addr_image<=pixel_index[17:0];//图像存储RAM的写入地址
addr<=gray;//将灰度值作为灰度直方图RAM的地址
//wea<=1'b0;//读出使能,读出对应灰度值的当前统计数量
state_hiostogram<=3'd1;
end
end
else begin
pixel_index<=32'd0;
addr<=8'd0;
state<=(state<<1);
end
end
3'd1:begin
if(!wen_hist)begin
pin_clken_hist<=1'b1;//让灰度直方图计算模块读入灰度值和已经统计的数量
wea_image<=1'b1;//写入图像存储RAM
state_hiostogram<=3'd0;
end
end
endcase
end
//建立映射关系
3'b010:begin
wen<=1'b1;//此时电路一直忙直至完成整个流程
case (state_reflct)
//根据已经赋值好的addr_1读出灰度分布RAM中的数值,并给出像素时钟给reflect模块
2'b00:begin
if(!bussy_reflect) begin
addr_trans<=addr;
pin_clken_reflect<=1'b1;
gray_in_reflect<=addr;
state_reflct<=2'b01;
end
end
//地址+1,准备读出下一个灰度值对应的像素数量
2'b01:begin
if(addr<255) begin
addr<=addr+1'b1;
end
else begin
addr<=8'd0;
state<=(state<<1);
end
pin_clken_reflect<=1'b0;
state_reflct<=2'b00;
end
endcase
end
//映射
3'b100:begin
case(state_grayout)
2'b00:begin
pout_clken<=1'b0;
//wea_image<=1'b0;//前面已经置0了
addr_image<=pixel_index;
state_grayout<=2'b01;
end
2'b01:begin
addr_trans<=dout_image;
state_grayout<=2'b10;
end
2'b10:begin
pout_clken<=1'b1;
state_grayout<=2'b00;
if(pixel_index<total_pixel) begin
pixel_index<=pixel_index+1'b1;
end
end
endcase
end
endcase
end
end
wire rst_hist=rst_n&state[0];//此处可以使state发生变化后一直复位此模块,输出wen_hist=0,wea=0,addr=0,din=0
histogram u_histogram (
.clk(clk),
.rst_n(rst_hist),
.gray(gray),//灰度图片输入
.pin_clken(pin_clken_hist),//对应的像素时钟
.pixel_count(dout),//读入RAM中对应灰度值的数量
.wen(wen_hist),//高电平表示此模块正忙,无法接受数据
.wea(wea),//写入RAM使能
.addr(),//写入地址
.pixel_count_r(din)//写入值
);
//定义一个位宽为32bit,位深为256的单口RAM,存储灰度直方图
blk_mem_gen_0 Ram32x256 (
.clka(clk), // input wire clka
.wea(wea), // input wire [0 : 0] wea,高电平时为写使能
.addra(addr), // input wire [7 : 0] addr
.dina(din), // input wire [31 : 0] din
.douta(dout) // output wire [31 : 0] dout
);
//图像存储的RAM 512x512
image_mem u_image_mem (
.clka(clk), // input wire clk
.wea(wea_image), // input wire [0 : 0] wea
.addra(addr_image), // input wire [17 : 0] addra
.dina(din_image), // input wire [7 : 0] dina
.douta(dout_image) // output wire [7 : 0] douta
);
wire rst_reflect=rst_n&state[1];
reflect u_reflect(
.clk(clk),
.rst_n(rst_reflect),
//.gray(gray_in_reflect),
.gray(addr),
.count(dout),
.pin_clken(pin_clken_reflect),
.gray_reflect(din_trans),
.pout_clken(wea_trans),
.bussy(bussy_reflect)
);
//存储灰度变换关系的RAM
transform256x8 u_transform256x8 (
.clka(clk), // input wire clka
.wea(wea_trans), // input wire [0 : 0] wea
.addra(addr_trans), // input wire [7 : 0] addra
.dina(din_trans), // input wire [7 : 0] dina
.douta(hist_gray) // output wire [7 : 0] douta
);
endmodule
直方图统计
这个模块实际上就是完成根据输入的灰度值,对Ram32x256里对应地址的内容加1。
`timescale 1ns / 1ps
//此模块为了统计灰度图像的直方图,输入为灰度图像,输出直方图到RAM中进行存储
module histogram (
input clk,
input rst_n,
input [7:0] gray,
input pin_clken,
input [31:0] pixel_count,//从外部输入,记录直方图已经统计的像素数量
output reg wen,//高电平表示此模块正忙,无法接受数据
output reg wea,//对外部的RAM读写使能,0为读,1为写
output reg [7:0] addr,//根据输入的灰度值输出地址
output reg [31:0] pixel_count_r//在已经统计的像素数量的基础上,加上1输出
);
reg [1:0] state;
always @(posedge clk or negedge rst_n)
begin
if(!rst_n)
begin
state<=2'd0;
wen<=1'b0;
wea<=1'b0;
addr<=8'd0;
pixel_count_r<=32'd0;
end
else begin
case(state)
2'd0: begin
wea<=1'b0;//读取直方图RAM的数据
if(pin_clken) begin
state<=2'd1;
addr<=gray;//准备直方图RAM写入的地址
wen<=1'b1;//处理当前输入的灰度,不接收其它数据
end
end
2'd1: begin
state<=2'd0;
pixel_count_r<=pixel_count+1'b1;//直方图已经统计的像素数量加1
wea<=1'b1;//写入数据
wen<=1'b0;//通知外部模块可以在下一个时钟周期输入
end
endcase
end
end
endmodule
映射关系计算
这个模块根据Ram32x256中的直方图内容计算概率密度,并且为了避免浮点计算,特意使用像素点数量代替对应的概率密度进行计算,但是整个过程的花费周期被延长了。
`timescale 1ns / 1ps
// 建立映射关系的模块,即映射从原图像素到输出像素的映射关系
module reflect(
input clk,
input rst_n,
input [7:0] gray,//输入像素灰度级
input [31:0] count,//输入像素的数量
input pin_clken,//输入像素时钟
output reg [7:0] gray_reflect,//输出的映射灰度级
output reg pout_clken,//输出像素时钟
output reg bussy//模块忙
);
//灰度概率密度函数p(gray)=num(gray)/total_pixel_num 其分布函数为s(i)=sum(p(gray)) k=0...i
//因为上面的计算中涉及大量的浮点计算,为了避免浮点计算,特意采用下面的方法进行改进
//改进后的概率密度函数p'(gray)=num(gray) 分布函数s'(i)=sum(p'(gray))) k=0...i
//原来的映射关系gray_r=s(i)*255(gray->gray_r)
//新的映射关系gray_r'*gray_level<=s'(i)*255<(gray_r+1)'*(gray_level)(gray->gray_r')
//gray_level=total_pixel_num/256(总的灰度级)
//使用数量代替概率密度避免浮点计算
reg[31:0] cdf;//记录灰度分布函数
reg [31:0] cdf_r=20'd0;//
parameter gray_level = 11'd1024;
parameter gray_level_half = 10'd512;
reg [7:0] gray_reflect_r;
reg [1:0] i;//状态机
always @(posedge clk or negedge rst_n) begin
if(!rst_n) begin
cdf<= 20'd0;
gray_reflect<=8'd0;
gray_reflect_r<=8'd0;
i<=3'd0;
pout_clken<=1'b0;
bussy<=1'b0;
end
else begin
case(i)
//第一步接收输入,累加分布函数
2'd0:begin
if(pin_clken)begin
//此处是为了加速一下那些不存在的像素,他们的映射将会是0
if(count!=32'd0) begin
cdf<=cdf+count;
bussy<=1'b1;
i<=2'd1;
end
else
pout_clken<=1'b1;
end
pout_clken<=1'b0;
end
//建立映射关系进行输出
2'd1:begin
if(cdf>=cdf_r&&cdf<cdf_r+gray_level_half)begin
gray_reflect<=gray_reflect_r;
pout_clken<=1'b1;
i<=1'b0;
bussy<=1'b0;
end
else if(cdf>=cdf_r+gray_level_half&&cdf<=cdf_r+gray_level) begin
if(gray_reflect_r<255) begin
gray_reflect<=gray_reflect_r+1'b1;
end
pout_clken<=1'b1;
i<=1'b0;
bussy<=1'b0;
end
else begin
if(gray_reflect_r<255) begin
gray_reflect_r<=gray_reflect_r+1'b1;
end
cdf_r<=cdf_r+gray_level;
end
end
endcase
end
end
endmodule
仿真文件
仿真文件中调用了Verilog的系统函数,读取最开始预处理得到的txt文件的内容输入到上述模块里面,再将模块输入存储到txt文件里面。
`timescale 1ns / 100ps
module tb_top_histogram();
integer file_hist;
integer file_gray;
initial begin
file_gray=$fopen("./image_gray.txt","r");
file_hist=$fopen("./histogram.txt","w");
end
reg clk=1'b0;
always #1 clk=~clk;
reg rst_n=1'b0;
initial begin
#2;
rst_n<=1'b1;
end
reg pin_clken;
reg [7:0] gray;
wire [7:0] gray_hist;
wire pout_clken;
wire wen;
reg [1:0] i=2'd0;
always @(posedge clk or negedge rst_n) begin
if(!rst_n) begin
gray<=8'd0;
pin_clken<=1'b0;
end
else begin
case (i)
2'd0:begin
$fscanf(file_gray,"%d\n",gray);
i<=i+1;
pin_clken<=1'b0;
end
2'd1:begin
if(!wen)begin
pin_clken<=1'b1;
i<=2'd0;
end
end
default:begin
i<=2'd0;
pin_clken<=1'b0;
end
endcase
end
end
top_histogram u_top_histogram(
.clk(clk),
.rst_n(rst_n),
.gray(gray),
.pin_clken(pin_clken),
.wen(wen),
.hist_gray(gray_hist),
.pout_clken(pout_clken)
);
parameter image_width=512;
parameter image_height=512;
parameter total_pixels=image_width*image_height;
integer j=0;
always @(posedge clk) begin
if(pout_clken)begin
$fwrite(file_hist,"%d\n",gray_hist);
j=j+1;
end
end
always @(posedge clk) begin
if(j==total_pixels)begin
$fclose(file_hist);
$fclose(file_gray);
$stop;
end
end
endmodule
数据后处理
仿真结束后会在Histogram_equalize\Histogram_equalize.sim\sim_1\behav\modelsim路径下生成histogram.txt的文件,将其拷贝到matlab脚本路径下,将其转为图像进行显示。
%hex数据转为灰度图进行显示
clear;clc;
image_width=512;
image_height=512;
file_path="./histogram.txt";
file=fopen(file_path,"r");
image_buffer=uint8(zeros(image_height,image_width));
for x=1:image_height
for y=1:image_width
image_buffer(x,y)=uint8(fscanf(file,"%d",1));
end
end
imshow(image_buffer),title('还原图片');
imwrite(image_buffer,"histogram.jpg");
fclose(file);
结果分析
用于Verilog与软件算法对比肯定会出现精度不足等情况,因此采用PNSR作为两者相似程度的衡量标准,PNSR值越大,表示两者越相似。
首先是采用matlab进行直方图均衡:
% 读取图像
img = imread('lena.png');
% 将图像转换为灰度图像
gray_img = rgb2gray(img);
% 计算灰度直方图
histogram = imhist(gray_img);
% 计算累积分布函数
cdf = cumsum(histogram) / numel(gray_img);
% 对图像进行直方图均衡化
equalized_img = cdf(gray_img + 1);
% 显示原始图像和均衡化后的图像
subplot(1, 2, 1);
imshow(gray_img),title('原始图像');
subplot(1, 2, 2);
imshow(equalized_img),title('均衡化后的图像');
imwrite(equalized_img,"histogram_matlab.jpg");
然后计算PNSR:
%计算两个图像的PNSR
clear;clc;
image_width=512;
image_height=512;
image1_buffer=imread("./histogram.jpg");
image2_buffer=imread("./histogram_matlab.jpg");
MSE=0.0;
for x=1:image_height
for y=1:image_width
MSE=MSE+(image1_buffer(x,y)-image2_buffer(x,y))^2;
end
end
MSE=single(MSE)/(image_width*image_height);
PNSR=20*log10(255/sqrt(MSE));
figure(1);
subplot(1,2,1);
imshow(image1_buffer),title('图片1');
subplot(1,2,2);
imshow(image2_buffer),title('图片2');
以下是两张图片的结果对比,一张是Verilog实现的算法的结果,一张是matlab实现的算法的结果:
可以看出两者的相似度是非常高的,证明了本次实现的硬件算法和matlab算法是等价的。
用于仿真的clk周期为2ns,所以如果在实际电路中采用50Mhz频率的时钟,那么处理一帧512x512的图像就需要47ms。
总结
本次实验实现了Verilog的直方图均衡化算法,并且和matlab实现的算法进行了对比,证明了两者的等价性。本次实验也对算法的实时性进行了估计,显然处理单张图像的实时性较差,但实际应用时应该结合流水线的模式,牺牲硬件资源换取速度的提升,但是需要注意的是,由于RAM空间对应的是实际的电路,并不像C/C++数组那样灵活,因此当图像分辨率发生改变时,需要重新生成RAM IP核以及修改硬件参数。项目工程