As you might have guessed, L****Y found this code. I have not verified it, so use at your own risk.
"""
# spherical_coords.py
Convert between (theta, phi) and (azimuth, elevation) coordinate systems.
Author: Danny Price
License: MIT
## Overview
Both (theta, phi) and (azimuth, elevation) are spherical coordinate systems.
The first is used commonly in antenna pattern measurements, while the second
is used extensively in ground-based astronomy. Conversion between the two
essentially requires moving the location of the pole from one axis to another.
In (theta, phi) coordinates, phi is the angle from the y-axis toward the z-axis,
as measured from the yz plane. Phi runs [0, 360) degrees. Theta runs [0, 180).
In (azimuth, elevation) coordinates, azimuth is the angle from the x-axis
toward y-ax, as measured from the xy plane. Azimuth runs [-180, 180) degrees.
Elevation runs [-90 to 90) degrees.
A more complete explanation of the coordinate systems can be found
in [1] and [2]. Coordinate transforms, as coded up, are based on the
formulas presented in [1].
## SUPER IMPORTANT NOTE TO READ BEFORE USE
When writing some unit tests, I noticed that a round-trip conversion does
not return the identical angle pair that you started with. My first thought
was that this is due to the interval bounds on inverse functions (arcsin etc),
but some things seem a little odd still. So I'm releasing this as a gist, so
people can decide whether it's a useful start, and hopefully someone will
point to the problem and fix it one day.
## References
[1] phitheta2azel MATLAB documentation, Mathworks.
url: http://www.mathworks.com/help/phased/ref/phitheta2azel.html
[2] COORDINATE SYSTEM PLOTTING FOR ANTENNA MEASUREMENTS, Gregory F. Masters
and Stuart F. Gregson. Tech. Report, Nearfield Systems, 2007.
url: http://ww2.nearfield.com/amta/AMTA07-0092-GFM_SFG.pdf
"""
import numpy as np
def azel_to_thetaphi(az, el):
""" Az-El to Theta-Phi conversion.
Args:
az (float or np.array): Azimuth angle, in radians
el (float or np.array): Elevation angle, in radians
Returns:
(theta, phi): Tuple of corresponding (theta, phi) angles, in radians
"""
cos_theta = np.cos(el) * np.cos(az)
tan_phi = np.tan(el) / np.sin(az)
theta = np.arccos(cos_theta)
phi = np.arctan2(np.tan(el), np.sin(az))
phi = (phi + 2 * np.pi) % (2 * np.pi)
return theta, phi
def thetaphi_to_azel(theta, phi):
""" Az-El to Theta-Phi conversion.
Args:
theta (float or np.array): Theta angle, in radians
phi (float or np.array): Phi angle, in radians
Returns:
(az, el): Tuple of corresponding (azimuth, elevation) angles, in radians
"""
sin_el = np.sin(phi) * np.sin(theta)
tan_az = np.cos(phi) * np.tan(theta)
el = np.arcsin(sin_el)
az = np.arctan(tan_az)
return az, el
def thetaphi_to_uv(theta, phi):
""" Theta-phi to U-V direction cosines
Args:
theta (float or np.array): Theta angle, in radians
phi (float or np.array): Phi angle, in radians
Returns:
(u, v): Tuple of corresponding (u, v) direction cosines
"""
u = np.sin(theta) * np.cos(phi)
v = np.sin(theta) * np.sin(phi)
w = np.cos(theta)
return u, v
def azel_to_uv(az, el):
""" Azimuth-elevation to U-V direction cosines
Args:
az (float or np.array): azimuth angle, in radians
el (float or np.array): elevation angle, in radians
Returns:
(u, v): Tuple of corresponding (u, v) direction cosines
"""
u = np.cos(el) * np.sin(az)
v = np.sin(el)
w = np.cos(az) * np.cos(el)
return u, v
def uv_to_thetaphi(u, v):
""" U-V direction cosines to theta-phi coordinates
Args:
u (float or np.array): U direction cosine
v (float or np.array): V direction cosine
Returns:
(theta, phi): Tuple of corresponding (theta, phi) angles, in radians
"""
theta = np.arcsin(np.sqrt(u**2 + v**2))
phi = np.arctan2(u, v)
phi = (phi + 2 * np.pi) % (2 * np.pi)
return theta, phi
def uv_to_azel(u, v):
""" U-V direction cosines to azimuth-elevation coordinates
Args:
u (float or np.array): U direction cosine
v (float or np.array): V direction cosine
Returns:
(az, el): Tuple of corresponding (azimuth, elevation) angles, in radians
"""
az = np.arctan2(u, np.sqrt(1 - u**2 - v**2))
el = np.arcsin(v)
return az, el
def plot_compare(a, b, ap, bp, title='', show=False, nf=True):
""" Helper function for debugging, plotting before & after """
if nf:
plt.figure()
plt.subplot(211)
plt.plot(a, ap, '.', label=title)
plt.legend()
plt.subplot(212)
plt.plot(b, bp, '.', label=title)
if show:
plt.show()
if __name__ == "__main__":
import pylab as plt
# Generate vectors covering intervals of theta & phi
theta = np.linspace(0, np.pi, 101)
phi = np.linspace(0, 2*np.pi, 101)
# Convert to other coordinates and back
u, v = thetaphi_to_uv(theta, phi)
theta_p, phi_p = uv_to_thetaphi(u, v)
plot_compare(theta, phi, theta_p, phi_p,
title='T-P -> U-V -> T-P')
az, el = thetaphi_to_azel(theta, phi)
theta_p, phi_p = azel_to_thetaphi(az, el)
plot_compare(theta, phi, theta_p, phi_p,
title='T-P -> A-E -> T-P', nf=False)
u, v = thetaphi_to_uv(theta, phi)
az, el = uv_to_azel(u, v)
theta_p, phi_p = azel_to_thetaphi(az, el)
plot_compare(theta, phi, theta_p, phi_p,
title='T-P -> U-V -> A-E -> T-P', nf=False)
u, v = thetaphi_to_uv(theta, phi)
az, el = uv_to_azel(u, v)
u, v = azel_to_uv(az, el)
theta_p, phi_p = uv_to_thetaphi(u, v)
plot_compare(theta, phi, theta_p, phi_p,
title='T-P -> U-V -> A-E -> U-V -> T-P', nf=False)
# Convert to Az-el, two methods
u, v = thetaphi_to_uv(theta, phi)
az, el = uv_to_azel(u, v)
plot_compare(theta, phi, az, el,
title='T-P -> U-V -> A-E')
az, el = thetaphi_to_azel(theta, phi)
plot_compare(theta, phi, az, el,
title='T-P -> A-E', nf=False)
# Generate vectors covering intervals of azimuth & elevation
az = np.linspace(-np.pi, np.pi, 101)
el = np.linspace(-np.pi/2, np.pi/2, 101)
# Convert to UV and back, compare
u, v = azel_to_uv(az, el)
az_p, el_p = uv_to_azel(u, v)
plot_compare(az, el, az_p, el_p,
title='A-E -> U-V -> A-E')
theta, phi = azel_to_thetaphi(az, el)
az_p, el_p = thetaphi_to_azel(theta, phi)
plot_compare(az, el, az_p, el_p,
title='A-E -> T-P -> A-E', nf=False)
u, v = azel_to_uv(az, el)
theta, phi = uv_to_thetaphi(u, v)
az_p, el_p = thetaphi_to_azel(theta, phi)
plot_compare(az, el, az_p, el_p,
title='A-E -> U-V -> T-P -> A-E', nf=False)
u, v = azel_to_uv(az, el)
theta, phi = uv_to_thetaphi(u, v)
u, v = thetaphi_to_uv(theta, phi)
az_p, el_p = uv_to_azel(u, v)
plot_compare(az, el, az_p, el_p,
title='A-E -> U-V -> T-P -> U-V -> A-E', nf=False)
# Convert to Theta-phi, two methods
u, v = thetaphi_to_uv(az, el)
theta, phi = uv_to_thetaphi(u, v)
plot_compare(az, el, theta, phi,
title='A-E -> U-V -> T-P')
theta, phi = azel_to_thetaphi(az, el)
plot_compare(az, el, theta, phi,
title='A-E -> T-P', nf=False)
plt.show()
