aboutsummaryrefslogtreecommitdiff
diff options
context:
space:
mode:
-rw-r--r--dist.cpp20
-rw-r--r--forceScheme.cpp45
-rw-r--r--glyph.cpp103
-rw-r--r--glyph.h42
-rw-r--r--lamp.cpp50
-rw-r--r--main.cpp48
-rw-r--r--main_view.qml40
-rw-r--r--mp.h11
-rw-r--r--pm.pro16
-rw-r--r--pm.qrc5
-rw-r--r--scatterplot.cpp60
-rw-r--r--scatterplot.h26
12 files changed, 466 insertions, 0 deletions
diff --git a/dist.cpp b/dist.cpp
new file mode 100644
index 0000000..aaa2167
--- /dev/null
+++ b/dist.cpp
@@ -0,0 +1,20 @@
+#include "mp.h"
+
+double mp::euclidean(const arma::rowvec &x1, const arma::rowvec &x2)
+{
+ return arma::norm(x1 - x2, 2);
+}
+
+arma::mat mp::dist(const arma::mat &X, double (*distCalc)(const arma::rowvec &, const arma::rowvec &))
+{
+ arma::uword n = X.n_rows;
+ arma::mat D(n, n, arma::fill::zeros);
+
+ for (arma::uword i = 0; i < n; i++)
+ for (arma::uword j = 0; j < i; j++) {
+ D(i, j) = distCalc(X.row(i), X.row(j));
+ D(j, i) = D(i, j);
+ }
+
+ return D;
+}
diff --git a/forceScheme.cpp b/forceScheme.cpp
new file mode 100644
index 0000000..3683c10
--- /dev/null
+++ b/forceScheme.cpp
@@ -0,0 +1,45 @@
+#include "mp.h"
+
+#include <algorithm>
+#include <vector>
+
+typedef arma::uvec V;
+
+arma::mat mp::forceScheme(const arma::mat &D,
+ arma::mat &Y,
+ size_t maxIter,
+ double tol,
+ double fraction)
+{
+
+ arma::uword n = (arma::uword) Y.n_rows;
+ V i(n), j(n);
+ for (arma::uword k = 0; k < n; k++)
+ i[k] = j[k] = k;
+
+ double prev_delta_sum = 1. / 0.;
+ for (size_t iter = 0; iter < maxIter; iter++) {
+ double delta_sum = 0;
+
+ arma::shuffle(i);
+ for (V::iterator a = i.begin(); a != i.end(); a++) {
+ arma::shuffle(j);
+ for (V::iterator b = j.begin(); b != j.end(); b++) {
+ if (*a == *b)
+ continue;
+
+ arma::rowvec direction(Y.row(*b) - Y.row(*a));
+ double d2 = std::max(arma::norm(direction, 2), mp::EPSILON);
+ double delta = (D(*a, *b) - d2) / fraction;
+ delta_sum += fabs(delta);
+ Y.row(*b) += delta * (direction / d2);
+ }
+ }
+
+ if (fabs(prev_delta_sum - delta_sum) < tol)
+ break;
+ prev_delta_sum = delta_sum;
+ }
+
+ return Y;
+}
diff --git a/glyph.cpp b/glyph.cpp
new file mode 100644
index 0000000..f143ee8
--- /dev/null
+++ b/glyph.cpp
@@ -0,0 +1,103 @@
+#include "glyph.h"
+
+#include <cmath>
+#include <QSGGeometry>
+#include <QSGGeometryNode>
+#include <QSGMaterial>
+#include <QSGFlatColorMaterial>
+
+const float PI = 3.1415f;
+
+Glyph::Glyph(QQuickItem *parent)
+ : QQuickItem(parent)
+ , m_size(0)
+{
+ setFlag(QQuickItem::ItemHasContents);
+}
+
+void Glyph::setSize(qreal size)
+{
+ if (size == m_size)
+ return;
+
+ m_size = size;
+ setWidth(size);
+ setHeight(size);
+ emit sizeChanged();
+ update();
+}
+
+void Glyph::setColor(const QColor &color)
+{
+ if (color == m_color)
+ return;
+
+ m_color = color;
+ emit colorChanged();
+ update();
+}
+
+int calculateCircleVertexCount(qreal radius)
+{
+ // 10 * sqrt(r) \approx 2*pi / acos(1 - 1 / (4*r))
+ return (int) (10.0 * sqrt(radius));
+}
+
+void updateCircleGeometry(QSGGeometry *geometry, const QRectF &bounds)
+{
+ int vertexCount = geometry->vertexCount();
+ float cy = bounds.center().x();
+ float cx = bounds.center().y();
+
+ float theta = 2 * PI / float(vertexCount);
+ float c = cosf(theta);
+ float s = sinf(theta);
+ float x = bounds.width() / 2;
+ float y = 0;
+
+ QSGGeometry::Point2D *vertexData = geometry->vertexDataAsPoint2D();
+ for (int i = 0; i < vertexCount; i++) {
+ vertexData[i].set(x + cx, y + cy);
+
+ float t = x;
+ x = c*x - s*y;
+ y = s*t + c*y;
+ }
+}
+
+QSGNode *Glyph::updatePaintNode(QSGNode *oldNode, UpdatePaintNodeData *)
+{
+ QSGGeometryNode *node = 0;
+ QSGGeometry *geometry = 0;
+ QSGFlatColorMaterial *material = 0;
+ int vertexCount = calculateCircleVertexCount(m_size / 2);
+
+ if (!oldNode) {
+ node = new QSGGeometryNode;
+
+ geometry = new QSGGeometry(QSGGeometry::defaultAttributes_Point2D(), vertexCount);
+ geometry->setDrawingMode(GL_POLYGON);
+ node->setGeometry(geometry);
+ node->setFlag(QSGNode::OwnsGeometry);
+
+ material = new QSGFlatColorMaterial;
+ material->setColor(m_color);
+ node->setMaterial(material);
+ node->setFlag(QSGNode::OwnsMaterial);
+ } else {
+ node = static_cast<QSGGeometryNode *>(oldNode);
+ geometry = node->geometry();
+ geometry->allocate(vertexCount);
+ }
+
+ updateCircleGeometry(geometry, boundingRect());
+ node->markDirty(QSGNode::DirtyGeometry);
+
+ material = static_cast<QSGFlatColorMaterial *>(node->material());
+ if (material->color() != m_color) {
+ material->setColor(m_color);
+ node->markDirty(QSGNode::DirtyMaterial);
+ }
+
+ return node;
+}
diff --git a/glyph.h b/glyph.h
new file mode 100644
index 0000000..7532d04
--- /dev/null
+++ b/glyph.h
@@ -0,0 +1,42 @@
+#ifndef GLYPH_H
+#define GLYPH_H
+
+#include <QColor>
+#include <QtQuick/QQuickItem>
+
+class Glyph : public QQuickItem
+{
+ Q_OBJECT
+ Q_PROPERTY(QColor color READ color WRITE setColor NOTIFY colorChanged)
+ Q_PROPERTY(qreal size READ size WRITE setSize NOTIFY sizeChanged)
+
+public:
+ enum GlyphType {
+ GLYPH_CIRCLE,
+ GLYPH_SQUARE,
+ GLYPH_STAR,
+ GLYPH_CROSS
+ };
+
+ Glyph(QQuickItem *parent = 0);
+
+ QSGNode *updatePaintNode(QSGNode *oldNode, UpdatePaintNodeData *);
+
+ QColor color() const { return m_color; }
+ void setColor(const QColor &color);
+
+ qreal size() const { return m_size; }
+ void setSize(qreal size);
+
+signals:
+ void colorChanged();
+ void sizeChanged();
+
+public slots:
+
+private:
+ QColor m_color;
+ qreal m_size;
+};
+
+#endif // GLYPH_H
diff --git a/lamp.cpp b/lamp.cpp
new file mode 100644
index 0000000..a9d7a8b
--- /dev/null
+++ b/lamp.cpp
@@ -0,0 +1,50 @@
+#include "mp.h"
+
+#include <algorithm>
+
+arma::mat mp::lamp(const arma::mat &X, const arma::uvec &sampleIndices, const arma::mat &Ys)
+{
+ arma::mat Xs = X.rows(sampleIndices);
+ arma::uword sampleSize = sampleIndices.n_elem;
+ arma::mat projection(X.n_rows, 2);
+
+ for (arma::uword i = 0; i < X.n_rows; i++) {
+ arma::rowvec point = X.row(i);
+
+ // calculate alphas
+ arma::rowvec alphas(sampleSize);
+ for (arma::uword j = 0; j < sampleSize; j++) {
+ double dist = arma::accu(arma::square(Xs.row(j) - point));
+ alphas[j] = 1. / std::max(dist, mp::EPSILON);
+ }
+
+ double alphas_sum = arma::accu(alphas);
+ arma::rowvec alphas_sqrt = arma::sqrt(alphas);
+
+ // calculate \tilde{X} and \tilde{Y}
+ arma::rowvec Xtil = arma::sum(alphas * Xs, 0) / alphas_sum;
+ arma::rowvec Ytil = arma::sum(alphas * Ys, 0) / alphas_sum;
+
+ // calculate \hat{X} and \hat{Y}
+ arma::mat Xhat = Xs;
+ Xhat.each_row() -= Xtil;
+ arma::mat Yhat = Ys;
+ Yhat.each_row() -= Ytil;
+
+ // calculate A and B
+ arma::mat At = Xhat.t();
+ At.each_row() %= alphas_sqrt;
+ arma::mat B = Yhat;
+ B.each_col() %= alphas_sqrt.t();
+
+ arma::mat U, V;
+ arma::vec s;
+ arma::svd(U, s, V, At * B);
+ arma::mat M = U.cols(0, 1) * V.t();
+
+ // the projection of point i
+ projection.row(i) = (point - Xtil) * M + Ytil;
+ }
+
+ return projection;
+}
diff --git a/main.cpp b/main.cpp
new file mode 100644
index 0000000..863a20d
--- /dev/null
+++ b/main.cpp
@@ -0,0 +1,48 @@
+#include <cmath>
+#include <QGuiApplication>
+#include <QtQuick/QQuickView>
+
+#include "mp.h"
+#include "glyph.h"
+#include "scatterplot.h"
+
+arma::uvec getSample(arma::uword n)
+{
+ return arma::randi<arma::uvec>((arma::uword) 3*sqrt(n), arma::distr_param(0, n-1));
+}
+
+arma::mat getProjection(const arma::mat &X)
+{
+ arma::uword n = X.n_rows;
+ arma::uvec sampleIndices = getSample(n);
+ arma::mat Ys = arma::randn(sampleIndices.n_elem, 2);
+ Ys = mp::forceScheme(mp::dist(X.rows(sampleIndices)), Ys);
+ return mp::lamp(X, sampleIndices, Ys);
+}
+
+int main(int argc, char **argv)
+{
+ QGuiApplication app(argc, argv);
+
+ qmlRegisterType<Scatterplot>("PM", 1, 0, "Scatterplot");
+ qmlRegisterType<Glyph>("PM", 1, 0, "Glyph");
+
+ QQuickView view;
+ QSurfaceFormat format = view.format();
+ format.setSamples(16);
+ view.setFormat(format);
+ view.setResizeMode(QQuickView::SizeRootObjectToView);
+ view.setSource(QUrl("qrc:///main_view.qml"));
+
+ if (argc > 1) {
+ arma::mat X;
+ X.load(argv[1], arma::raw_ascii);
+ arma::mat projection = getProjection(X);
+ Scatterplot *plot = view.rootObject()->findChild<Scatterplot *>("plot", Qt::FindDirectChildrenOnly);
+ plot->setData(projection);
+ }
+
+ view.show();
+
+ return app.exec();
+}
diff --git a/main_view.qml b/main_view.qml
new file mode 100644
index 0000000..bf143fc
--- /dev/null
+++ b/main_view.qml
@@ -0,0 +1,40 @@
+import QtQuick 2.0
+import PM 1.0
+
+Item {
+ width: 480
+ height: 480
+
+ Scatterplot {
+ id: plot
+ objectName: "plot"
+ anchors.fill: parent
+ }
+
+ Item {
+ id: messageBox
+ anchors.right: parent.right
+ anchors.left: parent.left
+ anchors.bottom: parent.bottom
+
+ Rectangle {
+ color: Qt.rgba(1, 1, 1, 0.7)
+ radius: 5
+ border.width: 1
+ border.color: "white"
+ anchors.fill: messageLabel
+ anchors.margins: -10
+ }
+
+ Text {
+ id: messageLabel
+ color: "black"
+ wrapMode: Text.WordWrap
+ text: "The background here is a squircle rendered with raw OpenGL using the 'beforeRender()' signal in QQuickWindow. This text label and its border is rendered using QML"
+ anchors.right: parent.right
+ anchors.left: parent.left
+ anchors.bottom: parent.bottom
+ anchors.margins: 20
+ }
+ }
+}
diff --git a/mp.h b/mp.h
new file mode 100644
index 0000000..b3f9c75
--- /dev/null
+++ b/mp.h
@@ -0,0 +1,11 @@
+#include <armadillo>
+
+namespace mp {
+
+static const double EPSILON = 1e-3;
+double euclidean(const arma::rowvec &x1, const arma::rowvec &x2);
+arma::mat dist(const arma::mat &X, double (*distCalc)(const arma::rowvec &, const arma::rowvec &) = euclidean);
+arma::mat lamp(const arma::mat &X, const arma::uvec &sampleIndices, const arma::mat &Ys);
+arma::mat forceScheme(const arma::mat &D, arma::mat &Y, size_t maxIter = 20, double tol = 1e-3, double fraction = 0.25);
+
+} // namespace mp
diff --git a/pm.pro b/pm.pro
new file mode 100644
index 0000000..a8d87d7
--- /dev/null
+++ b/pm.pro
@@ -0,0 +1,16 @@
+QT += qml quick
+
+QMAKE_LIBS += -larmadillo
+HEADERS += glyph.h \
+ scatterplot.h \
+ mp.h
+SOURCES += main.cpp \
+ glyph.cpp \
+ scatterplot.cpp \
+ lamp.cpp \
+ forceScheme.cpp \
+ dist.cpp
+RESOURCES += pm.qrc
+
+target.path = .
+INSTALLS += target
diff --git a/pm.qrc b/pm.qrc
new file mode 100644
index 0000000..a3059be
--- /dev/null
+++ b/pm.qrc
@@ -0,0 +1,5 @@
+<RCC>
+ <qresource prefix="/">
+ <file>main_view.qml</file>
+ </qresource>
+</RCC>
diff --git a/scatterplot.cpp b/scatterplot.cpp
new file mode 100644
index 0000000..e775661
--- /dev/null
+++ b/scatterplot.cpp
@@ -0,0 +1,60 @@
+#include "scatterplot.h"
+
+#include <cstdio>
+#include "glyph.h"
+#include <QSGNode>
+#include <QSGGeometryNode>
+
+const int GLYPH_SIZE = 10;
+
+Scatterplot::Scatterplot()
+{
+ setFlag(QQuickItem::ItemHasContents);
+}
+
+Scatterplot::~Scatterplot()
+{
+}
+
+void Scatterplot::setData(const arma::mat &data)
+{
+ if (data.n_cols != 2)
+ return;
+
+ m_data = data;
+ qreal xmin = m_data.col(0).min(),
+ xmax = m_data.col(0).max(),
+ ymin = m_data.col(1).min(),
+ ymax = m_data.col(1).max();
+
+ for (arma::uword i = 0; i < m_data.n_rows; i++) {
+ arma::rowvec row = m_data.row(i);
+
+ Glyph *glyph = new Glyph();
+
+ glyph->setSize(5);
+ glyph->setX((row[0] - xmin) / (xmax - xmin) * width());
+ glyph->setY((row[1] - ymin) / (ymax - ymin) * height());
+
+ glyph->setParent(this);
+ }
+
+ update();
+}
+
+QSGNode *Scatterplot::updatePaintNode(QSGNode *oldNode, UpdatePaintNodeData *)
+{
+ QSGNode *node = 0;
+
+ if (!oldNode) {
+ node = new QSGNode;
+ for (QObjectList::const_iterator it = children().begin(); it != children().end(); it++)
+ node->appendChildNode(static_cast<Glyph *>(*it)->updatePaintNode(0, 0));
+ } else {
+ node = static_cast<QSGNode *>(oldNode);
+ }
+
+ node->markDirty(QSGNode::DirtyGeometry);
+
+ return node;
+}
diff --git a/scatterplot.h b/scatterplot.h
new file mode 100644
index 0000000..f38ce78
--- /dev/null
+++ b/scatterplot.h
@@ -0,0 +1,26 @@
+#ifndef SCATTERPLOT_H
+#define SCATTERPLOT_H
+
+#include <armadillo>
+#include <vector>
+#include <QQuickItem>
+
+class Scatterplot : public QQuickItem
+{
+ Q_OBJECT
+public:
+ Scatterplot();
+ ~Scatterplot();
+
+ QSGNode *updatePaintNode(QSGNode *oldNode, UpdatePaintNodeData *);
+ void setData(const arma::mat &data);
+
+signals:
+
+public slots:
+
+private:
+ arma::mat m_data;
+};
+
+#endif // SCATTERPLOT_H