Try changing preallocations in get_spin_coords

JuliaHealth · Sep 2, 2024 · 00d7ea8 · 00d7ea8
1 parent f77f59e
commit 00d7ea8
Show file tree

Hide file tree

Showing 3 changed files with 43 additions and 34 deletions.
diff --git a/KomaMRIBase/src/motion/motionlist/MotionList.jl b/KomaMRIBase/src/motion/motionlist/MotionList.jl
@@ -109,27 +109,21 @@ function get_spin_coords(
 ) where {T<:Real}
     # Buffers for positions:
     xt, yt, zt = x .+ 0*t, y .+ 0*t, z .+ 0*t
-    # Buffers for displacements:
-    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)
         idx = get_idx(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)
-        displacement_z!(@view(uz[idx, :]), m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
-        xt .+= ux; yt .+= uy; zt .+= uz
-        ux .*= zero(T); uy .*= zero(T); uz .*= zero(T)
+        @view(xt[idx, :]) .+= displacement_x(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
+        @view(yt[idx, :]) .+= displacement_y(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
+        @view(zt[idx, :]) .+= displacement_z(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
     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)
         idx = get_idx(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)
-        xt .+= ux; yt .+= uy; zt .+= uz
-        ux .*= zero(T); uy .*= zero(T); uz .*= zero(T)
+        @view(xt[idx, :]) .+= displacement_x(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
+        @view(yt[idx, :]) .+= displacement_y(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
+        @view(zt[idx, :]) .+= displacement_z(m.action, @view(xt[idx, :]), @view(yt[idx, :]), @view(zt[idx, :]), t_unit)
     end
     return xt, yt, zt
 end

diff --git a/KomaMRIBase/src/motion/motionlist/actions/ArbitraryAction.jl b/KomaMRIBase/src/motion/motionlist/actions/ArbitraryAction.jl
@@ -98,5 +98,41 @@ function displacement_z!(
     return nothing
 end
 
+
+
+
+function displacement_x(
+    action::ArbitraryAction{T},
+    x::AbstractArray{T},
+    y::AbstractArray{T},
+    z::AbstractArray{T},
+    t::AbstractArray{T},
+) where {T<:Real}
+    itp = interpolate(action.dx, Gridded(Linear()), Val(size(action.dx,1)))
+    return resample(itp, t)
+end
+
+function displacement_y(
+    action::ArbitraryAction{T},
+    x::AbstractArray{T},
+    y::AbstractArray{T},
+    z::AbstractArray{T},
+    t::AbstractArray{T},
+) where {T<:Real}
+    itp = interpolate(action.dy, Gridded(Linear()), Val(size(action.dy,1)))
+    return resample(itp, t)
+end
+
+function displacement_z(
+    action::ArbitraryAction{T},
+    x::AbstractArray{T},
+    y::AbstractArray{T},
+    z::AbstractArray{T},
+    t::AbstractArray{T},
+) where {T<:Real}
+    itp = interpolate(action.dz, Gridded(Linear()), Val(size(action.dz,1)))
+    return resample(itp, t)
+end
+
 include("arbitraryactions/Path.jl")
 include("arbitraryactions/FlowPath.jl")
diff --git a/KomaMRICore/test/runtests.jl b/KomaMRICore/test/runtests.jl
@@ -469,33 +469,12 @@ end
     include("initialize_backend.jl")
     include(joinpath(@__DIR__, "test_files", "utils.jl"))
 
-    for i in 1:8
+    for i in 1:15
         sig_jemris = signal_brain_motion_jemris()
         seq = seq_epi_100x100_TE100_FOV230()
         sys = Scanner()
         obj = phantom_brain_arbitrary_motion()
 
-        # vx = 0.0f0
-        # vy = 0.1f0
-        # vz = 0.0f0
-        # t = collect(0:0.1:10)
-
-        # obj = obj |> f32 |> gpu 
-        # t   = t   |> f32 |> gpu
-
-        # x, y, z = get_spin_coords(obj.motion, obj.x, obj.y, obj.z, t')
-
-        # obj = obj |> cpu
-        # t   = t   |> cpu
-
-        # x = x |> cpu
-        # y = y |> cpu
-        # z = z |> cpu
-
-        # @test x ≈ (obj.x .+ vx .* t')
-        # @test y ≈ (obj.y .+ vy .* t')
-        # @test z ≈ (obj.z .+ vz .* t')
-
         sim_params = Dict{String, Any}(
             "gpu"=>USE_GPU,
             "sim_method"=>KomaMRICore.Bloch(),