Eigen+PCA投影处理 opencv绘图展示

#include <iostream>
#include <opencv2/opencv.hpp>
#include <opencv2/core/core.hpp>
#include <opencv2/imgproc/imgproc.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <eigen3/Eigen/Core>
#include <eigen3/Eigen/Dense>

Eigen::MatrixXf extractPixel(cv::Mat img,int k,int num)
{
    int rows=img.rows;
    int cols=img.cols;
    Eigen::MatrixXf res(num,3);
    int c=0;
    for (int i=0;i<rows;i++) {
        for (int j=0;j<cols;j++) {
            int m=img.at<int>(i,j);
            if(m==uchar(k))
            {
                res(c,0)=i;res(c,1)=j;res(c,2)=0;c++;
            }
        }
    }
    res.resize(c,3);
    return res;
}
Eigen::RowVectorXf featurenormalize(Eigen::MatrixXf &X)
{
    //计算每一维度均值
    Eigen::MatrixXf meanval = X.colwise().mean();//每一列的矩阵，列降维
    Eigen::RowVectorXf meanvecRow = meanval;
    //样本均值化为0
    X.rowwise() -= meanvecRow;
    return meanvecRow;
}

void computeCov(Eigen::MatrixXf &X, Eigen::MatrixXf &C)
{
    //计算协方差矩阵C = XTX / n-1;
    C = X.adjoint()*X ;
    C = C.array() / X.rows() - 1;
}
void computeEig(Eigen::MatrixXf &C, Eigen::MatrixXf &vec, Eigen::MatrixXf &val)
{
    //计算特征值和特征向量，使用selfadjont按照对阵矩阵的算法去计算，可以让产生的vec和val按照有序排列
    Eigen::SelfAdjointEigenSolver<Eigen::MatrixXf> eig(C);
    vec = eig.eigenvectors();
    val = eig.eigenvalues();
}
cv::Mat ColinearMat(cv::Mat input,int label)
{
    cv::Mat pcaMat=cv::Mat::zeros(300,500,CV_8UC1);
    cv::Mat labels, stats, centroids;
    int nccomps = connectedComponentsWithStats(input, labels, stats, centroids); //可以得到连通域的个数
    for (int comp = 1; comp < nccomps; ++comp)
    {
        //erase small block
        if(stats.at<int>(comp, 4)<50)
            continue;
        //extract compenet pixel
        Eigen::MatrixXf obj = extractPixel(labels,comp,stats.at<int>(comp, 4));
        //zero-mean
        Eigen::MatrixXf obj_mean=obj;
        Eigen::RowVectorXf meanvecRow=featurenormalize(obj_mean);
        //coveriance
        Eigen::MatrixXf cov(3,3);
        computeCov(obj_mean, cov);
        //calculate feature vec & val
        Eigen::MatrixXf featureVec,featureVal;
        computeEig(cov, featureVec, featureVal);
        //projectPCA
        Eigen::Vector3f direction=featureVec.col(2).transpose();
        direction.normalize();
        Eigen::MatrixXf point_proj=obj_mean*(direction*direction.transpose());
        point_proj.rowwise() += meanvecRow;
        //draw mat
        for (int i=0;i<point_proj.rows();i++)
        {
            int x=int(point_proj(i,0));
            int y=int(point_proj(i,1));
            pcaMat.at<uchar>(x,y)=label;
        }
    }
    return pcaMat;
}
int main()
{
    cv::Mat srcImg=cv::imread("/home/freja/component2.png",0);
    int label=2;
     cv::Mat pcaMat = ColinearMat(srcImg,label);
     cv::imshow("img", srcImg);
     cv::imshow("pcaMat", pcaMat);
    cv::waitKey(0);
//    */
    std::cout << "Hello World!" << std::endl;
    return 0;
}