Interp2 yields undesired results on center daigonal.

I am doing cross spectral matrices on some data where they look something like this:

Now when I interpolate using interp2 with 3 times as many points I get this:

[X,Y]   = meshgrid(18:33);
[Xq,Yq] = meshgrid(18:1/3:33);
imagesc(18:33,18:33,interp2(X,Y,coMatSimple,Xq,Yq))
set(gca,'YDir','Normal')
colorbar

While this is technically correct, the problem with this is that the diagonal from the bottom left to top right corners should be all 1's (by definition of what I'm doing) and it is for my original data, but when I interpolate it's only 1 at the original data points and not anywhere else on that diagonal.

Is there a way I can get around this problem, maybe another function I can use instead of interp2 or some parameter that I can change to force the center diagonal values to be one and then have the interpolation function take that into account when interpolating?

[Walter Roberson]

Could you perhaps just assign ones to the diagonal afterwards?

Answer 1

[Stephen23]

There is nothing wrong with interp2.

Here is a minimal working example of your situation:

>> interp2([0,1;1,0])
ans =
   0.00000   0.50000   1.00000
   0.50000   0.50000   0.50000
   1.00000   0.50000   0.00000
>>

The middle value is 0.5 (as it should be), because it interpolates between the zero values as well as the one values. There is no reason why the middle value should be one, because the ones do not have some magical higher priority that overrides the zeros, and so the original zeros also have an influence on the interpolated values.

Someone might offer better solutions, but here are two ideas to try:

1. rotate the data matrix by 45 degrees, interpolate, and then rotate back. You could use the image processing tools to do the rotation. Keep in mind that this will lose some data precision as well.
2. specify every output point along the antidiagonal and then use a scattered interpolant. This will be slightly more complex to set up (you will have to convert the data matrices into vectors and add the antidiagonal values), but has the advantages that there will be no loss of precision prior to interpolating, and you can specify the antidiagonal values exactly. Here is a demonstration of this:

>> M = [0,1;1,0] % input data
M =
     0     1
     1     0
>> [Ri,Ci] = ndgrid(0:1,0:1); % input locations
>> Rv = [Ri(:);0.5]; % add antidiagonal location
>> Cv = [Ci(:);0.5]; % add antidiagonal location
>> Mv = [ M(:);1];   % add antidiagonal data value
>> F = TriScatteredInterp(Rv,Cv,Mv);
>> [Ro,Co] = ndgrid(0:0.5:1,0:0.5:1); % output locations
>> F(Ro,Co)
ans =
            0          0.5            1
          0.5            1          0.5
            1          0.5            0
>>

Perfect!

[Jacob Ward]

Right, I realize that what interp2 is doing is correct (I'll edit my question to reflect that). I'm just saying for my specific situation I want all the center diagonal values to be 1's and I'm wondering if there's any way to accomplish that and have interp2 or whatever function anyone would recommend to take that into account.

[Jacob Ward]

I thought about maybe rotating it, looks like that may be the way I should go. Thanks for the help and quick response!

[Jacob Ward]

I'm intrigued by the second idea. I think I'll give both a shot!

[Stephen23]

@Jacob Ward: please make a comment and show an image with the interpolated data. It would be nice to see how it compares to the images in your original question.

Answer 2

[Cam Salzberger]

Hello Jacob,

interp2 uses a bilinear interpolation algorithm by default. Here is an image illustrating what that means:

To use your plot as an example, let's say we wanted to calculate the interpolated value at the corner of the two bottom left blue squares - around (18.5,18.5). MATLAB would first linearly interpolate the values at (18.5,18) and (18.5,19). This looks like they'll both be about 0.8ish. Then it'll linearly interpolate the value at (18.5,18.5), using the two previously calculated points. Since they're the same, the new value will also be about 0.8ish.

Were you to rotate your image such that the lines of similar value aligned with the axes of interpolation (basic x and y), then you might achieve more "accurate" results on the diagonal. However, you would be sacrificing "accuracy" in any areas where the trends already aligned with the original x and y axes. It's simply a limitation of the method.

You could try the 'cubic' or 'spline' interpolation methods, and see if they achieve more accurate results. I doubt you'll get perfect "1"s across the diagonal though, since they still depend on the x and y axes as their directions for interpolation.

If you could tell us what you are trying to do the interpolation for, we may be able to suggest alternatives. Otherwise, I hope this at least gives you a clearer picture of what's going on.

-Cam

[Jacob Ward]

Thanks for the explanation. I think I definitely have a better idea of what's going on now. I think the trends I'm looking for are actually all perpendicular or parallel to the diagonal so rotating my data by 45° before interpolating just might do the trick. Thanks for your help!