Spherical to Cartesian function

[daviddoria]

Is there a library that has a function that will take (theta,phi) and return (x,y,z)?

Thanks,

Dave

[bvdet]

I don't know of one. This code requires a radial coordinate.

class Pt(object):
    def __init__(self, x=0.0, y=0.0, z=0.0):
        self.x = x
        self.y = y
        self.z = z
        
    def __str__(self):
        return '(%0.4f, %0.4f, %0.4f)' % (self.x, self.y, self.z)

    def __repr__(self):
        return 'Pt(%f, %f, %f)' % (self.x, self.y, self.z)

    def __add__(self, other):
        return Point(self.x+other.x, self.y+other.y, self.z+other.z)

    def __sub__(self, other):
        return Point(self.x-other.x, self.y-other.y, self.z-other.z)
    
    def __mul__(self, f):
        return Point(self.x*f, self.y*f, self.z*f)

    def dist(self, other):
        p = self-other
        return (p.x**2 + p.y**2 + p.z**2)**0.5

    def toSpherical(self):
        r = mag(self)
        theta = atan2(sqrt(self.x**2+self.y**2), self.z)
        phi = atan2(self.y, self.x)
        return SphericalPt(r, theta, phi)

class SphericalPt(object):
    def __init__(self, r, theta, phi):
        # radial coordinate, zenith angle, azimuth angle
        self.r = r
        self.theta = theta
        self.phi = phi

    def __str__(self):
        return '(%0.4f, %0.4f, %0.4f)' % (self.r, self.theta, self.phi)

    def __repr__(self):
        return 'SphericalPt(%f, %f, %f)' % (self.r, self.theta, self.phi)

    def toCartesian(self):
        x = self.r*cos(self.phi)*sin(self.theta)
        y = self.r*sin(self.phi)*sin(self.theta)
        z = self.r*cos(self.theta)
        return Pt(x,y,z)

[daviddoria]

Thanks, I haven't used classes in python yet so this is a good example for me.. For future readers, there is a small typo. In the Pt class, the return statements say Point(...) when they should say Pt(...).

It's interesting this isn't in a library anywhere though, it seems like deg2rad, cart2sphere, etc type functions would be pretty handy.

Dave

[bvdet]

You are welcome, daviddoria. I added this routine to the geometry library I have been developing for about three years.

-BV

[daviddoria]

Cool - is it available for download?

[bvdet]

Cool - is it available for download?

No. The work I have done interfaces with SDS/2 software. SDS/2 has a built-in interpreter, and we run scripts inside the SDS/2 3D model to automate tasks, extract and manipulate information, etc. Some of the code can be viewed here, but I have not updated it in several months.

What type of geometry work are you interested in? Following are some of my geometry applications:
LineLineIntersect3D (two point pairs)
Plane3D (three points)
three point circle
rotate a point about an arbitrary axis
DistancePointLine3D
LinePlaneIntersect3D
DistancePointPlane3D
PlanePlaneIntersect3D
BasisTransToGlobal - translate a point in a defined orthonormal basis to a point in the standard basis set
BasisTransToLocal - translate a point in the standard basis set to a point in a defined orthonormal basis
CircleCircleIntersect3D

I have no special knowledge nor am I an expert in geometry. Everything I have developed was from research done at various sites on the internet, particularly here.

[remiverhoeven]

The only problem with that code is that if you are transfering r,phi,theta to x,y,z you're out of luck with either of them is equal to 0.
This happends - for example - if you are trying to find x,y,z coordinates on a plane XY while Z=0.

How would you solve this problem?

[TrustyTony]

I do not get the problem, only I got not found error from mag function found in code, so I replaced it with distance from origin. I did some renaming to make the program more "PEP8 understandable"

from __future__ import division
from math import *

class Point(object):
    def __init__(self, x=0.0, y=0.0, z=0.0):
        self.x = x
        self.y = y
        self.z = z
        
    def __str__(self):
        return '(%0.4f, %0.4f, %0.4f)' % (self.x, self.y, self.z)

    def __repr__(self):
        return 'Point(%f, %f, %f)' % (self.x, self.y, self.z)

    def __add__(self, other):
        return Point(self.x + other.x, self.y + other.y, self.z + other.z)

    def __sub__(self, other):
        return Point(self.x-other.x, self.y-other.y, self.z-other.z)
    
    def __mul__(self, f):
        return Point(self.x*f, self.y*f, self.z*f)

    def dist(self, other):
        p = self - other
        return (p.x ** 2  +  p.y ** 2  +  p.z ** 2) ** 0.5

    def to_spherical(self):
        r = self.dist(Point(0, 0, 0))
        theta = atan2(hypot(self.x, self.y), self.z)
        phi = atan2(self.y, self.x)
        return SphericalPoint(r, theta, phi)

