#include <memory>

#include "logcategories.h"
#include "subdivision.h"

Subdivision::Subdivision(QOpenGLFunctions_4_3_Core *f)
{
    this->f =f;
}

Subdivision::~Subdivision()
{
    delete neighborsShader;
    delete copyShader;
    delete edgeShader;
    delete vertexShader;
}

QString Subdivision::formatTimeMeasurement(int time){
   return QString("%1m %2s %3ms").arg( time / 60000         , 2, 10, QChar('0'))
                                  .arg((time % 60000) / 1000, 2, 10, QChar('0'))
                                   .arg((time % 1000)       , 3, 10, QChar('0'));
}

void Subdivision::init() {
    QString source = QLatin1String(":/subdivision-neighbors.compute");
    neighborsShader = initComputeShaderProgram(source);

    source = QLatin1String(":/subdivision-copy.compute");
    copyShader = initComputeShaderProgram(source);

    source = QLatin1String(":/subdivision-edge.compute");
    edgeShader = initComputeShaderProgram(source);

    source = QLatin1String(":/subdivision-vertex.compute");
    vertexShader = initComputeShaderProgram(source);
}

QOpenGLShaderProgram *Subdivision::initComputeShaderProgram(QString &source){
    qCDebug(log_subdiv) << "Compiling compute shader ...";
    QOpenGLShader *computeShader = new QOpenGLShader(QOpenGLShader::Compute);
    if(!computeShader->compileSourceFile(source)){
        qCCritical(log_subdiv) << "Compute Shader" << source << "failed" << computeShader->log();
        exit(5);
    }

    qCDebug(log_subdiv) << "Adding compute shader ...";
    QOpenGLShaderProgram *shader = new QOpenGLShaderProgram();
    shader->addShader(computeShader);
    qCDebug(log_subdiv) << "Linking compute shader ...";
    if(!shader->link()){
        qCCritical(log_subdiv) << "Linking compute shader failed:" << shader->log();
        exit(5);
    }
    qCDebug(log_subdiv) << "Linking compute shader done";

    return shader;
}

void Subdivision::subdivide(Mesh *mesh, int level) {
    // For now we only look at the first mesh entry
    Mesh::Node root = mesh->getRootNode();

    int first_mesh_index = -1;
    if (!root.getFirstMeshIndex(first_mesh_index)) {
        qCCritical(log_subdiv) << "No mesh found, aborting subdivision";
        return;
    }


    Mesh::MeshEntry *current_mesh = mesh->getMeshEntry(first_mesh_index);
    Input input;

    int currentMax = current_mesh->buffers.size()-1;
    //qCDebug(log_subdiv) << "Subdiv Level" << level << "Requested. Current Max:" << currentMax;
    if(level <= currentMax)
        return;

    std::shared_ptr<Mesh::SubdivEntry> entry = current_mesh->buffers[currentMax];

    input.vb_handle = entry->VB_handle;
    input.vertex_buffer = entry->vertices;
    input.index_irregular_buffer = entry->indices_irregular;
    input.index_regular_buffer = patchIBToTriangleIB(entry->indices_regular);

    if (input.index_irregular_buffer.isEmpty()) {
        //  current_mesh->update(input.vb_handle, input.vertex_buffer, input.index_buffer, entry->indices_regular);
    } else {
        Tables tables = precomputeTables(input);
        Result result = runShader(input, tables);

        QVector<Triangle> triangles;
        ibToTriangles(&result.vertex_buffer, tables.index_buffer, triangles);
        QVector<Triangle> all_triangles;
        ibToTriangles(&result.vertex_buffer, tables.extra_triangles, all_triangles);
        all_triangles += triangles;

        QMap<Triangle, Triangle::Neighbors> neighbors;
        buildNeighborsMap(result.vertex_buffer, all_triangles, neighbors);

        QVector<Triangle> regular;
        QVector<Triangle> irregular;
        findRegular(triangles, neighbors, regular, irregular);

        qCDebug(log_subdiv) << "Indices" << tables.index_buffer.length();
        qCDebug(log_subdiv) << "subdivide: regular" << regular.length();
        qCDebug(log_subdiv) << "subdivide: irregular" << irregular.length();

        QVector<unsigned int> irregular_ib;
        trianglesToIB(irregular, irregular_ib);

        updateIrregularForDraw(irregular, neighbors, result);

        QVector<unsigned int> patches;
        getPatchIndexBuffer(regular, neighbors, patches);
        qCDebug(log_subdiv) << "patches" << patches.length();

        current_mesh->update(result.vb_handle, result.vertex_buffer, result.vertex_buffer_irregular, irregular_ib, patches);
    }
}

struct Edge {
    unsigned int a;
    unsigned int b;

    friend bool operator<(const Edge &l, const Edge &r) {
        if (l.a < r.a) {
            return true;
        } else if (l.a > r.a) {
            return false;
        } else {
            return l.b < r.b;
        }
    }
};

/**
 * @brief Subdivision::precomputeTables
 * @param mesh
 * @param edgeIndices_base Wird gefüllt mit 4er-Arrays, die die indices in den alten Vertexbuffer enthalten, von denen der jeweilige edge-point abhängt
 * @param vertexIndices Wird gefüllt mit Folgen von ints, die die indices in den alten Vertexbuffer enthalten, von denen die neuen vertex-points abhängen
 * @param vertexIndicesOffsets enthält offsets in vertexIndices
 * @param newIndexBuffer Wird gefüllt mit neuem Indexbuffer. Der Indexbuffer geht davon aus, dass im neuen Vertex-Buffer
 * erst die Vertex-points (in der selben Reihenfolge wie im ursprünglichen Vertex-Buffer) und dann die Edgepoints kommen.
 *
 * Einträge in neuem Vertex-buffer: Verschobene Knoten aus altem Vertex-Buffer :: edge-vertices erstellt aus edgeIndices_base.
 */
