基于插值算法和Gardner定时误差检测的OOK信号定时同步的FPGA实现

  本文介绍如何用FPGA实现基于插值算法的OOK信号定时同步，Verilog代码参考杜勇《数字调制解调技术的MATLAB与FPGA实现》。我们的目标是用外部提供50MHz时钟的zynq7100芯片实现400MHz采样频率和100Mbps的OOK数字基带信号的定时同步。
  采用传统的锁相环技术实现定时同步时，本地时钟需要有较高的频率。当数据采样频率很高，并且本地时钟受到器件性能限制而不能远高于采样频率时，锁相环技术性能不佳。插值算法可以不改变采样时钟的频率和相位来实现位同步信号的调整，同时，插值算法可以根据采样值以及数控振荡器输出的采样时刻信号和误差信号获取最佳采样值。
  插值位同步算法的框图如下图所示（图片源于文献[1]）。主要模块为括插值滤波器（INTERPOLATOR）、定时误差检测器（TIMING ERROR DETECTOR）、环路滤波器（LOOP FILTER）和数控振荡器（CONTROLLER）。下面我们一一介绍并给出Verilog代码。

插值滤波器

  插值滤波器的功能是速率转换，上图显示了内插滤波器的原理。设采样周期为$T_s$，符号周期为$T$，插值周期为$T_i$，插值滤波器的脉冲响应为$h_I(t)$。由上图可知，插值滤波器输出的连续时间信号为 $$y(t)=\sum _{m}x(m)h_I(t-m{T}_s)$$ 按插值周期对$y(t)$进行采样得 $$y(k{T}_i)=\sum _{m}x(m)h_I(k{T}_i-m{T}_s)\text{\hspace{1em}(1)}$$ 上式中，$m$是采样信号的索引，下面定义一个滤波器的索引 $$i=int[\frac{k{T}_i}{T_s}]-m$$ 基准点索引 $$m_k=int[\frac{k{T}_i}{T_s}]$$ 分数间隔 $${\mu }_k=\frac{k{T}_i}{T_s}-m_k$$ 公式$(1)$可以重写为 $$y(k{T}_i)=y[(m_k+{\mu }_k)T_s]=\sum _{i}x[(m_k-i)T_s]h_I[(i+{\mu }_k)T_s]$$ 其中$k{T}_i$表示第$k$个插值点，$m_k{T}_s$表示第$k$个插值点前相邻的采样点， ${\mu }_k{T}_s$是插值点$k{T}_i$与采样点$m_k{T}_s$之间的时间间隔。下图展示了采样点与插值点之间的关系。  根据文献[2]，我们采用Farrow结构的插值滤波器，如下图所示。  插值滤波器有三条纵向支路和一条横向支路，它们的计算公式为 $$\begin{array}{rl}& f_1=0.5x(m)-0.5x(m-1)-0.5x(m-2)+0.5x(m-3)\\ & f_2=-0.5x(m)+1.5x(m-1)-0.5x(m-2)-0.5x(m-3)\\ & f_3=x(m-2)\\ & y(k)=f_1{\mu }_k^2+f_2{\mu }_k+f_3\end{array}$$
  下面给出插值滤波器的Verilog代码。其中，mult18_16模块是vivado的乘法器IP核，设置了2级流水线；$f_1、f_2、f_3$的放在always内用时序逻辑描述（书上用的是组合逻辑），否则时序仿真时会出现很多毛刺造成性能下降。

