Skip to content

Commit

Permalink
added interractivity for 1D inversion plots
Browse files Browse the repository at this point in the history
  • Loading branch information
aris committed Nov 6, 2024
1 parent e65109a commit 8d4d771
Show file tree
Hide file tree
Showing 5 changed files with 156 additions and 42 deletions.
4 changes: 2 additions & 2 deletions docs/src/functions.md
Original file line number Diff line number Diff line change
Expand Up @@ -127,8 +127,8 @@ We basically take the `plot()` function offered by GLMakie and extend it to type
For example, for 1D inversions, we have:

```@docs
plot(::NMRInversions.inv_out_1D)
plot!(::Union{Makie.Figure,Makie.GridPosition}, ::NMRInversions.inv_out_1D )
plot(res::inv_out_1D)
plot(res_mat::VecOrMat{inv_out_1D};kwargs...)
```

Expand Down
168 changes: 134 additions & 34 deletions ext/gui_ext.jl
Original file line number Diff line number Diff line change
Expand Up @@ -73,65 +73,165 @@ end
## Plots for 1D inversions

"""
plot(res::NMRInversions.inv_out_1D...; kwargs...)
plot(res_mat::VecOrMat{inv_out_1D};kwargs...)
Plot the results contained in a `inv_out_1D` structure.
This function can take any number of `inv_out_1D` structures as input.
Plot the results contained in a vector of `inv_out_1D` structures.
The keyword (optional) arguments are:
- `selections` : Whether to highlight the selections in the plots (default is `false`).
"""
function Makie.plot(res::NMRInversions.inv_out_1D...; kwargs...)
function Makie.plot(res_mat::AbstractVecOrMat{NMRInversions.inv_out_1D}; selections = false)

f = Figure(size=(800, 400))
plot!(f, res...; kwargs...)
fig = Figure(size=(700, 600))

