Hello Matplotlib Users.
I am trying to work my way out to interpolate a surface from polar coordinates (R,Theta,Z) to rectangular coordinates(X,Y,Z)
Basically I am looking for an equivalent for POLAR_SURFACE.pro routine in IDL using matplotlib/Scipy?
http://idlastro.gsfc.nasa.gov/idl_html_help/POLAR_SURFACE.html
Does anyone of you have done that before or could tell me what functions can I use to get that?
Regards
Bhargav Vaidya
(Quick Note: technically, you have cylindrical coordinates because you have Z, not Phi)
While I am sure there are some more optimal methods, the straightforward way is to do something like the following:
import numpy as np
from scipy.interpolate import griddata
orig_x = Rnp.cos(Theta)
orig_y = Rnp.sin(Theta)
orig_z = Z
grid_x, grid_y, grid_z = np.mgrid[orig_x.min():orig_x.max():10j,
orig_y.min():orig_y.max():10j,
orig_z.min():orig_z.max():10j]
new_vals = griddata((orig_x, orig_y, orig_z), orig_vals, (grid_x, grid_y, grid_z)
After running this, you should have new_vals which should have the same shape as grid_x. You can also specify the interpolation method used by griddata through the ‘method’ keyword argument. The ‘10j’ in the np.mgrid[] call merely specifies the number of points you want in that axis and can be changed to any other integer you want (just keep the ‘j’ there to make it work the way we want).
I hope that helps!
Ben Root
Hello Benjamin
Thanks for the reply… I was a little vague in my first question before…
Here is my specific problem
I have a 3D matrix say of density (rho [ r, theta , phi] )
where r, theta ,phi have (64,32,32) elements
theta and phi values are in radians.
I have to interpolate rho(r, theta) at some phi to rho(x,z) and rho(r, phi) at some theta to rho(x,y) …
However I used the a different code then the one you mentioned … I am using the scipy.ndimage.geometric_transform
I am attaching the figure obtained using a code below
THE PART OF THE CODE TO CONVERT THE POLAR TO CARTESIAN IS GIVEN BELOW
I have the following problems :

I am not happy with the interpolation. This is partly because I am using pcolormesh for plotting data as imshow is not working for my irregularly spaced grid. I have also tried higher order for the geometric transform but the difference is not that appreciable. Can I improve the interpolation is some way?

In the right panels of the plot above, I dont get the radial and vertical coordinate as the right cartesian values… I think this can be solved if I give the right x y and z values to pcolormesh … but I am not sure how to get the right values for them so that they are consistent with the transform.
import numpy as np
import scipy as S
import scipy.ndimage
def get_polar_plot(self,var,**kwargs):
def r2theta(outcoords, inputshape, origin):
xindex, yindex = outcoords
x0, y0 = origin
x = xindex  x0
y = yindex  y0
r = np.sqrt(x**2 + y**2)
theta = np.arctan2(np.piy, x)
theta_index = np.round((theta+ np.pi) * inputshape[1] / (1.015*np.pi))
return (r,theta_index)
def r2phi(outcoords, inputshape, origin):
xindex, yindex = outcoords
x0, y0 = origin
x = xindex  x0
y = yindex  y0
r = np.sqrt(x**2 + y**2)
phi = np.arctan2(y, x)
phi_index = np.round((phi + np.pi) * inputshape[1] / (2.1*np.pi))
return (r,phi_index)
x1 = kwargs.get('x1') # This is Radial coordinate
x2 = kwargs.get('x2') # This is Either Theta or Phi
data_arr2d = var
if kwargs.get('rtheta',False) == True:
var_cart = S.ndimage.geometric_transform(data_arr2d, r2theta,order=1,output_shape = (data_arr2d.shape[0]*2, data_arr2d.shape[0]),extra_keywords = {'inputshape':data_arr2d.shape,'origin':(data_arr2d.shape[0],0)})
if kwargs.get('rphi',False) == True:
var_cart = S.ndimage.geometric_transform(data_arr2d, r2phi, order=1,output_shape = (data_arr2d.shape[0]*2, data_arr2d.shape[0]*2),extra_keywords = {'inputshape':data_arr2d.shape,'origin':(data_arr2d.shape[0],data_arr2d.shape[0])})
return var_cart
Hope you will help me in this regard.
Regards
Bhargav
···
On Wed, Feb 23, 2011 at 2:09 PM, bhargav vaidya <coolastro@…287…> wrote: