矩的计算

double M[7] = {0};        //HU不变矩
bool HuMoment(Mat &image)
{
	int bmpWidth = image.cols;
	int bmpHeight = image.rows;
	int bmpStep = image.step;
	int bmpChannels = image.channels();

	uchar* pBmpBuf = (uchar*)image.data;

	double m00 = 0, m11 = 0, m20 = 0, m02 = 0, m30 = 0, m03 = 0, m12 = 0, m21 = 0; //中心矩
	double x0 = 0, y0 = 0; //计算中心矩时所使用的临时变量
	double u20 = 0, u02 = 0, u11 = 0, u30 = 0, u03 = 0, u12 = 0, u21 = 0; //规划化后的中心矩
	double t1 = 0, t2 = 0, t3 = 0, t4 =  0, t5 = 0; //临时变量


	//获得图像的区域的重心(普通矩)
	double s10 = 0, s01 = 0, s00 = 0; //0阶矩和1阶矩
	for (int j = 0; j <bmpHeight; j++) //y
	{
		for (int i = 0; i <bmpWidth; i++) //x
		{
			s10 += i * pBmpBuf[ j * bmpStep + i];
			s01 += j * pBmpBuf[ j * bmpStep + i];
			s00 += pBmpBuf[ j * bmpStep + i];
		}
	}
	int center_x = (int)( s10 / s00 + 0.5);
	int center_y = (int)( s01 / s00 + 0.5);

   //计算二阶 三阶矩(重心矩)
    m00 = s00;
   for (int j = 0; j <bmpHeight; j++) 
   {
	   for (int i = 0; i <bmpWidth; i++) 
	   {
		   x0 = i - center_x;
		   y0 = j - center_y;
		   
		   m11 += x0 * y0 * pBmpBuf[ j * bmpStep + i];
		   m20 += x0 * x0 * pBmpBuf[ j * bmpStep + i];
		   m02 += y0 * y0 * pBmpBuf[ j * bmpStep + i];
		   m03 += y0 * y0 * y0 *pBmpBuf[ j * bmpStep + i];
		   m30 += x0 * x0 * x0 * pBmpBuf[ j * bmpStep + i];
		   m12 += x0 * y0 * y0 * pBmpBuf[ j * bmpStep + i];
		   m21 += x0 * x0 * y0 * pBmpBuf[ j * bmpStep + i];
	   }
   }

   //计算规范化后的中心矩
   u20 = m20 / pow( m00, 2);
   u02 = m02 / pow( m00, 2);
   u11 = m11 / pow( m00, 2);
   u30 = m30 / pow( m00, 2.5);
   u03 = m03 / pow( m00, 2.5);
   u12 = m12 / pow( m00, 2.5);
   u21 = m21 / pow( m00, 2.5);

   //计算中间变量
   t1 = ( u20 - u02);
   t2 = ( u30 - 3 * u12);
   t3 = ( 3 * u21 - u03);
   t4 = ( u30 + u12);
   t5 = ( u21 + u03);

   //计算不变矩
   M[0] = u20 + u02;
   M[1] = t1 * t1 + 4 * u11 * u11; 
   M[2] = t2 * t2 + t3 * t3; 
   M[3] = t4 * t4 + t5 * t5;
   M[4] = t2 * t4 * (t4 * t4 - 3 * t5 * t5) + t3 * t5 * ( 3 * t4 * t4 - t5 * t5); 
   M[5] = t1 * ( t4 * t4 - t5 * t5) + 4 * u11 * t4 * t5;
   M[6] = t3 * t4 * ( t4 * t4 - 3 * t5 * t5) - t2 * t5 * ( 3 * t4 * t4 - t5 * t5);


   return true;

}

 

 

opencv API:

int main(int argc, char** argv) 
{ 
    Mat image = imread(argv[1]);  
    cvtColor(image, image, CV_BGR2GRAY); 
    Moments mts = moments(image); 
    double hu[7]; 
    HuMoments(mts, hu); 
    for (int i=0; i<7; i++) 
    { 
        cout << log(abs(hu[i])) <<endl; 
    } 
   return 0; 
}

opencv里对Hu矩的计算有直接的API，它分为了两个函数：moments()函数用于计算中心矩，HuMoments函数用于由中心矩计算Hu矩。

Moments moments(InputArray array, bool binaryImage=false )

参数说明
输入参数：array是一幅单通道，8-bits的图像，或一个二维浮点数组(Point of Point2f)。binaryImage用来指示输出图像是否为一幅二值图像，如果是二值图像，则图像中所有非0像素看作为1进行计算。
输出参数：moments是一个类：

class Moments
{
public:
    Moments();
    Moments(double m00, double m10, double m01, double m20, double m11,
            double m02, double m30, double m21, double m12, double m03 );
    Moments( const CvMoments& moments );
    operator CvMoments() const;
}

里面保存了图像的2阶与3阶中心矩的值。

void HuMoments(const Moments& moments, double* hu)
参数说明：
输入参数：moments即为上面一个函数计算得到的moments类型。
输出参数：hu是一个含有7个数的数组。