gfsk调制原理

1. FSK调制原理

我们知道FSK的调制的时域表达式为：

\[s(t)=cos(2\pi(f_c+I_nf_d)t) \tag{1.1} \]

其中\(I_n=\pm 1, \pm3, \cdots, \pm(M-1)\)为数字映射信息，\(f_d\)为载波调制频率偏移量。FSK调制通过\(I_nf_d\)来传达要传输的数字信息，但是在FSK调制连续数字信息时，频率的突发式切换带来了载波相位的不连续性，导致调制信号频谱旁瓣过大，频谱利用率低。所以后来提出了连续相位频移键控(Continuous Phase Frequency Shift Keying, CPFSK)。
为了表示CPFSK信号，以PAM信号：

\[d(t)=\sum_{n}I_{n}g(t-nT) \tag{1.2} \]

开始，式中\(I_{n}\)表示幅度序列，它是由信息序列\({a_{n}}\)来的\(k\)比特二进制数字组映射到幅度电平\(I_{n} = \pm 1, \pm 3, \cdots, \pm(M-1)\)得到的。而\(g(t)\)是一个幅度为 \(1/2T\) 且持续时间为\(T\)秒的矩形脉冲。信号\(d(t)\)用来对载波进行频率调制，调制之后的载波信号可以表示为：

\[s(t)=Acos[2\pi f_{c}t+\phi(t,I)+\phi_{0}] \tag{1.3} \]

其中，\(f_{c}\)表示载波频率，\(\phi_{0}\)表示载波的初始相位，\(\phi(t,I)\)为载波随时间变化所对应的相位，可以表示为：

\[\phi(t;I)=4\pi f_{d}T\int_{-\infty}^{t}d(\tau)d\tau=4\pi f_{d}T\int_{-\infty}^{t}[\sum_{n}I_{n}g(\tau-nT)]d\tau \tag{1.4} \]

式中\(f_{d}\)表示峰值频率偏移，\(d(t)\)信号是不连续的，但是其积分是连续的(注，式中的"\(4\)"单纯是为了后续的计算方便添加的，因为可以把矩形脉冲幅度里面的1/2抵消掉，然后恰好可以把"\(4\pi\)"变成"\(2\pi\)"，使之能够满足角频率与频率的转换关系)。所以载波的相位变化是连续的，调制后载波信号的相位在时间区间\([nT,(n+1)T]\) 可以表示为：

\[\begin{aligned} \phi(t;I)&=2\pi f_{d}T\sum_{k=-\infty}^{n-1}I_{k}+2\pi f_{d}q(t-nT)I_{n} \\ &=\theta_{n}+2\pi h I_{n}q(t-nT) \end{aligned} \tag{1.5} \]

式中，\(h\)表示调制指数，调制指数越大，调制后的两个载波频率差值越大；\(\theta_{n}\)是时间区间\([0,(n-1)T]\)内的累积相位，\(q(t)\)为相位响应函数，是频率成形滤波器\(g_{T}(t)\)的积分，它们可以表示为：

\[h=2f_{d}T \tag{1.6} \]

\[\theta_{n} = \pi h \sum_{k=-\infty}^{n-1} I_{k} \tag{1.7} \]

