X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=appli%2Fexamples%2Fexample_Image_Dijkstra_EndPointDetection.cxx;fp=appli%2Fexamples%2Fexample_Image_Dijkstra_EndPointDetection.cxx;h=b2d3e7c6c8c21529126619648a152a2f99e88176;hb=5e326afe442245572b6c3ec98ebeec8b45f9012f;hp=47585ad30327975306210ad1eb51fe2a5cd912c8;hpb=8480ad8f9b3e3b556aa8f3f25ff07daf1db6dc78;p=FrontAlgorithms.git

diff --git a/appli/examples/example_Image_Dijkstra_EndPointDetection.cxx b/appli/examples/example_Image_Dijkstra_EndPointDetection.cxx
index 47585ad..b2d3e7c 100644
--- a/appli/examples/example_Image_Dijkstra_EndPointDetection.cxx
+++ b/appli/examples/example_Image_Dijkstra_EndPointDetection.cxx
@@ -20,26 +20,11 @@
 #include <fpa/VTK/ImageMPR.h>
 #include <fpa/VTK/Image3DObserver.h>
 
-/*
-  #include <itkMinimumMaximumImageCalculator.h>
-  #include <itkInvertIntensityImageFilter.h>
-  #include <itkDanielssonDistanceMapImageFilter.h>
-  #include <itkImageFileWriter.h>
-
-  #include <vtkCamera.h>
-  #include <vtkImageActor.h>
-  #include <vtkInteractorStyleImage.h>
-  #include <vtkPointHandleRepresentation3D.h>
-  #include <vtkProperty.h>
-  #include <vtkRenderer.h>
-  #include <vtkRenderWindow.h>
-  #include <vtkRenderWindowInteractor.h>
-  #include <vtkSeedRepresentation.h>
-  #include <vtkSeedWidget.h>
-  #include <vtkSmartPointer.h>
-
-  #include <fpa/Image/Dijkstra.h>
-*/
+#include <vtkCellArray.h>
+#include <vtkImageMarchingCubes.h>
+#include <vtkPoints.h>
+#include <vtkPolyData.h>
+#include <vtkSmartPointer.h>
 
 // -------------------------------------------------------------------------
 const unsigned int Dim = 3;
@@ -54,6 +39,9 @@ typedef itk::ImageToVTKImageFilter< TImage > TVTKInputImage;
 template< class I >
 void ReadImage( typename I::Pointer& image, const std::string& filename );
 
