Spherical to Cartesian function

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)

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

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

-BV

Cool - is it available for download?

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.

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?

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)

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.

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))

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?

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.

There is an astropy where direct interconversion is possible

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