meshfactory.cpp

Go to the documentation of this file.

//=============================================================================================================
 
//=============================================================================================================
// INCLUDES
//=============================================================================================================
 
#include "meshfactory.h"
#include "renderable/brainsurface.h"
#include <Eigen/Geometry>
#include <QMap>
#include <QPair>
#include <cmath>
 
//=============================================================================================================
// STATIC HELPERS
//=============================================================================================================
 
namespace {
 
Eigen::Matrix3f extractRotation(const QMatrix4x4 &m)
{
    Eigen::Matrix3f rot;
    for (int r = 0; r < 3; ++r)
        for (int c = 0; c < 3; ++c)
            rot(r, c) = m(r, c);
    return rot;
}
 
} // anonymous namespace
 
//=============================================================================================================
// DEFINE MEMBER METHODS
//=============================================================================================================
 
void MeshFactory::buildIcosphere(QVector<Eigen::Vector3f> &vertices,
                                  QVector<Eigen::Vector3i> &faces,
                                  int subdivisions)
{
    // Golden ratio
    const float phi = (1.0f + std::sqrt(5.0f)) / 2.0f;
 
    vertices.clear();
    vertices.reserve(12 + 30 * subdivisions); // rough estimate
    vertices << Eigen::Vector3f(-1,  phi, 0) << Eigen::Vector3f( 1,  phi, 0)
             << Eigen::Vector3f(-1, -phi, 0) << Eigen::Vector3f( 1, -phi, 0)
             << Eigen::Vector3f( 0, -1,  phi) << Eigen::Vector3f( 0,  1,  phi)
             << Eigen::Vector3f( 0, -1, -phi) << Eigen::Vector3f( 0,  1, -phi)
             << Eigen::Vector3f( phi, 0, -1) << Eigen::Vector3f( phi, 0,  1)
             << Eigen::Vector3f(-phi, 0, -1) << Eigen::Vector3f(-phi, 0,  1);
    for (auto &v : vertices)
        v.normalize();
 
    faces.clear();
    faces.reserve(20);
    faces << Eigen::Vector3i(0,11,5)  << Eigen::Vector3i(0,5,1)
          << Eigen::Vector3i(0,1,7)   << Eigen::Vector3i(0,7,10)
          << Eigen::Vector3i(0,10,11) << Eigen::Vector3i(1,5,9)
          << Eigen::Vector3i(5,11,4)  << Eigen::Vector3i(11,10,2)
          << Eigen::Vector3i(10,7,6)  << Eigen::Vector3i(7,1,8)
          << Eigen::Vector3i(3,9,4)   << Eigen::Vector3i(3,4,2)
          << Eigen::Vector3i(3,2,6)   << Eigen::Vector3i(3,6,8)
          << Eigen::Vector3i(3,8,9)   << Eigen::Vector3i(4,9,5)
          << Eigen::Vector3i(2,4,11)  << Eigen::Vector3i(6,2,10)
          << Eigen::Vector3i(8,6,7)   << Eigen::Vector3i(9,8,1);
 
    // Subdivide
    for (int s = 0; s < subdivisions; ++s) {
        QMap<QPair<int,int>, int> cache;
        auto midpoint = [&](int a, int b) -> int {
            auto key = qMakePair(qMin(a, b), qMax(a, b));
            if (cache.contains(key))
                return cache[key];
            int idx = vertices.size();
            vertices.append((vertices[a] + vertices[b]).normalized());
            cache[key] = idx;
            return idx;
        };
 
        QVector<Eigen::Vector3i> newFaces;
        newFaces.reserve(faces.size() * 4);
        for (const auto &f : faces) {
            int ab = midpoint(f(0), f(1));
            int bc = midpoint(f(1), f(2));
            int ca = midpoint(f(2), f(0));
            newFaces << Eigen::Vector3i(f(0), ab, ca)
                     << Eigen::Vector3i(f(1), bc, ab)
                     << Eigen::Vector3i(f(2), ca, bc)
                     << Eigen::Vector3i(ab,   bc, ca);
        }
        faces = std::move(newFaces);
    }
}
 
//=============================================================================================================
 
int MeshFactory::sphereVertexCount(int subdivisions)
{
    QVector<Eigen::Vector3f> verts;
    QVector<Eigen::Vector3i> faces;
    buildIcosphere(verts, faces, subdivisions);
    return verts.size();
}
 
//=============================================================================================================
 
