[Insight-users] 3D segmentation of LUNGS ( Region growing algorithm ) Need help

wassim_belhadj at topnet.tn wassim_belhadj at topnet.tn
Wed Mar 4 06:27:58 EST 2009


Hi,

Thank you very much for your help.

I was able to do the segmentation on an "unsigned char" raw file ( 512 *
512 * 494 ).

The filter gave me a binary mask of the two lungs.

I tried to apply the same algorithm to separate 2D dicom (512 *512) that
form the Raw file. 

These Dicom file are encoded in 16 bit (signed short).

I changed only the lowerThreshold (-950) and the  upperThreshold (-550) and
SetReplaceValue to 1100.

I have not had the same result. Some images contain only the left side of
the lungs. 

That was not the case on the images coded on 8 bit.

I doubt that itk::Image does not support "signed short" type.


Can you help me please.

Thank you in advance

/***********************/
   typedef    signed short			InternalPixelType; 
   typedef    signed short			OutputPixelType; 
   const	     unsigned int			Dimension = 2 ; 

   typedef itk::Image< InternalPixelType, Dimension >  InternalImageType; 
   typedef itk::Image< OutputPixelType, Dimension > OutputImageType; 

    ....
/***********************/




On Mon, 23 Feb 2009 07:29:43 -0500, Luis Ibanez <luis.ibanez at kitware.com>
wrote:
> Hi Wassim,
> 
> 
> The error below simply means that the type of ImageFileWriter
> that you are instantiating doesn't match the output image type
> of the the connectedThreshold filter.
> 
> 
> Please change the declaration:
> 
>  >> >   typedef itk::ConnectedThresholdImageFilter< InternalImageType,
>  >> > InternalImageType > ConnectedFilterType;
> 
> 
> to
> 
> 
>  >> >   typedef itk::ConnectedThresholdImageFilter< InternalImageType,
>  >> > OutputImageType > ConnectedFilterType;
> 
> 
> So that the output type of the region growing filter match
> the expected input type of the writer.
> 
> 
> 
>     Regards,
> 
> 
>         Luis
> 
> 
> 
> --------------------------------
> wassim_belhadj at topnet.tn wrote:
>> Hi Luis,
>> 
>> You said me that I don't need a Cast filter. 
>> 
>> If I don't use a CastingFilterType::Pointer  I have  this Error message
>> at
>> line :
>> 
>> writer->SetInput( connectedThreshold->GetOutput() )
>> 
>> cannot convert parameter 1 from 'itk::Image<TPixel,VImageDimension> *'
to
>> 'const itk::Image<TPixel,VImageDimension> *'
>> 
>> I am using ImageJ and trial Volview version for looking at the output
>> image.I have got the same result. 
>> 
>> 
>> I located my seed points  on a 2D Slice ( lungs parts).
>> 
>> Thank you in advance
>> 
>> 
>> 
>> On Fri, 20 Feb 2009 12:30:19 -0500, Luis Ibanez
<luis.ibanez at kitware.com>
>> wrote:
>> 
>>>Hi Wassim,
>>>
>>>
>>>     Important Rule:   Start Simple
>>>
>>>
>>>You don't need the Smoothing filter (until proven otherwise),
>>>You don't need the Cast filter
>>>You don't need the ITK-VTK import export code.
>>>
>>>
>>>
>>>You code could be reduced to:
>>>
>>>
>>> > #include "itkConnectedThresholdImageFilter.h"
>>> > #include "itkImage.h"
>>> > #include "itkImageFileReader.h"
>>> > #include "itkImageFileWriter.h"
>>> >
>>> > int main( int argc, char *argv[])
>>> > {
>>> >
>>> >   typedef   signed short         InternalPixelType;
>>> >   typedef   unsigned char        OutputPixelType;
>>> >   const     unsigned int		   Dimension = 3;
>>> >   typedef itk::Image< InternalPixelType, Dimension >
>>>                                                InternalImageType;
>>> >   typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
>>> >
>>> >   typedef  itk::ImageFileReader< InternalImageType > ReaderType;
>>> >   typedef  itk::ImageFileWriter<  OutputImageType  > WriterType;
>>> >
>>> >   ReaderType::Pointer reader = ReaderType::New();
>>> >   WriterType::Pointer writer = WriterType::New();
>>> >
>>> > 	reader->SetFileName( "c:/VTKEdge_Base_Test/images/Input.mhd" );
>>> > 	writer->SetFileName( "c:/VTKEdge_Base_Test/images/Output.mhd" );
>>> >
>>> >
>>> >   typedef itk::ConnectedThresholdImageFilter< InternalImageType,
>>> > InternalImageType > ConnectedFilterType;
>>> >
>>> >   ConnectedFilterType::Pointer connectedThreshold =
>>> > ConnectedFilterType::New();
>>> >
>>> >   connectedThreshold->SetInput( reader->GetOutput() );
>>> >   writer->SetInput( connectedThreshold->GetOutput() );
>>> >
>>> >
>>> >   const InternalPixelType lowerThreshold =  -700  ;
>>> >   const InternalPixelType upperThreshold = 200 ;
>>> >
>>> >   connectedThreshold->SetLower(  lowerThreshold  );
>>> >   connectedThreshold->SetUpper(  upperThreshold  );
>>> >
>>> >   connectedThreshold->SetReplaceValue( 255 );
>>> >
>>> >   InternalImageType::IndexType  index;
>>> >   index[0] = 180;
>>> >   index[1] =  296 ;
>>> >
>>> >   connectedThreshold->SetSeed( index );
>>> >
>>> >   InternalImageType::IndexType index2;
>>> >   index2[0] = 153;
>>> >   index2[1] = 274;
>>> >
>>> >
>>> >   connectedThreshold->SetSeed( index2 );
>>> >
>>> >   try
>>> >     {
>>> >     writer->Update();
>>> >     }
>>> >   catch( itk::ExceptionObject & excep )
>>> >     {
>>> >     std::cerr << "Exception caught !" << std::endl;
>>> >     std::cerr << excep << std::endl;
>>> >     }
>>> >
>>> >   return 0;
>>> > }
>>>
>>>
>>>
>>>
>>>BTW: in your code below you were calling SetInput()
>>>      in the threshold filter twice, with two different
>>>      inputs. Only the second would take place. This is
>>>      fixed in the code above.
>>>
>>>
>>>Also, please let us know what visualization tool
>>>       you are using for looking at the output image.
>>>       That is: How are you arriving to the conclusion
>>>                that the output image is all set to
>>>                zero values.
>>>
>>>
>>>If the image is actually all set to zero values, then
>>>the most common reason is that your seed points are
>>>located in incorrect places (e.g. in places whose
>>>intensity is not in the range -700HU to 200HU.
>>>
>>>
>>>
>>>     Regards,
>>>
>>>
>>>
>>>        Luis
>>>
>>>
>>>
>>>------------------------------------
>>>wassim_belhadj at topnet.tn wrote:
>>>
>>>>Hi Luis,
>>>>
>>>>Thank you very much for your help.
>>>>
>>>>I just followed the steps you gave me. It's a bit difficult for a
>>>>beginner
>>>>ITK user.
>>>>
>>>>
>>>>The segmentation was not conclusive (A black Volume result).
>>>>
>>>>
>>>>can you help me ?
>>>>
>>>>
>>>>Below my Input.mhd :
>>>>
>>>>
>> 
>>
/***************************************************************************************************/
>> 
>>
/***************************************************************************************************/
>> 
>>>>NDims = 3
>>>>DimSize = 512 512 494
>>>>ElementSpacing =0.777344 0.777344 0.7
>>>>Position = 0 0 0
>>>>ElementByteOrderMSB = False 
>>>>ElementType = MET_SHORT
>>>>HeaderSize = -1
>>>>ElementDataFile = prodhom.raw
>>>>
>>>>
>> 
>>
/***************************************************************************************************/
>> 
>>
/***************************************************************************************************/
>> 
>>>>This is the code that I used :
>>>>
>>>>
>> 
>>
/***************************************************************************************************/
>> 
>>
/***************************************************************************************************/
>> 
>>>>#include "itkConnectedThresholdImageFilter.h"
>>>>#include "itkImage.h"
>>>>#include "itkCastImageFilter.h"
>>>>#include "itkCurvatureFlowImageFilter.h"
>>>>#include "itkImageFileReader.h"
>>>>#include "itkImageFileWriter.h"
>>>>#include "itkVTKImageExport.h"
>>>>#include "itkVTKImageImport.h"
>>>>
>>>>int main( int argc, char *argv[])
>>>>{
>>>>
>>>>  typedef   signed short         InternalPixelType;
>>>>  typedef	  signed short		   OutputPixelType;
>>>>  const     unsigned int		   Dimension = 3;
>>>>  typedef itk::Image< InternalPixelType, Dimension > 
InternalImageType;
>>>>  typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
>>>>  typedef itk::CastImageFilter< InternalImageType, OutputImageType
>>>>
>>>>
>>>>>CastingFilterType;
>>>>
>>>>  CastingFilterType::Pointer caster = CastingFilterType::New();        

>> 
>> 
>>>>       
>>>>  typedef  itk::ImageFileReader< InternalImageType > ReaderType;
>>>>  typedef  itk::ImageFileWriter<  OutputImageType  > WriterType;
>>>>
>>>>  ReaderType::Pointer reader = ReaderType::New();
>>>>  WriterType::Pointer writer = WriterType::New();
>>>>
>>>>	reader->SetFileName( "c:/VTKEdge_Base_Test/images/Input.mhd" );
>>>>	writer->SetFileName( "c:/VTKEdge_Base_Test/images/Output.mhd" );
>>>>
>>>>  typedef itk::CurvatureFlowImageFilter< InternalImageType,
>>>>InternalImageType > CurvatureFlowImageFilterType;
>>>>
>>>>  CurvatureFlowImageFilterType::Pointer smoothing =
>>>>CurvatureFlowImageFilterType::New();
>>>>
>>>>  typedef itk::ConnectedThresholdImageFilter< InternalImageType,
>>>>InternalImageType > ConnectedFilterType;
>>>>
>>>>  ConnectedFilterType::Pointer connectedThreshold =
>>>>ConnectedFilterType::New();
>>>>
>>>>  smoothing->SetInput( reader->GetOutput() );
>>>>  connectedThreshold->SetInput( smoothing->GetOutput() );
>>>>  connectedThreshold->SetInput( reader->GetOutput() );
>>>>  caster->SetInput( connectedThreshold->GetOutput() );
>>>>  writer->SetInput( caster->GetOutput() );
>>>>
>>>>  smoothing->SetNumberOfIterations( 5 );
>>>>  smoothing->SetTimeStep( 0.125 );
>>>>
>>>>  const InternalPixelType lowerThreshold =  -700  ;
>>>>  const InternalPixelType upperThreshold = 200 ;
>>>>
>>>>  connectedThreshold->SetLower(  lowerThreshold  );
>>>>  connectedThreshold->SetUpper(  upperThreshold  );
>>>>
>>>>  connectedThreshold->SetReplaceValue( 255 );
>>>>
>>>>  InternalImageType::IndexType  index;
>>>>  index[0] = 180;
>>>>  index[1] =  296 ;
>>>>
>>>>  connectedThreshold->SetSeed( index );
>>>>  
>>>>  InternalImageType::IndexType index2;
>>>>  index2[0] = 153;
>>>>  index2[1] = 274;
>>>>
>>>>
>>>>  connectedThreshold->SetSeed( index2 );
>>>>
>>>>  try
>>>>    {
>>>>    writer->Update();
>>>>    }
>>>>  catch( itk::ExceptionObject & excep )
>>>>    {
>>>>    std::cerr << "Exception caught !" << std::endl;
>>>>    std::cerr << excep << std::endl;
>>>>    }
>>>>
>>>>  return 0;
>>>>}
>>>>
>> 
>>
/***************************************************************************************************/
>> 
>>
/***************************************************************************************************/
>> 
>>
/***************************************************************************************************/
>> 
>>>>Best regards, Wassim
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>
>>>>On Wed, 18 Feb 2009 10:20:55 -0500, Luis Ibanez
>> 
>> <luis.ibanez at kitware.com>
>> 
>>>>wrote:
>>>>
>>>>
>>>>>Hi Wassim,
>>>>>
>>>>>Please do the following:
>>>>>
>>>>>0) Extract a region of interest that includes the Trackea,
>>>>>   and the full Lungs, but not the head of the patient
>>>>>   (e.g. cut below the head, somewhere in the neck).
>>>>>   You can use the RegionOfInterestImageFilter for this.
>>>>>
>>>>
>>>>
>>
http://www.itk.org/Insight/Doxygen/html/classitk_1_1RegionOfInterestImageFilter.html
>> 
>>>>>1) Use the ConnectedThresholdImageFilter
>>>>>
>>>>
>>>>
>>
http://www.itk.org/Insight/Doxygen/html/classitk_1_1ConnectedThresholdImageFilter.html
>> 
>>>>>2) Set seed points inside the Trackea
>>>>>   (one seed point should be enough
>>>>>
>>>>>3) Set the lower threshold to -700 Hounsfield Units
>>>>>
>>>>>4) Set the upper threshodl to -200 Hounsfield Units
>>>>>
>>>>>5) Run the filter
>>>>>
>>>>>
>>>>>This shoudld give you a binary mask of the two lungs.
>>>>>
>>>>>If you find many holes inside the mask, you could
>>>>>clean it up with the Morphological filters or the
>>>>>Voting Hole filling filters.
>>>>>
>>>>>
>>>>>        Start typing....   :-)
>>>>>
>>>>>
>>>>>Once you are done, please keep in mind that we are
>>>>>happy to get "case studies" as papers to the Insight
>>>>>Journal. http://www.insight-journal.org/
>>>>>
>>>>>These, more than "papers", are actually "technical
>>>>>reports". Yours could be entitled "Using ITK for
>>>>>Lung segmentation".
>>>>>
>>>>>Please let us know if you need any help on writing
>>>>>this report.
>>>>>
>>>>>
>>>>>
>>>>>    Regards,
>>>>>
>>>>>
>>>>>       Luis
>>>>>
>>>>>
>>>>>
>>>>>------------
>>>>>wassim_belhadj at topnet.tn wrote:
>>>>>
>>>>>
>>>>>>Hi, 
>>>>>>
>>>>>>I'am using VTK and Qt4 for volume rendering.
>>>>>>
>>>>>>I'm trying to define values of pixel that define LUNGS.
>>>>>>
>>>>>>I assigned opacity to each pixel by checking its intensity value 
with
>>>>>>Hounsfield Units values.
>>>>>>
>>>>>>I  obtained a volume that contains the lungs and skin/air. 
>>>>>>
>>>>>>
>>>>>>I varied the opacity values to eliminate the skin that surrounds the
>>>>>>lungs
>>>>>>but I have not succeeded.
>>>>>>
>>>>>>VTK users asked me to do a segmentation (region growing algorithm )
to
>>>>>>display ONLY the lungs. 
>>>>>>
>>>>>>http://www.vtk.org/pipermail/vtkusers/2009-February/099524.html
>>>>>>
>>>>>>Can you help me? 
>>>>>>
>>>>>>How do I proceed (ITK-VTK)?
>>>>>>
>>>>>>Have you an example of 3D segmentation (region growing algorithm)?
>>>>>>
>>>>>>Regards, Wassim
>>>>>>_____________________________________
>>>>>>Powered by www.kitware.com
>>>>>>
>>>>>>Visit other Kitware open-source projects at
>>>>>>http://www.kitware.com/opensource/opensource.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