While doing some tests, I discovered something crazy.
If you use the usual definitions for coordinates transformations (like these), take a single point of coordinates (1, 1, 1). Now transform it in spherical coordinates. And then back to Cartesian. Everything ok.
Now take (-1, 1, 1) or (1, -1, 1) or (-1, -1, 1). And have fun discovering that, depending on the choice for the azimuthal angle one of the coordinates (either x or y) will return with reversed sign.
Try it.
It took me some time to (re)discover that I need to change the phi definition according to the value of x and y.
Here's a simple python code on how it should be done:
import numpy as np
import sys
try:
x= float(sys.argv[1]); y=float(sys.argv[2]) ; z=float(sys.argv[3])
except:
print ("Error\nUsage:", sys.argv[0], "x y z"); sys.exit(1)
print ("you inserted\nx\ty\tz\n%f\t%f\t%f\n" %(x,y,z))
rho = np.sqrt((x*x)+(y*y)+(z*z))
theta= np.arccos(z/rho)
phi = np.arctan(y/x)
if y<0 and x>=0:
phi= np.arcsin(y/np.sqrt(x*x+y*y))
elif x<0 and y>=0:
phi=np.arccos(x/(np.sqrt(x*x+y*y)))
elif x<0 and y<0:
phi+=np.pi
print ("rho\ttheta\tphi\n%f\t%f\t%f\n" %(rho, theta, phi))
x1= rho*np.sin(theta)*np.cos(phi)
y1= rho*np.sin(theta)*np.sin(phi)
z1= rho*np.cos(theta)
print ("x1\ty1\tz1\n%f\t%f\t%f\n" %(x1, y1, z1))
The cool part is to add pi to the arctan definition when both x and y are negative, and to use the definitions with arcsin and arccos when only y or only x is negative. If you do like that, things never loose sign.
I don't know if it is me or what, but I find irritating that nobody ever told me that if you go like Cartesian->Spherical->Cartesian you have to be so careful.
Well, now it works.