std::shared_ptr<BrainSurface> MeshFactory::createSphere(const QVector3D &center,
                                                          float radius,
                                                          const QColor &color,
                                                          int subdivisions)
{
    QVector<Eigen::Vector3f> verts;
    QVector<Eigen::Vector3i> faces;
    buildIcosphere(verts, faces, subdivisions);
 
    const int nV = verts.size();
    const int nT = faces.size();
 
    Eigen::MatrixX3f rr(nV, 3), nn(nV, 3);
    for (int i = 0; i < nV; ++i) {
        nn(i, 0) = verts[i].x();
        nn(i, 1) = verts[i].y();
        nn(i, 2) = verts[i].z();
        rr(i, 0) = verts[i].x() * radius + center.x();
        rr(i, 1) = verts[i].y() * radius + center.y();
        rr(i, 2) = verts[i].z() * radius + center.z();
    }
 
    Eigen::MatrixX3i tris(nT, 3);
    for (int i = 0; i < nT; ++i)
        tris.row(i) = faces[i];
 
    auto surf = std::make_shared<BrainSurface>();
    surf->createFromData(rr, nn, tris, color);
    return surf;
}
 
//=============================================================================================================
 
std::shared_ptr<BrainSurface> MeshFactory::createPlate(const QVector3D &center,
                                                         const QMatrix4x4 &orientation,
                                                         const QColor &color,
                                                         float size)
{
    const float hw = size / 2.0f;
    const float hh = size / 2.0f;
    const float hd = size * 0.05f; // thin plate
 
    const Eigen::Matrix3f rot = extractRotation(orientation);
 
    Eigen::Vector3f corners[8] = {
        {-hw, -hh, -hd}, { hw, -hh, -hd}, { hw,  hh, -hd}, {-hw,  hh, -hd},
        {-hw, -hh,  hd}, { hw, -hh,  hd}, { hw,  hh,  hd}, {-hw,  hh,  hd}
    };
    for (auto &c : corners)
        c = rot * c;
 
    const int faceIndices[6][4] = {
        {4,5,6,7}, {1,0,3,2},
        {3,7,6,2}, {0,1,5,4},
        {1,2,6,5}, {0,4,7,3}
    };
 
    Eigen::MatrixX3f rr(24, 3), nn(24, 3);
    Eigen::MatrixX3i tris(12, 3);
 
    for (int f = 0; f < 6; ++f) {
        Eigen::Vector3f v0 = corners[faceIndices[f][0]];
        Eigen::Vector3f v1 = corners[faceIndices[f][1]];
        Eigen::Vector3f v2 = corners[faceIndices[f][2]];
        Eigen::Vector3f e1 = v1 - v0;
        Eigen::Vector3f e2 = v2 - v0;
        Eigen::Vector3f fn = e1.cross(e2).normalized();
        int base = f * 4;
        for (int k = 0; k < 4; ++k) {
            Eigen::Vector3f v = corners[faceIndices[f][k]];
            rr.row(base + k) << v.x() + center.x(), v.y() + center.y(), v.z() + center.z();
            nn.row(base + k) = fn;
        }
        tris.row(f * 2)     << base, base + 1, base + 2;
        tris.row(f * 2 + 1) << base, base + 2, base + 3;
    }
 
    auto surf = std::make_shared<BrainSurface>();
    surf->createFromData(rr, nn, tris, color);
    return surf;
}
 
//=============================================================================================================
 
