//
// This file is part of MorphoGraphX - http://www.MorphoGraphX.org
// Copyright (C) 2012-2016 Richard S. Smith and collaborators.
//
// If you use MorphoGraphX in your work, please cite:
//   http://dx.doi.org/10.7554/eLife.05864
//
// MorphoGraphX is free software, and is licensed under under the terms of the 
// GNU General (GPL) Public License version 2.0, http://www.gnu.org/licenses.
// 
#include <MeshProcessCellMesh.hpp>
#include <GraphUtils.hpp>
#include <Progress.hpp>

#include <queue>

namespace mgx
{
  using std::swap;

  bool FixBorderTriangles::run(Mesh* mesh)
  {
    vvGraph& S = mesh->graph();
    if(S.empty())
      throw(QString("Mesh empty."));

    // Only attempt fixing corners on a normal mesh
    if(mesh->meshType() != "MGXM")
      throw(QString("Mesh must be of type MGXM."));

    // Start the progress bar
    progressStart("Fixing border triangles", 0);

    // Fix Border triangles
    int count = fixBorderTris(S);

    // Update the GUI
    mesh->updateAll();

    SETSTATUS(QString("%1 border triangles fixed").arg(count));
    return true;
  }
  REGISTER_PROCESS(FixBorderTriangles);