return f
res = res_mat[1]
# Make axes
if res.seq in [NMRInversions.IR]
ax1 = Axis(fig[1:5, 1:5], xlabel="time (s)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[1:5, 6:10], xlabel="T₁ (s)", xscale=log10)
ax3 = Axis(fig[6:10,1:5], xlabel= "time (s)", ylabel="Residuals (a.u.)")

elseif res.seq in [NMRInversions.CPMG]
ax1 = Axis(fig[1:5, 1:5], xlabel="time (s)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[1:5, 6:10], xlabel="T₂ (s)", xscale=log10)
ax3 = Axis(fig[6:10,1:5], xlabel= "time (s)", ylabel="Residuals (a.u.)")

elseif res.seq in [NMRInversions.PFG]
ax1 = Axis(fig[1:5, 1:5], xlabel="b factor (s/m² e-9)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[1:5, 6:10], xlabel="D (m²/s)", xscale=log10)
ax3 = Axis(fig[6:10,1:5], xlabel= "b factor (s/m² e-9)", ylabel="Residuals (a.u.)")
end

linkxaxes!(ax1, ax3)

for r in res_mat
draw_on_axes(ax1, ax2,ax3, r, selections = selections)
end

return fig
end


"""
plot!(fig, res...; title)
plot(res::inv_out_1D)
Plot the results contained in a `inv_out_1D` structure on a figure or a grid position object.
This function can take any number of `inv_out_1D` structures as input.
Run the GUI to plot the 1D inversion results and select peaks you want to label.
"""
function Makie.plot(res::NMRInversions.inv_out_1D)

The arguments are:
fig = plot([res], selections = true)

- `fig` : The figure or grid position object.
- `res` : One or more `inv_out_1D` structures containing the fit results.
slider = IntervalSlider(fig[6,6:10], range = [1:length(res.X)...], startvalues = (1, length(res.X)))

Keyword (optional) arguments:
- `title` : Title of the plot.
"""
function Makie.plot!(fig::Union{Makie.Figure,Makie.GridPosition}, res::NMRInversions.inv_out_1D...;
title="")
int_low= lift(slider.interval) do i
return res.X[i[1]]
end

# Make axes
if res[1].seq in [NMRInversions.IR]
ax1 = Axis(fig[:, 1], xlabel="time (s)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[:, 2], xlabel="T₁ (s)", xscale=log10)
int_high= lift(slider.interval) do i
return res.X[i[2]]
end

elseif res[1].seq in [NMRInversions.CPMG]
ax1 = Axis(fig[:, 1], xlabel="time (s)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[:, 2], xlabel="T₂ (s)", xscale=log10)
weighted_average = lift(slider.interval) do i
return res.f[i[1]:i[2]]' * res.X[i[1]:i[2]] / sum(res.f[i[1]:i[2]])
end

elseif res[1].seq in [NMRInversions.PFG]
ax1 = Axis(fig[:, 1], xlabel="b factor (s/m² e-9)", ylabel="Signal (a.u.)")
ax2 = Axis(fig[:, 2], xlabel="D (m²/s)", xscale=log10)
labeltext = lift(slider.interval) do i
return "Selected range from $(round(res.X[i[1]], sigdigits=2)) to $(round(res.X[i[2]], sigdigits=2))
\n Weighted average = $(round(weighted_average[], sigdigits=2))"
end

Label(fig[7,6:10], labeltext, tellwidth=false)
vlines!(fig.content[2], int_low, color=:red)
vlines!(fig.content[2], int_high, color=:red)

button_label = Button(fig[8,6:10], label = "Save selection")
button_reset = Button(fig[9,6:10], label = "Reset selections")
button_save = Button(fig[10,6:10], label = "Save and exit")

on(button_label.clicks) do _
push!(res.selections, slider.interval[])
empty!(fig.content[1])
empty!(fig.content[2])
empty!(fig.content[3])
draw_on_axes(fig.content[1], fig.content[2], fig.content[3], res, selections = true)
vlines!(fig.content[2], int_low, color=:red)
vlines!(fig.content[2], int_high, color=:red)
end

for r in res
draw_on_axes(ax1, ax2, r)
on(button_reset.clicks) do _
empty!(res.selections)
empty!(fig.content[1])
empty!(fig.content[2])
empty!(fig.content[3])
draw_on_axes(fig.content[1], fig.content[2], fig.content[3], res)
vlines!(fig.content[2], int_low, color=:red)
vlines!(fig.content[2], int_high, color=:red)
end

on(button_save.clicks) do _
savedir = NMRInversions.save_file(res.title, filterlist = "png")
f = plot([res], selections = true)
save(savedir, f)
end

return fig
end

function draw_on_axes(ax1, ax2, res::NMRInversions.inv_out_1D)

scatter!(ax1, res.x, res.y)
lines!(ax1, res.xfit, res.yfit)
lines!(ax2, res.X, res.f)
function draw_on_axes(ax1, ax2, ax3, res::NMRInversions.inv_out_1D; selections = false)

c = length(ax2.scene.plots)

scatter!(ax1, res.x, res.y, colormap=:tab10, colorrange=(1, 10), color=c)
lines!(ax1, res.xfit, res.yfit, colormap=:tab10, colorrange=(1, 10), color=c)
lines!(ax2, res.X, res.f, colormap=:tab10, colorrange=(1, 10), color=c)
scatter!(ax3, res.x, res.r, colormap=:tab10, colorrange=(1, 10), color=c)
lines!(ax3, res.x, res.r, colormap=:tab10, colorrange=(1, 10), color=c)

if selections

if res.seq in [NMRInversions.IR]
Xlabel = ["<T₁> = ", " (s)"]
elseif res.seq in [NMRInversions.CPMG]
Xlabel = ["<T₂> = ", " (s)"]
elseif res.seq in [NMRInversions.PFG]
Xlabel = ["<D> = ", " (m²/s)"]
end
wa = weighted_averages(res)
for (i,s) in enumerate(res.selections)
clr = only(
Makie.numbers_to_colors(
[i],
Makie.to_colormap(:tab10),
identity,
Vec2(1, 10),
Makie.automatic,
Makie.automatic,
RGBAf(0, 0, 0, 0)
)
)

band!(
ax2,
res.X[s[1]:s[2]],
zeros(length(res.f[s[1]:s[2]])),
res.f[s[1]:s[2]],
color = clr, alpha = 0.5
)

height = (1 - (i-1) * 0.1) * maximum(res.f)
text!(
ax2, res.X[2], height,
text = Xlabel[1]*"$(round(wa[i], sigdigits = 2))"*Xlabel[2],
align = (:left, :top), color = clr
)

#=vlines!(ax2, res.X[s[1]], color = clr, alpha = 0.8)=#
#=vlines!(ax2, res.X[s[2]],color = clr, alpha = 0.8)=#
end
end
end


Expand Down
3 changes: 1 addition & 2 deletions src/inversion_functions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -82,9 +82,8 @@ function invert(seq::Type{<:pulse_sequence1D}, x::AbstractArray, y::Vector;
X .= X ./ 1e9
end

weighted_average = f'X / sum(f)

return inv_out_1D(seq, x, y, x_fit, y_fit, X, f, r, SNR, α, weighted_average,"")
return inv_out_1D(seq, x, y, x_fit, y_fit, X, f, r, SNR, α, [],"")

end

Expand Down
17 changes: 16 additions & 1 deletion src/misc.jl
Original file line number Diff line number Diff line change
Expand Up @@ -244,4 +244,19 @@ function compress_data(t_direct::AbstractVector, G::AbstractMatrix, bins::Int=64
usv1 = svd(sqrt(W0) * K1) #paper (13)
# surface(G, camera=(110, 25), xscale=:log10)

end
end


export weighted_averages
"""
weighted_averages(r::inv_out_1D)
Return a vector with the weighted averages for the selections in `r.selections`.
"""
function weighted_averages(r::inv_out_1D)

wa = Vector(undef, length(r.selections))
for (i,s) in enumerate(r.selections)
wa[i] = r.f[s[1]:s[2]]' * r.X[s[1]:s[2]] / sum(r.f[s[1]:s[2]])
end
return wa
end
6 changes: 3 additions & 3 deletions src/types.jl
Original file line number Diff line number Diff line change
Expand Up @@ -207,7 +207,7 @@ end

export inv_out_1D
"""
inv_out_1D(seq, x, y, xfit, yfit, X, f, r, SNR, α, wa)
inv_out_1D(seq, x, y, xfit, yfit, X, f, r, SNR, α, selections, title)
Output of the invert function for 1D pulse sequences.
A structure containing the following fields:
Expand All @@ -221,7 +221,7 @@ A structure containing the following fields:
- `r`, the residuals.
- `SNR`, the signal-to-noise ratio.
- `α`, the regularization parameter.
- `wa`, the weighted average of the inversion results (e.g. the mean relaxation time or diffusion coefficient).
- `selections`, a vector of tuples whose elements are indicate the first and last index of the selected peaks.
- `title`, a title describing the data.
"""
Expand All @@ -236,7 +236,7 @@ mutable struct inv_out_1D
r::Vector
SNR::Real
alpha::Real
wa::Real
selections::Vector{Tuple}
title::String
end

Expand Down

0 comments on commit 8d4d771

Please sign in to comment.