Subdivision::Tables Subdivision::precomputeTables(Input input) {
    Tables tables;

    QVector<unsigned int> ib = input.index_irregular_buffer;
    qCDebug(log_subdiv_trace) << "Index Buffer Irregular: " << ib;

    QVector<unsigned int> ib_regular = input.index_regular_buffer;
    qCDebug(log_subdiv_trace) << "Index Buffer Regular: " << ib_regular;

    QVector<Vertex> vb = input.vertex_buffer;
    qCDebug(log_subdiv_trace) << "Vertex Buffer: " << vb;

    QTime subTimer;
    subTimer.start();
    QVector<Triangle> triangles;
    ibToTriangles(&vb, ib, triangles);

    QVector<Triangle> triangles_regular;
    ibToTriangles(&vb, ib_regular, triangles_regular);

    qCDebug(log_timing) << "building Triangles:" << formatTimeMeasurement(subTimer.elapsed());

    QMap<Triangle, Triangle::Neighbors> neighbors;
    QVector<Triangle> all_triangles = triangles + triangles_regular;
    buildNeighborsMap(vb, all_triangles, neighbors);

    subTimer.restart();
    precomputeEdgeTable(tables, triangles, neighbors, (unsigned int) vb.length());
    qCDebug(log_timing) << "precomputeEdgeTable:" << formatTimeMeasurement(subTimer.elapsed());

    QMap<unsigned int, unsigned int> modified_vertices;
    subTimer.restart();
    precomputeVertexTable(tables, input, modified_vertices);
    qCDebug(log_timing) << "precomputeVertexTable:" << formatTimeMeasurement(subTimer.elapsed());

    subTimer.restart();
    updateIndexBuffer(tables.index_buffer, modified_vertices);
    updateIndexBuffer(tables.extra_triangles, modified_vertices);
    qCDebug(log_timing) << "updateIndexBuffer:" << formatTimeMeasurement(subTimer.elapsed());

    qCDebug(log_subdiv) << "Precompute Tables Done";
    return tables;
}

