diff --git a/modules/imgproc/src/bilateral_filter.dispatch.cpp b/modules/imgproc/src/bilateral_filter.dispatch.cpp index 7b992303b6..9ac499f685 100644 --- a/modules/imgproc/src/bilateral_filter.dispatch.cpp +++ b/modules/imgproc/src/bilateral_filter.dispatch.cpp @@ -13,6 +13,7 @@ // Copyright (C) 2000-2008, 2018, Intel Corporation, all rights reserved. // Copyright (C) 2009, Willow Garage Inc., all rights reserved. // Copyright (C) 2014-2015, Itseez Inc., all rights reserved. +// Copyright (C) 2025, Advanced Micro Devices, all rights reserved. // Third party copyrights are property of their respective owners. // // Redistribution and use in source and binary forms, with or without modification, @@ -174,8 +175,8 @@ bilateralFilter_8u( const Mat& src, Mat& dst, int d, return; } - double gauss_color_coeff = -0.5/(sigma_color*sigma_color); - double gauss_space_coeff = -0.5/(sigma_space*sigma_space); + float gauss_color_coeff = (float)(-0.5/(sigma_color*sigma_color)); + float gauss_space_coeff = (float)(-0.5/(sigma_space*sigma_space)); if( d <= 0 ) radius = cvRound(sigma_space*1.5); @@ -195,15 +196,31 @@ bilateralFilter_8u( const Mat& src, Mat& dst, int d, int* space_ofs = &_space_ofs[0]; // initialize color-related bilateral filter coefficients + i = 0; +#if (CV_SIMD || CV_SIMD_SCALABLE) + int nlanes = VTraits::vlanes(); + v_float32 v_gauss_color_coeff = vx_setall_f32(gauss_color_coeff); + float counter[16] = {0.0, 1., 2., 3., 4., 5., 6., 7., + 8., 9., 10., 11., 12., 13., 14., 15.}; + v_float32 v_i = vx_load(counter); + v_float32 v_inc = vx_setall_f32(float(nlanes)); - for( i = 0; i < 256*cn; i++ ) + for( ; i < (256*cn) - nlanes; i += nlanes ) + { + v_float32 v_color_weight = v_mul(v_mul(v_i,v_i), v_gauss_color_coeff); + v_store(color_weight + i, v_exp(v_color_weight)); + v_i = v_add(v_i, v_inc); + } +#endif + for(; i < 256*cn; i++ ) + { color_weight[i] = (float)std::exp(i*i*gauss_color_coeff); + } // initialize space-related bilateral filter coefficients for( i = -radius, maxk = 0; i <= radius; i++ ) { j = -radius; - for( ; j <= radius; j++ ) { double r = std::sqrt((double)i*i + (double)j*j); diff --git a/modules/imgproc/src/bilateral_filter.simd.hpp b/modules/imgproc/src/bilateral_filter.simd.hpp index 77e0328678..ab718eac91 100644 --- a/modules/imgproc/src/bilateral_filter.simd.hpp +++ b/modules/imgproc/src/bilateral_filter.simd.hpp @@ -13,6 +13,7 @@ // Copyright (C) 2000-2008, 2018, Intel Corporation, all rights reserved. // Copyright (C) 2009, Willow Garage Inc., all rights reserved. // Copyright (C) 2014-2015, Itseez Inc., all rights reserved. +// Copyright (C) 2025, Advanced Micro Devices, all rights reserved. // Third party copyrights are property of their respective owners. // // Redistribution and use in source and binary forms, with or without modification, @@ -73,13 +74,71 @@ public: { } +#if (CV_SIMD || CV_SIMD_SCALABLE) + static void expand(const v_uint8& v_input, v_uint32& v_out0, v_uint32& v_out1, v_uint32& v_out2, v_uint32& v_out3) + { + v_uint16 d0, d1; + v_expand(v_input, d0, d1); + v_expand(d0, v_out0, v_out1); + v_expand(d1, v_out2, v_out3); + } + + static void computeBilateral(const v_uint8& v_val_u8,const v_uint8& v_abs_diff, v_float32& kweight, const float* color_weight, + v_float32& v_wsum0, v_float32& v_sum0, v_float32& v_wsum1, v_float32& v_sum1, v_float32& v_wsum2, v_float32& v_sum2, v_float32& v_wsum3, v_float32& v_sum3) + { + v_uint32 d0, d1, d2, d3; + v_float32 w0, w1, w2, w3; + v_uint32 val0, val1, val2, val3; + expand(v_abs_diff, d0, d1, d2, d3); + expand(v_val_u8, val0, val1, val2, val3); + w0 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(d0))); + w1 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(d1))); + w2 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(d2))); + w3 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(d3))); + v_wsum0 = v_add(v_wsum0, w0); + v_sum0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum0); + v_wsum1 = v_add(v_wsum1, w1); + v_sum1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum1); + v_wsum2 = v_add(v_wsum2, w2); + v_sum2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum2); + v_wsum3 = v_add(v_wsum3, w3); + v_sum3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum3); + } + + static void computeBilateral(const v_uint8& v_val_u8,const v_uint8& v_abs_diff, const float* color_weight, + v_float32& v_wsum0, v_float32& v_sum0, v_float32& v_wsum1, v_float32& v_sum1, v_float32& v_wsum2, v_float32& v_sum2, v_float32& v_wsum3, v_float32& v_sum3) + { + v_uint32 d0, d1, d2, d3; + v_float32 w0, w1, w2, w3; + v_uint32 val0, val1, val2, val3; + expand(v_abs_diff, d0, d1, d2, d3); + expand(v_val_u8, val0, val1, val2, val3); + w0 = v_lut(color_weight, v_reinterpret_as_s32(d0)); + w1 = v_lut(color_weight, v_reinterpret_as_s32(d1)); + w2 = v_lut(color_weight, v_reinterpret_as_s32(d2)); + w3 = v_lut(color_weight, v_reinterpret_as_s32(d3)); + v_wsum0 = v_add(v_wsum0, w0); + v_sum0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum0); + v_wsum1 = v_add(v_wsum1, w1); + v_sum1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum1); + v_wsum2 = v_add(v_wsum2, w2); + v_sum2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum2); + v_wsum3 = v_add(v_wsum3, w3); + v_sum3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum3); + } +#endif + virtual void operator() (const Range& range) const CV_OVERRIDE { CV_INSTRUMENT_REGION(); int i, j, cn = dest->channels(), k; Size size = dest->size(); - +#if (CV_SIMD || CV_SIMD_SCALABLE) + int nlanes = VTraits::vlanes(); + int nlanes_2 = 2*nlanes; + int nlanes_4 = 4*nlanes; +#endif for( i = range.start; i < range.end; i++ ) { const uchar* sptr = temp->ptr(i+radius) + radius*cn; @@ -87,365 +146,233 @@ public: if( cn == 1 ) { - AutoBuffer buf(alignSize(size.width, CV_SIMD_WIDTH) + size.width + CV_SIMD_WIDTH - 1); - memset(buf.data(), 0, buf.size() * sizeof(float)); - float *sum = alignPtr(buf.data(), CV_SIMD_WIDTH); - float *wsum = sum + alignSize(size.width, CV_SIMD_WIDTH); - k = 0; - for(; k <= maxk-4; k+=4) + k = 0; j=0; + +#if (CV_SIMD || CV_SIMD_SCALABLE) + for ( ;j <= size.width - nlanes_4; j += nlanes_4) { - const uchar* ksptr0 = sptr + space_ofs[k]; - const uchar* ksptr1 = sptr + space_ofs[k+1]; - const uchar* ksptr2 = sptr + space_ofs[k+2]; - const uchar* ksptr3 = sptr + space_ofs[k+3]; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight0 = vx_setall_f32(space_weight[k]); - v_float32 kweight1 = vx_setall_f32(space_weight[k+1]); - v_float32 kweight2 = vx_setall_f32(space_weight[k+2]); - v_float32 kweight3 = vx_setall_f32(space_weight[k+3]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes()) + const uchar* sptr_j = sptr + j; + v_float32 v_wsum0 = vx_setzero_f32(); + v_float32 v_wsum1 = vx_setzero_f32(); + v_float32 v_wsum2 = vx_setzero_f32(); + v_float32 v_wsum3 = vx_setzero_f32(); + v_float32 v_sum0 = vx_setzero_f32(); + v_float32 v_sum1 = vx_setzero_f32(); + v_float32 v_sum2 = vx_setzero_f32(); + v_float32 v_sum3 = vx_setzero_f32(); + v_uint8 v_sptr8 = vx_load(sptr_j); + + k=0; + if(maxk==5) { - v_uint32 rval = vx_load_expand_q(sptr + j); + const uchar* ksptrline1 = sptr_j + space_ofs[0]; + const uchar* ksptrline2 = sptr_j + space_ofs[1]; + const uchar* ksptrline3 = sptr_j + space_ofs[4]; - v_uint32 val = vx_load_expand_q(ksptr0 + j); - v_float32 w = v_mul(kweight0, v_lut(color_weight, v_reinterpret_as_s32(v_absdiff(val, rval)))); - v_float32 v_wsum = v_add(vx_load_aligned(wsum + j), w); - v_float32 v_sum = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val)), w, vx_load_aligned(sum + j)); + v_float32 kweight = vx_setall_f32(space_weight[0]);//same weight for all, expect centre one - val = vx_load_expand_q(ksptr1 + j); - w = v_mul(kweight1, v_lut(color_weight, v_reinterpret_as_s32(v_absdiff(val, rval)))); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val)), w, v_sum); + v_uint8 v_val_u8_line1 = vx_load(ksptrline1); + v_uint8 v_val_u8_line2_0 = vx_load(ksptrline2); + v_uint8 v_val_u8_line2_1 = vx_load(ksptrline2 + 1); + v_uint8 v_val_u8_line2_2 = vx_load(ksptrline2 + 2); + v_uint8 v_val_u8_line3 = vx_load(ksptrline3); - val = vx_load_expand_q(ksptr2 + j); - w = v_mul(kweight2, v_lut(color_weight, v_reinterpret_as_s32(v_absdiff(val, rval)))); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val)), w, v_sum); + //compute abs diff + v_uint8 v_abs_diff_line1 = v_absdiff(v_val_u8_line1, v_sptr8); + v_uint8 v_abs_diff_line2_0 = v_absdiff(v_val_u8_line2_0, v_sptr8); + v_uint8 v_abs_diff_line2_1 = v_absdiff(v_val_u8_line2_1, v_sptr8); + v_uint8 v_abs_diff_line2_2 = v_absdiff(v_val_u8_line2_2, v_sptr8); + v_uint8 v_abs_diff_line3 = v_absdiff(v_val_u8_line3, v_sptr8); - val = vx_load_expand_q(ksptr3 + j); - w = v_mul(kweight3, v_lut(color_weight, v_reinterpret_as_s32(v_absdiff(val, rval)))); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val)), w, v_sum); + computeBilateral( v_val_u8_line1, v_abs_diff_line1, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line2_0, v_abs_diff_line2_0, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line2_1, v_abs_diff_line2_1, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line2_2, v_abs_diff_line2_2, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line3, v_abs_diff_line3, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); - v_store_aligned(wsum + j, v_wsum); - v_store_aligned(sum + j, v_sum); + k = maxk; } -#endif -#if CV_SIMD128 - v_float32x4 kweight4 = v_load(space_weight + k); -#endif - for (; j < size.width; j++) + else if(maxk==13) { -#if CV_SIMD128 - v_uint32x4 rval = v_setall_u32(sptr[j]); - v_uint32x4 val(ksptr0[j], ksptr1[j], ksptr2[j], ksptr3[j]); - v_float32x4 w = v_mul(kweight4, v_lut(this->color_weight, v_reinterpret_as_s32(v_absdiff(val, rval)))); - wsum[j] += v_reduce_sum(w); - sum[j] += v_reduce_sum(v_mul(v_cvt_f32(v_reinterpret_as_s32(val)), w)); -#else - int rval = sptr[j]; + const uchar* ksptrline1 = sptr_j + space_ofs[0]; + const uchar* ksptrline5 = sptr_j + space_ofs[12];//last element + const uchar* ksptrline2 = sptr_j + space_ofs[1]; + const uchar* ksptrline3 = sptr_j + space_ofs[4]; + const uchar* ksptrline4 = sptr_j + space_ofs[9]; - int val = ksptr0[j]; - float w = space_weight[k] * color_weight[std::abs(val - rval)]; - wsum[j] += w; - sum[j] += val * w; + v_float32 kweight = vx_setall_f32(space_weight[0]); - val = ksptr1[j]; - w = space_weight[k+1] * color_weight[std::abs(val - rval)]; - wsum[j] += w; - sum[j] += val * w; + //compute line 1 and 5 + v_uint8 v_val_u8_line1 = vx_load(ksptrline1); + v_uint8 v_val_u8_line5 = vx_load(ksptrline5); + v_uint8 v_abs_diff_line1 = v_absdiff(v_val_u8_line1, v_sptr8); + v_uint8 v_abs_diff_line5 = v_absdiff(v_val_u8_line5, v_sptr8); + computeBilateral( v_val_u8_line1, v_abs_diff_line1, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line5, v_abs_diff_line5, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); - val = ksptr2[j]; - w = space_weight[k+2] * color_weight[std::abs(val - rval)]; - wsum[j] += w; - sum[j] += val * w; + //compute line 2 and 4 + v_uint8 v_val_u8_line2_0 = vx_load(ksptrline2); + v_uint8 v_val_u8_line2_1 = vx_load(ksptrline2 + 1); + v_uint8 v_val_u8_line2_2 = vx_load(ksptrline2 + 2); - val = ksptr3[j]; - w = space_weight[k+3] * color_weight[std::abs(val - rval)]; - wsum[j] += w; - sum[j] += val * w; -#endif + v_uint8 v_val_u8_line4_0 = vx_load(ksptrline4); + v_uint8 v_val_u8_line4_1 = vx_load(ksptrline4 + 1); + v_uint8 v_val_u8_line4_2 = vx_load(ksptrline4 + 2); + + v_uint8 v_abs_diff_line2_0 = v_absdiff(v_val_u8_line2_0, v_sptr8); + v_uint8 v_abs_diff_line2_1 = v_absdiff(v_val_u8_line2_1, v_sptr8); + v_uint8 v_abs_diff_line2_2 = v_absdiff(v_val_u8_line2_2, v_sptr8); + v_uint8 v_abs_diff_line4_0 = v_absdiff(v_val_u8_line4_0, v_sptr8); + v_uint8 v_abs_diff_line4_1 = v_absdiff(v_val_u8_line4_1, v_sptr8); + v_uint8 v_abs_diff_line4_2 = v_absdiff(v_val_u8_line4_2, v_sptr8); + v_float32 kweight_1 = vx_setall_f32(space_weight[1]); + v_float32 kweight_2 = vx_setall_f32(space_weight[2]); + computeBilateral( v_val_u8_line2_0, v_abs_diff_line2_0, kweight_1, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line2_1, v_abs_diff_line2_1, kweight_2, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line2_2, v_abs_diff_line2_2, kweight_1, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + + computeBilateral( v_val_u8_line4_0, v_abs_diff_line4_0, kweight_1, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line4_1, v_abs_diff_line4_1, kweight_2, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line4_2, v_abs_diff_line4_2, kweight_1, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + + //compute line 3 + v_uint8 v_val_u8_line3_0 = vx_load(ksptrline3); + v_uint8 v_val_u8_line3_1 = vx_load(ksptrline3 + 1); + v_uint8 v_val_u8_line3_2 = vx_load(ksptrline3 + 2); + v_uint8 v_val_u8_line3_3 = vx_load(ksptrline3 + 3); + v_uint8 v_val_u8_line3_4 = vx_load(ksptrline3 + 4); + + v_uint8 v_abs_diff_line3_0 = v_absdiff(v_val_u8_line3_0, v_sptr8); + v_uint8 v_abs_diff_line3_1 = v_absdiff(v_val_u8_line3_1, v_sptr8); + v_uint8 v_abs_diff_line3_2 = v_absdiff(v_val_u8_line3_2, v_sptr8); + v_uint8 v_abs_diff_line3_3 = v_absdiff(v_val_u8_line3_3, v_sptr8); + v_uint8 v_abs_diff_line3_4 = v_absdiff(v_val_u8_line3_4, v_sptr8); + + computeBilateral( v_val_u8_line3_0, v_abs_diff_line3_0, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line3_1, v_abs_diff_line3_1, kweight_2, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line3_2, v_abs_diff_line3_2, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line3_3, v_abs_diff_line3_3, kweight_2, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + computeBilateral( v_val_u8_line3_4, v_abs_diff_line3_4, kweight, color_weight, + v_wsum0, v_sum0, v_wsum1, v_sum1, v_wsum2, v_sum2, v_wsum3, v_sum3); + + k = maxk; } + + for(; k < maxk; k++) + { + const uchar* ksptr = sptr_j + space_ofs[k]; + v_float32 kweight = vx_setall_f32(space_weight[k]); + + v_uint8 v_val_u8 = vx_load(ksptr); + v_uint8 v_abs_diff = v_absdiff(v_val_u8, v_sptr8); + + v_uint16 d0, d1; + v_uint32 diff0, diff1, diff2, diff3; + v_expand(v_abs_diff, d0, d1); + v_expand(d0, diff0, diff1); + v_expand(d1, diff2, diff3); + + v_uint16 v0, v1; + v_uint32 val0, val1, val2, val3; + v_expand(v_val_u8, v0, v1); + v_expand(v0, val0, val1); + v_expand(v1, val2, val3); + + v_float32 w0 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(diff0))); + v_float32 w1 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(diff1))); + v_float32 w2 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(diff2))); + v_float32 w3 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(diff3))); + + v_wsum0 = v_add(v_wsum0, w0); + v_sum0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum0); + v_wsum1 = v_add(v_wsum1, w1); + v_sum1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum1); + v_wsum2 = v_add(v_wsum2, w2); + v_sum2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum2); + v_wsum3 = v_add(v_wsum3, w3); + v_sum3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum3); + } + v_pack_u_store(dptr + j, v_pack(v_round(v_div(v_sum0, v_wsum0)), + v_round(v_div(v_sum1, v_wsum1)))); + + v_pack_u_store(dptr + j + nlanes_2, v_pack(v_round(v_div(v_sum2, v_wsum2)), + v_round(v_div(v_sum3, v_wsum3)))); } - for(; k < maxk; k++) - { - const uchar* ksptr = sptr + space_ofs[k]; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight = vx_setall_f32(space_weight[k]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes()) - { - v_uint32 val = vx_load_expand_q(ksptr + j); - v_float32 w = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(v_absdiff(val, vx_load_expand_q(sptr + j))))); - v_store_aligned(wsum + j, v_add(vx_load_aligned(wsum + j), w)); - v_store_aligned(sum + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(val)), w, vx_load_aligned(sum + j))); - } -#endif - for (; j < size.width; j++) - { - int val = ksptr[j]; - float w = space_weight[k] * color_weight[std::abs(val - sptr[j])]; - wsum[j] += w; - sum[j] += val * w; - } - } - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - for (; j <= size.width - 2*VTraits::vlanes(); j += 2*VTraits::vlanes()) - v_pack_u_store(dptr + j, v_pack(v_round(v_div(vx_load_aligned(sum + j), vx_load_aligned(wsum + j))), - v_round(v_div(vx_load_aligned(sum + j + VTraits::vlanes()), vx_load_aligned(wsum + j + VTraits::vlanes()))))); #endif for (; j < size.width; j++) { + uchar val0 = sptr[j]; + float wsumT = 0; + float sumT = 0; + for(k=0; k < maxk; k++) + { + const uchar* ksptr = sptr + space_ofs[k]; + uchar val = ksptr[j]; + float w = space_weight[k] * color_weight[std::abs(val - val0)]; + wsumT += w; + sumT += val * w; + } + // overflow is not possible here => there is no need to use cv::saturate_cast - CV_DbgAssert(fabs(wsum[j]) > 0); - dptr[j] = (uchar)cvRound(sum[j]/wsum[j]); + CV_DbgAssert(fabs(wsumT) > 0); + dptr[j] = (uchar)cvRound(sumT/wsumT); } + } else { CV_Assert( cn == 3 ); - AutoBuffer buf(alignSize(size.width, CV_SIMD_WIDTH)*3 + size.width + CV_SIMD_WIDTH - 1); - memset(buf.data(), 0, buf.size() * sizeof(float)); - float *sum_b = alignPtr(buf.data(), CV_SIMD_WIDTH); - float *sum_g = sum_b + alignSize(size.width, CV_SIMD_WIDTH); - float *sum_r = sum_g + alignSize(size.width, CV_SIMD_WIDTH); - float *wsum = sum_r + alignSize(size.width, CV_SIMD_WIDTH); - k = 0; - for(; k <= maxk-4; k+=4) - { - const uchar* ksptr0 = sptr + space_ofs[k]; - const uchar* ksptr1 = sptr + space_ofs[k+1]; - const uchar* ksptr2 = sptr + space_ofs[k+2]; - const uchar* ksptr3 = sptr + space_ofs[k+3]; - const uchar* rsptr = sptr; - j = 0; + j = 0; + const uchar* sptr_j = sptr; #if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight0 = vx_setall_f32(space_weight[k]); - v_float32 kweight1 = vx_setall_f32(space_weight[k+1]); - v_float32 kweight2 = vx_setall_f32(space_weight[k+2]); - v_float32 kweight3 = vx_setall_f32(space_weight[k+3]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), rsptr += 3*VTraits::vlanes(), - ksptr0 += 3*VTraits::vlanes(), ksptr1 += 3*VTraits::vlanes(), ksptr2 += 3*VTraits::vlanes(), ksptr3 += 3*VTraits::vlanes()) - { - v_uint8 kb, kg, kr, rb, rg, rr; - v_load_deinterleave(rsptr, rb, rg, rr); - - v_load_deinterleave(ksptr0, kb, kg, kr); - v_uint16 val0, val1, val2, val3, val4; - v_expand(v_absdiff(kb, rb), val0, val1); - v_expand(v_absdiff(kg, rg), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - v_expand(v_absdiff(kr, rr), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - - v_uint32 vall, valh; - v_expand(val0, vall, valh); - v_float32 w0 = v_mul(kweight0, v_lut(color_weight, v_reinterpret_as_s32(vall))); - v_float32 w1 = v_mul(kweight0, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j, v_add(w0, vx_load_aligned(wsum + j))); - v_store_aligned(wsum + j + VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + VTraits::vlanes()))); - v_expand(kb, val0, val2); - v_expand(val0, vall, valh); - v_store_aligned(sum_b + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_b + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + VTraits::vlanes()))); - v_expand(kg, val0, val3); - v_expand(val0, vall, valh); - v_store_aligned(sum_g + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_g + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + VTraits::vlanes()))); - v_expand(kr, val0, val4); - v_expand(val0, vall, valh); - v_store_aligned(sum_r + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j))); - v_store_aligned(sum_r + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + VTraits::vlanes()))); - - v_expand(val1, vall, valh); - w0 = v_mul(kweight0, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight0, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j + 2 * VTraits::vlanes(), v_add(w0, vx_load_aligned(wsum + j + 2 * VTraits::vlanes()))); - v_store_aligned(wsum + j + 3 * VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + 3 * VTraits::vlanes()))); - v_expand(val2, vall, valh); - v_store_aligned(sum_b + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_b + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + 3 * VTraits::vlanes()))); - v_expand(val3, vall, valh); - v_store_aligned(sum_g + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_g + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + 3 * VTraits::vlanes()))); - v_expand(val4, vall, valh); - v_store_aligned(sum_r + j + 2*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j + 2*VTraits::vlanes()))); - v_store_aligned(sum_r + j + 3*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + 3*VTraits::vlanes()))); - - v_load_deinterleave(ksptr1, kb, kg, kr); - v_expand(v_absdiff(kb, rb), val0, val1); - v_expand(v_absdiff(kg, rg), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - v_expand(v_absdiff(kr, rr), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - - v_expand(val0, vall, valh); - w0 = v_mul(kweight1, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight1, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j, v_add(w0, vx_load_aligned(wsum + j))); - v_store_aligned(wsum + j + VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + VTraits::vlanes()))); - v_expand(kb, val0, val2); - v_expand(val0, vall, valh); - v_store_aligned(sum_b + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_b + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + VTraits::vlanes()))); - v_expand(kg, val0, val3); - v_expand(val0, vall, valh); - v_store_aligned(sum_g + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_g + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + VTraits::vlanes()))); - v_expand(kr, val0, val4); - v_expand(val0, vall, valh); - v_store_aligned(sum_r + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j))); - v_store_aligned(sum_r + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + VTraits::vlanes()))); - - v_expand(val1, vall, valh); - w0 = v_mul(kweight1, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight1, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j + 2 * VTraits::vlanes(), v_add(w0, vx_load_aligned(wsum + j + 2 * VTraits::vlanes()))); - v_store_aligned(wsum + j + 3 * VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + 3 * VTraits::vlanes()))); - v_expand(val2, vall, valh); - v_store_aligned(sum_b + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_b + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + 3 * VTraits::vlanes()))); - v_expand(val3, vall, valh); - v_store_aligned(sum_g + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_g + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + 3 * VTraits::vlanes()))); - v_expand(val4, vall, valh); - v_store_aligned(sum_r + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_r + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + 3 * VTraits::vlanes()))); - - v_load_deinterleave(ksptr2, kb, kg, kr); - v_expand(v_absdiff(kb, rb), val0, val1); - v_expand(v_absdiff(kg, rg), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - v_expand(v_absdiff(kr, rr), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - - v_expand(val0, vall, valh); - w0 = v_mul(kweight2, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight2, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j, v_add(w0, vx_load_aligned(wsum + j))); - v_store_aligned(wsum + j + VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + VTraits::vlanes()))); - v_expand(kb, val0, val2); - v_expand(val0, vall, valh); - v_store_aligned(sum_b + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_b + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + VTraits::vlanes()))); - v_expand(kg, val0, val3); - v_expand(val0, vall, valh); - v_store_aligned(sum_g + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_g + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + VTraits::vlanes()))); - v_expand(kr, val0, val4); - v_expand(val0, vall, valh); - v_store_aligned(sum_r + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j))); - v_store_aligned(sum_r + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + VTraits::vlanes()))); - - v_expand(val1, vall, valh); - w0 = v_mul(kweight2, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight2, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j + 2 * VTraits::vlanes(), v_add(w0, vx_load_aligned(wsum + j + 2 * VTraits::vlanes()))); - v_store_aligned(wsum + j + 3 * VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + 3 * VTraits::vlanes()))); - v_expand(val2, vall, valh); - v_store_aligned(sum_b + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_b + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + 3 * VTraits::vlanes()))); - v_expand(val3, vall, valh); - v_store_aligned(sum_g + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_g + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + 3 * VTraits::vlanes()))); - v_expand(val4, vall, valh); - v_store_aligned(sum_r + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_r + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + 3 * VTraits::vlanes()))); - - v_load_deinterleave(ksptr3, kb, kg, kr); - v_expand(v_absdiff(kb, rb), val0, val1); - v_expand(v_absdiff(kg, rg), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - v_expand(v_absdiff(kr, rr), val2, val3); - val0 = v_add(val0, val2); val1 = v_add(val1, val3); - - v_expand(val0, vall, valh); - w0 = v_mul(kweight3, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight3, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j, v_add(w0, vx_load_aligned(wsum + j))); - v_store_aligned(wsum + j + VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + VTraits::vlanes()))); - v_expand(kb, val0, val2); - v_expand(val0, vall, valh); - v_store_aligned(sum_b + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_b + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + VTraits::vlanes()))); - v_expand(kg, val0, val3); - v_expand(val0, vall, valh); - v_store_aligned(sum_g + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_g + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + VTraits::vlanes()))); - v_expand(kr, val0, val4); - v_expand(val0, vall, valh); - v_store_aligned(sum_r + j, v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j))); - v_store_aligned(sum_r + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + VTraits::vlanes()))); - - v_expand(val1, vall, valh); - w0 = v_mul(kweight3, v_lut(color_weight, v_reinterpret_as_s32(vall))); - w1 = v_mul(kweight3, v_lut(color_weight, v_reinterpret_as_s32(valh))); - v_store_aligned(wsum + j + 2 * VTraits::vlanes(), v_add(w0, vx_load_aligned(wsum + j + 2 * VTraits::vlanes()))); - v_store_aligned(wsum + j + 3 * VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + 3 * VTraits::vlanes()))); - v_expand(val2, vall, valh); - v_store_aligned(sum_b + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_b + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_b + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_b + j + 3 * VTraits::vlanes()))); - v_expand(val3, vall, valh); - v_store_aligned(sum_g + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_g + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_g + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_g + j + 3 * VTraits::vlanes()))); - v_expand(val4, vall, valh); - v_store_aligned(sum_r + j + 2 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(vall)), w0, vx_load_aligned(sum_r + j + 2 * VTraits::vlanes()))); - v_store_aligned(sum_r + j + 3 * VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(valh)), w1, vx_load_aligned(sum_r + j + 3 * VTraits::vlanes()))); - } -#endif -#if CV_SIMD128 - v_float32x4 kweight4 = v_load(space_weight + k); -#endif - for(; j < size.width; j++, rsptr += 3, ksptr0 += 3, ksptr1 += 3, ksptr2 += 3, ksptr3 += 3) - { -#if CV_SIMD128 - v_uint32x4 rb = v_setall_u32(rsptr[0]); - v_uint32x4 rg = v_setall_u32(rsptr[1]); - v_uint32x4 rr = v_setall_u32(rsptr[2]); - v_uint32x4 b(ksptr0[0], ksptr1[0], ksptr2[0], ksptr3[0]); - v_uint32x4 g(ksptr0[1], ksptr1[1], ksptr2[1], ksptr3[1]); - v_uint32x4 r(ksptr0[2], ksptr1[2], ksptr2[2], ksptr3[2]); - v_float32x4 w = v_mul(kweight4, v_lut(this->color_weight, v_reinterpret_as_s32(v_add(v_add(v_absdiff(b, rb), v_absdiff(g, rg)), v_absdiff(r, rr))))); - wsum[j] += v_reduce_sum(w); - sum_b[j] += v_reduce_sum(v_mul(v_cvt_f32(v_reinterpret_as_s32(b)), w)); - sum_g[j] += v_reduce_sum(v_mul(v_cvt_f32(v_reinterpret_as_s32(g)), w)); - sum_r[j] += v_reduce_sum(v_mul(v_cvt_f32(v_reinterpret_as_s32(r)), w)); -#else - int rb = rsptr[0], rg = rsptr[1], rr = rsptr[2]; - - int b = ksptr0[0], g = ksptr0[1], r = ksptr0[2]; - float w = space_weight[k]*color_weight[std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)]; - wsum[j] += w; - sum_b[j] += b*w; sum_g[j] += g*w; sum_r[j] += r*w; - - b = ksptr1[0]; g = ksptr1[1]; r = ksptr1[2]; - w = space_weight[k+1] * color_weight[std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)]; - wsum[j] += w; - sum_b[j] += b*w; sum_g[j] += g*w; sum_r[j] += r*w; - - b = ksptr2[0]; g = ksptr2[1]; r = ksptr2[2]; - w = space_weight[k+2] * color_weight[std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)]; - wsum[j] += w; - sum_b[j] += b*w; sum_g[j] += g*w; sum_r[j] += r*w; - - b = ksptr3[0]; g = ksptr3[1]; r = ksptr3[2]; - w = space_weight[k+3] * color_weight[std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)]; - wsum[j] += w; - sum_b[j] += b*w; sum_g[j] += g*w; sum_r[j] += r*w; -#endif - } - } - for(; k < maxk; k++) + int n_8_lanes = VTraits::vlanes(); + for (; j <= size.width - n_8_lanes; j += n_8_lanes, sptr_j += 3*n_8_lanes, dptr += 3*n_8_lanes) { - const uchar* ksptr = sptr + space_ofs[k]; - const uchar* rsptr = sptr; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight = vx_setall_f32(space_weight[k]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), ksptr += 3*VTraits::vlanes(), rsptr += 3*VTraits::vlanes()) + const uchar* rsptr = sptr_j; + v_float32 v_wsum_0 = vx_setzero_f32(); + v_float32 v_wsum_1 = vx_setzero_f32(); + v_float32 v_wsum_2 = vx_setzero_f32(); + v_float32 v_wsum_3 = vx_setzero_f32(); + + v_float32 v_sum_b_0 = vx_setzero_f32(); + v_float32 v_sum_b_1 = vx_setzero_f32(); + v_float32 v_sum_b_2 = vx_setzero_f32(); + v_float32 v_sum_b_3 = vx_setzero_f32(); + + v_float32 v_sum_g_0 = vx_setzero_f32(); + v_float32 v_sum_g_1 = vx_setzero_f32(); + v_float32 v_sum_g_2 = vx_setzero_f32(); + v_float32 v_sum_g_3 = vx_setzero_f32(); + + v_float32 v_sum_r_0 = vx_setzero_f32(); + v_float32 v_sum_r_1 = vx_setzero_f32(); + v_float32 v_sum_r_2 = vx_setzero_f32(); + v_float32 v_sum_r_3 = vx_setzero_f32(); + + v_float32 v_one = vx_setall_f32(1.f); + + for(k=0; k < maxk; k++) { + const uchar* ksptr = sptr_j + space_ofs[k]; + v_float32 kweight = vx_setall_f32(space_weight[k]); + v_uint8 kb, kg, kr, rb, rg, rr; v_load_deinterleave(ksptr, kb, kg, kr); v_load_deinterleave(rsptr, rb, rg, rr); @@ -467,69 +394,71 @@ public: v_float32 w1 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(val1))); v_float32 w2 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(val2))); v_float32 w3 = v_mul(kweight, v_lut(color_weight, v_reinterpret_as_s32(val3))); - v_store_aligned(wsum + j , v_add(w0, vx_load_aligned(wsum + j))); - v_store_aligned(wsum + j + VTraits::vlanes(), v_add(w1, vx_load_aligned(wsum + j + VTraits::vlanes()))); - v_store_aligned(wsum + j + 2*VTraits::vlanes(), v_add(w2, vx_load_aligned(wsum + j + 2 * VTraits::vlanes()))); - v_store_aligned(wsum + j + 3*VTraits::vlanes(), v_add(w3, vx_load_aligned(wsum + j + 3 * VTraits::vlanes()))); + v_wsum_0 = v_add(w0, v_wsum_0); + v_wsum_1 = v_add(w1, v_wsum_1); + v_wsum_2 = v_add(w2, v_wsum_2); + v_wsum_3 = v_add(w3, v_wsum_3); + v_expand(b_l, val0, val1); v_expand(b_h, val2, val3); - v_store_aligned(sum_b + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_b + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, vx_load_aligned(sum_b + j + VTraits::vlanes()))); - v_store_aligned(sum_b + j + 2*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, vx_load_aligned(sum_b + j + 2*VTraits::vlanes()))); - v_store_aligned(sum_b + j + 3*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, vx_load_aligned(sum_b + j + 3*VTraits::vlanes()))); + v_sum_b_0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum_b_0); + v_sum_b_1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum_b_1); + v_sum_b_2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum_b_2); + v_sum_b_3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum_b_3); + v_expand(g_l, val0, val1); v_expand(g_h, val2, val3); - v_store_aligned(sum_g + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_g + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, vx_load_aligned(sum_g + j + VTraits::vlanes()))); - v_store_aligned(sum_g + j + 2*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, vx_load_aligned(sum_g + j + 2*VTraits::vlanes()))); - v_store_aligned(sum_g + j + 3*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, vx_load_aligned(sum_g + j + 3*VTraits::vlanes()))); + v_sum_g_0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum_g_0); + v_sum_g_1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum_g_1); + v_sum_g_2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum_g_2); + v_sum_g_3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum_g_3); + v_expand(r_l, val0, val1); v_expand(r_h, val2, val3); - v_store_aligned(sum_r + j , v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, vx_load_aligned(sum_r + j))); - v_store_aligned(sum_r + j + VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, vx_load_aligned(sum_r + j + VTraits::vlanes()))); - v_store_aligned(sum_r + j + 2*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, vx_load_aligned(sum_r + j + 2*VTraits::vlanes()))); - v_store_aligned(sum_r + j + 3*VTraits::vlanes(), v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, vx_load_aligned(sum_r + j + 3*VTraits::vlanes()))); + v_sum_r_0 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val0)), w0, v_sum_r_0); + v_sum_r_1 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val1)), w1, v_sum_r_1); + v_sum_r_2 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val2)), w2, v_sum_r_2); + v_sum_r_3 = v_muladd(v_cvt_f32(v_reinterpret_as_s32(val3)), w3, v_sum_r_3); } + v_float32 w0 = v_div(v_one, v_wsum_0); + v_float32 w1 = v_div(v_one, v_wsum_1); + v_float32 w2 = v_div(v_one, v_wsum_2); + v_float32 w3 = v_div(v_one, v_wsum_3); + + v_store_interleave(dptr, v_pack_u(v_pack(v_round(v_mul(w0, v_sum_b_0)), + v_round(v_mul(w1, v_sum_b_1))), + v_pack(v_round(v_mul(w2, v_sum_b_2)), + v_round(v_mul(w3, v_sum_b_3)))), + v_pack_u(v_pack(v_round(v_mul(w0, v_sum_g_0)), + v_round(v_mul(w1, v_sum_g_1))), + v_pack(v_round(v_mul(w2, v_sum_g_2)), + v_round(v_mul(w3, v_sum_g_3)))), + v_pack_u(v_pack(v_round(v_mul(w0, v_sum_r_0)), + v_round(v_mul(w1, v_sum_r_1))), + v_pack(v_round(v_mul(w2, v_sum_r_2)), + v_round(v_mul(w3, v_sum_r_3))))); + } #endif - for(; j < size.width; j++, ksptr += 3, rsptr += 3) + for(; j < size.width; j++, sptr_j += 3) + { + const uchar* rsptr = sptr_j; + float wsum = 0.f; + float sum_b = 0.f, sum_g = 0.f, sum_r = 0.f; + for(k=0; k < maxk; k++) { + const uchar* ksptr = sptr_j + space_ofs[k]; + int b = ksptr[0], g = ksptr[1], r = ksptr[2]; float w = space_weight[k]*color_weight[std::abs(b - rsptr[0]) + std::abs(g - rsptr[1]) + std::abs(r - rsptr[2])]; - wsum[j] += w; - sum_b[j] += b*w; sum_g[j] += g*w; sum_r[j] += r*w; + wsum += w; + sum_b += b*w; sum_g += g*w; sum_r += r*w; } - } - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 v_one = vx_setall_f32(1.f); - for(; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), dptr += 3*VTraits::vlanes()) - { - v_float32 w0 = v_div(v_one, vx_load_aligned(wsum + j)); - v_float32 w1 = v_div(v_one, vx_load_aligned(wsum + j + VTraits::vlanes())); - v_float32 w2 = v_div(v_one, vx_load_aligned(wsum + j + 2 * VTraits::vlanes())); - v_float32 w3 = v_div(v_one, vx_load_aligned(wsum + j + 3 * VTraits::vlanes())); - v_store_interleave(dptr, v_pack_u(v_pack(v_round(v_mul(w0, vx_load_aligned(sum_b + j))), - v_round(v_mul(w1, vx_load_aligned(sum_b + j + VTraits::vlanes())))), - v_pack(v_round(v_mul(w2, vx_load_aligned(sum_b + j + 2 * VTraits::vlanes()))), - v_round(v_mul(w3, vx_load_aligned(sum_b + j + 3 * VTraits::vlanes()))))), - v_pack_u(v_pack(v_round(v_mul(w0, vx_load_aligned(sum_g + j))), - v_round(v_mul(w1, vx_load_aligned(sum_g + j + VTraits::vlanes())))), - v_pack(v_round(v_mul(w2, vx_load_aligned(sum_g + j + 2 * VTraits::vlanes()))), - v_round(v_mul(w3, vx_load_aligned(sum_g + j + 3 * VTraits::vlanes()))))), - v_pack_u(v_pack(v_round(v_mul(w0, vx_load_aligned(sum_r + j))), - v_round(v_mul(w1, vx_load_aligned(sum_r + j + VTraits::vlanes())))), - v_pack(v_round(v_mul(w2, vx_load_aligned(sum_r + j + 2 * VTraits::vlanes()))), - v_round(v_mul(w3, vx_load_aligned(sum_r + j + 3 * VTraits::vlanes())))))); - } -#endif - for(; j < size.width; j++) - { - CV_DbgAssert(fabs(wsum[j]) > 0); - wsum[j] = 1.f/wsum[j]; - *(dptr++) = (uchar)cvRound(sum_b[j]*wsum[j]); - *(dptr++) = (uchar)cvRound(sum_g[j]*wsum[j]); - *(dptr++) = (uchar)cvRound(sum_r[j]*wsum[j]); + CV_DbgAssert(fabs(wsum) > 0); + wsum = 1.f/wsum; + *(dptr++) = (uchar)cvRound(sum_b*wsum); + *(dptr++) = (uchar)cvRound(sum_g*wsum); + *(dptr++) = (uchar)cvRound(sum_r*wsum); } } } @@ -552,6 +481,7 @@ void bilateralFilterInvoker_8u( int* space_ofs, float *space_weight, float *color_weight) { CV_INSTRUMENT_REGION(); + BilateralFilter_8u_Invoker body(dst, temp, radius, maxk, space_ofs, space_weight, color_weight); parallel_for_(Range(0, dst.rows), body, dst.total()/(double)(1<<16)); } @@ -577,7 +507,12 @@ public: int i, j, k; Size size = dest->size(); - +#if (CV_SIMD || CV_SIMD_SCALABLE) + int nlanes = VTraits::vlanes(); + int nlanes_2 = 2 * nlanes; + int nlanes_3 = 3 * nlanes; + int nlanes_4 = 4 * nlanes; +#endif for( i = range.start; i < range.end; i++ ) { const float* sptr = temp->ptr(i+radius) + radius*cn; @@ -585,363 +520,161 @@ public: if( cn == 1 ) { - AutoBuffer buf(alignSize(size.width, CV_SIMD_WIDTH) + size.width + CV_SIMD_WIDTH - 1); - memset(buf.data(), 0, buf.size() * sizeof(float)); - float *sum = alignPtr(buf.data(), CV_SIMD_WIDTH); - float *wsum = sum + alignSize(size.width, CV_SIMD_WIDTH); + j = 0; + const float* sptr_j = sptr; #if (CV_SIMD || CV_SIMD_SCALABLE) v_float32 v_one = vx_setall_f32(1.f); v_float32 sindex = vx_setall_f32(scale_index); -#endif - k = 0; - for(; k <= maxk - 4; k+=4) + + for(; j <= size.width - nlanes_4; j += nlanes_4, sptr_j += nlanes_4, dptr += nlanes_4) { - const float* ksptr0 = sptr + space_ofs[k]; - const float* ksptr1 = sptr + space_ofs[k + 1]; - const float* ksptr2 = sptr + space_ofs[k + 2]; - const float* ksptr3 = sptr + space_ofs[k + 3]; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight0 = vx_setall_f32(space_weight[k]); - v_float32 kweight1 = vx_setall_f32(space_weight[k+1]); - v_float32 kweight2 = vx_setall_f32(space_weight[k+2]); - v_float32 kweight3 = vx_setall_f32(space_weight[k+3]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes()) + v_float32 v_wsum0 = vx_setzero_f32(); + v_float32 v_wsum1 = vx_setzero_f32(); + v_float32 v_wsum2 = vx_setzero_f32(); + v_float32 v_wsum3 = vx_setzero_f32(); + v_float32 v_sum0 = vx_setzero_f32(); + v_float32 v_sum1 = vx_setzero_f32(); + v_float32 v_sum2 = vx_setzero_f32(); + v_float32 v_sum3 = vx_setzero_f32(); + + v_float32 rval0 = vx_load(sptr_j); + v_float32 rval1 = vx_load(sptr_j + nlanes); + v_float32 rval2 = vx_load(sptr_j + nlanes_2); + v_float32 rval3 = vx_load(sptr_j + nlanes_3); + for(k = 0; k < maxk; k++) { - v_float32 rval = vx_load(sptr + j); + const float* ksptr = sptr_j + space_ofs[k]; + v_float32 kweight = vx_setall_f32(space_weight[k]); - v_float32 val = vx_load(ksptr0 + j); - v_float32 knan = v_not_nan(val); - v_float32 alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex), v_not_nan(rval)), knan); - v_int32 idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - v_float32 w = v_and(v_mul(kweight0, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_float32 v_wsum = v_add(vx_load_aligned(wsum + j), w); - v_float32 v_sum = v_muladd(v_and(val, knan), w, vx_load_aligned(sum + j)); + //0th + v_float32 val0 = vx_load(ksptr); + v_float32 knan0 = v_not_nan(val0); + v_float32 alpha0 = v_and(v_and(v_mul(v_absdiff(val0, rval0), sindex), v_not_nan(rval0)), knan0); + v_int32 idx0 = v_trunc(alpha0); + alpha0 = v_sub(alpha0, v_cvt_f32(idx0)); + v_float32 w0 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx0), alpha0, v_mul(v_lut(this->expLUT, idx0), v_sub(v_one, alpha0)))), knan0); + v_wsum0 = v_add(v_wsum0, w0); + v_sum0 = v_muladd(v_and(val0, knan0), w0, v_sum0); - val = vx_load(ksptr1 + j); - knan = v_not_nan(val); - alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex), v_not_nan(rval)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight1, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_and(val, knan), w, v_sum); + //1st + v_float32 val1 = vx_load(ksptr + nlanes); + v_float32 knan1 = v_not_nan(val1); + v_float32 alpha1 = v_and(v_and(v_mul(v_absdiff(val1, rval1), sindex), v_not_nan(rval1)), knan1); + v_int32 idx1 = v_trunc(alpha1); + alpha1 = v_sub(alpha1, v_cvt_f32(idx1)); + v_float32 w1 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx1), alpha1, v_mul(v_lut(this->expLUT, idx1), v_sub(v_one, alpha1)))), knan1); + v_wsum1 = v_add(v_wsum1, w1); + v_sum1 = v_muladd(v_and(val1, knan1), w1, v_sum1); - val = vx_load(ksptr2 + j); - knan = v_not_nan(val); - alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex), v_not_nan(rval)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight2, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_and(val, knan), w, v_sum); + //2nd + v_float32 val2 = vx_load(ksptr + nlanes_2); + v_float32 knan2 = v_not_nan(val2); + v_float32 alpha2 = v_and(v_and(v_mul(v_absdiff(val2, rval2), sindex), v_not_nan(rval2)), knan2); + v_int32 idx2 = v_trunc(alpha2); + alpha2 = v_sub(alpha2, v_cvt_f32(idx2)); + v_float32 w2 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx2), alpha2, v_mul(v_lut(this->expLUT, idx2), v_sub(v_one, alpha2)))), knan2); + v_wsum2 = v_add(v_wsum2, w2); + v_sum2 = v_muladd(v_and(val2, knan2), w2, v_sum2); - val = vx_load(ksptr3 + j); - knan = v_not_nan(val); - alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex), v_not_nan(rval)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight3, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum = v_muladd(v_and(val, knan), w, v_sum); - - v_store_aligned(wsum + j, v_wsum); - v_store_aligned(sum + j, v_sum); + //3rd + v_float32 val3 = vx_load(ksptr + nlanes_3); + v_float32 knan3 = v_not_nan(val3); + v_float32 alpha3 = v_and(v_and(v_mul(v_absdiff(val3, rval3), sindex), v_not_nan(rval3)), knan3); + v_int32 idx3 = v_trunc(alpha3); + alpha3 = v_sub(alpha3, v_cvt_f32(idx3)); + v_float32 w3 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx3), alpha3, v_mul(v_lut(this->expLUT, idx3), v_sub(v_one, alpha3)))), knan3); + v_wsum3 = v_add(v_wsum3, w3); + v_sum3 = v_muladd(v_and(val3, knan3), w3, v_sum3); } -#endif -#if CV_SIMD128 - v_float32x4 v_one4 = v_setall_f32(1.f); - v_float32x4 sindex4 = v_setall_f32(scale_index); - v_float32x4 kweight4 = v_load(space_weight + k); -#endif - for (; j < size.width; j++) - { -#if CV_SIMD128 - v_float32x4 rval = v_setall_f32(sptr[j]); - v_float32x4 val(ksptr0[j], ksptr1[j], ksptr2[j], ksptr3[j]); - v_float32x4 knan = v_not_nan(val); - v_float32x4 alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex4), v_not_nan(rval)), knan); - v_int32x4 idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - v_float32x4 w = v_and(v_mul(kweight4, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one4, alpha)))), knan); - wsum[j] += v_reduce_sum(w); - sum[j] += v_reduce_sum(v_mul(v_and(val, knan), w)); -#else - float rval = sptr[j]; + v_store(dptr , v_div(v_add(v_sum0, v_and(rval0, v_not_nan(rval0))), v_add(v_wsum0, v_and(v_one, v_not_nan(rval0))))); + v_store(dptr + nlanes, v_div(v_add(v_sum1, v_and(rval1, v_not_nan(rval1))), v_add(v_wsum1, v_and(v_one, v_not_nan(rval1))))); + v_store(dptr + nlanes_2, v_div(v_add(v_sum2, v_and(rval2, v_not_nan(rval2))), v_add(v_wsum2, v_and(v_one, v_not_nan(rval2))))); + v_store(dptr + nlanes_3, v_div(v_add(v_sum3, v_and(rval3, v_not_nan(rval3))), v_add(v_wsum3, v_and(v_one, v_not_nan(rval3))))); + } + for (; j <= size.width - nlanes_2; j += nlanes_2, sptr_j += nlanes_2, dptr += nlanes_2) + { + v_float32 v_wsum0 = vx_setzero_f32(); + v_float32 v_wsum1 = vx_setzero_f32(); + v_float32 v_sum0 = vx_setzero_f32(); + v_float32 v_sum1 = vx_setzero_f32(); + v_float32 rval0 = vx_load(sptr_j); + v_float32 rval1 = vx_load(sptr_j + nlanes); + v_float32 rval0_not_nan = v_not_nan(rval0); + v_float32 rval1_not_nan = v_not_nan(rval1); - float val = ksptr0[j]; + for (k = 0; k < maxk; k++) + { + v_float32 kweight = vx_setall_f32(space_weight[k]); + const float* ksptr = sptr_j + space_ofs[k]; + + //0th + v_float32 val0 = vx_load(ksptr); + v_float32 knan0 = v_not_nan(val0); + v_float32 alpha0 = v_and(v_and(v_mul(v_absdiff(val0, rval0), sindex), rval0_not_nan), knan0); + v_int32 idx0 = v_trunc(alpha0); + alpha0 = v_sub(alpha0, v_cvt_f32(idx0)); + v_float32 w0 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx0), alpha0, v_mul(v_lut(this->expLUT, idx0), v_sub(v_one, alpha0)))), knan0); + v_wsum0 = v_add(v_wsum0, w0); + v_sum0 = v_muladd(v_and(val0, knan0), w0, v_sum0); + + //1st + v_float32 val1 = vx_load(ksptr + nlanes); + v_float32 knan1 = v_not_nan(val1); + v_float32 alpha1 = v_and(v_and(v_mul(v_absdiff(val1, rval1), sindex), rval1_not_nan), knan1); + v_int32 idx1 = v_trunc(alpha1); + alpha1 = v_sub(alpha1, v_cvt_f32(idx1)); + v_float32 w1 = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx1), alpha1, v_mul(v_lut(this->expLUT, idx1), v_sub(v_one, alpha1)))), knan1); + v_wsum1 = v_add(v_wsum1, w1); + v_sum1 = v_muladd(v_and(val1, knan1), w1, v_sum1); + } + v_store(dptr, v_div(v_add(v_sum0, v_and(rval0, rval0_not_nan)), v_add(v_wsum0, v_and(v_one, rval0_not_nan)))); + v_store(dptr + nlanes, v_div(v_add(v_sum1, v_and(rval1, rval1_not_nan)), v_add(v_wsum1, v_and(v_one, rval1_not_nan)))); + } +#endif + for (; j < size.width; j++, sptr_j++, dptr++) + { + float rval = *sptr_j; + float wsum = 0.f; + float sum = 0.f; + for (k = 0; k < maxk; k++) + { + const float* ksptr = sptr_j + space_ofs[k]; + float val = *ksptr; float alpha = std::abs(val - rval) * scale_index; int idx = cvFloor(alpha); alpha -= idx; if (!cvIsNaN(val)) { - float w = space_weight[k] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum[j] += val * w; - } - - val = ksptr1[j]; - alpha = std::abs(val - rval) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!cvIsNaN(val)) - { - float w = space_weight[k+1] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum[j] += val * w; - } - - val = ksptr2[j]; - alpha = std::abs(val - rval) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!cvIsNaN(val)) - { - float w = space_weight[k+2] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum[j] += val * w; - } - - val = ksptr3[j]; - alpha = std::abs(val - rval) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!cvIsNaN(val)) - { - float w = space_weight[k+3] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum[j] += val * w; - } -#endif - } - } - for(; k < maxk; k++) - { - const float* ksptr = sptr + space_ofs[k]; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight = vx_setall_f32(space_weight[k]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes()) - { - v_float32 val = vx_load(ksptr + j); - v_float32 rval = vx_load(sptr + j); - v_float32 knan = v_not_nan(val); - v_float32 alpha = v_and(v_and(v_mul(v_absdiff(val, rval), sindex), v_not_nan(rval)), knan); - v_int32 idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - - v_float32 w = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_store_aligned(wsum + j, v_add(vx_load_aligned(wsum + j), w)); - v_store_aligned(sum + j, v_muladd(v_and(val, knan), w, vx_load_aligned(sum + j))); - } -#endif - for (; j < size.width; j++) - { - float val = ksptr[j]; - float rval = sptr[j]; - float alpha = std::abs(val - rval) * scale_index; - int idx = cvFloor(alpha); - alpha -= idx; - if (!cvIsNaN(val)) - { - float w = space_weight[k] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum[j] += val * w; + float w = space_weight[k] * (cvIsNaN(rval) ? 1.f : (expLUT[idx] + alpha * (expLUT[idx + 1] - expLUT[idx]))); + wsum += w; + sum += val * w; } } - } - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes()) - { - v_float32 v_val = vx_load(sptr + j); - v_store(dptr + j, v_div(v_add(vx_load_aligned(sum + j), v_and(v_val, v_not_nan(v_val))), v_add(vx_load_aligned(wsum + j), v_and(v_one, v_not_nan(v_val))))); - } -#endif - for (; j < size.width; j++) - { - CV_DbgAssert(fabs(wsum[j]) >= 0); - dptr[j] = cvIsNaN(sptr[j]) ? sum[j] / wsum[j] : (sum[j] + sptr[j]) / (wsum[j] + 1.f); + CV_DbgAssert(fabs(wsum) >= 0); + *dptr = cvIsNaN(rval) ? sum / wsum : (sum + rval) / (wsum + 1.f); } } else { - CV_Assert( cn == 3 ); - AutoBuffer buf(alignSize(size.width, CV_SIMD_WIDTH)*3 + size.width + CV_SIMD_WIDTH - 1); - memset(buf.data(), 0, buf.size() * sizeof(float)); - float *sum_b = alignPtr(buf.data(), CV_SIMD_WIDTH); - float *sum_g = sum_b + alignSize(size.width, CV_SIMD_WIDTH); - float *sum_r = sum_g + alignSize(size.width, CV_SIMD_WIDTH); - float *wsum = sum_r + alignSize(size.width, CV_SIMD_WIDTH); + CV_Assert(cn == 3); + j = 0; + const float* sptr_j = sptr; #if (CV_SIMD || CV_SIMD_SCALABLE) v_float32 v_one = vx_setall_f32(1.f); v_float32 sindex = vx_setall_f32(scale_index); -#endif - k = 0; - for (; k <= maxk-4; k+=4) + + for (; j <= size.width - nlanes; j += nlanes, sptr_j += nlanes_3, dptr += nlanes_3) { - const float* ksptr0 = sptr + space_ofs[k]; - const float* ksptr1 = sptr + space_ofs[k+1]; - const float* ksptr2 = sptr + space_ofs[k+2]; - const float* ksptr3 = sptr + space_ofs[k+3]; - const float* rsptr = sptr; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight0 = vx_setall_f32(space_weight[k]); - v_float32 kweight1 = vx_setall_f32(space_weight[k+1]); - v_float32 kweight2 = vx_setall_f32(space_weight[k+2]); - v_float32 kweight3 = vx_setall_f32(space_weight[k+3]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), rsptr += 3 * VTraits::vlanes(), - ksptr0 += 3 * VTraits::vlanes(), ksptr1 += 3 * VTraits::vlanes(), ksptr2 += 3 * VTraits::vlanes(), ksptr3 += 3 * VTraits::vlanes()) + const float* rsptr = sptr_j; + v_float32 v_wsum = vx_setzero_f32(); + v_float32 v_sum_b = vx_setzero_f32(); + v_float32 v_sum_g = vx_setzero_f32(); + v_float32 v_sum_r = vx_setzero_f32(); + for (k = 0; k < maxk; k++) { - v_float32 kb, kg, kr, rb, rg, rr; - v_load_deinterleave(rsptr, rb, rg, rr); + const float* ksptr = sptr_j + space_ofs[k]; + v_float32 kweight = vx_setall_f32(space_weight[k]); - v_load_deinterleave(ksptr0, kb, kg, kr); - v_float32 knan = v_and(v_and(v_not_nan(kb), v_not_nan(kg)), v_not_nan(kr)); - v_float32 alpha = v_and(v_and(v_and(v_and(v_mul(v_add(v_add(v_absdiff(kb, rb), v_absdiff(kg, rg)), v_absdiff(kr, rr)), sindex), v_not_nan(rb)), v_not_nan(rg)), v_not_nan(rr)), knan); - v_int32 idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - v_float32 w = v_and(v_mul(kweight0, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_float32 v_wsum = v_add(vx_load_aligned(wsum + j), w); - v_float32 v_sum_b = v_muladd(v_and(kb, knan), w, vx_load_aligned(sum_b + j)); - v_float32 v_sum_g = v_muladd(v_and(kg, knan), w, vx_load_aligned(sum_g + j)); - v_float32 v_sum_r = v_muladd(v_and(kr, knan), w, vx_load_aligned(sum_r + j)); - - v_load_deinterleave(ksptr1, kb, kg, kr); - knan = v_and(v_and(v_not_nan(kb), v_not_nan(kg)), v_not_nan(kr)); - alpha = v_and(v_and(v_and(v_and(v_mul(v_add(v_add(v_absdiff(kb, rb), v_absdiff(kg, rg)), v_absdiff(kr, rr)), sindex), v_not_nan(rb)), v_not_nan(rg)), v_not_nan(rr)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight1, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum_b = v_muladd(v_and(kb, knan), w, v_sum_b); - v_sum_g = v_muladd(v_and(kg, knan), w, v_sum_g); - v_sum_r = v_muladd(v_and(kr, knan), w, v_sum_r); - - v_load_deinterleave(ksptr2, kb, kg, kr); - knan = v_and(v_and(v_not_nan(kb), v_not_nan(kg)), v_not_nan(kr)); - alpha = v_and(v_and(v_and(v_and(v_mul(v_add(v_add(v_absdiff(kb, rb), v_absdiff(kg, rg)), v_absdiff(kr, rr)), sindex), v_not_nan(rb)), v_not_nan(rg)), v_not_nan(rr)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight2, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum_b = v_muladd(v_and(kb, knan), w, v_sum_b); - v_sum_g = v_muladd(v_and(kg, knan), w, v_sum_g); - v_sum_r = v_muladd(v_and(kr, knan), w, v_sum_r); - - v_load_deinterleave(ksptr3, kb, kg, kr); - knan = v_and(v_and(v_not_nan(kb), v_not_nan(kg)), v_not_nan(kr)); - alpha = v_and(v_and(v_and(v_and(v_mul(v_add(v_add(v_absdiff(kb, rb), v_absdiff(kg, rg)), v_absdiff(kr, rr)), sindex), v_not_nan(rb)), v_not_nan(rg)), v_not_nan(rr)), knan); - idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - w = v_and(v_mul(kweight3, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_wsum = v_add(v_wsum, w); - v_sum_b = v_muladd(v_and(kb, knan), w, v_sum_b); - v_sum_g = v_muladd(v_and(kg, knan), w, v_sum_g); - v_sum_r = v_muladd(v_and(kr, knan), w, v_sum_r); - - v_store_aligned(wsum + j, v_wsum); - v_store_aligned(sum_b + j, v_sum_b); - v_store_aligned(sum_g + j, v_sum_g); - v_store_aligned(sum_r + j, v_sum_r); - } -#endif -#if CV_SIMD128 - v_float32x4 v_one4 = v_setall_f32(1.f); - v_float32x4 sindex4 = v_setall_f32(scale_index); - v_float32x4 kweight4 = v_load(space_weight + k); -#endif - for (; j < size.width; j++, rsptr += 3, ksptr0 += 3, ksptr1 += 3, ksptr2 += 3, ksptr3 += 3) - { -#if CV_SIMD128 - v_float32x4 rb = v_setall_f32(rsptr[0]); - v_float32x4 rg = v_setall_f32(rsptr[1]); - v_float32x4 rr = v_setall_f32(rsptr[2]); - v_float32x4 kb(ksptr0[0], ksptr1[0], ksptr2[0], ksptr3[0]); - v_float32x4 kg(ksptr0[1], ksptr1[1], ksptr2[1], ksptr3[1]); - v_float32x4 kr(ksptr0[2], ksptr1[2], ksptr2[2], ksptr3[2]); - v_float32x4 knan = v_and(v_and(v_not_nan(kb), v_not_nan(kg)), v_not_nan(kr)); - v_float32x4 alpha = v_and(v_and(v_and(v_and(v_mul(v_add(v_add(v_absdiff(kb, rb), v_absdiff(kg, rg)), v_absdiff(kr, rr)), sindex4), v_not_nan(rb)), v_not_nan(rg)), v_not_nan(rr)), knan); - v_int32x4 idx = v_trunc(alpha); - alpha = v_sub(alpha, v_cvt_f32(idx)); - v_float32x4 w = v_and(v_mul(kweight4, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one4, alpha)))), knan); - wsum[j] += v_reduce_sum(w); - sum_b[j] += v_reduce_sum(v_mul(v_and(kb, knan), w)); - sum_g[j] += v_reduce_sum(v_mul(v_and(kg, knan), w)); - sum_r[j] += v_reduce_sum(v_mul(v_and(kr, knan), w)); -#else - float rb = rsptr[0], rg = rsptr[1], rr = rsptr[2]; - bool r_NAN = cvIsNaN(rb) || cvIsNaN(rg) || cvIsNaN(rr); - - float b = ksptr0[0], g = ksptr0[1], r = ksptr0[2]; - bool v_NAN = cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r); - float alpha = (std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)) * scale_index; - int idx = cvFloor(alpha); - alpha -= idx; - if (!v_NAN) - { - float w = space_weight[k] * (r_NAN ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum_b[j] += b*w; - sum_g[j] += g*w; - sum_r[j] += r*w; - } - - b = ksptr1[0]; g = ksptr1[1]; r = ksptr1[2]; - v_NAN = cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r); - alpha = (std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!v_NAN) - { - float w = space_weight[k+1] * (r_NAN ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum_b[j] += b*w; - sum_g[j] += g*w; - sum_r[j] += r*w; - } - - b = ksptr2[0]; g = ksptr2[1]; r = ksptr2[2]; - v_NAN = cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r); - alpha = (std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!v_NAN) - { - float w = space_weight[k+2] * (r_NAN ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum_b[j] += b*w; - sum_g[j] += g*w; - sum_r[j] += r*w; - } - - b = ksptr3[0]; g = ksptr3[1]; r = ksptr3[2]; - v_NAN = cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r); - alpha = (std::abs(b - rb) + std::abs(g - rg) + std::abs(r - rr)) * scale_index; - idx = cvFloor(alpha); - alpha -= idx; - if (!v_NAN) - { - float w = space_weight[k+3] * (r_NAN ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum_b[j] += b*w; - sum_g[j] += g*w; - sum_r[j] += r*w; - } -#endif - } - } - for (; k < maxk; k++) - { - const float* ksptr = sptr + space_ofs[k]; - const float* rsptr = sptr; - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - v_float32 kweight = vx_setall_f32(space_weight[k]); - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), ksptr += 3*VTraits::vlanes(), rsptr += 3*VTraits::vlanes()) - { v_float32 kb, kg, kr, rb, rg, rr; v_load_deinterleave(ksptr, kb, kg, kr); v_load_deinterleave(rsptr, rb, rg, rr); @@ -952,14 +685,26 @@ public: alpha = v_sub(alpha, v_cvt_f32(idx)); v_float32 w = v_and(v_mul(kweight, v_muladd(v_lut(this->expLUT + 1, idx), alpha, v_mul(v_lut(this->expLUT, idx), v_sub(v_one, alpha)))), knan); - v_store_aligned(wsum + j, v_add(vx_load_aligned(wsum + j), w)); - v_store_aligned(sum_b + j, v_muladd(v_and(kb, knan), w, vx_load_aligned(sum_b + j))); - v_store_aligned(sum_g + j, v_muladd(v_and(kg, knan), w, vx_load_aligned(sum_g + j))); - v_store_aligned(sum_r + j, v_muladd(v_and(kr, knan), w, vx_load_aligned(sum_r + j))); + v_wsum = v_add(v_wsum, w); + v_sum_b = v_muladd(v_and(kb, knan), w, v_sum_b); + v_sum_g = v_muladd(v_and(kg, knan), w, v_sum_g); + v_sum_r = v_muladd(v_and(kr, knan), w, v_sum_r); } + + v_float32 b, g, r; + v_load_deinterleave(sptr_j, b, g, r); + v_float32 mask = v_and(v_and(v_not_nan(b), v_not_nan(g)), v_not_nan(r)); + v_float32 w = v_div(v_one, v_add(v_wsum, v_and(v_one, mask))); + v_store_interleave(dptr, v_mul(v_add(v_sum_b, v_and(b, mask)), w), v_mul(v_add(v_sum_g, v_and(g, mask)), w), v_mul(v_add(v_sum_r, v_and(r, mask)), w)); + } #endif - for (; j < size.width; j++, ksptr += 3, rsptr += 3) + for (; j < size.width; j++, sptr_j += 3) + { + const float* rsptr = sptr_j; + float wsum = 0.f, sum_b = 0.f, sum_g = 0.f, sum_r = 0.f; + for (k = 0; k < maxk; k++) { + const float* ksptr = sptr_j + space_ofs[k]; float b = ksptr[0], g = ksptr[1], r = ksptr[2]; bool v_NAN = cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r); float rb = rsptr[0], rg = rsptr[1], rr = rsptr[2]; @@ -969,44 +714,31 @@ public: alpha -= idx; if (!v_NAN) { - float w = space_weight[k] * (r_NAN ? 1.f : (expLUT[idx] + alpha*(expLUT[idx + 1] - expLUT[idx]))); - wsum[j] += w; - sum_b[j] += b*w; - sum_g[j] += g*w; - sum_r[j] += r*w; + float w = space_weight[k] * (r_NAN ? 1.f : (expLUT[idx] + alpha * (expLUT[idx + 1] - expLUT[idx]))); + wsum += w; + sum_b += b * w; + sum_g += g * w; + sum_r += r * w; } } - } - j = 0; -#if (CV_SIMD || CV_SIMD_SCALABLE) - for (; j <= size.width - VTraits::vlanes(); j += VTraits::vlanes(), sptr += 3*VTraits::vlanes(), dptr += 3*VTraits::vlanes()) - { - v_float32 b, g, r; - v_load_deinterleave(sptr, b, g, r); - v_float32 mask = v_and(v_and(v_not_nan(b), v_not_nan(g)), v_not_nan(r)); - v_float32 w = v_div(v_one, v_add(vx_load_aligned(wsum + j), v_and(v_one, mask))); - v_store_interleave(dptr, v_mul(v_add(vx_load_aligned(sum_b + j), v_and(b, mask)), w), v_mul(v_add(vx_load_aligned(sum_g + j), v_and(g, mask)), w), v_mul(v_add(vx_load_aligned(sum_r + j), v_and(r, mask)), w)); - } -#endif - for (; j < size.width; j++) - { - CV_DbgAssert(fabs(wsum[j]) >= 0); - float b = *(sptr++); - float g = *(sptr++); - float r = *(sptr++); + + CV_DbgAssert(fabs(wsum) >= 0); + float b = *(sptr_j); + float g = *(sptr_j+1); + float r = *(sptr_j+2); if (cvIsNaN(b) || cvIsNaN(g) || cvIsNaN(r)) { - wsum[j] = 1.f / wsum[j]; - *(dptr++) = sum_b[j] * wsum[j]; - *(dptr++) = sum_g[j] * wsum[j]; - *(dptr++) = sum_r[j] * wsum[j]; + wsum = 1.f / wsum; + *(dptr++) = sum_b * wsum; + *(dptr++) = sum_g * wsum; + *(dptr++) = sum_r * wsum; } else { - wsum[j] = 1.f / (wsum[j] + 1.f); - *(dptr++) = (sum_b[j] + b) * wsum[j]; - *(dptr++) = (sum_g[j] + g) * wsum[j]; - *(dptr++) = (sum_r[j] + r) * wsum[j]; + wsum = 1.f / (wsum + 1.f); + *(dptr++) = (sum_b + b) * wsum; + *(dptr++) = (sum_g + g) * wsum; + *(dptr++) = (sum_r + r) * wsum; } } }