CN108596885A - 基于cpu+fpga的快速sar图像变化检测方法 - Google Patents

Abstract

本发明公开了一种基于CPU+FPGA的快速SAR图像变化检测方法，本发明采用CPU+FPGA异构并行加速实现SAR图像的变化检测，与现有技术中SAR图像变化检测的方法在CPU下计算相比，FPGA的并行计算能力要远远高于CPU数据串行计算能力。并且采用FPGA与CPU结合工作的方式，FPGA负责将整个算法中数据的处理并行化，CPU负责在主机端控制任务的调度以及参数的大小。本发明采用OpenCL的并行编程标准对SAR图像变化检测的算法进行并行化编写，代码的可移植性比较高，在不同的异构平台都可以运行。

Description

基于CPU+FPGA的快速SAR图像变化检测方法

技术领域

本发明属于图像处理领域，具体涉及一种基于CPU+FPGA的快速SAR图像变化检测方法。

背景技术

变化检测问题属于图像处理领域，通常是指依据不同时间的多次观测来对一个物体的状态变化进行确定或者对某种现象的变化进行确定的过程。目前随着自然环境变化的加剧，城市也在迅速发展，各种自然灾害频频发生，图像变化检测技术也随之快速崛起。合成孔径雷达图像(SAR)是变换检测问题的主要数据来源。SAR图像具有全天候、全天时的工作能力，对地物有一定的穿透能力，且不受大气、气候等随机因素的影响，具有不可比拟的优点，所以对SAR图像变换检测的研究是很重要的。

目前，SAR图像变化检测的普遍的思路就是先对两幅同地区不同时间的图像做差异图，再对差异图进行分类比较。较为经典的SAR图像变化检测方法有图像比值法，图像差值法等基于简单代数运算的变化检测以及基于图像分类的变化检测。其中图像差值法算法简单，较容易理解，但是该方法主要通过分析光谱的变化产生差值图像，通常只是利用灰度的差异，容易产生不正确的信息变化。图像比值法在图像呈正态分布的假设下可以克服乘性噪声，但是实际上这个假设是不成立的，缺少适用性。

Yaoguo Zheng等在《Using Combined Difference Image and k-MeansClustering for SAR Image Change Detection》论文中提出一种基于图像融合和聚类的SAR图像变化检测的方法，用预处理后的图像生成差值差异图以及比值差异图，将差值差异图和比值差异图通过融合系数进行融合，对融合后的差异图进行聚类。此方法通过不同差异图的融合可以较好的保持图像区域连续性以及边缘信息。但是这种方法计算复杂度高，无法应用在SAR图像的快速变化检测的处理任务上。因此将算法进行并行实现来有效的降低处理时间是非常有必要的。

采用CPU+FPGA对算法进行加速，采用并行编程标准的开放运算语言OpenCL进行内核编程，OpenCL是第一个面向异构系统通用的并行编程的标准，也是一个统一的编程环境。用OpenCL编写的内核可以在FPGA、GPU和CPU等多种设备上执行，跨平台性良好。FPGA具有非常强大的数据并行计算能力，具有可重复编程的能力，设计灵活性很大，特别适合并行计算，可以处理密集的数据。此外，FPGA功耗较低，在诸如自动驾驶、智能导航和跟踪监控等对功耗有要求的领域，FPGA的低功耗特性就显得尤为重要。CPU可以完成参数的设定以及任务的调度等工作。

发明内容

本发明的目的在于克服上述不足，提供一种基于CPU+FPGA的快速SAR图像变化检测方法，通过对算法计算并行，优化访存方式等，减少耗时，满足实时性的要求。

为了达到上述目的，本发明包括以下步骤：

步骤一，初始化开放运算语言环境，选择异构计算平台以及计算设备，在计算平台上创建上下文，在上下文中创建命令队列，得到异构框架，创建并编译内核程序；

步骤二，将初始图像读入到主机内存中，再将图像数据传递到FPGA的全局内存中；

步骤三，初始化lee滤波，均值滤波；

步骤四，调用内核程序中的LeeFilter_3函数对两幅图像分别进行lee滤波，滤波后的图像数据分别保存在缓存对象lee1和缓存对象lee2中；

