<br>Hi Michael, <br><br>Thanks for adding the Print() statements<br>and posting the outcome to the mailing list.<br><br>The image data look ok.<br>and so it does the filter metadata.<br><br>It is good that you also added the code <br>
that you are using for generating the images.<br><br>A common suspect when you using the <br>import filter is to get wrong the ownership<br>of the image pixel buffer.<br><br>You are using <br><br> const unsigned int numberOfPixels = diImage-&gt;getSize();<br>

 dType* data = diImage-&gt;getData();<br>
 import-&gt;SetImportPointer( data, numberOfPixels, false );<br><br>I&#39;ll strongly suggest that you verify that the <br>&quot;numberOfPixels&quot; is correct here, and also<br>reconsider whether the flag in SetImportPointer()<br>
should be &quot;false&quot; or &quot;true&quot;.  <br><br>When set to true, the import filter will own the memory<br>buffer, which means that it will deallocate it when the<br>importer is destroyed.<br><br><br>The process of importing the image is a lot more<br>
suspicious than the resampling.  <br><br>Please do the following:<br><br>    1) Before the resampling code, instantiate two<br>        ImageFileWriters, one for the fixed image,<br>        and the other one for the moving image.<br>
       and write the images to a fileformat for which<br>       you have a trusty image viewer.<br><br>        (you could use Slicer, VV, ParaView, SNAP...)<br><br>   2) Run the code and save the images, no need <br>        to run the resampler.<br>
<br>    3) Open the images with the viewer and let us<br>        know if they are in good shape.<br><br>    4) I would expect that the images are corrupted <br>        at this point, but we will only know for  sure<br>        after you run this test.<br>
<br><br>Please let us now what you find.<br><br><br>      Thanks<br><br><br>           Luis<br><br><br><br>----------------------------------------------------------<br><div class="gmail_quote">On Wed, Jul 1, 2009 at 8:28 AM, Michael Schildt <span dir="ltr">&lt;<a href="mailto:michael.schildt@ifn-magdeburg.de">michael.schildt@ifn-magdeburg.de</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;"><div><div></div><div class="h5">Hello Luis,<br>
<br>
thank you for your suggestions.<br>
<br>
I inserted the Prints in my source code and here ist the output:<br>
--------------------------- fixedIm -------------------------------<br>
Image (00339668)<br>
 RTTI typeinfo:   class itk::Image&lt;float,3&gt;<br>
 Reference Count: 3<br>
 Modified Time: 20<br>
 Debug: Off<br>
 Observers:<br>
   none<br>
 Source: (none)<br>
 Source output index:  0<br>
 Release Data: Off<br>
 Data Released: False<br>
 Global Release Data: Off<br>
 PipelineMTime: 13<br>
 UpdateMTime: 19<br>
 LargestPossibleRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 BufferedRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 RequestedRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 Spacing: [1, 1, 1]<br>
 Origin: [0, 0, 0]<br>
 Direction:<br>
1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 IndexToPointMatrix:<br>
 1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 PointToIndexMatrix:<br>
 1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 PixelContainer:<br>
   ImportImageContainer (01390EA8)<br>
     RTTI typeinfo:   class itk::ImportImageContainer&lt;unsigned long,float&gt;<br>
     Reference Count: 1<br>
     Modified Time: 17<br>
     Debug: Off<br>
     Observers:<br>
       none<br>
     Pointer: 034B0020<br>
     Container manages memory: false<br>
     Size: 16777216<br>
     Capacity: 16777216<br>
<br>
--------------------------- movingIm ------------------------------<br>
Image (00337178)<br>
 RTTI typeinfo:   class itk::Image&lt;float,3&gt;<br>
 Reference Count: 4<br>
 Modified Time: 40<br>
 Debug: Off<br>
 Observers:<br>
   none<br>
 Source: (none)<br>
 Source output index:  0<br>
 Release Data: Off<br>
 Data Released: False<br>
 Global Release Data: Off<br>
 PipelineMTime: 33<br>
 UpdateMTime: 39<br>
 LargestPossibleRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 BufferedRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 RequestedRegion:<br>
   Dimension: 3<br>
   Index: [0, 0, 0]<br>
   Size: [256, 256, 256]<br>
 Spacing: [1, 1, 1]<br>
 Origin: [0, 0, 0]<br>
 Direction:<br>
1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 IndexToPointMatrix:<br>
 1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 PointToIndexMatrix:<br>
 1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 PixelContainer:<br>
   ImportImageContainer (0033BE98)<br>
     RTTI typeinfo:   class itk::ImportImageContainer&lt;unsigned long,float&gt;<br>
     Reference Count: 1<br>
     Modified Time: 37<br>
     Debug: Off<br>
     Observers:<br>
       none<br>
     Pointer: 074C0020<br>
     Container manages memory: false<br>
     Size: 16777216<br>
     Capacity: 16777216<br>
<br>
--------------------------- resampler -----------------------------<br>
ResampleImageFilter (00338A08)<br>
 RTTI typeinfo:   class itk::ResampleImageFilter&lt;class<br>
itk::Image&lt;float,3&gt;,cla<br>
s itk::Image&lt;float,3&gt;,double&gt;<br>
 Reference Count: 2<br>
 Modified Time: 64<br>
 Debug: Off<br>
 Observers:<br>
   none<br>
 Number Of Required Inputs: 1<br>
 Number Of Required Outputs: 1<br>
 Number Of Threads: 1<br>
 ReleaseDataFlag: Off<br>
 ReleaseDataBeforeUpdateFlag: Off<br>
 Input 0: (00337178)<br>
 Output 0: (00338B88)<br>
 AbortGenerateData: Off<br>
 Progress: 0<br>
 Multithreader:<br>
   RTTI typeinfo:   class itk::MultiThreader<br>
   Reference Count: 1<br>
   Modified Time: 49<br>
   Debug: Off<br>
   Observers:<br>
     none<br>
   Thread Count: 1<br>
   Global Maximum Number Of Threads: 128<br>
   Global Default Number Of Threads: 1<br>
 DefaultPixelValue: 0<br>
 Size: [256, 256, 256]<br>
 OutputStartIndex: [0, 0, 0]<br>
 OutputOrigin: [0, 0, 0]<br>
 OutputSpacing: [1, 1, 1]<br>
 OutputDirection: 1 0 0<br>
0 1 0<br>
0 0 1<br>
<br>
 Transform: 00337350<br>
 Interpolator: 00338860<br>
 UseReferenceImage: Off<br>
<br>
I already call Update on the Images but after the Set..- commands. I&#39;m<br>
not sure about the order of commands to import the images to itk. Here<br>
is my import method called to create fixedIm and movingIm, maybe this<br>
helps to localize the problem (MRIScalarLayer is a custom class just<br>
containing the image data and some infos):<br>
<br>
template&lt;typename dType&gt; typename itk::Image&lt;dType,3&gt;::Pointer<br>
dataInterface2itk(MRIScalarLayer&lt;dType&gt; *diImage)<br>
{<br>
 typedef itk::Image&lt;dType,3&gt; Image3DType;<br>
 typename Image3DType::Pointer image = Image3DType::New();<br>
 typedef itk::ImportImageFilter&lt;dType,3&gt; ImportFilterImage3DType;<br>
 typename ImportFilterImage3DType::Pointer import =<br>
ImportFilterImage3DType::New();<br>
 typename ImportFilterImage3DType::SizeType size;<br>
 size[0] = diImage-&gt;getSize((unsigned char)0); // size along X<br>
 size[1] = diImage-&gt;getSize((unsigned char)1); // size along Y<br>
 size[2] = diImage-&gt;getSize((unsigned char)2); // size along Z<br>
 typename ImportFilterImage3DType::IndexType start;<br>
 start.Fill( 0 );<br>
 typename ImportFilterImage3DType::RegionType region;<br>
 region.SetIndex( start );<br>
 region.SetSize( size );<br>
 import-&gt;SetRegion( region );<br>
 double spacing[ 3 ];<br>
 // hier sollte noch das gap mit beachtet werden, welches es nicht in<br>
itk gibt<br>
 spacing[0] = diImage-&gt;getSpacing((unsigned char)0); // along X direction<br>
 spacing[1] = diImage-&gt;getSpacing((unsigned char)1); // along Y direction<br>
 spacing[2] = diImage-&gt;getSpacing((unsigned char)2); // along Z direction<br>
 import-&gt;SetSpacing( spacing );<br>
 double origin[ 3 ];<br>
 origin[0] = diImage-&gt;getOrigin((unsigned char)0); // X coordinate<br>
 origin[1] = diImage-&gt;getOrigin((unsigned char)1); // Y coordinate<br>
 origin[2] = diImage-&gt;getOrigin((unsigned char)2); // Z coordinate<br>
 import-&gt;SetOrigin( origin );<br>
<br>
 const unsigned int numberOfPixels = diImage-&gt;getSize();<br>
 dType* data = diImage-&gt;getData();<br>
 import-&gt;SetImportPointer( data, numberOfPixels, false );<br>
 import-&gt;Update();<br>
 image = import-&gt;GetOutput();<br>
 image-&gt;Update();<br>
<br>
//  cout &lt;&lt; image;<br>
 return image;<br>
}<br>
<br>
Best requards,<br>
   Michael Schildt<br>
