[Insight-users] combining output from segmentation with label map attributes using Gaetan Lehmann's paper: "Label object representation and manipulation with ITK"
Luis Ibanez
luis.ibanez at kitware.com
Fri Dec 11 11:20:55 EST 2009
Hi John,
Have you saved the output of the connectedThreshold filter to
a file, and verified that indeed it contains the segmentation of
the Ventricles ?
Alos, you are instantiating the connectedThreshold filter by
using as output an image of pixel type Float.
typedef float InternalPixelType;
const unsigned int Dimension = 3;
typedef itk::Image< InternalPixelType, Dimension > InternalImageType;
typedef itk::ConnectedThresholdImageFilter<
InternalImageType,
InternalImageType > ConnectedFilterType;
Do you have any motivation for doing so ?
The output of a Connected Threshold filter is typically
a binary image, so you probably should do:
typedef unsigned char BinaryPixelType;
typedef itk::Image< BinaryPixelType, Dimension > BinaryImageType;
typedef itk::ConnectedThresholdImageFilter<
InternalImageType,
BinaryImageType > ConnectedFilterType;
and subsequently, you should instantiate the
BinaryImageToLabelMapFilter also using as
input this Binary image type.
Currently, you have:
typedef itk::BinaryImageToLabelMapFilter<
InternalImageType, LabelMapType > ConverterType;
.. and you should probably replace it with:
typedef itk::BinaryImageToLabelMapFilter<
BinaryImageType, LabelMapType > ConverterType;
Regards
Luis
----------------------------
On Tue, Dec 8, 2009 at 2:38 PM, John Drozd <john.drozd at gmail.com> wrote:
> Hello,
>
> I have segmented the ventricles from an image using a connected threshold
> image filter.
> I then tried to convert my outputted segmentation to a label map filter and
> then to get the label's attributes.
>
> I am using parts of the example code attribute_values.cxx that is described
> in section 9.4 of
> Gaetan Lehmann's paper: "Label object representation and manipulation with
> ITK"
>
> My problem is that my code is telling me that the number of labels is 0.
> How do I apply a label to my segmented ventricles?
>
> Below is my code:
>
> /*
> to run type:
> ./ConnectedThresholdImageFilter correctedsubject5.dcm outsubject5.dcm 103
> 142 95 17100 17300
> */
>
> #if defined(_MSC_VER)
> #pragma warning ( disable : 4786 )
> #endif
>
> #ifdef __BORLANDC__
> #define ITK_LEAN_AND_MEAN
> #endif
>
>
> #include "itkConnectedThresholdImageFilter.h"
>
> #include "itkImage.h"
> #include "itkCastImageFilter.h"
>
> #include "itkCurvatureFlowImageFilter.h"
>
> #include "itkImageFileReader.h"
> #include "itkImageFileWriter.h"
>
> #include "itkGDCMImageIO.h"
>
> #include "itkVersion.h"
>
> #include "itkOrientedImage.h"
> #include "itkMinimumMaximumImageFilter.h"
>
> #include "itkGDCMImageIO.h"
> #include "itkGDCMSeriesFileNames.h"
> #include "itkNumericSeriesFileNames.h"
>
> #include "itkImageSeriesReader.h"
> #include "itkImageSeriesWriter.h"
>
> #include "itkResampleImageFilter.h"
> #include "itkShiftScaleImageFilter.h"
>
> #include "itkIdentityTransform.h"
> #include "itkLinearInterpolateImageFunction.h"
>
> #include <itksys/SystemTools.hxx>
>
> #include "gdcm/src/gdcmFile.h"
> #include "gdcm/src/gdcmUtil.h"
>
> #include <string>
>
> //added per attribute_values.cxx
> #include "itkImageFileReader.h"
>
> #include "itkShapeLabelObject.h"
> //#include "itkStatisticsLabelObject.h"
>
> #include "itkLabelMap.h"
>
> #include "itkBinaryImageToShapeLabelMapFilter.h"
> //#include "itkBinaryImageToStatisticsLabelMapFilter.h"
>
> //end of added code
>
> int main( int argc, char *argv[])
> {
> if( argc < 7 )
> {
> std::cerr << "Missing Parameters " << std::endl;
> std::cerr << "Usage: " << argv[0];
> std::cerr << " inputImage outputImage seedX seedY seedZ lowerThreshold
> upperThreshold" << std::endl;
>
> return 1;
> }
>
> typedef float InternalPixelType;
>
> const unsigned int Dimension = 3;
>
> typedef itk::Image< InternalPixelType, Dimension > InternalImageType;
>
> typedef signed short OutputPixelType;
>
> typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
> typedef itk::Image< float, Dimension > OutputImageType2;
> typedef itk::CastImageFilter< InternalImageType, OutputImageType >
> CastingFilterType;
> CastingFilterType::Pointer caster = CastingFilterType::New();
>
> const unsigned int ImageDimension = 3;
> typedef signed short PixelType;
>
> typedef itk::Image< PixelType, ImageDimension > FixedImageType;
> typedef itk::Image< float, ImageDimension > FloatImageType;
>
> typedef itk::ImageFileReader< FixedImageType > ReaderType;
> typedef itk::ImageFileWriter< OutputImageType > WriterType;
> typedef itk::ImageFileWriter< FloatImageType > WriterType2;
>
> ReaderType::Pointer reader = ReaderType::New();
> WriterType::Pointer writer = WriterType::New();
> WriterType2::Pointer writer2 = WriterType2::New();
>
> typedef itk::GDCMImageIO ImageIOTypefixed;
> ImageIOTypefixed::Pointer gdcmImageIOfixed = ImageIOTypefixed::New();
> reader->SetImageIO( gdcmImageIOfixed );
>
> typedef itk::GDCMImageIO ImageIOTypefixed2;
> ImageIOTypefixed2::Pointer gdcmImageIOfixed2 = ImageIOTypefixed2::New();
>
> reader->SetFileName( argv[1] );
>
> reader->Update();
>
> typedef itk::CurvatureFlowImageFilter< InternalImageType,
> InternalImageType >
> CurvatureFlowImageFilterType;
>
> CurvatureFlowImageFilterType::Pointer smoothing =
> CurvatureFlowImageFilterType::New();
>
> typedef itk::ConnectedThresholdImageFilter< InternalImageType,
> InternalImageType > ConnectedFilterType;
>
> ConnectedFilterType::Pointer connectedThreshold =
> ConnectedFilterType::New();
>
> typedef signed short InputAPixelType;
> typedef float OutputBPixelType;
>
> typedef itk::Image< InputAPixelType, 3 > InputAImageType;
> typedef itk::Image< OutputBPixelType, 3 > OutputBImageType;
>
> typedef itk::CastImageFilter< InputAImageType, OutputBImageType >
> CastFilterType;
>
> CastFilterType::Pointer castFilter = CastFilterType::New();
>
>
> castFilter->SetInput( reader->GetOutput() );
>
>
> connectedThreshold->SetInput( castFilter->GetOutput() );
>
> caster->SetInput( connectedThreshold->GetOutput() );
>
>
> smoothing->SetNumberOfIterations( 20 ); //was 5
> smoothing->SetTimeStep( 0.125 );
>
> const InternalPixelType lowerThreshold = atof( argv[6] );
> const InternalPixelType upperThreshold = atof( argv[7] );
>
> connectedThreshold->SetLower( lowerThreshold );
> connectedThreshold->SetUpper( upperThreshold );
>
> connectedThreshold->SetReplaceValue( 255 );
>
> InternalImageType::IndexType index;
>
> index[0] = atoi( argv[3] );
> index[1] = atoi( argv[4] );
>
> //added
> index[2] = atoi( argv[5] );
>
> std::cout << index << std::endl;
>
> // Software Guide : BeginCodeSnippet
> connectedThreshold->SetSeed( index );
>
> //obtain a 5 x 5 bounding region of seeds
> int ii, jj, kk;
>
> ii = index[0];
> jj = index[1];
> kk = index[2];
>
> for (int i = ii; i < ii + 5; i++)
> for (int j = jj; j < jj + 5; j++)
> for (int k = kk; k < kk + 5; k++)
> {
>
> index[0] = i;
> index[1] = j;
> index[2] = k;
> connectedThreshold->AddSeed( index );
> }
>
> for (int i = ii; i > ii - 5; i--)
> for (int j = jj; j > jj - 5; j--)
> for (int k = kk; k > kk - 5; k--)
> {
>
> index[0] = i;
> index[1] = j;
> index[2] = k;
> connectedThreshold->AddSeed( index );
> }
>
> connectedThreshold->Print(std::cout,17100);
>
> typedef itk::MetaDataDictionary DictionaryType;
>
> DictionaryType inputdict = reader->GetMetaDataDictionary();
>
> writer->SetMetaDataDictionary( inputdict );
>
> writer->SetFileName( argv[2] );
>
> //added per attribute_values.cxx
>
> // define the object type. Here the ShapeLabelObject type
> // is chosen in order to read some attribute related to the shape
> // of the objects (by opposition to the content of the object, with
> // the StatisticsLabelObejct).
> typedef unsigned long LabelType;
>
> typedef itk::ShapeLabelObject< LabelType, 3 > LabelObjectType;
>
> typedef itk::LabelMap< LabelObjectType > LabelMapType;
>
> // convert the image in a collection of objects
> typedef itk::BinaryImageToLabelMapFilter< InternalImageType, LabelMapType
>> ConverterType;
>
> ConverterType::Pointer converter = ConverterType::New();
>
> converter->SetInput( connectedThreshold->GetOutput() );
>
> //converter->SetForegroundValue( atoi(argv[2]) );
> //converter->SetForegroundValue( 255 );
>
> // valuate the attributes with the dedicated filter: ShapeLabelMapFilter
>
> typedef itk::ShapeLabelMapFilter< LabelMapType > ShapeFilterType;
> ShapeFilterType::Pointer shape = ShapeFilterType::New();
> shape->SetInput( converter->GetOutput() );
>
> shape->Update();
>
>
> // then we can read the attribute values we're interested in. The
> BinaryImageToShapeLabelMapFilter
> // produce consecutive labels, so we can use a for loop and
> GetLabelObject() method to retrieve
> // the label objects. If the labels are not consecutive, the
> GetNthLabelObject() method must be
> // use instead of GetLabelObject(), or an iterator on the label object
> container of the label map.
> LabelMapType::Pointer labelMap = shape->GetOutput();
>
> std::cout << "Number of label objects = " <<
> labelMap->GetNumberOfLabelObjects() << std::endl;
>
>
> for( unsigned int label=1; label<=labelMap->GetNumberOfLabelObjects();
> label++ )
> {
> // we don't need a SmartPointer of the label object here, because the
> reference is kept in
> // in the label map.
> const LabelObjectType * labelObject = labelMap->GetLabelObject( label );
> std::cout << label << "\t" << labelObject->GetPhysicalSize() << "\t" <<
> labelObject->GetCentroid() << std::endl;
> }
>
> writer->SetInput( caster->GetOutput() );
>
> try
> {
> writer->Update();
> }
> catch( itk::ExceptionObject & excep )
> {
> std::cerr << "Exception caught !" << std::endl;
> std::cerr << excep << std::endl;
> }
>
> std::cout << "output from reader->GetOutput()->GetDirection()" <<
> std::endl;
> std::cout << reader->GetOutput()->GetDirection() << std::endl;
>
> std::cout << "output from castFilter->GetOutput()->GetDirection()" <<
> std::endl;
> std::cout << castFilter->GetOutput()->GetDirection() << std::endl;
>
> std::cout << "output from connectedThreshold->GetOutput()->GetDirection()"
> << std::endl;
> std::cout << connectedThreshold->GetOutput()->GetDirection() << std::endl;
>
> std::cout << "output from caster->GetOutput()->GetDirection()" <<
> std::endl;
> std::cout << caster->GetOutput()->GetDirection() << std::endl;
>
> return 0;
> }
>
>
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.html
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-users
>
>
More information about the Insight-users
mailing list