From cdb7da2186c852e47577d0c5066f3b198fc0fc30 Mon Sep 17 00:00:00 2001
From: Michiel Cottaar <michiel.cottaar@ndcn.ox.ac.uk>
Date: Sun, 28 Jan 2024 14:05:49 +0000
Subject: [PATCH] pretty print fixed gradient and pulses

---
 Project.toml                     |  1 +
 src/gradients/fixed_gradients.jl | 37 +++++++++++++++++++++++++++++---
 src/pulses/fixed_pulses.jl       | 25 +++++++++++++++++----
 3 files changed, 56 insertions(+), 7 deletions(-)

diff --git a/Project.toml b/Project.toml
index 9c5ae50..321acf7 100644
--- a/Project.toml
+++ b/Project.toml
@@ -7,6 +7,7 @@ version = "0.0.0"
 Ipopt = "b6b21f68-93f8-5de0-b562-5493be1d77c9"
 JuMP = "4076af6c-e467-56ae-b986-b466b2749572"
 Juniper = "2ddba703-00a4-53a7-87a5-e8b9971dde84"
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 Polynomials = "f27b6e38-b328-58d1-80ce-0feddd5e7a45"
 Printf = "de0858da-6303-5e67-8744-51eddeeeb8d7"
 QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
diff --git a/src/gradients/fixed_gradients.jl b/src/gradients/fixed_gradients.jl
index 27d904a..5ee85f8 100644
--- a/src/gradients/fixed_gradients.jl
+++ b/src/gradients/fixed_gradients.jl
@@ -1,7 +1,9 @@
 module FixedGradients
 
-import ...BuildingBlocks: GradientBlock, fixed, BuildingBlock
-import ...Variables: variables, duration, qval
+import Printf: @sprintf
+import LinearAlgebra: norm
+import ...BuildingBlocks: GradientBlock, fixed, BuildingBlock, BuildingBlockPrinter
+import ...Variables: variables, duration, qval, gradient_strength, slew_rate
 
 
 """
@@ -49,7 +51,36 @@ variables(::Type{<:FixedGradient}) = []
 
 duration(fg::FixedGradient) = maximum(fg.time)
 
-Base.show(io::IO, fb::FixedGradient) = print(io, "FixedGradient for $(duration(fb)) ms")
+gradient_strength(fg::FixedGradient) = max(maximum(abs.(fg.Gx)), maximum(abs.(fg.Gy)), maximum(abs.(fg.Gz)))
+
+slew_rate(fg::FixedGradient) = maximum(map((fg.Gx, fg.Gy, fg.Gz)) do gradient
+    return maximum(abs.(gradient[2:end] - gradient[1:end-1]) ./ (fg.time[2:end] - fg.time[1:end-1]))
+end)
+
+qvec(fg::FixedGradient) = map((fg.Gx, fg.Gy, fg.Gz)) do gradient
+    weights_double = fg.time[2:end] - fg.time[1:end-1]
+    weights = [weights_double[1] / 2, ((weights_double[1:end-1] + weights_double[2:end]) / 2)..., weights_double[end]/2]
+    return sum(weights .* gradient)
+end
+
+qval(fg::FixedGradient) = norm(qvec(fg))
+
+
+function Base.show(io::IO, printer::BuildingBlockPrinter{<:FixedGradient})
+    fp = printer.bb
+    pp(value::Number) = @sprintf("%.3g", value)
+
+    print(io, "FixedGradient(")
+    if isnothing(printer.start_time)
+        print(io, "duration=", pp(duration(fp)))
+    else
+        print(io, "t=", pp(printer.start_time), "-", pp(printer.start_time + duration(fp)))
+    end
+    for fn in (qval, gradient_strength, slew_rate)
+        print(io, ", $(nameof(fn))=", pp(fn(fp)))
+    end
+    print(io, ")")
+end
 
 """
     FixedInstantGradient(orientation, qval)
diff --git a/src/pulses/fixed_pulses.jl b/src/pulses/fixed_pulses.jl
index 49de3fe..bbc6332 100644
--- a/src/pulses/fixed_pulses.jl
+++ b/src/pulses/fixed_pulses.jl
@@ -1,7 +1,8 @@
 module FixedPulses
 
+import Printf: @sprintf
 import ...BuildingBlocks: RFPulseBlock, fixed, BuildingBlockPrinter
-import ...Variables: variables, duration, amplitude, effective_time
+import ...Variables: variables, duration, amplitude, effective_time, flip_angle
 
 
 """
@@ -22,6 +23,7 @@ struct FixedPulse <: RFPulseBlock
     function FixedPulse(time::AbstractVector{<:Number}, amplitude::AbstractVector{<:Number}, phase::AbstractVector{<:Number})
         @assert length(time) == length(amplitude)
         @assert length(time) == length(phase)
+        @assert all(time .>= 0)
         new(Float64.(time), Float64.(amplitude), Float64.(phase))
     end
 end
@@ -30,14 +32,29 @@ function FixedPulse(time::AbstractVector{<:Number}, amplitude::AbstractVector{<:
     return FixedPulse(time, amplitude, (time .- time[1]) .* (frequency * 360) .+ phase)
 end
 
-variables(::Type{<:FixedPulse}) = []
+variables(::Type{<:FixedPulse}) = [duration, flip_angle, effective_time]
 
 duration(fp::FixedPulse) = maximum(fp.time)
 amplitude(fp::FixedPulse) = maximum(abs.(fp.amplitude))
-effective_time(pulse::FixedPulse) = pulse.time[argmax(abs(pulse.amplitude))]
+effective_time(pulse::FixedPulse) = pulse.time[argmax(abs.(pulse.amplitude))]
+function flip_angle(pulse::FixedPulse)
+    weights_double = pulse.time[2:end] - pulse.time[1:end-1]
+    weights = [weights_double[1] / 2, ((weights_double[1:end-1] + weights_double[2:end]) / 2)..., weights_double[end]/2]
+    return sum(weights .* pulse.amplitude) * 360
+end
 
-Base.show(io::IO, fp::FixedPulse) = print(io, "FixedPulse for $(duration(fp)) ms (peak at $(effective_time(fp)) ms)")
+function Base.show(io::IO, printer::BuildingBlockPrinter{<:FixedPulse})
+    fp = printer.bb
+    pp(value::Number) = @sprintf("%.3g", value)
 
+    print(io, "FixedPulse(")
+    if isnothing(printer.start_time)
+        print(io, "duration=", pp(duration(fp)))
+    else
+        print(io, "t=", pp(printer.start_time), "-", pp(printer.start_time + duration(fp)))
+    end
+    print(io, ", flip_angle=", pp(flip_angle(fp)), ", effective_time=", pp(effective_time(fp)), ")")
+end
 
 
 """
-- 
GitLab