有关笛卡尔坐标于极坐标的转换函数void cartToPolar()

[king of code porter]

有关笛卡尔坐标于极坐标的转换函数void cartToPolar(const Mat& x, const Mat& y,
Mat& magnitude, Mat& angle,
bool angleInDegrees=false)
最近在学习OpenCV的时候遇到了这个函数，网上查了下，笛卡尔坐标转极坐标的公式是th = atan2(y,x);
r = hypot(x,y);r是半径，th是角度（弧度）。按照这个函数的意思根据这两个公式自己写了个函数来进行转换，发现效果不是很好。看了OpenCV里面的这个函数的源代码，好像这个函数调用了CPU的一些指令，没接触过，不知道OpenCV计算的时候有没有用到上面的两个函数？求高人指点，还有源代码中的那些语句不是很明白，不知道这个函数的核心函数是哪块地方，求高人指点和解释！附上这个函数的源代码

void cartToPolar( const Mat& X, const Mat& Y, Mat& Mag, Mat& Angle, bool angleInDegrees )
{
    float buf[2][MAX_BLOCK_SIZE];
    int i, j, type = X.type(), depth = X.depth(), cn = X.channels();

    CV_Assert( X.size() == Y.size() && type == Y.type() && (depth == CV_32F || depth == CV_64F));
    Mag.create( X.size(), type );
    Angle.create( X.size(), type );

    Size size = getContinuousSize( X, Y, Mag, Angle, cn );
    bool inplace = Mag.data == X.data || Mag.data == Y.data;

    if( depth == CV_32F )
    {
        const float *x = (const float*)X.data, *y = (const float*)Y.data;
        float *mag = (float*)Mag.data, *angle = (float*)Angle.data;
        size_t xstep = X.step/sizeof(x[0]), ystep = Y.step/sizeof(y[0]);
        size_t mstep = Mag.step/sizeof(mag[0]), astep = Angle.step/sizeof(angle[0]);

        for( ; size.height--; x += xstep, y += ystep, mag += mstep, angle += astep )
        {
            for( i = 0; i < size.width; i += MAX_BLOCK_SIZE )
            {
                int block_size = std::min(MAX_BLOCK_SIZE, size.width - i);
                Magnitude_32f( x + i, y + i, inplace ? buf[0] : mag + i, block_size );
                FastAtan2_32f( y + i, x + i, angle + i, block_size, angleInDegrees );
                if( inplace )
                    for( j = 0; j < block_size; j++ )
                        mag[i + j] = buf[0][j];
            }
        }
    }
    else
    {
        const double *x = (const double*)X.data, *y = (const double*)Y.data;
        double *mag = (double*)Mag.data, *angle = (double*)Angle.data;
        size_t xstep = X.step/sizeof(x[0]), ystep = Y.step/sizeof(y[0]);
        size_t mstep = Mag.step/sizeof(mag[0]), astep = Angle.step/sizeof(angle[0]);

        for( ; size.height--; x += xstep, y += ystep, mag += mstep, angle += astep )
        {
            for( i = 0; i < size.width; i += MAX_BLOCK_SIZE )
            {
                int block_size = std::min(MAX_BLOCK_SIZE, size.width - i);
                for( j = 0; j < block_size; j++ )
                {
                    buf[0][j] = (float)x[i + j];
                    buf[1][j] = (float)y[i + j];
                }
                FastAtan2_32f( buf[1], buf[0], buf[0], block_size, angleInDegrees );
                Magnitude_64f( x + i, y + i, mag + i, block_size );
				for( j = 0; j < block_size; j++ )
					angle[i + j] = buf[0][j];
            }
        }
    }
}

