基于FPGA的cordic算法的verilog初步实现

小灰灰 2022-05-23 08:39 236阅读 0赞

## https://www.cnblogs.com/sepeng/p/5819586.html    ##

##     ##

## [基于FPGA的cordic算法的verilog初步实现][FPGA_cordic_verilog] ##

最近在看cordic算法，由于还不会使用matlab，真是痛苦，一系列的笔算才大概明白了这个算法是怎么回事。于是尝试用verilog来实现。用verilog实现之前先参考软件的程序，于是先看了此博文[http://blog.csdn.net/liyuanbhu/article/details/8458769 ][http_blog.csdn.net_liyuanbhu_article_details_8458769]也不截图了，因为怕图形被其他博客网站检测到后屏蔽图片，造成此博文无法正常阅读。

阅读此博文，需要先阅读上面这个博文的内容。

这是此博文中的C代码。避免浮点运算，所以angle数组里面的角度值都扩大了256倍。此程序中计算点（10,20）与X轴之间的夹角。最终运行结果是16238。而16238/256=63.43

![复制代码][copycode.gif]

1 #include <stdio.h>
     2 #include <stdlib.h>
     3 
     4 //double my_atan2(double x, double y);
     5 
     6 
     7 double  my_atan2  (int x, int y)
     8 {
     9     const int angle[] = {
         11520, 6801, 3593, 1824, 916, 458, 229, 115, 57, 29, 14, 7, 4, 2, 1};
    10 
    11     int i = 0;
    12     int x_new, y_new;
    13     int angleSum = 0;
    14 
    15     x *= 1024;// ½« X Y ·Å´óÒ»Ð©£¬½á¹û»á¸ü×¼È·
    16     y *= 1024;
    17 
    18     printf ("org_x = %d, org_y=%d\n",x,y);
    19 
    20     for(i = 0; i < 15; i++)
    21     {
    22         if(y > 0)
    23         {
    24             x_new = x + (y >> i);
    25             y_new = y - (x >> i);
    26             x = x_new;
    27             y = y_new;
    28             angleSum += angle[i];
    29         }
    30         else
    31         {
    32             x_new = x - (y >> i);
    33             y_new = y + (x >> i);
    34             x = x_new;
    35             y = y_new;
    36             angleSum -= angle[i];
    37         }
    38         printf("Debug: i = %d x=%d, y =%d angleSum = %d, angle = %d\n", i,x,y ,angleSum, angle[i]);
    39     }
    40     return angleSum;
    41 }
    42 
    43 void main()
    44 {
    45     double z=0 ; 
    46     z = my_atan2(10.0, 20.0);
    47     printf("\n z = %lf \n", z);
    48 
    49 }

![复制代码][copycode.gif]

既然有了C 就很好进行verilog设计。

先谈架构。

1，先将角度数据放到一个rom中进行存储

2，取一个数，运算一个数。直到循环结束。

这也就是所谓的cordic算法的向量模式。用来求角度。先看顶层

![复制代码][copycode.gif]

1 //
     2 //
     3 //
     4 //
     5 
     6 module cordic_rotation_top (
     7                             clock ,
     8                             rst_n ,
     9                             x_crd,
    10                             y_crd,
    11                             ena ,
    12                             deg_sum
    13                             ); 
    14 input         clock ; 
    15 input         rst_n ; 
    16 input     [15:0]    x_crd ; 
    17 input     [15:0]    y_crd ; 
    18 input             ena ; 
    19 output     [15:0]     deg_sum ; 
    20 
    21 
    22 wire     [4:0] deg_addr ; 
    23 wire     [15:0] deg_data ; 
    24 
    25 alt_ip_rom_cordic u_rom (
    26                             .address  (deg_addr),
    27                             .clock    (clock),
    28                             .q        (deg_data)
    29                         );
    30 
    31 
    32 
    33 cordic_rotation  u_cord (
    34                         .clk      (clock),
    35                         .rst_n    (rst_n),
    36                         .ena      (ena),
    37                         .x_crd    (x_crd),
    38                         .y_crd    (y_crd),
    39                         .deg_data (deg_data),
    40                         .deg_addr (deg_addr),
    41                         .deg_sum  (deg_sum)
    42                         );
    43                         
    44 endmodule

![复制代码][copycode.gif]

rom的初始化文件为

![138275-20160829201138777-1035377068.png][]

再看运算单元。

![复制代码][copycode.gif]

1 module cordic_rotation (
     2                         clk ,
     3                         rst_n ,
     4                         ena ,
     5                         x_crd,
     6                         y_crd,
     7                         deg_data,
     8                         deg_addr,
     9                         deg_sum
    10                         );
    11                         
    12 input     clk ; 
    13 input     rst_n ; 
    14 input     ena ; 
    15 input     [15:0]    x_crd ; //x coordinate
    16 input     [15:0]    y_crd ; //y coordinate
    17 input     [15:0]     deg_data ; 
    18 output     [4:0]     deg_addr ; 
    19 output     reg[15:0]     deg_sum ; 
    20 
    21 // ------ rotation count  0 - 14   -------
    22 reg     [4:0] rot_cnt ; 
    23 reg     [4:0] rot_cnt_r ; 
    24 wire             opr_en ; 
    25 assign opr_en = ((rot_cnt_r<5'd15)&(ena)) ; 
    26 
    27 always @ (posedge clk or negedge rst_n)
    28     if (!rst_n) begin 
    29             rot_cnt <= 5'd0 ; 
    30             rot_cnt_r<= 5'd0; 
    31             end 
    32     else if (opr_en) begin 
    33             rot_cnt    <= rot_cnt + 5'd1 ; 
    34             rot_cnt_r<= rot_cnt ; 
    35             end 
    36     else if (!ena) begin 
    37             rot_cnt <= 5'd0 ; 
    38             rot_cnt_r<= rot_cnt ; 
    39             end 
    40 assign deg_addr = rot_cnt ; 
    41 
    42 //---------------------------------------
    43 reg                 cal_cnt ; 
    44 reg signed [16:0]     x_d,    y_d ; 
    45 always @ (posedge clk or negedge rst_n)
    46     if (!rst_n) begin 
    47             x_d     <= 17'd0 ; 
    48             y_d     <= 17'd0 ; 
    49             deg_sum <= 16'd0 ; 
    50             cal_cnt <= 1'd0 ; 
    51             end 
    52     else if (opr_en) begin 
    53             case (cal_cnt)
    54                 1'd0 : begin 
    55                         x_d     <= {
         1'd0,x_crd}; 
    56                         y_d     <= {
         1'd0,y_crd};
    57                         deg_sum <= 16'd0 ; 
    58                         cal_cnt <= 1'd1 ; 
    59                         end 
    60                 1'd1 : begin 
    61                         if ((y_d[16])|(y_d==16'd0)) begin 
    62                                 x_d     <= x_d - (y_d >>> rot_cnt_r); 
    63                                 y_d     <= y_d + (x_d >>> rot_cnt_r) ;
    64                                 deg_sum <= deg_sum - deg_data ; 
    65                                 end 
    66                         else begin 
    67                                 x_d     <= x_d + (y_d >>> rot_cnt_r); 
    68                                 y_d     <= y_d - (x_d >>> rot_cnt_r) ;
    69                                 deg_sum <= deg_sum + deg_data ; 
    70                                 end 
    71                         end 
    72             endcase 
    73             end 
    74     else begin 
    75             cal_cnt <= 1'd0 ; 
    76             end
    77     
    78 endmodule

![复制代码][copycode.gif]

rot\_cnt作为rom的地址。但是由于rom返回度数值需要一个周期，所以下面的运算单元需要等待一个周期才可以进行运算，于是有了rot\_cnt\_r这个参数。

于是问题又来了，上面的是向量模式。还有一个旋转模式。于是开始捣鼓。

稍微将上述的C代码进行修改，计算30度的cos30，以及sin30

![复制代码][copycode.gif]

1 #include <stdio.h>
     2 #include <stdlib.h>
     3 
     4 //double my_atan2(double x, double y);
     5 
     6 
     7 double  my_atan2  (int x, int y)
     8 {
     9     const int angle[] = {
         11520, 6801, 3593, 1824, 916, 458, 229, 115, 57, 29, 14, 7, 4, 2, 1};
    10 
    11     int i = 0;
    12     int x_new, y_new;
    13     int angleSum = 30*256;
    14 
    15     x *= 1024;// ½« X Y ·Å´óÒ»Ð©£¬½á¹û»á¸ü×¼È·
    16     y *= 1024;
    17 
    18     printf ("org_x = %d, org_y=%d\n",x,y);
    19 
    20     for(i = 0; i < 15; i++)
    21     {
    22         if(angleSum <= 0)
    23         {
    24             x_new = x + (y >> i);
    25             y_new = y - (x >> i);
    26             x = x_new;
    27             y = y_new;
    28             angleSum += angle[i];
    29         }
    30         else
    31         {
    32             x_new = x - (y >> i);
    33             y_new = y + (x >> i);
    34             x = x_new;
    35             y = y_new;
    36             angleSum -= angle[i];
    37         }
    38         printf("Debug: i = %d x=%d, y =%d angleSum = %d, angle = %d\n", i,x,y ,angleSum, angle[i]);
    39     }
    40     return angleSum;
    41 }
    42 
    43 void main()
    44 {
    45     double z=0 ; 
    46     z = my_atan2(1.0, 0.0);
    47     printf("\n z = %lf \n", z);
    48 
    49 }

![复制代码][copycode.gif]

结果是这样子的

![138275-20160829201711121-1132739946.png][]

先看已知条件：

1，cordic算法中的Kn极限值=1.6476 。 1/Kn=0.6073.

2，上述软件设置y=0. x=1.并且坐标值扩大了1024倍。从公式上，cosa= x sina= y 角度a=30度

3，运算结果x=1461. y=844

好，开始运算 cosa=1461/(1024\*1.6476)=0.8659。 cos30=0.8660.

sina=844/(1024\*1.6476) = 0.500 。 sin30 = 0.5

精度由loop次数相关

好，既然验证了这个软件程序是可以的。那么就将它转换成Verilog的程序

![复制代码][copycode.gif]

1 //cordic 算法的旋转模式 
     2 
     3 
     4 module cordic_rotation (
     5                         clk ,
     6                         rst_n ,
     7                         ena ,
     8                         x_crd,
     9                         y_crd,
    10                         deg_data,
    11                         deg_addr,
    12                         deg_sum
    13                         );
    14                         
    15 input     clk ; 
    16 input     rst_n ; 
    17 input     ena ; 
    18 input     [15:0]    x_crd ; //x coordinate
    19 input     [15:0]    y_crd ; //y coordinate
    20 input     [15:0]     deg_data ; 
    21 output     [4:0]     deg_addr ; 
    22 output     [15:0]     deg_sum ; 
    23 
    24 
    25 parameter   DEGREE = 30 ;
    26 parameter     DEG_PARA =  (DEGREE*256) ; 
    27 // ------ rotation count  0 - 14   -------
    28 reg     [4:0] rot_cnt ; 
    29 reg     [4:0] rot_cnt_r ; 
    30 wire             opr_en ; 
    31 assign opr_en = ((rot_cnt_r<5'd15)&(ena)) ; 
    32 
    33 always @ (posedge clk or negedge rst_n)
    34     if (!rst_n) begin 
    35             rot_cnt <= 5'd0 ; 
    36             rot_cnt_r<= 5'd0; 
    37             end 
    38     else if (opr_en) begin 
    39             rot_cnt    <= rot_cnt + 5'd1 ; 
    40             rot_cnt_r<= rot_cnt ; 
    41             end 
    42     else if (!ena) begin 
    43             rot_cnt <= 5'd0 ; 
    44             rot_cnt_r<= rot_cnt ; 
    45             end 
    46 assign deg_addr = rot_cnt ; 
    47 
    48 //---------------------------------------
    49 reg                 cal_cnt ; 
    50 reg signed [15:0]     deg_sum_r  ; 
    51 reg signed [16:0]     x_d,    y_d ; 
    52 always @ (posedge clk or negedge rst_n)
    53     if (!rst_n) begin 
    54             x_d     <= 17'd0 ; 
    55             y_d     <= 17'd0 ; 
    56             deg_sum_r <= 16'd0 ; 
    57             cal_cnt <= 1'd0 ; 
    58             end 
    59     else if (opr_en) begin 
    60             case (cal_cnt)
    61                 1'd0 : begin 
    62                         x_d     <= {
         1'd0,x_crd}; 
    63                         y_d     <= {
         1'd0,y_crd};
    64                         deg_sum_r <= DEG_PARA ; 
    65                         cal_cnt <= 1'd1 ; 
    66                         end 
    67                 1'd1 : begin 
    68                         if ((deg_sum_r[15])|(deg_sum_r==16'd0)) begin  //<=0
    69                                 x_d     <= x_d + (y_d >>> rot_cnt_r); 
    70                                 y_d     <= y_d - (x_d >>> rot_cnt_r) ;
    71                                 deg_sum_r <= deg_sum_r + deg_data ; 
    72                                 end 
    73                         else begin 
    74                                 x_d     <= x_d - (y_d >>> rot_cnt_r); 
    75                                 y_d     <= y_d + (x_d >>> rot_cnt_r) ;
    76                                 deg_sum_r <= deg_sum_r - deg_data ; 
    77                                 end 
    78                         end 
    79             endcase 
    80             end 
    81     else begin 
    82             cal_cnt <= 1'd0 ; 
    83             end
    84             
    85 assign deg_sum = deg_sum_r ; 
    86     
    87 endmodule

![复制代码][copycode.gif]

这个程序也是在上一个程序的基础上修改的。所以尽量少改动。而是修改了tb。所以需要看看TB是怎么弄的

![复制代码][copycode.gif]

1 ///
     2 //
     3 //
     4 //
     5 `timescale 1ns/100ps
     6 
     7 
     8 module cordic_rotation_top_tb ; 
     9 reg clock ,rst_n ; 
    10 reg [15:0]     x_crd ,y_crd ; 
    11 reg ena ; 
    12 
    13 wire [15:0] deg_sum ; 
    14 
    15 
    16 cordic_rotation_top u_top (
    17                             .clock   (clock),
    18                             .rst_n   (rst_n),
    19                             .x_crd   (x_crd),
    20                             .y_crd   (y_crd),
    21                             .ena     (ena),
    22                             .deg_sum (deg_sum)
    23                             ); 
    24                             
    25 always #10 clock = ~clock ; 
    26 initial begin 
    27         clock = 0 ;  rst_n =0 ; 
    28         x_crd = 1*1024 ; y_crd=0; 
    29         ena = 1 ; 
    30         #40 rst_n = 1 ; 
    31         #350 
    32 /*        
    33         ena = 0 ; 
    34         #40 
    35         x_crd = 10*1024 ; y_crd=10*1024; 
    36         ena = 1 ; 
    37         #350 
    38 
    39         ena = 0 ; 
    40         #40 
    41         x_crd = 20*1024 ; y_crd=10*1024; 
    42         ena = 1 ; 
    43         #350  */
    44         $stop ;
    45         end 
    46                             
    47 endmodule

![复制代码][copycode.gif]

将第28行替换为 x\_crd = 10\*1024 ; y\_crd=20\*1024; 就变成了上一个工程的测试代码。

此程序输出结果为：

![138275-20160829203434855-1238489975.png][]

和软件程序输出结果相同。

有没有发现问题，在计算这些东西的时候，很多网络博文说精度和loop次数相关。但是这个扩大256倍运算到第11次就已经定住了。后面的输出结果白循环了。

[FPGA_cordic_verilog]: https://www.cnblogs.com/sepeng/p/5819586.html
[http_blog.csdn.net_liyuanbhu_article_details_8458769]: http://blog.csdn.net/liyuanbhu/article/details/8458769
[copycode.gif]: /images/20220523/263b7b3b1be545b7a90ceff6ddbb8712.png
[138275-20160829201138777-1035377068.png]: /images/20220523/d465a7c6b0dc49919c7e27d49628d7d8.png
[138275-20160829201711121-1132739946.png]: /images/20220523/71af034351634dc1939edafea3f6ba74.png
[138275-20160829203434855-1238489975.png]: /images/20220523/a7ad2ecc688a489eb00ed50d95c95297.png