plm_gen

This routine computes the latitude dependent part $\lambda_{\ell m}$ of the spherical harmonics ( $Y_{\ell m}(\theta,\phi) = \lambda_{\ell m}(\theta) e^{i m \phi}$) of spin 0 and 2 (see HEALPix primer) used to synthetize or analyze HEALPix maps of temperature and polarisation. Since version 2.20, which introduced optimized algorithms for spherical harmonic transforms, it has become obsolete and should no longer be used.

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_{\mathrm{side}}$ value for which to compute the $\lambda_{\ell m}$.
nlmax I4B IN The maximum multipole order $\ell$ of the generated $\lambda_{\ell m}$.
nmmax I4B IN The maximum degree $m$ of the generated $\lambda_{\ell m}$.
plm(0:n_plm-1, 1:p) DP OUT The $\lambda_{\ell m}$ values, either for temperature only ($p=1$) or temperature and polarisation ($p=3$). The number of $\lambda_{\ell m}$ is n_plm = nsmax*(nmmax+1)*(2*nlmax-nmmax+2). They are stored in the order of increasing order $\ell$, 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 $N_{\mathrm{side}}=256$, up to 512 in $\ell$ 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_{\ell m}$ 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 $a_{\ell m}$ that can use precomputed $\lambda_{\ell m}$ generated by plm_gen
map2alm
routine analysing maps of temperature and polarisation that can use precomputed $\lambda_{\ell m}$ generated by plm_gen
plmgen
executable using plm_gen to compute the $\lambda_{\ell m}$ and writting them on disc

Version 3.82, 2022-07-28