EgtGeomKernel/SurfTriMeshUtilities.cpp

//----------------------------------------------------------------------------
//  EgalTech 2021-2021
//----------------------------------------------------------------------------
// File : SurfTriMeshUtilities.cpp     Data : 01.11.21       Versione : 2.3k1
// Contenuto : Implementazione funzioni di utilità di Superfici TriMesh.
//
//
//
// Modifiche : 25.10.21 LM Creazione modulo. 
//
//----------------------------------------------------------------------------

//--------------------------- Include ----------------------------------------
#include "stdafx.h"
#include "Triangulate.h"
#include "SurfTriMesh.h"
#include "DistPointLine.h"
#include <unordered_map>

using namespace std ;

//----------------------------------------------------------------------------
static bool
IsVertex( PNTULIST& PointList, PNTULIST::const_iterator itCurr)
{
  // recupero il punto precedente
   PNTULIST::const_iterator itPrev ;
   if ( itCurr == PointList.begin())
      itPrev = prev( PointList.end(), 2) ;
   else
      itPrev = prev( itCurr) ;
  // recupero il punto successivo
   auto itNext = next( itCurr) ;
   if ( itNext == PointList.end())
      itNext = next( PointList.begin()) ;
  // se cambia faccia adiacente tra prima e dopo, va bene
   if ( itPrev->second != itCurr->second)
      return true ;
  // se lati aperti e cambia direzione tra prima e dopo, va bene
   if ( itPrev->second == -1) {
      DistPointLine PointLineDistCalc( itCurr->first, itPrev->first, itNext->first) ;
      double dDist ;
      if ( PointLineDistCalc.GetDist( dDist)  &&  dDist > EPS_SMALL)
         return true ;
   }
  // altrimenti non va bene
   return false ;
}

//----------------------------------------------------------------------------
static bool
ChooseGoodStartPoint( PNTULIST& PointList)
{
  // se il punto iniziale è un vertice, non devo fare alcunché
   if ( IsVertex( PointList, PointList.begin()))
      return true ;
  // altrimenti cerco il vertice più vicino
   for ( auto it = next( PointList.begin()) ; it != PointList.end() ; ++it) {
      if ( IsVertex( PointList, it)) {
        // se ultimo punto non devo fare alcunché
         if ( next( it) == PointList.end())
            return false ;
        // cancello ultimo punto ( coincide con primo)
         PointList.pop_back() ;
        // sposto la parte iniziale dei punti alla fine
         PointList.splice( PointList.end(), PointList, PointList.begin(), it) ;
        // aggiungo punto finale come copia dell'iniziale
         PointList.push_back( PointList.front()) ;
        // ho finito
         return true ;
      }
   }
   return false ;
}

