+++ /dev/null
-#!/bin/bash
-
-input=$1
-dir=`dirname $input`
-base=`basename $1 .mhd`_florez-l
-output="$dir/$base"_segmentation.mhd
-mori="$dir/$base"_mori.mhd
-signal="$dir/$base"_mori_signal.txt
-labels="$dir/$base"_labels.mhd
-diff="$dir/$base"_diff.mhd
-seed="$dir"/seed_florez-l.txt
-
-(>&2 echo -n "Processing $input... ")
-echo "************************************************"
-./fpa_CTBronchi_Process \
- -in $input \
- -out $output \
- -seed `cat $seed` \
- -out_mori $mori \
- -out_signal $signal \
- -out_labels $labels \
- -out_diff $diff
-echo "************************************************"
-(>&2 echo "done.")
-
-## eof - $RCSfile$
_TStringArg labels( "l", "labels", "Input labels", true, "", "file" );
_TStringArg vesselness( "v", "vesselness", "Input vesselness", true, "", "file" );
_TStringArg out( "o", "output", "Output image", true, "", "file" );
- _TRealArg vThr( "a", "vesselness_thr", "Vesselness threshold", false, 0.05, "value (0.05)" );
+ _TRealArg vThr( "a", "vesselness_thr", "Vesselness threshold (%)", false, 5, "value (5)" );
_TPixelArg uThr( "u", "upper_thr", "Upper threshold", false, -400, "value (-400)" );
try
{
protected:
MoriLabelling( )
: Superclass( ),
- m_VesselnessThreshold( 0.05 ),
+ m_VesselnessThreshold( 5 ),
m_UpperThreshold( -650 )
{
ivqITKInputConfigureMacro( InputLabels, _TLabels );
minMax->SetImage( this->GetInputVesselness( ) );
minMax->Compute( );
this->m_MinVesselness =
- ( 1.0 - this->m_VesselnessThreshold ) *
+ ( 1.0 - ( this->m_VesselnessThreshold / double( 100 ) ) ) *
double( minMax->GetMaximum( ) );
}
+++ /dev/null
-// =========================================================================
-// @author Leonardo Florez Valencia
-// @email florez-l@javeriana.edu.co
-// =========================================================================
-
-#include <itkImage.h>
-#include <itkAbsoluteValueDifferenceImageFilter.h>
-#include <itkBinaryThresholdImageFilter.h>
-#include <CTBronchi/Functions.h>
-#include <CTBronchi/Mori.h>
-#include <CTBronchi/MoriLabelling.h>
-#include <CTBronchi/RandomWalker.h>
-
-// -------------------------------------------------------------------------
-const unsigned int Dim = 3;
-typedef short TInputPixel;
-typedef unsigned char TLabelPixel;
-typedef itk::Image< TInputPixel, Dim > TImage;
-typedef itk::Image< TLabelPixel, Dim > TLabelImage;
-
-// -------------------------------------------------------------------------
-int main( int argc, char* argv[] )
-{
- std::map< std::string, std::string > args;
- try
- {
- if( CTBronchi::ParseArgs( args, argc, argv ) )
- {
- // Parse seed
- TImage::PointType seed;
- char* str = new char[ args[ "seed" ].size( ) + 1 ];
- strcpy( str, args[ "seed" ].c_str( ) );
- seed[ 0 ] = std::atof( strtok( str, ";" ) );
- seed[ 1 ] = std::atof( strtok( NULL, ";" ) );
- seed[ 2 ] = std::atof( strtok( NULL, ";" ) );
-
- // Read input image
- TImage::Pointer input_image;
- CTBronchi::ReadImage( input_image, args[ "in" ] );
-
- // Mori segmentation
- TLabelImage::Pointer mori;
- TInputPixel opt_thr = CTBronchi::Mori( mori, input_image, seed, args );
-
- // Label image
- TLabelImage::Pointer labels;
- CTBronchi::Label( labels, input_image, mori, opt_thr, args );
-
- // Final labels
- TLabelImage::Pointer final_labels;
- CTBronchi::RandomWalker( final_labels, input_image, labels, args );
-
- // Extract label 1
- typedef itk::BinaryThresholdImageFilter< TLabelImage, TLabelImage > TBinThr;
- TBinThr::Pointer bin_thr = TBinThr::New( );
- bin_thr->SetInput( final_labels );
- bin_thr->SetLowerThreshold( 1 );
- bin_thr->SetUpperThreshold( 1 );
- bin_thr->SetInsideValue( 1 );
- bin_thr->SetOutsideValue( 0 );
- double t = CTBronchi::MeasureTime( bin_thr );
- std::cout << "Label extracted in " << t << " s" << std::endl;
- TLabelImage::Pointer segmentation = bin_thr->GetOutput( );
-
- // Save
- CTBronchi::WriteImage( segmentation, args[ "out" ] );
-
- // Save difference
- std::map< std::string, std::string >::const_iterator i =
- args.find( "out_diff" );
- if( i != args.end( ) )
- {
- typedef itk::AbsoluteValueDifferenceImageFilter< TLabelImage, TLabelImage, TLabelImage > TDiff;
- TDiff::Pointer diff = TDiff::New( );
- diff->SetInput1( mori );
- diff->SetInput2( segmentation );
- t = CTBronchi::MeasureTime( diff );
- std::cout << "Difference computed in " << t << " s" << std::endl;
- TLabelImage::Pointer difference = diff->GetOutput( );
- CTBronchi::WriteImage( difference, i->second );
-
- } // fi
- }
- else
- return( 1 );
- }
- catch( std::exception& err )
- {
- std::cerr
- << "===============================" << std::endl
- << "Error caught: " << std::endl
- << err.what( )
- << "===============================" << std::endl
- << std::endl;
- return( 1 );
-
- } // yrt
- return( 0 );
-}
-
-// eof - $RCSfile$
--- /dev/null
+#!/bin/bash
+
+## -- Command line options
+vesselness_sigma="0.5"
+vesselness_alpha1="0.5"
+vesselness_alpha2="2"
+mori_min_thr="-850"
+mori_kernel="20"
+mori_signal_thr="100"
+mori_signal_influence="0.5"
+mori_lower="-1024"
+mori_upper="0"
+mori_delta="1"
+labels_vesselness_thr="5"
+labels_upper_thr="-400"
+beta="2.5"
+epsilon="1e-5"
+while [[ "$#" -gt 0 ]]; do
+ key="$1"
+ case $key in
+ -input)
+ input=`realpath $2`
+ shift
+ ;;
+ -force)
+ force="1"
+ ;;
+ -vesselness_sigma)
+ vesselness_sigma="$2"
+ shift
+ ;;
+ -vesselness_alpha1)
+ vesselness_alpha1="$2"
+ shift
+ ;;
+ -vesselness_alpha2)
+ vesselness_alpha2="$2"
+ shift
+ ;;
+ -mori_min_thr)
+ mori_min_thr="$2"
+ shift
+ ;;
+ -mori_kernel)
+ mori_kernel="$2"
+ shift
+ ;;
+ -mori_signal_thr)
+ mori_signal_thr="$2"
+ shift
+ ;;
+ -mori_signal_influence)
+ mori_signal_influence="$2"
+ shift
+ ;;
+ -mori_lower)
+ mori_lower="$2"
+ shift
+ ;;
+ -mori_upper)
+ mori_upper="$2"
+ shift
+ ;;
+ -mori_delta)
+ mori_delta="$2"
+ shift
+ ;;
+ -labels_vesselness_thr)
+ labels_vesselness_thr="$2"
+ shift
+ ;;
+ -labels_upper_thr)
+ labels_upper_thr="$2"
+ shift
+ ;;
+ -beta)
+ beta="$2"
+ shift
+ ;;
+ -epsilon)
+ epsilon="$2"
+ shift
+ ;;
+ -seed)
+ seed="$2 $3 $4"
+ shift
+ shift
+ shift
+ ;;
+ *)
+ # Do nothing
+ ;;
+ esac
+ shift
+done
+
+## -- Check command line options
+if [ -z "$input" ] || [ -z "$seed" ] ; then
+ (>&2 echo "Usage: $0 -input <file> -seed <x y z> [-force]")
+ exit 1
+fi
+
+base_dir=`dirname $input | xargs realpath`
+base_name=$base_dir/`basename $input .mhd`
+
+vesselness=$base_name"_vesselness.mhd"
+mori=$base_name"_mori.mhd"
+mori_signal=$base_name"_mori_signal.txt"
+labels=$base_name"_labels.mhd"
+fastrw=$base_name"_fastrw.mhd"
+slicerw=$base_name"_slicerw.mhd"
+
+(>&2 echo "Processing $input... ")
+echo "************************************************"
+if [ ! -f $vesselness ] || [ -n "$force" ] ; then
+ ./fpa_CTBronchi_Vesselness \
+ -i $input -o $vesselness \
+ -s $vesselness_sigma -a $vesselness_alpha1 -b $vesselness_alpha2
+fi
+
+if [ ! -f $mori ] || [ -n "$force" ] ; then
+ ./fpa_CTBronchi_MoriSegmentation \
+ -i $input \
+ -o $mori \
+ -s $mori_signal \
+ -p "$seed" \
+ -t $mori_min_thr \
+ -k $mori_kernel \
+ -r $mori_signal_thr \
+ -f $mori_signal_influence \
+ -l $mori_lower \
+ -u $mori_upper \
+ -d $mori_delta
+fi
+
+if [ ! -f $labels ] || [ -n "$force" ] ; then
+ ./fpa_CTBronchi_MoriLabelling \
+ -i $input \
+ -l $mori \
+ -v $vesselness \
+ -o $labels \
+ -a $labels_vesselness_thr \
+ -u $labels_upper_thr
+fi
+
+if [ ! -f $fastrw ] || [ -n "$force" ] ; then
+ ./fpa_CTBronchi_FastRandomWalker \
+ -i $input \
+ -l $labels \
+ -o $fastrw \
+ -b $beta \
+ -e $epsilon
+fi
+
+if [ ! -f $slicerw ] || [ -n "$force" ] ; then
+ ./fpa_CTBronchi_SliceBySliceRandomWalker \
+ -i $input \
+ -l $mori \
+ -v $vesselness \
+ -o $slicerw \
+ -t $labels_vesselness_thr \
+ -b $beta \
+ -e $epsilon
+fi
+echo "************************************************"
+(>&2 echo "done.")
+
+## eof - $RCSfile$
_TStringArg out( "o", "output", "Output image", true, "", "file" );
_TRealArg beta( "b", "beta", "Beta", false, 2.5, "value (2.5)" );
_TRealArg eps( "e", "epsilon", "Epsilon", false, 1e-5, "value (1e-5)" );
+ _TRealArg vThr( "t", "vthreshold", "Vesselness thresnold (%)", false, 5, "value (5)" );
try
{
TCLAP::CmdLine cmd( "FastRandomWalker", ' ', "1.0.0" );
cmd.add( eps );
+ cmd.add( vThr );
cmd.add( beta );
cmd.add( out );
cmd.add( vesselness );
filter->SetInput( input_image );
filter->SetInputLabels( input_labels );
filter->SetInputVesselness( input_vesselness );
+ filter->SetBeta( beta.getValue( ) );
+ filter->SetVesselnessThreshold( vThr.getValue( ) );
+ filter->SetEpsilon( eps.getValue( ) );
+ double t = CTBronchi::MeasureTime( filter );
+ std::cout
+ << "SliceBySliceRandomWalker executed in "
+ << t << " s" << std::endl;
- // Prepare weight functor
- /* TODO
- typedef fpa::Functors::Dijkstra::Image::Gaussian< TImage, TScalar > TWeight;
- TWeight::Pointer weight = TWeight::New( );
- weight->SetBeta( beta.getValue( ) );
- weight->SetEpsilon( eps.getValue( ) );
-
- // Random walk
- typedef fpa::Filters::Image::RandomWalker< TImage, TLabels, TScalar > TFilter;
- TFilter::Pointer filter = TFilter::New( );
- filter->SetInputImage( input_image );
- filter->SetInputLabels( input_labels );
- filter->SetWeightFunction( weight );
- double t = CTBronchi::MeasureTime( filter );
- std::cout << "FastRandomWalk executed in " << t << " s" << std::endl;
-
- // Write result
- CTBronchi::WriteImage( filter->GetOutputLabels( ), out.getValue( ) );
- */
+ // Write result
+ CTBronchi::WriteImage( filter->GetOutput( ), out.getValue( ) );
}
catch( std::exception& err )
{
#include <itkImageFileReader.h>
#include <itkImageFileWriter.h>
-#include <fpa/Common/OriginalRandomWalker.h>
+#include <fpa/Common/RandomWalker.h>
#include <fpa/Filters/Image/RandomWalker.h>
#include <fpa/Functors/Dijkstra/Image/Gaussian.h>
fpa::Functors::Dijkstra::Image::Gaussian< _TImage, _TScalar >
_TFunction;
typedef
- fpa::Common::OriginalRandomWalker< _TImage, _TLabels, _TScalar >
+ fpa::Common::RandomWalker< _TImage, _TLabels, _TScalar >
_TFilter;
// Random walker function
#include <itkImage.h>
#include <itkImageFileReader.h>
#include <itkImageFileWriter.h>
-#include <fpa/Common/OriginalRandomWalker.h>
+#include <fpa/Common/RandomWalker.h>
#include <fpa/Functors/Dijkstra/Image/Gaussian.h>
// -------------------------------------------------------------------------
edge->TreatAsWeightOff( );
// Random walker
- typedef fpa::Common::OriginalRandomWalker< TImage, TLabels, TScalar > TFilter;
+ typedef fpa::Common::RandomWalker< TImage, TLabels, TScalar > TFilter;
TFilter::Pointer rw = TFilter::New( );
rw->SetInput( input_image_reader->GetOutput( ) );
rw->SetInputLabels( input_labels_reader->GetOutput( ) );
fpa::Common::SliceBySliceRandomWalker, itk::ImageToImageFilter
);
+ itkGetConstMacro( Epsilon, TScalar );
+ itkSetMacro( Epsilon, TScalar );
+
+ itkGetConstMacro( Beta, TScalar );
+ itkSetMacro( Beta, TScalar );
+
+ itkGetConstMacro( VesselnessThreshold, TScalar );
+ itkSetMacro( VesselnessThreshold, TScalar );
+
ivqITKInputMacro( InputLabels, TLabels );
ivqITKInputMacro( InputVesselness, TScalarImage );
void _Composite(
typename TScalarImage::Pointer& composite,
const TLabels* labels,
- const TScalarImage* vesselness
+ const TScalarImage* vesselness,
+ const TScalar& maxVess
);
template< class _TSlicePtr, class _TInput >
// Purposely not implemented
SliceBySliceRandomWalker( const Self& other );
Self& operator=( const Self& other );
+
+ protected:
+ TScalar m_Epsilon;
+ TScalar m_Beta;
+ TScalar m_VesselnessThreshold;
+
+ TScalar m_VesselnessValue;
};
} // ecapseman
#include <itkImageRegionConstIterator.h>
#include <itkImageRegionIterator.h>
#include <itkJoinSeriesImageFilter.h>
+#include <itkMaximumImageFilter.h>
#include <itkMinimumMaximumImageCalculator.h>
#include <fpa/Common/RandomWalker.h>
#include <fpa/Functors/Dijkstra/Image/Gaussian.h>
template< class _TImage, class _TLabels, class _TScalarImage >
fpa::Common::SliceBySliceRandomWalker< _TImage, _TLabels, _TScalarImage >::
SliceBySliceRandomWalker( )
- : Superclass( )
+ : Superclass( ),
+ m_Epsilon( TScalar( 1e-5 ) ),
+ m_Beta( TScalar( 1 ) ),
+ m_VesselnessThreshold( TScalar( 5 ) )
{
ivqITKInputConfigureMacro( InputLabels, TLabels );
ivqITKInputConfigureMacro( InputVesselness, TScalarImage );
const TScalarImage* vesselness = this->GetInputVesselness( );
this->AllocateOutputs( );
+ // Get vesselness range
+ typedef itk::MinimumMaximumImageCalculator< TScalarImage > _TMinMax;
+ typename _TMinMax::Pointer vMinMax = _TMinMax::New( );
+ vMinMax->SetImage( vesselness );
+ vMinMax->Compute( );
+ TScalar vMax = vMinMax->GetMaximum( );
+
// Some values
typename TImage::RegionType region = input->GetRequestedRegion( );
- typename TImage::SizeType regionIndex = region.GetIndex( );
+ typename TImage::IndexType regionIndex = region.GetIndex( );
typename TImage::SizeType sliceSize = region.GetSize( );
int numSlices = sliceSize[ 2 ];
sliceSize[ 2 ] = 0;
- this->m_VesselnessValue =
- this->m_VesselnessThreshold * this->m_VesselnessMax / TScalar( 100 );
+ this->m_VesselnessValue = TScalar( 0.01 );
+ this->m_VesselnessValue *= this->m_VesselnessThreshold * vMax;
// Composite image
typename TScalarImage::Pointer composite;
- this->_Composite( composite, labels, vesselness );
+ this->_Composite( composite, labels, vesselness, vMax );
// Extract slices
std::vector< typename _TSliceImage::Pointer > data3D( numSlices );
typename _TSliceScalarImage::Pointer composite_slice;
this->_Slice( input_slice, input, sliceRegion );
this->_Slice( labels_slice, labels, sliceRegion );
- this->_Slice( composite_slice, composite, sliceRegion );
+ this->_Slice( composite_slice, composite.GetPointer( ), sliceRegion );
// Update vectors with each image
data3D[ i ] = input_slice;
this->_RandomWalker( binaryTree, data3D, fusion, false );
// Join results
- typedef itk::JoinSeriesImageFilter< _TSliceImage, TImage > _TJoin;
+ typedef itk::JoinSeriesImageFilter< _TSliceLabels, TLabels > _TJoin;
typename _TJoin::Pointer join = _TJoin::New( );
for( int i = 0; i < numSlices; ++i )
join->SetInput( i, binaryTree[ i ] );
_Composite(
typename TScalarImage::Pointer& composite,
const TLabels* labels,
- const TScalarImage* vesselness
+ const TScalarImage* vesselness,
+ const TScalar& maxVess
)
{
- // Min-Max vesselness values
- typedef itk::MinimumMaximumImageCalculator< TScalarImage > _TMinMax;
- typename _TMinMax::Pointer minMax = _TMinMax::New( );
- minMax->SetImage( vesselness );
- minMax->Compute( );
- TScalar maxVess = minMax->GetMaximum( );
+ // Some values
typename TScalarImage::RegionType region = labels->GetRequestedRegion( );
// Composite image
typename _TInput::RegionType region
)
{
- typedef typename _TSlicePTr::ObjectType _TSlice;
+ typedef typename _TSlicePtr::ObjectType _TSlice;
typedef itk::ExtractImageFilter< _TInput, _TSlice > _TExtract;
typename _TExtract::Pointer extract = _TExtract::New( );
// Random walker function
typedef fpa::Functors::Dijkstra::Image::Gaussian< _TSliceData3D, _TSliceScalar > _TFunction;
- typename TFunction::Pointer edge = TFunction::New( );
+ typename _TFunction::Pointer edge = _TFunction::New( );
edge->SetBeta( this->m_Beta );
- edge->SetEpsilon( this->m_Eps );
+ edge->SetEpsilon( this->m_Epsilon );
edge->TreatAsWeightOff( );
// Real random walker
// Keep maximum
typedef itk::MaximumImageFilter< _TSliceBinary > _TMax;
- _TMax::Pointer max = _TMax::New();
+ typename _TMax::Pointer max = _TMax::New();
max->SetInput( 0, rw->GetOutput( ) );
max->SetInput( 1, next );
max->Update( );
git://git.creatis.insa-lyon.fr / FrontAlgorithms.git / commitdiff