Skip to content

Commit

Permalink
Merge pull request #210 from cncastillo/observables
Browse files Browse the repository at this point in the history
Improve UI Observables management
  • Loading branch information
cncastillo authored Nov 14, 2023
2 parents be764be + 16743c4 commit 90d345e
Show file tree
Hide file tree
Showing 123 changed files with 7,393 additions and 3,660 deletions.
5 changes: 5 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,11 @@ qte_timings.txt
*.nsys-rep
*.mem
docs/src/generated/
docs/src/assets/logo*
docs/src/assets/icon*
src/ui/assets/*
!src/ui/assets/.gitkeep
!src/ui/assets/logo-icon.png
examples/3.koma_paper/comparison_accuracy/**/*.mrd
examples/3.koma_paper/mrf/**/*.mrd
examples/3.koma_paper/comparison_accuracy/**/*.mat
Expand Down
129 changes: 74 additions & 55 deletions KomaMRICore/src/datatypes/Phantom.jl
Original file line number Diff line number Diff line change
@@ -1,21 +1,22 @@
"""
obj = Phantom(name, x, y, z, ρ, T1, T2, T2s, Δw, Dλ1, Dλ2, Dθ, ux, uy, uz)
The Phantom struct.
The Phantom struct. Most of its field names are vectors, with each element associated with
a property value representing a spin. This struct serves as an input for the simulation.
# Arguments
- `name`: (`::String`) name of the phantom
- `x`: (`::AbstractVector{T}`, `[m]`) vector of x-positions of the spins
- `y`: (`::AbstractVector{T}`, `[m]`) vector of y-positions of the spins
- `z`: (`::AbstractVector{T}`, `[m]`) vector of z-positions of the spins
- `ρ`: (`::AbstractVector{T}`) vector of proton density of the spins
- `T1`: (`::AbstractVector{T}`, `[s]`) vector of T1 parameters of the spins
- `T2`: (`::AbstractVector{T}`, `[s]`) vector of T2 parameters of the spins
- `T2s`: (`::AbstractVector{T}`, `[s]`) vector of T2s parameters of the spins
- `Δw`: (`::AbstractVector{T}`, `[rad/s]`) vector of off-resonance parameters of the spins
- `Dλ1`: (`::AbstractVector{T}`) vector of Dλ1 (diffusion) parameters of the spins
- `Dλ2`: (`::AbstractVector{T}`) vector of Dλ2 (diffusion) parameters of the spins
- `Dθ`: (`::AbstractVector{T}`) vector of Dθ (diffusion) parameters of the spins
- `name`: (`::String`) phantom name
- `x`: (`::AbstractVector{T<:Real}`, `[m]`) spin x-position vector
- `y`: (`::AbstractVector{T<:Real}`, `[m]`) spin y-position vector
- `z`: (`::AbstractVector{T<:Real}`, `[m]`) spin z-position vector
- `ρ`: (`::AbstractVector{T<:Real}`) spin proton density vector
- `T1`: (`::AbstractVector{T<:Real}`, `[s]`) spin T1 parameter vector
- `T2`: (`::AbstractVector{T<:Real}`, `[s]`) spin T2 parameter vector
- `T2s`: (`::AbstractVector{T<:Real}`, `[s]`) spin T2s parameter vector
- `Δw`: (`::AbstractVector{T<:Real}`, `[rad/s]`) spin off-resonance parameter vector
- `Dλ1`: (`::AbstractVector{T<:Real}`) spin Dλ1 (diffusion) parameter vector
- `Dλ2`: (`::AbstractVector{T<:Real}`) spin Dλ2 (diffusion) parameter vector
- `Dθ`: (`::AbstractVector{T<:Real}`) spin Dθ (diffusion) parameter vector
- `ux`: (`::Function`) displacement field in the x-axis
- `uy`: (`::Function`) displacement field in the y-axis
- `uz`: (`::Function`) displacement field in the z-axis
Expand All @@ -25,22 +26,22 @@ The Phantom struct.
# Examples
```julia-repl
julia> obj = Phantom()
julia> obj = Phantom(x=[0.0])
julia> obj.ρ
```
"""
@with_kw mutable struct Phantom{T<:Real}
name::String = "spins"
x::AbstractVector{T}
y::AbstractVector{T} = zeros(size(x))
z::AbstractVector{T} = zeros(size(x))
ρ::AbstractVector{T} = ones(size(x))
T1::AbstractVector{T} = ones(size(x)) * 1_000_000
T2::AbstractVector{T} = ones(size(x)) * 1_000_000
y::AbstractVector{T} = zeros(size(x))
z::AbstractVector{T} = zeros(size(x))
ρ::AbstractVector{T} = ones(size(x))
T1::AbstractVector{T} = ones(size(x)) * 1_000_000
T2::AbstractVector{T} = ones(size(x)) * 1_000_000
T2s::AbstractVector{T} = ones(size(x)) * 1_000_000
#Off-resonance related
Δw::AbstractVector{T} = zeros(size(x))
Δw::AbstractVector{T} = zeros(size(x))
#χ::Vector{SusceptibilityModel}
#Diffusion
Dλ1::AbstractVector{T} = zeros(size(x))
Expand Down Expand Up @@ -77,7 +78,9 @@ Base.isapprox(obj1::Phantom, obj2::Phantom) = begin
obj1. obj2.
end

