FINUFFT.jl Reference
For installation and basic usage, see the README at https://github.com/ludvigak/FINUFFT.jl
For documentation of the library functions that are being called, see the FINUFFT documentation at https://finufft.readthedocs.io
Index
FINUFFT.cufinufft_opts
FINUFFT.nufft_opts
FINUFFT._cufinufft_makeplan
FINUFFT._finufft_makeplan
FINUFFT.cufinufft_default_opts
FINUFFT.cufinufft_destroy!
FINUFFT.cufinufft_exec
FINUFFT.cufinufft_exec!
FINUFFT.cufinufft_makeplan
FINUFFT.cufinufft_setpts!
FINUFFT.finufft_default_opts
FINUFFT.finufft_destroy!
FINUFFT.finufft_exec
FINUFFT.finufft_exec!
FINUFFT.finufft_makeplan
FINUFFT.finufft_setpts!
FINUFFT.nufft1d1
FINUFFT.nufft1d1!
FINUFFT.nufft1d1!
FINUFFT.nufft1d2
FINUFFT.nufft1d2!
FINUFFT.nufft1d2!
FINUFFT.nufft1d3
FINUFFT.nufft1d3!
FINUFFT.nufft2d1
FINUFFT.nufft2d1!
FINUFFT.nufft2d1!
FINUFFT.nufft2d2
FINUFFT.nufft2d2!
FINUFFT.nufft2d2!
FINUFFT.nufft2d3
FINUFFT.nufft2d3!
FINUFFT.nufft3d1
FINUFFT.nufft3d1!
FINUFFT.nufft3d1!
FINUFFT.nufft3d2
FINUFFT.nufft3d2!
FINUFFT.nufft3d2!
FINUFFT.nufft3d3
FINUFFT.nufft3d3!
Types
FINUFFT.cufinufft_opts
— Typemutable struct cufinufft_opts
upsampfac :: Cdouble # upsampling ratio sigma, only 2.0 (standard) is implemented
# following options are for gpu #
gpu_method :: Cint # 1: nonuniform-pts driven, 2: shared mem (SM)
gpu_sort :: Cint # when NU-pts driven: 0: no sort (GM), 1: sort (GM-sort)
gpu_binsizex :: Cint # used for 2D, 3D subproblem method
gpu_binsizey :: Cint
gpu_binsizez :: Cint
gpu_obinsizex :: Cint # used for 3D spread block gather method
gpu_obinsizey :: Cint
gpu_obinsizez :: Cint
gpu_maxsubprobsize :: Cint
gpu_kerevalmeth :: Cint # 0: direct exp(sqrt()), 1: Horner ppval
gpu_spreadinterponly :: Cint # 0: NUFFT, 1: spread or interpolation only
gpu_maxbatchsize :: Cint
# multi-gpu support #
gpu_device_id :: Cint
gpu_stream :: Ptr{Cvoid}
modeord :: Cint # (type 1,2 only): 0 CMCL-style increasing mode order
# 1 FFT-style mode order
end
Options struct passed to cuFINUFFT, see C documentation.
FINUFFT.nufft_opts
— Typemutable struct nufft_opts
modeord :: Cint
chkbnds :: Cint
debug :: Cint
spread_debug :: Cint
showwarn :: Cint
nthreads :: Cint
fftw :: Cint
spread_sort :: Cint
spread_kerevalmeth :: Cint
spread_kerpad :: Cint
upsampfac :: Cdouble
spread_thread :: Cint
maxbatchsize :: Cint
spread_nthr_atomic :: Cint
spread_max_sp_size :: Cint
end
Options struct passed to the FINUFFT library.
Fields
This is a summary only; see FINUFFT documentation for full descriptions.
modeord :: Cint
0: CMCL-style increasing mode ordering (neg to pos), or
1: FFT-style mode ordering (affects type-1,2 only)
chkbnds :: Cint
0: don't check if input NU pts in [-3pi,3pi], 1: do
debug :: Cint
0: silent, 1: text basic timing output
spread_debug :: Cint
passed to spread_opts, 0 (no text) 1 (some) or 2 (lots)
showwarn :: Cint
Whether to print warnings to stderr. 0: silent, 1: print warnings
nthreads :: Cint
How many threads FINUFFT should use, or 0 (use max available in OMP)
fftw :: Cint
0:FFTW_ESTIMATE
, or 1:FFTW_MEASURE
(slow plan but faster FFTs)
spread_sort :: Cint
passed to spread_opts, 0 (don't sort) 1 (do) or 2 (heuristic)
spread_kerevalmeth :: Cint
passed to spread_opts, 0: exp(sqrt()), 1: Horner ppval (faster)
spread_kerpad :: Cint
passed to spread_opts, 0: don't pad to mult of 4, 1: do
upsampfac :: Cdouble
upsampling ratio sigma: 2.0 (standard), or 1.25 (small FFT), or
0.0 (auto).
spread_thread :: Cint
(for ntrans>1 only)
0: auto choice,
1: sequential multithreaded,
2: parallel singlethreaded spread.
maxbatchsize :: Cint
(for ntrans>1 only). max blocking size for vectorized, 0 for auto-set
spread_nthr_atomic :: Cint
if >=0, threads above which spreader OMP critical goes atomic
spread_max_sp_size :: Cint
if >0, overrides spreader (dir=1 only) max subproblem size
Functions
FINUFFT._cufinufft_makeplan
— Method_cufinufft_makeplan
Type-stable internal version of cufinufft_makeplan
FINUFFT._finufft_makeplan
— Method_finufft_makeplan
Type-stable internal version of finufft_makeplan
FINUFFT.cufinufft_default_opts
— Methodp = cufinufft_default_opts()
Return a FINUFFT.cufinufft_opts
struct with the default settings.
FINUFFT.cufinufft_destroy!
— Methodcufinufft_destroy!(plan::cufinufft_plan)
Destroy a cufinufft_plan
object, deallocating all memory used.
FINUFFT.cufinufft_exec!
— Methodcufinufft_exec!(plan::cufinufft_plan{T},
input::CUDA.CuArray{Complex{T}},
output::CUDA.CuArray{Complex{T}}
)
where T :: Float32 or Float64
Execute cuFINUFFT transform(s) with preallocated arrays on device.
FINUFFT.cufinufft_exec
— Methodcufinufft_exec(plan::cufinufft_plan{T},
input :: Array{Complex{T}} or CUDA.CuArray{Complex{T}}
) -> Array{Complex{T}} or CUDA.CuArray{Complex{T}}
where T :: Float32 or Float64
Execute cuFINFFT plan and return output in a newly allocated array.
output
type will match that of input
:
- If
input
isCUDA.CuArray
on device, thenoutput
is allocated on device. - If
input
isArray
on host, then it is copied to device before computation andoutput
is copied to host after computation.
FINUFFT.cufinufft_makeplan
— Methodcufinufft_makeplan(type::Integer,
n_modes_or_dim::Union{Array{Int64},Integer},
iflag::Integer,
ntrans::Integer,
eps::Real;
dtype=Float64,
kwargs...) -> plan::cufinufft_plan{dtype}
Create a cufinufft_plan
object. See finufft_makeplan
for arguments.
kwargs
(optional): Options set inFINUFFT.cufinufft_opts
.
FINUFFT.cufinufft_setpts!
— Methodcufinufft_setpts!(plan, xj [, yj[, zj[, s[, t[, u]]]]])
Input nonuniform points. See finufft_setpts!
for arguments.
Points can be either CUDA.CuArray
's on device or Array
's on host. The latter will be automatically copied to device before being passed to cuFINUFFT.
FINUFFT.finufft_default_opts
— Functionp = finufft_default_opts()
p = finufft_default_opts(dtype=Float32)
Return a nufft_opts
struct with the default FINUFFT settings. Set up the double precision variant by default.
See: https://finufft.readthedocs.io/en/latest/usage.html#options
FINUFFT.finufft_destroy!
— Methodstatus = finufft_destroy!(plan::finufft_plan{T}) where T <: finufftReal
This destroys a FINUFFT plan object: it deallocates all stored FFTW plans, nonuniform point sorting arrays, kernel Fourier transforms arrays, and any other allocations, and nulls the plan pointer.
An integer status code is returned, that is 0 if successful. If one attempts to destroy an already-destroyed plan, 1 is returned (see FINUFFT documentation for finufft_destroy).
FINUFFT.finufft_exec!
— Methodfinufft_exec!(plan::finufft_plan{T},
input::Array{Complex{T}},
output::Array{Complex{T}}) where T <: finufftReal
Execute single or many-vector FINUFFT transforms in a plan, with output written to preallocated array. See finufft_exec
for arguments.
FINUFFT.finufft_exec
— Methodoutput::Array{Complex{T}} = finufft_exec(plan::finufft_plan{T},
input::Array{Complex{T}}) where T <: finufftReal
Execute single or many-vector FINUFFT transforms in a plan.
output = finufft_exec(plan, input)
For plan
a previously created finufft_plan
object also containing all needed nonuniform point coordinates, do a single (or if ntrans>1
in the plan stage, multiple) NUFFT transform(s), with the strengths or Fourier coefficient inputs vector(s) from input
. The result of the transform(s) is returned as a (possibly multidimensional) array.
Inputs
- `plan` `finufft_plan` object, already planned and containing
nonuniform points.
- `input` strengths (types 1 or 3) or Fourier coefficients (type 2)
vector, matrix, or array of appropriate size. For type 1 and 3,
this is either a length-M vector (where M is the length of `xj`),
or an `(M,ntrans)` matrix when `ntrans>1`. For type 2, in 1D this is size `(ms,)`, in 2D size `(ms,mt)`, or in 3D size `(ms,mt,mu)`, or
each of these with an extra last dimension `ntrans` if `ntrans>1`.
Output
`output` vector of output strengths at targets (types 2 or 3), or array
of Fourier coefficients (type 1), or, if `ntrans>1`, a stack of
such vectors or arrays, of appropriate size.
Specifically, if `ntrans=1`, for type 1, in 1D
this is size `(ms,)`, in 2D size
`(ms,mt)`, or in 3D size `(ms,mt,mu)`; for types 2 and 3
it is a column vector of length `M` (the length of `xj` in type 2),
or `nk` (the length of `s` in type 3). If `ntrans>1` it is a stack
of such objects, ie, it has an extra last dimension `ntrans`.
FINUFFT.finufft_makeplan
— Methodfinufft_makeplan(type::Integer,
n_modes_or_dim::Union{Array{Int64},Integer},
iflag::Integer,
ntrans::Integer,
eps::Real;
dtype=Float64,
kwargs...) -> plan::finufft_plan{dtype}
Creates a finufft_plan
object for the guru interface to FINUFFT, of type 1, 2 or 3, and with given numbers of Fourier modes (unless type 3).
Inputs
type
transform type: 1, 2, or 3n_modes_or_dim
iftype
is 1 or 2, the number of Fourier modes in each dimension:ms
in 1D,[ms mt]
in 2D, or[ms mt mu]
in 3D. Its length thus sets the dimension, which must be 1, 2 or 3. Iftype
is 3, in contrast, its value fixes the dimension.iflag
if >=0, uses + sign in exponential, otherwise - sign.ntrans
number of transforms to compute simultaneouslyeps
relative precision requested (generally between 1e-15 and 1e-1), real, need not match type ofdtype
dtype
Float32
orFloat64
, plan for single precision or double precisionkwargs
(optional): for options, seenufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Returns
- finufft_plan struct
Examples
julia> p = finufft_makeplan(2,10,+1,1,1e-6);
creates a plan for a 1D type 2 Float64 transform with 10 Fourier modes and tolerance 1e-6.
julia> p = finufft_makeplan(1,[10 20],+1,1,1e-6);
creates a plan for a 2D type 1 Float64 transform with 10*20 Fourier modes.
julia> p = finufft_makeplan(3,2,+1,1,1e-4,dtype=Float32,nthreads=4);
creates a plan for a 2D type 3 Float32 transform with tolerance 1e-4, to use 4 threads.
FINUFFT.finufft_setpts!
— Methodfinufft_setpts!(plan, xj [, yj[, zj[, s[, t[, u]]]]])
Input nonuniform points for general FINUFFT transform(s).
Given an already-planned finufft_plan
, this reads in nonuniform point coordinate arrays xj
(and yj
if 2D or 3D, and zj
if 3D), and additionally in the type 3 case, nonuniform frequency target coordinate arrays s
(and t
if 2D or 3D, and u
if 3D). Empty arrays may be passed in the case of unused dimensions. For all types, sorting is done to internally store a reindexing of points, and for type 3 the spreading and FFTs are planned. These nonuniform points may then be used for multiple transforms.
Inputs
plan
afinufft_plan
object for one/many general nonuniform FFTsxj
Array{Float32} or Array{Float64}, vector of x-coords of all nonuniform pointsyj
empty (if dim<2), or vector of y-coords of all nonuniform pointszj
empty (if dim<3), or vector of z-coords of all nonuniform pointss
vector of x-coords of all nonuniform frequency targetst
empty (if dim<2), or vector of y-coords of all frequency targetsu
empty (if dim<3), or vector of z-coords of all frequency targets
FINUFFT.nufft1d1!
— Methodnufft1d1!(xj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-1 1D complex nonuniform FFT. Output written to fk
. See nufft1d1
.
FINUFFT.nufft1d1!
— Methodnufft1d1!(xj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft1d1
— Methodnufft1d1(xj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
ms :: Integer;
kwargs...
) -> Array{ComplexF64} or Array{ComplexF32}
Compute type-1 1D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f(k1) = SUM c[j] exp(+/-i k1 x(j)) for -ms/2 <= k1 <= (ms-1)/2
j=1
Inputs
xj
locations of nonuniform sources on interval [-3pi,3pi), length njcj
length-nj complex vector of source strengths. If length(cj)>nj, expects a stack of vectors (eg, a nj*ntrans matrix) each of which is transformed with the same source locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)ms
number of Fourier modes computed, may be even or odd; in either case, mode range is integers lying in [-ms/2, (ms-1)/2]- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- size
(ms,)
complex vector of Fourier coefficients f, or, ifntrans>1
, matrix of size(ms,ntrans)
.
FINUFFT.nufft1d2!
— Methodnufft1d2!(xj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-2 1D complex nonuniform FFT. Output written to cj
. See nufft1d2
.
FINUFFT.nufft1d2!
— Methodnufft1d2!(xj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft1d2
— Methodnufft1d2(xj :: Array{Float64} or Array{Float32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
) -> Array{ComplexF64}
Compute type-2 1D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
c[j] = SUM f[k1] exp(+/-i k1 x[j]) for j = 1,...,nj
k1
where sum is over -ms/2 <= k1 <= (ms-1)/2.
Input
xj
location of nonuniform targets on interval [-3pi,3pi), length njfk
complex Fourier coefficients. If a vector, length setsms
(with mode ordering given by opts.modeord). If a matrix, each column is transformed with the same nonuniform targets.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- complex
(nj,)
vector c of answers at targets, or, ifntrans>1
, matrix of size(nj,ntrans)
.
FINUFFT.nufft1d3!
— Methodnufft1d3!(xj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32},
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-3 1D complex nonuniform FFT. Output written to fk
. See nufft1d3
.
FINUFFT.nufft1d3
— Methodnufft1d3(xj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32};
kwargs...
) -> Array{ComplexF64}
Compute type-3 1D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f[k] = SUM c[j] exp(+-i s[k] x[j]), for k = 1, ..., nk
j=1
Inputs
xj
locations of nonuniform sources on R (real line), length-nj vector.cj
length-nj complex vector of source strengths. If length(cj)>nj, expects a size(nj,ntrans)
matrix each column of which is transformed with the same source and target locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)sk
frequency locations of nonuniform targets on R, length-nk vector.- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- complex vector f size '(nk,)
of values at targets, or, if
ntrans>1, a matrix of size
(nk,ntrans)`
FINUFFT.nufft2d1!
— Methodnufft2d1!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-1 2D complex nonuniform FFT. Output written to fk
. See nufft2d1
.
FINUFFT.nufft2d1!
— Methodnufft2d1!(xj :: CuArray{Float64} or CuArray{Float32},
yj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft2d1
— Methodnufft2d1(xj :: Array{Float64} or Array{Float32}
yj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
ms :: Integer,
mt :: Integer;
kwargs...
) -> Array{ComplexF64}
Compute type-1 2D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f[k1,k2] = SUM c[j] exp(+-i (k1 x[j] + k2 y[j]))
j=1
for -ms/2 <= k1 <= (ms-1)/2, -mt/2 <= k2 <= (mt-1)/2.
Inputs
xj
,yj
coordinates of nonuniform sources on the square [-3pi,3pi)^2, each a length-nj vectorcj
length-nj complex vector of source strengths. If length(cj)>nj, expects a stack of vectors (eg, a nj*ntrans matrix) each of which is transformed with the same source locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)ms
,mt
number of Fourier modes requested in x & y; each may be even or odd. In either case the mode range is integers lying in [-m/2, (m-1)/2]- kwargs (optional), see
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- size
(ms,mt)
complex matrix of Fourier coefficients f (ordering given by opts.modeord in each dimension;ms
fast,mt
slow), or, ifntrans>1
, a array of size(ms,mt,ntrans)
.
FINUFFT.nufft2d2!
— Methodnufft2d2!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-2 2D complex nonuniform FFT. Output written to cj
. See nufft2d2
.
FINUFFT.nufft2d2!
— Methodnufft2d2!(xj :: CuArray{Float64} or CuArray{Float32},
yj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft2d2
— Methodnufft2d2(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
) -> Array{ComplexF64}
Compute type-2 2D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
c[j] = SUM f[k1,k2] exp(+/-i (k1 x[j] + k2 y[j])) for j = 1,..,nj
k1,k2
where sum is over -ms/2 <= k1 <= (ms-1)/2, -mt/2 <= k2 <= (mt-1)/2,
Inputs
xj
,yj
coordinates of nonuniform targets on the square [-3pi,3pi)^2, each a vector of length njfk
complex Fourier coefficient matrix, whose size determines (ms,mt). (Mode ordering given by opts.modeord, in each dimension.) If a 3D array, 3rd dimension setsntrans
, and each ofntrans
matrices is transformed with the same nonuniform targets.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- complex size
(nj,)
vector c of answers at targets, or, ifntrans>1
, matrix of size(nj,ntrans)
.
FINUFFT.nufft2d3!
— Methodnufft2d3!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32},
tk :: Array{Float64} or Array{Float32},
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-3 2D complex nonuniform FFT. Output written to 'fk'. See nufft2d3
.
FINUFFT.nufft2d3
— Methodnufft2d3(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32},
tk :: Array{Float64} or Array{Float32};
kwargs...
) -> Array{ComplexF64}
Compute type-3 2D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f[k] = SUM c[j] exp(+-i (s[k] x[j] + t[k] y[j])), for k = 1, ..., nk
j=1
Inputs
xj
,yj
coordinates of nonuniform sources in R^2, each a length-nj vector.cj
complex vector(nj,)
of source strengths. If length(cj)>nj, expects a(nj,ntrans)
matrix, each column of which is transformed with the same source and target locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)sk
,tk
frequency coordinates of nonuniform targets in R^2, each a length-nk vector.- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- complex vector size
(nk,)
of values at targets, or, ifntrans>1
, a matrix of size(nk,ntrans)
FINUFFT.nufft3d1!
— Methodnufft3d1!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-1 3D complex nonuniform FFT. Output written to fk
. See nufft3d1
.
FINUFFT.nufft3d1!
— Methodnufft3d1!(xj :: CuArray{Float64} or CuArray{Float32},
yj :: CuArray{Float64} or CuArray{Float32},
zj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft3d1
— Methodnufft3d1(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
ms :: Integer,
mt :: Integer,
mu :: Integer;
kwargs...
) -> Array{ComplexF64}
Compute type-1 3D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f[k1,k2,k3] = SUM c[j] exp(+-i (k1 x[j] + k2 y[j] + k3 z[j]))
j=1
for -ms/2 <= k1 <= (ms-1)/2, -mt/2 <= k2 <= (mt-1)/2,
-mu/2 <= k3 <= (mu-1)/2.
Inputs
xj
,yj
,zj
coordinates of nonuniform sources on the cube [-3pi,3pi)^3, each a length-nj vectorcj
length-nj complex vector of source strengths. If length(cj)>nj, expects a stack of vectors (eg, a nj*ntrans matrix) each of which is transformed with the same source locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)ms
,mt
,mu
number of Fourier modes requested in x,y and z; each may be even or odd. In either case the mode range is integers lying in [-m/2, (m-1)/2]- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- size
(ms,mt,mu)
complex array of Fourier coefficients f (ordering given by opts.modeord in each dimension;ms
fastest,mu
slowest), or, ifntrans>1
, a 4D array of size(ms,mt,mu,ntrans)
.
FINUFFT.nufft3d2!
— Methodnufft3d2!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-2 3D complex nonuniform FFT. Output written to cj
. See nufft3d2
.
FINUFFT.nufft3d2!
— Methodnufft3d2!(xj :: CuArray{Float64} or CuArray{Float32},
yj :: CuArray{Float64} or CuArray{Float32},
zj :: CuArray{Float64} or CuArray{Float32},
cj :: CuArray{ComplexF64} or CuArray{ComplexF32},
iflag :: Integer,
eps :: Real,
fk :: CuArray{ComplexF64} or CuArray{ComplexF32};
kwargs...
)
CUDA version.
FINUFFT.nufft3d2
— Methodnufft3d2(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
iflag :: Integer,
eps :: Real,
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
) -> Array{ComplexF64}
Compute type-2 3D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
c[j] = SUM f[k1,k2,k3] exp(+/-i (k1 x[j] + k2 y[j] + k3 z[j]))
k1,k2,k3
for j = 1,..,nj
where sum is over -ms/2 <= k1 <= (ms-1)/2, -mt/2 <= k2 <= (mt-1)/2,
-mu/2 <= k3 <= (mu-1)/2.
Inputs
xj
,yj
,zj
coordinates of nonuniform targets on the cube [-3pi,3pi)^3, each a vector of length njfk
complex Fourier coefficient array, whose size sets(ms,mt,mu)
. (Mode ordering given by opts.modeord, in each dimension.) If a 4D array, 4th dimension setsntrans
, and each ofntrans
3D arrays is transformed with the same nonuniform targets.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- complex vector c of size
(nj,)
giving answers at targets, or, ifntrans>1
, matrix of size(nj,ntrans)
.
FINUFFT.nufft3d3!
— Methodnufft3d3!(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32},
tk :: Array{Float64} or Array{Float32},
uk :: Array{Float64} or Array{Float32},
fk :: Array{ComplexF64} or Array{ComplexF32};
kwargs...
)
Compute type-3 3D complex nonuniform FFT. Output written to fk
. See nufft3d3
.
FINUFFT.nufft3d3
— Methodnufft3d3(xj :: Array{Float64} or Array{Float32},
yj :: Array{Float64} or Array{Float32},
zj :: Array{Float64} or Array{Float32},
cj :: Array{ComplexF64} or Array{ComplexF32},
iflag :: Integer,
eps :: Real,
sk :: Array{Float64} or Array{Float32},
tk :: Array{Float64} or Array{Float32},
uk :: Array{Float64} or Array{Float32};
kwargs...
) -> Array{ComplexF64}
Compute type-3 3D complex nonuniform FFT. This computes, to relative precision eps, via a fast algorithm:
nj
f[k] = SUM c[j] exp(+-i (s[k] x[j] + t[k] y[j] + u[k] z[j])),
j=1
for k = 1, ..., nk
Inputs
xj
,yj
,zj
coordinates of nonuniform sources in R^3, each a length-nj vector.cj
complex(nj,)
vector of source strengths. If length(cj)>nj, expects a '(nj,ntrans)' matrix, each column of which is transformed with the same source and target locations.iflag
if >=0, uses + sign in exponential, otherwise - sign.eps
relative precision requested (generally between 1e-15 and 1e-1)sk
,tk,
uk` frequency coordinates of nonuniform targets in R^3, each a length-nk vector.- kwargs (optional). See
nufft_opts
and https://finufft.readthedocs.io/en/latest/opts.html
Output
- size
(nk,)
complex vector f values at targets, or, ifntrans>1
, a matrix of size(nk,ntrans)