使用OpenCV编写的LDA程序----C++ LDA代码

改写自OpenCV中的lda.cpp程序，通过改写的程序可以返回自己所需的信息（LDA算法过程中产生的我们感兴趣的中间值），实现算法的独立编译，也可以通过阅读程序，加深对LDA算法的理解。

// main.cpp : 定义控制台应用程序的入口点。
//

#include "stdafx.h"
#include <cxcore.hpp>
#include <vector>
#include <iostream>
#include "lda.h"
using namespace std;
using namespace cv;

int main(void)
{
	double data[6][2]={
  {0,1},{0,2},{1,4},{8,0},{8,2},{9,4}};
	Mat dmat=Mat(6,2,CV_64FC1,data);
	int labels[6]={0, 0, 0, 1, 1, 1};
	Mat lmat=Mat(1,6,CV_32SC1,labels);

	cout<<"--------------------------------"<<endl;

	MyLDA(dmat, lmat); 

	system("pause");
	return 0;
}

//lda.h

#ifndef  _MY_LDA_H
#define _MY_LDA_H

#include <cxcore.hpp>
#include <vector>
using namespace std;
using namespace cv;

/*********************
 _src:
 输入的采样数据，Mat类型，每行是一个样本数据
 _lbls:
 输入的类别标签，可以是矩阵或者向量

 还没有返回值，需要的朋友可以自己做进一步改写
/********************/
extern void MyLDA(InputArrayOfArrays _src, InputArray _lbls);

#endif 

//lda.cpp

#include "stdafx.h"
#include <cxcore.hpp>
#include <vector>
#include <set>
#include <iostream>
using namespace cv;
using namespace std;

// Removes duplicate elements in a given vector.
template<typename _Tp>
inline vector<_Tp> remove_dups(const vector<_Tp>& src) {
    typedef typename set<_Tp>::const_iterator constSetIterator;
    typedef typename vector<_Tp>::const_iterator constVecIterator;
    set<_Tp> set_elems;
    for (constVecIterator it = src.begin(); it != src.end(); ++it)
        set_elems.insert(*it);
    vector<_Tp> elems;
    for (constSetIterator it = set_elems.begin(); it != set_elems.end(); ++it)
        elems.push_back(*it);
    return elems;
}

static Mat argsort(InputArray _src, bool ascending=true)
{
    Mat src = _src.getMat();
    if (src.rows != 1 && src.cols != 1) {
        string error_message = "Wrong shape of input matrix! Expected a matrix with one row or column.";
        CV_Error(CV_StsBadArg, error_message);
    }
    int flags = CV_SORT_EVERY_ROW+(ascending ? CV_SORT_ASCENDING : CV_SORT_DESCENDING);
    Mat sorted_indices;
    sortIdx(src.reshape(1,1),sorted_indices,flags);
    return sorted_indices;
}

static Mat asRowMatrix(InputArrayOfArrays src, int rtype, double alpha=1, double beta=0) {
    // make sure the input data is a vector of matrices or vector of vector
    if(src.kind() != _InputArray::STD_VECTOR_MAT && src.kind() != _InputArray::STD_VECTOR_VECTOR) {
        string error_message = "The data is expected as InputArray::STD_VECTOR_MAT (a std::vector<Mat>) or _InputArray::STD_VECTOR_VECTOR (a std::vector< vector<...> >).";
        CV_Error(CV_StsBadArg, error_message);
    }
    // number of samples
    size_t n = src.total();
    // return empty matrix if no matrices given
    if(n == 0)
        return Mat();
    // dimensionality of (reshaped) samples
    size_t d = src.getMat(0).total();
    // create data matrix
    Mat data((int)n, (int)d, rtype);
    // now copy data
    for(int i = 0; i < (int)n; i++) {
        // make sure data can be reshaped, throw exception if not!
        if(src.getMat(i).total() != d) {
            string error_message = format("Wrong number of elements in matrix #%d! Expected %d was %d.", i, (int)d, (int)src.getMat(i).total());
            CV_Error(CV_StsBadArg, error_message);
        }
        // get a hold of the current row
        Mat xi = data.row(i);
        // make reshape happy by cloning for non-continuous matrices
        if(src.getMat(i).isContinuous()) {
            src.getMat(i).reshape(1, 1).convertTo(xi, rtype, alpha, beta);
        } else {
            src.getMat(i).clone().reshape(1, 1).convertTo(xi, rtype, alpha, beta);
        }
    }
    return data;
}