void Subdivision::buildNeighborsMap(QVector<Vertex> &vb, QVector<Triangle> &triangles, QMap<Triangle, Triangle::Neighbors> &neighbors) {
    QTime timer;
    timer.start();
    qCDebug(log_timing) << "buildNeighborsMap_gpu started";

    QTime subTimer;
    subTimer.start();

    GLuint vb_handle;
    f->glGenBuffers(1, &vb_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, vb_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, vb.size() * sizeof(Vertex), vb.data(), GL_DYNAMIC_DRAW);

    GLuint ib_handle;
    f->glGenBuffers(1, &ib_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, ib_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, triangles.size() * 3 * sizeof(GLuint), NULL, GL_DYNAMIC_DRAW);

    GLuint *ptr = (GLuint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
    for (int i = 0; i < triangles.size(); i++) {
        ptr[3 * i + 0] = triangles[i].u_idx();
        ptr[3 * i + 1] = triangles[i].v_idx();
        ptr[3 * i + 2] = triangles[i].w_idx();
    }
    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);

    GLuint offset_handle;
    f->glGenBuffers(1, &offset_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, offset_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, 2 * sizeof(GLuint), NULL, GL_DYNAMIC_DRAW);

    GLuint output_handle;
    f->glGenBuffers(1, &output_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, output_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, triangles.size() * 3 * sizeof(GLint), NULL, GL_DYNAMIC_DRAW);

    GLint *output_ptr = (GLint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
    memset(output_ptr, (GLint) -1, triangles.size() * 3 * sizeof(GLint));
    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);

    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, 0);

    qCDebug(log_timing) << "Copying to GPU:" << formatTimeMeasurement(subTimer.elapsed());

    subTimer.restart();
    neighborsShader->bind();

    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, vb_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, ib_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, offset_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, output_handle);

    GLuint offset_in = 0;
    GLuint offset_search = 0;
    while (offset_in <= (unsigned) triangles.size() || offset_search <= (unsigned) triangles.size()) {
        if (offset_search > (unsigned) triangles.size()) {
            offset_in += 128;
            offset_search = 0;
        }

        f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, offset_handle);
        GLuint *offset_ptr = (GLuint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
        offset_ptr[0] = offset_in;
        offset_ptr[1] = offset_search;
        f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);
        f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, 0);

        f->glDispatchCompute(128, 32, 1);

        offset_search += 128 * 32;
    }

    f->glMemoryBarrier(GL_SHADER_STORAGE_BARRIER_BIT);

    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, 0);

    neighborsShader->release();

    qCDebug(log_timing) << "Running shader:" << formatTimeMeasurement(subTimer.elapsed());

    subTimer.restart();
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, output_handle);
    GLint *out_ptr = (GLint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_READ_ONLY);
    for (int i = 0; i < triangles.length(); i++) {
        Triangle::Neighbors ns;

        GLint other = out_ptr[3 * i + 0];
        assert(other < triangles.length());
        if (other < 0) {
            ns.uv.triangle = NULL;
        } else {
            ns.uv.triangle = &(triangles[other]);
            if (out_ptr[3 * other + 0] == i) {
                ns.uv.edge.name = Triangle::Edge::Name::uv;
                ns.uv.edge.a = ns.uv.triangle->u_idx();
                ns.uv.edge.b = ns.uv.triangle->v_idx();
                ns.uv.edge.c = ns.uv.triangle->w_idx();
            } else if (out_ptr[3 * other + 1] == i) {
                ns.uv.edge.name = Triangle::Edge::Name::vw;
                ns.uv.edge.a = ns.uv.triangle->v_idx();
                ns.uv.edge.b = ns.uv.triangle->w_idx();
                ns.uv.edge.c = ns.uv.triangle->u_idx();
            } else if (out_ptr[3 * other + 2] == i) {
                ns.uv.edge.name = Triangle::Edge::Name::wu;
                ns.uv.edge.a = ns.uv.triangle->w_idx();
                ns.uv.edge.b = ns.uv.triangle->u_idx();
                ns.uv.edge.c = ns.uv.triangle->v_idx();
            }
        }

        other = out_ptr[3 * i + 1];
        assert(other < triangles.length());
        if (other < 0) {
            ns.vw.triangle = NULL;
        } else {
            ns.vw.triangle = &(triangles[other]);
            if (out_ptr[3 * other + 0] == i) {
                ns.vw.edge.name = Triangle::Edge::Name::uv;
                ns.vw.edge.a = ns.vw.triangle->u_idx();
                ns.vw.edge.b = ns.vw.triangle->v_idx();
                ns.vw.edge.c = ns.vw.triangle->w_idx();
            } else if (out_ptr[3 * other + 1] == i) {
                ns.vw.edge.name = Triangle::Edge::Name::vw;
                ns.vw.edge.a = ns.vw.triangle->v_idx();
                ns.vw.edge.b = ns.vw.triangle->w_idx();
                ns.vw.edge.c = ns.vw.triangle->u_idx();
            } else if (out_ptr[3 * other + 2] == i) {
                ns.vw.edge.name = Triangle::Edge::Name::wu;
                ns.vw.edge.a = ns.vw.triangle->w_idx();
                ns.vw.edge.b = ns.vw.triangle->u_idx();
                ns.vw.edge.c = ns.vw.triangle->v_idx();
            }
        }

        other = out_ptr[3 * i + 2];
        assert(other < triangles.length());
        if (other < 0) {
            ns.wu.triangle = NULL;
        } else {
            ns.wu.triangle = &(triangles[other]);
            if (out_ptr[3 * other + 0] == i) {
                ns.wu.edge.name = Triangle::Edge::Name::uv;
                ns.wu.edge.a = ns.wu.triangle->u_idx();
                ns.wu.edge.b = ns.wu.triangle->v_idx();
                ns.wu.edge.c = ns.wu.triangle->w_idx();
            } else if (out_ptr[3 * other + 1] == i) {
                ns.wu.edge.name = Triangle::Edge::Name::vw;
                ns.wu.edge.a = ns.wu.triangle->v_idx();
                ns.wu.edge.b = ns.wu.triangle->w_idx();
                ns.wu.edge.c = ns.wu.triangle->u_idx();
            } else if (out_ptr[3 * other + 2] == i) {
                ns.wu.edge.name = Triangle::Edge::Name::wu;
                ns.wu.edge.a = ns.wu.triangle->w_idx();
                ns.wu.edge.b = ns.wu.triangle->u_idx();
                ns.wu.edge.c = ns.wu.triangle->v_idx();
            }
        }

        neighbors.insert(triangles[i], ns);
    }
    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);

    qCDebug(log_timing) << "Building map:" << formatTimeMeasurement(subTimer.elapsed());

    f->glDeleteBuffers(1, &vb_handle);
    f->glDeleteBuffers(1, &ib_handle);
    f->glDeleteBuffers(1, &output_handle);

    qCDebug(log_timing) << "buildNeighborsMap_gpu done:" << formatTimeMeasurement(timer.elapsed());
}

/**
 * Returns triangles that use a vertex except for the starting triangle.
 *
 * @param triangle the starting triangle. Will not be returned.
 * @param current_neighbor a neighbor of the starting triangle.
 * @param neighbors neighbors of triangles
 * @return all neighbors of the first vertex of the neighbor except for the starting triangle.
 */
QVector<Triangle> Subdivision::vertexNeighbors(Triangle &triangle, Triangle::Neighbor current_neighbor, QMap<Triangle, Triangle::Neighbors> neighbors) {
    QVector<Triangle> ns;
    do {
        assert(current_neighbor.triangle != NULL);
        ns.push_back(*current_neighbor.triangle);
        current_neighbor = neighbors.value(*current_neighbor.triangle).get_neighbor(rotate_edge_name(current_neighbor.edge.name));
    } while (*current_neighbor.triangle != triangle);
    assert(current_neighbor.triangle != NULL);
    return ns;
}

