找回密码
 立即注册
首页 业界区 业界 高性能计算实践-OpenCV图像矩阵转置 transpose SIMD加速 ...

高性能计算实践-OpenCV图像矩阵转置 transpose SIMD加速(ippicv)复现

觞刈 2025-11-19 03:05:03
说明

矩阵转置是高性能计算中的经典问题。OpenCV 的 transpose 函数内部依赖 ippicv 库中的 ippiTranspose_8u_C1R 实现。本文将对该优化算法进行复现与分析。
与上一篇基于 cv::flip / ippiMirror 的图像翻转不同,矩阵转置不再是简单的行内倒序,而是将整幅图像在行列维度上重新映射。我们可以用块划分(tiling)的遍历方式来解决,同时加上各种优化技巧。
复现
  1. #ifdef _MSC_VER
  2. #define FORCE_INLINE __forceinline
  3. #elif defined(__GNUC__)
  4. #define FORCE_INLINE __attribute__((always_inline)) inline
  5. #else
  6. #define FORCE_INLINE inline
  7. #endif
  8. /**
  9. * 8x8 SSE 转置微核
  10. * 读取 8 行源数据 -> 转置 -> 连续写入 64 字节到 buffer
  11. */
  12. FORCE_INLINE void transpose_8x8_store_contiguous(const uint8_t* src0,
  13.                                                  const uint8_t* src1,
  14.                                                  const uint8_t* src2,
  15.                                                  const uint8_t* src3,
  16.                                                  const uint8_t* src4,
  17.                                                  const uint8_t* src5,
  18.                                                  const uint8_t* src6,
  19.                                                  const uint8_t* src7,
  20.                                                  uint8_t* pDst) {
  21.     __m128i r0 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src0));
  22.     __m128i r1 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src1));
  23.     __m128i r2 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src2));
  24.     __m128i r3 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src3));
  25.     __m128i r4 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src4));
  26.     __m128i r5 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src5));
  27.     __m128i r6 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src6));
  28.     __m128i r7 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(src7));
  29.     __m128i t0 = _mm_unpacklo_epi8(r0, r1);
  30.     __m128i t1 = _mm_unpacklo_epi8(r2, r3);
  31.     __m128i t2 = _mm_unpacklo_epi8(r4, r5);
  32.     __m128i t3 = _mm_unpacklo_epi8(r6, r7);
  33.     __m128i t4 = _mm_unpacklo_epi16(t0, t1);
  34.     __m128i t5 = _mm_unpacklo_epi16(t2, t3);
  35.     __m128i t6 = _mm_unpackhi_epi16(t0, t1);
  36.     __m128i t7 = _mm_unpackhi_epi16(t2, t3);
  37.     __m128i c0 = _mm_unpacklo_epi32(t4, t5);
  38.     __m128i c1 = _mm_unpackhi_epi32(t4, t5);
  39.     __m128i c2 = _mm_unpacklo_epi32(t6, t7);
  40.     __m128i c3 = _mm_unpackhi_epi32(t6, t7);
  41.     // 将转置后的 8x8 块 (64字节) 连续写入 buffer
  42.     // buffer 是 alignas(64) 的,始终对齐
  43.     _mm_store_si128(reinterpret_cast<__m128i*>(pDst + 0), c0);
  44.     _mm_store_si128(reinterpret_cast<__m128i*>(pDst + 16), c1);
  45.     _mm_store_si128(reinterpret_cast<__m128i*>(pDst + 32), c2);
  46.     _mm_store_si128(reinterpret_cast<__m128i*>(pDst + 48), c3);
  47. }
  48. /**
  49. * 64x64 Tile 转置核心优化
  50. * 使用 64x64 栈上缓存
  51. * 1. 读 Src (8x8 块),转置后线性写入 Tmp (Row-Major Block)
  52. * 2. 读 Tmp (Strided),合并后流式写入 Dst (Contiguous Rows)
  53. */
  54. template <bool UseStream>
  55. FORCE_INLINE void
  56. transpose_64x64_tile_impl(const uint8_t* pSrc, unsigned int srcStep, uint8_t* pDst, unsigned int dstStep) {
  57.     // 64x64 临时 Buffer
  58.     alignas(64) uint8_t tmp[64 * 64];
  59.     uint8_t* tmpPtr = tmp;
  60.     // 1. 读取源并填充 Buffer
  61.     // 策略:保持源图像的线性访问 (Y then X),这对性能至关重要
  62.     // 结果:tmp 中的块是按 "Row-Major Block" 顺序排列的
  63.     // 即:[B(0,0)] [B(0,1)] ... [B(0,7)] [B(1,0)] ...
  64.     // 预计算步长指针,减少循环内乘法
  65.     size_t srcStep8 = (size_t)srcStep * 8;
  66.     const uint8_t* s0 = pSrc;
  67.     for (int y = 0; y < 64; y += 8) {
  68.         const uint8_t* r0 = s0;
  69.         const uint8_t* r1 = s0 + srcStep;
  70.         const uint8_t* r2 = s0 + srcStep * 2;
  71.         const uint8_t* r3 = s0 + srcStep * 3;
  72.         const uint8_t* r4 = s0 + srcStep * 4;
  73.         const uint8_t* r5 = s0 + srcStep * 5;
  74.         const uint8_t* r6 = s0 + srcStep * 6;
  75.         const uint8_t* r7 = s0 + srcStep * 7;
  76.         for (int x = 0; x < 64; x += 8) {
  77.             transpose_8x8_store_contiguous(r0 + x, r1 + x, r2 + x, r3 + x, r4 + x, r5 + x, r6 + x, r7 + x, tmpPtr);
  78.             tmpPtr += 64; // buffer 线性写入
  79.         }
  80.         s0 += srcStep8;
  81.     }
  82.     // 2. 从 Buffer 读取并流式写入 Dst
  83.     // 目标:写入 Dst 的行
  84.     // Dst 的第 i 个条带 (由8行组成) 对应 Source 的第 i 个块列
  85.     // Source 的块列 i 包含块:B(0,i), B(1,i), ... B(7,i)
  86.     // 在 Row-Major 的 tmp 中,这些块的内存地址不是连续的,而是相隔 8个块 (8*64 = 512字节)
  87.     // 外层循环:遍历 8 个垂直条带 (Strip),对应 tmp 中的 Block Column 0..7
  88.     for (int colBlock = 0; colBlock < 8; ++colBlock) {
  89.         // 当前条带中 B(0, colBlock) 的起始地址
  90.         // 在 tmp 中,Block(row, col) 的索引是 row*8 + col
  91.         // Block(0, colBlock) 的偏移是 colBlock * 64
  92.         const uint8_t* bBase = tmp + colBlock * 64;
  93.         // 处理条带内的 8 行
  94.         for (int r = 0; r < 8; ++r) {
  95.             // 我们需要从 8 个垂直堆叠的块中,分别取出第 r 行
  96.             // Block stride = 512 bytes.
  97.             // Row offset inside block = r * 8 bytes.
  98.             int laneOffset = r * 8;
  99.             // 从 tmp 中以 512 字节 stride 读取
  100.             __m128i b0 = _mm_loadl_epi64((const __m128i*)(bBase + 0 * 512 + laneOffset));
  101.             __m128i b1 = _mm_loadl_epi64((const __m128i*)(bBase + 1 * 512 + laneOffset));
  102.             __m128i b2 = _mm_loadl_epi64((const __m128i*)(bBase + 2 * 512 + laneOffset));
  103.             __m128i b3 = _mm_loadl_epi64((const __m128i*)(bBase + 3 * 512 + laneOffset));
  104.             __m128i b4 = _mm_loadl_epi64((const __m128i*)(bBase + 4 * 512 + laneOffset));
  105.             __m128i b5 = _mm_loadl_epi64((const __m128i*)(bBase + 5 * 512 + laneOffset));
  106.             __m128i b6 = _mm_loadl_epi64((const __m128i*)(bBase + 6 * 512 + laneOffset));
  107.             __m128i b7 = _mm_loadl_epi64((const __m128i*)(bBase + 7 * 512 + laneOffset));
  108.             __m128i v0 = _mm_unpacklo_epi64(b0, b1);
  109.             __m128i v1 = _mm_unpacklo_epi64(b2, b3);
  110.             __m128i v2 = _mm_unpacklo_epi64(b4, b5);
  111.             __m128i v3 = _mm_unpacklo_epi64(b6, b7);
  112.             // 计算目标地址:
  113.             // 当前是第 colBlock 个条带,第 r 行 -> 全局行 colBlock*8 + r
  114.             uint8_t* dstRowPtr = pDst + (colBlock * 8 + r) * dstStep;
  115.             if (UseStream) { // 编译期优化,生成无分支代码
  116.                 // Stream 路径:要求 dstRowPtr 必须 16 字节对齐
  117.                 // 适用于 dstStep % 16 == 0 且 pDst 对齐的情况
  118.                 _mm_stream_si128(reinterpret_cast<__m128i*>(dstRowPtr + 0), v0);
  119.                 _mm_stream_si128(reinterpret_cast<__m128i*>(dstRowPtr + 16), v1);
  120.                 _mm_stream_si128(reinterpret_cast<__m128i*>(dstRowPtr + 32), v2);
  121.                 _mm_stream_si128(reinterpret_cast<__m128i*>(dstRowPtr + 48), v3);
  122.             } else {
  123.                 // StoreU 路径:安全处理任意对齐,且依然是 SIMD 向量化
  124.                 // 适用于 dstStep % 16 != 0 的情况
  125.                 _mm_storeu_si128(reinterpret_cast<__m128i*>(dstRowPtr + 0), v0);
  126.                 _mm_storeu_si128(reinterpret_cast<__m128i*>(dstRowPtr + 16), v1);
  127.                 _mm_storeu_si128(reinterpret_cast<__m128i*>(dstRowPtr + 32), v2);
  128.                 _mm_storeu_si128(reinterpret_cast<__m128i*>(dstRowPtr + 48), v3);
  129.             }
  130.         }
  131.     }
  132. }
  133. /**
  134. * 处理边缘的小块 (8x8 fallback)
  135. * 将 8x8 源块 (srcStep) 转置写入 8x8 目标块 (dstStep)
  136. */
  137. FORCE_INLINE void
  138. transpose_8x8_u8_to_strided(const uint8_t* pSrc, unsigned int srcStep, uint8_t* pDst, unsigned int dstStep) {
  139.     __m128i r0 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 0 * srcStep));
  140.     __m128i r1 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 1 * srcStep));
  141.     __m128i r2 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 2 * srcStep));
  142.     __m128i r3 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 3 * srcStep));
  143.     __m128i r4 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 4 * srcStep));
  144.     __m128i r5 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 5 * srcStep));
  145.     __m128i r6 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 6 * srcStep));
  146.     __m128i r7 = _mm_loadl_epi64(reinterpret_cast<const __m128i*>(pSrc + 7 * srcStep));
  147.     __m128i t0 = _mm_unpacklo_epi8(r0, r1);
  148.     __m128i t1 = _mm_unpacklo_epi8(r2, r3);
  149.     __m128i t2 = _mm_unpacklo_epi8(r4, r5);
  150.     __m128i t3 = _mm_unpacklo_epi8(r6, r7);
  151.     __m128i t4 = _mm_unpacklo_epi16(t0, t1);
  152.     __m128i t5 = _mm_unpacklo_epi16(t2, t3);
  153.     __m128i t6 = _mm_unpackhi_epi16(t0, t1);
  154.     __m128i t7 = _mm_unpackhi_epi16(t2, t3);
  155.     __m128i c0 = _mm_unpacklo_epi32(t4, t5);
  156.     __m128i c1 = _mm_unpackhi_epi32(t4, t5);
  157.     __m128i c2 = _mm_unpacklo_epi32(t6, t7);
  158.     __m128i c3 = _mm_unpackhi_epi32(t6, t7);
  159.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 0 * dstStep), c0);
  160.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 1 * dstStep), _mm_srli_si128(c0, 8));
  161.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 2 * dstStep), c1);
  162.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 3 * dstStep), _mm_srli_si128(c1, 8));
  163.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 4 * dstStep), c2);
  164.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 5 * dstStep), _mm_srli_si128(c2, 8));
  165.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 6 * dstStep), c3);
  166.     _mm_storel_epi64(reinterpret_cast<__m128i*>(pDst + 7 * dstStep), _mm_srli_si128(c3, 8));
  167. }
  168. /**
  169. * 核心转置内核,处理任意 WxH 块
  170. * 内部使用 64x64 Tile 优化,并处理 8x8 和 1x1 边缘
  171. */
  172. template <bool UseStream>
  173. int64_t icv_y8_owniTransposeWxH_8uC1_impl(const uint8_t* pSrc,
  174.                                           unsigned int srcStep,
  175.                                           uint8_t* pDst,
  176.                                           unsigned int dstStep,
  177.                                           int width,
  178.                                           int height) {
  179.     if (width <= 0 || height <= 0)
  180.         return 0;
  181.     constexpr int TILE = 64;
  182.     constexpr int MICRO = 8;
  183.     const int wMain = width & ~(TILE - 1);  // 64x 块主区域
  184.     const int hMain = height & ~(TILE - 1); // 64x 块主区域
  185.     // 1. 主循环 64x64 Tile (使用模板参数选择优化策略)
  186.     for (int y = 0; y < hMain; y += TILE) {
  187.         for (int x = 0; x < wMain; x += TILE) {
  188.             // Source Tile (x, y) 转置后写入 Dst Tile (y, x)
  189.             const uint8_t* srcTile = pSrc + y * srcStep + x;
  190.             uint8_t* dstTile = pDst + x * dstStep + y;
  191.             transpose_64x64_tile_impl<UseStream>(srcTile, srcStep, dstTile, dstStep);
  192.         }
  193.     }
  194.     // 2. 边缘处理 (通用代码,不依赖 UseStream,因为 storel 总是安全的)
  195.     // 高度为 hMain,宽度为 wTail
  196.     const int wTail = width - wMain;
  197.     if (wTail > 0) {
  198.         int wTailMain = wTail & ~(MICRO - 1); // 8x 块区域
  199.         int wTailTail = wTail - wTailMain;    // 1x 标量区域
  200.         for (int y = 0; y < hMain; y += MICRO) {
  201.             const uint8_t* srcRow = pSrc + y * srcStep + wMain;
  202.             uint8_t* dstCol = pDst + wMain * dstStep + y;
  203.             int xOff = 0;
  204.             // 8x8 块
  205.             for (; xOff < wTailMain; xOff += MICRO) {
  206.                 transpose_8x8_u8_to_strided(srcRow + xOff, srcStep, dstCol + xOff * dstStep, dstStep);
  207.             }
  208.             // 标量补齐
  209.             // (y, xOff) -> (xOff, y)
  210.             for (int k = 0; k < MICRO; ++k) { // 遍历 8 行
  211.                 for (int x = 0; x < wTailTail; ++x) {
  212.                     dstCol[(xOff + x) * dstStep + k] = srcRow[k * srcStep + (xOff + x)];
  213.                 }
  214.             }
  215.         }
  216.     }
  217.     // 3. 处理底部边缘 (Height non-64, 左侧部分)
  218.     // 高度为 hBottomTail,宽度为 wMain
  219.     const int hBottomTail = height - hMain;
  220.     if (hBottomTail > 0) {
  221.         int hBottomMain = hBottomTail & ~(MICRO - 1); // 8x 块区域
  222.         int hBottomTailTail = hBottomTail - hBottomMain; // 1x 标量区域
  223.         for (int x = 0; x < wMain; x += MICRO) {
  224.             const uint8_t* srcCol = pSrc + hMain * srcStep + x;
  225.             uint8_t* dstRow = pDst + x * dstStep + hMain;
  226.             int yOff = 0;
  227.             // 8x8 块
  228.             for (; yOff < hBottomMain; yOff += MICRO) {
  229.                 transpose_8x8_u8_to_strided(srcCol + yOff * srcStep, srcStep, dstRow + yOff, dstStep);
  230.             }
  231.             // 标量补齐
  232.             // (yOff, k) -> (k, yOff)
  233.             for (int k = 0; k < MICRO; ++k) { // 遍历 8 列
  234.                 for (int y = 0; y < hBottomTailTail; ++y) {
  235.                     dstRow[k * dstStep + (yOff + y)] = srcCol[(yOff + y) * srcStep + k];
  236.                 }
  237.             }
  238.         }
  239.     }
  240.     // 4. 处理右下角 (wTail x hBottomTail)
  241.     if (wTail > 0 && hBottomTail > 0) {
  242.         // C++ 标量实现
  243.         const uint8_t* srcCorner = pSrc + hMain * srcStep + wMain;
  244.         uint8_t* dstCorner = pDst + wMain * dstStep + hMain;
  245.         for (int y = 0; y < hBottomTail; ++y) {
  246.             for (int x = 0; x < wTail; ++x) {
  247.                 dstCorner[x * dstStep + y] = srcCorner[y * srcStep + x];
  248.             }
  249.         }
  250.     }
  251.     // 如果使用了 Stream (NT Store),需要 sfence 确保数据可见性
  252.     if (UseStream) {
  253.         _mm_sfence();
  254.     }
  255.     return 0;
  256. }
  257. /**
  258. * 核心转置内核 Dispatcher
  259. * 根据 dstStep 和 pDst 的对齐情况,分发到 Stream 版或 StoreU 版
  260. */
  261. int64_t icv_y8_owniTransposeWxH_8uC1(const uint8_t* pSrc,
  262.                                      unsigned int srcStep,
  263.                                      uint8_t* pDst,
  264.                                      unsigned int dstStep,
  265.                                      int width,
  266.                                      int height) {
  267.     // 检查对齐
  268.     // 1. pDst 地址必须 16 字节对齐
  269.     // 2. dstStep 必须是 16 的倍数
  270.     // 只有同时满足,才能在 64x64 块内部安全使用 stream 指令
  271.     bool isAligned = (((uintptr_t)pDst | (uintptr_t)dstStep) & 0xF) == 0;
  272.     if (isAligned) {
  273.         return icv_y8_owniTransposeWxH_8uC1_impl<true>(pSrc, srcStep, pDst, dstStep, width, height);
  274.     } else {
  275.         return icv_y8_owniTransposeWxH_8uC1_impl<false>(pSrc, srcStep, pDst, dstStep, width, height);
  276.     }
  277. }
  278. /**
  279. * 顶层转置函数:将整幅图像按 512x512 分块,调度到 icv_y8_owniTransposeWxH_8uC1
  280. */
  281. int64_t icv_transpose_8u_C1(const uint8_t* pSrc,
  282.                                          unsigned int srcStep,
  283.                                          uint8_t* pDst,
  284.                                          unsigned int dstStep,
  285.                                          int width,
  286.                                          int height) {
  287.     constexpr int TILE = 512;
  288.     if (width <= 0 || height <= 0) {
  289.         return 0;
  290.     }
  291.     const int h_main = height & ~(TILE - 1); // height - height % 512
  292.     const int w_main = width & ~(TILE - 1);  // width  - width  % 512
  293.     const int h_tail = height - h_main;      // height % 512
  294.     const int w_tail = width - w_main;       // width % 512
  295.     int64_t last_ret = 0; // 保存最后一次调用内核的返回值
  296.     // 1. 主 512x512 网格区域:0..h_main-1, 0..w_main-1
  297.     //    外层循环 Width (bj),内层展开 Height (bi)
  298.     //    这使得 Source 每次读取跳跃 512 行 (垂直),
  299.     //    而 Destination 每次写入跳跃 512 列 (水平,即连续内存),
  300.     //    这对写合并缓冲 (Write Combining) 非常友好
  301.     if (h_main > 0 && w_main > 0) {
  302.         const int blocksH = h_main / TILE; // 垂直方向块数
  303.         const int blocksW = w_main / TILE; // 水平方向块数
  304.         const int GROUP = 8;               // 8 个 512x512 块一组
  305.         // 外层遍历 Destination 的行 (即 Source 的列)
  306.         for (int bj = 0; bj < blocksW; ++bj) {
  307.             const int srcColOffset = bj * TILE;
  308.             const int dstRowOffset = bj * TILE * static_cast<int>(dstStep);
  309.             int bi = 0;
  310.             // 1a. 内层展开:处理 Source 的 8 个垂直块 (Vertical Blocks)
  311.             //     这会生成 Destination 的 8 个水平块 (Horizontal Blocks -> 连续写入)
  312.             for (; bi + GROUP - 1 < blocksH; bi += GROUP) {
  313.                 const int srcRowOffset = bi * TILE * static_cast<int>(srcStep);
  314.                 const int dstColOffset = bi * TILE;
  315.                 const uint8_t* srcBase = pSrc + srcRowOffset + srcColOffset;
  316.                 uint8_t* dstBase = pDst + dstRowOffset + dstColOffset;
  317.                 // Source 指针每次加 srcStep * TILE (垂直移动)
  318.                 // Dest 指针每次加 TILE (水平移动)
  319.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 0 * TILE * srcStep, srcStep,
  320.                                                         dstBase + 0 * TILE, dstStep, TILE, TILE);
  321.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 1 * TILE * srcStep, srcStep,
  322.                                                         dstBase + 1 * TILE, dstStep, TILE, TILE);
  323.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 2 * TILE * srcStep, srcStep,
  324.                                                         dstBase + 2 * TILE, dstStep, TILE, TILE);
  325.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 3 * TILE * srcStep, srcStep,
  326.                                                         dstBase + 3 * TILE, dstStep, TILE, TILE);
  327.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 4 * TILE * srcStep, srcStep,
  328.                                                         dstBase + 4 * TILE, dstStep, TILE, TILE);
  329.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 5 * TILE * srcStep, srcStep,
  330.                                                         dstBase + 5 * TILE, dstStep, TILE, TILE);
  331.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 6 * TILE * srcStep, srcStep,
  332.                                                         dstBase + 6 * TILE, dstStep, TILE, TILE);
  333.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 7 * TILE * srcStep, srcStep,
  334.                                                         dstBase + 7 * TILE, dstStep, TILE, TILE);
  335.             }
  336.             // 1b. 本行(列)剩余的 512x512 块(不足 8 个的一段)
  337.             for (; bi < blocksH; ++bi) {
  338.                 const int srcRowOffset = bi * TILE * static_cast<int>(srcStep);
  339.                 const int dstColOffset = bi * TILE;
  340.                 const uint8_t* srcBlock = pSrc + srcRowOffset + srcColOffset;
  341.                 uint8_t* dstBlock = pDst + dstRowOffset + dstColOffset;
  342.                 last_ret = icv_y8_owniTransposeWxH_8uC1(srcBlock, srcStep, dstBlock, dstStep, TILE, TILE);
  343.             }
  344.         }
  345.     }
  346.     // 2. 右侧边缘:宽度剩余 w_tail x 高度 h_main
  347.     //    这个区域的块尺寸是 w_tail x 512
  348.     if (w_tail > 0 && h_main > 0) {
  349.         const int blocksH = h_main / TILE;
  350.         for (int bi = 0; bi < blocksH; ++bi) {
  351.             const int srcRowOffset = bi * TILE * static_cast<int>(srcStep);
  352.             const int dstColOffset = bi * TILE;
  353.             const uint8_t* srcBlock = pSrc + srcRowOffset + w_main;
  354.             uint8_t* dstBlock = pDst + w_main * static_cast<int>(dstStep) + dstColOffset;
  355.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBlock, srcStep, dstBlock, dstStep, w_tail, TILE);
  356.         }
  357.     }
  358.     // 3. 底部边缘:宽度 w_main x 高度剩余 h_tail
  359.     //    区域被拆成若干 512x h_tail 的块,同样按宽度做 8x 展开
  360.     //    (这部分保持不变,因为高度 < 512,无法进行垂直展开)
  361.     if (h_tail > 0 && w_main > 0) {
  362.         const int blocksW = w_main / TILE;
  363.         const int GROUP = 8;
  364.         const int srcRowOffsetBase = h_main * static_cast<int>(srcStep);
  365.         const int dstColOffsetBase = h_main;
  366.         int bj = 0;
  367.         // 3a. 每次处理 8 个 512x h_tail 的块
  368.         for (; bj + GROUP - 1 < blocksW; bj += GROUP) {
  369.             const int srcColOffset = bj * TILE;
  370.             const int dstRowOffset = bj * TILE * static_cast<int>(dstStep);
  371.             const uint8_t* srcBase = pSrc + srcRowOffsetBase + srcColOffset;
  372.             uint8_t* dstBase = pDst + dstRowOffset + dstColOffsetBase;
  373.             // 水平展开
  374.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 0 * TILE,
  375.                                                     srcStep,
  376.                                                     dstBase + 0 * TILE * static_cast<int>(dstStep),
  377.                                                     dstStep,
  378.                                                     TILE,
  379.                                                     h_tail);
  380.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 1 * TILE,
  381.                                                     srcStep,
  382.                                                     dstBase + 1 * TILE * static_cast<int>(dstStep),
  383.                                                     dstStep,
  384.                                                     TILE,
  385.                                                     h_tail);
  386.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 2 * TILE,
  387.                                                     srcStep,
  388.                                                     dstBase + 2 * TILE * static_cast<int>(dstStep),
  389.                                                     dstStep,
  390.                                                     TILE,
  391.                                                     h_tail);
  392.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 3 * TILE,
  393.                                                     srcStep,
  394.                                                     dstBase + 3 * TILE * static_cast<int>(dstStep),
  395.                                                     dstStep,
  396.                                                     TILE,
  397.                                                     h_tail);
  398.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 4 * TILE,
  399.                                                     srcStep,
  400.                                                     dstBase + 4 * TILE * static_cast<int>(dstStep),
  401.                                                     dstStep,
  402.                                                     TILE,
  403.                                                     h_tail);
  404.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 5 * TILE,
  405.                                                     srcStep,
  406.                                                     dstBase + 5 * TILE * static_cast<int>(dstStep),
  407.                                                     dstStep,
  408.                                                     TILE,
  409.                                                     h_tail);
  410.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 6 * TILE,
  411.                                                     srcStep,
  412.                                                     dstBase + 6 * TILE * static_cast<int>(dstStep),
  413.                                                     dstStep,
  414.                                                     TILE,
  415.                                                     h_tail);
  416.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBase + 7 * TILE,
  417.                                                     srcStep,
  418.                                                     dstBase + 7 * TILE * static_cast<int>(dstStep),
  419.                                                     dstStep,
  420.                                                     TILE,
  421.                                                     h_tail);
  422.         }
  423.         // 3b. 本行剩余的 512x h_tail 块
  424.         for (; bj < blocksW; ++bj) {
  425.             const uint8_t* srcBlock = pSrc + srcRowOffsetBase + bj * TILE;
  426.             uint8_t* dstBlock = pDst + bj * TILE * static_cast<int>(dstStep) + dstColOffsetBase;
  427.             last_ret = icv_y8_owniTransposeWxH_8uC1(srcBlock, srcStep, dstBlock, dstStep, TILE, h_tail);
  428.         }
  429.     }
  430.     // 4. 右下角小块:w_tail x h_tail
  431.     if (h_tail > 0 && w_tail > 0) {
  432.         const uint8_t* srcBlock = pSrc + h_main * static_cast<int>(srcStep) + w_main;
  433.         uint8_t* dstBlock = pDst + w_main * static_cast<int>(dstStep) + h_main;
  434.         last_ret = icv_y8_owniTransposeWxH_8uC1(srcBlock, srcStep, dstBlock, dstStep, w_tail, h_tail);
  435.     }
  436.     return last_ret;
  437. }
复制代码
性能测试

普通c++基准代码如下,一般可以得到6倍加速
  1. template <typename T>
  2. void simple_transpose(
  3.         const T* pSrc, unsigned int srcStep,
  4.         T* pDst, unsigned int dstStep,
  5.         int width, int height)
  6. {
  7.     // 块大小根据 CPU L1 Cache 大小调整
  8.     // 对于 uint8_t,64x64 = 4KB,通常适合 L1 Cache
  9.     // 对于 float,可能需要减小到 32 或 16
  10.     constexpr int BLOCK_SIZE = 64;
  11.     // 以块为单位遍历 (i, j 指向块的左上角)
  12.     for (int i = 0; i < height; i += BLOCK_SIZE)
  13.     {
  14.         for (int j = 0; j < width; j += BLOCK_SIZE)
  15.         {
  16.             // 如果 i + BLOCK_SIZE 超过了 height,则只处理到 height 为止
  17.             // 解决了矩阵尺寸不能被 BLOCK_SIZE 整除的情况
  18.             const int i_max = std::min(i + BLOCK_SIZE, height);
  19.             const int j_max = std::min(j + BLOCK_SIZE, width);
  20.             for (int ii = i; ii < i_max; ++ii)
  21.             {
  22.                 for (int jj = j; jj < j_max; ++jj)
  23.                 {
  24.                     // Dst(row, col) = Src(col, row)
  25.                     pDst[jj * dstStep + ii] = pSrc[ii * srcStep + jj];
  26.                 }
  27.             }
  28.         }
  29.     }
  30. }
