Starts binocular polyfit

Uses data from both eyes for a single polynomial regression

EyeRecToo/EyeRecToo.pro

@@ -62,7 +62,8 @@ SOURCES +=\
     $${TOP}/src/globals.cpp \
     $${TOP}/src/AudioRecorder.cpp \
     $${TOP}/src/data/common.cpp \
-    $${TOP}/src/gaze-estimation/GazeEstimationMethod.cpp
+    $${TOP}/src/gaze-estimation/GazeEstimationMethod.cpp \
+    $${TOP}/src/gaze-estimation/BinocularPolyFit.cpp
 
 HEADERS  += \
     $${TOP}/src/MainWindow.h\
@@ -104,7 +105,8 @@ HEADERS  += \
     $${TOP}/src/data/FieldData.h \
     $${TOP}/src/data/InputData.h \
     $${TOP}/src/data/EyeData.h \
-    $${TOP}/src/data/Marker.h
+    $${TOP}/src/data/Marker.h \
+    $${TOP}/src/gaze-estimation/BinocularPolyFit.h
 
 FORMS    += \
     $${TOP}/src/MainWindow.ui \

EyeRecToo/src/GazeEstimation.cpp

@@ -13,12 +13,21 @@ GazeEstimation::GazeEstimation(QObject* parent)
     , lastOverlayIdx(0)
     , settings(nullptr)
 {
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY));
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY));
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY_XX_YY));
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY_XX_YY_XXYY));
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY_XX_YY_XYY_YXX));
-    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY));
+
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY));
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY));
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY_XX_YY));
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XXYY));
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX));
+    availableGazeEstimationMethods.push_back(make_shared<PolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY));
+
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY));
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY));
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY_XX_YY));
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XXYY));
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX));
+    availableGazeEstimationMethods.push_back(make_shared<BinocularPolyFit>(PolyFit::POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY));
+
     availableGazeEstimationMethods.push_back(make_shared<Homography>());
 }
 GazeEstimation::~GazeEstimation()

EyeRecToo/src/GazeEstimation.h

@@ -10,6 +10,7 @@
 
 #include "Reference.h"
 #include "data/DataTuple.h"
+#include "gaze-estimation/BinocularPolyFit.h"
 #include "gaze-estimation/GazeEstimationMethod.h"
 #include "gaze-estimation/Homography.h"
 #include "gaze-estimation/PolyFit.h"

EyeRecToo/src/gaze-estimation/BinocularPolyFit.cpp