std::shared_ptr<BrainSurface> MeshFactory::createBarbell(const QVector3D &center,
                                                           const QMatrix4x4 &orientation,
                                                           const QColor &color,
                                                           float size)
{
    const Eigen::Matrix3f rot = extractRotation(orientation);
 
    const float halfSpan = size * 0.6f;
    const float sphereR  = size * 0.2f;
    const float rodR     = size * 0.08f;
 
    QVector<Eigen::Vector3f> allVerts;
    QVector<Eigen::Vector3f> allNorms;
    QVector<Eigen::Vector3i> allTris;
 
    // ── Helper: append a sphere ────────────────────────────────────────
    auto appendSphere = [&](const QVector3D &pos, float radius) {
        QVector<Eigen::Vector3f> sv;
        QVector<Eigen::Vector3i> sf;
        buildIcosphere(sv, sf, 1); // single subdivision
 
        const int base = allVerts.size();
        allVerts.reserve(base + sv.size());
        allNorms.reserve(base + sv.size());
        allTris.reserve(allTris.size() + sf.size());
 
        for (const auto &v : sv) {
            Eigen::Vector3f vn = rot * v;
            Eigen::Vector3f vp = rot * (v * radius);
            allVerts.append(Eigen::Vector3f(vp.x() + pos.x(), vp.y() + pos.y(), vp.z() + pos.z()));
            allNorms.append(vn.normalized());
        }
        for (const auto &f : sf) {
            allTris.append(Eigen::Vector3i(base + f(0), base + f(1), base + f(2)));
        }
    };
 
    // ── Helper: append a box-rod ───────────────────────────────────────
    auto appendRod = [&]() {
        float hx = rodR, hy = halfSpan, hz = rodR;
        Eigen::Vector3f corners[8] = {
            {-hx, -hy, -hz}, { hx, -hy, -hz}, { hx,  hy, -hz}, {-hx,  hy, -hz},
            {-hx, -hy,  hz}, { hx, -hy,  hz}, { hx,  hy,  hz}, {-hx,  hy,  hz}
        };
        for (auto &c : corners) c = rot * c;
 
        const int faceIdx[6][4] = {
            {4,5,6,7}, {1,0,3,2}, {3,7,6,2}, {0,1,5,4}, {1,2,6,5}, {0,4,7,3}
        };
 
        const int base = allVerts.size();
        for (int f = 0; f < 6; ++f) {
            Eigen::Vector3f v0 = corners[faceIdx[f][0]];
            Eigen::Vector3f v1 = corners[faceIdx[f][1]];
            Eigen::Vector3f v2 = corners[faceIdx[f][2]];
            Eigen::Vector3f e1 = v1 - v0;
            Eigen::Vector3f e2 = v2 - v0;
            Eigen::Vector3f fn = e1.cross(e2).normalized();
            int faceBase = base + f * 4;
            for (int k = 0; k < 4; ++k) {
                Eigen::Vector3f v = corners[faceIdx[f][k]];
                allVerts.append(Eigen::Vector3f(v.x() + center.x(), v.y() + center.y(), v.z() + center.z()));
                allNorms.append(fn);
            }
            allTris.append(Eigen::Vector3i(faceBase, faceBase + 1, faceBase + 2));
            allTris.append(Eigen::Vector3i(faceBase, faceBase + 2, faceBase + 3));
        }
    };
 
    // Build the barbell
    QVector3D axis(rot(0,1), rot(1,1), rot(2,1));
    appendSphere(center + axis * halfSpan, sphereR);
    appendSphere(center - axis * halfSpan, sphereR);
    appendRod();
 
    // Pack into Eigen matrices
    Eigen::MatrixX3f rr(allVerts.size(), 3), nn(allNorms.size(), 3);
    Eigen::MatrixX3i tt(allTris.size(), 3);
    for (int i = 0; i < allVerts.size(); ++i) {
        rr.row(i) = allVerts[i].transpose();
        nn.row(i) = allNorms[i].transpose();
    }
    for (int i = 0; i < allTris.size(); ++i)
        tt.row(i) = allTris[i];
 
    auto surf = std::make_shared<BrainSurface>();
    surf->createFromData(rr, nn, tt, color);
    return surf;
}
 
//=============================================================================================================
 
std::shared_ptr<BrainSurface> MeshFactory::createCylinder(const QVector3D &from,
                                                          const QVector3D &to,
                                                          float radius,
                                                          const QColor &color,
                                                          int sides)
{
    if (sides < 3) sides = 3;
 
    const QVector3D axis = to - from;
    const float length = axis.length();
    if (length < 1e-8f)
        return std::make_shared<BrainSurface>();
 
    const QVector3D axisN = axis / length;
 
    // Build an orthonormal basis (axisN, u, v)
    QVector3D u = QVector3D::crossProduct(axisN, QVector3D(0, 0, 1));
    if (u.lengthSquared() < 1e-6f)
        u = QVector3D::crossProduct(axisN, QVector3D(1, 0, 0));
    u.normalize();
    const QVector3D v = QVector3D::crossProduct(axisN, u).normalized();
 
    // Two rings of vertices (bottom = from, top = to) + 2 center caps
    const int nVerts = 2 * sides + 2;
    const int nTris  = 4 * sides; // 2 * sides for barrel + sides for each cap
 
    Eigen::MatrixX3f verts(nVerts, 3);
    Eigen::MatrixX3f norms(nVerts, 3);
    Eigen::MatrixX3i tris(nTris, 3);
 
    constexpr float kTwoPi = 2.0f * 3.14159265358979f;
 
    // Bottom ring vertices [0..sides-1]
    for (int i = 0; i < sides; ++i) {
        const float angle = kTwoPi * static_cast<float>(i) / static_cast<float>(sides);
        const float ca = std::cos(angle);
        const float sa = std::sin(angle);
        const QVector3D normal = (u * ca + v * sa);
        const QVector3D pos = from + normal * radius;
        verts(i, 0) = pos.x();
        verts(i, 1) = pos.y();
        verts(i, 2) = pos.z();
        norms(i, 0) = normal.x();
        norms(i, 1) = normal.y();
        norms(i, 2) = normal.z();
    }
 
    // Top ring vertices [sides..2*sides-1]
    for (int i = 0; i < sides; ++i) {
        const float angle = kTwoPi * static_cast<float>(i) / static_cast<float>(sides);
        const float ca = std::cos(angle);
        const float sa = std::sin(angle);
        const QVector3D normal = (u * ca + v * sa);
        const QVector3D pos = to + normal * radius;
        verts(sides + i, 0) = pos.x();
        verts(sides + i, 1) = pos.y();
        verts(sides + i, 2) = pos.z();
        norms(sides + i, 0) = normal.x();
        norms(sides + i, 1) = normal.y();
        norms(sides + i, 2) = normal.z();
    }
 
    // Cap center vertices
    const int bottomCenter = 2 * sides;
    const int topCenter    = 2 * sides + 1;
    verts(bottomCenter, 0) = from.x(); verts(bottomCenter, 1) = from.y(); verts(bottomCenter, 2) = from.z();
    norms(bottomCenter, 0) = -axisN.x(); norms(bottomCenter, 1) = -axisN.y(); norms(bottomCenter, 2) = -axisN.z();
    verts(topCenter, 0) = to.x(); verts(topCenter, 1) = to.y(); verts(topCenter, 2) = to.z();
    norms(topCenter, 0) = axisN.x(); norms(topCenter, 1) = axisN.y(); norms(topCenter, 2) = axisN.z();
 
    // Barrel triangles (2 tris per side)
    int t = 0;
    for (int i = 0; i < sides; ++i) {
        const int i0 = i;
        const int i1 = (i + 1) % sides;
        const int i2 = sides + i;
        const int i3 = sides + (i + 1) % sides;
        tris(t, 0) = i0; tris(t, 1) = i2; tris(t, 2) = i1; ++t;
        tris(t, 0) = i1; tris(t, 1) = i2; tris(t, 2) = i3; ++t;
    }
 
    // Bottom cap
    for (int i = 0; i < sides; ++i) {
        tris(t, 0) = bottomCenter;
        tris(t, 1) = (i + 1) % sides;
        tris(t, 2) = i;
        ++t;
    }
 
    // Top cap
    for (int i = 0; i < sides; ++i) {
        tris(t, 0) = topCenter;
        tris(t, 1) = sides + i;
        tris(t, 2) = sides + (i + 1) % sides;
        ++t;
    }
 
    auto surf = std::make_shared<BrainSurface>();
    surf->createFromData(verts, norms, tris, color);
    return surf;
}
 