步骤五，根据缓存对象lee1和缓存对象lee2调用内核程序中的diff函数生成差值差异图和比值差异图，分别保存在缓存对象ds和缓存对象dl中；

步骤六，调用内核程序中的average函数对缓存对象ds进行均值滤波，以及median函数对缓存对象dl进行中值滤波；

步骤七，根据滤波后的缓存对象ds和缓存对象dl，以及图像融合系数a调用内核程序中的add函数进行差异图的融合，结果保存在缓存对象diff中；

步骤八，根据缓存对象diff以及循环迭代次数loop，在主机端循环调用内核程序kmeans中的函数，对diff进行聚类，将聚类后的结果保存在缓存对象cluster中；

步骤九，将缓存对象cluster通过PCI-E总线传输到CPU主机内存中，在主机端对变化和非变化区域分别标记为白色和黑色，将变化检测后的图像保存，完成图像变化的检测。

步骤一中，开放运算语言环境为OpenCL环境，异构计算平台采用Altera，用OpenCL提供的应用程序编程接口API创建上下文以及命令队列，得到异构框架，使用SDK创建并编译内核程序。

步骤二的具体方法如下，选取两张Ottawa地区不同时刻的灰度图作为原始图像，将原始图像读入到主机内存中，将其变换成行向量后，分别写到缓存对象I₁和I₂中再将数据传递到FPGA的全局内存中。

步骤四的具体方法如下：

第一步，根据步骤三中设置的lee滤波的窗口尺寸d1，对窗口区域内的像素求均值，每个工作项负责一个窗口内所有的数据计算；均值μ_Z：

第二步，根据第一步中的均值计算窗口区域内像素方差：

第三步，根据μ_Z和计算ENL以及SX2：

第四步，计算lee滤波后的像素值xcap：

其中ω为滤波窗口宽度，Z_ij为原始图像(i,j)处的像素值，Z为3×3窗口的中心点像素值。

步骤五的具体方法如下：

第一步，调用缓存对象lee1和lee2生成差值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象ds中，计算公式：

ds(x,y)＝|lee1(x,y)-lee2(x,y)|

第二步，调用缓存对象lee1和lee2生成比值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象dl中，计算公式：

其中lee1(x,y)为第一幅图lee滤波后在坐标(x,y)处的值，lee2(x,y)为第二幅图lee滤波后在坐标(x,y)处的值，ds(x,y)为差值差异图在(x,y)的值，dl(x,y)为比值差异图在(x,y)的值。

步骤六的具体方法如下：

第一步，对缓存对象ds，调用内核程序中的average函数进行均值滤波，滤波后保存在ds中；

第二步，对缓存对象dl，调用内核程序中的median函数进行中值滤波，滤波后保存在dl中。

步骤七中，调用内核程序中的add函数进行差异图融合，此操作分布在若干个工作组中，每个工作组中工作项的大小为64×64，所有工作项并行工作，每个工作项处理一个数据，差异图融合公式为：

diff＝α*ds+(1-α)*dl

其中ds为均值滤波后的差值差异图，dl为中值滤波后的比值差异图，diff为差值差异图和比值差异图融合后的差异图。

步骤八的，具体方法如下：

第一步，调用内核程序中的kmeans_1k函数将diff中的全局变量读入到工作组的私有变量中，计算缓存对象diff中每个点到聚类中心{μ₀,μ₁}的距离，确定每个点所属的类c⁽ⁱ⁾，将聚类结果保存到缓存对象cluster中，该操作分布在若干个工作组中，每个工作组中的工作项为256个，每个工作项处理一个数据点，计算如下：

其中diff⁽ⁱ⁾为差异图上第i个点；

第二步，得到每个点所属类别后，调用kmeans_2k统计每类样本的个数以及对应的diff中所有样本的特征和，使多个工作项同时进行此操作，每个工作项负责256个数据点的统计，计算公式如下：

为类别j中所有样本的特征和，

为类别j中的样本个数；

第三步，调用kmeans_3k函数对每类统计后的个数和对应的diff中所有样本的特征和求平均值，得到新的聚类中心μ_j：

第四步，用新的聚类中心重复上述操作，在主机端控制迭代次数，当loop大于8次后迭代停止，将缓存对象cluster读取到主机内存中，在主机端对变化和非变化像素点分别标记为白色和黑色，将变化检测后的图像保存。

