Merge pull request #321 from biometrics/eigenutil_cleanup

Eigenutil cleanup

Merge pull request #321 from biometrics/eigenutil_cleanup
Eigenutil cleanup
Scott Klum
2 parents 3122cac3 82822bb2
Showing 4 changed files with 70 additions and 154 deletions
openbr/core/eigenutils.cpp
openbr/core/eigenutils.h
openbr/plugins/eigen3.cpp
openbr/plugins/stasm4.cpp
@@ -4,80 +4,15 @@
 using namespace Eigen;
 using namespace cv;
-//Helper function to quickly write eigen matrix to disk. Not efficient.
-void writeEigen(MatrixXf X, QString filename) {
-    Mat m(X.rows(),X.cols(),CV_32FC1);
-    for (int i = 0; i < X.rows(); i++) {
-        for (int j = 0; j < X.cols(); j++) {
-            m.at<float>(i,j) = X(i,j);
-        }
-    }
-    QScopedPointer<br::Format> format(br::Factory<br::Format>::make(filename));
-    format->write(br::Template(m));
-}
-
-void writeEigen(MatrixXd X, QString filename) {
-    Mat m(X.rows(),X.cols(),CV_32FC1);
-    for (int i = 0; i < X.rows(); i++) {
-        for (int j = 0; j < X.cols(); j++) {
-            m.at<float>(i,j) = (float)X(i,j);
-        }
-    }
-    QScopedPointer<br::Format> format(br::Factory<br::Format>::make(filename));
-    format->write(br::Template(m));
-}
-
-void writeEigen(VectorXd X, QString filename) {
-    Mat m(X.size(),1,CV_32FC1);
-    for (int i = 0; i < X.rows(); i++) {
-        m.at<float>(i,0) = (float)X(i);
-    }
-    QScopedPointer<br::Format> format(br::Factory<br::Format>::make(filename));
-    format->write(br::Template(m));
-}
-
-void writeEigen(VectorXf X, QString filename) {
-    Mat m(X.size(),1,CV_32FC1);
-    for (int i = 0; i < X.rows(); i++) {
-        m.at<float>(i,0) = X(i);
-    }
-    QScopedPointer<br::Format> format(br::Factory<br::Format>::make(filename));
-    format->write(br::Template(m));
-}
-
-void printEigen(Eigen::MatrixXd X) {
-    for (int i = 0; i < X.rows(); i++) {
-        QString str;
-        for (int j = 0; j < X.cols(); j++) {
-            str.append(QString::number(X(i,j)) + " ");
-        }
-        qDebug() << str;
-    }
-}
-void printEigen(Eigen::MatrixXf X) {
-    for (int i = 0; i < X.rows(); i++) {
-        QString str;
-        for (int j = 0; j < X.cols(); j++) {
-            str.append(QString::number(X(i,j)) + " ");
-        }
-        qDebug() << str;
-    }
-}
-
-void printSize(Eigen::MatrixXf X) {
+void EigenUtils::printSize(Eigen::MatrixXf X) {
     qDebug() << "Rows=" << X.rows() << "\tCols=" << X.cols();
 }
-float eigMean(const Eigen::MatrixXf& x) {
-    return x.array().sum() / (x.rows() * x.cols());
+float EigenUtils::stddev(const Eigen::MatrixXf& x) {
+    return sqrt((x.array() - x.mean()).pow(2).sum() / (x.cols() * x.rows()));
 }
-float eigStd(const Eigen::MatrixXf& x) {
-    float mean = eigMean(x);
-    return sqrt((x.array() - mean).pow(2).sum() / (x.cols() * x.rows()));
-}
-
-MatrixXf removeRowCol(const MatrixXf X, int row, int col) {
+MatrixXf EigenUtils::removeRowCol(const MatrixXf X, int row, int col) {
     MatrixXf Y(X.rows() - 1,X.cols() - 1);
     for (int i1 = 0, i2 = 0; i1 < X.rows(); i1++) {
@@ -96,7 +31,7 @@ MatrixXf removeRowCol(const MatrixXf X, int row, int col) {
     return Y;
 }
-MatrixXf pointsToMatrix(const QList<QPointF> points, bool isAffine) {
+MatrixXf EigenUtils::pointsToMatrix(const QList<QPointF> points, bool isAffine) {
     MatrixXf P(points.size(), isAffine ? 3 : 2);
     for (int i = 0; i < points.size(); i++) {
         P(i, 0) = points[i].x();
@@ -107,7 +42,7 @@ MatrixXf pointsToMatrix(const QList&lt;QPointF&gt; points, bool isAffine) {
     return P;
 }
-QList<QPointF> matrixToPoints(const Eigen::MatrixXf P) {
+QList<QPointF> EigenUtils::matrixToPoints(const Eigen::MatrixXf P) {
     QList<QPointF> points;
     for (int i = 0; i < P.rows(); i++)
         points.append(QPointF(P(i, 0), P(i, 1)));
@@ -115,7 +50,7 @@ QList&lt;QPointF&gt; matrixToPoints(const Eigen::MatrixXf P) {
 }
 //Converts x y points in a single vector to two column matrix
-Eigen::MatrixXf vectorToMatrix(const Eigen::MatrixXf vector) {
+Eigen::MatrixXf EigenUtils::vectorToMatrix(const Eigen::MatrixXf vector) {
     int n = vector.rows();
     Eigen::MatrixXf matrix(n / 2, 2);
     for (int i = 0; i < n / 2; i++) {
@@ -126,7 +61,7 @@ Eigen::MatrixXf vectorToMatrix(const Eigen::MatrixXf vector) {
     return matrix;
 }
-Eigen::MatrixXf matrixToVector(const Eigen::MatrixXf matrix) {
+Eigen::MatrixXf EigenUtils::matrixToVector(const Eigen::MatrixXf matrix) {
     int n2 = matrix.rows();
     Eigen::MatrixXf vector(n2 * 2, 1);
     for (int i = 0; i < n2; i++) {
@@ -136,12 +71,3 @@ Eigen::MatrixXf matrixToVector(const Eigen::MatrixXf matrix) {
     }
     return vector;
 }
-
-Eigen::MatrixXf toEigen(const Mat m) {
-    if (m.type() != CV_32F)
-        qFatal("Mat to Eigen Converstation only supports CV_32F");
-
-    Eigen::MatrixXf data(m.rows, m.cols);
-    return Eigen::Map<const Eigen::MatrixXf>(m.ptr<float>(), m.rows, m.cols);
-}
-
@@ -18,32 +18,74 @@
 #define EIGENUTILS_H
 #include <QDataStream>
+#include <QDebug>
 #include <Eigen/Core>
 #include <assert.h>
+#include <opencv2/core/eigen.hpp>
 #include <opencv2/core/core.hpp>
-void writeEigen(Eigen::MatrixXf X, QString filename);
-void writeEigen(Eigen::MatrixXd X, QString filename);
-void writeEigen(Eigen::VectorXd X, QString filename);
-void writeEigen(Eigen::VectorXf X, QString filename);
-void printEigen(Eigen::MatrixXd X);
-void printEigen(Eigen::MatrixXf X);
-void printSize(Eigen::MatrixXf X);
+#include "openbr/core/qtutils.h"
-//Converts x y points in a single vector to two column matrix
-Eigen::MatrixXf vectorToMatrix(const Eigen::MatrixXf vector);
-Eigen::MatrixXf matrixToVector(const Eigen::MatrixXf matrix);
+namespace EigenUtils
+{
+    template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
+    QString matrixToString(const Eigen::Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &mat)
+    {
+        QString result;
+        if (mat.rows() > 1) result += "{ ";
+        for (int r=0; r<mat.rows(); r++) {
+            if ((mat.rows() > 1) && (r > 0)) result += "  ";
+            if (mat.cols() > 1) result += "[";
+            for (int c=0; c<mat.cols(); c++) {
+                result += QString::number(mat(r, c));
+                if (c < mat.cols() - 1) result += ", ";
+            }
+            if (mat.cols() > 1) result += "]";
+            if (r < mat.rows()-1) result += "\n";
+        }
+        if (mat.rows() > 1) result += " }";
+        return result;
+    }
+
+    template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
+    void writeMatrix(const Eigen::Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &mat, const QString &filename)
+    {
+        int r = mat.rows();
+        int c = mat.cols();
+
+        _Scalar *data = new _Scalar[r*c];
+        for (int i=0; i<r; i++)
+            for (int j=0; j<c; j++)
+                data[i*c+j] = mat(i, j);
+        int bytes = r*c*sizeof(_Scalar);
+        QByteArray byteArray((const char*)data,bytes);
+        QtUtils::writeFile(filename,byteArray);
+    }
+
+    void printSize(Eigen::MatrixXf X);
-//Remove row and column from the matrix:
-Eigen::MatrixXf removeRowCol(const Eigen::MatrixXf X, int row, int col);
+    // Converts x y points in a single vector to two column matrix
+    Eigen::MatrixXf vectorToMatrix(const Eigen::MatrixXf vector);
+    Eigen::MatrixXf matrixToVector(const Eigen::MatrixXf matrix);
-//Convert a point list into a matrix:
-Eigen::MatrixXf pointsToMatrix(const QList<QPointF> points, bool isAffine=false);
-QList<QPointF> matrixToPoints(const Eigen::MatrixXf P);
+    // Remove row and column from the matrix:
+    Eigen::MatrixXf removeRowCol(const Eigen::MatrixXf X, int row, int col);
-//Convert cv::Mat to Eigen
-Eigen::MatrixXf toEigen(const cv::Mat m);
+    // Convert a point list into a matrix:
+    Eigen::MatrixXf pointsToMatrix(const QList<QPointF> points, bool isAffine=false);
+    QList<QPointF> matrixToPoints(const Eigen::MatrixXf P);
+
+    // Compute the element-wise standard deviation
+    float stddev(const Eigen::MatrixXf& x);
+}
+
+template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
+inline QDebug operator<<(QDebug dbg, const Eigen::Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &mat)
+{
+    dbg.nospace() << EigenUtils::matrixToString(mat);
+    return dbg.space();
+}
 template<typename _Scalar, int _Rows, int _Cols, int _Options, int _MaxRows, int _MaxCols>
 inline QDataStream &operator<<(QDataStream &stream, const Eigen::Matrix< _Scalar, _Rows, _Cols, _Options, _MaxRows, _MaxCols > &mat)
@@ -83,56 +125,4 @@ inline QDataStream &amp;operator&gt;&gt;(QDataStream &amp;stream, Eigen::Matrix&lt; _Scalar, _Row
     return stream;
 }
-/*Compute the mean of the each column (dim == 1) or row (dim == 2)
-  of the matrix*/
-template<typename T>
-Eigen::MatrixBase<T> eigMean(const Eigen::MatrixBase<T>& x,int dim)
-{
-    if (dim == 1) {
-        Eigen::MatrixBase<T> y(1,x.cols());
-        for (int i = 0; i < x.cols(); i++)
-            y(i) = x.col(i).sum() / x.rows();
-        return y;
-    } else if (dim == 2) {
-        Eigen::MatrixBase<T> y(x.rows(),1);
-        for (int i = 0; i < x.rows(); i++)
-            y(i) = x.row(i).sum() / x.cols();
-        return y;
-    }
-    qFatal("A matrix can only have two dimensions");
-}
-
-/*Compute the element-wise mean*/
-float eigMean(const Eigen::MatrixXf& x);
-/*Compute the element-wise mean*/
-float eigStd(const Eigen::MatrixXf& x);
-
-/*Compute the std dev of the each column (dim == 1) or row (dim == 2)
-  of the matrix*/
-template<typename T>
-Eigen::MatrixBase<T> eigStd(const Eigen::MatrixBase<T>& x,int dim)
-{
-    Eigen::MatrixBase<T> mean = eigMean(x, dim);
-    if (dim == 1) {
-        Eigen::MatrixBase<T> y(1,x.cols());
-        for (int i = 0; i < x.cols(); i++) {
-            T value = 0;
-            for (int j = 0; j < x.rows(); j++)
-                value += pow(y(j, i) - mean(i), 2);
-            y(i) = sqrt(value / (x.rows() - 1));
-        }
-        return y;
-    } else if (dim == 2) {
-        Eigen::MatrixBase<T> y(x.rows(),1);
-        for (int i = 0; i < x.rows(); i++) {
-            T value = 0;
-            for (int j = 0; j < x.cols(); j++)
-                value += pow(y(i, j) - mean(j), 2);
-            y(i) = sqrt(value / (x.cols() - 1));
-        }
-        return y;
-    }
-    qFatal("A matrix can only have two dimensions");
-}
-
 #endif // EIGENUTILS_H
@@ -600,7 +600,7 @@ class SparseLDATransform : public Transform
         //Only works on binary class problems for now
         assert(ldaOrig.projection.cols() == 1);
-        float ldaStd = eigStd(ldaOrig.projection);
+        float ldaStd = EigenUtils::stddev(ldaOrig.projection);
         for (int i = 0; i < ldaOrig.projection.rows(); i++)
             if (abs(ldaOrig.projection(i)) > varThreshold * ldaStd)
                 selections.append(i);
@@ -180,7 +180,7 @@ private:
     void project(const Template &src, Template &dst) const
     {
         QList<float> paramList = src.file.getList<float>("affineParameters");
-        Eigen::MatrixXf points = pointsToMatrix(src.file.points(), true);
+        Eigen::MatrixXf points = EigenUtils::pointsToMatrix(src.file.points(), true);
         Eigen::MatrixXf affine = Eigen::MatrixXf::Zero(3, 3);
         for (int i = 0, cnt = 0; i < 2; i++)
             for (int j = 0; j < 3; j++, cnt++)
@@ -192,7 +192,7 @@ private:
         points =  affineInv * pointsT;
         dst = src;
         dst.file.clearPoints();
-        dst.file.setPoints(matrixToPoints(points.transpose()));
+        dst.file.setPoints(EigenUtils::matrixToPoints(points.transpose()));
     }
 };