2

I have to say that I'm terrified and surprised how little I know about basic math. Essentially what I have is an origin point (0,0,0), I know the radius of the circle (10), and I know both angles (theta and phi). Given this assumptions I want to calculate the projected point on the sphere. I've came up with the bottom code through reading the answers of https://stackoverflow.com/a/969880/1230358, https://stackoverflow.com/a/36369852/1230358, http://tutorial.math.lamar.edu/Classes/CalcII/SphericalCoords.aspx and https://en.wikipedia.org/wiki/Spherical_coordinate_system.

My current code:

#!/usr/bin/env python3
import math


PI = math.pi
PI_2 = PI / 2

def calc_sphere_coordinates(radius, phi, theta):
    # see: https://stackoverflow.com/a/969880/1230358
    # see: https://stackoverflow.com/q/19673067/1230358
    # see: http://mathinsight.org/spherical_coordinates
    # see: https://en.wikipedia.org/wiki/Spherical_coordinate_system
    # see: http://tutorial.math.lamar.edu/Classes/CalcII/SphericalCoords.aspx

    # φ phi is the polar angle, rotated down from the positive z-axis (slope)
    # θ theta is azimuthal angle, the angle of the rotation around the z-axis (aspect)

    # z        
    # |  x
    # | /
    # |/
    # +-------- y

    # both angles need to be in radians, not degrees!
    theta = theta * PI / 180
    phi = phi * PI / 180

    x = radius * math.sin(phi) * math.cos(theta)
    y = radius * math.sin(phi) * math.sin(theta)
    z = radius * math.cos(phi)
    return (x, y, z)

if __name__ == "__main__":

    # calculate point position in hemisphere by rotating down from positive z-axis
    for i in (10, 20, 30, 40, 50, 60, 70 , 80, 90):
        print(calc_sphere_coordinates(10, i, 0))

    print("-"*10)

    # calculate point position in hemisphere by rotating around the z axis
    for i in (10, 20, 30, 40, 50, 60, 70 , 80, 90):
        print(calc_sphere_coordinates(10, 0, i))

    print("-"*10)

    # calculate point position by rotating in both directions
    for i in (10, 20, 30, 40, 50, 60, 70 , 80, 90):
        print(calc_sphere_coordinates(10, i, 90-i))    

The output of the code is as follows:

(1.7364817766693033, 0.0, 9.84807753012208)
(3.420201433256687, 0.0, 9.396926207859085)
(4.999999999999999, 0.0, 8.660254037844387)
(6.4278760968653925, 0.0, 7.660444431189781)
(7.66044443118978, 0.0, 6.427876096865393)
(8.660254037844386, 0.0, 5.000000000000001)
(9.396926207859083, 0.0, 3.4202014332566884)
(9.84807753012208, 0.0, 1.7364817766693041)
(10.0, 0.0, 6.123233995736766e-16)
----------
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
(0.0, 0.0, 10.0)
----------
(0.30153689607045814, 1.7101007166283433, 9.84807753012208)
(1.16977778440511, 3.2139380484326963, 9.396926207859085)
(2.5, 4.330127018922192, 8.660254037844387)
(4.131759111665348, 4.92403876506104, 7.660444431189781)
(5.868240888334652, 4.92403876506104, 6.427876096865393)
(7.5, 4.330127018922192, 5.000000000000001)
(8.83022221559489, 3.2139380484326963, 3.4202014332566884)
(9.69846310392954, 1.7101007166283433, 1.7364817766693041)
(10.0, 0.0, 6.123233995736766e-16)

enter image description here

Shouldn't the line with (10.0, 0.0, 6.123233995736766e-16) have a z-coordinate of 0 and not 6.123233995736766e-16? Also rotating around the z-axis gives always the same result, no matter what angle is used (0.0, 0.0, 10.0).

Community
  • 1
  • 1
hetsch
  • 1,508
  • 2
  • 12
  • 27

1 Answers1

0

As far as I can see, your code is working just fine. Being honest, one has to admit that 6.123233995736766e-16 is pretty much 0 for all real applications, right?

Your question boils down to

why math.cos(math.pi / 2.0) does not equal zero

The reason is found in how floating point numbers are generated and stored in a computer. Try calculating 0.1+0.2 in python. Supprised?! In case you want to know more, just google floating point error or anything related.

ImportanceOfBeingErnest
  • 321,279
  • 53
  • 665
  • 712
  • All right, strangely I've never realized this problem. Indeed, I'm a little bit surprised :)! Still, I think one problem exists with the rotation around the z-axis. I've added a chart to better illustrate the problem. I don't think that my calculation equations are ok, because rotating around the z-axis (2nd chart) always returns the same value ```(0.0, 0.0, 10.0)``` resulting in all points being at the same position. Shouldn't they be scattered in a nice bow as in the first chart but arranged around the z-axis? – hetsch Sep 23 '16 at 08:03
  • @hetsch As I said, the equations shoudl be fine. What you're doing in the second chart is rotating a point about the z-axis which already sits on the this very axis, thus it'll stay where it is. Maybe you try `print [calc_sphere_coordinates(10, 90, i) for i in [15,30,45,60,75,90] ]` or something similar instead. – ImportanceOfBeingErnest Sep 29 '16 at 23:08