与现有技术相比，本发明采用CPU+FPGA异构并行加速实现SAR图像的变化检测，与现有技术中SAR图像变化检测的方法在CPU下计算相比，FPGA的并行计算能力要远远高于CPU数据串行计算能力，并且采用FPGA与CPU结合工作的方式，FPGA负责将整个算法中数据的处理并行化，CPU负责在主机端控制任务的调度以及参数的大小。本发明采用OpenCL的并行编程标准对SAR图像变化检测的算法进行并行化编写，代码的可移植性比较高，在不同的异构平台都可以运行。

附图说明

图1为本发明的实现流程框图；

图2为本发明对Ottawa地区的变化检测结果图；其中，(a)为Ottawa地区的A时间的图像；(b)为Ottawa地区的B时间的图像；(c)为用传统的图像融合方法进行变换检测后的结果图；(d)为用本发明方法进行变化检测后的结果图。

具体实施方式

下面结合附图和实施例对本发明做进一步说明。

实施例1

参照图1：本发明包括如下步骤：

步骤1：初始化开放运算语言OpenCL环境，选择异构计算平台以及计算设备，在计算平台上创建上下文，在上下文中创建命令队列，得到异构框架，创建并编译内核程序。

选取Altera异构计算平台，加速设备为FPGA，用OpenCL提供的应用程序编程接口API创建上下文以及命令队列，得到异构框架，使用SDK创建并编译内核程序。

步骤2：将初始图像读入到主机内存中，再将图像数据传递到FPGA的全局内存中。

选取两张Ottawa地区不同时刻的灰度图作为原始图像，将原始图像读入到主机内存中，将其变换成行向量后，写到缓存对象I₁和I₂中再将数据传递到FPGA的全局内存中。

步骤3：初始化lee滤波，均值滤波，中值滤波窗口尺寸分别为d1＝3，d2＝5，d3＝11，图像融合系数a＝0.2，以及kmeans聚类循环的次数loop＝8和两个聚类中心{μ₀,μ₁}＝{90,150}；

步骤4：调用内核程序中的LeeFilter_3函数对两幅图像分别进行lee滤波，滤波后的图像数据分别保存在缓存对象lee1和lee2中；

(4a)根据步骤3中设置的lee滤波的窗口尺寸d1＝3，对窗口区域内的像素求均值，每个工作项负责一个窗口内所有的数据计算。均值μ_Z：

(4b)根据(4a)中的均值计算窗口区域内像素方差：

(4c)根据μ_Z和计算ENL以及SX2：

(4d)计算lee滤波后的像素值xcap：

其中ω为滤波窗口宽度，Z_ij为原始图像(i,j)处的像素值，Z为3×3窗口的中心点像素值。

步骤5：根据缓存对象lee1和lee2调用内核程序中的diff函数生成差值差异图和比值差异图，分别保存在缓存对象ds和dl中；

(5a)调用缓存对象lee1和lee2生成差值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象ds中，计算公式：

ds(x,y)＝|lee1(x,y)-lee2(x,y)|

(5b)调用缓存对象lee1和lee2生成比值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象dl中，计算公式：

其中lee1(x,y)为第一幅图lee滤波后在坐标(x,y)处的值，lee2(x,y)为第二幅图lee滤波后在坐标(x,y)处的值，ds(x,y)为差值差异图在(x,y)的值，dl(x,y)为比值差异图在(x,y)的值。

步骤6：调用内核程序中的average函数对缓存对象ds进行均值滤波，以及median函数对缓存对象dl进行中值滤波；

(6a)对缓存对象ds，调用内核程序中的average函数进行均值滤波，滤波后保存在ds中。

(6b)对缓存对象dl，调用内核程序中的median函数进行中值滤波，滤波后保存在dl中。

步骤7：根据滤波后的缓存对象ds和dl，以及图像融合系数a调用内核程序中的add函数进行差异图的融合，结果保存在缓存对象diff中。此操作分布在若干个工作组中，每个工作组中工作项的大小为64×64，所有工作项并行工作，每个工作项处理一个数据，差异图融合公式为：

diff＝α*ds+(1-α)*dl

