potok/opencv
1
0
Fork 0
mirror of https://github.com/opencv/opencv.git synced 2026-07-21 19:33:03 +04:00
Code Packages Releases Wiki Activity

Merge pull request #29339 from amd:fast_mean_simd

Browse Source 
core: Vectorize meanStdDev #29339

- Add SIMD SumSqr_SIMD reductions (previously scalar) plus AVX-512 dispatch for the mean TU. 
- Uses universal intrinsics, so all SIMD backends benefit.

### Pull Request Readiness Checklist

See details at https://github.com/opencv/opencv/wiki/How_to_contribute#making-a-good-pull-request

- [x] I agree to contribute to the project under Apache 2 License.
- [x] To the best of my knowledge, the proposed patch is not based on a code under GPL or another license that is incompatible with OpenCV
- [x] The PR is proposed to the proper branch
- [ ] There is a reference to the original bug report and related work
- [x] There is accuracy test, performance test and test data in opencv_extra repository, if applicable
      Patch to opencv_extra has the same branch name.
- [x] The feature is well documented and sample code can be built with the project CMake
This commit is contained in:
Madan mohan Manokar
2026-06-23 15:43:15 +05:30
committed by GitHub
parent df0c9d7950
commit 4eb3f8ecb6
2 changed files with 295 additions and 3 deletions
Download Patch File Download Diff File Expand all files Collapse all files
modules/core/CMakeLists.txt
+1 -1
View File
@@ -8,7 +8,7 @@ ocv_add_dispatched_file(convert_scale SSE2 AVX2 LASX)
ocv_add_dispatched_file(count_non_zero SSE2 AVX2 LASX)
ocv_add_dispatched_file(has_non_zero SSE2 AVX2 LASX )
ocv_add_dispatched_file(matmul SSE2 SSE4_1 AVX2 AVX512_SKX NEON_DOTPROD LASX)
ocv_add_dispatched_file(mean SSE2 AVX2 LASX)
ocv_add_dispatched_file(mean SSE2 AVX2 AVX512_SKX AVX512_ICL LASX)
ocv_add_dispatched_file(merge SSE2 AVX2 LASX)
ocv_add_dispatched_file(split SSE2 AVX2 LASX)
ocv_add_dispatched_file(sum SSE2 AVX2 LASX)
modules/core/src/mean.simd.hpp
+294 -2
View File
@@ -2,6 +2,7 @@
// It is subject to the license terms in the LICENSE file found in the top-level directory
// of this distribution and at http://opencv.org/license.html
// Copyright (C) 2026, Advanced Micro Devices, Inc., all rights reserved.
#include "precomp.hpp"
#include "stat.hpp"
@@ -31,7 +32,35 @@ struct SumSqr_SIMD<uchar, int, int>
{
int operator () (const uchar * src0, const uchar * mask, int * sum, int * sqsum, int len, int cn) const
{
if (mask || (cn != 1 && cn != 2 && cn != 4))
if (mask)
return 0;
// cn==3: deinterleave so each lane belongs to one channel, then reduce per
// channel. sum/sqsum stay exact in s32 within the dispatcher's 1<<15 block.
if (cn == 3)
{
const int vl8 = VTraits<v_uint8>::vlanes();
v_int32 vs0 = vx_setzero_s32(), vs1 = vx_setzero_s32(), vs2 = vx_setzero_s32();
v_int32 vq0 = vx_setzero_s32(), vq1 = vx_setzero_s32(), vq2 = vx_setzero_s32();
auto acc = [](const v_uint8& ch, v_int32& vs, v_int32& vq)
{
v_uint16 lo, hi; v_expand(ch, lo, hi);
v_uint32 s0, s1, s2, s3; v_expand(lo, s0, s1); v_expand(hi, s2, s3);
vs = v_add(vs, v_reinterpret_as_s32(v_add(v_add(s0, s1), v_add(s2, s3))));
v_int16 d0 = v_reinterpret_as_s16(lo), d1 = v_reinterpret_as_s16(hi);
vq = v_add(vq, v_add(v_dotprod(d0, d0), v_dotprod(d1, d1)));
};
int e = 0;
for (; e <= len - vl8; e += vl8)
{
v_uint8 a, b, c; v_load_deinterleave(src0 + e * 3, a, b, c);
acc(a, vs0, vq0); acc(b, vs1, vq1); acc(c, vs2, vq2);
}
sum[0] += v_reduce_sum(vs0); sum[1] += v_reduce_sum(vs1); sum[2] += v_reduce_sum(vs2);
sqsum[0] += v_reduce_sum(vq0); sqsum[1] += v_reduce_sum(vq1); sqsum[2] += v_reduce_sum(vq2);
v_cleanup();
return e;
}
if (cn != 1 && cn != 2 && cn != 4)
return 0;
len *= cn;
@@ -98,7 +127,32 @@ struct SumSqr_SIMD<schar, int, int>
{
int operator () (const schar * src0, const uchar * mask, int * sum, int * sqsum, int len, int cn) const
{
if (mask || (cn != 1 && cn != 2 && cn != 4))
if (mask)
return 0;
if (cn == 3)
{
const int vl8 = VTraits<v_int8>::vlanes();
v_int32 vs0 = vx_setzero_s32(), vs1 = vx_setzero_s32(), vs2 = vx_setzero_s32();
v_int32 vq0 = vx_setzero_s32(), vq1 = vx_setzero_s32(), vq2 = vx_setzero_s32();
auto acc = [](const v_int8& ch, v_int32& vs, v_int32& vq)
{
v_int16 lo, hi; v_expand(ch, lo, hi);
v_int32 s0, s1, s2, s3; v_expand(lo, s0, s1); v_expand(hi, s2, s3);
vs = v_add(vs, v_add(v_add(s0, s1), v_add(s2, s3)));
vq = v_add(vq, v_add(v_dotprod(lo, lo), v_dotprod(hi, hi)));
};
int e = 0;
for (; e <= len - vl8; e += vl8)
{
v_int8 a, b, c; v_load_deinterleave(src0 + e * 3, a, b, c);
acc(a, vs0, vq0); acc(b, vs1, vq1); acc(c, vs2, vq2);
}
sum[0] += v_reduce_sum(vs0); sum[1] += v_reduce_sum(vs1); sum[2] += v_reduce_sum(vs2);
sqsum[0] += v_reduce_sum(vq0); sqsum[1] += v_reduce_sum(vq1); sqsum[2] += v_reduce_sum(vq2);
v_cleanup();
return e;
}
if (cn != 1 && cn != 2 && cn != 4)
return 0;
len *= cn;
@@ -160,6 +214,212 @@ struct SumSqr_SIMD<schar, int, int>
}
};
#if (CV_SIMD_64F || CV_SIMD_SCALABLE_64F)
// Expand a 16-bit vector into two s32 halves. The unsigned overload reinterprets
// after a widening expand (values 0..65535 stay non-negative in s32, so the
// subsequent v_cvt_f64 is exact); the signed overload expands directly. This is
// the only real difference between the ushort and short reductions below.
static inline void sumSqrExpandS32(const v_uint16& v, v_int32& lo, v_int32& hi)
{
v_uint32 a, b; v_expand(v, a, b);
lo = v_reinterpret_as_s32(a); hi = v_reinterpret_as_s32(b);
}
static inline void sumSqrExpandS32(const v_int16& v, v_int32& lo, v_int32& hi)
{
v_expand(v, lo, hi);
}
// Shared 16-bit -> {s32 sum, f64 sqsum} reduction for ushort/short. The register
// type is deduced from the pointer, so the body is identical for both depths.
template <typename T>
static inline int sumSqr16ToF64(const T* src0, int* sum, double* sqsum, int len, int cn)
{
typedef decltype(vx_load(src0)) VT; // v_uint16 or v_int16
if (cn == 3)
{
const int vl = VTraits<VT>::vlanes();
v_int32 vs0 = vx_setzero_s32(), vs1 = vx_setzero_s32(), vs2 = vx_setzero_s32();
v_float64 vq0 = vx_setzero_f64(), vq1 = vx_setzero_f64(), vq2 = vx_setzero_f64();
auto acc = [](const VT& ch, v_int32& vs, v_float64& vq)
{
v_int32 lo, hi; sumSqrExpandS32(ch, lo, hi);
vs = v_add(vs, v_add(lo, hi));
v_float64 f0 = v_cvt_f64(lo), f1 = v_cvt_f64_high(lo);
v_float64 f2 = v_cvt_f64(hi), f3 = v_cvt_f64_high(hi);
vq = v_fma(f0, f0, v_fma(f1, f1, v_fma(f2, f2, v_fma(f3, f3, vq))));
};
int e = 0;
for (; e <= len - vl; e += vl)
{
VT a, b, c; v_load_deinterleave(src0 + e * 3, a, b, c);
acc(a, vs0, vq0); acc(b, vs1, vq1); acc(c, vs2, vq2);
}
sum[0] += v_reduce_sum(vs0); sum[1] += v_reduce_sum(vs1); sum[2] += v_reduce_sum(vs2);
sqsum[0] += v_reduce_sum(vq0); sqsum[1] += v_reduce_sum(vq1); sqsum[2] += v_reduce_sum(vq2);
v_cleanup();
return e;
}
if (cn != 1 && cn != 2 && cn != 4)
return 0;
len *= cn;
const int vl16 = VTraits<VT>::vlanes();
const int vl32 = VTraits<v_int32>::vlanes();
const int vl64 = VTraits<v_float64>::vlanes();
// The lane->channel scatter below (i % cn) and the returned pixel count
// (x / cn) require the lane counts to be a multiple of cn. This always holds
// for fixed-width SIMD (vlanes is a power of two >= cn for cn in {1,2,4}),
// but scalable backends (e.g. RVV) may report a non-multiple; fall back to
// scalar there. vl16 == 2*vl32, so checking vl32 covers both scatter loops.
if (vl32 % cn != 0)
return 0;
v_int32 v_sum = vx_setzero_s32();
v_float64 v_sq0 = vx_setzero_f64(), v_sq1 = vx_setzero_f64();
v_float64 v_sq2 = vx_setzero_f64(), v_sq3 = vx_setzero_f64();
int x = 0;
for (; x <= len - vl16; x += vl16)
{
VT v_src = vx_load(src0 + x);
v_int32 lo, hi; sumSqrExpandS32(v_src, lo, hi);
v_sum = v_add(v_sum, v_add(lo, hi));
v_float64 f0 = v_cvt_f64(lo);
v_float64 f1 = v_cvt_f64_high(lo);
v_float64 f2 = v_cvt_f64(hi);
v_float64 f3 = v_cvt_f64_high(hi);
v_sq0 = v_fma(f0, f0, v_sq0);
v_sq1 = v_fma(f1, f1, v_sq1);
v_sq2 = v_fma(f2, f2, v_sq2);
v_sq3 = v_fma(f3, f3, v_sq3);
}
int CV_DECL_ALIGNED(CV_SIMD_WIDTH) sbuf[VTraits<v_int32>::max_nlanes];
double CV_DECL_ALIGNED(CV_SIMD_WIDTH) qbuf[VTraits<v_float64>::max_nlanes * 4];
v_store_aligned(sbuf, v_sum);
v_store_aligned(qbuf + vl64 * 0, v_sq0);
v_store_aligned(qbuf + vl64 * 1, v_sq1);
v_store_aligned(qbuf + vl64 * 2, v_sq2);
v_store_aligned(qbuf + vl64 * 3, v_sq3);
// sbuf lane j (j<vl32) holds samples j and j+vl32 -> channel j%cn (vl32%cn==0).
// qbuf index i (i<vl16) holds sample lane i -> channel i%cn (vl16%cn==0).
for (int i = 0; i < vl32; ++i)
sum[i % cn] += sbuf[i];
for (int i = 0; i < vl16; ++i)
sqsum[i % cn] += qbuf[i];
v_cleanup();
return x / cn;
}
// Shared 32-bit -> f64 reduction for int/float. Both sum and sqsum accumulate in
// f64 to match the scalar reference's double accumulation; values are widened to
// f64 before squaring so the square is computed in double (same as (double)v*v).
// The register type is deduced from the pointer, so int and float share one body.
template <typename T>
static inline int sumSqr32ToF64(const T* src0, double* sum, double* sqsum, int len, int cn)
{
typedef decltype(vx_load(src0)) VT; // v_int32 or v_float32
if (cn == 3)
{
const int vl = VTraits<VT>::vlanes();
v_float64 vs0 = vx_setzero_f64(), vs1 = vx_setzero_f64(), vs2 = vx_setzero_f64();
v_float64 vq0 = vx_setzero_f64(), vq1 = vx_setzero_f64(), vq2 = vx_setzero_f64();
auto acc = [](const VT& ch, v_float64& vs, v_float64& vq)
{
v_float64 f0 = v_cvt_f64(ch), f1 = v_cvt_f64_high(ch);
vs = v_add(vs, v_add(f0, f1));
vq = v_fma(f0, f0, v_fma(f1, f1, vq));
};
int e = 0;
for (; e <= len - vl; e += vl)
{
VT a, b, c; v_load_deinterleave(src0 + e * 3, a, b, c);
acc(a, vs0, vq0); acc(b, vs1, vq1); acc(c, vs2, vq2);
}
sum[0] += v_reduce_sum(vs0); sum[1] += v_reduce_sum(vs1); sum[2] += v_reduce_sum(vs2);
sqsum[0] += v_reduce_sum(vq0); sqsum[1] += v_reduce_sum(vq1); sqsum[2] += v_reduce_sum(vq2);
v_cleanup();
return e;
}
if (cn != 1 && cn != 2 && cn != 4)
return 0;
len *= cn;
const int vl32 = VTraits<VT>::vlanes();
const int vl64 = VTraits<v_float64>::vlanes();
// See note in sumSqr16ToF64: guard scalable backends whose lane count may
// not divide cn (the i % cn scatter / x / cn rely on it).
if (vl32 % cn != 0)
return 0;
v_float64 vs0 = vx_setzero_f64(), vs1 = vx_setzero_f64();
v_float64 vq0 = vx_setzero_f64(), vq1 = vx_setzero_f64();
int x = 0;
for (; x <= len - vl32; x += vl32)
{
VT v_src = vx_load(src0 + x);
v_float64 f0 = v_cvt_f64(v_src);
v_float64 f1 = v_cvt_f64_high(v_src);
vs0 = v_add(vs0, f0);
vs1 = v_add(vs1, f1);
vq0 = v_fma(f0, f0, vq0);
vq1 = v_fma(f1, f1, vq1);
}
double CV_DECL_ALIGNED(CV_SIMD_WIDTH) sbuf[VTraits<v_float64>::max_nlanes * 2];
double CV_DECL_ALIGNED(CV_SIMD_WIDTH) qbuf[VTraits<v_float64>::max_nlanes * 2];
v_store_aligned(sbuf, vs0); v_store_aligned(sbuf + vl64, vs1);
v_store_aligned(qbuf, vq0); v_store_aligned(qbuf + vl64, vq1);
for (int i = 0; i < vl32; ++i)
{
sum[i % cn] += sbuf[i];
sqsum[i % cn] += qbuf[i];
}
v_cleanup();
return x / cn;
}
template <>
struct SumSqr_SIMD<ushort, int, double>
{
int operator () (const ushort* src0, const uchar* mask, int* sum, double* sqsum, int len, int cn) const
{ return mask ? 0 : sumSqr16ToF64(src0, sum, sqsum, len, cn); }
};
template <>
struct SumSqr_SIMD<short, int, double>
{
int operator () (const short* src0, const uchar* mask, int* sum, double* sqsum, int len, int cn) const
{ return mask ? 0 : sumSqr16ToF64(src0, sum, sqsum, len, cn); }
};
template <>
struct SumSqr_SIMD<int, double, double>
{
int operator () (const int* src0, const uchar* mask, double* sum, double* sqsum, int len, int cn) const
{ return mask ? 0 : sumSqr32ToF64(src0, sum, sqsum, len, cn); }
};
template <>
struct SumSqr_SIMD<float, double, double>
{
int operator () (const float* src0, const uchar* mask, double* sum, double* sqsum, int len, int cn) const
{ return mask ? 0 : sumSqr32ToF64(src0, sum, sqsum, len, cn); }
};
#endif // CV_SIMD_64F
#endif
template<typename T, typename ST, typename SQT>
@@ -253,6 +513,21 @@ static int sumsqr_(const T* src0, const uchar* mask, ST* sum, SQT* sqsum, int le
sum[0] = s0;
sqsum[0] = sq0;
}
else if( cn == 2 )
{
ST s0 = sum[0], s1 = sum[1];
SQT sq0 = sqsum[0], sq1 = sqsum[1];
for( i = 0; i < len; i++, src += 2 )
if( mask[i] )
{
T v0 = src[0], v1 = src[1];
s0 += v0; sq0 += (SQT)v0*v0;
s1 += v1; sq1 += (SQT)v1*v1;
nzm++;
}
sum[0] = s0; sum[1] = s1;
sqsum[0] = sq0; sqsum[1] = sq1;
}
else if( cn == 3 )
{
ST s0 = sum[0], s1 = sum[1], s2 = sum[2];
@@ -269,6 +544,23 @@ static int sumsqr_(const T* src0, const uchar* mask, ST* sum, SQT* sqsum, int le
sum[0] = s0; sum[1] = s1; sum[2] = s2;
sqsum[0] = sq0; sqsum[1] = sq1; sqsum[2] = sq2;
}
else if( cn == 4 )
{
ST s0 = sum[0], s1 = sum[1], s2 = sum[2], s3 = sum[3];
SQT sq0 = sqsum[0], sq1 = sqsum[1], sq2 = sqsum[2], sq3 = sqsum[3];
for( i = 0; i < len; i++, src += 4 )
if( mask[i] )
{
T v0 = src[0], v1 = src[1], v2 = src[2], v3 = src[3];
s0 += v0; sq0 += (SQT)v0*v0;
s1 += v1; sq1 += (SQT)v1*v1;
s2 += v2; sq2 += (SQT)v2*v2;
s3 += v3; sq3 += (SQT)v3*v3;
nzm++;
}
sum[0] = s0; sum[1] = s1; sum[2] = s2; sum[3] = s3;
sqsum[0] = sq0; sqsum[1] = sq1; sqsum[2] = sq2; sqsum[3] = sq3;
}
else
{
for( i = 0; i < len; i++, src += cn )