[Insight-users] Re-sort a 3D Volume

Gaëtan Lehmann gaetan.lehmann at jouy.inra.fr
Fri May 18 05:25:33 EDT 2007


Le 18 mai 07 à 09:40, Andreas Keil a écrit :

> I am using SliceBySliceImageFilter and it works just fine :-)
>

Great :-)

Please try to write a short review, so it can be included in the  
toolkit :-)

Gaëtan


> Andreas
>
>
> -----Original Message-----
> From: insight-users-bounces+andreas.keil=cs.tum.edu at itk.org
> [mailto:insight-users-bounces+andreas.keil=cs.tum.edu at itk.org] On  
> Behalf
> Of Gaëtan Lehmann
> Sent: Wednesday, May 16, 2007 23:19
> To: Insight Users
> Cc: Luis Ibanez; Roger Bramon Feixas
> Subject: Re: [Insight-users] Re-sort a 3D Volume
>
>
> Hi all,
>
> Just a comment for the future users (a with the hope to get some
> reviews): the SliceBySliceImageFilter should be able to do what is
> done here, but in a much simple and lisible way.
> That code should be easier to modify to work with uncommon cases  
> though.
>
> Gaëtan
>
>
>
> Le 16 mai 07 à 18:18, d f a écrit :
>
>> Hi all,
>>
>> I recently needed something just like what Roger was nice enough to
>> share last month (http://public.kitware.com/pipermail/insight-users/
>> 2007-April/021796.html ).  I modified it a bit for my needs, and
>> decided to send this along in case it helps anyone else (although
>> I'm not sure how often other people would find the additional
>> functionality useful).
>>
>> This just uses a basic Gaussian filter for the processing, but I
>> was also able to implement a full level set segmentation, operating
>> on 2D slices at a time, by using this as a template.  (In case you
>> are curious, I was comparing some tradeoffs of implementing 3D vs
>> 2D segmentation on 3D images)
>>
>> drf
>>
>>
>> //  This programs reads in a 3D image volume, extracts 2D slices,
>> //  and performs some 2D filtering on each. The results are inserted
>> //  into a 3D output volume of identical dimensions to the input.  
>> This
>> //  is used to perform 2D processing on 3D images without having to
>> //  write large numbers of 2D images to disk,  simplifying the batch
>> //  files used to perform processing.
>> //
>> //  Any 2D processing pipeline can be used here; the object that
>> //  passes its output to the TileImageFilter must be updated inside
>> //  the loop (lines 167,168).  In this case, a discrete gaussian  
>> image
>> //  filter is applied to the 2D images before being passed to the
>> tiler.
>> //
>> //  Additionally, 2D slices along the x,y, or z dimensions can be
>> //  operated on.  The input parameter is collapseDim (0,1,or 2).
>> //  For example, a value of 2 "collapses" (or steps through) the
>> //  z-dimension, yielding many 2D images from the x-y plane.
>> Likewise,
>> //  collapseDim=0 collapses x, yielding images from the y-z plane.
>> The
>> //  result is permuted so the output has the same dimensions as the
>> original.
>> //
>> //  Thanks to Luis Ibanez for the help and to Roger Bramon Feixas,
>> //  who actually coded the loop that puts the images into the tiler:
>> //  http://public.kitware.com/pipermail/insight-users/2007-April/
>> 021796.html
>>
>>
>> #if defined(_MSC_VER)
>> #pragma warning ( disable : 4786 )
>> #endif
>>
>> #ifdef __BORLANDC__
>> #define ITK_LEAN_AND_MEAN
>> #endif
>>
>> #include "itkImageFileReader.h"
>> #include "itkImageFileWriter.h"
>> #include "itkExtractImageFilter.h"
>> #include "itkTileImageFilter.h"
>> #include "itkDiscreteGaussianImageFilter.h "
>> #include "itkRescaleIntensityImageFilter.h"
>> #include "itkImage.h"
>> #include "itkPermuteAxesImageFilter.h"
>>
>>
>> int main( int argc, char ** argv )
>> {
>>   //  Check input arguments
>>   if( argc < 3 )
>>     {
>>     std::cerr << "Usage: " << std::endl;
>>     std::cerr << argv[0] << " input3DImageFile  output3DImageFile
>> collapseDim" << std::endl;
>>     return EXIT_FAILURE;
>>     }
>>
>>
>>   //  Typedefs
>>   typedef float         InputPixelType;
>>   typedef float         MiddlePixelType;
>>   typedef float         OutputPixelType;
>>   typedef unsigned char WritePixelType;
>>
>>   typedef itk::Image< InputPixelType,  3 >    InputImageType;
>>   typedef itk::Image< MiddlePixelType, 2 >    MiddleImageType;
>>   typedef itk::Image< OutputPixelType, 3 >    OutputImageType;
>>   typedef itk::Image< WritePixelType,  3 >    WriteImageType;
>>
>>   typedef itk::ImageFileReader< InputImageType  >  ReaderType;
>>   typedef itk::ImageFileWriter< WriteImageType  >  WriterType;
>>
>>   typedef itk::ExtractImageFilter<
>>                InputImageType, MiddleImageType  >  ExtractFilterType;
>>   typedef itk::DiscreteGaussianImageFilter<
>>                MiddleImageType, MiddleImageType >  FilterType;
>>   typedef itk::TileImageFilter<
>>                MiddleImageType, OutputImageType >  TileFilterType;
>>   typedef itk::RescaleIntensityImageFilter<
>>                OutputImageType, WriteImageType  >  RescaleFilterType;
>>   typedef itk::PermuteAxesImageFilter<
>>                OutputImageType                  >  PermuteFilterType;
>>
>>
>>   //  Input/Output filenames
>>   const char * inputFilename  = argv[1];
>>   const char * outputFilename = argv[2];
>>
>>
>>   // Reader/Writer instantiations
>>   ReaderType::Pointer reader = ReaderType::New();
>>   WriterType::Pointer writer = WriterType::New();
>>   reader->SetFileName( inputFilename  );
>>   writer->SetFileName( outputFilename );
>>
>>
>>   //  Read input file
>>   try
>>     {
>>     reader->Update();
>>     }
>>   catch( itk::ExceptionObject & err )
>>     {
>>     std::cerr << "ExceptionObject caught !" << std::endl;
>>     std::cerr << err << std::endl;
>>     return EXIT_FAILURE;
>>     }
>>
>>
>>   //  Find input image size
>>   InputImageType::RegionType inputRegion =
>>            reader->GetOutput()->GetLargestPossibleRegion();
>>   InputImageType::SizeType size = inputRegion.GetSize();
>>   const unsigned int collapseDim  = atoi( argv[3] );
>>   InputImageType::SizeType::SizeValueType numSlices = size
>> [collapseDim];
>>   size[collapseDim] = 0;   // collapse dimension ( 3D->2D )
>>
>>
>>   //  Set region for 2D extraction
>>   InputImageType::IndexType start = inputRegion.GetIndex();
>>   InputImageType::RegionType desiredRegion;
>>   desiredRegion.SetSize(  size  );
>>
>>
>>   //  Extract filter instantiation
>>   ExtractFilterType::Pointer extractFilter = ExtractFilterType::New 
>> ();
>>
>>
>>   //  2D Gaussian filter instantiation
>>   FilterType::Pointer filter = FilterType::New();
>>   const double gaussianVariance = 1;
>>   const unsigned int maxKernelWidth = 20;
>>   filter->SetVariance( gaussianVariance );
>>   filter->SetMaximumKernelWidth( maxKernelWidth );
>>
>>
>>   //  Tiling filter for inserting 2D results into 3D output volume
>>   //  layout=[1,1,0] is code to automatically create necessary
>> image size
>>   //     for a given number of slices
>>   TileFilterType::Pointer tileFilter = TileFilterType::New();
>>   TileFilterType::LayoutArrayType layout;
>>   layout[0] = 1;
>>   layout[1] = 1;
>>   layout[2] = 0;
>>   tileFilter->SetLayout( layout );
>>
>>
>>   //  Rescale used to cast float->char for image writing
>>   RescaleFilterType::Pointer rescaler = RescaleFilterType::New();
>>   rescaler->SetOutputMinimum(   0 );
>>   rescaler->SetOutputMaximum( 255 );
>>
>>
>>   //  Set up pipeline
>>   extractFilter->SetInput( reader->GetOutput() );
>>   filter->SetInput( extractFilter->GetOutput() );
>>   tileFilter->SetInput( filter->GetOutput() );
>>
>>
>>   //  Vector of pointers to the many 2D extracted images
>>   std::vector< MiddleImageType::Pointer > extracts;
>>
>>
>>   //  Extract 2D slices, perform processing
>>   for(unsigned int sliceNumber = 0; sliceNumber<numSlices;
>> sliceNumber++)
>>   {
>>       start[collapseDim] = sliceNumber;
>>       desiredRegion.SetIndex( start );
>>
>>       extractFilter->SetExtractionRegion( desiredRegion );
>>
>>       filter->Update();
>>       extracts.push_back( filter->GetOutput() );
>>       extracts.back()->DisconnectPipeline();
>>
>>       if (sliceNumber != 0)
>>       {
>>           tileFilter->PushBackInput( extracts.back() );
>>       }
>>   }
>>
>>   //  Update tiling filter
>>   tileFilter->Update();
>>
>>
>>   //  Permute image axes to match input
>>   PermuteFilterType::Pointer permute = PermuteFilterType::New();
>>   unsigned int permuteAxes[3];
>>   switch (collapseDim)
>>   {
>>     case 0:
>>         {
>>         permuteAxes[0] = 2;
>>         permuteAxes[1] = 0;
>>         permuteAxes[2] = 1;
>>         break;
>>         }
>>     case 1:
>>         {
>>         permuteAxes[0] = 0;
>>         permuteAxes[1] = 2;
>>         permuteAxes[2] = 1;
>>         break;
>>         }
>>     case 2:
>>         {
>>         permuteAxes[0] = 0;
>>         permuteAxes[1] = 1;
>>         permuteAxes[2] = 2;
>>         break;
>>         }
>>     default:
>>         {
>>         permuteAxes[0] = 0;
>>         permuteAxes[1] = 1;
>>         permuteAxes[2] = 2;
>>         break;
>>         }
>>   }
>>   permute->SetOrder( permuteAxes );
>>   permute->SetInput( tileFilter->GetOutput() );
>>
>>
>>   //  Remainder of pipeline
>>   //  Cast float->char for image writing
>>   rescaler->SetInput( permute->GetOutput() );
>>   writer->SetInput( rescaler->GetOutput() );
>>
>>
>>   //  Update pipeline, write file
>>   try
>>     {
>>     writer->Update();
>>     }
>>   catch( itk::ExceptionObject & err )
>>     {
>>     std::cerr << "ExceptionObject caught !" << std::endl;
>>     std::cerr << err << std::endl;
>>     return EXIT_FAILURE;
>>     }
>>
>>   return EXIT_SUCCESS;
>> }
>>
>>
>>
>>
>> On 4/9/07, Roger Bramon Feixas <rogerbramon at gmail.com > wrote:
>>> Hi Luis,
>>>
>>> Problem solved. Now, I would like share with us the most important
>>> parts of  my algorithm, if somebody has a similar problem:
>>>
>>> Befor the loop:
>>> ----------------------
>>>
>>> ExtractImageType::Pointer extractFilter = ExtractImageType::New();
>>> extractFilter->SetInput( myItkImageType3D );
>>>
>>> TileFilterType::Pointer tileFilter = TileFilterType::New();
>>>
>>> TileFilterType::LayoutArrayType layout;
>>> layout[0] = 1;
>>> layout[1] = 1;
>>> layout[2] = 0;
>>>
>>> tileFilter->SetLayout( layout );
>>>
>>> std::vector< ItkImageType2D::Pointer > extracts;
>>>
>>> Into the loop:
>>> ------------------
>>>
>>> extractFilter->SetExtractionRegion( region );
>>> extractFilter->Update();
>>> extracts.push_back( extractFilter->GetOutput() );
>>> extracts.back()->DisconnectPipeline();
>>> tileFilter->PushBackInput( extracts.back() );
>>>
>>> After the loop:
>>> -------------------
>>> tileFilter->Update();
>>>
>>>
>>>
>>> Thanks for all,
>>>
>>> Roger
>>>
>>>
>>> El 08/04/2007, a las 0:33, Luis Ibanez escribió:
>>>
>>>>
>>>> Hi Roger,
>>>>
>>>> You will need to run the ExtractImageFilter first in a loop,
>>>> and at each iteration extract a different slice.
>>>>
>>>> You should put the output of the ExtractImageFilter into a
>>>> std::vector of SmartPointers, and after pushing it in the
>>>> vector you *MUST* call DisconnectPipeline() on that image.
>>>>
>>>> Otherwise, all your "slice" images are using the *same
>> SmartPointer*.
>>>> That is, the ExtractFilter output is reusing the same instance
>> of the
>>>> image object every time.
>>>>
>>>> I'm guessing that what you see as output of the TileImageFilter
>>>> is that all the tiles contain a replication of the "Last" slice
>>>> that you extracted with the ExtractImageFilter.
>>>>
>>>>
>>>> Adding the call to DisconnectPipeline() will force the Extract
>>>> ImageFilter to create a new instance of the image object for
>>>> storing its output at every iteration.
>>>>
>>>>
>>>>    Regards,
>>>>
>>>>
>>>>       Luis
>>>>
>>>>
>>>> ---------------------------------
>>>> Roger Bramon Feixas wrote:
>>>>> Hi,
>>>>> I need to re-order the slices of a 3D volume. I tried to extract
>>>>> 2D  slices using a "itkExtractImageFilter" and create a new
>> volume
>>>>> using  "itkTileImageFilter" but it doesn't work. The extraction
>>>>> works good  because I write the output in an image and it's
>>>>> correct. The problem  is in "itkTileImageFilter". I receive data
>>>>> empty error when I want to  use tileImageFilter output. That
>> is my
>>>>> algorism:
>>>>> before the loop:
>>>>>  - I instantiate the TileImageFilter
>>>>>  - I set the layout in 1,1,0
>>>>> into the loop over the images:
>>>>>  - I allocated a different ExtractImageFilter for each image
>>>>>  - I fed the TileImageFilter with the output of each different
>>>>> ExtractImageFilter: tileFilter->PushBackInput( extractFilter-
>>>>>> GetOutput() ). I also tried to do that using: tileFilter-
>>>>>> SetInput ( i, extractFilter->GetOutput() )
>>>>> after the loop:
>>>>>  TileImageFilter->Update();
>>>>> After the loop, the layout is : 1,1,X where X is the number of
>>>>> slices  that I pushed. Also, I read in a insight-user mail that
>>>>> they solved a  similar problem saving all of ExtractimageFilter
>>>>> objects in a vector  but in my way it doesn't work.
>>>>> I hope anybody know anything about my problem.
>>>>> Thanks,
>>>>> Roger
>>>>> _______________________________________________
>>>>> Insight-users mailing list
>>>>> Insight-users at itk.org
>>>>> http://www.itk.org/mailman/listinfo/insight-users
>>>
>>> _______________________________________________
>>> Insight-users mailing list
>>> Insight-users at itk.org
>>> http://www.itk.org/mailman/listinfo/insight-users
>>>
>>
>> _______________________________________________
>> Insight-users mailing list
>> Insight-users at itk.org
>> http://www.itk.org/mailman/listinfo/insight-users
>
> --
> 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
>
>
>
>
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users

--
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



-------------- next part --------------
A non-text attachment was scrubbed...
Name: PGP.sig
Type: application/pgp-signature
Size: 186 bytes
Desc: =?ISO-8859-1?Q?Ceci_est_une_signature_=E9lectronique_PGP?=
Url : http://public.kitware.com/pipermail/insight-users/attachments/20070518/94190611/PGP-0001.pgp


More information about the Insight-users mailing list