厳密な重心の求め方

移動体の追跡や粒子のカウントのチュートリアル中では輪郭構成点から近似的な中心座標を求める手法を紹介していましたが、今回は厳密な重心を求める手法を紹介します。まずは重心の定義から見ていきましょう。なお、同様のことをPythonでも記述可能です(こちら)。

幾何学的には、ある図形の、そのまわりでの一次モーメントが 0 であるような点のこと。図形 D （およびその周辺）の各点 r が密度 f(r) を持つなら、その重心 g とは、

$\int_D (\boldsymbol{g} - \boldsymbol{r})f(\boldsymbol{r})\,dV = \boldsymbol{0}$

を満たす点 g である（g が D 外の点であることもあり得る）。

密度が一定の場合、単体に限って言うなら、全頂点の各座標の値の算術平均をその座標の値として持つ点はその単体の重心となる。
[Wikipediaより]

moments関数について

ポリゴンまたはラスタライズされた形状の，3次までのモーメントを求める、moments関数を使用します。

関数momentsは，ベクタ形状またはラスタライズされた形状の、3次までのモーメントを求めます。

$\texttt{m} _{ji}= \sum _{x,y} \left ( \texttt{array} (x,y) \cdot x^j \cdot y^i \right ),$
また，中心モーメントは次のように求められます。

$\texttt{mu} _{ji}= \sum _{x,y} \left ( \texttt{array} (x,y) \cdot (x - \bar{x} )^j \cdot (y - \bar{y} )^i \right )$

ここで $$(\bar{x}, \bar{y})$$は重心を表します：

$\bar{x} = \frac{\texttt{m}_{10}}{\texttt{m}_{00}} , \; \bar{y} = \frac{\texttt{m}_{01}}{\texttt{m}_{00}}$ [OpenCV 2.2 documentation より]

では、実際に次の画像について重心を求めてみましょう。画像をもとに重心を求めるやり方と、輪郭をもとに重心を求めるやり方の2通りのコードを以下に示しています。画像をimread関数で読み込む部分は共通です。

いずれの場合にも次のような、重心の位置に円が描かれた画像が表示されます。

画像をもとに重心を求める

moments関数はグレースケールの画像をもとにモーメントを計算するのでimreadで画像を読み込むときに、IMREAD_GRAYSCALE オプションをつけます。moments(img, false)としていますが、この第2引数が真(true)のとき、すべての非0ピクセルが1として扱われ、均一な密度を持つ板の重心を求めることができます。逆に言えば、この第2引数が偽(false)のときには、密度の不均一な図形の重心を求めることができます。このとき、密度は画素値に対応し、白(255)ほど高密度、黒(0)ほど低密度となります。このプログラムは重心座標を次のように出力します。

x: 143.288574 y: 102.619743

#include <stdio.h>
#include "opencv/cv.h"
#include "opencv/highgui.h"
using namespace cv;

int main(){
    
    Mat img = imread("sample.jpg", IMREAD_GRAYSCALE);
    
    Moments mu = moments( img, false );
    Point2f mc = Point2f( mu.m10/mu.m00 , mu.m01/mu.m00 );
    circle( img, mc, 4, Scalar(100), 2, 4);
    
    printf("x: %f  y: %f", mc.x, mc.y);
    
    imshow("img",img);
    waitKey(0);
    return 0;
}

輪郭をもとに重心を求める

moments関数が輪郭をもとにモーメントを計算するときには、内部的にはグリーンの定理を用いて求めています。なお、findContours関数もグレースケールの画像をもとに輪郭を計算するのでimreadで画像を読み込むときに、IMREAD_GRAYSCALE オプションをつけます。下のコードでは最大の輪郭長を持つ輪郭について重心を求めています。輪郭をもとにモーメントを計算する場合にはmoments関数に第2引数を指定する必要はなく、指定しても無視されます。このプログラムは重心座標を次のように出力します。

x: 142.963867 y: 102.192329

#include <stdio.h>
#include "opencv/cv.h"
#include "opencv/highgui.h"
using namespace cv;

int main(){
    
    Mat img = imread("sample.jpg", IMREAD_GRAYSCALE);
    Mat img2 =img.clone();
    
    vector<vector<Point> > contours;
    findContours(img2, contours, CV_RETR_EXTERNAL, CV_CHAIN_APPROX_NONE);
    
    double max_size=0;
    int max_id=-1;
    
    for(int i=0; i<contours.size();i++){
        if(contours[i].size()>max_size){
            max_size=contours[i].size();
            max_id=i;
        }
    }
    
    Moments mu = moments( contours[max_id]);
    Point2f mc = Point2f( mu.m10/mu.m00 , mu.m01/mu.m00 );
    
    circle( img, mc, 4, Scalar(100), 2, 4);
    
    printf("x: %f  y: %f", mc.x, mc.y);
    
    imshow("img",img);
    waitKey(0);
    return 0;
}

実は、全く同じ画像を入力に用いても、画像から求めた重心と、輪郭から求めた重心は少しだけ異なります。これは、画像の解像度が有限であるためです。