vtkImageReslice Output Appears to be "Shifted Up and to the Left" by one Pixel

Support
#1

Hello,
I am using vtkImageReslice to resample a simple image. The output image is not as expected, though. It appears to be “shifted up and to the left”. Here is what I mean:

Let’s say the input to vtkImageReslice is a 5x5x1 image like this:
0 1 2 3 4
10 11 12 13 14
20 21 22 23 24
30 31 32 33 34
And Lets say that this input image spacing is 2x2x2, and the origin is (0, 0, 0).

If I set the reslice axes to be the identity matrix, set output extent to be 10x10x1, set the output origin to be (0, 0, 0), set output spacing to be 1x1x1, and set interpolation mode to be nearest neighbor, I expect the output image to look like this:
0 0 1 1 2 2 3 3 4 4
0 0 1 1 2 2 3 3 4 4
10 10 11 11 12 12 13 13 14 14
10 10 11 11 12 12 13 13 14 14
20 20 21 21 22 22 23 23 24 24
20 20 21 21 22 22 23 23 24 24
etc…

But the output looks like this:
0 1 1 2 2 3 3 4 4 4
10 11 11 12 12 13 13 14 14
10 11 11 12 12 13 13 14 14 14
20 21 21 22 22 23 23 24 24 24
20 21 21 22 22 23 23 24 24 24
etc…

Notice that it appears as though the the top row and left-most column were truncated during reslicing or interpolation. If this is a question whose answer can be found easily by searching, please forgive me. I tried, but perhaps I do not know the correct search terms to use.

I would be glad to post my test code if needed.

Thank you.

#2

@dgobbi are you available to help me with this? I hope I’ve worded my question correctly. Thank you!

#3

HI Michelle,

This has to with the way that pixel coordinates are defined.

For your input image,

  • point (0.0, 0.0, 0.0) has a value of 0
  • point (2.0, 0.0, 0.0) has a value of 1
  • etc

For the output image,

  • point (0.0, 0.0, 0.0) has a value of 0
  • point (1.0, 0.0, 0.0) is halfway and rounds up (as per usual rounding behavior)
  • point (2.0, 0.0, 0.0) has a value of 1

Let’s look at this another way. Consider the output for linear interpolation:

  • point (0.0, 0.0, 0.0) has a value of 0
  • point (1.0, 0.0, 0.0) has a value of 0.5
  • point (2.0, 0.0, 0.0) has a value of 1

The mistake that you made is that you assumed that the origin should be the same for the input and the output. Remember that the coords of an image pixel denote the position of the center of the pixel. So what happens when you divide a pixel in four? Consider the upper-left pixel in the image:

+-------+          +---+---+      
|       |          | + |   |      
|   +   |          |---+---|      
|       |          |   |   |      
+-------+          +---+---+

See? The origin in fact must change if you want the 2x2 block of subpixels to exactly match the bounds of the original pixel.

Edit: I should show what happens when the origin is kept the same. In this case, the new pixels positions are midway between the original positions. This is typical for data interpolation.

+-------+-------+                           
|       |       |            +---+---+---+      
|   +   |       |            | + |   |   |      
|       |       |            |---+---+---+      
+-------+-------+            |   |   |   |      
|       |       |            |---+---+---+      
|       |       |            |   |   |   |      
|       |       |            +---+---+---+            
+-------+-------+
#4

Thank you for helping me, David. Your explanation is very clear, and I appreciate the illustrations. Now instead of:
reslicer->SetOutputOrigin(0,0,0)
I do this:
reslicer->SetOutputOrigin(pixelWidth/2, pixelHeight/2, 0)
and it works.

But! It works as long as I leave the input image origin at (0,0, 0).
When I initially create the input image (vtkImageData), I leave it’s origin as the default (0, 0, 0). Then when I reslice, I set the ouput image’s origin as the center of the first pixel, as described above. And this works. But if I try to set the origin of the input image as the center of it’s first pixel as well, the output is not as expected. Why is this? Is it because vtkImageData origin is in world coordinates and vtkImageReslice outputOrigin is in image coordinates?

  • Michelle
How split image data with plane form DICOM images
#5

Both vtkImageData and vtkImageReslice use the same coords, which are VTK data coords (and data coords are the same as world coords if the Prop3D transform is identity, which is the default).

Going back to the figure, the point that you want to be invariant when you resample the data is the corner, p:

p-------+          p---+---+      
|       |          | o |   |      
|   o   |          |---+---|      
|       |          |   |   |      
+-------+          +---+---+

The corner p is displaced from the origin by half the spacing (in fact, by negative half the spacing). So to do your nearest-neighbor resampling, the relationship between the input origin and the output origin is as follows since you want p to be invariant:

input_origin - 0.5*input_spacing = output_origin - 0.5*output_spacing

So let’s say that the input spacing is 2 and the output spacing is 1. Then if your input origin is (0,0,0), the corner point is in fact at (-1,-1,0) and:

output_origin = input_origin - 0.5*input_spacing + 0.5*output_spacing
              = 0.0 - 1.0 + 0.5
              = -0.5

Let’s say you instead want p to be zero. Then, for input spacing 2 and output spacing 1, you would have:

input_origin = 1.0
output_origin = 0.5

So SetOutputOrigin(pixelWidth/2, pixelHeight/2, 0) should actually give incorrect results when the input origin is (0,0,0).

#6

Thank you. This is clear.