From 1c44eaf8bd404640382dd2ed5ae1660de73e6c91 Mon Sep 17 00:00:00 2001 From: Alexander Smorkalov Date: Sun, 24 May 2026 17:28:53 +0300 Subject: [PATCH] Move divSpectrums to core module. --- modules/core/include/opencv2/core.hpp | 16 ++ modules/core/misc/js/gen_dict.json | 4 +- modules/core/src/dxt.cpp | 196 ++++++++++++++ modules/core/test/test_dxt.cpp | 274 ++++++++++++++++++++ modules/imgproc/include/opencv2/imgproc.hpp | 16 -- modules/imgproc/misc/js/gen_dict.json | 1 - modules/imgproc/src/phasecorr.cpp | 198 -------------- modules/imgproc/test/test_pc.cpp | 274 -------------------- platforms/js/opencv_js.config.py | 3 +- 9 files changed, 489 insertions(+), 493 deletions(-) diff --git a/modules/core/include/opencv2/core.hpp b/modules/core/include/opencv2/core.hpp index f9c575cc20..d096625e08 100644 --- a/modules/core/include/opencv2/core.hpp +++ b/modules/core/include/opencv2/core.hpp @@ -2284,6 +2284,22 @@ or not (false). CV_EXPORTS_W void mulSpectrums(InputArray a, InputArray b, OutputArray c, int flags, bool conjB = false); +/** @brief Performs the per-element division of the first Fourier spectrum by the second Fourier spectrum. + * + * The function cv::divSpectrums performs the per-element division of the first array by the second array. + * The arrays are CCS-packed or complex matrices that are results of a real or complex Fourier transform. + * + * @param a first input array. + * @param b second input array of the same size and type as src1 . + * @param c output array of the same size and type as src1 . + * @param flags operation flags; currently, the only supported flag is cv::DFT_ROWS, which indicates that + * each row of src1 and src2 is an independent 1D Fourier spectrum. If you do not want to use this flag, then simply add a `0` as value. + * @param conjB optional flag that conjugates the second input array before the multiplication (true) + * or not (false). + */ +CV_EXPORTS_W void divSpectrums(InputArray a, InputArray b, OutputArray c, + int flags, bool conjB = false); + /** @brief Returns the optimal DFT size for a given vector size. DFT performance is not a monotonic function of a vector size. Therefore, when you calculate diff --git a/modules/core/misc/js/gen_dict.json b/modules/core/misc/js/gen_dict.json index f1e892d86b..723dee403c 100644 --- a/modules/core/misc/js/gen_dict.json +++ b/modules/core/misc/js/gen_dict.json @@ -3,8 +3,8 @@ { "": [ "absdiff", "add", "addWeighted", "bitwise_and", "bitwise_not", "bitwise_or", "bitwise_xor", "cartToPolar", - "compare", "convertScaleAbs", "copyMakeBorder", "countNonZero", "determinant", "dft", "divide", "eigen", - "exp", "flip", "getOptimalDFTSize","gemm", "hconcat", "inRange", "invert", "kmeans", "log", "magnitude", + "compare", "convertScaleAbs", "copyMakeBorder", "countNonZero", "determinant", "dft", "divide", "divSpectrums", + "eigen", "exp", "flip", "getOptimalDFTSize","gemm", "hconcat", "inRange", "invert", "kmeans", "log", "magnitude", "max", "mean", "meanStdDev", "merge", "min", "minMaxLoc", "mixChannels", "multiply", "norm", "normalize", "perspectiveTransform", "polarToCart", "pow", "randn", "randu", "reduce", "repeat", "rotate", "setIdentity", "setRNGSeed", "solve", "solvePoly", "split", "sqrt", "subtract", "trace", "transform", "transpose", "vconcat", diff --git a/modules/core/src/dxt.cpp b/modules/core/src/dxt.cpp index a4fad757a5..9d95f5d901 100644 --- a/modules/core/src/dxt.cpp +++ b/modules/core/src/dxt.cpp @@ -3773,6 +3773,202 @@ void cv::mulSpectrums( InputArray _srcA, InputArray _srcB, } } +void cv::divSpectrums( InputArray _srcA, InputArray _srcB, OutputArray _dst, int flags, bool conjB) +{ + Mat srcA = _srcA.getMat(), srcB = _srcB.getMat(); + int depth = srcA.depth(), cn = srcA.channels(), type = srcA.type(); + int rows = srcA.rows, cols = srcA.cols; + int j, k; + + CV_Assert( type == srcB.type() && srcA.size() == srcB.size() ); + CV_Assert( type == CV_32FC1 || type == CV_32FC2 || type == CV_64FC1 || type == CV_64FC2 ); + + _dst.create( srcA.rows, srcA.cols, type ); + Mat dst = _dst.getMat(); + + CV_Assert(dst.data != srcA.data); // non-inplace check + CV_Assert(dst.data != srcB.data); // non-inplace check + + bool is_1d = (flags & DFT_ROWS) || (rows == 1 || (cols == 1 && + srcA.isContinuous() && srcB.isContinuous() && dst.isContinuous())); + + if( is_1d && !(flags & DFT_ROWS) ) + cols = cols + rows - 1, rows = 1; + + int ncols = cols*cn; + int j0 = cn == 1; + int j1 = ncols - (cols % 2 == 0 && cn == 1); + + if( depth == CV_32F ) + { + const float* dataA = srcA.ptr(); + const float* dataB = srcB.ptr(); + float* dataC = dst.ptr(); + float eps = FLT_EPSILON; // prevent div0 problems + + size_t stepA = srcA.step/sizeof(dataA[0]); + size_t stepB = srcB.step/sizeof(dataB[0]); + size_t stepC = dst.step/sizeof(dataC[0]); + + if( !is_1d && cn == 1 ) + { + for( k = 0; k < (cols % 2 ? 1 : 2); k++ ) + { + if( k == 1 ) + dataA += cols - 1, dataB += cols - 1, dataC += cols - 1; + dataC[0] = dataA[0] / (dataB[0] + eps); + if( rows % 2 == 0 ) + dataC[(rows-1)*stepC] = dataA[(rows-1)*stepA] / (dataB[(rows-1)*stepB] + eps); + if( !conjB ) + for( j = 1; j <= rows - 2; j += 2 ) + { + double denom = (double)dataB[j*stepB]*dataB[j*stepB] + + (double)dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + (double)eps; + + double re = (double)dataA[j*stepA]*dataB[j*stepB] + + (double)dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; + + double im = (double)dataA[(j+1)*stepA]*dataB[j*stepB] - + (double)dataA[j*stepA]*dataB[(j+1)*stepB]; + + dataC[j*stepC] = (float)(re / denom); + dataC[(j+1)*stepC] = (float)(im / denom); + } + else + for( j = 1; j <= rows - 2; j += 2 ) + { + + double denom = (double)dataB[j*stepB]*dataB[j*stepB] + + (double)dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + (double)eps; + + double re = (double)dataA[j*stepA]*dataB[j*stepB] - + (double)dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; + + double im = (double)dataA[(j+1)*stepA]*dataB[j*stepB] + + (double)dataA[j*stepA]*dataB[(j+1)*stepB]; + + dataC[j*stepC] = (float)(re / denom); + dataC[(j+1)*stepC] = (float)(im / denom); + } + if( k == 1 ) + dataA -= cols - 1, dataB -= cols - 1, dataC -= cols - 1; + } + } + + for( ; rows--; dataA += stepA, dataB += stepB, dataC += stepC ) + { + if( is_1d && cn == 1 ) + { + dataC[0] = dataA[0] / (dataB[0] + eps); + if( cols % 2 == 0 ) + dataC[j1] = dataA[j1] / (dataB[j1] + eps); + } + + if( !conjB ) + for( j = j0; j < j1; j += 2 ) + { + double denom = (double)dataB[j]*dataB[j] + (double)dataB[j+1]*dataB[j+1] + (double)eps; + double re = (double)dataA[j]*dataB[j] + (double)dataA[j+1]*dataB[j+1]; + double im = (double)dataA[j+1]*dataB[j] - (double)dataA[j]*dataB[j+1]; + dataC[j] = (float)(re / denom); + dataC[j+1] = (float)(im / denom); + } + else + for( j = j0; j < j1; j += 2 ) + { + double denom = (double)dataB[j]*dataB[j] + (double)dataB[j+1]*dataB[j+1] + (double)eps; + double re = (double)dataA[j]*dataB[j] - (double)dataA[j+1]*dataB[j+1]; + double im = (double)dataA[j+1]*dataB[j] + (double)dataA[j]*dataB[j+1]; + dataC[j] = (float)(re / denom); + dataC[j+1] = (float)(im / denom); + } + } + } + else + { + const double* dataA = srcA.ptr(); + const double* dataB = srcB.ptr(); + double* dataC = dst.ptr(); + double eps = DBL_EPSILON; // prevent div0 problems + + size_t stepA = srcA.step/sizeof(dataA[0]); + size_t stepB = srcB.step/sizeof(dataB[0]); + size_t stepC = dst.step/sizeof(dataC[0]); + + if( !is_1d && cn == 1 ) + { + for( k = 0; k < (cols % 2 ? 1 : 2); k++ ) + { + if( k == 1 ) + dataA += cols - 1, dataB += cols - 1, dataC += cols - 1; + dataC[0] = dataA[0] / (dataB[0] + eps); + if( rows % 2 == 0 ) + dataC[(rows-1)*stepC] = dataA[(rows-1)*stepA] / (dataB[(rows-1)*stepB] + eps); + if( !conjB ) + for( j = 1; j <= rows - 2; j += 2 ) + { + double denom = dataB[j*stepB]*dataB[j*stepB] + + dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + eps; + + double re = dataA[j*stepA]*dataB[j*stepB] + + dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; + + double im = dataA[(j+1)*stepA]*dataB[j*stepB] - + dataA[j*stepA]*dataB[(j+1)*stepB]; + + dataC[j*stepC] = re / denom; + dataC[(j+1)*stepC] = im / denom; + } + else + for( j = 1; j <= rows - 2; j += 2 ) + { + double denom = dataB[j*stepB]*dataB[j*stepB] + + dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + eps; + + double re = dataA[j*stepA]*dataB[j*stepB] - + dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; + + double im = dataA[(j+1)*stepA]*dataB[j*stepB] + + dataA[j*stepA]*dataB[(j+1)*stepB]; + + dataC[j*stepC] = re / denom; + dataC[(j+1)*stepC] = im / denom; + } + if( k == 1 ) + dataA -= cols - 1, dataB -= cols - 1, dataC -= cols - 1; + } + } + + for( ; rows--; dataA += stepA, dataB += stepB, dataC += stepC ) + { + if( is_1d && cn == 1 ) + { + dataC[0] = dataA[0] / (dataB[0] + eps); + if( cols % 2 == 0 ) + dataC[j1] = dataA[j1] / (dataB[j1] + eps); + } + + if( !conjB ) + for( j = j0; j < j1; j += 2 ) + { + double denom = dataB[j]*dataB[j] + dataB[j+1]*dataB[j+1] + eps; + double re = dataA[j]*dataB[j] + dataA[j+1]*dataB[j+1]; + double im = dataA[j+1]*dataB[j] - dataA[j]*dataB[j+1]; + dataC[j] = re / denom; + dataC[j+1] = im / denom; + } + else + for( j = j0; j < j1; j += 2 ) + { + double denom = dataB[j]*dataB[j] + dataB[j+1]*dataB[j+1] + eps; + double re = dataA[j]*dataB[j] - dataA[j+1]*dataB[j+1]; + double im = dataA[j+1]*dataB[j] + dataA[j]*dataB[j+1]; + dataC[j] = re / denom; + dataC[j+1] = im / denom; + } + } + } +} /****************************************************************************************\ Discrete Cosine Transform diff --git a/modules/core/test/test_dxt.cpp b/modules/core/test/test_dxt.cpp index cc382588b0..f5112f06cd 100644 --- a/modules/core/test/test_dxt.cpp +++ b/modules/core/test/test_dxt.cpp @@ -879,9 +879,283 @@ void CxCore_MulSpectrumsTest::prepare_to_validation( int /*test_case_idx*/ ) } } +////////////////////// DivSpectrums //////////////////////// +class CV_DivSpectrumsTest : public cvtest::ArrayTest +{ +public: + CV_DivSpectrumsTest(); +protected: + void run_func(); + void get_test_array_types_and_sizes( int, vector >& sizes, vector >& types ); + void prepare_to_validation( int test_case_idx ); + int flags; +}; + + +CV_DivSpectrumsTest::CV_DivSpectrumsTest() : flags(0) +{ + // Allocate test matrices. + test_array[INPUT].push_back(NULL); // first input DFT as a CCS-packed array or complex matrix. + test_array[INPUT].push_back(NULL); // second input DFT as a CCS-packed array or complex matrix. + test_array[OUTPUT].push_back(NULL); // output DFT as a complex matrix. + test_array[REF_OUTPUT].push_back(NULL); // reference output DFT as a complex matrix. + test_array[TEMP].push_back(NULL); // first input DFT converted to a complex matrix. + test_array[TEMP].push_back(NULL); // second input DFT converted to a complex matrix. + test_array[TEMP].push_back(NULL); // output DFT as a CCV-packed array. +} + +void CV_DivSpectrumsTest::get_test_array_types_and_sizes( int test_case_idx, vector >& sizes, vector >& types ) +{ + cvtest::ArrayTest::get_test_array_types_and_sizes(test_case_idx, sizes, types); + RNG& rng = ts->get_rng(); + + // Get the flag of the input. + const int rand_int_flags = cvtest::randInt(rng); + flags = rand_int_flags & (CV_TEST_DXT_MUL_CONJ | DFT_ROWS); + + // Get input type. + const int rand_int_type = cvtest::randInt(rng); + int type; + + if (rand_int_type % 4) + { + type = CV_32FC1; + } + else if (rand_int_type % 4 == 1) + { + type = CV_32FC2; + } + else if (rand_int_type % 4 == 2) + { + type = CV_64FC1; + } + else + { + type = CV_64FC2; + } + + for( size_t i = 0; i < types.size(); i++ ) + { + for( size_t j = 0; j < types[i].size(); j++ ) + { + types[i][j] = type; + } + } + + // Inputs are CCS-packed arrays. Prepare outputs and temporary inputs as complex matrices. + if( type == CV_32FC1 || type == CV_64FC1 ) + { + types[OUTPUT][0] += CV_DEPTH_MAX; + types[REF_OUTPUT][0] += CV_DEPTH_MAX; + types[TEMP][0] += CV_DEPTH_MAX; + types[TEMP][1] += CV_DEPTH_MAX; + } +} + +/// Helper function to convert a ccs array of depth_t into a complex matrix. +template +static void convert_from_ccs_helper( const Mat& src0, const Mat& src1, Mat& dst ) +{ + const int cn = src0.channels(); + int srcstep = cn; + int dststep = 1; + + if( !dst.isContinuous() ) + dststep = (int)(dst.step/dst.elemSize()); + + if( !src0.isContinuous() ) + srcstep = (int)(src0.step/src0.elemSize1()); + + Complex *dst_data = dst.ptr >(); + const depth_t* src0_data = src0.ptr(); + const depth_t* src1_data = src1.ptr(); + dst_data->re = src0_data[0]; + dst_data->im = 0; + const int n = dst.cols + dst.rows - 1; + const int n2 = (n+1) >> 1; + + if( (n & 1) == 0 ) + { + dst_data[n2*dststep].re = src0_data[(cn == 1 ? n-1 : n2)*srcstep]; + dst_data[n2*dststep].im = 0; + } + + int delta0 = srcstep; + int delta1 = delta0 + (cn == 1 ? srcstep : 1); + + if( cn == 1 ) + srcstep *= 2; + + for( int i = 1; i < n2; i++, delta0 += srcstep, delta1 += srcstep ) + { + depth_t t0 = src0_data[delta0]; + depth_t t1 = src0_data[delta1]; + + dst_data[i*dststep].re = t0; + dst_data[i*dststep].im = t1; + + t0 = src1_data[delta0]; + t1 = -src1_data[delta1]; + + dst_data[(n-i)*dststep].re = t0; + dst_data[(n-i)*dststep].im = t1; + } +} + +/// Helper function to convert a ccs array into a complex matrix. +static void convert_from_ccs( const Mat& src0, const Mat& src1, Mat& dst, const int flags ) +{ + if( dst.rows > 1 && (dst.cols > 1 || (flags & DFT_ROWS)) ) + { + const int count = dst.rows; + const int len = dst.cols; + const bool is2d = (flags & DFT_ROWS) == 0; + for( int i = 0; i < count; i++ ) + { + const int j = !is2d || i == 0 ? i : count - i; + const Mat& src0row = src0.row(i); + const Mat& src1row = src1.row(j); + Mat dstrow = dst.row(i); + convert_from_ccs( src0row, src1row, dstrow, 0 ); + } + + if( is2d ) + { + const Mat& src0row = src0.col(0); + Mat dstrow = dst.col(0); + convert_from_ccs( src0row, src0row, dstrow, 0 ); + + if( (len & 1) == 0 ) + { + const Mat& src0row_even = src0.col(src0.cols - 1); + Mat dstrow_even = dst.col(len/2); + convert_from_ccs( src0row_even, src0row_even, dstrow_even, 0 ); + } + } + } + else + { + if( dst.depth() == CV_32F ) + { + convert_from_ccs_helper( src0, src1, dst ); + } + else + { + convert_from_ccs_helper( src0, src1, dst ); + } + } +} + +/// Helper function to compute complex number (nu_re + nu_im * i) / (de_re + de_im * i). +static std::pair divide_complex_numbers( const double nu_re, const double nu_im, + const double de_re, const double de_im, + const bool conj_de ) +{ + if ( conj_de ) + { + return divide_complex_numbers( nu_re, nu_im, de_re, -de_im, false /* conj_de */ ); + } + + const double result_de = de_re * de_re + de_im * de_im + DBL_EPSILON; + const double result_re = nu_re * de_re + nu_im * de_im; + const double result_im = nu_re * (-de_im) + nu_im * de_re; + return std::pair(result_re / result_de, result_im / result_de); +} + +/// Helper function to divide a DFT in src1 by a DFT in src2 with depths depth_t. The DFTs are +/// complex matrices. +template +static void div_complex_helper( const Mat& src1, const Mat& src2, Mat& dst, int flags ) +{ + CV_Assert( src1.size == src2.size && src1.type() == src2.type() ); + dst.create( src1.rows, src1.cols, src1.type() ); + const int cn = src1.channels(); + int cols = src1.cols * cn; + + for( int i = 0; i < dst.rows; i++ ) + { + const depth_t *src1_data = src1.ptr(i); + const depth_t *src2_data = src2.ptr(i); + depth_t *dst_data = dst.ptr(i); + for( int j = 0; j < cols; j += 2 ) + { + std::pair result = + divide_complex_numbers( src1_data[j], src1_data[j + 1], + src2_data[j], src2_data[j + 1], + (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); + dst_data[j] = (depth_t)result.first; + dst_data[j + 1] = (depth_t)result.second; + } + } +} + +/// Helper function to divide a DFT in src1 by a DFT in src2. The DFTs are complex matrices. +static void div_complex( const Mat& src1, const Mat& src2, Mat& dst, const int flags ) +{ + const int type = src1.type(); + CV_Assert( type == CV_32FC2 || type == CV_64FC2 ); + + if ( src1.depth() == CV_32F ) + { + return div_complex_helper( src1, src2, dst, flags ); + } + else + { + return div_complex_helper( src1, src2, dst, flags ); + } +} + +void CV_DivSpectrumsTest::prepare_to_validation( int /* test_case_idx */ ) +{ + Mat &src1 = test_mat[INPUT][0]; + Mat &src2 = test_mat[INPUT][1]; + Mat &ref_dst = test_mat[REF_OUTPUT][0]; + const int cn = src1.channels(); + // Inputs are CCS-packed arrays. Convert them to complex matrices and get the expected output + // as a complex matrix. + if( cn == 1 ) + { + Mat &converted_src1 = test_mat[TEMP][0]; + Mat &converted_src2 = test_mat[TEMP][1]; + convert_from_ccs( src1, src1, converted_src1, flags ); + convert_from_ccs( src2, src2, converted_src2, flags ); + div_complex( converted_src1, converted_src2, ref_dst, flags ); + } + // Inputs are complex matrices. Get the expected output as a complex matrix. + else + { + div_complex( src1, src2, ref_dst, flags ); + } +} + +void CV_DivSpectrumsTest::run_func() +{ + const Mat &src1 = test_mat[INPUT][0]; + const Mat &src2 = test_mat[INPUT][1]; + const int cn = src1.channels(); + + // Inputs are CCS-packed arrays. Get the output as a CCS-packed array and convert it to a + // complex matrix. + if ( cn == 1 ) + { + Mat &dst = test_mat[TEMP][2]; + cv::divSpectrums( src1, src2, dst, flags, (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); + Mat &converted_dst = test_mat[OUTPUT][0]; + convert_from_ccs( dst, dst, converted_dst, flags ); + } + // Inputs are complex matrices. Get the output as a complex matrix. + else + { + Mat &dst = test_mat[OUTPUT][0]; + cv::divSpectrums( src1, src2, dst, flags, (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); + } +} + TEST(Core_DCT, accuracy) { CxCore_DCTTest test; test.safe_run(); } TEST(Core_DFT, accuracy) { CxCore_DFTTest test; test.safe_run(); } TEST(Core_MulSpectrums, accuracy) { CxCore_MulSpectrumsTest test; test.safe_run(); } +TEST(Core_DivSpectrums, accuracy) { CV_DivSpectrumsTest test; test.safe_run(); } + class Core_DFTComplexOutputTest : public cvtest::BaseTest { diff --git a/modules/imgproc/include/opencv2/imgproc.hpp b/modules/imgproc/include/opencv2/imgproc.hpp index 98e5b619fb..682da3c500 100644 --- a/modules/imgproc/include/opencv2/imgproc.hpp +++ b/modules/imgproc/include/opencv2/imgproc.hpp @@ -2951,22 +2951,6 @@ An example is shown below: */ CV_EXPORTS_W void createHanningWindow(OutputArray dst, Size winSize, int type); -/** @brief Performs the per-element division of the first Fourier spectrum by the second Fourier spectrum. - -The function cv::divSpectrums performs the per-element division of the first array by the second array. -The arrays are CCS-packed or complex matrices that are results of a real or complex Fourier transform. - -@param a first input array. -@param b second input array of the same size and type as src1 . -@param c output array of the same size and type as src1 . -@param flags operation flags; currently, the only supported flag is cv::DFT_ROWS, which indicates that -each row of src1 and src2 is an independent 1D Fourier spectrum. If you do not want to use this flag, then simply add a `0` as value. -@param conjB optional flag that conjugates the second input array before the multiplication (true) -or not (false). -*/ -CV_EXPORTS_W void divSpectrums(InputArray a, InputArray b, OutputArray c, - int flags, bool conjB = false); - //! @} imgproc_motion //! @addtogroup imgproc_misc diff --git a/modules/imgproc/misc/js/gen_dict.json b/modules/imgproc/misc/js/gen_dict.json index eeb876cd79..c8221a5d55 100644 --- a/modules/imgproc/misc/js/gen_dict.json +++ b/modules/imgproc/misc/js/gen_dict.json @@ -35,7 +35,6 @@ "dilate", "distanceTransform", "distanceTransformWithLabels", - "divSpectrums", "drawContours", "drawMarker", "ellipse", diff --git a/modules/imgproc/src/phasecorr.cpp b/modules/imgproc/src/phasecorr.cpp index 3fbe62883d..8d6c3cd85d 100644 --- a/modules/imgproc/src/phasecorr.cpp +++ b/modules/imgproc/src/phasecorr.cpp @@ -159,204 +159,6 @@ static void magSpectrums( InputArray _src, OutputArray _dst) } } -void divSpectrums( InputArray _srcA, InputArray _srcB, OutputArray _dst, int flags, bool conjB) -{ - Mat srcA = _srcA.getMat(), srcB = _srcB.getMat(); - int depth = srcA.depth(), cn = srcA.channels(), type = srcA.type(); - int rows = srcA.rows, cols = srcA.cols; - int j, k; - - CV_Assert( type == srcB.type() && srcA.size() == srcB.size() ); - CV_Assert( type == CV_32FC1 || type == CV_32FC2 || type == CV_64FC1 || type == CV_64FC2 ); - - _dst.create( srcA.rows, srcA.cols, type ); - Mat dst = _dst.getMat(); - - CV_Assert(dst.data != srcA.data); // non-inplace check - CV_Assert(dst.data != srcB.data); // non-inplace check - - bool is_1d = (flags & DFT_ROWS) || (rows == 1 || (cols == 1 && - srcA.isContinuous() && srcB.isContinuous() && dst.isContinuous())); - - if( is_1d && !(flags & DFT_ROWS) ) - cols = cols + rows - 1, rows = 1; - - int ncols = cols*cn; - int j0 = cn == 1; - int j1 = ncols - (cols % 2 == 0 && cn == 1); - - if( depth == CV_32F ) - { - const float* dataA = srcA.ptr(); - const float* dataB = srcB.ptr(); - float* dataC = dst.ptr(); - float eps = FLT_EPSILON; // prevent div0 problems - - size_t stepA = srcA.step/sizeof(dataA[0]); - size_t stepB = srcB.step/sizeof(dataB[0]); - size_t stepC = dst.step/sizeof(dataC[0]); - - if( !is_1d && cn == 1 ) - { - for( k = 0; k < (cols % 2 ? 1 : 2); k++ ) - { - if( k == 1 ) - dataA += cols - 1, dataB += cols - 1, dataC += cols - 1; - dataC[0] = dataA[0] / (dataB[0] + eps); - if( rows % 2 == 0 ) - dataC[(rows-1)*stepC] = dataA[(rows-1)*stepA] / (dataB[(rows-1)*stepB] + eps); - if( !conjB ) - for( j = 1; j <= rows - 2; j += 2 ) - { - double denom = (double)dataB[j*stepB]*dataB[j*stepB] + - (double)dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + (double)eps; - - double re = (double)dataA[j*stepA]*dataB[j*stepB] + - (double)dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; - - double im = (double)dataA[(j+1)*stepA]*dataB[j*stepB] - - (double)dataA[j*stepA]*dataB[(j+1)*stepB]; - - dataC[j*stepC] = (float)(re / denom); - dataC[(j+1)*stepC] = (float)(im / denom); - } - else - for( j = 1; j <= rows - 2; j += 2 ) - { - - double denom = (double)dataB[j*stepB]*dataB[j*stepB] + - (double)dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + (double)eps; - - double re = (double)dataA[j*stepA]*dataB[j*stepB] - - (double)dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; - - double im = (double)dataA[(j+1)*stepA]*dataB[j*stepB] + - (double)dataA[j*stepA]*dataB[(j+1)*stepB]; - - dataC[j*stepC] = (float)(re / denom); - dataC[(j+1)*stepC] = (float)(im / denom); - } - if( k == 1 ) - dataA -= cols - 1, dataB -= cols - 1, dataC -= cols - 1; - } - } - - for( ; rows--; dataA += stepA, dataB += stepB, dataC += stepC ) - { - if( is_1d && cn == 1 ) - { - dataC[0] = dataA[0] / (dataB[0] + eps); - if( cols % 2 == 0 ) - dataC[j1] = dataA[j1] / (dataB[j1] + eps); - } - - if( !conjB ) - for( j = j0; j < j1; j += 2 ) - { - double denom = (double)dataB[j]*dataB[j] + (double)dataB[j+1]*dataB[j+1] + (double)eps; - double re = (double)dataA[j]*dataB[j] + (double)dataA[j+1]*dataB[j+1]; - double im = (double)dataA[j+1]*dataB[j] - (double)dataA[j]*dataB[j+1]; - dataC[j] = (float)(re / denom); - dataC[j+1] = (float)(im / denom); - } - else - for( j = j0; j < j1; j += 2 ) - { - double denom = (double)dataB[j]*dataB[j] + (double)dataB[j+1]*dataB[j+1] + (double)eps; - double re = (double)dataA[j]*dataB[j] - (double)dataA[j+1]*dataB[j+1]; - double im = (double)dataA[j+1]*dataB[j] + (double)dataA[j]*dataB[j+1]; - dataC[j] = (float)(re / denom); - dataC[j+1] = (float)(im / denom); - } - } - } - else - { - const double* dataA = srcA.ptr(); - const double* dataB = srcB.ptr(); - double* dataC = dst.ptr(); - double eps = DBL_EPSILON; // prevent div0 problems - - size_t stepA = srcA.step/sizeof(dataA[0]); - size_t stepB = srcB.step/sizeof(dataB[0]); - size_t stepC = dst.step/sizeof(dataC[0]); - - if( !is_1d && cn == 1 ) - { - for( k = 0; k < (cols % 2 ? 1 : 2); k++ ) - { - if( k == 1 ) - dataA += cols - 1, dataB += cols - 1, dataC += cols - 1; - dataC[0] = dataA[0] / (dataB[0] + eps); - if( rows % 2 == 0 ) - dataC[(rows-1)*stepC] = dataA[(rows-1)*stepA] / (dataB[(rows-1)*stepB] + eps); - if( !conjB ) - for( j = 1; j <= rows - 2; j += 2 ) - { - double denom = dataB[j*stepB]*dataB[j*stepB] + - dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + eps; - - double re = dataA[j*stepA]*dataB[j*stepB] + - dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; - - double im = dataA[(j+1)*stepA]*dataB[j*stepB] - - dataA[j*stepA]*dataB[(j+1)*stepB]; - - dataC[j*stepC] = re / denom; - dataC[(j+1)*stepC] = im / denom; - } - else - for( j = 1; j <= rows - 2; j += 2 ) - { - double denom = dataB[j*stepB]*dataB[j*stepB] + - dataB[(j+1)*stepB]*dataB[(j+1)*stepB] + eps; - - double re = dataA[j*stepA]*dataB[j*stepB] - - dataA[(j+1)*stepA]*dataB[(j+1)*stepB]; - - double im = dataA[(j+1)*stepA]*dataB[j*stepB] + - dataA[j*stepA]*dataB[(j+1)*stepB]; - - dataC[j*stepC] = re / denom; - dataC[(j+1)*stepC] = im / denom; - } - if( k == 1 ) - dataA -= cols - 1, dataB -= cols - 1, dataC -= cols - 1; - } - } - - for( ; rows--; dataA += stepA, dataB += stepB, dataC += stepC ) - { - if( is_1d && cn == 1 ) - { - dataC[0] = dataA[0] / (dataB[0] + eps); - if( cols % 2 == 0 ) - dataC[j1] = dataA[j1] / (dataB[j1] + eps); - } - - if( !conjB ) - for( j = j0; j < j1; j += 2 ) - { - double denom = dataB[j]*dataB[j] + dataB[j+1]*dataB[j+1] + eps; - double re = dataA[j]*dataB[j] + dataA[j+1]*dataB[j+1]; - double im = dataA[j+1]*dataB[j] - dataA[j]*dataB[j+1]; - dataC[j] = re / denom; - dataC[j+1] = im / denom; - } - else - for( j = j0; j < j1; j += 2 ) - { - double denom = dataB[j]*dataB[j] + dataB[j+1]*dataB[j+1] + eps; - double re = dataA[j]*dataB[j] - dataA[j+1]*dataB[j+1]; - double im = dataA[j+1]*dataB[j] + dataA[j]*dataB[j+1]; - dataC[j] = re / denom; - dataC[j+1] = im / denom; - } - } - } - -} - static void fftShift(InputOutputArray _out) { Mat out = _out.getMat(); diff --git a/modules/imgproc/test/test_pc.cpp b/modules/imgproc/test/test_pc.cpp index 40e290322a..d5d093b08b 100644 --- a/modules/imgproc/test/test_pc.cpp +++ b/modules/imgproc/test/test_pc.cpp @@ -151,278 +151,4 @@ TEST(Imgproc_PhaseCorrelatorTest, float32_overflow) { EXPECT_NEAR(std::abs(phaseShift.y), 0.0, 1.0); } -////////////////////// DivSpectrums //////////////////////// -class CV_DivSpectrumsTest : public cvtest::ArrayTest -{ -public: - CV_DivSpectrumsTest(); -protected: - void run_func(); - void get_test_array_types_and_sizes( int, vector >& sizes, vector >& types ); - void prepare_to_validation( int test_case_idx ); - int flags; -}; - - -CV_DivSpectrumsTest::CV_DivSpectrumsTest() : flags(0) -{ - // Allocate test matrices. - test_array[INPUT].push_back(NULL); // first input DFT as a CCS-packed array or complex matrix. - test_array[INPUT].push_back(NULL); // second input DFT as a CCS-packed array or complex matrix. - test_array[OUTPUT].push_back(NULL); // output DFT as a complex matrix. - test_array[REF_OUTPUT].push_back(NULL); // reference output DFT as a complex matrix. - test_array[TEMP].push_back(NULL); // first input DFT converted to a complex matrix. - test_array[TEMP].push_back(NULL); // second input DFT converted to a complex matrix. - test_array[TEMP].push_back(NULL); // output DFT as a CCV-packed array. -} - -void CV_DivSpectrumsTest::get_test_array_types_and_sizes( int test_case_idx, vector >& sizes, vector >& types ) -{ - cvtest::ArrayTest::get_test_array_types_and_sizes(test_case_idx, sizes, types); - RNG& rng = ts->get_rng(); - - // Get the flag of the input. - const int rand_int_flags = cvtest::randInt(rng); - flags = rand_int_flags & (CV_TEST_DXT_MUL_CONJ | DFT_ROWS); - - // Get input type. - const int rand_int_type = cvtest::randInt(rng); - int type; - - if (rand_int_type % 4) - { - type = CV_32FC1; - } - else if (rand_int_type % 4 == 1) - { - type = CV_32FC2; - } - else if (rand_int_type % 4 == 2) - { - type = CV_64FC1; - } - else - { - type = CV_64FC2; - } - - for( size_t i = 0; i < types.size(); i++ ) - { - for( size_t j = 0; j < types[i].size(); j++ ) - { - types[i][j] = type; - } - } - - // Inputs are CCS-packed arrays. Prepare outputs and temporary inputs as complex matrices. - if( type == CV_32FC1 || type == CV_64FC1 ) - { - types[OUTPUT][0] += CV_DEPTH_MAX; - types[REF_OUTPUT][0] += CV_DEPTH_MAX; - types[TEMP][0] += CV_DEPTH_MAX; - types[TEMP][1] += CV_DEPTH_MAX; - } -} - -/// Helper function to convert a ccs array of depth_t into a complex matrix. -template -static void convert_from_ccs_helper( const Mat& src0, const Mat& src1, Mat& dst ) -{ - const int cn = src0.channels(); - int srcstep = cn; - int dststep = 1; - - if( !dst.isContinuous() ) - dststep = (int)(dst.step/dst.elemSize()); - - if( !src0.isContinuous() ) - srcstep = (int)(src0.step/src0.elemSize1()); - - Complex *dst_data = dst.ptr >(); - const depth_t* src0_data = src0.ptr(); - const depth_t* src1_data = src1.ptr(); - dst_data->re = src0_data[0]; - dst_data->im = 0; - const int n = dst.cols + dst.rows - 1; - const int n2 = (n+1) >> 1; - - if( (n & 1) == 0 ) - { - dst_data[n2*dststep].re = src0_data[(cn == 1 ? n-1 : n2)*srcstep]; - dst_data[n2*dststep].im = 0; - } - - int delta0 = srcstep; - int delta1 = delta0 + (cn == 1 ? srcstep : 1); - - if( cn == 1 ) - srcstep *= 2; - - for( int i = 1; i < n2; i++, delta0 += srcstep, delta1 += srcstep ) - { - depth_t t0 = src0_data[delta0]; - depth_t t1 = src0_data[delta1]; - - dst_data[i*dststep].re = t0; - dst_data[i*dststep].im = t1; - - t0 = src1_data[delta0]; - t1 = -src1_data[delta1]; - - dst_data[(n-i)*dststep].re = t0; - dst_data[(n-i)*dststep].im = t1; - } -} - -/// Helper function to convert a ccs array into a complex matrix. -static void convert_from_ccs( const Mat& src0, const Mat& src1, Mat& dst, const int flags ) -{ - if( dst.rows > 1 && (dst.cols > 1 || (flags & DFT_ROWS)) ) - { - const int count = dst.rows; - const int len = dst.cols; - const bool is2d = (flags & DFT_ROWS) == 0; - for( int i = 0; i < count; i++ ) - { - const int j = !is2d || i == 0 ? i : count - i; - const Mat& src0row = src0.row(i); - const Mat& src1row = src1.row(j); - Mat dstrow = dst.row(i); - convert_from_ccs( src0row, src1row, dstrow, 0 ); - } - - if( is2d ) - { - const Mat& src0row = src0.col(0); - Mat dstrow = dst.col(0); - convert_from_ccs( src0row, src0row, dstrow, 0 ); - - if( (len & 1) == 0 ) - { - const Mat& src0row_even = src0.col(src0.cols - 1); - Mat dstrow_even = dst.col(len/2); - convert_from_ccs( src0row_even, src0row_even, dstrow_even, 0 ); - } - } - } - else - { - if( dst.depth() == CV_32F ) - { - convert_from_ccs_helper( src0, src1, dst ); - } - else - { - convert_from_ccs_helper( src0, src1, dst ); - } - } -} - -/// Helper function to compute complex number (nu_re + nu_im * i) / (de_re + de_im * i). -static std::pair divide_complex_numbers( const double nu_re, const double nu_im, - const double de_re, const double de_im, - const bool conj_de ) -{ - if ( conj_de ) - { - return divide_complex_numbers( nu_re, nu_im, de_re, -de_im, false /* conj_de */ ); - } - - const double result_de = de_re * de_re + de_im * de_im + DBL_EPSILON; - const double result_re = nu_re * de_re + nu_im * de_im; - const double result_im = nu_re * (-de_im) + nu_im * de_re; - return std::pair(result_re / result_de, result_im / result_de); -} - -/// Helper function to divide a DFT in src1 by a DFT in src2 with depths depth_t. The DFTs are -/// complex matrices. -template -static void div_complex_helper( const Mat& src1, const Mat& src2, Mat& dst, int flags ) -{ - CV_Assert( src1.size == src2.size && src1.type() == src2.type() ); - dst.create( src1.rows, src1.cols, src1.type() ); - const int cn = src1.channels(); - int cols = src1.cols * cn; - - for( int i = 0; i < dst.rows; i++ ) - { - const depth_t *src1_data = src1.ptr(i); - const depth_t *src2_data = src2.ptr(i); - depth_t *dst_data = dst.ptr(i); - for( int j = 0; j < cols; j += 2 ) - { - std::pair result = - divide_complex_numbers( src1_data[j], src1_data[j + 1], - src2_data[j], src2_data[j + 1], - (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); - dst_data[j] = (depth_t)result.first; - dst_data[j + 1] = (depth_t)result.second; - } - } -} - -/// Helper function to divide a DFT in src1 by a DFT in src2. The DFTs are complex matrices. -static void div_complex( const Mat& src1, const Mat& src2, Mat& dst, const int flags ) -{ - const int type = src1.type(); - CV_Assert( type == CV_32FC2 || type == CV_64FC2 ); - - if ( src1.depth() == CV_32F ) - { - return div_complex_helper( src1, src2, dst, flags ); - } - else - { - return div_complex_helper( src1, src2, dst, flags ); - } -} - -void CV_DivSpectrumsTest::prepare_to_validation( int /* test_case_idx */ ) -{ - Mat &src1 = test_mat[INPUT][0]; - Mat &src2 = test_mat[INPUT][1]; - Mat &ref_dst = test_mat[REF_OUTPUT][0]; - const int cn = src1.channels(); - // Inputs are CCS-packed arrays. Convert them to complex matrices and get the expected output - // as a complex matrix. - if( cn == 1 ) - { - Mat &converted_src1 = test_mat[TEMP][0]; - Mat &converted_src2 = test_mat[TEMP][1]; - convert_from_ccs( src1, src1, converted_src1, flags ); - convert_from_ccs( src2, src2, converted_src2, flags ); - div_complex( converted_src1, converted_src2, ref_dst, flags ); - } - // Inputs are complex matrices. Get the expected output as a complex matrix. - else - { - div_complex( src1, src2, ref_dst, flags ); - } -} - -void CV_DivSpectrumsTest::run_func() -{ - const Mat &src1 = test_mat[INPUT][0]; - const Mat &src2 = test_mat[INPUT][1]; - const int cn = src1.channels(); - - // Inputs are CCS-packed arrays. Get the output as a CCS-packed array and convert it to a - // complex matrix. - if ( cn == 1 ) - { - Mat &dst = test_mat[TEMP][2]; - cv::divSpectrums( src1, src2, dst, flags, (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); - Mat &converted_dst = test_mat[OUTPUT][0]; - convert_from_ccs( dst, dst, converted_dst, flags ); - } - // Inputs are complex matrices. Get the output as a complex matrix. - else - { - Mat &dst = test_mat[OUTPUT][0]; - cv::divSpectrums( src1, src2, dst, flags, (flags & CV_TEST_DXT_MUL_CONJ) != 0 ); - } -} - -TEST(Imgproc_DivSpectrums, accuracy) { CV_DivSpectrumsTest test; test.safe_run(); } - }} // namespace diff --git a/platforms/js/opencv_js.config.py b/platforms/js/opencv_js.config.py index 5108d1ed5b..ced62a19a0 100644 --- a/platforms/js/opencv_js.config.py +++ b/platforms/js/opencv_js.config.py @@ -3,7 +3,7 @@ core = { '': [ 'absdiff', 'add', 'addWeighted', 'bitwise_and', 'bitwise_not', 'bitwise_or', 'bitwise_xor', 'cartToPolar', - 'compare', 'convertScaleAbs', 'copyMakeBorder', 'countNonZero', 'determinant', 'dft', 'divide', 'eigen', + 'compare', 'convertScaleAbs', 'copyMakeBorder', 'countNonZero', 'determinant', 'dft', 'divide', 'divSpectrums', 'eigen', 'exp', 'flip', 'getOptimalDFTSize','gemm', 'hconcat', 'inRange', 'invert', 'kmeans', 'log', 'magnitude', 'max', 'mean', 'meanStdDev', 'merge', 'min', 'minMaxLoc', 'mixChannels', 'multiply', 'norm', 'normalize', 'perspectiveTransform', 'polarToCart', 'pow', 'randn', 'randu', 'reduce', 'repeat', 'rotate', 'setIdentity', 'setRNGSeed', @@ -49,7 +49,6 @@ imgproc = { 'dilate', 'distanceTransform', 'distanceTransformWithLabels', - 'divSpectrums', 'drawContours', 'drawMarker', 'ellipse',