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
SHTOOLS.PlmBar!
— FunctionPlmBar!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the 4π-normalized associated Legendre functions.
See also: PlmBar
, PlmBar_d1!
, PlBar!
, PlmON!
, PlmSchmidt!
, PLegendreA!
, PlmIndex
SHTOOLS.PlmBar
— Functionp = PlmBar(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the 4π-normalized associated Legendre functions.
See also: PlmBar!
, PlmBar_d1
, PlBar
, PlmON
, PlmSchmidt
, PLegendreA
, PlmIndex
SHTOOLS.PlmBar_d1!
— FunctionPlmBar_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the 4π-normalized associated Legendre functions and first derivatives.
See also: PlmBar_d1
, PlmBar_d1
, PlBar_d1!
, PlmON_d1!
, PlmSchmidt_d1!
, PLegendreA_d1!
, PlmIndex
SHTOOLS.PlmBar_d1
— Functionp, dp = PlmBar_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the 4π-normalized associated Legendre functions and first derivatives.
See also: PlmBar
, PlmBar_d1!
, PlBar_d1
, PlmON_d1
, PlmSchmidt_d1
, PLegendreA_d1
, PlmIndex
SHTOOLS.PlBar!
— FunctionPlBar!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the 4π-normalized Legendre polynomials.
See also: PlBar
, PlBar_d1!
, PlmBar!
, PlON!
, PlSchmidt!
, PLegendre!
SHTOOLS.PlBar
— Functionp = PlBar(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the 4π-normalized associated Legendre functions.
See also: PlBar!
, PlBar_d1
, PlmBar
, PlON
, PlSchmidt
, PLegendre
SHTOOLS.PlBar_d1!
— FunctionPlBar_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the 4π-normalized associated Legendre functions and first derivatives.
See also: PlBar_d1
, PlBar_d1
, PlmBar_d1!
, PlON_d1!
, PlSchmidt_d1!
, PLegendre_d1!
SHTOOLS.PlBar_d1
— Functionp, dp = PlBar_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the 4π-normalized associated Legendre functions and first derivatives.
See also: PlBar
, PlBar_d1!
, PlmBar_d1
, PlON_d1
, PlSchmidt_d1
, PLegendre_d1
Orthonormalized
SHTOOLS.PlmON!
— FunctionPlmON!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the orthonormalized associated Legendre functions.
See also: PlmON
, PlmON_d1!
, PlON!
, PlmBar!
, PlmSchmidt!
, PLegendreA!
, PlmIndex
SHTOOLS.PlmON
— Functionp = PlmON(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the orthonormalized associated Legendre functions.
See also: PlmON!
, PlmON_d1
, PlON
, PlmBar
, PlmSchmidt
, PLegendreA
, PlmIndex
SHTOOLS.PlmON_d1!
— FunctionPlmON_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the orthonormalized associated Legendre functions and first derivatives.
See also: PlmON_d1
, PlmON_d1
, PlON_d1!
, PlmBar_d1!
, PlmSchmidt_d1!
, PLegendreA_d1!
, PlmIndex
SHTOOLS.PlmON_d1
— Functionp, dp = PlmON_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the orthonormalized associated Legendre functions and first derivatives.
See also: PlmON
, PlmON_d1!
, PlON_d1
, PlmBar_d1
, PlmSchmidt_d1
, PLegendreA_d1
, PlmIndex
SHTOOLS.PlON!
— FunctionPlON!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the orthonormalized associated Legendre functions.
See also: PlON
, PlON_d1!
, PlmON!
, PlBar!
, PlSchmidt!
, PLegendre!
SHTOOLS.PlON
— Functionp = PlON(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the orthonormalized associated Legendre functions.
See also: PlON!
, PlON_d1
, PlmON
, PlBar
, PlSchmidt
, PLegendre
SHTOOLS.PlON_d1!
— FunctionPlON_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the orthonormalized associated Legendre functions and first derivatives.
See also: PlON_d1
, PlON_d1
, PlmON_d1!
, PlBar_d1!
, PlSchmidt_d1!
, PLegendre_d1!
SHTOOLS.PlON_d1
— Functionp, dp = PlON_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the orthonormalized associated Legendre functions and first derivatives.
See also: PlON
, PlON_d1!
, PlmON_d1
, PlBar_d1
, PlSchmidt_d1
, PLegendre_d1
Schmidt semi-normalized
SHTOOLS.PlmSchmidt!
— FunctionPlmSchmidt!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the Schmidt-normalized associated Legendre functions.
See also: PlmSchmidt
, PlmSchmidt_d1!
, PlSchmidt!
, PlmBar!
, PlmON!
, PLegendreA!
, PlmIndex
SHTOOLS.PlmSchmidt
— Functionp = PlmSchmidt(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the Schmidt-normalized associated Legendre functions.
See also: PlmSchmidt!
, PlmSchmidt_d1
, PlSchmidt
, PlmBar
, PlmON
, PLegendreA
, PlmIndex
SHTOOLS.PlmSchmidt_d1!
— FunctionPlmSchmidt_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the Schmidt-normalized associated Legendre functions and first derivatives.
See also: PlmSchmidt_d1
, PlmSchmidt_d1
, PlSchmidt_d1!
, PlmBar_d1!
, PlmON_d1!
, PLegendreA_d1!
, PlmIndex
SHTOOLS.PlmSchmidt_d1
— Functionp, dp = PlmSchmidt_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
cnorm::Integer=0,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the Schmidt-normalized associated Legendre functions and first derivatives.
See also: PlmSchmidt
, PlmSchmidt_d1!
, PlSchmidt_d1
, PlmBar_d1
, PlmON_d1
, PLegendreA_d1
, PlmIndex
SHTOOLS.PlSchmidt!
— FunctionPlSchmidt!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the Schmidt-normalized associated Legendre functions.
See also: PlSchmidt
, PlSchmidt_d1!
, PlmSchmidt!
, PlBar!
, PlON!
, PLegendre!
SHTOOLS.PlSchmidt
— Functionp = PlSchmidt(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the Schmidt-normalized associated Legendre functions.
See also: PlSchmidt!
, PlSchmidt_d1
, PlmSchmidt
, PlBar
, PlON
, PLegendre
SHTOOLS.PlSchmidt_d1!
— FunctionPlSchmidt_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the Schmidt-normalized associated Legendre functions and first derivatives.
See also: PlSchmidt_d1
, PlSchmidt_d1
, PlmSchmidt_d1!
, PlBar_d1!
, PlON_d1!
, PLegendre_d1!
SHTOOLS.PlSchmidt_d1
— Functionp, dp = PlSchmidt_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the Schmidt-normalized associated Legendre functions and first derivatives.
See also: PlSchmidt
, PlSchmidt_d1!
, PlmSchmidt_d1
, PlBar_d1
, PlON_d1
, PLegendre_d1
Unnormalized
SHTOOLS.PLegendreA!
— FunctionPLegendreA!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the unnormalized associated Legendre functions.
See also: PLegendreA
, PLegendreA_d1!
, PLegendre!
, PlmBar!
, PlmON!
, PlmSchmidt!
, PlmIndex
SHTOOLS.PLegendreA
— Functionp = PLegendreA(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the unnormalized associated Legendre functions.
See also: PLegendreA!
, PLegendreA_d1
, PLegendre
, PlmBar
, PlmON
, PlmSchmidt
, PlmIndex
SHTOOLS.PLegendreA_d1!
— FunctionPLegendreA_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the unnormalized associated Legendre functions and first derivatives.
See also: PLegendreA_d1
, PLegendreA_d1
, PLegendre_d1!
, PlmBar_d1!
, PlmON_d1!
, PlmSchmidt_d1!
, PlmIndex
SHTOOLS.PLegendreA_d1
— Functionp, dp = PLegendreA_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
csphase::Integer=1,
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the unnormalized associated Legendre functions and first derivatives.
See also: PLegendreA
, PLegendreA_d1!
, PLegendre_d1
, PlmBar_d1
, PlmON_d1
, PlmSchmidt_d1
, PlmIndex
SHTOOLS.PLegendre!
— FunctionPLegendre!(p::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the unnormalized associated Legendre functions.
See also: PLegendre
, PLegendre_d1!
, PLegendreA!
, PlBar!
, PlON!
, PlSchmidt!
SHTOOLS.PLegendre
— Functionp = PLegendre(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
Compute all the unnormalized associated Legendre functions.
See also: PLegendre!
, PLegendre_d1
, PLegendreA
, PlBar
, PlON
, PlSchmidt
SHTOOLS.PLegendre_d1!
— FunctionPLegendre_d1!(p::AbstractVector{Cdouble},
dp::AbstractVector{Cdouble},
lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
Compute all the unnormalized associated Legendre functions and first derivatives.
See also: PLegendre_d1
, PLegendre_d1
, PLegendreA_d1!
, PlBar_d1!
, PlON_d1!
, PlSchmidt_d1!
SHTOOLS.PLegendre_d1
— Functionp, dp = PLegendre_d1(lmax::Integer,
z::Union{AbstractFloat,Integer};
exitstatus::Union{Nothing,Ref{<:Integer}}=nothing)
p::Vector{Cdouble}
dp::Vector{Cdouble}
Compute all the unnormalized associated Legendre functions and first derivatives.
See also: PLegendre
, PLegendre_d1!
, PLegendre_d1!
, PlBar_d1
, PlON_d1
, PlSchmidt_d1
Utilities
SHTOOLS.PlmIndex
— Functionindex = PlmIndex(l::Integer,
m::Integer)
index::Int
Compute the index of an array of Legendre functions corresponding to degree l
and angular order m
.
Spherical harmonic transforms
Equally sampled (N×N) and equally spaced (N×2N) grids
SHTOOLS.SHExpandDH!
— Functioncilm, 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!
SHTOOLS.SHExpandDH
— Functioncilm, 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
SHTOOLS.MakeGridDH!
— FunctionMakeGridDH!(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!
SHTOOLS.MakeGridDH
— FunctionMakeGridDH(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
SHTOOLS.SHExpandDHC!
— FunctionSHExpandDHC!(...)
Alias for SHExpandDH!
SHTOOLS.SHExpandDHC
— FunctionSHExpandDHC(...)
Alias for SHExpandDH
SHTOOLS.MakeGridDHC!
— FunctionMakeGridDHC!(...)
Alias for MakeGridDH!
SHTOOLS.MakeGridDHC
— FunctionMakeGridDHC(...)
Alias for MakeGridDH
SHTOOLS.MakeGradientDH!
— FunctionMakeGradientDH!(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!
SHTOOLS.MakeGradientDH
— Functiontheta, 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
Gauss-Legendre quadrature grids
SHTOOLS.SHGLQ!
— Functionzero, 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!
SHTOOLS.SHGLQ
— Functionzero, 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
SHTOOLS.SHExpandGLQ!
— Functioncilm = 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!
SHTOOLS.SHExpandGLQ
— Functioncilm = 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
SHTOOLS.MakeGridGLQ!
— Functiongridglq = 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!
SHTOOLS.MakeGridGLQ
— Functiongridglq = 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
SHTOOLS.SHExpandGLQC!
— FunctionSHExpandGLQC!(...)
Alias for SHExpandGLQ!
SHTOOLS.SHExpandGLQC
— FunctionSHExpandGLQC(...)
Alias for SHExpandGLQ
SHTOOLS.MakeGridGLQC!
— FunctionMakeGridGLQC!(...)
Alias for SHExpandGLQ!
SHTOOLS.MakeGridGLQC
— FunctionMakeGridGLQC(...)
Alias for SHExpandGLQ
SHTOOLS.GLQGridCoord!
— Functionlatglq, 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!
SHTOOLS.GLQGridCoord
— Functionlatglq, 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
Other routines
SHTOOLS.SHExpandLSQ!
— Functioncilm, 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!
SHTOOLS.SHExpandLSQ
— Functioncilm, 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
SHTOOLS.MakeGrid2d!
— Functiongrid, 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
SHTOOLS.MakeGrid2d
— Functiongrid, 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
SHTOOLS.MakeGridPoint
— Functionvalue = 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
SHTOOLS.MakeGridPointC
— FunctionMakeGridPointC(...)
Alias for MakeGridPoint
.
SHTOOLS.SHMultiply!
— Functioncilmout = 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!
SHTOOLS.SHMultiply
— Functioncilmout = 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
Spherical harmonic I/O, storage, and conversions
Spherical harmonic storage
SHTOOLS.SHCilmToCindex!
— Functioncindex = 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!
SHTOOLS.SHCilmToCindex
— Functioncindex = 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
SHTOOLS.SHCindexToCilm!
— Functioncilm = 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!
SHTOOLS.SHCindexToCilm
— Functioncilm = 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
SHTOOLS.SHCilmToVector!
— FunctionSHCilmToVector!(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
SHTOOLS.SHCilmToVector
— Functionvector = 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!
SHTOOLS.SHVectorToCilm!
— FunctionSHVectorToCilm!(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
SHTOOLS.SHVectorToCilm
— Functioncilm = 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!
SHTOOLS.YlmIndexVector
— Functionindex = YlmIndexVector(i, l, m)
Determine the index of a 1D ordered vector of spherical harmonic coefficients corresponding to i
, l
, and m
.
See also: SHCilmToVector!
, SHCilmToVector
, SHVectorToCilm!
, SHVectorToCilm
Spherical harmonic conversions
Miscellaneous
SHTOOLS.Wigner3j
— Functionw3j, 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
SHTOOLS.Wigner3j!
— Functionw3j, 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.
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