\[q(t) = \left \{\begin{aligned} &0, t \le 0 \\ &t/2T, 0<t<T \\ &1/2, t>T \end{aligned} \right. \tag{1.8} \]

\[q(t)=\int_{0}^{1}g_{T}(\tau)d\tau \tag{1.9} \]

CPFSK调制是用将要发送的基带数据控制一个单一频率的载波，载波的相位变化是连续的，有效地抑制了频谱旁瓣，调制后信号的带宽变窄。

2. GFSK调制原理

为了进一步降低CPFSK调制信号的频谱旁瓣，首先要对发送的二进制信号映射的矩形脉冲进行高斯低通滤波，然后将滤波后的信号进行载波调制，这种调制技术称为GFSK。GFSK调制信号的包络恒幅、功率谱集中、抗干扰能力强，同时GFSK调制技术实现简单，在低速率设备以及低成本无线传感网络中，应用广泛。

高斯滤波器的带宽影响着GFSK调制信号的频谱，所以高斯滤波器的3dB带宽和二进制码元周期T的乘积BT是高斯滤波器的主要参数，BT取值越小，脉冲响应持续的时间越长，在实际设计中，通常将脉冲响应的持续时间根据需要进行截断，以便数字电路实现。

高斯滤波器的冲激响应表达式为：

\[h(t)=\frac{1}{\sqrt{2\pi}\sigma T}e^{\frac{-t^{2}}{2\sigma^{2}T^{2}}} \tag{2.1} \]

式中

\[\sigma = \sqrt{ln(2)}/2\pi BT \tag{2.2} \]

假设输入信号\(x(t)\)是要发送的二进制数据所对应的双极性方波，双极性方波经过高斯滤波器后的输出为：

\[g(t)=h(t)\otimes m(t;T) \tag{2.3} \]

上式中，矩形脉冲定义为：

\[m(t;T)= \left \{ \begin{aligned} & 1/T, |t|<t/2 \\ & 0, \text{otherwise} \end{aligned} \right. \tag{2.4} \]

经过推导得出表达式：

\[g(t)=\frac{1}{2T}[Q(2\pi BT\frac{t-T/2}{T\sqrt{ln(2)}})-Q(2\pi BT\frac{t+T/2}{T\sqrt{ln(2)}})] \tag{2.5} \]

式中，Q函数为：

\[Q(t) = \frac{1}{\sqrt{2\pi}} \int_{t}^{\infty}e^{-\frac{\tau^{2}}{2}}d\tau \tag{2.6} \]

假设序列x(n)是对输入信号x(t)的离散时间采样，是双极性不归零码，则GFSK调制后信号的相位表达式为：

\[\theta(t) = 2\pi h \int_{-\infty}^{t}\sum_{n=-\infty}^{\infty}x[n]g(\tau-nT)d\tau \tag{2.7} \]

GFSK调制信号在载波频率为\(f_{c}\)时的表达式为：

\[s(t) = \sqrt{\frac{2E_{b}}{T}}cos[2\pi f_{c}t+\theta(t)+\phi_{0}] \tag{2.8} \]

总结：其实GFSK是在CPFSK的基础上对矩形脉冲进行了高斯滤波。从最后的结果来看其实就是将卷积矩形脉冲修改为卷积“矩形脉冲经过高斯滤波器之后的脉冲”。

3. GFSK的MATLAB代码实现

首先实现“经过高斯滤波的矩形脉冲”，也即是式\((2.5)\)中的\(g(t)\)。

function [Gaussian_filter_q] = Gaussian_filter(B, symb_rate, delay, ts)
% INPUTS:
%   B         : the 3dB bandwidth of Gaussian filter
%   symb_rate : rate of the symbol, in ble it is 1e6 symbol/s
%   delay     : it controls how many symbols the Gaussian filter cover
%   ts        : the reciprocal of sample rate
  symb_T  = 1 / symb_rate;
  Gau_t   = -delay*symb_T:ts:delay*symb_T;
  K       = 2*pi*B / sqrt(log(2));
  Gau_det = K.*(Gau_t - symb_T / 2);
  Gau_sum = K.*(Gau_t + symb_T / 2);
  Gaussian_filter_q = 1/(2*symb_T)*[qfunc(Gau_det)-qfunc(Gau_sum)];
end

然后设定参数进行测试：

clear
clc

h = 0.5 ; % the modulation index h = 2f_{d}T
sample_rate = 26e6; % the sample rate
sample_time = 1/sample_rate;
symbol_rate = 1e6;
oversample = sample_rate / symbol_rate;
B = 0.5e6; % the 3db bandwidth is 500Mhz
delay = 1.5 % the gaussian filter covers 3 symbol
tx_bits = [1,0,0,1,0,0,1,0,1,1,1,0,1,1,0,0,1]; % the transmit bits

Gaussian_filter_coeff = Gaussian_filter(B, symbol_rate, delay, sample_time);
tx_data = 2*tx_bits - 1; % NRZ code
tx_data_sample = upsample(tx_data, oversample);
tx_data_filtered = 2*pi*h*conv(Gaussian_filter_coeff,tx_data_sample);

tx_phase = sample_time*cumsum(tx_data_filterd);
tx_cos = cos(tx_phase);
tx_sin = sin(tx_phase);

tx_gfsk = complex(tx_cos,tx_sin);

plot(tx_cos); 
hold on 
plot(tx_sin);
hold off

最终波形如下：

The gfsk wave.