static void sortMatrixColumnsByIndices(InputArray _src, InputArray _indices, OutputArray _dst) {
    if(_indices.getMat().type() != CV_32SC1) {
        CV_Error(CV_StsUnsupportedFormat, "cv::sortColumnsByIndices only works on integer indices!");
    }
    Mat src = _src.getMat();
    vector<int> indices = _indices.getMat();
    _dst.create(src.rows, src.cols, src.type());
    Mat dst = _dst.getMat();
    for(size_t idx = 0; idx < indices.size(); idx++) {
        Mat originalCol = src.col(indices[idx]);
        Mat sortedCol = dst.col((int)idx);
        originalCol.copyTo(sortedCol);
    }
}

static Mat sortMatrixColumnsByIndices(InputArray src, InputArray indices) {
    Mat dst;
    sortMatrixColumnsByIndices(src, indices, dst);
    return dst;
}


template<typename _Tp> static bool
isSymmetric_(InputArray src) {
    Mat _src = src.getMat();
    if(_src.cols != _src.rows)
        return false;
    for (int i = 0; i < _src.rows; i++) {
        for (int j = 0; j < _src.cols; j++) {
            _Tp a = _src.at<_Tp> (i, j);
            _Tp b = _src.at<_Tp> (j, i);
            if (a != b) {
                return false;
            }
        }
    }
    return true;
}

template<typename _Tp> static bool
isSymmetric_(InputArray src, double eps) {
    Mat _src = src.getMat();
    if(_src.cols != _src.rows)
        return false;
    for (int i = 0; i < _src.rows; i++) {
        for (int j = 0; j < _src.cols; j++) {
            _Tp a = _src.at<_Tp> (i, j);
            _Tp b = _src.at<_Tp> (j, i);
            if (std::abs(a - b) > eps) {
                return false;
            }
        }
    }
    return true;
}

static bool isSymmetric(InputArray src, double eps=1e-16)
{
    Mat m = src.getMat();
    switch (m.type()) {
        case CV_8SC1: return isSymmetric_<char>(m); break;
        case CV_8UC1:
            return isSymmetric_<unsigned char>(m); break;
        case CV_16SC1:
            return isSymmetric_<short>(m); break;
        case CV_16UC1:
            return isSymmetric_<unsigned short>(m); break;
        case CV_32SC1:
            return isSymmetric_<int>(m); break;
        case CV_32FC1:
            return isSymmetric_<float>(m, eps); break;
        case CV_64FC1:
            return isSymmetric_<double>(m, eps); break;
        default:
            break;
    }
    return false;
}


//------------------------------------------------------------------------------
// cv::subspaceProject
//------------------------------------------------------------------------------
Mat subspaceProject(InputArray _W, InputArray _mean, InputArray _src) {
    // get data matrices
    Mat W = _W.getMat();
    Mat mean = _mean.getMat();
    Mat src = _src.getMat();
    // get number of samples and dimension
    int n = src.rows;
    int d = src.cols;
    // make sure the data has the correct shape
    if(W.rows != d) {
        string error_message = format("Wrong shapes for given matrices. Was size(src) = (%d,%d), size(W) = (%d,%d).", src.rows, src.cols, W.rows, W.cols);
        CV_Error(CV_StsBadArg, error_message);
    }
    // make sure mean is correct if not empty
    if(!mean.empty() && (mean.total() != (size_t) d)) {
        string error_message = format("Wrong mean shape for the given data matrix. Expected %d, but was %d.", d, mean.total());
        CV_Error(CV_StsBadArg, error_message);
    }
    // create temporary matrices
    Mat X, Y;
    // make sure you operate on correct type
    src.convertTo(X, W.type());
    // safe to do, because of above assertion
    if(!mean.empty()) {
        for(int i