OpenMP并行计算在侧向散射激光雷达数据处理中的应用∗

刘俊见 陶宗明

（1.陆军军官学院研究生管理大队 合肥 230031）（2.陆军军官学院基础部 合肥 230031）

1 引言

大气气溶胶是指液体或固体微粒散布在大气中形成的相对稳定的悬浮体系，对人体健康环境和全球气候都有一定的影响［1］，是当前大气科学领域中的热点问题之一。在大气环境监测领域，准确了解气溶胶随时空的分布规律具有十分重要的意义。

激光雷达以它精细的时间、空间和光谱分辨率、大的垂直跨度和高的探测精度很快便成为探测大气气溶胶空间分布的有力工具，被广泛应用于环境与大气监测等领域［2～6］。基于CCD的侧向散射激光雷达是正在研究中的新技术，它将发射装置与接收装置分两处放置，解决了传统后向散射激光雷达中的过渡区和盲区问题，而且由于使用了侧向技术，其在近距离段具有极高的空间分辨率［7～10］。

为了让侧向散射激光雷达在白天也能工作，克服白天背景光过强的有效方法是减小曝光时间和多次曝光求信号的平均值。我们所用CCD的图片的像素是3352*2532，大小约为16M，多次曝光（如300次）提取信号的平均值需要占用较多的计算时间，这对数据的实时处理造成一定的困难。

随着多核处理器架构被业界广泛采用，通过并行计算缩短执行时间已成为高性能计算领域使用最广的方法［11］。OpenMP（Open Multi-Processing）是一种典型的共享内存平台下的编程API，支持跨平台共享内存方式的多线程并发，其通过提供对并行算法的高层抽象描述，降低了并行编程的复杂性和难度［12］。

本文通过分析侧向散射激光雷达数据处理流程的特点，找出迫切需要采用并行计算提高效率的关键任务。基于OpenMP标准，设计并行优化方案，实现侧向散射激光雷达数据处理计算并行化，通过实验分析并检验其并行效果，以期有效挖掘不同计算设备的潜能，为采用并行计算技术处理侧向散射激光雷达数据提供一种高效、经济的实施方案。

2 侧向散射激光雷达数据处理中的关键计算任务分析

激光雷达大气探测的监控周期一般为10分钟，即每间隔10分钟采集一组数据。对于基于CCD的侧向散射激光雷达，由于采集的图像数据占用存储空间庞大，以一张CCD图像16M、每组300对图像（每对为两张图片，用于减少背景光、干扰光）为例，白天一小时就要采集3600张（6*300*2），需要占用接近60G的硬盘空间，因此需要在一个周期内完成该组数据的采集、处理、保存并删除初始CCD图像数据，这对数据处理的效率提出了很高的要求。从侧向散射激光雷达数据处理的需求出发，计算主要集中在以下几个方面，包括信号提取、配准、反演、可视化显示等运算。

基于上述的侧向散射激光雷达数据处理中的几个关键任务，下面对它们的并行需求作一个简要分析。在进行反演计算时，一般采取分像元处理的策略，由于各位置（或光束方向上各像元）之间的计算是没有依赖关系的，因此，在多核环境下，可以同时启动多个进程，每个进程分别完成对应某部分像元的反演计算，既能充分利用设备的系统资源，又避免了设计对应的并行算法；可视化的主要目的是为了方便研究人员对原始数据做出正确解释，相对于其他的处理任务，其处理时间灵活，甚至可以在一天的观测任务结束后再绘制当天的气溶胶后向散射系数时空分布图；配准的计算量相对较少，主要是循环及方差计算，直接在已有的串行代码基础上嵌入指令语句将其并行化就可以取得较好的加速效果。因而信号提取是侧向散射激光雷达数据处理中迫切需要设计对应并行算法的关键。

针对夜间图像，信号提取主要指在减少背景光、干扰光的基础上，通过高斯拟合去噪并提取对应角宽度内的散射光子数。减少背景光、干扰光主要通过连续拍摄曝光时间相同的两张CCD相片［13］，如图1所示，（a）是包含所有光的CCD相片，（b）是缺少激光散射光束的CCD相片，将两次拍摄的CCD对应像元信号相减，获得的结果如（c）所示，这就减少了背景光和干扰光的影响；提取散射光子数主要通过高斯拟合逐一计算沿光束方向的每个像元对应的信号光子数与噪声光子数，从而得到光子数与位置（或光束方向上像元）的一一对应关系。

图1 减少背景光、干扰光

信号提取的难度主要在白天图像的处理上，相对于夜间，白天作为主要背景光的日光太强，曝光时间过长会导致图像过曝，因此需采取减少单张图像曝光时间、多张图像叠加的方式来获取CCD图像数据。为了减少局部误差、获得较高的信噪比，白天每次最少需要采集300对～500对的CCD图像数据，每对图像都需要减少背景光、干扰光，去噪并提取散射光子数，随后将提取的散射光子数逐对逐像元进行叠加求和，得到最终的结果。同时为了应对系统崩溃、采集与处理不同步等意外情况，数据处理算法必须要有快速处理之前积压的多组数据的能力，要能高效率地从上千张CCD图像中提取出每个时段的信号，这对数据处理的效率提出了更高的要求。

