[Insight-users] itkMergeLabelMapFilter.h
John Drozd
john.drozd at gmail.com
Thu Dec 17 13:12:35 EST 2009
Hi Gaetan,
My particular brain subject had the left and right brain ventricle as
separate regions, so I inputed a seed for the left ventricle and a seed for
the right ventricle, and created a separate label map for each ventricle.
(I could not figure out how to create a label map with two labelled
segmentations, so I am trying to merge these two label maps into a combined
label map.)
Below is my code and thank you:
/*
to run type:
./ConnectedThresholdImagechangedtoFuzzyConnectednessImageFilter
"/trumpet/downloads/ReadAtlasDicomSeriesAndReadSubjectDicomSeries_Plugin/ReadAtlasDicomSeriesAndReadSubjectDicomSeries/ReadAtlasDicomSeriesAndReadSubjectDicomSeries3DSlicerDec22009/normal002_S_0295.dcm"
outsubjectnormal4final.dcm outfuzzymapnormal.mhd 100 151 95 86 152 76 2.5 10
*/
#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 "itkLabelMap.h"
#include "itkBinaryImageToShapeLabelMapFilter.h"
//end of added code
//added per label.cxx
#include "itkLabelObject.h"
#include "itkLabelMap.h"
#include "itkBinaryImageToLabelMapFilter.h"
#include "itkLabelMapToLabelImageFilter.h"
//end of added code
#include "itkSimpleFilterWatcher.h"
//added per FuzzyConnectednessImageFilter.cxx
#include "itkSimpleFuzzyConnectednessScalarImageFilter.h"
#include "itkConfidenceConnectedImageFilter.h"
//end of added code
#include "itkMergeLabelMapFilter.h"
int main( int argc, char *argv[])
{
if( argc < 7 )
{
std::cerr << "Missing Parameters " << std::endl;
std::cerr << "Usage: " << argv[0];
std::cerr << " inputImage outputImage outputAffinityMap seedX seedY
seedZ seedX seedY seedZ multiplier tolerancedivisor" << std::endl;
return 1;
}
typedef unsigned char BinaryPixelType;
const unsigned int Dimension = 3;
typedef itk::Image< BinaryPixelType, Dimension > BinaryImageType;
typedef signed short InputPixelType;
typedef itk::Image< InputPixelType, Dimension > InputImageType;
typedef unsigned char LabelType;
const unsigned int dim = 3;
typedef itk::ShapeLabelObject< LabelType, dim > LabelObjectType;
typedef itk::LabelMap< LabelObjectType > LabelMapType;
typedef itk::LabelMap< LabelObjectType > LabelMapType2;
typedef itk::BinaryImageToShapeLabelMapFilter< BinaryImageType,
LabelMapType> ConverterType;
ConverterType::Pointer converter = ConverterType::New();
typedef itk::BinaryImageToShapeLabelMapFilter< BinaryImageType,
LabelMapType> ConverterType2;
ConverterType2::Pointer converter2 = ConverterType2::New();
typedef itk::ImageFileReader< InputImageType > ReaderType;
typedef itk::ImageFileWriter< BinaryImageType > WriterType;
ReaderType::Pointer reader = ReaderType::New();
WriterType::Pointer writer = WriterType::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::ConfidenceConnectedImageFilter<
InputImageType,
BinaryImageType
>
ConfidenceConnectedFilterType;
ConfidenceConnectedFilterType::Pointer confidenceConnectedFilter =
ConfidenceConnectedFilterType::New();
typedef itk::SimpleFuzzyConnectednessScalarImageFilter<
InputImageType,
BinaryImageType
>
FuzzySegmentationFilterType;
FuzzySegmentationFilterType::Pointer fuzzysegmenter =
FuzzySegmentationFilterType::New();
typedef FuzzySegmentationFilterType::FuzzySceneType FuzzySceneType;
typedef itk::ImageFileWriter< FuzzySceneType > FuzzyWriterType;
FuzzyWriterType::Pointer fwriter = FuzzyWriterType::New();
fwriter->SetFileName( argv[3] );
BinaryImageType::IndexType index;
index[0] = atoi( argv[4] );
index[1] = atoi( argv[5] );
index[2] = atoi( argv[6] );
const double varianceMultiplier = atof( argv[10] );
confidenceConnectedFilter->SetInput( reader->GetOutput() );
confidenceConnectedFilter->SetMultiplier( varianceMultiplier );
confidenceConnectedFilter->SetNumberOfIterations( 12 );
fuzzysegmenter->SetInput( reader->GetOutput() );
double varianceEstimation = 0;
double meanEstimation = 0;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation << std::endl;
std::cout << "Variance estimation = " << varianceEstimation << std::endl;
if( varianceEstimation > meanEstimation/atoi( argv[11] ) )
{
index[0] = index[0]-10;
index[1] = index[1]-10;
index[2] = index[2]-10;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
}
while( varianceEstimation > meanEstimation/atoi( argv[11] ) )
{
index[0] = index[0] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
if ( varianceEstimation > meanEstimation/atoi( argv[11] ))
{
index[1] = index[1] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
}
if ( varianceEstimation > meanEstimation/atoi( argv[11] ))
{
index[2] = index[2] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
}
}
std::cout << index << std::endl;
std::cout << "Mean estimation = " << meanEstimation << std::endl;
std::cout << "Variance estimation = " << varianceEstimation << std::endl;
fuzzysegmenter->SetObjectSeed( index );
fuzzysegmenter->SetMean( meanEstimation );
fuzzysegmenter->SetVariance( varianceEstimation );
fuzzysegmenter->SetThreshold( 0.5 );
fuzzysegmenter->Update();
//begin second label
FuzzySegmentationFilterType::Pointer fuzzysegmenter2 =
FuzzySegmentationFilterType::New();
index[0] = atoi( argv[7] );
index[1] = atoi( argv[8] );
index[2] = atoi( argv[9] );
confidenceConnectedFilter->SetInput( reader->GetOutput() );
confidenceConnectedFilter->SetMultiplier( varianceMultiplier );
confidenceConnectedFilter->SetNumberOfIterations( 12 );
fuzzysegmenter2->SetInput( reader->GetOutput() );
varianceEstimation = 0;
meanEstimation = 0;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation << std::endl;
std::cout << "Variance estimation = " << varianceEstimation << std::endl;
if( varianceEstimation > meanEstimation/atoi( argv[11] ) )
{
index[0] = index[0]-10;
index[1] = index[1]-10;
index[2] = index[2]-10;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
}
while( varianceEstimation > meanEstimation/atoi( argv[11] ) )
{
index[0] = index[0] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
if ( varianceEstimation > meanEstimation/atoi( argv[11] ))
{
index[1] = index[1] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
}
if ( varianceEstimation > meanEstimation/atoi( argv[11] ))
{
index[2] = index[2] + 1;
confidenceConnectedFilter->SetSeed( index );
confidenceConnectedFilter->Update();
meanEstimation = confidenceConnectedFilter->GetMean();
varianceEstimation = confidenceConnectedFilter->GetVariance();
std::cout << "Mean estimation = " << meanEstimation <<
std::endl;
std::cout << "Variance estimation = " << varianceEstimation <<
std::endl;
}
}
std::cout << index << std::endl;
std::cout << "Mean estimation = " << meanEstimation << std::endl;
std::cout << "Variance estimation = " << varianceEstimation << std::endl;
fuzzysegmenter2->SetObjectSeed( index );
fuzzysegmenter2->SetMean( meanEstimation );
fuzzysegmenter2->SetVariance( varianceEstimation );
fuzzysegmenter2->SetThreshold( 0.5 );
fuzzysegmenter2->Update();
//end second label
typedef itk::MetaDataDictionary DictionaryType;
DictionaryType inputdict = reader->GetMetaDataDictionary();
writer->SetMetaDataDictionary( inputdict );
writer->SetFileName( argv[2] );
converter->SetInput( fuzzysegmenter->GetOutput() );
converter->SetFullyConnected( true );
converter->SetInputForegroundValue( 255 );
converter->SetOutputBackgroundValue( 0 );
converter->Update();
converter2->SetInput( fuzzysegmenter2->GetOutput() );
converter2->SetFullyConnected( true );
converter2->SetInputForegroundValue( 255 );
converter2->SetOutputBackgroundValue( 0 );
converter2->Update();
LabelMapType::Pointer labelMap = converter->GetOutput();
LabelMapType2::Pointer labelMap2 = converter2->GetOutput();
std::cout.precision(20);
std::cout << "labelMap->GetNumberOfLabelObjects() = " <<
labelMap->GetNumberOfLabelObjects() << std::endl;
std::cout << "labelMap2->GetNumberOfLabelObjects() = " <<
labelMap2->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->GetSize() << "\t" <<
labelObject->GetPhysicalSize() << "\t" << labelObject->GetCentroid() <<
std::endl;
}
for( unsigned int label=1; label<=labelMap2->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 * labelObject2 = labelMap2->GetLabelObject( label
);
std::cout << label << "\t" << labelObject2->GetSize() << "\t" <<
labelObject2->GetPhysicalSize() << "\t" << labelObject2->GetCentroid() <<
std::endl;
}
//added to merge the two label maps
typedef itk::MergeLabelMapFilter< LabelMapType > ConverterMergeType;
ConverterMergeType::Pointer convertermerge = ConverterMergeType::New();
/*******************************
how do I merge labelMap and LabelMap2 ?
********************************/
//end
typedef itk::LabelMapToLabelImageFilter< LabelMapType, BinaryImageType >
L2IType;
L2IType::Pointer l2i = L2IType::New();
l2i->SetInput( converter->GetOutput() );
writer->SetInput( l2i->GetOutput() );
try
{
writer->Update();
}
catch( itk::ExceptionObject & excep )
{
std::cerr << "Exception caught !" << std::endl;
std::cerr << excep << std::endl;
}
fwriter->SetInput( fuzzysegmenter->GetFuzzyScene() );
fwriter->Update();
return 0;
}
2009/12/17 Gaëtan Lehmann <gaetan.lehmann at jouy.inra.fr>
>
> Le 17 déc. 09 à 18:22, John Drozd a écrit :
>
>
> Hi,
>>
>> Can anyone show me how to use the MergeLabelMapFilter to merge two label
>> maps?
>>
>
>
> Hi John,
>
> You can look at Testing/Code/Review/itkMergeLabelMapFilterTest1.cxx, as an
> example.
>
> Can you be more specific on what you want to do?
>
> Regards,
>
> Gaëtan
>
>
> --
> Gaëtan Lehmann
> Biologie du Développement et de la Reproduction
> INRA de Jouy-en-Josas (France)
> tel: +33 1 34 65 29 66 fax: 01 34 65 29 09
> http://voxel.jouy.inra.fr http://www.itk.org
> http://www.mandriva.org http://www.bepo.fr
>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20091217/198454dc/attachment-0001.htm>
More information about the Insight-users
mailing list