Skip to content
Snippets Groups Projects
Verified Commit 1b66a519 authored by Michiel Cottaar's avatar Michiel Cottaar
Browse files

Support new variables interface

parent c5cb1ed2
No related branches found
No related tags found
No related merge requests found
...@@ -129,7 +129,7 @@ for variable_func in keys(variables) ...@@ -129,7 +129,7 @@ for variable_func in keys(variables)
value = alt_var(bb) value = alt_var(bb)
catch e catch e
if e isa VariableNotAvailable if e isa VariableNotAvailable
throw(VariableNotAvailable(typeof(bb), variable_func)) throw(VariableNotAvailable(typeof(bb), Variables.$variable_func))
end end
rethrow() rethrow()
end end
...@@ -138,9 +138,9 @@ for variable_func in keys(variables) ...@@ -138,9 +138,9 @@ for variable_func in keys(variables)
elseif value isa AbstractArray{<:Number} elseif value isa AbstractArray{<:Number}
return backward.(value) return backward.(value)
end end
throw(VariableNotAvailable(typeof(bb), variable_func, alt_var)) throw(VariableNotAvailable(typeof(bb), Variables.$variable_func, alt_var))
end end
throw(VariableNotAvailable(typeof(bb), variable_func)) throw(VariableNotAvailable(typeof(bb), Variables.$variable_func))
end end
end end
...@@ -208,21 +208,21 @@ function Base.show(io::IO, printer::BuildingBlockPrinter) ...@@ -208,21 +208,21 @@ function Base.show(io::IO, printer::BuildingBlockPrinter)
end end
try try
numeric_value = _robust_value(fn(block)) numeric_value = _robust_value(fn(block))
if isnothing(numeric_value)
continue
end
if numeric_value isa AbstractVector
printed_value = "[" * join(map(v -> @sprintf("%.3g", v), numeric_value), ", ") * "]"
else
printed_value = @sprintf("%.3g", numeric_value)
end
print(io, "$(nameof(fn))=$(printed_value), ")
catch e catch e
if e isa VariableNotAvailable if e isa VariableNotAvailable
continue continue
end end
rethrow() rethrow()
end end
if isnothing(numeric_value)
continue
end
if numeric_value isa AbstractVector
printed_value = "[" * join(map(v -> @sprintf("%.3g", v), numeric_value), ", ") * "]"
else
printed_value = @sprintf("%.3g", numeric_value)
end
print(io, "$(nameof(fn))=$(printed_value), ")
end end
print(io, ")") print(io, ")")
end end
...@@ -239,8 +239,8 @@ If set to `:min` or `:max`, minimising or maximising this function will be added ...@@ -239,8 +239,8 @@ If set to `:min` or `:max`, minimising or maximising this function will be added
""" """
function set_simple_constraints!(block::BuildingBlock, kwargs) function set_simple_constraints!(block::BuildingBlock, kwargs)
for (key, value) in kwargs for (key, value) in kwargs
if key in keys(alternative_variables) if variables[key] in keys(alternative_variables)
alt_var, forward, backward, to_invert = alternative_variables[Variables.$variable_func] alt_var, forward, backward, to_invert = alternative_variables[variables[key]]
invert_value(value::VariableType) = forward(value) invert_value(value::VariableType) = forward(value)
invert_value(value::Symbol) = invert_value(Val(value)) invert_value(value::Symbol) = invert_value(Val(value))
invert_value(::Val{:min}) = to_invert ? Val(:max) : Val(:min) invert_value(::Val{:min}) = to_invert ? Val(:max) : Val(:min)
......
module Abstract module Abstract
import ...BuildingBlocks: BuildingBlock, ContainerBlock, RFPulseBlock, get_children_indices import ...BuildingBlocks: BuildingBlock, ContainerBlock, RFPulseBlock, get_children_indices, VariableNotAvailable
import ...Wait: WaitBlock import ...Wait: WaitBlock
import ...Readouts: InstantReadout import ...Readouts: InstantReadout
import ...Gradients: GradientBlock, split_gradient import ...Gradients: GradientBlock, split_gradient
...@@ -59,9 +59,9 @@ for (func_symbol, _) in [Dict(all_variables_symbols)[:pulse]..., :effective_time ...@@ -59,9 +59,9 @@ for (func_symbol, _) in [Dict(all_variables_symbols)[:pulse]..., :effective_time
@eval function $func_symbol(ao::AbstractOverlapping, args...; kwargs...) @eval function $func_symbol(ao::AbstractOverlapping, args...; kwargs...)
single_pulse = pulses(ao) single_pulse = pulses(ao)
if iszero(length(single_pulse)) if iszero(length(single_pulse))
error("Building block does not contain any RF pulse, so cannot compute $func_symbol.") throw(VariableNotAvailable(typeof(ao), $func_symbol))
elseif length(single_pulse) > 1 elseif length(single_pulse) > 1
error("Building block has multipl pulses, so cannot compute $func_symbol.") error("Building block has multipl pulses, so cannot compute ", $func_symbol, ".")
end end
$func_symbol(single_pulse[1], args...; kwargs...) $func_symbol(single_pulse[1], args...; kwargs...)
end end
......
module SingleLines module SingleLines
import ....BuildingBlocks: set_simple_constraints! import ....BuildingBlocks: set_simple_constraints!
import ....Readouts: ADC import ....Readouts: ADC
import ....Variables: dwell_time, fov, resolution, nsamples, effective_time, gradient_strength, slew_rate, voxel_size, rise_time, duration, variables import ....Variables: dwell_time, inverse_fov, resolution, nsamples, effective_time, gradient_strength, slew_rate, inverse_voxel_size, rise_time, duration, variables
import ...Abstract: AbstractOverlapping, waveform, interruptions import ...Abstract: AbstractOverlapping, waveform, interruptions
import ...GradientPulses: TrapezoidGradient import ...GradientPulses: TrapezoidGradient
...@@ -27,8 +27,8 @@ effective_time(sl::SingleLine) = rise_time(sl.grad) + effective_time(sl.adc) ...@@ -27,8 +27,8 @@ effective_time(sl::SingleLine) = rise_time(sl.grad) + effective_time(sl.adc)
dwell_time(sl::SingleLine) = dwell_time(sl.adc) dwell_time(sl::SingleLine) = dwell_time(sl.adc)
slew_rate(sl::SingleLine) = slew_rate(sl.grad)[1] slew_rate(sl::SingleLine) = slew_rate(sl.grad)[1]
gradient_strenth(sl::SingleLine) = slew_rate(sl) * rise_time(sl.grad) gradient_strenth(sl::SingleLine) = slew_rate(sl) * rise_time(sl.grad)
fov_inverse(sl::SingleLine) = 1e3 * dwell_time(sl) * gradient_strenth(sl) inverse_fov(sl::SingleLine) = 1e3 * dwell_time(sl) * gradient_strenth(sl)
voxel_size_inverse(sl::SingleLine) = 1e3 * duration(sl.adc) * gradient_strenth(sl) inverse_voxel_size(sl::SingleLine) = 1e3 * duration(sl.adc) * gradient_strenth(sl)
nsamples(sl::SingleLine) = nsamples(sl.adc) nsamples(sl::SingleLine) = nsamples(sl.adc)
resolution(sl::SingleLine) = resolution(sl.adc) resolution(sl::SingleLine) = resolution(sl.adc)
duration(sl::SingleLine) = duration(sl.grad) duration(sl::SingleLine) = duration(sl.grad)
......
module ADCs module ADCs
import JuMP: @constraint import JuMP: @constraint, value
import ...Variables: variables, dwell_time, duration, effective_time, get_free_variable, VariableType, nsamples, resolution, oversample import ...Variables: variables, dwell_time, duration, effective_time, get_free_variable, VariableType, nsamples, resolution, oversample
import ...BuildingBlocks: BuildingBlock, apply_simple_constraint!, set_simple_constraints!, fixed import ...BuildingBlocks: BuildingBlock, apply_simple_constraint!, set_simple_constraints!, fixed
import ...BuildSequences: global_model import ...BuildSequences: global_model
...@@ -53,9 +53,10 @@ resolution(adc::ADC) = adc.resolution ...@@ -53,9 +53,10 @@ resolution(adc::ADC) = adc.resolution
function fixed(adc::ADC) function fixed(adc::ADC)
# round nsamples during fixing # round nsamples during fixing
nsamples = Int(round(adc.nsamples)) r = Int(round(value(resolution(adc))))
dwell_time = duration(adc) / nsamples n = Int(round(value(nsamples(adc))))
return ADC(nsamples, dwell_time, time_to_center(adc)) oversample = n // r
return ADC(r, value(adc.dwell_time), value(adc.time_to_center), oversample)
end end
end end
\ No newline at end of file
...@@ -16,16 +16,19 @@ all_variables_symbols = [ ...@@ -16,16 +16,19 @@ all_variables_symbols = [
:amplitude => "The maximum amplitude of an RF pulse in kHz", :amplitude => "The maximum amplitude of an RF pulse in kHz",
:phase => "The angle of the phase of an RF pulse in KHz", :phase => "The angle of the phase of an RF pulse in KHz",
:frequency => "The off-resonance frequency of an RF pulse (relative to the Larmor frequency of water) in KHz", :frequency => "The off-resonance frequency of an RF pulse (relative to the Larmor frequency of water) in KHz",
:bandwidth => "Bandwidth of the RF pulse in kHz. To set constraints it is often better to use [`inverse_bandwidth`](@ref)", :bandwidth => "Bandwidth of the RF pulse in kHz. To set constraints it is often better to use [`inverse_bandwidth`](@ref).",
:inverse_bandwidth => "Inverse of bandwidth of the RF pulse in 1/kHz. Also see [`bandwidth`](@ref).",
:N_left => "The number of zero crossings of the RF pulse before the main peak", :N_left => "The number of zero crossings of the RF pulse before the main peak",
:N_right => "The number of zero crossings of the RF pulse after the main peak", :N_right => "The number of zero crossings of the RF pulse after the main peak",
:slice_thickness => "Slice thickness of an RF pulse that is active during a gradient in mm. To set constraints it is often better to use [`inverse_slice_thickness`](@ref).", :slice_thickness => "Slice thickness of an RF pulse that is active during a gradient in mm. To set constraints it is often better to use [`inverse_slice_thickness`](@ref).",
:inverse_slice_thickness => "Inverse of slice thickness of an RF pulse that is active during a gradient in 1/mm. Also, see [`slice_thickness`](@ref).",
], ],
# gradients # gradients
:gradient => [ :gradient => [
:qvec => "The spatial range and orientation on which the displacements can be detected due to this gradient in rad/um.", :qvec => "The spatial range and orientation on which the displacements can be detected due to this gradient in rad/um.",
:qval => "The spatial range on which the displacements can be detected due to this gradient in rad/um (i.e., norm of [`qvec`](@ref)).", :qval => "The spatial range on which the displacements can be detected due to this gradient in rad/um (i.e., norm of [`qvec`](@ref)). To set constraints it is often better to use [`qvec`](@ref) or [`qval_square`](@ref).",
:qval_square => "Square of [`qval`](@ref) in rad^2/um^2.",
:δ => "Effective duration of a gradient pulse ([`rise_time`](@ref) + [`flat_time`](@ref)) in ms.", :δ => "Effective duration of a gradient pulse ([`rise_time`](@ref) + [`flat_time`](@ref)) in ms.",
:rise_time => "Time for gradient pulse to reach its maximum value in ms.", :rise_time => "Time for gradient pulse to reach its maximum value in ms.",
:flat_time => "Time of gradient pulse at maximum value in ms.", :flat_time => "Time of gradient pulse at maximum value in ms.",
...@@ -35,8 +38,10 @@ all_variables_symbols = [ ...@@ -35,8 +38,10 @@ all_variables_symbols = [
:readout => [ :readout => [
:dwell_time => "Time between two samples in an `ADC` in ms.", :dwell_time => "Time between two samples in an `ADC` in ms.",
:nsamples => "Number of samples during a readout. During the optimisation this might produce non-integer values.", :nsamples => "Number of samples during a readout. During the optimisation this might produce non-integer values.",
:fov => "Size of the field of view in mm.", :fov => "Size of the field of view in mm. To set constraints it is often better to use [`inverse_fov`](@ref).",
:voxel_size => "Size of each voxel in mm.", :inverse_fov => "Inverse of size of the field of view in 1/mm. Also see [`fov`](@ref).",
:voxel_size => "Size of each voxel in mm. To set constraints it is often better to use [`inverse_voxel_size`](@ref).",
:inverse_voxel_size => "Inverse of voxel size in 1/mm. Also see [`voxel_size`](@ref).",
:resolution => "Number of voxels in the final readout. During the optimisation this might produce non-integer values, but this will be fixed after optimsation.", :resolution => "Number of voxels in the final readout. During the optimisation this might produce non-integer values, but this will be fixed after optimsation.",
:oversample => "How much to oversample with ([`nsamples`](@ref) / [`resolution`](@ref))", :oversample => "How much to oversample with ([`nsamples`](@ref) / [`resolution`](@ref))",
] ]
...@@ -60,54 +65,10 @@ end ...@@ -60,54 +65,10 @@ end
# helper functions # helper functions
"""
inverse_slice_thickness(pulse)
Inverse of [`slice_thickness`](@ref) in ms.
It is defined separately from `slice_thickness` to avoid unnecessary divisions in any constraint.
"""
function inverse_slice_thickness end
"""
inverse_bandwidth(pulse)
Inverse of [`bandwidth`](@ref) in 1/mm.
It is defined separately from `bandwidth` to avoid unnecessary divisions in any constraint.
"""
function inverse_bandwidth end
"""
inverse_fov(readout)
Inverse of [`fov`](@ref) in 1/mm.
It is defined separately from `fov` to avoid unnecessary divisions in any constraint.
"""
function inverse_fov end
"""
inverse_voxel_size(readout)
Inverse of [`voxel_size`](@ref) in 1/mm.
It is defined separately from `voxel_size` to avoid unnecessary divisions in any constraint.
"""
function inverse_voxel_size end
"""
qval_square(gradient)
Square of [`qval`](@ref) in rad/um.
It is defined separately from `qval` to avoid unnecessary square root in any constraint.
"""
function qval_square(bb, args...; kwargs...) function qval_square(bb, args...; kwargs...)
vec = qvec(bb, args...; kwargs...) vec = qvec(bb, args...; kwargs...)
return vec[1]^2 + vec[2]^2 + vec[3]^2 return vec[1]^2 + vec[2]^2 + vec[3]^2
end end
qval(bb; kwargs...) = sqrt(qval_square(bb))
alternative_variables = Dict( alternative_variables = Dict(
qval => (qval_square, n->n^2, sqrt, false), qval => (qval_square, n->n^2, sqrt, false),
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment