NASA - Jet Propulsion Laboratory
    + View the NASA Portal
Search JPL
Jet Propulsion Laboratory Home Earth Solar System Stars & Galaxies Technology
Introduction Background Software Links


plm_gen


This routine computes the latitude dependent part $\lambda_{lm}$ of the spherical harmonics ( $Y_{lm}(\theta,\phi) = \lambda_{lm}(\theta) e^{i m \phi}$) of spin 0 and 2 (see HEALPix primer) used to synthetize or analyze HEALPix maps of temperature and polarisation.


Location in HEALPix directory tree: src/f90/mod/alm_tools.f90


FORMAT

call plm_gen( nsmax, nlmax, nmmax, plm )


ARGUMENTS

name & dimensionality kind in/out description
       
nsmax I4B IN The ${N_{\rm side}}$ value for which to compute the $\lambda_{lm}$.
nlmax I4B IN The maximum multipole order l of the generated $\lambda_{lm}$.
nmmax I4B IN The maximum degree m of the generated $\lambda_{lm}$.
plm(0:n_plm-1, 1:p) DP OUT The $\lambda_{lm}$ values, either for temperature only (p=1) or temperature and polarisation (p=3). The number of $\lambda_{lm}$ is n_plm = nsmax*(nmmax+1)*(2*nlmax-nmmax+2). They are stored in the order of increasing order l, increasing degree m, for all the HEALPix ring colatitudes $\theta$ from North Pole to Equator, ie $\lambda_{00}(\theta_1)$, $\lambda_{10}(\theta_1)$, $\lambda_{20}(\theta_1)$, ..., $\lambda_{11}(\theta_1)$, $\lambda_{21}(\theta_1)$; ..., $\lambda_{00}(\theta_2)$ ...


EXAMPLE:

use healpix_types
use alm_tools, only : plm_gen
integer(I4B) :: nside, lmax, mmax, n_plm
real(DP), dimension(:,:), allocatable :: plm
...
nside=256 ; lmax=512 ; mmax=lmax
npix=nside2npix(nside)
n_plm=nside*(mmax+1)*(2*lmax-mmax+2)
allocate(plm(0:n_plm-1,1:3))
...
call plm_gen(nside, lmax, mmax, plm)

Computes the spherical harmonics for temperature and polarisation for Nside= 256, up to 512 in l and m.


MODULES & ROUTINES

This section lists the modules and routines used by plm_gen.

compute_lam_mm, get_pixel_layout,
gen_lamfac,gen_mfac, gen_normpol,
gen_recfac, init_rescale, l_min_ylm
Ancillary routines used for $\lambda_{lm}$ recursion
misc_utils
module, containing:
assert_alloc
routine to print error message, when an array can not be allocated properly


RELATED ROUTINES

This section lists the routines related to plm_gen


alm2map
routine generating maps of temperature and polarisation from their alm that can use precomputed $\lambda_{lm}$ generated by plm_gen
map2alm
routine analysing maps of temperature and polarisation that can use precomputed $\lambda_{lm}$ generated by plm_gen
plmgen
executable using plm_gen to compute the $\lambda_{lm}$ and writting them on disc

Eric Hivon 2010-06-18
Privacy / Copyright
FIRST GOV Contact: NASA Home Page Site Manager:
Webmaster:

CL 03-2650