<br>
Luis Ibanez schrieb:<br>
</div></div><blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;">
<br><div class="im">
Hi Michael,<br>
<br>
<br>
       This looks like a bug in your code      :-(<br>
<br>
<br>
<br>
Typical suspects are:<br>
<br>
<br>
A) The fixed image is not available yet at<br>
     the time you call GetLargestPossibleRegion(),<br>
     GetOrigin(), GetSpacing().<br>
<br>
      You should call  <br>
                    fixedIm-&gt;Update()<br>
<br>
        before you call any of the above methods.<br>
<br>
<br>
We will learn a lot about other possible causes if<br>
you add to your code the following:<br>
<br>
        fixedIm-&gt;Print( std::cout );<br>
        movingIm-&gt;Print( std::cout );<br>
        resampler-&gt;Print( std::cout );<br>
<br>
before the line<br>
<br>
         resampler-&gt;Update();<br>
<br>
and if you post the resulting output text to the<br>
mailing list.<br>
<br>
<br>
      Please let us know what you find.<br>
<br>
<br>
           Thanks<br>
<br>
<br>
                 Luis<br>
<br>
<br>
----------------------------------------------------------------<br></div><div><div></div><div class="h5">
On Fri, Jun 26, 2009 at 6:42 AM, Michael Schildt &lt;<a href="mailto:michael.schildt@ifn-magdeburg.de" target="_blank">michael.schildt@ifn-magdeburg.de</a> &lt;mailto:<a href="mailto:michael.schildt@ifn-magdeburg.de" target="_blank">michael.schildt@ifn-magdeburg.de</a>&gt;&gt; wrote:<br>

<br>
    Hello!<br>
<br>
    I have some trouble with ResampleImageFilter. I&#39;m not sure if it is a<br>
    bug in my program or in ITK. It works as expected in Debug<br>
    configuration<br>
    but crashes in Release mode. I traced to the source line in ITK<br>
    were it<br>
    appear in RelWithDebMode configuration. I attach some information and<br>
    source.<br>
<br>
    Best reguards,<br>
      Michael Schildt<br>
<br>
    The method based on the example in the itkSoftwareGuide:<br>
    ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////<br>
    template &lt;typename dataT, typename TrfType&gt;<br>
    typename itk::Image&lt;dataT,3&gt;::Pointer transformImage(typename<br>
    itk::Image&lt;dataT,3&gt;::Pointer fixedIm,<br>
                                                       typename<br>
    itk::Image&lt;dataT,3&gt;::Pointer movingIm,<br>
                                                       typename<br>
    TrfType::ConstPointer transform)<br>
    {<br>
     typedef itk::Image&lt;dataT, 3&gt; ImageType;<br>
     typedef TrfType TransformType;<br>
     typedef itk::LinearInterpolateImageFunction&lt;ImageType,double&gt;<br>
    InterpolatorType;<br>
     //typedef itk::WindowedSincInterpolateImageFunction&lt;ImageType,3&gt;<br>
    InterpolatorType;<br>
<br>
     typedef itk::ResampleImageFilter&lt;ImageType,ImageType&gt; ResamplerType;<br>
<br>
     typename InterpolatorType::Pointer interpolator =<br>
    InterpolatorType::New();<br>
<br>
     typename ResamplerType::Pointer resampler = ResamplerType::New();<br>
     resampler-&gt;SetInput( movingIm );<br>
<br>
     resampler-&gt;SetTransform( transform.GetPointer() );<br>
     resampler-&gt;SetInterpolator( interpolator.GetPointer() );<br>
     resampler-&gt;SetSize( fixedIm-&gt;GetLargestPossibleRegion().GetSize() );<br>
     resampler-&gt;SetOutputOrigin( fixedIm-&gt;GetOrigin() );<br>
     resampler-&gt;SetOutputSpacing( fixedIm-&gt;GetSpacing() );<br>
     resampler-&gt;SetDefaultPixelValue( 0 );<br>
<br>
     // resample the moving image<br>
     resampler-&gt;Update();<br>
<br>
     return resampler-&gt;GetOutput();<br>
    }<br>
    ////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////<br>
<br>
    It appears in the Update() method of the resampler. While debugging it<br>
    seems like a problem in itkDataObject.cxx at line marked with **<br>
<br>
    //----------------------------------------------------------------------------<br>
    void<br>
    DataObject<br>
    ::UpdateOutputData()<br>
    {<br>
     // If we need to update due to PipelineMTime, or the fact that our<br>
     // data was released, then propagate the UpdateOutputData to the<br>
    source<br>
     // if there is one.<br>
     if ( m_UpdateMTime &lt; m_PipelineMTime || m_DataReleased ||<br>
         this-&gt;RequestedRegionIsOutsideOfTheBufferedRegion() )<br>
      {<br>
      if ( m_Source )      // &lt;---m_Source is valid, but when i step over<br>
    here, it gets an invalid pointer<br>
        {<br>
    **      m_Source-&gt;UpdateOutputData(this);      // maybe call on<br>
    invalid<br>
    pointer?<br>
        }<br>
      }<br>
    }<br>
    //----------------------------------------------------------------------------<br>