void Subdivision::precomputeEdgeTable(Subdivision::Tables &tables, QVector<Triangle> &triangles, QMap<Triangle, Triangle::Neighbors> &neighbors, unsigned int offset) {
    //compute edge table
    //Format: first two entries: edge vertices. last two entries: distant vertices.

    // Maps edges to the index of the new vertex added on that edge. The keys are the indices of the two vertices defining the edge and must be in the order uv, vw or uw.
    QMap<Edge, unsigned int> edge_indices;

    QMap<Triangle, bool> irregular_triangles;
    QVectorIterator<Triangle> it(triangles);
    while (it.hasNext()) {
        irregular_triangles.insert(it.next(), true);
    }

    // Regular neighbors of irregular triangles
    // We have to generate edge vertices without putting them in the index buffer for rendering patches.
    QMap<Triangle, bool> extra_triangles;

    for (int i = 0; i < triangles.length(); i++){
        //schaue alle dreiecke an
        Triangle *triangle = &triangles[i];

        QVector<unsigned int> edge_indices_buffer;//push all edge indices into this, then check if enough edge indices were found. if yes, copy to Tables and add to new index buffer.

        /*
         * Find edge xy in all other triangles, add to edge_indices
         * Find edge yz
         * Find edge zx
         * - We assume that we have a closed surface, i.e. each triangle has exactly one adjacent triangle at each edge
         * - Each edge point will be computed twice. This allows us to implicitly know the adjacent points.
         */

        Triangle::Neighbors ns = neighbors.value(*triangle);
        QVector<Triangle> neighboring_triangles;
        neighboring_triangles += vertexNeighbors(*triangle, ns.uv, neighbors);
        neighboring_triangles += vertexNeighbors(*triangle, ns.vw, neighbors);
        neighboring_triangles += vertexNeighbors(*triangle, ns.wu, neighbors);
        QVectorIterator<Triangle> extra(neighboring_triangles);
        while (extra.hasNext()) {
            Triangle extra_triangle = extra.next();
            if (!irregular_triangles.contains(extra_triangle)) {
                extra_triangles.insert(extra_triangle, true);
            }
        }

        // indices of the three vertices added to the edges of this triangle
        unsigned int uv, vw, wu;

        Edge edge;

        edge = { triangle->u_idx(), triangle->v_idx() };
        if (edge_indices.contains(edge)) {
            uv = edge_indices.value(edge);
        } else {
            uv = offset++;
            edge_indices.insert(edge, uv);

            edge = { ns.uv.edge.a, ns.uv.edge.b };
            edge_indices.insert(edge, uv);

            tables.edge_indices.push_back(triangle->u_idx());
            tables.edge_indices.push_back(triangle->v_idx());
            tables.edge_indices.push_back(triangle->w_idx());
            tables.edge_indices.push_back(ns.uv.edge.c);
        }

        edge = { triangle->v_idx(), triangle->w_idx() };
        if (edge_indices.contains(edge)) {
            vw = edge_indices.value(edge);
        } else {
            vw = offset++;
            edge_indices.insert(edge, vw);

            edge = { ns.vw.edge.a, ns.vw.edge.b };
            edge_indices.insert(edge, vw);

            tables.edge_indices.push_back(triangle->v_idx());
            tables.edge_indices.push_back(triangle->w_idx());
            tables.edge_indices.push_back(triangle->u_idx());
            tables.edge_indices.push_back(ns.vw.edge.c);
        }

        edge = { triangle->w_idx(), triangle->u_idx() };
        if (edge_indices.contains(edge)) {
            wu = edge_indices.value(edge);
        } else {
            wu = offset++;
            edge_indices.insert(edge, wu);

            edge = { ns.wu.edge.a, ns.wu.edge.b };
            edge_indices.insert(edge, wu);

            tables.edge_indices.push_back(triangle->w_idx());
            tables.edge_indices.push_back(triangle->u_idx());
            tables.edge_indices.push_back(triangle->v_idx());
            tables.edge_indices.push_back(ns.wu.edge.c);
        }

        //add indices to new index buffer
        tables.index_buffer.push_back(triangle->u_idx());
        tables.index_buffer.push_back(uv);
        tables.index_buffer.push_back(wu);

        tables.index_buffer.push_back(triangle->v_idx());
        tables.index_buffer.push_back(vw);
        tables.index_buffer.push_back(uv);

        tables.index_buffer.push_back(triangle->w_idx());
        tables.index_buffer.push_back(wu);
        tables.index_buffer.push_back(vw);

        tables.index_buffer.push_back(uv);
        tables.index_buffer.push_back(vw);
        tables.index_buffer.push_back(wu);
    }

    QMapIterator<Triangle, bool> extra_it(extra_triangles);
    while (extra_it.hasNext()) {
        extra_it.next();
        Triangle triangle = extra_it.key();
        int uv = -1;
        int vw = -1;
        int wu = -1;
        Edge edge;

        Triangle::Neighbors ns = neighbors.value(triangle);

        if (ns.uv.triangle != NULL) {
            edge = { triangle.u_idx(), triangle.v_idx() };
            if (edge_indices.contains(edge)) {
                uv = edge_indices.value(edge);
            } else {
                uv = offset++;
                edge_indices.insert(edge, uv);

                edge = { ns.uv.edge.a, ns.uv.edge.b };
                edge_indices.insert(edge, uv);

                tables.edge_indices.push_back(triangle.u_idx());
                tables.edge_indices.push_back(triangle.v_idx());
                tables.edge_indices.push_back(triangle.w_idx());
                tables.edge_indices.push_back(ns.uv.edge.c);
            }
        }

        if (ns.vw.triangle != NULL) {
            edge = { triangle.v_idx(), triangle.w_idx() };
            if (edge_indices.contains(edge)) {
                vw = edge_indices.value(edge);
            } else {
                vw = offset++;
                edge_indices.insert(edge, vw);

                edge = { ns.vw.edge.a, ns.vw.edge.b };
                edge_indices.insert(edge, vw);

                tables.edge_indices.push_back(triangle.v_idx());
                tables.edge_indices.push_back(triangle.w_idx());
                tables.edge_indices.push_back(triangle.u_idx());
                tables.edge_indices.push_back(ns.vw.edge.c);
            }
        }

        if (ns.wu.triangle != NULL) {
            edge = { triangle.w_idx(), triangle.u_idx() };
            if (edge_indices.contains(edge)) {
                wu = edge_indices.value(edge);
            } else {
                wu = offset++;
                edge_indices.insert(edge, wu);

                edge = { ns.wu.edge.a, ns.wu.edge.b };
                edge_indices.insert(edge, wu);

                tables.edge_indices.push_back(triangle.w_idx());
                tables.edge_indices.push_back(triangle.u_idx());
                tables.edge_indices.push_back(triangle.v_idx());
                tables.edge_indices.push_back(ns.wu.edge.c);
            }
        }

        if (uv >= 0 && wu >= 0) {
            tables.extra_triangles.push_back(triangle.u_idx());
            tables.extra_triangles.push_back(uv);
            tables.extra_triangles.push_back(wu);
        }

        if (vw >= 0 && uv >= 0) {
            tables.extra_triangles.push_back(triangle.v_idx());
            tables.extra_triangles.push_back(vw);
            tables.extra_triangles.push_back(uv);
        }

        if (wu >= 0 && vw >= 0) {
            tables.extra_triangles.push_back(triangle.w_idx());
            tables.extra_triangles.push_back(wu);
            tables.extra_triangles.push_back(vw);
        }

        if (uv >= 0 && vw >= 0 && wu >= 0) {
            tables.extra_triangles.push_back(uv);
            tables.extra_triangles.push_back(vw);
            tables.extra_triangles.push_back(wu);
        }
    }

    qCDebug(log_subdiv) << "Done with edge table. " << tables.edge_indices.length();
    if (log_subdiv_trace().isDebugEnabled()) {
        qCDebug(log_subdiv) << "Eedges found. Table: " << tables.edge_indices;
        qCDebug(log_subdiv) << "Table (long version):";
        for (int i = 0; i < tables.edge_indices.length(); i++){
            qCDebug(log_subdiv) << "blub:" << tables.edge_indices[i];
        }
    }
}