"""Separate object spins in a sub-group."""
"""
Separate object spins in a sub-group
"""
Base.getindex(obj::Phantom, p::AbstractRange) = begin
Phantom(name=obj.name,
x=obj.x[p],
Expand Down Expand Up @@ -166,9 +169,9 @@ end
end

"""
phantom = brain_phantom2D(;axis="axial", ss=4)
obj = brain_phantom2D(; axis="axial", ss=4)
Creates a two-dimentional brain phantom struct.
Creates a two-dimensional brain Phantom struct.
# References
- B. Aubert-Broche, D.L. Collins, A.C. Evans: "A new improved version of the realistic
Expand All @@ -178,24 +181,27 @@ Creates a two-dimentional brain phantom struct.
- https://brainweb.bic.mni.mcgill.ca/brainweb
# Keywords
- `axis`: (`::String`, `="axial"`, opts=[`"axial"`]) orientation of the phantom
- `ss`: (`::Real`, `=4`) subsampling parameter in all axis
- `axis`: (`::String`, `="axial"`, opts=[`"axial"`, `"coronal"`, `"sagittal"`]) orientation of the phantom
- `ss`: (`::Integer`, `=4`) subsampling parameter in all axis
# Returns
- `phantom`: (`::Phantom`) 2D Phantom struct
- `obj`: (`::Phantom`) Phantom struct
# Examples
```julia-repl
julia> obj = brain_phantom2D()
julia> obj = brain_phantom2D(; axis="sagittal", ss=1)
julia> plot_phantom_map(obj, :ρ)
```
"""
function brain_phantom2D(;axis="axial", ss=4)
function brain_phantom2D(; axis="axial", ss=4)

# Get data from .mat file
path = @__DIR__
data = MAT.matread(path*"/phantom/brain2D.mat")

class = data[axis][1:ss:end,1:ss:end]

# Define spin position vectors
Δx = .5e-3*ss
M, N = size(class)
FOVx = (M-1)*Δx #[m]
Expand All @@ -204,6 +210,7 @@ function brain_phantom2D(;axis="axial", ss=4)
y = -FOVy/2:Δx:FOVy/2 #spin coordinates
x, y = x .+ y'*0, x*0 .+ y' #grid points

# Define spin property vectors
T2 = (class.==23)*329 .+ #CSF
(class.==46)*83 .+ #GM
(class.==70)*70 .+ #WM
Expand Down Expand Up @@ -255,22 +262,26 @@ function brain_phantom2D(;axis="axial", ss=4)
T1 = T1*1e-3
T2 = T2*1e-3
T2s = T2s*1e-3
phantom = Phantom{Float64}(name="brain2D_"*axis,
x=y[ρ.!=0],
y=x[ρ.!=0],
z=0*x[ρ.!=0],
ρ=ρ[ρ.!=0],
T1=T1[ρ.!=0],
T2=T2[ρ.!=0],
T2s=T2s[ρ.!=0],
Δw=Δw[ρ.!=0])
phantom

# Define and return the Phantom struct
obj = Phantom{Float64}(
name = "brain2D_"*axis,
x = y[ρ.!=0],
y = x[ρ.!=0],
z = 0*x[ρ.!=0],
ρ = ρ[ρ.!=0],
T1 = T1[ρ.!=0],
T2 = T2[ρ.!=0],
T2s = T2s[ρ.!=0],
Δw = Δw[ρ.!=0],
)
return obj
end

"""
phantom = brain_phantom3D(;ss=4)
obj = brain_phantom3D(; ss=4)
Creates a three-dimentional brain phantom struct.
Creates a three-dimentional brain Phantom struct.
# References
- B. Aubert-Broche, D.L. Collins, A.C. Evans: "A new improved version of the realistic
Expand All @@ -280,23 +291,26 @@ Creates a three-dimentional brain phantom struct.
- https://brainweb.bic.mni.mcgill.ca/brainweb
# Keywords
- `ss`: (`::Real`, `=4`) subsampling parameter in all axis
- `ss`: (`::Integer`, `=4`) subsampling parameter in all axes
# Returns
- `phantom`: (`::Phantom`) 3D Phantom struct
- `obj`: (`::Phantom`) 3D Phantom struct
# Examples
```julia-repl
julia> obj = brain_phantom3D()
julia> obj = brain_phantom3D(; ss=5)
julia> plot_phantom_map(obj, :ρ)
```
"""
function brain_phantom3D(;ss=4,start_end=[160, 200])

# Get data from .mat file
path = @__DIR__
data = MAT.matread(path*"/phantom/brain3D.mat")

class = data["data"][1:ss:end,1:ss:end,start_end[1]:ss:start_end[2]]

# Define spin position vectors
Δx = .5e-3*ss
M, N, Z = size(class)
FOVx = (M-1)*Δx #[m]
Expand All @@ -309,6 +323,7 @@ function brain_phantom3D(;ss=4,start_end=[160, 200])
y = 0*xx .+ 1*yy .+ 0*zz
z = 0*xx .+ 0*yy .+ 1*zz

# Define spin property vectors
T2 = (class.==23)*329 .+ #CSF
(class.==46)*83 .+ #GM
(class.==70)*70 .+ #WM
Expand Down Expand Up @@ -360,14 +375,18 @@ function brain_phantom3D(;ss=4,start_end=[160, 200])
T1 = T1*1e-3
T2 = T2*1e-3
T2s = T2s*1e-3
phantom = Phantom{Float64}(name ="brain3D",
x = y[ρ.!=0],
y = x[ρ.!=0],
z = z[ρ.!=0],
ρ = ρ[ρ.!=0],
T1 = T1[ρ.!=0],
T2 = T2[ρ.!=0],
T2s = T2s[ρ.!=0],
Δw = Δw[ρ.!=0])
phantom

# Define and return the Phantom struct
obj = Phantom{Float64}(
name = "brain3D",
x = y[ρ.!=0],
y = x[ρ.!=0],
z = z[ρ.!=0],
ρ = ρ[ρ.!=0],
T1 = T1[ρ.!=0],
T2 = T2[ρ.!=0],
T2s = T2s[ρ.!=0],
Δw = Δw[ρ.!=0],
)
return obj
end
5 changes: 3 additions & 2 deletions KomaMRICore/src/datatypes/Scanner.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,13 +2,14 @@
sys = Scanner(B0, B1, Gmax, Smax, ADC_Δt, seq_Δt, GR_Δt, RF_Δt,
RF_ring_down_T, RF_dead_time_T, ADC_dead_time_T)
The Scanner struct.
The Scanner struct. It contains hardware limitations of the MRI resonator. It is an input
for the simulation.
# Arguments
- `B0`: (`::Real`, `=1.5`, `[T]`) main magnetic field strength
- `B1`: (`::Real`, `=10e-6`, `[T]`) maximum RF amplitude
- `Gmax`: (`::Real`, `=60e-3`, `[T/m]`) maximum gradient amplitude
- `Smax`: (`::Real`, `=500`, `[mT/m/ms]`) gradient maximum slew-rate
- `Smax`: (`::Real`, `=500`, `[mT/m/ms]`) gradient's maximum slew-rate
- `ADC_Δt`: (`::Real`, `=2e-6`, `[s]`) ADC raster time
- `seq_Δt`: (`::Real`, `=1e-5`, `[s]`) sequence-block raster time
- `GR_Δt`: (`::Real`, `=1e-5`, `[s]`) gradient raster time
Expand Down
89 changes: 36 additions & 53 deletions KomaMRICore/src/datatypes/Sequence.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,16 +8,17 @@
seq = Sequence(GR::Array{Grad,1}, RF::Array{RF,1})
seq = Sequence(GR::Array{Grad,1}, RF::Array{RF,1}, A::ADC, DUR, DEF)
The Sequence struct.
The Sequence struct. It contains events of an MRI sequence. Most field names (except for the
DEF field) consist of matrices or vectors, where each column index represents a sequence
block. This struct serves as an input for the simulation.
# Arguments
- `GR`: (`::Matrix{Grad}`) gradient matrix, rows are for (x,y,z) and columns are for time
- `RF`: (`::Matrix{RF}`) RF matrix, the 1 row is for the coil and columns are for time
- `ADC`: (`::Vector{ADC}`) ADC vector in time
- `DUR`: (`::Vector{Float64}`, `[s]`) duration of each sequence-block, this enables
delays after RF pulses to satisfy ring-down times
- `GR`: (`::Matrix{Grad}`) gradient matrix. Rows for x-y-z amplitudes and columns are for blocks
- `RF`: (`::Matrix{RF}`) RF matrix. The 1 row is for the coil and columns are for blocks
- `ADC`: (`::Array{ADC,1}`) ADC block vector
- `DUR`: (`::Vector`, `[s]`) duration block vector
- `DEF`: (`::Dict{String, Any}`) dictionary with relevant information of the sequence.
The possible keys are [`"AdcRasterTime"`, `"GradientRasterTime"`, `"Name"`, `"Nz"`,
Possible keys could be [`"AdcRasterTime"`, `"GradientRasterTime"`, `"Name"`, `"Nz"`,
`"Num_Blocks"`, `"Nx"`, `"Ny"`, `"PulseqVersion"`, `"BlockDurationRaster"`,
`"FileName"`, `"RadiofrequencyRasterTime"`]
Expand All @@ -42,60 +43,42 @@ mutable struct Sequence
end
end

#MAIN CONSTRUCTORS
function Sequence(GR) #If no RF is defined, just use a zero amplitude pulse
M,N = size(GR)
Sequence([i <= M ? GR[i,j] : Grad(0, 0) for i=1:3, j=1:N],
reshape([RF(0, 0) for i = 1:N],1,:),
[ADC(0, 0) for i = 1:N],
GR.dur,
Dict()
)
# Main Constructors
function Sequence(GR)
M, N = size(GR)
gr = [i <= M ? GR[i,j] : Grad(0, 0) for i in 1:3, j in 1:N]
rf = reshape([RF(0, 0) for i in 1:N], 1, :)
adc = [ADC(0, 0) for _ = 1:N]
return Sequence(gr, rf, adc, GR.dur, Dict())
end
function Sequence(GR, RF) #If no ADC is DEFined, just use a ADC with 0 samples
@assert size(GR,2) .== size(RF,2) "The number of Gradient, and RF objects must be the same."
M,N = size(GR)
Sequence([i <= M ? GR[i,j] : Grad(0, 0) for i=1:3, j=1:N],
RF,
[ADC(0, 0) for i = 1:N],
maximum([GR.dur RF.dur],dims=2)[:],
Dict()
)
function Sequence(GR, RF)
@assert size(GR, 2) .== size(RF, 2) "The number of Gradient, and RF objects must be the same."
M, N = size(GR)
gr = [i <= M ? GR[i,j] : Grad(0, 0) for i in 1:3, j in 1:N]
adc = [ADC(0, 0) for _ in 1:N]
dur = maximum([GR.dur RF.dur], dims=2)[:]
return Sequence(gr, RF, adc, dur, Dict())
end
function Sequence(GR, RF, ADC)
@assert size(GR,2) .== size(RF,2) .== length(ADC) "The number of Gradient, RF, and ADC objects must be the same."
M,N = size(GR)
Sequence([i <= M ? GR[i,j] : Grad(0, 0) for i=1:3, j=1:N],
RF,
ADC,
maximum([GR.dur RF.dur ADC.dur],dims=2)[:],
Dict())
@assert size(GR, 2) .== size(RF, 2) .== length(ADC) "The number of Gradient, RF, and ADC objects must be the same."
M, N = size(GR)
gr = [i <= M ? GR[i,j] : Grad(0, 0) for i in 1:3, j in 1:N]
dur = maximum([GR.dur RF.dur ADC.dur], dims=2)[:]
return Sequence(gr, RF, ADC, dur, Dict())
end
function Sequence(GR, RF, ADC, DUR)
@assert size(GR,2) .== size(RF,2) .== length(ADC) .== length(DUR) "The number of Gradient, RF, ADC, and DUR objects must be the same."
M,N = size(GR)
Sequence([i <= M ? GR[i,j] : Grad(0, GR[1,j].T, GR[1,j].rise, GR[1,j].fall, GR[1,j].delay) for i=1:3, j=1:N],
RF,
ADC,
maximum([GR.dur RF.dur ADC.dur DUR],dims=2)[:],
Dict())
@assert size(GR, 2) .== size(RF, 2) .== length(ADC) .== length(DUR) "The number of Gradient, RF, ADC, and DUR objects must be the same."
M, N = size(GR)
gr = [i <= M ? GR[i,j] : Grad(0, GR[1,j].T, GR[1,j].rise, GR[1,j].fall, GR[1,j].delay) for i in 1:3, j in 1:N]
dur = maximum([GR.dur RF.dur ADC.dur DUR], dims=2)[:]
return Sequence(gr, RF, ADC, dur, Dict())
end

#OTHER CONSTRUCTORS
# Other constructors
Sequence(GR::Array{Grad,1}) = Sequence(reshape(GR,1,:))
Sequence(GR::Array{Grad,1}, RF::Array{RF,1}) = Sequence(
reshape(GR,:,1),
reshape(RF,1,:),
[ADC(0,0) for i = 1:size(GR,2)]
)
Sequence(GR::Array{Grad,1}, RF::Array{RF,1}, A::ADC, DUR, DEF) = Sequence(
reshape(GR,:,1),
reshape(RF,1,:),
[A],
[DUR],
DEF
)
Sequence() = Sequence([Grad(0,0)])
Sequence(GR::Array{Grad,1}, RF::Array{RF,1})= Sequence(reshape(GR, :, 1), reshape(RF, 1, :), [ADC(0, 0) for i in 1:size(GR, 2)])
Sequence(GR::Array{Grad,1}, RF::Array{RF,1}, A::ADC, DUR, DEF)= Sequence(reshape(GR, :, 1), reshape(RF, 1, :), [A], [DUR], DEF)
Sequence() = Sequence([Grad(0, 0)])

"""
str = show(io::IO, s::Sequence)
Expand Down
Loading

0 comments on commit 90d345e

Please sign in to comment.