c++版 mulSpectrums

但是没有测试效果：
//----------------------------
 #define VAL(buf, elem) (((T*)((char*)data ## buf + (step ## buf * (elem))))[0])
 #define MUL_SPECTRUMS_COL(A, B, C) \
 VAL(C, 0) = VAL(A, 0) * VAL(B, 0); \
 for (size_t j = 1; j <= rows - 2; j += 2) \
 { \
 double a_re = VAL(A, j), a_im = VAL(A, j + 1); \
 double b_re = VAL(B, j), b_im = VAL(B, j + 1); \
 if (conjB) b_im = -b_im; \
 double c_re = a_re * b_re - a_im * b_im; \
 double c_im = a_re * b_im + a_im * b_re; \
 VAL(C, j) = (T)c_re; VAL(C, j + 1) = (T)c_im; \
 } \
 if ((rows & 1) == 0) \
 VAL(C, rows-1) = VAL(A, rows-1) * VAL(B, rows-1)
 template <typename T, bool conjB> static inline
 void mulSpectrums_processCol_noinplace(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows)
 {
 MUL_SPECTRUMS_COL(A, B, C);
 }
 template <typename T, bool conjB> static inline
 void mulSpectrums_processCol_inplaceA(const T* dataB, T* dataAC, size_t stepB, size_t stepAC, size_t rows)
 {
 MUL_SPECTRUMS_COL(AC, B, AC);
 }
 template <typename T, bool conjB, bool inplaceA> static inline
 void mulSpectrums_processCol(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows)
 {
 if (inplaceA)
 mulSpectrums_processCol_inplaceA<T, conjB>(dataB, dataC, stepB, stepC, rows);
 else
 mulSpectrums_processCol_noinplace<T, conjB>(dataA, dataB, dataC, stepA, stepB, stepC, rows);
 }
 #undef MUL_SPECTRUMS_COL
 #undef VAL
 #define VAL(buf, elem) (data ## buf[(elem)])
 #define MUL_SPECTRUMS_ROW(A, B, C) \
 for (size_t j = j0; j < j1; j += 2) \
 { \
 double a_re = VAL(A, j), a_im = VAL(A, j + 1); \
 double b_re = VAL(B, j), b_im = VAL(B, j + 1); \
 if (conjB) b_im = -b_im; \
 double c_re = a_re * b_re - a_im * b_im; \
 double c_im = a_re * b_im + a_im * b_re; \
 VAL(C, j) = (T)c_re; VAL(C, j + 1) = (T)c_im; \
 }
 template <typename T, bool conjB> static inline
 void mulSpectrums_processRow_noinplace(const T* dataA, const T* dataB, T* dataC, size_t j0, size_t j1)
 {
 MUL_SPECTRUMS_ROW(A, B, C);
 }
 template <typename T, bool conjB> static inline
 void mulSpectrums_processRow_inplaceA(const T* dataB, T* dataAC, size_t j0, size_t j1)
 {
 MUL_SPECTRUMS_ROW(AC, B, AC);
 }
 template <typename T, bool conjB, bool inplaceA> static inline
 void mulSpectrums_processRow(const T* dataA, const T* dataB, T* dataC, size_t j0, size_t j1)
 {
 if (inplaceA)
 mulSpectrums_processRow_inplaceA<T, conjB>(dataB, dataC, j0, j1);
 else
 mulSpectrums_processRow_noinplace<T, conjB>(dataA, dataB, dataC, j0, j1);
 }
 #undef MUL_SPECTRUMS_ROW
 #undef VAL
 template <typename T, bool conjB, bool inplaceA> static inline
 void mulSpectrums_processCols(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows, size_t cols)
 {
 mulSpectrums_processCol<T, conjB, inplaceA>(dataA, dataB, dataC, stepA, stepB, stepC, rows);
 if ((cols & 1) == 0)
 {
 mulSpectrums_processCol<T, conjB, inplaceA>(dataA + cols - 1, dataB + cols - 1, dataC + cols - 1, stepA, stepB, stepC, rows);
 }
 }
 template <typename T, bool conjB, bool inplaceA> static inline
 void mulSpectrums_processRows(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows, size_t cols, size_t j0, size_t j1, bool is_1d_CN1)
 {
 while (rows-- > 0)
 {
 if (is_1d_CN1)
 dataC[0] = dataA[0]*dataB[0];
 mulSpectrums_processRow<T, conjB, inplaceA>(dataA, dataB, dataC, j0, j1);
 if (is_1d_CN1 && (cols & 1) == 0)
 dataC[j1] = dataA[j1]*dataB[j1];
 dataA = (const T*)(((char*)dataA) + stepA);
 dataB = (const T*)(((char*)dataB) + stepB);
 dataC = (T*)(((char*)dataC) + stepC);
 }
 }
 template <typename T, bool conjB, bool inplaceA> static inline
 void mulSpectrums_Impl_(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows, size_t cols, size_t j0, size_t j1, bool is_1d, bool isCN1)
 {
 if (!is_1d && isCN1)
 {
 mulSpectrums_processCols<T, conjB, inplaceA>(dataA, dataB, dataC, stepA, stepB, stepC, rows, cols);
 }
 mulSpectrums_processRows<T, conjB, inplaceA>(dataA, dataB, dataC, stepA, stepB, stepC, rows, cols, j0, j1, is_1d && isCN1);
 }
 template <typename T, bool conjB> static inline
 void mulSpectrums_Impl(const T* dataA, const T* dataB, T* dataC, size_t stepA, size_t stepB, size_t stepC, size_t rows, size_t cols, size_t j0, size_t j1, bool is_1d, bool isCN1)
 {
 if (dataA == dataC)
 mulSpectrums_Impl_<T, conjB, true>(dataA, dataB, dataC, stepA, stepB, stepC, rows, cols, j0, j1, is_1d, isCN1);
 else
 mulSpectrums_Impl_<T, conjB, false>(dataA, dataB, dataC, stepA, stepB, stepC, rows, cols, j0, j1, is_1d, isCN1);
 }
 void mulSpectrums( InputArray _srcA, InputArray _srcB,
 OutputArray _dst, int flags, bool conjB )
 {
 // CV_INSTRUMENT_REGION()
 // CV_OCL_RUN(_dst.isUMat() && _srcA.dims() <= 2 && _srcB.dims() <= 2,
 // ocl_mulSpectrums(_srcA, _srcB, _dst, flags, conjB))
 Mat srcA = _srcA.getMat(), srcB = _srcB.getMat();
 int depth = srcA.depth(), cn = srcA.channels(), type = srcA.type();
 size_t rows = srcA.rows, cols = srcA.cols;
 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();
 // correct inplace support
 // Case 'dst.data == srcA.data' is handled by implementation,
 // because it is used frequently (filter2D, matchTemplate)
 if (dst.data == srcB.data)
 srcB = srcB.clone(); // workaround for B only
 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;
 bool isCN1 = cn == 1;
 size_t j0 = isCN1 ? 1 : 0;
 size_t j1 = cols*cn - (((cols & 1) == 0 && cn == 1) ? 1 : 0);
 if (depth == CV_32F)
 {
 const float* dataA = srcA.ptr<float>();
 const float* dataB = srcB.ptr<float>();
 float* dataC = dst.ptr<float>();
 if (!conjB)
 mulSpectrums_Impl<float, false>(dataA, dataB, dataC, srcA.step, srcB.step, dst.step, rows, cols, j0, j1, is_1d, isCN1);
 else
 mulSpectrums_Impl<float, true>(dataA, dataB, dataC, srcA.step, srcB.step, dst.step, rows, cols, j0, j1, is_1d, isCN1);
 }
 else
 {
 const double* dataA = srcA.ptr<double>();
 const double* dataB = srcB.ptr<double>();
 double* dataC = dst.ptr<double>();
 if (!conjB)
 mulSpectrums_Impl<double, false>(dataA, dataB, dataC, srcA.step, srcB.step, dst.step, rows, cols, j0, j1, is_1d, isCN1);
 else
 mulSpectrums_Impl<double, true>(dataA, dataB, dataC, srcA.step, srcB.step, dst.step, rows, cols, j0, j1, is_1d, isCN1);
 }
 }
 //-------------------------
