Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Pass to plan_fft all optional keywords of NFFTPlan #25

Merged
merged 1 commit into from
Apr 20, 2017
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
52 changes: 30 additions & 22 deletions src/NFFT.jl
Original file line number Diff line number Diff line change
Expand Up @@ -52,15 +52,18 @@ end

@inline dim{D,DIM}(::NFFTPlan{D,DIM}) = DIM

@doc """
"""
NFFTPlan(x, N, ...) -> plan
Compute `D` dimensional NFFT plan for sampling locations `x` (a vector or a `D`-by-`M` matrix) that can be applied on arrays of size `N` (a tuple of length `D`).
The optional arguments control the accuracy.
"""->
It takes as optional keywords all the keywords supported by `plan_fft` function (like
`flags` and `timelimit`). See documentation of `plan_fft` for reference.
"""
function NFFTPlan{D,T}(x::AbstractMatrix{T}, N::NTuple{D,Int}, m=4, sigma=2.0,
window=:kaiser_bessel, K=2000)
window=:kaiser_bessel, K=2000; kwargs...)
if !isa(x, Matrix)
x = collect(x)
end
Expand All @@ -75,8 +78,8 @@ function NFFTPlan{D,T}(x::AbstractMatrix{T}, N::NTuple{D,Int}, m=4, sigma=2.0,

M = size(x,2)

FP = plan_fft!(tmpVec)
BP = plan_bfft!(tmpVec)
FP = plan_fft!(tmpVec; kwargs...)
BP = plan_bfft!(tmpVec; kwargs...)

# Create lookup table
win, win_hat = getWindow(window)
Expand All @@ -102,20 +105,24 @@ function NFFTPlan{D,T}(x::AbstractMatrix{T}, N::NTuple{D,Int}, m=4, sigma=2.0,
NFFTPlan{D,0,T}(N, M, x, m, sigma, n, K, windowLUT, windowHatInvLUT, FP, BP, tmpVec )
end

function NFFTPlan(x::AbstractVector, N::Integer, m=4, sigma=2.0, window=:kaiser_bessel, K=2000)
NFFTPlan(reshape(x,1,length(x)), (N,), m, sigma, window, K)
function NFFTPlan(x::AbstractVector, N::Integer, m=4, sigma=2.0, window=:kaiser_bessel,
K=2000; kwargs...)
NFFTPlan(reshape(x,1,length(x)), (N,), m, sigma, window, K; kwargs...)
end


# Directional NFFT
@doc """
"""
NFFTPlan(x, d, N, ...) -> plan
Compute *directional* NFFT plan:
A 1D plan that is applied along dimension `d` of a `D` dimensional array of size `N` with sampling locations `x` (a vector).
"""->
It takes as optional keywords all the keywords supported by `plan_fft` function (like
`flags` and `timelimit`). See documentation of `plan_fft` for reference.
"""
function NFFTPlan{D,T}(x::AbstractVector{T}, dim::Integer, N::NTuple{D,Int64}, m=4,
sigma=2.0, window=:kaiser_bessel, K=2000)
sigma=2.0, window=:kaiser_bessel, K=2000; kwargs...)
n = ntuple(d->round(Int, sigma*N[d]), D)

sz = [N...]
Expand All @@ -124,8 +131,8 @@ function NFFTPlan{D,T}(x::AbstractVector{T}, dim::Integer, N::NTuple{D,Int64}, m

M = length(x)

FP = plan_fft!(tmpVec, dim)
BP = plan_bfft!(tmpVec, dim)
FP = plan_fft!(tmpVec, dim; kwargs...)
BP = plan_bfft!(tmpVec, dim; kwargs...)

# Create lookup table
win, win_hat = getWindow(window)
Expand All @@ -147,12 +154,13 @@ function NFFTPlan{D,T}(x::AbstractVector{T}, dim::Integer, N::NTuple{D,Int64}, m
NFFTPlan{D,dim,T}(N, M, reshape(x,1,M), m, sigma, n, K, windowLUT, windowHatInvLUT, FP, BP, tmpVec)
end

function NFFTPlan{D,T}(x::Matrix{T}, dim::Integer, N::NTuple{D,Int}, m=4, sigma=2.0, window=:kaiser_bessel, K=2000)
function NFFTPlan{D,T}(x::Matrix{T}, dim::Integer, N::NTuple{D,Int}, m=4, sigma=2.0,
window=:kaiser_bessel, K=2000; kwargs...)
if size(x,1) != 1 && size(x,2) != 1
throw(DimensionMismatch())
end

NFFTPlan(vec(x), dim, N, m, sigma, window, K)
NFFTPlan(vec(x), dim, N, m, sigma, window, K; kwargs...)
end


Expand Down Expand Up @@ -182,13 +190,13 @@ end

### nfft functions ###

@doc """
"""
nfft!(p, f, fHat) -> fHat
Calculate the NFFT of `f` with plan `p` and store the result in `fHat`.
Both `f` and `fHat` must be complex arrays.
"""->
"""
function nfft!{T}(p::NFFTPlan, f::AbstractArray{T}, fHat::StridedArray{T})
consistencyCheck(p, f, fHat)

Expand All @@ -203,7 +211,7 @@ function nfft!{T}(p::NFFTPlan, f::AbstractArray{T}, fHat::StridedArray{T})
return fHat
end

@doc """
"""
nfft(p, f) -> fHat
For a **non**-directional `D` dimensional plan `p` this calculates the NFFT of a `D` dimensional array `f` of size `N`.
Expand All @@ -213,7 +221,7 @@ For a **non**-directional `D` dimensional plan `p` this calculates the NFFT of a
For a **directional** `D` dimensional plan `p` both `f` and `fHat` are `D`
dimensional arrays, and the dimension specified in the plan creation is
affected.
"""->
"""
function nfft{D,T}(p::NFFTPlan{D,0}, f::AbstractArray{T,D})
fHat = zeros(T, p.M)
nfft!(p, f, fHat)
Expand All @@ -234,13 +242,13 @@ function nfft{D,DIM,T}(p::NFFTPlan{D,DIM}, f::AbstractArray{T,D})
end


@doc """
"""
nfft_adjoint!(p, fHat, f) -> f
Calculate the adjoint NFFT of `fHat` and store the result in `f`.
Both `f` and `fHat` must be complex arrays.
"""->
"""
function nfft_adjoint!(p::NFFTPlan, fHat::AbstractArray, f::StridedArray)
consistencyCheck(p, f, fHat)

Expand All @@ -254,7 +262,7 @@ function nfft_adjoint!(p::NFFTPlan, fHat::AbstractArray, f::StridedArray)
return f
end

@doc """
"""
nfft_adjoint(p, f) -> fHat
For a **non**-directional `D` dimensional plan `p` this calculates the adjoint NFFT of a length `M` vector `fHat`
Expand All @@ -264,7 +272,7 @@ For a **non**-directional `D` dimensional plan `p` this calculates the adjoint N
For a **directional** `D` dimensional plan `p` both `f` and `fHat` are `D`
dimensional arrays, and the dimension specified in the plan creation is
affected.
"""->
"""
function nfft_adjoint{D,DIM,T}(p::NFFTPlan{D,DIM}, fHat::AbstractArray{T})
f = Array{T}(p.N)
nfft_adjoint!(p, fHat, f)
Expand Down
4 changes: 2 additions & 2 deletions test/test.jl
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@ sigma = 2.0

M = prod(N)
x = rand(D,M) - 0.5
p = NFFTPlan(x, N, m, sigma, window)
p = NFFTPlan(x, N, m, sigma, window, flags = FFTW.ESTIMATE)

fHat = rand(M) + rand(M)*im
f = ndft_adjoint(p, fHat)
Expand All @@ -34,7 +34,7 @@ end
@testset "Abstract sampling points" begin
M, N = rand(100:200, 2)
x = linspace(-0.4, 0.4, M)
p = NFFTPlan(x, N)
p = NFFTPlan(x, N, flags = FFTW.MEASURE)
end

@testset "Directional NFFT $D dim" for D in 2:3 begin
Expand Down