#ifndef __CPM__ALGORITHMS__QUADEDGE__EDGEDECIMATIONFILTER__HXX__
#define __CPM__ALGORITHMS__QUADEDGE__EDGEDECIMATIONFILTER__HXX__

#include <itkTriangleHelper.h>

// -------------------------------------------------------------------------
template< class I, class O, class C >
cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
EdgeDecimationFilter( )
  : Superclass( ),
    m_Relocate( true ),
    m_CheckOrientation( false ),
    m_Element( NULL )
{
  this->m_JoinVertexFunction = TOperator::New( );
  this->m_PriorityQueue = TPriorityQueue::New( );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
~EdgeDecimationFilter( )
{
  TOutPrimalEdge* edge;

  while( !( this->m_PriorityQueue->Empty( ) ) )
  {
    edge = this->m_PriorityQueue->Peek( )->m_Element;
    this->m_PriorityQueue->Pop( );

    typename TQueueMap::iterator it = this->m_QueueMapper.find( edge );
    delete it->second;
    this->m_QueueMapper.erase( it );

  } // elihw
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
FillPriorityQueue( )
{
  typename O::Pointer output = this->GetOutput( );

  this->m_JoinVertexFunction->SetMesh( output );

  /* TODO: fill queue with qedges
     typename O::CellsContainerIterator it = output->GetEdgeCells( )->Begin( );
     typename O::CellsContainerIterator end = output->GetEdgeCells( )->End( );

     OutputEdgeCellType* edge;

     // cache for use in MeasureEdge
     this->GetOutput( ) = this->GetOutput( );

     while( it != end )
     {
     edge = dynamic_cast< OutputEdgeCellType*  >( it.Value( ) );

     if( edge )
     {
     PushElement( edge->GetQEGeom( ) );
     }
     ++it;
     }
 */
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
PushElement( TOutPrimalEdge* iEdge )
{
  TOutPointIdentifier id_org = iEdge->GetOrigin( );
  TOutPointIdentifier id_dest = iEdge->GetDestination( );

  TOutPrimalEdge* temp =( id_org < id_dest ) ? iEdge : iEdge->GetSym( );
  TScalar   measure = MeasureEdge( temp );

  TPriorityQueueItem* qi = new TPriorityQueueItem( temp,
                                                         TPriority( false, measure ) );

  this->m_QueueMapper[ temp ] = qi;
  this->m_PriorityQueue->Push( qi );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
bool cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
IsEdgeOKToBeProcessed( TOutPrimalEdge* iEdge )
{
  return( true );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
Extract( )
{
  typename O::Pointer output = this->GetOutput( );

  do
  {
    this->m_Element = this->m_PriorityQueue->Peek( )->m_Element;
    this->m_Priority = this->m_PriorityQueue->Peek( )->m_Priority;

    this->m_PriorityQueue->Pop( );
    typename TQueueMap::iterator it = this->m_QueueMapper.find( this->m_Element );
    delete it->second;
    this->m_QueueMapper.erase( it );
  }
  while( !IsEdgeOKToBeProcessed( this->m_Element ) );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
DeleteElement( TOutPrimalEdge* iEdge )
{
  if( iEdge ) // this test can be removed
  {
    TOutPrimalEdge* temp =( iEdge->GetOrigin( ) < iEdge->GetDestination( ) ) ?
      iEdge : iEdge->GetSym( );

    typename TQueueMap::iterator map_it = this->m_QueueMapper.find( temp );
    if( map_it != this->m_QueueMapper.end( ) )
    {
      if( !map_it->second->m_Priority.first )
      {
        TPriorityQueueItem* e( map_it->second );
        this->m_PriorityQueue->DeleteElement( e );
        delete map_it->second;
        this->m_QueueMapper.erase( map_it );
      }
    }
  }
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
DeletePoint(
  const TOutPointIdentifier& iIdToBeDeleted,
  const TOutPointIdentifier& iRemaining
  )
{
  /* TODO
     ( void )iRemaining;
     this->GetOutput( )->DeletePoint( iIdToBeDeleted );
  */
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
PushOrUpdateElement( TOutPrimalEdge* iEdge )
{
  TOutPrimalEdge* temp = iEdge;

  if( temp->GetOrigin( ) > temp->GetDestination( ) )
  {
    temp = temp->GetSym( );
  }

  typename TQueueMap::iterator map_it = this->m_QueueMapper.find( temp );

  TScalar measure = MeasureEdge( temp );
  if( map_it != this->m_QueueMapper.end( ) )
  {
    if( !map_it->second->m_Priority.first )
    {
      map_it->second->m_Priority.second = measure;
      this->m_PriorityQueue->Update( map_it->second );
    }
  }
  else
  {
    TPriorityQueueItem* qi = new TPriorityQueueItem( temp,
                                                           TPriority( false, measure ) );
    this->m_QueueMapper[ temp ] = qi;
    this->m_PriorityQueue->Push( qi );
  }
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
JoinVertexFailed( )
{
  typename TOperator::TEdgeStatus
    status = this->m_JoinVertexFunction->GetEdgeStatus( );
  switch( status )
  {
  default:
  case TOperator::EDGE_NULL:
  case TOperator::MESH_NULL:
  case TOperator::FACE_ISOLATED:
    break;
  case TOperator::EDGE_ISOLATED:
    itkDebugMacro( "EDGE_ISOLATED, at iteration: " << this->m_Iteration );
    TagElementOut( this->m_Element );
    break;
    // more than 2 common vertices in 0-ring of org and dest respectively
  case TOperator::TOO_MANY_COMMON_VERTICES:
    itkDebugMacro( "TOO_MANY_COMMON_VERTICES, at iteration "
                  << this->m_Iteration );
    itkDebugMacro( << this->m_Element->GetOrigin( ) << " -> "
                   << this->m_Element->GetDestination( ) );
    this->TagElementOut( this->m_Element );
    break;
    // ******************************************************************
    // Tetrahedron case
  case TOperator::TETRAHEDRON_CONFIG:
    itkDebugMacro( "TETRAHEDRON_CONFIG, at iteration " << this->m_Iteration );

    this->TagElementOut( this->m_Element );
    this->TagElementOut( this->m_Element->GetOnext( ) );
    this->TagElementOut( this->m_Element->GetOprev( ) );
    this->TagElementOut( this->m_Element->GetSym( ) );
    this->TagElementOut( this->m_Element->GetSym( )->GetOnext( ) );
    this->TagElementOut( this->m_Element->GetSym( )->GetOprev( ) );
    this->TagElementOut( this->m_Element->GetOnext( )->GetLnext( ) );
    break;
    // ******************************************************************
    // Samosa case
  case TOperator::SAMOSA_CONFIG:
    itkDebugMacro( "SAMOSA_CONFIG, at iteration " << this->m_Iteration );
    this->RemoveSamosa( );
    break;
    // ******************************************************************
    // Eye case
  case TOperator::EYE_CONFIG:
    itkDebugMacro( "EYE_CONFIG, at iteration " << this->m_Iteration );
    this->RemoveEye( );
    break;
  case TOperator::EDGE_JOINING_DIFFERENT_BORDERS:
    itkDebugMacro( "EDGE_JOINING_DIFFERENT_BORDERS, at iteration "
                  << this->m_Iteration );
    this->TagElementOut( this->m_Element );
    break;
  }
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
bool cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
ProcessWithoutAnyTopologicalGuarantee( )
{
  TOutPoint   pt;

  TOutPointIdentifier id_org = this->m_Element->GetOrigin( );
  TOutPointIdentifier id_dest = this->m_Element->GetDestination( );
  TOutPointIdentifier idx =( id_org < id_dest ) ? id_org : id_dest;

  bool to_be_processed( true );

  if( this->m_Relocate )
  {
    pt = Relocate( this->m_Element );
  }
  else
  {
    pt = this->GetOutput( )->GetPoint( idx );
  }

///TODO use CheckOrientation!!!
//   if( this->m_CheckOrientation )
//     to_be_processed = CheckOrientation( this->m_Element, idx, pt );

    if( !to_be_processed )
    {
      return( false );
    }

    std::list< TOutPrimalEdge*  > list_qe_to_be_deleted;
    TOutPrimalEdge*               temp = this->m_Element->GetOnext( );

    while( temp != this->m_Element )
    {
      list_qe_to_be_deleted.push_back( temp );
      temp = temp->GetOnext( );
    }

    temp = this->m_Element->GetSym( )->GetOnext( );
    while( temp != this->m_Element->GetSym( ) )
    {
      list_qe_to_be_deleted.push_back( temp );
      temp = temp->GetOnext( );
    }

    typename std::list< TOutPrimalEdge*  >::iterator
      it = list_qe_to_be_deleted.begin( );

    while( it != list_qe_to_be_deleted.end( ) )
    {
      DeleteElement( *it );
      ++it;
    }

    if( !this->m_JoinVertexFunction->Evaluate( this->m_Element ) )
    {
      it = list_qe_to_be_deleted.begin( );

      while( it != list_qe_to_be_deleted.end( ) )
      {
        PushOrUpdateElement( *it );
        ++it;
      }

      JoinVertexFailed( );
    }
    else
    {
      TOutPointIdentifier old_id = this->m_JoinVertexFunction->GetOldPointID( );

      TOutPointIdentifier new_id =( old_id == id_dest ) ? id_org : id_dest;
      DeletePoint( old_id, new_id );

      TOutPrimalEdge* edge = this->GetOutput( )->FindEdge( new_id );
      if( edge == ITK_NULLPTR )
      {
        itkDebugMacro( "edge == 0, at iteration " << this->m_Iteration );
        return( false );
      }

      if( this->m_Relocate )
      {
        this->GetOutput( )->SetPoint( new_id, pt );
      }

      temp = edge;

      do
      {
        PushOrUpdateElement( temp );
        temp = temp->GetOnext( );
      }
      while( temp != edge );
    }
    return( false );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
unsigned int cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
CheckQEProcessingStatus( )
{
  TOutPrimalEdge* qe = this->m_Element;
  TOutPrimalEdge* qe_sym = qe->GetSym( );

  /* TODO
     bool LeftIsTriangle = qe->IsLnextOfTriangle( );
     bool RightIsTriangle = qe->GetSym( )->IsLnextOfTriangle( );
  */
  bool LeftIsTriangle = false;
  bool RightIsTriangle = false;

  if( LeftIsTriangle || RightIsTriangle )
  {
    if( LeftIsTriangle && RightIsTriangle )
    {
      // two triangles
      bool OriginOrderIsTwo;// TODO: =( qe->GetOrder( ) == 2 );
      bool DestinationOrderIsTwo;// TODO: =( qe_sym->GetOrder( ) == 2 );

      if( OriginOrderIsTwo || DestinationOrderIsTwo )
      {
        if( OriginOrderIsTwo && DestinationOrderIsTwo )
        {
          // isolated component made of two triangles
          // sharing same points but with opposite orientation
          // looks like a samosa
          itkDebugMacro( "RemoveSamosa" );
          return( 1 );
        } // end if( OriginOrderIsTwo && DestinationOrderIsTwo )
        else
        {
          // two triangles share three points and two edges
          // the last edge is duplicated = two edge cells
          // having the same points. It is a valid manifold case
          // but you have to decimate it the right way.
          // from the top the drawing of that case looks like an Eye
          itkDebugMacro( "RemoveEye" );
          return( 2 );
        } // end else if( OriginOrderIsTwo && DestinationOrderIsTwo )
      }   // end if( OriginOrderIsTwo || DestinationOrderIsTwo )
      else  // if( OriginOrderIsTwo || DestinationOrderIsTwo )
      {
        if( NumberOfCommonVerticesIn0Ring( ) > 2 )
        {
          // both points have more than 2 edges on their O-ring
          itkDebugMacro( "NumberOfCommonVerticesIn0Ring( ) > 2" );
          return( 3 );
        } //end if( NumberOfCommonVerticesIn0Ring( ) > 2 )
        else
        {
          return( 0 );
        }
      } // end else if( OriginOrderIsTwo || DestinationOrderIsTwo )
    }   // end if( LeftIsTriangle && RightIsTriangle )
    else  // if( LeftIsTriangle && RightIsTriangle )
    {
      if( NumberOfCommonVerticesIn0Ring( ) > 1 )
      {
        itkDebugMacro( "NumberOfCommonVerticesIn0Ring( ) > 1" );
        return( 4 );
      }  // end if( NumberOfCommonVerticesIn0Ring( ) > 1 )
      else // if( NumberOfCommonVerticesIn0Ring( ) > 1 )
      {
        if( RightIsTriangle )
        {
          return( 5 );
        }
        else
        {
          return( 6 );
        }
      } // end else if( NumberOfCommonVerticesIn0Ring( ) > 1 )
    }   // end else if( LeftIsTriangle && RightIsTriangle )
  }     // end if( LeftIsTriangle || RightIsTriangle )
  else    // if( LeftIsTriangle || RightIsTriangle )
  {
    if( NumberOfCommonVerticesIn0Ring( ) > 0 )
    {
      return( 7 );
    }
    else
    {
      return( 0 );
    }
  } // end if( LeftIsTriangle || RightIsTriangle )

  //   return( 0 );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
bool cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
ProcessWithTopologicalGuarantee( )
{
  if( this->m_Priority.first )
  {
    return( true );
  }

  ProcessWithoutAnyTopologicalGuarantee( );
  return( false );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
unsigned long cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
NumberOfCommonVerticesIn0Ring( ) const
{
  TOutPrimalEdge* qe = this->m_Element;
  TOutPrimalEdge* e_it  = qe->GetOnext( );

  std::list< TOutPointIdentifier > dir_list, sym_list, intersection_list;
  do
  {
    dir_list.push_back( e_it->GetDestination( ) );
    e_it = e_it->GetOnext( );
  }
  while( e_it != qe );

  qe = qe->GetSym( );
  e_it = qe;

  do
  {
    sym_list.push_back( e_it->GetDestination( ) );
    e_it = e_it->GetOnext( );
  }
  while( e_it != qe );

  dir_list.sort( );
  sym_list.sort( );

  std::set_intersection( dir_list.begin( ), dir_list.end( ),
                         sym_list.begin( ), sym_list.end( ),
                         std::back_inserter( intersection_list ) );

  return( static_cast< unsigned long >( intersection_list.size( ) ) );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
RemoveSamosa( )
{
  DeleteElement( this->m_Element->GetLnext( ) );
  DeleteElement( this->m_Element->GetLprev( ) );
  DeleteElement( this->m_Element->GetRnext( ) );
  DeleteElement( this->m_Element->GetRprev( ) );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
TagElementOut( TOutPrimalEdge* iEdge )
{
  typename TQueueMap::iterator map_it = this->m_QueueMapper.find( iEdge );

  if( map_it != this->m_QueueMapper.end( ) )
  {
    map_it->second->m_Priority.first = true;
    map_it->second->m_Priority.second = static_cast< TScalar >( 0. );
    this->m_PriorityQueue->Update( map_it->second );
  }
  else
  {
    TPriorityQueueItem* qi = new TPriorityQueueItem( iEdge,
                                                     TPriority( true, static_cast< TScalar >( 0. ) ) );

    this->m_QueueMapper[ iEdge ] = qi;
    this->m_PriorityQueue->Push( qi );
  }
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
void cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
RemoveEye( )
{
  TOutPrimalEdge* qe = this->m_Element;
  TOutPrimalEdge* qe_sym = this->m_Element->GetSym( );

  /* TODO
     if( qe->GetSym( )->GetOrder( ) == 2 )
     {
     qe = qe_sym;
     }
  */

  TagElementOut( qe );
  TagElementOut( qe->GetOnext( ) );
  TagElementOut( qe->GetSym( )->GetOnext( ) );
  TagElementOut( qe->GetSym( )->GetOprev( ) );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
bool cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
CheckOrientation( TOutPrimalEdge* iEdge,
                 const TOutPointIdentifier& iId,
                 const TOutPoint& iPt )
{
  typename O::Pointer           output = this->GetOutput( );
  typename O::CellsContainerPointer cells = output->GetCells( );

  std::list< typename O::CellIdentifier > r1, r2, elements_to_be_tested;
  TOutPrimalEdge*                     qe = iEdge;
  TOutPrimalEdge*                     qe_it = qe->GetOnext( );

  do
  {
    r1.push_back( qe_it->GetLeft( ) );
    qe_it = qe_it->GetOnext( );
  }
  while( qe_it != qe );

  qe = iEdge->GetSym( );
  qe_it = qe->GetOnext( );

  do
  {
    r2.push_back( qe_it->GetLeft( ) );
    qe_it = qe_it->GetOnext( );
  }
  while( qe_it != qe );

  r1.sort( );
  r2.sort( );

  std::set_symmetric_difference( r1.begin( ), r1.end( ),
                                 r2.begin( ), r2.end( ),
                                 std::back_inserter( elements_to_be_tested ) );

  typename std::list< typename O::CellIdentifier >::iterator
    it = elements_to_be_tested.begin( );

  typedef itk::TriangleHelper< TOutPoint > TTriangle;

  bool                  orientation_ok( true );
  typename O::CellIdentifier  c_id( 0 );

  /* TODO: iterate over polygons
  OutputPolygonType*    poly;
  TOutPointIdentifier p_id;

  int             k( 0 ), replace_k( 0 );
  TOutPoint pt[ 3 ];

  OutputVectorType n_bef, n_aft;

  while( ( it != elements_to_be_tested.end( ) ) && orientation_ok )
  {
    c_id =* it;
    poly = dynamic_cast< OutputPolygonType*  >( cells->GetElement( c_id ) );
    qe = poly->GetEdgeRingEntry( );
    qe_it = qe;
    k = 0;
    do
    {
      p_id = qe_it->GetOrigin( );
      if( p_id == iId )
      {
        replace_k = k;
      }
      pt[ k++ ] = output->GetPoint( p_id );
      qe_it = qe_it->GetLnext( );
    }
    while( qe_it != qe );

    n_bef = TTriangle::ComputeNormal( pt[ 0 ], pt[ 1 ], pt[ 2 ] );
    switch( replace_k )
    {
    default:
    case 0:
      n_aft = TTriangle::ComputeNormal( iPt, pt[ 1 ], pt[ 2 ] );
      break;
    case 1:
      n_aft = TTriangle::ComputeNormal( pt[ 0 ], iPt, pt[ 2 ] );
      break;
    case 2:
      n_aft = TTriangle::ComputeNormal( pt[ 0 ], pt[ 1 ], iPt );
      break;
    }

    orientation_ok =( n_bef * n_aft ) < 0.;
    ++it;
  }
  */
  return( orientation_ok );
}

// -------------------------------------------------------------------------
template< class I, class O, class C >
bool cpm::Algorithms::QuadEdge::EdgeDecimationFilter< I, O, C >::
IsCriterionSatisfied( )
{
  if( this->m_PriorityQueue->Empty( ) )
  {
    return( true );
  }
  else
  {
    return( this->m_Criterion->is_satisfied( this->GetOutput( ), 0, this->m_Priority.second ) );
  }
}

#endif // __CPM__ALGORITHMS__QUADEDGE__EDGEDECIMATIONFILTER__HXX__

// eof - $RCSfile$