其中ds为均值滤波后的差值差异图，dl为中值滤波后的比值差异图，diff为差值差异图和比值差异图融合后的差异图。

步骤8：根据缓存对象diff以及循环迭代次数loop，在主机端循环调用内核程序kmeans中的函数，对diff进行聚类，将聚类后的结果保存在缓存对象cluster中。

(8a)调用内核程序中的kmeans_1k函数将diff中的全局变量读入到工作组的私有变量中，计算缓存对象diff中每个点到聚类中心{μ₀,μ₁}的距离，确定每个点所属的类c⁽ⁱ⁾，将聚类结果保存到缓存对象cluster中。该操作分布在若干个工作组中，每个工作组中的工作项为256个，每个工作项处理一个数据点，计算如下：

其中diff⁽ⁱ⁾为差异图上第i个点。

(8b)得到每个点所属类别后，调用kmeans_2k统计每类样本的个数以及对应的diff中所有样本的特征和，使多个工作项同时进行此操作，每个工作项负责256个数据点的统计。计算公式如下：

为类别j中所有样本的特征和，

为类别j中的样本个数；

(8c)调用kmeans_3k函数对每类统计后的个数和对应的diff中所有样本的特征和求平均值，得到新的聚类中心μ_j：

(8d)用新的聚类中心重复上述操作，在主机端控制迭代次数，当loop大于8次后迭代停止，将缓存对象cluster读取到主机内存中。在主机端对变化和非变化像素点分别标记为白色和黑色，将变化检测后的图像保存。

以下结合实验，对本发明的效果作进一步的描述。

1、实验环境：

异构平台为Altera异构开发平台，其中主机端CPU为Intel Core i5-3470，内存16G，FPGA板卡为DE5-NET，软件环境为Visual Studio 2015和Altera SDK for Open CL16.0，操作系统为64-bit Windows10。

2、实验内容：

实验1，对输入的Ottawa地区的两幅不同时间的图像分别用传统基于图像融合的方法和本发明进行变化检测操作，结果如图2。

图2(a)和图2(b)分别为Ottawa地区的不同时间的两张图像；

图2(c)为用传统的图像融合方法进行变换检测后的结果图；

图2(d)为用本发明方法进行变化检测后的结果图；

由图2可以看出，本发明方法与传统的图像融合方法对SAR进行变化检测的检测结果在视觉上效果相当。

实验2，对输入不同分辨率大小的SAR图像分别用传统图像融合方法和本发明方法进行变化检测操作，两种方法耗时的对比，对比结果如表1。

表1本发明与传统SAR图像变化检测的耗时(单位：s)

原始图像大小	290*350	301*301	400*400	4096*4096
					传统方法耗时	7.4097	7.3521	10.314	249.729
本发明方法耗时	0.8616	0.8549	1.3134	33.081

由表1可见本发明提出的基于CPU+FPGA的SAR图像变化检测方法相比于传统的基于图像融合变化检测的方法，耗时明显较少，对于分辨率较小的图可以做到实时检测。

综上所述，本发明基于CPU+FPGA的快速SAR图像变化检测方法，有效的解决了传统变化检测方法进行SAR图像变化检测时速度慢的问题，在保证原有的检测准确率前提下，加快检测速度，实现对SAR图像的快速变化检测。

Claims (8)

1.基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，包括以下步骤：

步骤一，初始化开放运算语言环境，选择异构计算平台以及计算设备，在计算平台上创建上下文，在上下文中创建命令队列，得到异构框架，创建并编译内核程序；

步骤二，将初始图像读入到主机内存中，再将图像数据传递到FPGA的全局内存中；

步骤三，初始化lee滤波，均值滤波；

步骤四，调用内核程序中的LeeFilter_3函数对两幅图像分别进行lee滤波，滤波后的图像数据分别保存在缓存对象lee1和缓存对象lee2中；

步骤五，根据缓存对象lee1和缓存对象lee2调用内核程序中的diff函数生成差值差异图和比值差异图，分别保存在缓存对象ds和缓存对象dl中；

步骤六，调用内核程序中的average函数对缓存对象ds进行均值滤波，以及median函数对缓存对象dl进行中值滤波；

