X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=appli%2FCTBronchi%2FMoriSegmentation.cxx;h=90a1bf3333b608aa93392ab9210a82bde94310db;hb=9817c556a0b8b5e3b332d45f07faa84d91afb2d0;hp=39b9e0dcb3bcec9978554fa19e2c6bdec5a8f35d;hpb=e977d1bb57ca61759a1fb1a1d0d8937718ef8722;p=FrontAlgorithms.git diff --git a/appli/CTBronchi/MoriSegmentation.cxx b/appli/CTBronchi/MoriSegmentation.cxx index 39b9e0d..90a1bf3 100644 --- a/appli/CTBronchi/MoriSegmentation.cxx +++ b/appli/CTBronchi/MoriSegmentation.cxx @@ -1,137 +1,133 @@ -#include +// ========================================================================= +// @author Leonardo Florez Valencia +// @email florez-l@javeriana.edu.co +// ========================================================================= + +#include +#include +#include #include -#include -#include -#include +#include +#include "Functions.h" // ------------------------------------------------------------------------- const unsigned int Dim = 3; -typedef short TPixel; -typedef unsigned short TLabel; - -typedef itk::Image< TPixel, Dim > TInputImage; -typedef itk::Image< TLabel, Dim > TLabelImage; -typedef fpa::Image::Mori< TInputImage, TLabelImage > TFilter; +typedef short TPixel; +typedef unsigned char TLabel; +typedef itk::Image< TPixel, Dim > TImage; +typedef itk::Image< TLabel, Dim > TLabels; // ------------------------------------------------------------------------- int main( int argc, char* argv[] ) { - // Get arguments - if( argc < 16 ) - { - std::cerr - << "Usage: " << argv[ 0 ] << std::endl - << " input_image output_image output_signal" << std::endl - << " init_threshold end_threshold delta" << std::endl - << " minimum_threshold" << std::endl - << " inside_value outside_value" << std::endl - << " signal_kernel_size signal_threshold signal_influence" - << std::endl - << " seed_x seed_y seed_z" - << std::endl; - return( 1 ); - - } // fi - std::string input_image_filename = argv[ 1 ]; - std::string output_image_filename = argv[ 2 ]; - std::string output_signal_filename = argv[ 3 ]; - TPixel init_threshold = std::atoi( argv[ 4 ] ); - TPixel end_threshold = std::atoi( argv[ 5 ] ); - TPixel delta = std::atoi( argv[ 6 ] ); - TPixel minimum_threshold = std::atoi( argv[ 7 ] ); - TLabel inside_value = std::atoi( argv[ 8 ] ); - TLabel outside_value = std::atoi( argv[ 9 ] ); - unsigned long signal_kernel_size = std::atoi( argv[ 10 ] ); - double signal_threshold = std::atof( argv[ 11 ] ); - double signal_influence = std::atof( argv[ 12 ] ); - TInputImage::PointType pseed; - for( int i = 0; i < Dim; ++i ) - pseed[ i ] = std::atof( argv[ 13 + i ] ); - - // Read image - itk::ImageFileReader< TInputImage >::Pointer input_image_reader = - itk::ImageFileReader< TInputImage >::New( ); - input_image_reader->SetFileName( input_image_filename ); - - // Prepare filter - TFilter::Pointer filter = TFilter::New( ); - filter->SetInput( input_image_reader->GetOutput( ) ); - filter->SetSeed( pseed ); - filter->SetThresholds( init_threshold, end_threshold, delta ); - filter->SetMinimumThreshold( minimum_threshold ); - filter->SetInsideValue( inside_value ); - filter->SetOutsideValue( outside_value ); - filter->SetSignalKernelSize( signal_kernel_size ); - filter->SetSignalThreshold( signal_threshold ); - filter->SetSignalInfluence( signal_influence ); - - // Show some information - std::cout << "----------------------------------------------" << std::endl; - std::cout << "Image: " << input_image_filename << std::endl; - - // Execute pipeline - std::chrono::time_point< std::chrono::high_resolution_clock > ts, te; - std::chrono::duration< double > tr; + typedef TCLAP::ValueArg< std::string > _TStringArg; + typedef TCLAP::ValueArg< double > _TRealArg; + typedef TCLAP::ValueArg< unsigned long > _TUIntArg; + typedef TCLAP::ValueArg< TPixel > _TPixelArg; + + // Parse input line + _TStringArg in( "i", "input", "Input image", true, "", "file" ); + _TStringArg out( "o", "output", "Output image", true, "", "file" ); + _TStringArg signal( "s", "signal", "Output signal", false, "", "file" ); + _TStringArg aSeed( "p", "seed", "Input seed", true, "", "\"x y z\"" ); + _TPixelArg mThr( "t", "min_thr", "Minimum threshold", false, -850, "value (-850)" ); + _TUIntArg sKernel( "k", "kernel", "Kernel size", false, 20, "value (20)" ); + _TRealArg sThr( "r", "signal_thr", "Signal threshold", false, 100, "value (100)" ); + _TRealArg sInfluence( "f", "signal_influence", "Signal influence", false, 0.5, "value (0.5)" ); + _TPixelArg lower( "l", "lower", "Lower threshold", false, -1024, "value (-1024)" ); + _TPixelArg upper( "u", "upper", "Upper threshold", false, 0, "value (0)" ); + _TPixelArg delta( "d", "delta", "Delta threshold", false, 1, "value (1)" ); try { - ts = std::chrono::high_resolution_clock::now( ); - input_image_reader->Update( ); - te = std::chrono::high_resolution_clock::now( ); - tr = te - ts; - std::cout << "Read time: " << tr.count( ) << " s" << std::endl; - - ts = std::chrono::high_resolution_clock::now( ); - filter->Update( ); - te = std::chrono::high_resolution_clock::now( ); - tr = te - ts; - std::cout - << "Mori time: " << tr.count( ) << " s" << std::endl - << "Optimum threshold: " << filter->GetOptimumThreshold( ) << std::endl - << "Number of thresholds: " - << filter->GetNumberOfEvaluatedThresholds( ) << std::endl; + TCLAP::CmdLine cmd( "Mori segmentation", ' ', "1.0.0" ); + cmd.add( delta ); + cmd.add( upper ); + cmd.add( lower ); + cmd.add( sInfluence ); + cmd.add( sThr ); + cmd.add( sKernel ); + cmd.add( mThr ); + cmd.add( aSeed ); + cmd.add( signal ); + cmd.add( out ); + cmd.add( in ); + cmd.parse( argc, argv ); } - catch( std::exception& err ) + catch( TCLAP::ArgException& err ) { - std::cerr << "Error caught: " << err.what( ) << std::endl; + std::cerr + << "===============================" << std::endl + << "Error caught: " << std::endl + << err.error( ) << " " << err.argId( ) + << "===============================" << std::endl + << std::endl; return( 1 ); } // yrt - // Save output image - itk::ImageFileWriter< TLabelImage >::Pointer output_image_writer = - itk::ImageFileWriter< TLabelImage >::New( ); - output_image_writer->SetInput( filter->GetThresholdedOutput( ) ); - output_image_writer->SetFileName( output_image_filename ); + // Get seed + std::istringstream sSeed( aSeed.getValue( ) ); + TImage::PointType seed; + sSeed >> seed[ 0 ] >> seed[ 1 ] >> seed[ 2 ]; + try { - ts = std::chrono::high_resolution_clock::now( ); - output_image_writer->Update( ); - te = std::chrono::high_resolution_clock::now( ); - tr = te - ts; - std::cout << "Write time: " << tr.count( ) << " s" << std::endl; + // Read image + TImage::Pointer input_image; + CTBronchi::ReadImage( input_image, in.getValue( ) ); + + // Mori segmentation + typedef fpa::Filters::Image::Mori< TImage, TLabels > _TMori; + _TMori::Pointer mori = _TMori::New( ); + mori->SetInput( input_image ); + mori->SetSeed( seed ); + mori->SetInsideValue( 1 ); + mori->SetOutsideValue( 0 ); + mori->SetMinimumThreshold( mThr.getValue( ) ); + mori->SetSignalKernelSize( sKernel.getValue( ) ); + mori->SetSignalThreshold( sThr.getValue( ) ); + mori->SetSignalInfluence( sInfluence.getValue( ) ); + mori->SetThresholds( + lower.getValue( ), upper.getValue( ), delta.getValue( ) + ); + double t = CTBronchi::MeasureTime( mori ); + std::cout << "Mori executed in " << t << " s" << std::endl; + + // Write result + CTBronchi::WriteImage( mori->GetOutput( ), out.getValue( ) ); + + // Write signal + if( signal.getValue( ) != "" ) + { + std::stringstream signalStr; + unsigned long nthr = mori->GetNumberOfEvaluatedThresholds( ); + signalStr << "# nThr = " << nthr << std::endl; + signalStr << "# Opt = " << mori->GetOptimumThreshold( ) << std::endl; + for( unsigned long j = 0; j < nthr; ++j ) + { + typename _TMori::TPeak p; + double x, y; + mori->GetSignalValues( j, x, y, p ); + signalStr << x << " " << y << std::endl; + + } // rof + std::ofstream oStrSignals( signal.getValue( ).c_str( ) ); + oStrSignals << signalStr.str( ); + oStrSignals.close( ); + + } // fi } catch( std::exception& err ) { - std::cerr << "Error caught: " << err.what( ) << std::endl; + std::cerr + << "===============================" << std::endl + << "Error caught: " << std::endl + << err.what( ) + << "===============================" << std::endl + << std::endl; return( 1 ); } // yrt - - // Save output signal - std::ofstream osignal( output_signal_filename.c_str( ) ); - for( - unsigned long i = 0; i < filter->GetNumberOfEvaluatedThresholds( ); ++i - ) - { - double x, y; - TFilter::TPeak p; - filter->GetSignalValues( i, x, y, p ); - osignal << x << " " << y << std::endl; - - } // rof - osignal.close( ); - std::cout << "----------------------------------------------" << std::endl; - return( 0 ); }