option(fpa_BUILD_CTBronchi "Build bronchi analysis from CT images applications?" OFF)
if(fpa_BUILD_CTBronchi)
+ configure_file(
+ CTBronchi_process.sh
+ ${PROJECT_BINARY_DIR}/CTBronchi_process.sh
+ COPYONLY
+ )
include_directories(
${PROJECT_SOURCE_DIR}/lib
${PROJECT_BINARY_DIR}/lib
--- /dev/null
+#!/bin/bash
+
+## Command line arguments
+if [ "$#" -lt 5 ]; then
+ echo "Usage: $0 input_image [index/point] seed_x seed_y seed_z"
+ exit 1
+fi
+
+exec_dir=`dirname $0`
+mori_seg=`dirname $0`/fpa_CTBronchi_MoriSegmentation
+mori_lab=`dirname $0`/fpa_CTBronchi_MoriLabelling
+random_walker=`dirname $0`/fpa_CTBronchi_RandomWalker
+input_image=$1
+seed_type=$2
+seed_x=$3
+seed_y=$4
+seed_z=$5
+
+base_name=`dirname $input_image`/`basename $input_image .mhd`
+
+mori_output_image="$base_name"_mori.mhd
+mori_output_signal="$base_name"_mori_signal.txt
+mori_init_threshold=-1024
+mori_end_threshold=0
+mori_delta=1
+mori_minimum_threshold=-850
+mori_inside_value=255
+mori_outside_value=0
+mori_signal_kernel_size=20
+mori_signal_threshold=500
+mori_signal_influence=500
+
+labels_output_image="$base_name"_labels.mhd
+label_upper_threshold=-600
+label_inside=1
+label_outside=2
+
+random_walker_output_image="$base_name"_rw.mhd
+random_walker_alpha=0
+random_walker_beta=20
+
+$mori_seg \
+ $input_image $mori_output_image $mori_output_signal \
+ $mori_init_threshold \
+ $mori_end_threshold \
+ $mori_delta \
+ $mori_minimum_threshold \
+ $mori_inside_value \
+ $mori_outside_value \
+ $mori_signal_kernel_size \
+ $mori_signal_threshold \
+ $mori_signal_influence \
+ $seed_type \
+ $seed_x $seed_y $seed_z
+
+$mori_lab \
+ $input_image $mori_output_image $labels_output_image \
+ $label_upper_threshold \
+ $mori_inside_value \
+ $label_inside \
+ $label_outside
+
+$random_walker \
+ $input_image $labels_output_image $random_walker_output_image \
+ $label_inside \
+ $random_walker_alpha \
+ $random_walker_beta
+
+## eof - $RCSfile$
#ifndef __CTBronchi__MoriLabelling__h__
#define __CTBronchi__MoriLabelling__h__
-#include <itkImageToImageFilter.h>
-#include <itkHessianRecursiveGaussianImageFilter.h>
-#include <itkHessian3DToVesselnessMeasureImageFilter.h>
-
-#include <itkImageFileWriter.h>
-
+#include <fpa/Base/RegionGrow.h>
+#include <fpa/Base/MarksInterface.h>
+#include <fpa/Base/Functors/RegionGrow/BinaryThreshold.h>
+#include <fpa/Image/Algorithm.h>
+#include <fpa/Image/LabelledSeedsInterface.h>
namespace CTBronchi
{
*/
template< class _TInputImage, class _TLabelImage >
class MoriLabelling
- : public itk::ImageToImageFilter< _TLabelImage, _TLabelImage >
+ : public fpa::Base::RegionGrow< fpa::Image::Algorithm< _TInputImage, _TLabelImage, fpa::Base::MarksInterface< typename _TInputImage::IndexType >, fpa::Image::LabelledSeedsInterface< typename _TInputImage::IndexType, typename _TInputImage::PointType, typename _TInputImage::PixelType, typename _TLabelImage::PixelType, typename _TLabelImage::PixelType, typename _TInputImage::IndexType::LexicographicCompare > > >
{
public:
- typedef MoriLabelling Self;
- typedef itk::ImageToImageFilter< _TLabelImage, _TLabelImage > Superclass;
- typedef itk::SmartPointer< Self > Pointer;
- typedef itk::SmartPointer< const Self > ConstPointer;
-
typedef _TInputImage TInputImage;
typedef _TLabelImage TLabelImage;
- typedef typename TInputImage::PixelType TPixel;
- typedef typename TLabelImage::PixelType TLabel;
- typedef typename TLabelImage::RegionType TRegion;
+ typedef typename TInputImage::PixelType TInputValue;
+ typedef typename TInputImage::PointType TPoint;
+ typedef typename TInputImage::IndexType TVertex;
+ typedef typename TLabelImage::PixelType TOutputValue;
+ typedef typename TVertex::LexicographicCompare TVertexCompare;
- public:
- itkNewMacro( Self );
- itkTypeMacro( MoriLabelling, itk::ImageToImageFilter );
+ typedef fpa::Base::MarksInterface< TVertex > TMarksInterface;
+ typedef fpa::Image::LabelledSeedsInterface< TVertex, TPoint, TInputValue, TOutputValue, TOutputValue, TVertexCompare > TSeedsInterface;
+ typedef fpa::Image::Algorithm< TInputImage, TLabelImage, TMarksInterface, TSeedsInterface > TAlgorithm;
- itkGetConstMacro( UpperThreshold, TPixel );
- itkSetMacro( UpperThreshold, TPixel );
+ typedef MoriLabelling Self;
+ typedef fpa::Base::RegionGrow< TAlgorithm > Superclass;
+ typedef itk::SmartPointer< Self > Pointer;
+ typedef itk::SmartPointer< const Self > ConstPointer;
- itkGetConstMacro( InsideValue, TLabel );
- itkSetMacro( InsideValue, TLabel );
+ typedef typename TSeedsInterface::TNode TNode;
+ typedef typename TSeedsInterface::TNodes TNodes;
- itkGetConstMacro( InsideLabel, TLabel );
- itkSetMacro( InsideLabel, TLabel );
+ typedef fpa::Base::Functors::RegionGrow::BinaryThreshold< TInputValue > TThresholdFunction;
- itkGetConstMacro( OutsideLabel, TLabel );
- itkSetMacro( OutsideLabel, TLabel );
+ public:
+ itkNewMacro( Self );
+ itkTypeMacro( MoriLabelling, fpa::Base::RegionGrow );
+
+ itkGetConstMacro( InsideLabel, TOutputValue );
+ itkSetMacro( InsideLabel, TOutputValue );
+
+ itkGetConstMacro( OutsideLabel, TOutputValue );
+ itkSetMacro( OutsideLabel, TOutputValue );
public:
const TLabelImage* GetInputLabelImage( ) const;
const TInputImage* GetInputRawImage( ) const;
void SetInputRawImage( TInputImage* image );
+ TInputValue GetUpperThreshold( ) const;
+ void SetUpperThreshold( TInputValue t );
+
protected:
MoriLabelling( );
virtual ~MoriLabelling( );
- virtual void BeforeThreadedGenerateData( ) override
- {
- this->Superclass::BeforeThreadedGenerateData( );
-
- const TInputImage* raw = this->GetInputRawImage( );
- typename TInputImage::SpacingType spac = raw->GetSpacing( );
- double sigma = spac[ 0 ];
- for( unsigned int d = 1; d < TInputImage::ImageDimension; ++d )
- sigma = ( spac[ d ] < sigma )? spac[ d ]: sigma;
- sigma *= 1.5;
-
- typedef itk::HessianRecursiveGaussianImageFilter< TInputImage > _THessian;
- typename _THessian::Pointer hessian = _THessian::New( );
- hessian->SetInput( raw );
- hessian->SetSigma( sigma );
-
- typedef itk::Hessian3DToVesselnessMeasureImageFilter< double > _TVesselness;
- typename _TVesselness::Pointer vesselness = _TVesselness::New( );
- vesselness->SetInput( hessian->GetOutput( ) );
- vesselness->SetAlpha1( 0.5 );
- vesselness->SetAlpha2( 2.0 );
-
- typename itk::ImageFileWriter< typename _TVesselness::OutputImageType >::Pointer w =
- itk::ImageFileWriter< typename _TVesselness::OutputImageType >::New( );
- w->SetInput( vesselness->GetOutput( ) );
- w->SetFileName( "vessel.mhd" );
- w->Update( );
- }
-
- virtual void ThreadedGenerateData( const TRegion& region, itk::ThreadIdType threadId ) override;
+ virtual TNodes _UnifySeeds( ) override;
+ virtual void _UpdateOutputValue( TNode& n ) override;
private:
// Purposely not implemented.
Self& operator=( const Self& other );
protected:
- TPixel m_UpperThreshold;
- TLabel m_InsideValue;
- TLabel m_InsideLabel;
- TLabel m_OutsideLabel;
+ TOutputValue m_InsideLabel;
+ TOutputValue m_OutsideLabel;
};
} // ecapseman
#ifndef __CTBronchi__MoriLabelling__hxx__
#define __CTBronchi__MoriLabelling__hxx__
+#include <itkImageRegionConstIteratorWithIndex.h>
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+const typename CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+TLabelImage* CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+GetInputLabelImage( ) const
+{
+ return( this->GetLabels( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+void CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+SetInputLabelImage( TLabelImage* image )
+{
+ this->SetLabels( image );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+const typename CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+TInputImage* CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+GetInputRawImage( ) const
+{
+ return( this->GetInput( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+void CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+SetInputRawImage( TInputImage* image )
+{
+ this->SetInput( image );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+typename CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+TInputValue CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+GetUpperThreshold( ) const
+{
+ const TThresholdFunction* func =
+ dynamic_cast< const TThresholdFunction* >( this->GetValuePredicate( ) );
+ if( func != NULL )
+ return( func->GetUpper( ) );
+ else
+ return( TInputValue( 0 ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+void CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+SetUpperThreshold( TInputValue t )
+{
+ TThresholdFunction* func =
+ dynamic_cast< TThresholdFunction* >( this->GetValuePredicate( ) );
+ if( func != NULL )
+ func->SetUpper( t );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+MoriLabelling( )
+ : Superclass( ),
+ m_InsideLabel( TOutputValue( 0 ) ),
+ m_OutsideLabel( TOutputValue( 0 ) )
+{
+ typename TThresholdFunction::Pointer func = TThresholdFunction::New( );
+ this->SetPredicate( func );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+~MoriLabelling( )
+{
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+typename CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+TNodes CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+_UnifySeeds( )
+{
+ this->m_Seeds.clear( );
+ const TLabelImage* lbl = this->GetLabels( );
+ if( lbl == NULL )
+ {
+ std::ostringstream msg;
+ msg << "itk::ERROR: CTBronchi::MoriLabelling (" << this
+ << "): Labelled image not defined.";
+ ::itk::ExceptionObject e(
+ __FILE__, __LINE__, msg.str( ).c_str( ), ITK_LOCATION
+ );
+ throw e;
+
+ } // fi
+
+ // Iterate over labels
+ typename TLabelImage::RegionType reg = lbl->GetRequestedRegion( );
+ itk::ImageRegionConstIteratorWithIndex< TLabelImage > lIt( lbl, reg );
+ for( lIt.GoToBegin( ); !lIt.IsAtEnd( ); ++lIt )
+ {
+ if( lIt.Get( ) > 0 )
+ {
+ bool is_seed = false;
+ for( unsigned int d = 0; d < TLabelImage::ImageDimension; ++d )
+ {
+ for( int s = -1; s <= 1; s += 2 )
+ {
+ TVertex neigh = lIt.GetIndex( );
+ neigh[ d ] += s;
+ if( reg.IsInside( neigh ) )
+ is_seed |= ( lbl->GetPixel( neigh ) == 0 );
+
+ } // rof
+
+ } // rof
+
+ if( !is_seed )
+ {
+ typename TSeedsInterface::TNode node;
+ node.Vertex = lIt.GetIndex( );
+ node.Parent = lIt.GetIndex( );
+ node.FrontId = lIt.Get( );
+ node.Value = this->m_InsideLabel;
+ this->_Mark( node.Vertex, node.FrontId );
+ this->_UpdateOutputValue( node );
+ }
+ else
+ {
+ typename TSeedsInterface::TSeed seed;
+ seed.Vertex = lIt.GetIndex( );
+ seed.IsPoint = false;
+ seed.FrontId = lIt.Get( );
+ this->m_Seeds.push_back( seed );
+
+ } // fi
+
+ } // fi
+
+ } // rof
+
+ // Ok, finish initialization
+ return( this->Superclass::_UnifySeeds( ) );
+}
+
+// -------------------------------------------------------------------------
+template< class _TInputImage, class _TLabelImage >
+void CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
+_UpdateOutputValue( TNode& n )
+{
+ const TInputImage* input = this->GetInputRawImage( );
+ const TLabelImage* labels = this->GetInputLabelImage( );
+}
+
+
+/* TODO
+ TOutputValue m_InsideValue;
+ TOutputValue m_InsideLabel;
+ TOutputValue m_OutsideLabel;
+*/
+
+
+
+
+
+
+
+
+
+
+/* TODO
#include <itkImageRegionConstIterator.h>
#include <itkImageRegionIterator.h>
CTBronchi::MoriLabelling< _TInputImage, _TLabelImage >::
MoriLabelling( )
: Superclass( ),
- m_UpperThreshold( TLabel( 0 ) ),
- m_InsideValue( TLabel( 0 ) ),
- m_InsideLabel( TLabel( 0 ) ),
- m_OutsideLabel( TLabel( 0 ) )
+ m_UpperThreshold( TOutputValue( 0 ) ),
+ m_InsideValue( TOutputValue( 0 ) ),
+ m_InsideLabel( TOutputValue( 0 ) ),
+ m_OutsideLabel( TOutputValue( 0 ) )
{
this->SetNumberOfRequiredInputs( 2 );
}
if( this->m_UpperThreshold < iIt.Get( ) )
oIt.Set( this->m_OutsideLabel );
else
- oIt.Set( TLabel( 0 ) );
+ oIt.Set( TOutputValue( 0 ) );
}
else
oIt.Set( this->m_InsideLabel );
} // elihw
}
+*/
#endif // __CTBronchi__MoriLabelling__hxx__
#include <fpa/Image/RandomWalker.h>
#include <fpa/Image/Functors/Dijkstra/Gaussian.h>
+#include <ivq/ITK/ExtractLabelFunction.h>
+#include <ivq/ITK/ImageUnaryFunctionFilter.h>
// -------------------------------------------------------------------------
const unsigned int Dim = 3;
typedef short TPixel;
typedef unsigned short TLabel;
+typedef unsigned char TBinary;
typedef double TScalar;
typedef itk::Image< TPixel, Dim > TInputImage;
typedef itk::Image< TLabel, Dim > TLabelImage;
+typedef itk::Image< TBinary, Dim > TBinaryImage;
typedef fpa::Image::RandomWalker< TInputImage, TLabelImage, TScalar > TFilter;
typedef fpa::Image::Functors::Dijkstra::Gaussian< TInputImage, TScalar > TWeight;
+typedef ivq::ITK::ExtractLabelFunction< TLabel, TBinary > TLabelFunction;
+typedef ivq::ITK::ImageUnaryFunctionFilter< TLabelImage, TBinaryImage > TLabelExtractor;
+
// -------------------------------------------------------------------------
int main( int argc, char* argv[] )
{
// Get arguments
- if( argc < 6 )
+ if( argc < 7 )
{
std::cerr
<< "Usage: " << argv[ 0 ] << std::endl
<< " input_image label_image output_image" << std::endl
- << " alpha(0) beta(100)"
+ << " label alpha(0) beta(100)"
<< std::endl;
return( 1 );
std::string input_image_filename = argv[ 1 ];
std::string label_image_filename = argv[ 2 ];
std::string output_image_filename = argv[ 3 ];
- double alpha = std::atof( argv[ 4 ] );
- double beta = std::atof( argv[ 5 ] );
+ TLabel label = TLabel( std::atoi( argv[ 4 ] ) );
+ double alpha = std::atof( argv[ 5 ] );
+ double beta = std::atof( argv[ 6 ] );
// Read images
itk::ImageFileReader< TInputImage >::Pointer input_image_reader =
filter->SetLabels( label_image_reader->GetOutput( ) );
filter->SetWeightFunction( weight );
+ // Label extraction
+ TLabelFunction::Pointer label_function = TLabelFunction::New( );
+ label_function->SetLabel( label );
+ label_function->SetInsideValue( 255 );
+ label_function->SetOutsideValue( 0 );
+
+ TLabelExtractor::Pointer label_extractor = TLabelExtractor::New( );
+ label_extractor->SetInput( filter->GetMarks( ) );
+ label_extractor->SetFunction( label_function );
+
+
// Show some information
std::cout << "----------------------------------------------" << std::endl;
std::cout << "Image: " << input_image_filename << std::endl;
tr = te - ts;
std::cout
<< "Labelling time: " << tr.count( ) << " s" << std::endl;
+
+ ts = std::chrono::high_resolution_clock::now( );
+ label_extractor->Update( );
+ te = std::chrono::high_resolution_clock::now( );
+ tr = te - ts;
+ std::cout
+ << "Extract label time: " << tr.count( ) << " s" << std::endl;
}
catch( std::exception& err )
{
} // yrt
// Save output image
- itk::ImageFileWriter< TLabelImage >::Pointer output_image_writer =
- itk::ImageFileWriter< TLabelImage >::New( );
- output_image_writer->SetInput( filter->GetMarks( ) );
+ itk::ImageFileWriter< TBinaryImage >::Pointer output_image_writer =
+ itk::ImageFileWriter< TBinaryImage >::New( );
+ output_image_writer->SetInput( label_extractor->GetOutput( ) );
output_image_writer->SetFileName( output_image_filename );
try
{