并行Shi-Tomasi图像角点特征检测算法
作者： 李雨嫣 陈华
摘 要：在大规模图像的处理中，串行角点检测算法存在着运算量大、耗时长的问题。基于API、OpenMP及PPL多核CPU技术，提出了三种改进的并行角点检测算法。实验结果显示，三种并行算法对不同尺寸的图片均具有较好的加速效果。此外，实验采取不同的线程数量进行测试。基于API的并行检测算法加速比在四线程情况下平均可达3.02，在八线程情况下平均可达3.85。基于OpenMP的并行检测算法加速比在四线程情况下平均可达2.26，在八线程情况下平均可达3.30。基于API和OpenMP的并行检测算法均具有良好的多核可扩展性。三种并行算法中，基于API的并行检测算法具有最高最稳定的加速效果。
关键词：Shi-Tomasi; 角点检测; API; OpenMP; PPL
中图分类号：TP391 文献标识码：A
文章编号：1009-3044（2022）06-0100-05
开放科学（资源服务）标识码（OSID）：
图像中关键区域的形状由角点所决定，角点能够反应图像中重要的特征信号，因而成为图像的重要局部特征之一[1]。角点在影像重建[2-3]、图像配准[4]和目标跟踪[5-6]等多种计算机视觉处理任务中均有着广泛的应用。
对角点的定义一般可以分为三种，即图像边界曲线具有极大曲率值的点、图像中梯度值和梯度变化率很高的点和图像边界方向变化不连续的点[7]。针对不同的定义，目前存在不同的角点检测算法[8]。
1988年，Chris Harris和Mike Stephens两人在H.Moravec算法的基础上，应用邻近像素点灰度差值的概念对角点进行判断，提出了基于自相关矩阵的角点提取算法，称为Harris算法[9]。1993年，Jianbo Shi和Carlo Tomasi两人对Harris算法中的分数公式进行了改进，提出了Shi-Tomasi角点检测算法[10]。Shi-Tomasi角点检测算法是基于图像灰度的一种检测方法，该方法通过计算每个像素点的曲率以及梯度进行角点检测，可以有效地避免由图像边缘编码带来的分割和边缘提取问题，是目前常用的角点检测算法之一。
为了解决传统串行算法中针对大规模图像运算耗时长的问题，肖汉等[11]提出了一种基于多GPU的Harris角点检测并行算法。郭志全[12]提出基于数据并行和任务并行的两种Harris角点检测算法的并行设计思路，并通过SMT-PAAG仿真器进行实现。王俊超[13]采用对图像进行网格化分块的方法，使用OpenMP对Shi-Tomasi角点检测算法进行了并行化设计。朱超[14]等分别基于CPU和GPU对Harris角点检测算法进行了并行化设计。
本文研究的目的是采用API，OpenMP和PPL三种并行CPU加速技术，改善Shi-Tomasi角点检测算法耗时长的问题。
1 串行算法描述
1.1 Shi-Tomasi角点检测算法原理
将一个窗口置于图像之上，对其进行移动。设图像I（x，y）在点（x，y）处以偏移量（u，v）进行移动后，产生的灰度差值E（x，y，u，v）为：
[Ex，y，u，v=x，y∈Swx，yIx+u，y+v-Ix，y2] （1）
（1） 式中，[S]是移动窗口的区域，w（x，y）是高斯加权函数。使用泰勒展开公式对I（x+u，y+v）进行近似，得到下式：
[Ix+u，y+v=Ix，y+∂I∂xu+∂I∂yv+Ou2，v2] （2）
将（2）式代入（1）式当中可得：
[Ex，y，u，v=u，vx，y∈Swx，yI2xIxIyIxIyI2yuv] （3）
设[M]矩阵为：
[M=x，y∈Swx，yI2xIxIyIxIyI2y] （4）
代入可得：
[Ex，y，u，v=u，vMuv] （5）
（5） 式中，M矩阵是关于x和y的二阶函数，根据其特征值的不同，可以将像素点分为角点、直线和平面三种情况。若该点是角点，则因其移动窗口朝任意方向移动都会导致窗口内灰度产生较大的变化，故M矩阵的两个特征值都较大。
在Shi-Tomasi算法当中，采用一个角点响应值R作为判断角点质量的依据。当R值大于设定的阈值时，判断该像素点为一个角点[9]。
[R=minλ1，λ2] （6）
2 并行算法设计
算法1 基于Windows API的Shi-Tomasi并行算法
输入：图像
输出：标记出角点的特征图像
算法描述：采取数据并行的思想，将输入图像切割为n等分的子区域，每一个区域由一个线程进行处理。每一个线程完成所分配区域的初始化、梯度值计算、相关矩阵M的计算、角点响应值R的计算、极大值抑制、Shi-Tomasi角点响应的判断及角点标记。定义结构体picture存储子区域的相关信息，用于主线程向从线程传递参数，使用CreateThread指令创建从线程。以四线程为例，该方法的并行计算流程如图1所示。
伪代码如下：
读取图像;
图像灰度化;
分割原始图像和灰度图像;
定义picture结构体变量pi，并将分割子区域信息存储到结构体中;
启动计时器;
Shitomasi_demo（原始图像，灰度图像，输出窗口句柄，角点个数）;//Shi-Tomasi角点检测函数
结束计时器;
计算并输出串行算法时间;
启动计时器;
hThread[i] = CreateThread（0， 0， threadHarris， （LPVOID）&pi， 0， NULL）;//创建从线程
结束计时器;
计算并输出并行计算时间;
计算并输出加速比;
算法2 基于OpenMP的Shi-Tomasi并行算法
输入：图像
输出：标记出角点的特征图像
算法描述：本文在FOLK/JOIN模型的基础上进行设计。采取数据并行的思想，将输入图像切割为n等分的子区域，每一个区域由一个线程进行处理。每一个线程完成所分配区域的初始化，梯度值计算，相关矩阵M的计算，角点响应值R的计算，极大值抑制，Shi-Tomasi角点响应的判断以及角点标记。使用sections和section指令创建多线程，将每个子区域分配到不同的线程中并行执行。以四线程为例，该方法的并行计算流程如图2所示。
伪代码如下：
读取图像;
图像灰度化;
分割原始图像和灰度图像;
启动计时器;
Shitomasi_demo（原始图像， 灰度图像， 输出窗口句柄， 角点个数）;//Shi-Tomasi角点检测函数
结束计时器;
计算并输出串行算法时间;
启动计时器;
#pragma omp parallel sections
{
#pragma omp section
Shitomasi_demo（原图像分割区域1， 灰度图像分割区域1， 输出窗口句柄， 角点个数/4）;
#pragma omp section
Shitomasi_demo（原图像分割区域2， 灰度图像分割区域2， 输出窗口句柄， 角点个数/4）;
#pragma omp section
Shitomasi_demo（原图像分割区域3， 灰度图像分割区域3， 输出窗口句柄， 角点个数/4）;
#pragma omp section
Shitomasi_demo（原图像分割区域4， 灰度图像分割区域4， 输出窗口句柄， 角点个数/4）;
}//创建从线程，分配计算任务
结束计时器;
计算并输出并行计算时间;
计算并输出加速比;
算法3 基于PPL的Shi-Tomasi并行算法
输入：图像
输出：标记出角点的特征图像
算法描述：采取数据并行的思想，将输入图像切割为n等分的子区域，每一个区域由一个线程进行处理。每一个线程完成所分配区域的初始化、梯度值计算、相关矩阵M的计算、角点响应值R的计算、极大值抑制、Shi-Tomasi角点响应的判断以及角点标记。使用task_group和parallel_invoke指令[15]将计算任务分给不同的线程。以四线程为例，该方法的并行计算流程如图3所示。
伪代码如下：
读取图像;
图像灰度化;
分割原始图像和灰度图像;
启动计时器;
Shitomasi_demo（原始图像， 灰度图像， 输出窗口句柄， 角点个数）;//Shi-Tomasi角点检测函数
结束计时器;
计算并输出串行算法时间;
启动计时器;
task_group tasks;//构造task_group对象
parallel_invoke（
[&]{ tasks.run（[&]{ Harris_Demo（原图像分割区域1， 灰度图像分割区域1， 输出窗口句柄， 角点个数/4）; }）; }，
[&]{ tasks.run（[&]{ Harris_Demo（原图像分割区域1， 灰度图像分割区域1， 输出窗口句柄， 角点个数/4）; }）; }，
[&]{ tasks.run（[&]{ Harris_Demo（原图像分割区域1， 灰度图像分割区域1， 输出窗口句柄， 角点个数/4）; }）; }，
[&]{ tasks.run（[&]{ Harris_Demo（原图像分割区域1， 灰度图像分割区域1， 输出窗口句柄， 角点个数/4）; }）; }
）;//创建从线程，将任务添加到task_group对象中并等待
tasks.wait（）;//同步
结束计时器;
计算并输出并行计算时间;
计算并输出加速比;
3 实验测试与结果分析
3.1 实验环境
硬件：CPU（4核）：Intel（R） Core（TM） i5-8300H CPU @ 2.30GHz
软件：操作系统：Windows 10
使用工具：Visual Studio2013
3.2 实验测试用例
部分测试用例如图4、图5所示：