Slightly reworked rewinding; added values rewinding.

Needs a solution to the problem of which values must be displayed
and/or interpolated. Currently, whenever the user rewinds, the
current error measure is displayed, regardless of what was being
displayed before. This will probably be trivial to solve once we
have a nice way of changing the current measure.

* Also changed all OpenMP-powered for loops to use signed integers,
  requirements of OMP2.x (which is what MSVC supports currently)
* The above change comes with a new header for utility functions

11 files changed, 111 insertions, 50 deletions

diff --git a/barchart.cpp b/barchart.cpp
index 70c6004..3adf441 100644
--- a/barchart.cpp
+++ b/barchart.cpp
@@ -105,6 +105,7 @@ void BarChart::brushItem(int item)
         m_brushedItem = item;
         emit itemBrushed(m_brushedItem);
     } else {
+        // safe comparison: we just checked for negative values
         if (m_values.n_elem == 0 || item > m_values.n_elem - 1) {
             return;
         }
diff --git a/dist.cpp b/dist.cpp
index 0a241ec..f307564 100644
--- a/dist.cpp
+++ b/dist.cpp
@@ -1,5 +1,7 @@
 #include "mp.h"
 
+#include "utils.h"
+
 double mp::euclidean(const arma::rowvec &x1, const arma::rowvec &x2)
 {
     return arma::norm(x1 - x2, 2);
@@ -7,11 +9,11 @@ double mp::euclidean(const arma::rowvec &x1, const arma::rowvec &x2)
 
 arma::mat mp::dist(const arma::mat &X, mp::DistFunc dfunc)
 {
-    arma::uword n = X.n_rows;
+    int n = uintToInt<arma::uword, int>(X.n_rows);
     arma::mat D(n, n, arma::fill::zeros);
 
-    #pragma omp parallel for shared(X, D)
-    for (arma::uword i = 0; i < n; i++) {
+    #pragma omp parallel for shared(X, D, n)
+    for (int i = 0; i < n; i++) {
         for (arma::uword j = 0; j < i; j++) {
             D(i, j) = dfunc(X.row(i), X.row(j));
             D(j, i) = D(i, j);
diff --git a/lamp.cpp b/lamp.cpp
index 293ce30..77bbac4 100644
--- a/lamp.cpp
+++ b/lamp.cpp
@@ -2,6 +2,8 @@
 
 #include <algorithm>
 
+#include "utils.h"
+
 static const double EPSILON = 1e-3;
 
 arma::mat mp::lamp(const arma::mat &X, const arma::uvec &sampleIndices, const arma::mat &Ys)
@@ -13,11 +15,12 @@ arma::mat mp::lamp(const arma::mat &X, const arma::uvec &sampleIndices, const ar
 
 void mp::lamp(const arma::mat &X, const arma::uvec &sampleIndices, const arma::mat &Ys, arma::mat &Y)
 {
+    int n = uintToInt<arma::uword, int>(X.n_rows);
     const arma::mat &Xs = X.rows(sampleIndices);
     arma::uword sampleSize = sampleIndices.n_elem;
 
-    #pragma omp parallel for shared(X, Xs, Ys, Y)
-    for (arma::uword i = 0; i < X.n_rows; i++) {
+    #pragma omp parallel for shared(X, Xs, Ys, Y, n)
+    for (int i = 0; i < n; i++) {
         const arma::rowvec &point = X.row(i);
 
         // calculate alphas
diff --git a/main.cpp b/main.cpp
index b1a8e8b..c2d79df 100644
--- a/main.cpp
+++ b/main.cpp
@@ -174,8 +174,6 @@ int main(int argc, char **argv)
             m->rpPlot, SLOT(setXY(const arma::mat &)));
     QObject::connect(&manipulationHandler, SIGNAL(rpChanged(const arma::mat &)),
             m->splat, SLOT(setSites(const arma::mat &)));
-    QObject::connect(plotTC, SIGNAL(tChanged(double)),
-            &manipulationHandler, SLOT(setRewind(double)));
 
     // Keep both scatterplots and the splat scaled equally and relative to the
     // full plot
@@ -250,6 +248,23 @@ int main(int argc, char **argv)
     QObject::connect(&rpSelectionHandler, SIGNAL(selectionChanged(const std::vector<bool> &)),
             &projectionObserver, SLOT(setRPSelection(const std::vector<bool> &)));
 
+    // Connect visual components (but not barcharts) to rewinding mechanism
+    QObject::connect(plotTC, SIGNAL(tChanged(double)),
+            &manipulationHandler, SLOT(setRewind(double)));
+    QObject::connect(&manipulationHandler, SIGNAL(cpRewound(const arma::mat &)),
+            m->cpPlot, SLOT(setXY(const arma::mat &)));
+    QObject::connect(&manipulationHandler, SIGNAL(rpRewound(const arma::mat &)),
+            m->rpPlot, SLOT(setXY(const arma::mat &)));
+    QObject::connect(&manipulationHandler, SIGNAL(rpRewound(const arma::mat &)),
+            m->splat, SLOT(setSites(const arma::mat &)));
+
+    QObject::connect(plotTC, SIGNAL(tChanged(double)),
+            m->projectionObserver, SLOT(setRewind(double)));
+    QObject::connect(m->projectionObserver, SIGNAL(cpValuesRewound(const arma::vec &)),
+            m->cpPlot, SLOT(setColorData(const arma::vec &)));
+    QObject::connect(m->projectionObserver, SIGNAL(rpValuesRewound(const arma::vec &)),
+            m->splat, SLOT(setValues(const arma::vec &)));
+
     // General component set up
     plotTC->setAcceptedMouseButtons(Qt::MiddleButton);
     m->cpPlot->setDragEnabled(true);
diff --git a/manipulationhandler.cpp b/manipulationhandler.cpp
index 7d3e3c0..707d4c7 100644
--- a/manipulationhandler.cpp
+++ b/manipulationhandler.cpp
@@ -68,11 +68,7 @@ void ManipulationHandler::setRewind(double t)
     }
 
     arma::mat Y = m_Y * t + m_prevY * (1.0 - t);
-    emit cpChanged(Y.rows(m_cpIndices));
-    emit rpChanged(Y.rows(m_rpIndices));
-
-    // NOTE: this signal was supposed to be emitted, but since we don't want
-    // anything besides graphical objects to know the projection is being
-    // rewound, this is (for now) left out.
-    // emit mapChanged(Y);
+    emit cpRewound(Y.rows(m_cpIndices));
+    emit rpRewound(Y.rows(m_rpIndices));
+    emit mapRewound(Y);
 }
diff --git a/manipulationhandler.h b/manipulationhandler.h
index 1e17446..e099da8 100644
--- a/manipulationhandler.h
+++ b/manipulationhandler.h
@@ -24,6 +24,9 @@ signals:
     void cpChanged(const arma::mat &cpY) const;
     void rpChanged(const arma::mat &rpY) const;
     void mapChanged(const arma::mat &Y) const;
+    void cpRewound(const arma::mat &cpY) const;
+    void rpRewound(const arma::mat &rpY) const;
+    void mapRewound(const arma::mat &Y) const;
 
 public slots:
     void setCP(const arma::mat &Ys);
diff --git a/measures.cpp b/measures.cpp
index daf3c78..2a3908f 100644
--- a/measures.cpp
+++ b/measures.cpp
@@ -4,15 +4,19 @@
 #include <cmath>
 #include <algorithm>
 
+#include "utils.h"
+
+static const float EPSILON = 1e-6f;
+
 arma::vec mp::neighborhoodPreservation(const arma::mat &distA,
                                        const arma::mat &distB,
                                        arma::uword k)
 {
-    arma::uword n = distA.n_rows;
+    int n = uintToInt<arma::uword, int>(distA.n_rows);
     arma::vec np(n);
 
     #pragma omp parallel for shared(np, n)
-    for (arma::uword i = 0; i < n; i++) {
+    for (int i = 0; i < n; i++) {
         arma::uvec nnA(k);
         arma::uvec nnB(k);
         arma::vec dist(k);
@@ -39,6 +43,32 @@ arma::vec mp::silhouette(const arma::mat &distA,
     return arma::vec(distA.n_rows, arma::fill::zeros);
 }
 
+void mp::aggregatedError(const arma::mat &distX,
+                         const arma::mat &distY,
+                         arma::vec &v)
+{
+    int n = uintToInt<arma::uword, int>(v.n_elem);
+    double maxX = distX.max();
+    double maxY = distY.max();
+
+    #pragma omp parallel for shared(maxX, maxY, distX, distY, v, n)
+    for (int i = 0; i < n; i++) {
+        v[i] = 0;
+        for (int j = 0; j < n; j++) {
+            if (i == j) {
+                continue;
+            }
+
+            double diff = fabs(distY(i, j) / maxY - distX(i, j) / maxX);
+            if (diff < EPSILON) {
+                continue;
+            }
+
+            v[i] += diff;
+        }
+    }
+}
+
 /*
 double mp::stress(const arma::mat &Dp, const arma::mat &Dq)
 {
@@ -66,16 +96,16 @@ arma::vec mp::klDivergence(const arma::mat &P, const arma::mat &Q)
     return diver;
 }
 
-void mp::klDivergence(const arma::mat &P, const arma::mat &Q, arma::vec &diver)
+void mp::klDivergence(const arma::mat &P, const arma::mat &Q, arma::vec &diverg)
 {
     assert(P.n_rows == P.n_cols);
     assert(Q.n_rows == Q.n_cols);
     assert(P.n_rows == Q.n_cols);
-    assert(diver.n_elem == P.n_rows);
+    assert(diverg.n_elem == P.n_rows);
 
     arma::uword n = P.n_rows;
     for (arma::uword i = 0; i < n; i++)
-        diver(i) = klDivergence(P.row(i), Q.row(i));
+        diverg(i) = klDivergence(P.row(i), Q.row(i));
 }
 
 double mp::klDivergence(const arma::rowvec &pi, const arma::rowvec &qi)
diff --git a/mp.h b/mp.h
index 2f32c8e..8f555a2 100644
--- a/mp.h
+++ b/mp.h
@@ -13,9 +13,9 @@ arma::mat dist(const arma::mat &X, DistFunc dfunc = euclidean);
 void knn(const arma::mat &dmat, arma::uword i, arma::uword k, arma::uvec &nn, arma::vec &dist);
 
 // Evaluation measures
-typedef arma::vec (*MeasureFunc)(const arma::mat &distA, const arma::mat &distB);
 arma::vec neighborhoodPreservation(const arma::mat &distA, const arma::mat &distB, arma::uword k = 10);
 arma::vec silhouette(const arma::mat &distA, const arma::mat &distB, const arma::vec &labels);
+void aggregatedError(const arma::mat &distX, const arma::mat &distY, arma::vec &v);
 
 // Techniques
 arma::mat lamp(const arma::mat &X, const arma::uvec &sampleIndices, const arma::mat &Ys);
diff --git a/projectionobserver.cpp b/projectionobserver.cpp
index 66c140e..bf73ded 100644
--- a/projectionobserver.cpp
+++ b/projectionobserver.cpp
@@ -7,31 +7,6 @@
 #include "mp.h"
 #include "numericrange.h"
 
-static const float EPSILON = 1e-6f;
-
-static void aggregatedError(const arma::mat &distX, const arma::mat &distY, arma::vec &v)
-{
-    double maxX = distX.max();
-    double maxY = distY.max();
-
-    #pragma omp parallel for shared(maxX, maxY, distX, distY, v)
-    for (arma::uword i = 0; i < v.n_elem; i++) {
-        v[i] = 0;
-        for (arma::uword j = 0; j < v.n_elem; j++) {
-            if (i == j) {
-                continue;
-            }
-
-            float diff = fabs(distY(i, j) / maxY - distX(i, j) / maxX);
-            if (diff < EPSILON) {
-                continue;
-            }
-
-            v[i] += diff;
-        }
-    }
-}
-
 ProjectionObserver::ProjectionObserver(const arma::mat &X,
                                        const arma::uvec &cpIndices)
     : m_type(OBSERVER_CURRENT)
@@ -82,10 +57,10 @@ void ProjectionObserver::setMap(const arma::mat &Y)
 
     m_Y = Y;
     m_distY = mp::dist(Y);
-    aggregatedError(m_distX, m_distY, m_values);
+    mp::aggregatedError(m_distX, m_distY, m_values);
 
     // method called for the first time; set original Y
-    if (m_origY.n_elem == 0) {
+    if (m_origY.n_elem != m_Y.n_elem) {
         m_origY = m_Y;
         m_origDistY = m_distY;
         m_origValues = m_values;
@@ -172,7 +147,7 @@ bool ProjectionObserver::emitValuesChanged() const
         emit valuesChanged(m_values);
         return true;
     case OBSERVER_DIFF_PREVIOUS:
-        if (m_prevValues.n_elem > 0) {
+        if (m_prevValues.n_elem == m_values.n_elem) {
             arma::vec diff = m_values - m_prevValues;
             emit rpValuesChanged(diff(m_rpIndices));
             emit valuesChanged(diff);
@@ -180,7 +155,7 @@ bool ProjectionObserver::emitValuesChanged() const
         }
         return false;
     case OBSERVER_DIFF_ORIGINAL:
-        if (m_origValues.n_elem > 0) {
+        if (m_origValues.n_elem == m_values.n_elem) {
             arma::vec diff = m_values - m_origValues;
             emit rpValuesChanged(diff(m_rpIndices));
             emit valuesChanged(diff);
@@ -191,3 +166,16 @@ bool ProjectionObserver::emitValuesChanged() const
         return false;
     }
 }
+
+void ProjectionObserver::setRewind(double t)
+{
+    if (m_prevValues.n_elem != m_values.n_elem) {
+        return;
+    }
+
+    arma::vec values = m_values * t + m_prevValues * (1.0 - t);
+
+    emit cpValuesRewound(values(m_cpIndices));
+    emit rpValuesRewound(values(m_rpIndices));
+    emit valuesRewound(values);
+}
diff --git a/projectionobserver.h b/projectionobserver.h
index badca4d..db68c73 100644
--- a/projectionobserver.h
+++ b/projectionobserver.h
@@ -19,12 +19,16 @@ signals:
     void valuesChanged(const arma::vec &values) const;
     void cpValuesChanged(const arma::vec &values) const;
     void rpValuesChanged(const arma::vec &values) const;
+    void valuesRewound(const arma::vec &values) const;
+    void cpValuesRewound(const arma::vec &values) const;
+    void rpValuesRewound(const arma::vec &values) const;
 
 public slots:
     void setMap(const arma::mat &Y);
     bool setType(int type);
     void setCPSelection(const std::vector<bool> &cpSelection);
     void setRPSelection(const std::vector<bool> &rpSelection);
+    void setRewind(double t);
 
 private:
     bool emitValuesChanged() const;
@@ -37,7 +41,7 @@ private:
     bool m_cpSelectionEmpty, m_rpSelectionEmpty;
     std::vector<int> m_cpSelection, m_rpSelection;
 
-    // alpha(i, j): the influence CP j has on RP i 
+    // alpha(i, j): the influence CP j has on RP i
     void computeAlphas();
     arma::mat m_alphas, m_influences;
 
diff --git a/utils.h b/utils.h
new file mode 100644
index 0000000..ba13f27
--- /dev/null
+++ b/utils.h
@@ -0,0 +1,19 @@
+#include <limits>
+#include <stdexcept>
+
+/*
+ * Credits to:
+ * http://stackoverflow.com/questions/13150449/efficient-unsigned-to-signed-cast-avoiding-implementation-defined-behavior
+ */
+template<typename Uint, typename Int>
+Int uintToInt(Uint x)
+{
+    if (x <= std::numeric_limits<Int>::max())
+        return static_cast<Int>(x);
+
+    if (x >= std::numeric_limits<Int>::min())
+        return static_cast<Int>(x - std::numeric_limits<Int>::min())
+             + std::numeric_limits<Int>::min();
+
+    throw std::overflow_error("given value does not fit integer type");
+}