Rearrange data in two-dimensional array according to transformation from polar to Cartesian coordinates

user_3767239 0
Viewed 1716

I have a two-dimensional array that represents function values at positions in a polar coordinate system. For example:

import numpy as np

radius = np.linspace(0, 1, 50)
angle = np.linspace(0, 2*np.pi, radius.size)
r_grid, a_grid = np.meshgrid(radius, angle)
data = np.sqrt((r_grid/radius.max())**2
               + (a_grid/angle.max())**2)

Here the data is arranged in a rectangular grid corresponding to the polar coordinates. I want to rearrange the data in the array such that the axes represent the corresponding Cartesian coordinate system. The old versus new layout can be visualized as follows:

import matplotlib.pyplot as plt

fig, (ax1, ax2) = plt.subplots(ncols=2, figsize=plt.figaspect(0.5))
ax1.set(title='Polar coordinates', xlabel='Radius', ylabel='Angle')
ax1.pcolormesh(r_grid, a_grid, data)
ax2.set(title='Cartesian coordinates', xlabel='X', ylabel='Y')
x_grid = r_grid * np.cos(a_grid)
y_grid = r_grid * np.sin(a_grid)
ax2.pcolormesh(x_grid, y_grid, data)

Example

Here the coordinates are explicitly given and the plot is adjusted accordingly. I want the data to be rearranged in the data array itself instead. It should contain all values, optionally filling with zeros to fit the shape (similar to scipy.ndimage.rotate(..., reshape=True)).

If I manually loop over the polar arrays to compute the Cartesian coordinates, the result contains empty regions which ideally should be filled as well:

new = np.zeros_like(data)
visits = np.zeros_like(new)
for r, a, d in np.nditer((r_grid, a_grid, data)):
    i = 0.5 * (1 + r * np.sin(a)) * new.shape[0]
    j = 0.5 * (1 + r * np.cos(a)) * new.shape[1]
    i = min(int(i), new.shape[0] - 1)
    j = min(int(j), new.shape[1] - 1)
    new[i, j] += d
    visits[i, j] += 1
new /= np.maximum(visits, 1)
ax2.imshow(new, origin='lower')

Example attempt

Is there a way to achieve the transformation while avoiding empty regions in the resulting data array?

polar-coordinates coordinate-transformation scipy numpy python
3 Answers

tl;dr: No, not without changing some conditions of your problem.

The artefact you are seeing is a property of the transformation. It is not due to the fixed resolution in angle for all radii. Hence, it is not due to a wrong or bad implementation of the transformation. The Cartesian grid simply implies a higher special resolution at these areas as there are resolved points from the polar map.

  • The only "clean" way (that I can think of right now) to handle this is to have an adjustable resolution in the polar coordinates to account for the 1/r scaling. (If you input data allows it)

  • A somewhat cheating way of visualizing this without the gaps would to randomly distribute them over the gaps. The argument here is, that you do not have the resolution to decide in which bin they were to begin with. So you could just randomly throw them in one which might have been a possible origin and not throw them all in same one (as you are doing right now). However, I would like discourage this stronlgy. It just gives you a prettier plot. Note, that this is somewhat equivalent to the behaviour of the upper right plot in your question.

user_60493080

This doesn't really give the expected result, but maybe will help you in achieving a solution after some needed correction...


import numpy as np

radius = np.linspace(0, 1, 50)
angle = np.linspace(0, 2*np.pi, radius.size)
r_grid, a_grid = np.meshgrid(radius, angle)
data = np.sqrt((r_grid/radius.max())**2
               + (a_grid/angle.max())**2)


def polar_to_cartesian(data):
    new = np.zeros_like(data) * np.nan
    x = np.linspace(-1, 1, new.shape[1])
    y = np.linspace(-1, 1, new.shape[0])
    for i in range(new.shape[0]):
        for j in range(new.shape[1]):
            x0, y0 = x[j], y[i]
            r, a = np.sqrt(x0**2 + y0**2), np.arctan2(y0, x0)
            data_i = np.argmin(np.abs(a_grid[:, 0] - a))
            data_j = np.argmin(np.abs(r_grid[0, :] - r))
            val = data[data_i, data_j]

            if r <= 1:
                new[i, j] = val

    return new

new = polar_to_cartesian(data)
fig, ax = plt.subplots()
ax.imshow(new, origin='lower')

enter image description here

EDIT: Modified using np.arctan2 according to the suggestions of OP.

user_78129120

You could loop over the Cartesian array, transforming each grid point to polar coordinates and approximating the function value by interpolation from your polar grid data. You may still want to leave the corner regions blank though, for lack of close enough data.

I don't think there is a better way, unless of course you have access to the original function.

user_130141720
Related