From 09017f0bc0e9ccc29fb09593497bbc79b4e5ff55 Mon Sep 17 00:00:00 2001 From: "y.qiu" Date: Fri, 20 Sep 2024 22:50:25 +0800 Subject: [PATCH] let integral as standalone --- CMakeLists.txt | 2 + grayMatch.cpp | 175 +----------------------------------------------- integral.cpp | 176 +++++++++++++++++++++++++++++++++++++++++++++++++ integral.h | 8 +++ 4 files changed, 188 insertions(+), 173 deletions(-) create mode 100644 integral.cpp create mode 100644 integral.h diff --git a/CMakeLists.txt b/CMakeLists.txt index ff0780d..3257693 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -24,6 +24,8 @@ add_library(algo SHARED serialize.cpp privateType.h apiExport.h + integral.h + integral.cpp ) target_include_directories(algo PRIVATE ${OpenCV_INCLUDE_DIRS}) diff --git a/grayMatch.cpp b/grayMatch.cpp index 8a72d24..cd5468c 100644 --- a/grayMatch.cpp +++ b/grayMatch.cpp @@ -211,6 +211,8 @@ cv::Size computeRotationSize(const cv::Size &dstSize, const cv::Size &templateSi } #ifdef CV_SIMD +#include "integral.h" + void matchTemplateSimd(const cv::Mat &src, const cv::Mat &templateImg, cv::Mat &result) { result = cv::Mat(src.size() - templateImg.size() + cv::Size(1, 1), CV_32FC1); @@ -254,179 +256,6 @@ void matchTemplateSimd(const cv::Mat &src, const cv::Mat &templateImg, cv::Mat & } } } - -inline void expand(const cv::v_int32 &src, cv::v_float64 &low, cv::v_float64 &high) { - low = cv::v_cvt_f64(src); - high = cv::v_cvt_f64_high(src); -} - -inline void integralSum(const cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { - auto sum = src + cv::v_rotate_left<1>(src); - sum += cv::v_rotate_left<2>(sum); - sum += cv::v_rotate_left<4>(sum); - - cv::v_uint32 v1; - cv::v_uint32 v2; - cv::v_expand(sum, v1, v2); - v1 += pre; - v2 += pre; - pre = cv::v_setall_u32(cv::v_extract_n(v2)); - - cv::v_float64 v3; - cv::v_float64 v4; - expand(cv::v_reinterpret_as_s32(v1), v3, v4); - cv::v_store(dst, v3 + cv::v_load(prevDst)); - cv::v_store(dst + cv::v_float64::nlanes, v4 + cv::v_load(prevDst + cv::v_float64::nlanes)); - - expand(cv::v_reinterpret_as_s32(v2), v3, v4); - cv::v_store(dst + cv::v_float64::nlanes * 2, - v3 + cv::v_load(prevDst + cv::v_float64::nlanes * 2)); - cv::v_store(dst + cv::v_float64::nlanes * 3, - v4 + cv::v_load(prevDst + cv::v_float64::nlanes * 3)); -} - -inline void integralSqSum(cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { - cv::v_uint32 v1; - cv::v_uint32 v2; - cv::v_expand(src, v1, v2); - - { - auto shift1 = cv::v_rotate_left<1>(src); - cv::v_uint32 v3; - cv::v_uint32 v4; - cv::v_expand(shift1, v3, v4); - - v1 += v3; - v2 += v4; - - v4 = cv::v_extract<2>(v1, v2); - v2 += v4; - - v3 = cv::v_rotate_left<2>(v1); - v1 += v3; - - v1 += pre; - v2 += v1; - - pre = cv::v_setall_u32(cv::v_extract_n(v2)); - } - - cv::v_float64 v3; - cv::v_float64 v4; - expand(cv::v_reinterpret_as_s32(v1), v3, v4); - cv::v_store(dst, v3 + cv::v_load(prevDst)); - cv::v_store(dst + cv::v_float64::nlanes, v4 + cv::v_load(prevDst + cv::v_float64::nlanes)); - - expand(cv::v_reinterpret_as_s32(v2), v3, v4); - cv::v_store(dst + cv::v_float64::nlanes * 2, - v3 + cv::v_load(prevDst + cv::v_float64::nlanes * 2)); - cv::v_store(dst + cv::v_float64::nlanes * 3, - v4 + cv::v_load(prevDst + cv::v_float64::nlanes * 3)); -} - -/* -inline void integralSqSum(cv::v_uint32 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { - src += cv::v_rotate_left<1>(src); - src += cv::v_rotate_left<2>(src); - src += pre; - pre = cv::v_setall_u32(cv::v_extract_n(src)); - - cv::v_float64 v1; - cv::v_float64 v2; - expand(cv::v_reinterpret_as_s32(src), v1, v2); - - cv::v_store(dst, v1 + cv::v_load(prevDst)); - cv::v_store(dst + cv::v_float64::nlanes, v2 + cv::v_load(prevDst + cv::v_float64::nlanes)); -} - -inline void integralSqSum(cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { - cv::v_uint32 v1; - cv::v_uint32 v2; - cv::v_expand(src, v1, v2); - integralSqSum(v1, dst, prevDst, pre); - integralSqSum(v2, dst + cv::v_uint32::nlanes, prevDst + cv::v_uint32::nlanes, pre); -} -*/ - -inline void integralSum(cv::v_uint16 &v1, cv::v_uint16 &v2, double *dst, double *prevDst, - cv::v_uint32 &pre) { - integralSum(v1, dst, prevDst, pre); - integralSum(v2, dst + cv::v_uint16::nlanes, prevDst + cv::v_uint16::nlanes, pre); -} - -inline void integralSqSum(cv::v_uint16 &v1, cv::v_uint16 &v2, double *dst, double *prevDst, - cv::v_uint32 &pre) { - v1 = cv::v_mul_wrap(v1, v1); - v2 = cv::v_mul_wrap(v2, v2); - - integralSqSum(v1, dst, prevDst, pre); - integralSqSum(v2, dst + cv::v_uint16::nlanes, prevDst + cv::v_uint16::nlanes, pre); -} - -void integralSimd(const cv::Mat &src, cv::Mat &sum, cv::Mat &sqSum) { - const auto size = src.size() + cv::Size(1, 1); - sum.create(size, CV_64FC1); - sqSum.create(size, CV_64FC1); - memset(sum.data, 0, sum.step[ 0 ]); - memset(sqSum.data, 0, sqSum.step[ 0 ]); - - const auto *srcStart = src.ptr(); - const auto srcStep = src.step[ 0 ]; - auto *sumStart = sum.ptr(1) + 1; - const auto sumStep = sum.step[ 0 ] / sum.step[ 1 ]; - auto *sqSumStart = sqSum.ptr(1) + 1; - const auto sqSumStep = sqSum.step[ 0 ] / sqSum.step[ 1 ]; - const auto end = src.cols - cv::v_uint8::nlanes; - for (int y = 0; y < src.rows; y++) { - auto *srcPtr = srcStart + srcStep * y; - auto *sumPtr = sumStart + sumStep * y; - auto *preSumPtr = sumStart + sumStep * (y - 1); - sumPtr[ -1 ] = 0; - - cv::v_uint32 prevSum = cv::vx_setzero_u32(); - int x = 0; - for (; x < end; x += cv::v_uint8::nlanes) { - cv::v_uint16 v1; - cv::v_uint16 v2; - cv::v_expand(cv::v_load(srcPtr + x), v1, v2); - - integralSum(v1, v2, sumPtr + x, preSumPtr + x, prevSum); - } - } - - for (int y = 0; y < src.rows; y++) { - auto *srcPtr = srcStart + srcStep * y; - auto *sqSumPtr = sqSumStart + sqSumStep * y; - auto *preSqSumPtr = sqSumStart + sqSumStep * (y - 1); - sqSumPtr[ -1 ] = 0; - - cv::v_uint32 prevSqSum = cv::vx_setzero_u32(); - int x = 0; - for (; x < end; x += cv::v_uint8::nlanes) { - cv::v_uint16 v1; - cv::v_uint16 v2; - cv::v_expand(cv::v_load(srcPtr + x), v1, v2); - - integralSqSum(v1, v2, sqSumPtr + x, preSqSumPtr + x, prevSqSum); - } - } - - const auto start = src.cols - src.cols % cv::v_uint8::nlanes; - for (int y = 0; y < src.rows; y++) { - auto *srcPtr = srcStart + srcStep * y; - auto *sumPtr = sumStart + sumStep * y; - auto *sqSumPtr = sqSumStart + sqSumStep * y; - auto *preSumPtr = sumStart + sumStep * (y - 1); - auto *preSqSumPtr = sqSumStart + sqSumStep * (y - 1); - for (int x = start; x < src.cols; x++) { - auto val = srcPtr[ x ]; - auto sqVal = val * val; - - sumPtr[ x ] = sumPtr[ x - 1 ] + val + preSumPtr[ x ] - preSumPtr[ x - 1 ]; - sqSumPtr[ x ] = sqSumPtr[ x - 1 ] + sqVal + preSqSumPtr[ x ] - preSqSumPtr[ x - 1 ]; - } - } -} #endif void ccoeffDenominator(const cv::Mat &src, const cv::Size &templateSize, cv::Mat &result, diff --git a/integral.cpp b/integral.cpp new file mode 100644 index 0000000..6bea837 --- /dev/null +++ b/integral.cpp @@ -0,0 +1,176 @@ +#include "integral.h" + +#include + +inline void expand(const cv::v_int32 &src, cv::v_float64 &low, cv::v_float64 &high) { + low = cv::v_cvt_f64(src); + high = cv::v_cvt_f64_high(src); +} + +inline void integralSum(const cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { + auto sum = src + cv::v_rotate_left<1>(src); + sum += cv::v_rotate_left<2>(sum); + sum += cv::v_rotate_left<4>(sum); + + cv::v_uint32 v1; + cv::v_uint32 v2; + cv::v_expand(sum, v1, v2); + v1 += pre; + v2 += pre; + pre = cv::v_setall_u32(cv::v_extract_n(v2)); + + cv::v_float64 v3; + cv::v_float64 v4; + expand(cv::v_reinterpret_as_s32(v1), v3, v4); + cv::v_store(dst, v3 + cv::v_load(prevDst)); + cv::v_store(dst + cv::v_float64::nlanes, v4 + cv::v_load(prevDst + cv::v_float64::nlanes)); + + expand(cv::v_reinterpret_as_s32(v2), v3, v4); + cv::v_store(dst + cv::v_float64::nlanes * 2, + v3 + cv::v_load(prevDst + cv::v_float64::nlanes * 2)); + cv::v_store(dst + cv::v_float64::nlanes * 3, + v4 + cv::v_load(prevDst + cv::v_float64::nlanes * 3)); +} + +inline void integralSqSum(cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { + cv::v_uint32 v1; + cv::v_uint32 v2; + cv::v_expand(src, v1, v2); + + { + auto shift1 = cv::v_rotate_left<1>(src); + cv::v_uint32 v3; + cv::v_uint32 v4; + cv::v_expand(shift1, v3, v4); + + v1 += v3; + v2 += v4; + + v4 = cv::v_extract<2>(v1, v2); + v2 += v4; + + v3 = cv::v_rotate_left<2>(v1); + v1 += v3; + + v1 += pre; + v2 += v1; + + pre = cv::v_setall_u32(cv::v_extract_n(v2)); + } + + cv::v_float64 v3; + cv::v_float64 v4; + expand(cv::v_reinterpret_as_s32(v1), v3, v4); + cv::v_store(dst, v3 + cv::v_load(prevDst)); + cv::v_store(dst + cv::v_float64::nlanes, v4 + cv::v_load(prevDst + cv::v_float64::nlanes)); + + expand(cv::v_reinterpret_as_s32(v2), v3, v4); + cv::v_store(dst + cv::v_float64::nlanes * 2, + v3 + cv::v_load(prevDst + cv::v_float64::nlanes * 2)); + cv::v_store(dst + cv::v_float64::nlanes * 3, + v4 + cv::v_load(prevDst + cv::v_float64::nlanes * 3)); +} + +/* +inline void integralSqSum(cv::v_uint32 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { + src += cv::v_rotate_left<1>(src); + src += cv::v_rotate_left<2>(src); + src += pre; + pre = cv::v_setall_u32(cv::v_extract_n(src)); + + cv::v_float64 v1; + cv::v_float64 v2; + expand(cv::v_reinterpret_as_s32(src), v1, v2); + + cv::v_store(dst, v1 + cv::v_load(prevDst)); + cv::v_store(dst + cv::v_float64::nlanes, v2 + cv::v_load(prevDst + cv::v_float64::nlanes)); +} + +inline void integralSqSum(cv::v_uint16 &src, double *dst, double *prevDst, cv::v_uint32 &pre) { + cv::v_uint32 v1; + cv::v_uint32 v2; + cv::v_expand(src, v1, v2); + integralSqSum(v1, dst, prevDst, pre); + integralSqSum(v2, dst + cv::v_uint32::nlanes, prevDst + cv::v_uint32::nlanes, pre); +} +*/ + +inline void integralSum(cv::v_uint16 &v1, cv::v_uint16 &v2, double *dst, double *prevDst, + cv::v_uint32 &pre) { + integralSum(v1, dst, prevDst, pre); + integralSum(v2, dst + cv::v_uint16::nlanes, prevDst + cv::v_uint16::nlanes, pre); +} + +inline void integralSqSum(cv::v_uint16 &v1, cv::v_uint16 &v2, double *dst, double *prevDst, + cv::v_uint32 &pre) { + v1 = cv::v_mul_wrap(v1, v1); + v2 = cv::v_mul_wrap(v2, v2); + + integralSqSum(v1, dst, prevDst, pre); + integralSqSum(v2, dst + cv::v_uint16::nlanes, prevDst + cv::v_uint16::nlanes, pre); +} + +void integralSimd(const cv::Mat &src, cv::Mat &sum, cv::Mat &sqSum) { + const auto size = src.size() + cv::Size(1, 1); + sum.create(size, CV_64FC1); + sqSum.create(size, CV_64FC1); + memset(sum.data, 0, sum.step[ 0 ]); + memset(sqSum.data, 0, sqSum.step[ 0 ]); + + const auto *srcStart = src.ptr(); + const auto srcStep = src.step[ 0 ]; + auto *sumStart = sum.ptr(1) + 1; + const auto sumStep = sum.step[ 0 ] / sum.step[ 1 ]; + auto *sqSumStart = sqSum.ptr(1) + 1; + const auto sqSumStep = sqSum.step[ 0 ] / sqSum.step[ 1 ]; + const auto end = src.cols - cv::v_uint8::nlanes; + for (int y = 0; y < src.rows; y++) { + auto *srcPtr = srcStart + srcStep * y; + auto *sumPtr = sumStart + sumStep * y; + auto *preSumPtr = sumStart + sumStep * (y - 1); + sumPtr[ -1 ] = 0; + + cv::v_uint32 prevSum = cv::vx_setzero_u32(); + int x = 0; + for (; x < end; x += cv::v_uint8::nlanes) { + cv::v_uint16 v1; + cv::v_uint16 v2; + cv::v_expand(cv::v_load(srcPtr + x), v1, v2); + + integralSum(v1, v2, sumPtr + x, preSumPtr + x, prevSum); + } + } + + for (int y = 0; y < src.rows; y++) { + auto *srcPtr = srcStart + srcStep * y; + auto *sqSumPtr = sqSumStart + sqSumStep * y; + auto *preSqSumPtr = sqSumStart + sqSumStep * (y - 1); + sqSumPtr[ -1 ] = 0; + + cv::v_uint32 prevSqSum = cv::vx_setzero_u32(); + int x = 0; + for (; x < end; x += cv::v_uint8::nlanes) { + cv::v_uint16 v1; + cv::v_uint16 v2; + cv::v_expand(cv::v_load(srcPtr + x), v1, v2); + + integralSqSum(v1, v2, sqSumPtr + x, preSqSumPtr + x, prevSqSum); + } + } + + const auto start = src.cols - src.cols % cv::v_uint8::nlanes; + for (int y = 0; y < src.rows; y++) { + auto *srcPtr = srcStart + srcStep * y; + auto *sumPtr = sumStart + sumStep * y; + auto *sqSumPtr = sqSumStart + sqSumStep * y; + auto *preSumPtr = sumStart + sumStep * (y - 1); + auto *preSqSumPtr = sqSumStart + sqSumStep * (y - 1); + for (int x = start; x < src.cols; x++) { + auto val = srcPtr[ x ]; + auto sqVal = val * val; + + sumPtr[ x ] = sumPtr[ x - 1 ] + val + preSumPtr[ x ] - preSumPtr[ x - 1 ]; + sqSumPtr[ x ] = sqSumPtr[ x - 1 ] + sqVal + preSqSumPtr[ x ] - preSqSumPtr[ x - 1 ]; + } + } +} diff --git a/integral.h b/integral.h new file mode 100644 index 0000000..ad5fd5b --- /dev/null +++ b/integral.h @@ -0,0 +1,8 @@ +#ifndef INTEGRAL_H +#define INTEGRAL_H + +#include + +void integralSimd(const cv::Mat &src, cv::Mat &sum, cv::Mat &sqSum); + +#endif // INTEGRAL_H