class SphericalPoint(object):
    def __init__(self, r, theta, phi):
        # radial coordinate, zenith angle, azimuth angle
        self.r = r
        self.theta = theta
        self.phi = phi

    def degrees(self):
        return 'SphericalPoint(%.4f, %.4f deg, %.4f deg)' % (self.r, degrees(self.theta) % 360, degrees(self.phi) % 360)
    
    def __str__(self):
        return '(%0.4f, %0.4f, %0.4f)' % (self.r, self.theta, self.phi)

    def __repr__(self):
        return 'SphericalPoint(%f, %f, %f)' % (self.r, self.theta, self.phi)

    def to_cartesian(self):
        x = self.r * cos(self.phi) * sin(self.theta)
        y = self.r * sin(self.phi) * sin(self.theta)
        z = self.r * cos(self.theta)
        return Point(x, y, z)

if __name__=='__main__':
    for p in (Point(x, y, 0).to_spherical() for x in range(-10, 10) for y in range(-10, 10)):
        print(p)

[remiverhoeven]

imagine any input where theta is 0, but phi is not. This point could be anywhere in the ZY plane. However, the formula: y=r*sin(phi)*sin(theta) would always return 0 (sin(0)=0, hence ...*0 = 0!). Therefore, it is not describing the points on the ZY axis.

[TrustyTony]

OK, this shows the problem then:

if __name__=='__main__':
    for p in (SphericalPoint(r, 0, phi) for r in range(1, 4) for phi in range(-2, 2)):
        print(p, p.to_cartesian(), p.to_cartesian().to_spherical())

Output:

(SphericalPoint(1.000000, 0.000000, -2.000000), Point(-0.000000, -0.000000, 1.000000), SphericalPoint(1.000000, 0.000000, -3.141593))
(SphericalPoint(1.000000, 0.000000, -1.000000), Point(0.000000, -0.000000, 1.000000), SphericalPoint(1.000000, 0.000000, -0.000000))
(SphericalPoint(1.000000, 0.000000, 0.000000), Point(0.000000, 0.000000, 1.000000), SphericalPoint(1.000000, 0.000000, 0.000000))
(SphericalPoint(1.000000, 0.000000, 1.000000), Point(0.000000, 0.000000, 1.000000), SphericalPoint(1.000000, 0.000000, 0.000000))
(SphericalPoint(2.000000, 0.000000, -2.000000), Point(-0.000000, -0.000000, 2.000000), SphericalPoint(2.000000, 0.000000, -3.141593))
(SphericalPoint(2.000000, 0.000000, -1.000000), Point(0.000000, -0.000000, 2.000000), SphericalPoint(2.000000, 0.000000, -0.000000))
(SphericalPoint(2.000000, 0.000000, 0.000000), Point(0.000000, 0.000000, 2.000000), SphericalPoint(2.000000, 0.000000, 0.000000))
(SphericalPoint(2.000000, 0.000000, 1.000000), Point(0.000000, 0.000000, 2.000000), SphericalPoint(2.000000, 0.000000, 0.000000))
(SphericalPoint(3.000000, 0.000000, -2.000000), Point(-0.000000, -0.000000, 3.000000), SphericalPoint(3.000000, 0.000000, -3.141593))
(SphericalPoint(3.000000, 0.000000, -1.000000), Point(0.000000, -0.000000, 3.000000), SphericalPoint(3.000000, 0.000000, -0.000000))
(SphericalPoint(3.000000, 0.000000, 0.000000), Point(0.000000, 0.000000, 3.000000), SphericalPoint(3.000000, 0.000000, 0.000000))
(SphericalPoint(3.000000, 0.000000, 1.000000), Point(0.000000, 0.000000, 3.000000), SphericalPoint(3.000000, 0.000000, 0.000000))

[bvdet]

We use these conversions to calculate the vertical slope angle and horiz slope angle of a vector in R3. Using these slopes as defaults, we recalculate new points in R3 using the slopes entered into a dialog box. So, we are starting with Point objects, converting to a SphericalPoint object, and back to Point objects. We should not encounter the problem you are describing. Am I wrong?

[remiverhoeven]

We use these conversions to calculate the vertical slope angle and horiz slope angle of a vector in R3. Using these slopes as defaults, we recalculate new points in R3 using the slopes entered into a dialog box. So, we are starting with Point objects, converting to a SphericalPoint object, and back to Point objects. We should not encounter the problem you are describing. Am I wrong?

Actually it might. The problem with this function is the calculation of the spherical coordinates is well defined. Atan2 is contiuous between -pi/2 and +pi/2 so will not cause any particular problems. When both x and y are 0 however, the atan is not defined; but since you do not need the theta or phi for that calculation (you purely need the R); thats okay too.

When you convert a theta,phi,r base back to x,y,z you still run into the problem of having a multiplication error between r*sin*sin and r*cos*sin that causes unwanted 0's in your outcome. Except if your vectors are limited in some way from obtaining unwanted (0-producing) values.

[Ratan_2]

There is an astropy where direct interconversion is possible

[Ratan_2]

There is one problem I am facing what should be the range of radial coordinate given radius of earth is 6378.14 km. I want to convert spherical cordinates ot cartesian in relation to earth