void Subdivision::precomputeVertexTable(Subdivision::Tables &tables, Input &input, QMap<unsigned int, unsigned int> &modified_vertices) {
    //compute vertex table
    //Format: First entry: original vertex. Other entries: adjacent vertices.
    QVector<QVector<unsigned int> > duplicates;//for debugging
    QVector<unsigned int> ib_combined = input.index_irregular_buffer + input.index_regular_buffer;
    unsigned int offset = 0;
    unsigned int output_offset = input.vertex_buffer.size() + tables.edge_indices.size() / 4;
    for (int i = 0; i < input.vertex_buffer.length(); i++){
        //hat evtl viel redundanz

        QVector<Vertex> adj_v;//helfer
        tables.vertex_offsets.push_back(offset);
        tables.vertex_indices.push_back(i);
        modified_vertices.insert(i, output_offset + i);
        offset++;

        Vertex originalVertex = input.vertex_buffer[i];

        QVector<unsigned int> d;

        for (int j = 0; j < ib_combined.length(); j++){

            //suche verweise auf vertex im indexbuffer
            if (input.vertex_buffer[ib_combined[j]].samePos(originalVertex)){
                d.push_back(j);
                unsigned int i1, i2; //indices for neighbour vertices
                switch (j % 3) {
                case 0:
                    i1 = ib_combined[j+1];
                    i2 = ib_combined[j+2];
                    break;
                case 1:
                    i1 = ib_combined[j-1];
                    i2 = ib_combined[j+1];
                    break;
                case 2:
                    i1 = ib_combined[j-2];
                    i2 = ib_combined[j-1];
                    break;
                }

                Vertex v1,v2; //neighbour vertices;
                v1 = input.vertex_buffer[i1];
                v2 = input.vertex_buffer[i2];
                //check if vertices where already used to influence this vertex
                //Note: Can't use contain. Contain uses equal (would have to implement equal, but we're only interested in the position here.
                bool found1 = false;
                bool found2 = false;
                for (int k = 0; k < adj_v.length(); k++){
                    if (adj_v[k].samePos(v1)) found1 = true;
                    if (adj_v[k].samePos(v2)) found2 = true;
                }
                if (!found1){
                    adj_v.push_back(v1);
                    tables.vertex_indices.push_back(i1);
                    offset++;
                }
                if (!found2){
                    adj_v.push_back(v2);
                    tables.vertex_indices.push_back(i2);
                    offset++;
                }
                //if one adjacent vertex is referenced by multiple different indices, only its first index will appear in the adjacent_indices list.1
            }
        }
        if (!duplicates.contains(d)) duplicates.push_back(d);
    }

    tables.vertex_offsets.push_back(offset);

    qCDebug(log_subdiv) << "Done with vertex index table. ";
    qCDebug(log_subdiv) << "Duplicates: "<<duplicates.length();
    for (int i = 0; i < duplicates.length(); i++){
        qCDebug(log_subdiv_trace) << duplicates[i];
    }
}

void Subdivision::updateIndexBuffer(QVector<unsigned int> &index_buffer, QMap<unsigned int, unsigned int> map) {
    for (int i = 0; i < index_buffer.length(); i++) {
        unsigned int old = index_buffer[i];
        index_buffer[i] = map.value(old, old);
    }
}

