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