1 // =========================================================================
2 // @author Leonardo Florez Valencia
3 // @email florez-l@javeriana.edu.co
4 // =========================================================================
9 #include <itkImageFileReader.h>
10 #include <itkImageFileWriter.h>
11 #include <fpa/Filters/Image/Mori.h>
13 // -------------------------------------------------------------------------
14 const unsigned int Dim = 3;
16 typedef unsigned char TLabel;
17 typedef itk::Image< TPixel, Dim > TImage;
18 typedef itk::Image< TLabel, Dim > TLabels;
20 // -------------------------------------------------------------------------
21 double MeasureTime( itk::ProcessObject* f )
23 std::chrono::time_point< std::chrono::high_resolution_clock > s, e;
24 std::chrono::duration< double > t;
25 s = std::chrono::high_resolution_clock::now( );
27 e = std::chrono::high_resolution_clock::now( );
32 // -------------------------------------------------------------------------
33 template< class _TImagePtr >
34 void ReadImage( _TImagePtr& image, const std::string& fname )
36 typedef typename _TImagePtr::ObjectType _TImage;
37 typedef itk::ImageFileReader< _TImage > _TReader;
38 typename _TReader::Pointer reader = _TReader::New( );
39 reader->SetFileName( fname );
40 double t = MeasureTime( reader );
41 std::cout << "Read " << fname << " in " << t << " s" << std::endl;
42 image = reader->GetOutput( );
43 image->DisconnectPipeline( );
46 // -------------------------------------------------------------------------
47 template< class _TImage >
48 void WriteImage( const _TImage*image, const std::string& fname )
50 typedef itk::ImageFileWriter< _TImage > _TWriter;
51 typename _TWriter::Pointer writer = _TWriter::New( );
52 writer->SetFileName( fname );
53 writer->SetInput( image );
54 double t = MeasureTime( writer );
55 std::cout << "Wrote " << fname << " in " << t << " s" << std::endl;
58 // -------------------------------------------------------------------------
60 std::map< std::string, std::string >& args, int argc, char* argv[]
63 std::set< std::string > mandatory;
64 mandatory.insert( "in" );
65 mandatory.insert( "out" );
66 mandatory.insert( "labels" );
68 args[ "upper_threshold" ] = "-600";
73 std::string cmd = argv[ i ] + 1;
74 args[ cmd ] = argv[ i + 1 ];
80 for( std::string t: mandatory )
81 complete &= ( args.find( t ) != args.end( ) );
86 << "Usage: " << argv[ 0 ] << std::endl
87 << "\t-in filename" << std::endl
88 << "\t-out filename" << std::endl
89 << "\t-labels filename" << std::endl
90 << "\t[-upper_threshold value]" << std::endl;
97 // -------------------------------------------------------------------------
98 int main( int argc, char* argv[] )
100 std::map< std::string, std::string > args;
103 if( ParseArgs( args, argc, argv ) )
106 TImage::Pointer input_image;
107 ReadImage( input_image, args[ "in" ] );
110 TLabels::Pointer input_labels;
111 ReadImage( input_labels, args[ "labels" ] );
115 typedef fpa::Filters::Image::Mori< TImage, TLabels > TMori;
116 TMori::Pointer mori = TMori::New( );
117 mori->SetInput( input_image );
118 mori->SetSeed( seed );
119 mori->SetInsideValue( 1 );
120 mori->SetOutsideValue( 0 );
121 mori->SetMinimumThreshold(
122 TPixel( std::atof( args[ "minimum_threshold" ].c_str( ) ) )
124 mori->SetSignalKernelSize(
125 std::atoi( args[ "signal_kernel_size" ].c_str( ) )
127 mori->SetSignalThreshold(
128 std::atof( args[ "signal_threshold" ].c_str( ) )
130 mori->SetSignalInfluence(
131 std::atof( args[ "signal_influence" ].c_str( ) )
134 TPixel( std::atof( args[ "lower_threshold" ].c_str( ) ) ),
135 TPixel( std::atof( args[ "upper_threshold" ].c_str( ) ) ),
136 TPixel( std::atof( args[ "delta_threshold" ].c_str( ) ) )
138 double t = MeasureTime( mori );
139 std::cout << "Mori executed in " << t << " s" << std::endl;
140 WriteImage( mori->GetOutput( ), args[ "out" ] );
142 std::map< std::string, std::string >::const_iterator i =
143 args.find( "out_signal" );
144 if( i != args.end( ) )
146 std::stringstream signal;
147 unsigned long nthr = mori->GetNumberOfEvaluatedThresholds( );
148 signal << "# nThr = " << nthr << std::endl;
149 signal << "# Opt = " << mori->GetOptimumThreshold( ) << std::endl;
150 for( unsigned long j = 0; j < nthr; ++j )
152 typename TMori::TPeak p;
154 mori->GetSignalValues( j, x, y, p );
155 signal << x << " " << y << std::endl;
158 std::ofstream signals_str( i->second.c_str( ) );
159 signals_str << signal.str( );
160 signals_str.close( );
168 catch( std::exception& err )
171 << "===============================" << std::endl
172 << "Error caught: " << std::endl
174 << "===============================" << std::endl