[Insight-users] Attempt to rotate image resulting in black image.

[Luis Ibanez]

Hi Michael,


We define as "Origin" of an image the physical coordinates of its
pixel with index [0,0,0].

Rotations, are performed by default, with respect to the origin
of the coordinate system (the point with physical coordinates
(0.0, 0.0, 0.0), which doesn't have to coincide with the first
pixel of the image, nor with the central pixel of the image.

In ITK we assume that your image is a representation of some
physical reality (for example a patient's body), and therefore
all transformations should be performed in the context of the
physical coordinates, not the grid of pixels.

The purpose of the "Origin" is that if you use your physical
image acquisition device (microscope, telescope, CT scanner,
MRI scanner, ultrasound device, lidar...) and you acquire
images of different sections of an object that has a physical
manifestation in the real world, then you should be able to
recreate a mosaic of such real object by simply placing the
images in a common coordinate system according to the coordinates
of their origin.


----


In order to understand the behavior of the transformations,
you should think of your image in the context of the physical
coordinate system.

If you want your image to rotate around the one of the pixel
in one of the image corners the you should do the following:


   A) Compute the Index of that pixel

   B) Call image->TransformIndexToPhysicalPoint( index, point )
      to compute the physical coordinates corresponding to this
      pixel. This transformation will take into account the
      image origin, pixel spacing and image orientation.

   C) Call transform->SetCenter( point )
      to tell the transform that you want to use this point
      as the center of rotation.

   D) Connect the transform to the resample filter.

   E) call Update() in the resample filter.



Please let us know if you still find any problem,


     Thanks,


        Luis



---------------------------
Michael Jackson wrote:
> Just to follow up with this, If I force my origins to be (0,0) then  the 
> rotation seems to work. So I guess I am not understanding what  exactly 
> the origins are? I have reread through several sections of the  ITK 
> guide and I thought I had it figured out. The origins of where the  
> image is taken is actually stored in our data files so I just thought  I 
> would use that. The values are in Microns and are typically around  the 
> 48,000 to 50,000 range.
> 
> Is there an implicit units (like mm) in vtkImage that I don't know  
> about? I just can not get this figured out.
> 
>  Explanations are truly appreciated at this point.
> 
> _________________________________________________________
> Mike Jackson                  mike.jackson at bluequartz.net
> BlueQuartz Software                    www.bluequartz.net
> Principal Software Engineer                  Dayton, Ohio
> 
> 
> 
> On Jun 3, 2009, at 5:50 PM, Michael Jackson wrote:
> 
>> I am trying to rotate an image about the "upper left" corner. I had  
>> to write my own ImageSource derived class to read the data from and  
>> HDF5 file. I am pretty sure the "Hdf5 reader" works as I can create  
>> the ImageReader instance from it, read an image and write the image  
>> to disk as a tiff file.
>>
>>  Where my code seems to be failing is in the application of an  
>> AffineTransformFilter. I am basically taking the code straight from  
>> one of the examples adapting it to use my ImageSource derived reader  
>> and running it. And I get the dreaded "black image". Here is the  code 
>> and the output from the program. Oddly, if I switch the reader  back 
>> to the usual ImageFileReader and read up a tiff image I get the  
>> correct result, so maybe my custom reader class is not quite correct?
>>
>>
>> #define USE_MXA_FILE 1
>> int main(int argc, char **argv) {
>>  int exampleAction =0;
>>  const     unsigned int   Dimension = 2;
>>  typedef   unsigned char  InputPixelType;
>>  typedef   unsigned char  OutputPixelType;
>>
>>  typedef itk::Image< InputPixelType,  Dimension >   InputImageType;
>>  typedef itk::Image< OutputPixelType, Dimension >   OutputImageType;
>>
>> #if USE_MXA_FILE
>>  typedef itkR3DImageIO  ReaderType;
>> #else
>>  typedef itk::ImageFileReader< InputImageType >  ReaderType;
>> #endif
>>  typedef itk::ImageFileWriter< OutputImageType >  WriterType;
>>
>>  ReaderType::Pointer reader = ReaderType::New();
>>  WriterType::Pointer writer = WriterType::New();
>>
>>  reader->SetFileName( argv[1] );
>> #if USE_MXA_FILE
>>  reader->SetDatasetPath(argv[2]);
>> #endif
>>  reader->Update();
>>  writer->SetFileName( argv[3]);
>>  typedef itk::ResampleImageFilter<InputImageType, OutputImageType >   
>> FilterType;
>>  FilterType::Pointer filter = FilterType::New();
>>  typedef itk::AffineTransform< double, Dimension >  TransformType;
>>  TransformType::Pointer transform = TransformType::New();
>>
>>  typedef  itk::NearestNeighborInterpolateImageFunction<InputImageType, 
>> double  >  InterpolatorType;
>>  InterpolatorType::Pointer interpolator = InterpolatorType::New();
>>  filter->SetInterpolator( interpolator );
>>  filter->SetDefaultPixelValue( 0 );
>>  const InputImageType::SpacingType& spacing = reader->GetOutput()- 
>> >GetSpacing();
>>  const InputImageType::PointType& origin = reader->GetOutput()- 
>> >GetOrigin();
>>  InputImageType::SizeType size = reader->GetOutput()- 
>> >GetLargestPossibleRegion().GetSize();
>>  filter->SetOutputOrigin( origin );
>>  filter->SetOutputSpacing( spacing );
>>  filter->SetSize( size );
>>
>>  filter->SetInput( reader->GetOutput() );
>>  writer->SetInput( filter->GetOutput() );
>>
>>  std::cout << "Input Size: " << size[0] << "," << size[1] <<  std::endl;
>>  std::cout << "Input Origin: " << origin[0] << "," << origin[1] <<  
>> std::endl;
>>  std::cout << "Input Spacing: " << spacing[0] << "," << spacing[1]  << 
>> std::endl;
>>
>>  TransformType::OutputVectorType translation1;
>>  translation1[0] =   -origin[0];
>>  translation1[1] =   -origin[1];
>>  transform->Translate( translation1 );
>>
>>  const double degreesToRadians = atan(1.0) / 45.0;
>>  transform->Rotate2D( -10.0 * degreesToRadians, false );
>>
>>  TransformType::OutputVectorType translation2;
>>  translation2[0] =   origin[0];
>>  translation2[1] =   origin[1];
>>  transform->Translate( translation2, false );
>>  filter->SetTransform( transform );
>>  try
>>    {
>>    writer->Update();
>>    }
>>  catch( itk::ExceptionObject & excep )
>>    {
>>    std::cerr << "Exception caught !" << std::endl;
>>    std::cerr << excep << std::endl;
>>    }
>>
>>  return EXIT_SUCCESS;
>> }
>>
>>
>> output:
>> GenerateData()
>> GenerateData()
>> Input Size: 1292,968
>> Input Origin: 46978.9,48347.8
>> Input Spacing: 0.207987,0.207987
>> GenerateData()
>>
>> Any help would be great.
>> _________________________________________________________
>> Mike Jackson                  mike.jackson at bluequartz.net
>> BlueQuartz Software                    www.bluequartz.net
>> Principal Software Engineer                  Dayton, Ohio
>>
>>
>>
> 
> _____________________________________
> Powered by www.kitware.com
> 
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.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
>