Magnitude_64f(const double* x, const double* y, double* mag, int len)
{
    int i = 0;
    
#if CV_SSE2   
    if( checkHardwareSupport(CV_CPU_SSE2) )
    {
        for( ; i <= len - 4; i += 4 )
        {
            __m128d x0 = _mm_loadu_pd(x + i), x1 = _mm_loadu_pd(x + i + 2);
            __m128d y0 = _mm_loadu_pd(y + i), y1 = _mm_loadu_pd(y + i + 2);
            x0 = _mm_add_pd(_mm_mul_pd(x0, x0), _mm_mul_pd(y0, y0));
            x1 = _mm_add_pd(_mm_mul_pd(x1, x1), _mm_mul_pd(y1, y1));
            x0 = _mm_sqrt_pd(x0); x1 = _mm_sqrt_pd(x1);
            _mm_storeu_pd(mag + i, x0); _mm_storeu_pd(mag + i + 2, x1);
        }
    }
#endif
    
    for( ; i < len; i++ )
    {
        double x0 = x[i], y0 = y[i];
        mag[i] = std::sqrt(x0*x0 + y0*y0);
    }
    
    return CV_OK;
}

static CvStatus CV_STDCALL FastAtan2_32f(const float *Y, const float *X, float *angle, int len, bool angleInDegrees=true )
{
    if( !Y || !X || !angle || len < 0 )
        return CV_BADFACTOR_ERR;
	
	int i = 0;
	float scale = angleInDegrees ? (float)(180/CV_PI) : 1.f;

#if CV_SSE2
    if( checkHardwareSupport(CV_CPU_SSE) )
    {
        Cv32suf iabsmask; iabsmask.i = 0x7fffffff;
        __m128 eps = _mm_set1_ps((float)DBL_EPSILON), absmask = _mm_set1_ps(iabsmask.f);
        __m128 _90 = _mm_set1_ps((float)(CV_PI*0.5)), _180 = _mm_set1_ps((float)CV_PI), _360 = _mm_set1_ps((float)(CV_PI*2));
        __m128 zero = _mm_setzero_ps(), _0_28 = _mm_set1_ps(0.28f), scale4 = _mm_set1_ps(scale);
        
        for( ; i <= len - 4; i += 4 )
        {
            __m128 x4 = _mm_loadu_ps(X + i), y4 = _mm_loadu_ps(Y + i);
            __m128 xq4 = _mm_mul_ps(x4, x4), yq4 = _mm_mul_ps(y4, y4);
            __m128 xly = _mm_cmplt_ps(xq4, yq4);
            __m128 z4 = _mm_div_ps(_mm_mul_ps(x4, y4), _mm_add_ps(_mm_add_ps(_mm_max_ps(xq4, yq4),
                _mm_mul_ps(_mm_min_ps(xq4, yq4), _0_28)), eps));
            
            // a4 <- x < y ? 90 : 0;
            __m128 a4 = _mm_and_ps(xly, _90);
            // a4 <- (y < 0 ? 360 - a4 : a4) == ((x < y ? y < 0 ? 270 : 90) : (y < 0 ? 360 : 0))
            __m128 mask = _mm_cmplt_ps(y4, zero);
            a4 = _mm_or_ps(_mm_and_ps(_mm_sub_ps(_360, a4), mask), _mm_andnot_ps(mask, a4));
            // a4 <- (x < 0 && !(x < y) ? 180 : a4)
            mask = _mm_andnot_ps(xly, _mm_cmplt_ps(x4, zero));
            a4 = _mm_or_ps(_mm_and_ps(_180, mask), _mm_andnot_ps(mask, a4));
            
            // a4 <- (x < y ? a4 - z4 : a4 + z4)
            a4 = _mm_mul_ps(_mm_add_ps(_mm_xor_ps(z4, _mm_andnot_ps(absmask, xly)), a4), scale4);
            _mm_storeu_ps(angle + i, a4);
        }
    }
#endif
	
	for( ; i < len; i++ )
	{
		float x = X[i], y = Y[i];
		float a, x2 = x*x, y2 = y*y;
		if( y2 <= x2 )
			a = x*y/(x2 + 0.28f*y2 + (float)DBL_EPSILON) + (float)(x < 0 ? CV_PI : y >= 0 ? 0 : CV_PI*2);
		else
			a = (float)(y >= 0 ? CV_PI*0.5 : CV_PI*1.5) - x*y/(y2 + 0.28f*x2 + (float)DBL_EPSILON);
		angle[i] = a*scale;
	}

    return CV_OK;
}

在线求解

[king of code porter]

没人回复，自己先顶一个

[应澜lst]

楼主知道这个函数怎么用了吗？求指导