表示和描述（3）

接上篇 表示和描述（2）

主分量描述

主分量描述适用于边界和区域。用在区域（图像）上可以抽取方差最大的分量（主分量），用在边界上可以对其做缩放、平移和旋转的归一化。

主分量（PCA）一般用于数据降维，因为大特征值对应图像细节（高频）。

主分量的计算过程可归纳为（具体分析参考）：由多组特征向量计算均值向量m，由此得到协方差矩阵，计算该矩阵（实对称）的特征值，按从大到小排列，选择前n个特征值对应的特征向量，按行组成变换矩阵A，这样，对原来的每一条向量x，其映射得到的向量y=A(x-m)。

1）使用主分量描述图像

考虑一系列M*N大小的6波段遥感图，将其数据组织为(M*N)*6的二维矩阵，每一行即为一个特征向量。

// 准备数据
    int i, j, k;
    Mat src[6];
    src[0] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(a)(WashingtonDC_Band1_512).tif", 0);
    src[1] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(b)(WashingtonDC_Band2_512).tif", 0);
    src[2] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(c)(WashingtonDC_Band3_512).tif", 0);
    src[3] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(d)(WashingtonDC_Band4_512).tif", 0);
    src[4] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(e)(WashingtonDC_Band5_512).tif", 0);
    src[5] = imread("../DIP_CODE+IMAGE/IMAGE/dipum_images_ch11/Fig1125(f)(WashingtonDC_Band6_512).tif", 0);

    Mat srcData(src[0].rows*src[0].cols, 6, CV_32FC1);
    for (i = 0; i < src[0].rows; i++)
    {
        uchar* ptr[6];
        for (k = 0; k < 6; k++)
            ptr[k] = src[k].ptr<uchar>(i);
        
        for (j = 0; j < src[0].cols; j++)
        {
            float* ptrD = srcData.ptr<float>(i*src[0].cols + j);
            for (k = 0; k < 6; k++)
                ptrD[k] = ptr[k][j] / 255.0;
        }
    }

    // PCA计算
    // 取前4组特征值
    PCA pca(srcData, Mat(), CV_PCA_DATA_AS_ROW, 4);

    cout << "eigenvalues:" << endl;
    for (i = 0; i < pca.eigenvalues.rows; i++)
        cout <<"            "<< pca.eigenvalues.at<float>(i, 0) << endl;

    cout << "eigenvectors:" << endl;
    for (i = 0; i < pca.eigenvectors.rows; i++)
    {
        float* ptr = pca.eigenvectors.ptr<float>(i);
        for (j = 0; j < pca.eigenvectors.cols; j++)
            cout << ptr[j] << "  ";
        cout << endl;
    }

    //
    for (i = 0; i < srcData.rows; i++)
    {
        Mat proj = pca.project(srcData.row(i));
        Mat bProj = pca.backProject(proj);

        float* ptrD = srcData.ptr<float>(i);
       
        for (j = 0; j < 6; j++)
            ptrD[j] = bProj.at<float>(0, j);
    }
    
    Mat dst[6];
    for (i = 0; i < 6; i++)
        dst[i] = Mat(src[0].size(), CV_8UC1);
    for (i = 0; i < src[0].rows; i++)
    {
        uchar* ptr[6];
        for (k = 0; k < 6; k++)
            ptr[k] = dst[k].ptr<uchar>(i);
        for (j = 0; j < src[0].cols; j++)
        {
            float* ptrD = srcData.ptr<float>(i*src[0].cols + j);
            for (k = 0; k < 6; k++)
                ptr[k][j] = saturate_cast<uchar>(ptrD[k] * 255);
        }
    }

    namedWindow("band1_");
    namedWindow("band2_");
    namedWindow("band3_");
    namedWindow("band4_");
    namedWindow("band5_");
    namedWindow("band6_");
    imshow("band1_", dst[0]);
    imshow("band2_", dst[1]);
    imshow("band3_", dst[2]);
    imshow("band4_", dst[3]);
    imshow("band5_", dst[4]);
    imshow("band6_", dst[5]);

    waitKey(0);

图3-1为取前4个特征值恢复得到的第6波段图像，可以看到，原始图比较模糊，但恢复的图像很清晰，因为模糊的部分（对应小特征值）已被丢弃

图3-1. 左为原图，右为PCA后的图

2）使用主分量归一化边界

假设共N个点，则组织数据为N*2

对所有特征向量，应用y=A(x-m)即获得归一化描述

Mat srcData = (Mat_<float>(4, 2) << 1, 1, 2, 4, 4, 2, 5, 5);

    PCA pca(srcData, Mat(), CV_PCA_DATA_AS_ROW);

    int i, j;
    Mat normali(srcData.size(), srcData.type());
    for (i = 0; i < srcData.rows; i++)
    {
        Mat nor = normali.row(i);
        Mat src = srcData.row(i);
        nor = (pca.eigenvectors * ((src - pca.mean)).t()).t();
    }

    cout << "normalization:" << endl;
    for (i = 0; i < normali.rows; i++)
    {
        float* ptr = normali.ptr<float>(i);
        for (j = 0; j < normali.cols; j++)
            cout << ptr[j] << "    ";
        cout << endl;
    }

    //
    for (i = 0; i < srcData.rows; i++)
    {
        Mat nor = normali.row(i);
        Mat src = srcData.row(i);
        pca.project(src, nor);
    }

    cout << "_normalization:" << endl;
    for (i = 0; i < normali.rows; i++)
    {
        float* ptr = normali.ptr<float>(i);
        for (j = 0; j < normali.cols; j++)
            cout << ptr[j] << "    ";
        cout << endl;
    }

从代码效果来看，OpenCV提供的pca.project就是实现的y=A(x-m)

关系描述

关系描述旨在用一套模式来描述边界或区域。

关系描述首先定义一些基本描述元，如用a表示向右一个像素，b表示向下一个像素；然后定义一套规则，如(1) S->aA (2)A->bS (3) A->b。则规则组合(1)(2)(1)(3)表示的边界为（图4-1）

图4-1. 关系描述（《数字图像处理》）

另一种是用树的形式描述每个边界或区域的从属关系。