Peter M. Groen / openbr

Commit dd3004354f3c861bd0d6cf2323a80497e2b7f62d

Authored by Brendan Klare
1 parent a6db81c8

Training of new procrustes method finished

Inline Side-by-side
Showing 1 changed file with 159 additions and 0 deletions
openbr/plugins/landmarks.cpp
  View file @dd30043
... ... @@ -5,9 +5,11 @@
5 5 #include "openbr/core/eigenutils.h"
6 6 #include <QString>
7 7 #include <Eigen/SVD>
  8 +#include <Eigen/Dense>
8 9  
9 10 using namespace std;
10 11 using namespace cv;
  12 +using namespace Eigen;
11 13  
12 14 namespace br
13 15 {
... ... @@ -138,6 +140,163 @@ class ProcrustesTransform : public Transform
138 140 BR_REGISTER(Transform, ProcrustesTransform)
139 141  
140 142 /*!
  143 + */
  144 +class ProcrustesAlignTransform : public Transform
  145 +{
  146 + Q_OBJECT
  147 +
  148 + Q_PROPERTY(float width READ get_width WRITE set_width RESET reset_width STORED false)
  149 + BR_PROPERTY(float, width, true)
  150 +
  151 + Eigen::MatrixXf referenceShape;
  152 +
  153 + MatrixXf getRotation(MatrixXf ref, MatrixXf sample) {
  154 + MatrixXf R = ref.transpose() * sample;
  155 + JacobiSVD<MatrixXf> svd(R, ComputeFullU | ComputeFullV);
  156 + R = svd.matrixU() * svd.matrixV();
  157 +
  158 + return R;
  159 + }
  160 +
  161 + //Converts x y points in a single vector to two column matrix
  162 + MatrixXf vectorToMatrix(MatrixXf vector) {
  163 + int n = vector.rows();
  164 + MatrixXf matrix(n / 2, 2);
  165 + for (int i = 0; i < n / 2; i++) {
  166 + for (int j = 0; j < 2; j++) {
  167 + matrix(i, j) = vector(i * 2 + j);
  168 + }
  169 + }
  170 + return matrix;
  171 + }
  172 +
  173 + MatrixXf matrixToVector(MatrixXf matrix) {
  174 + int n2 = matrix.rows();
  175 + MatrixXf vector(n2 * 2, 1);
  176 + for (int i = 0; i < n2; i++) {
  177 + for (int j = 0; j < 2; j++) {
  178 + vector(i * 2 + j) = matrix(i, j);
  179 + }
  180 + }
  181 + return vector;
  182 + }
  183 +
  184 + void train(const TemplateList &data)
  185 + {
  186 + MatrixXf points(data[0].file.points().size() * 2, data.size());
  187 +
  188 + // Normalize all sets of points
  189 + for (int j = 0; j < data.size(); j++) {
  190 + QList<QPointF> imagePoints = data[j].file.points();
  191 +
  192 + float meanX = 0,
  193 + meanY = 0;
  194 + for (int i = 0; i < imagePoints.size(); i++) {
  195 + points(i * 2, j) = imagePoints[i].x();
  196 + points(i * 2 + 1, j) = imagePoints[i].y();
  197 + meanX += imagePoints[i].x();
  198 + meanY += imagePoints[i].y();
  199 + }
  200 + meanX /= imagePoints.size();
  201 + meanY /= imagePoints.size();
  202 +
  203 + for (int i = 0; i < imagePoints.size(); i++) {
  204 + points(i * 2, j) -= meanX;
  205 + points(i * 2 + 1, j) -= meanY;
  206 + }
  207 + }
  208 +
  209 +
  210 + //normalize scale
  211 + for (int i = 0; i < points.cols(); i++)
  212 + points.col(i) = points.col(i) / points.col(i).norm();
  213 +
  214 + //Normalize rotation
  215 + MatrixXf refPrev;
  216 + referenceShape = vectorToMatrix(points.rowwise().sum() / points.cols());
  217 + float diff = FLT_MAX;
  218 +
  219 + for (int j = 0; j < points.cols(); j++) {
  220 + MatrixXf p = vectorToMatrix(points.col(j));
  221 + MatrixXf R = getRotation(referenceShape, p);
  222 + p = p * R.transpose();
  223 + points.col(j) = matrixToVector(p);
  224 + }
  225 +
  226 + referenceShape = vectorToMatrix(points.rowwise().sum() / points.cols());
  227 + }
  228 +
  229 + void procustesAlgin(const Template &src, Template &dst, MatrixXf p) {
  230 +
  231 + QList<QPointF> imagePoints = src.file.points();
  232 + MatrixXf p(imagePoints.size() * 2, 1);
  233 + for (int i = 0; i < imagePoints.size(); i++) {
  234 +
  235 + }
  236 + }
  237 +
  238 + void project(const Template &src, Template &dst) const
  239 + {
  240 + MatrixXf p(src.file.points().size() * 2, 1);
  241 +
  242 +
  243 + // Normalize all sets of points
  244 + for (int j = 0; j < data.size(); j++) {
  245 + QList<QPointF> imagePoints = data[j].file.points();
  246 +
  247 + float meanX = 0,
  248 + meanY = 0;
  249 + for (int i = 0; i < imagepoints.size(); i++) {
  250 + points(i * 2, j) = imagepoints[i].x();
  251 + points(i * 2 + 1, j) = imagepoints[i].y();
  252 + meanX += imagePoints[i].x();
  253 + meanY += imagePoints[i].y();
  254 + }
  255 + meanX /= imagePoints.size();
  256 + meanY /= imagePoints.size();
  257 +
  258 + for (int i = 0; i < imagePoints.size(); i++) {
  259 + points(i * 2, j) -= meanX;
  260 + points(i * 2 + 1, j) -= meanY;
  261 + }
  262 + }
  263 +
  264 +
  265 + //normalize scale
  266 + for (int i = 0; i < points.cols(); i++)
  267 + points.col(i) = points.col(i) / points.col(i).norm();
  268 +
  269 + //Normalize rotation
  270 + MatrixXf refPrev;
  271 + referenceShape = vectorToMatrix(points.rowwise().sum() / points.cols());
  272 + float diff = FLT_MAX;
  273 +
  274 + for (int j = 0; j < points.cols(); j++) {
  275 + MatrixXf p = vectorToMatrix(points.col(j));
  276 + MatrixXf R = getRotation(referenceShape, p);
  277 + p = p * R.transpose();
  278 + points.col(j) = matrixToVector(p);
  279 + }
  280 +
  281 + referenceShape = vectorToMatrix(points.rowwise().sum() / points.cols());
  282 +
  283 + }
  284 +
  285 + void store(QDataStream &stream) const
  286 + {
  287 + stream << referenceShape;
  288 + }
  289 +
  290 + void load(QDataStream &stream)
  291 + {
  292 + stream >> referenceShape;
  293 + }
  294 +
  295 +};
  296 +
  297 +BR_REGISTER(Transform, ProcrustesAlignTransform)
  298 +
  299 +/*!
141 300 * \ingroup transforms
142 301 * \brief Creates a Delaunay triangulation based on a set of points
143 302 * \author Scott Klum \cite sklum
... ...