基于CUDA架构的AMTI算法实现

金 莉，王仁志，宋万杰

(西安电子科技大学 雷达信号处理国家重点实验室，陕西 西安710071)

雷达所处的大环境中存在各种杂波，与目标信号相比，杂波往往具有更大的功率，所以雷达信号处理的关键在于从杂波背景中检测出目标。雷达信号处理系统中，自适应动目标显式(AMTI)因为对静杂波和动杂波都有抑制能力，成为常用的杂波抑制方式。实现AMTI技术的方法主要有两种：一种是两级动目标显示级联形式实现，让第一级凹口位于多普勒零频处，抑制静止杂波；第二级设计凹口抑制动杂波；另一种是设计双凹口杂波抑制滤波器[1]。

本文采用NVIDIA的图形处理单元(Graphics Processing Unit ,GPU)来实现AMTI算法。自GPU概念提出以来，GPU峰值性能一直以超过摩尔定律的速度增加，平均每6个月翻一番[2-3]。相对CPU而言，GPU具有浮点计算能力强、性价比高、能耗低等优点，被广泛用于通用计算中。特别是统一设备架构(Compute Unified Device Architecture ,CUDA)的发布，将图形硬件和应用接口程序进行封装，使开发人员可以直接采用C语言进行编程，方便易学[4]。以GPU为平台，采用CUDA技术，在单指令多数据(Single Instruction Multiple Data，SIMD)编程模型下实现AMTI算法，发挥GPU强大的计算能力，提高运算效率。

1 AMTI算法

为了滤除杂波，需要在杂波信号位置形成凹口。因为动杂波不位于零频附近，所以可以在动杂波谱位置形成凹口，也可以将动杂波移动至凹口附近来滤除动杂波。即有两种方法实现AMTI技术：一种是将杂波信号谱中心移到零频凹口处，然后采用MTI滤波器进行滤波处理；另一种是将抑制杂波的凹口移动至杂波谱中心直接实现滤波，这种方法需要事先设计凹口位于不同频率点的滤波器，并且将滤波器权系数存储在一个库里，称为权系数库法。在工程中常使用这种方法，首先设计出凹口位于零频处的MTI滤波器，然后在频率轴平移得到滤波器权系数库，当估计到杂波频谱中心后，直接调用对应的滤波器进行对消处理[5-6]。本文采用基于权系数库实现的AMTI算法，图1为该算法的原理框图。

图1 基于权库法的AMTI原理框图

基于权系数库的AMTI算法的思想是：先估计动杂波的多普勒范围f∈(fdmin,fdmax)。然后，在这个估计的多普勒范围内，等间隔的设计一组MTI杂波抑制滤波器，最后将该杂波抑制滤波器的滤波器系数直接放入存储设备中备用[7-8]。具体实现过程中，首先需要确定调用的MTI滤波器，具体方法为：先估计杂波的多普勒中心频率f，然后依据f调用对应的MTI杂波抑制滤波器滤波。因此，基于权系数库的AMTI滤波器的设计重点在于准确估计动杂波的多普勒中心频率f。在实际的工程实现中，动杂波多普勒中心频率f通过计算相邻脉冲间相位差得到。

考虑到雷达回波信号x(t)的表达式为

x(t)=A(t)ej(wdt+φ0)+n(t)

(1)

式中，A(t)为回波信号的幅度；wd为杂波的多普勒频率，并有wd=2πfd；w0表示初回波信号的初始相位；加性噪声n(t)在不同脉冲间互不相关。假设Tr为脉冲的重复周期，那么雷达接收到的下一个Tr的发射脉冲回波为

x(t-Tr)=A(t-Tr)ej(wd(t-Tr)+φ0)+n(t-Tr)

(2)

那么，两次回波信号的相关函数可表示为

Rxx=E[x(t)x*(t-Tr)]=E[A(t)A(t-Tr)]ejwdTr

(3)

式(3)中A(t)为窄带信号，因此，可以得出式E[A(t)A(t-Tr)]=E[|A(t)|2]的结果为实数。因此，杂波的多普勒频率wd为

(4)

进一步

(5)

式(5)中，i表示不同脉冲回波信号序列。采用相邻两个脉冲间相位差的方法，得到f表达式为

(6)

为了准确估计杂波多普勒中心频率，实际工程中常采用质心法，质心法就是将回波信息中的幅度值也考虑在内，多次估计的杂波多普勒中心频率，最后进行平滑处理[9]。进一步提高f估计的准确性。

2 CUDA编程模型

CUDA编程模型将CPU作为主机(Host),GPU作为协处理器(Co-processor)或者设备(Device)。在一个系统中可以存在若干个主机和设备[10]。在CUDA编程模中，CPU和GPU协同工作，CPU主要负责整个串行运算的执行，控制算法执行流程，初始化数据将其传递到GPU，GPU接收CPU传递过来的数据，进行并行计算，计算完成后，将计算结果传递回主机[11-12]。另外，GPU和CPU拥有各自独立的存储空间，CPU端为主机内存，GPU端为显存，GPU对显存的操作一般与CPU对内存操作一样，需要调用CUDA API中的内存管理函数。

