<div dir="ltr">Hi Simon,<div><br></div><div>Thanks for your very clear description of the problem you are facing.</div><div><br></div><div>The message that you are receiving is quite common, and it is the </div><div>result of the fact that when ITK generates a 2D image from a 3D </div>
<div>image, the direction matrix the 3D image (which is initially a 3x3 </div><div>matrix) must be reduced to a 2D direction matrix which is a 2x2 </div><div>matrix).</div><div><br></div><div>There are multiple strategies for performing this reduction, none</div>
<div>of which is "optimal" nor "correct", and hence, the choice is left</div><div>to the developer..... so... here you are :-)</div><div><br></div><div>An example of when this happens,</div><div>is if you have a 3D direction matrix such as</div>
<div><br></div><div>1 0 0 </div><div>0 0 -1</div><div>0 1 0</div><div><br></div><div>and with the simple strategy of keeping the first two dimensions,</div><div>we would end up in 2D with a matrix</div><div><br></div><div>
1 0 </div><div>0 0</div><div><br></div><div>which is singular.... </div><div>and hence you get the exception that you are observing now.</div><div><br></div><div><br></div><div>Printing out the values of your direction matrix </div>
<div>will be enlightening.</div><div><br></div><div>std::cout << reader->GetOutput()->GetDirection() << std::endl;</div><div><br></div><div><br></div><div><br></div><div>When you opt for:</div><div><br></div>
<div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">filter-></span><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">SetDirectionCollapseToIndentit</span><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">y();</span><br>
</div><div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px"><br></span></div><div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">you are forcing the filter to produce the 2D matrix </span></div>
<div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px"><br></span></div><div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">1 0 </span></div><div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">0 1</span></div>
<div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px"><br></span></div><div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px">which brings you one step closer to happiness.</span></div>
<div><span style="font-family:arial,sans-serif;font-size:12.666666984558105px"><br></span></div><div><br></div><div>What I see out of order in your code is the line</div><div><br></div><div>filter->InPlaceOn();<br></div>
<div><br></div><div>which doesn't work well in the Extract image</div><div>filter when you are extracting a 2D image out</div><div>of a 3D image.</div><div><br></div><div>The InPlaceOn() method is intended for saving</div>
<div>memory in the cases where an output image is</div><div>of the same type of an input image.</div><div><br></div><div>In your case, the two images are of different</div><div>dimension, and hence the output can not be</div>
<div>stored in the input image.</div><div><br></div><div><br></div><div><br></div><div> Hope this helps.</div><div><br></div><div><br></div><div>BTW, could you please elaborate more, on what</div><div>you mean by: "the slice spacing information is lost" ?</div>
<div><br></div><div>In principle, the reduction from 3D to 2D fully </div><div>collapses a third dimension, so the spacing is </div><div>indeed meaningless after the extraction.</div><div><br></div><div>One alternative is to Extract, without collapsing</div>
<div>dimension, and in this way the output image is</div><div>a 3D (not 2D) image, of a single slice.</div><div><br></div><div><br></div><div> Thanks</div><div><br></div><div> Luis</div><div><br></div></div><div class="gmail_extra">
<br><br><div class="gmail_quote">On Sat, Mar 8, 2014 at 3:01 PM, Simon Wimmer <span dir="ltr"><<a href="mailto:wimmersimon@gmail.com" target="_blank">wimmersimon@gmail.com</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div dir="ltr"><div class="gmail_quote"><div dir="ltr">Hello ITK community,<div><br></div><div>I've started to use ITK two weeks ago for an university project and have been quite satisfied with what I achieved with it so far.</div>
<div>But today, I stumbled across a problem that I could not resolve by consulting the existing resources on the web and old mailing list threads.</div>
<div>For our project we have to work with 3D MRI data that is provided as an .mhd file.</div><div>We need to extract slices orthogonal to every axis from that data.</div><div>For that matter I consulted the docs and found the example code in ImageReadExtractWrite.cxx which makes use of the ExtractImageFilter. This method works perfectly for extracting a slice in z direction as the example code does.</div>
<div>So I tried to make the code extract a slice in x direction by simply changing the lines</div><div><br></div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
size[2] = 0;<br>InputImageType::IndexType start = inputRegion.GetIndex();<br>const unsigned int sliceNumber = atoi( argv[3] );<br>start[2] = sliceNumber;</blockquote><div><br></div><div>to</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
size[0] = 0;<br>InputImageType::IndexType start = inputRegion.GetIndex();<br>const unsigned int sliceNumber = atoi( argv[3] );<br>start[0] = sliceNumber; </blockquote><div><br></div><div>But this code throws an exception:</div>
</div><div><br></div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Description: itk::ERROR: ExtractImageFilter(037F46F0): Invalid submatrix extracted for collapsed direction.</blockquote>
<div><br></div><div>If I change</div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">filter->SetDirectionCollapseToSubmatrix();</blockquote>
<div>to </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"> filter->SetDirectionCollapseToIndentity();</blockquote>
<div>the code runs cleanly, but produces an awkward result image because seemingly the slice spacing information is lost.</div></div><div><br></div><div>Do you have any hints on how I could resolve the issue/where I have gone wrong or can you point out different approaches I could try to pursue?</div>
<div><br></div><div>The whole code:</div><div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">#include "itkImageFileReader.h"<br>
#include "itkImageFileWriter.h"<br>#include "itkExtractImageFilter.h"<br>#include "itkImage.h"<br><br>int main( int argc, char ** argv )<br>{<br> // Verify the number of parameters in the command line<br>
if( argc < 3 )<br> {<br> std::cerr << "Usage: " << std::endl;<br> std::cerr << argv[0] << " input3DImageFile output2DImageFile " << std::endl;<br> std::cerr << " sliceNumber " << std::endl;<br>
return EXIT_FAILURE;<br> }<br> typedef signed short InputPixelType;<br> typedef signed short OutputPixelType;<br> typedef itk::Image< InputPixelType, 3 > InputImageType;<br> typedef itk::Image< OutputPixelType, 2 > OutputImageType;<br>
typedef itk::ImageFileReader< InputImageType > ReaderType;<br> typedef itk::ImageFileWriter< OutputImageType > WriterType;<br> // Here we recover the file names from the command line arguments<br> //<br>
const char * inputFilename = argv[1];<br> const char * outputFilename = argv[2];<br> ReaderType::Pointer reader = ReaderType::New();<br> WriterType::Pointer writer = WriterType::New();<br> reader->SetFileName( inputFilename );<br>
writer->SetFileName( outputFilename );<br> typedef itk::ExtractImageFilter< InputImageType,<br> OutputImageType > FilterType;<br> FilterType::Pointer filter = FilterType::New();<br>
filter->InPlaceOn();<br> filter->SetDirectionCollapseToSubmatrix();<br> reader->UpdateOutputInformation();<br> InputImageType::RegionType inputRegion =<br> reader->GetOutput()->GetLargestPossibleRegion();<br>
InputImageType::SizeType size = inputRegion.GetSize();<br> size[0] = 0;<br> InputImageType::IndexType start = inputRegion.GetIndex();<br> const unsigned int sliceNumber = atoi( argv[3] );<br> start[0] = sliceNumber;<br>
InputImageType::RegionType desiredRegion;<br> desiredRegion.SetSize( size );<br> desiredRegion.SetIndex( start );<br> filter->SetExtractionRegion( desiredRegion );<br> filter->SetInput( reader->GetOutput() );<br>
writer->SetInput( filter->GetOutput() );<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>}</blockquote></div><div><br></div><div>Regards,</div><div><br></div><div>Simon</div><div><br></div></div>
</div><br></div>
<br>_______________________________________________<br>
Community mailing list<br>
<a href="mailto:Community@itk.org">Community@itk.org</a><br>
<a href="http://public.kitware.com/cgi-bin/mailman/listinfo/community" target="_blank">http://public.kitware.com/cgi-bin/mailman/listinfo/community</a><br>
<br></blockquote></div><br></div>