`timescale 1ns / 1ps
module interpolate_filter(
    input resetn,
    input clk,
    input signed [14:0] data_in,   // 采样数据，采样频率同clk频率
    input signed [15:0] uk,        //分数间隔
    output signed [17:0] data_out  //插值滤波输出，输出速率同clk频率
    );
reg  signed [14:0] din_1,din_2,din_3,din_4,din_5,din_6;
reg  signed [15:0] u_1,u_2;
wire signed [33:0] f1_u,f2_u;
wire signed [33:0] f1_u2;
reg  signed [33:0] f2_u_1,f2_u_2;    
reg signed [17:0] f1;
reg signed [17:0] f2;
reg signed [17:0] f3;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        din_1 <= 15'd0;
        din_2 <= 15'd0;
        din_3 <= 15'd0;
        din_4 <= 15'd0;
        din_5 <= 15'd0;
        din_6 <= 15'd0;
        u_1 <= 16'd0;
        u_2 <= 16'd0;
        f2_u_1 <= 32'd0;
        f2_u_2 <= 32'd0;
        f1 <= 18'd0;
        f2 <= 18'd0;
        f3 <= 18'd0;
    end
    else begin
        din_1 <= data_in;
        din_2 <= din_1;
        din_3 <= din_2;
        din_4 <= din_3;
        din_5 <= din_4;
        din_6 <= din_5;
        u_1 <= uk;
        u_2 <= u_1;
        f2_u_1 <= f2_u;
        f2_u_2 <= f2_u_1;
        f1 <= {{4{data_in[14]}},data_in[14:1]}-{{4{din_1[14]}},din_1[14:1]}-{{4{din_2[14]}},din_2[14:1]}+{{4{din_3[14]}},din_3[14:1]};
        f2 <= {{3{din_1[14]}},din_1}+{{4{din_1[14]}},din_1[14:1]}-{{4{data_in[14]}},data_in[14:1]}-{{4{din_2[14]}},din_2[14:1]}-{{4{din_3[14]}},din_3[14:1]};
        f3 <= {{3{din_6[14]}},din_6}; //一个乘法器延时2个clk
    end
end
//乘法器2级流水线
mult18_16 u1 (
  .CLK(clk),  // input wire CLK
  .A(f1),      // input wire [17 : 0] A
  .B(uk),      // input wire [15 : 0] B
  .P(f1_u)      // output wire [33 : 0] P
);
//乘法器2级流水线
mult18_16 u2 (
  .CLK(clk),  // input wire CLK
  .A(f2),      // input wire [17 : 0] A
  .B(uk),      // input wire [15 : 0] B
  .P(f2_u)      // output wire [33 : 0] P
);
//乘法器2级流水线
mult18_16 u3 (
  .CLK(clk),  // input wire CLK
  .A(f1_u[32:15]),  // input wire [17 : 0] A //matlab仿真得f1_u有效位宽不超过31bit
  .B(u_2),      // input wire [15 : 0] B
  .P(f1_u2)      // output wire [33 : 0] P
);
wire signed [18:0] dt;
assign dt = (!resetn)? 19'd0:(f2_u_2[33:15]+f1_u2[33:15]+{f3,1'b0});//小数位对齐
reg signed [17:0] dt_1;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        dt_1 <= 18'd0; 
    end
    else begin
        dt_1 <= dt[17:0];
    end
end
assign data_out = dt_1;
endmodule

Gardner定时误差检测算法

  Gardner定时误差检测算法根据插值滤波器输出的插值数据进行定时误差检测。一个符号周期内只需要两个插值点，一个插值点出现在数据的峰值时刻，另一个插值点出现在两个数据峰值的中间时刻。
  Gardner定时误差检测算法的公式为 $${\mu }_t(k)=I(k-\frac{1}{2})[I(k)-I(k-1)]$$ 其中，$I(k)$表示第$k$个码元数据选通时的插值，$I(k-\frac{1}{2})$表示第$k$和$k-1$个码元中间时刻的插值。我们还可以用插值信号的正负号代替插值的实际值，虽然有一定的性能损失，但是提高了追踪能力，而且系统实现简单，公式为μt(k)=I(k−12){sgn[I(k)]−sgn[I(k−1)]}\mu_t(k)=I(k-\frac{1}{2})\{sgn[I(k)]-sgn[I(k-1)]\}μt​(k)=I(k−21​){sgn[I(k)]−sgn[I(k−1)]}
  Gardner定时误差检测和环路滤波器的代码写在一个模块内，后面一起给出代码。

环路滤波器和数字振荡器

环路滤波器

  采用二阶环路滤波器，滤波器框图如下图所示。
  环路滤波器的传递函数为 $$H(z)=C_1+\frac{C_2{z}^{-1}}{1-z^{-1}}$$ $C_1、C_2$可以下式计算 $$\begin{array}{rl}& C_1=\frac{8B_L{T}_s}{3}\\ & C_2=\frac{32(B_L{T}_s{)}^2}{9}\end{array}$$其中，$B_L{T}_s$为单边噪声带宽与采样周期的乘积，通常要求$B_L{T}_s\ll 0.1$。为了在FPGA中用移位运算代替乘除运算并简化设计，我们取$C_1=2^{-8}，C_2=0$。
  下面给出定时误差检测和环路滤波器的Verilog代码。环路滤波器输出的误差信号量化为16bit定点数，其中包括1bit符号位和15bit的小数位。

module ted_loop_filter #(
    parameter SPS_2 = 2   // 上采样率（采样速率与数据速率之比）的一半，
    )(
    input resetn,
    input clk,
    input strobe,
    input signed [17:0] data_in,   //插值滤波器输出得插值数据
    output signed [17:0] data_out, //最佳采样时刻得插值数据，用于判决0、1
    output signed [15:0] wk,       //环路滤波器输出定时误差信号，15 bit小数位
    output sync             //位同步信号
    );
reg [3:0] strobe_cnt;
reg sk;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        strobe_cnt <= 4'd0;
        sk <= 1'b0;
    end
    else begin
        if(strobe) begin
            strobe_cnt <= (strobe_cnt >= SPS_2-1'b1)?4'd0:(strobe_cnt+1'b1);
            sk <= (strobe_cnt >= SPS_2-1'b1)?(1'b1):1'b0;  // sk翻转周期位符号周期，作为位定时时钟输出
        end
    end
end
assign sync = sk;
reg signed [17:0] din_1,din_2,din_3;
reg signed [17:0] dout;
reg signed [17:0] err,err_1;
reg signed [15:0] w;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        din_1 <= 18'd0;
        din_2 <= 18'd0;
        din_3 <= 18'd0;
        dout <= 18'd0;
        err <= 18'd0;
        err_1 <= 18'd0;
        w <= 16'b0100000000000000;
    end
    else begin
        if((strobe_cnt==0 || strobe_cnt==SPS_2-1) && strobe) begin
            din_1 <= data_in;
            din_2 <= din_1;
            din_3 <= din_2;
            err <= (!din_1[17] && din_3[17])?{din_2[17:1],1'b0}:((din_1[17] && !din_3[17])?(-{din_2[17:1],1'b0}):18'd0);
            if(sk) begin
                dout <= din_1;
                err_1 <= err;
                //w(ms+1)=w(ms)+c1*(err(ms)-err(ms-1))+c2*err(ms), c1 = 2^(-8)， c2≈0
                w <= w+{{6{err[17]}},err[17:8]}-{{6{err_1[17]}},err_1[17:8]};
            end
        end
    end
end
assign wk = w;
assign data_out = dout;  
endmodule

数控振荡器(NCO)

  数控振荡器是一个相位递减器，它的差分方程为 $$\eta (m+1)=[\eta (m)-w(m)]mod1$$ 其中$\eta (m+1)$是nco寄存器的值，$w(m)$是环路滤波器输出的定时误差值。

  nco的原理如上图所示。根据相似三角形的性质可以得出下式$$\frac{{\mu }_k{T}_s}{\eta (m_k)}=\frac{(1-{\mu }_k)T_s}{1-\eta (m_k+1)}$$ 分数间隔为$${\mu }_k=\frac{\eta (m_k)}{1-\eta (m_k+1)+\eta (m_k)}=\frac{\eta (m_k)}{w(m_k)}\text{\hspace{1em}(2)}$$
  当环路稳定时，nco每隔$T_s$时间就递减$w(m)$，即平均nco每$\frac{1}{w(m)}$个时钟就下溢一次，那么插值周期与采样周期有近似关系$$T_i=\frac{T_s}{w(m)}$$Gardner定时误差算法每个符号内有两个插值点参与运算，故$T\approx 2T_i$。Farrow插值滤波器每个插值需要4个采样点，故$T\approx 4T_s$。所以$w(m)\approx 0.5$。公式$(2)$可以改写为$${\mu }_k\approx 2\eta (m_k)$$或者一阶修正公式$${\mu }_k\approx 2\eta (m_k)[2-2w(m_k-1)]$$
  下面给出定数控振荡器的Verilog代码。理论分析中的小数（分数）都采用16bit定点数量化。理论上，当nco寄存器的值下溢时进行一次插值。在FPGA中，插值是在不断进行的，但是只有nco寄存器的值下溢时获得的插值有用。

module nco(
    input resetn,
    input clk,
    input signed [15:0] wk,    //环路滤波器输出定时误差信号，15 bit小数位
    output signed [15:0] uk,   //NCO输出的插值间隔小数，15 bit小数位
    output strobe       //NCO输出的插值计算选通信号，高电平有效
    );
reg signed [16:0] nkt;
reg signed [16:0] ut;
reg str;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        nkt <= 17'b00110000000000000;   //0.75
        ut <= 17'b00100000000000000;    //0.5
        str <= 1'b0;
    end
    else begin
        if(nkt < {wk[15],wk}) begin // 负值+1，相当于mod(1);
            nkt <= 17'b01000000000000000+nkt-{wk[15],wk};
            ut <= {nkt[14:0],1'b0}; //取出nkt减去wk之前的值，乘以2作为u值输出
            str <= 1'b1;
        end
        else begin
            nkt <= nkt-{wk[15],wk};
            str <= 1'b0;
        end
    end
end
assign uk = ut;
assign strobe = str;
endmodule

定时同步顶层文件

  将前面几个模块整合起来，获得完整的定时同步模块。

(* dont_touch = "yes" *) module bit_sync #(
    parameter SPS_2 = 2   // 上采样率（采样速率与数据速率之比）的一半，
    )
    (
    input resetn,
    input clk,
    input signed [14:0] data_in,  //采样数据
    output data_out,    //位同步后0、1bit
    output sync         //位同步脉冲
    );
(* dont_touch = "yes" *) reg signed [14:0] din;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        din <= 15'd0;
    end
    else begin
        din <= data_in;
    end
end
(* dont_touch = "yes" *) wire signed [15:0] uk,wk;
(* dont_touch = "yes" *) wire signed [17:0] data_interpolate;
(* dont_touch = "yes" *) wire strobe;
(* dont_touch = "yes" *) wire [17:0] dout;
(* dont_touch = "yes" *) wire bit_sync;
(* dont_touch = "yes" *) interpolate_filter u1(
    .resetn(resetn),
    .clk(clk),
    .data_in(din),   // 采样数据，采样频率同clk频率
    .uk(uk),           //分数间隔
    .data_out(data_interpolate)      //插值滤波输出，输出速率同clk频率
    );
(* dont_touch = "yes" *) ted_loop_filter #(
    .SPS_2(SPS_2)   // 上采样率（采样速率与数据速率之比）的一半，
    )
    u2 (
    .resetn(resetn),
    .clk(clk),
    .strobe(strobe),
    .data_in(data_interpolate),   //插值滤波器输出得插值数据
    .data_out(dout), //最佳采样时刻得插值数据，用于判决0、1
    .wk(wk),       //环路滤波器输出定时误差信号，15 bit小数位
    .sync(bit_sync)      //位同步信号
    );
(* dont_touch = "yes" *) nco u3(
    .resetn(resetn),
    .clk(clk),
    .wk(wk),    //环路滤波器输出定时误差信号，15 bit小数位
    .uk(uk),   //NCO输出的插值间隔小数，15 bit小数位
    .strobe(strobe)       //NCO输出的插值计算选通信号，高电平有效
    );
(* dont_touch = "yes" *) reg bit_dout;
(* dont_touch = "yes" *) reg bit_sync_1;
(* dont_touch = "yes" *) reg bit_sync_2;
always @(posedge clk or negedge resetn) begin
    if(!resetn) begin
        bit_dout <= 1'b0;
        bit_sync_1 <= 1'b0;
        bit_sync_2 <= 1'b0;
    end
    else begin
        bit_sync_1 <= bit_sync;
        if(bit_sync & (!bit_sync_1)) begin
            bit_dout <= !dout[17];
            bit_sync_2 <= 1'b1;
        end
        else begin
            bit_dout <= bit_dout;
            bit_sync_2 <= 1'b0;
        end
    end
end
assign data_out = bit_dout;
assign sync = bit_sync_2; //上升沿
endmodule

  我们要实现400MHz采样频率，同步100Mbps的OOK信号，而FPGA的外部晶振只有50MHz，所以还要再放一个PLL IP核将时钟倍频至400MHz。将PLL IP核和位同步模块封装在一起。

module test_top(
    input resetn,
    input clk,       // 50MHz
    input signed [14:0] data_in,  //采样数据
    output clk400M,
    output data_out,    //位同步后0、1bit
    output sync         //位同步脉冲
    );
wire locked;
clk_wiz_0 u1
   (
    // Clock out ports
    .clk_out1(clk400M),     // output clk_out1
    // Status and control signals
    .reset(~resetn), // input reset
    .locked(locked),       // output locked
   // Clock in ports
    .clk_in1(clk));      // input clk_in1
bit_sync #(
    .SPS_2(2)   // 上采样率（采样速率与数据速率之比）的一半，
    )
    u2 (
    .resetn(resetn&locked),
    .clk(clk400M),
    .data_in(data_in),  //采样数据
    .data_out(data_out),    //位同步后0、1bit
    .sync(sync)         //位同步脉冲
    );

endmodule

时序仿真

  创建vivado工程，选择zynq7100芯片，分析，综合，编写testbench文件，然后进行综合后时序仿真。
testbench文件如下。adc_data.txt是用MATLAB生成的存放波形数据的文件，定时同步结果输出到bit_sync.txt用于计算误码率。

module tb_test_top();
parameter LEN = 4915200;
reg resetn;
reg clk;
reg [14:0] data [LEN-1:0];
reg [14:0] data_in;
wire clk400M;
wire data_out;
wire sync;
integer k;
integer file;
initial begin
    resetn = 1'b0;
    clk = 1'b1;
    data_in <= 15'd0;
    k = 0;
    file = $fopen("bit_sync.txt");
    $readmemh("adc_data.txt",data);
    #100
    resetn = 1'b1;
    #1600
    for(k = 0; k < LEN; k = k+1) begin
        #2.5 //400Mhz
        data_in <= data[k];
    end
    #2.5
    data_in <= 15'd0;
    #100 $finish;
end
always #10 clk = ~clk;  // 50Mhz
always @(posedge sync)
	   $fdisplay(file,"%d",data_out);
test_top u1 (
    .resetn(resetn),
    .clk(clk),
    .data_in(data_in),  //采样数据
    .clk400M(clk400M),
    .data_out(data_out),    //位同步后0、1bit
    .sync(sync)         //位同步脉冲
    );
endmodule

  等一会儿，时序仿真很慢，然后能看到仿真波形。

参考文献

[1] Gardner F M. Interpolation in Digital Modems-Part I: Fundamentals[J]. IEEE Trans. Commun, 1993, 41(3):501-507.

[2] Erup L, Gardner F M, Harris R A. Interpolation in digital modems-Part II : Implementation and performance[J]. IEEE Trans. Commun, 1993, 41(6):998 - 1008.

[3] Gardner F M. A BPSK/QPSK Timing-Error Detector for Sampled Receivers[J]. IEEE Trans. Commun, 1986, 34(5):423-429.

[4] 杜勇. 数字调制解调技术的MATLAB与FPGA实现[M]. 北京:电子工业出版社, 2014.