### Tags: algorithm, amimplementing, complex, following, interesting, interp2, interpolation, matlab, numbers, programming, requires, running, stepsfimage

# Interpolation over complex numbers

On Programmer » Matlab

10,168 words with 9 Comments; publish: Tue, 29 Apr 2008 14:55:00 GMT; (20046.88, « »)

Hi,

I'm running into an interesting problem with interp2. I am

implementing an algorithm that requires the following steps:

fimage = FFT2(image)

POLARTRANSFORM(fimage)

...

CARTESIANTRANSFORM(fimage)

image2 = IFFT2(fimage)

What I've discovered is the process of polar/cartesian transforming

introduces serious errors into the phase component of the FFT'd image

but NOT into the magnitude component. I do the polar/cartesian polar

transforms using interp2 as follows:

function [pim, radius, theta] = c2p(im, N, rlower, rupper)

[rows, cols] = size(im);

cx = cols/2.0 + 0.5; % Add 0.5 because index starts at 1

cy = rows/2.0 + 0.5;

[theta, radius] = meshgrid(linspace(0, 2*pi, N),

logspace(log10(rlower), log10(rupper), N));

xi = radius.*cos(theta) + cx; % Locations in image to interpolate

data

yi = radius.*sin(theta) + cy; % from.

[x,y] = meshgrid([1:cols],[1:rows]);

% Older versions of interp2 may not handle complex numbers correctly!

% Need to interpolate the magnitude and phase separately to avoid

problems.

pim = interp2(x, y, abs(im), xi, yi, 'cubic').*exp(i*interp2(x, y,

angle(im), xi, yi, 'cubic'));

mask = isnan(pim);

pim(mask) = 0;

Is there a reason I am getting junk of the phase interpolation but

not for the magnitude?

*http://matlab.todaysummary.com/q_matlab_32000.html*

All Comments

Leave a comment...

- 9 Comments
- "Elvis Dieguez" <edieguez.matlab.todaysummary.com.ieeeREMOVEME.org> wrote in message
news:ef4cfe8.-1.matlab.todaysummary.com.webcrossing.raydaftYaTP...

> Hi,

> I'm running into an interesting problem with interp2. I am

> implementing an algorithm that requires the following steps:

> fimage = FFT2(image)

> POLARTRANSFORM(fimage)

> ...

> CARTESIANTRANSFORM(fimage)

> image2 = IFFT2(fimage)

> What I've discovered is the process of polar/cartesian transforming

> introduces serious errors into the phase component of the FFT'd image

> but NOT into the magnitude component. I do the polar/cartesian polar

> transforms using interp2 as follows:

> function [pim, radius, theta] = c2p(im, N, rlower, rupper)

> [rows, cols] = size(im);

> cx = cols/2.0 + 0.5; % Add 0.5 because index starts at 1

> cy = rows/2.0 + 0.5;

> [theta, radius] = meshgrid(linspace(0, 2*pi, N),

> logspace(log10(rlower), log10(rupper), N));

> xi = radius.*cos(theta) + cx; % Locations in image to interpolate

> data

> yi = radius.*sin(theta) + cy; % from.

> [x,y] = meshgrid([1:cols],[1:rows]);

> % Older versions of interp2 may not handle complex numbers correctly!

> % Need to interpolate the magnitude and phase separately to avoid

> problems.

> pim = interp2(x, y, abs(im), xi, yi, 'cubic').*exp(i*interp2(x, y,

> angle(im), xi, yi, 'cubic'));

> mask = isnan(pim);

> pim(mask) = 0;

> Is there a reason I am getting junk of the phase interpolation but

> not for the magnitude?

Without thinking too hard or looking at your code, my guess is a problem

with phase unwrapping.

#1; Tue, 29 Apr 2008 14:56:00 GMT

- "Elvis Dieguez" <edieguez.matlab.todaysummary.com.ieeeREMOVEME.org> wrote in message
- Elvis Dieguez wrote:
> I'm running into an interesting problem with interp2. I am

> implementing an algorithm that requires the following steps:

> fimage = FFT2(image)

> POLARTRANSFORM(fimage)

> ...

> CARTESIANTRANSFORM(fimage)

> image2 = IFFT2(fimage)

> What I've discovered is the process of polar/cartesian transforming

> introduces serious errors into the phase component of the FFT'd

> image

> but NOT into the magnitude component. I do the polar/cartesian

> polar

> transforms using interp2 as follows:

As Ken pointed out you run into a phase unwrapping problem. Phase

information has a discontinuity at 2*pi which is not handled by

imresize. For example: if two neighboring pixels have the angles 0

and 2*pi imresize may use the average value for the output pixel.

This obviously introduces a huge error. Instead of using magnitude

and phase images as input to imresize you could use real and

imaginary part images. These images can be resized savely.

-Herbert

#2; Tue, 29 Apr 2008 14:57:00 GMT

- Elvis Dieguez wrote:
- Herbert Ramoser wrote:
>

> Elvis Dieguez wrote:

> transforming

> As Ken pointed out you run into a phase unwrapping problem. Phase

> information has a discontinuity at 2*pi which is not handled by

> imresize. For example: if two neighboring pixels have the angles 0