  bool FixMeshCorners::run(Mesh* mesh, bool autoSegment, bool selectBadVertices, int maxRun)
  {
    vvGraph& S = mesh->graph();
    if(S.empty())
      throw(QString("Mesh empty."));

    // Only attempt fixing corners on a normal mesh
    if(mesh->meshType() != "MGXM")
      throw(QString("Mesh must be of type MGXM."));

    // select unlabeled/troublesome vertices and put them in a list during the first
    // run of fix corner. The subsequent re-segmentation will run only on them,
    // so it goes faster. Afterwards we take them out of the list if they are labeled.
    // new vertices resulting from subdivision are also in the list.
    set_vector<vertex> vCheckSwap;

    // Important for consistency among the mesh processes
    const std::vector<vertex>& avs = mesh->activeVertices();
    set_vector<vertex> vCheck;

    // Starts by checking border and unlabeled vertices
    forall(vertex v, avs)
        if(v->label <= 0)
        vCheck.insert(v);

    progressStart(QString("Fixing border triangles for mesh %1").arg(mesh->userId()), 0, false);

    int tricount = 0, prevtricount = INT_MAX, totaltricount = 0;
    int vdelcnt = 0;
    int nRun;

    for(nRun = 0; not vCheck.empty() and prevtricount != 0 and nRun < maxRun; ++nRun) {
      if(nRun == 0)
        prevtricount = 0;
      std::vector<Point3ul> tris;
      // In first run we check all vertices, after re-segmentation we check only vertices which caused trouble in
      // previous runs.
      forall(vertex v, vCheck) {
        // Check for triangles with all vertices on the margin
        if(v->label == -1 and v->margin) {
          forall(vertex n, S.neighbors(v)) {
            vertex m = S.nextTo(v, n);
            if(!S.uniqueTri(v, n, m) or !n->margin or n->label != -1 or !m->margin or m->label != -1)
              continue;
            tris.push_back(Point3ul(v.id(), n.id(), m.id()));
          }
        }
      }

      // Margin triangles must have a vertex with valence 2, delete it
      vvGraph N;
      forall(const Point3ul& tri, tris) {
        for(int i = 0; i < 2; i++) {
          vertex v(tri[i]);
          if(S.valence(v) <= 2) {
            // Add the neighbor, who might now be on the margin
            forall(const vertex& n, S.neighbors(v))
                N.insert(n);
            vdelcnt++;
            S.erase(v);
            // Remove in neighbor set if there
            if(N.contains(v))
              N.erase(v);
          }
        }
      }
      // Check if neighbors of deleted vertices now on the margin
      markMargins(N, S);

      tris.clear();
      forall(vertex v, vCheck) {
        if(v->label == -1) {
          forall(vertex n, S.neighbors(v)) {         // Check for triangles with all vertices on the border
            vertex m = S.nextTo(v, n);
            if(!S.uniqueTri(v, n, m) or n->label != -1 or m->label != -1)
              continue;
            if(S.edge(n, m))
              tris.push_back(Point3ul(v.id(), n.id(), m.id()));
          }
        } else if(v->label == 0) {       // Check for unlabeled vertices that have only border neighbors
          bool allborder = true;
          forall(vertex n, S.neighbors(v))
              if(n->label != -1)
              allborder = false;

          if(allborder) {
            forall(vertex n, S.neighbors(v)) {
              vertex m = S.nextTo(v, n);
              if(S.edge(n, m))
                tris.push_back(Point3ul(v.id(), n.id(), m.id()));
            }
          }
        }

        // Subdivide internal borders ("islands") within skinny cells
        if(v->label == -1) {
          forall(vertex n, S.neighbors(v)) {
            if(n->label != -1)
              continue;
            if(S.valence(v) <= 1)
              Information::out << "Problem:with valence" << S.valence(v) << endl;

            vertex nn = S.nextTo(v, n);
            vertex pn = S.prevTo(v, n);

            if(nn->label > 0 and nn->label == pn->label and S.edge(n, nn))
              tris.push_back(Point3ul(v.id(), n.id(), nn.id()));
          }
        }
      }

      // Subdivide triangles and accumulate new vertices to check on next iteration
      tricount = 0;
      forall(const Point3ul& tri, tris) {
        vertex v(tri.x()), n(tri.y()), m(tri.z());
        v->label = m->label = n->label = 0;
        std::vector<vertex> newCheck;
        newCheck.push_back(v);
        newCheck.push_back(m);
        newCheck.push_back(n);
        if(m == S.nextTo(v, n)) {
          tricount++;
          if(S.edge(v, n) and S.edge(n, m) and S.edge(m, v))
            subdivideBisect(S, v, n, m, &newCheck);
          else
            Information::out << "Error:FixMeshCorners:Attempting to subdivide bad triangle" << endl;
        }
        vCheckSwap.insert(newCheck.begin(), newCheck.end());
      }
      markMargins(S);

      // re-segment mesh only for troublesome vertices (unlabeled, subdivided etc.)
      if(autoSegment and not vCheckSwap.empty()) {
        SegmentMesh sm(*this);
        bool result;
        result = sm.run(mesh, 0, vCheckSwap.vector());
        if(!result)
          return false;
      }

      totaltricount += tricount;

      Information::out << "nRun " << nRun << " tricount: " << tricount << " prevtricount: "
                       << prevtricount << " totaltricount: " << totaltricount
                       << " size vCheck: " << vCheck.size() << endl;

      prevtricount = tricount;

      // Put in vCheck the list of vertices to check on the next round
      swap(vCheck, vCheckSwap);
      vCheckSwap.clear();
    }

    mesh->updateAll();

    SETSTATUS(vdelcnt << " margin vertices deleted and " << totaltricount
              << " border triangles subdivided for Mesh "
              << mesh->userId() << " in " << nRun << " runs.");

    // un-select fixed vertices after last run, keep unlabeled and bad normals selected
    forall(vertex v, S)
        if(v->label == 0)
        vCheck.insert(v);
    if(not vCheck.empty()) {
      QString are_selected;
      if(selectBadVertices) {
        forall(vertex v, S)
            v->selected = false;
        forall(vertex v, vCheck)
            v->selected = true;
        mesh->setMeshView("Selected");
        are_selected = " All bad vertices are selected.";
      }
      QString errMsg = "Remaining bad vertices after fix corners. "
                       "%1 unlabeled vertices found.%2";
      return setErrorMessage(errMsg.arg(vCheck.size()).arg(are_selected));
    }

    return true;
  }
  REGISTER_PROCESS(FixMeshCorners);
  
  bool ConvertToMgxm::run(Mesh* mesh)
  {
    if(!mesh->tissue().toMgxm())
      throw(QString("Unexpected error converting to normal (MGXM) mesh"));

    SETSTATUS(QString("Created normal mesh (MGXM) in mesh %1 with %2 vertices.")
              .arg(mesh->userId()).arg(mesh->tissue().S.size()));
    // Update the GUI
    mesh->updateAll();

    return true;
  }
  REGISTER_PROCESS(ConvertToMgxm);

