Spherical Harmonic Transformations

SHExpandDH(grid, lmax_calc, degmax=None, norm=4, sampling=2, csphase=1)

This routine will expand a grid containing N samples in both longitude and latitude (or N x 2N, see sampling) into spherical harmonics. This routine makes use of the sampling theorem presented in Driscoll and Healy (1994) and employs FFTs when calculating the sine and cosine terms. The number of samples, N, must be EVEN for this routine to work, and the spherical harmonic expansion is exact if the function is band limited to degree N/2-1. Legendre functions are computed on the fly using the scaling methodology presented in Holmes and Featherston (2002). When norm is 1,2 or 4, these are accurate to about degree 2800. When norm is 3, the routine is only stable to about degree 15.

Parameters:

• grid (numpy array) – the 2D grid to decompose into SH coefficients. If sampling == 1: (N,N) and if sampling == 2: (N,2N)
• lmax_calc (int) – maximum spherical harmonic degree for which the coefficients will be calculated
• degmax (integer or None) – maximum degree represented in the cilm array (if degmax=None, maximum degree in the array is lmax_calc). If degmax > lmax_calc then the coefficients for degrees greater than lmax_calc are set to zero
• norm (int) – spherical harmonic normalization
• sampling (int) – the shape of the data array that will be decomposed in SH coefficients
• csphase (int) – Condon-Shortley phase factor

Returns:

coefficients array of shape (2, l+1, l+1) where l=max(degmax, lmax_calc)

Return type:

numpy array

SHExpandLSQ(d, lat, lon, lmax, norm=4, csphase=1)

This routine will expand a set of discrete data points into spherical harmonics using a least squares inversion. When there are more data points than spherical harmonic coefficients (len(d) > (lmax+1)**2) the solution of the overdetermined system will be determined by least squares. If there are more coefficients than data points, then the solution of the underdetermined system will be determined by minimizing the solution norm. (See LAPACK DGELS documentation).

Note that this routine takes lots of memory (~ 8*nmax*(lmax+1)**2 bytes) and is very slow for large lmax.

Parameters:

• d (numpy array) – 1D array of raw data points
• lat (numpy array) – 1D array of corresponding latitudes
• lon (numpy array) – 1D array of corresponding longitudes
• lmax (int) – maximum spherical harmonic degree for which the coefficients will be calculated
• norm (int) – spherical harmonic normalization
• csphase (int) – Condon-Shortley phase factor

Returns:

coefficients array of shape (2, lmax+1, lmax+1)

Return type:

numpy array

Returns:

residual sum of squares misfit for an overdetermined inversion

Return type:

real

MakeGrid2D(cilm, ny, nx, norm=4, csphase=1, north=90.0, south=-90.0, east=360.0, west=0.0)

Given the Spherical Harmonic coefficients cilm, this routine will compute a 2D grid with equal latitude and longitude spacings. The spacing is determined based on the bounding coordinate values represented by north, south, east, west and the number of latitude and longitude points - ny, nx respectively. If the mix of these values does not produce equal spacing in both latitude and longitude, then the program will raise an error.

The values in the returned grid are in the order of latitude going from 90 to -90 and longitudes going from 0 to 360 (both these longitudes present in the returned grid). If the optional parameters north, south, east and west are specified, the upper-left and lower right coordinates of the output grid are (north, west) and (south, east), respectively.

Note that since this routined does not use FFTs, this routine is therefore slower than MakeGridDH. However MakeGridDH does not allow for the specification of arbitray dimensions of the returned grid.

Parameters:

• cilm (numpy array) – SH coefficients array of shape (2, lmax+1, lmax+1), where lmax represents the maximum degree for which coefficients are present in the array
• ny (int) – number of latitude points in the output grid
• nx (int) – number of longitude points in the output grid
• norm (int) – spherical harmonic normalization
• csphase (int) – Condon-Shortley phase factor
• north (real) – The north bounding latitude in the returned grid
• south (real) – The south bounding latitude in the returned grid
• east (real) – The east bounding latitude in the returned grid
• west (real) – The west bounding latitude in the returned grid

Returns:

2D grid that has been constructed from the coefficients. If sampling==1, the shape of the array is (N,N) and for sampling==2 the shape is (N,2N), where N=2*l+1.

Return type:

numpy array

MakeGridDH(cilm, l, norm=4, csphase=1, sampling=2)

Parameters:

• cilm (numpy array) – SH coefficients array of shape (2, lmax+1, lmax+1), where lmax represents the maximum degree for which coefficients are present in the array
• norm (int) – spherical harmonic normalization
• csphase (int) – Condon-Shortley phase factor
• sampling (int) – the shape of the data array that will be decomposed in SH coefficients

Returns:

2D grid that has been constructed from the coefficients. If sampling==1, the shape of the array is (N,N) and for sampling==2 the shape is (N,2N), where N=2*l+1.

Return type:

numpy array