//=============================================================================================================
 
std::shared_ptr<BrainSurface> MeshFactory::createBatchedSpheres(const QVector<QVector3D> &positions,
                                                                  float radius,
                                                                  const QColor &color,
                                                                  int subdivisions)
{
    if (positions.isEmpty())
        return std::make_shared<BrainSurface>();
 
    // Build template sphere once
    QVector<Eigen::Vector3f> tmplVerts;
    QVector<Eigen::Vector3i> tmplFaces;
    buildIcosphere(tmplVerts, tmplFaces, subdivisions);
 
    const int nVPerPt = tmplVerts.size();
    const int nTPerPt = tmplFaces.size();
    const int nPts = positions.size();
 
    // Scale template
    Eigen::MatrixX3f templateV(nVPerPt, 3), templateN(nVPerPt, 3);
    for (int i = 0; i < nVPerPt; ++i) {
        templateN(i, 0) = tmplVerts[i].x();
        templateN(i, 1) = tmplVerts[i].y();
        templateN(i, 2) = tmplVerts[i].z();
        templateV(i, 0) = tmplVerts[i].x() * radius;
        templateV(i, 1) = tmplVerts[i].y() * radius;
        templateV(i, 2) = tmplVerts[i].z() * radius;
    }
    Eigen::MatrixX3i templateT(nTPerPt, 3);
    for (int i = 0; i < nTPerPt; ++i)
        templateT.row(i) = tmplFaces[i];
 
    // Merge all spheres into single mesh
    Eigen::MatrixX3f allVerts(nPts * nVPerPt, 3);
    Eigen::MatrixX3f allNorms(nPts * nVPerPt, 3);
    Eigen::MatrixX3i allTris(nPts * nTPerPt, 3);
 
    for (int p = 0; p < nPts; ++p) {
        const QVector3D &pos = positions[p];
        const int vOff = p * nVPerPt;
        const int tOff = p * nTPerPt;
        for (int iv = 0; iv < nVPerPt; ++iv) {
            allVerts(vOff + iv, 0) = templateV(iv, 0) + pos.x();
            allVerts(vOff + iv, 1) = templateV(iv, 1) + pos.y();
            allVerts(vOff + iv, 2) = templateV(iv, 2) + pos.z();
            allNorms.row(vOff + iv) = templateN.row(iv);
        }
        for (int it = 0; it < nTPerPt; ++it) {
            allTris(tOff + it, 0) = templateT(it, 0) + vOff;
            allTris(tOff + it, 1) = templateT(it, 1) + vOff;
            allTris(tOff + it, 2) = templateT(it, 2) + vOff;
        }
    }
 
    auto surf = std::make_shared<BrainSurface>();
    surf->createFromData(allVerts, allNorms, allTris, color);
    return surf;
}