一种基于cuda的多体制信号频率快速捕获实现方法
【专利摘要】一种基于CUDA的多体制信号频率快速捕获实现方法，其特点在于：基于高性能PC，利用CUDA架构，以并行的方式实现多体制信号频率的快速捕获；实现步骤包括：设置多普勒及多普勒变化率分段，确定各GPU内核函数的线程块数目和每个线程块中的线程数目，将信号数据读入GPU卡全局存储空间，启动下变频和积分模块，进行并行FFT运算，以并行归约方式查找最大值谱线及对应的多普勒变化率分段和多普勒分段，计算载波频偏和多普勒变化率估值完成捕获。本发明使信号频率快速捕获算法符合CUDA架构和GPU工作特点，能够充分发挥GPU并行计算能力，同时适合多种体制的信号，有效提高捕获速度且易于实现。
【专利说明】—种基于CUDA的多体制信号频率快速捕获实现方法
【技术领域】
[0001]本发明涉及一种信号频率快速捕获实现方法，特别是一种基于CUDA的多体制信号频率快速捕获方法。
【背景技术】
[0002]信号的快速捕获是测控通信的关键技术之一，随着深空探测技术的发展，信号捕获技术的难度和复杂度大大增加，当信号极其微弱且多普勒动态范围大、变化率高时，信号捕获一般采用二维FFT捕获方法，即首先对多普勒分段，再在每个多普勒频点对多普勒变化率分段，并在每个分段上进行相关、积分和FFT运算，查找最大谱线，这使得信号捕获过程的运算量呈级数增长。现有的研究中一般采用高性能FPGA来实现这种二维FFT信号频率快速捕获算法，这种方式需要开发专门的硬件板卡，实现复杂，且FPGA的资源有限，一般一次只能实现一种体制信号的捕获，通用性差、可扩展性差。
[0003]近年来，随着商用高性能计算机系统的普及，基于现代高性能PC或服务器平台的信号处理技术成为新的研究热点，并且，统一计算设备架构(Computer Unified DeviceArchitecture, CUDA)开启了使用GPU强大计算能力做通用计算的大门，GPU具有强大的并行计算能力，如果能够通过CUDA架构在GPU上实现复杂的二维FFT信号快速捕获算法，就能够大大提高信号捕获的速度，同时适应多种体制，实现手段简单，可扩展性强。
【发明内容】
[0004]本发明的技术解决问题是:克服现有技术的不足之处，提供了种基于CUDA的多体制信号频率快速捕获方法，以高性能PC为平台，基于CUDA架构，以并行的方式实现多体制信号频率的快速捕获；本发明使信号频率快速捕获算法符合CUDA架构和GPU工作特点，能够充分发挥GPU并行计算能力，同时适合多种体制的信号，有效提高了捕获速度且易于实现。
[0005]本发明的技术解决方案是:一种基于CUDA的多体制信号频率快速捕获方法，步骤如下:
[0006](I) CPU通过CUDA架构选择GPU协处理器并初始化CPU和GPU的运行环境，同时对运行参数进行初始化设置，所述运行参数包括:信号数据文件路径、信号码速率、信号采样频率、信号载波频率初始值、信号多普勒范围、信号多普勒变化率范围、信号的FFT分析点数和信号体制；
[0007](2) CPU根据步骤(I)中设置的信号数据文件路径将信号数据从数据文件读入CPU内存中；
[0008](3) CPU根据步骤(I)中设置的信号码速率确定信号的FFT采样频率和信号的FFT分析带宽，然后根据信号的FFT采样频率计算对步骤(I)中设置的信号多普勒范围进行分段的步进值，进而确定多普勒分段数，根据信号的FFT分析带宽计算对步骤(I)中设置的信号多普勒变化率范围进行分段的步进值，进而确定多普勒变化率分段数，最后根据步骤(I)中设置的信号采样频率和本步骤获得的信号FFT采样频率计算信号的积分长度；
[0009](4) CPU分配GPU全局存储空间，并初始化已分配的GPU全局存储空间；
[0010](5)CPU根据步骤(3)中获得的信号的积分长度设置GPU中在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目；
[0011](6)CPU根据步骤(3)中获得的多普勒变化率分段数和步骤(I)中设置的信号FFT分析点数设置并行一维FFT运算的参数，所述并行一维FFT运算的参数包括并行的一维FFT运算的个数和每个一维FFT运算的FFT分析点数；
[0012](7)CPU根据步骤(3)获得的多普勒变化率分段数和步骤(I)中设置的信号FFT分析点数设置GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，所述二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组，对所有多普勒分段进行二维最大值查找的结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数；
[0013](S)CPU根据步骤(3)获得的多普勒变化率分段数和多普勒分段数设置GPU对步骤(7)中所述的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，所述二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组；
[0014](9) CPU根据步骤(3)获得的多普勒分段数设置GPU对步骤⑶所述的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目；
[0015](1)CPU将步骤⑵中读入到CPU内存中的信号数据从CPU内存拷贝至步骤(4)确定的GPU全局存储空间中；
[0016](Il)GPU利用步骤(5)中确定的在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，以并行的方式完成对步骤(10)读入到GPU全局存储空间中的信号数据的下变频处理，然后以归约算法进行并行求和，并对并行求和后的结果进行鉴相处理，将鉴相处理结果写回步骤(4)确定的GPU全局存储空间中，所述每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数；
[0017](12)GPU利用步骤(6)中确定的并行一维FFT运算的参数，对步骤(11)中的处理结果进行并行一维FFT运算，并将结果写回步骤(4)确定的GPU全局存储空间，所述并行一维FFT运算结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数；
[0018](13)GPU根据步骤(7)中确定的对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(12)的并行一维FFT运算结果分别进行取模获得FFT模值，并分别对FFT模值以归约算法并行查找最大值和最大值对应的FFT点号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组；
[0019](14)GPU循环进行步骤(11)至步骤(13)，循环次数为多普勒分段数，获得对所有多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果，所述结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数；
[0020](15)GPU根据步骤(8)中确定的对步骤(14)获得的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(14)获得的二维数组分别以归约算法并行查找FFT模值在每一个多普勒分段下所有多普勒变化率范围内的最大值和对应的多普勒变化率分段号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述对步骤(14)获得的二维数组进行二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组；
[0021](16)GPU根据步骤(9)中确定的对步骤(15)获得的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程并行的方式完成对步骤(15)获得的一维数组以归约算法并行查找FFT模值在所有多普勒分段中的最大值和对应的最大多普勒分段号，并将结果写回步骤(4)确定的GPU全局存储空间；
[0022](17)CPU将步骤(16)中获得的FFT模值在所有多普勒分段中的最大值及其对应的多普勒分段号、所述FFT模值在所有多普勒分段中的最大值在步骤(14)获得的对应的FFT点号和在步骤(15)获得的对应的多普勒变化率分段号从GPU全局存储空间中拷贝回CPU内存，并计算载波频偏和多普勒变化率。
[0023]所述信号包括BPSK信号和标准TT&C信号。
[0024]所述步骤(3)中，CPU根据步骤⑴中设置的信号码速率确定信号的FFT采样频率和信号的FFT分析带宽，然后根据信号的FFT采样频率计算对步骤(I)中设置的信号多普勒范围进行分段的步进值，进而确定多普勒分段数，根据信号的FFT分析带宽计算对步骤
(I)中设置的信号多普勒变化率范围进行分段的步进值，进而确定多普勒变化率分段数，最后根据步骤(I)中设置的信号采样频率和本步骤获得的信号FFT采样频率计算信号的积分长度，具体为:
[0025]当信号体制为二相编码体制，即BPSK时:
[0026]Fs = CodeRate ；
[0027]RBff = Fs/FFTSize ；
[0028]PerDoppler = Fs/4 ；
[0029]PerDopplerRate = 3 X RBff X RBff ；
[0030]Seg = 2 X Doppler/PerDoppler+1 ；
[0031]Segl = 2 X DopplerRate/PerDoppleRate+1 ；
[0032]SumSize = SampleRate/Fs+0.5 ；
[0033]当信号体制为统一测控体制，即标准TT&C时，RBff根据信号的情况和经验选取:
[0034]Fs = RBffXFFTSize ；
[0035]PerDoppler = Fs/4 ；
[0036]PerDopplerRate = 3 X RBff X RBff ；
[0037]Seg = 2 X Doppler/PerDoppler+1 ；
[0038]Segl = 2 X DopplerRate/PerDoppleRate+1 ；
[0039]SumSize = SampleRate/Fs+0.5 ；
[0040]式中，CodeRate为信号码速率，Fs为FFT采样频率，RBff为信号的FFT分析带宽，Doppler为信号多普勒范围，DopplerRate为多普勒变化率范围，PerDoppler为多普勒范围分段的步进值，Seg为多普勒分段数，PerDopplerRate为多普勒变化率范围分段的步进值，Segl为多普勒变化率分段数，FFTSize为信号的FFT分析点数，SumSize为积分长度，SampleRate为信号采样频率。
[0041]所述步骤(5)中，CPU根据步骤(3)中获得的信号的积分长度设置GPU中在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，具体为:
[0042](I)每个线程块的线程数目由信号数据的积分长度决定，当积分长度小于等于1024时，线程数目等于积分长度，当积分长度大于1024时，线程数目等于1024 ；
[0043](2)令线程块数目为(X，y，z)，其中X为多普勒变化率的分段数，y为FFT运算的分析点数，z值为I。
[0044]所述步骤(7)中，CPU根据步骤(3)获得的多普勒变化率分段数和步骤⑴中设置的信号FFT分析点数设置GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0045](I)每个线程块中的线程数目由四分之一 FFT分析点数的值决定，当四分之一 FFT分析点数的值小于等于1024时，线程数目等于四分之一FFT分析点数的值，当四分之一FFT分析点数的值大于1024时，线程数目等于1024 ；
[0046](2)线程块数目等于多普勒变化率的分段数。
[0047]所述步骤(8)中，CPU根据步骤(3)获得的多普勒变化率分段数和多普勒分段数设置GPU对步骤(7)中所述的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0048](I)每个线程块中的线程数目由多普勒变化率分段数决定，当多普勒变化率分段数小于等于1024时，线程数目等于多普勒变化率分段数，当多普勒变化率分段数大于1024时，线程数目等于1024;
[0049](2)线程块数目等于多普勒的分段数。
[0050]所述步骤(9)中，CPU根据步骤(3)获得的多普勒分段数设置GPU对步骤(8)所述的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0051](I)每个线程块的线程数目由多普勒分段数决定，当多普勒分段数小于等于1024时，线程数目等于多普勒分段数，当多普勒分段数大于1024时，线程数目等于1024 ；
[0052](2)线程块数目等于I。
[0053]所述步骤(17)中计算载波频偏和多普勒变化率，具体为:
[0054]载波频偏由公式:
[0055]EstCarrFreq = (MaxDopplerID-Seg/2)XPerDoppIer+MaxFFTIDX RBff/2 ；
[0056]给出，式中，MaxDopplerID为FFT模值在所有多普勒分段中最大值对应的多普勒分段号，MaxFFTID为FFT模值在所有多普勒分段中最大值对应的FFT点号；
[0057]多普勒变化率由公式:
[0058]EstDopplerRate = (MaxDopplerRateID-Segl/2)XPerDopplerRate
[0059]给出，式中MaxDopplerRateID为FFT模值在所有多普勒分段中最大值对应的多普勒变化率分段号。
[0060]本发明与现有技术相比的有益效果是:
[0061](I)易于实现:基于FPGA的信号频率快速捕获方法需要开发专门的硬件板卡和固件程序，实现难度大，实现周期长，且对硬件的调试工具有限，手段有限，而本发明基于CPU+GPU的开发平台易于获取，且利用CUDA架构，可以在C语言环境下完成对GPU的调用，代码编写容易，调试过程简单。
[0062](2)运行速度快:本发明将信号频率快速捕获算法拆分为多个GPU处理过程，每个处理过程可以通过线程块和线程两层并行进行加速，并且本发明对信号快速捕获算法的每个处理过程进行了合理的并行优化设计，给出了每个处理过程线程块数目和线程数目的确定方法，充分发挥了 GPU的并行计算能力，保证了本发明方法的运行速度。
[0063](3)适用范围广，基于CPU+GPU的开发平台资源充足，本发明通过对信号频率快速捕获算法的处理过程进行结构优化，将BPSK体制和标准TT&C体制信号的频率快速捕获过程进行了融合，并且本发明针对此两种体制分别给出了运行参数的确定方法，使本发明方法能够同时满足两种体制信号频率快速捕获的要求。
[0064](4)可扩展性强:基于FPGA的信号频率快速捕获方法因为受到硬件自身特点的约束，资源受到限制，当信号的多普勒范围和多普勒变化率范围增大，多普勒分段数和多普勒变化率分段数增多，信号的FFT分析点数增大时都可能导致FPGA的资源不够，从而需要更换FPGA器件并且重新设计硬件板卡，可扩展性差，而本发明采用的CPU+GPU的开发平台和CUDA架构均具有很强的可扩展性，当信号的多普勒范围和多普勒变化率范围增大，多普勒分段数和多普勒变化率分段数增多，信号的FFT分析点数增大时只需要根据本发明中给出的各处理过程线程块数目和线程数目的确定方法重新确定各处理过程的线程块数目和线程数目即可。
【专利附图】
【附图说明】
[0065]图1是本发明采用的CUDA架构的示意图；
[0066]图2是本发明流程示意图。
【具体实施方式】
[0067]本发明采用自动运行的方式实施。实施例使用CUDA架构提供的扩展C语言作为开发工具。
[0068]CUDA架构示意如图1所示，实施例使用高性能PC和新一代基于Ifepler架构的GPU协处理器实施和运行。CUDA架构中，高性能PC的CPU与GPU协同工作，CPU负责进行逻辑性强的事务处理和串行计算，包括启动CUDA并行处理过程前进行数据准备和设备初始化的工作，以及在各CUDA并行处理过程间进行串行计算，GPU则专注于执行高度线程化的并行处理任务，运行在GPU上的CUDA并行处理过程称为kernel (内核处理过程)，一个kernel存在两个层次的并行，即Grid (线程网格)中的block (线程块)间并行和block (线程块)中的thread (线程)间并行。
[0069]实施例中将信号码速率设为CodeRate、信号采样率设为SampleRate、信号载波频率初始值设为CarrFreq、信号多普勒范围设为Doppler、信号多普勒变化率范围设为DopplerRate、信号的FFT分析点数设为FFTSize、信号的FFT采样频率设为Fs、信号的FFT分析带宽设为RBW、多普勒范围分段的步进值设为PerDoppler、多普勒分段数设为Seg、多普勒变化率范围分段的步进值设为PerDopplerRate、多普勒变化率分段数设为Segl、信号的积分长度设为SumSize、最大谱线对应的多普勒分段号设为MaxDopplerlD，一次FFT分析所有点中最大值对应的FFT点号设为MaxFFTID、最大谱线对应的多普勒变化率分段号设为MaxDoppIerRateIDο
[0070]如图2所示，实施例包括以下步骤:
[0071](I) CPU通过调用CUDA的API函数cudaSetDevice来选择使用O号GPU协处理器并初始化运行环境，同时从配置文件读取信号数据文件路径、信号码速率(CodeRate)、信号采样率(SampleRate)、信号载波频率初始值(CarrFreq)、信号多普勒范围(Doppler)、信号多普勒变化率范围(DopplerRate)、信号的FFT分析点数(FFTSize)、信号体制等运行参数；
[0072](2) CPU根据步骤(I)中设置的信号数据文件路径将信号数据从数据文件读入CPU内存中；
[0073](3) CPU根据步骤(I)中设置的信号码速率(CodeRate)确定信号的FFT采样频率(Fs)和信号的FFT分析带宽(RBW)，然后根据信号的FFT采样频率(Fs)计算对步骤(I)中设置的信号多普勒范围(Doppler)进行分段的步进值(PerDoppler)，进而确定多普勒分段数(Seg)，根据信号的FFT分析带宽(RBW)计算对步骤(I)中设置的信号多普勒变化率范围(DopplerRate)进行分段的步进值(PerDopplerRate),进而确定多普勒变化率分段数(Segl)，最后根据步骤(I)中设置的信号采样频率(SampleRate)和本步骤获得的信号FFT采样频率(Fs)计算信号的积分长度(SumSize)，具体为:
[0074]当信号体制为二相编码体制(BPSK)时:
[0075]Fs = CodeRate ；
[0076]RBW = Fs/FFTSize ；
[0077]PerDoppler = Fs/4 ；
[0078]PerDopplerRate = 3 X RBWX RBW ；
[0079]Seg = 2 X Doppler/PerDoppler+1 ；
[0080]Segl = 2 X DopplerRate/PerDoppleRate+1 ；
[0081]SumSize = SampleRate/Fs+0.5 ；
[0082]当信号体制为统一测控体制(标准TT&C)时，RBW根据信号的情况和经验选取，本实施例中RBW选取为10:
[0083]Fs = RBWXFFTSize ；
[0084]PerDoppler = Fs/4 ；
[0085]PerDopplerRate = 3 X RBWX RBW ；
[0086]Seg = 2 X Doppler/PerDoppler+1 ；
[0087]Segl = 2 X DopplerRate/PerDoppleRate+1 ；
[0088]SumSize = SampleRate/Fs+0.5 ；
[0089](4)CPU根据(3)的计算结果，分配GPU全局存储空间，并初始化已分配的GPU全局存储空间，具体为:
[0090]首先CPU调用CUDA的API函数cudaMalloc在GPU全局存储空间中分配存储空间，用以存储步骤(10)将会用到的从CPU内存中拷贝至GPU全局存储空间的信号数据，该存储空间的大小等于信号的积分长度(SumSize)和信号的FFT分析点数(FFTSize)的乘积，然后调用CUDA的API函数cudaMemset对已分配好的GPU全局存储空间进行初始化；
[0091]其二 CPU调用CUDA的API函数cudaMalloc在GPU全局存储空间中分配存储空间，用以存储步骤(11)的处理结果，该存储空间的大小等于信号的多普勒变化率的分段数(Segl)与信号的FFT分析点数(FFTSize)的乘积，然后调用CUDA的API函数cudaMemset对已分配好的GPU全局存储空间进行初始化；
[0092]其三CPU调用CUDA的API函数cudaMalloc在GPU全局存储空间中分配存储空间，用以存储步骤(13)的处理结果，该存储空间的大小等于信号的多普勒的分段数(Seg)与信号的多普勒变化率的分段数(Segl)的乘积，然后调用CUDA的API函数cudaMemset对已分配好的GPU全局存储空间进行初始化；
[0093]其四CPU调用CUDA的API函数cudaMalloc在GPU全局存储空间中分配存储空间，用以存储步骤(15)的处理结果，该存储空间的大小等于信号的多普勒的分段数(Seg)，然后调用CUDA的API函数cudaMemset对已分配好的GPU全局存储空间进行初始化；
[0094]最后CPU调用CUDA的API函数cudaMalloc在GPU全局存储空间中分配存储空间，用以存储步骤(16)的处理结果，该存储空间为一个包含4个单精度浮点数的结构体，然后调用CUDA的API函数cudaMemset对已分配好的GPU全局存储空间进行初始化；
[0095](5)CPU根据步骤(3)中获得的信号的积分长度设置GPU中在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，具体为:
[0096]每个线程块的线程数目由信号数据的积分长度决定，当积分长度小于等于1024时，线程数目等于积分长度，当积分长度大于1024时，线程数目等于1024 ；
[0097]令线程块数目为(X，y, z)，其中X为多普勒变化率的分段数，y为FFT运算的分析点数，z值为I ；
[0098](6) CPU根据步骤(3)中获得的多普勒变化率分段数和步骤(I)中设置的信号FFT分析点数设置并行一维FFT运算的参数，所述并行一维FFT运算的参数包括并行的一维FFT运算的个数和每个一维FFT运算的FFT分析点数，本实施例中具体为:
[0099]当信号体制为二相编码体制(BPSK)时:
[0100]并行的一维FFT运算的个数等于多普勒变化率分段数；
[0101]每个一维FFT运算的FFT分析点数等于4096 ；
[0102]当信号体制为统一测控体制(标准TT&C)时:
[0103]并行的一维FFT运算的个数等于多普勒变化率分段数；
[0104]每个一维FFT运算的FFT分析点数等于16384 ；
[0105](7)CPU根据步骤(3)获得的多普勒变化率分段数和步骤⑴中设置的信号FFT分析点数设置GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0106]每个线程块中的线程数目由四分之一 FFT分析点数的值决定，当四分之一 FFT分析点数的值小于等于1024时，线程数目等于四分之一 FFT分析点数的值，当四分之一 FFT分析点数的值大于1024时，线程数目等于1024 ；
[0107]线程块数目等于多普勒变化率的分段数；
[0108]本步骤所述二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组，对所有多普勒分段进行本步骤所述二维最大值查找的结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数；
[0109](S)CPU根据步骤(3)获得的多普勒变化率分段数和多普勒分段数设置GPU对步骤(7)中所述的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0110]每个线程块中的线程数目由多普勒变化率分段数决定，当多普勒变化率分段数小于等于1024时，线程数目等于多普勒变化率分段数，当多普勒变化率分段数大于1024时，线程数目等于1024;
[0111]线程块数目等于多普勒的分段数；
[0112]本步骤所述二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组；
[0113](9)CPU根据步骤(3)获得的多普勒分段数设置GPU对步骤⑶所述的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，具体为:
[0114]每个线程块的线程数目由多普勒分段数决定，当多普勒分段数小于等于1024时，线程数目等于多普勒分段数，当多普勒分段数大于1024时，线程数目等于1024 ；
[0115]线程块数目等于I;
[0116](1)CPU将步骤⑵中读入到CPU内存中的信号数据从CPU内存拷贝至步骤(4)确定的GPU全局存储空间中；
[0117]具体实施时，CPU调用CUDA的API函数cudaMemcpy完成；
[0118](Il)GPU利用步骤(5)中确定的在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，以并行的方式完成对步骤(10)读入到GPU全局存储空间中的信号数据的下变频处理，然后以归约算法进行并行求和，并对并行求和后的结果进行鉴相处理，将鉴相处理结果写回步骤(4)确定的GPU全局存储空间中，所述每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数；
[0119]具体实施时当步骤(5)确定的GPU中在每个多普勒分段下所有多普勒变化率分段下对信号进行下变频和积分处理的每个线程块中的线程数目等于步骤(3)计算得到的信号的积分长度时，每个线程块中的每一个线程对应一个信号采样点的下变频过程和用归约算法实现的积分过程，当每个线程块中的线程数目小于步骤(3)计算得到的信号的积分长度时，每个线程块中的每一个线程对应多个采样点的下变频过程和用归约算法实现的积分过程；
[0120]具体实施时，当信号体制为二相编码体制(BPSK)时，下变频和积分内核函数中保留鉴相处理过程，当信号体制为统一测控体制(标准TT&C)时，下变频和积分内核函数中屏蔽鉴相处理过程。
[0121](12)GPU利用步骤(6)中确定的并行一维FFT运算的参数，对步骤(11)中的处理结果进行并行一维FFT运算，并将结果写回写回步骤(4)确定的GPU全局存储空间，所述并行一维FFT运算结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数；
[0122]具体实施时CPU使用GPU厂商NVIDIA提供的基于GPU的通用函数库CUFFT，首先根据步骤(6)中确定的并行一维FFT运算的参数,调用CUFFT函数库中的cufftPlanld函数创建一个一维CUFFT句柄，然后调用CUFFT函数库中的cufTtExecC2C函数对步骤(11)的结果启动GPU进行并行一批多个一维的FFT运算，运算结果写回存放步骤(11)结果的同一片GPU全局存储空间；
[0123](13)GPU根据步骤(7)中确定的对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(12)的并行一维FFT运算结果分别进行取模获得FFT模值，并分别对FFT模值以归约算法并行查找最大值和最大值对应的FFT点号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组；
[0124]具体实施时当步骤(7)确定的GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的每个线程块中的线程数目等于四分之一 FFT分析点数时，每个线程块中的每一个线程对应一个FFT运算结果的取模和归约求最大值过程，当每个线程块中的线程数目小于四分之一 FFT分析点数时，每个线程块中的每一个线程对应多个FFT运算结果的取模和归约求最大值过程；
[0125](14)GPU循环进行步骤(11)至步骤(13)，循环次数为多普勒分段数，获得对所有多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果，所述结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数；
[0126](15)GPU根据步骤(8)中确定的对步骤(14)获得的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(14)获得的二维数组分别以归约算法并行查找FFT模值在每一个多普勒分段下所有多普勒变化率范围内的最大值和对应的多普勒变化率分段号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述对步骤(14)获得的二维数组进行二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组；
[0127]具体实施时当步骤(8)确定的对步骤(14)获得的二维数组进行二维最大值查找的每个线程块中的线程数目等于多普勒变化率的分段数时，每个线程块中的每一个线程对应一个多普勒变化率分段FFT模值的归约求最大值过程，当每个线程块中的线程数目小于多普勒变化率的分段数时，每个线程块中的每一个线程对应多个多普勒变化率分段FFT模值的归约求最大值过程；
[0128](16)GPU根据步骤(9)中确定的对步骤(15)获得的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程并行的方式完成对步骤(15)获得的一维数组以归约算法并行查找FFT模值在所有多普勒分段中的最大值和对应的最大多普勒分段号，并将结果写回GPU全局存储空间；
[0129]具体实施时当步骤(9)确定的对步骤(15)获得的维度与多普勒分段数相同的一维数组进行一维最大值查找的每个线程块中的线程数目等于多普勒分段数时，每个线程块中的每一个线程对应一个多普勒分段最大谱线值的归约求最大值过程，当每个线程块中的线程数目小于多普勒分段数时，每个线程块中的每一个线程对应多个多普勒分段最大谱线值的归约求最大值过程；
[0130](17)CPU将步骤(16)中获得的FFT模值在所有多普勒分段中的最大值及其对应的多普勒分段号、所述FFT模值在所有多普勒分段中的最大值在步骤(14)获得的对应的FFT点号和在步骤(15)获得的对应的多普勒变化率分段号从GPU全局存储空间中拷贝回CPU内存，并计算载波频偏和多普勒变化率，具体为:
[0131]载波频偏由公式:
[0132]EstCarrFreq = (MaxDopplerID-Seg/2)XPerDoppIer+MaxFFTIDX RBff/2 ；
[0133]给出，式中，MaxDopplerID为FFT模值在所有多普勒分段中最大值对应的多普勒分段号，MaxFFTID为FFT模值在所有多普勒分段中最大值对应的FFT点号；
[0134]多普勒变化率由公式:
[0135]EstDopplerRate = (MaxDopplerRateID-Segl/2)XPerDopplerRate
[0136]给出，式中MaxDopplerRateID为FFT模值在所有多普勒分段中最大值对应的多普勒变化率分段号。
[0137]信号的快速捕获是测控通信的关键技术之一，随着深空探测技术的发展，信号捕获技术的难度和复杂度大大增加，当信号极其微弱且多普勒动态范围大、变化率高时，信号捕获一般采用二维FFT捕获方法，即首先对多普勒分段，再在每个多普勒频点对多普勒变化率分段，并在每个分段上进行相关、积分和FFT运算，查找最大谱线，这使得信号捕获过程的运算量呈级数增长。现有的研究中一般采用高性能FPGA来实现这种二维FFT信号频率快速捕获算法，这种方式需要开发专门的硬件板卡，实现复杂，且FPGA的资源有限，一般一次只能实现一种体制信号的捕获，通用性差、可扩展性差。
[0138]本发明充分李永乐CUDA架构和GPU的工作特点，能够充分发挥GPU并行计算能力，同时适合多种体制的信号，有效提高了捕获速度且易于实现，适用于深空探测中测控信号的快速捕获，尤其是信号极其微弱且多普勒动态范围大、变化率高时，可以有效地实现信号频率的快速捕获。
[0139]本发明说明书中未作详细描述的内容属于本领域技术人员的公知技术。
【权利要求】
1.一种基于CUDA的多体制信号频率快速捕获方法，其特征在于步骤如下: (1)CPU通过CUDA架构选择GPU协处理器并初始化CPU和GPU的运行环境，同时对运行参数进行初始化设置，所述运行参数包括:信号数据文件路径、信号码速率、信号采样频率、信号载波频率初始值、信号多普勒范围、信号多普勒变化率范围、信号的FFT分析点数和信号体制； (2)CPU根据步骤(I)中设置的信号数据文件路径将信号数据从数据文件读入CPU内存中； (3)CPU根据步骤(I)中设置的信号码速率确定信号的FFT采样频率和信号的FFT分析带宽，然后根据信号的FFT采样频率计算对步骤(I)中设置的信号多普勒范围进行分段的步进值，进而确定多普勒分段数，根据信号的FFT分析带宽计算对步骤(I)中设置的信号多普勒变化率范围进行分段的步进值，进而确定多普勒变化率分段数，最后根据步骤(I)中设置的信号采样频率和本步骤获得的信号FFT采样频率计算信号的积分长度； (4)CPU分配GPU全局存储空间，并初始化已分配的GPU全局存储空间； (5)CPU根据步骤(3)中获得的信号的积分长度设置GPU中在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目； (6)CPU根据步骤(3)中获得的多普勒变化率分段数和步骤⑴中设置的信号FFT分析点数设置并行一维FFT运算的参数，所述并行一维FFT运算的参数包括并行的一维FFT运算的个数和每个一维FFT运算的FFT分析点数； (7)CPU根据步骤(3)获得的多普勒变化率分段数和步骤(I)中设置的信号FFT分析点数设置GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，所述二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组，对所有多普勒分段进行二维最大值查找的结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数； (8)CPU根据步骤(3)获得的多普勒变化率分段数和多普勒分段数设置GPU对步骤(7)中所述的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，所述二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组； (9)CPU根据步骤(3)获得的多普勒分段数设置GPU对步骤⑶所述的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目； (10)CPU将步骤(2)中读入到CPU内存中的信号数据从CPU内存拷贝至步骤(4)确定的GPU全局存储空间中； (11)GPU利用步骤(5)中确定的在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，以并行的方式完成对步骤(10)读入到GPU全局存储空间中的信号数据的下变频处理，然后以归约算法进行并行求和，并对并行求和后的结果进行鉴相处理，将鉴相处理结果写回步骤(4)确定的GPU全局存储空间中，所述每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数； (12)GPU利用步骤(6)中确定的并行一维FFT运算的参数，对步骤(11)中的处理结果进行并行一维FFT运算，并将结果写回步骤(4)确定的GPU全局存储空间，所述并行一维FFT运算结果为一个二维数组，该数组的一个维度为多普勒变化率分段数，另一个维度为信号的FFT分析点数； (13)GPU根据步骤(7)中确定的对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(12)的并行一维FFT运算结果分别进行取模获得FFT模值，并分别对FFT模值以归约算法并行查找最大值和最大值对应的FFT点号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果为一个维度与多普勒变化率分段数相同的一维数组； (14)GPU循环进行步骤(11)至步骤(13)，循环次数为多普勒分段数，获得对所有多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的结果，所述结果为一个二维数组，该数组的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数； (15)GPU根据步骤(8)中确定的对步骤(14)获得的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程块及多个线程并行的方式完成对步骤(14)获得的二维数组分别以归约算法并行查找FFT模值在每一个多普勒分段下所有多普勒变化率范围内的最大值和对应的多普勒变化率分段号，并将获得的结果写回步骤(4)确定的GPU全局存储空间，所述对步骤(14)获得的二维数组进行二维最大值查找的结果为一个维度与多普勒分段数相同的一维数组； (16)GPU根据步骤(9)中确定的对步骤(15)获得的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，以多个线程并行的方式完成对步骤(15)获得的一维数组以归约算法并行查找FFT模值在所有多普勒分段中的最大值和对应的最大多普勒分段号，并将结果写回步骤(4)确定的GPU全局存储空间； (17)CPU将步骤(16)中获得的FFT模值在所有多普勒分段中的最大值及其对应的多普勒分段号、所述FFT模值在所有多普勒分段中的最大值在步骤(14)获得的对应的FFT点号和在步骤(15)获得的对应的多普勒变化率分段号从GPU全局存储空间中拷贝回CPU内存，并计算载波频偏和多普勒变化率。
2.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述信号包括BPSK信号和标准TT&C信号。
3.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(3)中，CPU根据步骤(I)中设置的信号码速率确定信号的FFT采样频率和信号的FFT分析带宽，然后根据信号的FFT采样频率计算对步骤(I)中设置的信号多普勒范围进行分段的步进值，进而确定多普勒分段数，根据信号的FFT分析带宽计算对步骤(I)中设置的信号多普勒变化率范围进行分段的步进值，进而确定多普勒变化率分段数，最后根据步骤(I)中设置的信号采样频率和本步骤获得的信号FFT采样频率计算信号的积分长度，具体为: 当信号体制为二相编码体制，即BPSK时:
Ps — CodeR&te ；
RBW = Fs/FFTSize ；
PerDoppler = Fs/4 ；
PerDopplerRate = 3 X RBff X RBff ；
Seg = 2 X Doppler/PerDoppler+1 ；
Segl = 2 X DopplerRate/PerDoppleRate+1 ；
SumSize = SampleRate/Fs+0.5 ； 当信号体制为统一测控体制，即标准TT&C时，RBW根据信号的情况和经验选取:
Fs = RBffXFFTSize ；
PerDoppler = Fs/4 ；
PerDopplerRate = 3 X RBff X RBff ；
Seg = 2 X Doppler/PerDoppler+1 ；
Segl = 2 X DopplerRate/PerDoppleRate+1 ；
SumSize = SampleRate/Fs+0.5 ； 式中，CodeRate为信号码速率，Fs为FFT采样频率，RBW为信号的FFT分析带宽，Doppler为信号多普勒范围，DopplerRate为多普勒变化率范围，PerDoppler为多普勒范围分段的步进值，Seg为多普勒分段数，PerDopplerRate为多普勒变化率范围分段的步进值，Segl为多普勒变化率分段数，FFTSize为信号的FFT分析点数，SumSize为积分长度，SampleRate为信号采样频率。
4.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(5)中，CPU根据步骤(3)中获得的信号的积分长度设置GPU中在每个多普勒分段下所有多普勒变化率分段对信号进行下变频和积分处理的线程块数目和每个线程块中的线程数目，具体为: (1)每个线程块的线程数目由信号数据的积分长度决定，当积分长度小于等于1024时，线程数目等于积分长度，当积分长度大于1024时，线程数目等于1024 ； (2)令线程块数目为(x，y，z)，其中X为多普勒变化率的分段数，y为FFT运算的分析点数，z值为I。
5.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(7)中，CPU根据步骤(3)获得的多普勒变化率分段数和步骤(I)中设置的信号FFT分析点数设置GPU对每个多普勒分段下所有多普勒变化率分段的FFT模值进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为: (1)每个线程块中的线程数目由四分之一FFT分析点数的值决定，当四分之一FFT分析点数的值小于等于1024时，线程数目等于四分之一 FFT分析点数的值，当四分之一 FFT分析点数的值大于1024时，线程数目等于1024 ； (2)线程块数目等于多普勒变化率的分段数。
6.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(8)中，CPU根据步骤(3)获得的多普勒变化率分段数和多普勒分段数设置GPU对步骤(7)中所述的一个维度为多普勒分段数，另一个维度为多普勒变化率分段数的二维数组进行二维最大值查找的线程块数目和每个线程块中的线程数目，具体为: (I)每个线程块中的线程数目由多普勒变化率分段数决定，当多普勒变化率分段数小于等于1024时，线程数目等于多普勒变化率分段数，当多普勒变化率分段数大于1024时，线程数目等于1024; (2)线程块数目等于多普勒的分段数。
7.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(9)中，CPU根据步骤(3)获得的多普勒分段数设置GPU对步骤(8)所述的维度与多普勒分段数相同的一维数组进行一维最大值查找的线程块数目和每个线程块中的线程数目，具体为: (1)每个线程块的线程数目由多普勒分段数决定，当多普勒分段数小于等于1024时，线程数目等于多普勒分段数，当多普勒分段数大于1024时，线程数目等于1024 ； (2)线程块数目等于I。
8.根据权利要求1所述的一种基于CUDA的多体制信号频率快速捕获实现方法，其特征在于:所述步骤(17)中计算载波频偏和多普勒变化率，具体为: 载波频偏由公式:
EstCarrFreq = (MaxDopplerID-Seg/2)XPerDoppIer+MaxFFTIDX RBff/2 ； 给出，式中，MaxDopplerID为FFT模值在所有多普勒分段中最大值对应的多普勒分段号，MaxFFTID为FFT模值在所有多普勒分段中最大值对应的FFT点号； 多普勒变化率由公式:
EstDopplerRate = (MaxDopplerRateID-Segl/2)XPerDopplerRate 给出，式中MaxDoppIerRateID为FFT模值在所有多普勒分段中最大值对应的多普勒变化率分段号。
【文档编号】H04L12/26GK104202209SQ201410381743
【公开日】2014年12月10日 申请日期:2014年8月5日 优先权日:2014年8月5日
【发明者】陈蓉, 王静温 申请人:北京遥测技术研究所, 航天长征火箭技术有限公司