Skip to content

Commit

Permalink
Update docs
Browse files Browse the repository at this point in the history
  • Loading branch information
JeffFessler committed Dec 4, 2022
1 parent 13b2fda commit 3b607af
Show file tree
Hide file tree
Showing 2 changed files with 44 additions and 4 deletions.
11 changes: 10 additions & 1 deletion docs/lit/examples/1-mirt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,4 +2,13 @@
# # [MIRT overview](@id 1-mirt)
#---------------------------------------------------------

# This is a placeholder
#=
The Julia package
[`MIRT.jl`](https://github.com/JeffFessler/MIRTjim.jl)
provides tools for performing image reconstruction
and solving related inverse problems.
For more tools and more examples,
see
[JuliaImageRecon](https://github.com/JuliaImageRecon)
=#
37 changes: 34 additions & 3 deletions docs/lit/examples/2-nufft.jl
Original file line number Diff line number Diff line change
Expand Up @@ -3,8 +3,9 @@
#---------------------------------------------------------

#=
Examples illustrating the `nufft` option in the Julia package
[`MIRTjim`](https://github.com/JeffFessler/MIRTjim.jl).
Examples illustrating the `nufft` methods
in the Julia package
[`MIRT`](https://github.com/JeffFessler/MIRT.jl).
The `nufft` functions in this package
are wrappers around
Expand Down Expand Up @@ -38,7 +39,7 @@ This page was generated from a single Julia file:

using Plots; default(markerstrokecolor = :auto, label="")
using MIRTjim: prompt
using MIRT: nufft_errors
using MIRT: nufft_init, nufft_errors, dtft
using InteractiveUtils: versioninfo

# The following line is helpful when running this file as a script;
Expand All @@ -47,6 +48,36 @@ using InteractiveUtils: versioninfo
isinteractive() && prompt(:prompt);


#=
## 1D example
This code illustrates how to call the NUFFT.
=#
Ω = range(-1,1,301) * π # digital frequencies (radians/sample)
N = 16 # # of samples
(nufft, nufft_adjoint, Anufft) = nufft_init(Ω, N)
x = randn(ComplexF32, N) # random 1D signal
y1 = nufft(x) # one way to compute NUFFT
y2 = Anufft * x # another way to compute NUFFT
@assert y1 == y2
yd = dtft(Ω, x) # exact (slow) nonuniform discrete Fourier transform
maximum(abs, yd - y1) # worst error for this 1D signal

# This plot shows that the NUFFT is very close to the exact nonuniform DFT.
xticks = ([-π,0,π], ["", "0", "π"]); xlims = (-π,π); xlabel="Ω"
p1 = plot(ylabel="Real part"; xlabel, xlims, xticks)
plot!(Ω, real(y1), label="real(NUFFT)", line=:dash)
plot!(Ω, real(yd), label="real(DTFT)", line=:dot)
p2 = plot(ylabel="Imag part"; xlabel, xlims, xticks)
plot!(Ω, imag(y1), label="imag(NUFFT)", line=:dash)
plot!(Ω, imag(yd), label="imag(NUFFT)", line=:dot)
p3 = plot(Ω, real(yd - y1), label="real(Error)"; xlabel, xlims, xticks)
p4 = plot(Ω, imag(yd - y1), label="imag(Error)"; xlabel, xlims, xticks)
plot(p1, p2, p3, p4)

#
prompt()

#=
## Plot worst-case errors vs ``N``
Expand Down

0 comments on commit 3b607af

Please sign in to comment.