Edit on GitHub

Extract k-space

Description

This example described how to extract kspace from cartesian datasets, undersample the data and create back and AcquisitionData structure that can be used for reconstruction.

This possibility can be easily combined with the BartIO package to reconstruct the kspace.

Setup

using CairoMakie
using ImageUtils: shepp_logan
using MRIReco, MRISimulation
using InteractiveUtils: versioninfo

function plot_im2D(im2D;title::String="")
    f = Figure()
    ax = Axis(f[1, 1],aspect = DataAspect(), yreversed = true, title = title)
    image!(ax, im2D')
    hidedecorations!(ax, grid = false)
    f
end
plot_im2D (generic function with 1 method)

Let's create a non-square AcquisitionData structure and perform a standard reconstruction

N = 128
N2 = 96
96

image

x = shepp_logan(N)
128×128 Matrix{Float64}:
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 ⋮                        ⋮              ⋱                      ⋮         
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  …  0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0
 0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0     0.0  0.0  0.0  0.0  0.0  0.0  0.0

simulation

params = Dict{Symbol, Any}()
params[:simulation] = "fast"
params[:trajName] = "Cartesian"
params[:numSamplingPerProfile] = N
params[:numProfiles] = floor(Int64, N2)

acqData = simulation(x[1:N,1:N2], params)

params[:reco] = "direct"
params[:reconSize] = (N,N2)

x_approx = reconstruction(acqData, params);

let's show the results for first bin

plot_im2D(abs.(x_approx[:,:,1,1,1,1]),title = "Reco from AcquisitionData")

We can extract the cartesian kspace (works for 2D or 3D) with the function kDataCart

kspace = kDataCart(acqData)
size(kspace)
(128, 96, 1, 1, 1, 1)

Dimensions of the kspace are : kx, ky, kz, Channels, Echoes, Repetitions

Let's see a standard reconstruction with the ifft function works :

x_approx2 = ifftshift(ifft(ifftshift(kspace)))
plot_im2D(abs.(x_approx2[:,:,1,1,1,1]),title = "Reco from extracted kspace")

We can create create and apply a mask on the kspace

mask = ones(eltype(kspace),size(kspace))
mask[1:2:end,:,:,:,:,:] .= 0
acq_u = AcquisitionData(kspace .* mask)
params[:reco] = "direct"
params[:reconSize] = (N,N2)

x_u = reconstruction(acq_u, params)
plot_im2D(abs.(x_u[:,:,1,1,1,1]),title = "Reco from undersampled AcquisitionData")

Reproducibility

This page was generated with the following version of Julia:

using InteractiveUtils
io = IOBuffer();
versioninfo(io);
split(String(take!(io)), '\n')
12-element Vector{SubString{String}}:
 "Julia Version 1.10.2"
 "Commit bd47eca2c8a (2024-03-01 10:14 UTC)"
 "Build Info:"
 "  Official https://julialang.org/ release"
 "Platform Info:"
 "  OS: Linux (x86_64-linux-gnu)"
 "  CPU: 4 × AMD EPYC 7763 64-Core Processor"
 "  WORD_SIZE: 64"
 "  LIBM: libopenlibm"
 "  LLVM: libLLVM-15.0.7 (ORCJIT, znver3)"
 "Threads: 1 default, 0 interactive, 1 GC (on 4 virtual cores)"
 ""

And with the following package versions

import Pkg; Pkg.status()
Status `~/work/MRIReco.jl/MRIReco.jl/docs/Project.toml`
  [13f3f980] CairoMakie v0.11.10
⌃ [e30172f5] Documenter v0.25.2
  [8ad4436d] ImageUtils v0.2.11
  [98b081ad] Literate v2.16.1
  [f7771a9a] MRIBase v0.4.2 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRIBase#master`
  [c57eb701] MRICoilSensitivities v0.1.3 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRICoilSensitivities#master`
  [5a6f062f] MRIFiles v0.3.0 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRIFiles#master`
  [fb1137e3] MRIOperators v0.2.1 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRIOperators#master`
  [bdf86e05] MRIReco v0.9.0 `~/work/MRIReco.jl/MRIReco.jl`
  [9be66c26] MRISampling v0.1.2 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRISampling#master`
  [8988da37] MRISimulation v0.1.2 `/home/runner/work/MRIReco.jl/MRIReco.jl:MRISimulation#master`
  [b77e0a4c] InteractiveUtils
Info Packages marked with ⌃ have new versions available and may be upgradable.

This page was generated using Literate.jl.