JuliaHealth/KomaMRI.jl

New motion approach to combine SimpleMotion and ArbitraryMotion into the same phantom

Opened this issue · 3 comments

Idea for new motion related types and structs redistribution, in order to be able to combine them and make them more alike:

abstract type Motion end

abstract type SimpleMotion <: Motion end
abstract type ArbitraryMotion <: Motion end

Specific motion structures:

Translation <: SimpleMotion
Rotation <: SimpleMotion
...
PeriodicHeartBeat <: SimpleMotion
Trajectory <: ArbitraryMotion
PeriodicTrajectory <: ArbitraryMotion
FlowTrajectory <: ArbitraryMotion
PeriodicFlowTrajectory <: ArbitraryMotion

That way, phantom struct will look like this:

@with_kw mutable struct Phantom{T<:Real}
    ...
    motion::Vector{<:Motion{T}} = []
    ...
end

and motion will be initialized, for example, like this:

motion = [Translation(...), Rotation(...), Trajectory(...), HeartBeat(...)]

This is a very good idea! We just need to think about how this will impact:

  • obj1 + obj2: obj1 has motion1 and obj2 has motion2 (for example, one ball going up, and another one going down). For ArbitraryMotion's (now Trajectory) concatenating with zeros will work, for SimpleMotion's I don't know.
  • obj.phantom: phantom file read/write

Another possible approach regarding periodicity:

  • KomaMRIBase/src/timing/Periodicity.jl:
abstract type TimeScale{T<:Real} end

struct Periodic{T<:Real} <: TimeScale{T}
   period::T
   asymmetry::T
end

struct TimeRange{T<:Real} <: TimeScale{T} 
   t_start::T
   t_end::T
end
  • KomaMRIBase/src/datatypes/phantom/motion/simplemotion/Translation.jl:
struct Translation{T<:Real, TS<:TimeScale{T}} <: SimpleMotionType{T}
   dx::T
   dy::T
   dz::T
   times::TS
end

This way we could avoid re-defining structs and methods for periodic types.
A SimpleMotiontype, for example, could be instantiated as follows:
tr = Translation(dx, dy, dz, TimeRange(t_start, t_end))

Regarding flow:

function reset_magnetization!(M::Mag{T}, Mxy::AbstractArray{T}, motion<:Vector{Motion{T}})
   for m in motion
      reset_magnetization!(M, Mxy, m)
   end
   return nothing
end

function reset_magnetization!(M::Mag{T}, Mxy::AbstractArray{T}, motion<:Motion{T})
   return nothing
end

function reset_magnetization!(M::Mag{T}, Mxy::AbstractArray{T}, motion::FlowTrajectory{T})
   # do things with flow (create interpolation functions, resolve them, etc):
   Mxy .= ...
   M.z .= ...
   return nothing
end

This will be called from simulation functions:

function run_spin_precession!(...)
   ...
   #Mxy precession and relaxation, and Mz relaxation
   tp = cumsum(seq.Δt) # t' = t - t0
   dur = sum(seq.Δt)   # Total length, used for signal relaxation
   Mxy = [M.xy M.xy .* exp.(1im .* ϕ .- tp' ./ p.T2)] #This assumes Δw and T2 are constant in time
   M.z .= M.z .* exp.(-dur ./ p.T1) .+ p.ρ .* (1 .- exp.(-dur ./ p.T1))
   reset_magnetization!(M, Mxy, phantom.motion)
   M.xy .= Mxy[:, end]
   #Acquired signal
   sig .= transpose(sum(Mxy[:, findall(seq.ADC)]; dims=1)) #<--- TODO: add coil sensitivities
return nothing