void Subdivision::splitRegular(Mesh *mesh) {
    QTime totalTimer;
    totalTimer.start();
    qCDebug(log_timing) << "splitRegular(mesh) started";
    Mesh::Node root = mesh->getRootNode();

    int first_mesh_index = -1;
    if (!root.getFirstMeshIndex(first_mesh_index)) {
        qCCritical(log_subdiv) << "No mesh found, aborting subdivision";
        return;
    }

    Mesh::MeshEntry *current_mesh = mesh->getMeshEntry(first_mesh_index);

    QVector<Triangle> all_triangles;
    ibToTriangles(&current_mesh->buffers[0]->vertices, current_mesh->buffers[0]->indices_irregular, all_triangles);
    QVector<unsigned int> old_regular_ib = patchIBToTriangleIB(current_mesh->buffers[0]->indices_regular);
    ibToTriangles(&current_mesh->buffers[0]->vertices, old_regular_ib, all_triangles);

    QMap<Triangle, Triangle::Neighbors> neighbors;
    buildNeighborsMap(current_mesh->buffers[0]->vertices, all_triangles, neighbors);

    QVector<Triangle> regular;
    QVector<Triangle> irregular;
    findRegular(all_triangles, neighbors, regular, irregular);

    QVector<unsigned int> irregular_ib;
    trianglesToIB(irregular, irregular_ib);

    QVector<unsigned int> patches;
    getPatchIndexBuffer(regular, neighbors, patches);

    current_mesh->buffers[0]->indices_regular = patches;
    current_mesh->buffers[0]->indices_irregular = irregular_ib;
    current_mesh->buffers[0]->updateIndices();


    qCDebug(log_timing) << "splitRegular(mesh): time: " <<formatTimeMeasurement(totalTimer.elapsed());
    qCDebug(log_subdiv) << "splitRegular(mesh): regular: " << regular.length();
    qCDebug(log_subdiv) << "splitRegular(mesh): irregular: " << irregular.length();
}

void Subdivision::getPatchIndexBuffer(QVector<Triangle> &triangles, QMap<Triangle, Triangle::Neighbors> &neighbors, QVector<unsigned int> &patch_ib) {
    QTime totalTimer;
    totalTimer.start();

    QVector<unsigned int> patch;
    patch.resize(12);

    QVectorIterator<Triangle> it(triangles);
    while (it.hasNext()) {
        Triangle triangle = it.next();

        // These three vertices come from the triangle we're looking at
        patch[3] = triangle.u_idx();
        patch[6] = triangle.v_idx();
        patch[7] = triangle.w_idx();

        // These three vertices are those completing the three neighboring triangles
        Triangle::Neighbors ns = neighbors.value(triangle);
        patch[2] = ns.uv.edge.c;
        patch[10] = ns.vw.edge.c;
        patch[4] = ns.wu.edge.c;

        // TODO: there's a lot of repetition in the rest of this loop but for now it's good enough

        // To get the two vertices at the top we go from the left neighbor to its top neighbor and the the top neighbor of that triangle
        Triangle::Neighbors left_neighbors = neighbors.value(*ns.uv.triangle);
        Triangle::Neighbor top_left = left_neighbors.get_neighbor(rotate_edge_name(ns.uv.edge.name));

        Triangle::Neighbors top_left_neighbors = neighbors.value(*top_left.triangle);
        Triangle *top = top_left_neighbors.get_neighbor(rotate_edge_name(top_left.edge.name)).triangle;

        // Then we find the one vertex shared with the main triangle and use the other two vertices
        if (triangle.u().samePos(top->u())) {
            patch[0] = top->w_idx();
            patch[1] = top->v_idx();
        } else if (triangle.u().samePos(top->v())) {
            patch[0] = top->u_idx();
            patch[1] = top->w_idx();
        } else if (triangle.u().samePos(top->w())) {
            patch[0] = top->v_idx();
            patch[1] = top->u_idx();
        } else {
            assert(false);
        }

        // The four vertices at the bottom are similar to those at the top, but we start with the bottom neighbor and go left/right

        Triangle::Neighbors bottom_neighbors = neighbors.value(*ns.vw.triangle);
        Triangle::Neighbor bottom_middle_left = bottom_neighbors.get_neighbor(rotate_edge_name(ns.vw.edge.name));
        Triangle::Neighbors bottom_middle_left_neighbors = neighbors.value(*bottom_middle_left.triangle);
        Triangle *bottom_left = bottom_middle_left_neighbors.get_neighbor(rotate_edge_name(bottom_middle_left.edge.name)).triangle;

        if (triangle.v().samePos(bottom_left->u())) {
            patch[5] = bottom_left->v_idx();
            patch[9] = bottom_left->w_idx();
        } else if (triangle.v().samePos(bottom_left->v())) {
            patch[5] = bottom_left->w_idx();
            patch[9] = bottom_left->u_idx();
        } else if (triangle.v().samePos(bottom_left->w())) {
            patch[5] = bottom_left->u_idx();
            patch[9] = bottom_left->v_idx();
        } else {
            assert(false);
        }

        Triangle::Neighbor bottom_middle_right = bottom_neighbors.get_neighbor(rotate_edge_name(rotate_edge_name(ns.vw.edge.name)));
        Triangle::Neighbors bottom_middle_right_neighbors = neighbors.value(*bottom_middle_right.triangle);
        Triangle *bottom_right = bottom_middle_right_neighbors.get_neighbor(rotate_edge_name(rotate_edge_name(bottom_middle_right.edge.name))).triangle;

        if (triangle.w().samePos(bottom_right->u())) {
            patch[8] = bottom_right->w_idx();
            patch[11] = bottom_right->v_idx();
        } else if (triangle.w().samePos(bottom_right->v())) {
            patch[8] = bottom_right->u_idx();
            patch[11] = bottom_right->w_idx();
        } else if (triangle.w().samePos(bottom_right->w())) {
            patch[8] = bottom_right->v_idx();
            patch[11] = bottom_right->u_idx();
        } else {
            assert(false);
        }

        patch_ib += patch;
    }

    qCDebug(log_timing) << "PatchIndexBuffer done:" << formatTimeMeasurement(totalTimer.elapsed());
}

