[Insight-users] Image index limits

[Bill Lorensen]

This code is incorrect:
 for (int i=0; i<10; i++)
  {
    pixelIndex[0] +=   pixelIndex[0];
    ImageType::PixelType   pixelValue = image->GetPixel( pixelIndex );
    std::cout << "pixelValue = " << pixelValue << std::endl;
  }

pixelIndex[0] rapidly gets out of range.


On Tue, Jan 5, 2010 at 4:16 PM, Oleksandr Dzyubak <adzyubak at gmail.com> wrote:
> Hi Luis,
>
> The problem of empty slices when processing large images
> still bothers me and this time I decided to make direct tests.
> I wrote a tiny program which allocates and fills a buffer (1340x1300x720)
> with
> constant values I=100. A code snippet is attached below.
>
> 1) This program does not read an external file
> so we do not use my "corrupted" images anymore plus
> we exclude influence of the ImageRead class.
>
> 2) Next, I directly read intensities from a buffer before and after
> a "magic boundary" (slice #104) and print the values to a console thus
> we can see whether or not this failure is caused by the allocation
> procedure.
>
> 3) In addition to 2), I read the buffer content and write it to a file and
> measure
> the intensities again to see whether the ImageWrite class causes a problem.
>
> 4) Finally, I did massive tests for all possible pixel types.
>
> Believe it or not, but it works perfectly well and produces
> correct intensities over the entire volume 1340x1300x720.
>
> So I conclude that at least this part works correctly.
>
> Next step. Now I use the example from the ITK library
>
> InsightToolkit-3.16.0/Examples/IO/ImageReadRegionOfInterestWrite.cxx
>
> and I try to extract couple slices before and after the "magic boundary".
>
> Wow! After the "magic boundary", I get all the slices with zeros.
>
> What bothers me the most is that I actually need the classes which use
> the Gaussians and Gaussian derivatives. Surprisingly those algorithms
> have the same problem: after the "magic boundary", all the slices
> they produce are with zeros.
>
> As to me, it looks like they all use the same technique to fly over indices
> which results in a sort of "partial volume effect"
> discarding intensities with large index values.
>
> Unless I solve this problem, I cannot go any further with
> using Gaussians and Gaussian derivatives for processing large images.
>
> Luis, I would appreciate any suggestion or any crazy hint.
>
> Regards,
>
> Alex
>
>
> ******** BeginCodeSnippet *********
>
> #include "itkImage.h"
> #include "itkImageFileWriter.h"
>
> int main(int argc, char * argv[])
> {
>  if( argc < 2 )
>   {
>   std::cerr << "Usage: " << std::endl;
>   std::cerr << argv[0] << "  outputImageFile" << std::endl;
>   return 1;
>   }
>
> //  typedef unsigned char   PixelType;
> //  typedef unsigned int    PixelType;
> //  typedef short          PixelType;
>  typedef float          PixelType;
>
>  const unsigned int      Dimension = 3;
>  typedef itk::Image< PixelType, Dimension > ImageType;
>
>  ImageType::Pointer image = ImageType::New();
>
>  // The image region should be initialized
>  ImageType::IndexType start;
>  ImageType::SizeType  size;
>
>  start[0] =   0;  // first index on X
>  start[1] =   0;  // first index on Y
>  start[2] =   0;  // first index on Z
>
>  size[0]  = 1340;  // size along X
>  size[1]  = 1300;  // size along Y
>  size[2]  =  720;  // size along Z
>
>  ImageType::RegionType region;
>  region.SetSize( size );
>  region.SetIndex( start );
>
>  // Pixel data is allocated
>  image->SetRegions( region );
>  image->Allocate();
>
>  // The image buffer is initialized to a particular value
>  ImageType::PixelType  initialValue = 100;
>  image->FillBuffer( initialValue );
>
>  ImageType::IndexType pixelIndex;
>
>  pixelIndex[0] = 670;   // x position
>  pixelIndex[1] = 650;   // y position
>  pixelIndex[2] =  10;   // z position
>
>  for (int i=0; i<10; i++)
>   {
>     pixelIndex[0] +=   pixelIndex[0];
>     ImageType::PixelType   pixelValue = image->GetPixel( pixelIndex );
>     std::cout << "pixelValue = " << pixelValue << std::endl;
>   }
>
>  pixelIndex[0] = 670;   // x position
>  pixelIndex[1] = 650;   // y position
>  pixelIndex[2] = 360;   // z position
>
>  for (int i=0; i<10; i++)
>   {
>     pixelIndex[0] +=   pixelIndex[0];
>     ImageType::PixelType   pixelValue = image->GetPixel( pixelIndex );
>     std::cout << "pixelValue = " << pixelValue << std::endl;
>   }
>
> //  image->SetPixel(   pixelIndex,   pixelValue+1  );
>
>  typedef itk::ImageFileWriter< ImageType > WriterType;
>  WriterType::Pointer writer = WriterType::New();
>
>  writer->SetFileName( argv[1] );
>  writer->SetInput( image );
>
>  try
>   {
>   writer->Update();
>   }
>  catch( itk::ExceptionObject & exp )
>   {
>   std::cerr << "Exception caught !" << std::endl;
>   std::cerr << exp << std::endl;
>   }
>
>  return 0;
> }
>
> ******** EndCodeSnippet ********
>
>
> Oleksandr Dzyubak wrote:
>>
>> Hi Luis,
>>
>> Yes, indeed.
>>
>> This time I tried both pixel types, "float" and "signed short"
>> and compared. The results are the same. It goes well up to
>> the slice # 104. From the slice # 105 and on, it's all zeros.
>>
>> It is some sort of "saturation" or something specific to my machine?
>> Did you try on your machine? What's your result?
>>
>> Regards,
>>
>> Alex
>>
>> Luis Ibanez wrote:
>>>
>>> Hi Alex,
>>>
>>>  Did you changed to "float" the pixel type of
>>>
>>>      ImageReadRegionOfInterestWrite.cxx    ?
>>>
>>>    by default the pixel type is "signed short".
>>>
>>>   and your image is of pixel type "float".
>>>
>>>   ITK will perform silent C-Style casting
>>>   at the moment of reading the image.
>>>
>>>
>>> Please let us know,
>>>
>>>     Thanks,
>>>
>>>          Luis
>>>
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.html
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://www.itk.org/mailman/listinfo/insight-users
>