  bool ConvertToMgxc::run(Mesh* mesh, double wallMax)
  {

    if(!mesh->tissue().toMgxc(wallMax))
      throw(QString("Unexpected error converting to cell (MGXC) mesh"));

    mesh->updateCentersNormals();

    SETSTATUS(QString("Created cell mesh (MGXC) in mesh %1 with %2 vertices, %3 cells.")
              .arg(mesh->userId()).arg(mesh->tissue().S.size()).arg(mesh->labelCenter().size()));

    // Update the GUI
    mesh->updateAll();
    /*
    correctNormals(mesh->graph());

    forall(const vertex& v, mesh->graph()) {
      v->selected = false;
      if(v->nrml.z() < 0){
        std::cout << "c " << v->nrml << "//" << v->label << std::endl;
      } else {
        v->selected = true;
        std::cout << ".";
      }
    }
*/

    return true;
  }
  REGISTER_PROCESS(ConvertToMgxc);

  bool ConvertToMgx2d::run(Mesh* mesh, double wallMax)
  {
    // Convert the MGX2D
    if(!mesh->tissue().toMgx2d(wallMax))
      throw(QString("Unexpected error converting to 2D cell tissue (MGX2D)"));

    SETSTATUS(QString("Created cell tissue in mesh %1 with %2 cells.")
              .arg(mesh->userId()).arg(mesh->tissue().C.size()));
    // Update the GUI
    mesh->updateAll();

    return true;
  }
  REGISTER_PROCESS(ConvertToMgx2d);

  bool ConvertToMgx3d::run(Mesh* mesh, double tolerance, double neighborMinArea)
  {
    // Convert the MGXC
    if(!mesh->tissue().toMgx3d(tolerance, neighborMinArea))
      throw(QString("Unexpected error converting to 3D cell tissue"));

    SETSTATUS(QString("Created 3D cell tissue in mesh %1 with %2 cells.")
              .arg(mesh->userId()).arg(mesh->tissue().C.size()));
    // Update the GUI
    mesh->updateAll();

    return true;
  }
  REGISTER_PROCESS(ConvertToMgx3d);

  bool SurfaceFromMGX3D::run(Mesh* mesh1, Mesh* mesh2)
  {
    if(mesh1->meshType() != "MGX3D")
      throw(QString("Mesh type (%1) doesn't have shared triangles, mesh type must be (MGX3D)").arg(mesh1->meshType()));
    const cellGraph& Cells = mesh1->cells();

    // look which vertices belong to more than one cell (shared)
    // and grab vertex labels from cell labels
    std::map<vertex,int> vertLabel;
    std::map<vertex,int>::iterator itV;
    forall(const cell& c, Cells){
      forall(const vertex& v, c->S){
        itV = vertLabel.find(v);
        // insert new vertex in map
        if(itV == vertLabel.end())
          vertLabel[v] = c->label;
        else {
          // if vertex already in map, it should shared between cells
          if(itV->second != c->label)
            vertLabel[v] = -1;
          else
            Information::out << "duplicate vertex in cell nr " << c->label << endl;
        }
      }
    }

    // list of all vertices and triangles that will be in final mesh
    std::vector<vertex> vertices;
    std::vector<Point3i> triangles;

    forall(const cell& c, Cells){
      forall(const vertex& v, c->S){
        forall(const vertex& n, c->S.neighbors(v)) {
          const vertex& m = c->S.nextTo(v, n);
          if(!c->S.uniqueTri(v, n, m))
            continue;
          // do not consider triangles with only shared vertices
          if((vertLabel[v] == -1) and (vertLabel[n] == -1) and (vertLabel[m] == -1))
            continue;
          // put vertices into vector, use vertex position in vector
          // as vertex index (to identify triangles)
          std::vector<vertex>::iterator it;
          int vId, nId, mId;
          // insert vertex v
          it = std::find(vertices.begin(), vertices.end(), v);
          if(it!= vertices.end()){
            // v already in list
            vId = it - vertices.begin();
          } else {
            // add v to list
            vertices.push_back(v);
            vId = vertices.size() -1;
          }
          // insert vertex n
          it = std::find(vertices.begin(), vertices.end(), n);
          if(it!= vertices.end())
            nId = it - vertices.begin();
          else {
            vertices.push_back(n);
            nId = vertices.size() -1;
          }
          // insert vertex m
          it = std::find(vertices.begin(), vertices.end(), m);
          if(it!= vertices.end())
            mId = it - vertices.begin();
          else {
            vertices.push_back(m);
            mId = vertices.size() -1;
          }
          triangles.push_back(Point3i(vId, nId, mId));
        }
      }
    }

    //make list of new vertices to create mesh2, copy position from mesh1.
    std::vector<vertex> newVertices;
    for (std::vector<vertex>::iterator it = vertices.begin(); it != vertices.end(); ++it){
      Point3d pos = (*it)->pos;
      vertex newV;
      newV->label = vertLabel[*it];
      newV->pos = pos;
      newVertices.push_back(newV);
    }
    // create new graph
    vvGraph newS;
    Mesh::meshFromTriangles(newS, newVertices, triangles, true);

    // Put new graph in mesh2
    mesh2->reset();
    mesh2->setMeshType("MGXM");
    mesh2->graph().swap(newS);
    mesh2->setShowLabel("Label");
    mesh2->updateAll();

    return true;
  }
  REGISTER_PROCESS(SurfaceFromMGX3D);