+#include "BinocularPolyFit.h"
+
+using namespace std;
+using namespace cv;
+
+BinocularPolyFit::BinocularPolyFit(const PlType& plType)
+    : PolyFit(plType)
+    , binocularUnknowns(0)
+    , binocularCalibrated(false)
+{
+    binocularUnknowns = 2 * (PolyFit::unknowns - 1) + 1;
+}
+
+string BinocularPolyFit::description() const
+{
+    return "BI" + PolyFit::description();
+}
+
+bool BinocularPolyFit::calibrate(vector<CollectionTuple>& calibrationTuples, QString& errorMsg)
+{
+    if (!PolyFit::calibrate(calibrationTuples, errorMsg))
+        return false;
+
+    PointVector leftPupil;
+    PointVector rightPupil;
+    PointVector gaze;
+
+    for (const auto& tuple : calibrationTuples) {
+        leftPupil.insert(PolyFit::normalize(tuple.lEye.pupil.center, tuple.lEye.input));
+        rightPupil.insert(PolyFit::normalize(tuple.rEye.pupil.center, tuple.rEye.input));
+        gaze.insert(PolyFit::normalize(tuple.field.collectionMarker.center, tuple.field.input));
+    }
+
+    if (gaze.size() < binocularUnknowns) {
+        errorMsg = QString("Not enough calibration points (%1/%2).").arg(gaze.size()).arg(binocularUnknowns);
+        return false;
+    }
+
+    binocularCalibrated = calibrate(plType, leftPupil.x, leftPupil.y, rightPupil.x, rightPupil.y, gaze.x, bcx, errorMsg)
+        && calibrate(plType, leftPupil.x, leftPupil.y, rightPupil.x, rightPupil.y, gaze.y, bcy, errorMsg);
+
+    return binocularCalibrated;
+}
+
+bool BinocularPolyFit::calibrate(const PlType& type, const Mat& lx, const Mat& ly, const Mat& rx, const Mat& ry, const Mat& z, Mat1d& c, QString& errorMsg)
+{
+
+    Mat ones = Mat::ones(lx.rows, 1, DataType<double>::type);
+    Mat lpl = buildMat(type, lx, ly);
+    Mat rpl = buildMat(type, rx, ry);
+    Mat tmp;
+    Mat A;
+    hconcat(ones, lpl, tmp);
+    hconcat(tmp, rpl, A);
+
+    z.convertTo(z, A.type());
+
+    bool result = true;
+    try {
+        result &= solve(A, z, c, DECOMP_SVD);
+    } catch (...) {
+        errorMsg = "cv::solve exception";
+        result = false;
+    }
+
+    if (result) {
+        // TODO: is this still necessary?
+        for (int i = 0; i < c.rows; i++)
+            if (std::isnan(c(i))) {
+                errorMsg = "found NaN coefficients.";
+                result = false;
+            }
+    }
+
+    return result;
+}
+
+void BinocularPolyFit::estimateBinocular2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular)
+{
+    (void)tuple;
+    bool leftValid = left.pupilConfidence >= minPupilConfidence;
+    bool rightValid = right.pupilConfidence >= minPupilConfidence;
+
+    auto evaluate = [&](DataTuple tuple, const Mat1d& c) {
+        Point2d nlp = PolyFit::normalize(tuple.lEye.pupil.center, tuple.lEye.input);
+        Point2d nrp = PolyFit::normalize(tuple.rEye.pupil.center, tuple.rEye.input);
+        double res = c.at<double>(0); // constant contribution
+        res += PolyFit::evaluateVariable(plType, nlp.x, nlp.y, c.rowRange(1, PolyFit::unknowns)); // left eye contribution
+        res += PolyFit::evaluateVariable(plType, nrp.x, nrp.y, c.rowRange(PolyFit::unknowns, binocularUnknowns)); // right eye contribution
+        return res;
+    };
+
+    if (leftValid && rightValid && binocularCalibrated) {
+        binocular.gp = { static_cast<float>(tuple.field.width * evaluate(tuple, bcx)), static_cast<float>(tuple.field.height * evaluate(tuple, bcy)) };
+        binocular.valid = true;
+        binocular.pupilConfidence = 0.5 * (left.pupilConfidence + right.pupilConfidence);
+    } else {
+        if (leftValid)
+            binocular = left;
+        if (rightValid)
+            binocular = right;
+    }
+}

EyeRecToo/src/gaze-estimation/BinocularPolyFit.h

+#ifndef BINOCULARPOLYFIT_H
+#define BINOCULARPOLYFIT_H
+
+#include "PolyFit.h"
+
+class BinocularPolyFit : public PolyFit {
+public:
+    BinocularPolyFit(const PlType& plType);
+    bool calibrate(std::vector<CollectionTuple>& calibrationTuples, QString& errorMsg) override;
+    std::string description() const override;
+
+private:
+    void estimateBinocular2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular) override;
+    int binocularUnknowns;
+    bool binocularCalibrated;
+
+    cv::Mat1d bcx;
+    cv::Mat1d bcy;
+
+    bool calibrate(const PlType& type, const cv::Mat& lx, const cv::Mat& ly, const cv::Mat& rx, const cv::Mat& ry, const cv::Mat& z, cv::Mat1d& c, QString& errorMsg);
+};
+
+#endif // BINOCULARPOLYFIT_H

EyeRecToo/src/gaze-estimation/GazeEstimationMethod.cpp

 #include "gaze-estimation/GazeEstimationMethod.h"
 
-void GazeEstimationMethod::estimateBinocular2d(GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular)
+void GazeEstimationMethod::estimateBinocular2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular)
 {
+    (void)tuple;
     bool leftValid = left.valid && (left.pupilConfidence >= minPupilConfidence);
     bool rightValid = right.valid && (right.pupilConfidence >= minPupilConfidence);
 
@@ -25,14 +26,15 @@ cv::Point2f GazeEstimationMethod::estimate2d(const DataTuple& tuple, GazeEstimat
     right.pupilConfidence = tuple.rEye.pupil.confidence;
 
     implEstimate2d(tuple, left, right);
-    estimateBinocular2d(left, right, binocular);
+    estimateBinocular2d(tuple, left, right, binocular);
     // TODO: unproject to 3d
     return binocular.gp;
 }
 
-void GazeEstimationMethod::estimateBinocular3d(GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular)
+void GazeEstimationMethod::estimateBinocular3d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular)
 {
     // TODO
+    (void)tuple;
     (void)left;
     (void)right;
     (void)binocular;

EyeRecToo/src/gaze-estimation/GazeEstimationMethod.h

@@ -41,7 +41,7 @@ public:
         MONO_RIGHT = 3,
     };
 