步骤七，根据滤波后的缓存对象ds和缓存对象dl，以及图像融合系数a调用内核程序中的add函数进行差异图的融合，结果保存在缓存对象diff中；

步骤八，根据缓存对象diff以及循环迭代次数loop，在主机端循环调用内核程序kmeans中的函数，对diff进行聚类，将聚类后的结果保存在缓存对象cluster中；

步骤九，将缓存对象cluster通过PCI-E总线传输到CPU主机内存中，在主机端对变化和非变化区域分别标记为白色和黑色，将变化检测后的图像保存，完成图像变化的检测。

2.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤一中，开放运算语言环境为OpenCL环境，异构计算平台采用Altera，用OpenCL提供的应用程序编程接口API创建上下文以及命令队列，得到异构框架，使用SDK创建并编译内核程序。

3.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤二的具体方法如下，选取两张Ottawa地区不同时刻的灰度图作为原始图像，将原始图像读入到主机内存中，将其变换成行向量后，分别写到缓存对象I₁和I₂中再将数据传递到FPGA的全局内存中。

4.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤四的具体方法如下：

第一步，根据步骤三中设置的lee滤波的窗口尺寸d1，对窗口区域内的像素求均值，每个工作项负责一个窗口内所有的数据计算；均值μ_Z：

第二步，根据第一步中的均值计算窗口区域内像素方差：

第三步，根据μ_Z和计算ENL以及SX2：

第四步，计算lee滤波后的像素值xcap：

其中ω为滤波窗口宽度，Z_ij为原始图像(i,j)处的像素值，Z为3×3窗口的中心点像素值。

5.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤五的具体方法如下：

第一步，调用缓存对象lee1和lee2生成差值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象ds中，计算公式：

ds(x,y)＝|lee1(x,y)-lee2(x,y)|

第二步，调用缓存对象lee1和lee2生成比值差异图，此算法操作分布在若干个工作项上，每个工作组内工作项的大小为64×64，每个工作项只负责一个数据的操作，将结果保存到缓存对象dl中，计算公式：

其中lee1(x,y)为第一幅图lee滤波后在坐标(x,y)处的值，lee2(x,y)为第二幅图lee滤波后在坐标(x,y)处的值，ds(x,y)为差值差异图在(x,y)的值，dl(x,y)为比值差异图在(x,y)的值。

6.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤六的具体方法如下：

第一步，对缓存对象ds，调用内核程序中的average函数进行均值滤波，滤波后保存在ds中；

第二步，对缓存对象dl，调用内核程序中的median函数进行中值滤波，滤波后保存在dl中。

7.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤七中，调用内核程序中的add函数进行差异图融合，此操作分布在若干个工作组中，每个工作组中工作项的大小为64×64，所有工作项并行工作，每个工作项处理一个数据，差异图融合公式为：

diff＝α*ds+(1-α)*dl

其中ds为均值滤波后的差值差异图，dl为中值滤波后的比值差异图，diff为差值差异图和比值差异图融合后的差异图。

8.根据权利要求1所述的一种基于CPU+FPGA的快速SAR图像变化检测方法，其特征在于，步骤八的，具体方法如下：

第一步，调用内核程序中的kmeans_1k函数将diff中的全局变量读入到工作组的私有变量中，计算缓存对象diff中每个点到聚类中心{μ₀,μ₁}的距离，确定每个点所属的类c⁽ⁱ⁾，将聚类结果保存到缓存对象cluster中，该操作分布在若干个工作组中，每个工作组中的工作项为256个，每个工作项处理一个数据点，计算如下：

其中diff⁽ⁱ⁾为差异图上第i个点；

第二步，得到每个点所属类别后，调用kmeans_2k统计每类样本的个数以及对应的diff中所有样本的特征和，使多个工作项同时进行此操作，每个工作项负责256个数据点的统计，计算公式如下：

为类别j中所有样本的特征和，

为类别j中的样本个数；

第三步，调用kmeans_3k函数对每类统计后的个数和对应的diff中所有样本的特征和求平均值，得到新的聚类中心μ_j：

第四步，用新的聚类中心重复上述操作，在主机端控制迭代次数，当loop大于8次后迭代停止，将缓存对象cluster读取到主机内存中，在主机端对变化和非变化像素点分别标记为白色和黑色，将变化检测后的图像保存。