3 动态设置并行循环的线程数量

在侧向散射激光雷达数据处理的实际情况中，由于不同情况下对数据精度及成本的双重考虑，数据处理程序将运行在不同的设备环境里，这些设备可能是双核、4核甚至更多，同时由于超线程技术的应用，某些型号的CPU可以同时执行多重线程。同时对于不同的数据处理任务，由于需要处理的数据量及循环次数不同，执行该任务时需要的线程数也不同。因此需要根据设备性能的不同来自动设置合适的线程数量，以满足不同情况下的需求，使数据处理程序具有较强的可扩展性，当运算量发生变化时，仍能通过自动设置使线程数量满足要求。

通常认为，在具体计算需要使用多少线程时，主要需要考虑以下两点［14］：

1）当循环次数比较少时，如果分成过多数量的线程来执行，可能会使得总运行时间高于较少线程或一个线程执行的情况，并且会增加能耗。

2）如果设置的线程数量远大于CPU核数的话，那么存在着大量的任务切换和调度等开销，也会降低整体效率。

由于不同数据处理任务的循环次数差距很大，夜间进行信号提取时，每组数据为300～500对CCD图片，那么对每组图片进行信号提取需要进行300～500次循环；而进行配准计算时的循环次数只有10～30次，那么如何根据循环的次数和CPU线程数来动态地设置线程的数量呢？通过分析侧向散射激光雷达数据处理任务的特点，确定了动态设置线程数的需求如下：

1）多个线程运行时，由于现有CPU线程数一般不超过16，因此每个线程运行的循环次数应不低于总循环次数的十六分之一，同时不少于4次。

2）总的运行线程数不超过CPU总线程数，以避免任务切换开销。

基于上文确定的需求，动态设置线程数的程序流程图如图3所示。其中Num为CPU总线程数，Run_num为设置的执行线程数。

图2 动态设置线程数程序流程图

实际编程中，将其编为一个独立的功能函数，其返回值即是针对此次循环应设置的线程数，每次并行化循环开始前调用该函数来获取合适的线程数量，并使用num_threads子句设置线程数。

4 并行方案及并行方法

并行算法设计常见的模式有数据分解模式、分治模式、流水线模式、任务并行模式、任务图调度模式、动态任务调度模式等，针对夜间信号提取这一关键任务，从其数据处理流程特点与需求出发，选取了两种最适合的模式分别进行了并行算法的设计。

该任务的串行算法A可以简化为三层循环，如下所示，内层循环完成一对CCD图像沿光束方向上各像元信号的提取，output即为提取出的序号为j的像元的散射光子数；第二层循环完成数量为amount的一组CCD图像的信号提取，每对图像提取出的散射光子数进行累加求和并存储在数组a（j）中；外层循环是针对因为意外堆积的多组数据，find函数采用递归策略扫描文件夹内所有的CCD图像路径，返回值find（）为图像对数。

4.1 方法一

数据分解模式即是将数据分解为N个独立的数据子块，每个线程处理其中的一个或多个数据块。由于每对CCD图像之间都是相互独立的，符合数据分解模式的要求，因此选择其作为第一种方法完成并行算法的设计，此方法记为算法B。

由于算法A为三层循环结构，且最外层循环次数一般较少，若对外层循环并行化，实际创建的线程数很可能小于CPU的总线程数，不能很好地发挥CPU性能，会大大降低并行算法的效率；同时很难通过对外层循环进行调度使内层循环达到负载均衡。针对上述问题，本文在进行并行算法设计时，首先将外层循环和第二层循环进行合并，合并后的两层循环结构如下所示。

其中number为采集完成的图像组数，计算流程如图3所示，随后在内层循环通过判断图像名的标志位，将提取出的散射光子数按组别分别累加。

图3 计算图像组数流程图

针对合并后的两层循环，若对内层循环进行并行化，就必须要解决a（j）等数组的同时操作问题。内层循环中，每个线程完成信号提取后，都需要和a（j）累加，因此每个线程都要对a（j）进行读写操作。对这类问题，OpenMP提供了临界区、数据规约等方式来进行线程之间的同步，但是不论哪种方法都会在一定程度上增加系统开销。本文选择只对外层循环进行并行化，避开内层循环的锁竞争，经过比较，该方法的加速性能与采用数据规约基本一致，并优于采用临界区的方法。

虽然需要处理的图像对数M难以被线程数N整除，CCD图像数据无法被划分为计算量相同的N个数据子块，即不使用size参数时不同线程实际分配的迭代次数会相差1，但相对于每个线程分配的任务总量而言，某一个或某几个线程多处理一对图像几乎不会对总体加速性能产生影响，因此本文采用静态调度方法实现并行算法的负载均衡，关键代码如下。

4.2 方法二

