[vtkusers] odd vtkImageData interpolation behavior

[David Feng]

I tracked down the problem in vtkImageData's FindCell.  There are 
actually two:

1 - Even when a dimension is of size 1 (the Y dimension in my case), 
there is still a voxel cell in that dimension.  The indexing code 
automatically sets all slice indexing to 0 when it multiplies the z 
voxel index by (dims[1]-1).  This is fixed with a simple check on the 
first two dimensions.  If they are of size one, set them to size two.  
Obviously it doesn't matter if the dimensions are used in order.

2 - The interpolation weights that a computed before returning from 
FromCell assume that the dimensions will be used in order.  The values 
need to be shifted if dimensions are to be ignored.  The shifting is a 
bit ugly. 

Both of these problems are fixed by inserting the following code after 
the InterpolateFunctions function call inside of FindCell, somewhere 
around line 769 of vtkimagedata.cxx:

 //
  // dimensions are used for cell indexing later.  if a dimension has 
only one
  // voxel, we need to make sure it still has one cell.  Also,
  // the coordinate weights need to be shifted to ignore that dimension.
  // This is ugly.  The order of these ifs makes a difference.
  //
 
  //
  // is y trivial?  if so, keep all even blocks of two weights
  //
  if (dims[1] == 1)
    {
    dims[1] = 2;
    weights[2] = weights[4];
    weights[3] = weights[5];
    }

  //
  // is x trivial? if so, keep all even indexed weights
  //
  if (dims[0] == 1)
    {
    dims[0] = 2;
    weights[1] = weights[2];
    weights[2] = weights[4];
    weights[3] = weights[6];
    }

I'll try to make a patch and submit it.

David

David Feng wrote:
> I was indeed missing something.
> I thought that VTK was using cell-centered data, which means that each 
> voxel contains a single data point at its center.  In fact VTK puts a 
> data point at every vertex of the voxel.  This means that there is one 
> fewer voxel in each dimension than there are data points, so the cell 
> indexing was in fact correct.
> I still don't know why my XZ image is not interpolating correctly, 
> though.
>
> David
>
> David Feng wrote:
>> I noticed some odd behavior when interpolating a vtkImageData 
>> object.  Rather than make a 2D image and rotate it to the XZ plane, I 
>> made it a 3D image with dimensions (x,1,y).  When I call FindCell 
>> simply using the world-space pixel coordinates, the interpolated 
>> values were not what I expected.  Taking a closer look at the 
>> function, I noticed the following:
>>
>> return loc[2] * (dims[0]-1)*(dims[1]-1) +
>>          loc[1] * (dims[0]-1) + loc[0];
>>
>> ...where loc is the interpolated image-space pixel coordinate and 
>> dims are the dimensions of the image.  The return value is supposed 
>> to the 1D array index within the image.  With this formula, a pixel 
>> location of (0,0,5) and (0,0,0) have the same 1D index for an image 
>> such as mine (XZ instead of XY, with dimensions (x,1,y)).   I haven't 
>> looked more closely, but indexing into a 3D should actually look like 
>> this, if I'm not mistaken:
>>
>> return loc[2] * dims[0] * dims[1] +
>>          loc[1] * dims[0] + loc[0];
>>
>> I'm not sure why the dimensions are being decremented in the original 
>> code.  It would seem to make the indexing off by several pixels for 
>> each row and slice.  When I made the change, I started getting the 
>> values I expected.  A simple way to test this is to use a 
>> vtkProbeFilter to probe an XZ image using itself as the input sample 
>> points -- without the change, my input is different from my output, 
>> which should not be the case.   Am I missing something?
>>
>> Thanks,
>> David
>>
>>
>>
>
>