//----------------------------------------------------------------------------
static bool
AdjustLoop( PNTULIST& PointList, double dMaxEdgeLen, bool& bModif)
{
  // Ciclo sui punti del loop
   auto itLast = PointList.begin() ;
   for ( auto it = next( itLast) ; it != PointList.end() ; ++ it) {
   
     // Se dal punto corrente inizia un segmento adiacente a un'altra faccia
      if ( itLast->second != it->second) {
        // Elimino i punti interni
         auto itNextToLast = next( itLast) ;
         for ( auto itInn = itNextToLast ; itInn != it ; ) {
            itInn = PointList.erase( itInn) ;
            bModif = true ;
         }
        // Se la lunghezza del segmento supera il limite imposto
         double dSegLen = Dist( it->first, itLast->first) ;
         if ( dSegLen > dMaxEdgeLen) {
           // determino il numero di step
            double dRatio = dSegLen / dMaxEdgeLen ;
            int nStepCount = int( dRatio) + 1 ;
           // inserisco i punti
            auto itAdd = it ;
            for ( int nP = 1 ; nP < nStepCount ; ++ nP) {
               double dCoeff = double( nP) / nStepCount ;
               itAdd = PointList.insert( itAdd, POINTU( Media( itLast->first, it->first, 1 - dCoeff), itLast->second)) ;
               bModif = true ;
            }
         }
        // Nuovo punto di riferimento
         itLast = it ;
      }

     // Se sono due segmenti liberi non allineati
      else if ( itLast->second == - 1) {
        // Calcolo se i punti compresi fra gli estremi sono allineati
         bool bAreAligned = true ;
         auto itNextToLast = next( itLast) ;
         for ( auto itInn = itNextToLast ; itInn != it  &&  bAreAligned ; ++ itInn) {
            DistPointLine PointLineDistCalc( itInn->first, itLast->first, it->first) ;
            double dDist ;
            if ( PointLineDistCalc.GetDist( dDist))
               bAreAligned = ( dDist < EPS_SMALL) ;
         }
        // Se i punti sono allineati
         if ( bAreAligned) {
           // Verifico se il successivo punto non è più allineato
            auto itNextToCurr = next( it) ;
            if ( itNextToCurr != PointList.end()) {
               for ( auto itInn = itNextToLast ; itInn != itNextToCurr  &&  bAreAligned ; ++ itInn) {
                  DistPointLine PointLineDistCalc( itInn->first, itLast->first, itNextToCurr->first) ;
                  double dDist ;
                  if ( PointLineDistCalc.GetDist( dDist))
                     bAreAligned = ( dDist < EPS_SMALL) ;
               }
            }
           // Se ho trovato un insieme massimale di punti allineati li processo
            if ( ! bAreAligned  ||  itNextToCurr == PointList.end()) {
              // Elimino i punti interni
               for ( auto itInn = itNextToLast ; itInn != it ; ) {
                  itInn = PointList.erase( itInn) ;
                  bModif = true ;
               }
              // Se la lunghezza del segmento supera il limite imposto
               double dSegLen = Dist( it->first, itLast->first) ;
               if ( dSegLen > dMaxEdgeLen) {
                 // determino il numero di step
                  double dRatio = dSegLen / dMaxEdgeLen ;
                  int nStepCount = int( dRatio) + 1 ;
                 // inserisco i punti
                  auto itAdd = it ;
                  for ( int nP = 1 ; nP < nStepCount ; ++ nP) {
                     double dCoeff = double( nP) / nStepCount ;
                     itAdd = PointList.insert( itAdd, POINTU( Media( itLast->first, it->first, 1 - dCoeff), itLast->second)) ;
                     bModif = true ;
                  }
               }
              // Nuovo punto di riferimento
               itLast = it ;
            }
         }
      }
   }

   return true ;
}

