SHTOOLS.jl

The SHTOOLS.jl package wraps SHTOOLS, the Spherical Harmonic Tools.

SHTOOLS is an archive of Fortran 95 software that can be used to perform spherical harmonic transforms, multitaper spectral analyses, expansions of functions into Slepian bases, and standard operations on global gravitational and magnetic field data.

Legendre functions

4π normalized

Orthonormalized

Schmidt semi-normalized

Unnormalized

Utilities

SHTOOLS.PlmIndexFunction
index = PlmIndex(l::Integer,
                 m::Integer)
index::Int

Compute the index of an array of Legendre functions corresponding to degree l and angular order m.

source

Spherical harmonic transforms

Equally sampled (N×N) and equally spaced (N×2N) grids

SHTOOLS.SHExpandDH!Function
cilm, lmax = SHExpandDH!(cilm::AbstractArray{Cdouble,3},
                         griddh::AbstractArray{Cdouble,2},
                         n::Integer;
                         norm::Integer=1,
                         sampling::Integer=1,
                         csphase::Integer=1,
                         lmax_calc::Union{Nothing,Integer}=nothing,
                         exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::AbstractArray{Cdouble,3}
lmax:Int

cilm, lmax = SHExpandDH!(cilm::AbstractArray{Complex{Cdouble},3},
                         griddh::AbstractArray{Complex{Cdouble},2},
                         n::Integer;
                         norm::Integer=1,
                         sampling::Integer=1,
                         csphase::Integer=1,
                         lmax_calc::Union{Nothing,Integer}=nothing,
                         exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::AbstractArray{Complex{Cdouble},3}
lmax:Int

Expand an equally sampled or equally spaced real or complex map into real or complex spherical harmonics using Driscoll and Healy’s (1994) sampling theorem.

See also: SHExpandDH, MakeGridDH!, MakeGradientDH!

source
SHTOOLS.SHExpandDHFunction
cilm, lmax = SHExpandDH(griddh::AbstractArray{Cdouble,2},
                        n::Integer;
                        norm::Integer=1,
                        sampling::Integer=1,
                        csphase::Integer=1,
                        lmax_calc::Union{Nothing,Integer}=nothing,
                        exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::Array{Cdouble,3}
lmax:Int

cilm, lmax = SHExpandDH(griddh::AbstractArray{Complex{Cdouble},2},
                        n::Integer;
                        norm::Integer=1,
                        sampling::Integer=1,
                        csphase::Integer=1,
                        lmax_calc::Union{Nothing,Integer}=nothing,
                        exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::Array{Complex{Cdouble},3}
lmax:Int

Expand an equally sampled or equally spaced real or complex map into real or complex spherical harmonics using Driscoll and Healy’s (1994) sampling theorem.

See also: SHExpandDH!, MakeGridDH, MakeGradientDH

source
SHTOOLS.MakeGridDH!Function
MakeGridDH!(griddh::AbstractArray{Cdouble,2},
            cilm::AbstractArray{Cdouble,3},
            lmax::Integer;
            norm::Integer=1,
            sampling::Integer=1,
            csphase::Integer=1,
            lmax_calc::Union{Nothing,Integer}=nothing,
            extend::Integer=0,
            exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
griddh::AbstractArray{Cdouble,2}
n::Int

MakeGridDH!(griddh::AbstractArray{Complex{Cdouble},2},
            cilm::AbstractArray{Complex{Cdouble},3},
            lmax::Integer;
            norm::Integer=1,
            sampling::Integer=1,
            csphase::Integer=1,
            lmax_calc::Union{Nothing,Integer}=nothing,
            extend::Integer=0,
            exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
griddh::AbstractArray{Complex{Cdouble},2}
n::Int

Create a real or complex 2D map from a set of real or complex spherical harmonic coefficients that conforms with Driscoll and Healy’s (1994) sampling theorem.

See also: SHExpandDH!, MakeGridDH, MakeGradientDH!

source
SHTOOLS.MakeGridDHFunction
MakeGridDH(cilm::AbstractArray{Cdouble,3},
           lmax::Integer;
           norm::Integer=1,
           sampling::Integer=1,
           csphase::Integer=1,
           lmax_calc::Union{Nothing,Integer}=nothing,
           extend::Integer=0,
           exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
griddh::Array{Cdouble,2}
n::Int

MakeGridDH(cilm::AbstractArray{Complex{Cdouble},3},
           lmax::Integer;
           norm::Integer=1,
           sampling::Integer=1,
           csphase::Integer=1,
           lmax_calc::Union{Nothing,Integer}=nothing,
           extend::Integer=0,
           exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
griddh::Array{Complex{Cdouble},2}
n::Int

Create a real or complex 2D map from a set of real or complex spherical harmonic coefficients that conforms with Driscoll and Healy’s (1994) sampling theorem.

See also: SHExpandDH, MakeGridDH!, MakeGradientDH

source
SHTOOLS.MakeGradientDH!Function
MakeGradientDH!(theta::AbstractArray{Cdouble,2},
                phi::AbstractArray{Cdouble,2},
                cilm::AbstractArray{Cdouble,3},
                lmax::Integer;
                norm::Integer=1,
                sampling::Integer=1,
                csphase::Integer=1,
                lmax_calc::Union{Nothing,Integer}=nothing,
                extend::Integer=0,
                radius::Union{Nothing,AbstractFloat,Integer}=nothing,
                exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
theta::AbstractArray{Cdouble,2}
phi::AbstractArray{Cdouble,2}
n::Int

Compute the gradient of a scalar function and return grids of the two horizontal components that conform with Driscoll and Healy’s (1994) sampling theorem.

See also: MakeGradientDH, SHExpandDH!, MakeGridDH!

source
SHTOOLS.MakeGradientDHFunction
theta, phi, n = MakeGradientDH(cilm::AbstractArray{Cdouble,3},
                               lmax::Integer;
                               norm::Integer=1,
                               sampling::Integer=1,
                               csphase::Integer=1,
                               lmax_calc::Union{Nothing,Integer}=nothing,
                               extend::Integer=0,
                               exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
theta::Array{Cdouble,2}
phi::Array{Cdouble,2}
n::Int

Compute the gradient of a scalar function and return grids of the two horizontal components that conform with Driscoll and Healy’s (1994) sampling theorem.

See also: MakeGradientDH!, SHExpandDH, MakeGridDH

source

Gauss-Legendre quadrature grids

SHTOOLS.SHGLQ!Function
zero, w = SHGLQ!(zero::AbstractVector{Cdouble},
                 w::AbstractVector{Cdouble},
                 plx::Union{Nothing,AbstractArray{Cdouble,2}},
                 lmax::Integer;
                 norm::Integer=1,
                 csphase::Integer=1,
                 cnorm::Integer=0,
                 exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
zero::AbstractVector{Cdouble}
w::AbstractVector{Cdouble}

Precompute weights, nodes, and associated Legendre functions used in the GLQ-based spherical harmonics routines.

See also: SHGLQ, SHExpandGLQ!, MakeGridGLQ!, GLQGridCoord!

source
SHTOOLS.SHGLQFunction
zero, w = SHGLQ(plx::Union{Nothing,AbstractArray{Cdouble,2}},
                lmax::Integer;
                norm::Integer=1,
                csphase::Integer=1,
                cnorm::Integer=0,
                exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
zero::Vector{Cdouble}
w::Vector{Cdouble}

Precompute weights, nodes, and associated Legendre functions used in the GLQ-based spherical harmonics routines.

See also: SHGLQ!, SHExpandGLQ, MakeGridGLQ, GLQGridCoord

source
SHTOOLS.SHExpandGLQ!Function
cilm = SHExpandGLQ!(cilm::AbstractArray{Cdouble,3},
                    lmax::Integer,
                    gridglq::AbstractArray{Cdouble,2},
                    w::AbstractVector{Cdouble},
                    plx::Union{Nothing,AbstractArray{Cdouble,2}},
                    zero::Union{Nothing,AbstractVector{Cdouble}};
                    norm::Integer=1,
                    csphase::Integer=1,
                    lmax_calc::Union{Nothing,Integer}=nothing,
                    exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::AbstractArray{Cdouble,3},

cilm = SHExpandGLQ!(cilm::AbstractArray{Complex{Cdouble},3},
                    lmax::Integer,
                    gridglq::AbstractArray{Complex{Cdouble},2},
                    w::AbstractVector{Cdouble},
                    plx::Union{Nothing,AbstractArray{Cdouble,2}},
                    zero::Union{Nothing,AbstractVector{Cdouble}};
                    norm::Integer=1,
                    csphase::Integer=1,
                    lmax_calc::Union{Nothing,Integer}=nothing,
                    exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::AbstractArray{Complex{Cdouble},3},

Expand a 2D real or complex map sampled on the Gauss-Legendre quadrature nodes into real or complex spherical harmonics.

See also: SHExpandGLQ, SHGLQ!, MakeGridGLQ!, GLQGridCoord!

source
SHTOOLS.SHExpandGLQFunction
cilm = SHExpandGLQ(lmax::Integer,
                   gridglq::AbstractArray{Cdouble,2},
                   w::AbstractVector{Cdouble},
                   plx::Union{Nothing,AbstractArray{Cdouble,2}},
                   zero::Union{Nothing,AbstractVector{Cdouble}};
                   norm::Integer=1,
                   csphase::Integer=1,
                   lmax_calc::Union{Nothing,Integer}=nothing,
                   exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::Array{Cdouble,3},

cilm = SHExpandGLQ(lmax::Integer,
                   gridglq::AbstractArray{Complex{Cdouble},2},
                   w::AbstractVector{Cdouble},
                   plx::Union{Nothing,AbstractArray{Cdouble,2}},
                   zero::Union{Nothing,AbstractVector{Cdouble}};
                   norm::Integer=1,
                   csphase::Integer=1,
                   lmax_calc::Union{Nothing,Integer}=nothing,
                   exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::Array{Complex{Cdouble},3},

Expand a 2D real or complex map sampled on the Gauss-Legendre quadrature nodes into real or complex spherical harmonics.

See also: SHExpandGLQ!, SHGLQ, MakeGridGLQ, GLQGridCoord

source
SHTOOLS.MakeGridGLQ!Function
gridglq = MakeGridGLQ!(gridglq::AbstractArray{Cdouble,2},
                       cilm::AbstractArray{Cdouble,3},
                       lmax::Integer,
                       plx::Union{Nothing,AbstractArray{Cdouble,2}},
                       zero::Union{Nothing,AbstractVector{Cdouble}};
                       norm::Integer=1,
                       csphase::Integer=1,
                       lmax_calc::Union{Nothing,Integer}=nothing,
                       extend::Integer=0,
                       exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
gridqlq::AbstractArray{Cdouble,2}

gridglq = MakeGridGLQ!(gridglq::AbstractArray{Complex{Cdouble},2},
                       cilm::AbstractArray{Complex{Cdouble},3},
                       lmax::Integer,
                       plx::Union{Nothing,AbstractArray{Cdouble,2}},
                       zero::Union{Nothing,AbstractVector{Cdouble}};
                       norm::Integer=1,
                       csphase::Integer=1,
                       lmax_calc::Union{Nothing,Integer}=nothing,
                       extend::Integer=0,
                       exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
gridqlq::AbstractArray{Complex{Cdouble},2}

Create a 2D real or complex map from a set of real or complex spherical harmonic coefficients sampled on a the Gauss-Legendre quadrature nodes.

See also: MakeGridGLQ, SHGLQ!, SHExpandGLQ!, GLQGridCoord!

source
SHTOOLS.MakeGridGLQFunction
gridglq = MakeGridGLQ(cilm::AbstractArray{Cdouble,3},
                      lmax::Integer,
                      plx::Union{Nothing,AbstractArray{Cdouble,2}},
                      zero::Union{Nothing,AbstractVector{Cdouble}};
                      norm::Integer=1,
                      csphase::Integer=1,
                      lmax_calc::Union{Nothing,Integer}=nothing,
                      extend::Integer=0,
                      exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
gridqlq::Array{Cdouble,2}

gridglq = MakeGridGLQ(cilm::AbstractArray{Complex{Cdouble},3},
                      lmax::Integer,
                      plx::Union{Nothing,AbstractArray{Cdouble,2}},
                      zero::Union{Nothing,AbstractVector{Cdouble}};
                      norm::Integer=1,
                      csphase::Integer=1,
                      lmax_calc::Union{Nothing,Integer}=nothing,
                      extend::Integer=0,
                      exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
gridqlq::Array{Complex{Cdouble},2}

Create a 2D real or complex map from a set of real or complex spherical harmonic coefficients sampled on a the Gauss-Legendre quadrature nodes.

See also: MakeGridGLQ!, SHGLQ, SHExpandGLQ, GLQGridCoord

source
SHTOOLS.GLQGridCoord!Function
latglq, longlq, nlat, nlong =
    GLQGridCoord!(latglq::AbstractVector{Cdouble},
                  longlq::AbstractVector{Cdouble},
                  lmax::Integer;
                  extend::Integer=0,
                  exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
latglq::AbstractVector{Cdouble},
longlq::AbstractVector{Cdouble},
nlat::Int
nlong::Int

Compute the latitude and longitude coordinates used in Gauss-Legendre quadrature grids.

See also: GLQGridCoord, SHExpandGLQ!, MakeGridGLQ!

source
SHTOOLS.GLQGridCoordFunction
latglq, longlq, nlat, nlong =
    GLQGridCoord(lmax::Integer;
                 extend::Integer=0,
                 exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
latglq::Vector{Cdouble},
longlq::Vector{Cdouble},
nlat::Int
nlong::Int

Compute the latitude and longitude coordinates used in Gauss-Legendre quadrature grids.

See also: GLQGridCoord!, SHExpandGLQ, MakeGridGLQ

source

Other routines

SHTOOLS.SHExpandLSQ!Function
cilm, chi = SHExpandLSQ!(cilm::AbstractArray{Cdouble,3},
                         d::AbstractVector{Cdouble},
                         lat::AbstractVector{Cdouble},
                         lon::AbstractVector{Cdouble},
                         nmax::Integer,
                         lmax::Integer;
                         norm::Integer=1,
                         csphase::Integer=1,
                         weights::Vector{Cdouble},
                         exitstatus::Optional{Ref{<:Integer}}=nothing)
cilm::AbstractArray{Cdouble,3},
chi::Float64

Expand a set of irregularly sampled data points into spherical harmonics using a (weighted) least squares inversion.

See also: SHExpandLSQ, MakeGridPoint, SHExpandDH!, SHExpandGLQ!

source
SHTOOLS.SHExpandLSQFunction
cilm, chi = SHExpandLSQ(d::AbstractVector{Cdouble},
                        lat::AbstractVector{Cdouble},
                        lon::AbstractVector{Cdouble},
                        nmax::Integer,
                        lmax::Integer;
                        norm::Integer=1,
                        csphase::Integer=1,
                        weights::Vector{Cdouble},
                        exitstatus::Optional{Ref{<:Integer}}=nothing)
cilm::Array{Cdouble,3},
chi::Cdouble

Expand a set of irregularly sampled data points into spherical harmonics using a (weighted) least squares inversion.

See also: SHExpandLSQ!, MakeGridPoint, SHExpandDH, SHExpandGLQ

source
SHTOOLS.MakeGrid2d!Function
grid, nlat, nlong =
    MakeGrid2d!(grid::AbstractArray{Cdouble,2},
                cilm::AbstractArray{Cdouble,3},
                lmax::Integer,
                interval::Union{AbstractFloat,Integer};
                norm::Integer=1,
                csphase::Integer=1,
                f::Optional{Union{AbstractFloat,Integer}}=nothing,
                a::Optional{Union{AbstractFloat,Integer}}=nothing,
                north::Union{AbstractFloat,Integer}=90.0,
                south::Union{AbstractFloat,Integer}=-90.0,
                east::Union{AbstractFloat,Integer}=360.0,
                west::Union{AbstractFloat,Integer}=0.0,
                dealloc::Bool=false,
                exitstatus::Optional{Ref{<:Integer}}=nothing)
grid::AbstractArray{Cdouble,2}
nlat::Int
nlong::Int

Create a 2D cylindrical map with arbitrary grid spacing from a set of spherical harmonic coefficients.

See also: MakeGrid2d, MakeGridPoint

source
SHTOOLS.MakeGrid2dFunction
grid, nlat, nlong =
    MakeGrid2d(cilm::AbstractArray{Cdouble,3},
               lmax::Integer,
               interval::Union{AbstractFloat,Integer};
               norm::Integer=1,
               csphase::Integer=1,
               f::Optional{Union{AbstractFloat,Integer}}=nothing,
               a::Optional{Union{AbstractFloat,Integer}}=nothing,
               north::Union{AbstractFloat,Integer}=90.0,
               south::Union{AbstractFloat,Integer}=-90.0,
               east::Union{AbstractFloat,Integer}=360.0,
               west::Union{AbstractFloat,Integer}=0.0,
               dealloc::Bool=false,
               exitstatus::Optional{Ref{<:Integer}}=nothing)
grid::Array{Cdouble,2}
nlat::Int
nlong::Int

Create a 2D cylindrical map with arbitrary grid spacing from a set of spherical harmonic coefficients.

See also: MakeGrid2d!, MakeGridPoint

source
SHTOOLS.MakeGridPointFunction
value = MakeGridPoint(cilm::AbstractArray{Cdouble,3},
                      lmax::Integer,
                      lat::Union{AbstractFloat,Integer},
                      long::Union{AbstractFloat,Integer};
                      norm::Integer=1,
                      csphase::Integer=1,
                      dealloc::Bool=false)
value::Float64

value = MakeGridPoint(cilm::AbstractArray{Complex{Cdouble},3},
                      lmax::Integer,
                      lat::Union{AbstractFloat,Integer},
                      long::Union{AbstractFloat,Integer};
                      norm::Integer=1,
                      csphase::Integer=1,
                      dealloc::Bool=false)
value::Complex{Float64}

Evaluate a real or complex function expressed in real or complex spherical harmonics at a single point.

See also: SHExpandLSQ!, SHExpandLSQ, MakeGrid2d!, MakeGrid2d

source
SHTOOLS.SHMultiply!Function
cilmout = SHMultiply!(cilmout::AbstractArray{Cdouble,3},
                      cilm1::AbstractArray{Cdouble,3},
                      lmax1::Integer,
                      cilm2::AbstractArray{Cdouble,3},
                      lmax2::Integer;
                      precomp::Bool=false,
                      norm::Integer=1,
                      csphase::Integer=1,
                      exitstatus::Optional{Ref{<:Integer}}=nothing)
cilmout::AbstractArray{Cdouble,3}

Multiply two functions and determine the spherical harmonic coefficients of the resulting function.

See also: SHMultiply, SHExpandDH!, SHExpandGLQ!, SHExpandLSQ!

source
SHTOOLS.SHMultiplyFunction
cilmout = SHMultiply(cilm1::AbstractArray{Cdouble,3},
                     lmax1::Integer,
                     cilm2::AbstractArray{Cdouble,3},
                     lmax2::Integer;
                     precomp::Bool=false,
                     norm::Integer=1,
                     csphase::Integer=1,
                     exitstatus::Optional{Ref{<:Integer}}=nothing)
cilmout::Array{Cdouble,3}

Multiply two functions and determine the spherical harmonic coefficients of the resulting function.

See also: SHMultiply!, SHExpandDH, SHExpandGLQ, SHExpandLSQ

source

Spherical harmonic I/O, storage, and conversions

Spherical harmonic storage

SHTOOLS.SHCilmToCindex!Function
cindex = SHCilmToCindex!(cindex::AbstractArray{Cdouble,2},
                         cilm::AbstractArray{Cdouble,3};
                         degmax::Optional{Cint}=nothing,
                         exitstatus::Optional{Ref{<:Integer}}=nothing)
cindex::AbstractArray{Cdouble,2}

Convert a three-dimensional array of spherical harmonic coefficients to a two-dimensional indexed array.

See also: SHCilmToCindex, SHCindexToCilm!, SHCilmToVector!

source
SHTOOLS.SHCilmToCindexFunction
cindex = SHCilmToCindex(cilm::AbstractArray{Cdouble,3};
                        degmax::Optional{Cint}=nothing,
                        exitstatus::Optional{Ref{<:Integer}}=nothing)
cindex::Array{Cdouble,2}

Convert a three-dimensional array of spherical harmonic coefficients to a two-dimensional indexed array.

See also: SHCilmToCindex!, SHCindexToCilm, SHCilmToVector

source
SHTOOLS.SHCindexToCilm!Function
cilm = SHCindexToCilm!(cilm::AbstractArray{Cdouble,3},
                       cindex::AbstractArray{Cdouble,2};
                       degmax::Optional{Cint}=nothing,
                       exitstatus::Optional{Ref{<:Integer}}=nothing)
cilm::AbstractArray{Cdouble,3},

Convert a two-dimensional indexed array of spherical harmonic coefficients to a three-dimensional array.

See also: SHCindexToCilm, SHCilmToCindex!, SHVectorToCilm!

source
SHTOOLS.SHCindexToCilmFunction
cilm = SHCindexToCilm(cindex::AbstractArray{Cdouble,2};
                      degmax::Optional{Cint}=nothing,
                      exitstatus::Optional{Ref{<:Integer}}=nothing)
cilm::AbstractArray{Cdouble,3},

Convert a two-dimensional indexed array of spherical harmonic coefficients to a three-dimensional array.

See also: SHCindexToCilm!, SHCilmToCindex, SHVectorToCilm

source
SHTOOLS.SHCilmToVector!Function
SHCilmToVector!(vector::AbstractVector{Cdouble},
                cilm::AbstractArray{Cdouble,3},
                lmax::Integer;
                exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)

Convert a three-dimensional array of real spherical harmonic coefficients to a 1-dimensional indexed vector.

See also: SHCilmToVector

source
SHTOOLS.SHCilmToVectorFunction
vector = SHCilmToVector(cilm::AbstractArray{Cdouble,3},
                        lmax::Integer;
                        exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
vector::Vector{Cdouble}

Convert a three-dimensional array of real spherical harmonic coefficients to a 1-dimensional indexed array.

See also: SHCilmToVector!

source
SHTOOLS.SHVectorToCilm!Function
SHVectorToCilm!(cilm::AbstractArray{Cdouble,3},
                vector::AbstractVector{Cdouble},
                lmax::Integer;
                exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)

Convert a 1-dimensional indexed vector of real spherical harmonic coefficients to a 3-dimensional array.

See also: SHVectorToCilm

source
SHTOOLS.SHVectorToCilmFunction
cilm = SHVectorToCilm(cilm::AbstractArray{Cdouble,3},
                      vector::AbstractVector{Cdouble},
                      lmax::Integer;
                      exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
cilm::Array{Cdouble,3}

Convert a 1-dimensional indexed vector of real spherical harmonic coefficients to a 3-dimensional array.

See also: SHVectorToCilm!

source

Spherical harmonic conversions

Miscellaneous

SHTOOLS.Wigner3jFunction
w3j, jmin, jmax = Wigner3j(j2::Integer, j3::Integer,
                            m1::Integer, m2::Integer, m3::Integer,
                            exitstatus::Union{Nothing,Ref{<:Integer}} = nothing)

Compute the Wigner 3j symbol

\[\left(\begin{array}{ccc} j & j_{2} & j_{3}\\ m_{1} & m_{2} & m_{3} \end{array}\right)\]

for all integral values of j that satisfy max(abs(m1), abs(j2 - j3)) <= j <= j2 + j3. The values of jmin and jmax are set to max(abs(m1), abs(j2 - j3)) and j2 + j3 respectively. The vector w3j contains the computed 3j symbols. Its first index corresponds to j = jmin, and last index corresponds to j = jmax.

See also: Wigner3j!

SHTOOLS page: Wigner3j

source
SHTOOLS.Wigner3j!Function
w3j, jmin, jmax = Wigner3j!(w3j::AbstractVector{Cdouble}, j2::Integer, j3::Integer,
                            m1::Integer, m2::Integer, m3::Integer,
                            exitstatus::Union{Nothing,Ref{<:Integer}} = nothing)

Compute the Wigner 3j symbol

\[\left(\begin{array}{ccc} j & j_{2} & j_{3}\\ m_{1} & m_{2} & m_{3} \end{array}\right)\]

for all integral values of j that satisfy max(abs(m1), abs(j2 - j3)) <= j <= j2 + j3. The values of jmin and jmax are set to max(abs(m1), abs(j2 - j3)) and j2 + j3 respectively. The vector w3j is filled with the computed 3j symbols. Its first index corresponds to j = jmin, and it has its first jnum = jmax - jmin + 1 indices populated.

Note

w3j must be sufficiently long to contain at least j2+j3+1 elements. Pre-existing values in w3j[jnum:j2+j3+1] may be overwritten with zero.

See also: Wigner3j

SHTOOLS page: Wigner3j

source