<br>
    Stack backtrace:<br>
<br>
       BrainGUI.exe!std::_Vector_iterator&lt;unsigned<br>
    int,std::allocator&lt;unsigned int&gt; &gt;::operator+()  + 0x17 Bytes    C++<br>
<br>
          BrainGUI.exe!std::vector&lt;unsigned<br>
        int,std::allocator&lt;unsigned int&gt; ::resize(unsigned int<br>
        _Newsize=1, unsigned int _Val=0)  Zeile 721 +<br>
<br>
    0x36 Bytes    C++<br>
       BrainGUI.exe!std::vector&lt;bool,std::allocator&lt;bool&gt;<br>
<br>
        ::_Insert_x(std::_Vb_const_iterator&lt;unsigned<br>
<br>
    int,int,std::vector&lt;bool,std::allocator&lt;bool&gt; &gt; &gt; _Where=..., unsigned<br>
    int _Count=1)  Zeile 2312    C++<br>
       BrainGUI.exe!std::vector&lt;bool,std::allocator&lt;bool&gt;<br>
<br>
        ::_Insert_n(std::_Vb_const_iterator&lt;unsigned<br>
<br>
    int,int,std::vector&lt;bool,std::allocator&lt;bool&gt; &gt; &gt; _Where=..., unsigned<br>
    int _Count=1, bool _Val=false)  Zeile 2292    C++<br>
       BrainGUI.exe!std::vector&lt;bool,std::allocator&lt;bool&gt;<br>
<br>
        ::resize(unsigned int _Newsize=1, bool _Val=false)  Zeile 2044<br>
           C++<br>
<br>
       BrainGUI.exe!itk::ProcessObject::CacheInputReleaseDataFlags()<br>
    Zeile 1055    C++<br>
       BrainGUI.exe!itk::ProcessObject::UpdateOutputData(itk::DataObject *<br>
    __formal=0x02079270)  Zeile 964    C++<br>
       BrainGUI.exe!itk::DataObject::UpdateOutputData()  Zeile 420 + 0x8<br>
    Bytes    C++<br>
<br>
    BrainGUI.exe!transformImage&lt;float,itk::VersorRigid3DTransform&lt;double&gt;<br>
<br>
        (itk::SmartPointer&lt;itk::Image&lt;float,3&gt; &gt; fixedIm={...},<br>
<br>
    itk::SmartPointer&lt;itk::Image&lt;float,3&gt; &gt; movingIm={...},<br>
    itk::SmartPointer&lt;itk::VersorRigid3DTransform&lt;double&gt; const &gt;<br>
    transform={...})  Zeile 123    C++<br>
       BrainGUI.exe!VersorRigid3DRegistration::doRegistration(bool<br>
    newVolume=true, void (const double &amp;, const<br>
    std::basic_string&lt;char,std::char_traits&lt;char&gt;,std::allocator&lt;char&gt;<br>
    &gt; &amp;)*<br>
    progress=0x004029e1)  Zeile 254 + 0x81 Bytes    C++<br>
       BrainGUI.exe!RegThread::Entry()  Zeile 273 + 0xf Bytes    C++<br>
       BrainGUI.exe!wxThreadInternal::DoThreadStart()  + 0xa5 Bytes    C++<br>
       BrainGUI.exe!wxThreadInternal::WinThreadStart()  + 0x6a Bytes       C++<br>
       BrainGUI.exe!_callthreadstartex()  Zeile 348 + 0x6 Bytes    C<br>
       BrainGUI.exe!_threadstartex(void * ptd=0x04ef1118)  Zeile 326 + 0x5<br>
    Bytes    C<br>
<br>
    System:<br>
      - Windows XP Prof SP3<br>
      - Visual Studio 2008 Express Edition<br>
      - ITK 3.14<br>
      - CMake 2.6.4<br>
      - wxWidgets 2.8.9 (wxPack Installer)<br>
      - Pentium 4<br>
      - 2GB Ram<br>
<br>
<br>
    _____________________________________<br></div></div>
    Powered by <a href="http://www.kitware.com" target="_blank">www.kitware.com</a> &lt;<a href="http://www.kitware.com" target="_blank">http://www.kitware.com</a>&gt;<div class="im"><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>
    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>
<br>
<br>
</div></blockquote><div><div></div><div class="h5">
<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>
Please keep messages on-topic and check the ITK FAQ at: <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>
</div></div></blockquote></div><br>