//----------------------------------------------------------------------------
bool
SurfTriMesh::SimplifyFacets( double dMaxEdgeLen, bool bForced)
{
  // La trimesh deve essere valida
   if ( ! IsValid())
      return false ;
  // Se la lunghezza massima del lato del triangolo sul bordo della faccia è nulla, non devo fare alcunché
   if ( dMaxEdgeLen < EPS_SMALL)
      return true ;

  // Recupero il numero delle facce (esegue anche una verifica delle stesse)
   int nFacetCnt = GetFacetCount() ;

  // Ciclo sulle facce della mesh per trovare quelle da ritriangolare
   unordered_map< int, pair< PNTVECTOR, INTVECTOR>> FacetMap ;
   for ( int nF = 0 ; nF < nFacetCnt ; ++ nF) {
      
     // Recupero i loop della faccia (il parametro indica la faccia adiacente)
      POLYLINEVECTOR LoopVec ;
      GetFacetLoops( nF, LoopVec) ;

     // Ciclo sui loop della faccia
      bool bToRetriangulate = bForced ;
      for ( int nL = 0 ; nL < int( LoopVec.size()) ; ++ nL) {

        // Lista dei punti del loop
         PNTULIST& PointList = LoopVec[nL].GetUPointList() ;

        // Mi assicuro che il punto iniziale/finale non sia all'interno di un possibile segmento
         if ( ! ChooseGoodStartPoint( PointList))
            continue ;

        // Sistemo il loop
         bool bModif = false ;
         if ( ! AdjustLoop( PointList, dMaxEdgeLen, bModif))
            return false ;
         if ( bModif)
            bToRetriangulate = true ;
      }

     // Se non richiesta ritriangolazione dai bordi, verifico se ci sono vertici di triangoli interni
      if ( false  &&  ! bToRetriangulate) {
        // numero dei triangoli nella faccia
         INTVECTOR vFacetTria ;
         GetAllTriaInFacet( nF, vFacetTria) ;
         int nTriaCnt = int( vFacetTria.size()) ;
        // numero dei lati di contorno della faccia
         int nSideCnt = 0 ;
         for ( int nL = 0 ; nL < int( LoopVec.size()) ; ++ nL)
            nSideCnt += LoopVec[nL].GetLineNbr() ;
        // numero dei buchi della faccia
         int nHoleCnt = int( LoopVec.size()) - 1 ;
        // dalla formula di Eulero adattata al caso ( nTriaCnt = nSideCnt + 2 * ( nHoleCnt -1))
         if ( nTriaCnt != nSideCnt + 2 * ( nHoleCnt - 1))
            bToRetriangulate = true ;
      }

     // Se da ritriangolare
      if ( bToRetriangulate) {
        // Eseguo la ritriangolazione della faccia
         PNTVECTOR vPt ;
         INTVECTOR vTr ;
         if ( Triangulate().Make( LoopVec, vPt, vTr)) {
            FacetMap.emplace( nF, make_pair( vPt, vTr)) ;
         }
        // Se non riesco a triangolare anche solo questa faccia, interrompo tutto
         else
            return false ;
      }
   }

  // Ciclo sulle facce da ritriangolare per eliminare i triangoli (nel contempo salvo flag colore)
   INTVECTOR vDelTria ;
   unordered_map< int, int> ColorMap ;
   for ( auto itF = FacetMap.begin() ; itF != FacetMap.end() ; ++ itF) {
     // Recupero i triangoli della faccia
      INTVECTOR vFacetTria ;
      GetAllTriaInFacet( itF->first, vFacetTria) ;
      vDelTria.insert( vDelTria.end(), vFacetTria.begin(), vFacetTria.end()) ;
     // Salvo il colore della faccia da flag di un suo triangolo
      ColorMap.emplace( itF->first, m_vTria[m_vFacet[itF->first]].nTFlag) ;
   }
  // Cancello i triangoli
   for ( int nT : vDelTria)
      RemoveTriangle( nT) ;
   
  // Applico le nuove triangolazioni delle facce
   for ( auto itF = FacetMap.begin() ; itF != FacetMap.end() ; ++ itF) {
      const PNTVECTOR& vPt = itF->second.first ;
      const INTVECTOR& vTr = itF->second.second ;
     // Inserisco i nuovi triangoli
      for ( int n = 0 ; n < int( vTr.size()) - 2 ; n += 3) {
         int nNewId[3] = { AddVertex( vPt[vTr[n]]),
                           AddVertex( vPt[vTr[n + 1]]),
                           AddVertex( vPt[vTr[n + 2]])} ;
         auto itCol = ColorMap.find( itF->first) ;
         int nTFlag = ( itCol != ColorMap.end() ? itCol->second : 0) ;
         int nNewTriaId = AddTriangle( nNewId, nTFlag) ;
      }
   }
   
  // dichiaro necessità ricalcolo della grafica e di hashgrids3d
   m_OGrMgr.Reset() ;
   ResetHashGrids3d() ;

  // Eseguo aggiustamenti
   return ( AdjustVertices()  &&  DoCompacting()) ;
}

//----------------------------------------------------------------------------
bool
SurfTriMesh::AddChainToChain( const Chain& ChainToAdd, PNTVECTOR& OrigChain)
{
  // Se la catena da aggiungere è vuota, non devo fare alcunchè
   if ( ChainToAdd.size() == 0)
      return true ;
  // Se la catena originale è vuota, non è possibile aggiungere nulla
   if ( OrigChain.size() == 0)
      return false ;
  // Se la catena originale è chiusa non posso aggiungere nulla
   int nLastOrig = max( int( OrigChain.size()) - 1, 0) ;
   if ( AreSamePointApprox( OrigChain[0], OrigChain[nLastOrig]))
      return false ;
   int nLastToAdd = max( int( ChainToAdd.size()) - 1, 0) ;  
   if ( AreSamePointApprox( OrigChain[nLastOrig], ChainToAdd[0].ptSt)) {
      for ( int nPt = 1 ; nPt <= nLastToAdd ; ++ nPt) {
         if ( nPt == nLastToAdd) {
            if ( ! AreSamePointApprox(OrigChain[0], ChainToAdd[nPt].ptSt))
               OrigChain.emplace_back( ChainToAdd[nPt].ptSt) ;
         }
         else if ( nPt == 1) {
            if ( ! AreSamePointApprox( OrigChain[nLastOrig], ChainToAdd[nPt].ptSt))
               OrigChain.emplace_back( ChainToAdd[nPt].ptSt) ;
         }
         else
            OrigChain.emplace_back( ChainToAdd[nPt].ptSt) ;
      }
      return true ;
   }
   else
      return false ;
}

