硬件平台
基于XLINX公司生产的AX7035开发板，具有HDMI输出输出，可以满足在没有示波器条件下输入输出回环测试。项目中仅使用了ROM ip核用来存储查找表计算根号、对数、cos、sin，可以移植到其他任意开发中，但HDMI输出波形可能无法观测到，只能通过示波器显示。
设计内容
设计内容主要分为两部分：高斯分布序列产生和HDMI显示。该项目侧重点是高斯白噪声产生，我主要介绍LFSR序列发生器和Box Muller转换设计思路。
LFSR伪随机数生成
该模块产生32位均匀分布序列，循环周期是2^64 = 1.8*10^19 。
利用64位斐波那契型LFSR，反馈多项式为 x^64 + x^63 + x^61 + x^60 + 1，可通过调用函数改变初始种子seed参数以产生不同的随机序列。
同时修改参数中的MIN和MAX大小设置输出随机序列的范围，理论上可以设置为-2,147,483,648——+2,147,483,648，也就是±2^31。
module LFSR_fibonacci #(
parameter seed = 64'd18446744070000000000,
parameter MIN = 1,
parameter MAX = 100
)(
input wire clk,
input wire rst_n,
input wire enable,
output wire [31:0] rand_out
);
wire signed [63:0] rand64;
wire signed [31:0] rand_out1;
wire signed [31:0] rand_out2;
reg [63:0] lfsr_reg;
always @(posedge clk or posedge rst_n) begin
if (!rst_n) begin
lfsr_reg <= seed; // 初始值可以根据需要设置为其他值
end else if (enable) begin
// LFSR 反馈多项式为 x^64 + x^63 + x^61 + x^60 + 1
lfsr_reg <= {lfsr_reg[62:0], lfsr_reg[63] ^ lfsr_reg[60] ^ lfsr_reg[61] ^ lfsr_reg[63]};
end
end
assign rand64 = lfsr_reg;
assign rand_out1 = rand64[47:16] + {rand64[15:0],rand64[63:48]};
assign rand_out2 = rand64[31:0] + rand64[63:32];
assign rand_out = (rand_out1 + rand_out2) % (MAX - MIN+1'b1) + MIN;
endmodule
注意，我将LFSR设置为64位，而输出是经过64位随机序列变换得到，因此产生的随机序列周期是64位线性移位寄存器产生的伪随机数周期，但输出位数为32位。如果有对于资源要求比较高可以减小线性移位寄存器位数，通过查找响应位数特征多项式（不同位数有不同的特征多项式，这决定了产生随机序列的周期，使用特征多项式可以达到最大周期2……N-1）。
对于输出变换的解释，由LFSR直接产生的随机序列在一段时间内小数特别多，不符合我们对于均匀随机序列的要求，具体可编写测试文件观测到。因此，我对输出64位随机序列进行变换，这里我用两个rand_out1和rand_out2相加的结果取模，rand_out1和rand_out2分别是rand64取不同位数得到。当然，这里的方法不唯一，感兴趣的可以自行设计。
同时，需要注意的是，在输出级采用MAX和MIN对输出范围做了限制，输出的随机序列最大最小值均是有符号数！想要获得任意范围定点数可以自行另外设计输出逻辑。
更详细的LFSR设计原理和结构可以参照
https://cloud.tencent.com/developer/article/2287083?areaId=106001
Box-Muller转换产生高斯分布序列
Box-Muller transformation算法有两种形式，一种称之为基本形式：
其中U1，U2是均匀分布的随机数，Z1和Z2是满足高斯分布的随机数。可以看出，要得到高斯分布序列需要两个均匀分布序列，同时经过算数运算得到。
另一种表达形式是极坐标形式我没不使用，原因是极坐标形式Box-Muller转换虽然不涉及cos、sin运算，但引入了除法，这在FPGA设计中需要占用较多逻辑资源，我就使用基本形式对其变换。
利用使用ROM查找表，设计ln、根号、正余弦计算器。根据输入变量大小实现地址映射，考虑到设计精度等因素，对输入输出进行一定比例放大。对于cos只需要存储1/4周期的数据，其余数据均可通过地址映射得到这里不展开细讲。
module cos(
input clk,
input [15:0] x_in,
input enable,
output [15:0] result_out
);
wire [9:0] addr_ROM;
wire [13:0] ROM_out;
cos_ROM cos_ROM_0(
.clka (~clk),
.ena (enable),
.addra (addr_ROM),
.douta (ROM_out)
);
wire [15:0] angle_3600;
assign angle_3600 = (x_in[15] == 1'b1) ? ((~x_in[15:0]+1'b1) % 16'd3600) : (x_in[15:0] % 16'd3600);
wire sign;
wire [15:0] addr_angle;
assign addr_angle = (angle_3600 <= 16'd900) ? (angle_3600 ) :
(angle_3600 < 16'd1800) ? (16'd1800 - angle_3600) :
(angle_3600 < 16'd2700) ? (angle_3600 - 16'd1800) :
(angle_3600 < 16'd3600) ? (16'd3600 - angle_3600) : 16'd901;
assign sign = (angle_3600 <= 16'd900) ? (1'b0 ) :
(angle_3600 < 16'd1800) ? (1'b1 ) :
(angle_3600 < 16'd2700) ? (1'b1 ) :
(angle_3600 < 16'd3600) ? (1'b0 ) : 1'b0;
assign addr_ROM = addr_angle[9:0];
reg sign_delay0;
always @(posedge clk) begin
sign_delay0 <= sign;
end
assign result_out = sign_delay0 ? (~{2'b0,ROM_out} + 1'b1): {2'b0,ROM_out}; //
endmodule
输入是角度制，同时需要将角度放大10倍，也就是3600是一个周期。输出16位有符号数，为是结果放大10000倍。
另外，是利用流水线的计算模式，ln和cos是第一级流水线，ln计算的结果乘以-2后再取根号，因此需要把cos计算的结果延迟一个节拍，再相乘输出。
Box Muller 转换流水线代码如下
module Box_Muller_transform(
input wire clk_data,
input wire clk_cal,
input wire enable,
input wire signed [15:0] data_in_1,
input wire signed [15:0] data_in_2,
output reg signed [31:0] result_out_0,
output reg signed [31:0] result_out_1
);
wire signed [15:0] cos_sin_in;
wire signed [15:0] sqrt_in;
wire signed [15:0] sqrt_out, ln_out, cos_out, sin_out;
reg signed [15:0] cos_out_delay, sin_out_delay;
reg signed [15:0] x_in_1, x_in_2;
always @(posedge clk_data) begin
if (enable) begin
x_in_1 <= data_in_1;
x_in_2 <= data_in_2;
end
end
// instantiate sqrt module
assign sqrt_in = ln_out * (-16'd2);
sqrt sqrt_inst (
.clk(clk_cal),
.enable(enable),
.x_in(sqrt_in),
.result_out(sqrt_out)
);
// instantiate ln module
ln ln_inst (
.clk(clk_cal),
.enable(enable),
.x_in(x_in_1),
.result_out(ln_out)
);
// instantiate cos module
assign cos_sin_in = x_in_2 * 16'd36 / 16'd10;
cos cos_inst (
.clk(clk_cal),
.enable(enable),
.x_in(cos_sin_in),
.result_out(cos_out)
);
sin sin_inst (
.clk(clk_cal),
.enable(enable),
.x_in(cos_sin_in),
.result_out(sin_out)
);
// delay cos result
always @(posedge clk_cal) begin
if (enable) begin
cos_out_delay <= cos_out;
sin_out_delay <= sin_out;
end
end
// calculate final result using your formula
always @(posedge clk_data) begin
if (enable) begin
result_out_0 <= cos_out_delay * sqrt_out;
result_out_1 <= sin_out_delay * sqrt_out;
end
end
endmodule
高斯序列产生
该模块例化了两个LFSR随机数生成模块（利用不同的种子）和一个Box Muller转换模块，将两个均匀分布序列转换为高斯分布序列。
同时，可以在模块参数中修改均值和方差，以达到产生任意分布正态分布的作用。默认为标准正态分布。
/默认产生32位有符号数，结果除以1M为实际值
module Guassian_generate #(
parameter average = 32'd0, //均值
parameter variance = 32'd0 //方差
)(
//input sys_clk,
input clk_data,
input clk_tra,
input rst_n,
input enable,
debug
//output wire [31:0] rand_out_1,
//output wire [31:0] rand_out_2,
output wire signed [31:0] rand_out1,
output wire signed [31:0] rand_out2
);
wire [31:0] rand_out_1, rand_out_2;
LFSR_fibonacci #(
.seed(64'd18446744070000000000),
.MIN(32'd1),
.MAX(32'd1000)
) Uniform_Distribution_generator_0(
.clk(clk_data),
.rst_n(rst_n),
.enable(1'b1),
.rand_out(rand_out_1)
);
LFSR_fibonacci #(
.seed(64'd14742147901846518232),
.MIN(32'd1),
.MAX(32'd1000)
) Uniform_Distribution_generator_1(
.clk(clk_data),
.rst_n(rst_n),
.enable(1'b1),
.rand_out(rand_out_2)
);
//将均匀分布输出数据转化为16位输入
wire [15:0] tra_data_in_1, tra_data_in_2;
assign tra_data_in_1 = rand_out_1[15:0];
assign tra_data_in_2 = rand_out_2[15:0];
wire signed [31:0] Guassian_0;
wire signed [31:0] Guassian_1;
Box_Muller_transform Box_Muller_transform(
.clk_data(clk_data),
.clk_cal(clk_tra),
.enable(1'b1),
.data_in_1(tra_data_in_1),
.data_in_2(tra_data_in_2),
.result_out_0(Guassian_0),
.result_out_1(Guassian_1)
);
assign rand_out1 = enable ? ((Guassian_0[31] == 1'b1) ? (~((~Guassian_0 + 1'b1) / variance ) + 1'b1) + average : Guassian_0 / variance + average) : 32'd0;
assign rand_out2 = enable ? ((Guassian_1[31] == 1'b1) ? (~((~Guassian_1 + 1'b1) / variance ) + 1'b1) + average : Guassian_1 / variance + average) : 32'd0;
endmodule
注意该模块默认产生32位有符号数，结果除以1000000为实际值。这里实际值是标准正态分布的结果，因为根据Box Muller转换，输入均匀分布的范围是0-1，我利用LFSR产生1-1000的随机数，相当与输入放大了1000倍。
仿真测试
测试LFSR伪随机数生成均匀分布序列
使用vivado内置仿真环境
matlab仿真
延长仿真时间生成10M个样本
LFSR仿真测试代码和MATLB验证代码同下正态高斯测试一样，只需稍加修改即可。
高斯序列产生测试Verilog代码
`timescale 1ns / 1ps
//
// Company:
// Engineer:
//
// Create Date: 2024/01/14 12:37:28
// Design Name:
// Module Name: Guassian_generate_tb
// Project Name:
// Target Devices:
// Tool Versions:
// Description:
//
// Dependencies:
//
// Revision:
// Revision 0.01 - File Created
// Additional Comments:
//
//
module Guassian_generate_tb;
reg clk_data;
reg clk_tra;
reg rst_n;
reg enable;
wire signed [31:0] rand_out1;
Guassian_generate#(
.average(32'd0),
.variance(32'd0)
)Guassian_generate_inst (
.clk_data(clk_data),
.clk_tra(clk_tra),
.rst_n(rst_n),
.enable(enable),
.rand_out1(rand_out1)
);
wire [7:0] rand;
assign rand = rand_out1[7:0];
// Clock generation
initial begin
clk_tra = 0;
forever #3 clk_tra = ~clk_tra;
end
initial begin
clk_data = 0;
#7;
forever #10 clk_data = ~clk_data;
end
integer Guss0;
//integer Guss1;
initial begin
Guss0=$fopen("randGuss0.txt"); //打开所创建的文件
// Guss1=$fopen("randGuss1.txt"); //打开所创建的文件
// rand0=$fopen("rand0.txt"); //打开所创建的文件
// rand1=$fopen("rand1.txt"); //打开所创建的文件
end
always @(negedge clk_data) begin
if(enable) begin
$fdisplay(Guss0,"%d",$signed(rand_out1)); //$fdisplay(dout_file,"%d",$signed(num)); //保存有符号数据
// $fdisplay(Guss1,"%d",$signed(Z_standard1)); //$fdisplay(dout_file,"%d",$signed(num)); //保存有符号数据
// $fdisplay(rand0,"%d",$signed(rand_out_1)); //$fdisplay(dout_file,"%d",$signed(num)); //保存有符号数据
// $fdisplay(rand1,"%d",$signed(rand_out_2)); //$fdisplay(dout_file,"%d",$signed(num)); //保存有符号数据
end
end
// 在仿真结束时关闭文件
initial begin
#10000000; // 确保仿真足够长的时间
$fclose(Guss0);
// $fclose(Guss1);
// $fclose(rand0);
// $fclose(rand1);
$stop;
end
// Testbench stimulus
initial begin
rst_n = 0;
enable = 0;
#100
rst_n =1;
#100
enable = 1;
// Finish simulation
#10000000;
end
endmodule
将测试文件生成的随机数存储在randGuss0.txt文件中，再通过matlab读取进行分析。
可以看出产生的随机序列为一个标准的正态分布！
matlab测试代码如下。注：这里没有使用绝对路径查看txt文件，'randGuss0.txt'使用相对路径保存的文件再sim目录下，请自行查找。
%生成rand64向量
vivado_data_Guss0 = importdata('randGuss0.txt'); %read files
% 统计数据分布
minValue = min(vivado_data_Guss0);
maxValue = max(vivado_data_Guss0);
meanValue = mean(vivado_data_Guss0);
medianValue = median(vivado_data_Guss0);
stdDeviation = std(vivado_data_Guss0);
% 显示统计结果
fprintf('最小值：%f
', minValue);
fprintf('最大值：%f
', maxValue);
fprintf('平均值：%f
', meanValue);
fprintf('中位数：%f
', medianValue);
fprintf('标准差：%f
', stdDeviation);
% 绘制直方图
figure(1)
histogram(vivado_data_Guss0, 100, 'Normalization', 'probability');
title('guss1数据分布直方图');
xlabel('数值');
ylabel('频率');
DC输出高斯分布序列AD采集，HDMI显示
