...

[FrontAlgorithms.git] / examples / BronchiiInitialSegmentationWithMori.cxx

diff --git a/examples/BronchiiInitialSegmentationWithMori.cxx b/examples/BronchiiInitialSegmentationWithMori.cxx
index ed0c294132e264d2a11558f45616a870ac89bdd9..fa940d5a1ce5af997fcbf20bbbdfa430f6eb2254 100644 (file)
--- a/examples/BronchiiInitialSegmentationWithMori.cxx
+++ b/examples/BronchiiInitialSegmentationWithMori.cxx
@@ -52,6 +52,7 @@ int main( int argc, char* argv[] )
   } // fi
   std::string input_image_filename = GetFullPath( argv[ 1 ] );
   std::string output_image_filename = argv[ 2 ];
+  std::string output_auxiliary_image_filename = output_image_filename + "_aux.mhd";
 
   // Try to guess initial seed
   std::string seed_filename = GetFullPathToDirectory( argv[ 1 ] ) + "seed.txt";
@@ -66,6 +67,7 @@ int main( int argc, char* argv[] )
   } // fi
   TImage::IndexType seed;
   seed_str >> seed[ 0 ] >> seed[ 1 ] >> seed[ 2 ];
+  seed_str.close( );
 
   // Read image
   typedef itk::ImageFileReader< TImage > TReader;
@@ -77,7 +79,7 @@ int main( int argc, char* argv[] )
   TFilter::Pointer filter = TFilter::New( );
   filter->SetInput( reader->GetOutput( ) );
   filter->SetSeed( seed );
-  filter->SetThresholdRange( -1024, 0 );
+  filter->SetThresholdRange( -200, 100, 1 );
   filter->SetInsideValue( 1 );
   filter->SetOutsideValue( 0 );
 
@@ -98,6 +100,57 @@ int main( int argc, char* argv[] )
     return( 1 );
 
   } // yrt
+
+  // Write auxiliary image
+  typedef itk::ImageFileWriter< TFilter::TAuxImage > TAuxWriter;
+  TAuxWriter::Pointer aux_writer = TAuxWriter::New( );
+  aux_writer->SetInput( filter->GetAuxiliaryImage( ) );
+  aux_writer->SetFileName( output_auxiliary_image_filename );
+
+  // Execute pipeline
+  try
+  {
+    aux_writer->Update( );
+  }
+  catch( std::exception& err )
+  {
+    std::cerr << "Error caught: " << err.what( ) << std::endl;
+    return( 1 );
+
+  } // yrt
+
+  // Write result signal
+  std::string output_signal_filename = output_image_filename + ".csv";
+  std::ofstream output_signal_str( output_signal_filename.c_str( ) );
+  TFilter::TCurve curve = filter->GetCurve( );
+  unsigned long max_count = 0;
+  for( TFilter::TCurve::value_type d : curve )
+  {
+    output_signal_str << d.first << " " << d.second << std::endl;
+    if( max_count < d.second )
+      max_count = d.second;
+
+  } // rof
+  output_signal_str.close( );
+
+  // Write gnuplot script
+  std::string output_gnuplot_filename = output_image_filename + ".gnuplot";
+  std::ofstream output_gnuplot_str( output_gnuplot_filename.c_str( ) );
+  unsigned int thr_pos = filter->GetOptimumThreshold( );
+  int thr = curve[ thr_pos ].first;
+  output_gnuplot_str
+    << "set term png" << std::endl
+    << "set output \"" << output_image_filename << ".png\"" << std::endl
+    << "set arrow 1 from " << thr
+    << ",0 to " << thr << "," << max_count
+    << std::endl
+    << "show arrow 1" << std::endl
+    << "plot \"" << output_signal_filename
+    << "\" using 1:2 with linespoints title \"Evolution ("
+    << thr << "," << thr_pos << ") " << "\""
+    << std::endl;
+  output_gnuplot_str.close( );
+
   return( 0 );
 }