#include <streambuf>
#include <tclap/CmdLine.h>
+#include <itkAndImageFilter.h>
+#include <itkBinaryThresholdImageFilter.h>
#include <itkMinimumMaximumImageCalculator.h>
#include <itkInvertIntensityImageFilter.h>
#include <itkHessianRecursiveGaussianImageFilter.h>
#include <itkHessian3DToVesselnessMeasureImageFilter.h>
+#include <fpa/Common/SliceBySliceRandomWalker.h>
+#include <fpa/Filters/Image/ExtractSkeleton.h>
#include <fpa/Filters/Image/Mori.h>
+#include <fpa/Filters/Image/RandomWalker.h>
+#include <fpa/Functors/Dijkstra/Image/Gaussian.h>
#include "MoriLabelling.h"
#include "Process.h"
m_InsideLabel( 1 ),
m_OutsideLabel( 2 )
{
- this->m_StrArg[ "input" ] = TArgument( "", true );
- this->m_StrArg[ "seed_file" ] = TArgument( "", false );
- this->m_StrArg[ "seed" ] = TArgument( "", false );
- this->m_StrArg[ "seed_type" ] = TArgument( "point", false );
-
- this->m_DblArg[ "beta" ] = TArgument( "2.5", false );
- this->m_DblArg[ "epsilon" ] = TArgument( "1e-5", false );
- this->m_DblArg[ "vesselness_sigma" ] = TArgument( "0.5", false );
- this->m_DblArg[ "vesselness_alpha1" ] = TArgument( "0.5", false );
- this->m_DblArg[ "vesselness_alpha2" ] = TArgument( "2", false );
- this->m_DblArg[ "mori_min_thr" ] = TArgument( "-850", false );
- this->m_DblArg[ "mori_kernel" ] = TArgument( "20", false );
- this->m_DblArg[ "mori_signal_thr" ] = TArgument( "100", false );
- this->m_DblArg[ "mori_signal_influence" ] = TArgument( "0.5", false );
- this->m_DblArg[ "mori_lower" ] = TArgument( "-1024", false );
- this->m_DblArg[ "mori_upper" ] = TArgument( "0", false );
- this->m_DblArg[ "mori_delta" ] = TArgument( "1", false );
- this->m_DblArg[ "slicerw_thr" ] = TArgument( "5", false );
- this->m_DblArg[ "fastrw_thr" ] = TArgument( "65", false );
- this->m_DblArg[ "labels_upper_thr" ] = TArgument( "-400", false );
+ this->m_StrArgs[ "input" ] = TStrArg( "", "I", true );
+ this->m_StrArgs[ "vesselness" ] = TStrArg( "", "V", false );
+ this->m_StrArgs[ "mori" ] = TStrArg( "", "M", false );
+ this->m_StrArgs[ "labels" ] = TStrArg( "", "L", false );
+ this->m_StrArgs[ "fastrw" ] = TStrArg( "", "F", false );
+ this->m_StrArgs[ "slicerw" ] = TStrArg( "", "R", false );
+ this->m_StrArgs[ "andrw" ] = TStrArg( "", "A", false );
+ this->m_StrArgs[ "fastrw_skeleton" ] = TStrArg( "", "G", false );
+ this->m_StrArgs[ "slicerw_skeleton" ] = TStrArg( "", "U", false );
+ this->m_StrArgs[ "andrw_skeleton" ] = TStrArg( "", "W", false );
+ this->m_StrArgs[ "fastrw_points" ] = TStrArg( "", "B", false );
+ this->m_StrArgs[ "slicerw_points" ] = TStrArg( "", "C", false );
+ this->m_StrArgs[ "andrw_points" ] = TStrArg( "", "D", false );
+ this->m_StrArgs[ "seed_file" ] = TStrArg( "", "P", false );
+ this->m_StrArgs[ "seed" ] = TStrArg( "", "S", false );
+ this->m_StrArgs[ "seed_type" ] = TStrArg( "point", "T", false );
+
+ this->m_DblArgs[ "beta" ] = TDblArg( 2.5, "b", false );
+ this->m_DblArgs[ "epsilon" ] = TDblArg( 1e-5, "e", false );
+ this->m_DblArgs[ "vesselness_sigma" ] = TDblArg( 0.5, "a", false );
+ this->m_DblArgs[ "vesselness_alpha1" ] = TDblArg( 0.5, "c", false );
+ this->m_DblArgs[ "vesselness_alpha2" ] = TDblArg( 2, "d", false );
+ this->m_DblArgs[ "mori_min_thr" ] = TDblArg( -850, "f", false );
+ this->m_DblArgs[ "mori_kernel" ] = TDblArg( 20, "g", false );
+ this->m_DblArgs[ "mori_signal_thr" ] = TDblArg( 100, "i", false );
+ this->m_DblArgs[ "mori_signal_influence" ] = TDblArg( 0.5, "j", false );
+ this->m_DblArgs[ "mori_lower" ] = TDblArg( -1024, "k", false );
+ this->m_DblArgs[ "mori_upper" ] = TDblArg( 0, "l", false );
+ this->m_DblArgs[ "mori_delta" ] = TDblArg( 1, "m", false );
+ this->m_DblArgs[ "slicerw_thr" ] = TDblArg( 5, "n", false );
+ this->m_DblArgs[ "fastrw_thr" ] = TDblArg( 65, "o", false );
+ this->m_DblArgs[ "labels_upper_thr" ] = TDblArg( -400, "p", false );
}
// -------------------------------------------------------------------------
void CTBronchi::Process::
ParseArguments( int argc, char* argv[] )
{
- typedef TCLAP::ValueArg< std::string > _TStrArg;
- typedef TCLAP::ValueArg< TScalar > _TDblArg;
- std::map< std::string, _TStrArg* > strArg;
- std::map< std::string, _TDblArg* > dblArg;
- unsigned char flag = 'a';
-
- TArguments::iterator sIt = this->m_StrArg.begin( );
- for( ; sIt != this->m_StrArg.end( ); ++sIt )
- {
- std::stringstream fStr;
- fStr << flag;
- flag++;
- if( flag == 'h' )
- flag++;
-
- _TStrArg* a = new _TStrArg(
- fStr.str( ), sIt->first, sIt->first,
+ typedef TCLAP::ValueArg< TString > _TStrArg;
+ typedef TCLAP::ValueArg< TScalar > _TDblArg;
+ std::map< TString, _TStrArg* > strArgs;
+ std::map< TString, _TDblArg* > dblArgs;
+
+ // Load arguments
+ TStrArgs::iterator sIt = this->m_StrArgs.begin( );
+ for( ; sIt != this->m_StrArgs.end( ); ++sIt )
+ strArgs[ sIt->first ] = new _TStrArg(
std::get< 1 >( sIt->second ),
+ sIt->first, sIt->first,
+ std::get< 2 >( sIt->second ),
std::get< 0 >( sIt->second ),
- "string"
+ TString( "string (" ) + std::get< 0 >( sIt->second ) + ")"
);
- strArg[ sIt->first ] = a;
-
- } // rof
- TArguments::iterator dIt = this->m_DblArg.begin( );
- for( ; dIt != this->m_DblArg.end( ); ++dIt )
+ TDblArgs::iterator dIt = this->m_DblArgs.begin( );
+ for( ; dIt != this->m_DblArgs.end( ); ++dIt )
{
- std::stringstream fStr;
- fStr << flag;
- flag++;
- if( flag == 'h' )
- flag++;
-
- std::istringstream vStr( std::get< 0 >( dIt->second ) );
- TScalar v;
- vStr >> v;
-
- _TDblArg* a = new _TDblArg(
- fStr.str( ), dIt->first, dIt->first,
+ std::stringstream v;
+ v << "value (" << std::get< 0 >( dIt->second ) << ")";
+ dblArgs[ dIt->first ] = new _TDblArg(
std::get< 1 >( dIt->second ),
- v, "value"
+ dIt->first, dIt->first,
+ std::get< 2 >( dIt->second ),
+ std::get< 0 >( dIt->second ),
+ v.str( )
);
- dblArg[ dIt->first ] = a;
} // rof
- std::map< std::string, _TStrArg* >::iterator saIt;
- std::map< std::string, _TDblArg* >::iterator daIt;
+ // Parse arguments
+ std::map< TString, _TStrArg* >::iterator saIt;
+ std::map< TString, _TDblArg* >::iterator daIt;
try
{
TCLAP::CmdLine cmd( "CTBronchi pipeline", ' ', "1.0.0" );
- for( saIt = strArg.begin( ); saIt != strArg.end( ); ++saIt )
+ for( saIt = strArgs.begin( ); saIt != strArgs.end( ); ++saIt )
cmd.add( *( saIt->second ) );
- for( daIt = dblArg.begin( ); daIt != dblArg.end( ); ++daIt )
+ for( daIt = dblArgs.begin( ); daIt != dblArgs.end( ); ++daIt )
cmd.add( *( daIt->second ) );
cmd.parse( argc, argv );
}
} // yrt
- for( saIt = strArg.begin( ); saIt != strArg.end( ); ++saIt )
+ // Get values
+ for( saIt = strArgs.begin( ); saIt != strArgs.end( ); ++saIt )
{
- std::get< 0 >( this->m_StrArg[ saIt->first ] ) = saIt->second->getValue( );
+ std::get< 0 >( this->m_StrArgs[ saIt->first ] ) =
+ saIt->second->getValue( );
delete saIt->second;
} // rof
- for( daIt = dblArg.begin( ); daIt != dblArg.end( ); ++daIt )
+ for( daIt = dblArgs.begin( ); daIt != dblArgs.end( ); ++daIt )
{
- std::stringstream vStr;
- vStr << daIt->second->getValue( );
- std::get< 0 >( this->m_DblArg[ daIt->first ] ) = vStr.str( );
+ std::get< 0 >( this->m_DblArgs[ daIt->first ] ) =
+ daIt->second->getValue( );
delete daIt->second;
} // rof
- // Compute base filename
- std::string iname = std::get< 0 >( this->m_StrArg[ "input" ] );
- this->m_BaseFileName = iname.substr( 0, iname.find_last_of( "." ) );
+ // Compute base filenames
+ TString bname = std::get< 0 >( this->m_StrArgs[ "input" ] );
+ bname = bname.substr( 0, bname.find_last_of( "." ) );
+ const unsigned int N = 6;
+ const unsigned int M = 6;
+ TString names[ N + M ] =
+ {
+ "vesselness",
+ "mori",
+ "labels",
+ "fastrw",
+ "slicerw",
+ "andrw",
+ "fastrw_skeleton",
+ "slicerw_skeleton",
+ "andrw_skeleton",
+ "fastrw_points",
+ "slicerw_points",
+ "andrw_points"
+ };
+ for( unsigned int i = 0; i < N; ++i )
+ if( std::get< 0 >( this->m_StrArgs[ names[ i ] ] ) == "" )
+ std::get< 0 >( this->m_StrArgs[ names[ i ] ] ) =
+ bname + "_" + names[ i ] + ".mha";
+ for( unsigned int i = 0; i < M; ++i )
+ if( std::get< 0 >( this->m_StrArgs[ names[ i + N ] ] ) == "" )
+ std::get< 0 >( this->m_StrArgs[ names[ i + N ] ] ) =
+ bname + "_" + names[ i + N ] + ".txt";
}
// -------------------------------------------------------------------------
void CTBronchi::Process::
Update( )
{
- this->_Input( );
- this->_Vesselness( );
- this->_Mori( );
- this->_MoriLabelling( );
+ this->_Input( this->m_Input );
+ this->_Vesselness( this->m_Input, this->m_Vesselness );
+ this->_Mori( this->m_Input, this->m_Mori, this->m_Seed );
+ this->_MoriLabelling(
+ this->m_Input, this->m_Mori, this->m_Vesselness, this->m_Labels
+ );
+ this->_FastRW( this->m_Input, this->m_Labels, this->m_FastRW );
+ this->_SliceRW(
+ this->m_Input, this->m_Mori, this->m_Vesselness, this->m_SliceRW
+ );
+ this->_AndImages( this->m_FastRW, this->m_SliceRW, this->m_AndRW );
+ this->_Skeleton(
+ this->m_SliceRW, this->m_SliceRWSkeleton,
+ std::get< 0 >( this->m_StrArgs[ "slicerw_skeleton" ] )
+ );
+ this->_Skeleton(
+ this->m_AndRW, this->m_AndRWSkeleton,
+ std::get< 0 >( this->m_StrArgs[ "andrw_skeleton" ] )
+ );
+ this->_Skeleton(
+ this->m_FastRW, this->m_FastRWSkeleton,
+ std::get< 0 >( this->m_StrArgs[ "fastrw_skeleton" ] )
+ );
}
// -------------------------------------------------------------------------
+template< class _TImage >
void CTBronchi::Process::
-_Input( )
+_Input( _TImage& input )
{
- std::string iname = std::get< 0 >( this->m_StrArg[ "input" ] );
- std::string sname = std::get< 0 >( this->m_StrArg[ "seed_file" ] );
- std::string seed = std::get< 0 >( this->m_StrArg[ "seed" ] );
- std::string stype = std::get< 0 >( this->m_StrArg[ "seed_type" ] );
+ TString sname = std::get< 0 >( this->m_StrArgs[ "seed_file" ] );
+ TString seed = std::get< 0 >( this->m_StrArgs[ "seed" ] );
+ TString stype = std::get< 0 >( this->m_StrArgs[ "seed_type" ] );
if( sname == "" && seed == "" )
throw std::runtime_error( "No seed given." );
if( sname != "" )
sSeed >> idx[ 0 ] >> idx[ 1 ] >> idx[ 2 ];
// Read image
- double t = this->m_Input.Load( iname );
+ TString iname = std::get< 0 >( this->m_StrArgs[ "input" ] );
+ double t = input.Load( iname );
if( t < 0 )
- std::runtime_error( std::string( "Could not read \"" ) + iname + "\"" );
+ std::runtime_error( TString( "Could not read \"" ) + iname + "\"" );
std::cout << "\"" << iname << "\" read in " << t << " s." << std::endl;
// Synch seed
if( stype == "point" )
- this->m_Input.Get( )->TransformPhysicalPointToIndex( pnt, idx );
+ input.Get( )->TransformPhysicalPointToIndex( pnt, idx );
else
- this->m_Input.Get( )->TransformIndexToPhysicalPoint( idx, pnt );
+ input.Get( )->TransformIndexToPhysicalPoint( idx, pnt );
this->m_Seed = TSeed( idx, pnt );
}
// -------------------------------------------------------------------------
+template< class _TInput, class _TVesselness >
void CTBronchi::Process::
-_Vesselness( )
+_Vesselness( _TInput& input, _TVesselness& vesselness )
{
- std::string vname = this->m_BaseFileName + "_vesselness.mha";
- double t = this->m_Vesselness.Load( vname );
+ TString iname = std::get< 0 >( this->m_StrArgs[ "vesselness" ] );
+ double t = vesselness.Load( iname );
if( t < 0 )
{
t = 0;
- // Arguments
- std::stringstream v;
- v << std::get< 0 >( this->m_DblArg[ "vesselness_sigma" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "vesselness_alpha1" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "vesselness_alpha2" ] );
- std::istringstream w( v.str( ) );
- TScalar s, a1, a2;
- w >> s >> a1 >> a2;
-
// Min-max
- typedef itk::MinimumMaximumImageCalculator< TPixelImage::TImage > _TMinMax;
- _TMinMax::Pointer minMax = _TMinMax::New( );
- minMax->SetImage( this->m_Input.Get( ) );
+ typedef itk::MinimumMaximumImageCalculator< typename _TInput::TImage > _TMinMax;
+ typename _TMinMax::Pointer minMax = _TMinMax::New( );
+ minMax->SetImage( input.Get( ) );
minMax->Compute( );
// Invert intensities
- typedef CTBronchi::Filter< itk::InvertIntensityImageFilter< TPixelImage::TImage > > _TInverter;
+ typedef CTBronchi::Filter< itk::InvertIntensityImageFilter< typename _TInput::TImage > > _TInverter;
_TInverter inverter;
- inverter.Get( )->SetInput( this->m_Input.Get( ) );
+ inverter.Get( )->SetInput( input.Get( ) );
inverter.Get( )->SetMaximum( minMax->GetMaximum( ) );
double t1 = inverter.Update( );
t += t1;
std::cout << "Inversion executed in " << t1 << " s" << std::endl;
// Compute hessian image
- typedef CTBronchi::Filter< itk::HessianRecursiveGaussianImageFilter< TPixelImage::TImage > > _THessian;
+ typedef CTBronchi::Filter< itk::HessianRecursiveGaussianImageFilter< typename _TInput::TImage > > _THessian;
_THessian hessian;
hessian.Get( )->SetInput( inverter.Get( )->GetOutput( ) );
- hessian.Get( )->SetSigma( s );
+ hessian.Get( )->SetSigma(
+ std::get< 0 >( this->m_DblArgs[ "vesselness_sigma" ] )
+ );
t1 = hessian.Update( );
t += t1;
std::cout << "Hessian executed in " << t1 << " s" << std::endl;
// Vesselness
- typedef CTBronchi::Filter< itk::Hessian3DToVesselnessMeasureImageFilter< TScalar > > _TVesselness;
- _TVesselness vesselness;
- vesselness.Get( )->SetInput( hessian.Get( )->GetOutput( ) );
- vesselness.Get( )->SetAlpha1( a1 );
- vesselness.Get( )->SetAlpha2( a2 );
- t1 = vesselness.Update( );
+ typedef CTBronchi::Filter< itk::Hessian3DToVesselnessMeasureImageFilter< typename _TVesselness::TImage::PixelType > > _TVesselnessFilter;
+ _TVesselnessFilter vFilter;
+ vFilter.Get( )->SetInput( hessian.Get( )->GetOutput( ) );
+ vFilter.Get( )->SetAlpha1(
+ std::get< 0 >( this->m_DblArgs[ "vesselness_alpha1" ] )
+ );
+ vFilter.Get( )->SetAlpha2(
+ std::get< 0 >( this->m_DblArgs[ "vesselness_alpha2" ] )
+ );
+ t1 = vFilter.Update( );
t += t1;
std::cout << "Hessian measure computed in " << t1 << " s." << std::endl;
- this->m_Vesselness.Set( vesselness.Get( )->GetOutput( ) );
- t1 = this->m_Vesselness.Save( vname );
+ vesselness.Set( vFilter.Get( )->GetOutput( ) );
+ t1 = vesselness.Save( iname );
t += t1;
- std::cout << "\"" << vname << "\" saved in " << t1 << " s." << std::endl;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
} // fi
std::cout << "Vesselness computed in " << t << " s." << std::endl;
}
// -------------------------------------------------------------------------
+template< class _TInput, class _TMori, class _TSeed >
void CTBronchi::Process::
-_Mori( )
+_Mori( _TInput& input, _TMori& mori, _TSeed& seed )
{
- std::string mname = this->m_BaseFileName + "_mori.mha";
- double t = this->m_Mori.Load( mname );
+ TString iname = std::get< 0 >( this->m_StrArgs[ "mori" ] );
+ double t = mori.Load( iname );
if( t < 0 )
{
t = 0;
- // Arguments
- std::stringstream v;
- v << std::get< 0 >( this->m_DblArg[ "mori_min_thr" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_kernel" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_signal_thr" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_signal_influence" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_lower" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_upper" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "mori_delta" ] );
- std::istringstream w( v.str( ) );
- TScalar mThr, sKernel, sThr, sInfluence, lower, upper, delta;
- w >> mThr >> sKernel >> sThr >> sInfluence >> lower >> upper >> delta;
-
// Mori segmentation
- typedef CTBronchi::Filter< fpa::Filters::Image::Mori< TPixelImage::TImage, TLabelImage::TImage > > _TMori;
- _TMori mori;
- mori.Get( )->SetInput( this->m_Input.Get( ) );
- mori.Get( )->SetSeed( this->m_Seed.second );
- mori.Get( )->SetInsideValue( this->m_InsideLabel );
- mori.Get( )->SetOutsideValue( this->m_UndefinedLabel );
- mori.Get( )->SetMinimumThreshold( mThr );
- mori.Get( )->SetSignalKernelSize( sKernel );
- mori.Get( )->SetSignalThreshold( sThr );
- mori.Get( )->SetSignalInfluence( sInfluence );
- mori.Get( )->SetThresholds( lower, upper, delta );
- double t1 = mori.Update( );
+ typedef CTBronchi::Filter< fpa::Filters::Image::Mori< TPixelImage::TImage, TLabelImage::TImage > > _TMoriFilter;
+ _TMoriFilter mFilter;
+ mFilter.Get( )->SetInput( input.Get( ) );
+ mFilter.Get( )->SetSeed( seed.second );
+ mFilter.Get( )->SetInsideValue( this->m_InsideLabel );
+ mFilter.Get( )->SetOutsideValue( this->m_UndefinedLabel );
+ mFilter.Get( )->SetMinimumThreshold(
+ std::get< 0 >( this->m_DblArgs[ "mori_min_thr" ] )
+ );
+ mFilter.Get( )->SetSignalKernelSize(
+ std::get< 0 >( this->m_DblArgs[ "mori_kernel" ] )
+ );
+ mFilter.Get( )->SetSignalThreshold(
+ std::get< 0 >( this->m_DblArgs[ "mori_signal_thr" ] )
+ );
+ mFilter.Get( )->SetSignalInfluence(
+ std::get< 0 >( this->m_DblArgs[ "mori_signal_influence" ] )
+ );
+ mFilter.Get( )->SetThresholds(
+ std::get< 0 >( this->m_DblArgs[ "mori_lower" ] ),
+ std::get< 0 >( this->m_DblArgs[ "mori_upper" ] ),
+ std::get< 0 >( this->m_DblArgs[ "mori_delta" ] )
+ );
+ double t1 = mFilter.Update( );
t += t1;
std::cout << "Segmentation computed in " << t1 << " s." << std::endl;
- this->m_Mori.Set( mori.Get( )->GetOutput( ) );
- t1 = this->m_Mori.Save( mname );
+ mori.Set( mFilter.Get( )->GetOutput( ) );
+ t1 = mori.Save( iname );
t += t1;
- std::cout << "\"" << mname << "\" saved in " << t1 << " s." << std::endl;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
} // fi
std::cout << "Mori segmentation computed in " << t << " s." << std::endl;
}
// -------------------------------------------------------------------------
+template< class _TInput, class _TMori, class _TVesselness, class _TLabels >
void CTBronchi::Process::
-_MoriLabelling( )
+_MoriLabelling(
+ _TInput& input, _TMori& mori, _TVesselness& vesselness, _TLabels& labels
+ )
{
- std::string mname = this->m_BaseFileName + "_labels.mha";
- double t = this->m_Labels.Load( mname );
+ TString iname = std::get< 0 >( this->m_StrArgs[ "labels" ] );
+ double t = labels.Load( iname );
if( t < 0 )
{
t = 0;
- // Arguments
- std::stringstream v;
- v << std::get< 0 >( this->m_DblArg[ "fastrw_thr" ] ) << " "
- << std::get< 0 >( this->m_DblArg[ "labels_upper_thr" ] );
- std::istringstream w( v.str( ) );
- TScalar vThr, uThr;
- w >> vThr >> uThr;
-
// Mori labelling
- typedef CTBronchi::Filter< CTBronchi::MoriLabelling< TPixelImage::TImage, TLabelImage::TImage, TScalarImage::TImage > > _TLabelling;
+ typedef CTBronchi::Filter< CTBronchi::MoriLabelling< typename _TInput::TImage, typename _TLabels::TImage, typename _TVesselness::TImage > > _TLabelling;
_TLabelling labelling;
- labelling.Get( )->SetInput( this->m_Input.Get( ) );
- labelling.Get( )->SetInputLabels( this->m_Mori.Get( ) );
- labelling.Get( )->SetInputVesselness( this->m_Vesselness.Get( ) );
- labelling.Get( )->SetVesselnessThreshold( vThr );
- labelling.Get( )->SetUpperThreshold( uThr );
+ labelling.Get( )->SetInput( input.Get( ) );
+ labelling.Get( )->SetInputLabels( mori.Get( ) );
+ labelling.Get( )->SetInputVesselness( vesselness.Get( ) );
+ labelling.Get( )->SetVesselnessThreshold( std::get< 0 >( this->m_DblArgs[ "fastrw_thr" ] ) );
+ labelling.Get( )->SetUpperThreshold( std::get< 0 >( this->m_DblArgs[ "labels_upper_thr" ] ) );
labelling.Get( )->SetInsideValue( this->m_InsideLabel );
labelling.Get( )->SetOutsideValue( this->m_OutsideLabel );
labelling.Get( )->SetFillValue( this->m_OutsideLabel );
t += t1;
std::cout << "Labelling computed in " << t1 << " s." << std::endl;
- this->m_Labels.Set( labelling.Get( )->GetOutput( ) );
- t1 = this->m_Labels.Save( mname );
+ labels.Set( labelling.Get( )->GetOutput( ) );
+ t1 = labels.Save( iname );
t += t1;
- std::cout << "\"" << mname << "\" saved in " << t1 << " s." << std::endl;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
} // fi
std::cout << "Mori labelling computed in " << t << " s." << std::endl;
}
+// -------------------------------------------------------------------------
+template< class _TInput, class _TLabels, class _TFastRW >
+void CTBronchi::Process::
+_FastRW( _TInput& input, _TLabels& labels, _TFastRW& fastrw )
+{
+ TString iname = std::get< 0 >( this->m_StrArgs[ "fastrw" ] );
+ double t = fastrw.Load( iname );
+ if( t < 0 )
+ {
+ t = 0;
+
+ // Prepare weight functor
+ typedef fpa::Functors::Dijkstra::Image::Gaussian< typename _TInput::TImage, typename _TFastRW::TImage::PixelType > TWeight;
+ typename TWeight::Pointer weight = TWeight::New( );
+ weight->SetBeta( std::get< 0 >( this->m_DblArgs[ "beta" ] ) );
+ weight->SetEpsilon( std::get< 0 >( this->m_DblArgs[ "epsilon" ] ) );
+
+ // Random walk
+ typedef CTBronchi::Filter< fpa::Filters::Image::RandomWalker< typename _TInput::TImage, typename _TLabels::TImage, typename _TFastRW::TImage::PixelType > > TFilter;
+ TFilter filter;
+ filter.Get( )->SetInputImage( input.Get( ) );
+ filter.Get( )->SetInputLabels( labels.Get( ) );
+ filter.Get( )->SetWeightFunction( weight );
+ double t1 = filter.Update( );
+ t += t1;
+ std::cout << "Fast random walker executed in " << t1 << " s" << std::endl;
+
+ // Extract label
+ typedef CTBronchi::Filter< itk::BinaryThresholdImageFilter< typename _TLabels::TImage, typename _TLabels::TImage > > _TExtract;
+ _TExtract extract;
+ extract.Get( )->SetInput( filter.Get( )->GetOutputLabels( ) );
+ extract.Get( )->SetInsideValue( this->m_InsideLabel );
+ extract.Get( )->SetOutsideValue( this->m_UndefinedLabel );
+ extract.Get( )->SetLowerThreshold( this->m_InsideLabel );
+ extract.Get( )->SetUpperThreshold( this->m_InsideLabel );
+ t1 = extract.Update( );
+ t += t1;
+ std::cout << "Extract labels executed in " << t1 << " s" << std::endl;
+
+ fastrw.Set( extract.Get( )->GetOutput( ) );
+ t1 = fastrw.Save( iname );
+ t += t1;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
+
+ } // fi
+ std::cout << "Fast random walker computed in " << t << " s." << std::endl;
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TLabels, class _TVesselness, class _TSliceRW >
+void CTBronchi::Process::
+_SliceRW(
+ _TInput& input, _TLabels& labels, _TVesselness& vesselness, _TSliceRW& slicerw
+ )
+{
+ TString iname = std::get< 0 >( this->m_StrArgs[ "slicerw" ] );
+ double t = slicerw.Load( iname );
+ if( t < 0 )
+ {
+ t = 0;
+
+ // Random walk
+ typedef CTBronchi::Filter< fpa::Common::SliceBySliceRandomWalker< typename _TInput::TImage, typename _TLabels::TImage, typename _TVesselness::TImage > > TFilter;
+ TFilter filter;
+ filter.Get( )->SetInput( input.Get( ) );
+ filter.Get( )->SetInputLabels( labels.Get( ) );
+ filter.Get( )->SetInputVesselness( vesselness.Get( ) );
+ filter.Get( )->SetBeta( std::get< 0 >( this->m_DblArgs[ "beta" ] ) );
+ filter.Get( )->SetVesselnessThreshold( std::get< 0 >( this->m_DblArgs[ "slicerw_thr" ] ) );
+ filter.Get( )->SetEpsilon( std::get< 0 >( this->m_DblArgs[ "epsilon" ] ) );
+ double t1 = filter.Update( );
+ t += t1;
+ std::cout << "Extract labels executed in " << t1 << " s" << std::endl;
+
+ slicerw.Set( filter.Get( )->GetOutput( ) );
+ t1 = slicerw.Save( iname );
+ t += t1;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
+
+ } // fi
+ std::cout << "Slice by slice random walker computed in " << t << " s." << std::endl;
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput >
+void CTBronchi::Process::
+_AndImages( _TInput& a, _TInput& b, _TInput& c )
+{
+ TString iname = std::get< 0 >( this->m_StrArgs[ "andrw" ] );
+ double t = c.Load( iname );
+ if( t < 0 )
+ {
+ t = 0;
+
+ // Random walk
+ typedef CTBronchi::Filter< itk::AndImageFilter< typename _TInput::TImage > > TFilter;
+ TFilter filter;
+ filter.Get( )->SetInput( 0, a.Get( ) );
+ filter.Get( )->SetInput( 1, b.Get( ) );
+ double t1 = filter.Update( );
+ t += t1;
+ std::cout << "And segmentations executed in " << t1 << " s" << std::endl;
+
+ c.Set( filter.Get( )->GetOutput( ) );
+ t1 = c.Save( iname );
+ t += t1;
+ std::cout << "\"" << iname << "\" saved in " << t1 << " s." << std::endl;
+
+ } // fi
+ std::cout << "And segmentations computed in " << t << " s." << std::endl;
+}
+
+// -------------------------------------------------------------------------
+template< class _TInput, class _TSkeleton >
+void CTBronchi::Process::
+_Skeleton( _TInput& a, _TSkeleton& s, const TString& fname )
+{
+ double t = s.Load( fname );
+ if( t < 0 )
+ {
+ t = 0;
+
+ typedef CTBronchi::Filter< fpa::Filters::Image::ExtractSkeleton< typename _TInput::TImage > > TFilter;
+ TFilter filter;
+ filter.Get( )->SetInput( a.Get( ) );
+ filter.Get( )->SeedFromMaximumDistanceOff( );
+ filter.Get( )->SetSeed( this->m_Seed.first );
+ filter.Get( )->GetDistanceMap( )->InsideIsPositiveOn( );
+ filter.Get( )->GetDistanceMap( )->SquaredDistanceOff( );
+ filter.Get( )->GetDistanceMap( )->UseImageSpacingOn( );
+ double t1 = filter.Update( );
+ t += t1;
+ std::cout << "Skeleton executed in " << t1 << " s" << std::endl;
+
+ s.Set( filter.Get( )->GetOutput( ) );
+ t1 = s.Save( fname );
+ t += t1;
+ std::cout << "\"" << fname << "\" saved in " << t1 << " s." << std::endl;
+
+ } // fi
+ std::cout << "Skeleton computed in " << t << " s." << std::endl;
+}
+
// eof - $RCSfile$
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / blobdiff