  bool RelabelFragments::run(Mesh *mesh1, int relabelThreshold)
  { 
    vvGraph& S = mesh1->graph();
    std::queue<vertex> vtxQ;
    std::set<vertex> visited, inQ;
    std::map<int, std::vector<std::set<vertex> > > labelPiecesMap;
    typedef std::pair<int, std::vector<std::set<vertex> > > PieceVec;

    bool keepParents = true;

    int maxLabel = -1;

    // go through all vtx
    forall(const vertex& v, S){
      if(visited.find(v) != visited.end() or v->label <= 0) continue;

      if(v->label > maxLabel) maxLabel = v->label;

      inQ.insert(v); // if not visited -> queue for visit
      
      std::set<vertex> currentPiece;
      currentPiece.insert(v);

      while(inQ.size() > 0){
        const vertex& currentVtx = *(inQ.begin());
        inQ.erase(currentVtx);
        visited.insert(currentVtx);
        // emtpy the queue -> visit all neighbors
        forall(const vertex& n, S.neighbors(currentVtx)){
          if(visited.find(n) != visited.end() or n->label <= 0) continue;
          inQ.insert(n);
          currentPiece.insert(n);
        }

      }
      labelPiecesMap[v->label].push_back(currentPiece);

    }

    forall(PieceVec p, labelPiecesMap){
      uint sizeMax = 0;
      uint idxMax = 0;
      // find the largest piece
      for(uint i = 0; i < p.second.size(); i++){
        if(sizeMax < p.second[i].size()){
          sizeMax = p.second[i].size();
          idxMax = i;
        }
      }

      // all other pieces are unlabelled
      for(uint i = 0; i < p.second.size(); i++){
        if(i == idxMax) continue;
        int newLabel = 0;
        if(p.second[i].size() > relabelThreshold and relabelThreshold != -1) newLabel = ++maxLabel;
        if(keepParents) mesh1->parents()[newLabel] = mesh1->parents()[p.first];
        forall(const vertex& v, p.second[i]) v->label = newLabel;
      }
    }

    mesh1->updateAll();
    return true;
  }
  REGISTER_PROCESS(RelabelFragments);


  bool ExtendCells::run(Mesh *m)
    {
  
    vvGraph& S = m->graph();
    markMargins(S);
    std::set<vertex> vtxToCheck;

    forall(const vertex& v, S){
      if(v->label != 0) continue;
      vtxToCheck.insert(v);
    }

    bool abort = false;

    while(!vtxToCheck.empty() and !abort){
      std::set<vertex> vtxToCheckNew;
      //std::cout << "check " << vtxToCheck.size() << std::endl;

      forall(const vertex& v, vtxToCheck){
        bool found = false;
        //int n = 0;
        forall(const vertex& n, S.neighbors(v)){
          if(n->label > 0){
            found = true;
            v->label = n->label;
            break;
          } 

          std::set<int> neighborL;

          forall(const vertex& k, S.neighbors(n)){
            if(k->label < 1) continue;

            neighborL.insert(k->label);
            //std::cout << "i " << k->label << std::endl;
          }
          //std::cout << "n " << neighborL.size() << std::endl;
          if(neighborL.size() != 1) continue;
          n->label = *(neighborL.begin());
          found = true;
          //std::cout << "found " << n->label << std::endl;
          break;

        }
        if(!found){
          vtxToCheckNew.insert(v);
        } 
      }
      if(vtxToCheck.size() == vtxToCheckNew.size()) abort = true;

      vtxToCheck = vtxToCheckNew;
      markMargins(S);
    }

    m->updateAll();

    return true;
    }
    REGISTER_PROCESS(ExtendCells);


}