复制代码
为了节约篇幅就不展示功能测试了。
性能测试代码如下
  1. class TransposeFixture : public benchmark::Fixture
  2. {
  3. public:
  4.     void SetUp(const ::benchmark::State& state) override
  5.     {
  6.         width = state.range(0);
  7.         height = state.range(1);
  8.         srcStep = width;
  9.         src = std::make_unique<uint8_t[]>(static_cast<size_t>(height) * srcStep);
  10.         dst = std::make_unique<uint8_t[]>(static_cast<size_t>(width) * height);
  11.         for (int i = 0; i < height; ++i) {
  12.             for (int j = 0; j < width; ++j) {
  13.                 src[i * srcStep + j] = static_cast<uint8_t>((i + j) % 256);
  14.             }
  15.         }
  16.         srcMat = cv::Mat(height, width, CV_8UC1, src.get(), srcStep);
  17.         dstMat = cv::Mat(width, height, CV_8UC1, dst.get(), height);
  18.         bytes_per_iteration = static_cast<int64_t>(width) * height * 2;
  19.     }
  20.     void TearDown(const ::benchmark::State&) override
  21.     {
  22.         src.reset();
  23.         dst.reset();
  24.     }
  25. protected:
  26.     int width{};
  27.     int height{};
  28.     int srcStep{};
  29.     int64_t bytes_per_iteration{};
  30.     std::unique_ptr<uint8_t[]> src;
  31.     std::unique_ptr<uint8_t[]> dst;
  32.     cv::Mat srcMat;  // 预创建的cv::Mat对象
  33.     cv::Mat dstMat;  // 预创建的cv::Mat对象
  34. };
  35. BENCHMARK_DEFINE_F(TransposeFixture, Optimized)(benchmark::State& state)
  36. {
  37.     for (auto _ : state) {
  38.         icv_transpose_8u_C1(
  39.                 src.get(), srcStep,
  40.                 dst.get(), height,
  41.                 width, height
  42.         );
  43.         benchmark::DoNotOptimize(dst.get());
  44.         benchmark::ClobberMemory();
  45.     }
  46.     state.SetBytesProcessed(state.iterations() * bytes_per_iteration);
  47. }
  48. BENCHMARK_DEFINE_F(TransposeFixture, OpenCV)(benchmark::State& state)
  49. {
  50.     for (auto _ : state) {
  51.         cv::transpose(srcMat, dstMat);
  52.         benchmark::DoNotOptimize(dstMat.data);
  53.         benchmark::ClobberMemory();
  54.     }
  55.     state.SetBytesProcessed(state.iterations() * bytes_per_iteration);
  56. }
复制代码
性能测试结果如下
  1. TransposeFixture/Optimized/4096/4096    1538876 ns      1537400 ns          498 bytes_per_second=20.3265Gi/s
  2. TransposeFixture/OpenCV/4096/4096       6065054 ns      5998884 ns          112 bytes_per_second=5.2093Gi/s
  3. TransposeFixture/Optimized/2050/1920     285198 ns       288771 ns         2489 bytes_per_second=25.3882Gi/s
  4. TransposeFixture/OpenCV/2050/1920        279878 ns       284630 ns         2635 bytes_per_second=25.7576Gi/s
复制代码
可以看到在大图上能够完全超越,在普通图像上性能接近。

来源:程序园用户自行投稿发布,如果侵权,请联系站长删除
免责声明:如果侵犯了您的权益,请联系站长,我们会及时删除侵权内容,谢谢合作!

相关推荐

2025-11-27 06:33:53

举报

您需要登录后才可以回帖 登录 | 立即注册