/**
 * Turns a patch index buffer into an index buffer of triangles.
 */
QVector<unsigned int> Subdivision::patchIBToTriangleIB(QVector<unsigned int> ib) {
    QVector<unsigned int> new_ib;
    for (int i = 0; i < ib.length(); i += 12) {
        new_ib.push_back(ib[i+3]);
        new_ib.push_back(ib[i+6]);
        new_ib.push_back(ib[i+7]);
    }
    return new_ib;
}

bool isVertexRegular(Triangle &triangle, Triangle::Neighbor current_neighbor, QMap<Triangle, Triangle::Neighbors> neighbors) {
    unsigned int count = 1;
    assert(current_neighbor.triangle != NULL);
    while (count < 6 && *current_neighbor.triangle != triangle) {
        count++;
        current_neighbor = neighbors.value(*current_neighbor.triangle).get_neighbor(rotate_edge_name(current_neighbor.edge.name));
        assert(current_neighbor.triangle != NULL);
    }
    return count == 6 && *current_neighbor.triangle == triangle;
}

void Subdivision::findRegular(QVector<Triangle> triangles, QMap<Triangle, Triangle::Neighbors> neighbors, QVector<Triangle> &regular, QVector<Triangle> &irregular) {
    QTime timer;
    timer.start();
    qCDebug(log_timing) << "findRegular started";

    QVectorIterator<Triangle> it(triangles);
    while (it.hasNext()) {
        Triangle triangle = it.next();
        Triangle::Neighbors ns = neighbors.value(triangle);

        bool is_regular = true;
        if (!isVertexRegular(triangle, ns.uv, neighbors)) {
            is_regular = false;
        }
        if (is_regular && !isVertexRegular(triangle, ns.vw, neighbors)) {
            is_regular = false;
        }
        if (is_regular && !isVertexRegular(triangle, ns.wu, neighbors)) {
            is_regular = false;
        }

        if (is_regular) {
            regular.push_back(triangle);
        } else {
            irregular.push_back(triangle);
        }
    }

    qCDebug(log_timing) << "findRegular done:" << formatTimeMeasurement(timer.elapsed());
}

Subdivision::Result Subdivision::runShader(Input input, Tables &tables) {
    qCDebug(log_subdiv) << "Running compute shader";

    Result result;

    // Create an input buffer for the vertex indices
    GLuint vertex_indices_handle;
    f->glGenBuffers(1, &vertex_indices_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, vertex_indices_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, tables.vertex_indices.size() * sizeof(GLuint), tables.vertex_indices.data(), GL_DYNAMIC_DRAW);

    // Create an input buffer for the vertex offsets
    GLuint vertex_offsets_handle;
    f->glGenBuffers(1, &vertex_offsets_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, vertex_offsets_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, tables.vertex_offsets.size() * sizeof(GLuint), tables.vertex_offsets.data(), GL_DYNAMIC_DRAW);

    // Create an input buffer for the edge indices
    GLuint edge_indices_handle;
    f->glGenBuffers(1, &edge_indices_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, edge_indices_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, tables.edge_indices.size() * sizeof(GLuint), tables.edge_indices.data(), GL_DYNAMIC_DRAW);

    // Create an output buffer
    // the new vertex buffer contains the old vertices + the new edge vertices + the modified non-edge vertices
    unsigned int output_size = input.vertex_buffer.size() + tables.edge_indices.size() / 4 + (tables.vertex_offsets.size() - 1);
    GLuint output_handle;
    f->glGenBuffers(1, &output_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, output_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, output_size * sizeof(Vertex), NULL, GL_DYNAMIC_DRAW);
    result.vb_handle = output_handle;

    // Create a buffer for holding the offset into the output
    GLuint offset_handle;
    f->glGenBuffers(1, &offset_handle);
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, offset_handle);
    f->glBufferData(GL_SHADER_STORAGE_BUFFER, sizeof(GLuint), NULL, GL_DYNAMIC_DRAW);

    // copy old vertices
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, 0);
    QTime timer;
    timer.start();
    runCopyShader(input.vertex_buffer.size(), input.vb_handle, output_handle);
    int copyTime = timer.elapsed();

    // calculate positions of new vertices
    timer.restart();
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, offset_handle);
    GLuint *offset_buffer = (GLuint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
    *offset_buffer = input.vertex_buffer.size();
    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);
    runEdgeShader(tables.edge_indices.size() / 4, input.vb_handle, edge_indices_handle, output_handle, offset_handle);
    int edgeTime = timer.elapsed();

    // move old vertices
    timer.restart();
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, offset_handle);
    offset_buffer = (GLuint *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
    *offset_buffer = input.vertex_buffer.size() + tables.edge_indices.size() / 4;
    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);
    runVertexShader(tables.vertex_offsets.size() - 1, input.vb_handle, vertex_indices_handle, vertex_offsets_handle, output_handle, offset_handle);
    int vertexTime = timer.elapsed();

    timer.restart();
    // Map the output buffer so we can read the results
    f->glBindBuffer(GL_SHADER_STORAGE_BUFFER, output_handle);
    Vertex *ptr;
    ptr = (Vertex *) f->glMapBuffer(GL_SHADER_STORAGE_BUFFER, GL_WRITE_ONLY);
    for (unsigned int i = 0; i < output_size; i++) {
        result.vertex_buffer.push_back(ptr[i]);
    }
    int bufferbildTime = timer.elapsed();

    f->glUnmapBuffer(GL_SHADER_STORAGE_BUFFER);

    // Delete the buffers the free the resources
    f->glDeleteBuffers(1, &vertex_indices_handle);
    f->glDeleteBuffers(1, &vertex_offsets_handle);
    f->glDeleteBuffers(1, &edge_indices_handle);
    f->glDeleteBuffers(1, &offset_handle);

    qCDebug(log_timing)<<"Compute Shader done. Copyshader"<<copyTime<<"ms; Vertexshader"<<vertexTime<<"ms; EdgeShader"<<edgeTime<<" BufferToCPU"<<bufferbildTime<<"ms";
    return result;
}

