Sunday, August 15, 2010

Spherical to Cartesian Coords - Thanks Python!

Game and graphics programming requires the modeling of three-dimensional geometries in software.  This stuff is relatively simple in theory but easy to screw up and difficult to debug. 

The math used in 3-D geometrical programming has its limitations and special cases.  One example is the gimbal lock problem that occurs when rotating a body using matrices and Euler angles.

Another problem that I recently stumbled across occurred when trying to find the intersection point of a line and plane.  If the line is pointing away from the plane, the intersection point ends up being behind the line's start point.  Essentially, if we shoot an arrow directly away from a target, it will still hit the target, but with the back of the arrow.

The problem I describe above with the intersection point, line and plane is not really a problem at all.  The special case is easily detected while computing the result and the code can respond as required.

The real issue is that we must have a complete understanding of the mathematical techniques we use in our software.

There are two ways to understand math, one is to use analysis, the other is intuition/visualization.  Analysis is the stronger of the two methods and includes things such as proofs.  Intuition, visualization and experimentation is the approach that I take.

My latest interest is converting spherical coordinates into Cartesian coordinates.  I can easily plug and chug using the standard formulas, but after learning a few hard lessons, I know I had better gain some understanding.

When it comes to three-dimensional programming, it is not as easy to visualize 3-D space as one would imagine.  I needed some help to visualize this stuff.

Locations on earth are identified using spherical coordinates, using the intersection of two angles, an angle of inclination (latitude) and an angle of azimuth (longitude).  In my investigation, I would like to consider inclination angles from 0 to 180 degrees, azimuths from 0 to 360 degrees and a spherical radius of one.



The plot simplifies the problem by reducing three-dimensions to a two-dimensional display with circles of common inclination.  This is essentially a polar-plot. Each circle of common inclination (latitude) represents two inclinations, the inclination (northern hemisphere) and 180 - inclination (southern hemisphere).

As a two-dimensional polar plot, the problem is solved using the azimuth angle, and the distance from the origin to the proper circle of common inclination.

It appears from this work, that spherical to cartesian coordinate conversions is a simple and safe process (with the exception of when the inclination angle is 0 or 180 degrees, of course).


import numpy as np
import matplotlib.pylab as plt

def x(inclin, azi):
""" compute the x cartesian coord """
   return np.sin(inclin) * np.cos(azi)

def y(inclin, azi):
""" compute the y cartesian coord """
   return np.sin(inclin) * np.sin(azi)

def z(inclin):
""" compute the z cartesian coord """
   return np.cos(inclin)



inclin = []
inclin.append([np.pi/2 for i in range(0,20)])
inclin.append([np.pi/3 for i in range(0,20)])
inclin.append([np.pi/4 for i in range(0,20)])
inclin.append([np.pi/6 for i in range(0,20)])
inclin.append([np.pi/8 for i in range(0,20)])
inclin.append([np.pi/12 for i in range(0,20)])
inclin.append([np.pi/24 for i in range(0,20)])
inclin.append([0 for i in range(0,20)])
azi = np.arange(-np.pi, np.pi, (2 * np.pi)/float(20))

plt.plot(
   x(inclin[0],azi), y(inclin[0],azi), "g-o",
   x(inclin[1],azi), y(inclin[1],azi), "r-o",
   x(inclin[2],azi), y(inclin[2],azi), "c-o",
   x(inclin[3],azi), y(inclin[3],azi), "m-o",
   x(inclin[4],azi), y(inclin[4],azi), "y-o",
   x(inclin[5],azi), y(inclin[5],azi), "k-o",
   x(inclin[6],azi), y(inclin[6],azi), "b-o",
   x(inclin[7],azi), y(inclin[7],azi), "w-o")
plt.legend(
   ('90', '60, 120', '45, 135',
   '30, 150', '23, 157',
   '15, 165', '8, 172', '0, 180'))
plt.grid(True)
plt.ylabel('Y Cartesian Coord')
plt.title('Cartesian from Spherical')
plt.xlabel('X Cartesian Coord')
plt.show()


No comments: