<br>Hi Lars,<br><br>It looks like this is a case in which:<br><br><br>1) The documentation incorrectly recommended a misuse of the class.<br>2) The misuse was later detected as a buggy behavior<br>3) The buggy behavior has been identified and fixed.<br>
<br><br>The "feature" of reducing 3D to 2D dimensions only<br>works in this filter due to a happy coincidence resulting <br>from the way in which the image iterators are implemented.<br><br>Note for example, that the extracted slice can only be the <br>
one corresponding to the slowest changing index.<br><br><br>Instead of using the CastImageFilter, you should use the <br><br> ExtractImageFilter<br><br>which is explicitly designed for dealing with dimensionality<br>
reduction between the input image and the output image. <br><br>For example, extracting one 2D slice from a 3D image,<br>as illustrated in the Example:<br><br> Insight/Examples/IO/ImageReadExtractWrite.cxx<br><br>Note that the ExtractImageFilter will do both the pixel <br>
casting and the dimensionality reduction.<br><br>In the meantime we have fixed the documentation <br>of the itkCastImageFilter.<br><a href="http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&view=log">http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&view=log</a><br>
<br>with the following diff:<br><a href="http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&r1=1.16&r2=1.17">http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&r1=1.16&r2=1.17</a><br>
<br><br><br> Regards,<br><br><br> Luis<br><br><br>----------------------------------------------------------------------------------------------------<br><div class="gmail_quote">On Sun, Apr 18, 2010 at 10:52 AM, Lars Friedrich Lars <span dir="ltr"><<a href="mailto:lars-friedrich@gmx.net">lars-friedrich@gmx.net</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">Hello,<br>
<br>
as can be retrieved from itk::CastImageFilter's Doxygen-entry, this filter can obviously be used for casting a 3D-volume to a 2D-slice (dimension-reduction).<br>
<br>
"... If you attempt to cast an image to a lower dimension, the first "slice" (or line or volume) will be extracted. ..."<br>
<br>
I've been using this functionality for a while in an application. However, today I changed something in the code (at completely different code-fragment) and from this point my application started to throw segfaults when executing it.<br>
I traced the bug (in current ITK-CVS-version) and found out that casting from a 3D-slice to a 2D-image causes the segfault. Inspecting the code of itk::CastImageFilter and its superclass itk::UnaryFunctorImageFilter showed me that the code lines 97-112 of itkUnaryFunctorImageFilter.txx seem to cause the segfault:<br>
<br>
for (i=0; i < Superclass::InputImageDimension; ++i)<br>
{<br>
outputSpacing[i] = inputSpacing[i];<br>
outputOrigin[i] = inputOrigin[i];<br>
for (j=0; j < Superclass::OutputImageDimension; j++)<br>
{<br>
if (j < Superclass::InputImageDimension)<br>
{<br>
outputDirection[j][i] = inputDirection[j][i];<br>
}<br>
else<br>
{<br>
outputDirection[j][i] = 0.0;<br>
}<br>
}<br>
}<br>
<br>
When casting from 3D to 2D the InputImageDimension=3, but the OutputImageDimension=2. Writing to the 3rd elements of outputSpacing, outputOrigin and outputDirection causes the memory-violation.<br>
<br>
Please correct me if I am wrong.<br>
<br>
However, it is scary and funny at the same time that my application has worked without throwing segfaults and without any obvious errors for a few weeks ...<br>
<br>
The example program below shows the bug (if it is a bug).<br>
<br>
regards,<br>
<br>
lars<br>
<br>
<br>
// TestReductionCasting.cxx<br>
<br>
#include <iostream><br>
#include <stdlib.h><br>
#include <time.h><br>
<br>
#include <itkCastImageFilter.h><br>
#include <itkImage.h><br>
#include <itkImageRegionIterator.h><br>
<br>
int main(int argc, char *argv[])<br>
{<br>
std::cerr << "Constructing 3D image ...";<br>
typedef itk::Image<float, 3> Image3DType;<br>
typedef itk::ImageRegionIterator<Image3DType> IteratorType;<br>
<br>
Image3DType::RegionType region;<br>
Image3DType::SizeType size;<br>
Image3DType::IndexType start;<br>
Image3DType::SpacingType spacing;<br>
Image3DType::PointType origin;<br>
Image3DType::DirectionType orientation;<br>
<br>
size[0] = 100; size[1] = 200; size[2] = 1;<br>
start.Fill(0);<br>
region.SetIndex(start);<br>
region.SetSize(size);<br>
spacing[0] = 0.25; spacing[1] = 0.5; spacing[2] = 1;<br>
origin[0] = 100; origin[1] = 200; origin[2] = -50;<br>
orientation.SetIdentity();<br>
Image3DType::Pointer image = Image3DType::New();<br>
image->SetRegions(region);<br>
image->SetSpacing(spacing);<br>
image->SetOrigin(origin);<br>
image->SetDirection(orientation);<br>
image->Allocate();<br>
<br>
srand(time(NULL));<br>
IteratorType it(image, region);<br>
for (it.GoToBegin(); !it.IsAtEnd(); ++it)<br>
{<br>
it.Set(static_cast<float>(rand() % 100));<br>
}<br>
std::cerr << " DONE\n";<br>
<br>
std::cerr << "Casting to 2D image (cast image filter) ...";<br>
typedef itk::Image<unsigned char, 2> Image2DType;<br>
typedef itk::CastImageFilter<Image3DType, Image2DType> CasterType;<br>
<br>
CasterType::Pointer caster = CasterType::New();<br>
caster->SetInput(image);<br>
try<br>
{<br>
caster->Update(); // CAUSES SEGFAULT!!!<br>
std::cerr << " DONE\n";<br>
}<br>
catch (itk::ExceptionObject &e)<br>
{<br>
std::cerr << " ERROR\n";<br>
}<br>
<br>
return EXIT_SUCCESS;<br>
}<br>
<br>
--<br>
GRATIS für alle GMX-Mitglieder: Die maxdome Movie-FLAT!<br>
Jetzt freischalten unter <a href="http://portal.gmx.net/de/go/maxdome01" target="_blank">http://portal.gmx.net/de/go/maxdome01</a><br>
_____________________________________<br>
Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a><br>
<br>
Visit other Kitware open-source projects at<br>
<a href="http://www.kitware.com/opensource/opensource.html" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
<br>
Kitware offers ITK Training Courses, for more information visit:<br>
<a href="http://www.kitware.com/products/protraining.html" target="_blank">http://www.kitware.com/products/protraining.html</a><br>
<br>
Please keep messages on-topic and check the ITK FAQ at:<br>
<a href="http://www.itk.org/Wiki/ITK_FAQ" target="_blank">http://www.itk.org/Wiki/ITK_FAQ</a><br>
<br>
Follow this link to subscribe/unsubscribe:<br>
<a href="http://www.itk.org/mailman/listinfo/insight-users" target="_blank">http://www.itk.org/mailman/listinfo/insight-users</a><br>
</blockquote></div><br>