小小博客
0%

SLIC超像素分割算法的并行计算实现

发表于

SLIC 超像素分割算法的并行计算实现

在图像处理领域，超像素分割作为一种重要的预处理技术，能够将图像划分为具有语义意义的像素簇，为后续的目标检测、图像分割等任务奠定基础。

SLIC 算法由 Radhakrishna Achanta 等人于 2012 年提出，其核心思想是通过线性迭代聚类将图像分割为具有相似颜色和空间位置的超像素。与传统算法相比，SLIC 具有计算效率高、参数少（近乎零参数）的特点，尤其适合实时性要求高的场景。

算法的主要步骤包括：

  1. 颜色空间转换：将 RGB 图像转换到 LAB 颜色空间，LAB 空间的欧氏距离更符合人眼感知
  2. 种子点初始化：在图像中均匀采样种子点，并基于边缘强度调整位置
  3. 迭代聚类：以种子点为中心，计算邻域像素的颜色距离和空间距离，更新像素所属聚类
  4. 连通性增强：消除孤立小区域，保证超像素的连通性

算法优化策略与实现细节

编译优化：开启 O3 选项

编译优化是最直接的性能提升手段。通过 GCC 的-O3选项，编译器会执行一系列优化操作，包括指令调度、循环展开、常量折叠等。

1
g++ -std=c++11 -O3 -fopenmp -mavx -mavx2 -march=native bestgai.cpp -o slictest

AVX 向量化

AVX（Advanced Vector Extensions）指令集支持单指令多数据操作，能同时处理多个数据元素。在 SLIC 算法中，颜色空间转换部分是向量化优化的重点：

1
2
3
4
5
6
7
8
9
10
11
12
// RGB到LAB转换的AVX向量化实现
void SLIC::RGB2XYZ(const int& sR, const int& sG, const int& sB, double& X, double& Y, double& Z) {
    // 使用AVX指令同时处理4个数据元素
    __m256d rgb = _mm256_set_pd(0.0, sB, sG, sR);
    __m256d mmm = _mm256_set1_pd(255.0);
    __m256d result = _mm256_div_pd(rgb, mmm);
    double resultArray[4];
    _mm256_storeu_pd(resultArray, result);
    double R = resultArray[0], G = resultArray[1], B = resultArray[2];
    
    // 后续计算...
}

向量化的核心在于将标量操作转换为向量操作，充分利用 CPU 的 SIMD 单元，减少指令周期数。

OpenMP 并行：多线程分工加速迭代计算

SLIC 算法中的迭代聚类过程具有天然的并行性，通过 OpenMP 实现多线程并行计算：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
// 并行化RGB到LAB转换
void SLIC::DoRGBtoLABConversion(const unsigned int*& ubuff, double*& lvec, double*& avec, double*& bvec) {
    int sz = m_width * m_height;
    lvec = new double[sz];
    avec = new double[sz];
    bvec = new double[sz];
    
    // 使用parallel for指令并行化循环
    #pragma omp parallel for
    for (int j = 0; j < sz; j++) {
        int r = (ubuff[j] >> 16) & 0xFF;
        int g = (ubuff[j] >> 8) & 0xFF;
        int b = (ubuff[j]) & 0xFF;
        RGB2LAB(r, g, b, lvec[j], avec[j], bvec[j]);
    }
}

更复杂的并行优化体现在迭代聚类的分块处理 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
void SLIC::PerformSuperpixelSegmentation_VariableSandM(...) {
    int xiancheng = omp_get_max_threads();
    int group = sz / xiancheng;
    
    // 分块处理图像行
    #pragma omp parallel
    {
        int qq = omp_get_thread_num();
        int start = qq * group;
        int endle = start + group;
        if (qq == xiancheng - 1) endle = m_height;
        
        #pragma omp parallel for
        for (int y = start; y < endle; y++) {
            for (int n = 0; n < numk; n++) {
                // 计算当前种子点的影响区域
                int wq1 = (int)(kseedsy[n] - offset);
                int wq2 = (int)(kseedsy[n] + offset);
                if (y >= wq1 && y < wq2) {
                    int x1 = max(0, (int)(kseedsx[n] - offset));
                    int x2 = min(m_width, (int)(kseedsx[n] + offset));
                    
                    for (int x = x1; x < x2; x++) {
                        // 计算颜色距离和空间距离
                        int i = y * m_width + x;
                        double l = m_lvec[i], a = m_avec[i], b = m_bvec[i];
                        double distlab = (l - kseedsl[n])*(l - kseedsl[n]) + 
                                        (a - kseedsa[n])*(a - kseedsa[n]) + 
                                        (b - kseedsb[n])*(b - kseedsb[n]);
                        double distxy = (x - kseedsx[n])*(x - kseedsx[n]) + 
                                       (y - kseedsy[n])*(y - kseedsy[n]);
                        double dist = distlab / maxlab[n] + distxy * invxywt;
                        
                        if (dist < distvec[i]) {
                            distvec[i] = dist;
                            klabels[i] = n;
                        }
                    }
                }
            }
        }
    }
}

循环调整

调整循环顺序、提取公共计算 

1
2
3
4
5
6
7
8
9
10
11
for (column = 0; column < 100; column++) {
    for (row = 0; row < 5; row++) {
        sum += table[row][column];
    }
}

// 优化后：按列遍历，提升缓存命中率
for (row = 0; row < 5; row++) {
    for (column = 0; column < 100; column++) {
        sum += table[row][column];
    }

将图像遍历从按行改为按列，并结合分块并行，进一步提升了数据局部性： 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
int group1 = m_height / xiancheng;
#pragma omp parallel
{
    int qq = omp_get_thread_num();
    int start = qq * group1;
    int endle = start + group1;
    if (qq == xiancheng - 1) endle = m_height;
    
    #pragma omp parallel for
    for (int y = start; y < endle; y++) {
        for (int n = 0; n < numk; n++) {
            // 种子点邻域计算...
        }
    }
}

按需计算

将使用频率低的数组改为按需计算，避免全程存储：

1
2
3
4
5
6
7
8
9
10
11
vector<double> distxy(sz, DBL_MAX);
vector<double> distlab(sz, DBL_MAX);
vector<double> distvec(sz, DBL_MAX);

// 优化后：按需计算，仅存储必要结果
double distlab = (l - kseedsl[n])*(l - kseedsl[n]) + 
                (a - kseedsa[n])*(a - kseedsa[n]) + 
                (b - kseedsb[n])*(b - kseedsb[n]);
double distxy = (x - kseedsx[n])*(x - kseedsx[n]) + 
               (y - kseedsy[n])*(y - kseedsy[n]);
double dist = distlab / maxlab[n] + distxy * invxywt;

性能优化结果

以下是各优化步骤的性能对比：

优化策略 计算时间（ms） 加速比（%）
基准版本 25120 0.00
编译优化 - O3 8600 293.49
AVX 向量化 8539 294.71
OpenMP 并行 1560 1610.25
按需计算 1436 1749.30
循环层次调整 388 6474.22
最终优化版本 350 -