由于信号提取首先需要从硬盘中读取大量的数据，很大程度上存在着IO瓶颈，本文选择流水线模式作为第二种方法对并行算法进行优化，此算法记为算法C。流水线模式借鉴工业生产中的流水线思想，将一个数据处理流程分解为几个步骤，每个步骤采用一个或多个线程来进行处理，在处理具有IO瓶颈的数据处理流程有较大优势。

针对信号提取这一流程，我们可以将其分解为如图4所示的三个步骤。

图4 对信号提取过程进行流水线分解

串行算法采用一个线程依次执行上述的三个步骤；而采用流水线模式进行并行实现时，分别使用一个或多个线程来读取CCD图像数据，减少背景光、干扰光，以及提取散射光子数。这样就使得三个步骤的处理并行进行，提高了处理效率，如图5所示。在并行实现的过程中，当读取数据线程获取了一部分CCD图像数据后，后续的线程就可以开始进行数据处理了。

图5 对信号提取过程进行并行处理

但是从上述的分析中也可以看出，运用流水线模式对信号提取过程进行并行实现时，数据从上一个步骤传递到下一个步骤时，这两个步骤是并行执行的，因而存在数据竞争，需要使用同步操作对数据进行保护，针对此问题，本文主要采用临界区的方法来实现同步操作。

并行实现时，为避免线程资源的浪费以及非线性读取对读入速度的影响，笔者首先指定了一个线程来读取数据。在实现后续两个步骤的并行化时，笔者发现很难将这两个步骤分开实现并行化，主要问题在于以下两点：

1）程序设计时被要求能在不同性能的设备上高效运行，不同设备可供使用的线程资源差距明显，很难对其进行线程分配。

2）上文中提到了数据保护的问题，并行时每增加一个步骤，实际上也增加了用于数据同步的系统开销。

在进行了多次试验后，笔者决定将后续两个步骤合并，共用数据处理线程来进行并行化，结果显示其加速性能要优于分开并行化的方案。实际上，方法二主要是在方法一的基础上对数据读取进行了一定的优化，数据处理部分的并行化与方法一相同。

5 并行测试及结果分析

按照第4节中的两种算法，对针对夜间信号提取这一关键与难点进行并行化，并运用型号为AlienwareM14xR1的笔记本电脑进行测试，其CPU为主频2.20G的i7-2670QM，内存16G，操作系统为Windows7，编程环境为VS2015。从300对CCD图像中提取散射光子数，算法B和算法C的总运行时间见如图6所示，其中由于算法C无法单线程执行，其线程数为1时直接执行算法A（即串行算法）。

图6 算法B和算法C运行时间比较

图6的横轴为采用的线程数，分别从1到8；纵轴为多线程模式下的该流程运行时间。可以看出，采用线程数较少时，算法B运行时间要比算法C少，而随着线程数增加算法C与算法B运行时间差距逐渐缩小，线程数为4时，两种算法的运行时间基本相同，当线程数大于4后，算法C的运行时间开始少于算法B。

为了更好地对两种并行算法结果进行分析，我们引入了并行加速比和并行效率对并行算法的性能进行评价，其定义如下［15］：

其中，p指使用的线程数；T1为串行算法运行时间，Tp为p个线程并行执行的时间，Sp即为并行加速比，Ep为并行效率。因此，Sp的理论值为p，而Ep的理论值为1，这里给出两种算法的并行加速比和并行效率，如图7。图7（a）是加速比，虚线是其理论值p。随着采用的线程数增多，导致系统开销增大，两种算法的加速比逐渐都远离理论值，但很明显随着线程数增加，算法C开始逐渐优于算法B。如图7（b）所示，并行效率的变化趋势也是如此。当线程数为8时，算法B的并行效率约为0.49，相当于计算速度提高到串行算法的4倍；算法C的并行效率约为0.6，相当于计算速度提高到串行算法的5倍。

图7 两种算法的并行加速比和并行效率比较

算法C由于数据读取单独占用一个线程，以及临界区导致的系统开销，使其在线程数较少时的加速性能并不理想，但是随着线程数增加，算法B的加速性能受限于IO瓶颈，算法C的优势也随之体现出来，当线程数大于6以后，再增加线程数，算法B加速性能的提升微乎其微，而算法C的加速性能仍能获得一个相当可观的提升。因而在实际应用中，可以将两种算法结合起来，无论采用多少线程，都可以获得一个相对可观的并行效率。

6 结语

本文对基于CCD的侧向散射激光雷达数据处理中的并行计算问题进行了研究。针对并行计算的重点与难点，设计了两种并行算法，并采用OpenMP对其进行并行化，阐述了并行化过程的技术难点，并对并行结果进行数值测试与分析。分析结果表明，文中提出的两种并行算法合理可行，结合起来后能充分利用不同性能的设备，提高侧向激光雷达数据处理的运算速度，为激光雷达大气探测提供了计算性能的保障。同时，侧向激光雷达数据处理速度某种程度上仍受限于IO瓶颈，在后续的研究中，笔者将利用聚合IO技术或其他并行IO技术对数据处理性能进行进一步的提升，为进一步利用基于CCD的侧向散射激光雷达探测大气气溶胶奠定基础。