#include <chrono>
#include <iomanip>
#include <itkImage.h>
#include <itkImageFileReader.h>

#include <fpa/Base/SkeletonWriter.h>
#include <fpa/Image/SkeletonFilter.h>

// -------------------------------------------------------------------------
static const unsigned int VDim = 3;
typedef unsigned char                                   TPixel;
typedef double                                          TScalar;
typedef itk::Image< TPixel, VDim >                      TImage;
typedef itk::ImageFileReader< TImage >                  TReader;
typedef fpa::Image::SkeletonFilter< TImage, TScalar >   TFilter;
typedef fpa::Base::SkeletonWriter< TFilter::TSkeleton > TWriter;

// -------------------------------------------------------------------------
int main( int argc, char* argv[] )
{
  // Get arguments
  if( argc != 3 && argc != 3 + VDim )
  {
    std::cerr
      << "Usage: " << argv[ 0 ]
      << " input_image output_skeleton [";
    for( unsigned int i = 0; i < VDim; ++i )
      std::cerr << " s_" << i;
    std::cerr << " ]" << std::endl;
    return( 1 );

  } // fi
  std::string input_image_filename = argv[ 1 ];
  std::string output_skeleton_filename = argv[ 2 ];

  TReader::Pointer reader = TReader::New( );
  reader->SetFileName( input_image_filename );
  try
  {
    reader->Update( );
  }
  catch( std::exception& err )
  {
    std::cerr << "ERROR: " << err.what( ) << std::endl;
    return( 1 );

  } // yrt

  TFilter::Pointer filter = TFilter::New( );
  filter->SetInput( reader->GetOutput( ) );

  if( argc > 3 )
  {
    filter->SeedFromMaximumDistanceOff( );
    TImage::PointType pnt;
    for( int i = 0; i < VDim; ++i )
      pnt[ i ] = std::atof( argv[ i + 3 ] );
    TImage::IndexType seed;
    if( !( reader->GetOutput( )->TransformPhysicalPointToIndex( pnt, seed ) ) )
    {
      std::cerr << "ERROR: seed outside image." << std::endl;
      return( 1 );

    } // fi
    filter->AddSeed( seed );
  }
  else
    filter->SeedFromMaximumDistanceOn( );

  std::chrono::time_point< std::chrono::system_clock > start, end;
  start = std::chrono::system_clock::now( );
  filter->Update( );
  end = std::chrono::system_clock::now( );
  std::chrono::duration< double > elapsed_seconds = end - start;

  TWriter::Pointer writer = TWriter::New( );
  writer->SetInput( filter->GetSkeleton( ) );
  writer->SetFileName( output_skeleton_filename );
  try
  {
    writer->Update( );
  }
  catch( std::exception& err )
  {
    std::cerr << "ERROR: " << err.what( ) << std::endl;
    return( 1 );

  } // yrt

  // Show data
  std::cout << "Time: " << elapsed_seconds.count( ) << "s" << std::endl;
  return( 0 );
}

// eof - $RCSfile$