m基于stft方法的跳频信号解跳和GFSK信号解调,并输出跳频图案

1.算法描述

   扩频通信即扩展信号频带范围的一种通信方式，著名的香农定理告诉我们，当信道容量一定时，信号带宽与信噪比可以互换，也就是说理论上完全有可能在极低的信噪比环境下通过增大信号带宽来保证通信质量，这便是扩频通信的理论基础。跳频通信作为扩频通信的一种方式，具有极强的抗干扰、抗截获的能力。在军事领域中，由于载波频率不断跳变使得敌方想要通过侦查搜索到当前的频率位置是十分困难的；在民用领域中，通信频带内具有很强的噪声信号，由于跳频信号跳变速度极快，频率集很大，很大一部分时间可以避开干扰频带继续通信，相比于传统的定频通信，跳频信号受这种强干扰的影响大大降低，极大的提高了通信质量。

   为了减小低信噪比下干扰和噪声对跳频信号检测的影响,提出一种基于时频分析的多跳频信号盲检测算法。针对跳频信号、定频信号、高斯白噪声具有的不同时频分布特点,该算法利用短时傅里叶变换。

   STFT 由于满足叠加原理属于线性时频表示，它只能将信号分解为基本的分量，从本质上来说还是属于傅里叶变换，不能描述信号的瞬时功率谱密度即信号的能量分布，而双线性时频表示且很好的弥补了这个缺点。双线性时频表示也叫做二次型时频表示，常用的二次型时频表示是维格纳-威尔分布（Wigner-Ville Distribution，WVD），它可以看成在某个时刻信号的自相关函数的傅里叶变换，因此可以表示信号能量在时频域中的分布。

   STFT描述非平稳信号的主要思想就是假设在所加的时间窗内的信号是平稳的，对信号在窗内的傅里叶变换后，将窗沿时间轴平移的过程便是 STFT。但是傅里叶变换仅仅是在频域上描述一个信号，而不是表示信号的能量，所以使用功率谱来分析信号。平稳随机信号的功率谱密度是信号自相关函数的傅里叶变换，在实际运用中我们只能采集到平稳随机过程中的有限个数据，需要通过这有限个数据对其功率谱密度进行估计。通常我们使用周期图法来估计信号的功率谱密度。

   跳频信号属于一种时间上有明显特征的一种信号，从单次加窗的功率谱中是无法识别出调频信号的，因此需要以一定的时间间隔为单位，获得若干次功率谱数据，并对若干次功率谱数据的联合检测来识别跳频信号。因此我们以时间进程来对信号进行多次采集和功率谱估计。获取一个信号切片分为两步进行，分别是信号采集和功率谱估计，如此重复获取多个信号切片，以每一个信号切片得到的功率谱数据来代表这个信号切片时间内信号的功率谱数据。

  “高斯频移键控GFSK (Gauss frequency Shift Keying),是在调制之前通过 一个高斯低通滤波器来限制信号的频谱宽度,以减小两个不同频率的载波切换时 的跳变能量,使得在相同的数据传输速率时频道间距可以变得更紧密。它是一种 连续相位频移键控调制技术,起源于FSK(Frequency- shift keying)。

2.仿真效果预览
matlab2022a仿真结果如下：

3.MATLAB核心程序

txBlockFFT  = zeros(frameNumber,Nfft);
%对信号进行短时傅里叶变换，得到跳频图案
for i = 0:frameNumber-1
    start                         = i*Nfft;
    txBlockFFT(i+1,:)             = fftshift(fft(x(start+1:start+Nfft)));
    %找出跳频图案中横轴上每个时间上对应的最大峰值，从时间开始到结束写个循环就行，
    [maxValue(i+1),maxIndex(i+1)] = max(abs(txBlockFFT(i+1,:)));
end
%计算峰值均值
Level = 10*mean(maxValue);
 
delta_f     = fs/Nfft;
f           = delta_f: delta_f: fs;
f           = f - fs/2;
 
%detect vaalid signal
validIndexCount             = 1;
validIndex                  = zeros(1,1);
validIndex(validIndexCount) = 0;
validFrameCount             = 0;
 
figure;
for i = 1:frameNumber-1
    i
    %峰值搜索
    if(max(abs(txBlockFFT(i+1,:))> Level))
        
        validFrameCount             = validFrameCount + 1;
        validFrame(validFrameCount) = i+1;
        if(abs(maxIndex(i+1) - validIndex(validIndexCount)) > 40)
           validIndexCount             = validIndexCount +1;
           validIndex(validIndexCount) = maxIndex(i+1) ;
           detectFHFc                  = validIndex(2:end)*fs/Nfft- fs/2;
        end
        
        detectFHResult(i+1) = maxIndex(i+1)*fs/Nfft- fs/2+fc;
        t                   = (0:length(detectFHResult)-1)*Nfft/fs;
        plot(t,detectFHResult,'c*');
        hold on
    end
    axis([0,0.1983,5.74e9,5.86e9]);
    pause(0.001);
end
hold on;
y_label = fc*ones(1,length(detectFHResult));
plot(t,y_label,'r');
ylim([fc-50e6 fc+50e6]);
xlabel('时间s')
ylabel('频率Hz');
grid on;
 
 
FT = zeros(1000,length(t));
%归一化频率
indx=find(detectFHResult==0);
detectFHResult1=detectFHResult;
detectFHResult1(indx)=[];
 
detectFHResult2 = detectFHResult-min(detectFHResult1);
indx=find(detectFHResult2<0);
detectFHResult2(indx)=0;
detectFHResult3 = floor(1000*detectFHResult2/max(detectFHResult2));
for i = 1:length(t)
    if detectFHResult3(i)>0
       FT(detectFHResult3(i),i) = 1;
    end
end