diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml
index 682d274a0..0a79ef058 100644
--- a/.github/workflows/CI.yml
+++ b/.github/workflows/CI.yml
@@ -19,7 +19,7 @@ jobs:
version:
- '1.10' # Replace this with the minimum Julia version that your package supports. E.g. if your package requires Julia 1.5 or higher, change this to '1.5'.
- '1' # Leave this line unchanged. '1' will automatically expand to the latest stable 1.x release of Julia.
- os: [ubuntu-latest, windows-latest, macos-12] # macos-latest] <- M1 Mac was generating problems #386, commented for now
+ os: [ubuntu-latest, windows-latest, macos-latest]
arch: [x64]
include:
- os: ubuntu-latest
diff --git a/KomaMRIBase/src/KomaMRIBase.jl b/KomaMRIBase/src/KomaMRIBase.jl
index 2f2744b0e..2f819bf89 100644
--- a/KomaMRIBase/src/KomaMRIBase.jl
+++ b/KomaMRIBase/src/KomaMRIBase.jl
@@ -26,8 +26,8 @@ include("timing/KeyValuesCalculation.jl")
include("datatypes/Sequence.jl")
include("datatypes/sequence/Delay.jl")
# Motion
-include("motion/motionlist/MotionList.jl")
-include("motion/nomotion/NoMotion.jl")
+include("motion/MotionList.jl")
+include("motion/NoMotion.jl")
# Phantom
include("datatypes/Phantom.jl")
# Simulator
@@ -53,7 +53,7 @@ export MotionList, NoMotion, Motion
export Translate, TranslateX, TranslateY, TranslateZ
export Rotate, RotateX, RotateY, RotateZ
export HeartBeat, Path, FlowPath
-export TimeRange, Periodic
+export TimeRange, Periodic, TimeCurve
export SpinRange, AllSpins
export get_spin_coords
# Secondary
diff --git a/KomaMRIBase/src/datatypes/Phantom.jl b/KomaMRIBase/src/datatypes/Phantom.jl
index 94ee8fa3f..a00d90a90 100644
--- a/KomaMRIBase/src/datatypes/Phantom.jl
+++ b/KomaMRIBase/src/datatypes/Phantom.jl
@@ -17,7 +17,7 @@ a property value representing a spin. This struct serves as an input for the sim
- `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
-- `motion`: (`::Union{NoMotion, MotionList{T<:Real}}`) motion
+- `motion`: (`::Union{NoMotion, Motion{T<:Real} MotionList{T<:Real}}`) motion
# Returns
- `obj`: (`::Phantom`) Phantom struct
@@ -47,7 +47,7 @@ julia> obj.ρ
Dθ::AbstractVector{T} = zeros(eltype(x), size(x))
#Diff::Vector{DiffusionModel} #Diffusion map
#Motion
- motion::Union{NoMotion, MotionList{T}} = NoMotion()
+ motion::Union{NoMotion, Motion{T}, MotionList{T}} = NoMotion()
end
const NON_STRING_PHANTOM_FIELDS = Iterators.filter(x -> fieldtype(Phantom, x) != String, fieldnames(Phantom))
diff --git a/KomaMRIBase/src/datatypes/simulation/DiscreteSequence.jl b/KomaMRIBase/src/datatypes/simulation/DiscreteSequence.jl
index fe8467db5..cbf6a2be5 100644
--- a/KomaMRIBase/src/datatypes/simulation/DiscreteSequence.jl
+++ b/KomaMRIBase/src/datatypes/simulation/DiscreteSequence.jl
@@ -85,8 +85,8 @@ based on simulation parameters.
# Returns
- `seqd`: (`::DiscreteSequence`) DiscreteSequence struct
"""
-function discretize(seq::Sequence; sampling_params=default_sampling_params())
- t, Δt = get_variable_times(seq; Δt=sampling_params["Δt"], Δt_rf=sampling_params["Δt_rf"])
+function discretize(seq::Sequence; sampling_params=default_sampling_params(), motion=NoMotion())
+ t, Δt = get_variable_times(seq; Δt=sampling_params["Δt"], Δt_rf=sampling_params["Δt_rf"], motion=motion)
B1, Δf = get_rfs(seq, t)
Gx, Gy, Gz = get_grads(seq, t)
tadc = get_adc_sampling_times(seq)
diff --git a/KomaMRIBase/src/motion/motionlist/Action.jl b/KomaMRIBase/src/motion/Action.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/Action.jl
rename to KomaMRIBase/src/motion/Action.jl
diff --git a/KomaMRIBase/src/motion/motionlist/actions/ArbitraryAction.jl b/KomaMRIBase/src/motion/Interpolation.jl
similarity index 59%
rename from KomaMRIBase/src/motion/motionlist/actions/ArbitraryAction.jl
rename to KomaMRIBase/src/motion/Interpolation.jl
index c62a7b4bd..09b699867 100644
--- a/KomaMRIBase/src/motion/motionlist/actions/ArbitraryAction.jl
+++ b/KomaMRIBase/src/motion/Interpolation.jl
@@ -29,19 +29,6 @@ const Interpolator2D = Interpolations.GriddedInterpolation{
Itp<:Interpolations.Gridded,
K<:Tuple{AbstractVector{TNodes}, AbstractVector{TNodes}},
}
-
-abstract type ArbitraryAction{T<:Real} <: AbstractAction{T} end
-
-function Base.getindex(action::ArbitraryAction, p)
- return typeof(action)([getfield(action, d)[p,:] for d in fieldnames(typeof(action))]...)
-end
-function Base.view(action::ArbitraryAction, p)
- return typeof(action)([@view(getfield(action, d)[p,:]) for d in fieldnames(typeof(action))]...)
-end
-
-Base.:(==)(m1::ArbitraryAction, m2::ArbitraryAction) = (typeof(m1) == typeof(m2)) & reduce(&, [getfield(m1, field) == getfield(m2, field) for field in fieldnames(typeof(m1))])
-Base.:(≈)(m1::ArbitraryAction, m2::ArbitraryAction) = (typeof(m1) == typeof(m2)) & reduce(&, [getfield(m1, field) ≈ getfield(m2, field) for field in fieldnames(typeof(m1))])
-
function GriddedInterpolation(nodes, A, ITP)
return Interpolations.GriddedInterpolation{eltype(A), length(nodes), typeof(A), typeof(ITP), typeof(nodes)}(nodes, A, ITP)
end
@@ -70,26 +57,10 @@ function resample(itp::Interpolator2D, t)
return itp.(id, t)
end
-function displacement_x!(ux, action::ArbitraryAction, x, y, z, t)
- itp = interpolate(action.dx, Gridded(Linear()), Val(size(action.dx,1)), t)
- ux .= resample(itp, t)
- return nothing
-end
-
-function displacement_y!(uy, action::ArbitraryAction, x, y, z, t)
- itp = interpolate(action.dy, Gridded(Linear()), Val(size(action.dy,1)), t)
- uy .= resample(itp, t)
- return nothing
-end
-
-function displacement_z!(uz, action::ArbitraryAction, x, y, z, t)
- itp = interpolate(action.dz, Gridded(Linear()), Val(size(action.dz,1)), t)
- uz .= resample(itp, t)
- return nothing
+function interpolate_times(t, t_unit, periodic, tq)
+ itp = GriddedInterpolation((t, ), t_unit, Gridded(Linear()))
+ return extrapolate(itp, periodic ? Interpolations.Periodic() : Flat()).(tq)
end
_similar(a, N) = similar(a, N)
-_similar(a::Real, N) = zeros(typeof(a), N)
-
-include("arbitraryactions/Path.jl")
-include("arbitraryactions/FlowPath.jl")
\ No newline at end of file
+_similar(a::Real, N) = zeros(typeof(a), N)
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/motionlist/Motion.jl b/KomaMRIBase/src/motion/Motion.jl
similarity index 68%
rename from KomaMRIBase/src/motion/motionlist/Motion.jl
rename to KomaMRIBase/src/motion/Motion.jl
index 869d0e06c..6eed9d3c7 100644
--- a/KomaMRIBase/src/motion/motionlist/Motion.jl
+++ b/KomaMRIBase/src/motion/Motion.jl
@@ -11,7 +11,7 @@ that are affected by that motion.
# Arguments
- `action`: (`::AbstractAction{T<:Real}`) action, such as [`Translate`](@ref) or [`Rotate`](@ref)
-- `time`: (`::AbstractTimeSpan{T<:Real}`, `=TimeRange(0.0)`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`, `=TimeRange(0.0)`) time information about the motion
- `spins`: (`::AbstractSpinSpan`, `=AllSpins()`) spin indexes affected by the motion
# Returns
@@ -28,22 +28,22 @@ julia> motion = Motion(
"""
@with_kw mutable struct Motion{T<:Real}
action::AbstractAction{T}
- time ::AbstractTimeSpan{T} = TimeRange(zero(typeof(action).parameters[1]))
- spins ::AbstractSpinSpan = AllSpins()
+ time ::TimeCurve{T} = TimeRange(t_start=zero(typeof(action).parameters[1]), t_end=eps(typeof(action).parameters[1]))
+ spins ::AbstractSpinSpan = AllSpins()
end
# Main constructors
function Motion(action)
T = first(typeof(action).parameters)
- return Motion(action, TimeRange(zero(T)), AllSpins())
+ return Motion(action, TimeRange(t_start=zero(T), t_end=eps(T)), AllSpins())
end
-function Motion(action, time::AbstractTimeSpan)
+function Motion(action, time::TimeCurve)
T = first(typeof(action).parameters)
return Motion(action, time, AllSpins())
end
function Motion(action, spins::AbstractSpinSpan)
T = first(typeof(action).parameters)
- return Motion(action, TimeRange(zero(T)), spins)
+ return Motion(action, TimeRange(t_start=zero(T), t_end=eps(T)), spins)
end
# Custom constructors
@@ -54,7 +54,7 @@ end
- `dx`: (`::Real`, `[m]`) translation in x
- `dy`: (`::Real`, `[m]`) translation in y
- `dz`: (`::Real`, `[m]`) translation in z
-- `time`: (`::AbstractTimeSpan{T<:Real}`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`) time information about the motion
- `spins`: (`::AbstractSpinSpan`) spin indexes affected by the motion
# Returns
@@ -65,7 +65,7 @@ end
julia> translate = Translate(0.01, 0.02, 0.03, TimeRange(0.0, 1.0), SpinRange(1:10))
```
"""
-function Translate(dx, dy, dz, time=TimeRange(zero(eltype(dx))), spins=AllSpins())
+function Translate(dx, dy, dz, time=TimeRange(t_start=zero(eltype(dx)), t_end=eps(eltype(dx))), spins=AllSpins())
return Motion(Translate(dx, dy, dz), time, spins)
end
@@ -76,7 +76,7 @@ end
- `pitch`: (`::Real`, `[º]`) rotation in x
- `roll`: (`::Real`, `[º]`) rotation in y
- `yaw`: (`::Real`, `[º]`) rotation in z
-- `time`: (`::AbstractTimeSpan{T<:Real}`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`) time information about the motion
- `spins`: (`::AbstractSpinSpan`) spin indexes affected by the motion
# Returns
@@ -87,7 +87,7 @@ end
julia> rotate = Rotate(15.0, 0.0, 20.0, TimeRange(0.0, 1.0), SpinRange(1:10))
```
"""
-function Rotate(pitch, roll, yaw, time=TimeRange(zero(eltype(pitch))), spins=AllSpins())
+function Rotate(pitch, roll, yaw, time=TimeRange(t_start=zero(eltype(pitch)), t_end=eps(eltype(pitch))), spins=AllSpins())
return Motion(Rotate(pitch, roll, yaw), time, spins)
end
@@ -98,7 +98,7 @@ end
- `circumferential_strain`: (`::Real`) contraction parameter
- `radial_strain`: (`::Real`) contraction parameter
- `longitudinal_strain`: (`::Real`) contraction parameter
-- `time`: (`::AbstractTimeSpan{T<:Real}`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`) time information about the motion
- `spins`: (`::AbstractSpinSpan`) spin indexes affected by the motion
# Returns
@@ -109,7 +109,7 @@ end
julia> heartbeat = HeartBeat(-0.3, -0.2, 0.0, TimeRange(0.0, 1.0), SpinRange(1:10))
```
"""
-function HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, time=TimeRange(zero(eltype(circumferential_strain))), spins=AllSpins())
+function HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, time=TimeRange(t_start=zero(eltype(circumferential_strain)), t_end=eps(eltype(circumferential_strain))), spins=AllSpins())
return Motion(HeartBeat(circumferential_strain, radial_strain, longitudinal_strain), time, spins)
end
@@ -120,7 +120,7 @@ end
- `dx`: (`::AbstractArray{T<:Real}`, `[m]`) displacements in x
- `dy`: (`::AbstractArray{T<:Real}`, `[m]`) displacements in y
- `dz`: (`::AbstractArray{T<:Real}`, `[m]`) displacements in z
-- `time`: (`::AbstractTimeSpan{T<:Real}`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`) time information about the motion
- `spins`: (`::AbstractSpinSpan`) spin indexes affected by the motion
# Returns
@@ -137,7 +137,7 @@ julia> path = Path(
)
```
"""
-function Path(dx, dy, dz, time=TimeRange(zero(eltype(dx))), spins=AllSpins())
+function Path(dx, dy, dz, time=TimeRange(t_start=zero(eltype(dx)), t_end=eps(eltype(dx))), spins=AllSpins())
return Motion(Path(dx, dy, dz), time, spins)
end
@@ -149,7 +149,7 @@ end
- `dy`: (`::AbstractArray{T<:Real}`, `[m]`) displacements in y
- `dz`: (`::AbstractArray{T<:Real}`, `[m]`) displacements in z
- `spin_reset`: (`::AbstractArray{Bool}`) reset spin state flags
-- `time`: (`::AbstractTimeSpan{T<:Real}`) time information about the motion
+- `time`: (`::TimeCurve{T<:Real}`) time information about the motion
- `spins`: (`::AbstractSpinSpan`) spin indexes affected by the motion
# Returns
@@ -167,7 +167,7 @@ julia> flowpath = FlowPath(
)
```
"""
-function FlowPath(dx, dy, dz, spin_reset, time=TimeRange(zero(eltype(dx))), spins=AllSpins())
+function FlowPath(dx, dy, dz, spin_reset, time=TimeRange(t_start=zero(eltype(dx)), t_end=eps(eltype(dx))), spins=AllSpins())
return Motion(FlowPath(dx, dy, dz, spin_reset), time, spins)
end
@@ -178,13 +178,30 @@ Base.:(≈)(m1::Motion, m2::Motion) = (typeof(m1) == typeof(m2)) & reduce(&, [g
""" Motion sub-group """
function Base.getindex(m::Motion, p)
idx, spin_range = m.spins[p]
- return idx !== nothing ? Motion(m.action[idx], m.time, spin_range) : nothing
+ return idx !== nothing ? Motion(m.action[idx], m.time, spin_range) : NoMotion()
end
function Base.view(m::Motion, p)
idx, spin_range = @view(m.spins[p])
- return idx !== nothing ? Motion(@view(m.action[idx]), m.time, spin_range) : nothing
+ return idx !== nothing ? Motion(@view(m.action[idx]), m.time, spin_range) : NoMotion()
+end
+
+"""
+ x, y, z = get_spin_coords(motion, x, y, z, t)
+"""
+function get_spin_coords(
+ m::Motion{T}, x::AbstractVector{T}, y::AbstractVector{T}, z::AbstractVector{T}, t
+) where {T<:Real}
+ ux, uy, uz = x .* (0*t), y .* (0*t), z .* (0*t) # Buffers for displacements
+ t_unit = unit_time(t, m.time.t, m.time.t_unit, m.time.periodic, m.time.periods)
+ idx = get_indexing_range(m.spins)
+ displacement_x!(@view(ux[idx, :]), m.action, @view(x[idx]), @view(y[idx]), @view(z[idx]), t_unit)
+ displacement_y!(@view(uy[idx, :]), m.action, @view(x[idx]), @view(y[idx]), @view(z[idx]), t_unit)
+ displacement_z!(@view(uz[idx, :]), m.action, @view(x[idx]), @view(y[idx]), @view(z[idx]), t_unit)
+ return x .+ ux, y .+ uy, z .+ uz
end
# Auxiliary functions
-times(m::Motion) = times(m.time)
-is_composable(m::Motion) = is_composable(m.action)
\ No newline at end of file
+times(m::Motion) = times(m.time.t, m.time.periods)
+is_composable(m::Motion) = is_composable(m.action)
+add_jump_times!(t, m::Motion) = add_jump_times!(t, m.action, m.time)
+add_jump_times!(t, ::AbstractAction, ::TimeCurve) = nothing
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/motionlist/MotionList.jl b/KomaMRIBase/src/motion/MotionList.jl
similarity index 69%
rename from KomaMRIBase/src/motion/motionlist/MotionList.jl
rename to KomaMRIBase/src/motion/MotionList.jl
index 953ad394f..d03fb6b78 100644
--- a/KomaMRIBase/src/motion/motionlist/MotionList.jl
+++ b/KomaMRIBase/src/motion/MotionList.jl
@@ -1,6 +1,7 @@
-include("Action.jl")
+include("Interpolation.jl")
include("SpinSpan.jl")
-include("TimeSpan.jl")
+include("TimeCurve.jl")
+include("Action.jl")
include("Motion.jl")
"""
@@ -36,27 +37,24 @@ struct MotionList{T<:Real}
motions::Vector{<:Motion{T}}
end
+# NOTE: this constructor must be simplified once the Vector{<:Motion} approach is accomplished:
+# https://github.com/JuliaHealth/KomaMRI.jl/issues/480
""" Constructors """
-MotionList(motions::Motion...) = length([motions]) > 0 ? MotionList([motions...]) : @error "You must provide at least one motion as input argument. If you do not want to define motion, use `NoMotion{T}()`"
-
-""" MotionList sub-group """
-function Base.getindex(mv::MotionList{T}, p) where {T<:Real}
- motion_array_aux = Motion{T}[]
- for m in mv.motions
- m[p] !== nothing ? push!(motion_array_aux, m[p]) : nothing
- end
- return length(motion_array_aux) > 0 ? MotionList(motion_array_aux) : NoMotion()
-end
-function Base.view(mv::MotionList{T}, p) where {T<:Real}
- motion_array_aux = Motion{T}[]
- for m in mv.motions
- @view(m[p]) !== nothing ? push!(motion_array_aux, @view(m[p])) : nothing
- end
- return length(motion_array_aux) > 0 ? MotionList(motion_array_aux) : NoMotion()
+function MotionList(motions::Motion...)
+ if length(motions) == 0
+ return NoMotion()
+ elseif length(motions) == 1
+ return motions[1]
+ else
+ return MotionList([motions...])
+ end
end
-""" Addition of MotionLists """
-function Base.vcat(m1::MotionList{T}, m2::MotionList{T}, Ns1::Int, Ns2::Int) where {T<:Real}
+# NOTE: these vcat methods must be simplified once the Vector{<:Motion} approach is accomplished:
+# https://github.com/JuliaHealth/KomaMRI.jl/issues/480
+""" Addition of MotionLists """
+# MotionList + MotionList
+function Base.vcat(m1::MotionList{T}, m2::MotionList{T}, Ns1, Ns2) where {T<:Real}
mv_aux = Motion{T}[]
for m in m1.motions
m_aux = copy(m)
@@ -69,7 +67,50 @@ function Base.vcat(m1::MotionList{T}, m2::MotionList{T}, Ns1::Int, Ns2::Int) whe
m_aux.spins = SpinRange(m_aux.spins.range .+ Ns1)
push!(mv_aux, m_aux)
end
- return MotionList(mv_aux)
+ return MotionList(mv_aux...)
+end
+# Motion + Motion
+function Base.vcat(m1::Motion{T}, m2::Motion{T}, Ns1, Ns2) where {T<:Real}
+ mv_aux = Motion{T}[]
+ m_aux = copy(m1)
+ m_aux.spins = expand(m_aux.spins, Ns1)
+ push!(mv_aux, m_aux)
+ m_aux = copy(m2)
+ m_aux.spins = expand(m_aux.spins, Ns2)
+ m_aux.spins = SpinRange(m_aux.spins.range .+ Ns1)
+ push!(mv_aux, m_aux)
+ return MotionList(mv_aux...)
+end
+# Motion + MotionList
+Base.vcat(m1::MotionList{T}, m2::Motion{T}, Ns1, Ns2) where {T<:Real} = vcat(m2, m1, Ns2, Ns1)
+function Base.vcat(m1::Motion{T}, m2::MotionList{T}, Ns1, Ns2) where {T<:Real}
+ mv_aux = Motion{T}[]
+ m_aux = copy(m1)
+ m_aux.spins = expand(m_aux.spins, Ns1)
+ push!(mv_aux, m_aux)
+ for m in m2.motions
+ m_aux = copy(m)
+ m_aux.spins = expand(m_aux.spins, Ns2)
+ m_aux.spins = SpinRange(m_aux.spins.range .+ Ns1)
+ push!(mv_aux, m_aux)
+ end
+ return MotionList(mv_aux...)
+end
+
+""" MotionList sub-group """
+function Base.getindex(mv::MotionList{T}, p) where {T<:Real}
+ motion_array_aux = Motion{T}[]
+ for m in mv.motions
+ m[p] isa NoMotion ? nothing : push!(motion_array_aux, m[p])
+ end
+ return MotionList(motion_array_aux...)
+end
+function Base.view(mv::MotionList{T}, p) where {T<:Real}
+ motion_array_aux = Motion{T}[]
+ for m in mv.motions
+ @view(m[p]) isa NoMotion ? nothing : push!(motion_array_aux, @view(m[p]))
+ end
+ return MotionList(motion_array_aux...)
end
""" Compare two MotionLists """
@@ -116,7 +157,7 @@ function get_spin_coords(
ux, uy, uz = xt .* zero(T), yt .* zero(T), zt .* zero(T)
# Composable motions: they need to be run sequentially. Note that they depend on xt, yt, and zt
for m in Iterators.filter(is_composable, ml.motions)
- t_unit = unit_time(t, m.time)
+ t_unit = unit_time(t, m.time.t, m.time.t_unit, m.time.periodic, m.time.periods)
idx = get_indexing_range(m.spins)
displacement_x!(@view(ux[idx, :]), m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
displacement_y!(@view(uy[idx, :]), m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
@@ -126,7 +167,7 @@ function get_spin_coords(
end
# Additive motions: these motions can be run in parallel
for m in Iterators.filter(!is_composable, ml.motions)
- t_unit = unit_time(t, m.time)
+ t_unit = unit_time(t, m.time.t, m.time.t_unit, m.time.periodic, m.time.periods)
idx = get_indexing_range(m.spins)
displacement_x!(@view(ux[idx, :]), m.action, @view(x[idx]), @view(y[idx]), @view(z[idx]), t_unit)
displacement_y!(@view(uy[idx, :]), m.action, @view(x[idx]), @view(y[idx]), @view(z[idx]), t_unit)
@@ -159,6 +200,12 @@ If `motionset::MotionList`, this function sorts its motions.
- `nothing`
"""
function sort_motions!(m::MotionList)
- sort!(m.motions; by=m -> times(m)[1])
+ sort!(m.motions; by=m -> m.time.t_start)
return nothing
end
+
+function add_jump_times!(t, ml::MotionList)
+ for m in ml.motions
+ add_jump_times!(t, m)
+ end
+end
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/nomotion/NoMotion.jl b/KomaMRIBase/src/motion/NoMotion.jl
similarity index 73%
rename from KomaMRIBase/src/motion/nomotion/NoMotion.jl
rename to KomaMRIBase/src/motion/NoMotion.jl
index 43f27cb08..233d00a35 100644
--- a/KomaMRIBase/src/motion/nomotion/NoMotion.jl
+++ b/KomaMRIBase/src/motion/NoMotion.jl
@@ -17,7 +17,9 @@ Base.getindex(mv::NoMotion, p) = mv
Base.view(mv::NoMotion, p) = mv
""" Addition of NoMotions """
+# NoMotion + NoMotion
Base.vcat(m1::NoMotion, m2::NoMotion, Ns1, Ns2) = m1
+# NoMotion + MotionList
Base.vcat(m1::MotionList, m2::NoMotion, Ns1, Ns2) = vcat(m2, m1, 0, Ns1)
function Base.vcat(m1::NoMotion, m2::MotionList{T}, Ns1, Ns2) where {T}
mv_aux = Motion{T}[]
@@ -29,6 +31,14 @@ function Base.vcat(m1::NoMotion, m2::MotionList{T}, Ns1, Ns2) where {T}
end
return MotionList(mv_aux)
end
+# NoMotion + Motion
+Base.vcat(m1::Motion, m2::NoMotion, Ns1, Ns2) = vcat(m2, m1, 0, Ns1)
+function Base.vcat(m1::NoMotion, m2::Motion{T}, Ns1, Ns2) where {T}
+ m_aux = copy(m2)
+ m_aux.spins = expand(m_aux.spins, Ns2)
+ m_aux.spins = SpinRange(m_aux.spins.range .+ Ns1)
+ return m_aux
+end
""" Compare two NoMotions """
Base.:(==)(m1::NoMotion, m2::NoMotion) = true
@@ -39,3 +49,4 @@ function get_spin_coords(
) where {T<:Real}
return x, y, z
end
+add_jump_times!(t, ::NoMotion) = nothing
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/motionlist/SpinSpan.jl b/KomaMRIBase/src/motion/SpinSpan.jl
similarity index 92%
rename from KomaMRIBase/src/motion/motionlist/SpinSpan.jl
rename to KomaMRIBase/src/motion/SpinSpan.jl
index 21d1a7396..fab1b9787 100644
--- a/KomaMRIBase/src/motion/motionlist/SpinSpan.jl
+++ b/KomaMRIBase/src/motion/SpinSpan.jl
@@ -69,3 +69,5 @@ Base.length(sr::SpinRange) = length(sr.range)
get_indexing_range(spins::SpinRange) = spins.range
expand(sr::SpinRange, Ns::Int) = sr
intersect_idx(a, b) = findall(x -> x in a, b)
+intersect_idx(a, b::BitVector) = findall(x -> x in a, findall(x->x==true, b))
+intersect_idx(a::BitVector, b) = findall(x -> x in findall(x->x==true, a), b)
diff --git a/KomaMRIBase/src/motion/TimeCurve.jl b/KomaMRIBase/src/motion/TimeCurve.jl
new file mode 100644
index 000000000..235f6c9ca
--- /dev/null
+++ b/KomaMRIBase/src/motion/TimeCurve.jl
@@ -0,0 +1,136 @@
+@doc raw"""
+ timecurve = TimeCurve(t, t_unit, periodic, periods)
+
+TimeCurve struct. It is a specialized type that defines a time curve, which represents
+the temporal behavior of motion. This curve is defined by two vectors:
+`t` and `t_unit`, which represent the horizontal (x-axis) and vertical (y-axis) axes
+of the curve, respectively. To some extent, this curve can be associated with animation curves,
+commonly found in software for video editing, 3D scene creation, or video game development.
+
+Additionally, the TimeCurve struct contains two more fields, independent of each other:
+`periodic` is a Boolean that indicates whether the time curve should be repeated periodically.
+`periods` contains as many elements as repetitions are desired in the time curve.
+Each element specifies the scaling factor for that repetition.
+
+# Arguments
+- `t`: (`::AbstractVector{<:Real}`, `[s]`) time vector
+- `t_unit`: (`::AbstractVector{<:Real}`) y vector, it needs to be scaled between 0 and 1
+- `periodic`: (`::Bool`, `=false`) indicates whether the time curve should be periodically repeated
+- `periods`: (`::Union{<:Real,AbstractVector{<:Real}}`, `=1.0`): represents the relative duration
+ of each period with respect to the baseline duration defined by `t[end] - t[1]`.
+ In other words, it acts as a scaling factor to lengthen or shorten specific periods.
+ This allows for the creation of patterns such as arrhythmias or other variations in periodicity.
+
+# Returns
+- `timecurve`: (`::TimeCurve`) TimeCurve struct
+
+# Examples
+1\. Non-periodic motion with a single repetition:
+```julia-repl
+julia> timecurve = TimeCurve(t=[0.0, 0.2, 0.4, 0.6], t_unit=[0.0, 0.2, 0.5, 1.0])
+```
+
+
+2\. Periodic motion with a single repetition:
+```julia-repl
+julia> timecurve = TimeCurve(t=[0.0, 0.2, 0.4, 0.6], t_unit=[0.0, 1.0, 1.0, 0.0], periodic=true)
+```
+
+
+3\. Non-periodic motion with multiple repetitions:
+```julia-repl
+julia> timecurve = TimeCurve(t=[0.0, 0.2, 0.4, 0.6], t_unit=[0.0, 1.0, 1.0, 0.0], periods=[1.0, 0.5, 1.5])
+```
+
+
+4\. Periodic motion with multiple repetitions:
+```julia-repl
+julia> timecurve = TimeCurve(t=[0.0, 0.2, 0.4, 0.6], t_unit=[0.0, 1.0, 1.0, 0.0], periods=[1.0, 0.5, 1.5], periodic=true)
+```
+
+"""
+@with_kw struct TimeCurve{T<:Real}
+ t::AbstractVector{T}
+ t_unit::AbstractVector{T}
+ periodic::Bool = false
+ periods::Union{T,AbstractVector{T}} = oneunit(eltype(t))
+ t_start::T = t[1]
+ t_end::T = t[end]
+ @assert check_unique(t) "Vector t=$(t) contains duplicate elements. Please ensure all elements in t are unique and try again"
+end
+
+check_unique(t) = true
+check_unique(t::Vector) = length(t) == length(unique(t))
+
+# Main Constructors
+TimeCurve(t, t_unit, periodic, periods) = TimeCurve(t=t, t_unit=t_unit, periodic=periodic, periods=periods)
+TimeCurve(t, t_unit) = TimeCurve(t=t, t_unit=t_unit)
+# Custom constructors
+"""
+ timerange = TimeRange(t_start, t_end)
+
+The `TimeRange` function is a custom constructor for the `TimeCurve` struct.
+It allows defining a simple time interval, with start and end times.
+
+# Arguments
+- `t_start`: (`::Real`, `[s]`, `=0.0`) start time
+- `t_end`: (`::Real`, `[s]`, `=1.0`) end time
+
+# Returns
+- `timerange`: (`::TimeCurve`) TimeCurve struct
+
+# Examples
+```julia-repl
+julia> timerange = TimeRange(t_start=0.6, t_end=1.4)
+```
+
+"""
+TimeRange(t_start::T, t_end::T) where T = TimeCurve(t=[t_start, t_end], t_unit=[zero(T), oneunit(T)])
+TimeRange(; t_start=0.0, t_end=1.0) = TimeRange(t_start, t_end)
+"""
+ periodic = Periodic(period, asymmetry)
+
+The `Periodic` function is a custom constructor for the `TimeCurve` struct.
+It allows defining time intervals that repeat periodically with a triangular period.
+It includes a measure of asymmetry in order to recreate a asymmetric period.
+
+# Arguments
+- `period`: (`::Real`, `[s]`, `=1.0`) period duration
+- `asymmetry`: (`::Real`, `=0.5`) temporal asymmetry factor. Between 0 and 1.
+
+# Returns
+- `periodic`: (`::TimeCurve`) TimeCurve struct
+
+# Examples
+```julia-repl
+julia> periodic = Periodic(period=1.0, asymmetry=0.2)
+```
+
+"""
+Periodic(period::T, asymmetry::T) where T = TimeCurve(t=[zero(T), period*asymmetry, period], t_unit=[zero(T), oneunit(T), zero(T)])
+Periodic(; period=1.0, asymmetry=0.5) = Periodic(period, asymmetry)
+
+""" Compare two TimeCurves """
+Base.:(==)(t1::TimeCurve, t2::TimeCurve) = reduce(&, [getfield(t1, field) == getfield(t2, field) for field in fieldnames(typeof(t1))])
+Base.:(≈)(t1::TimeCurve, t2::TimeCurve) = reduce(&, [getfield(t1, field) ≈ getfield(t2, field) for field in fieldnames(typeof(t1))])
+
+""" times & unit_time """
+# Although the implementation of these two functions when `per` is a vector is valid
+# for all cases, it performs unnecessary and costly operations when `per` is a scalar.
+# Therefore, it has been decided to use method dispatch between these two cases.
+function times(t, per::Real)
+ return per .* t
+end
+function times(t, per::AbstractVector)
+ tr = repeat(t, length(per))
+ scale = repeat(per, inner=[length(t)])
+ offsets = repeat(vcat(0, cumsum(per)[1:end-1]), inner=[length(t)])
+ tr .= (tr .* scale) .+ offsets
+ return tr
+end
+function unit_time(tq, t, t_unit, periodic, per::Real)
+ return interpolate_times(t .* per, t_unit, periodic, tq)
+end
+function unit_time(tq, t, t_unit, periodic, per::AbstractVector)
+ return interpolate_times(times(t, per), repeat(t_unit, length(per)), periodic, tq)
+end
diff --git a/KomaMRIBase/src/motion/actions/ArbitraryAction.jl b/KomaMRIBase/src/motion/actions/ArbitraryAction.jl
new file mode 100644
index 000000000..5507d24f8
--- /dev/null
+++ b/KomaMRIBase/src/motion/actions/ArbitraryAction.jl
@@ -0,0 +1,32 @@
+abstract type ArbitraryAction{T<:Real} <: AbstractAction{T} end
+
+function Base.getindex(action::ArbitraryAction, p)
+ return typeof(action)([getfield(action, d)[p,:] for d in fieldnames(typeof(action))]...)
+end
+function Base.view(action::ArbitraryAction, p)
+ return typeof(action)([@view(getfield(action, d)[p,:]) for d in fieldnames(typeof(action))]...)
+end
+
+Base.:(==)(m1::ArbitraryAction, m2::ArbitraryAction) = (typeof(m1) == typeof(m2)) & reduce(&, [getfield(m1, field) == getfield(m2, field) for field in fieldnames(typeof(m1))])
+Base.:(≈)(m1::ArbitraryAction, m2::ArbitraryAction) = (typeof(m1) == typeof(m2)) & reduce(&, [getfield(m1, field) ≈ getfield(m2, field) for field in fieldnames(typeof(m1))])
+
+function displacement_x!(ux, action::ArbitraryAction, x, y, z, t)
+ itp = interpolate(action.dx, Gridded(Linear()), Val(size(action.dx,1)), t)
+ ux .= resample(itp, t)
+ return nothing
+end
+
+function displacement_y!(uy, action::ArbitraryAction, x, y, z, t)
+ itp = interpolate(action.dy, Gridded(Linear()), Val(size(action.dy,1)), t)
+ uy .= resample(itp, t)
+ return nothing
+end
+
+function displacement_z!(uz, action::ArbitraryAction, x, y, z, t)
+ itp = interpolate(action.dz, Gridded(Linear()), Val(size(action.dz,1)), t)
+ uz .= resample(itp, t)
+ return nothing
+end
+
+include("arbitraryactions/Path.jl")
+include("arbitraryactions/FlowPath.jl")
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/motionlist/actions/SimpleAction.jl b/KomaMRIBase/src/motion/actions/SimpleAction.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/actions/SimpleAction.jl
rename to KomaMRIBase/src/motion/actions/SimpleAction.jl
diff --git a/KomaMRIBase/src/motion/motionlist/actions/arbitraryactions/FlowPath.jl b/KomaMRIBase/src/motion/actions/arbitraryactions/FlowPath.jl
similarity index 84%
rename from KomaMRIBase/src/motion/motionlist/actions/arbitraryactions/FlowPath.jl
rename to KomaMRIBase/src/motion/actions/arbitraryactions/FlowPath.jl
index 5faec8bdb..5046b49a6 100644
--- a/KomaMRIBase/src/motion/motionlist/actions/arbitraryactions/FlowPath.jl
+++ b/KomaMRIBase/src/motion/actions/arbitraryactions/FlowPath.jl
@@ -34,4 +34,9 @@ julia> flowpath = FlowPath(
dy::AbstractArray{T}
dz::AbstractArray{T}
spin_reset::AbstractArray{Bool}
+end
+
+function add_jump_times!(t, a::FlowPath, tc::TimeCurve)
+ jump_times = (tc.t_end - tc.t_start)/(size(a.spin_reset)[2]-1) * (getindex.(findall(a.spin_reset .== 1), 2) .- 1) .- 1e-6
+ append!(t, jump_times)
end
\ No newline at end of file
diff --git a/KomaMRIBase/src/motion/motionlist/actions/arbitraryactions/Path.jl b/KomaMRIBase/src/motion/actions/arbitraryactions/Path.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/actions/arbitraryactions/Path.jl
rename to KomaMRIBase/src/motion/actions/arbitraryactions/Path.jl
diff --git a/KomaMRIBase/src/motion/motionlist/actions/simpleactions/HeartBeat.jl b/KomaMRIBase/src/motion/actions/simpleactions/HeartBeat.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/actions/simpleactions/HeartBeat.jl
rename to KomaMRIBase/src/motion/actions/simpleactions/HeartBeat.jl
diff --git a/KomaMRIBase/src/motion/motionlist/actions/simpleactions/Rotate.jl b/KomaMRIBase/src/motion/actions/simpleactions/Rotate.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/actions/simpleactions/Rotate.jl
rename to KomaMRIBase/src/motion/actions/simpleactions/Rotate.jl
diff --git a/KomaMRIBase/src/motion/motionlist/actions/simpleactions/Translate.jl b/KomaMRIBase/src/motion/actions/simpleactions/Translate.jl
similarity index 100%
rename from KomaMRIBase/src/motion/motionlist/actions/simpleactions/Translate.jl
rename to KomaMRIBase/src/motion/actions/simpleactions/Translate.jl
diff --git a/KomaMRIBase/src/motion/motionlist/TimeSpan.jl b/KomaMRIBase/src/motion/motionlist/TimeSpan.jl
deleted file mode 100644
index 551e4965a..000000000
--- a/KomaMRIBase/src/motion/motionlist/TimeSpan.jl
+++ /dev/null
@@ -1,143 +0,0 @@
-abstract type AbstractTimeSpan{T<:Real} end
-
-"""
- timerange = TimeRange(t_start, t_end)
-
-TimeRange struct. It is a specialized type that inherits from AbstractTimeSpan and
-defines a time interval, with start and end times.
-
-# Arguments
-- `t_start`: (`::Real`, `[s]`) start time
-- `t_end`: (`::Real`, `[s]`) end time
-
-# Returns
-- `timerange`: (`::TimeRange`) TimeRange struct
-
-# Examples
-```julia-repl
-julia> timerange = TimeRange(0.0, 1.0)
-```
-"""
-@with_kw struct TimeRange{T<:Real} <: AbstractTimeSpan{T}
- t_start ::T
- t_end ::T
- @assert t_end >= t_start "t_end must be greater or equal than t_start"
-end
-
-""" Constructors """
-TimeRange(t_start) = TimeRange(t_start, t_start)
-
-""" times """
-times(ts::TimeRange) = [ts.t_start, ts.t_end]
-
-"""
- t_unit = unit_time(t, time_range)
-
-The `unit_time` function normalizes a given array of time values t
-to a unit interval [0, 1] based on a specified start time `t_start` and end time `t_end`.
-This function is used for non-periodic motions, where each element of t is transformed
-to fit within the range [0, 1] based on the provided start and end times.
-
-
-
-# Arguments
-- `t`: (`::AbstractArray{T<:Real}`, `[s]`) array of time values to be normalized
-- `time_range`: (`::TimeRange{T<:Real}`, `[s]`) time interval (defined by `t_start` and `t_end`) over which we want to normalise
-
-# Returns
-- `t_unit`: (`::AbstractArray{T<:Real}`, `[s]`) array of normalized time values
-
-# Examples
-```julia-repl
-julia> t_unit = KomaMRIBase.unit_time([0.0, 1.0, 2.0, 3.0, 4.0, 5.0], TimeRange(1.0, 4.0))
-6-element Vector{Float64}:
- 0.0
- 0.0
- 0.333
- 0.666
- 1.0
- 1.0
-```
-"""
-function unit_time(t, ts::TimeRange{T}) where {T<:Real}
- if ts.t_start == ts.t_end
- return (t .>= ts.t_start) .* oneunit(T)
- else
- tmp = max.((t .- ts.t_start) ./ (ts.t_end - ts.t_start), zero(T))
- return min.(tmp, oneunit(T))
- end
-end
-
-
-"""
- periodic = Periodic(period, asymmetry)
-
-Periodic struct. It is a specialized type that inherits from AbstractTimeSpan,
-designed to work with time intervals that repeat periodically. It includes a measure of
-asymmetry in order to recreate a asymmetric period.
-
-# Arguments
-- `period`: (`::Real`, `[s]`) period duration
-- `asymmetry`: (`::Real`, `=0.5`) temporal asymmetry factor. Between 0 and 1.
-
-# Returns
-- `periodic`: (`::Periodic`) Periodic struct
-
-# Examples
-```julia-repl
-julia> periodic = Periodic(1.0, 0.2)
-```
-"""
-@with_kw struct Periodic{T<:Real} <: AbstractTimeSpan{T}
- period::T
- asymmetry::T = typeof(period)(0.5)
-end
-
-""" Constructors """
-Periodic(period) = Periodic(period, typeof(period)(0.5))
-
-""" times """
-times(ts::Periodic{T}) where {T<:Real} = [zero(T), ts.period * ts.asymmetry, ts.period]
-
-"""
- t_unit = unit_time(t, periodic)
-
-The `unit_time` function normalizes a given array
-of time values t to a unit interval [0, 1] for periodic motions,
-based on a specified period and an asymmetry factor.
-This function is useful for creating triangular waveforms
-or normalizing time values in periodic processes.
-
-
-
-# Arguments
-- `t`: (`::AbstractArray{T<:Real}`, `[s]`) array of time values to be normalized
-- `periodic`: (`::Periodic{T<:Real}`, `[s]`) information about the `period` and the temporal `asymmetry`
-# Returns
-- `t_unit`: (`::AbstractArray{T<:Real}`, `[s]`) array of normalized time values
-
-# Examples
-```julia-repl
-julia> t_unit = KomaMRIBase.unit_time([0.0, 1.0, 2.0, 3.0, 4.0, 5.0], Periodic(4.0, 0.5))
-6-element Vector{Float64}:
- 0.0
- 0.5
- 1.0
- 0.5
- 0.0
- 0.5
-```
-"""
-function unit_time(t, ts::Periodic{T}) where {T<:Real}
- t_rise = ts.period * ts.asymmetry
- t_fall = ts.period * (oneunit(T) - ts.asymmetry)
- t_relative = mod.(t, ts.period)
- if t_rise == 0
- t_unit = ifelse.(t_relative .< t_rise, zero(T), oneunit(T) .- t_relative ./ t_fall)
- elseif t_fall == 0
- t_unit = ifelse.(t_relative .< t_rise, t_relative ./ t_rise, oneunit(T))
- else
- t_unit = ifelse.( t_relative .< t_rise, t_relative ./ t_rise, oneunit(T) .- (t_relative .- t_rise) ./ t_fall)
- end
- return t_unit
-end
\ No newline at end of file
diff --git a/KomaMRIBase/src/timing/TimeStepCalculation.jl b/KomaMRIBase/src/timing/TimeStepCalculation.jl
index e1b76ff27..d77e383bb 100644
--- a/KomaMRIBase/src/timing/TimeStepCalculation.jl
+++ b/KomaMRIBase/src/timing/TimeStepCalculation.jl
@@ -81,13 +81,15 @@ This function returns non-uniform time points that are relevant in the sequence
samples for RF excitation (by nominal we mean that the time separation should be at most
`Δt_rf` when the samples are regarded by [`KomaMRI.is_RF_on`](@ref), otherwise the time
points are not necessary and the separation will be bigger)
+- `motion`: (`::Union{NoMotion, Motion, MotionList}`, `=NoMotion()`) phantom motion,
+ from which it may be necessary to extract key time points relevant for the simulation
# Returns
- `t`: (`::Vector{Float64}`, `[s]`) time array with non-uniform time values
- `Δt`: (`::Vector{Float64}`, `[s]`) delta time array with the separation between two
adjacent time points of the `t` time array
"""
-function get_variable_times(seq; Δt=1e-3, Δt_rf=1e-5)
+function get_variable_times(seq; Δt=1e-3, Δt_rf=1e-5, motion=NoMotion())
t = Float64[]
ϵ = MIN_RISE_TIME # Small Float64
T0 = get_block_start_times(seq)
@@ -120,6 +122,7 @@ function get_variable_times(seq; Δt=1e-3, Δt_rf=1e-5)
end
append!(t, t_block)
end
+ add_jump_times!(t, motion)
# Removing repeated points
sort!(unique!(t))
# Fixes a problem with ADC at the start and end of the seq
diff --git a/KomaMRIBase/test/runtests.jl b/KomaMRIBase/test/runtests.jl
index 238e74b88..2fddb1abb 100644
--- a/KomaMRIBase/test/runtests.jl
+++ b/KomaMRIBase/test/runtests.jl
@@ -389,7 +389,129 @@ end
@test size(obj1) == size(ρ)
@test length(obj1) == length(ρ)
- # Test obtaining spin positions
+ simplemotion = MotionList(
+ Translate(0.05, 0.05, 0.0, Periodic(period=0.5, asymmetry=0.5)),
+ Rotate(0.0, 0.0, 90.0, TimeRange(t_start=0.05, t_end=0.5), SpinRange(1:3))
+ )
+
+ Ns = length(obj1)
+ Nt = 3
+ t_start = 0.0
+ t_end = 1.0
+ arbitrarymotion = Path(0.01 .* rand(Ns, Nt), 0.01 .* rand(Ns, Nt), 0.01 .* rand(Ns, Nt), TimeRange(t_start, t_end), SpinRange(2:2:4))
+
+ obs1 = Phantom(
+ name,
+ x,
+ y,
+ z,
+ ρ,
+ T1,
+ T2,
+ T2s,
+ Δw,
+ Dλ1,
+ Dλ2,
+ Dθ,
+ simplemotion
+ )
+
+ # Test phantom subset (simple range)
+ rng = 1:2:5
+ obs2 = Phantom(
+ name,
+ x[rng],
+ y[rng],
+ z[rng],
+ ρ[rng],
+ T1[rng],
+ T2[rng],
+ T2s[rng],
+ Δw[rng],
+ Dλ1[rng],
+ Dλ2[rng],
+ Dθ[rng],
+ simplemotion[rng],
+ )
+ @test obs1[rng] == obs2
+ @test @view(obs1[rng]) == obs2
+
+ obs1.motion = arbitrarymotion
+ obs2.motion = arbitrarymotion[rng]
+ @test obs1[rng] == obs2
+
+ # Test addition of phantoms
+ oba = Phantom(
+ name,
+ [x; x[rng]],
+ [y; y[rng]],
+ [z; z[rng]],
+ [ρ; ρ[rng]],
+ [T1; T1[rng]],
+ [T2; T2[rng]],
+ [T2s; T2s[rng]],
+ [Δw; Δw[rng]],
+ [Dλ1; Dλ1[rng]],
+ [Dλ2; Dλ2[rng]],
+ [Dθ; Dθ[rng]],
+ vcat(obs1.motion, obs2.motion, length(obs1), length(obs2))
+ )
+ @test obs1 + obs2 == oba
+
+ # Test phantom subset (BitVector range)
+ obs3 = copy(obs1)
+ obs4 = copy(obs1)
+ rng = obs3.x .> 0
+ obs3.motion = Translate(5e-4, 6e-4, 7e-4, TimeRange(0.0, 1.0), SpinRange(rng))
+ obs4.motion = Translate(5e-4, 6e-4, 7e-4, TimeRange(0.0, 1.0), SpinRange(1:length(obs4)))
+ @test obs3[rng] == obs4[rng]
+ @test obs3[rng].motion == obs4.motion[rng]
+
+ # Test scalar multiplication of a phantom
+ c = 7
+ obc = Phantom(name=name, x=x, y=y, z=z, ρ=c*ρ, T1=T1, T2=T2, T2s=T2s, Δw=Δw, Dλ1=Dλ1, Dλ2=Dλ2, Dθ=Dθ)
+ @test c * obj1 == obc
+
+ #Test brain phantom 2D
+ ph = brain_phantom2D()
+ @test ph.name == "brain2D_axial"
+ @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
+
+ #Test brain phantom 3D
+ ph = brain_phantom3D()
+ @test ph.name == "brain3D"
+ @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 1]
+
+ #Test pelvis phantom 2D
+ ph = pelvis_phantom2D()
+ @test ph.name == "pelvis2D"
+ @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
+
+ #Test heart phantom
+ ph = heart_phantom()
+ @test ph.name == "LeftVentricle"
+ @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
+end
+
+@testitem "Motion" tags=[:base] begin
+ # Test Motion constructors
+ @testset "Constructors" begin
+ action = Rotate(10.0, 20.0, 40.0)
+ spins = AllSpins()
+
+ # TimeCurve constructors
+ time = TimeRange(t_start=0.0, t_end=1.0)
+ time = Periodic(period=1.0, asymmetry=0.5)
+ time = TimeCurve([0.0, eps()], [0.0, 1.0])
+
+ m = Motion(action, time, spins)
+
+ @test Motion(action) == m
+ @test Motion(action, time) == m
+ @test Motion(action, spins) == m
+ end
+
+ # Tests obtaining spin positions
@testset "NoMotion" begin
ph = Phantom(x=[1.0, 2.0], y=[1.0, 2.0])
t_start=0.0; t_end=1.0
@@ -405,19 +527,11 @@ end
t = collect(range(t_start, t_end, 11))
dx, dy, dz = [1.0, 0.0, 0.0]
vx, vy, vz = [dx, dy, dz] ./ (t_end - t_start)
- translation = MotionList(Translate(dx, dy, dz, TimeRange(t_start, t_end)))
+ translation = Translate(dx, dy, dz, TimeRange(t_start, t_end))
xt, yt, zt = get_spin_coords(translation, ph.x, ph.y, ph.z, t')
@test xt == ph.x .+ vx.*t'
@test yt == ph.y .+ vy.*t'
@test zt == ph.z .+ vz.*t'
- # ----- t_start = t_end --------
- t_start = t_end = 0.0
- t = [-0.5, -0.25, 0.0, 0.25, 0.5]
- translation = MotionList(Translate(dx, dy, dz, TimeRange(t_start, t_end)))
- xt, yt, zt = get_spin_coords(translation, ph.x, ph.y, ph.z, t')
- @test xt == ph.x .+ dx*[0, 0, 1, 1, 1]'
- @test yt == ph.y .+ dy*[0, 0, 1, 1, 1]'
- @test zt == ph.z .+ dz*[0, 0, 1, 1, 1]'
end
@testset "PeriodicTranslate" begin
ph = Phantom(x=[1.0], y=[1.0])
@@ -427,7 +541,7 @@ end
asymmetry = 0.5
dx, dy, dz = [1.0, 0.0, 0.0]
vx, vy, vz = [dx, dy, dz] ./ (t_end - t_start)
- periodictranslation = MotionList(Translate(dx, dy, dz, Periodic(period, asymmetry)))
+ periodictranslation = Translate(dx, dy, dz, Periodic(period=period, asymmetry=asymmetry))
xt, yt, zt = get_spin_coords(periodictranslation, ph.x, ph.y, ph.z, t')
@test xt == ph.x .+ vx.*t'
@test yt == ph.y .+ vy.*t'
@@ -440,7 +554,7 @@ end
pitch = 45.0
roll = 0.0
yaw = 45.0
- rotation = MotionList(Rotate(pitch, roll, yaw, TimeRange(t_start, t_end)))
+ rotation = Rotate(pitch, roll, yaw, TimeRange(t_start, t_end))
xt, yt, zt = get_spin_coords(rotation, ph.x, ph.y, ph.z, t')
r = vcat(ph.x, ph.y, ph.z)
R = rotz(π*yaw/180) * roty(π*roll/180) * rotx(π*pitch/180)
@@ -448,14 +562,6 @@ end
@test xt[end ,end] ≈ rot_x
@test yt[end ,end] ≈ rot_y
@test zt[end ,end] ≈ rot_z
- # ----- t_start = t_end --------
- t_start = t_end = 0.0
- t = [-0.5, -0.25, 0.0, 0.25, 0.5]
- rotation = MotionList(Rotate(pitch, roll, yaw, TimeRange(t_start, t_end)))
- xt, yt, zt = get_spin_coords(rotation, ph.x, ph.y, ph.z, t')
- @test xt ≈ [ph.x ph.x rot_x rot_x rot_x]
- @test yt ≈ [ph.y ph.y rot_y rot_y rot_y]
- @test zt ≈ [ph.z ph.z rot_z rot_z rot_z]
end
@testset "PeriodicRotation" begin
ph = Phantom(x=[1.0], y=[1.0])
@@ -466,7 +572,7 @@ end
pitch = 45.0
roll = 0.0
yaw = 45.0
- periodicrotation = MotionList(Rotate(pitch, roll, yaw, Periodic(period, asymmetry)))
+ periodicrotation = Rotate(pitch, roll, yaw, Periodic(period=period, asymmetry=asymmetry))
xt, yt, zt = get_spin_coords(periodicrotation, ph.x, ph.y, ph.z, t')
r = vcat(ph.x, ph.y, ph.z)
R = rotz(π*yaw/180) * roty(π*roll/180) * rotx(π*pitch/180)
@@ -482,26 +588,13 @@ end
circumferential_strain = -0.1
radial_strain = 0.0
longitudinal_strain = -0.1
- heartbeat = MotionList(HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, TimeRange(t_start, t_end)))
+ heartbeat = HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, TimeRange(t_start, t_end))
xt, yt, zt = get_spin_coords(heartbeat, ph.x, ph.y, ph.z, t')
r = sqrt.(ph.x .^ 2 + ph.y .^ 2)
θ = atan.(ph.y, ph.x)
@test xt[:,end] == ph.x .* (1 .+ circumferential_strain * maximum(r) .* cos.(θ))
@test yt[:,end] == ph.y .* (1 .+ circumferential_strain * maximum(r) .* sin.(θ))
@test zt[:,end] == ph.z .* (1 .+ longitudinal_strain)
- # ----- t_start = t_end --------
- t_start = t_end = 0.0
- t = [-0.5, -0.25, 0.0, 0.25, 0.5]
- heartbeat = MotionList(HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, TimeRange(t_start, t_end)))
- xt, yt, zt = get_spin_coords(heartbeat, ph.x, ph.y, ph.z, t')
- r = sqrt.(ph.x .^ 2 + ph.y .^ 2)
- θ = atan.(ph.y, ph.x)
- dx = (1 .+ circumferential_strain * maximum(r) .* cos.(θ))
- dy = (1 .+ circumferential_strain * maximum(r) .* sin.(θ))
- dz = (1 .+ longitudinal_strain)
- @test xt == [ph.x ph.x ph.x .* dx ph.x .* dx ph.x .* dx]
- @test yt == [ph.y ph.y ph.y .* dy ph.y .* dy ph.y .* dy]
- @test zt == [ph.z ph.z ph.z .* dz ph.z .* dz ph.z .* dz]
end
@testset "PeriodicHeartBeat" begin
ph = Phantom(x=[1.0], y=[1.0])
@@ -512,7 +605,7 @@ end
circumferential_strain = -0.1
radial_strain = 0.0
longitudinal_strain = -0.1
- periodicheartbeat = MotionList(HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, Periodic(period, asymmetry)))
+ periodicheartbeat = HeartBeat(circumferential_strain, radial_strain, longitudinal_strain, Periodic(period=period, asymmetry=asymmetry))
xt, yt, zt = get_spin_coords(periodicheartbeat, ph.x, ph.y, ph.z, t')
r = sqrt.(ph.x .^ 2 + ph.y .^ 2)
θ = atan.(ph.y, ph.x)
@@ -530,7 +623,7 @@ end
dx = rand(Ns, Nt)
dy = rand(Ns, Nt)
dz = rand(Ns, Nt)
- arbitrarymotion = MotionList(Path(dx, dy, dz, TimeRange(t_start, t_end)))
+ arbitrarymotion = Path(dx, dy, dz, TimeRange(t_start, t_end))
t = range(t_start, t_end, Nt)
xt, yt, zt = get_spin_coords(arbitrarymotion, ph.x, ph.y, ph.z, t')
@test xt == ph.x .+ dx
@@ -545,7 +638,7 @@ end
dx = rand(Ns, Nt)
dy = rand(Ns, Nt)
dz = rand(Ns, Nt)
- arbitrarymotion = MotionList(Path(dx, dy, dz, TimeRange(t_start, t_end)))
+ arbitrarymotion = Path(dx, dy, dz, TimeRange(t_start, t_end))
t = range(t_start, t_end, Nt)
xt, yt, zt = get_spin_coords(arbitrarymotion, ph.x, ph.y, ph.z, t')
@test xt == ph.x .+ dx
@@ -553,98 +646,38 @@ end
@test zt == ph.z .+ dz
end
- simplemotion = MotionList(
- Translate(0.05, 0.05, 0.0, Periodic(period=0.5, asymmetry=0.5)),
- Rotate(0.0, 0.0, 90.0, TimeRange(t_start=0.05, t_end=0.5), SpinRange(1:3))
- )
-
- Ns = length(obj1)
- Nt = 3
- t_start = 0.0
- t_end = 1.0
- arbitrarymotion = MotionList(Path(0.01 .* rand(Ns, Nt), 0.01 .* rand(Ns, Nt), 0.01 .* rand(Ns, Nt), TimeRange(t_start, t_end), SpinRange(2:2:4)))
-
- # Test phantom subset
- obs1 = Phantom(
- name,
- x,
- y,
- z,
- ρ,
- T1,
- T2,
- T2s,
- Δw,
- Dλ1,
- Dλ2,
- Dθ,
- simplemotion
- )
- rng = 1:2:5
- obs2 = Phantom(
- name,
- x[rng],
- y[rng],
- z[rng],
- ρ[rng],
- T1[rng],
- T2[rng],
- T2s[rng],
- Δw[rng],
- Dλ1[rng],
- Dλ2[rng],
- Dθ[rng],
- simplemotion[rng],
- )
- @test obs1[rng] == obs2
- @test @view(obs1[rng]) == obs2
-
- obs1.motion = arbitrarymotion
- obs2.motion = arbitrarymotion[rng]
- @test obs1[rng] == obs2
-
- # Test addition of phantoms
- oba = Phantom(
- name,
- [x; x[rng]],
- [y; y[rng]],
- [z; z[rng]],
- [ρ; ρ[rng]],
- [T1; T1[rng]],
- [T2; T2[rng]],
- [T2s; T2s[rng]],
- [Δw; Δw[rng]],
- [Dλ1; Dλ1[rng]],
- [Dλ2; Dλ2[rng]],
- [Dθ; Dθ[rng]],
- vcat(obs1.motion, obs2.motion, length(obs1), length(obs2))
- )
- @test obs1 + obs2 == oba
-
- # Test scalar multiplication of a phantom
- c = 7
- obc = Phantom(name=name, x=x, y=y, z=z, ρ=c*ρ, T1=T1, T2=T2, T2s=T2s, Δw=Δw, Dλ1=Dλ1, Dλ2=Dλ2, Dθ=Dθ)
- @test c * obj1 == obc
-
- #Test brain phantom 2D
- ph = brain_phantom2D()
- @test ph.name == "brain2D_axial"
- @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
-
- #Test brain phantom 3D
- ph = brain_phantom3D()
- @test ph.name == "brain3D"
- @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 1]
-
- #Test pelvis phantom 2D
- ph = pelvis_phantom2D()
- @test ph.name == "pelvis2D"
- @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
-
- #Test heart phantom
- ph = heart_phantom()
- @test ph.name == "LeftVentricle"
- @test KomaMRIBase.get_dims(ph) == Bool[1, 1, 0]
+ @testset "FlowPath" begin
+ # 1 spin
+ ph = Phantom(x=[1.0], y=[1.0])
+ Ns = length(ph)
+ t_start = 0.0
+ t_end = 1.0
+ Nt = 10
+ dx = rand(Ns, Nt)
+ dy = rand(Ns, Nt)
+ dz = rand(Ns, Nt)
+ arbitrarymotion = FlowPath(dx, dy, dz, collect(Bool.(zeros(Ns, Nt))), TimeRange(t_start, t_end))
+ t = range(t_start, t_end, Nt)
+ xt, yt, zt = get_spin_coords(arbitrarymotion, ph.x, ph.y, ph.z, t')
+ @test xt == ph.x .+ dx
+ @test yt == ph.y .+ dy
+ @test zt == ph.z .+ dz
+ # More than 1 spin
+ ph = Phantom(x=[1.0, 2.0], y=[1.0, 2.0])
+ Ns = length(ph)
+ t_start = 0.0
+ t_end = 1.0
+ Nt = 10
+ dx = rand(Ns, Nt)
+ dy = rand(Ns, Nt)
+ dz = rand(Ns, Nt)
+ arbitrarymotion = FlowPath(dx, dy, dz, collect(Bool.(zeros(Ns, Nt))), TimeRange(t_start, t_end))
+ t = range(t_start, t_end, Nt)
+ xt, yt, zt = get_spin_coords(arbitrarymotion, ph.x, ph.y, ph.z, t')
+ @test xt == ph.x .+ dx
+ @test yt == ph.y .+ dy
+ @test zt == ph.z .+ dz
+ end
end
@testitem "Scanner" tags=[:base] begin
diff --git a/KomaMRICore/Project.toml b/KomaMRICore/Project.toml
index 8be0ccb09..82fcb819d 100644
--- a/KomaMRICore/Project.toml
+++ b/KomaMRICore/Project.toml
@@ -28,7 +28,7 @@ KomaoneAPIExt = "oneAPI"
AMDGPU = "0.9, 1"
Adapt = "3, 4"
CUDA = "3, 4, 5"
-Functors = "0.4"
+Functors = "0.4, 0.5"
KernelAbstractions = "0.9"
KomaMRIBase = "0.9"
Metal = "1.2"
diff --git a/KomaMRICore/src/rawdata/ISMRMRD.jl b/KomaMRICore/src/rawdata/ISMRMRD.jl
index 1ee1ae048..369fd05fe 100644
--- a/KomaMRICore/src/rawdata/ISMRMRD.jl
+++ b/KomaMRICore/src/rawdata/ISMRMRD.jl
@@ -154,7 +154,7 @@ function signal_to_raw_data(
profiles = Profile[]
t_acq = get_adc_sampling_times(seq)
Nadcs = sum(is_ADC_on.(seq))
- NadcsPerImage = floor(Int, Nadcs / Nz)
+ NadcsPerImage = max(floor(Int, Nadcs / Nz), 1)
scan_counter = 0
nz = 0
current = 1
diff --git a/KomaMRICore/src/simulation/Flow.jl b/KomaMRICore/src/simulation/Flow.jl
index 55f4be49b..a9727b95d 100644
--- a/KomaMRICore/src/simulation/Flow.jl
+++ b/KomaMRICore/src/simulation/Flow.jl
@@ -81,13 +81,7 @@ function replace_view(replace_by, idx)
return replace_by
end
-function get_mask(spin_reset, t::Real)
+function get_mask(spin_reset, t)
itp = KomaMRIBase.interpolate(spin_reset, KomaMRIBase.Gridded(KomaMRIBase.Constant{KomaMRIBase.Next}()), Val(size(spin_reset, 1)), t)
return KomaMRIBase.resample(itp, t)
-end
-function get_mask(spin_reset, t::AbstractArray)
- itp = KomaMRIBase.interpolate(spin_reset, KomaMRIBase.Gridded(KomaMRIBase.Constant{KomaMRIBase.Next}()), Val(size(spin_reset, 1)), t)
- mask = KomaMRIBase.resample(itp, t)
- mask .= (cumsum(mask; dims=2) .== 1)
- return mask
end
\ No newline at end of file
diff --git a/KomaMRICore/src/simulation/Functors.jl b/KomaMRICore/src/simulation/Functors.jl
index 8434818c7..f24b57681 100644
--- a/KomaMRICore/src/simulation/Functors.jl
+++ b/KomaMRICore/src/simulation/Functors.jl
@@ -110,8 +110,7 @@ adapt_storage(T::Type{<:Real}, xs::MotionList) = MotionList(paramtype.(T, xs.mot
@functor HeartBeat
@functor Path
@functor FlowPath
-@functor TimeRange
-@functor Periodic
+@functor TimeCurve
# Spinor
@functor Spinor
# DiscreteSequence
diff --git a/KomaMRICore/src/simulation/SimulatorCore.jl b/KomaMRICore/src/simulation/SimulatorCore.jl
index 6ba42e583..aff8a56f7 100644
--- a/KomaMRICore/src/simulation/SimulatorCore.jl
+++ b/KomaMRICore/src/simulation/SimulatorCore.jl
@@ -336,7 +336,7 @@ function simulate(
""" maxlog=1
end
# Simulation init
- seqd = discretize(seq; sampling_params=sim_params) # Sampling of Sequence waveforms
+ seqd = discretize(seq; sampling_params=sim_params, motion=obj.motion) # Sampling of Sequence waveforms
parts, excitation_bool = get_sim_ranges(seqd; Nblocks=sim_params["Nblocks"]) # Generating simulation blocks
t_sim_parts = [seqd.t[p[1]] for p in parts]
append!(t_sim_parts, seqd.t[end])
diff --git a/KomaMRICore/test/runtests.jl b/KomaMRICore/test/runtests.jl
index be10f384f..022b281cd 100644
--- a/KomaMRICore/test/runtests.jl
+++ b/KomaMRICore/test/runtests.jl
@@ -441,7 +441,7 @@ end
motions = [
Translate(0.1, 0.1, 0.0, TimeRange(0.0, 1.0)),
Rotate(0.0, 0.0, 45.0, TimeRange(0.0, 1.0)),
- HeartBeat(-0.6, 0.0, 0.0, Periodic(1.0)),
+ HeartBeat(-0.6, 0.0, 0.0, Periodic(period=1.0)),
Path([0.0 0.0], [0.0 1.0], [0.0 0.0], TimeRange(0.0, 10.0)),
FlowPath([0.0 0.0], [0.0 1.0], [0.0 0.0], [false false], TimeRange(0.0, 10.0))
]
diff --git a/KomaMRIFiles/src/Phantom/Phantom.jl b/KomaMRIFiles/src/Phantom/Phantom.jl
index 601e46e5c..fe0c796aa 100644
--- a/KomaMRIFiles/src/Phantom/Phantom.jl
+++ b/KomaMRIFiles/src/Phantom/Phantom.jl
@@ -47,7 +47,7 @@ function import_motion!(phantom_fields::Array, motion_group::HDF5.Group)
end
push!(motion_array, Motion(; motion_fields...))
end
- push!(phantom_fields, (:motion, MotionList(motion_array)))
+ push!(phantom_fields, (:motion, MotionList(motion_array...)))
end
function import_motion_field!(motion_fields::Array, motion::HDF5.Group, name::String)
@@ -57,19 +57,17 @@ function import_motion_field!(motion_fields::Array, motion::HDF5.Group, name::St
get_subtypes(t::Type) = reduce(vcat,(subtypes(t)))
get_subtype_strings(t::Type) = last.(split.(string.(get_subtypes(t::Type)), "."))
- subtype_strings = reduce(vcat, get_subtype_strings.([
+ subtype_strings = [reduce(vcat, get_subtype_strings.([
KomaMRIBase.SimpleAction,
KomaMRIBase.ArbitraryAction,
- KomaMRIBase.AbstractTimeSpan,
KomaMRIBase.AbstractSpinSpan
- ]))
+ ])); "TimeCurve"]
- subtype_vector = reduce(vcat, get_subtypes.([
+ subtype_vector = [reduce(vcat, get_subtypes.([
KomaMRIBase.SimpleAction,
KomaMRIBase.ArbitraryAction,
- KomaMRIBase.AbstractTimeSpan,
KomaMRIBase.AbstractSpinSpan
- ]))
+ ])); KomaMRIBase.TimeCurve]
motion_subfields = []
for (i, subtype_string) in enumerate(subtype_strings)
@@ -89,6 +87,9 @@ function import_motion_subfield!(motion_subfields::Array, subfield_value::Union{
return nothing
end
function import_motion_subfield!(motion_subfields::Array, subfield_value::String, key::String)
+ if subfield_value in ["true", "false"]
+ return push!(motion_subfields, subfield_value == "true" ? true : false)
+ end
endpoints = parse.(Int, split(subfield_value, ":"))
range = length(endpoints) == 3 ? (endpoints[1]:endpoints[2]:endpoints[3]) : (endpoints[1]:endpoints[2])
push!(motion_subfields, range)
@@ -127,13 +128,14 @@ function write_phantom(
contrast[String(x)] = getfield(obj, x)
end
# Motion
- if (typeof(obj.motion) <: MotionList) & store_motion
+ if !(obj.motion isa NoMotion) & store_motion
motion_group = create_group(fid, "motion")
export_motion!(motion_group, obj.motion)
end
return close(fid)
end
+export_motion!(motion_group::HDF5.Group, motion::Motion) = export_motion!(motion_group, MotionList([motion]))
function export_motion!(motion_group::HDF5.Group, motion_list::MotionList)
KomaMRIBase.sort_motions!(motion_list)
for (counter, m) in enumerate(motion_list.motions)
@@ -165,4 +167,7 @@ function export_motion_subfield!(field_group::HDF5.Group, subfield::Array, subna
end
function export_motion_subfield!(field_group::HDF5.Group, subfield::BitMatrix, subname::String)
field_group[subname] = Int.(subfield)
+end
+function export_motion_subfield!(field_group::HDF5.Group, subfield::Bool, subname::String)
+ field_group[subname] = string(subfield)
end
\ No newline at end of file
diff --git a/KomaMRIPlots/Project.toml b/KomaMRIPlots/Project.toml
index eddc20610..4dc4aa2b9 100644
--- a/KomaMRIPlots/Project.toml
+++ b/KomaMRIPlots/Project.toml
@@ -9,6 +9,7 @@ Kaleido_jll = "f7e6163d-2fa5-5f23-b69c-1db539e41963"
KomaMRIBase = "d0bc0b20-b151-4d03-b2a4-6ca51751cb9c"
MAT = "23992714-dd62-5051-b70f-ba57cb901cac"
PlotlyJS = "f0f68f2c-4968-5e81-91da-67840de0976a"
+QMRIColors = "2bec176e-9e8c-4764-a62f-295118d1ec05"
Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
[weakdeps]
@@ -23,7 +24,8 @@ Kaleido_jll = "0.1"
KomaMRIBase = "0.9"
MAT = "0.10"
PlotlyJS = "0.18"
-PlutoPlotly = "0.4, 0.5"
+PlutoPlotly = "0.4, 0.5, 0.6"
+QMRIColors = "1"
Reexport = "1"
julia = "1.9"
diff --git a/KomaMRIPlots/src/KomaMRIPlots.jl b/KomaMRIPlots/src/KomaMRIPlots.jl
index ad46a4c24..5d8190f93 100644
--- a/KomaMRIPlots/src/KomaMRIPlots.jl
+++ b/KomaMRIPlots/src/KomaMRIPlots.jl
@@ -2,6 +2,7 @@ module KomaMRIPlots
using KomaMRIBase
using MAT, Interpolations, PlotlyJS
+using QMRIColors
include("ui/PlotBackends.jl")
include("ui/DisplayFunctions.jl")
diff --git a/KomaMRIPlots/src/ui/DisplayFunctions.jl b/KomaMRIPlots/src/ui/DisplayFunctions.jl
index 696418669..7ed1e63c4 100644
--- a/KomaMRIPlots/src/ui/DisplayFunctions.jl
+++ b/KomaMRIPlots/src/ui/DisplayFunctions.jl
@@ -980,7 +980,7 @@ function plot_kspace(seq::Sequence; width=nothing, height=nothing, darkmode=fals
"resetCameraLastSave3d",
"orbitRotation",
"resetCameraDefault3d",
- ],
+ ],
)
return plot_koma(p, l; config)
end
@@ -1003,7 +1003,6 @@ Plots a phantom map for a specific spin parameter given by `key`.
- `colorbar`: (`::Bool`, `=true`) boolean to indicate whether to display a colorbar
- `max_spins`:(`::Int`, `=100_000`) maximum number of displayed spins
- `time_samples`:(`::Int`, `=0`) intermediate time samples between motion `t_start` and `t_end`
-- `max_time_samples`:(`::Int`, `=100`) maximum number of time samples
- `frame_duration_ms`:(`::Int`, `=250`) time in miliseconds between two frames
# Returns
@@ -1031,27 +1030,37 @@ function plot_phantom_map(
view_2d=sum(KomaMRIBase.get_dims(obj)) < 3,
colorbar=true,
max_spins=20_000,
- time_samples=0,
- max_time_samples=100,
+ time_samples= obj.motion isa NoMotion ? 0 : length(times(obj.motion)),
kwargs...,
)
function process_times(::NoMotion)
return [zero(eltype(obj.x))]
end
- function process_times(motion::MotionList)
- KomaMRIBase.sort_motions!(motion)
- t = times(motion)
- if length(t)>1
- # Interpolate time points (as many as indicated by time_samples)
- itp = interpolate((1:(time_samples + 1):(length(t) + time_samples * (length(t) - 1)), ), t, Gridded(Linear()))
- t = itp.(1:(length(t) + time_samples * (length(t) - 1)))
+ function process_times(motion)
+ time_nodes = times(motion)
+ if length(time_nodes) > 1
+ time_min, time_max = minimum(time_nodes), maximum(time_nodes)
+ t = collect(range(time_min, time_max, length=time_samples))
+ merged_times = sort(unique([t; time_nodes]))
+ if length(merged_times) > time_samples
+ to_keep = time_nodes
+ remaining = setdiff(merged_times, to_keep)
+ num_needed = time_samples - length(to_keep)
+ extra = remaining[round.(Int, range(1, length(remaining), length=num_needed))]
+ merged_times = sort([to_keep; extra])
+ elseif length(merged_times) < time_samples
+ extra_points = setdiff(t, merged_times)
+ extra_needed = time_samples - length(merged_times)
+ additional = sort(extra_points)[1:extra_needed]
+ merged_times = sort([merged_times; additional])
+ end
+ return merged_times
+ else
+ return time_nodes
end
- # Decimate time points so their number is smaller than max_time_samples
- ss = length(t) > max_time_samples ? length(t) ÷ max_time_samples : 1
- return t[1:ss:end]
end
-
+
function decimate_uniform_phantom(obj, num_points::Int)
ss = Int(ceil(length(obj) / num_points))
return obj[1:ss:end]
@@ -1071,38 +1080,20 @@ function plot_phantom_map(
unit = " ms"
if key == :T1
cmax_key = 2500 / factor
- colors = MAT.matread(path * "/assets/T1cm.mat")["T1colormap"][1:70:end, :]
- N, _ = size(colors)
- idx = range(0, 1; length=N) #range(0,T,N) works in Julia 1.7
- colormap = [
- (
- idx[n],
- string("rgb(",
- floor(Int, colors[n,1] * 255), ",",
- floor(Int, colors[n,2] * 255), ",",
- floor(Int, colors[n,3] * 255), ")"
- )
- )
- for n in 1:N
- ]
+ colors =
+ replace.(string.(relaxationColorMap("T1") .* 255), "RGB{Float64}" => "rgb")
+ N = length(colors)
+ indices = range(0.0; stop=1.0, length=N)
+ colormap = [(idx, color) for (idx, color) in zip(indices, colors)]
elseif key == :T2 || key == :T2s
if key == :T2
cmax_key = 250 / factor
end
- colors = MAT.matread(path * "/assets/T2cm.mat")["T2colormap"][1:70:end, :]
- N, _ = size(colors)
- idx = range(0, 1; length=N) #range(0,T,N) works in Julia 1.7
- colormap = [
- (
- idx[n],
- string("rgb(",
- floor(Int, colors[n,1] * 255), ",",
- floor(Int, colors[n,2] * 255), ",",
- floor(Int, colors[n,3] * 255), ")"
- )
- )
- for n in 1:N
- ]
+ colors =
+ replace.(string.(relaxationColorMap("T2") .* 255), "RGB{Float64}" => "rgb")
+ N = length(colors)
+ indices = range(0.0; stop=1.0, length=N)
+ colormap = [(idx, color) for (idx, color) in zip(indices, colors)]
end
elseif key == :x || key == :y || key == :z
factor = 1e2
diff --git a/docs/make.jl b/docs/make.jl
index fae44c450..b096cf139 100644
--- a/docs/make.jl
+++ b/docs/make.jl
@@ -63,7 +63,7 @@ makedocs(;
prettyurls=true,
sidebar_sitename=false,
collapselevel=1,
- assets=["assets/hide-documenter-example-output.css"],
+ assets=["assets/hide-documenter-example-output.css","assets/center-images.css"],
),
clean=false,
)
diff --git a/docs/src/assets/center-images.css b/docs/src/assets/center-images.css
new file mode 100644
index 000000000..e8940f978
--- /dev/null
+++ b/docs/src/assets/center-images.css
@@ -0,0 +1,5 @@
+img {
+ display: block; /* Ensures the image behaves as a block element */
+ margin-left: auto;
+ margin-right: auto;
+}
\ No newline at end of file
diff --git a/docs/src/assets/create-remote-step1.png b/docs/src/assets/create-remote-step1.png
old mode 100644
new mode 100755
index b8c18e218..da3917a0c
Binary files a/docs/src/assets/create-remote-step1.png and b/docs/src/assets/create-remote-step1.png differ
diff --git a/docs/src/assets/create-remote-step2.png b/docs/src/assets/create-remote-step2.png
old mode 100644
new mode 100755
index 722ce1331..7e1b7564d
Binary files a/docs/src/assets/create-remote-step2.png and b/docs/src/assets/create-remote-step2.png differ
diff --git a/docs/src/assets/periodic.svg b/docs/src/assets/periodic.svg
new file mode 100644
index 000000000..f84409d10
--- /dev/null
+++ b/docs/src/assets/periodic.svg
@@ -0,0 +1,303 @@
+
+
+
+
diff --git a/docs/src/assets/time-curve-1.svg b/docs/src/assets/time-curve-1.svg
new file mode 100644
index 000000000..a637ec4db
--- /dev/null
+++ b/docs/src/assets/time-curve-1.svg
@@ -0,0 +1,327 @@
+
+
+
+
diff --git a/docs/src/assets/time-curve-2.svg b/docs/src/assets/time-curve-2.svg
new file mode 100644
index 000000000..101da15a5
--- /dev/null
+++ b/docs/src/assets/time-curve-2.svg
@@ -0,0 +1,425 @@
+
+
+
+
diff --git a/docs/src/assets/time-curve-3.svg b/docs/src/assets/time-curve-3.svg
new file mode 100644
index 000000000..d69b1d68a
--- /dev/null
+++ b/docs/src/assets/time-curve-3.svg
@@ -0,0 +1,405 @@
+
+
+
+
diff --git a/docs/src/assets/time-curve-4.svg b/docs/src/assets/time-curve-4.svg
new file mode 100644
index 000000000..ab64f7145
--- /dev/null
+++ b/docs/src/assets/time-curve-4.svg
@@ -0,0 +1,425 @@
+
+
+
+
diff --git a/docs/src/assets/time-range.svg b/docs/src/assets/time-range.svg
new file mode 100644
index 000000000..a80537e25
--- /dev/null
+++ b/docs/src/assets/time-range.svg
@@ -0,0 +1,299 @@
+
+
+
+
diff --git a/docs/src/assets/unit-time-triangular.svg b/docs/src/assets/unit-time-triangular.svg
deleted file mode 100644
index 5a91bd907..000000000
--- a/docs/src/assets/unit-time-triangular.svg
+++ /dev/null
@@ -1,355 +0,0 @@
-
-
-
-
diff --git a/docs/src/assets/unit-time.svg b/docs/src/assets/unit-time.svg
deleted file mode 100644
index c7cc6d3f9..000000000
--- a/docs/src/assets/unit-time.svg
+++ /dev/null
@@ -1,405 +0,0 @@
-
-
-
-
diff --git a/docs/src/how-to/5-contribute-to-koma.md b/docs/src/how-to/5-contribute-to-koma.md
index b8726d6bf..7507d0ad8 100644
--- a/docs/src/how-to/5-contribute-to-koma.md
+++ b/docs/src/how-to/5-contribute-to-koma.md
@@ -13,7 +13,7 @@ If you're interested in contributing to Koma, this document will guide you throu
### 1. Clone KomaMRI repository
To install the dev version of Koma, we will use the Julia REPL:
-```julia-repl
+```julia
pkg> dev KomaMRI
```
This command will clone KomaMRI.jl's repository (`dev` version) to your `~/.julia/dev/KomaMRI/` directory if you are in a MacOS or Linux operative system, or `C:\Users\\.julia\dev\KomaMRI\` if you are using Windows, where `` should be replaced with your Windows user.
@@ -36,21 +36,20 @@ Now, you need to ensure that your GitHub account is connected to VSCode. This al
- Sign in to your GitHub account if you're not already signed in.
>💡You can also check if your `git` credentials are correctlly added to your machine by writing in the VScode terminal:
-```shell
-git config --global user.name
-git config --global user.email
-```
+>```shell
+>git config --global user.name
+>git config --global user.email
+>```
### 4. Open your forked repository in VSCode
In VSCode, click on **File** -> **Open Folder...** and select your `~/.julia/dev/KomaMRI/` directory (`C:\Users\\.julia\dev\KomaMRI\` if you are using Windows).
-
-
Now add the fork URL by clicking **Source Control** -> **...** -> **Remote** -> **Add Remote...**
-
-
+```@raw html
+ 
+```
This will create the option to provide a repository URL. Here is where you will paste your fork URL and give it the name `my-fork`.
@@ -62,13 +61,29 @@ The Julia extension should automatically detect the `KomaMRI` environment. To ch
### 5. KomaMRI monorepo setup
-As KomaMRI.jl contains multiple packages in one GitHub repository, you need to specify that you want to use your local copies (instead of the ones available on the Julia registries) with:
-```shell
-(KomaMRI) pkg> dev ./KomaMRIBase ./KomaMRICore ./KomaMRIFiles ./KomaMRIPlots
+As KomaMRI.jl contains multiple packages in one GitHub repository, you need to specify that you want to use your local copies (instead of the ones available on the Julia registries) and using the `instantiate` command to install all the required packages (specified in `Project.toml`) with the following script:
+
+```julia
+using Pkg
+# Koma sub-packages dev setup
+koma_subpkgs = ["KomaMRICore", "KomaMRIFiles", "KomaMRIPlots"]
+for pkg in koma_subpkgs
+ Pkg.activate(pkg)
+ Pkg.develop(path = "./KomaMRIBase")
+end
+# Main package (KomaMRI) dev setup
+Pkg.activate(".")
+for pkg in koma_subpkgs
+ Pkg.develop(path = "./$pkg")
+end
+Pkg.instantiate()
```
-Finally, use the `instantiate` command to install all the required packages (specified in the `Project.toml`):
-```shell
-(KomaMRI) pkg> instantiate
+In case you want to contribute specifically in documentation, you will need to use the `docs` enviroment with the following script:
+
+```julia
+Pkg.activate("docs")
+Pkg.develop(path = ".")
+Pkg.instantiate()
```
This will also include all the specific package versions into the `Manifest.toml`. The `Manifest.toml` should not be updated to the repo when making a commit or pull request. Thus, it is present in the `.gitignore`.
@@ -79,8 +94,9 @@ If you did correctly follow the previous steps you will have correctly created y
To create this new branch, go to **Source Control** -> **...** -> **Branch** -> **Create Branch form...**
-
-
+```@raw html
+ 
+```
This will open a menu to select an starting point for your branch. Select `my-fork/master` as your starting point, and give it the name `my-new-feature`.
@@ -96,8 +112,9 @@ To do this, in VScode go to the Source Control panel in the Activity Bar.
Assuming you are currently in your `my-new-feature` branch, the Source Control panel should show your changes to the project and the option to create a commit message.
-
-
+```@raw html
+ 
+```
If you hove over the `Changes` tab, it should show a `+` icon. Press it to stage all changes in the project.
Write down a message that describes your changes you are stageing to the project, and press the Commit button.
@@ -112,14 +129,17 @@ If you want to send your commited new version of the repository, you can create
To create this pull request, in VScode, go to the `GitHub Pull Request` panel in the Activity Bar and hove over the `Pull request` tab. This should show a Create pull request icon to press.
-
+```@raw html
+ 
+```
In the `Create` tab that appears, select `JuliaHealth/master` as the base and the branch you are working with to merge.
-To finish your pull request, give your pull request a description that explains the issue or feature you are addresing in your branch, and press the Create button.
-
-
+To finish your pull request, give it a name with a clear mention of the subject of the contribution you made, and a description that explains the issue or feature you are addresing in your branch, and press the Create button.
+```@raw html
+ 
+```
>💡 **Tips for a successful Pull Request:**
> - Try to address one issue or feature per pull request to make it easier for reviewers.
> - Provide all the context necesary, including all the information of the related issue or added feature.
diff --git a/docs/src/reference/2-koma-base.md b/docs/src/reference/2-koma-base.md
index 2fc69cf23..5cc479ff5 100644
--- a/docs/src/reference/2-koma-base.md
+++ b/docs/src/reference/2-koma-base.md
@@ -44,12 +44,12 @@ FlowPath
FlowPath(dx, dy, dz, spin_reset, time, spins)
```
-### `AbstractTimeSpan` types and related functions
+### `TimeCurve` types and related functions
```@docs
+TimeCurve
TimeRange
-Periodic
-unit_time
+Periodic
```
### `AbstractSpinSapn` types