//----------------------------------------------------------------------------
// Una faccia di una trimesh ha una sola componente connessa, con un loop esterno e possibili loop interni
bool
SurfTriMesh::DistPointFacet( const Point3d& ptP, const POLYLINEVECTOR& vPolyVec, double& dPointFacetDist)
{
  // Verifico la presenza del loop esterno
   if ( vPolyVec.size() < 1)
      return false ;
  // Proietto il punto sul piano della faccia, utilizzando il loop esterno
   Plane3d plPlane ;
   double dArea ;
   if ( ! vPolyVec[0].IsClosedAndFlat( plPlane, dArea))
      return false ;
   double dDistPtPl = DistPointPlane( ptP, plPlane) ;
   Point3d ptProjP = ptP + dDistPtPl * plPlane.GetVersN() ;
  // Verifico se il punto proiettato è esterno al loop esterno
   int nPtOut = -1 ;
   if ( ! IsPointInsidePolyLine( ptProjP, vPolyVec[0], EPS_SMALL))
      nPtOut = 0 ;
  // Verifico se il punto proiettato è interno ai loop interni (quindi esterno alla faccia)
   for ( int nLoop = 1 ; nLoop < int( vPolyVec.size()) && nPtOut < 0 ; ++ nLoop) {
      Plane3d plPlane ;
      double dArea ;
      if ( ! vPolyVec[nLoop].IsClosedAndFlat( plPlane, dArea))
         return false ;
      if ( IsPointInsidePolyLine( ptProjP, vPolyVec[nLoop], EPS_SMALL))
         nPtOut = nLoop ;
   }
  // Se il punto si proietta sulla faccia, la distanza dalla faccia coincide con quella dal piano
   if ( nPtOut < 0) {
      dPointFacetDist = abs( dDistPtPl) ;
      return true ;
   }
  // Altrimenti calcolo la minima distanza del punto dalla polilinea del contorno a cui è esterno
   double dDist ;
   if ( DistPointPolyLine( ptP, vPolyVec[nPtOut], dDist)) {
      dPointFacetDist = dDist ; 
      return true ;
   }
   return false ;
}

