• 用FPGA CORDIC IP核实现信号的相位检测,计算相位角


    用FPGA CORDIC IP核实现信号的相位检测

    1.matlab仿真

    波形仿真代码:

    代码功能:生成一个点频信号s,求出s的实部和虚部;并且结算相位角atan2。画出图形,并且将Q和I数据写入文件中。

    %代码功能:生成一个点频信号s,求出s的实部和虚部;并且结算相位角atan2。画出图形,并且将Q和I数据写入
    clc; clear;close all;
    
    F1=1; %信号频率
    Fs=65536; %采样频率					Fs=N,能采一个周期
    P1=45; %信号初始相位(单位:°),90-cos函数
    N=65536; %采样点数
    t=[0:1/Fs:(N-1)/Fs]; %采样时刻
    A=2^15-1; %信号幅度
    
    %生成点频信号
    s=A*exp(1i *(2*pi*F1*t + pi*P1/180));
    
    % IQ分解,分别提取实部和虚部
    I = real(s);
    Q = imag(s);
    
    % 计算相位角
    phase = atan2(Q, I);
    
    % 提取初始相位(t=0处的相位值)
    initial_phase = phase(1);
    %显示初始相位值,将其转换为π的倍数进行显示 
    disp(['初始相位值:', num2str(initial_phase/pi),'pi']);
    
    
    % 绘制结果图,画出signal信号的实部和虚部
    figure;
    subplot(3,1,1);
    plot(t, real(s), 'b');
    title('In-phase Component (I)');
    
    subplot(3,1,2);
    plot(t, imag(s), 'r');
    title('Quadrature Component (Q)');
    xlabel('Time (s)');
    
    % 计算相位角
    subplot(3,1,3);
    plot(t, phase);
    xlabel('时间');
    ylabel('相位');
    title('相位随时间变化');
    %将纵坐标转化为Π的倍数
    yticks([-1*pi, -0.5*pi, 0, 0.5*pi, pi]);
    yticklabels({'-π', '-0.5π', '0', '0.5π', 'π'});
    
    
    %创建 coe 文件  Idata
     fild = fopen('Idata_65536x15bit.coe','wt');
     %写入 coe 文件头
     %固定写法,表示写入的数据是 10 进制表示
     fprintf(fild, '%s\n','memory_initialization_radix=10;');
     %固定写法,下面开始写入数据
     fprintf(fild, '%s\n\n','memory_initialization_vector ='); 
     for i = 1:N
     s2(i) = round(I(i)); %对小数四舍五入以取整
     fprintf(fild, '%d',s2(i)); %数据写入
     if i==N
     fprintf(fild, '%s\n',';'); %最后一个数据用;
     else
     fprintf(fild,',\n'); % 其他数据用,
     end
     end
     fclose(fild); % 写完了,关闭文件
    
     %创建 coe 文件  Qdata
     fild = fopen('Qdata_65536x15bit.coe','wt');
     %写入 coe 文件头
     %固定写法,表示写入的数据是 10 进制表示
     fprintf(fild, '%s\n','memory_initialization_radix=10;');
     %固定写法,下面开始写入数据
     fprintf(fild, '%s\n\n','memory_initialization_vector ='); 
     for i = 1:N
     s3(i) = round(Q(i)); %对小数四舍五入以取整
     fprintf(fild, '%d',s3(i)); %数据写入
     if i==N
     fprintf(fild, '%s\n',';'); %最后一个数据用;
     else
     fprintf(fild,',\n'); % 其他数据用,
     end
     end
     fclose(fild); % 写完了,关闭文件
    
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25
    • 26
    • 27
    • 28
    • 29
    • 30
    • 31
    • 32
    • 33
    • 34
    • 35
    • 36
    • 37
    • 38
    • 39
    • 40
    • 41
    • 42
    • 43
    • 44
    • 45
    • 46
    • 47
    • 48
    • 49
    • 50
    • 51
    • 52
    • 53
    • 54
    • 55
    • 56
    • 57
    • 58
    • 59
    • 60
    • 61
    • 62
    • 63
    • 64
    • 65
    • 66
    • 67
    • 68
    • 69
    • 70
    • 71
    • 72
    • 73
    • 74
    • 75
    • 76
    • 77
    • 78
    • 79
    • 80
    • 81
    • 82
    • 83
    • 84

    ​ 上面的代码除了生成下面的波形结果外,还将数据写入文件“Idata_65536x15bit.coe” 和 “Qdata_65536x15bit.coe”

    可以在电脑win图标旁边直接搜索这两个文件,默认是在MATLAB文件中的一个文件夹中。

    波形结果

    在这里插入图片描述

    2.FPGA实现

    生成的点频信号signal 以及I,Q的分解,可以用ROM来输入。

    用FPGA实现的关键以及难点是计算相位角phase = atan2(Q, I);

    ROM

    在这里插入图片描述

    在这里插入图片描述

    在这里插入图片描述

    创建两个ROM,分别将Idata_65536x15bit 和 Qdata_65536x15bit 写入。

    数学知识补充:atan和atan2

    参考文章:atan2函数和atan函数 - 知乎 (zhihu.com)

    在MATLAB中,atan和atan2函数都用于计算角度,但它们之间有一些重要的区别。

    1. atan函数:
    • atan函数是计算反正切值的标准函数,其语法为y = atan(x)。
    • atan函数返回的角度范围是[-π/2, π/2],即从-90度到90度之间。
    • atan函数只能接受一个参数,即y = atan(x),其中x是输入的实数值。
    1. atan2函数:
    • atan2函数是计算反正切值的扩展函数,其语法为y = atan2(y, x)。
    • atan2函数返回的角度范围是[-π, π],即从-180度到180度之间。
    • atan2函数可以处理所有四个象限的情况,避免了由于分母为零而导致的错误。
    • atan2函数接受两个参数,即y = atan2(y, x),其中y是输入的虚部,x是输入的实部。

    ​ 主要区别在于atan函数只能处理一个参数且返回值范围是[-π/2, π/2],而atan2函数可以处理两个参数且返回值范围是[-π, π],并且能够处理所有四个象限的情况。在处理复数的相位角度时,通常会使用atan2函数来确保得到正确的结果

    这两个函数的转换关系

    atan函数转换为atan2函数:(图片来自知乎)

    在这里插入图片描述

    • x > 0, 是一、四象限的点,等价atan;
    • x < 0,是二、三象限的点,根据y的范围来确定:
      • y > 0,是第二象限的点,先用atan得到第四象限的弧度制,再加上Π就转换到了第二象限;
      • y < 0.是第三象限的点,先用atan得到第一象限的弧度制,再减去Π就转换到了第二象限;
    • x = 0,根据y的正负来确定
      • y > 0, 值为Π/2
      • y < 0,值为-Π/2
    • 注意,当x=0,y=0时,atan2(0,0)=0 (matlab中这样规定)

    FPGA代码思路

    VIVADO的CORDIC IP核中的Arc Tan ,可以直接计算atan2(自动将atan转换为atan2)。

    CORDIC IP核使用

    参考视频:FPGA IP之CORDIC使用与仿真_哔哩哔哩_bilibili

    参考文章:FPGA数字信号处理(十四)Vivado Cordic IP核计算arctan_fpga arctan-CSDN博客

    FPGA 代码

    DDS_IQ模块,是DDS波形产生模块,读取ROM中存入的数据,生成两个波形 I_data 和 Q_data

    module DDS_IQ
    (
        input 				clk,    	//系统时钟
        input 				rst_n,
    	
    	output signed [15:0] I_data,
    	output signed [15:0] Q_data
    
    );
    
    reg [15:0]r_Fword;		//频率控制字寄存器
    reg [1:0]r_Pword;		//相位控制字寄存器
    reg [31:0] Fcnt;		//累加寄存器
    wire [15:0] I_rom_addr; 	//ROM地址,宽度:16位
    wire [15:0] Q_rom_addr; 	//ROM地址,宽度:16位
    
    //将值存入寄存器
    always@(posedge clk or negedge rst_n)begin
    	if(!rst_n)begin
    		r_Fword <= 16'd0;
    		r_Pword <= 2'd0;
    	end
    	else begin
    		r_Fword <= 16'd1000;	
    		r_Pword <= 2'd0;	
    	end
    end
    
    //累加
    always@(posedge clk or negedge rst_n)begin
    	if(!rst_n)
    		Fcnt <= 32'd0;
    	else
    		Fcnt <= Fcnt + r_Fword;
    end
    
    //相位调制器
    assign I_rom_addr = Fcnt[31:15] + r_Pword;	//截取高位,并加上相位累加器的值
    assign Q_rom_addr = Fcnt[31:15] + r_Pword;
    
    blk_mem_gen_I I_value(
      .clka(clk),    // input wire clka
      .ena(1'b1),      // input wire ena
      .addra(I_rom_addr),  // input wire [15 : 0] addra
      .douta(I_data)  // output wire [15 : 0] douta
    );
    blk_mem_gen_Q Q_value (
      .clka(clk),    // input wire clka
      .ena(1'b1),      // input wire ena
      .addra(Q_rom_addr),  // input wire [15 : 0] addra
      .douta(Q_data)  // output wire [15 : 0] douta
    );
    
    endmodule
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25
    • 26
    • 27
    • 28
    • 29
    • 30
    • 31
    • 32
    • 33
    • 34
    • 35
    • 36
    • 37
    • 38
    • 39
    • 40
    • 41
    • 42
    • 43
    • 44
    • 45
    • 46
    • 47
    • 48
    • 49
    • 50
    • 51
    • 52
    • 53
    • 54

    顶层 atan_top 模块,计算atan(Q/I)

    module atan_top(
    	input 				clk,    	//系统时钟
    	input 				rst_n,
    	input 				[15:0] I_data,
    	input 				[15:0] Q_data,
    
        output signed 		out_valid,   //输出有效信号
        output signed [15:0]theta 		//arctan计算结果
    	
        );
    
    //例化DDS_IQ模块,将两个信号引入	
    DDS_IQ u_DDS_IQ(
    	.clk		(clk),    	
    	.rst_n		(rst_n),
    	.I_data		(I_data),
    	.Q_data		(Q_data)
    );
    
    //输入I和Q,out=arctan(Q/I);
    //tdata端口,虚部Q在前,实部I在后
    cordic_1 u_cordic_1 (
      .aclk(clk),                                        // input wire aclk
      .aresetn(rst_n),                                  // input wire aresetn
      .s_axis_cartesian_tvalid(1'b1),  // input wire s_axis_cartesian_tvalid
      .s_axis_cartesian_tdata({Q_data[15],Q_data[15:1],I_data[15],I_data[15:1]}),    // input wire [31 : 0] s_axis_cartesian_tdata
      .m_axis_dout_tvalid(out_valid),            // output wire m_axis_dout_tvalid
      .m_axis_dout_tdata(theta)              // output wire [15 : 0] m_axis_dout_tdata
    );
    //输入的32位中,[15:0]为实部I,[16:31]为虚部Q
    	
    endmodule
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25
    • 26
    • 27
    • 28
    • 29
    • 30
    • 31
    • 32
    TB文件
    module cordic_tb_top();
    
    //接口声明
      reg clk;
      reg rst_n;
      
      wire signed out_valid;
      wire signed [15:0]theta;
      wire signed[15:0]I_data;
      wire signed[15:0]Q_data;
    
    //initial handle = $fopen("F:/MY_WORK/3U_phase_discrimination/cordic.txt");//打开文件
    
    /*
    //对被测的设计进行例化
    DDS_IQ u_DDS_IQ(
    	.clk		(clk),    	
    	.rst_n		(rst_n),
    	.I_data		(I_data),
    	.Q_data		(Q_data)
    );
    */
    DDS_IQ u_DDS_IQ(
    	.clk		(clk),    	
    	.rst_n		(rst_n),
    	.I_data		(I_data),
    	.Q_data		(Q_data)
    );
    
    atan_top u_atan_top(
    	.clk		(clk),    	
    	.rst_n		(rst_n),
    	.I_data		(I_data),
    	.Q_data		(Q_data),
    	.out_valid	(out_valid),  
    	.theta 		(theta)
    );
    
    
    //产生时钟 50MHZ
    initial clk = 1;
    always #10 clk = ~clk;
    
    //测试激励产生
    initial begin
    	rst_n = 0;
    
    	#200;
    	rst_n = 1;
    
    end
    
    /*
    always@(posedge clk) begin
      if(out_valid)
          $fdisplay(handle,"%b",theta);//写数据
    end
    */
    
    endmodule
    
    • 1
    • 2
    • 3
    • 4
    • 5
    • 6
    • 7
    • 8
    • 9
    • 10
    • 11
    • 12
    • 13
    • 14
    • 15
    • 16
    • 17
    • 18
    • 19
    • 20
    • 21
    • 22
    • 23
    • 24
    • 25
    • 26
    • 27
    • 28
    • 29
    • 30
    • 31
    • 32
    • 33
    • 34
    • 35
    • 36
    • 37
    • 38
    • 39
    • 40
    • 41
    • 42
    • 43
    • 44
    • 45
    • 46
    • 47
    • 48
    • 49
    • 50
    • 51
    • 52
    • 53
    • 54
    • 55
    • 56
    • 57
    • 58
    • 59
    • 60
    仿真结果

    注意这三个波形为模拟波形,有符号数

    在这里插入图片描述

    结果分析

    如何判断自己的仿真结果是正确的?

    就以上面的仿真图片和下面的结果来验证一下

    在这里插入图片描述

    数据1:Q = 0 ; I = -32767 ; atan2 = 25736

    数据2:Q = -25279; I = -20848; atan2 = -18517

    这些数据如何转化为二进制,可以查看CORDIC IP核中的这个界面:

    在这里插入图片描述

    Q格式数据可以用Fix格式数据表示。

    对于有符号数,表示为Fix(1+X+N)_N,X表示整数位数,N表示小数位数。

    但是在这里的结果验证中,不需要转换数据格式就可以


    • 数据1:属于《数学知识补充:atan和atan2》 的第二种情况,x<0, y>=0

    ​ atan2(Q,I) = atan(Q/I) + Π = 0 + Π = Π

    ​ 而根据输出波形,atan2 = 25736,将它转换:25736/(2^13) = 3.1416

    • 数据2:属于《数学知识补充:atan和atan2》 的第三种情况,x<0, y<0

    ​ atan2(Q,I) = atan(Q/I) - Π

    ​ 其中,atan(Q/I) = atan(-25279/-20848) = atan(1.2125) = 0.88115 Radians

    ​ ∴ atan2(Q,I) = 0.88115 - Π = -2.26

    而根据输出波形,atan2 = -18517,将它转换:-18517/(2^13) = -2.26


    仿真结果符合事实

  • 相关阅读:
    H5的基础
    普冉PY32系列(八) GPIO模拟和硬件SPI方式驱动无线收发芯片XN297LBW
    pytest(10)-常用执行参数说明
    大家都能看得懂的源码 - ahooks useSet 和 useMap
    ARM开发初级-STM32F4 USART串口的应用-学习笔记06
    JavaScript DOM基础
    计算机网络——网络安全
    云呐|固定资产管理系统功能包括哪些?
    第一章:最新版零基础学习 PYTHON 教程(第十一节 - Python 语句中的 – Python For 循环)
    E900V22C刷入CoreELEC、挂载云盘
  • 原文地址:https://blog.csdn.net/www_haha__/article/details/136446319