-    virtual std::string description() = 0;
+    virtual std::string description() const = 0;
     virtual bool calibrate(std::vector<CollectionTuple>& calibrationTuples, QString& error) = 0;
     void estimate(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
     cv::Point2f estimate2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
@@ -66,8 +66,8 @@ private:
     virtual void implEstimate2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right) = 0;
     virtual void implEstimate3d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right) = 0;
 
-    virtual void estimateBinocular2d(GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
-    virtual void estimateBinocular3d(GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
+    virtual void estimateBinocular2d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
+    virtual void estimateBinocular3d(const DataTuple& tuple, GazeEstimate& left, GazeEstimate& right, GazeEstimate& binocular);
 };
 Q_DECLARE_METATYPE(enum GazeEstimationMethod::InputType)
 

EyeRecToo/src/gaze-estimation/Homography.h

@@ -14,7 +14,7 @@ class Homography : public GazeEstimationMethod {
 public:
     Homography();
     bool calibrate(std::vector<CollectionTuple>& calibrationTuples, QString& errorMsg) override;
-    std::string description() override { return "Homography"; }
+    std::string description() const override { return "Homography"; }
 
 private:
     bool has3d() override { return false; }

EyeRecToo/src/gaze-estimation/PolyFit.cpp

@@ -4,49 +4,48 @@
 using namespace std;
 using namespace cv;
 
-PolyFit::PolyFit(PlType plType)
-    : leftCalibrated(false)
+PolyFit::PolyFit(const PlType& plType)
+    : plType(plType)
+    , leftCalibrated(false)
     , rightCalibrated(false)
 {
-    this->plType = plType;
-
     switch (plType) {
-    case POLY_1_X_Y_XY:
+    case POLY_X_Y_XY:
         unknowns = 4;
         break;
-    case POLY_1_X_Y_XY_XX_YY:
+    case POLY_X_Y_XY_XX_YY:
         unknowns = 6;
         break;
-    case POLY_1_X_Y_XY_XX_YY_XXYY:
+    case POLY_X_Y_XY_XX_YY_XXYY:
         unknowns = 7;
         break;
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX:
         unknowns = 8;
         break;
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY:
         unknowns = 9;
         break;
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
         unknowns = 10;
     default:
         break;
     }
 }
 
-string PolyFit::description()
+string PolyFit::description() const
 {
     switch (plType) {
-    case POLY_1_X_Y_XY:
+    case POLY_X_Y_XY:
         return "POLY_X_Y_XY";
-    case POLY_1_X_Y_XY_XX_YY:
+    case POLY_X_Y_XY_XX_YY:
         return "POLY_X_Y_XY_XX_YY";
-    case POLY_1_X_Y_XY_XX_YY_XXYY:
+    case POLY_X_Y_XY_XX_YY_XXYY:
         return "POLY_X_Y_XY_XX_YY_XXYY";
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX:
         return "POLY_X_Y_XY_XX_YY_XYY_YXX";
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY:
         return "POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY";
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
         return "POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY";
     default:
         return "Unknown";
@@ -94,6 +93,10 @@ void PolyFit::implEstimate2d(const DataTuple& tuple, GazeEstimate& left, GazeEst
 }
 
 double PolyFit::evaluate(const PlType& type, const double x, const double y, const Mat1d& c)
+{
+    return evaluateVariable(type, x, y, c.rowRange(1, unknowns)) + c.at<double>(0);
+}
+double PolyFit::evaluateVariable(const PlType& type, const double x, const double y, const Mat1d& c)
 {
     double xx = x * x;
     double yy = y * y;
@@ -107,28 +110,28 @@ double PolyFit::evaluate(const PlType& type, const double x, const double y, con
     double estimate;
 
     switch (type) {
-    case POLY_1_X_Y_XY:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3);
+    case POLY_X_Y_XY:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2);
         break;
 
-    case POLY_1_X_Y_XY_XX_YY:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3) + xx * c.at<double>(4) + yy * c.at<double>(5);
+    case POLY_X_Y_XY_XX_YY:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2) + xx * c.at<double>(3) + yy * c.at<double>(4);
         break;
 
-    case POLY_1_X_Y_XY_XX_YY_XXYY:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3) + xx * c.at<double>(4) + yy * c.at<double>(5) + xxyy * c.at<double>(6);
+    case POLY_X_Y_XY_XX_YY_XXYY:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2) + xx * c.at<double>(3) + yy * c.at<double>(4) + xxyy * c.at<double>(5);
         break;
 
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3) + xx * c.at<double>(4) + yy * c.at<double>(5) + xyy * c.at<double>(6) + yxx * c.at<double>(7);
+    case POLY_X_Y_XY_XX_YY_XYY_YXX:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2) + xx * c.at<double>(3) + yy * c.at<double>(4) + xyy * c.at<double>(5) + yxx * c.at<double>(6);
         break;
 
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3) + xx * c.at<double>(4) + yy * c.at<double>(5) + xyy * c.at<double>(6) + yxx * c.at<double>(7) + xxyy * c.at<double>(8);
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2) + xx * c.at<double>(3) + yy * c.at<double>(4) + xyy * c.at<double>(5) + yxx * c.at<double>(6) + xxyy * c.at<double>(7);
         break;
 
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
-        estimate = c.at<double>(0) + x * c.at<double>(1) + y * c.at<double>(2) + xy * c.at<double>(3) + xx * c.at<double>(4) + yy * c.at<double>(5) + xyy * c.at<double>(6) + yxx * c.at<double>(7) + xxx * c.at<double>(8) + yyy * c.at<double>(9);
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
+        estimate = x * c.at<double>(0) + y * c.at<double>(1) + xy * c.at<double>(2) + xx * c.at<double>(3) + yy * c.at<double>(4) + xyy * c.at<double>(5) + yxx * c.at<double>(6) + xxx * c.at<double>(7) + yyy * c.at<double>(8);
         break;
 
     default:
@@ -138,7 +141,7 @@ double PolyFit::evaluate(const PlType& type, const double x, const double y, con
     return estimate;
 }
 
-bool PolyFit::calibrate(const PlType& type, const cv::Mat& x, const cv::Mat& y, const cv::Mat& z, cv::Mat1d& c, QString& errorMsg)
+cv::Mat PolyFit::buildMat(const PlType& type, const cv::Mat& x, const cv::Mat& y)
 {
     // Possible terms; not using more than third degree atm
     Mat xx, yy, xy, xyy, yxx, xxyy, xxx, yyy, A;
@@ -151,75 +154,77 @@ bool PolyFit::calibrate(const PlType& type, const cv::Mat& x, const cv::Mat& y,
     multiply(xx, yy, xxyy);
     multiply(x, xx, xxx);
     multiply(y, yy, yyy);
-    A = Mat::ones(z.rows, unknowns, CV_64F);
+    A = Mat::zeros(x.rows, unknowns - 1, DataType<double>::type);
 
     switch (type) {
-    case POLY_1_X_Y_XY:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        break;
-
-    case POLY_1_X_Y_XY_XX_YY:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        xx.copyTo(A.col(4));
-        yy.copyTo(A.col(5));
-        break;
-
-    case POLY_1_X_Y_XY_XX_YY_XXYY:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        xx.copyTo(A.col(4));
-        yy.copyTo(A.col(5));
-        xxyy.copyTo(A.col(6));
-        break;
-
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        xx.copyTo(A.col(4));
-        yy.copyTo(A.col(5));
-        xyy.copyTo(A.col(6));
-        yxx.copyTo(A.col(7));
-        break;
-
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        xx.copyTo(A.col(4));
-        yy.copyTo(A.col(5));
-        xyy.copyTo(A.col(6));
-        yxx.copyTo(A.col(7));
-        xxyy.copyTo(A.col(8));
-        break;
-
-    case POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
-        // A.col(0) was initialized with one
-        x.copyTo(A.col(1));
-        y.copyTo(A.col(2));
-        xy.copyTo(A.col(3));
-        xx.copyTo(A.col(4));
-        yy.copyTo(A.col(5));
-        xyy.copyTo(A.col(6));
-        yxx.copyTo(A.col(7));
-        xxx.copyTo(A.col(8));
-        yyy.copyTo(A.col(9));
-        break;
-
-    default:
+    case POLY_X_Y_XY:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        break;
+
+    case POLY_X_Y_XY_XX_YY:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        xx.copyTo(A.col(3));
+        yy.copyTo(A.col(4));
+        break;
+
+    case POLY_X_Y_XY_XX_YY_XXYY:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        xx.copyTo(A.col(3));
+        yy.copyTo(A.col(4));
+        xxyy.copyTo(A.col(5));
+        break;
+
+    case POLY_X_Y_XY_XX_YY_XYY_YXX:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        xx.copyTo(A.col(3));
+        yy.copyTo(A.col(4));
+        xyy.copyTo(A.col(5));
+        yxx.copyTo(A.col(6));
+        break;
+
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        xx.copyTo(A.col(3));
+        yy.copyTo(A.col(4));
+        xyy.copyTo(A.col(5));
+        yxx.copyTo(A.col(6));
+        xxyy.copyTo(A.col(7));
+        break;
+
+    case POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY:
+        x.copyTo(A.col(0));
+        y.copyTo(A.col(1));
+        xy.copyTo(A.col(2));
+        xx.copyTo(A.col(3));
+        yy.copyTo(A.col(4));
+        xyy.copyTo(A.col(5));
+        yxx.copyTo(A.col(6));
+        xxx.copyTo(A.col(7));
+        yyy.copyTo(A.col(8));
         break;
     }
 
+    return A;
+}
+
+bool PolyFit::calibrate(const PlType& type, const cv::Mat& x, const cv::Mat& y, const cv::Mat& z, cv::Mat1d& c, QString& errorMsg)
+{
+
+    Mat ones = Mat::ones(x.rows, 1, DataType<double>::type);
+    Mat pl = buildMat(type, x, y);
+    Mat A;
+    hconcat(ones, pl, A);
+
     z.convertTo(z, A.type());
 
     bool result = true;

EyeRecToo/src/gaze-estimation/PolyFit.h

@@ -23,24 +23,24 @@ public:
         y.push_back(p.y);
         z.push_back(p.y);
     }
-    size_t size() { return static_cast<size_t>(x.rows); }
+    size_t size() const { return static_cast<size_t>(x.rows); }
     cv::Mat x, y, z;
 };
 
 class PolyFit : public GazeEstimationMethod {
 public:
     enum PlType {
-        POLY_1_X_Y_XY,
-        POLY_1_X_Y_XY_XX_YY,
-        POLY_1_X_Y_XY_XX_YY_XXYY,
-        POLY_1_X_Y_XY_XX_YY_XYY_YXX,
-        POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXYY,
-        POLY_1_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY
+        POLY_X_Y_XY,
+        POLY_X_Y_XY_XX_YY,
+        POLY_X_Y_XY_XX_YY_XXYY,
+        POLY_X_Y_XY_XX_YY_XYY_YXX,
+        POLY_X_Y_XY_XX_YY_XYY_YXX_XXYY,
+        POLY_X_Y_XY_XX_YY_XYY_YXX_XXX_YYY
     };
-    PolyFit(PlType plType);
+    PolyFit(const PlType& plType);
 
     bool calibrate(std::vector<CollectionTuple>& calibrationTuples, QString& errorMsg) override;
-    std::string description() override;
+    std::string description() const override;
 
 private:
     bool has3d() override { return false; }
@@ -52,6 +52,7 @@ private:
         (void)right;
     }
 
+protected:
     // Implementation related
     bool leftCalibrated;
     bool rightCalibrated;
@@ -64,8 +65,10 @@ private:
 
     unsigned int unknowns;
 
+    cv::Mat buildMat(const PlType& type, const cv::Mat& x, const cv::Mat& y);
     bool calibrate(const PlType& type, const cv::Mat& x, const cv::Mat& y, const cv::Mat& z, cv::Mat1d& c, QString& errorMsg);
     double evaluate(const PlType& type, const double x, const double y, const cv::Mat1d& c);
+    double evaluateVariable(const PlType& type, const double x, const double y, const cv::Mat1d& c);
 
     template <typename T>
     cv::Point2d normalize(const T& p, const cv::Mat& m) { return { p.x / m.cols, p.y / m.rows }; }