运行在GPU上的函数称为核函数(Kernel),核函数并不是一个完整的程序，而是整个CUDA程序中一个可以被并行执行的步骤[13-15]。如图2所示，一个完整的CUDA程序是由主机端的串行代码和设备端的并行函数Kernel(·)共同组成，各个并行函数之间为串行关系。

图2 CUDA并行编程模型

3 AMTI算法实现流程

为满足雷达系统要求，本算法采用5参差AMTI，滤波器间隔是5 Hz，从0～90 Hz共19组AMTI双凹口滤波器，图3为AMTI算法实现流程图。

图3 基于权系数库的AMTI算法实现流程图

基于CUDA架构的AMTI算法实现，主要流程如下：

步骤1在CPU上进行数据初始化，内存以及显存分配，将准备就绪的数据通过CUDA API中的拷贝函数拷贝至GPU；

步骤2首先在GPU上完成杂波谱中心f的估计，并判断频率值的正负值，将该值存放在单独数组中；

步骤3接下来在GPU上进行AMTI处理，具体处理流程如图3所示，最终实现抑制杂波的目的；

步骤4再次调用CUDA API中的拷贝函数就计算结果从设备拷贝回主机；

步骤5在CPU上释放分配的内存以及显存。

图4为采用CUDA编程架构的AMTI算法实现流程图。

GPU上实现AMT处理，主要是抽取AMTI中的并行部分，分配到各线程块中的线程进行执行。就杂波谱中心估计来说，包含脉冲数和数据量两层循环，但GPU中只要合理分配线程及线程块，让每个线程执行一个脉冲一组数据的运算，就可以大幅提高运算效率。本实现中，分别用x维线程控制数据量，y维线程控制脉冲数，运算完成后批量增加循环过程中的累加值，代码如下：

int tidx=threadIdx.x+blockIdx.x*blockDim.x;∥数据量

int tidy=threadIdx.y+blockIdx.y*blockDim.y;∥脉冲数

while(tidx

while(tidy

{

∥计算得到滤波器系数相关地址及符号

tidy+=blockDim.y*gridDim.y; 距离单元/个

}

∥计算系数偏移地址去除系数完成AMTI处理

tidx+=blockDim.x*gridDim.x; ∥批量并行运算

}

需要说明的是，在AMTI实现过程中，需要计算反正切值。反正切值是通过查表法得到。事先计算好得到每个角度对应的正切值，存放在内存中以备使用。为了提高正切值的精度，可以采用先乘数再取整的方法。具体实现方法为，将求得的值首先乘以一个倍数，然后进行取整，放入表格中，根据具体雷达系统对数据精度的要求，确定该倍数大小。

4 仿真结果及分析

仿真采用的操作系统为Windows 7 SP1，显卡是计算能力为2.0的NVIDIA Tesla C2050，该显卡支持双精度浮点运算，满足雷达系统对精度的要求。为了提高运算精度，数据处理均采用双精度浮点数，图5为AMTI前脉压数据。

图5 脉压原始数据

实验分别在CPU和GPU上实现AMTI算法，其结果分别如图6和图7所示。

图6 CPU平台实现AMTI

从图6和图7可以看到滤除了地杂波，目标幅度达到103数量级，可以认为杂波抑制效果良好且CPU和GPU处理结果一致，但两者的处理速度有很大差别。

图7 GPU平台实现AMT

表1 CPU与GPU所运行时间对比/ms

从表1可以看出， GPU比CPU的数据运算效率更高。需要说明的是，表格中括号内数字表示用于主机和设备之间数据传输的时间，例如AMTI算法中，总时间0.645 ms中仅0.02 ms用于数据处理，括号中为用于数据传输的时间，即GPU有很大一部分时间都是用在数据传输上，所以在GPU运算过程中，应当尽量减数据传输时间。

5 结束语

实验表明，通过在GPU上实现AMTI算法，不仅达到了和CPU上一样的效果，而且加快了运算速度，提高了运算效率。但是在CUDA程序中，主机与设备之间的数据传输会耗费很多时间，所以可以采用以下3种方法减少传输时间：(1)由于在设备和主机之间进行数据传输会耗费相当多的时间，所以应该在两次数据传输之间在GPU上进行尽可能多的运算，或者采用流的处理方式，将数据进行分块处理，传输的同时进行运算；(2)因为纹理存储器的存储方式为只读，处理速度相对其他存储器速度较快，可以采用纹理存储器来保存一些反三角函数表，滤波器系数等，以优化数据处理效率；(3)合理安排线程块以及线程数量，例如当使用二维线程块以及线程索引时，可以利用dim3这个变量进行定义，最好以一个线程束的倍数定义线程和线程块数量来最大限度提高运算速率。总的来说，利用CUDA处理数据相比CPU而言，较好地满足了目前信息处理对高精度、实时性要求。