> and 2*pi imresize may use the average value for the output pixel.

> This obviously introduces a huge error. Instead of using magnitude

> and phase images as input to imresize you could use real and

> imaginary part images. These images can be resized savely.

Don't know much about 2D, but in 1D you get unwanted effects if you

interpolate the real and imaginary parts rather than the magnitude

and phase.

#3; Tue, 29 Apr 2008 14:58:00 GMT

- Herbert Ramoser wrote:
- On Feb 7, 4:07 am, "Steve Amphlett"
<Firstname.Lastname.matlab.todaysummary.com.where_I_work.com> wrote:

> Don't know much about 2D, but in 1D you get unwanted effects if you

> interpolate the real and imaginary parts rather than the magnitude

> and phase.

>

Really?

I don't believe it.

But if you're right, then all the work by hundreds of researchers on

developing global tide models using data from TOPEX/Poseidon

oceanographic satellite is down the drain.

Can you demonstrate a single case where interpolation on real and

imaginary parts gives "unwanted" effects.

Seems to me interpolation on real and imag parts is the ONLY valid way

to interpolate in 1D or 2D without getting screwed up with phase, as

mentioned by others.

#4; Tue, 29 Apr 2008 14:59:00 GMT

- On Feb 7, 4:07 am, "Steve Amphlett"
- NZTideMan wrote:
>

> On Feb 7, 4:07 am, "Steve Amphlett"

> <Firstname.Lastname.matlab.todaysummary.com.where_I_work.com> wrote:

> you

> magnitude

> Really?

> I don't believe it.

> But if you're right, then all the work by hundreds of researchers

> on

> developing global tide models using data from TOPEX/Poseidon

> oceanographic satellite is down the drain.

> Can you demonstrate a single case where interpolation on real and

> imaginary parts gives "unwanted" effects.

> Seems to me interpolation on real and imag parts is the ONLY valid

> way

> to interpolate in 1D or 2D without getting screwed up with phase,

> as

> mentioned by others.

Ok, let's say you have two sine waves at very similar frequencies and

amplitudes but significantly different phases and want to smoothly

move from one to the other. If you interpolate the real and

imaginary parts independently the amplitude will dip as you move from

one point to the other, whereas you'd probably want the amplitude and

phase to gradually adjust each in a monotonic way. What you really

want in the complex plane is a spiral from one point to the other,

not a straight line.

#5; Tue, 29 Apr 2008 15:00:00 GMT

- NZTideMan wrote:
- Steve Amphlett wrote:
>

> Ok, let's say you have two sine waves at very similar frequencies

> and

> amplitudes but significantly different phases and want to smoothly

> move from one to the other. If you interpolate the real and

> imaginary parts independently the amplitude will dip as you move

> from

> one point to the other, whereas you'd probably want the amplitude

> and

> phase to gradually adjust each in a monotonic way. What you really

> want in the complex plane is a spiral from one point to the other,

> not a straight line.

... I'm not advocating averaging or any other DSP calcs in general

using amplitude and phase. Just highlighting a case where unusual

interpolation requirements may lean toward it. As in other areas,

one curve that's mathematically in the middle of two others may not

be what's desired.

#6; Tue, 29 Apr 2008 15:01:00 GMT

- Steve Amphlett wrote:
- On Feb 7, 9:48 am, "Steve Amphlett" <m....matlab.todaysummary.com.home.now> wrote:
> Steve Amphlett wrote:

>

> ... I'm not advocating averaging or any other DSP calcs in general

> using amplitude and phase. Just highlighting a case where unusual

> interpolation requirements may lean toward it. As in other areas,

> one curve that's mathematically in the middle of two others may not

> be what's desired.

In your (pathological) case, I'd interpolate on real and imag on one

signal, then the other and superpose.

#7; Tue, 29 Apr 2008 15:02:00 GMT

- On Feb 7, 9:48 am, "Steve Amphlett" <m....matlab.todaysummary.com.home.now> wrote:
- NZTideMan wrote:
>

> In your (pathological) case, I'd interpolate on real and imag on

> one

> signal, then the other and superpose.

Pathological? It comes up all the time in the simulation world. If

you have noise predictions at a set of finite (engine) speeds and

want to synthesise a slow transient through them, this is exactly

what you have.

#8; Tue, 29 Apr 2008 15:03:00 GMT

- NZTideMan wrote:
- Herbert Ramoser wrote:
I've tried doing the interpolation over real and imaginary parts and

the errors are different but just as bad. I believe you are correct

in diagnosing the problem (i.e. the phase discontinuity). I've tried

using 'unwrap' but that doesn't solve the problem. Could it be that

'unwrap' works column or row wise thereby still leaving

discontinuities on the 2D surface?

> As Ken pointed out you run into a phase unwrapping problem. Phase

> information has a discontinuity at 2*pi which is not handled by

> imresize. For example: if two neighboring pixels have the angles 0

> and 2*pi imresize may use the average value for the output pixel.

> This obviously introduces a huge error. Instead of using magnitude

> and phase images as input to imresize you could use real and

> imaginary part images. These images can be resized savely.

> -Herbert

#9; Tue, 29 Apr 2008 15:04:00 GMT

- Herbert Ramoser wrote: