Planetary maps


Planetary mapping routines.

phi = 0 faces toward the observer phi = pi thus faces away from the observer theta=pi/2 is the z-axis or ‘north pole’ theta=-pi/2 is the ‘south pole’

maps.errseries(param, phi, theta, rangle, meas, err, retmodel=False)[source]

Give the chi-squared to an input timeseries to get the wedge coefficients and inclination angle.

param – [inclination, phi0, wedge coefficients]

phi,theta – from MAKEGRID

rangle – rotation angles at which flux was measured

meas – measured flux

err – error on measurement (one-sigma)

retmodel – whether to return the model, instead of the chi-squared.

from scipy import optimize optimize.fmin(maps.errseries, guess, args=(...))
maps.makegrid(nphi, ntheta)[source]

Make grids of phi and theta values with the specified number of points in each direction. Phi ranges from 0 to 2pi, and theta ranges from -pi/2 to pi/2.

Returns (phi, theta)

maps.makespot(spotlat, spotlon, spotrad, phi, theta)[source]
spotlat : scalar

Latitude of spot center, in radians, from 0 to pi

spotlon : scalar

Longitude of spot center, in radians, from 0 to 2pi

spotrad : scalar

Radius of spot, in radians.

phi, theta : 2D NumPy arrays

output from makegrid(). Theta ranges from -pi/2 to +pi/2.

import maps
nlat, nlon = 60, 30
phi, theta = maps.makegrid(nlat, nlon)
# Make a small spot centered near, but not at, the equator:
equator_spot = maps.makespot(0, 0, 0.4, phi, theta)
# Make a larger spot centered near, but not at, the pole:
pole_spot = maps.makespot(1.2, 0, 0.7, phi, theta)
import maps
nlat, nlon = 60, 30
map =, nlon, i=0., deltaphi=0.)
phi = map.corners_latlon.mean(2)[:,1].reshape(nlon, nlat)
theta = map.corners_latlon.mean(2)[:,0].reshape(nlon, nlat) - np.pi/2.
# Make a small spot centered near, but not at, the equator:
equator_spot = maps.makespot(0, 0, 0.4, phi, theta)
# Make a larger spot centered near, but not at, the pole:
pole_spot = maps.makespot(1.2, 0, 0.7, phi, theta)
maps.makespot_old(phi, theta, da=None, vmap=None, inc=None, rot=None, long=0, lat=0, siz=0.7853981633974483, plotalot=False)[source]

Make a spot-map.

EXAMPLE: import maps inc, rot = pi/2., 0. phi, theta = maps.makegrid(120,60) phi2, theta2 = maps.rotcoord(phi,theta,inc,rot) vmap = maps.visiblemap(phi,theta,inc,rot)

maps.makespot(phi,theta,inc=pi/2,rot=0.,lat=pi/4,long=1.*pi) maps.makespot(phi2,theta2,vmap=vmap,lat=pi/4,long=1.*pi)

class, nlat=10, type='latlon', deltaphi=0, i=0)[source]

Very handy spherical mapping object. :INPUTS:

nlon, nlat : scalars
If mod==’latlon’, these inputs specify the number of grid cells across map, in latitude and longitude.
i : scalar
the inclination, is in units of radians. Zero means we see the object equator-on; pi/2 means we see it pole-on.
type : str
For now, only valid entry is ‘latlon’. Eventually, this could be expanded to allow other types of projection grids.
deltaphi : scalar
Rotation of map, specified in radians.
OUTPUT:A map-class object with various useful fields. Most of these fields refer to the coordinates (either Cartesian or spherical polar) or the projected radial velocities at the corners of specified grid cells, or the approximate projected areas of these grid cells.
NOTES:I have not been as careful as I should be in this code – my original goal was speed rather than exactitude. This means that some values are returned as ‘nan’, and the projected areas are only roughly correct. There’s plenty of room for improvement!

Compute velocity profile for normalized velocity values v.

v : NumPy array

Velocity normalized by the maximum rotation velocity observed; i.e., to convert v to true velocities, multiply by 2piR/P.

maps.polyarea(x, y)[source]

Compute the area of a polygon whose vertices are at the points (x,y).

x, y : 1D sequences

Cartesian coordinates of the (non-intersecting) polygon.


maps.projarea(phi, theta)[source]

Return a 2D map of the projected area. Note that you need to determine for yourself which parts of the areal map are actually visible to the observer.

Assumes phi and theta are grids straight out of meshgrid.

Assumes dA = cos(t) dt df (theta=t, phi=f)
and thus da = cos(f) cos(t)**2 dt df
import maps, pylab phi, theta = maps.makegrid(300,240) vis = maps.visiblemap(phi,theta,pi/2,0) da = maps.projarea(phi,theta) pylab.imshow(da*vis) pylab.title(‘visible projected area sums to: %f (pi)’ % (da*vis).sum())
maps.rotcoord(phi, theta, iangle, rangle)[source]

rotate coordinate system from local (planetary) coordinates into observer-oriented coordinages. rangle and iangle are the rotation and inclination angles of the planet in radians.

returns (phi2,theta2)

phi2 will be in the range (0,2pi) and theta2 will be in the range (-pi/2,pi/2)

maps.visiblemap(phi, theta, iangle, rangle)[source]

Return a 2D boolean map that’s True for the planetary latitude, longitude values visible from an observer. For rangle (rotation angle) zero and iangle (inclination angle) equal to pi/2, phi=pi is toward the observer.

maps.wedgebasis(phi, theta, iangle, rangle, nwedge, phi0, fwedge=None)[source]

Return a set of basis functions for the flux from each wedge.

phi,theta are the observer-centered grids from MAPS.MAKEGRID iangle is the inclination of the system (0 = pole-on) rangle (seq.) is a sequence of (0,2pi) rotation values at which

the wedge-based flux is evaluated.

nwedge is the number of wedges phi0 defines the center of wedge zero.

fwedge – an optional sequence to set individual flux values for each wedge.

maps.wedgemap(phi, theta, phi0, phi1)[source]

Return a 2D boolean map that’s True for the latitude, longitude values in a given longitudinal planet ‘wedge’.

phi0,phi1 should be in the range (0,2*pi)

maps.wedgestack(phi, theta, nwedge, phi0)[source]

Return a stack of 2D boolean maps via wedgemap.

phi0 defines the center of wedge zero.

Previous topic

Planetary phase curve routines

Next topic

Infrared Time-series Photometry

This Page