+template< class I >
+void SaveImage( const I* image, const std::string& filename );
+
 template< class I, class O >
 void DistanceMap(
   const typename I::Pointer& input, typename O::Pointer& output
@@ -62,16 +50,19 @@ void DistanceMap(
 // -------------------------------------------------------------------------
 int main( int argc, char* argv[] )
 {
-  if( argc < 2 )
+  if( argc < 5 )
   {
     std::cerr
       << "Usage: " << argv[ 0 ]
-      << " input_image"
+      << " input_image distancemap output_costmap output_labels"
       << std::endl;
     return( 1 );
 
   } // fi
   std::string input_image_fn = argv[ 1 ];
+  std::string distancemap_fn = argv[ 2 ];
+  std::string output_costmap_fn = argv[ 3 ];
+  std::string output_labels_fn = argv[ 4 ];
 
   // Read image
   TImage::Pointer input_image;
@@ -88,6 +79,7 @@ int main( int argc, char* argv[] )
     return( 1 );
 
   } // yrt
+
   TVTKInputImage::Pointer vtk_input_image = TVTKInputImage::New( );
   vtk_input_image->SetInput( input_image );
   vtk_input_image->Update( );
@@ -98,6 +90,13 @@ int main( int argc, char* argv[] )
   view.SetSize( 800, 800 );
   view.SetImage( vtk_input_image->GetOutput( ) );
 
+  vtkSmartPointer< vtkImageMarchingCubes > mc =
+    vtkSmartPointer< vtkImageMarchingCubes >::New( );
+  mc->SetInputData( vtk_input_image->GetOutput( ) );
+  mc->SetValue( 0, 1e-1 );
+  mc->Update( );
+  view.AddPolyData( mc->GetOutput( ), 1, 1, 1, 0.4 );
+
   // Allow some interaction and wait for, at least, one seed
   view.Render( );
   while( view.GetNumberOfSeeds( ) == 0 )
@@ -148,25 +147,142 @@ int main( int argc, char* argv[] )
     << std::difftime( end, start )
     << " s." << std::endl;
 
-  /* TODO
-  // Save final total cost map
-  itk::ImageFileWriter< TOutputImage >::Pointer output_image_writer =
-    itk::ImageFileWriter< TOutputImage >::New( );
-  output_image_writer->SetFileName( output_image_fn );
-  output_image_writer->SetInput( filter->GetOutput( ) );
-  try
+  // Minimum spanning tree
+  const TFilter::TMinimumSpanningTree* mst =
+    filter->GetMinimumSpanningTree( );
+
+  // Build branches from endpoints to seed
+  vtkSmartPointer< vtkPoints > endpoints_points =
+    vtkSmartPointer< vtkPoints >::New( );
+  vtkSmartPointer< vtkCellArray > endpoints_cells =
+    vtkSmartPointer< vtkCellArray >::New( );
+  vtkSmartPointer< vtkCellArray > endpoints_vertices =
+    vtkSmartPointer< vtkCellArray >::New( );
+
+  const TFilter::TUniqueVertices& endpoints = filter->GetEndPoints( );
+  TFilter::TUniqueVertices::const_iterator eIt = endpoints.begin( );
+  for( ; eIt != endpoints.end( ); ++eIt )
   {
-    output_image_writer->Update( );
-  }
-  catch( itk::ExceptionObject& err )
+    TFilter::TVertices path;
+    mst->GetPath( path, *eIt, filter->GetSeed( 0 ) );
+
+    std::cout << *eIt << std::endl;
+    TFilter::TVertices::const_iterator pIt = path.begin( );
+    for( ; pIt != path.end( ); ++pIt )
+    {
+      TImage::PointType pnt;
+      input_image->TransformIndexToPhysicalPoint( *pIt, pnt );
+      endpoints_points->InsertNextPoint( pnt[ 0 ], pnt[ 1 ], pnt[ 2 ] );
+
+      if( pIt != path.begin( ) )
+      {
+        endpoints_cells->InsertNextCell( 2 );
+        endpoints_cells->InsertCellPoint( endpoints_points->GetNumberOfPoints( ) - 2 );
+        endpoints_cells->InsertCellPoint( endpoints_points->GetNumberOfPoints( ) - 1 );
+      }
+      else
+      {
+        std::cout << "ok" << std::endl;
+        endpoints_vertices->InsertNextCell( 1 );
+        endpoints_vertices->InsertCellPoint( endpoints_points->GetNumberOfPoints( ) - 1 );
+
+      } // fi
+
+
+
+
+    } // rof
+
+  } // rof
+  std::cout << endpoints.size( ) << std::endl;
+      
+  vtkSmartPointer< vtkPolyData > endpoints_polydata =
+    vtkSmartPointer< vtkPolyData >::New( );
+  endpoints_polydata->SetPoints( endpoints_points );
+  endpoints_polydata->SetLines( endpoints_cells );
+  endpoints_polydata->SetVerts( endpoints_vertices );
+
+
+
+
+
+
+
+
+
+
+
+  // Bifurcations
+  vtkSmartPointer< vtkPoints > bifurcations_points =
+    vtkSmartPointer< vtkPoints >::New( );
+  vtkSmartPointer< vtkCellArray > bifurcations_vertices =
+    vtkSmartPointer< vtkCellArray >::New( );
+
+  const TFilter::TUniqueVertices& bifurcations = filter->GetBifurcationPoints( );
+  TFilter::TUniqueVertices::const_iterator bIt = bifurcations.begin( );
+  for( ; bIt != bifurcations.end( ); ++bIt )
   {
-    std::cerr
-      << "Error while writing image to " << output_image_fn << ": "
-      << err << std::endl;
-    return( 1 );
+    std::cout << *bIt << std::endl;
 
-  } // yrt
-  */
+    TImage::PointType pnt;
+    input_image->TransformIndexToPhysicalPoint( *bIt, pnt );
+    bifurcations_points->InsertNextPoint( pnt[ 0 ], pnt[ 1 ], pnt[ 2 ] );
+    bifurcations_vertices->InsertNextCell( 1 );
+    bifurcations_vertices->InsertCellPoint( bifurcations_points->GetNumberOfPoints( ) - 1 );
+
+  } // rof
+
+  vtkSmartPointer< vtkPolyData > bifurcations_polydata =
+    vtkSmartPointer< vtkPolyData >::New( );
+  bifurcations_polydata->SetPoints( bifurcations_points );
+  bifurcations_polydata->SetVerts( bifurcations_vertices );
+
+
+
+
+
+  // Build branches from endpoints to seed
+  vtkSmartPointer< vtkPoints > branches_points =
+    vtkSmartPointer< vtkPoints >::New( );
+  vtkSmartPointer< vtkCellArray > branches_cells =
+    vtkSmartPointer< vtkCellArray >::New( );
+
+  const TFilter::TBranches& branches = filter->GetBranches( );
+  TFilter::TBranches::const_iterator brIt = branches.begin( );
+  for( ; brIt != branches.end( ); ++brIt )
+  {
+    TFilter::TBranch::const_iterator brsIt = brIt->second.begin( );
+    for( ; brsIt != brIt->second.end( ); ++brsIt )
+    {
+      TImage::PointType pnt0, pnt1;
+      input_image->TransformIndexToPhysicalPoint( brIt->first, pnt0 );
+      input_image->TransformIndexToPhysicalPoint( brsIt->first, pnt1 );
+      branches_points->InsertNextPoint( pnt0[ 0 ], pnt0[ 1 ], pnt0[ 2 ] );
+      branches_points->InsertNextPoint( pnt1[ 0 ], pnt1[ 1 ], pnt1[ 2 ] );
+
+      branches_cells->InsertNextCell( 2 );
+      branches_cells->InsertCellPoint( branches_points->GetNumberOfPoints( ) - 2 );
+      branches_cells->InsertCellPoint( branches_points->GetNumberOfPoints( ) - 1 );
+
+    } // rof
+
+  } // rof
+
+  vtkSmartPointer< vtkPolyData > branches_polydata =
+    vtkSmartPointer< vtkPolyData >::New( );
+  branches_polydata->SetPoints( branches_points );
+  branches_polydata->SetLines( branches_cells );
+
+  view.AddPolyData( endpoints_polydata, 1, 0, 0, 1 );
+  view.AddPolyData( bifurcations_polydata, 0, 0, 1, 1 );
+  view.AddPolyData( branches_polydata, 0, 1, 0, 1 );
+  view.Render( );
+  view.Start( );
+
+  // Save final total cost map
+  SaveImage( filter->GetOutput( ), output_costmap_fn );
+  SaveImage( dmap.GetPointer( ), distancemap_fn );
+  SaveImage( filter->GetLabelImage( ), output_labels_fn );
   return( 0 );
 }
 
@@ -197,6 +313,27 @@ void ReadImage( typename I::Pointer& image, const std::string& filename )
   image->DisconnectPipeline( );
 }
 
+// -------------------------------------------------------------------------
+template< class I >
+void SaveImage( const I* image, const std::string& filename )
+{
+  typename itk::ImageFileWriter< I >::Pointer writer =
+    itk::ImageFileWriter< I >::New( );
+  writer->SetInput( image );
+  writer->SetFileName( filename );
+  try
+  {
+    writer->Update( );
+  }
+  catch( itk::ExceptionObject& err )
+  {
+    std::cerr
+      << "Error saving \"" << filename << "\": " << err
+      << std::endl;
+
+  } // yrt
+}
+
 // -------------------------------------------------------------------------
 template< class I, class O >
 void DistanceMap(
@@ -225,3 +362,13 @@ void DistanceMap(
 }
 
 // eof - $RCSfile$
+
+
+
+
+
+
+
+
+
+