void Subdivision::runCopyShader(GLuint size, GLuint vb_in, GLuint vb_out) {
    copyShader->bind();

    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, vb_in);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, vb_out);

    // Run the shader
    f->glDispatchCompute(size, 1, 1);

    // Wait for the shader to complete and the data to be written back to the global memory
    f->glMemoryBarrier(GL_SHADER_STORAGE_BARRIER_BIT);

    // Unbind the buffers from the shader
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, 0);

    copyShader->release();
}

void Subdivision::runVertexShader(GLuint size, GLuint vb_handle, GLuint vertex_indices_handle, GLuint vertex_offsets_handle, GLuint output_handle, GLuint offset_handle) {
    vertexShader->bind();

    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, vb_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, vertex_indices_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, vertex_offsets_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, output_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 4, offset_handle);

    // Run the shader
    f->glDispatchCompute(size, 1, 1);

    // Wait for the shader to complete and the data to be written back to the global memory
    f->glMemoryBarrier(GL_SHADER_STORAGE_BARRIER_BIT);

    // Unbind the buffers from the shader
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 4, 0);

    vertexShader->release();
}

void Subdivision::runEdgeShader(GLuint size, GLuint vb_handle, GLuint edge_indices_handle, GLuint output_handle, GLuint offset_handle) {
    edgeShader->bind();

    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, vb_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, edge_indices_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, output_handle);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, offset_handle);

    // Run the shader
    f->glDispatchCompute(size, 1, 1);

    // Wait for the shader to complete and the data to be written back to the global memory
    f->glMemoryBarrier(GL_SHADER_STORAGE_BARRIER_BIT);

    // Unbind the buffers from the shader
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 0, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 1, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 2, 0);
    f->glBindBufferBase(GL_SHADER_STORAGE_BUFFER, 3, 0);

    edgeShader->release();
}

void Subdivision::updateIrregularForDraw(const QVector<Triangle> &triangles,QMap<Triangle, Triangle::Neighbors> &neighbors, Result &result){
    result.vertex_buffer_irregular.resize(result.vertex_buffer.size());

    QVectorIterator<Triangle> it(triangles);
    while (it.hasNext()) {
        Triangle triangle = it.next();

        result.vertex_buffer_irregular[triangle.u_idx()] =  result.vertex_buffer.at(triangle.u_idx());
        result.vertex_buffer_irregular[triangle.v_idx()] =  result.vertex_buffer.at(triangle.v_idx());
        result.vertex_buffer_irregular[triangle.w_idx()] =  result.vertex_buffer.at(triangle.w_idx());

        Triangle::Neighbors ns = neighbors.value(triangle);

        result.vertex_buffer_irregular[triangle.u_idx()].pos = updateIrregularVertexForDraw(triangle.u(), vertexNeighbors(triangle, ns.uv, neighbors));
        result.vertex_buffer_irregular[triangle.v_idx()].pos = updateIrregularVertexForDraw(triangle.v(), vertexNeighbors(triangle, ns.vw, neighbors));
        result.vertex_buffer_irregular[triangle.w_idx()].pos = updateIrregularVertexForDraw(triangle.w(), vertexNeighbors(triangle, ns.wu, neighbors));
    }
}

bool Subdivision::containsVertex(QVector<Vertex> &vertices, Vertex vertex) {
    for (int i = 0; i < vertices.length(); i++) {
        if (vertex.samePos(vertices[i])) {
            return true;
        }
    }
    return false;
}

QVector3D Subdivision::updateIrregularVertexForDraw(Vertex currentCorner, QVector<Triangle> neighboring_triangles) {
    QVector<Vertex> surroundingVertex;

    QVectorIterator<Triangle> extra(neighboring_triangles);
    while (extra.hasNext()) {
        Triangle extra_triangle = extra.next();

        if (!currentCorner.samePos(extra_triangle.u()) && !containsVertex(surroundingVertex, extra_triangle.u()))
            surroundingVertex.push_back(extra_triangle.u());
        if (!currentCorner.samePos(extra_triangle.v()) && !containsVertex(surroundingVertex, extra_triangle.v()))
            surroundingVertex.push_back(extra_triangle.v());
        if (!currentCorner.samePos(extra_triangle.w()) && !containsVertex(surroundingVertex, extra_triangle.w()))
            surroundingVertex.push_back(extra_triangle.w());
    }

    QVector3D newPos(0.0, 0.0, 0.0);

    double n = surroundingVertex.length();
    double a = 5.0/8.0 - pow((3.0 / 8.0 + 2.0/8.0 * cos(2.0 * M_PI / n)), 2.0);
    double omega = 3.0 * n / (8.0 * a);
    for (int i = 0; i < surroundingVertex.length(); i++) {
        newPos += surroundingVertex[i].pos;
    }
    newPos *= 1.0 / (omega + n);
    newPos += omega / (omega + n) * currentCorner.pos;

    return newPos;
}
