<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 &quot;feature&quot; 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&amp;view=log">http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&amp;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&amp;r1=1.16&amp;r2=1.17">http://public.kitware.com/cgi-bin/viewcvs.cgi/Code/BasicFilters/itkCastImageFilter.h?root=Insight&amp;r1=1.16&amp;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">&lt;<a href="mailto:lars-friedrich@gmx.net">lars-friedrich@gmx.net</a>&gt;</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&#39;s Doxygen-entry, this filter can obviously be used for casting a 3D-volume to a 2D-slice (dimension-reduction).<br>
<br>
&quot;... If you attempt to cast an image to a lower dimension, the first &quot;slice&quot; (or line or volume) will be extracted. ...&quot;<br>
<br>
I&#39;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 &lt; Superclass::InputImageDimension; ++i)<br>
      {<br>
      outputSpacing[i] = inputSpacing[i];<br>
      outputOrigin[i] = inputOrigin[i];<br>
      for (j=0; j &lt; Superclass::OutputImageDimension; j++)<br>
        {<br>
        if (j &lt; 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 &lt;iostream&gt;<br>
#include &lt;stdlib.h&gt;<br>
#include &lt;time.h&gt;<br>
<br>
#include &lt;itkCastImageFilter.h&gt;<br>
#include &lt;itkImage.h&gt;<br>
#include &lt;itkImageRegionIterator.h&gt;<br>
<br>
int main(int argc, char *argv[])<br>
{<br>
  std::cerr &lt;&lt; &quot;Constructing 3D image ...&quot;;<br>
  typedef itk::Image&lt;float, 3&gt; Image3DType;<br>
  typedef itk::ImageRegionIterator&lt;Image3DType&gt; 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-&gt;SetRegions(region);<br>
  image-&gt;SetSpacing(spacing);<br>
  image-&gt;SetOrigin(origin);<br>
  image-&gt;SetDirection(orientation);<br>
  image-&gt;Allocate();<br>
<br>
  srand(time(NULL));<br>
  IteratorType it(image, region);<br>
  for (it.GoToBegin(); !it.IsAtEnd(); ++it)<br>
  {<br>
    it.Set(static_cast&lt;float&gt;(rand() % 100));<br>
  }<br>
  std::cerr &lt;&lt; &quot; DONE\n&quot;;<br>
<br>
  std::cerr &lt;&lt; &quot;Casting to 2D image (cast image filter) ...&quot;;<br>
  typedef itk::Image&lt;unsigned char, 2&gt; Image2DType;<br>
  typedef itk::CastImageFilter&lt;Image3DType, Image2DType&gt; CasterType;<br>
<br>
  CasterType::Pointer caster = CasterType::New();<br>
  caster-&gt;SetInput(image);<br>
  try<br>
  {<br>
    caster-&gt;Update(); // CAUSES SEGFAULT!!!<br>
    std::cerr &lt;&lt; &quot; DONE\n&quot;;<br>
  }<br>
  catch (itk::ExceptionObject &amp;e)<br>
  {<br>
    std::cerr &lt;&lt; &quot; ERROR\n&quot;;<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>