int main()
{
 bool HOG = true;
 bool FIXEDWINDOW = false;
 bool MULTISCALE = true;
 int w=10;
 int h=10;
 fftwf_complex * r   =
 (fftwf_complex*) fftwf_malloc(sizeof (fftwf_complex) *w*h );
 for (size_t i = 0; i < h; i++) {
 for (size_t j = 0; j < w; j++) {
 r[i*w+j][0]= i*j*0.1+1;
 r[i*w+j][1]= i*j*0.2+1;
 }
 }
 fftwf_complex * r2   =
 (fftwf_complex*) fftwf_malloc(sizeof (fftwf_complex) *w*h );
 for (size_t i = 0; i < h; i++) {
 for (size_t j = 0; j < w; j++) {
 r2[i*w+j][0]= i*j*0.5+2;
 r2[i*w+j][1]= i*j*0.9+2;
 }
 }
 SimpleTracker  tracker(HOG, FIXEDWINDOW, MULTISCALE );
fftwf_complex * res=tracker.mulSpectrums(r,r2,w,h,0);
for (size_t i = 0; i < h; i++) {
 for (size_t j = 0; j < w; j++) {
 printf("%f+%fi ", r2[i*w+j][0],r2[i*w+j][1]);
 }
 printf("\n");
}

文章来源: blog.csdn.net，作者：网奇，版权归原作者所有，如需转载，请联系作者。

原文链接：blog.csdn.net/jacke121/article/details/61195152