//----------------------------------------------------------------------------
bool
SurfTriMesh::ChangeStart( const Point3d& ptNewStart, PNTVECTOR& Loop)
{
  // Cerco il tratto del loop chiuso più vicino al punto
   int nMinSeg = - 1 ;  
   double dMinSqDinst = DBL_MAX ;
   for ( int nPt = 0 ; nPt < int( Loop.size()) ; ++ nPt) {
     // Estremi del segmento corrente del loop
      Point3d ptSegSt = Loop[nPt] ;
      Point3d ptSegEn = Loop[( nPt + 1) % int( Loop.size())] ;
     // Distanza del punto dal segmento del loop
      DistPointLine dDistCalc( ptNewStart, ptSegSt, ptSegEn) ;
      double dSqDist ;
      dDistCalc.GetSqDist( dSqDist) ;
      if ( dSqDist < dMinSqDinst) {
         dMinSqDinst = dSqDist ;
         nMinSeg = nPt ;
      }
   }
  // Se il punto non sta sul loop, errore
   if ( dMinSqDinst > SQ_EPS_SMALL)
      return false ;
  // Verifico che il punto stia su un vertice, in tal caso non devo fare nulla
   bool bOnStart = AreSamePointApprox( Loop[nMinSeg], ptNewStart) ;
   bool bOnEnd = AreSamePointApprox( Loop[( nMinSeg + 1) % int( Loop.size())], ptNewStart) ;
   if ( bOnStart  ||  bOnEnd) {
      if ( bOnEnd) {
         ++ nMinSeg ;
         if ( nMinSeg % int( Loop.size()) == 0)
            return true ;
      }
      PNTVECTOR vTempVec ;
      for ( int nPt = 0 ; nPt < nMinSeg ; ++ nPt)
         vTempVec.emplace_back( Loop[nPt]) ;
      int nSize = int( Loop.size()) ;
      for ( int nPt = 0 ; nPt < nSize - nMinSeg ; ++ nPt) {
         Loop[nPt] = Loop[nPt + nMinSeg] ;
      }
      for ( int nPt = 0 ; nPt < int( vTempVec.size()) ; ++ nPt) {
         Loop[nPt + nSize - nMinSeg] = vTempVec[nPt] ;
      }
      return true ;
   }
  // Ridimensiono il loop
   Loop.resize( Loop.size() + 1) ;
  // Copio i primi punti
   PNTVECTOR LoopTemp ;
   for ( int nPt = 0 ; nPt <= nMinSeg ; ++ nPt)
      LoopTemp.emplace_back( Loop[nPt]) ;
  // Aggiungo il nuovo punto all'inizio
   Loop[0] = ptNewStart ;
  // Sposto gli ultimi in testa
  int nLastPointNum = int( Loop.size()) - 1 - nMinSeg ;
   for ( int nPt = 1 ; nPt <= nLastPointNum ; ++ nPt) {
      Loop[nPt] = Loop[nPt + nMinSeg] ;
   }
  // Porto i primi in fondo 
   for ( int nPt = 0 ; nPt < int( LoopTemp.size()) ; ++ nPt) {
      Loop[nPt + nLastPointNum] = LoopTemp[nPt] ;
   }
   return true ;
}

//----------------------------------------------------------------------------
bool
SurfTriMesh::SplitAtPoint( const Point3d& ptStop, const PNTVECTOR& Loop, PNTVECTOR& Loop1, PNTVECTOR& Loop2)
{
  // Cerco il tratto del loop chiuso più vicino al punto
   int nMinSeg = -1 ;
   double dMinSqDinst = DBL_MAX ;
   for ( int nPt = 0 ; nPt < int( Loop.size()) ; ++ nPt) {
     // Estremi del segmento corrente del loop
      Point3d ptSegSt = Loop[nPt] ;
      Point3d ptSegEn = Loop[( nPt + 1) % int(Loop.size())] ;
     // Distanza del punto dal segmento del loop
      DistPointLine dDistCalc( ptStop, ptSegSt, ptSegEn) ;
      double dSqDist ;
      dDistCalc.GetSqDist( dSqDist) ;
      if ( dSqDist < dMinSqDinst) {
         dMinSqDinst = dSqDist ;
         nMinSeg = nPt ;
      }
   }
  // Se il punto non sta sul loop, errore
   if ( dMinSqDinst > SQ_EPS_SMALL)
      return false ;
  // Verifico che il punto stia su un vertice, in tal caso non devo aggiungerlo
   bool bFirst = AreSamePointApprox( Loop[nMinSeg], ptStop) ;
   bool bLast = AreSamePointApprox( Loop[( nMinSeg + 1) % int( Loop.size())], ptStop) ;
  // Se il punto è sul vertice finale del segmento, aggiungo il vertice alla lista da inglobare al primo loop
   if ( bLast)
      ++ nMinSeg ;
  // Inglobo fino a nSeg nel primo loop
   for ( int nPt = 0 ; nPt <= nMinSeg  &&  nPt < int( Loop.size()) ; ++ nPt)
      Loop1.emplace_back( Loop[nPt]) ;
  // Se il punto è interno al segmento, lo inglobo in entrambi i loop
   if ( ! ( bFirst  ||  bLast)) {
      Loop1.emplace_back( ptStop) ;
      Loop2.emplace_back( ptStop) ;
   }
   else {
      if ( nMinSeg != int( Loop.size()) )
         Loop2.emplace_back( Loop[nMinSeg]) ;
   }
  // Inglobo gli ultimi vertici in Loop2
   for ( int nPt = nMinSeg + 1 ; nPt < int( Loop.size()) ; ++ nPt)
      Loop2.emplace_back( Loop[nPt]) ;
   Loop2.emplace_back( Loop[0]) ;

   return true ;
}